Skip to content

AFNI/NIfTI Server

Sections
Personal tools
You are here: Home » AFNI » Documentation

Doxygen Source Code Documentation


Main Page   Alphabetical List   Data Structures   File List   Data Fields   Globals   Search  

plug_stavg.c File Reference

#include "afni.h"

Go to the source code of this file.


Defines

#define TEST_VOXEL   6177
#define TEST_TIME   0
#define RMB_DEBUG   0
#define _STAVG_NUM_METHODS   (sizeof(method_strings)/sizeof(char *))
#define _STAVG_METH_MEAN   0
#define _STAVG_METH_SIGMA   1
#define FREEUP(x)   if((x) != NULL){free((x)); (x)=NULL;}
#define FREE_WORKSPACE

Functions

char * STAVG_main (PLUGIN_interface *)
float ** avg_epochs (THD_3dim_dataset *dset, float *ref, int user_maxlength, int no1, int meth, PLUGIN_interface *plint)
MRI_IMARRdset_to_mri (THD_3dim_dataset *dset)
DEFINE_PLUGIN_PROTOTYPE PLUGIN_interface * PLUGIN_init (int ncall)

Variables

char helpstring []
char * yes_no_strings [] = { "No" , "Yes" }
char * method_strings [] = { "Mean" , "Sigma" }
PLUGIN_interface * global_plint = NULL
int M_maxlength

Define Documentation

#define _STAVG_METH_MEAN   0
 

Definition at line 44 of file plug_stavg.c.

Referenced by avg_epochs(), and PLUGIN_init().

#define _STAVG_METH_SIGMA   1
 

Definition at line 45 of file plug_stavg.c.

Referenced by avg_epochs().

#define _STAVG_NUM_METHODS   (sizeof(method_strings)/sizeof(char *))
 

Definition at line 43 of file plug_stavg.c.

Referenced by PLUGIN_init(), and STAVG_main().

#define FREE_WORKSPACE
 

Value:

do{ FREEUP(bptr) ; FREEUP(sptr) ; FREEUP(fptr) ;  \
      FREEUP(fout) ; \
      FREEUP(fxar) ; \
    } while(0)

Definition at line 205 of file plug_stavg.c.

Referenced by STAVG_main().

#define FREEUP      if((x) != NULL){free((x)); (x)=NULL;}
 

Definition at line 203 of file plug_stavg.c.

#define RMB_DEBUG   0
 

Definition at line 15 of file plug_stavg.c.

#define TEST_TIME   0
 

Definition at line 14 of file plug_stavg.c.

#define TEST_VOXEL   6177
 

Definition at line 13 of file plug_stavg.c.


Function Documentation

float ** avg_epochs THD_3dim_dataset   dset,
float *    ref,
int    user_maxlength,
int    no1,
int    meth,
PLUGIN_interface *    plint
 

Definition at line 479 of file plug_stavg.c.

References _STAVG_METH_MEAN, _STAVG_METH_SIGMA, THD_3dim_dataset::daxes, DSET_BRICK_FACTOR, DSET_BRICK_TYPE, DSET_load, DSET_NUM_TIMES, dset_to_mri(), DSET_unload, free, FREE_IMARR, IMAGE_IN_IMARR, M_maxlength, malloc, MRI_FLOAT_PTR, mri_to_float(), THD_dataxes::nxx, THD_dataxes::nyy, nz, THD_dataxes::nzz, PLUTO_popup_meter(), PLUTO_set_meter(), and ref.

Referenced by STAVG_main().

