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 "mrilib.h"
00008 
00009 /*----------------
00010   Another quickie.
00011 ------------------*/
00012 
00013 int main( int argc , char * argv[] )
00014 {
00015    int narg , nvox , ii,jj,kk,vv , mcount , iv , mc , ndset,ndval , i,j,k ;
00016    THD_3dim_dataset *mask_dset=NULL , **input_dset=NULL ;
00017    float mask_bot=666.0 , mask_top=-666.0 ;
00018    byte * mmm   = NULL ;
00019    char * oname = NULL , * obuf , * otemp ;
00020    FILE * ofile ;
00021    MRI_IMAGE * flim ;
00022    float * flar ;
00023    int noijk=0 , yes_xyz=0 ;
00024    int yes_index=0 ;                    /*-- 09 May 2003 [rickr] --*/
00025    byte *cmask=NULL ; int ncmask=0 ;
00026    int verb=1 ;
00027 
00028 #define BOXLEN   7
00029 #define BOX_XYZ  1
00030 #define BOX_DIC  2
00031 #define BOX_NEU  3
00032 #define BOX_IJK  4
00033    int box_num=0 ; float *box_dat=NULL ;   /* 09 May 2003 - RWCox */
00034    int nx,ny,nz,nxy,nxyz ;
00035 
00036    if( argc < 2 || strcmp(argv[1],"-help") == 0 ){
00037       printf("Usage: 3dmaskdump [options] dataset dataset ...\n"
00038              "Writes to an ASCII file values from the input datasets\n"
00039              "which satisfy the mask criteria given in the options.\n"
00040              "If no options are given, then all voxels are included.\n"
00041              "This might result in a GIGANTIC output file.\n"
00042              "Options:\n"
00043              "  -mask mset   Means to use the dataset 'mset' as a mask:\n"
00044              "                 Only voxels with nonzero values in 'mset'\n"
00045              "                 will be printed from 'dataset'.  Note\n"
00046              "                 that the mask dataset and the input dataset\n"
00047              "                 must have the same number of voxels.\n"
00048              "  -mrange a b  Means to further restrict the voxels from\n"
00049              "                 'mset' so that only those mask values\n"
00050              "                 between 'a' and 'b' (inclusive) will\n"
00051              "                 be used.  If this option is not given,\n"
00052              "                 all nonzero values from 'mset' are used.\n"
00053              "                 Note that if a voxel is zero in 'mset', then\n"
00054              "                 it won't be included, even if a < 0 < b.\n"
00055                              /*-- 09 May 2003: add -index option [rickr] */
00056              "  -index       Means to write out the dataset index values.\n"
00057              "  -noijk       Means not to write out the i,j,k values.\n"
00058              "  -xyz         Means to write the x,y,z coordinates from\n"
00059              "                 the 1st input dataset at the start of each\n"
00060              "                 output line.  These coordinates are in\n"
00061              "                 the 'RAI' order.\n"
00062              "  -o fname     Means to write output to file 'fname'.\n"
00063              "                 [default = stdout, which you won't like]\n"
00064              "\n"
00065              "  -cmask 'opts' Means to execute the options enclosed in single\n"
00066              "                  quotes as a 3dcalc-like program, and produce\n"
00067              "                  produce a mask from the resulting 3D brick.\n"
00068              "       Examples:\n"
00069              "        -cmask '-a fred+orig[7] -b zork+orig[3] -expr step(a-b)'\n"
00070              "                  produces a mask that is nonzero only where\n"
00071              "                  the 7th sub-brick of fred+orig is larger than\n"
00072              "                  the 3rd sub-brick of zork+orig.\n"
00073              "        -cmask '-a fred+orig -expr 1-bool(k-7)'\n"
00074              "                  produces a mask that is nonzero only in the\n"
00075              "                  7th slice (k=7); combined with -mask, you\n"
00076              "                  could use this to extract just selected voxels\n"
00077              "                  from particular slice(s).\n"
00078              "       Notes: * You can use both -mask and -cmask in the same\n"
00079              "                  run - in this case, only voxels present in\n"
00080              "                  both masks will be dumped.\n"
00081              "              * Only single sub-brick calculations can be\n"
00082              "                  used in the 3dcalc-like calculations -\n"
00083              "                  if you input a multi-brick dataset here,\n"
00084              "                  without using a sub-brick index, then only\n"
00085              "                  its 0th sub-brick will be used.\n"
00086              "              * Do not use quotes inside the 'opts' string!\n"
00087              "\n"
00088              "  -xbox x y z   Means to put a 'mask' down at the dataset (not DICOM)\n"
00089              "                  coordinates of 'x y z' mm.  By default, this box is\n"
00090              "                  1 voxel wide in each direction.  You can specify\n"
00091              "                  instead a range of coordinates using a colon ':'\n"
00092              "                  after the coordinates; for example:\n"
00093              "                    -xbox 22:27 31:33 44\n"
00094              "                  means a box from (x,y,z)=(22,31,44) to (27,33,44).\n"
00095              "\n"
00096              "  -dbox x y z   Means the same as -xbox, but the coordinates are in\n"
00097              "                  DICOM order (+x=Left, +y=Posterior, +z=Superior).\n"
00098              "                  These coordinates correspond to those you'd enter\n"
00099              "                  into the 'Jump to (xyz)' control in AFNI, and to\n"
00100              "                  those output by default from 3dclust.\n"
00101              "  -nbox x y z   Means the same as -xbot, but the coordinates are in\n"
00102              "                  'neuroscience' order (+x=Right, +y=Anterior, +z=Superior)\n"
00103              "\n"
00104              "  -ibox i j k   Means to put a 'mask' down at the voxel indexes\n"
00105              "                  given by 'i j k'.  By default, this picks out\n"
00106              "                  just 1 voxel.  Again, you can use a ':' to specify\n"
00107              "                  a range (now in voxels) of locations.\n"
00108              "       Notes: * Boxes are cumulative; that is, if you specify more\n"
00109              "                  than 1 box, you'll get more than one region.\n"
00110              "              * If a -mask and/or -cmask option is used, then\n"
00111              "                  the intersection of the boxes with these masks\n"
00112              "                  determines which voxels are output; that is,\n"
00113              "                  a voxel must be inside some box AND inside the\n"
00114              "                  mask in order to be selected for output.\n"
00115              "              * If boxes select more than 1 voxel, the output lines\n"
00116              "                  are NOT necessarily in the order of the options on\n"
00117              "                  the command line.\n"
00118              "              * Coordinates (for -xbox, -dbox, and -nbox) are relative\n"
00119              "                  to the first dataset on the command line.\n"
00120              "\n"
00121              "  -quiet        Means not to print progress messages to stderr.\n"
00122              "\n"
00123              "Inputs after the last option are datasets whose values you\n"
00124              "want to be dumped out.  These datasets (and the mask) can\n"
00125              "use the sub-brick selection mechanism (described in the\n"
00126              "output of '3dcalc -help') to choose which values you get.\n"
00127              "\n"
00128              "Each selected voxel gets one line of output:\n"
00129              "  i j k val val val ....\n"
00130              "where (i,j,k) = 3D index of voxel in the dataset arrays,\n"
00131              "and val = the actual voxel value.  Note that if you want\n"
00132              "the mask value to be output, you have to include that\n"
00133              "dataset in the dataset input list again, after you use\n"
00134              "it in the '-mask' option.\n"
00135              "\n"
00136              "N.B.: This program doesn't work with complex-valued datasets!\n"
00137             ) ;
00138 
00139       printf("\n" MASTER_SHORTHELP_STRING ) ;
00140 
00141       exit(0) ;
00142    }
00143 
00144    /*-- 20 Apr 2001: addto the arglist, if user wants to [RWCox] --*/
00145 
00146    mainENTRY("3dmaskdump main"); machdep() ; PRINT_VERSION("3dmaskdump") ;
00147 
00148    { int new_argc ; char ** new_argv ;
00149      addto_args( argc , argv , &new_argc , &new_argv ) ;
00150      if( new_argv != NULL ){ argc = new_argc ; argv = new_argv ; }
00151    }
00152 
00153    AFNI_logger("3dmaskdump",argc,argv) ;
00154 
00155    /* scan argument list */
00156 
00157    narg = 1 ;
00158    while( narg < argc && argv[narg][0] == '-' ){
00159 
00160       if( strcmp(argv[narg],"-quiet") == 0 ){                /* 09 May 2003 - RWC */
00161         verb = 0 ; narg++ ; continue ;
00162       }
00163 
00164       if( strcmp(argv[narg]+2,"box") == 0 ){                 /* 09 May 2003 - RWC */
00165         float xbot,xtop , ybot,ytop , zbot,ztop , btyp ;
00166         int nn ;
00167         char code = *(argv[narg]+1) ;   /* should be 'x', 'd' , 'n', or 'i' */
00168         switch( code ){
00169           case 'x': btyp = BOX_XYZ ; break ;
00170           case 'd': btyp = BOX_DIC ; break ;
00171           case 'n': btyp = BOX_NEU ; break ;
00172           case 'i': btyp = BOX_IJK ; break ;
00173           default:
00174             fprintf(stderr,"** Unknown 'box' option %s\n",argv[narg]); exit(1);
00175         }
00176         if( narg+3 >= argc ){
00177           fprintf(stderr,"** need 3 arguments after %s\n",argv[narg]); exit(1);
00178         }
00179         nn = sscanf( argv[narg+1] , "%f:%f" , &xbot , &xtop ) ;
00180         if( nn < 1 ){
00181           fprintf(stderr,"** Can't decode %s after %s\n",argv[narg+1],argv[narg]); exit(1);
00182         } else if( nn == 1 ){
00183           xtop=xbot ;
00184         }
00185         nn = sscanf( argv[narg+2] , "%f:%f" , &ybot , &ytop ) ;
00186         if( nn < 1 ){
00187           fprintf(stderr,"** Can't decode %s after %s\n",argv[narg+2],argv[narg]); exit(1);
00188         } else if( nn == 1 ){
00189           ytop=ybot ;
00190         }
00191         nn = sscanf( argv[narg+3] , "%f:%f" , &zbot , &ztop ) ;
00192         if( nn < 1 ){
00193           fprintf(stderr,"** Can't decode %s after %s\n",argv[narg+3],argv[narg]); exit(1);
00194         } else if( nn == 1 ){
00195           ztop=zbot ;
00196         }
00197         box_dat = (float *) realloc( box_dat , sizeof(float)*BOXLEN*(box_num+1) ) ;
00198         box_dat[0+BOXLEN*box_num] = xbot ;
00199         box_dat[1+BOXLEN*box_num] = xtop ;
00200         box_dat[2+BOXLEN*box_num] = ybot ;
00201         box_dat[3+BOXLEN*box_num] = ytop ;
00202         box_dat[4+BOXLEN*box_num] = zbot ;
00203         box_dat[5+BOXLEN*box_num] = ztop ;
00204         box_dat[6+BOXLEN*box_num] = btyp ; box_num++ ;
00205         narg += 4 ; continue ;
00206       }
00207 
00208       /*-- 09 May 2003: option to output index value (for Mike B) [rickr] --*/
00209       if( strcmp(argv[narg],"-index") == 0 ){
00210          yes_index = 1 ;
00211          narg++ ; continue ;
00212       }
00213 
00214       if( strcmp(argv[narg],"-noijk") == 0 ){
00215          noijk = 1 ;
00216          narg++ ; continue ;
00217       }
00218 
00219       if( strcmp(argv[narg],"-xyz") == 0 ){   /* 23 Mar 2003 */
00220          yes_xyz = 1 ;
00221          narg++ ; continue ;
00222       }
00223 
00224       if( strcmp(argv[narg],"-cmask") == 0 ){  /* 16 Mar 2000 */
00225          if( narg+1 >= argc ){
00226             fprintf(stderr,"** -cmask option requires a following argument!\n");
00227             exit(1) ;
00228          }
00229          cmask = EDT_calcmask( argv[++narg] , &ncmask ) ;
00230          if( cmask == NULL ){
00231             fprintf(stderr,"** Can't compute -cmask!\n"); exit(1);
00232          }
00233          narg++ ; continue ;
00234       }
00235 
00236       if( strncmp(argv[narg],"-mask",5) == 0 ){
00237          if( mask_dset != NULL ){
00238             fprintf(stderr,"** Cannot have two -mask options!\n") ; exit(1) ;
00239          }
00240          if( narg+1 >= argc ){
00241             fprintf(stderr,"** -mask option requires a following argument!\n");
00242             exit(1) ;
00243          }
00244          mask_dset = THD_open_dataset( argv[++narg] ) ;
00245          if( mask_dset == NULL ){
00246             fprintf(stderr,"** Cannot open mask dataset!\n") ; exit(1) ;
00247          }
00248          if( DSET_BRICK_TYPE(mask_dset,0) == MRI_complex ){
00249             fprintf(stderr,"** Cannot deal with complex-valued mask dataset!\n");
00250             exit(1) ;
00251          }
00252          narg++ ; continue ;
00253       }
00254 
00255       if( strncmp(argv[narg],"-mrange",5) == 0 ){
00256          if( narg+2 >= argc ){
00257             fprintf(stderr,"** -mrange option requires 2 following arguments!\n");
00258              exit(1) ;
00259          }
00260          mask_bot = strtod( argv[++narg] , NULL ) ;
00261          mask_top = strtod( argv[++narg] , NULL ) ;
00262          if( mask_top < mask_top ){
00263             fprintf(stderr,"** -mrange inputs are illegal!\n") ; exit(1) ;
00264          }
00265          narg++ ; continue ;
00266       }
00267 
00268       if( strcmp(argv[narg],"-o") == 0 ){
00269          if( narg+1 >= argc ){
00270             fprintf(stderr,"** -o needs an argument after it!\n"); exit(1);
00271          }
00272          oname = argv[++narg] ;
00273          if( ! THD_filename_ok(oname) ){
00274             fprintf(stderr,"** name after -o is illegal!\n"); exit(1) ;
00275          }
00276          if( THD_is_file(oname) ){
00277             fprintf(stderr,"** file %s already exists!\n",oname); exit(1);
00278          }
00279          narg++ ; continue ;
00280       }
00281 
00282       fprintf(stderr,"** Unknown option: %s\n",argv[narg]) ; exit(1) ;
00283    }
00284 
00285    /* should have at least one more argument */
00286 
00287    ndset = argc - narg ;
00288    if( ndset <= 0 ){
00289       fprintf(stderr,"** No input dataset!?\n") ; exit(1) ;
00290    }
00291 
00292    /* open all input datasets */
00293 
00294    input_dset = (THD_3dim_dataset **) malloc( sizeof(THD_3dim_dataset *)*ndset ) ;
00295    for( ndval=ii=0 ; ii < ndset ; ii++ ){
00296       input_dset[ii] = THD_open_dataset( argv[narg+ii] ) ;
00297       if( input_dset[ii] == NULL ){
00298          fprintf(stderr,"** Can't open dataset %s!\n",argv[narg+ii]) ;
00299          exit(1) ;
00300       }
00301 
00302       if( DSET_BRICK_TYPE(input_dset[ii],0) == MRI_complex ){
00303          fprintf(stderr,"** Cannot deal with complex-valued input dataset!\n");
00304          exit(1) ;
00305       }
00306 
00307       if( ii == 0 ){
00308          nvox = DSET_NVOX(input_dset[0]) ;
00309          nx   = DSET_NX(input_dset[0]) ;
00310          ny   = DSET_NY(input_dset[0]) ; nxy  = nx*ny ;
00311          nz   = DSET_NZ(input_dset[0]) ; nxyz = nxy*nz ;
00312       } else {
00313          if( DSET_NVOX(input_dset[ii]) != nvox ){
00314             fprintf(stderr,"** Dataset %s does not match in size!\n",argv[narg+ii]);
00315             exit(1) ;
00316          }
00317       }
00318 
00319       ndval += DSET_NVALS(input_dset[ii]) ;
00320    }
00321 
00322    /* make a byte mask from mask dataset */
00323 
00324    if( mask_dset == NULL ){
00325       mmm = NULL ;
00326       if( verb ) fprintf(stderr,"++ %d voxels in the entire dataset (no mask)\n",nvox) ;
00327    } else {
00328       if( DSET_NVOX(mask_dset) != nvox ){
00329          fprintf(stderr,"** Input and mask datasets are not same dimensions!\n");
00330          exit(1) ;
00331       }
00332       mmm = THD_makemask( mask_dset , 0 , mask_bot,mask_top ) ;
00333       mcount = THD_countmask( nvox , mmm ) ;
00334       if( mcount <= 0 ){
00335          fprintf(stderr,"** No voxels in the mask!\n") ; exit(1) ;
00336       }
00337       if( verb ) fprintf(stderr,"++ %d voxels in the mask\n",mcount) ;
00338       DSET_delete(mask_dset) ;
00339    }
00340 
00341    /* 16 Mar 2000: deal with the cmask */
00342 
00343    if( cmask != NULL ){
00344       if( ncmask != nvox ){
00345          fprintf(stderr,"** Input and cmask datasets are not same dimensions!\n");
00346          exit(1) ;
00347       }
00348       if( mmm != NULL ){
00349          for( ii=0 ; ii < nvox ; ii++ )
00350             mmm[ii] = (mmm[ii] && cmask[ii]) ;
00351          free(cmask) ;
00352          mcount = THD_countmask( nvox , mmm ) ;
00353          if( mcount <= 0 ){
00354             fprintf(stderr,"** No voxels in the mask+cmask!\n") ; exit(1) ;
00355          }
00356          if( verb ) fprintf(stderr,"++ %d voxels in the mask+cmask\n",mcount) ;
00357       } else {
00358          mmm = cmask ;
00359          mcount = THD_countmask( nvox , mmm ) ;
00360          if( mcount <= 0 ){
00361             fprintf(stderr,"** No voxels in the cmask!\n") ; exit(1) ;
00362          }
00363          if( verb ) fprintf(stderr,"++ %d voxels in the cmask\n",mcount) ;
00364       }
00365    }
00366 
00367    /* 09 May 2003: make a mask corresponding to the boxen - RWC */
00368 
00369    if( box_num > 0 ){
00370      byte *bmask = calloc(1,nvox) ;
00371      int bb, ibot,itop, jbot,jtop, kbot,ktop , btyp , ii,jj,kk ;
00372      float   xbot,xtop, ybot,ytop, zbot,ztop ;
00373      THD_fvec3 dv,xv ;
00374      THD_3dim_dataset *dset = input_dset[0] ;
00375      float xmin=dset->daxes->xxmin , xmax=dset->daxes->xxmax ;
00376      float ymin=dset->daxes->yymin , ymax=dset->daxes->yymax ;
00377      float zmin=dset->daxes->zzmin , zmax=dset->daxes->zzmax ;
00378 
00379      for( bb=0 ; bb < box_num ; bb++ ){
00380        xbot = box_dat[0+BOXLEN*bb]; xtop = box_dat[1+BOXLEN*bb];
00381        ybot = box_dat[2+BOXLEN*bb]; ytop = box_dat[3+BOXLEN*bb];
00382        zbot = box_dat[4+BOXLEN*bb]; ztop = box_dat[5+BOXLEN*bb];
00383        btyp = box_dat[6+BOXLEN*bb];
00384 
00385        if( btyp != BOX_IJK ){                      /* convert coords to indexes */
00386 
00387          if( btyp == BOX_NEU ){                    /* coords from Neuroscience to DICOM */
00388            xbot = -xbot; xtop = -xtop; ybot = -ybot; ytop = -ytop; btyp = BOX_DIC;
00389          }
00390          if( btyp == BOX_DIC ){                    /* coords from DICOM to dataset */
00391            LOAD_FVEC3(dv,xbot,ybot,zbot) ;
00392            xv = THD_dicomm_to_3dmm( dset , dv ) ;
00393            UNLOAD_FVEC3(xv,xbot,ybot,zbot) ;
00394            LOAD_FVEC3(dv,xtop,ytop,ztop) ;
00395            xv = THD_dicomm_to_3dmm( dset , dv ) ;
00396            UNLOAD_FVEC3(xv,xtop,ytop,ztop) ;
00397          }
00398          if( xbot < xmin && xtop < xmin ) continue ; /* skip box if outside dataset */
00399          if( xbot > xmax && xtop > xmax ) continue ;
00400          if( ybot < ymin && ytop < ymin ) continue ;
00401          if( ybot > ymax && ytop > ymax ) continue ;
00402          if( zbot < zmin && ztop < zmin ) continue ;
00403          if( zbot > zmax && ztop > zmax ) continue ;
00404          LOAD_FVEC3(dv,xbot,ybot,zbot) ;
00405          xv = THD_3dmm_to_3dfind( dset , dv ) ;   /* coords from dataset to index */
00406          UNLOAD_FVEC3(xv,xbot,ybot,zbot) ;
00407          LOAD_FVEC3(dv,xtop,ytop,ztop) ;
00408          xv = THD_3dmm_to_3dfind( dset , dv ) ;
00409          UNLOAD_FVEC3(xv,xtop,ytop,ztop) ;
00410        }
00411        ibot = rint(xbot) ; jbot = rint(ybot) ; kbot = rint(zbot) ;  /* round */
00412        itop = rint(xtop) ; jtop = rint(ytop) ; ktop = rint(ztop) ;
00413        if( ibot > itop ){ btyp = ibot; ibot = itop; itop = btyp; }  /* flip? */
00414        if( jbot > jtop ){ btyp = jbot; jbot = jtop; jtop = btyp; }
00415        if( kbot > ktop ){ btyp = kbot; kbot = ktop; ktop = btyp; }
00416 
00417        /* skip box if outside dataset */
00418        if ( itop < 0 || ibot >= nx-1 ) continue;
00419        if ( jtop < 0 || jbot >= ny-1 ) continue;
00420        if ( ktop < 0 || kbot >= nz-1 ) continue;
00421 
00422        /* constrain values to dataset dimensions */
00423        if ( ibot < 0 ) ibot = 0;  if ( itop >= nx-1 ) itop = nx-1;
00424        if ( jbot < 0 ) jbot = 0;  if ( jtop >= ny-1 ) jtop = ny-1;
00425        if ( kbot < 0 ) kbot = 0;  if ( ktop >= nz-1 ) ktop = nz-1;
00426 
00427        for( kk=kbot ; kk <= ktop ; kk++ )
00428         for( jj=jbot ; jj <= jtop ; jj++ )
00429          for( ii=ibot ; ii <= itop ; ii++ )
00430            bmask[ii+jj*nx+kk*nxy] = 1 ;
00431      }
00432 
00433      mcount = THD_countmask( nvox , bmask ) ;
00434      if( verb ) fprintf(stderr,"++ %d voxels in the boxes\n",mcount) ;
00435      if( mcount == 0 ){
00436        fprintf(stderr,"** Can't continue with no voxels insides boxes!\n"); exit(1);
00437      }
00438      if( mmm != NULL ){
00439        for( ii=0 ; ii < nvox ; ii++ )
00440           mmm[ii] = (mmm[ii] && bmask[ii]) ;
00441        free(bmask) ;
00442        mcount = THD_countmask( nvox , mmm ) ;
00443        if( mcount <= 0 ){
00444           fprintf(stderr,"** No voxels in the mask+boxes!\n") ; exit(1) ;
00445        }
00446        if( verb ) fprintf(stderr,"++ %d voxels in the mask+boxes\n",mcount) ;
00447      } else {
00448        mmm = bmask ;
00449      }
00450    } /* end of box mask */
00451 
00452    /* read in input dataset bricks */
00453 
00454    for( ii=0 ; ii < ndset ; ii++ ){
00455       DSET_load(input_dset[ii]) ;
00456       if( !DSET_LOADED(input_dset[ii]) ){
00457          fprintf(stderr,"*** Can't load dataset %s\n",argv[narg+ii]); exit(1);
00458       }
00459    }
00460 
00461    /* open the output file */
00462 
00463    if( oname != NULL ){
00464       ofile = fopen( oname , "w" ) ;
00465       if( ofile == NULL ){
00466          fprintf(stderr,"*** Can't open output file %s\n",oname); exit(1);
00467       }
00468    } else {
00469       ofile = stdout ;
00470    }
00471 
00472    /* output string buffers */
00473 
00474    /*-- 09 May 2003: add room for the index (9 -> 10)    [rickr] --*/
00475    obuf = (char *) malloc( sizeof(char) * (ndval+10) * 16 ) ;
00476 
00477    /* loop over voxels */
00478 
00479    for( ii=0 ; ii < nvox ; ii++ ){
00480       if( mmm != NULL && mmm[ii] == 0 ) continue ;  /* skip this voxel */
00481 
00482       obuf[0] = '\0' ;
00483 
00484       /*-- 09 May 2003: optionally display index    [rickr] --*/
00485       if ( yes_index ){
00486          otemp = MV_format_fval((float)ii);
00487          strcat(obuf,otemp); strcat(obuf," ");
00488       }
00489 
00490       if( ! noijk ){
00491          i = DSET_index_to_ix( input_dset[0] , ii ) ;  /* voxel indexes */
00492          j = DSET_index_to_jy( input_dset[0] , ii ) ;
00493          k = DSET_index_to_kz( input_dset[0] , ii ) ;
00494 
00495          otemp = MV_format_fval((float)i); strcat(obuf,otemp); strcat(obuf," ");
00496          otemp = MV_format_fval((float)j); strcat(obuf,otemp); strcat(obuf," ");
00497          otemp = MV_format_fval((float)k); strcat(obuf,otemp); strcat(obuf," ");
00498       }
00499 
00500       if( yes_xyz ){                                  /* 23 Mar 2003 */
00501         THD_ivec3 ind ; THD_fvec3 vec ;
00502         i = DSET_index_to_ix( input_dset[0] , ii ) ;
00503         j = DSET_index_to_jy( input_dset[0] , ii ) ;
00504         k = DSET_index_to_kz( input_dset[0] , ii ) ;
00505         LOAD_IVEC3(ind,i,j,k) ;
00506         vec = THD_3dind_to_3dmm ( input_dset[0] , ind ) ;
00507         vec = THD_3dmm_to_dicomm( input_dset[0] , vec ) ;
00508         otemp = MV_format_fval(vec.xyz[0]); strcat(obuf,otemp); strcat(obuf," ");
00509         otemp = MV_format_fval(vec.xyz[1]); strcat(obuf,otemp); strcat(obuf," ");
00510         otemp = MV_format_fval(vec.xyz[2]); strcat(obuf,otemp); strcat(obuf," ");
00511       }
00512 
00513       for( kk=0 ; kk < ndset ; kk++ ){
00514          flim = THD_extract_series( ii , input_dset[kk] , 0 ) ;
00515          flar = MRI_FLOAT_PTR(flim) ;
00516          for( vv=0 ; vv < flim->nx ; vv++ ){
00517             otemp = MV_format_fval(flar[vv]) ;
00518             strcat(obuf,otemp) ; strcat(obuf," ") ;
00519          }
00520          mri_free(flim) ;
00521       }
00522 
00523       jj = strlen(obuf) ; obuf[jj-1] = '\0' ; /* kill last blank */
00524 
00525       fprintf(ofile,"%s\n",obuf) ;
00526    }
00527 
00528    fflush(ofile) ; fsync(fileno(ofile)) ; fclose(ofile) ;
00529    exit(0) ;
00530 }
00531