#! /usr/bin/env perl
#
# Surface Feature Extractor
#
# Oliver Lyttelton oliver@bic.mni.mcgill.ca
#
# Copyright Alan C. Evans
# Professor of Neurology
# McGill University
#

use strict;
use warnings "all";
use Getopt::Tabular;
use File::Basename;
use File::Temp qw/ tempdir /;


my($Help, $Usage, $me);
my(@opt_table, $source, $target, $outxfm, @args, $tmpdir);

$me = &basename($0);

my $feature_choice = 0;
my $measurement_choice = 0;
my $extraction_transform = undef;
my $blurfwhm = undef;
my $alignment_transform = undef;
my $model_mid_surface = undef;
my $verbose = 0;
my $clobber = 0;
my $template_paint = undef;
my $template_labels = undef;
my $flip_transform = undef;
my $volume_feature = undef;
my $volume_feature_intersect = undef;

   
$Help = <<HELP;
|    $me extracts a surface feature set with options to 
|    choose the feature,  the type of measurement (roi vs vbc)
|    vary the extraction space (native or tal) the blurring kernel,
|    the registration space (tal or nl tal),
| 
| Problems or comments should be sent to: oliver\@bic.mni.mcgill.ca
HELP

$Usage = "Usage: $me [options] white.obj gray.obj output.txt\n".
         "       $me -help to list options\n\n";

@opt_table = (
        ["Feature Choices", "section"],
        
        ["-thickness","const","0",\$feature_choice,
        "Cortical Thickness (default)"] ,
        ["-surface_area","const","1",\$feature_choice,
        "Surface Area"] ,
        ["-complexity","const","2",\$feature_choice,
        "Cortical Complexity"] ,
          ["-cortical_volume","const","3",\$feature_choice,
        "Cortical Volume"] ,
          ["-callosal_distance","const","4",\$feature_choice,
        "Distance to Corpus Callosum"] ,
          ["-volume_feature","const","5",\$feature_choice,
        "Volumetric feature copied onto surface"] ,
         ["-volume_feature_mnc", "string", 1, \$volume_feature,
        "Volume from which to copy feature"],
          ["-volume_feature_intersect", "string", 1, \$volume_feature_intersect,
        "0 = White ... 2 = Mid ... 4= Gray"], 
           
       
                                                                  
        
        ["Measurement Type", "section"],
        
        ["-vbc","const","0",\$measurement_choice,
        "Vertex_based Corticometry (default)"] ,
         ["-asymm_vbc", "const", "2", \$measurement_choice,
        "flip and diff the vbc measurements for evaluation of asymmetries (nl registration recommended)"],           
            ["-roi","const","1",\$measurement_choice,
        "Region of Interest Analysis (must specify template paint and labels)"] ,      
        ["-asymm_vbc_transform","string",1,\$flip_transform,
        "Flip_Transform for asymmetry analysis"],
        ["-roi_template_paint", "string", 1, \$template_paint,
        "Node file with label paint. (NB. Template paint should correspond to alignment space)"], 
        
        ["Extraction Space", "section"],   
        
        ["-extraction_transform","string",1,\$extraction_transform,
        "Extraction_Transform for native or 7P extraction (otherwise 9P talairach space will be used"],
        
        ["Alignment Space", "section"],
        
        ["-alignment_transform","string",1,\$alignment_transform,
        "Alignment_Transform for Non-Linear Alignment (otherwise the unregistered node correspondence will be used). For Surface Area and  Cortical Volume vbc, and for all roi analysis it must be a target->source mapping. For anything else it must be a source ->target mapping"] ,
       
        ["Blurring","section"],
        ["-blurfwhm","string",1,\$blurfwhm,
        "Blurring Kernel (default 0. NB No Blurring with roi analysis)"],
        ["-model_mid_surface","string",1,\$model_mid_surface,
          "model object for area and volume blurring (for *highly experimental* area and volume corticometry :-> )"],
        
        
        ["Control", "section"], 
        
        ["-verbose","const","1", \$verbose,
        "be verbose" ],
        ["-clobber", "const","1", \$clobber,
        "clobber existing check features" ],  

   );

