#!/bin/tcsh -f
# midefacer - sources
if(-e $FREESURFER_HOME/sources.csh) then
  source $FREESURFER_HOME/sources.csh
endif

set VERSION = '$Id$';
set scriptname = `basename $0`

set fmt = mgz
set invol = ();
set outvol = ();
set facemask = ()
set ForceUpdate = 0
set DoPostHeadSurf = 0
set threads = 1
set xmask = ()
set MaskSegToHead = 1;
set DoSamsegXMask = 1;
set nSegDil = 1;
set SamsegFast = 1
set XEyeBalls = 1; #Samseg only
set DoSynthSegXMask = 0;
set DoBBR = 0
set BBRDoF = 9
set defacingsubject = defacing-atlas
set outdir = ()
set FillType = 1;
set FillConstIn = ()
set FillConstOut = ()
set fhi = ()
set DoEars = 1;
set DoForehead = 0;
set DoBackOfHead = 0
set DoPics = 0
set SubjCode = ()
set PicAnnot = ()
set DoRipple = 0
set RippleAmp = 1.5
set RipplePeriod = 30
set dwatermark = 5;
set atlasdir = ()
set imconvert = /usr/bin/convert
set XOptsFile = ()
set SamsegJson = ()
set UseOrMask = 0
set InitReg = ()
set tmpdir = ();
set cleanup = 1;
set statframe = 0
set LF = ();

set UseXvfb = 1
which Xvfb > /dev/null
if($status) set UseXvfb = 0

