#!/usr/bin/perl 
# Script to automate iterative runs of CRITICA

(@ARGV>=3) || usage();

$fraction=0.8;
$iterations=3;
$addlong=0;
$startiteration=0;
$dicodonscores="";
$initscores="";
$options="";
$fixedoptions="";
process_options();

$tripletsfile=pop(@ARGV);
$seqfile=pop(@ARGV);
$outputprefix=pop(@ARGV);

if(!-e($seqfile)) {
  die("iterate-critica: No such file as $seqfile!\n");
}
if(!-e($tripletsfile)) {
  die("iterate-critica: No such file as $tripletsfile!\n");
}

$totalinits=$outputprefix.".init";


for($i=$startiteration;$i<=$iterations;$i++) {
  $options.=$fixedoptions;
  $scratch=$outputprefix.".tmp";
  $scratch2=$outputprefix.".rmp";
  $criticaoutput=$outputprefix.$i;
  $crit=$criticaoutput.".crit";
  $cdsfile=$criticaoutput.".cds";
  $rejectfile=$criticaoutput.".reject";
  $sdscores=$criticaoutput.".sdscores";
  $dicodonscores=$criticaoutput.".dicod";
  $initscores=$criticaoutput.".init";
  system("critica $seqfile $tripletsfile $options > $crit");
  system("best-inits $crit > $scratch");
  if (($i==0) && ($addlong>0)) {
    system("addlongorfs -orf-aa-length=$addlong $scratch $seqfile >> $scratch");
  }
  system("remove-overlaps 20 $scratch > $scratch2 2>$rejectfile");
  system("remove-overlaps 20 $scratch2 > $cdsfile 2>>$rejectfile");
  if ($i<$iterations) {
    system("best-sd $crit > $scratch");
    system("compute-SDscores $crit $scratch > $sdscores");
    if ($i>0) {
      system("dicodontable $cdsfile $seqfile > $dicodonscores");
    }
    else {
     system("dicodontable -fraction-coding=$fraction $cdsfile $seqfile > $dicodonscores");
   }
    system("list-inits $cdsfile $seqfile > $scratch");
    if (!-e($totalinits)) {
      system("list-total-inits $seqfile > $totalinits");
    }
    system("compute-init-scores $totalinits $scratch > $initscores");
    $options="-dicodon-scores=$dicodonscores -init-scores=$initscores -sd-scores=$sdscores";
  }
}
$criticaoutput=$outputprefix.$i;
$cdsfile=$criticaoutput.".cds";
unlink $scratch, $scratch2;

sub usage {
  printf("iterate-critica: [options] output-name seq-file triplets-scores\n");
  printf("Valid iterate-critica options:\n");
  printf("\t-iterations=number (default: 3)\n");
  printf("\t-scoring-matrix=file\n");
  printf("\t-initial-dicod=file\n");
  printf("\t-initial-init=file\n");
  printf("\t-threshold=value\n");
  printf("\t-fraction-coding=value\n");
  printf("\t-add-longorfs=length\n");
  printf("\t-genetic-code=number\n");
  printf("\t-strict-threshold\n"); 
  die("\t-initial-sd=file\n");
}

sub process_options {
  for($i=0;$i<@ARGV;$i++) {
    if (substr($ARGV[$i],0,1) eq "-") {
      $arg=$ARGV[$i];
      if (index($arg,"-iterations=")>-1) {
        ($dummy,$iterations)=split("=",$arg);
      }
      elsif (index($arg,"-scoring-matrix=")>-1) {
        $fixedoptions.=" $arg";
      }
      elsif (index($arg,"-initial-dicod=")>-1) {
        ($dummy,$dicodonscores)=split("=",$arg);
        $startiteration=1;
        $options.="-dicodon-scores=$dicodonscores";
      }
      elsif (index($arg,"-initial-init=")>-1) {
        ($dummy,$initscores)=split("=",$arg);
        $startiteration=1;
        $options.=" -init-scores=$initscores";
      }
      elsif (index($arg,"-threshold=")>-1) {
        $fixedoptions.=" $arg";
      }
      elsif (index($arg,"-genetic-code=")>-1) {
        $fixedoptions.=" $arg";
      }
      elsif (index($arg,"-alpha=")>-1) {
        $fixedoptions.=" $arg";
      }
      elsif (index($arg,"-strict-threshold")>-1) {
        $fixedoptions.=" -strict-threshold";
      }
      elsif (index($arg,"-add-longorfs=")>-1) {
        ($dummy,$addlong)=split("=",$arg);
      }
      elsif (index($arg,"-initial-sd=")>-1) {
        ($dummy,$sdscores)=split("=",$arg);
        $startiteration=1;
        $options.=" -sd-scores=$sdscores";
      }
      elsif (index($arg,"-fraction-coding=")>-1) {
        ($dummy,$fraction)=split("=",$arg);
      }
      else {
        die("iterate-critica: option $arg not valid!\n");
      }
    }
  }
}
