#! /bin/tcsh -f

#
# synth-tr-sess -- synthesize a session based on task-rest data
#
# Original Author: Doug Greve
# CVS Revision Info:
#    $Author: nicks $
#    $Date: 2007/01/09 22:41:19 $
#    $Revision: 1.2 $
#
# Copyright (C) 2002-2007,
# The General Hospital Corporation (Boston, MA). 
# All rights reserved.
#
# Distribution, usage and copying of this software is covered under the
# terms found in the License Agreement file named 'COPYING' found in the
# FreeSurfer source code root directory, and duplicated here:
# https://surfer.nmr.mgh.harvard.edu/fswiki/FreeSurferOpenSourceLicense
#
# General inquiries: freesurfer@nmr.mgh.harvard.edu
# Bug reports: analysis-bugs@nmr.mgh.harvard.edu
#


set VERSION = '$Id: synth-tr-sess,v 1.2 2007/01/09 22:41:19 nicks Exp $';

set inputargs = ($argv);
set DateStr = "`date '+%y%m%d%H%M'`"

set analysis = ();
set maskcon = ();
set maskmap = ();
set maskthresh = 2;

set uberbasestem = ();

set srcsess = ();
set trgsess = ();

set restrlf = ();
set restbasestem = ();
set synthnoise = 0;
set nonoise = 0;

set outbasestem = ();

set monly = 0;
set QOE = 1;
set MLF = ();

set logunique = 0;
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

set PWD = `getpwdcmd`;
if($status) exit 1;

goto parse_args;
parse_args_return:

goto check_params;
check_params_return:

##### Create a log file ######
set logdir = `pwd`/log;
mkdir -p $logdir
if(! -e $logdir) then
  echo "WARNING: could not create $logdir"
  set LF = /dev/null
else
  if($logunique) then
    set LF = $logdir/synth-tr-sess-$DateStr.log
  else
    set LF = $logdir/synth-tr-sess.log
  endif
  if(-e $LF) mv $LF $LF.old
endif
echo "--------------------------------------------------------------"
echo "synth-tr-sess logfile is $LF"
echo "--------------------------------------------------------------"

echo "synth-tr-sess log file" >> $LF
echo $VERSION >> $LF
pwd           >> $LF
echo $0         >> $LF
echo $inputargs >> $LF
uname -a        >> $LF
date            >> $LF

set StartTime = `date`;

# get full path for cfg and info files #
pushd $analysis > /dev/null;
set analysisdir = `pwd`;
popd > /dev/null;
set cfgfile  = $analysisdir/analysis.cfg
set infofile = $analysisdir/analysis.info
if(! -e $infofile) then
  echo "ERROR: cannot find $infofile"
  exit 1;
endif

## Get parameters from the info file ##
set fsd       = `cat $infofile | awk '{if($1 == "fsd") print $2}'`;
set taskbasestem  = `cat $infofile | awk '{if($1 == "funcstem") print $2}'`;
set parname   = `cat $infofile | awk '{if($1 == "parname") print $2}'`;
set taskrlf  = `cat $infofile | awk '{if($1 == "runlistfile") print $2}'`;