set inputargs = ($argv);
set PrintHelp = 0;
if($#argv == 0) goto usage_exit;
set n = `echo $argv | grep -e -help | wc -l` 
if($n != 0) then
  set PrintHelp = 1;
  goto usage_exit;
endif
set n = `echo $argv | grep -e -version | wc -l` 
if($n != 0) then
  echo $VERSION
  exit 0;
endif
goto parse_args;
parse_args_return:
goto check_params;
check_params_return:

set ctab = $atlasdir/deface.ctab
set templatesurf = $atlasdir/average.head.dec44.us7.remesh3.ripple.surf
set templatevol = $atlasdir/average.head.mgz
set atlas2samseglta = $atlasdir/reg.samseg-default.lta
set templateface = $atlasdir/average.face.dec44.us7.remesh3.label
set upperfacenoeyes = $atlasdir/average.face.upper.no-eyes.dec44.us7.remesh3.label
set upperfacenoeyesnose = $atlasdir/average.face.upper.no-eyes.no-nose.dec44.us7.remesh3.label
set nose = $atlasdir/average.face.nose.dec44.us7.remesh3.label
set lowerface = $atlasdir/average.face.lower.dec44.us7.remesh3.label
set watermark = $atlasdir/watermark.dec44.us7.remesh3.label
set watermark = $atlasdir/mideface.watermark.dec44.us7.remesh3.label
set earleft  = $atlasdir/average.face.ear.left.dec44.us7.remesh3.label
set earright = $atlasdir/average.face.ear.right.dec44.us7.remesh3.label
set eyes  = $atlasdir/average.face.eyes.dec44.us7.remesh3.label
set forehead  = $atlasdir/average.face.forehead.dec44.us7.remesh3.label
set backofhead  = $atlasdir/average.back-of-head.dec44.us7.remesh3.label
set bbrface = $atlasdir/average.face.upper.dec44.us7.remesh3.label

set labellist = ($lowerface $upperfacenoeyesnose $nose $eyes)

set StartTime = `date`;
set tSecStart = `date '+%s'`;
set year  = `date +%Y`
set month = `date +%m`
set day   = `date +%d`
set hour   = `date +%H`
set min    = `date +%M`

if($#outdir == 0) set outdir = `dirname $outvol`
mkdir -p $outdir
pushd $outdir > /dev/null
set outdir = `pwd`;
popd > /dev/null

if($#tmpdir == 0) then
  if(-dw /scratch)   set tmpdir = /scratch/tmpdir.mideface.$$
  if(! -dw /scratch) set tmpdir = $outdir/tmpdir.mideface.$$
endif
mkdir -p $tmpdir/log
set invol  = `getfullpath $invol`
set outvol = `getfullpath $outvol`
set tmpdir = `getfullpath $tmpdir`

# Set up log file
set stem = `fname2stem $outvol`
if($#LF == 0) set LF = $stem.mideface.log
if($LF != /dev/null) rm -f $LF
echo "Log file for mideface" >> $LF
date  | tee -a $LF
echo "" | tee -a $LF
echo "setenv SUBJECTS_DIR $SUBJECTS_DIR" | tee -a $LF
echo "cd `pwd`"  | tee -a $LF
echo $0 $inputargs | tee -a $LF
ls -l $0  | tee -a $LF
echo "" | tee -a $LF
cat $FREESURFER_HOME/build-stamp.txt | tee -a $LF
echo $VERSION | tee -a $LF
uname -a  | tee -a $LF
echo "pid $$" | tee -a $LF
if($?PBS_JOBID) then
  echo "pbsjob $PBS_JOBID"  >> $LF
endif
if($?SLURM_JOB_ID) then
  echo SLURM_JOB_ID $SLURM_JOB_ID >> $LF
endif

echo "UseXvfb $UseXvfb" | tee -a $LF
set runxvfb = fsxvfb
if(! $UseXvfb) set runxvfb = ()

if($#XOptsFile) cp $XOptsFile $tmpdir
if($#InitReg)   cp $InitReg   $tmpdir

#========================================================
# Register to samseg atlas, maybe run samseg for xmask
set samsegdir = $tmpdir/samseg
set samseglta = $samsegdir/template.lta
if($DoSamsegXMask == 0) then
  # Run samseg but only to get the reg
  set ud = (`UpdateNeeded $samseglta $invol $InitReg`)
  if($ud || $ForceUpdate) then
    set cmd = (run_samseg -i $invol -o $tmpdir/samseg --reg-only --threads $threads)
    if($#InitReg) set cmd0 = ($cmd0 --init-reg $InitReg)
    set xopts = `fsr-getxopts run_samseg $XOptsFile`;
    set cmd = ($cmd $xopts);
    echo $cmd | tee -a $LF
    fs_time $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
else
  # Running full samseg may require 14GB of mem and add ~20min to run time
  # In theory, samseg could be used instead of seghead
  if($#SamsegJson) cp $SamsegJson $tmpdir
  set samseg = $samsegdir/seg.mgz
  set ud = (`UpdateNeeded $samseg $invol $InitReg`)
  if($ud || $ForceUpdate) then
    set cmd = (samseg --i $invol --o $samsegdir --threads $threads)
    if($SamsegFast)  set cmd = ($cmd --block-coordinate-descent)
    if($#SamsegJson) set cmd = ($cmd --options $SamsegJson)
    if($#InitReg) set cmd = ($cmd --initlta $InitReg);
    set xopts = `fsr-getxopts samseg $XOptsFile`;
    set cmd = ($cmd $xopts);
    echo $cmd | tee -a $LF
    fs_time $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
  set seg = $samseg
  set segidlist = (2 3 7 8 16 24 41 42 46 47 165) # 165=skull
  if($XEyeBalls) set segidlist = ($segidlist  259) # these are the eyeballs
  # Add subcort because it can help in the 2nd mkheadsurf
  set segidlist = ($segidlist  4  5 10 11 12 13 17 18 26 28 30 31 14 15 16 77 85)
  set segidlist = ($segidlist 43 44 49 50 51 52 53 54 58 60 62 63)
endif

if($DoSynthSegXMask) then
  # Segment using synthseg to create exclusion mask
  # Note: synthseg may need as much as 35GB of RAM and adds about 5min single-threaded
  set synthseg = $tmpdir/synthseg.mgz
  set ud = (`UpdateNeeded $synthseg $invol`)
  if($ud || $ForceUpdate) then
    set cmd = (synthseg --i $invol --o $synthseg --threads $threads)
    set xopts = `fsr-getxopts synthseg $XOptsFile`;
    set cmd = ($cmd $xopts);
    echo $cmd | tee -a $LF
    fs_time $cmd |& tee -a $LF
    if($status) goto error_exit
    set cmd = (mri_vol2vol --mov $synthseg --targ $invol --regheader \
     --interp nearest --o $synthseg)
    echo $cmd | tee -a $LF
    $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
  set seg = $synthseg
  set segidlist = (2 3 7 8 16 41 42 46 47)
  # Add subcort because it can help in the 2nd mkheadsurf
  set segidlist = ($segidlist  4  5 10 11 12 13 17 18 26 28 30 31 14 15 16 77 85)
  set segidlist = ($segidlist 43 44 49 50 51 52 53 54 58 60 62 63)
endif

# Compute reg from the defacing atlas to the subject
set atlas2inputlta = $tmpdir/atlas2input.lta
set ud = (`UpdateNeeded $atlas2inputlta $samseglta`)
if($ud || $ForceUpdate) then
  set cmd = (mri_concatenate_lta $atlas2samseglta $samseglta $atlas2inputlta)
  echo $cmd | tee -a $LF
  $cmd |& tee -a $LF
  if($status) goto error_exit
endif

# Compute the inverse of the reg
set input2atlaslta = $tmpdir/input2atlas.lta
set ud = (`UpdateNeeded $input2atlaslta $atlas2inputlta`)
if($ud || $ForceUpdate) then
  set cmd = (lta_convert --subject $defacingsubject --invert \
    --inlta $atlas2inputlta --outlta $input2atlaslta)
  echo $cmd | tee -a $LF
  $cmd |& tee -a $LF
  if($status) goto error_exit
endif

if($DoBBR) then
  # Run BBR to refine the registration to just the face
  # Careful using here. Faces are highly variable and the surface can get mangled.
  set bbrlta = $tmpdir/input2atlas.bbr.lta
  set bbrinvlta = $tmpdir/atlas2input2.bbr.lta
  set ud = (`UpdateNeeded $bbrlta $input2atlaslta`)
  if($ud || $ForceUpdate) then
    # Always use t1 since this is tissue vs background. Including zero
    # voxels is important to keep it from changing the scale to push
    # surface outside the FoV
    set cmd = (bbregister --mov $invol --t1 --wm-proj-abs 2  --gm-proj-abs 2 \
      --no-cortex-label --label $bbrface --init-reg $input2atlaslta --tolf 1e-5 \
      --reg $bbrlta --no-pass1 --$BBRDoF --lh-only  --no-brute2 --include-zero-voxels \
      --surf $templatesurf --target-volume $templatevol)
    set xopts = `fsr-getxopts bbregister $XOptsFile`;
    set cmd = ($cmd $xopts);
    echo $cmd | tee -a $LF
    $cmd |& tee -a $LF
    if($status) goto error_exit
    set cmd = (lta_convert --subject $defacingsubject --invert --inlta $bbrlta --outlta $bbrinvlta)
    echo $cmd | tee -a $LF
    fs_time $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
  set atlas2inputlta = $bbrinvlta
endif

# Map the atlas head/face surface to the subject
set templatesurfsubj = $tmpdir/template.surf
set ud = (`UpdateNeeded $templatesurfsubj $templatesurf $atlas2inputlta`)
if($ud || $ForceUpdate) then
  set cmd = (mris_apply_reg --lta $templatesurf $atlas2inputlta $templatesurfsubj)
  echo $cmd | tee -a $LF
  $cmd |& tee -a $LF
  if($status) goto error_exit
endif

# Create a head mask and surface of this subject
set headsurf = $tmpdir/head.surf
set headmask = $tmpdir/head.mgz
set ud = (`UpdateNeeded $headsurf $invol`)
if($ud || $ForceUpdate) then
  set cmd = (mkheadsurf -i $invol -o $headmask -surf $headsurf \
    -rescale -no-inflate -fill-holes-islands -mc -no-curv)
  # Adding samseg here makes the whole thing much more robust
  # because samseg is bascially provides a head mask
  if($UseOrMask && $DoSamsegXMask) set cmd = ($cmd -or-mask $samsegdir/seg.mgz)
  if($#fhi) set cmd = ($cmd -fhi $fhi);
  set xopts = `fsr-getxopts mkheadsurf $XOptsFile`;
  set cmd = ($cmd $xopts);
  echo $cmd | tee -a $LF
  fs_time $cmd |& tee -a $LF
  if($status) goto error_exit
endif

# Create statistics of the head volume for QA
set headstats = $tmpdir/headvol.stats
set ud = (`UpdateNeeded $headstats $headmask`)
if($ud || $ForceUpdate) then
  set cmd = (mri_segstats --seg $headmask --id 1 --sum $headstats)
  echo $cmd | tee -a $LF
  fs_time $cmd |& tee -a $LF
  if($status) goto error_exit
endif

if($DoSamsegXMask || $DoSynthSegXMask) then
  set segdil = $tmpdir/seg.bin.dil.mgz
  # Binarize, and maybe dilate, the samseg/synthseg output. A subset
  # of cerebral segids are selected to assure that only cerebral
  # structures are excluded. At the titme of writing, synthseg only
  # produces cerebral structures, so it is not an issue, but at some
  # point it might start segmenting extracerebral.  Probably the only
  # important seg is cortex, maybe skull and eyes.
  set ud = (`UpdateNeeded $segdil $seg`)
  if($ud || $ForceUpdate) then
    set cmd = (mri_binarize --i $seg --match $segidlist --o $segdil)
    if($nSegDil > 0) set cmd = ($cmd  --dilate $nSegDil)
    # This will keep stray skull voxels from interfering
    if($MaskSegToHead) set cmd = ($cmd --mask $headmask)
    echo $cmd | tee -a $LF
    $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
  set xmask = $segdil
endif

# Finally, perform the defacing
set minsurf = $tmpdir/min.surf
set maxsurf = $tmpdir/max.surf
set distdat = $tmpdir/dist.samples
set distlab = $tmpdir/dist.bounds.dat
set stats   = $tmpdir/stats.dat
#set facemask = $tmpdir/face.mask.mgz
set ud = (`UpdateNeeded $outvol $templatesurfsubj $invol $headmask $xmask $labellist $watermark`)
if($ud || $ForceUpdate) then
  set cmd = (mri_defacer --statframe $statframe --i $invol --ts $templatesurfsubj \
   --min $minsurf --max $maxsurf --m $facemask --o $outvol \
   --hm $headmask --distdat $distdat --stats $stats \
   --w $watermark $dwatermark) # --distbounds $distlab 
  if($FillType == 2) set cmd = ($cmd --fill-const $FillConstIn $FillConstOut)
  if($DoRipple) set cmd = ($cmd --ripple $RippleAmp $RipplePeriod)
  foreach label ($labellist)
    set cmd = ($cmd --l $label)
  end
  if($DoEars)       set cmd = ($cmd  --l $earleft  2 3 --l $earright 2 3)
  if($DoForehead)   set cmd = ($cmd  --l $forehead 2 3);
  if($DoBackOfHead) set cmd = ($cmd  --l $backofhead 2 3);
  if($#xmask) set cmd = ($cmd --xmask $xmask)
  set xopts = `fsr-getxopts mri_defacer $XOptsFile`;
  set cmd = ($cmd $xopts);
  echo $cmd | tee -a $LF
  fs_time $cmd |& tee -a $LF
  if($status) goto error_exit
  # Add the color table
  set cmd = (mri_convert --ctab $ctab $facemask $facemask)
  echo $cmd | tee -a $LF
  fs_time $cmd |& tee -a $LF
  if($status) goto error_exit
endif

# For QA, map the facemask onto the head surf
set fmoverlay = $tmpdir/facemask-on-headsurf.mgz
set ud = (`UpdateNeeded $fmoverlay $facemask $headsurf`)
if($ud || $ForceUpdate) then
  set cmd = (mri_vol2surf --vol2surf $facemask $headsurf 0 0 0 regheader novsm 0 $fmoverlay)
  echo $cmd | tee -a $LF
  $cmd |& tee -a $LF
  if($status) goto error_exit
endif

# For QA, generate a head surface after defacing
set posthead = $tmpdir/head.defaced.mgz
set postheadsurf = $tmpdir/head.defaced.surf
if($DoPostHeadSurf) then
  set ud = (`UpdateNeeded $posthead $outvol`)
  if($ud || $ForceUpdate) then
    set xopts = `fsr-getxopts mkheadsurf $XOptsFile`;
    set cmd = (mkheadsurf -i $outvol -o $posthead -surf $postheadsurf \
      -rescale -no-inflate -fill-holes-islands -mc -no-curv $xopts)
    if($#fhi) set cmd = ($cmd -fhi $fhi);
    # can't do samseg here, duh, but can use xmask if there
    if($UseOrMask && $#xmask) set cmd = ($cmd -or-mask $xmask)
    echo $cmd | tee -a $LF
    fs_time $cmd |& tee -a $LF
    if($status) goto error_exit
  endif
  # For QA, compute distance from true to new face
  set dist = $tmpdir/dist-to-face.mgz
  set distthresh = $tmpdir/dist-to-face.thresh.mgz
  set ud = `UpdateNeeded $distthresh $headsurf $postheadsurf`
  if($ud) then
    # Compute the distance between the true and refaced surfaces
    set cmd = (mris_diff --min-dist $headsurf $postheadsurf 1 $dist) #1=exact
    echo $cmd | tee -a $LF
    $cmd | tee -a $LF
    if($status) goto error_exit 
    # Threshold to exclude anything beyond 1mm
    set cmd = (mri_binarize --min 0.01 --max 1 --i $dist --o $distthresh)
    echo $cmd | tee -a $LF
    $cmd | tee -a $LF
    if($status) goto error_exit 
    set cmd = (mri_mask $dist $distthresh $dist)
    echo $cmd | tee -a $LF
    $cmd | tee -a $LF
    if($status) goto error_exit 
  endif
endif

# Take some pictures for QA
if($#SubjCode == 0) then
  set facepic = $tmpdir/face-before.png
  set postfacepic = $tmpdir/face-after.png
  set gif = $tmpdir/face-before+after.gif
else
  # Putting SubjCode in the name makes it easier for a slide show
  set facepic = $tmpdir/$SubjCode.face-before.png
  set postfacepic = $tmpdir/$SubjCode.face-after.png
  set gif = $tmpdir/$SubjCode.face-before+after.gif
endif
if($DoPics) then
  set ud = `UpdateNeeded $facepic $headsurf $fmoverlay`
  if($ud || $ForceUpdate) then
    # :overlay=${fmoverlay}:"overlay_threshold=.01,.011,.5"
    set cmd = (freeview -f ${headsurf}:affinexfm=${input2atlaslta}\
       --hide-3d-slices -viewport 3d -cam Azimuth -90 -ss $facepic --quit)
    set xopts = `fsr-getxopts freeview-before $XOptsFile`;
    set cmd = ($cmd $xopts);
    echo $cmd | tee -a $LF
    $runxvfb $cmd |& tee -a $LF
    if($status) goto error_exit
    if(-e $imconvert) then
      set cmd = ($imconvert -crop 398x557+224+38! $PicAnnot $facepic $facepic)
      $cmd |& tee -a $LF
      if($status) goto error_exit
    endif
  endif
  if(-e $postheadsurf) then
    set ud = `UpdateNeeded $postfacepic $postheadsurf`
    if($ud || $ForceUpdate) then
      set cmd = (freeview -f ${postheadsurf}:affinexfm=$input2atlaslta --hide-3d-slices \
        -viewport 3d -cam Azimuth -90 -ss $postfacepic --quit)
      set xopts = `fsr-getxopts freeview-after $XOptsFile`;
      set cmd = ($cmd $xopts);
      echo $cmd | tee -a $LF
      $runxvfb $cmd |& tee -a $LF
      if($status) goto error_exit
      if(-e $imconvert) then
        set cmd = ($imconvert -crop 398x557+224+38! $postfacepic $PicAnnot $postfacepic)
        $cmd |& tee -a $LF
        if($status) goto error_exit
      endif
    endif
  endif
  if(-e $facepic && -e $postfacepic && -e $imconvert) then
    # This is make an animation
    set ud = `UpdateNeeded $gif $facepic $postfacepic`
    if($ud || $ForceUpdate) then
      set cmd = ($imconvert -loop 0 -delay 60 $facepic $postfacepic $gif)
      echo $cmd | tee -a $LF
      $cmd |& tee -a $LF
      if($status) goto error_exit
    endif
  endif
endif

# Make a command line for viewing
set cmd = (fsvglrun freeview $invol $outvol --hide-3d-slices -viewport 3d \
  -f ${headsurf}:edgecolor=yellow)
#  -f ${minsurf}:edgecolor=blue:visible=0:lock=1 ${maxsurf}:edgecolor=cyan:visible=0:lock=1
if($DoPostHeadSurf) then
  set cmd = ($cmd -f ${postheadsurf}:edgecolor=red:visible=0)
else
  set cmd = ($cmd -f ${templatesurfsubj}:edgecolor=magenta:visible=0:lock=1)
endif

echo ""  |tee -a $LF
echo $cmd  |tee -a $LF
echo ""  |tee -a $LF
echo $cmd > $tmpdir/viewcmd

#========================================================

touch $donefile

# Cleanup
if($cleanup) rm -rf $tmpdir

# Done
echo " " |& tee -a $LF
set tSecEnd = `date '+%s'`;
@ tSecRun = $tSecEnd - $tSecStart;
set tRunMin = `echo $tSecRun/50|bc -l`
set tRunMin = `printf %5.2f $tRunMin`
set tRunHours = `echo $tSecRun/3600|bc -l`
set tRunHours = `printf %5.2f $tRunHours`
echo "Started at $StartTime " |& tee -a $LF
echo "Ended   at `date`" |& tee -a $LF
echo "Mideface-Run-Time-Sec $tSecRun" |& tee -a $LF
echo "Mideface-Run-Time-Min $tRunMin" |& tee -a $LF
echo "Mideface-Run-Time-Hours $tRunHours" |& tee -a $LF
echo " " |& tee -a $LF
echo "mideface Done" |& tee -a $LF
exit 0

###############################################

############--------------##################
error_exit:
echo "ERROR:"
echo $cmd | tee $errfile

exit 1;
###############################################

############--------------##################
parse_args:
set cmdline = ($argv);
while( $#argv != 0 )

  set flag = $argv[1]; shift;
  
  switch($flag)

    case "--o":
      if($#argv < 1) goto arg1err;
      set outvol = $argv[1]; shift;
      set stem = `fname2stem $outvol`
      if($status) then
        echo "$outvol is not a recognized format, treating as an output folder"
        set outdir = $outvol;
        set outvol = ()
        set DoPostHeadSurf = 1
      endif
      breaksw

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

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

    case "--fill-const":
      if($#argv < 2) goto arg2err;
      set FillConstIn  = $argv[1]; shift;
      set FillConstOut = $argv[1]; shift;
      set FillType = 2
      breaksw

    case "--fill-zero":
    case "--zero":
      set FillConstIn  = 0;
      set FillConstOut = 0;
      set FillType = 2
      breaksw

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

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

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

    case "--imconvert":
      if($#argv < 1) goto arg1err;
      set imconvert = $argv[1]; shift;
      if(! -e $imconvert) then
        echo "ERROR: cannot find $imconvert"
        exit 1
      endif
      breaksw

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

    case "--synthseg":
    case "--xmask-synthseg":
      if($#argv < 1) goto arg1err;
      set nSegDil = $argv[1]; shift;
      set DoSynthSegXMask = 1
      breaksw

    case "--samseg":
    case "--xmask-samseg":
      if($#argv < 1) goto arg1err;
      set nSegDil = $argv[1]; shift;
      set DoSamsegXMask = 1
      breaksw

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

    case "--xmask-eyeballs":
      set XEyeBalls = 1; #Samseg only
      breaksw
    case "--no-xmask-eyeballs":
      set XEyeBalls = 0;
      breaksw

    case "--samseg-json":
      if( $#argv < 1) goto arg1err;
      set SamsegJson = $argv[1]; shift;
      breaksw
    case "--no-samseg-json":
      set SamsegJson = ()
      breaksw
    case "--samseg-fast":
      set DoSamsegXMask = 1
      set nSegDil = 1
      set SamsegFast = 1;
      breaksw
    case "--no-samseg-fast":
      set SamsegFast = 0;
      breaksw
    case "--no-samseg":
      set SamsegFast = 0;
      set DoSamsegXMask = 0
      breaksw

    case "--atlas":
      if($#argv < 1) goto arg1err;
      set atlasdir = $argv[1]; shift;
      if(! -e $atlasdir) then
        echo "ERROR: cannot find $atlasdir"
        exit 1;
      endif
      breaksw

    case "--force":
     set ForceUpdate = 1
     breaksw
    case "--no-force":
     set ForceUpdate = 0
     breaksw

    case "--post":
     set DoPostHeadSurf = 1
     breaksw
    case "--no-post":
     set DoPostHeadSurf = 0
     breaksw

    case "--ears":
      set DoEars = 1
      breaksw
    case "--no-ears":
      set DoEars = 0
      breaksw

    case "--forehead":
      set DoForehead = 1
      breaksw
    case "--no-forehead":
      set DoForehead = 0
      breaksw

    case "--back-of-head":
      set DoBackOfHead = 1
      breaksw
    case "--no-back-of-head":
      set DoBackOfHead = 0
      breaksw

    case "--pics":
      set DoPics = 1
      breaksw
    case "--no-pics":
      set DoPics = 0
      breaksw

    case "--bbr":
      set DoBBR = 1
      breaksw
    case "--no-bbr":
      set DoBBR = 0
      breaksw
    case "--bbr-dof-6":
      set DoBBR = 1
      set BBRDoF = 6
      breaksw
    case "--bbr-dof-9":
      set DoBBR = 1
      set BBRDoF = 9
      breaksw
    case "--bbr-dof-12":
      set DoBBR = 1
      set BBRDoF = 12
      breaksw

    case "--mask-seg-to-head":
      set MaskSegToHead = 1;
      breaksw
    case "--no-mask-seg-to-head":
      set MaskSegToHead = 0;
      breaksw

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

    case "--no-xvfb":
      set UseXvfb = 0
      breaksw

    case "--code":
      if($#argv < 1) goto arg1err;
      set SubjCode = $argv[1]; shift;
      set PicAnnot = ("-fill white -pointsize 15 -annotate +15+15 $SubjCode")
      breaksw

    case "--ripple":
      if($#argv < 2) goto arg2err;
      set DoRipple = 1
      set RippleAmp = $argv[1]; shift;
      set RipplePeriod = $argv[1]; shift;
      breaksw

    case "--init-reg":
      if($#argv < 1) goto arg1err;
      set InitReg = $argv[1]; shift;
      if(! -e $InitReg) then
        echo "ERROR: cannot find $InitReg"
        exit 1
      endif
      breaksw

    case "--check":
      if($#argv < 1) goto arg1err;
      set vol = $argv[1]; shift;
      set outfile = ()
      if($#argv) then 
        set outfile = $argv[1]; shift;
        rm -f $outfile
      endif
      set cmd = (mri_defacer --check-code $vol $outfile)
      $cmd
      exit $status
      breaksw

    case "--apply":
      if($#argv < 4) goto arg4err;
      set vol = $argv[1]; shift;
      set facemask = $argv[1]; shift;
      set reg = $argv[1]; shift;
      set output = $argv[1]; shift;
      set udreg = $reg
      if($reg == regheader) set udreg = ()
      set ud = `UpdateNeeded $output $vol $facemask $udreg`
      if($ud) then
        set cmd = (mri_defacer --statframe $statframe --apply $vol $facemask $reg $output)
        echo $cmd 
        $cmd
        exit $status
      else
        echo "$output does not need to be updated"
        echo "If you want to force rerunning, then delete $output"
        exit 0
      endif
      breaksw

    case "--no-ripple":
      set DoRipple = 0
      breaksw

    case "--expert":
      if( $#argv < 1) goto arg1err;
      set XOptsFile = $argv[1]; shift;
      fsr-checkxopts $XOptsFile
      if($status) exit 1
      set XOptsFile = `getfullpath $XOptsFile`
      breaksw

    case "--or-mask"
      set UseOrMask = 1;
      breaksw
    case "--no-or-mask"
      set UseOrMask = 0;
      breaksw

    case "--mgz"
      set fmt = mgz;
      breaksw
    case "--nii"
      set fmt = nii
      breaksw
    case "--nii.gz"
      set fmt = nii.gz
      breaksw

    case "--display":
      if($#argv < 1) goto arg1err;
      setenv FSXVFB_START_D $argv[1]; shift;
      breaksw

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

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

    case "--tmp":
    case "--tmpdir":
      if($#argv < 1) goto arg1err;
      set tmpdir = $argv[1]; shift;
      set cleanup = 0;
      breaksw

    case "--nocleanup":
      set cleanup = 0;
      breaksw

    case "--cleanup":
      set cleanup = 1;
      breaksw

    case "--debug":
      set verbose = 1;
      set echo = 1;
      breaksw

    default:
      echo ERROR: Flag $flag unrecognized. 
      echo $cmdline
      exit 1
      breaksw
  endsw

end

goto parse_args_return;
############--------------##################

############--------------##################
check_params:

if($#invol == 0) then
  echo "ERROR: must spec input"
  exit 1;
endif
if(! -e $invol) then
  echo "ERROR: cannot find $invol"
  exit 1;
endif

if($#atlasdir == 0) set atlasdir = $FREESURFER/average/mideface-atlas
if(! -e $atlasdir) then
  echo "ERROR: cannot find atlasdir $atlasdir"
  exit 1
endif

if($#outdir) then
  if($#outvol == 0) set outvol = $outdir/defaced.$fmt
  if($#tmpdir == 0) set tmpdir = $outdir
  if($#facemask == 0) set facemask = $outdir/face.mask.mgz
  if($#LF == 0) then
    set LF = $outdir/log/mideface.log
    if(-e $LF) mv $LF $outdir/log/mideface.$$.log
  endif
  set cleanup = 0
endif

if($#outvol == 0) then
  echo "ERROR: must spec output"
  exit 1;
endif
if($#facemask == 0) then
  set facemask = $tmpdir/face.mask.mgz
endif

if($#xmask && $DoSynthSegXMask) then
  echo "ERROR: cannot use --xmask and --synthseg"
  exit 1
endif
if($DoSynthSegXMask && $DoSamsegXMask) then
  echo "ERROR: cannot use both --xmask-synthseg and --xmask-samseg"
  exit 1
endif

if($SamsegFast && $#SamsegJson) then
  echo "ERROR: cannot spec --samseg-fast and --samseg-json"
  exit 1
endif
if($SamsegFast) set SamsegJson = $atlasdir/samseg-options.json

set errfile = $tmpdir/log/mideface.error
rm -f $errfile

set donefile = $tmpdir/log/mideface.done
rm -f $donefile

goto check_params_return;
############--------------##################

############--------------##################
arg1err:
  echo "ERROR: flag $flag requires one argument"
  exit 1
############--------------##################
arg2err:
  echo "ERROR: flag $flag requires two arguments"
  exit 1
############--------------##################
arg4err:
  echo "ERROR: flag $flag requires four arguments"
  exit 1
############--------------##################

############--------------##################
usage_exit:
  echo ""
  echo "mideface : minimally invasive defacing tool"
  echo ""
  echo " --i volume to deface"
  echo " --o defaced input"
  echo " --facemask facemask"
  echo " --odir outputdir (turns on PostHeadSurf)"
  echo " --xmask xmask (exclusion mask)"
  echo " --xmask-samseg ndilations : segment input using samseg (14GB, +~20-40min)"
  echo " --samseg-json json: configure samseg"
  echo " --samseg-fast : configure samseg to run quickly; sets ndil=1 (default)"
  echo " --no-samseg-fast : do NOT configure samseg to run quickly"
  echo " --init-reg reg.lta : initialize samseg with reg (good in case samseg reg fails)"
  echo " --xmask-synthseg ndilations : segment input using synthseg (35GB, +~20min)"
  echo " --no-mask-seg-to-head : include stray skull voxels outside of the head mask"
  echo " --fill-const constIn constOut"
  echo " --fill-zero"
  echo " --fhi fhi : fhi value to use in MRIchangeType(); default is to use that in mri_seghead"
  echo " --no-ears : do not include ears in the defacing"
  echo " --back-of-head : include back of head in the defacing"
  echo " --forehead : include forehead in the defacing (risks removing brain)"
  echo " --pics : take pics (--no-pics)"
  echo " --code codename : embed codename in pics"
  echo " --imconvert /path/to/convert : path to imagemagik convert binary (for pics)"
  echo " --no-post : do not make a head surface after defacing"
  echo " --threads nthreads"
  echo " --statframe statframe : 0-based frame to get stats from (default is $statframe); use BEFORE --apply"
  echo " --force : force reprocessing (not applicable if --odir has not been used)"
  echo " --nii : use nifti format as output (only when output files are not specified)"
  echo " --nii.gz : use compressed nifti format as output (only when output files are not specified)"
  echo " --mgz : use compressed mgh format as output (default)"
  echo " --atlas atlasdir"
  echo " --expert xopts"
  echo " --display DisplayNo : set Xvfb display no for taking pics"
  echo ""
  echo "To apply midface output to a second volume"
  echo " --apply vol facemask reg output : apply to another volume (use regheader if no reg needed)"
  echo ""
  echo "To check whther a volume has been defaced"
  echo " --check vol <outfile>"
  echo "  This will print a 1 to stdout if this is an output of mideface or a 0 otherwise"
  echo "  outfile is an optional argument. If present, then the 0/1 will be printed in ascii to this file"
  echo ""

  if(! $PrintHelp) exit 1;
  echo $VERSION
  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

See https://surfer.nmr.mgh.harvard.edu/fswiki/MiDeFace for more info