00483 {
00484 
00485    int     numepochs, lastindex;
00486    int     nvox, numims, nx, ny, nz;
00487    int     kim, ne, te, tinc, nia;
00488    int     ii, kk;
00489    int     maxlength, minlength;
00490    int     datum;
00491    float   fac;       /* scaling factor for current subbrick */
00492    double  scaledval; /* temp var for scaled value to be used repeatedly */
00493    float ** fxar;
00494    float ** outar;    /* output averaged time-series */
00495    float * sumsq = NULL; /* sum of squared voxel values */
00496    float * pNumAvg;  /* array for number of pts to avg at each time*/
00497    int   * pTimeIndex; /* array of time markers (1st img of each epoch) */
00498    int   * pEpochLength; /* array of epoch lengths */
00499    float ** tempar;
00500    MRI_IMARR *inimar;
00501    
00502    
00503    nx = dset->daxes->nxx;
00504    ny = dset->daxes->nyy;
00505    nz = dset->daxes->nzz;
00506    nvox = dset->daxes->nxx * dset->daxes->nyy * dset->daxes->nzz;
00507    numims = DSET_NUM_TIMES(dset);
00508    datum = DSET_BRICK_TYPE( dset , 0 ) ; /* get old dataset datum type */
00509    
00510    PLUTO_popup_meter(plint) ;
00511    
00512    DSET_load(dset);
00513    
00514    inimar =  dset_to_mri(dset);
00515    if( inimar == NULL ) return NULL ;
00516 
00517    fxar = (float **) malloc( sizeof( float *) * numims);
00518    if( datum == MRI_float){
00519       for( kk=0; kk<numims; kk++){
00520          fxar[kk] = MRI_FLOAT_PTR(IMAGE_IN_IMARR(inimar,kk));
00521       }
00522    }
00523    else{
00524       for( kk=0; kk<numims; kk++){
00525          fxar[kk] = MRI_FLOAT_PTR(mri_to_float(IMAGE_IN_IMARR(inimar,kk)));
00526       }
00527    }
00528 
00529    nia = 0;    /* number of images (timepoints) averaged  where num epochs > 1*/
00530    
00531    if( RMB_DEBUG ) fprintf(stderr, "Start stavg\n");
00532    
00533    /* determine number of epochs to average */
00534    if( RMB_DEBUG ) fprintf(stderr, "Determining number of epochs...");
00535    numepochs = 1;
00536    for( kim=0; kim < numims; kim++ ){
00537       if( ref[kim] > 0) numepochs++;
00538    }
00539    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00540 
00541    /* set array of epoch lengths */
00542    if( RMB_DEBUG ) fprintf(stderr, "Set array of epoch lengths...");
00543    pEpochLength = (int *)malloc(sizeof(int) * numepochs);
00544    for( ne=0; ne < numepochs; ne++) pEpochLength[ne] = 0;
00545    ne = 0;
00546    for( kim=0; kim < numims; kim++ ){
00547       if( ref[kim] > 0) ne++;
00548       pEpochLength[ne]++;
00549    }
00550    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00551 
00552    /* set array of time markers (1st img of each epoch) */
00553    if( RMB_DEBUG ) fprintf(stderr, "Set array of time markers...");
00554    pTimeIndex = (int *)malloc(sizeof(int) * (numepochs - 1));
00555    lastindex = 0;
00556    minlength = numims;
00557    maxlength = 0;
00558    for( ne=0; ne < (numepochs-1); ne++){
00559       pTimeIndex[ne] = lastindex + pEpochLength[ne];
00560       lastindex = pTimeIndex[ne];
00561       if(pEpochLength[ne+1] > maxlength) maxlength = pEpochLength[ne+1];
00562       if(pEpochLength[ne+1] < minlength) minlength = pEpochLength[ne+1];
00563    }
00564    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00565 
00566    if(maxlength > user_maxlength) maxlength = user_maxlength;
00567 
00568    if( RMB_DEBUG ) fprintf(stderr, "init...");
00569    pNumAvg = (float *) malloc( sizeof(float) * maxlength);
00570    outar = (float **) malloc( sizeof(float *) * maxlength);
00571    for( te=0; te < maxlength; te++){
00572       outar[te] = (float *) malloc( sizeof(float) * nvox);
00573       for( ii=0; ii<nvox; ii++){
00574          outar[te][ii] = 0.0;
00575       }
00576    }
00577    if (meth == _STAVG_METH_SIGMA) {
00578        sumsq = (float *) malloc( sizeof(float) * nvox);
00579    }
00580    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00581 
00582    if( RMB_DEBUG ) fprintf(stderr, "Start averaging...");
00583    for( te=0; te < maxlength; te++){
00584       pNumAvg[te] = 0.0;   
00585 
00586      switch(meth) {
00587      default: case _STAVG_METH_MEAN:
00588 
00589       for( ne=0; ne < (numepochs-1); ne++){
00590          tinc = pTimeIndex[ne] + te;
00591          if( te < pEpochLength[ne+1] ){
00592             fac = DSET_BRICK_FACTOR(dset, tinc);
00593             if (fac == 0.0 || fac == 1.0) {
00594                 for( ii=0; ii<nvox; ii++){
00595                     outar[te][ii] += fxar[tinc][ii];
00596                 }
00597             } else {
00598                 for( ii=0; ii<nvox; ii++){
00599                     outar[te][ii] += fxar[tinc][ii] * fac;
00600                 }
00601             }
00602             pNumAvg[te]++;
00603          }
00604       }
00605       for( ii=0; ii<nvox; ii++){
00606          outar[te][ii] = outar[te][ii]/pNumAvg[te];
00607       }
00608       break;
00609 
00610      case _STAVG_METH_SIGMA:
00611 
00612       /* sigma statistic is actually the unbiased standard deviation
00613          calculated with this computational formula:
00614            sqrt( (sum(x^2) - sum(x)^2 / n) / (n-1) )  */
00615       for( ii=0; ii<nvox; ii++){
00616           sumsq[ii] = 0.0;
00617       }
00618       for( ne=0; ne < (numepochs-1); ne++){
00619          tinc = pTimeIndex[ne] + te;
00620          if( te < pEpochLength[ne+1] ){
00621             fac = DSET_BRICK_FACTOR(dset, tinc);
00622             if (fac == 0.0 || fac == 1.0) {
00623                 for( ii=0; ii<nvox; ii++){
00624                     outar[te][ii] += fxar[tinc][ii];
00625                     sumsq[ii] += fxar[tinc][ii] * fxar[tinc][ii];
00626                 }
00627             } else {
00628                 for( ii=0; ii<nvox; ii++){
00629                     scaledval = fxar[tinc][ii] * fac;
00630                     outar[te][ii] += scaledval;
00631                     sumsq[ii] += scaledval * scaledval;
00632                 }
00633             }
00634             pNumAvg[te]++;
00635          }
00636       }
00637       for( ii=0; ii<nvox; ii++){
00638           outar[te][ii] = sqrt( (sumsq[ii] -
00639                                  outar[te][ii] * outar[te][ii] /
00640                                  pNumAvg[te]) /
00641                                 (pNumAvg[te] - 1) );
00642       }
00643       break;
00644 
00645      }
00646 
00647       if( pNumAvg[te] > 1) nia ++;
00648       PLUTO_set_meter(plint, (100*(te+1))/maxlength ) ;
00649    }
00650    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00651    if ( no1 ){     /* ignore images with only one average */
00652       if( nia < maxlength) maxlength = nia;
00653    }
00654    
00655    M_maxlength = maxlength;
00656    
00657    
00658    if( RMB_DEBUG ) fprintf(stderr, "malloc output...");
00659    tempar = (float **) malloc(sizeof(float *) * maxlength);
00660    for( te=0 ; te < maxlength ; te++ ){
00661       tempar[te] = (float *) malloc( sizeof(float) * nvox);
00662    }
00663    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00664    
00665    if( RMB_DEBUG ) fprintf(stderr, "convert to output...");
00666    for( te=0; te < maxlength; te++){
00667       for( ii=0; ii<nvox; ii++){
00668          tempar[te][ii] = outar[te][ii];
00669       }
00670    }
00671    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00672 
00673    /* toss arrays */
00674    if( RMB_DEBUG ) fprintf(stderr, "free mem...");
00675    FREE_IMARR( inimar );
00676    if (meth == _STAVG_METH_SIGMA) free(sumsq);
00677    free( outar );
00678    free( pNumAvg );
00679    free( pTimeIndex );
00680    free( pEpochLength );
00681    if( RMB_DEBUG ) fprintf(stderr, "done\n");
00682    DSET_unload(dset);
00683    
00684    return(tempar);
00685 }

