1 /* ----------------------------- MNI Header -----------------------------------
2 @NAME : mincmakescalar
3 @INPUT : argc, argv - command line arguments
4 @OUTPUT : (none)
5 @RETURNS : status
6 @DESCRIPTION: Program to make a scalar minc file from a vector minc file.
7 @METHOD :
8 @GLOBALS :
9 @CALLS :
10 @CREATED : August 7, 1997 (Peter Neelin)
11 @MODIFIED :
12 * $Log: mincmakescalar.c,v $
13 * Revision 6.10 2008-01-17 02:33:02 rotor
14 * * removed all rcsids
15 * * removed a bunch of ^L's that somehow crept in
16 * * removed old (and outdated) BUGS file
17 *
18 * Revision 6.9 2008/01/12 19:08:15 stever
19 * Add __attribute__ ((unused)) to all rcsid variables.
20 *
21 * Revision 6.8 2007/12/11 12:43:01 rotor
22 * * added static to all global variables in main programs to avoid linking
23 * problems with libraries (compress in mincconvert and libz for example)
24 *
25 * Revision 6.7 2005/08/26 21:07:17 bert
26 * Use #if rather than #ifdef with MINC2 symbol, and be sure to include config.h whereever MINC2 is used
27 *
28 * Revision 6.6 2005/01/28 19:34:09 bert
29 * Warn user if vector_dimension is not the last dimension in the image
30 *
31 * Revision 6.5 2004/11/01 22:38:38 bert
32 * Eliminate all references to minc_def.h
33 *
34 * Revision 6.4 2004/04/27 15:32:33 bert
35 * Added -2 option
36 *
37 * Revision 6.3 2001/04/24 13:38:43 neelin
38 * Replaced NC_NAT with MI_ORIGINAL_TYPE.
39 *
40 * Revision 6.2 2001/04/17 18:40:20 neelin
41 * Modifications to work with NetCDF 3.x
42 * In particular, changed NC_LONG to NC_INT (and corresponding longs to ints).
43 * Changed NC_UNSPECIFIED to NC_NAT.
44 * A few fixes to the configure script.
45 *
46 * Revision 6.1 1999/10/19 14:45:25 neelin
47 * Fixed Log subsitutions for CVS
48 *
49 * Revision 6.0 1997/09/12 13:24:20 neelin
50 * Release of minc version 0.6
51 *
52 * Revision 5.0 1997/08/21 13:25:19 neelin
53 * Release of minc version 0.5
54 *
55 * Revision 1.2 1997/08/07 16:23:19 neelin
56 * Change -length to -magnitude.
57 *
58 * Revision 1.1 1997/08/07 16:09:07 neelin
59 * Initial revision
60 *
61 ---------------------------------------------------------------------------- */
62
63 #if HAVE_CONFIG_H
64 #include "config.h"
65 #endif
66
67 #include <stdlib.h>
68 #include <stdio.h>
69 #include <string.h>
70 #include <ctype.h>
71 #include <math.h>
72 #include <float.h>
73 #include <minc.h>
74 #include <ParseArgv.h>
75 #include <time_stamp.h>
76 #include <voxel_loop.h>
77
78 #ifndef TRUE
79 # define TRUE 1
80 # define FALSE 0
81 #endif
82
83 #define DEFAULT_RANGE DBL_MAX
84 #define NCOPTS_DEFAULT NC_VERBOSE | NC_FATAL
85 #define INVALID_DATA (-DBL_MAX)
86
87 /* Types */
88 typedef enum
89 {CONV_DEFAULT, CONV_AVERAGE, CONV_MAGNITUDE, CONV_GREY, CONV_LINEAR}
90 Conversion_Type;
91
92 /* Double_Array structure */
93 typedef struct {
94 int numvalues;
95 double *values;
96 } Double_Array;
97
98 /* Structure for program info */
99 typedef struct {
100 Conversion_Type conversion_type;
101 int num_coefficients;
102 double *linear_coefficients;
103 } Program_Data;
104
105 /* Function prototypes */
106 static void do_makescalar(void *caller_data, long num_voxels,
107 int input_num_buffers, int input_vector_length,
108 double *input_data[],
109 int output_num_buffers, int output_vector_length,
110 double *output_data[], Loop_Info *loop_info);
111 static int get_double_list(char *dst, char *key, char *nextarg);
112 static long get_vector_length(int mincid);
113
114 /* Argument variables */
115 static int clobber = FALSE;
116 #if MINC2
117 static int v2format = FALSE;
118 #endif /* MINC2 */
119 static int verbose = TRUE;
120 static nc_type datatype = MI_ORIGINAL_TYPE;
121 static int is_signed = FALSE;
122 static double valid_range[2] = {0.0, 0.0};
123 static int buffer_size = 10 * 1024;
124 static Conversion_Type conversion_type = CONV_DEFAULT;
125 static Double_Array linear_coefficients = {0, NULL};
126
127 /* Argument table */
128 static ArgvInfo argTable[] = {
129 #if MINC2
130 {"-2", ARGV_CONSTANT, (char *) TRUE, (char *) &v2format,
131 "Produce a MINC 2.0 format output file."},
132 #endif /* MINC2 */
133 {"-clobber", ARGV_CONSTANT, (char *) TRUE, (char *) &clobber,
134 "Overwrite existing file."},
135 {"-noclobber", ARGV_CONSTANT, (char *) FALSE, (char *) &clobber,
136 "Don't overwrite existing file (default)."},
137 {"-verbose", ARGV_CONSTANT, (char *) TRUE, (char *) &verbose,
138 "Print out log messages (default)."},
139 {"-quiet", ARGV_CONSTANT, (char *) FALSE, (char *) &verbose,
140 "Do not print out log messages."},
141 {"-buffer_size", ARGV_INT, (char *) 1, (char *) &buffer_size,
142 "Set the internal buffer size (in kb)."},
143 {"-filetype", ARGV_CONSTANT, (char *) MI_ORIGINAL_TYPE, (char *) &datatype,
144 "Use data type of first file (default)."},
145 {"-byte", ARGV_CONSTANT, (char *) NC_BYTE, (char *) &datatype,
146 "Write out byte data."},
147 {"-short", ARGV_CONSTANT, (char *) NC_SHORT, (char *) &datatype,
148 "Write out short integer data."},
149 {"-int", ARGV_CONSTANT, (char *) NC_INT, (char *) &datatype,
150 "Write out 32-bit integer data."},
151 {"-long", ARGV_CONSTANT, (char *) NC_INT, (char *) &datatype,
152 "Superseded by -int."},
153 {"-float", ARGV_CONSTANT, (char *) NC_FLOAT, (char *) &datatype,
154 "Write out single-precision floating-point data."},
155 {"-double", ARGV_CONSTANT, (char *) NC_DOUBLE, (char *) &datatype,
156 "Write out double-precision floating-point data."},
157 {"-signed", ARGV_CONSTANT, (char *) TRUE, (char *) &is_signed,
158 "Write signed integer data."},
159 {"-unsigned", ARGV_CONSTANT, (char *) FALSE, (char *) &is_signed,
160 "Write unsigned integer data (default if type specified)."},
161 {"-valid_range", ARGV_FLOAT, (char *) 2, (char *) valid_range,
162 "Valid range for output data."},
163 {"-magnitude", ARGV_CONSTANT, (char *) CONV_MAGNITUDE,
164 (char *) &conversion_type,
165 "Compute magnitude of vectors (default)."},
166 {"-average", ARGV_CONSTANT, (char *) CONV_AVERAGE,
167 (char *) &conversion_type,
168 "Average components of vectors."},
169 {"-rgbtogrey", ARGV_CONSTANT, (char *) CONV_GREY, (char *) &conversion_type,
170 "Convert RGB to greyscale."},
171 {"-rgbtogray", ARGV_CONSTANT, (char *) CONV_GREY, (char *) &conversion_type,
172 "Synonym for rgbtogrey."},
173 {"-linear", ARGV_FUNC, (char *) get_double_list,
174 (char *) &linear_coefficients,
175 "Specify comma-separated list of coefficients for linear combination."},
176 {NULL, ARGV_END, NULL, NULL, NULL}
177 };
178
179 static const char str_wrong_dimension_order[] = {
180 "Your input file contains an image with a vector dimension, but the\n"
181 "vector dimension isn't the last (i.e. fastest-varying) dimension in\n"
182 "the image. Please restructure the file using mincreshape before\n"
183 "attempting to use mincmakescalar.\n"
184 };
185
186 /* Main program */
187
main(int argc,char * argv[])188 int main(int argc, char *argv[])
189 {
190 char *infile, *outfile;
191 char *arg_string;
192 int inmincid;
193 int input_vector_length;
194 int ivalue;
195 Loop_Options *loop_options;
196 Program_Data program_data;
197
198 /* Save time stamp and args */
199 arg_string = time_stamp(argc, argv);
200
201 /* Get arguments */
202 if (ParseArgv(&argc, argv, argTable, 0) || (argc != 3)) {
203 (void) fprintf(stderr, "\nUsage: %s [options] <in.mnc> <out.mnc>\n",
204 argv[0]);
205 (void) fprintf(stderr, " %s -help\n\n",
206 argv[0]);
207 exit(EXIT_FAILURE);
208 }
209 infile = argv[1];
210 outfile = argv[2];
211
212 /* Check for conflicting options */
213 if ((conversion_type != CONV_DEFAULT) &&
214 (linear_coefficients.numvalues > 0)) {
215 (void) fprintf(stderr,
216 "Do not specify -linear with other conversion options\n");
217 exit(EXIT_FAILURE);
218 }
219
220 /* Set up conversion information */
221 if (conversion_type == CONV_DEFAULT) conversion_type = CONV_MAGNITUDE;
222 program_data.conversion_type = conversion_type;
223 program_data.num_coefficients = 0;
224 program_data.linear_coefficients = NULL;
225
226 /* Check for coefficients for linear combination */
227 if (linear_coefficients.numvalues > 0) {
228 conversion_type = CONV_LINEAR;
229 program_data.conversion_type = conversion_type;
230 program_data.num_coefficients = linear_coefficients.numvalues;
231 program_data.linear_coefficients =
232 malloc(linear_coefficients.numvalues *
233 sizeof(*program_data.linear_coefficients));
234 for (ivalue=0; ivalue < linear_coefficients.numvalues; ivalue++) {
235 program_data.linear_coefficients[ivalue] =
236 linear_coefficients.values[ivalue];
237 }
238 }
239
240 /* Open the input file and get the vector length */
241 inmincid = miopen(infile, NC_NOWRITE);
242 input_vector_length = get_vector_length(inmincid);
243 if (input_vector_length < 0) {
244 fprintf(stderr, str_wrong_dimension_order);
245 exit(EXIT_FAILURE);
246 }
247 if (input_vector_length < 1) input_vector_length = 1;
248
249 /* Check that this length is okay */
250 if ((conversion_type == CONV_GREY) &&
251 (input_vector_length != 3) && (input_vector_length > 1)) {
252 (void) fprintf(stderr, "Input file does not contain RGB data\n");
253 exit(EXIT_FAILURE);
254 }
255 if ((conversion_type == CONV_LINEAR) &&
256 (input_vector_length != program_data.num_coefficients) &&
257 (input_vector_length > 1)) {
258 (void) fprintf(stderr,
259 "Input vector length does not match number of linear coefficients\n");
260 exit(EXIT_FAILURE);
261 }
262
263 /* Set up looping options */
264 loop_options = create_loop_options();
265 set_loop_clobber(loop_options, clobber);
266 set_loop_verbose(loop_options, verbose);
267 #if MINC2
268 set_loop_v2format(loop_options, v2format);
269 #endif /* MINC2 */
270 set_loop_datatype(loop_options, datatype, is_signed,
271 valid_range[0], valid_range[1]);
272 set_loop_output_vector_size(loop_options, 1);
273 set_loop_buffer_size(loop_options, (long) buffer_size * 1024);
274 set_loop_first_input_mincid(loop_options, inmincid);
275
276 /* Do loop */
277 voxel_loop(1, &infile, 1, &outfile, arg_string, loop_options,
278 do_makescalar, (void *) &program_data);
279
280 /* Free stuff */
281 if (program_data.linear_coefficients != NULL) {
282 free(program_data.linear_coefficients);
283 }
284 if (linear_coefficients.values != NULL) {
285 free(linear_coefficients.values);
286 }
287
288 exit(EXIT_SUCCESS);
289 }
290
291 /* ----------------------------- MNI Header -----------------------------------
292 @NAME : do_makescalar
293 @INPUT : caller_data - pointer to structure containing lookup info
294 num_voxels - number of voxels to work on
295 input_num_buffers - number of input buffers
296 input_vector_length - length of input vector dimension
297 input_data - vector of pointers to input buffer data
298 output_num_buffers - number of output buffers
299 output_vector_length - length of output vector dimension
300 start - vector specifying start of hyperslab (not used)
301 count - vector specifying count of hyperslab (not used)
302 @OUTPUT : output_data - vector of pointers to output buffer data
303 @RETURNS : (nothing)
304 @DESCRIPTION: Routine to loop through an array of vector values and convert
305 them to scalar.
306 @METHOD :
307 @GLOBALS :
308 @CALLS :
309 @CREATED : August 7, 1997 (Peter Neelin)
310 @MODIFIED :
311 ---------------------------------------------------------------------------- */
do_makescalar(void * caller_data,long num_voxels,int input_num_buffers,int input_vector_length,double * input_data[],int output_num_buffers,int output_vector_length,double * output_data[],Loop_Info * loop_info)312 static void do_makescalar(void *caller_data, long num_voxels,
313 int input_num_buffers, int input_vector_length,
314 double *input_data[],
315 int output_num_buffers, int output_vector_length,
316 double *output_data[], Loop_Info *loop_info)
317 /* ARGSUSED */
318 {
319 Program_Data *program_data;
320 long ivoxel;
321 double *input_vector;
322 double value, result;
323 int ivalue;
324 static double grey_coefficients[3] = {0.299, 0.587, 0.114};
325
326 /* Get pointer to lookup info */
327 program_data = (Program_Data *) caller_data;
328
329 /* Check that values correspond */
330 if ((input_num_buffers != 1) || (output_num_buffers != 1) ||
331 (output_vector_length != 1)) {
332 (void) fprintf(stderr, "Bad internal values.\n");
333 exit(EXIT_FAILURE);
334 }
335 if ((program_data->conversion_type == CONV_GREY) &&
336 (input_vector_length != 3)) {
337 (void) fprintf(stderr,
338 "Input must have three components for -grey conversion.\n");
339 exit(EXIT_FAILURE);
340 }
341
342 /* Loop through the voxels */
343 for (ivoxel=0; ivoxel < num_voxels; ivoxel++) {
344
345 /* Handle the special case of scalar input */
346 if (input_vector_length <= 1) {
347 output_data[0][ivoxel] = input_data[0][ivoxel];
348 continue;
349 }
350
351 /* Get location of input value */
352 input_vector = &input_data[0][ivoxel*input_vector_length];
353
354 /* Loop over components of vector */
355 result = 0.0;
356 for (ivalue=0; ivalue < input_vector_length; ivalue++) {
357
358 /* Get value and check for invalid values */
359 value = input_vector[ivalue];
360 if (value == INVALID_DATA) {
361 result = INVALID_DATA;
362 break;
363 }
364
365 /* Sum things up according to scheme */
366 switch (program_data->conversion_type) {
367 case CONV_AVERAGE:
368 result += value;
369 break;
370 case CONV_MAGNITUDE:
371 result += value * value;
372 break;
373 case CONV_GREY:
374 result += grey_coefficients[ivalue] * value;
375 break;
376 case CONV_LINEAR:
377 result += program_data->linear_coefficients[ivalue] * value;
378 break;
379 }
380
381 } /* End of loop over components */
382
383 /* Finish up the calculation */
384 if (result != INVALID_DATA) {
385 switch (program_data->conversion_type) {
386 case CONV_AVERAGE:
387 result /= (double) input_vector_length;
388 break;
389 case CONV_MAGNITUDE:
390 result = sqrt(result);
391 break;
392 }
393 }
394
395 /* Save the result */
396 output_data[0][ivoxel] = result;
397 }
398
399 return;
400 }
401
402 /* ----------------------------- MNI Header -----------------------------------
403 @NAME : get_double_list
404 @INPUT : dst - client data passed by ParseArgv
405 key - matching key in argv
406 nextarg - argument following key in argv
407 @OUTPUT : (none)
408 @RETURNS : TRUE since nextarg is used.
409 @DESCRIPTION: Gets a list (array) of double values.
410 @METHOD :
411 @GLOBALS :
412 @CALLS :
413 @CREATED : March 8, 1995 (Peter Neelin)
414 @MODIFIED :
415 ---------------------------------------------------------------------------- */
get_double_list(char * dst,char * key,char * nextarg)416 static int get_double_list(char *dst, char *key, char *nextarg)
417 {
418 #define VECTOR_SEPARATOR ','
419
420 int num_elements;
421 int num_alloc;
422 double *double_list;
423 double dvalue;
424 char *cur, *end, *prev;
425 Double_Array *double_array;
426
427 /* Check for a following argument */
428 if (nextarg == NULL) {
429 (void) fprintf(stderr,
430 "\"%s\" option requires an additional argument\n",
431 key);
432 exit(EXIT_FAILURE);
433 }
434
435 /* Get pointers to array variables */
436 double_array = (Double_Array *) dst;
437
438 /* Set up pointers to end of string and first non-space character */
439 end = nextarg + strlen(nextarg);
440 cur = nextarg;
441 while (isspace(*cur)) cur++;
442 num_elements = 0;
443 num_alloc = 0;
444 double_list = NULL;
445
446 /* Loop through string looking for doubles */
447 while (cur!=end) {
448
449 /* Get double */
450 prev = cur;
451 dvalue = strtod(prev, &cur);
452 if (cur == prev) {
453 (void) fprintf(stderr,
454 "expected vector of doubles for \"%s\", but got \"%s\"\n",
455 key, nextarg);
456 exit(EXIT_FAILURE);
457 }
458
459 /* Add the value to the list */
460 num_elements++;
461 if (num_elements > num_alloc) {
462 num_alloc += 20;
463 if (double_list == NULL) {
464 double_list =
465 malloc(num_alloc * sizeof(*double_list));
466 }
467 else {
468 double_list =
469 realloc(double_list, num_alloc * sizeof(*double_list));
470 }
471 }
472 double_list[num_elements-1] = dvalue;
473
474 /* Skip any spaces */
475 while (isspace(*cur)) cur++;
476
477 /* Skip an optional comma */
478 if (*cur == VECTOR_SEPARATOR) cur++;
479
480 }
481
482 /* Update the global variables */
483 double_array->numvalues = num_elements;
484 if (double_array->values != NULL) {
485 free(double_array->values);
486 }
487 double_array->values = double_list;
488
489 return TRUE;
490 }
491
492 /* ----------------------------- MNI Header -----------------------------------
493 @NAME : get_vector_length
494 @INPUT : mincid - minc file id
495 @OUTPUT : (none)
496 @RETURNS : Length of vector dimension or zero if no such dimension.
497 @DESCRIPTION: Routine to get the length of the vector dimension in a minc file.
498 @METHOD :
499 @GLOBALS :
500 @CALLS :
501 @CREATED : November 30, 1994 (Peter Neelin)
502 @MODIFIED :
503 ---------------------------------------------------------------------------- */
get_vector_length(int mincid)504 static long get_vector_length(int mincid)
505 {
506 int imgid;
507 int ndims;
508 int dim[MAX_VAR_DIMS];
509 char dimname[MAX_NC_NAME];
510 long vector_length;
511 int i;
512
513 /* Get image variable id */
514 imgid = ncvarid(mincid, MIimage);
515
516 /* Get the image dimension info */
517 (void) ncvarinq(mincid, imgid, NULL, NULL, &ndims, dim, NULL);
518
519 /* Check for vector dimension */
520 (void) ncdiminq(mincid, dim[ndims-1], dimname, &vector_length);
521 if ((strcmp(dimname, MIvector_dimension) != 0) || (ndims <= 2)) {
522 vector_length = 0;
523
524 /* New deal - check for an actual vector_dimension anywhere in the
525 * file.
526 */
527 for (i = 0; i < ndims; i++) {
528 ncdiminq(mincid, dim[i], dimname, NULL);
529 if (!strcmp(dimname, MIvector_dimension)) {
530 vector_length = -1;
531 break;
532 }
533 }
534 }
535
536 return vector_length;
537 }
538
539