#!/bin/tcsh

@global_parse `basename $0` "$*" ; if ($status) exit 0

# ----------------------------------------------------------------------
# Create correlation volumes (correlate with average of local spheres)
# to look for scanner (or maybe motion) artifacts in EPI data.
#
# Optionally look for large clusters of high correlations.
#
#    inputs: EPI datasets
#    output: a directory containing correlation volumes
#
# This program was originally written to look for artifacts from coils
# (in a multi-channel system) where the signal would drop by 10-40% for
# significant periods of time.
# ----------------------------------------------------------------------


# ----------------------------------------------------------------------
# main input variables ($run is the typical one to set)
set nfirst = 3
set polort = 2          # polort for trend removal
set percentile = 80
set cthresh = 0.9       # use 0 to turn off
set sphere_rad = 20     # use 0 to turn off
set frac_limit = 0.02
set min_thr = 0.45
set rdir = corr_test.results
set din_list = ( )
set do_corr = yes       # use no for testing only
set do_clust = no       # use no to just create correlation volumes
set do_clean = no       # remove everything but the correlation datasets?
set clean_dir = .clean
set mask_in = ''

# 3dmerge/localstat options
set use_3dmerge = yes   # new default, as of 9 May, 2019
set merge_frad = 0.0    # 0 means to apply a full gaussian
                        #
                        # might go out to 0.5 of the half with at half max
                        # (make this distance correspond to $sphere_rad, so
                        # radius = $sphere_rad / $merge_frad)
set corr_mask = no      # should we automask before correlation blurring?


set prog = @radial_correlate

set version = "1.0, 23 February, 2024"