MRI_IMARR * dset_to_mri THD_3dim_dataset   dset
 

Definition at line 689 of file plug_stavg.c.

References ADDTO_IMARR, THD_3dim_dataset::daxes, DSET_ARRAY, DSET_BRICK_TYPE, DSET_NUM_TIMES, fout, IMARR_SUBIMAGE, INIT_IMARR, malloc, mri_fix_data_pointer(), mri_new_vol_empty(), THD_dataxes::nxx, THD_dataxes::nyy, nz, and THD_dataxes::nzz.

Referenced by avg_epochs().

00691 {
00692 
00693    int ii, kk, ntime, datum;
00694    int nvox, nx, ny, nz;
00695    int use_fac;
00696    
00697    MRI_IMARR * ims_in;
00698    MRI_IMAGE * im, *temp_im;
00699    
00700 
00701    byte   ** bptr  = NULL ;  /* one of these will be the array of */
00702    short  ** sptr  = NULL ;  /* pointers to input dataset sub-bricks */
00703    float  ** fptr  = NULL ;  /* (depending on input datum type) */
00704    
00705    float * fac  = NULL ;  /* array of brick scaling factors */
00706    
00707    float * fout;
00708    
00709 
00710    ntime = DSET_NUM_TIMES(dset) ;
00711    nx = dset->daxes->nxx;
00712    ny = dset->daxes->nyy;
00713    nz = dset->daxes->nzz;
00714    nvox = dset->daxes->nxx * dset->daxes->nyy * dset->daxes->nzz ;
00715    datum = DSET_BRICK_TYPE( dset , 0 ) ; /* get dataset datum type */
00716 
00717    switch( datum ){  /* pointer type depends on input datum type */
00718 
00719       default:
00720          return NULL  ;
00721 
00722       /** create array of pointers into old dataset sub-bricks **/
00723 
00724       /*--------- input is bytes ----------*/
00725       /* voxel #i at time #k is bptr[k][i] */
00726       /* for i=0..nvox-1 and k=0..ntime-1.  */
00727 
00728       case MRI_byte:
00729          bptr = (byte **) malloc( sizeof(byte *) * ntime ) ;
00730          if( bptr == NULL ) return NULL ;
00731          for( kk=0 ; kk < ntime ; kk++ )
00732             bptr[kk] = (byte *) DSET_ARRAY(dset,kk) ;
00733       break ;
00734 
00735       /*--------- input is shorts ---------*/
00736       /* voxel #i at time #k is sptr[k][i] */
00737       /* for i=0..nvox-1 and k=0..ntime-1.  */
00738 
00739       case MRI_short:
00740          sptr = (short **) malloc( sizeof(short *) * ntime ) ;
00741          if( sptr == NULL ) return NULL ;
00742          for( kk=0 ; kk < ntime; kk++ )
00743             sptr[kk] = (short *) DSET_ARRAY(dset,kk) ;
00744       break ;
00745 
00746       /*--------- input is floats ---------*/
00747       /* voxel #i at time #k is fptr[k][i] */
00748       /* for i=0..nvox-1 and k=0..ntime-1.  */
00749 
00750       case MRI_float:
00751          fptr = (float **) malloc( sizeof(float *) * ntime) ;
00752          if( fptr == NULL ) return NULL ;
00753          for( kk=0 ; kk < ntime; kk++ )
00754             fptr[kk] = (float *) DSET_ARRAY(dset,kk) ;
00755       break ;
00756 
00757    } /* end of switch on input type */
00758    
00759    INIT_IMARR(ims_in) ;
00760    for( kk=0 ; kk < ntime ; kk++ ){
00761       im = mri_new_vol_empty( nx , ny , nz , datum ) ;
00762       ADDTO_IMARR(ims_in,im) ;
00763    }
00764    
00765    for( kk=0 ; kk < ntime ; kk++ ){
00766       im = IMARR_SUBIMAGE(ims_in,kk) ;
00767       
00768       switch( datum ){
00769          case MRI_byte:  mri_fix_data_pointer( bptr[kk], im ) ; break ;
00770          case MRI_short: mri_fix_data_pointer( sptr[kk], im ) ; break ;
00771          case MRI_float: mri_fix_data_pointer( fptr[kk], im ) ; break ;
00772       }
00773    }
00774 
00775 
00776    
00777    return(ims_in);
00778 }

