1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
2 * filename: m-stack.c *
3 * *
4 * UTIL C-source: Medical Image Conversion Utility *
5 * *
6 * purpose : stack files as specified *
7 * *
8 * project : (X)MedCon by Erik Nolf *
9 * *
10 * Functions : MdcGetNormSliceSpacing() - Get spacing between slices *
11 * MdcStackSlices() - Stack single slice image files*
12 * MdcStackFrames() - Stack multi slice volume files*
13 * MdcStackFiles() - Main stack routine *
14 * *
15 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
16 /*
17 */
18
19 /*
20 Copyright (C) 1997-2021 by Erik Nolf
21
22 This program is free software; you can redistribute it and/or modify it
23 under the terms of the GNU General Public License as published by the
24 Free Software Foundation; either version 2, or (at your option) any later
25 version.
26
27 This program is distributed in the hope that it will be useful, but
28 WITHOUT ANY WARRANTY; without even the implied warranty of
29 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
30 Public License for more details.
31
32 You should have received a copy of the GNU General Public License along
33 with this program; if not, write to the Free Software Foundation, Inc.,
34 59 Place - Suite 330, Boston, MA 02111-1307, USA. */
35
36 /****************************************************************************
37 H E A D E R S
38 ****************************************************************************/
39
40 #include <stdio.h>
41 #include <math.h>
42
43 #include "medcon.h"
44
45 /****************************************************************************
46 D E F I N E S
47 ****************************************************************************/
48
49 static FILEINFO infi, outfi;
50
51 static int mdc_nrstack=0;
52
53 /****************************************************************************
54 F U N C T I O N S
55 ****************************************************************************/
56
MdcGetNormSliceSpacing(IMG_DATA * id1,IMG_DATA * id2)57 float MdcGetNormSliceSpacing(IMG_DATA *id1, IMG_DATA *id2)
58 {
59 /* slice_spacing = sqrt( (x2-x1)^2 + (y2-y1)^2 + (z2-z1)^2 ) */
60 float value, slice_spacing;
61 float dx, dy, dz;
62
63 slice_spacing = id1->slice_spacing;
64
65 dx = id1->image_pos_pat[0] - id2->image_pos_pat[0];
66 dy = id1->image_pos_pat[1] - id2->image_pos_pat[1];
67 dz = id1->image_pos_pat[2] - id2->image_pos_pat[2];
68
69 value = (float)sqrt((double)(dx*dx + dy*dy + dz*dz));
70
71 if (fabs(slice_spacing - value) <= MDC_FLT_EPSILON) {
72 /* insignificant difference, use original header value */
73 slice_spacing = id1->slice_spacing;
74 }else{
75 /* significant, use calculated value from image_pos_pat[] info */
76 slice_spacing = (float)value;
77 }
78
79 return(slice_spacing);
80
81 }
82
83 /* tomo : stack single slice image files into one 3D volume file (2D+ -> 3D)*/
84 /* planar: stack planar slice image files into one single file */
MdcStackSlices(void)85 char *MdcStackSlices(void)
86 {
87 FILEINFO *ifi, *ofi;
88 IMG_DATA *id1, *id2;
89 DYNAMIC_DATA *dd1, *dd2;
90 Uint32 d, nr_of_images;
91 int HAS_DYNAMIC_DATA = MDC_NO;
92 char *msg=NULL;
93 int *total = mdc_arg_total; /* total arguments of files & conversions */
94 int *convs = mdc_arg_convs; /* counter for each conversion format */
95 char **files = mdc_arg_files; /* array of pointers to input filenames */
96 int i, convert, c;
97 float time_frame_duration=0.;
98
99 ifi = &infi; ofi= &outfi;
100
101 /* initialize output FILEINFO */
102 MdcInitFI(ofi,"stack3d");
103 nr_of_images = total[MDC_FILES];
104
105 if ((ifi->dynnr > 0) && (ifi->dyndata != NULL)) HAS_DYNAMIC_DATA = MDC_YES;
106
107 /* read and stack the several single slice files */
108 for (i=0; i<total[MDC_FILES]; i++) {
109
110 /* open file */
111 if (MdcOpenFile(ifi,files[i]) != MDC_OK) {
112 MdcCleanUpFI(ofi); return("stack slices : Failure to open file");
113 }
114
115 /* read the file */
116 if (MdcReadFile(ifi,i,NULL) != MDC_OK) {
117 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
118 return("stack slices : Failure to read file");
119 }
120
121 if (i == 0) {
122 /* copy FILEINFO stuff from 1st image */
123 MdcCopyFI(ofi,ifi,MDC_NO,MDC_NO);
124
125 /* set some specific parameters */
126 ofi->dim[0] = 3;
127 ofi->dim[1] = ifi->dim[1];
128 ofi->dim[2] = ifi->dim[2];
129 ofi->dim[3] = nr_of_images;
130 ofi->pixdim[0] = 3.;
131 ofi->pixdim[1] = ifi->pixdim[1];
132 ofi->pixdim[2] = ifi->pixdim[2];
133
134 if (ofi->planar == MDC_NO) ofi->acquisition_type = MDC_ACQUISITION_TOMO;
135
136 if (!MdcGetStructDD(ofi,1)) {
137 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
138 return("stack slices : Couldn't alloc output DYNAMIC_DATA structs");
139 }else{
140 ofi->dyndata[0].nr_of_slices = nr_of_images;
141 }
142
143 if (!MdcGetStructID(ofi,nr_of_images)) {
144 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
145 return("stack slices : Couldn't alloc output ING_DATA structs");
146 }
147
148 /* remember time_frame_duration */
149 if (HAS_DYNAMIC_DATA == MDC_YES) {
150 time_frame_duration = ifi->dyndata[0].time_frame_duration;
151 }
152
153 }else{
154
155 if (HAS_DYNAMIC_DATA == MDC_YES) {
156 dd1 = &ofi->dyndata[0];
157 dd2 = &ifi->dyndata[0];
158
159 /* check time_frame_duration differences */
160 if (time_frame_duration != dd2->time_frame_duration) {
161 MdcPrntWarn("stack slices : Different image durations found");
162 }
163
164 /* planar = increment total time_frame_duration */
165 if (ofi->planar == MDC_YES) {
166 dd1->time_frame_duration += dd2->time_frame_duration;
167 }
168
169 }
170
171 }
172
173 /* sanity checks */
174 for (d=3; d < MDC_MAX_DIMS; d++) if (ifi->dim[d] > 1 ) {
175 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
176 return("stack slices : Only single slice (one image) files supported");
177 }
178
179 if (ifi->dim[3] == 0) {
180 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
181 return("stack slices : File without image found");
182 }
183
184 /* copy IMG_DATA info */
185 msg = MdcCopyID(&ofi->image[i],&ifi->image[0],MDC_YES);
186 if (msg != NULL) {
187 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
188 sprintf(mdcbufr,"stack slices : %s",msg); return(mdcbufr);
189 }
190
191 /* small checks for file integrity */
192 if (i > 0) {
193 id1 = &ifi->image[0]; /* current slice */
194 id2 = &ofi->image[i-1]; /* previous slice */
195
196 if (ifi->pat_slice_orient != ofi->pat_slice_orient) {
197 MdcPrntWarn("stack slices : Different 'patient_slice_orient' found");
198 }
199
200 if ((id1->width != id2->width) || (id1->height != id2->height)) {
201 MdcPrntWarn("stack slices : Different image dimensions found");
202 }
203 if (id1->slice_width != id2->slice_width) {
204 MdcPrntWarn("stack slices : Different slice thickness found");
205 }
206 if (id1->slice_spacing != id2->slice_spacing) {
207 MdcPrntWarn("stack slices : Different slice spacing found");
208 }
209 if (id1->type != id2->type) {
210 MdcPrntWarn("stack slices : Different pixel type found");
211 }
212 }
213
214 MdcCleanUpFI(ifi);
215
216 }
217
218 /* check all the images */
219 msg = MdcImagesPixelFiddle(ofi);
220 if (msg != NULL) {
221 MdcCleanUpFI(ofi);
222 sprintf(mdcbufr,"stack slices : %s",msg); return(mdcbufr);
223 }
224
225 if (ofi->planar == MDC_NO) {
226 /* check for orthogonal slices */
227 switch (ofi->pat_slice_orient) {
228 case MDC_SUPINE_HEADFIRST_SAGITTAL :
229 case MDC_SUPINE_FEETFIRST_SAGITTAL :
230 case MDC_PRONE_HEADFIRST_SAGITTAL :
231 case MDC_PRONE_FEETFIRST_SAGITTAL :
232 case MDC_SUPINE_HEADFIRST_CORONAL :
233 case MDC_SUPINE_FEETFIRST_CORONAL :
234 case MDC_PRONE_HEADFIRST_CORONAL :
235 case MDC_PRONE_FEETFIRST_CORONAL :
236 case MDC_SUPINE_HEADFIRST_TRANSAXIAL :
237 case MDC_SUPINE_FEETFIRST_TRANSAXIAL :
238 case MDC_PRONE_HEADFIRST_TRANSAXIAL :
239 case MDC_PRONE_FEETFIRST_TRANSAXIAL :
240 case MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL :
241 case MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL :
242 case MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL :
243 case MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL :
244 case MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL :
245 case MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL :
246 case MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL :
247 case MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL :
248 case MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL :
249 case MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL :
250 case MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL :
251 case MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL : break;
252 default:
253 MdcPrntWarn("stack slices : Probably file with Non-Orthogonal slices");
254 }
255 }
256
257 /* correct slice_spacing */
258 for (i=1; i<nr_of_images; i++) {
259 id1 = &ofi->image[i];
260 id2 = &ofi->image[i-1];
261 id1->slice_spacing=MdcGetNormSliceSpacing(id1,id2);
262 }
263 /* and also for the first image */
264 if (nr_of_images > 1) {
265 ofi->image[0].slice_spacing = ofi->image[1].slice_spacing;
266 }
267
268 /* apply read options */
269 msg = MdcApplyReadOptions(ofi);
270 if (msg != NULL) {
271 MdcCleanUpFI(ofi);
272 sprintf(mdcbufr,"stack slices : %s",msg); return(mdcbufr);
273 }
274
275 /* if requested, reverse slices */
276 if (MDC_SORT_REVERSE == MDC_YES) {
277 msg = MdcSortReverse(ofi);
278 if (msg != NULL) {
279 MdcCleanUpFI(ofi);
280 sprintf(mdcbufr,"stack slices : %s",msg); return(mdcbufr);
281 }
282 }
283
284 /* write the file */
285 if (total[MDC_CONVS] > 0) {
286 /* go through conversion formats */
287 for (c=1; c<MDC_MAX_FRMTS; c++) {
288 convert = convs[c];
289 /* write output format when selected */
290 while (convert -- > 0) {
291 if (MdcWriteFile(ofi, c, mdc_nrstack++, NULL) != MDC_OK) {
292 MdcCleanUpFI(ofi);
293 return("stack slices : Failure to write file");
294 }
295 }
296 }
297 }
298
299 MdcCleanUpFI(ofi);
300
301 return(NULL);
302
303 }
304
305 /* tomo : stack volumes at different time frames into one 4D file (3D+ -> 4D)*/
306 /* planar: stack planar dynamic files into one planar dynamic file */
MdcStackFrames(void)307 char *MdcStackFrames(void)
308 {
309
310 FILEINFO *ifi, *ofi;
311 Uint32 d, nr_of_frames, nr_of_images=0;
312 char *msg = NULL;
313 int *total = mdc_arg_total; /* total arguments of files & conversions */
314 int *convs = mdc_arg_convs; /* counter for each conversion format */
315 char **files = mdc_arg_files; /* array of pointers to input filenames */
316 int i, j, f, convert, c;
317
318 ifi = &infi; ofi = &outfi;
319
320 /* initialize output FILEINFO */
321 MdcInitFI(ofi,"stack4d");
322 nr_of_frames = total[MDC_FILES];
323
324 for (i=0, j=0, f=0; f < total[MDC_FILES]; f++) {
325
326 /* open file */
327 if (MdcOpenFile(ifi,files[f]) != MDC_OK) {
328 MdcCleanUpFI(ofi);
329 return("stack frames : Failure to open file");
330 }
331
332 /* read the file */
333 if (MdcReadFile(ifi,f,NULL) != MDC_OK) {
334 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
335 return("stack frames : Failure to read file");
336 }
337 MdcCloseFile(ifi->ifp); /* no further need */
338
339 /* sanity checks */
340 for (d=4; d<MDC_MAX_DIMS; d++) if (ifi->dim[d] > 1) {
341 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
342 return("stack frames : Only tomo volumes or planar dynamic supported");
343 }
344 if ((ifi->dim[3] == 1) && (ifi->planar == MDC_NO)) {
345 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
346 return("stack frames : Use option '-stacks' for single slice files");
347 }
348 if (ifi->dim[3] == 0) {
349 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
350 return("stack frames : File without images found");
351 }
352
353 if (f == 0) {
354 /* copy FILEINFO stuff from 1st file */
355 MdcCopyFI(ofi,ifi,MDC_NO,MDC_NO);
356
357 /* 4D -> dynamic */
358 ofi->acquisition_type = MDC_ACQUISITION_DYNAMIC;
359
360 /* get appropriate structs */
361 if (!MdcGetStructDD(ofi,nr_of_frames)) {
362 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
363 return("stack frames : Couldn't alloc output DYNAMIC_DATA structs");
364 }
365
366 /* set some specific parameters */
367 if (ofi->planar == MDC_YES) {
368 /* planar: asymmetric */
369 nr_of_images= ifi->number; /* increment */
370 ofi->dim[0] = 3;
371 ofi->dim[1] = ifi->dim[1];
372 ofi->dim[2] = ifi->dim[2];
373 ofi->dim[3] = nr_of_images;
374 ofi->pixdim[0] = 3.;
375 ofi->pixdim[1] = ifi->pixdim[1];
376 ofi->pixdim[2] = ifi->pixdim[2];
377 ofi->pixdim[3] = ifi->pixdim[3];
378 }else{
379 /* tomo : symmectric */
380 nr_of_images= ifi->number * nr_of_frames;
381 ofi->dim[0] = 4;
382 ofi->dim[1] = ifi->dim[1];
383 ofi->dim[2] = ifi->dim[2];
384 ofi->dim[3] = ifi->dim[3];
385 ofi->dim[4] = nr_of_frames;
386 ofi->pixdim[0] = 4.;
387 ofi->pixdim[1] = ifi->pixdim[1];
388 ofi->pixdim[2] = ifi->pixdim[2];
389 ofi->pixdim[3] = ifi->pixdim[3];
390 ofi->pixdim[4] = ofi->dyndata[0].time_frame_duration; /* tomo */
391 }
392
393 /* malloc all IMG_DATA structs */
394 if (!MdcGetStructID(ofi,nr_of_images)) {
395 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
396 return("stack frames : Couldn't alloc output IMG_DATA structs");
397 }
398
399 }else{
400
401 if (ofi->planar == MDC_YES) {
402 nr_of_images += ifi->number;
403 ofi->dim[3] = nr_of_images;
404 /* malloc IMG_DATA structs current frame */
405 if (!MdcGetStructID(ofi,nr_of_images)) {
406 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
407 return("stack frames : Couldn't alloc planar IMG_DATA structs");
408 }
409 }
410
411 /* copy DYNAMIC_DATA struct when available */
412 /* f=0 already copied at initial MdcCopyFI() */
413 if ((ifi->dynnr > 0) && (ifi->dyndata != NULL)) {
414 MdcCopyDD(&ofi->dyndata[f],&ifi->dyndata[0]);
415 }
416
417 /* suspectable differences */
418 if (ifi->pat_slice_orient != ofi->pat_slice_orient) {
419 MdcPrntWarn("stack frames : Different 'patient_slice_orient' found");
420 }
421 if (ifi->planar != ofi->planar) {
422 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
423 return("stack frames : wrongful mixture of tomo and planar frames");
424 }
425 }
426
427 for (i=0; i<ifi->dim[3]; i++, j++) {
428 /* copy IMG_DATA info */
429 msg = MdcCopyID(&ofi->image[j],&ifi->image[i],MDC_YES);
430 if (msg != NULL) {
431 MdcCleanUpFI(ofi); MdcCleanUpFI(ifi);
432 sprintf(mdcbufr,"stack frames : %s",msg); return(mdcbufr);
433 }
434 }
435
436 MdcCleanUpFI(ifi);
437
438 }
439
440 /* check all the images */
441 msg = MdcImagesPixelFiddle(ofi);
442 if (msg != NULL) {
443 MdcCleanUpFI(ofi);
444 sprintf(mdcbufr,"stack frames : %s",msg); return(mdcbufr);
445 }
446
447 if (ofi->planar == MDC_NO) {
448 /* check for orthogonal slices */
449 switch (ofi->pat_slice_orient) {
450 case MDC_SUPINE_HEADFIRST_SAGITTAL :
451 case MDC_SUPINE_FEETFIRST_SAGITTAL :
452 case MDC_PRONE_HEADFIRST_SAGITTAL :
453 case MDC_PRONE_FEETFIRST_SAGITTAL :
454 case MDC_SUPINE_HEADFIRST_CORONAL :
455 case MDC_SUPINE_FEETFIRST_CORONAL :
456 case MDC_PRONE_HEADFIRST_CORONAL :
457 case MDC_PRONE_FEETFIRST_CORONAL :
458 case MDC_SUPINE_HEADFIRST_TRANSAXIAL :
459 case MDC_SUPINE_FEETFIRST_TRANSAXIAL :
460 case MDC_PRONE_HEADFIRST_TRANSAXIAL :
461 case MDC_PRONE_FEETFIRST_TRANSAXIAL :
462 case MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL :
463 case MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL :
464 case MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL :
465 case MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL :
466 case MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL :
467 case MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL :
468 case MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL :
469 case MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL :
470 case MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL :
471 case MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL :
472 case MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL :
473 case MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL : break;
474 default:
475 MdcPrntWarn("stack frames : Probably file with Non-Orthogonal slices");
476 }
477 }
478
479 /* apply read options */
480 msg = MdcApplyReadOptions(ofi);
481 if (msg != NULL) {
482 MdcCleanUpFI(ofi);
483 sprintf(mdcbufr,"stack frames : %s",msg); return(mdcbufr);
484 }
485
486 /* write the file */
487 if (total[MDC_CONVS] > 0) {
488 /* go through conversion formats */
489 for (c=1; c<MDC_MAX_FRMTS; c++) {
490 convert = convs[c];
491 /* write output format when selected */
492 while (convert -- > 0) {
493 if (MdcWriteFile(ofi, c, mdc_nrstack++, NULL) != MDC_OK) {
494 MdcCleanUpFI(ofi);
495 return("stack frames : Failure to write file");
496 }
497 }
498 }
499 }
500
501 MdcCleanUpFI(ofi);
502
503 return(NULL);
504
505 }
506
MdcStackFiles(Int8 stack)507 char *MdcStackFiles(Int8 stack)
508 {
509 char *msg=NULL;
510
511 if (MDC_CONVERT != MDC_YES)
512 return("In order to stack specify an output format");
513
514 if (mdc_arg_total[MDC_FILES] == 1)
515 return("In order to stack at least two files are required");
516
517 switch (stack) {
518 case MDC_STACK_SLICES: msg = MdcStackSlices();
519 break;
520 case MDC_STACK_FRAMES: msg = MdcStackFrames();
521 break;
522 }
523
524 return(msg);
525 }
526