#! /bin/tcsh -f

#
# trac-paths
#
# Tractography for a single subject
#
# Original Author: Anastasia Yendiki
#
# Copyright © 2021 The General Hospital Corporation (Boston, MA) "MGH"
#
# Terms and conditions for use, reproduction, distribution and contribution
# are found in the 'FreeSurfer Software License Agreement' contained
# in the file 'LICENSE' found in the FreeSurfer distribution, and here:
#
# https://surfer.nmr.mgh.harvard.edu/fswiki/FreeSurferSoftwareLicense
#
# Reporting: freesurfer@nmr.mgh.harvard.edu
#
#

umask 002

set VERSION = 'trac-paths 7.3.2';
set ProgName = `basename $0`
set inputargs = ($argv)

#------------ Set default options -----------------------------------------#

set PrintHelp = 0	# Print help and exit

set debug = 0		# Generate more output

set DoTime = 0		# Time main commands
set fs_time = ()

set rcfile = ()		# Name of input run command file
set LF = ()		# Name of log file
set CF = ()		# Name of command file

set DoIsRunning = 1	# Create a lock file while processing continues

set RunIt = 1		# If 0, do everything but run commands (for debugging)

#------------ Parse input arguments ---------------------------------------#

if ($#argv == 0) goto usage_exit

set n = `echo $argv | egrep -e -help | wc -l`
if ($n != 0) then
  set PrintHelp = 1
  goto usage_exit
endif

set n = `echo $argv | egrep -e -version | wc -l`
if ($n != 0) then
  echo $VERSION
  exit 0
endif

source $FREESURFER_HOME/sources.csh

goto parse_args
parse_args_return:
goto check_params
check_params_return:

#------------ Read configuration file -------------------------------------#

source $rcfile

if ($intrareg == 1 || $intrareg == 2) then
  set reg = flt
else if ($intrareg == 3) then
  set reg = bbr
endif

if ($interreg == 1 || $interreg == 2) then
  set xspace = mni
else if ($interreg == 3) then
  set xspace = rob
else if ($interreg == 4) then
  set xspace = cvs
  set cvstempdir = `dirname $intertrg`          # Assuming ../mri/norm.mgz
  set cvstempdir = `dirname $cvstempdir`
  set cvstemp    = `basename $cvstempdir`
  set cvstempdir = `dirname $cvstempdir`
  set cvswarp = final_CVSmorph_to$cvstemp
else if ($interreg == 5) then
  set xspace = syn
  setenv MY_MORPHS_DO_NOT_CONFORM_DEAL_WITH_IT
else if ($interreg == 6) then
  set xspace = fnt
endif

if ($DoTime) then
  fs_time ls >& /dev/null
  if (! $status) set fs_time = fs_time
endif

printf '\n\n' >> $LF
echo "New invocation of $ProgName"  >> $LF
printf '\n\n' >> $LF
whoami >> $LF
hostname >> $LF
uname -a >> $LF
limit >> $LF
if (-e /usr/bin/free) then
  echo "" >> $LF
  /usr/bin/free >> $LF
  echo "" >> $LF
endif
if ("`uname -s`" == "Darwin") then
  echo "" >> $LF
  /usr/bin/top -l 1 | grep PhysMem >> $LF
  echo "" >> $LF
endif

#---------- Is this a longitudinal analysis? ----------------------------#

set dolong = 0

if ($?tplist) then
  if ($#tplist > 0) set dolong = 1
endif

#---------- Create output directories -----------------------------------#

set dtdir = $dtroot/$subj
set fsdir = $SUBJECTS_DIR/$subj
set dwidir = $dtdir/dmri
set xfmdir = $dwidir/xfms
set labdir = $dtdir/dlabel
set touchdir = $dtdir/touch

mkdir -p $dwidir $xfmdir $labdir

#---------- Handle run status files -------------------------------------#

# Delete the error file. This is created when error_exit is run.
set ErrorFile = $dtdir/scripts/trac-all.error
rm -f $ErrorFile

# Delete the done file. This is created when the program exits normally.
set DoneFile = $dtdir/scripts/$ProgName.done
rm -f $DoneFile

if ($DoIsRunning) then
  set IsRunningFile = $dtdir/scripts/IsRunning.trac

  if (-e $IsRunningFile) then
    echo ""
    echo "ERROR: it appears that an analysis is already running"
    echo "for $subj based on the presence of $IsRunningFile. It could"
    echo "also be that an analysis was running at one point but"
    echo "died in an unexpected way. If it is the case that there"
    echo "is a process running, you can kill it and start over or"
    echo "just let it run. If the process has died, you should type:"
    echo ""
    echo "rm $IsRunningFile"
    echo ""
    echo "and re-run. Or you can add -no-isrunning to the trac-all"
    echo "command-line. The contents of this file are:"
    echo "----------------------------------------------------------"
    cat  $IsRunningFile
    echo "----------------------------------------------------------"
    exit 1
  endif

  echo "------------------------------"	>  $IsRunningFile 
  echo "SUBJECT $subj"			>> $IsRunningFile 
  echo "DATE `date`"			>> $IsRunningFile 
  echo "USER `whoami`"			>> $IsRunningFile 
  echo "HOST `hostname`"		>> $IsRunningFile 
  echo "PROCESSID $$"			>> $IsRunningFile 
  echo "PROCESSOR `uname -m`"		>> $IsRunningFile 
  echo "OS `uname -s`"			>> $IsRunningFile 
  uname -a				>> $IsRunningFile 
  echo "PROGRAM $ProgName"		>> $IsRunningFile
  echo $VERSION				>> $IsRunningFile 
else
  set IsRunningFile = ()
endif

# Put a copy of myself (this script) in the scripts dir
cp $0 $dtdir/scripts/$ProgName.local-copy

if (! $RunIt) then
  echo "INFO: -dontrun flag is in effect so subsequent commands" |& tee -a $LF
  echo "may not have accurate arguments!" |& tee -a $LF
endif

echo "#-------------------------------------" |& tee -a $LF
echo "$0 $argv" |& tee -a $LF |& tee -a $CF

#------------ Select DWI set to use for path reconstruction ---------------#

echo "#-------------------------------------" |& tee -a $LF
echo "#@# DWI set selection `date`" |& tee -a $LF

if ($dolong) then
  set dirlist = `printf "$dtroot/%s/dmri " $tplist`
else
  set dirlist = $dwidir
endif

foreach dwidir ($dirlist)
  if ($#bmax) then			# Use all DWIs up to a maximum b-value
    set dwiset = dwi.bmax$bmax

    if (! -e $dwidir/$dwiset.nii.gz) then
      set cmd = dmri_bset
      set cmd = ($cmd --in   $dwidir/dwi.nii.gz)
      set cmd = ($cmd --bmax $bmax)
      set cmd = ($cmd --out  $dwidir/$dwiset.nii.gz)
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif
    endif
  else if ($#bshell) then		# Use DWIs from one or more shells
    set dwiset = dwi`printf '.b%s' $bshell`

    if (! -e $dwidir/$dwiset.nii.gz) then
      set cmd = dmri_bset
      set cmd = ($cmd --in  $dwidir/dwi.nii.gz)
      set cmd = ($cmd `printf '--b %s ' $bshell`)
      set cmd = ($cmd --out $dwidir/$dwiset.nii.gz)
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif
    endif
  else					# Use all DWIs
    set dwiset = dwi
  endif

  set cmd = (ln -sf $dwiset.nii.gz $dwidir/data.nii.gz)
  echo $cmd |& tee -a $LF |& tee -a $CF
  if ($RunIt) then
    $cmd |& tee -a $LF
    if ($status) goto error_exit
  endif

  set cmd = (ln -sf $dwiset.bvecs $dwidir/bvecs)
  echo $cmd |& tee -a $LF |& tee -a $CF
  if ($RunIt) then
    $cmd |& tee -a $LF
    if ($status) goto error_exit
  endif

  set cmd = (ln -sf $dwiset.bvals $dwidir/bvals)
  echo $cmd |& tee -a $LF |& tee -a $CF
  if ($RunIt) then
    $cmd |& tee -a $LF
    if ($status) goto error_exit
  endif
end

if (! $dolong) then
  echo "#-------------------------------------" |& tee -a $LF
  echo "#@# Path reconstruction `date`" |& tee -a $LF

  foreach ntrain ($ntrainlist)
    set avgmode = $avgname${ntrain}_${xspace}_$reg

    set outdir = $dtdir/dpath

    if ($dopathsubdirs) then
      set outdir = $outdir/${nsample}samp

      set ptype = ()

      if ($dosegprior) then
        set ptype = seg14
      endif
      if ($dotangprior) then
        if ($#ptype) then
          set ptype = $ptype.tang
        else
          set ptype = tang
        endif
      endif
      if ($doxyzprior) then
        if ($#ptype) then
          set ptype = $ptype.xyz
        else
          set ptype = xyz
        endif
      endif
      if (! $#ptype) then
        set ptype = none
      endif

      set outdir = $outdir/$ptype
    endif

    mkdir -p $outdir

    set outlist = `printf "$outdir/%s_$avgmode " $pathlist`

    set roi1list = `printf \
      "$labdir/$xspace/%s_${avgmode}_end1_dil.nii.gz " $pathlist`
    set roi2list = `printf \
      "$labdir/$xspace/%s_${avgmode}_end2_dil.nii.gz " $pathlist`

    set initlist = ()
    set sdplist = ()
    @ ipath = 1
    while ($ipath <= $#pathlist)
      set pname = $pathlist[$ipath]
      set nc    = $ncpts[$ipath]

      set initlist = ($initlist \
        $labdir/diff/${pname}_${avgmode}_cpts_$nc.txt)
      set sdplist = ($sdplist \
        $labdir/diff/${pname}_${avgmode}_cpts_${nc}_std.txt)

      @ ipath = $ipath + 1
    end

    set priorlist = ()
    foreach pathname ($pathlist)
      set priorlist = ($priorlist \
        $labdir/$xspace/${pathname}_${avgmode}_logprior_str_0.nii.gz \
        $labdir/$xspace/${pathname}_${avgmode}_logprior_str_1.nii.gz)
    end

    set npriorlist = ()
    foreach pathname ($pathlist)
      set npriorlist = ($npriorlist \
        $labdir/$xspace/${pathname}_${avgmode}_fsnnprior \
        $labdir/$xspace/${pathname}_${avgmode}_fsnnids)
    end

    if ($usetrunc) then
      set npriorlist = `printf "%s_all " $npriorlist`
    endif

    set lpriorlist = ()
    foreach pathname ($pathlist)
      set lpriorlist = ($lpriorlist \
        $labdir/$xspace/${pathname}_${avgmode}_fsprior \
        $labdir/$xspace/${pathname}_${avgmode}_fsids)
    end

    if ($usetrunc) then
      set lpriorlist = `printf "%s_all " $lpriorlist`
    endif

    set tpriorlist = `printf \
      "$labdir/$xspace/%s_${avgmode}_tangprior.txt " $pathlist`
    set cpriorlist = `printf \
      "$labdir/$xspace/%s_${avgmode}_curvprior.txt " $pathlist`

    if ($usetrunc) then
      set tpriorlist = `echo $tpriorlist | sed 's/\.txt/_all\.txt/'`
      set cpriorlist = `echo $cpriorlist | sed 's/\.txt/_all\.txt/'`
    endif

    if ($overwrite) then		# Clean up pre-existing directories
      foreach dname ($outlist)
        set cmd = (rm -rf $dname)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif
      end
    else				# Rename pre-existing directories
      foreach dname ($outlist)
        if (-e $dname) then
          @ vno = 0
          while (-e $dname.v$vno)
            @ vno = $vno + 1
          end
          set cmd = (mv $dname $dname.v$vno)
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif
        endif
      end
    endif

    if ($usemaskanat) then
      set brainmask = $labdir/diff/${segname}_mask.$reg.nii.gz
    else
      set brainmask = $labdir/diff/lowb_brain_mask.nii.gz
    endif

    @ irep = 0

    while ($irep <= 5)
      if ($irep) then
        cp $rcfile $rcfile.reinit
        echo "set dopriors = 1"			>> $rcfile.reinit
        echo "set reinit = 1"			>> $rcfile.reinit
        echo "set pathlist = ($rep_pathlist)"	>> $rcfile.reinit
        echo "set gmids = ($rep_gmids)"		>> $rcfile.reinit
        echo "set ncpts = ($rep_ncpts)"		>> $rcfile.reinit

        set cmd = trac-preproc
        set cmd = ($cmd -c $rcfile.reinit)
        set cmd = ($cmd -log $LF)
        set cmd = ($cmd -cmd $CF)
        set cmd = ($cmd -no-isrunning)
        echo $cmd
        $cmd
        if ($status) exit 1
      endif

      set cmd = $trcdir/dmri_paths
      set cmd = ($cmd --outdir $outlist)
      set cmd = ($cmd --dwi  $dwidir/data.nii.gz)
      set cmd = ($cmd --grad $dwidir/bvecs)
      set cmd = ($cmd --bval $dwidir/bvals)
      set cmd = ($cmd --mask $brainmask)
      set cmd = ($cmd --bpdir $dtdir/dmri.bedpostX)
      set cmd = ($cmd --ntr  $nstick)
      set cmd = ($cmd --fmin $fmin)
      set cmd = ($cmd --roi1 $roi1list)
      set cmd = ($cmd --roi2 $roi2list)
      if ($xspace == mni || $xspace == rob) then
        set cmd = ($cmd --reg $xfmdir/diff2$xspace.$reg.lta)
      else if ($xspace == cvs) then
        set cmd = ($cmd --reg $xfmdir/diff2anatorig.$reg.lta)
        if ($subj != $cvstemp) then
          set cmd = ($cmd --regnl $xfmdir/cvs/$cvswarp.m3z)
        endif
      else if ($xspace == syn) then
        set cmd = ($cmd --reg   $xfmdir/diff2syn.lta)
        set cmd = ($cmd --regnl $xfmdir/syn_warp.m3z)
      else if ($xspace == fnt) then
        set cmd = ($cmd --regnl $xfmdir/diff2fsl_warp.m3z)
      endif
      set cmd = ($cmd --init $initlist)
      if ($dosegprior) then
        set cmd = ($cmd --nprior $npriorlist)
        set cmd = ($cmd --lprior $lpriorlist)
        if ($usethalnuc) then
          set cmd = ($cmd --seg  $labdir/$xspace/${segname}+thalnuc.nii.gz)
        else
          set cmd = ($cmd --seg  $labdir/$xspace/$segname.nii.gz)
        endif
      endif
      if ($dotangprior) then
        set cmd = ($cmd --tprior $tpriorlist)
        set cmd = ($cmd --cprior $cpriorlist)
      endif
      if ($doxyzprior) then
        set cmd = ($cmd --prior $priorlist)
      endif
      set cmd = ($cmd --nb $nburnin --ns $nsample --nu $nupdate --nk $nkeep)
      if ($doinitprop) then
        set cmd = ($cmd --sdp $sdplist)
      endif
      if ($debug) then
        set cmd = ($cmd --debug)
      endif
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $fs_time $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif

      # Check if any paths need to be reinitialized
      set rep_pathlist = ()
      set rep_gmids = ()
      set rep_ncpts = ()

      @ ipath = 1
      while ($ipath <= $#pathlist)
        set pathname = $pathlist[$ipath]
        set pathdir = $outdir/${pathname}_$avgmode

        if (-e $pathdir/path.pd.nii.gz) then
          set pmax = `fslstats $pathdir/path.pd.nii.gz -R | awk '{print $2}'`
          set pthr = `echo "$pmax * .2" | bc -l`

          set cmd = fslmaths
          set cmd = ($cmd $pathdir/path.pd.nii.gz -thr $pthr)
          set cmd = ($cmd $pathdir/tmp.nii.gz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $fs_time $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          if (`fslstats $pathdir/tmp.nii.gz -V | awk '{print $2}'` \
           == `fslstats $pathdir/path.map.nii.gz -V | awk '{print $2}'`) then
            set rep_pathlist = ($rep_pathlist $pathname)
            set rep_gmids    = ($rep_gmids    $gmids[$ipath])
            set rep_ncpts    = ($rep_ncpts    $ncpts[$ipath])
          else
            set pathpre = ${pathname}_$avgmode

            set outlist    = `printf '%s\n' $outlist    | grep -v $pathpre`
            set roi1list   = `printf '%s\n' $roi1list   | grep -v $pathpre`
            set roi2list   = `printf '%s\n' $roi2list   | grep -v $pathpre`
            set initlist   = `printf '%s\n' $initlist   | grep -v $pathpre`
            set sdplist    = `printf '%s\n' $sdplist    | grep -v $pathpre`
            set priorlist  = `printf '%s\n' $priorlist  | grep -v $pathpre`
            set npriorlist = `printf '%s\n' $npriorlist | grep -v $pathpre`
            set lpriorlist = `printf '%s\n' $lpriorlist | grep -v $pathpre`
            set tpriorlist = `printf '%s\n' $tpriorlist | grep -v $pathpre`
            set cpriorlist = `printf '%s\n' $cpriorlist | grep -v $pathpre`
          endif

          set cmd = (rm -rf $pathdir/tmp.nii.gz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif
        endif

        @ ipath = $ipath + 1
      end

      if ($#rep_pathlist == 0) break;

      @ irep = $irep + 1
    end

    # Create reference paths for along-tract analysis
    set refdir = `dirname $trainsubjlist[1]`/$xspace
    if (-d $refdir) then
      set srclist = `printf "$refdir/%s.mean.trk " $pathlist`
      set trglist = `printf "$outdir/%s_$avgmode/path.ref.txt " $pathlist`

      # Map mean paths from template to individual
      set cmd = $trcdir/dmri_trk2trk
      set cmd = ($cmd --in     $srclist)
      set cmd = ($cmd --inref  $intertrg)
      if ($xspace == mni || $xspace == rob) then
        set cmd = ($cmd --reg $xfmdir/diff2$xspace.$reg.lta)
      else if ($xspace == cvs) then
        set cmd = ($cmd --reg $xfmdir/diff2anatorig.$reg.lta)
      if ($subj != $cvstemp) then
          set cmd = ($cmd --regnl $xfmdir/cvs/$cvswarp.m3z)
        endif
      else if ($xspace == syn) then
        set cmd = ($cmd --reg   $xfmdir/diff2syn.lta)
        set cmd = ($cmd --regnl $xfmdir/syn_warp.m3z)
      else if ($xspace == fnt) then
        set cmd = ($cmd --regnl $xfmdir/diff2fsl_warp.m3z)
      endif
      set cmd = ($cmd --inv)
      set cmd = ($cmd --smooth)
      set cmd = ($cmd --outasc $trglist)
      set cmd = ($cmd --outref $dwidir/dtifit_FA.nii.gz)
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $fs_time $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif
    endif

    foreach pathname ($pathlist)
      set pathpre = ${pathname}_$avgmode

      if (! -e $outdir/$pathpre/path.map.txt)	continue

      set pname = `echo $pathname | awk -v FS=_ '{print $1}'`

      # Get whole-path and along-path statistics
      set cmd = $trcdir/dmri_pathstats
      set cmd = ($cmd --intrc  $outdir/$pathpre)
      set cmd = ($cmd --dtbase $dwidir/dtifit)
      set cmd = ($cmd --path   $pname)
      set cmd = ($cmd --subj   $subj)
      if (-e $outdir/$pathpre/path.ref.txt) then
       	set cmd = ($cmd --invox path.ref.txt)
      endif
      set cmd = ($cmd --out    $outdir/$pathpre/pathstats.overall.txt)
      set cmd = ($cmd --outvox $outdir/$pathpre/pathstats.byvoxel.txt)
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $fs_time $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif

      # Project path endpoints onto cortical surface
      # Skip cerebellum pathways
      if ($pname == mcp)  continue

      # Is it an inter-hemispheric pathway?
      if ($pname =~ cc.* || $pname == acomm || $pname == mcp) then
        set isinter = 1
      else
        set isinter = 0
      endif

      foreach pt (1 2)
        # Skip non-cortical endpoints
        if ( ($pname =~ *.ar  && $pt == 1) || \
             ($pname =~ *.atr && $pt == 2) || \
             ($pname =~ *.cst && $pt == 2) || \
             ($pname =~ *.fx  && $pt == 1) || \
             ($pname =~ *.or  && $pt == 1) )	continue

        if ( $pname =~ lh.* || ($isinter == 1 && $pt == 2) ) then
          set hemi = lh
        else
          set hemi = rh
        endif

        if ($isinter) then
          set labelname = $hemi.$pname
        else
          set labelname = $pname
        endif

        # Project from volume in DWI space to surface in structural space
        set cmd = mri_vol2surf
        set cmd = ($cmd --mov $outdir/$pathpre/endpt$pt.pd.nii.gz)
        set cmd = ($cmd --reg $xfmdir/diff2anatorig.$reg.lta)
        set cmd = ($cmd --srcsubject $subj)
        set cmd = ($cmd --hemi $hemi)
        set cmd = ($cmd --projdist-avg $projmin $projmax $dproj)
        set cmd = ($cmd --trgsubject $subj)
        set cmd = ($cmd --o $outdir/$pathpre/endpt$pt.surf.mgz)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        set nproj = `echo "($projmax - ($projmin))/$dproj + 1" | bc -l`

        set cmd = fscalc
        set cmd = ($cmd     $outdir/$pathpre/endpt$pt.surf.mgz)
        set cmd = ($cmd mul $nproj)
        set cmd = ($cmd --o $outdir/$pathpre/endpt$pt.surf.mgz)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        set cmd = mri_binarize
        set cmd = ($cmd --min 1)
        set cmd = ($cmd --i $outdir/$pathpre/endpt$pt.surf.mgz)
        set cmd = ($cmd --o $outdir/$pathpre/endpt$pt.surf.bin.mgz)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        set nvert = `$cmd | grep "in range" | awk '{print $2}'`
        if ($nvert == 0)	continue

        set cmd = mri_cor2label
        set cmd = ($cmd --i $outdir/$pathpre/endpt$pt.surf.bin.mgz)
        set cmd = ($cmd --id 1)
        set cmd = ($cmd --l $outdir/$pathpre/endpt$pt.surf.label)
        set cmd = ($cmd --surf $subj $hemi white)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        set cmd = mri_label2label
        set cmd = ($cmd --srclabel $outdir/$pathpre/endpt$pt.surf.label)
        set cmd = ($cmd --s $subj --hemi $hemi --regmethod surface)
        set cmd = ($cmd --dilate 5 --erode 5)
        set cmd = ($cmd --trglabel $outdir/$pathpre/endpt$pt.surf.label)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        set cmd = mris_anatomical_stats
        set cmd = ($cmd -l $outdir/$pathpre/endpt$pt.surf.label)
        set cmd = ($cmd -f $outdir/$pathpre/endpt$pt.surf.stats)
        set cmd = ($cmd $subj $hemi)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif
      end
    end

    # Merge all existing pathways into a 4D volume
    set mergelist = ()
    foreach dname ($outdir/*_$avgmode)
      if (-e $dname/path.pd.nii.gz) then
        set mergelist = ($mergelist `basename $dname`/path.pd.nii.gz)
      endif
    end

    if ($#mergelist) then
      set cmd = $trcdir/dmri_mergepaths
      set cmd = ($cmd --indir $outdir)
      set cmd = ($cmd --in    $mergelist)
      set cmd = ($cmd --out   $outdir/merged_$avgmode.mgz)
      set cmd = ($cmd --ctab  $FREESURFER_HOME/FreeSurferColorLUT.txt)
      set cmd = ($cmd --thresh $pmin)
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif
    endif
  end
else
  echo "#-------------------------------------" |& tee -a $LF
  echo "#@# Path reconstruction (longitudinal) `date`" |& tee -a $LF

  foreach ntrain ($ntrainlist)
    set avgmode = $avgname${ntrain}_${xspace}_$reg

    set outdir = dpath

    if ($dopathsubdirs) then
      set outdir = $outdir/${nsample}samp

      if ($dosegprior && $dotangprior) then
        set outdir = $outdir/all
      else if ($dosegprior) then
        set outdir = $outdir/seg14
      else if ($doxyzprior) then
        set outdir = $outdir/xyz
      else if ($dotangprior) then
        set outdir = $outdir/tang
      else
        set outdir = $outdir/none
      endif
    endif

    set outlist = `printf "$outdir/%s_$avgmode " $pathlist`

    set roi1list = `printf \
      "$labdir/$xspace/%s_${avgmode}_end1_dil.nii.gz " $pathlist`
    set roi2list = `printf \
      "$labdir/$xspace/%s_${avgmode}_end2_dil.nii.gz " $pathlist`

    set initlist = ()
    set sdplist = ()
    @ ipath = 1
    while ($ipath <= $#pathlist)
      set pname = $pathlist[$ipath]
      set nc    = $ncpts[$ipath]

      set initlist = ($initlist \
        $labdir/anatorig/${pname}_${avgmode}_cpts_$nc.txt)
      set sdplist = ($sdplist \
        $labdir/anatorig/${pname}_${avgmode}_cpts_${nc}_std.txt)

      @ ipath = $ipath + 1
    end

    set priorlist = ()
    foreach pathname ($pathlist)
      set priorlist = ($priorlist \
        $labdir/$xspace/${pathname}_${avgmode}_logprior_str_0.nii.gz \
        $labdir/$xspace/${pathname}_${avgmode}_logprior_str_1.nii.gz)
    end

    set npriorlist = ()
    foreach pathname ($pathlist)
      set npriorlist = ($npriorlist \
        $labdir/$xspace/${pathname}_${avgmode}_fsnnprior \
        $labdir/$xspace/${pathname}_${avgmode}_fsnnids)
    end

    if ($usetrunc) then
      set npriorlist = `printf "%s_all " $npriorlist`
    endif

    set lpriorlist = ()
    foreach pathname ($pathlist)
      set lpriorlist = ($lpriorlist \
        $labdir/$xspace/${pathname}_${avgmode}_fsprior \
        $labdir/$xspace/${pathname}_${avgmode}_fsids)
    end

    if ($usetrunc) then
      set lpriorlist = `printf "%s_all " $lpriorlist`
    endif

    set tpriorlist = `printf \
      "$labdir/$xspace/%s_${avgmode}_tangprior.txt " $pathlist`
    set cpriorlist = `printf \
      "$labdir/$xspace/%s_${avgmode}_curvprior.txt " $pathlist`

    if ($usetrunc) then
      set tpriorlist = `echo $tpriorlist | sed 's/\.txt/_all\.txt/'`
      set cpriorlist = `echo $cpriorlist | sed 's/\.txt/_all\.txt/'`
    endif

    foreach subj_t ($tplist)
      if ($overwrite) then	# Clean up pre-existing directories
        foreach dname (`printf "$dtroot/$subj_t/%s " $outlist`)
          set cmd = (rm -rf $dname)
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif
        end
      else			# Rename pre-existing directories
        foreach dname (`printf "$dtroot/$subj_t/%s " $outlist`)
          if (-e $dname) then
            @ vno = 0
            while (-e $dname.v$vno)
              @ vno = $vno + 1
            end
            set cmd = (mv $dname $dname.v$vno)
            if ($RunIt) then
              $cmd |& tee -a $LF
              if ($status) goto error_exit
            endif
          endif
        end
      endif
    end

    if ($usemaskanat) then
      set brainmask = ${segname}_mask.$reg
      set basemask  = ${segname}_mask
    else
      set brainmask = lowb_brain_mask
      set basemask  = lowb_brain_mask.$reg
    endif

    set inlist = `printf "$dtroot/%s " $tplist`

    if ($usethalnuc) then
      set seglist = `printf \
        "$dtroot/%s/dlabel/$xspace/${segname}+thalnuc.nii.gz " $tplist`
    else
      set seglist = `printf \
        "$dtroot/%s/dlabel/$xspace/$segname.nii.gz " $tplist`
    endif

    @ irep = 0

    while ($irep <= 5)
      if ($irep) then
        cp $rcfile $rcfile.reinit
        echo "set dopriors = 1"                     >> $rcfile.reinit
        echo "set reinit = 1"                       >> $rcfile.reinit
        echo "set pathlist = ($rep_pathlist)"       >> $rcfile.reinit
        echo "set gmids = ($rep_gmids)"             >> $rcfile.reinit
        echo "set ncpts = ($rep_ncpts)"             >> $rcfile.reinit

        set cmd = trac-preproc
        set cmd = ($cmd -c $rcfile.reinit)
        set cmd = ($cmd -log $LF)
        set cmd = ($cmd -cmd $CF)
        set cmd = ($cmd -no-isrunning)
        echo $cmd
        $cmd
        if ($status) exit 1
      endif

      set cmd = $trcdir/dmri_paths
      set cmd = ($cmd --indir     $inlist)
      set cmd = ($cmd --outdir    $outlist)
      set cmd = ($cmd --dwi       dmri/data.nii.gz)
      set cmd = ($cmd --grad      dmri/bvecs)
      set cmd = ($cmd --bval      dmri/bvals)
      set cmd = ($cmd --mask      dlabel/diff/$brainmask.nii.gz)
      set cmd = ($cmd --bpdir     dmri.bedpostX)
      set cmd = ($cmd --ntr       $nstick)
      set cmd = ($cmd --fmin      $fmin)
      set cmd = ($cmd --basereg   dmri/xfms/anatorig2diff.$reg.lta)
      set cmd = ($cmd --basemask  $labdir/anatorig/$basemask.nii.gz)
      set cmd = ($cmd --roi1      $roi1list)
      set cmd = ($cmd --roi2      $roi2list)
      if ($xspace == mni || $xspace == rob) then
        set cmd = ($cmd --reg     $xfmdir/anatorig2$xspace.lta)
      else if ($xspace == cvs) then
        # Hack: Using one time point for now!
        if ($tplist[1] != $cvstemp) then
          set xfmdir_t = $dtroot/$tplist[1]/dmri/xfms
          set cmd = ($cmd --regnl $xfmdir_t/cvs/$cvswarp.m3z)
        endif
      else if ($xspace == syn) then
        set cmd = ($cmd --reg   $xfmdir/anatorig2syn.$reg.lta)
        set cmd = ($cmd --regnl $xfmdir/syn_warp.m3z)
      else if ($xspace == fnt) then
        set cmd = ($cmd --reg   $xfmdir/anatorig2diff.$reg.lta)
        set cmd = ($cmd --regnl $xfmdir/diff2fsl_warp.m3z)
      endif
      set cmd = ($cmd --init      $initlist)
      if ($dosegprior) then
        set cmd = ($cmd --nprior  $npriorlist)
        set cmd = ($cmd --lprior  $lpriorlist)
        set cmd = ($cmd --seg     $seglist)
      else if ($doxyzprior) then
        set cmd = ($cmd --prior   $priorlist)
      else if ($dotangprior) then
        set cmd = ($cmd --tprior  $tpriorlist)
        set cmd = ($cmd --cprior  $cpriorlist)
      endif
      set cmd = ($cmd --nb $nburnin --ns $nsample --nu $nupdate --nk $nkeep)
      if ($doinitprop) then
        set cmd = ($cmd --sdp     $sdplist)
      endif
      if ($debug) then
        set cmd = ($cmd --debug)
      endif
      echo $cmd |& tee -a $LF |& tee -a $CF
      if ($RunIt) then
        $fs_time $cmd |& tee -a $LF
        if ($status) goto error_exit
      endif

      # Check if any paths need to be reinitialized
      set rep_pathlist = ()
      set rep_gmids = ()
      set rep_ncpts = ()

      set outdir_t = $dtroot/$tplist[1]/$outdir

      @ ipath = 1
      while ($ipath <= $#pathlist)
        set pathname = $pathlist[$ipath]
        set pathdir = $outdir_t/${pathname}_$avgmode

        if (-e $pathdir/path.pd.nii.gz) then
          set pmax = `fslstats $pathdir/path.pd.nii.gz -R | awk '{print $2}'`
          set pthr = `echo "$pmax * .2" | bc -l`

          set cmd = fslmaths
          set cmd = ($cmd $pathdir/path.pd.nii.gz -thr $pthr)
          set cmd = ($cmd $pathdir/tmp.nii.gz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $fs_time $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          if (`fslstats $pathdir/tmp.nii.gz -V | awk '{print $2}'` \
           == `fslstats $pathdir/path.map.nii.gz -V | awk '{print $2}'`) then
            set rep_pathlist = ($rep_pathlist $pathname)
            set rep_gmids    = ($rep_gmids    $gmids[$ipath])
            set rep_ncpts    = ($rep_ncpts    $ncpts[$ipath])
          else
            set pathpre = ${pathname}_$avgmode

            set outlist    = `printf '%s\n' $outlist    | grep -v $pathpre`
            set roi1list   = `printf '%s\n' $roi1list   | grep -v $pathpre`
            set roi2list   = `printf '%s\n' $roi2list   | grep -v $pathpre`
            set initlist   = `printf '%s\n' $initlist   | grep -v $pathpre`
            set sdplist    = `printf '%s\n' $sdplist    | grep -v $pathpre`
            set priorlist  = `printf '%s\n' $priorlist  | grep -v $pathpre`
            set npriorlist = `printf '%s\n' $npriorlist | grep -v $pathpre`
            set lpriorlist = `printf '%s\n' $lpriorlist | grep -v $pathpre`
            set tpriorlist = `printf '%s\n' $tpriorlist | grep -v $pathpre`
            set cpriorlist = `printf '%s\n' $cpriorlist | grep -v $pathpre`
          endif

          set cmd = (rm -rf $pathdir/tmp.nii.gz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif
        endif

        @ ipath = $ipath + 1
      end

      if ($#rep_pathlist == 0) break;

      @ irep = $irep + 1
    end

    foreach subj_t ($tplist)
      set dwidir_t = $dtroot/$subj_t/dmri
      set xfmdir_t = $dwidir_t/xfms
      set outdir_t = $dtroot/$subj_t/$outdir

      # Create reference paths for along-tract analysis
      set refdir = `dirname $trainsubjlist[1]`/$xspace
      if (-d $refdir) then
        set srclist = `printf "$refdir/%s.mean.trk " $pathlist`
        set trglist = `printf "$outdir_t/%s_$avgmode/path.ref.txt " $pathlist`

        # Map mean paths from template to individual
        set cmd = $trcdir/dmri_trk2trk
        set cmd = ($cmd --in     $srclist)
        set cmd = ($cmd --inref  $intertrg)
        if ($xspace == mni || $xspace == rob) then
          set cmd = ($cmd --reg $xfmdir_t/diff2$xspace.$reg.lta)
        else if ($xspace == cvs) then
          set cmd = ($cmd --reg $xfmdir_t/diff2anatorig.$reg.lta)
        if ($subj != $cvstemp) then
            set cmd = ($cmd --regnl $xfmdir_t/cvs/$cvswarp.m3z)
          endif
        else if ($xspace == syn) then
          set cmd = ($cmd --reg   $xfmdir_t/diff2syn.lta)
          set cmd = ($cmd --regnl $xfmdir_t/syn_warp.m3z)
        else if ($xspace == fnt) then
          set cmd = ($cmd --regnl $xfmdir_t/diff2fsl_warp.m3z)
        endif
        set cmd = ($cmd --inv)
        set cmd = ($cmd --smooth)
        set cmd = ($cmd --outasc $trglist)
        set cmd = ($cmd --outref $dwidir_t/dtifit_FA.nii.gz)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $fs_time $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif
      endif

      foreach pathname ($pathlist)
        set pathdir = $outdir_t/${pathname}_$avgmode

        if (! -e $pathdir/path.map.txt)	continue

        set pname = `echo $pathname | awk -v FS=_ '{print $1}'`

        # Get whole-path and along-path statistics
        set cmd = $trcdir/dmri_pathstats
        set cmd = ($cmd --intrc  $pathdir)
        set cmd = ($cmd --dtbase $dwidir_t/dtifit)
        set cmd = ($cmd --path   $pname)
        set cmd = ($cmd --subj   $subj_t)
        if (-e $pathdir/path.ref.txt) then
       	  set cmd = ($cmd --invox path.ref.txt)
        endif
        set cmd = ($cmd --out    $pathdir/pathstats.overall.txt)
        set cmd = ($cmd --outvox $pathdir/pathstats.byvoxel.txt)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $fs_time $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif

        # Project path endpoints onto cortical surface
        # Skip cerebellum pathways
        if ($pname == mcp)  continue

        # Is it an inter-hemispheric pathway?
        if ($pname =~ cc.* || $pname == acomm || $pname == mcp) then
          set isinter = 1
        else
          set isinter = 0
        endif

        foreach pt (1 2)
          # Skip non-cortical endpoints
          if ( ($pname =~ *.ar  && $pt == 1) || \
               ($pname =~ *.atr && $pt == 2) || \
               ($pname =~ *.cst && $pt == 2) || \
               ($pname =~ *.fx  && $pt == 1) || \
               ($pname =~ *.or  && $pt == 1) )	continue

          if ( $pname =~ lh.* || ($isinter == 1 && $pt == 2) ) then
            set hemi = lh
          else
            set hemi = rh
          endif

          if ($isinter) then
            set labelname = $hemi.$pname
          else
            set labelname = $pname
          endif

          # Project from volume in DWI space to surface in structural space
          set cmd = mri_vol2surf
          set cmd = ($cmd --mov $pathdir/endpt$pt.pd.nii.gz)
          set cmd = ($cmd --reg $xfmdir_t/diff2anatorig.$reg.lta)
          set cmd = ($cmd --srcsubject $subj_t)
          set cmd = ($cmd --hemi $hemi)
          set cmd = ($cmd --projdist-avg $projmin $projmax $dproj)
          set cmd = ($cmd --trgsubject $subj_t)
          set cmd = ($cmd --o $pathdir/endpt$pt.surf.mgz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          set nproj = `echo "($projmax - ($projmin))/$dproj + 1" | bc -l`

          set cmd = fscalc
          set cmd = ($cmd     $pathdir/endpt$pt.surf.mgz)
          set cmd = ($cmd mul $nproj)
          set cmd = ($cmd --o $pathdir/endpt$pt.surf.mgz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          set cmd = mri_binarize
          set cmd = ($cmd --min 1)
          set cmd = ($cmd --i $pathdir/endpt$pt.surf.mgz)
          set cmd = ($cmd --o $pathdir/endpt$pt.surf.bin.mgz)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          set nvert = `$cmd | grep "in range" | awk '{print $2}'`
          if ($nvert == 0)	continue

          set cmd = mri_cor2label
          set cmd = ($cmd --i $pathdir/endpt$pt.surf.bin.mgz)
          set cmd = ($cmd --id 1)
          set cmd = ($cmd --l $pathdir/endpt$pt.surf.label)
          set cmd = ($cmd --surf $subj_t $hemi white)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          set cmd = mri_label2label
          set cmd = ($cmd --srclabel $pathdir/endpt$pt.surf.label)
          set cmd = ($cmd --s $subj_t --hemi $hemi --regmethod surface)
          set cmd = ($cmd --dilate 5 --erode 5)
          set cmd = ($cmd --trglabel $pathdir/endpt$pt.surf.label)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif

          set cmd = mris_anatomical_stats 
          set cmd = ($cmd -l $pathdir/endpt$pt.surf.label)
          set cmd = ($cmd -f $pathdir/endpt$pt.surf.stats)
          set cmd = ($cmd $subj $hemi)
          echo $cmd |& tee -a $LF |& tee -a $CF
          if ($RunIt) then
            $cmd |& tee -a $LF
            if ($status) goto error_exit
          endif
        end
      end

      # Merge all existing pathways into a 4D volume
      set mergelist = ()
      foreach dname ($outdir_t/*_$avgmode)
        if (-e $dname/path.pd.nii.gz) then
          set mergelist = ($mergelist `basename $dname`/path.pd.nii.gz)
        endif
      end

      if ($#mergelist) then
        set cmd = $trcdir/dmri_mergepaths
        set cmd = ($cmd --indir $outdir_t)
        set cmd = ($cmd --in    $mergelist)
        set cmd = ($cmd --out   $outdir_t/merged_$avgmode.mgz)
        set cmd = ($cmd --ctab  $FREESURFER_HOME/FreeSurferColorLUT.txt)
        set cmd = ($cmd --thresh $pmin)
        echo $cmd |& tee -a $LF |& tee -a $CF
        if ($RunIt) then
          $cmd |& tee -a $LF
          if ($status) goto error_exit
        endif
      endif
    end
  end
endif

# Remove the error file
if ($DoIsRunning && $#IsRunningFile) rm -f $IsRunningFile

# Create the done file
echo "------------------------------" > $DoneFile
echo "SUBJECT $subj" >> $DoneFile
echo "DATE `date`"     >> $DoneFile
echo "USER `whoami`"      >> $DoneFile
echo "HOST `hostname`" >> $DoneFile
echo "PROCESSOR `uname -m`" >> $DoneFile
echo "OS `uname -s`"       >> $DoneFile
uname -a         >> $DoneFile
echo $VERSION    >> $DoneFile
echo $0          >> $DoneFile

echo "#-------------------------------------" |& tee -a $LF
echo "$ProgName finished without error at `date`" |& tee -a $LF

exit 0
#############------------------------------------#######################
##################>>>>>>>>>>>>>.<<<<<<<<<<<<<<<<<#######################
#############------------------------------------#######################

############--------------##################
error_exit:
  uname -a | tee -a $LF
  echo "" |& tee -a $LF
  echo "$ProgName exited with ERRORS at `date`" |& tee -a $LF
  echo "" |& tee -a $LF

  # Remove IsRunningFile
  if($DoIsRunning && $#IsRunningFile) rm -f $IsRunningFile

  # Create an error file with date, cmd, etc of error
  if ($?ErrorFile) then
    echo "------------------------------" > $ErrorFile
    echo "SUBJECT $subj" >> $ErrorFile
    echo "DATE `date`"     >> $ErrorFile
    echo "USER `whoami`"      >> $ErrorFile
    echo "HOST `hostname`" >> $ErrorFile
    echo "PROCESSOR `uname -m`" >> $ErrorFile
    echo "OS `uname -s`"       >> $ErrorFile
    uname -a         >> $ErrorFile
    echo $VERSION    >> $ErrorFile
    echo $0          >> $ErrorFile
    echo "PWD `pwd`" >> $ErrorFile
    echo "CMD $cmd"  >> $ErrorFile
  endif

  # Finally exit
  exit 1

############--------------##################
parse_args:
set cmdline = ($argv)

while( $#argv != 0 )
  set flag = $argv[1]; shift;

  if ("$flag" == ";") break;

  switch($flag)
    case "-c":
      if ( $#argv < 1) goto arg1err;
      set rcfile = "$argv[1]"; shift;
      if (! -e "$rcfile") then
        echo "ERROR: cannot find $rcfile"
        goto error_exit
      endif
      if (! -r "$rcfile") then
        echo "ERROR: $rcfile exists but is not readable"
        goto error_exit
      endif
      breaksw

    case "-time":
      set DoTime = 1
      breaksw

    case "-notime":
      set DoTime = 0
      breaksw

    case "-log":
      if ( $#argv < 1) goto arg1err;
      set LF = $argv[1]; shift;
      breaksw

    case "-nolog":
      set LF = /dev/null
      breaksw

    case "-cmd":
      if ( $#argv < 1) goto arg1err;
      set CF = $argv[1]; shift;
      breaksw

    case "-nocmd":
      set CF = /dev/null
      breaksw

    case "-no-isrunning":
      set DoIsRunning = 0
      breaksw

    case "-umask":
      if ( $#argv < 1) goto arg1err;
      umask $1; shift;
      breaksw

    case "-grp":
      if ( $#argv < 1) goto arg1err;
      set grp = $argv[1];
      set curgrp = `id -gn`;
      if($grp != $curgrp) then
        echo "ERROR: current group $curgrp and specified group $grp differ"
        goto error_exit;
      endif
      breaksw

    case "-allowcoredump":
      limit coredumpsize unlimited
      breaksw

    case "-debug":
      set debug = 1
      breaksw

    case "-dontrun":
      set RunIt = 0
      breaksw

    case "-onlyversions":
      set DoVersionsOnly = 1
      breaksw

    default:
      echo "ERROR: flag $flag unrecognized"
      echo $cmdline
      goto error_exit
      breaksw
  endsw
end

goto parse_args_return;

############--------------##################
arg1err:
  echo "ERROR: flag $flag requires one argument"
  goto error_exit

############--------------##################
check_params:
  if (! $#rcfile) then
    echo "ERROR: a configuration (dmrirc) file must be specified"
    goto error_exit
  endif

  if (! $#LF) set LF = `dirname $rcfile`/trac-all.log
  if (! $#CF) set CF = `dirname $rcfile`/trac-all.cmd

goto check_params_return;

############--------------##################
usage_exit:
  echo ""
  echo "USAGE: $ProgName"
  echo ""
  echo "Required arguments:"
  echo "  -c <file>      : dmrirc file (see dmrirc.example)"
  echo ""
  echo "Other arguments:"
  echo "  -log <file>    : default is trac-all.log in the same dir as dmrirc"
  echo "  -nolog         : do not save a log file"
  echo "  -cmd <file>    : default is trac-all.cmd in the same dir as dmrirc"
  echo "  -nocmd         : do not save a cmd file"
  echo "  -no-isrunning  : do not check whether this subject is currently being processed"
  echo "  -umask umask   : set unix file permission mask (default 002)"
  echo "  -grp groupid   : check that current group is alpha groupid "
  echo "  -allowcoredump : set coredump limit to unlimited"
  echo "  -debug         : generate much more output"
  echo "  -dontrun       : do everything but execute each command"
  echo "  -version       : print version of this script and exit"
  echo "  -help          : print full contents of help"
  echo ""

  if(! $PrintHelp) exit 1

  echo $VERSION
  echo ""

  cat $0 | awk 'BEGIN{prt=0}{if(prt) print $0;if($1 == "BEGINHELP") prt=1}'

  exit 1

#---- Everything below here is printed out as part of help -----#
BEGINHELP

Tractography for a single subject.

This script is called by trac-all. Trac-all makes sure that a proper
configuration file is written locally (scripts/dmrirc.local) and passed
as an argument to this script.

SEE ALSO: trac-all, dmrirc.example