DEFINE_PLUGIN_PROTOTYPE PLUGIN_interface* PLUGIN_init int    ncall
 

Definition at line 79 of file plug_stavg.c.

References _STAVG_METH_MEAN, _STAVG_NUM_METHODS, ANAT_ALL_MASK, FUNC_FIM_MASK, global_plint, helpstring, method_strings, PLUTO_add_hint(), PLUTO_set_sequence(), STAVG_main(), and yes_no_strings.

00080 {
00081    PLUGIN_interface * plint ;     /* will be the output of this routine */
00082 
00083    if( ncall > 0 ) return NULL ;  /* one interfaces */
00084 
00085    /*---------------- set titles and call point ----------------*/
00086 
00087    plint = PLUTO_new_interface( "SingleTrial Avg" ,
00088                                 "Averaging of epochs in Single Trial data" ,
00089                                 helpstring ,
00090                                 PLUGIN_CALL_VIA_MENU , STAVG_main  ) ;
00091 
00092    PLUTO_add_hint( plint , "Averaging of epochs in Single Trial data" ) ;
00093 
00094    global_plint = plint ;  /* make global copy */
00095 
00096    PLUTO_set_sequence( plint , "z:Birn" ) ;
00097 
00098    /*--------- 1st line ---------*/
00099 
00100    PLUTO_add_option( plint ,
00101                      "Datasets" ,  /* label at left of input line */
00102                      "Datasets" ,  /* tag to return to plugin */
00103                      TRUE          /* is this mandatory? */
00104                    ) ;
00105 
00106    PLUTO_add_dataset(  plint ,
00107                        "Input" ,          /* label next to button   */
00108                        ANAT_ALL_MASK ,    /* take any anat datasets */
00109                        FUNC_FIM_MASK ,    /* only allow fim funcs   */
00110                        DIMEN_4D_MASK |    /* need 3D+time datasets  */
00111                        BRICK_ALLREAL_MASK /* need real-valued datasets */
00112                     ) ;
00113    PLUTO_add_hint( plint , "Input 3d+t dataset" ) ;
00114 
00115    PLUTO_add_string( plint ,
00116                      "Output" ,  /* label next to textfield */
00117                      0,NULL ,    /* no fixed strings to choose among */
00118                      19          /* 19 spaces for typing in value */
00119                    ) ;
00120    PLUTO_add_hint( plint , "Name of output dataset" ) ;
00121 
00122    /*---------- 2nd line --------*/
00123 
00124    PLUTO_add_option( plint ,
00125                      "Timing" ,
00126                      "Timing" ,
00127                      TRUE
00128                    ) ;
00129 
00130 
00131    PLUTO_add_timeseries(plint, "Stim. Timing");
00132    PLUTO_add_hint( plint , "Stimulus Timing (0 = no task, 1 = task)" ) ;
00133 
00134    PLUTO_add_number( plint ,
00135                      "delta" ,   
00136                      -1000 ,    
00137                      1000 ,  
00138                      0 ,    
00139                      0 ,   
00140                      TRUE
00141                    ) ;
00142    PLUTO_add_hint( plint , "Shift data timecourse by delta before splitting and averaging" ) ;
00143 
00144    /*---------- 3rd line: computation ----------*/
00145 
00146    PLUTO_add_option( plint ,
00147                      "Compute" ,  /* label at left of input line */
00148                      "Compute" ,  /* tag to return to plugin */
00149                      TRUE         /* is this mandatory? */
00150                    ) ;
00151 
00152    PLUTO_add_string( plint ,
00153                      "Method" ,           /* label next to chooser button */
00154                      _STAVG_NUM_METHODS,  /* number of strings in list */
00155                      method_strings ,     /* list of strings to choose among */
00156                      _STAVG_METH_MEAN     /* index of default string */
00157                    ) ;
00158 
00159    PLUTO_add_hint( plint , "Choose statistic to compute" ) ;
00160 
00161    /*---------- 4th line --------*/
00162 
00163    PLUTO_add_option( plint ,
00164                      "Parameters" ,  /* label at left of input line */
00165                      "Parameters" ,  /* tag to return to plugin */
00166                      FALSE            /* is this mandatory? */
00167                    ) ;
00168 
00169    PLUTO_add_number( plint ,
00170                      "maxlength" ,    /* label next to chooser */
00171                      0 ,         /* smallest possible value */
00172                      1000 ,        /* largest possible value */
00173                      0 ,         /* decimal shift (none in this case) */
00174                      15 ,         /* default value */
00175                      TRUE       /* allow user to edit value? */
00176                    ) ;
00177    PLUTO_add_hint( plint , "maximum # of timepoints of output dataset" ) ;
00178 
00179    PLUTO_add_string( plint ,
00180                      "no1?" ,               /* label next to chooser button */
00181                      2  ,               /* number of strings to choose among */
00182                      yes_no_strings ,  /* list of strings to choose among */
00183                      1                  /* index of default string */
00184                    ) ;
00185 
00186    PLUTO_add_hint( plint , "ignore timepoints where only one image is in average" ) ;
00187 
00188 
00189    /*--------- done with interface setup ---------*/
00190 
00191    return plint ;
00192 }