# Check arguments
&Getopt::Tabular::SetHelp($Help, $Usage);
&GetOptions (\@opt_table, \@ARGV) || exit 1;
die $Usage if($#ARGV != 2);
my $white = shift(@ARGV);
my $gray = shift(@ARGV);
my $output = shift(@ARGV);


# check for features
die "$me: Couldn't find input feature: $white\n" if (!-e $white);
die "$me: Couldn't find input feature: $gray\n" if (!-e $gray);
if(-e $output && !$clobber){
   die "$me: $output exists, -clobber to overwrite\n";
   };
   
if ($measurement_choice == 1 && $blurfwhm !=0){die ("no cakie and eaty. roi=no blur!");  };
if ($measurement_choice == 0 && !$alignment_transform &&($feature_choice==1 || $feature_choice==4)){die ("no cakie and eaty. vbc for surface_area and cortical_volume must use model-alignment.");  };
# make tmpdir
$tmpdir = &tempdir( "$me-XXXXXXXX", TMPDIR => 1, CLEANUP => 1 );



                
#set up internal variables
my $mid_surface = "${tmpdir}/mid.obj";
my $native_white = "${tmpdir}/native_white.obj";
my $native_gray  = "${tmpdir}/native_gray.obj";
my $raw_feature = "${tmpdir}/raw_feature.txt";
my $blurred_feature = "${tmpdir}/blurred_feature.txt"; 
 my $resampled_feature = "${tmpdir}/resample.txt";
my $flip_feature = "${tmpdir}/flip.txt";
my $diff_feature = "${tmpdir}/diff.txt";
my  $paint_on_target = "${tmpdir}/paint_on_target.txt";
my $node_paint = "${tmpdir}/node_paint.txt";

my $source_feature = $raw_feature;

my $source_white =$white;
my $source_gray = $gray;
 
#Optional pre-blurring transformation
  	if ($extraction_transform) {
       #perform_pre_blurring_transformation step
        &do_cmd("transform_objects",$source_white,$extraction_transform,$native_white);
        &do_cmd("transform_objects",$source_gray,$extraction_transform,$native_gray);
        
        $source_white = $native_white;
        $source_gray = $native_gray;
        $source_feature = $raw_feature;
	   
	}
    
#Extract measurement

  &do_cmd("average_surfaces",$mid_surface,"none","none",1,$source_white,$source_gray);

  if ($feature_choice ==0){
          #thickness
          &do_cmd("dump_rms",$source_white,$source_gray,$source_feature);
  }
 if ($feature_choice ==1){
          #surface_area                               
               &do_cmd("/data/aces/aces1/oliver/cvs/pipelines/civet/progs/cortex_area","-surface", $mid_surface,"-output",$source_feature);
  }
   if ($feature_choice ==2){
          #cortical_complexity
          &do_cmd("cortical_complexity",$mid_surface,$source_feature) ;
  }
     if ($feature_choice ==3){
          #cortical volume
           &do_cmd("cortical_volume",$source_white,$source_gray,$source_feature);
          
 }
 if ($feature_choice ==4){
          #volume_feature
           &do_cmd("volume_feature_on_surface",$source_white,$source_gray,$volume_feature,$volume_feature_intersect,$source_feature);
          
 }
   if ($feature_choice ==5){
          #callosal distance
          &do_cmd("callosal_distance",$source_white,$source_feature) ;
  }
 
 
 #ROI Analysis
 if ($measurement_choice==1) {
       #then we are doing ROI analysis 
       if ($alignment_transform){   
                 &do_cmd("surface-resample2","-clobber",$alignment_transform, $template_paint,$paint_on_target);
                 $template_paint =$paint_on_target ;
       }
          &do_cmd("roi_process","-clobber",$template_paint,$source_feature,$output);    
	}
   
 
 #VBC Cortical complexity, callosal distance, and thickness	
 if (($measurement_choice!=1)&&!($feature_choice==1 || $feature_choice==4)){ 
         if ($blurfwhm){
                &do_cmd( "diffuse", "-kernel", $blurfwhm, "-parametric", "0", $mid_surface, $source_feature, $blurred_feature);
                $source_feature = $blurred_feature; 
         }
          if ($alignment_transform){  
                  &do_cmd("surface-resample2","-clobber",$alignment_transform, $source_feature,$resampled_feature ); 
          }
 }
                                                            
#VBC Surface area and volume   
    if (($measurement_choice!=1) &&($feature_choice==1 || $feature_choice==4)) {               
           #create a node paint feature
           open(myFile,">${node_paint}");
           my $i = 1;
           while($i<40963){
           print myFile "${i}\n";
           $i=$i+1;
           }
           close myFile;
          &do_cmd("surface-resample2","-clobber",$alignment_transform, $node_paint,$paint_on_target ); 
          &do_cmd("zone_to_vertex_areas",$paint_on_target,$source_feature,$resampled_feature);
          $source_feature = $resampled_feature;
          if ($blurfwhm){
          &do_cmd( "diffuse", "-kernel", $blurfwhm, "-parametric", "0", $model_mid_surface, $source_feature, $blurred_feature);
         $source_feature = $blurred_feature;               
           }                            
}

#Optional Asymmetry stages                                                                      
  	if ($measurement_choice==2){
	    #perform_symmetry step
		&do_cmd("surface_resample2",$flip_transform,$source_feature, $flip_feature);
        #perform_symmetry step
        &do_cmd("surfmath","-clobber","-sub",$source_feature,$flip_feature,$diff_feature);
		$source_feature = $diff_feature;	
  }
  
  #copy to output
  &do_cmd("cp", $source_feature ,$output);
  


sub do_cmd { 
   print STDOUT "@_\n" if ($verbose);
   system(@_) == 0 or die;
}