1 #define MAIN
2 
3 #include "SUMA_suma.h"
4 #include "thd_segtools_fNM.h"
5 #include "SUMA_SegOpts.h"
6 #include "SUMA_SegFunc.h"
7 #include "matrix.h"
8 
9 
10 static int vn=0 ;
11 
12 
13 
14 static char shelp_Hist[] = {
15 "3dHist computes histograms using functions for generating priors.\n"
16 "If you are not sure you need this particular program, use 3dhistog instead.\n"
17 "\n"
18 "Example:\n"
19 "3dHist      -input sigs+orig \\n"
20 "               \n"
21 "\n"
22 "Options:\n"
23 "   -input DSET: Dset providing values for histogram. Exact 0s are not counted\n"
24 "   -dind SB: Use sub-brick SB from the input rather than 0\n"
25 "   -mask MSET: Provide mask dataset to select subset of input.\n"
26 "   -mask_range BOT TOP: Specify the range of values to consider from MSET.\n"
27 "                        Default is anything non-zero\n"
28 "   -cmask CMASK: Provide cmask expression. Voxels where expression is 0\n"
29 "                 are excluded from computations. For example:\n"
30 "            -cmask '-a T1.div.r+orig -b T1.uni.r+orig -expr step(a/b-10)'\n"
31 "   -thishist HIST.niml.hist: Read this previously created histogram instead\n"
32 "                             of forming one from DSET.\n"
33 "                             Obviously, DSET, or -mask options are not needed\n"
34 "   -prefix PREF: Write histogram to niml file called PREF.niml.hist \n"
35 "   -equalized PREF: Write a histogram equalized version of the input dataset\n"
36 "   Histogram Creation Parameters:\n"
37 "     By default, the program will select bin number, bin width, \n"
38 "     and range automatically. You can also set the parameters manually with \n"
39 "     the following options.\n"
40 "   -nbin K: Use K bins.\n"
41 "   -min MIN: Minimum intensity.\n"
42 "   -max MAX: Maximum intensity.\n"
43 "   -binwidth BW: Bin width\n"
44 "   -ignore_out: Do not count samples outside the user specified range.\n"
45 "   -rhist RHIST.niml.hist: Use previously created histogram to set range\n"
46 "                           and binwidth parameters.\n"
47 "\n"
48 "   -showhist: Display histogram to stdout\n"
49 "              You can also graph it with: 1dRplot HistOut.niml.hist\n"
50 "\n"
51 "   Histogram Queries:\n"
52 "   -at VAL: Set the value at which you want histogram values\n"
53 "   -get 'PAR1,PAR2,PAR3..': Return the following PAR* properties at VAL\n"
54 "                            Choose from:\n"
55 "                            freq: Frequency (normalized count)\n"
56 "                            count: Count\n"
57 "                            bin: Continuous bin location estimate\n"
58 "                            cdf: Cumulative count\n"
59 "                            rcdf: Reverse cumulative count (from the top)\n"
60 "                            ncdf: The normalized version of cdf\n"
61 "                            nrcdf: The reverse version of ncdf\n"
62 "                            outl: 1.0-(2*smallest tail area)\n"
63 "                               0 means VAL splits area in the middle\n"
64 "                               1 means VAL is at either end of the histogram\n"
65 "                            ALL: All the above.\n"
66 "                   You can select multiple ones with something like:\n"
67 "                        -get 'freq, count, bin' \n"
68 "\n"
69 "                   You can also set one of the PAR* to 'upvol' to get \n"
70 "                   the volume (liters) of voxels with values exceeding VAL\n"
71 "                   The use of upvol usually requires option -voxvol too.\n"
72 "  -voxvol VOL_MM3: A voxel's volume in mm^3. To be used with upvol if\n"
73 "                   no dataset is available or if you want to override\n"
74 "                   it.\n"
75 "  -val_at PAR PARVAL: Return the value (magnitude) where histogram property\n"
76 "                      PAR is equal to PARVAL\n"
77 "                      PAR can only be one of: cdf, rcdf, ncdf, nrcdf, upvol\n"
78 "                      For upvol, PARVAL is in Liters\n"
79 "  -quiet: Return a concise output to simplify parsing. For the moment, this\n"
80 "          option only affects output of option -val_at\n"
81 "\n"
82 "  Examples:\n"
83 "       #A histogram a la 3dhistog:\n"
84 "       3dHist -input T1+orig.\n"
85 "       #Getting parameters from previously created histogram:\n"
86 "       3dHist -thishist HistOut.niml.hist -at 144.142700 \n"
87 "       #Or the reverse query:\n"
88 "       3dHist -thishist HistOut.niml.hist -val_at ncdf 0.132564\n"
89 "       #Compute histogram and find dataset threshold (approximate)\n"
90 "       #such that 1.5 liters of voxels remain above it.\n"
91 "       3dHist -prefix toy -input flair_axial.nii.gz -val_at upvol 1.5 \n"
92 "\n"
93 /* Untested here
94 "   -mrange M0 M1: Consider MASK only for values between M0 and M1, inclusive\n"
95 "   -debug DBG: Set debug level\n"
96 */
97 "\n"
98 };
99 
100 
Hist_usage(int detail)101 void Hist_usage(int detail)
102 {
103    int i = 0;
104 
105    ENTRY("Hist_usage");
106 
107 
108    printf( "%s", shelp_Hist );
109    PRINT_COMPILE_DATE ;
110    exit(0);
111 }
112 
Hist_Default(char * argv[],int argc)113 SEG_OPTS *Hist_Default(char *argv[], int argc)
114 {
115    SEG_OPTS *Opt=NULL;
116 
117    ENTRY("Hist_Default");
118 
119    Opt = SegOpt_Struct();
120 
121    Opt->helpfunc = &Hist_usage;
122    Opt->ps = SUMA_Parse_IO_Args(argc, argv, "-talk;");
123    Opt->mset_name = NULL;
124    Opt->sig_name = NULL;
125    Opt->ndist_name = NULL;
126    Opt->uid[0] = '\0';
127    Opt->prefix = NULL;
128    Opt->crefix = NULL;
129    Opt->aset = NULL;
130    Opt->mset = NULL;
131    Opt->gset = NULL;
132    Opt->sig = NULL;
133    Opt->FDV = NULL;
134    Opt->pset = NULL;
135    Opt->cset = NULL;
136    Opt->debug = 0;
137    Opt->verbose = 1;
138    Opt->idbg = Opt->kdbg = Opt->jdbg = -1;
139    Opt->feats=Opt->clss=NULL;
140    Opt->feat_exp=NULL; Opt->featexpmeth=0; Opt->featsfam=NULL;
141    Opt->keys = NULL;
142    Opt->mixfrac=NULL;
143    Opt->logp = 1;
144    Opt->VoxDbg = -1;
145    Opt->VoxDbgOut = NULL;
146    Opt->rescale_p = 1;
147    Opt->openmp = 0;
148    Opt->labeltable_name = NULL;
149    Opt->pweight = 1;
150    Opt->cmask = NULL;
151    Opt->dimcmask = 0;
152    Opt->cmask_count=0;
153    Opt->mask_bot = 1.0;
154    Opt->mask_top = -1.0;
155    Opt->DO_p = TRUE;
156    Opt->DO_c = TRUE;
157    Opt->Writepcg_G_au = FALSE;
158    Opt->DO_r = FALSE;
159    Opt->group_classes = NULL;
160    Opt->group_keys = NULL;
161    Opt->fitmeth = SEG_LSQFIT;
162    Opt->proot = NULL;
163    Opt->cs = NULL;
164    Opt->Gcs = NULL;
165    Opt->fast = 1;
166    Opt->ShowThisDist = NULL;
167    Opt->wA = -123456999.0;
168    Opt->wL = -123456999.0;
169    Opt->this_xset_name=NULL;
170    Opt->smode = 0;
171    Opt->N_biasgroups = 0;
172    Opt->binwidth = 0;
173    Opt->UseTmp = 0;
174    Opt->range[0] = -9876543.0;
175    Opt->range[1] = -9876543.0;
176    Opt->bias_meth = NULL;
177    SUMA_RETURN(Opt);
178 }
179 
Hist_CheckOpts(SEG_OPTS * Opt)180 int Hist_CheckOpts(SEG_OPTS *Opt)
181 {
182    static char FuncName[]={"Hist_CheckOpts"};
183 
184    SUMA_ENTRY;
185 
186    if (  !Opt->sig_name   && !Opt->proot) {
187       SUMA_S_Err("Missing input");
188       SUMA_RETURN(0);
189    }
190 
191    SUMA_RETURN(1);
192 }
193 char ALL_GETS[]={"bin, count, freq, cdf, ncdf, rcdf, nrcdf, outl"};
194 
Hist_ParseInput(SEG_OPTS * Opt,char * argv[],int argc)195 SEG_OPTS *Hist_ParseInput (SEG_OPTS *Opt, char *argv[], int argc)
196 {
197    static char FuncName[]={"Hist_ParseInput"};
198    int kar, i, ind, exists;
199    char *outname, cview[10];
200    int brk = 0;
201    SUMA_GENERIC_ARGV_PARSE *ps=NULL;
202 
203    SUMA_ENTRY;
204 
205    brk = 0;
206    Opt->f1 = -1.0;
207    kar = 1;
208 	while (kar < argc) { /* loop accross command ine options */
209 		/*fprintf(stdout, "%s verbose: Parsing command line...\n", FuncName);*/
210 		if (strcmp(argv[kar], "-h") == 0 || strcmp(argv[kar], "-help") == 0) {
211 			 Opt->helpfunc(0);
212           exit (0);
213 		}
214 
215  		SUMA_SKIP_COMMON_OPTIONS(brk, kar);
216 
217       #ifdef USE_TRACING
218             if( strncmp(argv[kar],"-trace",5) == 0 ){
219                DBG_trace = 1 ;
220                brk = 1 ;
221             }
222             if( strncmp(argv[kar],"-TRACE",5) == 0 ){
223                DBG_trace = 2 ;
224                brk = 1 ;
225             }
226       #endif
227 
228       if (!brk && (strcmp(argv[kar], "-debug") == 0)) {
229          kar ++;
230 			if (kar >= argc)  {
231 		  		fprintf (stderr, "need argument after -debug \n");
232 				exit (1);
233 			}
234 			Opt->debug = atoi(argv[kar]);
235          brk = 1;
236 		}
237 
238       if (!brk && (strcmp(argv[kar], "-quiet") == 0)) {
239 			Opt->verbose = 0;
240          brk = 1;
241 		}
242 
243 
244       if (!brk && (strcmp(argv[kar], "-cmask") == 0)) {
245          kar ++;
246 			if (kar >= argc)  {
247 		  		ERROR_exit("-cmask option requires a following argument!\n");
248 			}
249 			Opt->cmask = EDT_calcmask( argv[kar] , &(Opt->dimcmask), 0 ) ;
250          if( Opt->cmask == NULL ) ERROR_exit("Can't compute -cmask!\n");
251          brk = 1;
252 		}
253 
254       if (!brk && (strcmp(argv[kar], "-mask") == 0)) {
255          kar ++;
256 			if (kar >= argc)  {
257 		  		fprintf (stderr, "need argument after -mask \n");
258 				exit (1);
259 			}
260 			Opt->mset_name = argv[kar];
261          brk = 1;
262       }
263 
264       if( !brk && (strncmp(argv[kar],"-mrange",5) == 0) ){
265          if( kar+2 >= argc )
266            ERROR_exit("-mrange option requires 2 following arguments!\n");
267          Opt->mask_bot = strtod( argv[++kar] , NULL ) ;
268          Opt->mask_top = strtod( argv[++kar] , NULL ) ;
269          if( Opt->mask_top < Opt->mask_bot )
270            ERROR_exit("-mrange inputs are illegal!\n") ;
271          brk = 1;
272       }
273 
274       if (!brk && (strcmp(argv[kar], "-input") == 0)) {
275          kar ++;
276 			if (kar >= argc)  {
277 		  		fprintf (stderr, "need argument after -input \n");
278 				exit (1);
279 			}
280 			while (kar < argc && argv[kar][0] != '-') {
281             Opt->sig_name =
282                SUMA_append_replace_string(Opt->sig_name, argv[kar], " ", 1);
283             ++kar;
284          }
285          if (kar < argc && argv[kar][0] == '-') --kar; /* unwind */
286          brk = 1;
287 		}
288 
289       if (!brk && (strcmp(argv[kar], "-rhist") == 0)) {
290          kar ++;
291 			if (kar >= argc)  {
292 		  		fprintf (stderr, "need argument after -rhist \n");
293 				exit (1);
294 			}
295 			Opt->ndist_name = argv[kar];
296          brk = 1;
297 		}
298 
299       if (!brk && (strcmp(argv[kar], "-nbin") == 0)) {
300          kar ++;
301 			if (kar >= argc)  {
302 		  		fprintf (stderr, "need argument after -nbin \n");
303 				exit (1);
304 			}
305 			Opt->UseTmp = atoi(argv[kar]);
306          brk = 1;
307 		}
308 
309       if (!brk && (strcmp(argv[kar], "-binwidth") == 0)) {
310          kar ++;
311 			if (kar >= argc)  {
312 		  		fprintf (stderr, "need value after -binwidth \n");
313 				exit (1);
314 			}
315 			Opt->binwidth = atof(argv[kar]);
316          brk = 1;
317 		}
318 
319       if (!brk && (strcmp(argv[kar], "-min") == 0)) {
320          kar ++;
321 			if (kar >= argc)  {
322 		  		fprintf (stderr, "need value after -min \n");
323 				exit (1);
324 			}
325 			Opt->range[0] = atof(argv[kar]);
326          brk = 1;
327 		}
328 
329       if (!brk && (strcmp(argv[kar], "-voxvol") == 0)) {
330          kar ++;
331 			if (kar >= argc)  {
332 		  		fprintf (stderr, "need value after -voxvol in mm^3\n");
333 				exit (1);
334 			}
335 			Opt->f1 = atof(argv[kar]);
336          brk = 1;
337 		}
338 
339       if (!brk && (strcmp(argv[kar], "-max") == 0)) {
340          kar ++;
341 			if (kar >= argc)  {
342 		  		fprintf (stderr, "need value after -max \n");
343 				exit (1);
344 			}
345 			Opt->range[1] = atof(argv[kar]);
346          brk = 1;
347 		}
348 
349       if( !brk && (strncmp(argv[kar],"-range",5) == 0) ){
350          if( kar+2 >= argc )
351            ERROR_exit("-range option requires 2 following arguments!\n");
352          Opt->range[0] = strtod( argv[++kar] , NULL ) ;
353          Opt->range[1] = strtod( argv[++kar] , NULL ) ;
354          if( Opt->range[0] < Opt->range[1] )
355            ERROR_exit("-range inputs are illegal!\n") ;
356          brk = 1;
357       }
358 
359       if (!brk && (strcmp(argv[kar], "-at") == 0)) {
360          kar ++;
361 			if (kar >= argc)  {
362 		  		fprintf (stderr, "need value after -at \n");
363 				exit (1);
364 			}
365 			Opt->wA = atof(argv[kar]);
366          brk = 1;
367 		}
368 
369       if (!brk && (strcmp(argv[kar], "-dind") == 0)) {
370          kar ++;
371 			if (kar >= argc)  {
372 		  		fprintf (stderr, "need integer after -dind \n");
373 				exit (1);
374 			}
375 			Opt->N_biasgroups = atoi(argv[kar]);
376          brk = 1;
377 		}
378 
379       if (!brk && (strcmp(argv[kar], "-get") == 0)) {
380          kar ++;
381 			if (kar >= argc)  {
382 		  		fprintf (stderr, "need string after -get \n");
383 				exit (1);
384 			}
385          if (strstr(argv[kar],"ALL")) {
386             Opt->feats = NI_strict_decode_string_list(ALL_GETS, ";, ");
387          } else {
388             if (!(Opt->feats = NI_strict_decode_string_list(argv[kar] ,";, "))){
389                SUMA_S_Errv("Bad option %s after -get", argv[kar]);
390                exit(1);
391             }
392          }
393          brk = 1;
394       }
395 
396       if (!brk && (strcmp(argv[kar], "-val_at") == 0)) {
397          kar ++;
398 			if (kar+1 >= argc)  {
399 		  		fprintf (stderr, "need string and value after -val_at \n");
400 				exit (1);
401 			}
402          Opt->this_xset_name = argv[kar];
403          if (strcmp(argv[kar],"cdf") &&
404              strcmp(argv[kar],"rcdf") &&
405              strcmp(argv[kar],"ncdf") &&
406              strcmp(argv[kar],"nrcdf") &&
407              strcmp(argv[kar],"upvol")   ) {
408             SUMA_S_Errv("String %s not allowed after -val_at\n",
409                         argv[kar]);
410             exit(1);
411          }
412          kar ++;
413          Opt->wL = atof(argv[kar]);
414          brk = 1;
415       }
416 
417       if (!brk && (strcmp(argv[kar], "-showhist") == 0)) {
418 			Opt->DO_r = TRUE;
419          brk = 1;
420 		}
421       if (!brk && (strcmp(argv[kar], "-ignore_out") == 0)) {
422 			Opt->smode = 1;
423          brk = 1;
424 		}
425 
426       if (!brk && (strcmp(argv[kar], "-thishist") == 0)) {
427          kar ++;
428 			if (kar >= argc)  {
429 		  		fprintf (stderr, "need argument after -thishist \n");
430 				exit (1);
431 			}
432 			Opt->proot = argv[kar];
433          brk = 1;
434 		}
435 
436       if (!brk && (strcmp(argv[kar], "-prefix") == 0)) {
437          kar ++;
438 			if (kar >= argc)  {
439 		  		fprintf (stderr, "need argument after -prefix \n");
440 				exit (1);
441 			}
442          Opt->prefix = (char*)calloc(strlen(argv[kar])+20, sizeof(char));
443          strcpy(Opt->prefix, argv[kar]); Opt->prefix[strlen(argv[kar])]='\0';
444          brk = 1;
445 		}
446 
447       if (!brk && (strcmp(argv[kar], "-equalized") == 0)) {
448          kar ++;
449 			if (kar >= argc)  {
450 		  		fprintf (stderr, "need argument after -equalized \n");
451 				exit (1);
452 			}
453          Opt->crefix = (char*)calloc(strlen(argv[kar])+20, sizeof(char));
454          strcpy(Opt->crefix, argv[kar]); Opt->crefix[strlen(argv[kar])]='\0';
455          brk = 1;
456 		}
457 
458 
459       if (!brk) {
460 			fprintf (stderr,"Option %s not understood. \n"
461                          "Try -help for usage\n", argv[kar]);
462 			suggest_best_prog_option(argv[0], argv[kar]);
463          exit (1);
464 		} else {
465 			brk = 0;
466 			kar ++;
467 		}
468 
469    }
470 
471    if (!Opt->prefix) Opt->prefix = strdup("./HistOut");
472    if (Opt->uid[0]=='\0') UNIQ_idcode_fill(Opt->uid);
473    if (Opt->VoxDbg > -1 && !Opt->VoxDbgOut) {
474       char stmp[256];
475       sprintf(stmp,"%d.GP.dbg", Opt->VoxDbg);
476       Opt->VoxDbgOut = fopen(stmp,"w");
477    }
478 
479 
480    SUMA_RETURN(Opt);
481 }
482 
483 
main(int argc,char ** argv)484 int main(int argc, char **argv)
485 {
486    static char FuncName[]={"3dHist"};
487    SEG_OPTS *Opt=NULL;
488    int ii=0;
489    SUMA_HIST *hh=NULL, *rhh=NULL;
490    SUMA_Boolean LocalHead = NOPE;
491 
492    SUMA_STANDALONE_INIT;
493 	SUMA_mainENTRY;
494 
495    SUMAg_DOv = SUMA_Alloc_DisplayObject_Struct (SUMA_MAX_DISPLAYABLE_OBJECTS);
496    Opt = Hist_Default(argv,  argc);
497    Opt = Hist_ParseInput (Opt, argv,  argc);
498 
499    if (!Hist_CheckOpts(Opt)) {
500       ERROR_exit("Failed on option check");
501    }
502 
503 
504    if (Opt->proot) { /* get this hist */
505       hh = SUMA_read_hist(Opt->proot);
506    } else {
507       if (Opt->ndist_name) { /* reference hist? */
508          rhh = SUMA_read_hist(Opt->ndist_name);
509          Opt->bias_meth = "NONE";
510       }
511       /* load the input data */
512       if (Opt->sig_name && !(Opt->sig = Seg_load_dset( Opt->sig_name ))) {
513          SUMA_S_Errv("Failed to load input %s\n", Opt->sig_name);
514          exit(1);
515       }
516 
517       if (Opt->mset_name) {
518          if (!(Opt->mset = Seg_load_dset( Opt->mset_name ))) {
519             exit(1);
520          }
521       }
522 
523       if (Opt->sig) {
524          Opt->cmask = MaskSetup(Opt, Opt->sig, 1,
525                    &(Opt->mset), &(Opt->cmask), Opt->dimcmask,
526                    Opt->mask_bot, Opt->mask_top, &(Opt->cmask_count));
527       }
528 
529       if (   Opt->binwidth != 0.0 || Opt->UseTmp != 0 ||
530              (int)Opt->range[0] != -9876543 || (int)Opt->range[1] != -9876543 ) {
531          Opt->bias_meth = "NONE";
532          SUMA_S_Note("Manual mode");
533          /* manual mode */
534          if (rhh) {
535             SUMA_S_Note("Overriding some reference histogram parameters");
536             if (Opt->UseTmp > 0) rhh->K = Opt->UseTmp;
537             if (Opt->binwidth != 0.0) rhh->W = Opt->binwidth;
538             if ((int)Opt->range[0] != -9876543) rhh->min = Opt->range[0];
539             if ((int)Opt->range[1] != -9876543) rhh->max = Opt->range[1];
540          }  else {
541             rhh = (SUMA_HIST *)SUMA_calloc(1, sizeof(SUMA_HIST));
542             rhh->K = Opt->UseTmp;
543             rhh->W = Opt->binwidth;
544             if ((int)Opt->range[0] != -9876543) rhh->min = Opt->range[0];
545                else rhh->min = 0.0;
546             if ((int)Opt->range[1] != -9876543) rhh->max = Opt->range[1];
547                else rhh->max = 0.0;
548          }
549          if (Opt->range[0] > Opt->range[1]) {
550             SUMA_S_Err ("Bad manual range values,"
551                         "both min and max must be specified");
552             exit(1);
553          }
554       }
555       /* Create the hist */
556       if (!(hh = SUMA_dset_hist(Opt->sig, Opt->N_biasgroups,
557                            Opt->cmask, DSET_PREFIX(Opt->sig), rhh,
558                            Opt->smode, 0.0, Opt->bias_meth))) {
559          SUMA_S_Errv("Failed to create histogram from %s\n",
560                      DSET_PREFIX(Opt->sig));
561          exit(1);
562       }
563       if (Opt->prefix) {
564          if (!SUMA_write_hist(hh, Opt->prefix)) {
565             SUMA_S_Errv("Failed to write hist to %s\n",
566                         Opt->prefix);
567             exit(1);
568          }
569 
570       }
571    }
572 
573    if (!hh) {
574       SUMA_S_Err("No Hist, Can't Travel.\n");
575       exit(1);
576    }
577 
578    if (Opt->crefix) {
579       THD_3dim_dataset *dout=NULL;
580       SUMA_LH("Equalizing");
581       if (!(dout = SUMA_dset_hist_equalize(Opt->sig, 0, Opt->cmask, hh))) {
582          SUMA_S_Err("Failed to equalize");
583          exit(1);
584       }
585       EDIT_dset_items( dout , ADN_prefix,  Opt->crefix,  ADN_none ) ;
586       tross_Copy_History( Opt->sig , dout ) ;
587       tross_Make_History( "3dHist" , argc,argv , dout ) ;
588 
589       DSET_write(dout);
590       DSET_delete(dout); dout=NULL;
591    }
592 
593 
594    if (Opt->DO_r) {
595       SUMA_Show_hist(hh, 0, NULL);
596    }
597 
598 
599    if (Opt->wA != -123456999.0) {
600       if (!Opt->feats) {
601          Opt->feats = NI_strict_decode_string_list(ALL_GETS, ";, ");
602       }
603       for (ii=0; ii<Opt->feats->num; ++ii) {
604          if (!strcmp(Opt->feats->str[ii], "upvol")) {
605             if (Opt->f1 < 0) {
606                if (Opt->sig) {
607                   Opt->f1 = SUMA_ABS(DSET_DX(Opt->sig)*
608                                     DSET_DY(Opt->sig)*DSET_DZ(Opt->sig));
609                } else {
610                   SUMA_S_Warn( "Have no -voxvol in mm^3 and no dataset "
611                               "from which to compute it, assuming 1mm^3");
612                   Opt->f1 = 1.0; // [PT: Feb 22, 2021] fix '='
613                }
614 
615             }
616             if (Opt->f1 == 0.0f) Opt->f1 = 1.0;
617          }
618       }
619       fprintf(stdout,"At value %f:\n", Opt->wA);
620       for (ii=0; ii<Opt->feats->num; ++ii) {
621          if (!strcmp(Opt->feats->str[ii], "upvol")) {
622             float vv;
623             if (Opt->f1 < 0) {
624                SUMA_S_Err( "Bad news seeting voxvol");
625                exit(1);
626             }
627             vv = SUMA_hist_value(hh, Opt->wA, "rcdf") / 1.0e6 * Opt->f1;
628             fprintf(stdout,"  %8s: %f\n",
629                Opt->feats->str[ii], vv);
630          } else {
631             fprintf(stdout,"  %8s: %f\n",
632                Opt->feats->str[ii],
633                SUMA_hist_value(hh, Opt->wA, Opt->feats->str[ii]));
634          }
635       }
636    }
637 
638    if (Opt->wL != -123456999.0) {
639       if (!Opt->this_xset_name || !strcmp(Opt->this_xset_name, "rcdf")) {
640          if (Opt->verbose)
641             fprintf(stdout,"Val: %f at %s: %f\n",
642                      SUMA_val_at_count(hh, Opt->wL, 0, 1),
643                      Opt->this_xset_name, Opt->wL);
644          else
645             fprintf(stdout,"%f\n", SUMA_val_at_count(hh, Opt->wL, 0, 1));
646       } else if (!strcmp(Opt->this_xset_name, "cdf")) {
647          if (Opt->verbose)
648             fprintf(stdout,"Val: %f at %s: %f\n",
649                      SUMA_val_at_count(hh, Opt->wL, 0, 0),
650                      Opt->this_xset_name, Opt->wL);
651          else
652             fprintf(stdout,"%f\n", SUMA_val_at_count(hh, Opt->wL, 0, 0));
653       } else if (!strcmp(Opt->this_xset_name, "ncdf")) {
654          if (Opt->verbose)
655             fprintf(stdout,"Val: %f at %s: %f\n",
656                      SUMA_val_at_count(hh, Opt->wL, 1, 0),
657                      Opt->this_xset_name, Opt->wL);
658          else
659             fprintf(stdout,"%f\n", SUMA_val_at_count(hh, Opt->wL, 1, 0));
660       } else if (!strcmp(Opt->this_xset_name, "nrcdf")) {
661          if (Opt->verbose)
662             fprintf(stdout,"Val: %f at %s: %f\n",
663                      SUMA_val_at_count(hh, Opt->wL, 1, 1),
664                      Opt->this_xset_name, Opt->wL);
665          else
666             fprintf(stdout,"%f\n", SUMA_val_at_count(hh, Opt->wL, 1, 1));
667       } else if (!strcmp(Opt->this_xset_name, "upvol")) {
668          float voxthr;
669          if (Opt->f1 < 0) {
670             if (Opt->sig) {
671                Opt->f1 = SUMA_ABS(DSET_DX(Opt->sig)*
672                                  DSET_DY(Opt->sig)*DSET_DZ(Opt->sig));
673             } else {
674                SUMA_S_Err( "Have no -voxvol in mm^3 and no dataset "
675                            "from which to compute it");
676                exit(1);
677             }
678          }
679          if (Opt->f1 == 0.0f) Opt->f1 = 1.0;
680          if (Opt->f1 < 0) {
681             SUMA_S_Err( "Bad news seeting voxvol");
682             exit(1);
683          }
684          voxthr = SUMA_val_at_count(hh, Opt->wL*1.0e6/Opt->f1, 0, 1);
685          if (Opt->verbose)
686             fprintf(stdout,"Val: %f at %s: %f\n",
687                      voxthr,
688                      Opt->this_xset_name, Opt->wL);
689          else
690             fprintf(stdout,"%f\n", voxthr);
691       }
692    }
693 
694    /* all done, free */
695    Opt = free_SegOpts(Opt);
696    if (rhh) SUMA_Free_hist(rhh);
697    if (hh) SUMA_Free_hist(hh);
698 
699    exit(0);
700 }
701