#! /bin/sh
# This is the LHEA perl script: /cvmfs/extras-fp7.egi.eu/extras/heasoft/hitomi/x86_64-unknown-linux-gnu-libc2.19-0/bin/hxiplot.pl
# The purpose of this special block is to make this script work with
# the user's local perl, regardless of where that perl is installed.
# The variable LHEAPERL is set by the initialization script to
# point to the local perl installation.
#-------------------------------------------------------------------------------
eval '
if [ "x$LHEAPERL" = x ]; then
  echo "Please run standard LHEA initialization before attempting to run /cvmfs/extras-fp7.egi.eu/extras/heasoft/hitomi/x86_64-unknown-linux-gnu-libc2.19-0/bin/hxiplot.pl."
  exit 3
elif [ "$LHEAPERL" = noperl ]; then
  echo "During LHEA initialization, no acceptable version of Perl was found."
  echo "Cannot execute script /cvmfs/extras-fp7.egi.eu/extras/heasoft/hitomi/x86_64-unknown-linux-gnu-libc2.19-0/bin/hxiplot.pl."
  exit 3
elif [ `$LHEAPERL -v < /dev/null 2> /dev/null | grep -ic "perl"` -eq 0 ]; then
  echo "LHEAPERL variable does not point to a usable perl."
  exit 3
else
  # Force Perl into 32-bit mode (to match the binaries) if necessary:
  if [ "x$HD_BUILD_ARCH_32_BIT" = xyes ]; then
    if [ `$LHEAPERL -V 2> /dev/null | grep -ic "USE_64_BIT"` -ne 0 ]; then
      VERSIONER_PERL_PREFER_32_BIT=yes
      export VERSIONER_PERL_PREFER_32_BIT
    fi
  fi
  exec $LHEAPERL -x $0 ${1+"$@"}
fi
'
if(0);
# Do not delete anything above this comment from an installed LHEA script!
#-------------------------------------------------------------------------------
#!/usr/bin/perl
#
# File name: ahplot_hxi.pl
# Author: M. S. Dutka NASA GSFC
# $Date: 2016/05/23 19:42:20 $
# Version: 0
#
# Create a region file and run the extractor tool to create an SXS
# image, lightcurve and spectrum file
#
# Tool Dependencies:
# 
# Library Dependencies:
#   gen/lib/perl/ahlog
#   gen/lib/perl/ahapp
#   gen/lib/perl/ahgen
#
# Modification History:
#

# Set up
use strict;
use warnings;

use ahlog ;
use ahapp ;
use ahgen qw (:ALL) ;
use File::Copy;

# turn on AUTOFLUSH
$| = 1;

#########################
# Variable Definitions 
#########################

my $ehkfile    =""; # ehkfile file
my $infile1    =""; # filtered event file
my $infile2    =""; # unfiltered event file
my $outfile    =""; # Output file
my $binsize    =0;  # Bin size in seconds for LC
my $offset     =""; # Offset X axis to 0?
my $maxpts     =""; # Maximum points per graph
my $pltcmd     =""; # Any legal PLT command
my $clobber    =""; # Overwrite existing files

my $telescop   ="HITOMI";

my $ahploterror = 0;  # ahplot exit status

#########################
#  Main Code Block 
#########################

ahapp::startup () ;

ahapp::begin_processing () ;

$ahploterror = get_parameters () ;
unless ( $ahploterror == 0 ) {
  ah_info "HIGH", ahapp::write_parameters () ;
  ahlog::ah_err "get_parameters" ;
  ahapp::end_processing($ahploterror);
}

$ahploterror = initialize () ;
unless ( $ahploterror == 0 ) {
  ah_info "HIGH", ahapp::write_parameters () ;
  ahlog::ah_err "initialize" ;
  ahapp::end_processing($ahploterror);
}

$ahploterror = do_work () ;
unless ( $ahploterror == 0 ) {
  ahlog::ah_err "do_work" ;
  ahapp::end_processing($ahploterror);
}

$ahploterror = finalize () ;
unless ( $ahploterror == 0 ) {
  ahlog::ah_err "finalize" ;
  ahapp::end_processing($ahploterror);
}

# We're done.
ahapp::end_processing($ahploterror);


#########################
# Subroutines
#########################

sub get_parameters {

  $ehkfile        = ahapp::query_parameter("ehkfile");
  $infile1     = ahapp::query_parameter("infile1");
  $infile2    = ahapp::query_parameter("infile2");
  $binsize    = ahapp::query_parameter("binsize");
  $outfile    = ahapp::query_parameter("outfile");
  $offset     = ahapp::query_parameter("offset");
  $maxpts     = ahapp::query_parameter("maxpts");
  $pltcmd     = ahapp::query_parameter("pltcmd");
  $clobber    = $ahapp::clobber ? "yes" : "no";

  return 0;

}

# ------------------------------------------------------------------------------

sub initialize {

  my $status = 0;
  
  # Input file checking
  if ($infile1 =~ /^@/) {
    if(isRequiredFileNotFound(substr($infile1,1))) { return 1; } 
  }
  # Load the input HXI normal event files into an array
  # and make sure each file exists
  my @infiles1 = readInputFileList($infile1);
  unless(@infiles1) {
    ahgen::ah_err "No input HXI normal event files.";
    return 1;
  }
  if (isRequiredFileNotFound(@infiles1)) { return 1; }

  my @infiles2 = readInputFileList($infile2);
  unless(@infiles2) {
    ahgen::ah_err "No input HXI normal event files.";
    return 1;
  }
  if (isRequiredFileNotFound(@infiles2)) { return 1; }



  if (removeOutputFileIfClobbering($outfile,$ahapp::clobber) ) { return 1; }

  

  # Check if the pltcmd is a text file
  # If so, make a copy of it to a temporary file
  # and add plot and quit commands.
  #  
  # These are failsafes in case the user did not put
  # either plot or quit into their pco file
  if ( $pltcmd =~ /^@/ ) {
    my $pltfile    ="ahplot.pco";
    ahapp::add_temp_file($pltfile);
    if(isRequiredFileNotFound(substr($pltcmd,1))) { return 1; }
    File::Copy::copy(substr($pltcmd,1),$pltfile);
    open PLT,">>$pltfile";
    print PLT "plot\n"; # force pgplot to plot the file
    print PLT "quit\n"; # force pgplot to quit
    close PLT;
    $pltcmd = "\@$pltfile";
  }

  


  # Write all parameters to this script to the log file.
  ah_info "HIGH", ahapp::write_parameters () ;                 

  return 0;

}

# ------------------------------------------------------------------------------

sub do_work {

  my $status = 0;
  
  # Create event file for layer 0 and layer 4 from unfiltered event data
  my $filename_uf_layer0 = substr($infile2,0,-4) . "_ufa_layer0.evt";
  my $filename_uf_layer4 = substr($infile2,0,-4) . "_ufa_layer4.evt";

  # Clean up existing files if clobber is set
  if ($clobber) {
    if (-e $filename_uf_layer0) {
      unlink $filename_uf_layer0;
    }
    if (-e $filename_uf_layer4) {
      unlink $filename_uf_layer4;
    }
  }

  $status = ahgen::run_ftool("ftselect","infile=$infile2","outfile=$filename_uf_layer0","expression=LAYER==0");
  if ( $status ) {
    ahlog::ah_err "ftselect failed.";
    return $status;
  }

  $status = ahgen::run_ftool("ftselect","infile=$infile2","outfile=$filename_uf_layer4","expression=LAYER==4");
  if ( $status ) {
    ahlog::ah_err "ftselect failed.";
    return $status;
  }

  # run extractor to create lightcurves for filtered events and unfilters layer 0 and layer 4
  my $filename = $infile1 . "[PI=     0:  2047]";
  my $lcnamefilter = substr($infile1,0,-4) . ".lc";
  if ($clobber) {
    if (-e $lcnamefilter) {
      unlink $lcnamefilter;
    }
  }
  $status = ahgen::run_ftool("extractor","filename=$filename", "eventsout=NONE",
                   "regionfile=NONE", "qdpfile=NONE", "fitsbinlc=$lcnamefilter",
                   "lcthresh=0", "lctzero=yes", "xronwn=NONE", "unbinlc=NONE", 
                   "phafile=NONE", "imgfile=NONE", "timefile=NONE", 
                   "adjustgti=no", "gstring=NONE", "timeorder=no", "xcolf=X",
                   "ycolf=Y", "xcolh=NONE", "ycolh=NONE", "xfkey=TLMAX",  
                   "yfkey=TLMAX", "xhkey=TLMAX", "yhkey=TLMAX", "phamax=TLMAX",
                   "specbin=1", "binh=1", "binf=1", "binlc=$binsize", 
                   "tcol=TIME", "ecol=PI", "ccol=NONE", "gcol=NONE", "gti=GTI",
                   "events=EVENTS", "gtitxt=NONE", "timeref=40000.00", 
                   "wtmapb=no", "gtinam=GTI", "exitnow=no");
  if ( $status ) {
    ahlog::ah_err "extractor failed.";
    return $status;
  }

  $filename = $filename_uf_layer0 . "[PI=     0:  2047]";
  my $lcnamelayer0 = substr($filename_uf_layer0,0,-4) . ".lc";
  if ($clobber) {
    if (-e $lcnamelayer0) {
      unlink $lcnamelayer0;
    }
  }
  $status = ahgen::run_ftool("extractor","filename=$filename", "eventsout=NONE",
                   "regionfile=NONE", "qdpfile=NONE", "fitsbinlc=$lcnamelayer0",
                   "lcthresh=0", "lctzero=yes", "xronwn=NONE", "unbinlc=NONE", 
                   "phafile=NONE", "imgfile=NONE", "timefile=NONE", 
                   "adjustgti=no", "gstring=NONE", "timeorder=no", "xcolf=X",
                   "ycolf=Y", "xcolh=NONE", "ycolh=NONE", "xfkey=TLMAX",  
                   "yfkey=TLMAX", "xhkey=TLMAX", "yhkey=TLMAX", "phamax=TLMAX",
                   "specbin=1", "binh=1", "binf=1", "binlc=$binsize", 
                   "tcol=TIME", "ecol=PI", "ccol=NONE", "gcol=NONE", "gti=GTI",
                   "events=EVENTS", "gtitxt=NONE", "timeref=40000.00", 
                   "wtmapb=no", "gtinam=GTI", "exitnow=no");
  if ( $status ) {
    ahlog::ah_err "extractor failed.";
    return $status;
  }


  $filename = $filename_uf_layer4 . "[PI=     0:  2047]";
  my $lcnamelayer4 = substr($filename_uf_layer4,0,-4) . ".lc";
  if ($clobber) {
    if (-e $lcnamelayer4) {
      unlink $lcnamelayer4;
    }
  }
  $status = ahgen::run_ftool("extractor","filename=$filename", "eventsout=NONE",
                   "regionfile=NONE", "qdpfile=NONE", "fitsbinlc=$lcnamelayer4", 
                   "lcthresh=0", "lctzero=yes", "xronwn=NONE", "unbinlc=NONE", 
                   "phafile=NONE", "imgfile=NONE", "timefile=NONE", 
                   "adjustgti=no", "gstring=NONE", "timeorder=no", "xcolf=X",
                   "ycolf=Y", "xcolh=NONE", "ycolh=NONE", "xfkey=TLMAX",  
                   "yfkey=TLMAX", "xhkey=TLMAX", "yhkey=TLMAX", "phamax=TLMAX",
                   "specbin=1", "binh=1", "binf=1", "binlc=$binsize", 
                   "tcol=TIME", "ecol=PI", "ccol=NONE", "gcol=NONE", "gti=GTI",
                   "events=EVENTS", "gtitxt=NONE", "timeref=40000.00", 
                   "wtmapb=no", "gtinam=GTI", "exitnow=no");
  if ( $status ) {
    ahlog::ah_err "extractor failed.";
    return $status;
  }
  
  # Run lcurve to put all the lightcurves (filtered events, layer 0 unfilterd and layer 4 unfilterd) 
  # into one file that is readable by fplot
  # Get the column number for the EHK column SAA_SXI
  my $ehkcol = ahgen::get_column_num($ehkfile,"EHK","SAA_SXI");
  unless ( $ehkcol ) { 
    ahgen::ah_err "Could not find column SAA_SXI in file $ehkfile\[EHK]";
    return 1;
  }

  # Check the EHK file for the TIMEDEL keyword
  my $ehk_timedel = ahgen::get_keyword($ehkfile,"EHK","TIMEDEL");
  if( ahgen::get_error_flag ) { 
    ahlog::ah_err "TIMEDEL keyword not defined in $ehkfile"; 
    return ahgen::get_error_flag;
  }
  $status = ahgen::run_ftool("lcurve","nser=4","cfile4=$ehkfile vy$ehkcol", 
                             "cfile3=$lcnamefilter",
                             "cfile2=$lcnamelayer4", 
                             "cfile1=$lcnamelayer0",
                             "window=-","dtnb=$binsize","nbint=900","plot=N",
                             "plotdev=/XW","plotdnum=4","outfile=out","clobber=$clobber");
  if ( $status ) {
    ahlog::ah_err "lcurve failed";
  }
  
  my $fplot_infile = "out.flc[RATE][col TIME;Layer0=RATE1;Layer4=RATE2;Filtered=RATE3;SAA=RATE4;]";
   
  # Run FPLOT to create the postscript file
  $status = ahgen::run_ftool("fplot",
                             "infile=$fplot_infile",
                             "maxpts=$maxpts",
                             "offset=$offset",
                             "xparm=TIME",
                             "yparm=Layer0 Layer4 Filtered SAA",
                             "rows=-",
                             "device=$outfile/ps",
                             "pltcmd=$pltcmd",
                             "binmode=DEFAULT",
                             "sensecase=no"
                           );

  if ( $status ) {
    ahlog::ah_err "fplot failed.";
    return $status;
  }

  return 0;

}

# ------------------------------------------------------------------------------

sub finalize {

  my $status = 0;

  return 0;
  
}

# ------------------------------------------------------------------------------