1 #include "afni.h"
2 #ifndef ALLOW_PLUGINS
ICOR_init(char * lab)3 PLUGIN_interface * ICOR_init(char *lab)
4 {
5 MCW_popup_message( THE_TOPSHELL ,
6 " \n"
7 " InstaCorr not available\n"
8 " since this copy of AFNI\n"
9 " was compiled without\n"
10 " support for plugins!\n " , MCW_USER_KILL ) ;
11 return NULL ;
12 }
13
GICOR_init(char * lab)14 PLUGIN_interface * GICOR_init(char *lab)
15 {
16 MCW_popup_message( THE_TOPSHELL ,
17 " \n"
18 " GrpInCorr not available\n"
19 " since this copy of AFNI\n"
20 " was compiled without\n"
21 " support for plugins!\n " , MCW_USER_KILL ) ;
22 return NULL ;
23 }
24 #else
25
26 /***********************************************************************
27 Pseudo-plugin to setup InstaCorr operations
28 ************************************************************************/
29
30 static int ncall=0 ;
31
32 static unsigned int called_before[26] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0} ;
33
34 /*--------------------- string to 'help' the user --------------------*/
35
36 static char i_helpstring[] =
37 "Purpose: control AFNI InstaCorr operations\n"
38 "\n"
39 "=================== The Two Steps to Using InstaCorr ===================\n"
40 "\n"
41 "(1) Use the Setup controller to prepare for the correlation computations.\n"
42 " * A detailed description of the controls is given farther below.\n"
43 "\n"
44 "(2) In an image viewer window:\n"
45 " * Right-click (or Ctrl-Left-click) to get a popup menu.\n"
46 "\n"
47 " * The top item is 'InstaCorr Set'.\n"
48 " ++ This item will only be enabled when Setup processing is complete.\n"
49 "\n"
50 " * Choosing this item will cause the current crosshair voxel to\n"
51 " be the seed time series that is correlated with all the others.\n"
52 "\n"
53 " * A new functional overlay dataset will be created to show the result,\n"
54 " which comprises the collection of correlation coefficients of each\n"
55 " processed voxel time series with the (processed) seed voxel time series.\n"
56 " ++ In AFNI controller 'A', this dataset will be named 'A_ICOR', etc.\n"
57 " ++ This correlation dataset will NOT be saved to disk automatically.\n"
58 " ++ You can save this dataset to disk with the 'Datamode -> Write OLay'\n"
59 " button in the main AFNI controller window.\n"
60 " ++ Pressing this button more than once will result in this dataset\n"
61 " being over-written by the latest version!\n"
62 "\n"
63 " * Each time you do 'InstaCorr Set', the functional overlay\n"
64 " will be updated to reflect the new correlation map.\n"
65 "\n"
66 " ** Alternative to using 'InstaCorr Set' on the popup menu:\n"
67 " ++ Hold down BOTH the Shift and Control keys on the keyboard\n"
68 " while clicking down the Left mouse button\n"
69 " (i.e., Shift-Ctrl-Left-click).\n"
70 " ++ This combination will jump the crosshairs to the selected point\n"
71 " AND then run 'InstaCorr Set' using this new voxel as the seed.\n"
72 " ++ You can do Shift-Ctrl-Left-click-and-hold, then drag the mouse\n"
73 " cursor around to change the seed location smoothly. The crosshairs\n"
74 " do not move with this smooth re-seeding until you finally release\n"
75 " the mouse button.\n"
76 "\n"
77 "CONTROLS\n"
78 "========\n"
79 "* TimeSeries:\n"
80 " Dataset = time series dataset to auto-correlate\n"
81 " [this dataset does NOT have to be the Underlay]\n"
82 " Start,End= indexes of start and stop times\n"
83 " Examples (assume 100 sub-bricks total):\n"
84 " 0,50 = correlate with first 51 time points\n"
85 " 10 = correlate with time points 10..99\n"
86 " 10+30 = correlate with 30 time points, 10..39\n"
87 " [if End is missing or 0, then it is the last sub-brick]\n"
88 " [N.B.: 'Start,End' replaces 'Ignore', as of Nov 2012]\n"
89 " -->>** NEW FEATURE: MULTIPLE SECTIONS [Oct 2014] **\n"
90 " You can input this field in the following format:\n"
91 " Start@Len,Num,Step\n"
92 " where Start = first index in the dataset to use\n"
93 " Len = length of each section to correlate\n"
94 " Num = number of sections to use (0 == maximum)\n"
95 " Step = step size between sections (between 1 and Len)\n"
96 " [The '@' character is the indicator that ]\n"
97 " [the multiple section information follows]\n"
98 " If you do this, you will get a multiple volume result,\n"
99 " where each volume is the correlation from a limited section\n"
100 " of data. Of course, the program will be slower (less Insta).\n"
101 " ** At present, only works with Pearson correlation :-(\n"
102 " Blur = FWHM in mm of Gaussian blurring to perform\n"
103 " [if a Mask is used, blurring is only inside the mask]\n"
104 "\n"
105 "* Mask:\n"
106 " Automask = Yes to compute an automask from the time series dataset\n"
107 " No to skip this automask step\n"
108 " Dataset = Dataset from which to draw a mask\n"
109 " [this dataset will be ignored if Automask is Yes]\n"
110 " Index = Sub-brick index to use for dataset-derived mask\n"
111 "\n"
112 "* Bandpass:\n"
113 " Lower = Smallest frequency to allow (in Hz) [can be 0]\n"
114 " Upper = Largest frequency to allow (in Hz) [must be > Lower]\n"
115 " [Even if Bandpass is turned off, each voxel time series]\n"
116 " [is detrended against a quadratic polynomial and then ]\n"
117 " [has the 0 and Nyquist frequencies removed; cf. Polort ]\n"
118 " Despike = If this is YES, then the time series have large spikes\n"
119 " filtered out BEFORE the detrending and bandpass operations.\n"
120 " This option is here to let you process datasets that have\n"
121 " a few large spikes, which would otherwise totally dominate\n"
122 " the correlation results.\n"
123 "\n"
124 "* Global Orts:\n"
125 " 1D file = Extra time series to remove from each voxel before\n"
126 " computing the correlations\n"
127 " [These are also bandpassed to avoid re-introducing any]\n"
128 " [of the frequency components rejected by Bandpass. ]\n"
129 " * If Start > 0, and if the 1D file is the same length (or more)\n"
130 " as the input dataset, then the first Start points of this 1D\n"
131 " file will also be ignored.\n"
132 " * If Start > 0, but this file is shorter than the input dataset,\n"
133 " then no initial points in this 1D file will be ignored.\n"
134 " * If the Global Ort file is too short, it will be ignored in toto.\n"
135 " #PC = Number of principal components of the collective input\n"
136 " time series vectors to compute and use as global orts.\n"
137 "\n"
138 "* Misc Opts:\n"
139 " IF environment variable AFNI_INSTACORR_SEEDBLUR is YES\n"
140 " SeedBlur = Extra radius about which to Gaussian blur when extracting\n"
141 " the seed voxel time series.\n"
142 " IF environment variable AFNI_INSTACORR_SEEDBLUR is NO or isn't set\n"
143 " SeedRad = Radius of sphere over which to average when extracting\n"
144 " the seed voxel time series.\n"
145 " These extra averages/blurs are done only inside the mask (if any), and\n"
146 " are applied AFTER the dataset is blurred (if Blur > 0).\n"
147 "\n"
148 " Polort = polynomial detrending level prior to Bandpass (if that's on);\n"
149 " -1 == no detrending [for Ziad]\n"
150 " 0 == mean removal [if you don't like mean things, I guess]\n"
151 " 1 == linear trend removal [for nonlinear thinkers?]\n"
152 " 2 == quadratic trend removal [the default]\n"
153 " You should only change this option from '2' if you understand\n"
154 " what you are doing.\n"
155 " ** If Polort <= 0 AND if Bandpass is turned off, then the FFT\n"
156 " filtering of the time series will not be done at all.\n"
157 " ** There is essentially no reason I can think of to set Polort = -1\n"
158 " AND to turn Bandpass on -- and the results won't be pleasant.\n"
159 " ** If Polort > 0, then FFT filtering will be done even if Bandpass\n"
160 " is off -- but only the 0 and Nyquist frequencies will be removed.\n"
161 " ** To get NO preprocessing of the time series before correlation,\n"
162 " set Polort = -1, turn Bandpass off, and don't input Global Orts.\n"
163 "\n"
164 " Method = Pearson (moment), Spearman (rank), or Quadrant Correlation.\n"
165 " ** For short time series, all methods will be 'Insta'. For longer\n"
166 " datasets, Spearman and Quadrant correlation will be perceptibly\n"
167 " slower than Pearson.\n"
168 "\n"
169 "* Iterate = This row allows you to specify that the results from the usual\n"
170 " InstaCorr process will be iterated -- that is, all the voxels\n"
171 " whose correlation is above a threshold will be averaged, and\n"
172 " then that average will be used as the seed, et cetera.\n"
173 " NOTE: voxels whose correlation is below minus Thresh will be\n"
174 " averaged into the new seed negatively. This choice is\n"
175 " to avoid cancellation.\n"
176 "\n"
177 " Count = Number of iterations, from 2..6.\n"
178 "\n"
179 " Thresh = Threshold at which to choose voxels for the averaging process.\n"
180 " * For 'fun', notice that if you set Thresh to 0.01 (the smallest\n"
181 " value allowed), and set Count to 6, then the resulting map\n"
182 " does not depend much on the starting seed location. That is,\n"
183 " the result converges to something similar to an eigenfunction\n"
184 " of the correlation operator.\n"
185 "\n"
186 "* ExtraSet = If this is chosen, then the seed voxel will be extracted\n"
187 " from Dataset, but its correlations will be done with the\n"
188 " voxels from Extraset.\n"
189 " * The number of time points and the grid spacing of\n"
190 " Extraset must match Dataset.\n"
191 " * The 2 time series dataset will be pre-processed identically.\n"
192 " * This option is a present for Ziad.\n"
193 " * If Iterate is used, the seed averaging at each step still\n"
194 " comes from Dataset, and the correlations are between the\n"
195 " Dataset-derived seed and ExtraSet.\n"
196 "\n"
197 "OPERATION\n"
198 "=========\n"
199 "* Once you have set the controls the way you want, press one of the ''Setup'\n"
200 " buttons, and the program will process (filter & blur) the data time series.\n"
201 "\n"
202 "* When this processing is finished, you will be ready to use 'InstaCorr Set'\n"
203 " (or Shift-Ctrl-Left-click) and have some InstaCorr fun!\n"
204 "\n"
205 "* The 'InstaCorr SeedJump' popup menu item will jump the crosshairs back\n"
206 " to the voxel that is the currently used InstaCorr seed.\n"
207 "\n"
208 "* The current seed time series is saved in the 1D 'library', and can be\n"
209 " plotted in a graph viewer using the 'FIM->Pick Ideal' menu item.\n"
210 "\n"
211 "* If you switch session directories, or switch views (e.g., +orig to +tlrc),\n"
212 " InstaCorr will be disabled and you'll have to use the 'Setup ICorr'\n"
213 " button again to re-initialize the computations.\n"
214 "\n"
215 "Author -- RW Cox -- May 2009\n"
216 ;
217
218 /*----------------- prototypes for internal routines -----------------*/
219
220 static char * ICOR_main( PLUGIN_interface * ) ;
221
222 /***********************************************************************
223 Set up the interface to the user
224 ************************************************************************/
225
ICOR_init(char * lab)226 PLUGIN_interface * ICOR_init( char *lab )
227 {
228 PLUGIN_interface *plint ; /* will be the output of this routine */
229 static char *yn[2] = { "No" , "Yes" } ;
230 /* only the first 4 of these methods is available to a "normal" user */
231 static char *meth_string[10] = { "Pearson" , "Spearman" ,
232 "Quadrant", "Ken Tau_b", "TicTacToe" ,
233 "BCpearson" , "VCpearson", "Euclidian",
234 "CityBlock" , "Quantile:9" } ;
235 char sk[32] , sc[32] ;
236 int gblur = AFNI_yesenv("AFNI_INSTACORR_SEEDBLUR") ;
237
238 if( lab == NULL ) lab = "\0" ;
239
240 /*---------------- set titles and call point ----------------*/
241
242 sprintf(sk,"%sSetup InstaCorr",lab) ;
243 plint = PLUTO_new_interface( "InstaCorr" ,
244 sk ,
245 i_helpstring ,
246 PLUGIN_CALL_VIA_MENU ,
247 (char *(*)())ICOR_main ) ;
248
249 sprintf(sk,"%sSetup+Keep",lab) ; sprintf(sc,"%sSetup+Quit",lab) ;
250 PLUTO_set_runlabels( plint , sk , sc ) ;
251
252 /*--------- make interface lines -----------*/
253
254 PLUTO_add_option ( plint , "TimeSeries" , "TimeSeries" , TRUE ) ;
255 PLUTO_add_dataset( plint , "Dataset" ,
256 ANAT_ALL_MASK , FUNC_ALL_MASK , DIMEN_4D_MASK | BRICK_ALLREAL_MASK ) ;
257 PLUTO_add_string ( plint , "Start,End" , 0,NULL,10 ) ;
258 PLUTO_add_number ( plint , "Blur" , 0,10,0,0,TRUE ) ;
259
260 PLUTO_add_option ( plint , "Mask" , "Mask" , TRUE ) ;
261 PLUTO_add_string ( plint , "Automask" , 2 , yn , 1 ) ;
262 PLUTO_add_dataset( plint , "Dataset" ,
263 ANAT_ALL_MASK , FUNC_ALL_MASK , DIMEN_ALL_MASK | BRICK_ALLREAL_MASK ) ;
264 PLUTO_add_number ( plint , "Index" , 0,ICOR_MAX_FTOP,0,0,TRUE ) ;
265
266 PLUTO_add_option( plint , "Bandpass(Hz)" , "Bandpass" , MAYBE ) ;
267 PLUTO_add_number( plint , "Lower" , 0, 9000,3, 10 , TRUE ) ;
268 PLUTO_add_number( plint , "Upper" , 0,10000,3,100 , TRUE ) ;
269 PLUTO_add_string( plint , "Despike" , 2 , yn , 0 ) ;
270
271 PLUTO_add_option ( plint , "Global Orts" , "GlobalOrts" , FALSE ) ;
272 PLUTO_add_timeseries( plint , "1D file" ) ;
273 PLUTO_add_number ( plint , "#PC" , 0 , 9 , 0 , 0 , TRUE ) ;
274
275 if( AFNI_yesenv("TCSV_EXPERIMENT") ) /* just for fun */
276 PLUTO_add_tcsv( plint , "[tc]sv file" ) ;
277
278 #if 0
279 PLUTO_add_option ( plint , "Slice Orts" , "SliceOrts" , FALSE ) ;
280 PLUTO_add_timeseries( plint , "1D file" ) ;
281 #endif
282
283 PLUTO_add_option( plint , "Misc Opts" , "MiscOpts" , FALSE ) ;
284 if( gblur ) PLUTO_add_number( plint , "SeedBlur" , 0,10,0,0,TRUE ) ;
285 else PLUTO_add_number( plint , "SeedRad" , -10,10,0,0,TRUE ) ;
286 PLUTO_add_number( plint , "Polort" , -1,2,0,2 , FALSE ) ;
287 { char *un = tross_username() ;
288 PLUTO_add_string( plint , "Method" ,
289 ( un != NULL && /* here is where we */
290 (strstr(un,"cox") != NULL || /* choose how many */
291 strstr(un,"ziad") != NULL || /* methods are seen */
292 AFNI_yesenv("AFNI_ICORR_UBER_USER") ) )
293 ? 10 : 4 ,
294 meth_string , 0 ) ;
295 }
296
297 PLUTO_add_option ( plint , "Iterate" , "Iterate" , FALSE ) ; /* 05 Feb 2015 */
298 PLUTO_add_number ( plint , "Count" , 2,6,0,2 , FALSE ) ;
299 PLUTO_add_number ( plint , "Thresh" , 1,9999,2,50, TRUE ) ;
300
301 PLUTO_add_option ( plint , "ExtraSet" , "ExtraSet" , FALSE ) ;
302 PLUTO_add_dataset( plint , "Extraset" ,
303 ANAT_ALL_MASK , FUNC_ALL_MASK , DIMEN_4D_MASK | BRICK_ALLREAL_MASK ) ;
304
305 return plint ;
306 }
307
308 #define ICOR_BIG (9.999f*ICOR_MAX_FTOP)
309
310 /***************************************************************************
311 Main routine for this plugin (will be called from AFNI).
312 If the return string is not NULL, some error transpired, and
313 AFNI will popup the return string in a message box.
314 ****************************************************************************/
315
ICOR_main(PLUGIN_interface * plint)316 static char * ICOR_main( PLUGIN_interface *plint )
317 {
318 char *tag ;
319 float fbot=-1.0f , ftop=ICOR_BIG ;
320 MRI_IMAGE *gortim=NULL ;
321 int gortnpc=0 ; /* 04 Nov 2016 */
322 THD_3dim_dataset *dset=NULL , *mset=NULL , *eset=NULL ;
323 int start=0,end=0 , mindex=0 , automask=0 , qq ; float blur=0.0f , sblur=0.0f ;
324 ICOR_setup *iset ; char *cpt ;
325 Three_D_View *im3d = plint->im3d ;
326 double etim ;
327 int polort = 2 ; /* 26 Feb 2010 */
328 int cmeth = NBISTAT_PEARSON_CORR ;
329 int despike = 0 ;
330 int clen=0,cnum=0,cstep=0 ;
331 int iter_count=0 ; float iter_thresh=0.0f ;
332
333 /*** ncall = 0 ; ***/
334
335 /* check rationality */
336
337 if( !IM3D_OPEN(im3d) || im3d->vwid->func->options_vedit_av->ival != VEDIT_INSTACORR ){
338 XtUnmapWidget(plint->wid->shell); return NULL;
339 }
340
341 /*--- loop over input option lines ---*/
342
343 while(1){
344 tag = PLUTO_get_optiontag(plint) ; /* which line? */
345 if( tag == NULL ) break ; /* none ==> done */
346
347 /** TimeSeries **/
348
349 if( strcmp(tag,"TimeSeries") == 0 ){
350 MCW_idcode *idc ; char *stend ;
351 idc = PLUTO_get_idcode(plint) ;
352 dset = PLUTO_find_dset(idc) ;
353 if( dset == NULL ) ERROR_message("Can't find TimeSeries dataset") ;
354 stend = PLUTO_get_string(plint) ;
355 blur = PLUTO_get_number(plint) ;
356
357 start = end = 0 ; clen = cnum = cstep = 0 ;
358 if( stend != NULL && *stend != '\0' ){
359 char *cpt ;
360 start = (int)strtod(stend,&cpt) ;
361 if( start < 0 ) start = 0 ;
362 while( isspace(*cpt) ) cpt++ ;
363 if( *cpt == ',' || *cpt == '+' ){
364 char qc = *cpt ;
365 if( !isdigit(*cpt) ) cpt++ ;
366 end = (int)strtod(cpt,NULL) ;
367 if( qc == '+' && end > 4 ) end = start + end-1 ;
368 } else if( *cpt == '@' ){ /* '@' stuff from 07 Oct 2014 */
369 clen = (int)strtod(cpt+1,&cpt) ;
370 if( clen < 5 ){
371 WARNING_message("Bad 'Start,End' @ string '%s'",stend) ; clen = 0 ;
372 } else if( *cpt == ',' ){
373 cnum = (int)strtod(cpt+1,&cpt) ;
374 if( *cpt == ',' )
375 cstep = (int)strtod(cpt+1,&cpt) ;
376 }
377 } else if( *cpt != '\0' ){
378 WARNING_message("Don't understand 'Start,End' string '%s'",stend) ;
379 }
380 }
381 continue ;
382 }
383
384 if( strcmp(tag,"Iterate") == 0 ){ /* 05 Feb 2015 */
385 iter_count = PLUTO_get_number(plint) ;
386 iter_thresh = PLUTO_get_number(plint) ;
387 continue ;
388 }
389
390 if( strcmp(tag,"ExtraSet") == 0 ){
391 MCW_idcode *idc ;
392 idc = PLUTO_get_idcode(plint) ;
393 eset = PLUTO_find_dset(idc) ;
394 if( eset == NULL ) ERROR_message("Can't find ExtraSet dataset") ;
395 continue ;
396 }
397
398 /** Mask **/
399
400 if( strcmp(tag,"Mask") == 0 ){
401 MCW_idcode *idc ; char *am ;
402 am = PLUTO_get_string(plint) ; automask = (am[0] == 'Y') ;
403 idc = PLUTO_get_idcode(plint) ; mset = PLUTO_find_dset(idc) ;
404 mindex = PLUTO_get_number(plint) ;
405 if( !automask && mset == NULL )
406 WARNING_message("No Masking selected?!") ;
407 else if( mset != NULL && automask )
408 WARNING_message("Mask dataset disabled when Automask is Yes") ;
409 continue ;
410 }
411
412 /** GlobalOrts **/
413
414 if( strcmp(tag,"GlobalOrts") == 0 ){
415 MRI_IMAGE *qim = PLUTO_get_timeseries(plint) ;
416 NI_element *nel ;
417 if( qim == NULL ) ERROR_message("Ignoring NULL 'Global Orts' time series") ;
418 else gortim = mri_copy(qim) ;
419 gortnpc = (int)PLUTO_get_number(plint) ;
420
421 if( AFNI_yesenv("TCSV_EXPERIMENT") ){ /* just for fun */
422 NI_element *nel = PLUTO_get_tcsv(plint) ;
423 char *vstring = NI_preview_string(nel,9,"--- Just For Fun ---") ;
424 if( vstring != NULL ){
425 (void) MCW_popup_message( THE_TOPSHELL , vstring , MCW_USER_KILL ) ;
426 free( vstring ) ;
427 }
428 }
429
430 continue ;
431 }
432
433 /** Bandpass **/
434
435 if( strcmp(tag,"Bandpass") == 0 ){
436 char *ds ;
437 fbot = PLUTO_get_number(plint) ;
438 ftop = PLUTO_get_number(plint) ;
439 if( fbot >= ftop ) ERROR_message("Ignoring disordered Bandpass frequencies") ;
440 ds = PLUTO_get_string(plint) ; despike = (ds[0] == 'Y') ;
441 continue ;
442 }
443
444 /** MiscOpts **/
445
446 if( strcmp(tag,"MiscOpts") == 0 ){
447 char *cm ;
448 sblur = PLUTO_get_number(plint) ;
449 polort = PLUTO_get_number(plint) ; /* 26 Feb 2010 */
450 cm = PLUTO_get_string(plint) ;
451 switch( cm[0] ){
452 default: cmeth = NBISTAT_PEARSON_CORR ; break ;
453 case 'S': cmeth = NBISTAT_SPEARMAN_CORR ; break ;
454 case 'K': cmeth = NBISTAT_KENDALL_TAUB ; break ;
455 case 'B': cmeth = NBISTAT_BC_PEARSON_M ; break ; /* 07 Mar 2011 */
456 case 'V': cmeth = NBISTAT_BC_PEARSON_V ; break ; /* 07 Mar 2011 */
457 case 'T': cmeth = NBISTAT_TICTACTOE_CORR; break ; /* 30 Mar 2011 */
458 case 'E': cmeth = NBISTAT_EUCLIDIAN_DIST; break ; /* 04 May 2012, ZSS*/
459 case 'C': cmeth = NBISTAT_CITYBLOCK_DIST; break ; /* 04 May 2012, ZSS*/
460 case 'Q':
461 if( cm[3] == 'n' ) cmeth = NBISTAT_QUANTILE_CORR ;
462 else cmeth = NBISTAT_QUADRANT_CORR ;
463 break ;
464 }
465 continue ;
466 }
467
468 /** should never transpire **/
469
470 return "** ICOR_main: table corruption! **" ;
471 }
472
473 /*** check inputs for stoopiditeeze ***/
474
475 if( dset == NULL )
476 return "** No TimeSeries dataset? **" ;
477 if( start >= DSET_NVALS(dset)-2 )
478 return "** 'Start' value is too large **" ;
479
480 if( end <= 0 || end <= start || end >= DSET_NVALS(dset) ) end = DSET_NVALS(dset)-1 ;
481
482 if( clen > 0 ){
483 int ss , nn , nv=DSET_NVALS(dset) ;
484 if( start+clen >= nv ) return "** length is too long **" ;
485 ss = (cstep <= 0 || cstep > clen) ? clen : cstep ;
486 nn = 1 + (nv-start-clen) / ss ;
487 if( cnum > nn || cnum == 0 ) cnum = nn ;
488 INFO_message("Section length=%d number=%d step=%d",clen,cnum,cstep) ;
489 if( cmeth != NBISTAT_PEARSON_CORR ){
490 ININFO_message("section analyses ('@') requires Pearson method") ;
491 cmeth = NBISTAT_PEARSON_CORR ;
492 }
493 }
494
495 if( end-start+1 < 9 ) {
496 WARNING_message("**************************\n"
497 " Too few samples in time series!\n"
498 " I hope you know what you are doing.\n");
499 if (polort >= 0) {
500 /* object even if we can get away with less. Otherwise
501 the < 9 condition has to be amended in
502 thd_bandpass.c's THD_bandpass_vectors() */
503 return "** TimeSeries dataset is way too short for InstaCorr **" ;
504 } else { /* allow it to proceed if series is not extremely short */
505 if ( end-start+1 < 3) {/* too much too little! */
506 return "** TimeSeries dataset is WAY too short for InstaCorr **" ;
507 }
508 }
509 }
510 if( eset != NULL &&
511 ( DSET_NVALS(dset) != DSET_NVALS(eset) || DSET_NVOX(dset) != DSET_NVOX(eset) ) )
512 return "** TimeSeries Dataset and Extraset don't match **" ;
513 if( !automask && mset != NULL && DSET_NVOX(mset) != DSET_NVOX(dset) )
514 return "** Mask dataset doesn't match up with TimeSeries dataset **" ;
515 if( !automask && mset != NULL && mindex >= DSET_NVALS(mset) )
516 return "** Mask dataset index is out of range **" ;
517 if( gortim != NULL && gortim->nx < end-start+1 )
518 return "** Global Orts file is too short for TimeSeries dataset **" ;
519
520 if( fbot >= ftop ){ fbot = 0.0f ; ftop = ICOR_BIG ; }
521 if( fbot < 0.0f ) fbot = 0.0f ;
522
523 if( polort == -1 && (fbot > 0.0f || ftop < ICOR_MAX_FTOP) ) /* 26 Feb 2010 */
524 WARNING_message("Combining Polort=-1 and Bandpass may give peculiar results!") ;
525
526 /** check if only thing changed is sblur -- don't need to re-prepare in that case **/
527 if( im3d->iset != NULL &&
528 im3d->iset->mv != NULL &&
529 im3d->iset->dset == dset &&
530 im3d->iset->eset == eset &&
531 im3d->iset->mset == mset &&
532 im3d->iset->gortim == gortim &&
533 im3d->iset->gortnpc == gortnpc &&
534 im3d->iset->start == start &&
535 im3d->iset->end == end &&
536 im3d->iset->clen == clen &&
537 im3d->iset->cnum == cnum &&
538 im3d->iset->cstep == cstep &&
539 im3d->iset->automask == automask &&
540 im3d->iset->mindex == mindex &&
541 im3d->iset->fbot == fbot &&
542 im3d->iset->ftop == ftop &&
543 im3d->iset->blur == blur &&
544 im3d->iset->despike == despike &&
545 im3d->iset->polort == polort &&
546 THD_instacorr_cmeth_needs_norm(im3d->iset->cmeth) == THD_instacorr_cmeth_needs_norm(cmeth) ){
547
548 INFO_message("InstaCorr setup: minor changes accepted") ;
549 im3d->iset->sblur = sblur ; im3d->iset->cmeth = cmeth ;
550 im3d->iset->iter_count = iter_count; im3d->iset->iter_thresh = iter_thresh;
551 im3d->iset->change = 1; return NULL ;
552 }
553
554 /** (re)create InstaCorr setup **/
555
556 DESTROY_ICOR_setup(im3d->iset) ;
557 INIT_ICOR_setup(iset) ;
558
559 iset->dset = dset ;
560 iset->eset = eset ;
561 iset->mset = (automask) ? NULL : mset ;
562 iset->gortim = gortim ;
563 iset->gortnpc = gortnpc ;
564 iset->start = start ;
565 iset->end = end ;
566 iset->automask = automask ;
567 iset->mindex = mindex ;
568 iset->fbot = fbot ;
569 iset->ftop = ftop ;
570 iset->despike = despike ; /* 14 Oct 2010 */
571 iset->blur = blur ;
572 iset->sblur = sblur ;
573 iset->polort = polort ; /* 26 Feb 2010 */
574 iset->cmeth = cmeth ; /* 01 Mar 2010 */
575 iset->prefix = (char *)malloc(sizeof(char)*16) ;
576 iset->change = 2; /* 07 May 2012 ZSS */
577 iset->clen = clen ; /* 07 Oct 2014 */
578 iset->cnum = cnum ;
579 iset->cstep = cstep ;
580 iset->iter_count = iter_count ; /* 05 Feb 2015 */
581 iset->iter_thresh = iter_thresh ;
582
583 cpt = AFNI_controller_label(im3d); sprintf(iset->prefix,"%c_ICOR",cpt[1]);
584
585 etim = PLUTO_elapsed_time() ;
586
587 /*** prepare the data for InstaCorr Set ***/
588
589 INSTACORR_LABEL_OFF(im3d) ;
590 SHOW_AFNI_PAUSE ;
591 /**************/ qq = THD_instacorr_prepare( iset ) ; /**************/
592 SHOW_AFNI_READY ;
593 if( qq == 0 ){
594 DESTROY_ICOR_setup(iset) ; return "** Error in InstaCorr setup!? **" ;
595 }
596 INSTACORR_LABEL_ON(im3d) ;
597
598 etim = PLUTO_elapsed_time() - etim ;
599
600 INFO_message("InstaCorr setup: %d voxels ready for work: %.2f sec",qq,etim) ;
601 ININFO_message("..... Use 'InstaCorr Set' to pick a seed voxel .....") ;
602 ININFO_message("..... (Mouse-right-click menu in image viewer) .....") ;
603 ININFO_message("..... (or Shift+Ctrl+left-click at seed point) .....") ;
604 ININFO_message("..... (Shift+Ctrl+left-click-and-drag for fun) .....") ;
605
606 im3d->iset = iset ;
607
608 ENABLE_INSTACORR(im3d) ; /* manage the widgets */
609 return NULL ;
610 }
611 #endif /* ALLOW_PLUGINS */
612
613 /*-------------------------------------------------------------------------*/
614 /* Seed location at crosshairs */
615
AFNI_icor_setref(Three_D_View * im3d)616 int AFNI_icor_setref( Three_D_View *im3d )
617 {
618 int ijk ;
619
620 ENTRY("AFNI_icor_setref") ;
621
622 if( !IM3D_OPEN(im3d) ) RETURN(-1) ;
623
624 ijk = AFNI_icor_setref_xyz( im3d , im3d->vinfo->xi ,
625 im3d->vinfo->yj , im3d->vinfo->zk ) ;
626 RETURN(ijk) ;
627 }
628
629 /*-------------------------------------------------------------------------*/
630 /* Seed location at anat dataset i,j,k voxel indexes */
631
AFNI_icor_setref_anatijk(Three_D_View * im3d,int ii,int jj,int kk)632 int AFNI_icor_setref_anatijk( Three_D_View *im3d , int ii,int jj,int kk )
633 {
634 int ijk ;
635 THD_ivec3 iv ; THD_fvec3 fv ;
636
637 ENTRY("AFNI_icor_setref_anatijk") ;
638
639 if( !IM3D_OPEN(im3d) ) RETURN(-1) ;
640
641 /* convert (i,j,k) to DICOM */
642
643 LOAD_IVEC3(iv,ii,jj,kk) ;
644 fv = THD_3dind_to_3dmm ( im3d->anat_now , iv ) ;
645 fv = THD_3dmm_to_dicomm( im3d->anat_now , fv ) ;
646
647 ijk = AFNI_icor_setref_xyz( im3d , fv.xyz[0] , fv.xyz[1] , fv.xyz[2] ) ;
648 RETURN(ijk) ;
649 }
650
651 /*----------------------------------------------------------------------------*/
652 /* Find the best index triple in the odset to match the uxyz
653 coordinates in the udset [07 May 2015].
654 *//*--------------------------------------------------------------------------*/
655
THD_find_closest_roundtrip(THD_3dim_dataset * odset,THD_3dim_dataset * udset,THD_fvec3 uxyz)656 THD_ivec3 THD_find_closest_roundtrip( THD_3dim_dataset *odset,
657 THD_3dim_dataset *udset, THD_fvec3 uxyz )
658 {
659 THD_ivec3 iv,jv,kv , ivbest ; THD_fvec3 xv,yv,zv ;
660 int di,dj,dk ; float dist , dbest=666666.6f ;
661 int nxo=DSET_NX(odset) ;
662 int nyo=DSET_NY(odset) ;
663 int nzo=DSET_NZ(odset) ;
664
665 xv = THD_dicomm_to_3dmm ( odset, uxyz ) ;
666 iv = THD_3dmm_to_3dind_no_wod ( odset, xv ) ; ivbest = iv ;
667 for( di=-1 ; di <= 1 ; di++ ){
668 for( dj=-1 ; dj <= 1 ; dj++ ){
669 for( dk=-1 ; dk <= 1 ; dk++ ){
670 jv.ijk[0] = iv.ijk[0] + di ; if( jv.ijk[0] < 0 || jv.ijk[0] >= nxo ) continue ;
671 jv.ijk[1] = iv.ijk[1] + dj ; if( jv.ijk[1] < 0 || jv.ijk[1] >= nyo ) continue ;
672 jv.ijk[2] = iv.ijk[2] + dk ; if( jv.ijk[2] < 0 || jv.ijk[2] >= nzo ) continue ;
673 yv = THD_3dind_to_dicomm_no_wod( odset, jv ) ;
674 xv = THD_dicomm_to_3dmm ( udset, yv ) ;
675 kv = THD_3dmm_to_3dind ( udset, xv ) ;
676 xv = THD_3dind_to_3dmm ( udset, kv ) ;
677 yv = THD_3dmm_to_dicomm ( udset, xv ) ;
678
679 dist = fabsf(uxyz.xyz[0]-yv.xyz[0])+fabsf(uxyz.xyz[1]-yv.xyz[1])+fabsf(uxyz.xyz[2]-yv.xyz[2]) ;
680 if( dist < dbest ){
681 ivbest = jv ; dbest = dist ;
682 }
683 #if 0
684 INFO_message("roundtrip: input xyz=%f %f %f output xyz=%f %f %f dist=%f dijk=%d %d %d %s",
685 uxyz.xyz[0] , uxyz.xyz[1] , uxyz.xyz[2] ,
686 yv.xyz[0] , yv.xyz[1] , yv.xyz[2] ,
687 fabsf(uxyz.xyz[0]-yv.xyz[0])+fabsf(uxyz.xyz[1]-yv.xyz[1])+fabsf(uxyz.xyz[2]-yv.xyz[2]) ,
688 di,dj,dk , (dist==dbest) ? "*" : "\0" ) ;
689 #endif
690 }}}
691
692 #if 0
693 INFO_message("iv nominal=%d %d %d ivbest=%d %d %d",
694 iv.ijk[0] , iv.ijk[1] , iv.ijk[2] , ivbest.ijk[0] , ivbest.ijk[1] , ivbest.ijk[2] ) ;
695 #endif
696
697 return ivbest ;
698 }
699
700 /*----------------------------------------------------------------------------*/
701 /* Seed location at DICOM x,y,z location */
702
AFNI_icor_setref_xyz(Three_D_View * im3d,float xx,float yy,float zz)703 int AFNI_icor_setref_xyz( Three_D_View *im3d , float xx,float yy,float zz )
704 {
705 MRI_IMAGE *iim=NULL; float *iar, rng; THD_fvec3 iv,jv; THD_ivec3 kv; int ijk,ic ;
706 THD_3dim_dataset *icoset ; THD_slist_find slf ; int nds=0 ;
707 double etim ;
708 MRI_IMARR *iimar=NULL; int nim=0 , qim ;
709
710 ENTRY("AFNI_icor_setref_xyz") ;
711
712 if( !IM3D_OPEN(im3d) ) RETURN(-1) ;
713
714 ic = AFNI_controller_index(im3d) ; if( ic < 0 || ic > 25 ) RETURN(-1) ;
715
716 /**** divert to Group InstaCorr? ****/
717
718 if( im3d->giset != NULL && im3d->giset->ready ){
719 ijk = AFNI_gicor_setref_xyz(im3d,xx,yy,zz) ; RETURN(ijk) ;
720 }
721
722 if( !ISVALID_ICOR_setup(im3d->iset) ) RETURN(-1) ;
723
724 /* find where we are working from, in dataset coordinates */
725
726 LOAD_FVEC3( iv , xx,yy,zz ) ;
727 jv = THD_dicomm_to_3dmm( im3d->iset->dset, iv ) ;
728
729 /* test if this point is possibly OK */
730
731 if( jv.xyz[0] < im3d->iset->dset->daxes->xxmin ||
732 jv.xyz[0] > im3d->iset->dset->daxes->xxmax ||
733 jv.xyz[1] < im3d->iset->dset->daxes->yymin ||
734 jv.xyz[1] > im3d->iset->dset->daxes->yymax ||
735 jv.xyz[2] < im3d->iset->dset->daxes->zzmin ||
736 jv.xyz[2] > im3d->iset->dset->daxes->zzmax ){
737
738 WARNING_message("InstaCorr set point (%.1f,%.1f,%.1f) outside dataset box",
739 xx,yy,zz ) ;
740 RETURN(-1) ;
741 }
742
743 /* convert to index in the InstaCorr dataset */
744
745 #if 0
746 kv = THD_3dmm_to_3dind_no_wod( im3d->iset->dset, jv ) ;
747 #else
748 kv = THD_find_closest_roundtrip(im3d->iset->dset,im3d->anat_now,iv) ;
749 #endif
750 ijk = DSET_ixyz_to_index( im3d->iset->dset, kv.ijk[0],kv.ijk[1],kv.ijk[2] ) ;
751
752 #if 0
753 INFO_message("AFNI_icor_setref_xyz: iv=%.3f,%.3f,%.3f kv=%d,%d,%d ijk=%d",
754 iv.xyz[0],iv.xyz[1],iv.xyz[2] ,
755 kv.ijk[0],kv.ijk[1],kv.ijk[2] , ijk ) ;
756 #endif
757
758 /* do the real work: ijk = voxel index */
759
760 etim = PLUTO_elapsed_time() ;
761
762 if( im3d->iset->clen <= 0 ){
763 iim = THD_instacorr( im3d->iset , ijk ) ;
764 if( iim == NULL ) RETURN(-1) ; /* did it fail? */
765 INIT_IMARR(iimar) ; ADDTO_IMARR(iimar,iim) ;
766 } else {
767 if( im3d->iset->cnum > 4 ) SHOW_AFNI_PAUSE ;
768 iimar = THD_instacorr_collection( im3d->iset , ijk ) ;
769 if( im3d->iset->cnum > 4 ) SHOW_AFNI_READY ;
770 if( iimar == NULL ) RETURN(-1) ; /* did it fail? */
771 }
772 nim = IMARR_COUNT(iimar) ;
773
774 if( ncall <= 1 )
775 ININFO_message(" InstaCorr elapsed time = %.2f sec: correlations",
776 PLUTO_elapsed_time()-etim) ;
777
778 /* 17 Mar 2010: save seed location we just did in the im3d struct */
779
780 im3d->vinfo->xi_icor = xx ; /* input DICOM coords */
781 im3d->vinfo->yj_icor = yy ;
782 im3d->vinfo->zk_icor = zz ;
783
784 /* saving 3D indexes in anat space also */
785 kv = THD_3dmm_to_3dind ( im3d->anat_now , jv ) ;
786 UNLOAD_IVEC3( kv , im3d->vinfo->i1_icor ,
787 im3d->vinfo->j2_icor , im3d->vinfo->k3_icor ) ;
788
789 /* find the output dataset */
790
791 slf = THD_dset_in_session( FIND_PREFIX , im3d->iset->prefix , im3d->ss_now ) ;
792
793 /* if it doesn't exist, or is not the right grid, create it now */
794
795 if( !ISVALID_DSET (slf.dset) ||
796 !EQUIV_DATAXES(slf.dset->daxes,im3d->iset->dset->daxes) ||
797 DSET_NVALS(slf.dset) != nim ){
798
799 icoset = EDIT_empty_copy( im3d->iset->dset ) ; /* make new dataset */
800 EDIT_dset_items( icoset ,
801 ADN_prefix , im3d->iset->prefix ,
802 ADN_nvals , nim ,
803 ADN_ntt , 0 ,
804 ADN_func_type , FUNC_BUCK_TYPE ,
805 ADN_type , HEAD_FUNC_TYPE ,
806 ADN_datum_all , MRI_float ,
807 ADN_none ) ;
808 DSET_superlock(icoset) ;
809
810 if( slf.dset != NULL ){ /* exists, but isn't right for us */
811
812 MCW_idcode old_idc = slf.dset->idcode ;
813 THD_delete_3dim_dataset(slf.dset,True) ; /* destroy the guts */
814 *slf.dset = *icoset ; /* copy the guts, keep the pointer */
815 slf.dset->idcode = old_idc ; /* and keep the idcode */
816 nds = slf.dset_index ;
817 INFO_message("trashed and re-used old dataset %s",im3d->iset->prefix) ;
818
819 } else { /* add to the session */
820 int vv = icoset->view_type ;
821 nds = im3d->ss_now->num_dsset ;
822 SET_SESSION_DSET(icoset, im3d->ss_now, nds, vv);
823 im3d->ss_now->num_dsset++ ;
824 AFNI_force_adoption( im3d->ss_now , False ) ;
825 AFNI_make_descendants( GLOBAL_library.sslist ) ;
826 INFO_message("created new dataset %s",im3d->iset->prefix) ;
827 }
828
829 /* just need to use existing dataset that matches */
830
831 } else {
832 icoset = slf.dset ; nds = slf.dset_index ;
833 DSET_mallocize(icoset); /* make sure not mmap-ed file ZSS */
834 }
835 icoset->dblk->diskptr->allow_directwrite = 1 ;
836
837 /* save the result into the output dataset */
838
839 for( qim=0 ; qim < nim ; qim++ ){
840 iim = IMARR_SUBIM(iimar,qim) ;
841 if( iim != NULL ){
842 iar = MRI_FLOAT_PTR(iim) ;
843 EDIT_substitute_brick( icoset , qim , MRI_float , iar ) ;
844 mri_clear_data_pointer(iim) ;
845 } else {
846 EDIT_substitute_brick( icoset , qim , MRI_float , NULL ) ;
847 }
848 }
849 DESTROY_IMARR(iimar) ; iim = NULL ;
850 DSET_KILL_STATS(icoset) ; THD_load_statistics(icoset) ;
851
852 /* 03 May 2010: add some attributes that say where this comes from */
853
854 { char buf[64] ;
855 THD_set_string_atr( icoset->dblk , "INSTACORR_PARENT" , DSET_HEADNAME(im3d->iset->dset) ) ;
856 sprintf(buf,"%d,%d,%d",im3d->vinfo->i1_icor,im3d->vinfo->j2_icor,im3d->vinfo->k3_icor) ;
857 THD_set_string_atr( icoset->dblk , "INSTACORR_SEEDIJK" , buf ) ;
858 if( im3d->iset->eset != NULL )
859 THD_set_string_atr( icoset->dblk, "INSTACORR_EXTRASET", DSET_HEADNAME(im3d->iset->eset) );
860 }
861
862 switch (im3d->iset->cmeth) {
863 case NBISTAT_EUCLIDIAN_DIST:
864 EDIT_BRICK_LABEL (icoset,0,"Inv.Euc.Dist") ;
865 rng = 10.0;
866 break;
867 case NBISTAT_CITYBLOCK_DIST:
868 EDIT_BRICK_LABEL (icoset,0,"Inv.City.Dist") ;
869 rng = 10.0;
870 break;
871 default:
872 EDIT_BRICK_LABEL (icoset,0,"Correlation") ;
873 rng = 0.7;
874 break;
875 }
876
877 for( qim=0 ; qim < nim ; qim++ )
878 EDIT_BRICK_TO_FICO(icoset,qim,im3d->iset->mv->nvals,1,im3d->iset->ndet) ;
879
880 DSET_BRICK_FDRCURVE_ALLKILL(icoset) ;
881 DSET_BRICK_MDFCURVE_ALLKILL(icoset) ;
882 flush_3Dview_sort(im3d,"T"); /* ZSS April 27 2012: Reset sorted threshold */
883
884 if( ncall <= 1 )
885 ININFO_message(" InstaCorr elapsed time = %.2f sec: dataset ops" ,
886 PLUTO_elapsed_time()-etim ) ;
887
888 if( AFNI_yesenv("AFNI_INSTACORR_FDR") ){
889 THD_create_all_fdrcurves(icoset) ;
890 if( ncall <= 1 )
891 ININFO_message(" InstaCorr elapsed time = %.2f sec: FDR curve" ,
892 PLUTO_elapsed_time()-etim ) ;
893 }
894
895 /* 10 May 2009: save seed timeseries into timeseries library */
896
897 if( im3d->iset->tseed != NULL ){
898 MRI_IMAGE *tsim ; float *tsar ;
899 tsim = mri_new( im3d->iset->mv->nvals + im3d->iset->start,1,MRI_float ) ;
900 tsar = MRI_FLOAT_PTR(tsim) ;
901 memcpy( tsar + im3d->iset->start , im3d->iset->tseed ,
902 sizeof(float)*im3d->iset->mv->nvals ) ;
903 tsim->name = (char *)malloc(sizeof(char)*16) ;
904 strcpy(tsim->name,im3d->iset->prefix) ; strcat(tsim->name,"_seed") ;
905 AFNI_replace_timeseries(tsim) ;
906 }
907
908 /* redisplay overlay */
909
910 if( called_before[ic] ) AFNI_ignore_pbar_top(1) ; /* 03 Jun 2014 */
911 if( im3d->fim_now != icoset || im3d->iset->change ){ /* switch to this dataset */
912 MCW_choose_cbs cbs ; char cmd[32] , *cpt=AFNI_controller_label(im3d) ;
913 cbs.ival = nds ;
914
915 AFNI_finalize_dataset_CB( im3d->vwid->view->choose_func_pb ,
916 (XtPointer)im3d , &cbs ) ;
917 AFNI_set_fim_index(im3d,0) ;
918 AFNI_set_thr_index(im3d,0) ;
919
920 if( !called_before[ic] ){
921 sprintf(cmd,"SET_FUNC_RANGE %c.%.2f",cpt[1], rng) ;
922 AFNI_driver(cmd) ;
923 }
924 }
925 if( MCW_val_bbox(im3d->vwid->func->range_bbox) ){
926 char cmd[32] , *cpt=AFNI_controller_label(im3d) ;
927 AFNI_ignore_pbar_top(0) ;
928 sprintf(cmd,"SET_FUNC_RANGE %c.%f",cpt[1],im3d->vinfo->fim_autorange) ;
929 AFNI_driver(cmd) ;
930 AFNI_ignore_pbar_top(1) ;
931 }
932 AFNI_reset_func_range(im3d) ; called_before[ic]++ ;
933 AFNI_ignore_pbar_top(0) ;
934
935 IM3D_CLEAR_TMASK(im3d) ; /* Mar 2013 */
936 IM3D_CLEAR_THRSTAT(im3d) ; /* 12 Jun 2014 */
937 if( VEDIT_good(im3d->vedset) ) im3d->vedset.flags = 1 ; /* 18 Jun 2014 */
938 if( MCW_val_bbox(im3d->vwid->view->see_func_bbox) == 0 ){ /* overlay is off */
939 char cmd[32] , *cpt=AFNI_controller_label(im3d) ;
940 sprintf(cmd,"SEE_OVERLAY %c.+",cpt[1]) ;
941 AFNI_driver(cmd) ;
942 } else { /* overlay is on */
943 AFNI_redisplay_func(im3d) ;
944 }
945 AFNI_set_thr_pval(im3d) ; AFNI_process_drawnotice(im3d) ;
946
947 im3d->iset->change = 0; /* reset change flag */
948
949 if( ncall <= 1 )
950 ININFO_message(" InstaCorr elapsed time = %.2f sec: redisplay" ,
951 PLUTO_elapsed_time()-etim ) ;
952 ncall++ ; RETURN(1) ;
953 }
954
955 /*---------------------------------------------------------------------------*/
956 /* Set seeds at same x,y,z location as the 'master' viewer.
957 This locking is only for individual dataset InstaCorr,
958 since only controller A can have Group InstaCorr running.
959 *//*-------------------------------------------------------------------------*/
960
AFNI_icor_setref_locked(Three_D_View * im3d)961 void AFNI_icor_setref_locked( Three_D_View *im3d )
962 {
963 Three_D_View *qq3d ; int ii,cc , glock ; static int busy=0 ;
964
965 ENTRY("AFNI_icor_setref_locked") ;
966
967 if( !IM3D_OPEN(im3d) ) EXRETURN ;
968 if( im3d->giset != NULL && im3d->giset->ready ) EXRETURN ;
969
970 glock = GLOBAL_library.controller_lock ;
971
972 if( busy ) EXRETURN ; /* routine already busy */
973 if( glock == 0 ) EXRETURN ; /* nothing to do */
974 if( GLOBAL_library.ignore_lock ) EXRETURN ; /* ordered not to do anything */
975
976 ii = AFNI_controller_index(im3d) ; /* which one am I? */
977 if( ii < 0 ) EXRETURN ; /* bad input: shouldn't happen */
978 if( ((1<<ii) & glock) == 0 ) EXRETURN ; /* input not locked */
979
980 busy = 1 ;
981
982 for( cc=0 ; cc < MAX_CONTROLLERS ; cc++ ){
983 qq3d = GLOBAL_library.controllers[cc] ; /* controller */
984 if( IM3D_OPEN(qq3d) && qq3d != im3d && ((1<<cc) & glock) != 0 )
985 (void)AFNI_icor_setref_xyz( qq3d ,
986 im3d->vinfo->xi_icor,
987 im3d->vinfo->yj_icor, im3d->vinfo->zk_icor);
988 }
989
990 busy = 0 ; EXRETURN ;
991 }
992
993 /******************************************************************************/
994 /*****--- Group InstaCorr stuff below here! ------------------------------*****/
995 /******************************************************************************/
996
997 #define GIQUIT \
998 do { free(im3d->giset); im3d->giset = NULL; EXRETURN; } while(0)
999
1000 static PLUGIN_interface *GICOR_plint = NULL ; /* 13 May 2010 */
1001
1002 #define GICOR_MAX_NSTAT 5 /* 14 May 2010 */
1003 #define GICOR_BASE_STAT 2
1004
1005 static void GICOR_refit_stat_menus(void) ;
1006
1007 #undef GIC_ALLOW_TTEST
1008
1009 #undef GIC_ALLOW_CLUST /* 19 Oct 2010 -- remove clustering */
1010
1011 /*---------------------------------------------------------------------------*/
1012 /*-- Called from afni_niml.c when 3dGroupInCorr sends a setup NIML element --*/
1013
GICOR_setup_func(NI_stream nsg,NI_element * nel)1014 void GICOR_setup_func( NI_stream nsg , NI_element *nel )
1015 {
1016 GICOR_setup *giset ;
1017 char *atr , *pre ;
1018 Three_D_View *im3d = A_CONTROLLER ; /* global variable */
1019 THD_slist_find sf ;
1020 THD_session *ss = im3d->ss_now ; int qs = ss->num_dsset , vv,qq ;
1021 THD_3dim_dataset *dset ; int nvals=2 ;
1022 static char *blab[6] = { "GIC_Delta" , "GIC_Zscore" ,
1023 "AAA_Delta" , "AAA_Zscore" ,
1024 "BBB_Delta" , "BBB_Zscore" } ;
1025 NI_str_array *labar=NULL ;
1026
1027 ENTRY("GICOR_setup_func") ;
1028
1029 /* if we are ready, nothing more to do (don't setup twice) */
1030
1031 if( im3d->giset != NULL && im3d->giset->ready ) EXRETURN ;
1032
1033 /* if can't shoehorn in another dataset, give up in disgust */
1034
1035 if( qs >= THD_MAX_SESSION_SIZE ){
1036 (void) MCW_popup_message( im3d->vwid->picture ,
1037 " \n"
1038 " ******* AFNI: ******* \n"
1039 " Can't use GrpInCorr \n"
1040 " because dataset table \n"
1041 " is completely full! \n " ,
1042 MCW_USER_KILL | MCW_TIMER_KILL ) ;
1043 NI_stream_closenow(nsg) ;
1044 DISABLE_GRPINCORR(im3d) ;
1045 SENSITIZE_INSTACORR(im3d,False) ;
1046 DESTROY_GICOR_setup(im3d->giset) ;
1047 }
1048
1049 /* create or clear out Group InstaCorr setup struct */
1050
1051 if( im3d->giset == NULL ){
1052 im3d->giset = (GICOR_setup *)calloc(1,sizeof(GICOR_setup)) ;
1053 } else {
1054 memset(im3d->giset,0,sizeof(GICOR_setup)) ; /* fixed an oopsie */
1055 }
1056 giset = im3d->giset ;
1057
1058 giset->ns = nsg ; /* save socket for I/O back to 3dGroupInCorr */
1059 giset->ready = 0 ; /* not ready yet */
1060 giset->busy = 0 ; /* not busy yet, either [18 Mar 2010] */
1061
1062 giset->apair = 0 ; /* Apr 2013 */
1063
1064 /* set various parameters from the NIML header */
1065
1066 atr = NI_get_attribute( nel , "ndset_A" ) ; if( atr == NULL ) GIQUIT;
1067 giset->ndset_A = (int)strtod(atr,NULL) ; if( giset->ndset_A < 2 ) GIQUIT;
1068
1069 atr = NI_get_attribute( nel , "ndset_B" ) ; if( atr == NULL ) GIQUIT;
1070 giset->ndset_B = (int)strtod(atr,NULL) ;
1071
1072 atr = NI_get_attribute( nel , "nvec" ) ; if( atr == NULL ) GIQUIT;
1073 giset->nvec = (int)strtod(atr,NULL) ; if( giset->nvec < 2 ) GIQUIT;
1074
1075 atr = NI_get_attribute( nel , "seedrad" ) ;
1076 if( atr != NULL ) giset->seedrad = (float)strtod(atr,NULL) ;
1077
1078 #ifdef GIC_ALLOW_TTEST
1079 atr = NI_get_attribute( nel , "ttest_opcode" ) ;
1080 if( atr != NULL ) giset->ttest_opcode = (int)strtod(atr,NULL) ;
1081 #endif
1082
1083 atr = NI_get_attribute( nel , "apair") ; /* Apr 2013 */
1084 if( YESSISH(atr) ) GICOR_set_apair_allow_bit(giset) ;
1085
1086 /* create output dataset, to be filled in from 3dGroupInCorr data later */
1087
1088 atr = NI_get_attribute( nel , "geometry_string" ); if( atr == NULL ) GIQUIT;
1089 pre = NI_get_attribute( nel , "target_name" ) ;
1090 if( pre == NULL || *pre == '\0' ) pre = "A_GRP_ICORR" ;
1091 dset = giset->dset = EDIT_geometry_constructor( atr , pre ) ;
1092 if( dset == NULL ) GIQUIT;
1093
1094 atr = NI_get_attribute( nel , "target_nvals" ) ;
1095 if( atr != NULL ){ nvals = (int)strtod(atr,NULL); nvals = MAX(1,nvals); }
1096 vv = AFNI_yesenv("AFNI_GROUPINCORR_ORIG") ;
1097 EDIT_dset_items( dset , ADN_nvals , nvals ,
1098 ADN_view_type , (vv) ? VIEW_ORIGINAL_TYPE
1099 : VIEW_TALAIRACH_TYPE ,
1100 ADN_brick_fac , NULL ,
1101 ADN_none ) ;
1102
1103 atr = NI_get_attribute( nel , "target_labels" ) ;
1104 if( atr != NULL )
1105 labar = NI_decode_string_list( atr , ";" ) ;
1106
1107 /* for each sub-brick in the dataset-to-be */
1108
1109 for( vv=0 ; vv < nvals ; vv++ ){
1110 EDIT_substitute_brick( dset, vv, MRI_float, NULL ) ; /* calloc sub-brick */
1111 if( labar != NULL && vv < labar->num ) /* and label-ize it */
1112 EDIT_BRICK_LABEL( dset , vv , labar->str[vv] ) ;
1113 else if( vv < 6 )
1114 EDIT_BRICK_LABEL( dset , vv , blab[vv] ) ;
1115 if( strstr( DSET_BRICK_LAB(dset,vv) , "_Zsc" ) != NULL )
1116 EDIT_BRICK_TO_FIZT(dset,vv) ; /* mark as a Z score */
1117 }
1118 DSET_superlock( dset ) ;
1119 giset->nvox = DSET_NVOX(dset) ;
1120
1121 if( labar != NULL ) NI_delete_str_array(labar) ; /* 14 May 2010 */
1122
1123 /* 14 May 2010: set various labels */
1124
1125 giset->label_AAA = giset->label_BBB = giset->toplabel = NULL ;
1126
1127 atr = NI_get_attribute( nel , "label_AAA") ;
1128 if( atr != NULL ) giset->label_AAA = strdup(atr) ;
1129
1130 atr = NI_get_attribute( nel , "label_BBB") ;
1131 if( atr != NULL ) giset->label_BBB = strdup(atr) ;
1132
1133 if( giset->label_AAA != NULL ){
1134 char *tlab = giset->toplabel ;
1135 if( tlab == NULL )
1136 giset->toplabel = tlab = (char *)malloc(sizeof(char)*256) ;
1137 if( giset->label_BBB == NULL )
1138 sprintf( tlab , "GrpInCorr: set AAA=%s" ,
1139 giset->label_AAA ) ;
1140 else
1141 sprintf( tlab , "GrpInCorr: set AAA=%s set BBB=%s" ,
1142 giset->label_AAA,giset->label_BBB ) ;
1143 PLUTO_set_toplabel( GICOR_plint , tlab ) ;
1144 }
1145
1146 giset->num_stat_available = 0 ;
1147 giset->lab_stat_available = NULL ;
1148 #if 0
1149 atr = NI_get_attribute( nel , "stats_available" ) ;
1150 labar = NI_decode_string_list( atr , ";" ) ;
1151 if( labar != NULL && labar->num > 0 ){
1152 giset->num_stat_available = labar->num ;
1153 giset->lab_stat_available = (char **)malloc(sizeof(char *)*labar->num) ;
1154 for( vv=0 ; vv < labar->num ; vv++ )
1155 giset->lab_stat_available[vv] = strdup(labar->str[vv]) ;
1156 NI_delete_str_array(labar) ;
1157 }
1158 #endif
1159
1160 /* add dataset to current session (change name if necessary) */
1161
1162 sf = THD_dset_in_session( FIND_PREFIX , pre , ss ) ;
1163 if( sf.dset != NULL ){
1164 int jj,nn ; char *npre ;
1165 nn = strlen(pre); npre = malloc(sizeof(char)*(nn+8)); strcpy(npre,pre);
1166 for( jj=1 ; ; jj++ ){ /* loop until we find something that works */
1167 sprintf(npre+nn,"_%d",jj) ;
1168 sf = THD_dset_in_session( FIND_PREFIX , npre , ss ) ;
1169 if( sf.dset == NULL ) break ;
1170 }
1171 EDIT_dset_items( dset , ADN_prefix , npre , ADN_none ) ;
1172 free(npre) ;
1173 }
1174
1175 POPDOWN_strlist_chooser ;
1176
1177 vv = dset->view_type ;
1178 if( vv < FIRST_VIEW_TYPE || vv > LAST_VIEW_TYPE ) vv = FIRST_VIEW_TYPE;
1179 /* null all other session dataset entries for this dataset */
1180 SET_SESSION_DSET(dset, ss, qs, vv);
1181 for( qq=FIRST_VIEW_TYPE ; qq <= LAST_VIEW_TYPE ; qq++ )
1182 if(qq!=vv)
1183 SET_SESSION_DSET(NULL, ss, qs, qq);
1184 ss->num_dsset++ ; giset->nds = qs ;
1185 giset->session = ss ;
1186
1187 UNDUMMYIZE ;
1188
1189 /* list of voxels to expect from each 3dGroupInCorr data */
1190
1191 if( nel->vec_len == 0 || nel->vec_num == 0 || nel->vec == NULL ){ /* all */
1192 giset->ivec = NULL ; giset->nivec = 0 ;
1193 /* INFO_message("DEBUG: GICOR_setup_func has ivec=NULL") ; */
1194 } else { /* make index list of voxels */
1195 int ii , nn , *iv=(int *)nel->vec[0] ;
1196 giset->ivec = (int *)calloc(sizeof(int),giset->nvec) ;
1197 nn = MIN(giset->nvec,nel->vec_len) ; giset->nivec = nn ;
1198 for( ii=0 ; ii < nn ; ii++ ) giset->ivec[ii] = iv[ii] ;
1199 /* INFO_message("DEBUG: GICOR_setup_func has ivec=int[%d]",nn) ; */
1200 }
1201
1202 /* 23 May 2012: extra string attributes to set? */
1203
1204 { ATR_string *aatr ; int nn ;
1205 char aaname[THD_MAX_NAME], *aastr, *nnatr, nnam[128], *cpt ;
1206
1207 for( nn=0 ; ; nn++ ){
1208 sprintf(nnam,"string_attribute_%06d",nn) ;
1209 nnatr = NI_get_attribute( nel , nnam ) ;
1210 if( nnatr == NULL || *nnatr == '\0' ) break ;
1211 cpt = strstr(nnatr," ==> ") ;
1212 if( cpt == NULL || cpt == nnatr || cpt-nnatr > 256 ) continue ;
1213 strncpy(aaname,nnatr,cpt-nnatr) ; aaname[cpt-nnatr] = '\0' ;
1214 cpt += 5 ; if( *cpt == '\0' ) continue ;
1215 aatr = (ATR_string *)XtMalloc(sizeof(ATR_string)) ;
1216 aatr->type = ATR_STRING_TYPE ;
1217 aatr->name = XtNewString(aaname) ;
1218 aatr->nch = strlen(cpt+1) ;
1219 aatr->ch = (char *)XtMalloc( sizeof(char) * aatr->nch ) ;
1220 memcpy( aatr->ch , cpt , sizeof(char) * aatr->nch ) ;
1221 THD_insert_atr( dset->dblk , (ATR_any *)aatr ) ;
1222 }
1223 }
1224
1225 giset->ready = 1 ; /* that is, ready to ROCK AND ROLL */
1226 GRPINCORR_LABEL_ON(im3d) ;
1227 SENSITIZE_INSTACORR(im3d,True) ;
1228
1229 #undef PUTENV
1230 #define PUTENV(nm,val) do{ if( getenv((nm)) == NULL ){ \
1231 char *str = (char *)malloc(256) ; \
1232 strcpy(str,(nm)); strcat(str,"="); \
1233 strcat(str,val); putenv(str); \
1234 }} while(0)
1235
1236 PUTENV("AFNI_THRESH_LOCK","VALUE") ;
1237 AFNI_set_all_thrlock_bboxes(NULL, -1);
1238 PUTENV("AFNI_RANGE_LOCK" ,"YES" ) ;
1239
1240 /* some messages to the screen */
1241
1242 INFO_message("Added 3dGroupInCorr dataset '%s' to controller %s",
1243 DSET_FILECODE(dset), AFNI_controller_label(im3d) ) ;
1244 ININFO_message("%d datasets in set A",giset->ndset_A) ;
1245 if( giset->ndset_B > 0 )
1246 ININFO_message("%d datasets in set B",giset->ndset_B) ;
1247
1248 ININFO_message("----------------------------------------------------") ;
1249 ININFO_message("----- AFNI is now connnected to 3dGroupInCorr! -----") ;
1250 ININFO_message("..... Use 'InstaCorr Set' to pick a seed voxel .....") ;
1251 ININFO_message("..... (Mouse-right-click menu in image viewer) .....") ;
1252 ININFO_message("----------------------------------------------------") ;
1253
1254 IM3D_CLEAR_TMASK(im3d) ; /* Mar 2013 */
1255 IM3D_CLEAR_THRSTAT(im3d) ; /* 12 Jun 2014 */
1256 GICOR_refit_stat_menus() ; /* 14 May 2010 */
1257
1258 if( im3d->vwid->func->options_vedit_av->ival != VEDIT_GRINCORR ){
1259 AV_assign_ival( im3d->vwid->func->options_vedit_av , VEDIT_GRINCORR ) ;
1260 AFNI_vedit_CB( im3d->vwid->func->options_vedit_av , im3d ) ;
1261 }
1262
1263 /* message for newbie users! [26 Apr 2016] */
1264
1265 MCW_popup_message( im3d->vwid->imag->topper ,
1266 "3dGroupInCorr is ready!\n"
1267 "* Use InstaCorr Set to\n"
1268 " choose a seed voxel\n"
1269 " (Mouse-right-click )\n"
1270 " (in an image viewer)\n"
1271 "* Or press Ctrl+Shift+\n"
1272 " Mouse-left-click in\n"
1273 " an image viewer. "
1274 ,
1275 MCW_USER_KILL | MCW_TIMER_KILL ) ;
1276
1277 EXRETURN ;
1278 }
1279
1280 /*----------------------------------------------------------------------------*/
1281
GICOR_process_message(NI_element * nel)1282 void GICOR_process_message( NI_element *nel ) /* Apr 2013 */
1283 {
1284 Three_D_View *im3d = A_CONTROLLER ;
1285
1286 ENTRY("GICOR_process_message") ;
1287 if( !IM3D_OPEN(im3d) || im3d->giset == NULL ) EXRETURN ;
1288 im3d->giset->busy = 0 ;
1289 process_NIML_textmessage(nel) ;
1290 EXRETURN ;
1291 }
1292
1293 /*----------------------------------------------------------------------------*/
1294 /* Called from afni_niml.c to process the received dataset from 3dGroupInCorr */
1295
GICOR_process_dataset(NI_element * nel,int ct_start)1296 void GICOR_process_dataset( NI_element *nel , int ct_start )
1297 {
1298 Three_D_View *im3d = A_CONTROLLER , *qq3d ;
1299 GICOR_setup *giset = im3d->giset ;
1300 float *nelar , *dsdar ;
1301 int nvec,nn,vv , vmul , ic ; float thr ;
1302
1303 int verb=0 ;
1304 if( AFNI_yesenv("AFNI_GIC_DEBUG") ) verb = 9 ; /* 07 Apr 2016 */
1305
1306 ENTRY("GICOR_process_dataset") ;
1307
1308 if( nel == NULL || nel->vec_num < 2 ){ /* should never happen */
1309 ERROR_message("AFNI GIC badly formatted dataset from 3dGroupInCorr!") ;
1310 EXRETURN ;
1311 }
1312
1313 ic = AFNI_controller_index(im3d) ; if( ic < 0 || ic > 25 ) EXRETURN ;
1314
1315 if( giset->toplabel != NULL )
1316 PLUTO_set_toplabel( GICOR_plint , giset->toplabel ) ;
1317
1318 nvec = nel->vec_len ; /* how many values in each column transmitted */
1319
1320 if( !IM3D_OPEN(im3d) ||
1321 giset == NULL ||
1322 !giset->ready ){ /* should never happen */
1323
1324 GRPINCORR_LABEL_OFF(im3d) ; SENSITIZE_INSTACORR(im3d,False) ;
1325 if( giset != NULL ) giset->ready = giset->busy = 0 ;
1326 AFNI_misc_CB(im3d->vwid->func->gicor_pb,(XtPointer)im3d,NULL) ;
1327 (void) MCW_popup_message( im3d->vwid->picture ,
1328 " \n"
1329 " ******* AFNI: *********\n"
1330 " 3dGrpInCorr sent data \n"
1331 " but setup isn't ready!\n " ,
1332 MCW_USER_KILL | MCW_TIMER_KILL ) ;
1333 EXRETURN ;
1334 }
1335
1336 /* copy NIML data into dataset */
1337
1338 if( verb > 8 )
1339 INFO_message("AFNI GIC: received %d vectors, length=%d",nel->vec_num,nvec) ;
1340
1341 for( vv=0 ; vv < DSET_NVALS(giset->dset) ; vv++ ){
1342 nelar = (float *)nel->vec[vv] ; /* NIML array */
1343 dsdar = (float *)DSET_ARRAY(giset->dset,vv) ; /* dataset array */
1344
1345 if( verb > 8 ){
1346 float mm,ss ; int nf ;
1347 nf = thd_floatscan( nvec , nelar ) ;
1348 meansigma_float( nvec , nelar , &mm,&ss ) ;
1349 ININFO_message(" vec#%02d nf=%d mean=%g sigma=%g",vv,nf,mm,ss) ;
1350 }
1351
1352 if( giset->ivec == NULL ){ /* all voxels */
1353 nn = MIN( giset->nvox , nvec ) ;
1354 if( verb > 8 ) ININFO_message(" memcpy-ing %d values into dataset",nn) ;
1355 memcpy(dsdar,nelar,sizeof(float)*nn) ;
1356 } else { /* some voxels */
1357 int *ivec=giset->ivec , kk ;
1358 nn = MIN( giset->nivec , nvec ) ;
1359 if( verb > 8 ){
1360 ININFO_message(" copying %d values into dataset",nn) ;
1361 for( kk=0 ; kk < nn ; kk++ ){
1362 ININFO_message(" dsdar[%d] = %g",ivec[kk],nelar[kk]) ;
1363 dsdar[ivec[kk]] = nelar[kk] ;
1364 }
1365 } else {
1366 for( kk=0 ; kk < nn ; kk++ ) dsdar[ivec[kk]] = nelar[kk] ;
1367 }
1368 }
1369 }
1370
1371 /* allow dset to be written out [ZSS Jan 2010] */
1372
1373 giset->dset->dblk->diskptr->allow_directwrite = 1 ;
1374
1375 /* switch to this dataset as overlay */
1376
1377 if( !EQUIV_DSETS(im3d->fim_now,giset->dset) ){
1378 MCW_choose_cbs cbs ; char cmd[32] , *cpt=AFNI_controller_label(im3d) ;
1379 if( verb > 8 )
1380 ININFO_message(" switching controller %c to GIC dataset",cpt[1]) ;
1381 if( verb > 8 )
1382 ININFO_message(" fim_now=%s;%s and giset=%s;%s",
1383 DSET_PREFIX(im3d->fim_now) , DSET_IDCODE_STR(im3d->fim_now) ,
1384 DSET_PREFIX(giset->dset) , DSET_IDCODE_STR(giset->dset) ) ;
1385 cbs.ival = giset->nds ;
1386 AFNI_finalize_dataset_CB( im3d->vwid->view->choose_func_pb ,
1387 (XtPointer)im3d , &cbs ) ;
1388 AFNI_set_fim_index(im3d,0) ;
1389 AFNI_set_thr_index(im3d,1) ;
1390 #if 1
1391 sprintf(cmd,"SET_FUNC_RANGE %c 0.4" , cpt[1]) ; AFNI_driver(cmd) ;
1392 sprintf(cmd,"SET_THRESHNEW %c 0.0" , cpt[1]) ; AFNI_driver(cmd) ;
1393 #endif
1394 } else if( verb > 8 ){
1395 ININFO_message(" not switching: fim_now=%s;%s and giset=%s;%s",
1396 DSET_PREFIX(im3d->fim_now) , DSET_IDCODE_STR(im3d->fim_now) ,
1397 DSET_PREFIX(giset->dset) , DSET_IDCODE_STR(giset->dset) ) ;
1398 }
1399
1400 /* self-threshold and clusterize? */
1401
1402 #undef THBOT
1403 #undef THTOP
1404 #undef THBIG
1405 #define THBIG 1.e+9f
1406 #define THBOT(t) ((thrsign==0 || thrsign==2) ? (-(t)) : (-THBIG))
1407 #define THTOP(t) ((thrsign==0 || thrsign==1) ? (t) : (THBIG))
1408
1409 vmul = giset->vmul ;
1410 thr = get_3Dview_func_thresh(im3d,1) ;
1411 #ifdef GIC_ALLOW_CLUST
1412 if( vmul > 0 && thr > 0.0f ){
1413 MRI_IMAGE *dsim , *tsim , *clim ;
1414 int thrsign=im3d->vinfo->thr_sign , pfun=im3d->vinfo->use_posfunc ;
1415 float thb,tht ;
1416
1417 thb = THBOT(thr) ; tht = THTOP(thr) ;
1418
1419 dsim = DSET_BRICK(giset->dset,im3d->vinfo->fim_index) ;
1420 tsim = DSET_BRICK(giset->dset,im3d->vinfo->thr_index) ;
1421 clim = mri_clusterize( 0.0f , vmul , dsim , thb,tht,tsim , pfun , NULL ) ;
1422 if( clim != NULL ){
1423 float *csar = MRI_FLOAT_PTR(clim) ;
1424 memcpy( MRI_FLOAT_PTR(dsim) , csar , sizeof(float)*clim->nvox ) ;
1425 mri_free(clim) ;
1426 }
1427 }
1428 #endif
1429 DSET_KILL_STATS(giset->dset) ; THD_load_statistics(giset->dset) ;
1430 AFNI_reset_func_range(im3d) ;
1431 flush_3Dview_sort(im3d,"T"); /* ZSS April 27 2012: Reset sorted threshold */
1432
1433 /* redisplay overlay */
1434
1435 if( verb > 8 ) ININFO_message(" redisplay functional overlay") ;
1436
1437 if( called_before[ic] ) AFNI_ignore_pbar_top(1) ; /* 03 Jun 2014 */
1438 IM3D_CLEAR_TMASK(im3d) ; /* Mar 2013 */
1439 IM3D_CLEAR_THRSTAT(im3d) ; /* 12 Jun 2014 */
1440 if( VEDIT_good(im3d->vedset) ) im3d->vedset.flags = 1 ; /* 18 Jun 2014 */
1441 if( MCW_val_bbox(im3d->vwid->view->see_func_bbox) == 0 ){ /* overlay = off */
1442 char cmd[32] , *cpt=AFNI_controller_label(im3d) ;
1443 sprintf(cmd,"SEE_OVERLAY %c.+",cpt[1]) ;
1444 AFNI_driver(cmd) ;
1445 } else { /* overlay = on */
1446 AFNI_redisplay_func(im3d) ;
1447 }
1448 AFNI_set_thr_pval(im3d) ; AFNI_process_drawnotice(im3d) ; called_before[ic]++ ;
1449
1450 for( vv=1 ; vv < MAX_CONTROLLERS ; vv++ ){ /* other controllers need redisplay? */
1451 qq3d = GLOBAL_library.controllers[vv] ;
1452 if( !IM3D_OPEN(qq3d) && qq3d != im3d ) continue ;
1453 if( qq3d->fim_now == giset->dset && MCW_val_bbox(qq3d->vwid->view->see_func_bbox) ){
1454 if( VEDIT_good(qq3d->vedset) ) qq3d->vedset.flags = 1 ; /* 18 Jun 2014 */
1455 AFNI_reset_func_range(qq3d) ; AFNI_redisplay_func(qq3d) ;
1456 }
1457 }
1458 AFNI_ignore_pbar_top(0) ;
1459
1460 if( im3d->vwid->func->options_vedit_av->ival != VEDIT_GRINCORR ){
1461 AV_assign_ival( im3d->vwid->func->options_vedit_av , VEDIT_GRINCORR ) ;
1462 AFNI_vedit_CB( im3d->vwid->func->options_vedit_av , im3d ) ;
1463 }
1464
1465 if( verb > 8 ) ININFO_message(" DONE with this data from 3dGIC") ;
1466
1467 giset->busy = 0 ; /* Not busy waiting anymore [18 Mar 2010] */
1468 GRPINCORR_LABEL_ON(im3d) ; /* 07 Apr 2010 */
1469 EXRETURN ;
1470 }
1471
1472 /*--------------------------------------------------------------------------*/
1473
AFNI_gicor_setapair_xyz(Three_D_View * im3d,float xx,float yy,float zz)1474 void AFNI_gicor_setapair_xyz( Three_D_View *im3d , float xx,float yy,float zz )
1475 {
1476 NI_element *nel ;
1477 char buf[256] ;
1478 GICOR_setup *giset = im3d->giset ;
1479 THD_fvec3 iv,jv; THD_ivec3 kv; int ijk,ii ;
1480
1481 ENTRY("AFNI_gicor_setapair_xyz") ;
1482
1483 if( !IM3D_OPEN(im3d) || giset == NULL || !giset->ready ) EXRETURN ; /* bad */
1484 if( !GICOR_apair_allow_bit(giset) ) EXRETURN ; /* not even a good idea */
1485 if( giset->busy ) EXRETURN ; /* it's busy over there */
1486 if( NI_stream_goodcheck(giset->ns,1) < 1 ) EXRETURN ; /* socket not good */
1487
1488 LOAD_FVEC3( iv , xx,yy,zz ) ;
1489 jv = THD_dicomm_to_3dmm( giset->dset, iv ) ;
1490
1491 if( jv.xyz[0] < giset->dset->daxes->xxmin ||
1492 jv.xyz[0] > giset->dset->daxes->xxmax ||
1493 jv.xyz[1] < giset->dset->daxes->yymin ||
1494 jv.xyz[1] > giset->dset->daxes->yymax ||
1495 jv.xyz[2] < giset->dset->daxes->zzmin ||
1496 jv.xyz[2] > giset->dset->daxes->zzmax ){
1497
1498 ERROR_message("GrpInCorr set Apair point outside dataset box -- ignored") ;
1499 EXRETURN ;
1500 }
1501
1502 kv = THD_3dmm_to_3dind_no_wod( giset->dset, jv ) ;
1503 ijk = DSET_ixyz_to_index( giset->dset, kv.ijk[0],kv.ijk[1],kv.ijk[2] ) ;
1504
1505 if( giset->ivec != NULL ){
1506 ii = bsearch_int( ijk , giset->nvec , giset->ivec ) ;
1507 if( ii < 0 ){
1508 ERROR_message("AFNI: GrpInCorr set Apair point not in mask from 3dGroupInCorr -- ignored") ;
1509 EXRETURN ;
1510 }
1511 }
1512
1513 nel = NI_new_data_element( "SETAPAIR_ijk" , 0 ) ;
1514
1515 sprintf( buf , "%d" , ijk ) ;
1516 NI_set_attribute( nel , "index" , buf ) ;
1517
1518 ii = NI_write_element( giset->ns , nel , NI_TEXT_MODE ) ;
1519 NI_free_element( nel ) ;
1520
1521 GICOR_set_apair_ready_bit(giset) ; /* apair is ready for beeswax */
1522 SENSITIZE_INSTACORR_GROUP(im3d,1) ;
1523 EXRETURN ;
1524 }
1525
1526 /*--------------------------------------------------------------------------*/
1527 /* Called to set the reference voxel,
1528 and to start the calculations over in 3dGroupInCorr.
1529 *//*------------------------------------------------------------------------*/
1530
AFNI_gicor_setref_xyz(Three_D_View * im3d,float xx,float yy,float zz)1531 int AFNI_gicor_setref_xyz( Three_D_View *im3d , float xx,float yy,float zz )
1532 {
1533 NI_element *nel ;
1534 char buf[256] ;
1535 GICOR_setup *giset = im3d->giset ;
1536 THD_fvec3 iv,jv; THD_ivec3 kv; int ijk,ii ;
1537
1538 int verb ;
1539 if( AFNI_yesenv("AFNI_GIC_DEBUG") ) verb = 9 ; /* 07 Apr 2016 */
1540
1541 ENTRY("AFNI_gicor_setref_xyz") ;
1542
1543 if( !IM3D_OPEN(im3d) ||
1544 giset == NULL ||
1545 !giset->ready ){ /* should not happen */
1546
1547 STATUS("something wrong in the setup") ;
1548 if( verb > 8 ) WARNING_message("AFNI GIC: something wrong in the setup") ;
1549
1550 GRPINCORR_LABEL_OFF(im3d) ; SENSITIZE_INSTACORR(im3d,False) ;
1551 if( giset != NULL ) giset->ready = giset->busy = 0 ;
1552 AFNI_misc_CB(im3d->vwid->func->gicor_pb,(XtPointer)im3d,NULL) ;
1553 RETURN(-1) ;
1554 }
1555
1556 /* if socket has gone bad, we're done */
1557
1558 STATUS("check socket to 3dGroupInCorr") ;
1559
1560 if( NI_stream_goodcheck(giset->ns,1) < 1 ){
1561 GRPINCORR_LABEL_OFF(im3d) ; SENSITIZE_INSTACORR(im3d,False) ;
1562 if( giset != NULL ) giset->ready = giset->busy = 0 ;
1563 AFNI_misc_CB(im3d->vwid->func->gicor_pb,(XtPointer)im3d,NULL) ;
1564 if( verb > 8 ) WARNING_message("AFNI GIC: connection to 3dGroupInCorr is broken") ;
1565 RETURN(-1) ;
1566 }
1567
1568 STATUS("check if already busy") ;
1569
1570 if( giset->busy ){ /* Already waiting? [18 Mar 2010] */
1571 #if 0
1572 MCW_flash_widget( 2 , im3d->vwid->func->gicor_label ) ; /* 07 Apr 2010 */
1573 #endif
1574 if( verb > 8 ) WARNING_message("AFNI GIC: already waiting for 3dGroupInCorr results") ;
1575 RETURN(0) ;
1576 }
1577
1578 if( GICOR_apair_allow_bit(giset) && /* Apr 2013 */
1579 !GICOR_apair_ready_bit(giset) &&
1580 !GICOR_apair_mirror_bit(giset) ){
1581
1582 ERROR_message("AFNI GIC: can't set InstaCorr seed before Set Apair") ;
1583 RETURN(0) ;
1584 }
1585
1586 /* find where we are working from, in dataset coordinates */
1587
1588 STATUS("check coordinates") ;
1589
1590 LOAD_FVEC3( iv , xx,yy,zz ) ;
1591 jv = THD_dicomm_to_3dmm( giset->dset, iv ) ;
1592
1593 if( jv.xyz[0] < giset->dset->daxes->xxmin ||
1594 jv.xyz[0] > giset->dset->daxes->xxmax ||
1595 jv.xyz[1] < giset->dset->daxes->yymin ||
1596 jv.xyz[1] > giset->dset->daxes->yymax ||
1597 jv.xyz[2] < giset->dset->daxes->zzmin ||
1598 jv.xyz[2] > giset->dset->daxes->zzmax ){
1599
1600 WARNING_message("AFNI GIC: seed point is outside dataset box") ;
1601 RETURN(-1) ;
1602 }
1603
1604 STATUS("transform coords to index" ) ;
1605
1606 #if 0
1607 kv = THD_3dmm_to_3dind_no_wod( giset->dset, jv ) ;
1608 #else
1609 kv = THD_find_closest_roundtrip(giset->dset,im3d->anat_now,iv) ;
1610 #endif
1611 ijk = DSET_ixyz_to_index( giset->dset, kv.ijk[0],kv.ijk[1],kv.ijk[2] ) ;
1612
1613 if( verb > 8 )
1614 INFO_message("AFNI GIC: AFNI_gicor_setref called --> ijk=%d",ijk) ;
1615
1616 if( giset->ivec != NULL ){
1617 STATUS("search for index in mask") ;
1618 ii = bsearch_int( ijk , giset->nvec , giset->ivec ) ;
1619 if( ii < 0 ){
1620 WARNING_message("AFNI GIC: seed point not in mask from 3dGroupInCorr") ;
1621 RETURN(-1) ;
1622 }
1623 }
1624
1625 /** Apr 2013: before sending the SETREF_ijk element,
1626 send the mirror point, if relevant to our lives **/
1627
1628 if( GICOR_apair_allow_bit(giset) && GICOR_apair_mirror_bit(giset) )
1629 AFNI_gicor_setapair_xyz( im3d , -xx, yy, zz ) ;
1630
1631 /* NOW send ijk node index to 3dGroupInCorr,
1632 which starts the glorious chain of calculations */
1633
1634 STATUS("create NIML element to send info") ;
1635
1636 nel = NI_new_data_element( "SETREF_ijk" , 0 ) ;
1637
1638 sprintf( buf , "%d" , ijk ) ;
1639 NI_set_attribute( nel , "index" , buf ) ;
1640
1641 sprintf( buf , "%g" , giset->seedrad ) ;
1642 NI_set_attribute( nel , "seedrad" , buf ) ;
1643
1644 #ifdef GIC_ALLOW_TTEST
1645 sprintf( buf , "%d" , giset->ttest_opcode ) ;
1646 NI_set_attribute( nel , "ttest_opcode" , buf ) ;
1647 #endif
1648
1649 STATUS("send NIML element") ;
1650
1651 if( verb > 8 ){
1652 INFO_message("AFNI GIC: sending command to 3dGroupInCorr") ;
1653 NI_sleep(1) ;
1654 }
1655
1656 ii = NI_write_element( giset->ns , nel , NI_TEXT_MODE ) ;
1657 NI_free_element( nel ) ;
1658 if( ii <= 0 ){
1659 ERROR_message("AFNI GIC: 3dGroupInCorr connection has failed :(") ;
1660 RETURN(-1) ;
1661 }
1662
1663 STATUS("This is my only chance at building a disreputable past.") ;
1664 /* 12 Apr 2010: save this seed location ZSS */
1665 LOAD_FVEC3( iv , xx,yy,zz ) ; /* reload to be safe, use of anat_now */
1666 im3d->vinfo->xi_icor = xx ; /* should do the trick. */
1667 im3d->vinfo->yj_icor = yy ;
1668 im3d->vinfo->zk_icor = zz ;
1669 jv = THD_dicomm_to_3dmm( im3d->anat_now , iv ) ;
1670 kv = THD_3dmm_to_3dind ( im3d->anat_now , jv ) ;
1671 UNLOAD_IVEC3( kv , im3d->vinfo->i1_icor ,
1672 im3d->vinfo->j2_icor , im3d->vinfo->k3_icor ) ;
1673 if( verb > 8 )
1674 INFO_message("AFNI GIC: called from xyz = %.2f %.2f %.2f (DICOM), "
1675 "seed radius %g",
1676 xx, yy, zz, giset->seedrad) ;
1677
1678 STATUS("mark that we're busy for now") ;
1679
1680 giset->busy = 1 ; /* Mark that we're busy right now [18 Mar 2010] */
1681 GRPINCORR_LABEL_WORKING(im3d) ; /* 07 Apr 2010 */
1682 RETURN(0) ;
1683 }
1684
1685 /****************************************************************************
1686 Set up the Group InstaCorr interface to the user
1687 *****************************************************************************/
1688
1689 static char g_helpstring[] =
1690 "Purpose: control AFNI Group InstaCorr operations\n"
1691 "\n"
1692 "Author -- RW Cox -- Dec 2009\n"
1693 ;
1694
1695 static char *topts[3] = { "pooled" , "unpooled" , "paired" } ;
1696
1697 static char * GICOR_main( PLUGIN_interface * ) ;
1698
GICOR_init(char * lab)1699 PLUGIN_interface * GICOR_init( char *lab )
1700 {
1701 int ntops , vv ;
1702 PLUGIN_interface *plint ; /* will be the output of this routine */
1703 char sk[32] ;
1704 GICOR_setup *giset = A_CONTROLLER->giset ;
1705
1706 ENTRY("GICOR_init") ;
1707
1708 if( lab == NULL ) lab = "\0" ;
1709
1710 if( !IM3D_OPEN(A_CONTROLLER) || giset == NULL || !giset->ready )
1711 RETURN(NULL) ;
1712
1713 /*---------------- set titles and call point ----------------*/
1714
1715 sprintf(sk,"%sSetup Group InstaCorr",lab) ;
1716 plint = PLUTO_new_interface( "GrpInCorr" ,
1717 sk ,
1718 g_helpstring ,
1719 PLUGIN_CALL_VIA_MENU ,
1720 (char *(*)())GICOR_main ) ;
1721 GICOR_plint = plint ;
1722
1723 PLUTO_set_runlabels( plint , "Set+Keep" , "Set+Quit" ) ;
1724
1725 /*--------- make interface lines -----------*/
1726
1727 PLUTO_add_option ( plint , "Parameters" , "Params" , TRUE ) ;
1728 PLUTO_add_number ( plint , "SeedRadius" , 0,16,0,(int)rint(giset->seedrad), TRUE ) ;
1729
1730 #ifdef GIC_ALLOW_TTEST
1731 ntops = (giset->ndset_A == giset->ndset_B) ? 3 : 2 ;
1732 PLUTO_add_string ( plint , "t-test" , ntops , topts , giset->ttest_opcode ) ;
1733 #endif
1734
1735 #ifdef GIC_ALLOW_CLUST
1736 PLUTO_add_option ( plint , "Clustering" , "Cluster" , TRUE ) ;
1737 PLUTO_add_number ( plint , "Voxels Min" , 0,9999,0 , 0,TRUE ) ;
1738 #endif
1739
1740 if( giset->toplabel != NULL )
1741 PLUTO_set_toplabel( GICOR_plint , giset->toplabel ) ;
1742
1743 if( giset->num_stat_available > 0 ){
1744 for( vv=0 ; vv < GICOR_MAX_NSTAT ; vv++ ){
1745 sprintf(sk,"Stat#%d",vv+1) ;
1746 PLUTO_add_option( plint , sk , sk , (vv==0) ? TRUE : FALSE ) ;
1747 PLUTO_add_string( plint , "What" ,
1748 giset->num_stat_available ,
1749 giset->lab_stat_available ,
1750 (vv < giset->num_stat_available) ? vv : 0 ) ;
1751 }
1752 for( vv=giset->num_stat_available ; vv < GICOR_MAX_NSTAT ; vv++ ){
1753 XtSetSensitive( plint->wid->opwid[vv+GICOR_BASE_STAT]->toggle, False ) ;
1754 }
1755 }
1756
1757 RETURN(plint) ;
1758 }
1759
1760 /*---------------------------------------------------------------------------*/
1761
GICOR_refit_stat_menus(void)1762 static void GICOR_refit_stat_menus(void)
1763 {
1764 GICOR_setup *giset = A_CONTROLLER->giset ;
1765 PLUGIN_interface *plint = GICOR_plint ;
1766 PLUGIN_option_widgets *opwid ;
1767 MCW_arrowval *av ;
1768 int vv , zz ;
1769
1770 ENTRY("GICOR_refit_stat_menus") ;
1771
1772 if( plint == NULL || plint->option_count <= GICOR_BASE_STAT ) EXRETURN ;
1773
1774 for( vv=0 ; vv < GICOR_MAX_NSTAT ; vv++ ){
1775 opwid = plint->wid->opwid[GICOR_BASE_STAT+vv] ;
1776 av = (MCW_arrowval *)opwid->chooser[0] ;
1777 zz = av->ival ; if( zz >= giset->num_stat_available ) zz = 0 ;
1778 refit_MCW_optmenu( av ,
1779 0 , giset->num_stat_available - 1 , zz , 0 ,
1780 MCW_av_substring_CB , giset->lab_stat_available ) ;
1781 if( vv >= giset->num_stat_available ){
1782 XmToggleButtonSetState( opwid->toggle, False,True ) ;
1783 XtSetSensitive ( opwid->toggle, False ) ;
1784 } else {
1785 XtSetSensitive ( opwid->toggle, True ) ;
1786 }
1787 }
1788
1789 EXRETURN ;
1790 }
1791
1792 /***************************************************************************
1793 Main routine for this plugin (will be called from AFNI).
1794 If the return string is not NULL, some error transpired, and
1795 AFNI will popup the return string in a message box.
1796 ****************************************************************************/
1797
GICOR_main(PLUGIN_interface * plint)1798 static char * GICOR_main( PLUGIN_interface *plint )
1799 {
1800 Three_D_View *im3d = plint->im3d ;
1801 float srad=0.0f ; int topcod=0 , vmul=0 ; char *tch ;
1802 GICOR_setup *giset ;
1803
1804 ENTRY("GICOR_main") ;
1805
1806 if( !IM3D_OPEN(im3d) ||
1807 im3d->giset == NULL ||
1808 !im3d->giset->ready ){ /* should not happen */
1809
1810 GRPINCORR_LABEL_OFF(im3d) ; SENSITIZE_INSTACORR(im3d,False) ;
1811 if( im3d->giset != NULL ) im3d->giset->ready = im3d->giset->busy = 0 ;
1812 XtUnmapWidget(plint->wid->shell) ;
1813 RETURN(" ************ AFNI: ************ \n 3dGroupInCorr is no longer enabled!? \n ") ;
1814 }
1815
1816 giset = im3d->giset ;
1817
1818 /* if socket has gone bad, we're done */
1819
1820 if( NI_stream_goodcheck(giset->ns,1) < 1 ){
1821 GRPINCORR_LABEL_OFF(im3d) ; SENSITIZE_INSTACORR(im3d,False) ;
1822 if( giset != NULL ) giset->ready = giset->busy = 0 ;
1823 XtUnmapWidget(plint->wid->shell) ;
1824 RETURN(" ************ AFNI: ************ \n 3dGroupInCorr is no longer connected! \n ") ;
1825 }
1826
1827 if( giset->toplabel != NULL )
1828 PLUTO_set_toplabel( GICOR_plint , giset->toplabel ) ;
1829
1830 PLUTO_next_option(plint) ;
1831 srad = PLUTO_get_number(plint) ;
1832 #ifdef GIC_ALLOW_TTEST
1833 tch = PLUTO_get_string(plint) ;
1834 topcod = PLUTO_string_index( tch , 3 , topts ) ;
1835 #endif
1836
1837 #ifdef GIC_ALLOW_CLUST
1838 PLUTO_next_option(plint) ;
1839 vmul = (int)PLUTO_get_number(plint) ;
1840 #endif
1841
1842 /* do something with these changes */
1843
1844 giset->seedrad = srad ;
1845 giset->vmul = vmul ;
1846 #ifdef GIC_ALLOW_TTEST
1847 giset->ttest_opcode = topcod ;
1848 #endif
1849
1850 RETURN(NULL) ;
1851 }
1852