#!/usr/bin/perl -w

use DBI;
use Getopt::Long;
use Term::ReadKey;
use Geo::Shapelib qw/:all/;
use Math::Trig;

$dbase    = "HiRISE";
$user     = "";
$pass     = "";
$host     = "pirldb.lpl.arizona.edu";
$ogrph    = 0;
$outrt    = "";
$help     = 0;

$opt = GetOptions ("dbase|d=s"         => \$dbase,
                   "host=s"            => \$host,
                   "user|u=s"          => \$user,
                   "pass|p=s"          => \$pass,
		   "graphic|g"         => \$ogrph,
		   "output|o=s"        => \$outrt,
                   "help|h"            => \$help);

if (!$opt || $help) {
 print "\n";
 print "obs_geom.pl -dbase -host -user -pass -help -output \n";
 print "\n";
 print "  dbase|d    database to connect to (default is HiRISE)\n";
 print "  host       hostname where the database is stored (default is pirldb.lpl.arizona.edu)\n";
 print "  user|u     username to connect to database (not username associated with suggestions)\n";
 print "  pass|p     password to connect to database (program will prompt you if not provided here)\n";
 print "  graphic|g  Export shapefile in 'ographic latitudes\n";
 print "  output|o   Root name for output e.g. 'xx' will produce files xx_yyyymmdd.shp etc...\n";
 print "  help|h     print this help message\n";
 print "\n";
 print "This program produces an ESRI shapefile containing HiRISE observations.\n";
 print "It does this by reading the HiCAT database.";
 print "Questions? Contact: Shane Byrne - shane\@lpl.arizona.edu\n";
 die("\n");
}

if ($user eq '') {
 print "MySQL username: ";
 $user = <STDIN>;
 chomp $user;
}

if ($pass eq '') {
 print "MySQL password: ";
 ReadMode('noecho');
 $pass = ReadLine(0);
 chomp $pass;
 ReadMode('normal');
}

if ($outrt eq '') {$outrt = "obs_geom"};

$dbh = DBI->connect("dbi:mysql:database=".$dbase.";host=".$host, $user, $pass);
$dbh->{RaiseError} = 1;
$dbh->{AutoCommit} = 1;

$sqlcmd1 = "select  
 Observation_Geometry.OBSERVATION_ID
,Observation_Geometry.BIN1_SCALED_PIXEL_WIDTH    
,Observation_Geometry.EMISSION_ANGLE             
,Observation_Geometry.INCIDENCE_ANGLE            
,Observation_Geometry.PHASE_ANGLE                
,Observation_Geometry.IMAGE_CENTER_LATITUDE      
,Observation_Geometry.IMAGE_CENTER_LONGITUDE     
,Observation_Geometry.SLANT_DISTANCE             
,Observation_Geometry.NORTH_AZIMUTH              
,Observation_Geometry.SUB_SOLAR_AZIMUTH          
,Observation_Geometry.SUB_SOLAR_LATITUDE         
,Observation_Geometry.SUB_SOLAR_LONGITUDE        
,Observation_Geometry.SUB_SPACECRAFT_LATITUDE    
,Observation_Geometry.SUB_SPACECRAFT_LONGITUDE   
,Observation_Geometry.SOLAR_LONGITUDE            
,Observation_Geometry.LOCAL_TIME                 
,Observation_Geometry.RED_LOWER_LEFT_LONGITUDE   
,Observation_Geometry.RED_LOWER_LEFT_LATITUDE    
,Observation_Geometry.RED_LOWER_RIGHT_LONGITUDE  
,Observation_Geometry.RED_LOWER_RIGHT_LATITUDE   
,Observation_Geometry.RED_UPPER_RIGHT_LONGITUDE  
,Observation_Geometry.RED_UPPER_RIGHT_LATITUDE   
,Observation_Geometry.RED_UPPER_LEFT_LONGITUDE   
,Observation_Geometry.RED_UPPER_LEFT_LATITUDE    
,Observation_Geometry.COLOR_LOWER_LEFT_LONGITUDE 
,Observation_Geometry.COLOR_LOWER_LEFT_LATITUDE  
,Observation_Geometry.COLOR_LOWER_RIGHT_LONGITUDE
,Observation_Geometry.COLOR_LOWER_RIGHT_LATITUDE 
,Observation_Geometry.COLOR_UPPER_RIGHT_LONGITUDE
,Observation_Geometry.COLOR_UPPER_RIGHT_LATITUDE 
,Observation_Geometry.COLOR_UPPER_LEFT_LONGITUDE 
,Observation_Geometry.COLOR_UPPER_LEFT_LATITUDE  
FROM Observation_Geometry";