char * STAVG_main PLUGIN_interface *   
 

Definition at line 213 of file plug_stavg.c.

References _STAVG_NUM_METHODS, abs, ADN_datum_all, ADN_malloc_type, ADN_none, ADN_ntt, ADN_nvals, ADN_prefix, avg_epochs(), DATABLOCK_MEM_MALLOC, THD_3dim_dataset::daxes, DSET_BRICK_TYPE, DSET_NUM_TIMES, DSET_NVALS_PER_TIME, DSET_TIMESTEP, DSET_TIMEUNITS, EDIT_coerce_scale_type(), EDIT_dset_items(), EDIT_empty_copy(), EDIT_substitute_brick(), fout, free, FREE_WORKSPACE, M_maxlength, malloc, method_strings, MRI_FLOAT_PTR, MRI_IMAGE::nx, THD_dataxes::nxx, THD_dataxes::nyy, nz, THD_dataxes::nzz, PLUTO_add_dset(), PLUTO_commandstring(), PLUTO_find_dset(), PLUTO_popup_meter(), PLUTO_prefix_ok(), PLUTO_set_meter(), PLUTO_string_index(), THD_delete_3dim_dataset(), tross_Append_History(), tross_Copy_History(), and yes_no_strings.

Referenced by PLUGIN_init().

