1 #if HAVE_CONFIG_H
2 #include "config.h"
3 #endif
4
5 #include <limits.h>
6 #include <float.h>
7 #include <minc.h>
8 #include <ParseArgv.h>
9 #include <volume_io.h>
10
11 #undef X /* These are used in nifti1_io */
12 #undef Y
13 #undef Z
14
15 #include <time_stamp.h>
16 #include "nifti1_io.h"
17 #include "analyze75.h"
18 #include "nifti1_local.h"
19
20
test_xform(mat44 m,int i,int j,int k)21 void test_xform(mat44 m, int i, int j, int k)
22 {
23 double x, y, z;
24
25 x = m.m[DIM_X][DIM_I] * i + m.m[DIM_X][DIM_J] * j + m.m[DIM_X][DIM_K] * k
26 + m.m[DIM_X][3];
27 y = m.m[DIM_Y][DIM_I] * i + m.m[DIM_Y][DIM_J] * j + m.m[DIM_Y][DIM_K] * k
28 + m.m[DIM_Y][3];
29 z = m.m[DIM_Z][DIM_I] * i + m.m[DIM_Z][DIM_J] * j + m.m[DIM_Z][DIM_K] * k
30 + m.m[DIM_Z][3];
31
32 printf("%d %d %d => ", i, j, k);
33 printf("%f %f %f\n", x, y, z);
34 }
35
usage(void)36 static int usage(void)
37 {
38 static const char msg[] = {
39 "nii2mnc: Convert NIfTI-1 files to MINC format\n"
40 "usage: nii2mnc [options] filename.nii [filename.mnc]\n"
41 };
42 fprintf(stderr, "%s", msg);
43 return (-1);
44 }
45
find_data_range(int datatype,long nvox,void * data,double range[2])46 static void find_data_range(int datatype,
47 long nvox,
48 void *data,
49 double range[2])
50 {
51 int i;
52
53 range[0] = DBL_MAX;
54 range[1] = -DBL_MAX;
55
56 for (i = 0; i < nvox; i++) {
57 double tmp;
58
59 switch (datatype) {
60 case DT_INT8:
61 tmp = (double) ((char *)data)[i];
62 break;
63 case DT_UINT8:
64 tmp = (double) ((unsigned char *)data)[i];
65 break;
66 case DT_INT16:
67 tmp = (double) ((short *)data)[i];
68 break;
69 case DT_UINT16:
70 tmp = (double) ((unsigned short *)data)[i];
71 break;
72 case DT_INT32:
73 tmp = (double) ((int *)data)[i];
74 break;
75 case DT_UINT32:
76 tmp = (double) ((unsigned int *)data)[i];
77 break;
78 case DT_FLOAT32:
79 tmp = (double) ((float *)data)[i];
80 break;
81 case DT_FLOAT64:
82 tmp = (double) ((double *)data)[i];
83 break;
84 default:
85 fprintf(stderr, "Data type %d not handled\n", datatype);
86 break;
87 }
88 if (tmp < range[0]) {
89 range[0] = tmp;
90 }
91 if (tmp > range[1]) {
92 range[1] = tmp;
93 }
94 }
95 }
96
97 int
main(int argc,char ** argv)98 main(int argc, char **argv)
99 {
100 /* NIFTI stuff */
101 nifti_image *nii_ptr;
102
103 /* MINC stuff */
104 int mnc_fd; /* MINC file descriptor */
105 nc_type mnc_mtype; /* MINC memory data type */
106 int mnc_msign; /* MINC !0 if signed data */
107 static nc_type mnc_vtype; /* MINC voxel data type */
108 static int mnc_vsign; /* MINC !0 if signed data */
109 int mnc_ndims; /* MINC image dimension count */
110 int mnc_dimids[MAX_VAR_DIMS]; /* MINC image dimension identifiers */
111 int mnc_iid; /* MINC Image variable ID */
112 long mnc_start[MAX_VAR_DIMS]; /* MINC data starts */
113 long mnc_count[MAX_VAR_DIMS]; /* MINC data counts */
114 char *mnc_hist; /* MINC history */
115 double mnc_vrange[2]; /* MINC valid min/max */
116 double mnc_srange[2]; /* MINC image min/max */
117 double mnc_time_step;
118 double mnc_time_start;
119 int mnc_spatial_axes[MAX_NII_DIMS];
120 double mnc_starts[MAX_SPACE_DIMS];
121 double mnc_steps[MAX_SPACE_DIMS];
122 double mnc_dircos[MAX_SPACE_DIMS][MAX_SPACE_DIMS];
123 Transform mnc_xform;
124 General_transform mnc_linear_xform;
125 int mnc_vid; /* Dimension variable id */
126 struct analyze75_hdr ana_hdr;
127
128 /* Other stuff */
129 char out_str[1024]; /* Big string for filename */
130 int i; /* Generic loop counter the first */
131 int j; /* Generic loop counter the second */
132 char *str_ptr; /* Generic ASCIZ string pointer */
133 int r; /* Result code. */
134 static int qflag = 0; /* Quiet flag (default is non-quiet) */
135 static int rflag = 1; /* Scan range flag */
136 static int oflag = DIMORDER_ZYX;
137 static int flip[MAX_SPACE_DIMS];
138
139 static char *mnc_ordered_dim_names[MAX_SPACE_DIMS];
140
141 static ArgvInfo argTable[] = {
142 {"-byte", ARGV_CONSTANT, (char *) NC_BYTE, (char *)&mnc_vtype,
143 "Write voxel data in 8-bit integer format."},
144 {"-short", ARGV_CONSTANT, (char *) NC_SHORT, (char *)&mnc_vtype,
145 "Write voxel data in 16-bit integer format."},
146 {"-int", ARGV_CONSTANT, (char *) NC_INT, (char *)&mnc_vtype,
147 "Write voxel data in 32-bit integer format."},
148 {"-float", ARGV_CONSTANT, (char *) NC_FLOAT, (char *)&mnc_vtype,
149 "Write voxel data in 32-bit floating point format."},
150 {"-double", ARGV_CONSTANT, (char *) NC_DOUBLE, (char *)&mnc_vtype,
151 "Write voxel data in 64-bit floating point format."},
152 {"-signed", ARGV_CONSTANT, (char *) 1, (char *)&mnc_vsign,
153 "Write integer voxel data in signed format."},
154 {"-unsigned", ARGV_CONSTANT, (char *) 0, (char *)&mnc_vsign,
155 "Write integer voxel data in unsigned format."},
156 {"-noscanrange", ARGV_CONSTANT, (char *) 0, (char *)&rflag,
157 "Do not scan data range."},
158 {"-quiet", ARGV_CONSTANT, (char *) 0, (char *)&qflag,
159 "Quiet operation"},
160 {"-transverse", ARGV_CONSTANT, (char *) DIMORDER_ZYX, (char *) &oflag,
161 "Assume transverse (ZYX) ordering for spatial dimensions"},
162 {"-sagittal", ARGV_CONSTANT, (char *) DIMORDER_XZY, (char *) &oflag,
163 "Assume sagittal (XZY) ordering for spatial dimensions"},
164 {"-coronal", ARGV_CONSTANT, (char *) DIMORDER_YZX, (char *) &oflag,
165 "Assume coronal (YZX) ordering for spatial dimensions"},
166 {"-xyz", ARGV_CONSTANT, (char *) DIMORDER_XYZ, (char *) &oflag,
167 "Assume XYZ ordering for spatial dimensions"},
168 {"-zxy", ARGV_CONSTANT, (char *) DIMORDER_ZXY, (char *) &oflag,
169 "Assume ZXY ordering for spatial dimensions"},
170 {"-yxz", ARGV_CONSTANT, (char *) DIMORDER_YXZ, (char *) &oflag,
171 "Assume YXZ ordering for spatial dimensions"},
172 {"-flipx", ARGV_CONSTANT, (char *) 1, (char *)&flip[DIM_X],
173 "Invert direction of X (left-right) axis"},
174 {"-flipy", ARGV_CONSTANT, (char *) 1, (char *)&flip[DIM_Y],
175 "Invert direction of Y (posterior-anterior) axis"},
176 {"-flipz", ARGV_CONSTANT, (char *) 1, (char *)&flip[DIM_Z],
177 "Invert direction of Z (inferior-superior) axis"},
178 {NULL, ARGV_END, NULL, NULL, NULL}
179 };
180
181
182 ncopts = 0; /* Clear global netCDF error reporting flag */
183
184 mnc_hist = time_stamp(argc, argv);
185 mnc_vtype = NC_NAT;
186
187 if (ParseArgv(&argc, argv, argTable, 0) || (argc < 2)) {
188 fprintf(stderr, "Too few arguments\n");
189 return usage();
190 }
191
192 if (argc == 2) {
193 strcpy(out_str, argv[1]);
194 str_ptr = strrchr(out_str, '.');
195 if (str_ptr != NULL) {
196 if (!strcmp(str_ptr, ".nii") || !strcmp(str_ptr, ".hdr")) {
197 *str_ptr = '\0';
198 strcat(out_str, ".mnc");
199 }
200 }
201 }
202 else if (argc == 3) {
203 strcpy(out_str, argv[2]);
204 }
205 else {
206 fprintf(stderr, "Filename argument required\n");
207 return usage();
208 }
209
210 /* Read in the entire NIfTI file. */
211 nii_ptr = nifti_image_read(argv[1], 1);
212
213 if (nii_ptr->nifti_type == 0) { /* Analyze file!!! */
214 FILE *fp;
215 int ss;
216 int must_swap;
217
218 fp = fopen(argv[1], "rb");
219 if (fp != NULL) {
220 fread(&ana_hdr, sizeof (ana_hdr), 1, fp);
221 fclose(fp);
222
223 must_swap = 0;
224 ss = ana_hdr.dime.dim[0];
225 if (ss != 0) {
226 if (ss < 0 || ss > 7) {
227 nifti_swap_2bytes(1, &(ss));
228 if (ss < 0 || ss > 7) {
229 /* We should never get here!! */
230 fprintf(stderr, "Bad dimension count!!\n");
231 }
232 else {
233 must_swap = 1;
234 }
235 }
236 }
237 else {
238 ss = ana_hdr.hk.sizeof_hdr;
239 if (ss != sizeof(ana_hdr)) {
240 nifti_swap_4bytes(1, &(ss));
241 if (ss != sizeof(ana_hdr)) {
242 /* We should never get here!! */
243 fprintf(stderr, "Bad header size!!\n");
244 }
245 else {
246 must_swap = 1;
247 }
248 }
249 }
250
251 if (must_swap) {
252 nifti_swap_4bytes(1, &(ana_hdr.hk.sizeof_hdr));
253 nifti_swap_4bytes(1, &(ana_hdr.hk.extents));
254 nifti_swap_2bytes(1, &(ana_hdr.hk.session_error));
255
256 nifti_swap_4bytes(1, &(ana_hdr.dime.compressed));
257 nifti_swap_4bytes(1, &(ana_hdr.dime.verified));
258 nifti_swap_4bytes(1, &(ana_hdr.dime.glmax));
259 nifti_swap_4bytes(1, &(ana_hdr.dime.glmin));
260 nifti_swap_2bytes(8, ana_hdr.dime.dim);
261 nifti_swap_4bytes(8, ana_hdr.dime.pixdim);
262 nifti_swap_2bytes(1, &(ana_hdr.dime.datatype));
263 nifti_swap_2bytes(1, &(ana_hdr.dime.bitpix));
264 nifti_swap_4bytes(1, &(ana_hdr.dime.vox_offset));
265 nifti_swap_4bytes(1, &(ana_hdr.dime.cal_max));
266 nifti_swap_4bytes(1, &(ana_hdr.dime.cal_min));
267 }
268
269 if (!qflag) {
270 printf("orient = %d\n", ana_hdr.hist.orient);
271 }
272 }
273 }
274
275 if (!qflag) {
276 nifti_image_infodump(nii_ptr);
277 }
278
279 /* Set up the spatial axis correspondence for the call to
280 * convert_transform_to_starts_and_steps()
281 */
282 switch (oflag) {
283 default:
284 case DIMORDER_ZYX:
285 mnc_ordered_dim_names[DIM_X] = MIxspace;
286 mnc_ordered_dim_names[DIM_Y] = MIyspace;
287 mnc_ordered_dim_names[DIM_Z] = MIzspace;
288 mnc_spatial_axes[DIM_X] = DIM_X;
289 mnc_spatial_axes[DIM_Y] = DIM_Y;
290 mnc_spatial_axes[DIM_Z] = DIM_Z;
291 break;
292 case DIMORDER_ZXY:
293 mnc_ordered_dim_names[DIM_X] = MIyspace;
294 mnc_ordered_dim_names[DIM_Y] = MIxspace;
295 mnc_ordered_dim_names[DIM_Z] = MIzspace;
296 mnc_spatial_axes[DIM_X] = DIM_Y;
297 mnc_spatial_axes[DIM_Y] = DIM_X;
298 mnc_spatial_axes[DIM_Z] = DIM_Z;
299 break;
300 case DIMORDER_XYZ:
301 mnc_ordered_dim_names[DIM_X] = MIzspace;
302 mnc_ordered_dim_names[DIM_Y] = MIyspace;
303 mnc_ordered_dim_names[DIM_Z] = MIxspace;
304 mnc_spatial_axes[DIM_X] = DIM_Z;
305 mnc_spatial_axes[DIM_Y] = DIM_Y;
306 mnc_spatial_axes[DIM_Z] = DIM_X;
307 break;
308 case DIMORDER_XZY:
309 mnc_ordered_dim_names[DIM_X] = MIyspace;
310 mnc_ordered_dim_names[DIM_Y] = MIzspace;
311 mnc_ordered_dim_names[DIM_Z] = MIxspace;
312 mnc_spatial_axes[DIM_X] = DIM_Y;
313 mnc_spatial_axes[DIM_Y] = DIM_Z;
314 mnc_spatial_axes[DIM_Z] = DIM_X;
315 break;
316 case DIMORDER_YZX:
317 mnc_ordered_dim_names[DIM_X] = MIxspace;
318 mnc_ordered_dim_names[DIM_Y] = MIzspace;
319 mnc_ordered_dim_names[DIM_Z] = MIyspace;
320 mnc_spatial_axes[DIM_X] = DIM_X;
321 mnc_spatial_axes[DIM_Y] = DIM_Z;
322 mnc_spatial_axes[DIM_Z] = DIM_Y;
323 break;
324 case DIMORDER_YXZ:
325 mnc_ordered_dim_names[DIM_X] = MIzspace;
326 mnc_ordered_dim_names[DIM_Y] = MIxspace;
327 mnc_ordered_dim_names[DIM_Z] = MIyspace;
328 mnc_spatial_axes[DIM_X] = DIM_Z;
329 mnc_spatial_axes[DIM_Y] = DIM_X;
330 mnc_spatial_axes[DIM_Z] = DIM_Y;
331 break;
332 }
333
334 switch (nii_ptr->datatype) {
335 case DT_INT8:
336 mnc_msign = 1;
337 mnc_mtype = NC_BYTE;
338 mnc_vrange[0] = CHAR_MIN;
339 mnc_vrange[1] = CHAR_MAX;
340 break;
341 case DT_UINT8:
342 mnc_msign = 0;
343 mnc_mtype = NC_BYTE;
344 mnc_vrange[0] = 0;
345 mnc_vrange[1] = UCHAR_MAX;
346 break;
347 case DT_INT16:
348 mnc_msign = 1;
349 mnc_mtype = NC_SHORT;
350 mnc_vrange[0] = SHRT_MIN;
351 mnc_vrange[1] = SHRT_MAX;
352 break;
353 case DT_UINT16:
354 mnc_msign = 0;
355 mnc_mtype = NC_SHORT;
356 mnc_vrange[0] = 0;
357 mnc_vrange[1] = USHRT_MAX;
358 break;
359 case DT_INT32:
360 mnc_msign = 1;
361 mnc_mtype = NC_INT;
362 mnc_vrange[0] = INT_MIN;
363 mnc_vrange[1] = INT_MAX;
364 break;
365 case DT_UINT32:
366 mnc_msign = 0;
367 mnc_mtype = NC_INT;
368 mnc_vrange[0] = 0;
369 mnc_vrange[1] = UINT_MAX;
370 break;
371 case DT_FLOAT32:
372 mnc_msign = 1;
373 mnc_mtype = NC_FLOAT;
374 mnc_vrange[0] = -FLT_MAX;
375 mnc_vrange[1] = FLT_MAX;
376 break;
377 case DT_FLOAT64:
378 mnc_msign = 1;
379 mnc_mtype = NC_DOUBLE;
380 mnc_vrange[0] = -DBL_MAX;
381 mnc_vrange[1] = DBL_MAX;
382 break;
383 default:
384 fprintf(stderr, "Data type %d not handled\n", nii_ptr->datatype);
385 break;
386 }
387
388 if (mnc_vtype == NC_NAT) {
389 mnc_vsign = mnc_msign;
390 mnc_vtype = mnc_mtype;
391 }
392
393 /* Open the MINC file. It should not already exist.
394 */
395 mnc_fd = micreate(out_str, NC_NOCLOBBER);
396 if (mnc_fd < 0) {
397 fprintf(stderr, "Can't create output file '%s'\n", out_str);
398 return (-1);
399 }
400
401 /* Create the necessary dimensions in the minc file
402 */
403
404 mnc_ndims = 0;
405
406 if (nii_ptr->nt > 1) {
407 mnc_dimids[mnc_ndims] = ncdimdef(mnc_fd, MItime, nii_ptr->nt);
408 mnc_start[mnc_ndims] = 0;
409 mnc_count[mnc_ndims] = nii_ptr->nt;
410 mnc_ndims++;
411
412 r = micreate_std_variable(mnc_fd, MItime, NC_INT, 0, NULL);
413 switch (nii_ptr->time_units) {
414 case NIFTI_UNITS_UNKNOWN:
415 case NIFTI_UNITS_SEC:
416 mnc_time_step = nii_ptr->dt;
417 mnc_time_start = nii_ptr->toffset;
418 break;
419 case NIFTI_UNITS_MSEC:
420 mnc_time_step = nii_ptr->dt / 1000;
421 mnc_time_start = nii_ptr->toffset / 1000;
422 break;
423 case NIFTI_UNITS_USEC:
424 mnc_time_step = nii_ptr->dt / 1000000;
425 mnc_time_start = nii_ptr->toffset / 1000000;
426 break;
427 default:
428 fprintf(stderr, "Unknown time units value %d\n",
429 nii_ptr->time_units);
430 break;
431 }
432 miattputdbl(mnc_fd, r, MIstart, mnc_time_start);
433 miattputdbl(mnc_fd, r, MIstep, mnc_time_step);
434 miattputstr(mnc_fd, r, MIunits, "s");
435 }
436
437 if (nii_ptr->nz > 1) {
438 mnc_dimids[mnc_ndims] = ncdimdef(mnc_fd, mnc_ordered_dim_names[DIM_Z],
439 nii_ptr->nz);
440 mnc_start[mnc_ndims] = 0;
441 mnc_count[mnc_ndims] = nii_ptr->nz;
442 mnc_ndims++;
443
444 r = micreate_std_variable(mnc_fd, mnc_ordered_dim_names[DIM_Z], NC_INT,
445 0, NULL);
446 miattputdbl(mnc_fd, r, MIstep, nii_ptr->dz);
447 miattputstr(mnc_fd, r, MIunits, "mm");
448 }
449
450 if (nii_ptr->ny > 1) {
451 mnc_dimids[mnc_ndims] = ncdimdef(mnc_fd, mnc_ordered_dim_names[DIM_Y],
452 nii_ptr->ny);
453 mnc_start[mnc_ndims] = 0;
454 mnc_count[mnc_ndims] = nii_ptr->ny;
455 mnc_ndims++;
456
457 r = micreate_std_variable(mnc_fd, mnc_ordered_dim_names[DIM_Y], NC_INT,
458 0, NULL);
459 miattputdbl(mnc_fd, r, MIstep, nii_ptr->dy);
460 miattputstr(mnc_fd, r, MIunits, "mm");
461 }
462
463 if (nii_ptr->nx > 1) {
464 mnc_dimids[mnc_ndims] = ncdimdef(mnc_fd, mnc_ordered_dim_names[DIM_X],
465 nii_ptr->nx);
466 mnc_start[mnc_ndims] = 0;
467 mnc_count[mnc_ndims] = nii_ptr->nx;
468 mnc_ndims++;
469
470 r = micreate_std_variable(mnc_fd, mnc_ordered_dim_names[DIM_X], NC_INT,
471 0, NULL);
472 miattputdbl(mnc_fd, r, MIstep, nii_ptr->dx);
473 miattputstr(mnc_fd, r, MIunits, "mm");
474 }
475
476 if (nii_ptr->nu > 1) {
477 mnc_dimids[mnc_ndims] = ncdimdef(mnc_fd, MIvector_dimension,
478 nii_ptr->nu);
479 mnc_start[mnc_ndims] = 0;
480 mnc_count[mnc_ndims] = nii_ptr->nu;
481 mnc_ndims++;
482 }
483
484 /* Create scalar image-min and image-max variables.
485 */
486 micreate_std_variable(mnc_fd, MIimagemax, NC_DOUBLE, 0, NULL);
487 micreate_std_variable(mnc_fd, MIimagemin, NC_DOUBLE, 0, NULL);
488
489 /* Create the group variables.
490 */
491 micreate_std_variable(mnc_fd, MIstudy, NC_INT, 0, NULL);
492 if (strlen(nii_ptr->descrip) > 0 && strlen(nii_ptr->descrip) < 79 ) {
493 int varid = micreate_std_variable(mnc_fd, MIpatient, NC_INT, 0, NULL);
494 (void) miattputstr(mnc_fd, varid, MIfull_name,
495 nii_ptr->descrip);
496 } else {
497 micreate_std_variable(mnc_fd, MIpatient, NC_INT, 0, NULL);
498 }
499 micreate_std_variable(mnc_fd, MIacquisition, NC_INT, 0, NULL);
500
501 /* Create the MINC image variable. If we can't, there is no
502 * further processing possible...
503 */
504 mnc_iid = micreate_std_variable(mnc_fd, MIimage, mnc_vtype, mnc_ndims,
505 mnc_dimids);
506 if (mnc_iid < 0) {
507 fprintf(stderr, "Can't create the image variable\n");
508 return (-1);
509 }
510
511 miattputstr(mnc_fd, mnc_iid, MIsigntype,
512 (mnc_vsign) ? MI_SIGNED : MI_UNSIGNED);
513
514
515 /* Calculate the starts, steps, and direction cosines. This only
516 * be done properly if the file is NIfTI-1 file. If it is an Analyze
517 * file we have to resort to other methods...
518 */
519
520 if (nii_ptr->nifti_type != 0 &&
521 (nii_ptr->sform_code != NIFTI_XFORM_UNKNOWN ||
522 nii_ptr->qform_code != NIFTI_XFORM_UNKNOWN)) {
523
524 make_identity_transform(&mnc_xform);
525
526 if (nii_ptr->sform_code != NIFTI_XFORM_UNKNOWN) {
527 if (!qflag) {
528 printf("Using s-form transform:\n");
529 }
530 for (i = 0; i < 4; i++) {
531 for (j = 0; j < 4; j++) {
532 Transform_elem(mnc_xform, i, j) = nii_ptr->sto_xyz.m[i][j];
533 if (!qflag) {
534 printf("%8.4f, ", nii_ptr->sto_xyz.m[i][j]);
535 }
536 }
537 if (!qflag) {
538 printf("\n");
539 }
540 }
541 }
542 else {
543 if (!qflag) {
544 printf("Using q-form transform:\n");
545 }
546 for (i = 0; i < 4; i++) {
547 for (j = 0; j < 4; j++) {
548 Transform_elem(mnc_xform, i, j) = nii_ptr->qto_xyz.m[i][j];
549 if (!qflag) {
550 printf("%8.4f, ", nii_ptr->qto_xyz.m[i][j]);
551 }
552 }
553 if (!qflag) {
554 printf("\n");
555 }
556 }
557 }
558
559 create_linear_transform(&mnc_linear_xform, &mnc_xform);
560
561 convert_transform_to_starts_and_steps(&mnc_linear_xform,
562 MAX_SPACE_DIMS,
563 NULL,
564 mnc_spatial_axes,
565 mnc_starts,
566 mnc_steps,
567 mnc_dircos);
568
569 }
570 else {
571 /* No official transform was found (possibly this is an Analyze
572 * file). Just use some reasonable defaults.
573 */
574 mnc_steps[mnc_spatial_axes[DIM_X]] = nii_ptr->dx;
575 mnc_steps[mnc_spatial_axes[DIM_Y]] = nii_ptr->dy;
576 mnc_steps[mnc_spatial_axes[DIM_Z]] = nii_ptr->dz;
577 mnc_starts[mnc_spatial_axes[DIM_X]] = -(nii_ptr->dx * nii_ptr->nx) / 2;
578 mnc_starts[mnc_spatial_axes[DIM_Y]] = -(nii_ptr->dy * nii_ptr->ny) / 2;
579 mnc_starts[mnc_spatial_axes[DIM_Z]] = -(nii_ptr->dz * nii_ptr->nz) / 2;
580
581 /* Unlike the starts and steps, the direction cosines do NOT change
582 * based upon the data orientation.
583 */
584 for (i = 0; i < MAX_SPACE_DIMS; i++) {
585 for (j = 0; j < MAX_SPACE_DIMS; j++) {
586 mnc_dircos[i][j] = (i == j) ? 1.0 : 0.0;
587 }
588 }
589 }
590
591 switch (nii_ptr->xyz_units) {
592 case NIFTI_UNITS_METER:
593 for (i = 0; i < MAX_SPACE_DIMS; i++) {
594 mnc_starts[i] *= 1000;
595 mnc_steps[i] *= 1000;
596 }
597 break;
598 case NIFTI_UNITS_MM:
599 break;
600 case NIFTI_UNITS_MICRON:
601 for (i = 0; i < MAX_SPACE_DIMS; i++) {
602 mnc_starts[i] /= 1000;
603 mnc_steps[i] /= 1000;
604 }
605 break;
606 default:
607 fprintf(stderr, "Unknown XYZ units %d\n", nii_ptr->xyz_units);
608 break;
609 }
610
611 /* Now we write the spatial axis information to the file. The starts,
612 * steps, and cosines have to be associated with the correct spatial
613 * axes. Also, we perform any orientation flipping that was requested.
614 */
615 for (i = 0; i < MAX_SPACE_DIMS; i++) {
616 if (!qflag) {
617 printf("%s start: %8.4f step: %8.4f cosines: %8.4f %8.4f %8.4f\n",
618 mnc_spatial_names[i],
619 mnc_starts[i],
620 mnc_steps[i],
621 mnc_dircos[i][DIM_X],
622 mnc_dircos[i][DIM_Y],
623 mnc_dircos[i][DIM_Z]);
624 }
625 mnc_vid = ncvarid(mnc_fd, mnc_spatial_names[i]);
626
627 /* If we selected "flipping" of the appropriate axis, do it here
628 */
629 if (flip[i]) {
630 miattputdbl(mnc_fd, mnc_vid, MIstart,
631 mnc_starts[i]+((mnc_count[i]-1)*mnc_steps[i]));
632 miattputdbl(mnc_fd, mnc_vid, MIstep, -mnc_steps[i]);
633 }
634 else {
635 miattputdbl(mnc_fd, mnc_vid, MIstart, mnc_starts[i]);
636 miattputdbl(mnc_fd, mnc_vid, MIstep, mnc_steps[i]);
637 }
638 ncattput(mnc_fd, mnc_vid, MIdirection_cosines, NC_DOUBLE,
639 MAX_SPACE_DIMS, mnc_dircos[i]);
640 }
641
642 /* Find the valid minimum and maximum of the data, in order to set the
643 * global image minimum and image maximum properly.
644 */
645 if (rflag) {
646 find_data_range(nii_ptr->datatype,
647 nii_ptr->nvox,
648 nii_ptr->data,
649 mnc_vrange);
650 }
651
652 if (nii_ptr->scl_slope != 0.0) {
653 /* Convert slope/offset to min/max
654 */
655 mnc_srange[0] = (mnc_vrange[0] * nii_ptr->scl_slope) + nii_ptr->scl_inter;
656 mnc_srange[1] = (mnc_vrange[1] * nii_ptr->scl_slope) + nii_ptr->scl_inter;
657 }
658 else {
659 mnc_srange[0] = mnc_vrange[0];
660 mnc_srange[1] = mnc_vrange[1];
661 }
662
663 ncattput(mnc_fd, mnc_iid, MIvalid_range, NC_DOUBLE, 2, mnc_vrange);
664 miattputstr(mnc_fd, NC_GLOBAL, MIhistory, mnc_hist);
665
666 /* Switch out of definition mode.
667 */
668 ncendef(mnc_fd);
669
670 /* Finally, write the values of the image-min, image-max, and image
671 * variables.
672 */
673 mivarput1(mnc_fd, ncvarid(mnc_fd, MIimagemin), mnc_start, NC_DOUBLE,
674 MI_SIGNED, &mnc_srange[0]);
675
676 mivarput1(mnc_fd, ncvarid(mnc_fd, MIimagemax), mnc_start, NC_DOUBLE,
677 MI_SIGNED, &mnc_srange[1]);
678
679 mivarput(mnc_fd, mnc_iid, mnc_start, mnc_count, mnc_mtype,
680 (mnc_msign) ? MI_SIGNED : MI_UNSIGNED, nii_ptr->data);
681
682 miclose(mnc_fd);
683
684 return (0);
685 }
686