$selection = " where ID != 54 AND ID != 59 ";

@fldnames = qw/FILTER OBS_ID BIN1_RES EMISSION INCIDENC PHASE CEN_LAT CEN_LON DISTANCE NORTH_AZ SOLAR_AZ SOL_LAT SOL_LON SPC_LAT SPC_LON LSUBS LOC_TIME ARCDATA_PT/;
@fldtypes = qw/String String Double Double Double Double Double Double Double Double Double Double Double Double Double Double Double String/;

$ff   = ((3396190.0/3376200.0)*(3396190.0/3376200.0));    # Flattening to convert ocentric lats to ographic
$r2d  = 57.295779513082323;                               # Convert radians to degrees (180/pi)
@filters = ("RED","COLOR");

($tm1, $tm2, $tm3) = (localtime)[3,4,5];
$tm2 = $tm2 + 1;
if (length($tm1) == 1) { $tm1 = "0".$tm1 };
if (length($tm2) == 1) { $tm2 = "0".$tm2 };
$nm = $outrt."_".($tm3+1900).$tm2.$tm1;

$shpfile = new Geo::Shapelib {
   Name => $nm,
   Shapetype => POLYGON,
   FieldNames => [@fldnames],
   FieldTypes => [@fldtypes]
};

print "Performing MySQL queries...\n";
@res = @{  $dbh->selectall_arrayref($sqlcmd1.$selection.";") };

foreach $fltr (@filters) {
#$fltr = "RED";

print "Processing $fltr records...\n";
if ($ogrph) {print "Converting latitudes to planetographic\n"};

foreach $j (@res) { 
 @res2 = @{$j};

 if ($fltr eq "RED")   {@nroi = @res2[16..23]};
 if ($fltr eq "COLOR") {@nroi = @res2[24..31]};
 $ele  = ($#nroi + 1)/2.0;

 if ($ele >= 3) { 
  if ($ogrph) {
   for ($i = 1; $i <= $ele; $i++) { $nroi[2.0*($i-1)+1] = atan2(($ff*tan($nroi[2.0*($i-1)+1]/$r2d)),1.0)*$r2d };
  }
  
  for ($i = 1; $i <= $ele; $i++) { 
   if (($nroi[2.0*($i-1)] - 180.0) > 0) { 
    $nroi[2.0*($i-1)] = $nroi[2.0*($i-1)]-360.0;
   }
  }
  
 if ($nroi[1] - $nroi[7]) { 
   @newroi = ( $nroi[2.0*($ele-1)], $nroi[2.0*($ele-1)+1] );
   for ($i = $ele-1; $i > 0; $i--) { push @newroi,( $nroi[2.0*($i-1)], $nroi[2.0*($i-1)+1] ) };  
   @nroi = @newroi;
 }
  
  @vert=[( $nroi[0], $nroi[1], 0.0, 0.0)];
  for ($i = 2; $i <= $ele; $i++) { push @vert,[( $nroi[2.0*($i-1)], $nroi[2.0*($i-1)+1], 0.0, 0.0)] };
  push @vert,[( $nroi[0], $nroi[1], 0.0, 0.0)];

  foreach $i (@res2) { 
   if ((defined $i) eq "") { $i  = 0.0 };
  }  

  $pole ="south";
  if ($res2[5] > 0) {$pole = "north"};
  push @{$shpfile->{ShapeRecords}}, [$fltr,@res2[0..15],"arcdata/".$pole."/hirise_browse/".$res2[0]."_RED.browse.jpg"];
  push @{$shpfile->{Shapes}}, {Vertices=>[@vert]};
 }

}
}

$fout = open PRJ, "> $nm.prj";
print PRJ 'GEOGCS["GCS_Mars_2000",DATUM["D_Mars_2000",SPHEROID["Mars_2000_IAU_IAG",3396190.0,169.8944472236118]],PRIMEM["Reference_Meridian",0.0],UNIT["Degree",0.0174532925199433]]\n';

print "Closing files and database connections...\n";
if ($fout) {close PRJ};
$dbh->disconnect;
$shpfile->save();