if ( $#argv < 1 ) goto SHOW_HELP

set ac = 1
while ( $ac <= $#argv )
   if ( "$argv[$ac]" == "-cthresh" ) then
      @ ac += 1
      set cthresh = $argv[$ac]
   else if ( "$argv[$ac]" == "-corr_mask" ) then
      @ ac += 1
      set corr_mask = $argv[$ac]
      if ( "$corr_mask" != 'yes' && "$corr_mask" != 'no' ) then
         echo "** -corr_mask requires either 'yes' or 'no'"
         exit 1
      endif
   else if ( "$argv[$ac]" == "-do_clean" ) then
      @ ac += 1
      set do_clean = $argv[$ac]
      if ( "$do_clean" != 'yes' && "$do_clean" != 'no' ) then
         echo "** -do_clean requires either 'yes' or 'no'"
         exit 1
      endif
   else if ( "$argv[$ac]" == "-do_clust" ) then
      @ ac += 1
      set do_clust = $argv[$ac]
      if ( "$do_clust" != 'yes' && "$do_clust" != 'no' ) then
         echo "** -do_clust requires either 'yes' or 'no'"
         exit 1
      endif
   else if ( "$argv[$ac]" == "-do_corr" ) then
      @ ac += 1
      set do_corr = $argv[$ac]
      if ( "$do_corr" != 'yes' && "$do_corr" != 'no' ) then
         echo "** -do_corr requires either 'yes' or 'no'"
         exit 1
      endif
   else if ( "$argv[$ac]" == "-frac_limit" ) then
      @ ac += 1
      set frac_limit = $argv[$ac]
   else if ( "$argv[$ac]" == "-nfirst" ) then
      @ ac += 1
      set nfirst = $argv[$ac]
   else if ( "$argv[$ac]" == "-merge_frad" ) then
      @ ac += 1
      set merge_frad = $argv[$ac]
   else if ( "$argv[$ac]" == "-min_thr" ) then
      @ ac += 1
      set min_thr = $argv[$ac]
   else if ( "$argv[$ac]" == "-mask" ) then
      @ ac += 1
      set mask_in = $argv[$ac]
      set corr_mask = yes
   else if ( "$argv[$ac]" == "-percentile" ) then
      @ ac += 1
      set percentile = $argv[$ac]
   else if ( "$argv[$ac]" == "-polort" ) then
      @ ac += 1
      set polort = $argv[$ac]
   else if ( "$argv[$ac]" == "-rdir" ) then
      @ ac += 1
      set rdir = $argv[$ac]
   else if ( "$argv[$ac]" == "-sphere_rad" ) then
      @ ac += 1
      set sphere_rad = $argv[$ac]
   else if ( "$argv[$ac]" == "-use_3dmerge" ) then
      @ ac += 1
      set use_3dmerge = $argv[$ac]
      if ( "$use_3dmerge" != 'yes' && "$use_3dmerge" != 'no' ) then
         echo "** -use_3dmerge requires either 'yes' or 'no'"
         exit 1
      endif
      # use a mask with 3dLocalstat
      if ( "$use_3dmerge" != 'no' ) then
         set corr_mask = yes
      endif
   else if ( "$argv[$ac]" == "-help" ) then
      goto SHOW_HELP
   else if ( "$argv[$ac]" == "-hist" ) then
      goto SHOW_HIST
   else if ( "$argv[$ac]" == "-ver" ) then
      echo "version $version"
      exit
   else if ( "$argv[$ac]" == "-verb" ) then
      set echo
   else
      # the rest are assumed to be input datasets
      set din_list = ( $argv[$ac-] )
      break
   endif
   @ ac += 1
end

if ( $do_clean == yes && $do_corr != yes ) then
   echo '** "-do_clean yes" requires "-do_corr yes"'
   exit 1
endif

# ----------------------------------------------------------------------
# create dset_list from din_list, and make sure input datasets have
# enough TRs for computation (require at least $nfirst + 2)
@ nneeded = $nfirst + 2
set dset_list = ()
set prefix_list = ()    # in case we name output based on input
foreach dset ( $din_list )
   set nvals = `3dinfo -nt $dset`
   if ( $status || "$nvals" == "NO-DSET" ) then
      echo ""
      echo "** failed to get NT from dset $dset, skipping..."
      echo ""
      continue
   else if ( $nvals < $nneeded ) then
      echo ""
      echo "** require $nneeded TRs, but have $nvals in $dset, skipping..."
      echo ""
      continue
   endif

   set prefix = `3dinfo -prefix_noext $dset`

   # otherwise, dataset looks good, add to list
   set dset_list = ( $dset_list $dset )
   set prefix_list = ( $prefix_list $prefix )
end

# if there are no remaining datasets, bail
if ( $#dset_list < 1 ) then
   echo "** missing input datasets"
   echo ""
   exit 1
endif

# note view
set view = `3dinfo -av_space $dset_list[1]`
if ( $status ) then
   echo "** cannot determine view in $dset_list[1]"
   exit 1
endif

# ----------------------------------------------------------------------
# make results dir, enter it and remove old results

if ( $do_corr == 'yes' ) then
   if ( -d $rdir ) then
      echo "** error, results dir already exists: $rdir"
      exit 1
   endif

   mkdir $rdir
   if ( $status ) then
      echo "** failed to create results directory, $rdir"
      exit 1
   endif

   # if cleaning, make a temporary dir for storage
   if ( $do_clean == yes ) then
      mkdir $rdir/$clean_dir
   endif

else if ( ! -d $rdir ) then
   echo "** error, missing results dir : $rdir"
   exit 1
endif


# ----------------------------------------------------------------------
# copy inputs (minus pre-ss trs) to results dir

if ( $do_corr == 'yes' ) then
   set rind = 1
   set dind = 1
   foreach index ( `count_afni -digits 1 1 $#dset_list` )
      set inset = $dset_list[$index]

      # special known datasets will get 'd'
      if ( $inset =~ errts* || $inset =~ all_runs* ) then
         set ichr  = d
         set ind02 = `ccalc -form '%02d' $dind`
         @ dind += 1
      else
         set ichr  = r     # r or d, for run or dset
         set ind02 = `ccalc -form '%02d' $rind`
         @ rind += 1
      endif

      3dTcat -prefix $rdir/epi.$ichr$ind02 $inset\[$nfirst..\$]

      if ( $status ) then
          echo "** failed to copy EPI dataset $inset, exiting..."
          exit 1
      endif
   end

   if ( $mask_in != "" ) then
      3dcopy $mask_in $rdir/
      if ( $status ) then
          echo "** failed to copy -mask dataset $mask_in, exiting..."
          exit 1
      endif
   endif
endif

# ----------------------------------------------------------------------
# enter results dir
cd $rdir

set res_list = ()
set thr_list = ()
set frac_list = ()

# make a file the has name correspondence
set name_file = dataset_names.txt
touch $name_file

# ----------------------------------------------------------------------
# main work: compute correlation datasets and surviving threshold fractions
set rind = 1
set dind = 1
foreach index ( `count_afni -digits 1 1 $#dset_list` )

   # resulting variables

   # ** cases of special datasets: errts, all_runs
   #    (any more and we make a list)
   set dset = $dset_list[$index]
   if ( $dset =~ errts* ) then
      set ichr = d
      set ind02 = `ccalc -form '%02d' $dind`
      set prefix = radcor.$sphere_rad.$ichr$ind02.errts
      @ dind += 1
   else if ( $dset =~ all_runs* ) then
      set ichr = d
      set ind02 = `ccalc -form '%02d' $dind`
      set prefix = radcor.$sphere_rad.$ichr$ind02.all_runs
      @ dind += 1
   else
      set ichr  = r     # r or d, for run or dset
      set ind02 = `ccalc -form '%02d' $rind`
      set prefix = radcor.$sphere_rad.$ichr$ind02
      @ rind += 1
   endif

   set eset = epi.$ichr$ind02
   set ulay = epi.ulay.$ichr$ind02      # save underlay dataset if clean?
   set mpre = $prefix.automask
   set mset = $mpre$view
   set sset = $prefix.scaled
   set cset = $prefix.corr

   set eset_orig = $eset                # store original

   if ( $mask_in != "" ) then
      set mset = $mask_in
   endif

   echo "$prefix    <---    $dset" >> $name_file

   # ----------------------------------------
   # remove pre-SS TRs and generate a mask
   if ( $do_corr == 'yes' && $mask_in == "" ) then
      3dAutomask -q -prefix $mpre $eset$view
   endif
   set nmask = `3dBrickStat -non-zero -count $mset`

   if ( $do_corr == 'yes' && -f tt.scale1$view.HEAD ) then
      \rm -f tt.*
   endif

   # ----------------------------------------
   # possibly detrend (after any automask)
   if ( $polort >= 0 ) then
      set detset = det.$ichr$ind02
      echo "-- detrend -polort $polort, new eset = $detset"
      3dTproject -quiet -polort $polort -prefix $detset -input $eset$view
      set eset = $detset
   endif

   # ----------------------------------------
   # compute the correlation dataset, either with the average scaled
   # time series or within the given sphere radius
   # (in either case, result is $cset$view)

   if ( $do_corr == 'yes' ) then
      echo ""
      echo "-- running correlation on dataset $dset_list[$index] ..."
      echo ""

      if ( $sphere_rad != 0 ) then
         set rad = $sphere_rad
         if ( $use_3dmerge == yes ) then
            # mask before blur?
            if ( $corr_mask == yes ) then
               set enew = rm.emasked.$ichr$ind02
               3dcalc -a $eset$view -b $mset -expr "a*b" -prefix $enew
               set eset = $enew
            endif

            # if $merge_frad is non-zero, use truncation (cubical extents?)
            if ( `ccalc -n "ispositive($merge_frad)"` ) then
               # scale rad by frac number of HWHM, and to FWHM
               set rad2 = `ccalc -n "2*$sphere_rad/$merge_frad"`
               # FIRFAC = 1.17741 when nhwhm=1, so scale
               set firfac = `ccalc "1.17741*$merge_frad"`
               echo "++ merge blur: rad $sphere_rad, rad2 $rad2, firfac $firfac"

               3dmerge -DAFNI_BLUR_FIRFAC=$firfac -1blur_fwhm $rad2 -doall \
                       -prefix sphere.mean.$rad.$ichr$ind02 $eset$view
            else
               # just convert to FWHM
               set rad2 = `ccalc -n "2*$sphere_rad"`
               echo "++ merge blur: rad $sphere_rad, FWHM $rad2"
               3dmerge -1blur_fwhm $rad2 -doall                         \
                       -prefix sphere.mean.$rad.$ichr$ind02 $eset$view
            endif
         else
            3dLocalstat -quiet -nbhd "SPHERE($rad)" -mask $mset \
               -stat mean -prefix sphere.mean.$rad.$ichr$ind02 $eset$view
         endif
         # correlate blur/ave dset with base dset
         3dTcorrelate -prefix $cset sphere.mean.$rad.$ichr$ind02$view $eset$view
      else
         # get a demeaned and scaled dataset (so sumsq = 1)
         # (only needed to use simple dot project, so forget)

         # 3dTstat -mean -datum float -prefix $prefix.mean $eset_orig$view
         # 3dcalc -a $eset$view -b $prefix.mean$view -c $mset \
         #        -expr 'c*(a-b)/b' -datum float -prefix tt.scale1
         # 3dTstat -sos -datum float -prefix tt.sos tt.scale1$view
         # 3dcalc -a tt.scale1$view -b tt.sos$view -c $mset   \
         #         -expr 'c*a/sqrt(b)' -prefix $sset

         3dmaskave -quiet -mask $mset $eset$view > $prefix.mean_signal.1D

         3dTcorr1D -mask $mset -pearson -prefix $cset \
                 $eset$view $prefix.mean_signal.1D
      endif

      # verify that we have a result (since no -e, as we grep)
      set testfile = $cset$view.HEAD
      if ( ! -f $testfile ) then
         echo "** missing correlation result $testfile"
         exit 1
      endif

      # extract ulay volume
      3dbucket -prefix $ulay $eset_orig$view'[0]'

      # if cleaning, hide what we want to keep
      # (to easily delete the rest)
      if ( $do_clean == yes ) then
          \mv $ulay* $clean_dir
          \cp $cset$view* $clean_dir
      endif
   endif

   # if we are not clustering, continue on to next dataset
   if ( $do_clust != 'yes' ) then
      echo "---- done\n"
      continue
   endif

   # ----------------------------------------
   # now look for the largest cluster of voxels thresholded at 80-percentile

   # threshold correlations via either cthresh or at the percentile
   set thr = $cthresh
   if ( $thr == 0 ) then
      set thr = `3dBrickStat -mask $mset -positive         \
                -percentile $percentile 1 $percentile           \
                $cset$view | awk '{print $2}'`
   endif
   set thr_list = ( $thr_list $thr )

   if ( $do_corr == 'yes' ) then
      3dclust -quiet -nosum -2thresh -10 $thr -dxyz=1 1.01 10 \
              $cset$view > $prefix.clust.txt
   endif

   set maxsize = `head -n 1 $prefix.clust.txt    \
                  | grep -v 'NO CLUSTERS' | awk '{print $1}'`

   if ( $maxsize == "" ) then
      echo "-- max clust @ thresh $thr in $cset$view : no clusters"
      set res_list = ( $res_list 0 )
      set frac_list = ( $frac_list 0 )
      continue
   endif

   # compute mask volume fraction
   set frac = `ccalc $maxsize.0/$nmask`

   echo "------------------------------------------------------------"
   echo "-- max clust @ thresh $thr in $cset$view : $maxsize"
   echo "   fraction of mask volume ($nmask) : $frac"

   # check for min threshold
   set fail = `ccalc -i "step($min_thr-$thr)"`
   if ( $fail ) then
      echo "== dataset $dset ($eset) looks okay"
      echo "   $percentile % threshold $thr below cutoff $min_thr, pass..."
      set res_list = ( $res_list 1 )
      set frac_list = ( $frac_list $frac )
      echo "------------------------------------------------------------"
      continue
   endif

   # check for min volume fraction
   set bad_level = `ccalc -i "1+2.0*$frac/$frac_limit"`
   set frac_list = ( $frac_list $frac )
   if ( $bad_level == 1 ) then
       set res_list = ( $res_list $bad_level )
       echo "== dataset $dset ($eset) looks okay"
       echo "   (max clust fraction $frac below limit $frac_limit)"
   else if ( $bad_level == 2 ) then
       set res_list = ( $res_list $bad_level )
       echo "== warning: dataset $dset ($eset) might be problematic"
       echo "   (max clust fraction $frac near limit $frac_limit)"
   else
       set res_list = ( $res_list $bad_level )
       echo "== error: dataset $dset ($eset) looks problematic"
       echo "   (max clust fraction $frac above limit $frac_limit)"
   endif
   echo "------------------------------------------------------------"

end  # foreach index

# if we did clustering, tabulate the results
if ( $do_clust == 'yes' ) then
   echo "============================================================"
   if ( $sphere_rad == 0 ) then
      echo "correlations are against spheres of radius $sphere_rad mm"
   else
      echo "correlations are against average time series over automask"
   endif
   echo ""
   echo "result    dataset (corr thresh / cluster fraction)"
   echo "------    --------------------------"

   set thr = $cthresh

   foreach index ( `count_afni -digits 1 1 $#dset_list` )
      set inset = $dset_list[$index]
      set res = $res_list[$index]

      if ( $res == 0 ) then
         echo -n "pass      "
      else if ( $res == 1 ) then
         echo -n "pass      "
      else if ( $res == 2 ) then
         echo -n "warning   "
      else
         echo -n "FAILURE   "
         ### [PT: Jun 15, 2021] make this a multiline if condition, to
         ### avoid possible badness on older tcsh versions; also, fix
         ### a likely mistake: 'set $thr = ...' -> 'set thr = ...'
         #if ( $thr == 0 ) set $thr = $thr_list[$index]
         if ( $thr == 0 ) then
            set thr = $thr_list[$index]
         endif
      endif

      echo "$inset ($thr_list[$index] / $frac_list[$index])"
   end
   echo "============================================================"
endif  # $do_clust

# ----------------------------------------------------------------------
# if do_clean, clean up and run away
if ( $do_clean == yes ) then
   echo "++ have do_clean, cleaning up..."
   echo ""

   \rm -f *
   \mv $clean_dir/* .
   \rmdir $clean_dir

   exit 0
endif


# ----------------------------------------------------------------------
# maybe suggest how the user could look at the results

# be sure to provide some threshold
set thr = 0.5

cat << EOF

   ___________________________________________________________________

   - to review, consider looking at correlations over epi time series
        run command: afni $rdir
        then set:    Underlay  = epi.SOMETHING
                     Overlay   = radcor.SOMETHING.corr
                     maybe threshold = $thr, maybe clusterize

   -  processed $#dset_list of $#din_list dataset(s)
   -  for name correspondence, see: $rdir/$name_file
   ___________________________________________________________________

EOF

# ----------------------------------------------------------------------
# DONE

exit


SHOW_HELP:
# ----------------------------------------------------------------------
# help string

echo "-----------------------------------------------------------------"
echo "$prog       - check datasets for correlation artifact"
echo ""
echo "    usage : $prog [options] datasets ..."
echo ""
echo "This program computes the correlation at each voxel with the average"
echo "time series in a 20 mm radius (by default).  If there is basically"
echo "one high-correlation cluster, it is suggestive of a coil artifact."
echo ""
echo "Note that significant motion can also cause such an effect.  But"
echo "while motion correlations will tend to follow the edge of the brain,"
echo "coil artifacts will tend to appear in large, dense clusters."
echo ""
echo "If people really care, I may add an option to see how large a sphere"
echo "might fit within the biggest cluster.  A big sphere would be more"
echo "suggestive of a coil artifact, rather than motion.  But adding such"
echo "an option sounds suspiciously like work."
echo ""
echo "  inputs: a list of EPI datasets (after any options)"
echo "  output: a directory containing correlation volumes (and more)"
echo ""
echo "-----------------------------------------------------------------"
echo ""
echo "Common examples (note that datasets are always passed last):"
echo ""
echo "  1a. Run default operation on a list of EPI datasets (so just create"
echo "      the correlation volumes)."
echo ""
echo "          $prog pb00.FT.*.HEAD"
echo ""
echo "  1b. Similar to 1a, but specify a results directory for correlations."
echo ""
echo "          $prog -rdir new.results pb00.FT.*.HEAD"
echo ""
echo "  2.  Do a cluster test on existing correlation volumes.  Note that"
echo "      this still uses the results directory variable, rdir."
echo ""
echo "          $prog -do_corr no -do_clust yes pb00.FT.*.HEAD"
echo ""
echo "  3.  Run a complete test, both creating the correlation volumes, and"
echo "      then looking for large clusters of high correlations."
echo "      Specify a mask."
echo ""
echo "          $prog -do_clust yes -mask full_mask.FT+orig pb00.FT.*.HEAD"
echo ""
echo "  4.  Run a complete test, but alter some clustering options."
echo "        - threshold at 0.7 (instead of the default 0.9)"
echo "        - increase the minimum cluster size (frac of mask) to 0.05"
echo "        - decrease the correlation sphere radius (from 20 mm) to 10 mm"
echo ""
echo "          $prog -do_clust yes                   \\"
echo "              -cthresh 0.7 -frac_limit 0.05 -sphere_rad 10  \\"
echo "              pb00.FT.*.HEAD"
echo ""
echo "-----------------------------------------------------------------"
echo ""
echo "Overview of processing steps: "
echo ""
echo "0. The first 3 TRs are removed from the input (see -nfirst),"
echo "   and an automask is created (limiting all future computations)."
echo "   Any -mask overrides the automask operation."
echo "   If -do_corr is 'no', this is skipped."
echo ""
echo "   (see -do_corr)"
echo ""
echo "1. The correlation dataset is created (unless -do_corr is 'no')."
echo ""
echo "   (see -sphere_rad, -do_corr, -do_clust)"
echo ""
echo "   At each voxel, compute the correlation either within a sphere"
echo "   or with the average masked time series."
echo ""
echo "   a. within a sphere (if -sphere_rad is not 0)"
echo ""
echo "      At each voxel, compute the average time series within a"
echo "      sphere of radius 20 mm (see -sphere_rad), and correlate the"
echo "      time series with this averaged result."
echo ""
echo "   b. with the average masked time series (if -sphere_rad is 0)"
echo ""
echo "      The demeaned data is scaled to have unit length (sumsq=1)."
echo "      Then compute the mean time series over the automask ROI"
echo "      (so across the expected brain)."
echo "      Correlate each voxel time series with the mean time series."
echo ""
echo "   If -do_clust is 'no', this is the last step."
echo ""
echo "2. Threshold the result (if -do_clust is 'yes')."
echo ""
echo "   (see -cthresh, -percentile, -do_clust)"
echo ""
echo "   Threshold the correlations either at a static value (see -cthresh),"
echo "   or at a certain percentile (see -percentile)."
echo ""
echo "   a. at r=cthresh (if -cthresh is not 0)"
echo ""
echo "      Simply threshold the correlations at this value, maybe 0.9."
echo ""
echo "      (see -cthresh)"
echo ""
echo "   b. at r=percentile (if -cthresh is 0)"
echo ""
echo "      Compute the given percentile (maybe 80), and threshold at"
echo "      that value, whatever it turns out to be."
echo ""
echo "      Note that when using an 80-percent threshold, for example,"
echo "      then 20-percent of the voxels should survive the cutoff."
echo "      Later, the question will be how they cluster."
echo ""
echo "      (see -percentile)"
echo ""
echo "3. if the percentile threshold is too small, considered the data okay"
echo ""
echo "   (see -min_thr)"
echo ""
echo "   In the case of -percentile above (meaning -cthresh is 0), if"
echo "   the resulting threshold is not large enough, then we do not"
echo "   expect the data to have a problem."
echo ""
echo "4. compare largest cluster to mask volume"
echo ""
echo "   (see -frac_limit)"
echo ""
echo "   Compute the size of the largest correlation cluster above the"
echo "   previous threshold (either -cthresh or via -percentile).  Then"
echo "   compute the fraction of the mask volume that this cluster"
echo "   occupies."
echo ""
echo "   If the largest cluster is a large fraction of the mask, then"
echo "   we expect there might be a problem (because most of the high"
echo "   correlation voxels are in one cluster)."
echo ""
echo "   Otherwise, if the high-correlation voxels are scattered about"
echo "   the volume, we do not expect any problem."
echo ""
echo "   For example, if the largest surviving cluster is more than 5%"
echo "   of the mask, the data is consider to FAIL (see -frac_limit)."
echo ""
echo "-----------------------------------------------------------------"
echo ""
echo "    usage : $prog [options] datasets ..."
echo ""
echo "---------------------------------------------"
echo ""
echo "general options:"
echo ""
echo "   -help             : show this help"
echo ""
echo "   -hist             : show modification history"
echo ""
echo "   -do_clean yes/no  : clean up at end, leaving only correlations"
echo "                       default = $do_clean"
echo ""
echo "                       In the case of computing correlations, this"
echo "                       option can be used to remove everything but those"
echo "                       correlation datasets, to save disk space."
echo ""
echo "   -do_clust yes/no  : clust correlation volumes? (yes or no)"
echo "                       default = $do_clust"
echo ""
echo "                       If 'no', only create the correlation volumes."
echo "                       Otherwise, run clustering and look for large"
echo "                       artifacts from bad coil channels."
echo ""
echo "   -do_corr yes/no   : create correlation volumes (yes or no)"
echo "                       default = $do_corr"
echo ""
echo "                       If 'yes', create the correlation volumes."
echo "                       If 'no', simply assume they already exist."
echo "                       This is for re-testing a previous execution."
echo ""
echo "   -polort POLORT    : detrend time series with given poly degree"
echo "                       default = $polort"
echo ""
echo "   -rdir RESULTS_DIR : directory to do computations in"
echo "                       default = $rdir"
echo ""
echo "   -use_3dmerge yes/no: use 3dmerge rather than 3dLocalstat"
echo "                       default = $use_3dmerge"
echo ""
echo "                       For computing a local average, 3dmerge can do"
echo "                       basically the same operation as 3dLocalstat, but"
echo "                       250 times as fast (divided by OpenMP speedup)."
echo ""
echo "                       One can make -merge_frad smaller to make the"
echo "                       results more similar, if desirable."
echo ""
echo "   -ver              : show version number"
echo ""
echo "   -verb             : make verbose: set echo"
echo ""
echo "---------------------------------------------"
echo ""
echo "computational options:"
echo ""
echo "   -cthesh THRESH    : threshold on correlation values"
echo "                       (if 0, use percentile, else use this)"
echo "                       default = $cthresh"
echo ""
echo "   -corr_mask yes/no : mask time series before corrlation blurring"
echo "                       default = $corr_mask"
echo ""
echo "                       This defines whether 3dmerge blurring is applied"
echo "                       to a masked dataset."
echo ""
echo "   -mask MASK_DSET   : specify a mask dataset to replace automask"
echo ""
echo "   -frac_limit LIMIT : min mask fraction surviving cluster"
echo "                       default = $frac_limit"
echo ""
echo "   -mask MASK_DSET   : specify a mask dataset to replace automask"
echo "                       default = automask"
echo "                       This mask is expected to cover the brain."
echo ""
echo "   -merge_frad FRAD  : specify a radius fraction for 3dmerge blurring"
echo "                       default = $merge_frad"
echo ""
echo "                       If FRAD is 1, the Gaussian blur kernel will"
echo "                       be applied with a shape out to the normal HWHM"
echo "                       (half width at half max).  That is to say the"
echo "                       the farthest neighbors would contribute 0.5 (half"
echo "                       max) of that of the central voxel."
echo ""
echo "                       FRAC is an inverse scalar on the blur size,"
echo "                       and a proportional scalar fraction on the size"
echo "                       where the blurring ends.  So sphere_rad is always"
echo "                       the applied blur size."
echo ""
echo "                       A smaller fraction will yield a flatter curve."
echo "                       For example, FRAD=0.5 yields a 0.84 relative"
echo "                       contribution at the radial distance, while."
echo "                       doubling the requested blur."
echo ""
echo "                    ** This leads to a cubical region of averaging,"
echo "                       rather than an intended spherical one.  It is not"
echo "                       a big deal, but is worth noting."
echo ""
echo "                       Use FRAD=0.0 to apply a full Gaussian, rather than"
echo "                       the truncated form."
echo ""
echo "   -nfirst NFIRST    : number of initial TRs to remove"
echo "                       default = $nfirst"
echo ""
echo "   -min_thr THR      : min percentile threshold to be considered"
echo "                       default = $min_thr"
echo ""
echo "   -percentile PERC  : percentile to use as threshold"
echo "                       default = $percentile"
echo ""
echo "   -sphere_rad RAD   : generate correlations within voxel spheres"
echo "                       (or Gaussian weighted versions)"
echo "                       (if 0, go against average time series)"
echo "                       default = $sphere_rad"
echo ""
echo "R Reynolds, Aug, 2011"
echo "------------------------------------------------------------"

exit

SHOW_HIST:

cat << EOF

-----------------------------------------------------------------
$prog modification history:

   0.1  : Aug 16, 2011: initial version (was @check_corr_artifact)
   0.2  : Sep  1, 2011: now $prog
          - new default: just create correlation volumes
          - add to AFNI distribution
          - add -do_corr and -do_clust
   0.3  : Aug 12, 2015: add -mask option
   0.4  : May 10, 2019: 3dLocalstat -> 3dmerge (much faster)
          - by default, use 3dmerge for local averaging
          - add -use_3dmerge, -corr_mask and -merge_frad
          - add -do_clean and -verb
          - on clean, save epi.ulay datasets
   0.5  : May 15, 2019:
          - handle special cases of errts and all_runs
          - do this by partitioning results into r* and d* dataset names
   0.6  : May 30, 2019: full Gaussian by default, not truncated
   0.7  : Mar  7, 2022: add -polort; minor cleanup
          - new default: detrend with polort $polort
   0.8  : Mar 17, 2022: change saved ulay to be from orig EPI
   0.9  : Apr  1, 2022: extract ulay also in non-clean case
   1.0  : Feb 23, 2024: fail if no corr dset is made (no -e in script)

EOF

# *** match version, above

exit