#!/usr/bin/perl -w

#
#  Calculate monthly mound credits from Harlan County Reservoir to Guide Rock
#  This program reads the annual SFI files and produces a table in HTML
#
#  Willem A. Schreuder
#  August 4, 2006
#

use strict;

(@ARGV==2) || die "Usage: $0 <first year> <last year>\n";
my $yr0 = shift @ARGV;
my $yr1 = shift @ARGV;

my %name = (1=>'Jan',2=>'Feb',3=>'Mar',4=>'Apr',5=>'May',6=>'Jun',7=>'Jul'=>8=>'Aug',9=>'Sep',10=>'Oct',11=>'Nov',12=>'Dec' , ANN=>'Annual');
my @per = (1 .. 12 , 'ANN');

#
#  Read FORTRAN record
#
sub get
{
   my ($fd) = @_;
   my ($buf,$tail);
   (read($fd,$buf,4) == 4) || die "Error reading record header\n";
   my $len = unpack('L',$buf);
   (read($fd,$buf,$len) == $len) || die "Error reading record (length=$len)\n";
   (read($fd,$tail,4) == 4) || die "Error reading record tail\n";
   (unpack('L',$tail)==$len) || die "Tail length mismatch $len\n";
   return $buf;
}

#  Read stream flow files for every year
my %vol;  #  Monthly volumes
foreach my $yr ($yr0 .. $yr1)
{
   #  Set gage IDs
   my $SIGR = ($yr<2009) ? 'SI253006RRGuideRock' : 'SI251003RRGuideRock';
   my $SIBH = 'SI234002RRBloHarlan';
   #  Read historical and mound data files
   foreach my $run ('','d')
   {
      #  Set sign to subtract mound run from historical run
      my $sign = $run ? -1 : +1;
      #  Open stream flow file
      open(SFI , "<../$yr/$yr$run.sfi") || die "Cannot open file ../$yr/$yr$run.sfi\n";
      binmode(SFI);
      #  Read number of variables
      my $buf = get(*SFI);
      my ($n) = unpack('L*',$buf);
      #  Read variable names
      $buf = get(*SFI);
      my %var;
      foreach my $k (0 .. $n-1)
      {
         my $var = substr($buf,20*$k+4,20);
         $var =~ s/ *$//;
         $var{$var} = $k;
      }
      #  Gage variable offsets
      exists($var{$SIGR}) || die "Missing gage $SIGR in../$yr/$yr$run.sfi\n";
      my $GR = $var{$SIGR};
      exists($var{$SIBH}) || die "Missing gage $SIBH in../$yr/$yr$run.sfi\n";
      my $BH = $var{$SIBH};
      #  Read initial record for time only
      $buf = get(*SFI);
      my $t0 = unpack('d*',$buf);
      #  Read data for all the months in the year
      foreach my $mo (1 .. 12)
      {
         #  Two time steps per month
         foreach my $k (1,2)
         {
            #  Read entire record
            $buf = get(*SFI);
            my ($t1,@q) = unpack('d*',$buf);
            #  Integrate volume and accumulate
            my $vol = $sign*($t1-$t0)*($q[$GR]-$q[$BH])/43560;
            $vol{$yr}{$mo} += $vol;
            $vol{$yr}{ANN} += $vol;
            #  Update start time
            $t0 = $t1;
         }
      }
   }
}
close(SFI);

#
#  Output volumes to table
#
open(HTM , '>hc2gr.htm') || die "Cannot open file hc2gr.htm\n";
print HTM "<html><table border>\n";
#  Header
print HTM "<tr><th colspan=14>Harlan County to Guide Rock Monthly Mound Credits (acre-feet)</th></tr>\n";
print HTM "<tr><th>Year</th>".join('',map {"<th>$name{$_}</th>"} @per)."</tr>\n";
#  Loop over years
foreach my $yr (sort keys %vol)
{
   print HTM "<tr align=right><td align=center>$yr</td>".join('',map {sprintf '<td>%.1f</td>' , $vol{$yr}{$_}} @per)."</tr>\n";
}
print HTM "</table></html>\n";
close(HTM);