00214 {
00215    MCW_idcode * idc ;                          /* input dataset idcode */
00216    THD_3dim_dataset * old_dset , * new_dset ;  /* input and output datasets */
00217    char * new_prefix , * str , * str2;         /* strings from user */
00218    int   meth;                                 /* chosen computation method */
00219    int   new_datum ,                           /* control parameters */
00220          old_datum , ntime ;
00221 
00222    int   te, ne, tinc, kim, nia;
00223    int   numepochs, minlength, maxlength, lastindex, navgpts;
00224    int   nvox , perc , new_units, old_units ;
00225    int   ii, ibot,itop , kk, jj; 
00226    int   no1, user_maxlength, delta;
00227    int   *pEpochLength, *pTimeIndex;
00228    int   nx, ny, nz, npix;
00229    float *pNumAvg;
00230    float old_dtime;
00231 
00232    MRI_IMAGE * stimim;
00233    
00234    MRI_IMARR *avgimar;
00235 
00236    byte   ** bptr  = NULL ;  /* one of these will be the array of */
00237    short  ** sptr  = NULL ;  /* pointers to input dataset sub-bricks */
00238    float  ** fptr  = NULL ;  /* (depending on input datum type) */
00239 
00240    float   * fxar  = NULL ;  /* array loaded from input dataset */
00241    float   * stimar = NULL ;
00242    float  ** fout  = NULL ;  /* will be array of output floats */
00243 
00244    float   * tar   = NULL ;  /* will be array of taper coefficients */
00245 
00246    float   * nstimar;
00247 
00248    /*--------------------------------------------------------------------*/
00249    /*----- Check inputs from AFNI to see if they are reasonable-ish -----*/
00250 
00251    /*--------- go to first input line ---------*/
00252 
00253    PLUTO_next_option(plint) ;
00254 
00255    idc      = PLUTO_get_idcode(plint) ;   /* get dataset item */
00256    old_dset = PLUTO_find_dset(idc) ;      /* get ptr to dataset */
00257    if( old_dset == NULL )
00258       return "*************************\n"
00259              "Cannot find Input Dataset\n"
00260              "*************************"  ;
00261 
00262    ntime = DSET_NUM_TIMES(old_dset) ;
00263    if( ntime < 2 )
00264       return "*****************************\n"
00265              "Dataset has only 1 time point\n"
00266              "*****************************"  ;
00267 
00268    ii = DSET_NVALS_PER_TIME(old_dset) ;
00269    if( ii > 1 )
00270       return "************************************\n"
00271              "Dataset has > 1 value per time point\n"
00272              "************************************"  ;
00273    
00274    old_datum = DSET_BRICK_TYPE( old_dset , 0 ) ; /* get old dataset datum type */
00275    new_datum = old_datum;
00276    old_dtime = DSET_TIMESTEP(old_dset);
00277    old_units = DSET_TIMEUNITS(old_dset);
00278    
00279    nvox = old_dset->daxes->nxx * old_dset->daxes->nyy * old_dset->daxes->nzz;
00280    npix = old_dset->daxes->nxx * old_dset->daxes->nyy;
00281    nx = old_dset->daxes->nxx;
00282 
00283 
00284    new_prefix = PLUTO_get_string(plint) ;   /* get string item (the output prefix) */
00285    if( ! PLUTO_prefix_ok(new_prefix) )      /* check if it is OK */
00286       return "************************\n"
00287              "Output Prefix is illegal\n"
00288              "************************"  ;
00289 
00290    /*--------- go to next input line ---------*/
00291 
00292    PLUTO_next_option(plint);
00293 
00294    stimim = PLUTO_get_timeseries(plint);
00295    if( stimim == NULL ) return "Please specify stimulus timing";
00296 
00297    if( stimim->nx < ntime ){
00298       return "**************************************\n"
00299              "Not enough pts in stimulus time-series\n"
00300              "**************************************";
00301    }
00302 
00303    stimar = MRI_FLOAT_PTR(stimim);
00304 
00305 
00306    delta = PLUTO_get_number(plint);
00307 
00308    if( abs(delta) > ntime ){
00309       return "************************\n"
00310              "Delta shift is too large\n"
00311              "************************";
00312    }
00313   
00314    /*initialize variables if not user specified */
00315    user_maxlength = ntime;
00316    no1 = 0;
00317 
00318    /*--------- go to next input line ---------*/
00319 
00320    PLUTO_next_option(plint);
00321 
00322    str  = PLUTO_get_string(plint) ;      /* get string item (the method) */
00323    meth = PLUTO_string_index( str ,      /* find it in list it is from */
00324                               _STAVG_NUM_METHODS ,
00325                               method_strings ) ;
00326 
00327    /*--------- see if the 4th option line is present --------*/
00328 
00329    str = PLUTO_get_optiontag( plint ) ;
00330    if( str != NULL ){
00331       user_maxlength = (int) PLUTO_get_number(plint) ;
00332       str2  = PLUTO_get_string(plint) ;      /* get string item (the method) */
00333       no1   = PLUTO_string_index( str2 ,      /* find it in list it is from */
00334                                  2 ,
00335                                  yes_no_strings) ;
00336    }
00337    
00338 
00339    /*------------------------------------------------------*/
00340    /*---------- At this point, the inputs are OK ----------*/
00341 
00342    PLUTO_popup_meter( plint ) ;  /* popup a progress meter */
00343 
00344    /*________________[ Main Code ]_________________________*/
00345   
00346    fout = avg_epochs( old_dset, stimar, user_maxlength, 1, meth, plint );
00347    if( fout == NULL ) return " \nError in avg_epochs() function!\n " ;
00348    
00349    if( RMB_DEBUG ) fprintf(stderr, "Done with avg_epochs\n");
00350    maxlength = M_maxlength;
00351    
00352    
00353    /*______________________________________________________*/
00354 
00355    
00356    new_dset = EDIT_empty_copy( old_dset ) ; /* start with copy of old one */
00357 
00358    { char * his = PLUTO_commandstring(plint) ;
00359      tross_Copy_History( old_dset , new_dset ) ;
00360      tross_Append_History( new_dset , his ) ; free( his ) ;
00361    }
00362    
00363    /*-- edit some of its internal parameters --*/
00364    ii = EDIT_dset_items(
00365            new_dset ,
00366               ADN_prefix      , new_prefix ,           /* filename prefix */
00367               ADN_malloc_type , DATABLOCK_MEM_MALLOC , /* store in memory */
00368               ADN_datum_all   , new_datum ,            /* atomic datum */
00369               ADN_nvals       , maxlength ,            /* # sub-bricks */
00370               ADN_ntt         , maxlength ,            /* # time points */
00371            /*   ADN_ttorg       , old_dtime ,  */              /* time origin */
00372            /*   ADN_ttdel       , old_dtime ,  */            /* time step */
00373            /*   ADN_ttdur       , old_dtime ,  */            /* time duration */
00374            /*   ADN_nsl         , 0 ,          */        /* z-axis time slicing */
00375            /*   ADN_tunits      , old_units ,  */        /* time units */
00376            ADN_none ) ;
00377 
00378    if( ii != 0 ){
00379       THD_delete_3dim_dataset( new_dset , False ) ;
00380       FREE_WORKSPACE ;
00381       return "***********************************\n"
00382              "Error while creating output dataset\n"
00383              "***********************************"  ;
00384    }
00385 
00386 
00387    /*------------------------------------------------------------*/
00388    /*------- The output is now in fout[kk][ii],
00389              for kk=0..maxlength-1 , ii=0..nvox-1.
00390              We must now put this into the output dataset -------*/
00391 
00392    switch( new_datum ){
00393 
00394       /*** output is floats is the simplest:
00395            we just have to attach the fout bricks to the dataset ***/
00396 
00397       case MRI_float:
00398          for( kk=0 ; kk < maxlength ; kk++ )
00399             EDIT_substitute_brick( new_dset , kk , MRI_float , fout[kk] ) ;
00400       break ;
00401 
00402       /*** output is shorts:
00403            we have to create a scaled sub-brick from fout ***/
00404 
00405       case MRI_short:{
00406          short * bout ;
00407          float fac ; 
00408 
00409          for( kk=0 ; kk < maxlength ; kk++ ){  /* loop over sub-bricks */
00410 
00411             /*-- get output sub-brick --*/
00412             bout = (short *) malloc( sizeof(short) * nvox ) ;
00413             if( bout == NULL ){
00414                fprintf(stderr,"\nFinal malloc error in plug_stavg!\n\a") ;
00415                return("Final malloc error in plug_stavg!"); ;
00416                /*  exit(1) ;*/
00417             }
00418 
00419             /*-- find scaling and then scale --*/
00420             /*fac = MCW_vol_amax( nvox,1,1 , MRI_float , fout[kk] ) ;*/
00421             fac = 1.0;
00422             EDIT_coerce_scale_type( nvox,fac ,
00423                                     MRI_float,fout[kk] , MRI_short,bout ) ;
00424             free( fout[kk] ) ;  /* don't need this anymore */
00425 
00426             /*-- put output brick into dataset, and store scale factor --*/
00427             EDIT_substitute_brick( new_dset , kk , MRI_short , bout ) ;
00428          }
00429       }
00430       break ;
00431 
00432       /*** output is bytes (byte = unsigned char)
00433            we have to create a scaled sub-brick from fout ***/
00434 
00435       case MRI_byte:{
00436          byte * bout ;
00437          float fac ;
00438 
00439          for( kk=0 ; kk < maxlength ; kk++ ){  /* loop over sub-bricks */
00440 
00441             /*-- get output sub-brick --*/
00442 
00443             bout = (byte *) malloc( sizeof(byte) * nvox ) ;
00444             if( bout == NULL ){
00445                fprintf(stderr,"\nFinal malloc error in plug_stavg!\n\a") ;
00446                return("Final malloc error in plug_stavg!"); ;
00447                /*               exit(1) ;*/
00448             }
00449 
00450             /*-- find scaling and then scale --*/
00451 
00452             fac = 1.0;
00453             EDIT_coerce_scale_type( nvox,fac ,
00454                                     MRI_float,fout[kk] , MRI_byte,bout ) ;
00455 
00456             free( fout[kk] ) ;  /* don't need this anymore */
00457 
00458             /*-- put output brick into dataset, and store scale factor --*/
00459 
00460             EDIT_substitute_brick( new_dset , kk , MRI_byte , bout ) ;
00461          }
00462 
00463       }
00464       break ;
00465 
00466    } /* end of switch on output data type */
00467 
00468    /*-------------- Cleanup and go home ----------------*/
00469 
00470    PLUTO_set_meter( plint , 100 ) ;  /* set progress meter to 100% */
00471 
00472    PLUTO_add_dset( plint , new_dset , DSET_ACTION_MAKE_CURRENT ) ;
00473 
00474    FREE_WORKSPACE ;
00475    return NULL ;  /* null string returned means all was OK */
00476 }