if($#outbasestem == 0)  set outbasestem  = $taskbasestem;
if($#restbasestem == 0) set restbasestem = $taskbasestem;

echo "----------------------------------" >> $LF
echo "src          $srcsess" >> $LF
echo "trg          $trgsess" >> $LF
echo "analysis     $analysis" >> $LF
echo "maskcontrast $maskcon" >> $LF
echo "maskmap      $maskmap" >> $LF
echo "maskthresh   $maskthresh" >> $LF
echo "ubermask     $uberbasestem" >> $LF
echo "synthnoise   $synthnoise" >> $LF
echo "nonoise      $nonoise" >> $LF
echo "restrlf      $restrlf" >> $LF
echo "reststem     $restbasestem" >> $LF
echo "----------------------------------" >> $LF

# Get full path to target session #
mkdir -p $trgsess/$fsd
pushd $trgsess > /dev/null
set trgsess = `pwd`;
popd > /dev/null

# Go to the source session FSD #
cd $srcsess/$fsd

if(! -e $analysis) then
  echo "ERROR: cannot find $analysis in source"
  exit 1;
endif

if(! -e $taskrlf) then
  echo "ERROR: could not find run list file $taskrlf" |& tee -a $LF
  exit 1;
endif
set TaskRunList = (`cat $taskrlf`);

if(! $nonoise) then
  if(! -e $restrlf) then
    echo "ERROR: could not find run list file $restrlf" |& tee -a $LF
    exit 1;
  endif
  set RestRunList = (`cat $restrlf`);

  if($#TaskRunList != $#RestRunList) then
    echo "ERROR: Rest and Task RLFs have different number of runs "|& tee -a $LF
    exit 1;
  endif
endif
set nRuns = $#TaskRunList;

# Copy some basic info to target #
cp $taskrlf $trgsess/$fsd
cp seq.info $trgsess/$fsd
if(-e register.dat) cp register.dat $trgsess/$fsd
if(-e ../subjectname ) cp ../subjectname $trgsess

mkdir -p $trgsess/$fsd/masks
set signalmaskstem = $trgsess/$fsd/masks/signal
if(-e masks/brain_000.hdr) then
  echo "Copying brain mask"
  cp masks/brain_???.* $trgsess/$fsd/masks
  cp masks/brain.bhdr $trgsess/$fsd/masks
  cp masks/brain.bhdr $trgsess/$fsd/masks/signal.bhdr
endif

if($#MLF == 0) set MLF = $trgsess/$fsd/synth_tr_sess-$$.m
rm -f $MLF
echo "% matlab file for synth-tr-sess%" >> $MLF
echo "% `date` %" >> $MLF
echo "cd `pwd`;"  >> $MLF
echo "QuitOnError = $QOE;" >> $MLF
echo "" >> $MLF

set nthRun = 1;
while ($nthRun <= $nRuns)
  set TaskRun = $TaskRunList[$nthRun];

  set outdir = $trgsess/$fsd/$TaskRun;
  mkdir -p $outdir

  set signalstem = `printf $analysis/signal/s%03d $nthRun`;
  set tmp = `getbfvdim $signalstem`;
  if($status) then
    echo "ERROR: $signalstem"
    exit 1;
  endif
  set nslices = $tmp[1];
  
  set taskstem = $TaskRun/$taskbasestem
  getbfvdim $taskstem > /dev/null
  if($status) then
    echo "ERROR: $taskstem"
    exit 1;
  endif

  set maskstem = ();
  if($#maskcon != 0) then
    set maskstem = $analysis/$maskcon/$maskmap
    getbfvdim $maskstem > /dev/null
    if($status) then
      echo "ERROR: $maskstem"
      exit 1;
    endif
  endif

  set uberstem = ();
  if($#uberbasestem != 0) then
    set uberstem = masks/$uberbasestem;
    getbfvdim $uberstem > /dev/null
    if($status) then
      echo "ERROR: $uberstem"
      exit 1;
    endif
  endif

  set reststem = ();
  if(! $nonoise) then
    set RestRun = $RestRunList[$nthRun];
    set reststem = $RestRun/$restbasestem
    getbfvdim $reststem > /dev/null
    if($status) then
      echo "ERROR: $reststem"
      exit 1;
    endif
  endif

  set parfile = $TaskRun/$parname
  if(! -e $parfile ) then
    echo "ERROR: cannot find $parfile"
    exit 1;
  endif
  cp $parfile $outdir

  set outstem = $outdir/$outbasestem
  if(-e $taskstem.bhdr)   cp $taskstem.bhdr $outstem.bhdr


#------------- This adds to the matlab file -------------#
tee -a $MLF > /dev/null <<EOF
  %--------------------------------%
  % Run $nthRun
  nthRun = $nthRun;
  signalstem      = '$signalstem';
  signalmaskstem  = '$signalmaskstem';
  uberstem    = '$uberstem';
  maskstem    = '$maskstem';
  maskthresh  = [$maskthresh];
  reststem    = '$reststem';
  outstem     = '$outstem';
  nonoise     = $nonoise;
  synthnoise  = $synthnoise;
  nslices     = $nslices;

  fprintf('\nnthRun = %d ----------------------------\n',$nthRun);

  for slice = 0:nslices-1
    fprintf('slice = %d --------------------\n',slice);

    signal = fast_ldbslice(signalstem,slice);
    if(isempty(signal))
      fprintf('ERROR: loading %s slice %d\n',signal,slice);
      return;
    end
    [nrows ncols ntrs] = size(signal);
    nvoxs = nrows*ncols;
    signal = reshape(signal, [nvoxs ntrs])'; %'
    out = signal;

    mask = [];
    if(~isempty(maskstem))
      mask = fast_ldbslice(maskstem,slice);
      if(isempty(mask))
        fprintf('ERROR: loading %s slice %d\n',maskstem,slice);
        return;
      end
      indzero = find(abs(mask) < maskthresh);
      out(:,indzero) = 0;
      if(nthRun==1) fprintf('mask: nzero = %d\n',length(indzero)); end
    end

    uber = [];
    if(~isempty(uberstem))
      uber = fast_ldbslice(uberstem,slice);
      if(isempty(uber))
        fprintf('ERROR: loading %s slice %d\n',uberstem,slice);
        return;
      end
      indzero = find(abs(uber) < 0.5);
      out(:,indzero) = 0;
      if(nthRun==1) fprintf('uber: nzero = %d\n',length(indzero)); end
    end

    if(~isempty(mask) & ~isempty(uber))
      signalmask = (abs(mask)>maskthresh) & uber;
    elseif(~isempty(mask))
      signalmask = (abs(mask)>maskthresh);
    elseif(~isempty(uber))
      signalmask = uber;
    else
      signalmask = [];
    end
    if(~isempty(signalmask)) 
      nsignalmask = length(find(signalmask));
      if(nthRun==1) fprintf('final: %d\n',nsignalmask); end
      fname = sprintf('%s_%03d.bshort',signalmaskstem,slice);
      fmri_svbfile(signalmask,fname);
    end

    if(~ nonoise)
      rest = fast_ldbslice(reststem,slice);
      if(isempty(rest))
        fprintf('ERROR: loading %s slice %d\n',reststem,slice);
        return;
      end
      rest = reshape(rest, [nvoxs ntrs])'; %'
      if(synthnoise)
        restmean = mean(rest);
        reststd = std(rest);
        rest = randn(size(rest));
        rest = rest .* repmat(reststd,[ntrs 1]);
        rest = rest + repmat(restmean,[ntrs 1]);
      end
      out = out + rest;
    end

    out = reshape(out', [nrows ncols ntrs]); %'
    outfile = sprintf('%s_%03d.bfloat',outstem,slice);
    %fprintf('Saving out to %s\n',outfile);
    fmri_svbfile(out,outfile);
  end

EOF
#----------------------------------------------------------------#
  echo "" >> $MLF

  @ nthRun = $nthRun + 1;
end

echo "---------------------------------------------------" >> $LF
cat $MLF >> $LF
echo "---------------------------------------------------" >> $LF

if(! $monly) then
  cat $MLF | matlab -nosplash -display iconic|& tee -a $LF
  rm -f $MLF 
endif

set EndTime = `date`;

echo "" |& tee -a $LF
echo "Started at $StartTime" |& tee -a $LF
echo "Ended   at $EndTime"   |& tee -a $LF

cp $LF $trgsess

echo "synth-tr-sess done"

echo ""

exit 0

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

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

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

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

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

    case "-analysis":
    case "-a":
      if ( $#argv == 0) goto arg1err;
      set analysis = $argv[1]; shift;
      breaksw
 
    case "-maskcon":
    case "-c":
      if ( $#argv == 0) goto arg1err;
      set maskcon = $argv[1]; shift;
      breaksw
 
    case "-maskmap":
    case "-m":
      if ( $#argv == 0) goto arg1err;
      set maskmap = $argv[1]; shift;
      breaksw

    case "-maskthresh":
    case "-thresh":
    case "-t":
      if ( $#argv == 0) goto arg1err;
      set maskthresh = $argv[1]; shift;
      breaksw

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

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

    case "-synth":
      set synthnoise = 1;
      breaksw

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

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

    case "-monly":
      if ( $#argv == 0) goto arg1err;
      set MLF = $argv[1]; shift;
      set monly = 1;
      set QOE = 0;
      breaksw

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

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

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

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

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

end

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

############--------------##################
check_params:
 
  if($#analysis == 0) then
    echo "ERROR: must specify analysis"
    exit 1;
  endif

  if($#srcsess == 0) then
    echo "ERROR: must specify source session"
    exit 1;
  endif

  if(! -e $srcsess) then
    echo "ERROR: cannot find source session $srcsess"
    exit 1;
  endif

  if($#trgsess == 0) then
    echo "ERROR: must specify target session"
    exit 1;
  endif

  if(! $nonoise && $#restrlf == 0) then
    echo "ERROR: must specify rest run list file"
    exit 1;
  endif

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

############--------------##################
arg1err:
  echo "ERROR: flag $flag requires one argument"
  exit 1
############--------------##################

############--------------##################
usage_exit:
  echo ""
  echo "USAGE: synth-tr-sess"
  echo ""
  echo "Arguments:";
  echo "   -src sessid : source session"
  echo "   -a   analysis : anlysis in source"
  echo "   -c   mask constrast : mask source signal"
  echo "   -m   mask map       : mask source signal"
  echo "   -uber maskstem : from fsd/mask"
  echo "   -rlf restrlf  : run list file for rest runs"
  echo "   -synth : use white noise (std based on rest)"
  echo "   -nonoise "
  echo "   -trg sessid : target session"
  echo ""
  echo "Yet even more options"
  echo "   -help  : hum, what's this?"
  echo "   -debug : print out ridiculous amounts of verbiage"
  echo "   -monly mlf : only create matlab file"
  echo "   -umask   mask : unix file permission mask"
  echo "   -logunique : make a unique log file name"
  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

This will synthesize a session based on task and rest runs.