00001 /*****************************************************************************
00002    Major portions of this software are copyrighted by the Medical College
00003    of Wisconsin, 1994-2000, and are released under the Gnu General Public
00004    License, Version 2.  See the file README.Copyright for details.
00005 ******************************************************************************/
00006    
00007 #include "afni.h"
00008 #include "afni_plugin.h"
00009 
00010 #ifndef ALLOW_PLUGINS
00011 #  error "Plugins not properly set up -- see machdep.h"
00012 #endif
00013 
00014 /***********************************************************************
00015   Plugin to extract 3D brick data 
00016 ************************************************************************/
00017 typedef struct 
00018         {
00019                   int nxx;                      /* number of voxels in the x direction */
00020                   int nyy;                      /* number of voxels in the y direction */
00021                   int nzz;                      /* number of voxels in the z direction */
00022                   char *dsetname; /* prefix of data set analyzed */
00023                   int errcode;          /* error code number returned from function */
00024                   int out;                      /* flag for writing to a file */
00025                   int format;
00026                   int iloc;
00027                   int xloc;
00028                   int yloc;
00029                   int zloc;
00030                   int fimonly;          /* set to 1 if no threshold is available */
00031                   int DoInt;
00032                   int DoThres;
00033                   int DoInd;
00034                   int intind; 
00035                   int   thrind;
00036                   int Nsub;
00037                   int isanat;
00038                   int isfunc;
00039                   float mini, maxi, minth, maxth;
00040                   char * strout;
00041                   FILE * outfile;
00042                   FILE * outlogfile;
00043                   char outname[PLUGIN_MAX_STRING_RANGE]; /* output data file name */
00044         }extract_data;
00045 
00046 
00047 static char helpstring[] =
00048   "                   3Ddump98 Plugin\n"
00049   "This plugin is used to write to an ascii file the data present in AFNI bricks.\n"
00050   "You can apply intenstity or threshold masks to the voxel data that is being extracted.\n\n" 
00051   "Plugin Inputs:\n\n"
00052   "   1- Dataset :\n"
00053   "      3D brick  -> 3D AFNI brick of the type :\n"
00054   "                    fim, fith, fico, fbuc, etc...\n"
00055   "                    spgr, epan.\n\n"
00056   "   2- SubBrick info : (Optional) \n"
00057   "      Intensity -> Index of the subbrick to be used\n"
00058   "                   as an intensity subbrick.\n"
00059   "      Threshold -> Index of the subbrick to be used \n"
00060   "                   as a threshold subbrick.\n" 
00061   "   While the subbrick indices are obvious when dealing with most bricks\n"
00062   "   You might need to specified them for bricks of the type bucket\n\n"
00063   "   3- Intensity Mask : (optional) \n"
00064   "      Minimum   -> Minimum boundary for intensity value (inclusive)\n"
00065   "      Maximum   -> Maximum boundary for intensity value (inclusive)\n"
00066   "   Data from voxels with intensity between Minimum and Maximum \n"
00067   "   are written to 'Filename'.\n\n"
00068   "   4- Threshold Mask : (optional) \n"
00069   "      Minimum   -> Minimum boundary for threshold value (inclusive)\n"
00070   "      Maximum   -> Maximum boundary for threshold value (inclusive)\n"
00071   "   Data from voxels with threshold value between Minimum and Maximum \n"
00072   "   are written to 'Filename'.\n\n"
00073   "   5- Output : \n"
00074   "      Filename  -> Name of ascii output file. \n"
00075   "                   If no name is specified, the default is\n"
00076   "                   the prefix of the inbut brick with the \n"
00077   "                   extension '.3Ddump' appended at the end.\n"
00078   "                   A LOG file, 'Filename.log' is also written to disk.\n"
00079   "                   The log file contains all the parameters settings used\n"
00080   "                   for generating 'Filename'.\n"
00081   "                   The format of 'Filename' is as follows :\n"
00082   "                   1- Voxel Index (VI) : Each voxel in an AFNI brick has a unique index.\n"
00083   "                                         Indices map directly to XYZ coordinates.\n"
00084   "                                         See AFNI plugin documentations for more info.\n"
00085   "                   2..4- Voxel coordinates (X Y Z) : Those are the voxel slice coordinates.\n"
00086   "                                                     You can see these coordinates in the upper left side\n"
00087   "                                                     of the AFNI window. To do so, you must first switch the\n"
00088   "                                                     voxel coordinate units from mm to slice coordinates.\n"
00089   "                                                     Define Datamode -> Misc -> Voxel Coords ?\n"
00090   "                                                     PS: The coords that show up in the graph window\n"
00091   "                                                     could be different from those in the upper left side \n"
00092   "                                                     of AFNI's main window.\n"
00093   "                   5..n- Subbrick values (Sb1 Sb2 ... Sbn) : Voxel values at each subbrick.\n\n"
00094   "If you have/find questions/comments/bugs about the plugin, \n"
00095   "send me an E-mail: ziad@image.bien.mu.edu\n\n"
00096   "                    Ziad Saad   Nov. 9 97, latest update Aug. 26 99.\n\n"
00097 ;
00098 
00099 /* Significant update Aug. 26 99 */ 
00100 /* The indexing into Storear has been swapped (columns became rows and rows columns.
00101 That's because each subbrick was not stored in a vector on N elements, rather in N vectors
00102 of 1 element each. That made the allocation process very slow, especially now that Bob has
00103 malloc and calloc going through his own macros. */
00104 
00105 /*-------- strings for output format  and some definitions -----------*/
00106 
00107 static char * yn_strings[] = { "n" , "y" }; 
00108 
00109 #define NUM_YN_STRINGS (sizeof(yn_strings)/sizeof(char *))
00110 
00111 #define YUP  1
00112 #define NOPE 0
00113 
00114 #define ERROR_FILEWRITE         2
00115 #define ERROR_OPTIONS           3
00116 
00117 /*---------- prototypes for internal routines ----------*/
00118 static int filexists (char *);
00119 
00120 static char * DUMP_main( PLUGIN_interface * ) ;
00121 
00122 static int Dumpit( extract_data* , THD_3dim_dataset * ) ;
00123 
00124 static void write_ud (extract_data*);
00125 
00126 static char **allocate2D (int rows,int cols,int element_size);
00127 
00128 static void free2D(char **a,int rows);
00129 
00130 static int equal_strings (char *s1,char *s2);
00131 
00132 
00133 /***********************************************************************
00134    Set up the interface to the user
00135 ************************************************************************/
00136 
00137 DEFINE_PLUGIN_PROTOTYPE
00138 
00139 PLUGIN_interface * PLUGIN_init( int ncall )
00140 {
00141    PLUGIN_interface * plint ;
00142 
00143    if( ncall > 0 ) return NULL ;  /* only one interface */
00144 
00145    /*-- set titles and call point --*/
00146 
00147    plint = PLUTO_new_interface( "3D Dump98" , "Ascii dump of 3D Dataset" , helpstring ,
00148                                  PLUGIN_CALL_VIA_MENU , DUMP_main  ) ;
00149 
00150    PLUTO_set_runlabels( plint , "Dump+Keep" , "Dump+Close" ) ;  /* 04 Nov 2003 */
00151 
00152    /*-- first line of input: Dataset --*/
00153 
00154    PLUTO_add_option( plint , "Dataset" , "Dataset" , TRUE ) ;
00155    PLUTO_add_dataset(plint , "3D brick" ,
00156                                     ANAT_ALL_MASK  ,  FUNC_ALL_MASK ,
00157                                     SESSION_ALL_MASK |
00158                                     DIMEN_3D_MASK    | BRICK_ALLREAL_MASK ) ;
00159    
00160         /*-- second line of input: intensity and threshold brick indices --*/
00161    
00162    
00163         PLUTO_add_option( plint ,
00164                      "SubBrik info" ,  /* label at left of input line */
00165                      "Index" ,  /* tag to return to plugin */
00166                      FALSE       /* is this mandatory? */
00167                    ) ;
00168         PLUTO_add_number( plint ,
00169                     "Intensity" ,  /* label next to chooser */
00170                     1 ,         /* smallest possible value */
00171                     10000 ,        /* largest possible value (inactivated for now)*/
00172                     0 ,         /* decimal shift (none in this case) */
00173                     0 ,         /* default value */
00174                     FALSE       /* allow user to edit value? */
00175                   ) ;
00176                                 
00177         PLUTO_add_number( plint ,
00178                     "Threshold" ,  /* label next to chooser */
00179                     1 ,         /* smallest possible value */
00180                     10000 ,        /* largest possible value (inactivated for now)*/
00181                     0 ,         /* decimal shift (none in this case) */
00182                     2 ,         /* default value */
00183                     FALSE       /* allow user to edit value? */
00184                   ) ;   
00185         /*-- Third line of input: intensity mask --*/
00186    
00187    PLUTO_add_option( plint ,
00188                      "Intensity Mask" ,  /* label at left of input line */
00189                      "Intensity" ,  /* tag to return to plugin */
00190                      FALSE       /* is this mandatory? */
00191                    ) ;
00192    
00193    PLUTO_add_number( plint ,
00194                     "Minimum" ,  /* label next to chooser */
00195                     -100000 ,         /* smallest possible value */
00196                     100000 ,        /* largest possible value (inactivated for now)*/
00197                     0 ,         /* decimal shift (none in this case) */
00198                     0 ,         /* default value */
00199                     TRUE       /* allow user to edit value? */
00200                   ) ;
00201                   
00202         PLUTO_add_number( plint ,
00203                     "Maximum" ,  /* label next to chooser */
00204                     -10000 ,         /* smallest possible value */
00205                     10000 ,        /* largest possible value (inactivated for now)*/
00206                     0 ,         /* decimal shift (none in this case) */
00207                     0 ,         /* default value */
00208                     TRUE       /* allow user to edit value? */
00209                   ) ;
00210    
00211    /*-- Fourth line of input: threshold mask --*/
00212    
00213    PLUTO_add_option( plint ,
00214                      "Threshold Mask" ,  /* label at left of input line */
00215                      "Threshold" ,  /* tag to return to plugin */
00216                      FALSE       /* is this mandatory? */
00217                    ) ;
00218    
00219    PLUTO_add_number( plint ,
00220                     "Minimum" ,  /* label next to chooser */
00221                     -10000 ,         /* smallest possible value */
00222                     10000 ,        /* largest possible value (inactivated for now)*/
00223                     0 ,         /* decimal shift (none in this case) */
00224                     0.5 ,         /* default value */
00225                     TRUE       /* allow user to edit value? */
00226                   ) ;
00227                   
00228         PLUTO_add_number( plint ,
00229                     "Maximum" ,  /* label next to chooser */
00230                     -10000 ,         /* smallest possible value */
00231                     10000 ,        /* largest possible value (inactivated for now)*/
00232                     0 ,         /* decimal shift (none in this case) */
00233                     1 ,         /* default value */
00234                     TRUE       /* allow user to edit value? */
00235                   ) ;
00236    
00237    /*---------- 5th line: output stuff ----------*/
00238 
00239    PLUTO_add_option( plint ,
00240                      "Output" ,  /* label at left of input line */
00241                      "Output" ,  /* tag to return to plugin */
00242                      TRUE        /* is this mandatory? */
00243                    ) ;
00244 
00245    PLUTO_add_string( plint ,
00246                      "Filename" ,  /* label next to textfield */
00247                      0,NULL ,    /* no fixed strings to choose among */
00248                      19          /* 19 spaces for typing in value */
00249                    ) ;
00250                   
00251    return plint ;
00252 }
00253 
00254 /***************************************************************************
00255   Main routine for this plugin (will be called from AFNI).
00256 ****************************************************************************/
00257 
00258 static char * DUMP_main( PLUGIN_interface * plint )
00259 {
00260    extract_data uda,*ud;
00261    MCW_idcode * idc ;
00262    THD_3dim_dataset * xset , * yset ;
00263    char * tag ;
00264    int demean ,ndmp,nprf;
00265    char *str, *nprfxstr, *mssg;
00266    float minx , maxx , minthr , maxthr ;
00267         
00268         str = (char *) calloc (PLUGIN_MAX_STRING_RANGE+10,sizeof(char));
00269         nprfxstr         = (char *) calloc (PLUGIN_MAX_STRING_RANGE+20,sizeof(char));
00270         /* Do not allocate more space for mssg, because AFNI would choke on it*/
00271         mssg = (char *) calloc (PLUGIN_MAX_STRING_RANGE,sizeof(char));
00272         
00273         if (str == NULL || nprfxstr == NULL || mssg == NULL ) 
00274                                                                           return "********************\n"
00275                                                                                                 "Could not Allocate\n"
00276                                                                                                 "a teeni weeni bit of\n"
00277                                                                                                 "Memory ! \n"
00278                                                                                                 "********************\n";
00279 
00280         ud = &uda;              /* ud now points to an allocated space */
00281    
00282    /*--------------------------------------------------------------------*/
00283    /*----- Check inputs from AFNI to see if they are reasonable-ish -----*/
00284         tag = PLUTO_get_optiontag(plint) ;
00285    
00286    if (tag == NULL)
00287         {
00288                 return "************************\n"
00289              "Bad 1st line option \n"
00290              "************************"  ;
00291         }       
00292         
00293    
00294    idc  = PLUTO_get_idcode(plint) ;     /* get 1st dataset item */
00295    xset = PLUTO_find_dset(idc) ;                   /* get ptr to dataset */
00296    if( xset == NULL )
00297       return "**********************\n"
00298              "Cannot find Dataset #1\n"
00299              "**********************"  ;
00300         
00301         ud->dsetname = DSET_FILECODE (xset);
00302         ud->Nsub = DSET_NVALS (xset);
00303         
00304         /*--------- loop over ramining options ---------*/
00305         ud->DoInt = NOPE;
00306         ud->DoThres = NOPE;
00307         ud->DoInd = NOPE;
00308         
00309         ud->intind = 1; /* The first subbrick is numbered 1 here, it makes more sense to the user.*/
00310         ud->thrind = 2; 
00311         
00312         do
00313                 {
00314                         tag = PLUTO_get_optiontag(plint) ;
00315                         if (tag == NULL) break;
00316                         
00317                         if (equal_strings (tag, "Index") == 1)
00318                                 {
00319                                         ud->DoInd = YUP;
00320                                         ud->intind = PLUTO_get_number(plint) ; 
00321                                         ud->thrind = PLUTO_get_number(plint) ;
00322                                         continue; 
00323                                 }
00324                         
00325                         if (equal_strings (tag, "Intensity") == 1)
00326                                 {
00327                                         ud->DoInt = YUP;
00328                                         ud->mini = PLUTO_get_number(plint) ; 
00329                                         ud->maxi = PLUTO_get_number(plint) ;
00330                                         continue;
00331                                 }
00332                         
00333                         if (equal_strings (tag, "Threshold") == 1)
00334                                 {
00335                                         ud->DoThres = YUP;
00336                                         ud->minth = PLUTO_get_number(plint) ; 
00337                                         ud->maxth = PLUTO_get_number(plint) ;
00338                                         continue;
00339                                 }
00340                         
00341                         if (equal_strings (tag, "Output") == 1)
00342                                         {
00343                                                 ud->strout = PLUTO_get_string(plint) ; 
00344                                                 continue;
00345                                         }
00346                         
00347                         
00348                 } while (1);
00349                 
00350         if ( ISFUNC(xset) ) ud->isfunc = YUP;
00351                 else ud->isfunc = NOPE;
00352         
00353         if ( ISANAT(xset) ) ud->isanat = YUP;
00354                 else ud->isanat = NOPE;
00355         
00356         if (xset->func_type == FUNC_FIM_TYPE)
00357                         ud->fimonly = 1;
00358                 else
00359                         ud->fimonly = 0;
00360         
00361         if (ud->isanat && (ud->DoThres== YUP || ud->DoInd == YUP))
00362                 {
00363                         return "*************************************\n"
00364                 "Can't use threshold or index options \n" 
00365                 "for fim or ANAT type bricks !\n"
00366                 "*************************************"  ;
00367                 }
00368         
00369 
00370         if (ud->DoInd == YUP && (ud->fimonly == YUP && ud->isfunc == YUP))
00371                 {
00372                         return "*******************************\n"
00373                 "Can't specify Indices for fim\n" 
00374                 "type bricks, they only have one !\n"
00375                 "*******************************"  ;
00376                 }
00377                         
00378 
00379         /* Check for plausibility of input parameters */
00380         if ((ud->DoInt && (ud->maxi < ud->mini)) || (ud->DoThres && (ud->maxth < ud->minth)))
00381                 {
00382                 return "**********************\n"
00383              "Something's wrong with\n" 
00384              "min and max  values.\n"
00385              "**********************"  ;
00386                 }
00387         if (ud->DoInd && (ud->intind > ud->Nsub || ud->thrind > ud->Nsub))
00388                 {
00389                 
00390                 return "**********************\n"
00391              "One or both of the indices\n" 
00392              "is larger than the maximum\n"
00393                                  "number of sub-bricks\n"
00394              "**********************"  ;
00395                 }
00396 
00397         /*------------------------------------------------------*/
00398    /*----- Open the output file for writing operation -----*/
00399    
00400         /*if strout is null or of length 0, use the a default name */
00401         if (ud->strout == NULL)
00402         nprf = 0;
00403    else
00404         nprf = 1;
00405         
00406         
00407    if (nprf == 1 && (int)strlen(ud->strout) == 0)
00408         nprf = 0;
00409                 
00410         
00411    if (nprf == 0)
00412         {
00413                         sprintf (nprfxstr,"%s.3Ddump",DSET_PREFIX(xset));
00414                         ud->strout = nprfxstr;
00415         }
00416 
00417    
00418         sprintf (str,"%s.log",ud->strout);      
00419                 
00420    if ((filexists(ud->strout) == 1) || (filexists(str) == 1))
00421         {
00422                 return "**************************************\n"
00423                        "Output file(s) exists, can't overwrite\n"
00424                        "**************************************\n";
00425         }
00426         
00427         ud->outfile = fopen (ud->strout,"w");
00428         ud->outlogfile = fopen (str,"w");
00429         
00430         if ((ud->outfile == NULL) || (ud->outlogfile == NULL))
00431                 {
00432                         return "*****************************************\n"
00433                        "Could not open Output file(s) for writing\n"
00434                        "Check permissions.\n"
00435                        "*****************************************\n";
00436                 }
00437         
00438    /*------------------------------------------------------*/
00439    /*---------- At this point, the inputs are OK ----------*/
00440 
00441    /*-- do the actual work --*/
00442    ndmp = Dumpit( ud ,  xset ) ;
00443 
00444    if( ndmp < -1.0 )
00445       {
00446                         switch (ndmp)
00447                                 {
00448                                         case -2:
00449                                                 return  "********************************************\n"
00450                                                 "Fatal Error: Could not allocate memory\n"
00451                                                                         "(this is a message brought to you by Dumpit)\n"
00452                                                 "********************************************\n"  ;
00453                                                 break;
00454                                         default:
00455                                                 return  "*********************************\n"
00456                                                 "Error while dumping data\n"
00457                                                 "*********************************"  ;
00458 
00459                                                 break;
00460                         
00461                                 }
00462                 }
00463 
00464    /*-- put the output to the screen --*/
00465 
00466    /* That output message was too long. AFNI was crashing, must ask Bob to allow more verbose messages */ 
00467          
00468    /*sprintf(mssg , "            Dataset %s was dumped.\n"
00469                  "%d voxels (%5f %% of total) met the boundary conditions.\n"
00470                   , DSET_FILECODE(xset) , ndmp , (float)ndmp/(float)(ud->nxx * ud->nyy * ud-> nzz)*100.0) ;*/
00471                                                 
00472    /* That's shorter */
00473         sprintf(mssg , "%d voxels (%5f %% of total) were dumped.\n" 
00474                          , ndmp , (float)ndmp/(float)(ud->nxx * ud->nyy * ud-> nzz)*100.0) ;
00475 
00476         PLUTO_popup_message( plint , mssg ) ;
00477         
00478         
00479         
00480         fclose (ud->outfile);
00481         fclose (ud->outlogfile);
00482         free (nprfxstr);        
00483         free (str);
00484         free (mssg); 
00485          
00486    return NULL ;  /* null string returned means all was OK */
00487 }
00488 
00489 
00490 static int Dumpit( extract_data* ud, THD_3dim_dataset * xset)
00491 {
00492    void  *  xar  ;
00493    void  * thar ;
00494    float * fxar  ;
00495         float ** Storear;
00496    int ii , jj, nxyz , fxar_new = 0  ,xpos,ypos,zpos , ndmp, pass;
00497 
00498         ndmp = -1;
00499         
00500    /*-- get datasets sizes --*/
00501 
00502         ud->nxx = xset->daxes->nxx;
00503         ud->nyy = xset->daxes->nyy;
00504         ud->nzz = xset->daxes->nzz;
00505         
00506    nxyz = xset->daxes->nxx * xset->daxes->nyy * xset->daxes->nzz ;
00507 
00508    /* Allocate space for data storage */
00509         fxar = (float *) malloc( sizeof(float) * nxyz ) ; fxar_new = 1 ;
00510    Storear = (float **) allocate2D (ud->Nsub ,nxyz ,sizeof(float)); /* changed from :  allocate2D (nxyz  ,ud->Nsub,sizeof(float)) */
00511         
00512         if (fxar == NULL || Storear == NULL)
00513                 {
00514                         return -2;
00515                 }
00516         
00517         /*-- load the first dataset into memory --*/
00518 
00519                 DSET_load( xset ) ;
00520                 
00521         /* Store Bricks in 2D array */  
00522    for (ii = 0;ii < ud->Nsub; ++ii)
00523                 {
00524                 xar   = DSET_ARRAY(xset,ii) ;        /* get the array */
00525                 EDIT_coerce_scale_type (nxyz,DSET_BRICK_FACTOR(xset,ii),
00526                                                         DSET_BRICK_TYPE(xset,ii), xar,  MRI_float,fxar ) ;   
00527                         
00528                         /* Store the iith sub-brick in Storear array */
00529                         for (jj = 0; jj < nxyz; ++jj)
00530                                         Storear[ii][jj] = fxar[jj]; /* changed from : Storear[jj][ii] */
00531                 }
00532                 
00533    DSET_unload( xset ) ;  /* don't need this in memory anymore */
00534         
00535         /* Dump the input info data to the log file */
00536         write_ud (ud);
00537         
00538 
00539    /* Now dump the data that meets the threshold */
00540 
00541    if( 1 ){
00542       for( ii=0 ; ii < nxyz ; ii++ ){
00543         pass = YUP;
00544         if (pass && ud->DoInt)
00545                 {
00546                         if (Storear[ud->intind-1][ii] < ud->mini || Storear[ud->intind-1][ii] > ud->maxi) /* changed both from : Storear[ii][ud->intind-1]*/
00547                                 pass = NOPE;
00548                 }
00549         
00550         if (pass && ud->DoThres)   
00551                 {
00552                         if (Storear[ud->thrind-1][ii] < ud->minth || Storear[ud->thrind-1][ii] > ud->maxth)/* changed both from : Storear[ii][ud->intind-1]*/
00553                                 pass = NOPE;
00554                 }
00555         
00556         if (pass)
00557         {
00558                 zpos = (int)ii / (int)(xset->daxes->nxx * xset->daxes->nyy);
00559                                 ypos = (int)(ii - zpos * xset->daxes->nxx * xset->daxes->nyy) / xset->daxes->nxx;
00560                                 xpos = ii - ( ypos * xset->daxes->nxx ) - ( zpos * xset->daxes->nxx * xset->daxes->nyy ) ;
00561                 
00562                 fprintf (ud->outfile,"%d\t%d\t%d\t%d\t",ii,xpos,ypos,zpos);
00563                         
00564                                 for (jj = 0; jj < ud->Nsub; ++jj)
00565                                         fprintf (ud->outfile," %f\t",Storear[jj][ii]); /* changed from: Storear[ii][jj] */
00566                                         
00567                                 fprintf (ud->outfile,"\n");
00568                                                 
00569                                 ++ndmp;
00570          }
00571         }
00572                 
00573                 /* increment by one because I want it to start at 0 */
00574         ++ndmp;
00575    }
00576         
00577         
00578         fprintf (ud->outlogfile,"\n%d voxel points met the threshold conditions\n",ndmp);
00579         
00580    /*-- free up arrays --*/
00581 
00582         
00583    if( fxar_new ) 
00584                 {
00585                 free(fxar) ;
00586                 }
00587          else 
00588                 {
00589                         DSET_unload(xset) ;
00590                 }
00591         
00592         
00593         free2D ((char **)Storear,ud->Nsub);     /* changed from free2D ((char **)Storear,nxyz); */
00594 
00595         
00596    return ndmp ;
00597 }
00598 
00599 
00600 /* ************************************************************ */ 
00601 /* function to check for file existence       */
00602 /* ************************************************************ */ 
00603         
00604 static int filexists (char *f_name)
00605 {/*filexists*/
00606         FILE *outfile;
00607         
00608         outfile = fopen (f_name,"r");
00609         if (outfile == NULL)
00610                 return (0);
00611         else 
00612                 fclose (outfile);
00613                 return (1);
00614                 
00615 }/*filexists*/
00616 
00617 /* ************************************************************ */ 
00618 /* function to create log file       */
00619 /* ************************************************************ */ 
00620 
00621 void write_ud (extract_data* ud)
00622         {
00623                 fprintf (ud->outlogfile,"\n\nUser Data Values \n");
00624                 fprintf (ud->outlogfile,"Input data set file name= %s\n",ud->dsetname);
00625                 fprintf (ud->outlogfile,"Is the file fim only type ? = %d\n",ud->fimonly);
00626                 fprintf (ud->outlogfile,"Number of Subbricks : %d\n",ud->Nsub);
00627                 fprintf (ud->outlogfile,"output file name = %s\n",ud->strout);
00628                 fprintf (ud->outlogfile,"Number of voxels in X direction = %d\n",ud->nxx);
00629                 fprintf (ud->outlogfile,"Number of voxels in Y direction = %d\n",ud->nyy);
00630                 fprintf (ud->outlogfile,"Number of voxels in Z direction = %d\n",ud->nzz);
00631                 fprintf (ud->outlogfile,"Do intensity mask ? = %d\n",ud->DoInt);
00632                 fprintf (ud->outlogfile,"Minimum intensity = %f\n",ud->mini);
00633                 fprintf (ud->outlogfile,"Maximum intensity = %f\n",ud->maxi);
00634                 fprintf (ud->outlogfile,"Do threshold mask ? = %d\n",ud->DoThres);
00635                 fprintf (ud->outlogfile,"Minimum threshold = %f\n",ud->minth);
00636                 fprintf (ud->outlogfile,"Maximum threshold = %f\n",ud->maxth);
00637                 fprintf (ud->outlogfile,"Select Intensity and Threshold indices = %d\n",ud->DoInd);
00638                 fprintf (ud->outlogfile,"Intensity index = %d\n",ud->intind);
00639                 fprintf (ud->outlogfile,"Threshold index = %d\n",ud->thrind);
00640                 fprintf (ud->outlogfile,"\nThe format for the output file is the following:\n");
00641            fprintf (ud->outlogfile,"VI\tX\tY\tZ\tSb1\tSb2\t... Sbn\n\n");
00642                 
00643                 
00644                 return;
00645         }
00646 
00647 /* ************************************************************ */ 
00648 /* function to allocate 2D arrays       */
00649 /* ************************************************************ */ 
00650 
00651 static char **allocate2D (int rows,int cols,int element_size)
00652 
00653 {
00654     int i, j;
00655     char **A;
00656 
00657 /* try to allocate the request */
00658     switch(element_size) {
00659         case sizeof(short): {    /* integer matrix */
00660             short **int_matrix;
00661             int_matrix = (short **)calloc(rows,sizeof(short *));
00662             if(!int_matrix) {
00663                 printf("\nError making pointers in %dx%d int matrix\n"
00664                             ,rows,cols);
00665                 return(NULL);
00666                 /*                exit(1);*/
00667             }
00668             for(i = 0 ; i < rows ; i++) {
00669                 int_matrix[i] = (short *)calloc(cols,sizeof(short));
00670                 if(!int_matrix[i]) {
00671                     printf("\nError making row %d in %dx%d int matrix\n"
00672                             ,i,rows,cols);
00673                     for(j=0;j<=i;j++) {
00674                       if(int_matrix[j])
00675                         free(int_matrix[j]);
00676                     }
00677                     free(int_matrix);
00678                     return(NULL);
00679                     /*                    exit(1);*/
00680                 }
00681             }
00682             A = (char **)int_matrix;
00683             break;
00684         }
00685         case sizeof(float): {    /* float matrix */
00686             float **float_matrix;
00687             float_matrix = (float **)calloc(rows,sizeof(float *));
00688             if(!float_matrix) {
00689                 printf("\nError making pointers in %dx%d float matrix\n"
00690                             ,rows,cols);
00691                 return(NULL);
00692                 /*                exit(1);*/
00693             }
00694             for(i = 0 ; i < rows ; i++) {
00695                 float_matrix[i] = (float *)calloc(cols,sizeof(float));
00696                 if(!float_matrix[i]) {
00697                     printf("\nError making row %d in %dx%d float matrix\n"
00698                             ,i,rows,cols);
00699                     for(j=0;j<=i;j++) {
00700                       if(float_matrix[j])
00701                         free(float_matrix[j]);
00702                     }
00703                     free(float_matrix);
00704                     return(NULL);
00705                     /*                    exit(1);*/
00706                 }
00707             }
00708             A = (char **)float_matrix;
00709             break;
00710         }
00711         case sizeof(double): {   /* double matrix */
00712             double **double_matrix;
00713             double_matrix = (double **)calloc(rows,sizeof(double *));
00714             if(!double_matrix) {
00715                 printf("\nError making pointers in %dx%d double matrix\n"
00716                             ,rows,cols);
00717                 return(NULL);
00718                 /*                exit(1);*/
00719             }
00720             for(i = 0 ; i < rows ; i++) {
00721                 double_matrix[i] = (double *)calloc(cols,sizeof(double));
00722                 if(!double_matrix[i]) {
00723                     printf("\nError making row %d in %dx%d double matrix\n"
00724                             ,i,rows,cols);
00725                     for(j=0;j<=i;j++) {
00726                       if(double_matrix[j])
00727                         free(double_matrix[j]);
00728                     }
00729                     free(double_matrix);
00730                     return(NULL);
00731 
00732                     /*                    exit(1);*/
00733                 }
00734             }
00735             A = (char **)double_matrix;
00736             break;
00737         }
00738         default:
00739             printf("\nERROR in matrix_allocate: unsupported type\n");
00740             return(NULL);
00741             /*            exit(1);*/
00742     }
00743     return(A);
00744 }
00745 
00746 /* ************************************************************ */ 
00747 /* function to free 2D arrays       */
00748 /* ************************************************************ */ 
00749 static void free2D(char **a,int rows)
00750     
00751 {
00752     int i;
00753     
00754 /* free each row of data */
00755     for(i = 0 ; i < rows ; i++) free(a[i]);
00756 
00757 /* free each row pointer */
00758     free((char *)a);
00759     a = NULL;           /* set to null for error */
00760     
00761         return;
00762 }
00763 
00764 /* ************************************************************ */ 
00765 /* function to Check if strings are equal       */
00766 /* ************************************************************ */ 
00767 
00768 
00769 int equal_strings (char *s1,char *s2)
00770 
00771  {
00772    int i=0;
00773    
00774    if (s1 == NULL && s2 == NULL) return (-2);
00775    
00776    if ((s1 == NULL && s2 != NULL) || (s1 != NULL && s2 == NULL)) return (-1);
00777    
00778    while (s1[i] == s2[i] 
00779                         && s1[i] != '\0' && s2[i] != '\0') ++i;
00780                         
00781         if (s1[i] == '\0' && s2[i] == '\0') return (1);
00782          else return (0);
00783  
00784  }