Variable Documentation

PLUGIN_interface* global_plint = NULL [static]
 

Definition at line 60 of file plug_stavg.c.

Referenced by PLUGIN_init().

char helpstring[] [static]
 

Initial value:

  "Purpose: Averaging epochs of single trial data\n"
  "\n"
  "Input items to this plugin are:\n"
  "   Datasets:   Input  = 3D+time dataset to process\n"
  "                      = reference time-series\n"
  "               Output = Prefix for new dataset\n"
  "   Additional Parameters\n"
  "               delta     = shift timeseries by delta before splitting and averaging\n"
  "               method    = type of statistic to calculate\n"
  "               maxlength = maximum avg ts length\n"
  "               no1?      = images w/ only one img in avg ignored\n"
  "Author -- RM Birn"

Definition at line 23 of file plug_stavg.c.

Referenced by PLUGIN_init().

int M_maxlength
 

Definition at line 61 of file plug_stavg.c.

Referenced by avg_epochs(), and STAVG_main().

char* method_strings[] = { "Mean" , "Sigma" } [static]
 

Definition at line 41 of file plug_stavg.c.

Referenced by PLUGIN_init(), and STAVG_main().

char* yes_no_strings[] = { "No" , "Yes" } [static]
 

Definition at line 40 of file plug_stavg.c.

Referenced by PLUGIN_init(), and STAVG_main().

 

Powered by Plone

This site conforms to the following standards: