#!/usr/bin/perl
#
# Make plots showing the properties of the ensemble photometric solution
#
$debug = 1;

$aper = "";

# subtract this from all JD values when plotting
$jd_base = 2455000.0;

$filter = "I";
$filter_lc = lc($filter);

$solvepht_out = "pg_solve" . $aper . "_" . $filter_lc . ".out";
$solvepht_sig = "pg_solve" . $aper . "_" . $filter_lc . ".sig";
$solvepht_img = "pg_solve" . $aper . "_" . $filter_lc . ".img";

# use this if the solvepht output doesn't start at 0.0 for the brightest stars
#$mag_offset = 40.491;
$mag_offset = 0.00;

$title = "'PG1633 field, " . $filter . " " . "3-pix aper diff mag, 12-inch RIT Obs Aug 27 2011 UT'";

# make the Postscript and PNG graphs of ensemble magsig vs. mag
sig_plot("pg_magsig" . $aper . "_" . $filter_lc . ".ps", "postscript color");

# make Postscript and PNG graphs of image adjustment vs. JD
em_plot("pg_emplot" . $aper . "_" . $filter_lc . ".ps", "postscript color");

exit 0;



############################################################################
# Make a plot showing "sigma" (stdev from mean mag) vs "mag" in solution.
#
# create a file with GNUPLOT commands, then execute GNUPLOT to read
#   from that file.
#
# usage:  sig_plot   outfile_file_name  set_term_options
#
sub sig_plot {
  
  my($output_file, $term_options);

  # get the arguments
  $output_file = $_[0];
  $term_options = $_[1];

  $cmdfile = "gnuplot.in";

  # determine the limits of the graph
  $minmag = -0.5;
  $maxmag = 8.0;
  #$minmag = 10;
  #$maxmag = 16.0;
  $minsig = -0.02;
  $maxsig = 0.30;
  
  open (CMDFILE, ">$cmdfile") || die("can't open $cmdfile for writing");
  printf CMDFILE "set output '$output_file' \n";
  printf CMDFILE "set term $term_options \n";
  printf CMDFILE "set grid \n";
  printf CMDFILE "set nokey \n";
  printf CMDFILE "set xlabel 'differential magnitude' \n";
  printf CMDFILE "set ylabel 'stdev from mean mag' \n";
  printf CMDFILE "set title $title\n";
  
  $cmd = "plot [$minmag:$maxmag][$minsig:$maxsig] '$solvepht_sig' using (\$2-$mag_offset):3 pt 6 ps 1.0 lw 2";
  printf CMDFILE "$cmd \n";
  printf CMDFILE "set output \n";
  printf CMDFILE "quit \n";
  close(CMDFILE) ;
  
  $retval = `gnuplot < $cmdfile`;
  printf "retval is ..%s..\n", $retval;

  if ($term_options =~ /postscript/) {
    $psfile = $output_file;
    $giffile = $psfile;
    $giffile =~ s/.ps/.png/;
    $cmd = "convert -rotate 90 $psfile $giffile";
    printf "cmd is ..$cmd.. \n";
    $ret = `$cmd`;
  }

  
}
  
############################################################################
# Make a plot showing "em" (image adjustment) vs Juliand Date
#
# create a file with GNUPLOT commands, then execute GNUPLOT to read
#   from that file.
#
# usage:  em_plot   outfile_file_name  set_term_options
#
sub em_plot {
  
  my($output_file, $term_options);

  # get the arguments
  $output_file = $_[0];
  $term_options = $_[1];

  $cmdfile = "gnuplot.in";

  # determine the limits of the graph
  $solveimg_file = "$solvepht_img";
  ($minjd, $maxjd) = get_minmax($solveimg_file, 3);
  $minjd -= $jd_base;
  $maxjd -= $jd_base;
  if ($debug > 0) {
    printf "em_plot: min JD %12.7f   max JD %12.7f \n", $minjd, $maxjd;
  }

  ($minem, $maxem) = get_minmax($solveimg_file, 4);
  if ($debug > 0) {
    printf "em_plot: min em %7.3f    max em %7.3f \n", $minem, $maxem;
  }
  # place absolute limits on the values of the adjustments, to avoid
  #   very distant outliers from throwing off the entire scale of the plot
if (0 == 1) {
  if ($minem < -10) {
    $minem = -10;
  }
  if ($maxem > 10) {
    $maxem = 10;
  }
}

  $minx = $minjd - ($maxjd - $minjd)*0.05;
  $maxx = $maxjd + ($maxjd - $minjd)*0.05;
  $miny = $minem - ($maxem - $minem)*0.05;
  $maxy = $maxem + ($maxem - $minem)*0.05;
  
  open (CMDFILE, ">$cmdfile") || die("can't open $cmdfile for writing");
  printf CMDFILE "set output '$output_file' \n";
  printf CMDFILE "set term $term_options \n";
  printf CMDFILE "set grid \n";
  printf CMDFILE "set nokey \n";
  printf CMDFILE "set xlabel 'Julian Date - $jd_base ' \n";
  printf CMDFILE "set ylabel 'image adjustment to solution' \n";
  printf CMDFILE "set title $title\n";
  
  $cmd = "plot [$minx:$maxx][$miny:$maxy] '$solvepht_img' ";
  $cmd = $cmd . " using (\$4-$jd_base):5 pt 6 ps 1.5 lw 2";
  printf CMDFILE "$cmd \n";
  printf CMDFILE "set output \n";
  printf CMDFILE "quit \n";
  close(CMDFILE) ;
  
  $retval = `gnuplot < $cmdfile`;
  printf "retval is ..%s..\n", $retval;

  if ($term_options =~ /postscript/) {
    $psfile = $output_file;
    $giffile = $psfile;
    $giffile =~ s/.ps/.png/;
    $cmd = "convert -rotate 90 $psfile $giffile";
    printf "cmd is ..$cmd.. \n";
    $ret = `$cmd`;
  }

  
}
  

#########################################################################
# PROCEDURE: get_minmax
#
# DESCRIPTION: open the given file, read numbers from the given 
#              column, find the min and max values, return them
#              in a two-element list
#
# RETURNS:
#      a list (min, max)
#
sub get_minmax {
  
  my($filename, 
     $column,
	  $value,
	  $sum,
	  $num,
	  $mean,
	  $i);

  $filename = $_[0];
  $column = $_[1];

  $num = 0;
  open(INFILE, $filename) || die ("get_mean: can't open file $filename");
  while (<INFILE>) {
    @words = split(/\s+/);
    $value = $words[$column];
	 if ($debug > 1) {
	   printf "  get_mean: next value is ..%s..\n", $value;
	 }
    if ($num == 0) {
	    $min = $value;
		 $max = $value;
	 } else {
	    if ($value < $min) {
		   $min = $value;
		 }
		 if ($value > $max) {
		   $max = $value;
		 }
    }
	 $num++;
  }
  close(INFILE);

  return($min, $max);
}