#!/usr/bin/perl
use strict;
use warnings;

# 
# Makes images for event files supplied by user over entire FOV. 
# If a region file is also given, an image will also be made for the region. 
#
#
# Requirements:
#    0) must be in same directory as input event files. 
#    1) text file which lists the events files to process, one per line
#
# Optional: 
#    0) region file which lists the regions of interest, one per line, in phyical units. 
#       Must be in ds9 region file format. SAS can accept only circles, ellipses, boxes, 
#       and annuli.
#


my $currDir;
chomp ($currDir = `pwd`);
my $unittype = 'N';
my $currLine;
my @files;
my $file;
my $m1_evt;
my $m2_evt;
my $pn_evt;
my $region_file=$currDir."/"."filler";
my $status=0;
my $process_pn = 0;
my $process_m1 = 0;
my $process_m2 = 0;
my $j;
my $inst;
my $region_flag=0;
my @extraction_region;
my $source;

####################################
# Set the environment variables.
####################################

my @res1;
my $baseDir;

if ($currDir =~ m/proc/) {
   @res1 = split(/proc/, $currDir);
   $baseDir = $res1[0];
}

if (-d $baseDir."ODF") { 
   $ENV{SAS_ODFPATH} = $baseDir."ODF";
   my $sumsas = glob ($baseDir."ODF/*SUM.SAS");
   $ENV{SAS_ODF} = $sumsas;
   $ENV{SAS_CCF} = $baseDir."ODF/ccf.cif";
} 
if (-d $baseDir."odf") {
   $ENV{SAS_ODFPATH} = $baseDir."odf";
   my $sumsas = glob ($baseDir."odf/*SUM.SAS");
   $ENV{SAS_ODF} = $sumsas;
   $ENV{SAS_CCF} = $baseDir."odf/ccf.cif";
}

####################################
# Get the event files.
####################################


print "Enter the file listing the event files: ";
chomp (my $event_files = <> );
$event_files = $currDir."/".$event_files;

print "Should I use a region file? If so, enter the file name; otherwise, enter N. ";
chomp (my $ans = <>) ;
if (($ans eq "N") or ($ans eq "n")) {
   $region_flag = 0 ;
} else {
   $region_flag = 1;
   $region_file = $currDir."/".$ans;
}

if (-e $event_files) {
   open(EVENTS, "<", $event_files);
   while ($currLine = <EVENTS>) {
      $currLine = $currDir."/".$currLine;
      chomp ($currLine);
      @files = split (/\s/, $currLine);
      foreach $file (@files) {  # Check the headers to find out which instruments we need.
         $inst = $currDir."/"."instru.list";
         `fkeyprint $file instrume outfile=$inst`;
         open (INSTRU, "<", $inst);
         while ($currLine = <INSTRU>) {
            if ($currLine =~ m/INSTRUME= 'EPN     '/ ) {
               $pn_evt = $file;
               $process_pn = 1;
            }
            if ($currLine =~ m/INSTRUME= 'EMOS1   '/ ) {
               $m1_evt = $file;
               $process_m1 = 1;
            }
            if ($currLine =~ m/INSTRUME= 'EMOS2   '/ ) {
               $m2_evt = $file;
               $process_m2 = 1;
            }
         }
         close (INSTRU);
         `rm $inst`;
      }
   }
   close (EVENTS);
}
else {
   print "No events file found; stopping here. \n" ;
   exit;
}

#####################################################
# Does the region file exist and have proper format?
# If so, grab the region info.
######################################################

my $count_regions=0;
my $test_ds9 = 0;
my $test_physical = 0;

if ($region_flag) { 
   if (-e $region_file) {
      print "Reading region file \n";
      open(REGION, "<", $region_file);
      while ($currLine = <REGION>) {
         if ($currLine =~ m/ DS9 /i) {
            $test_ds9 = 1;
         }
         if ($currLine =~ m'physical\n') {
            $test_physical = 1;
         }
      }
      close (REGION);
   }
   else {
      print "No region file found. \n" ;
      exit;
   }

   if ($test_ds9 == 0) {
      print "Region file is not in DS9 format. \n";
      exit;
   }
   if ($test_physical == 0) {
      print "Region file is not in correct units. \n";
      exit;
   }

   open(REGION, "<", $region_file);
   for ($j=0; $j<3; $j++) {   # Skip the header
      readline(REGION);
   }
   while ($currLine = <REGION>) {
      if (($currLine =~ m/circle/i) or ($currLine =~ m/ellipse/i) or ($currLine =~ m/box/i) or ($currLine =~ m/polygon/i) or ($currLine =~ m/annulus/i) or ($currLine =~ m/panda/i) or ($currLine =~ m/epanda/i) or ($currLine =~ m/bpanda/i)) {
         chomp($currLine);
         push (@extraction_region, $currLine);
         $count_regions ++;
      }
      else {
         print "Invalid shape in the region file! \n";
      }
   }
   close (REGION);
} 

#####################################
# Make the images of the full FOV.
#####################################

if (! -d $currDir.'/image') {
   `mkdir image`;
}

my $m1_evt_cp;
my $m2_evt_cp;
my $pn_evt_cp;

if ($process_m1 eq 1) {
    $m1_evt_cp = "image/m1_evt_cp";
    `cp $m1_evt $m1_evt_cp`;
}
if ($process_m2 eq 1) {
   $m2_evt_cp = "image/m2_evt_cp";
   `cp $m2_evt $m2_evt_cp`;
}
if ($process_pn eq 1) {
   $pn_evt_cp = "image/pn_evt_cp";
   `cp $pn_evt $pn_evt_cp`;
}

chdir $currDir.'/image';

if ($process_m1) {
   print "Making the M1 FOV image... \n";
   `evselect table=m1_evt_cp withimageset=yes imageset=m1_img_fov.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32`;
}
if ($process_m2) {
   print "Making the M2 FOV image... \n";
   `evselect table=m2_evt_cp withimageset=yes imageset=m2_img_fov.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32`;
}
if ($process_pn) {
   print "Making the PN FOV image... \n";
   `evselect table=pn_evt_cp withimageset=yes imageset=pn_img_fov.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32`;
}

#############################################################
# Make the images of the sources, if there is a region file.
##############################################################

if ($region_flag) {
   for ($j=0; $j<$count_regions; $j++) {
      $source = shift(@extraction_region);
      print "Making images of the source(s)... \n";
      if ($process_m1 == 1) {
         `evselect table=$m1_evt withimageset=yes imageset=m1_img_$j.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32 filtertype='expression' expression='((X,Y) in $source)'` ;
      }
      if ($process_m2 == 1) {
         `evselect table=$m2_evt withimageset=yes imageset=m2_img_$j.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32 filtertype='expression' expression='((X,Y) in $source)'`;
      }
      if ($process_pn == 1) {
         `evselect table=$pn_evt withimageset=yes imageset=pn_img_$j.fits xcolumn=X ycolumn=Y imagebinning=binSize ximagebinsize=32 yimagebinsize=32 filtertype='expression' expression='((X,Y) in $source)'`;
      }
   }
}

# Next to godliness.

if ($process_m1) {
   `rm m1_evt_cp`;
}
if ($process_m2) {
   `rm m2_evt_cp`;
}
if ($process_pn) {
   `rm pn_evt_cp`;
}

print "Done. \n";