#!/usr/bin/perl

# Routine to reformat the PDS IMG files in the Stardust directory into FITS 
# files with detached PDS labels in a new directory.

#$SRCDIR = "../TEST";
#$NEWDIR = "../data";
#$FITHDIR = "../TEST/fithdir";

###  Need to run twice, once for each directory ******
#$SRCDIR = "../WORK/ENC_IMG_NEW";
$SRCDIR = "../WORK/NAV_IMG_NEW";
$FITHDIR = "../WORK/fithdir";
$NEWDIR = "../data";

$EXPFILE = "../WORK/exposure_corrections.tab";
$SPICEFILE1 = "../WORK/Wild2_geometry1.dat";
$SPICEFILE2 = "../WORK/Wild2_geometry2.dat";
$WINFILE = "../WORK/window.txt";

$TMPLBL = "TMP.LBL";
$TMPHDR = "TMP.HDR";
$TMPDAT = "TMP.DAT";

opendir(SRC,$SRCDIR);
foreach $file (readdir(SRC))
  { next if ($file !~ /\.IMG$/);

    run_label($file);

#    system "spltimg < $SRCDIR/$file | swpinv | fitspad -n > $TMPDAT";

    $fits = "$file";
    $fits =~ tr/A-Z/a-z/;
    $fits =~ s/\.img$/.fith/;
    $fits = "$FITHDIR/$fits";

    system "cat $TMPHDR  > $fits";
    unlink $TMPHDR;
 }
closedir(SRC);


#===========================================================================

sub run_label

  # Open image file, read it to create new PDS label and FITS header 
  #  for the FITS file.

{ local ($file) = $_[0];
  my ($img,$label,$not_done,$window);
  my ($fitsfile,$junk);

  $img = "$SRCDIR/$file";
  $label = "$file";
  $label =~ tr/A-Z/a-z/;
  $label =~ s/\.img/.lbl/;
  $label = "$NEWDIR/$label";
  $fitsfile = $file;
  $fitsfile =~ tr/A-Z/a-z/;
  $fitsfile =~ s/img/fit/;

    $sfile = "$file";
    $sfile =~ tr/A-Z/a-z/;
    $sfile =~ s/\.img$//;

  # Open the various files:

  open(IMG,$img) || die "Could not open source image pipe ($!)";
  open(LBL,"| fixlen -c > $label") || die "Could not open label pipe ($!)";
  open(HDR,"| fixlen -r 81 | mkstrm | fitspad >$TMPHDR") || 
    die "Could not open header pipe ($!)";

  # read the new exposure info from the relevant file
  open(EXPC,$EXPFILE) || die "Could not open source file";
  $not_done = 1;
  while (($line=<EXPC>) && $not_done)
    {if ($line =~ $sfile) 
        {@expc = split(/\s+/,$line);
         $expcmnd = $expc[1];
         $exposure = $expc[2];
	 $xpos = 1024-$expc[4];
	 $ypos = $expc[3];
         $not_done=0;
         } 
    }
  close (EXPC);

  # read the new spice info from the relevant file
  open(SPICE1,$SPICEFILE1) || die "Could not open source file";
  $not_done = 1;
  while (($line=<SPICE1>) && $not_done)
    {if ($line =~ $sfile) 
        {@spice1 = split(/\s+/,$line);
         $cat = $spice1[1];
         $ra = $spice1[2];
         $dec = $spice1[3];
         $twist = $spice1[4];
         $cnorth = $spice1[5];
         $enorth = $spice1[6];
         $sunvec = $spice1[7];
         $pixscl = $spice1[8];
         $phase = $spice1[9];
         $otdist = $spice1[10];
         $osdist = $spice1[11];
         $tsdist = $spice1[12];
         $not_done=0;
         } 
    }

  close (SPICE1);

  # read the new spice info from the relevant file
  open(SPICE2,$SPICEFILE2) || die "Could not open source file";
  $not_done = 1;
  while (($line=<SPICE2>) && $not_done)
    {if ($line =~ $sfile) 
        {@spice2 = split(/\s+/,$line);
         $otx  = $spice2[1];
         $oty  = $spice2[2];
         $otz  = $spice2[3];
         $otvx = $spice2[4];
         $otvy = $spice2[5];
         $otvz = $spice2[6];
         $osx  = $spice2[7];
         $osy  = $spice2[8];
         $osz  = $spice2[9];
         $osvx = $spice2[10];
         $osvy = $spice2[11];
         $osvz = $spice2[12];
         $ra1 = $spice2[13];
         $ra2 = $spice2[14];
         $ra3 = $spice2[15];
         $ra4 = $spice2[16];
         $dec1 = $spice2[17];
         $dec2 = $spice2[18];
         $dec3 = $spice2[19];
         $dec4 = $spice2[20];
         $not_done=0;
         } 
    }
  close (SPICE2);

  $nwin = 0;
  $win = (' ',' ',' ',' ',' ',' ');
# read the windows file to see if the image was windowed.  If so add the info
### note that the positions of the Windows are computed to have removed the 
### overscan regions 
  open(WINDOWS,$WINFILE) || die "Could not open source file";
  $not_done = 1;
   while (($line=<WINDOWS>) && $not_done)
    {if ($line =~ $sfile) 
     {@wins = split(/ /,$line); 
     if ($wins[1] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[1];} 
     if ($wins[2] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[2];} 
     if ($wins[3] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[3];} 
     if ($wins[4] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[4];} 
     if ($wins[5] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[5];} 
     if ($wins[6] =~ /\[/) { $nwin = $nwin+1; $win[$nwin] = $wins[6];} 
    }}
  close(WINDOWS);

  # Digging through the label file:

  $not_done = 1;
  $window = 0;
  while (($line=<IMG>) && $not_done)
    { next if ($line =~ /LABEL_RECORDS/);
      next if ($line =~ /CHECKSUM/);
      next if ($line =~ /MAXIMUM/);
      next if ($line =~ /MINIMUM/); 
      next if ($line =~ /LINE_PREFIX_BYTES/); 
      next if ($line =~ /LINE_SUFFIX_BYTES/); 
      next if ($line =~ /^ *MEAN *=/);
      next if ($line =~ /SATURATED_PIXEL_COUNT/);
      next if ($line =~ /STANDARD_DEVIATION/);
      next if ($line =~ /TWIST_ANGLE_TYPE/);
      next if ($line =~ /SAMPLE_BIT_MASK/);
      next if ($line =~ /\^IMAGE_HISTOGRAM/);
      next if ($line =~ /QUATERNION/); 

      if ($line =~ /RECORD_BYTES\s*=\s*(.*)\s/)
         {if ($1 == 1048) {$bytes = 2;}
          else {$bytes = 4;}
         $records = int(1024*1044*$bytes/2880)+2 ;
         $line =~ s/=.*$/= 2880/;}

      elsif ($line =~ /FILE_RECORDS/)
        { $line =~ s/=.*$/= $records/; }

      elsif ($line =~ /DATA_SET_ID/)
        { $line =~ s/V1\.0/V2\.0/; }

      elsif ($line =~ /PRODUCT_ID/)
        { $line =~ s/IMG/FIT/; }

      elsif ($line =~ /PRODUCT_CREATION_TIME/)
        { ($sec,$min,$hr,$mday,$mon,$year) = localtime();
	  printf LBL "PRODUCT_CREATION_TIME = %4.4i-%2.2i-%2.2iT%2.2i:%2.2i:%2.2i \n",$year+1900,$mon+1,$mday,$hr,$min,$sec;
          next;}

      elsif ($line =~ /\^IMAGE\s+ /)
        { printf LBL "^FITS_HEADER                   = (\"$fitsfile\",1)\n";
          $line =~ s/[0-9]+/("$fitsfile",2)/;
        }

      elsif ($line =~ /SAMPLE_TYPE/)
        { printf LBL "  AXIS_ORDER_TYPE              = \"FIRST_INDEX_FASTEST\"\n";
          printf LBL "  UNIT                         = \"W/m^2/sR\"\n";
          $line =~ s/PC_REAL/IEEE_REAL/;
        } 

      elsif ($line =~ /LINES /)
        { $line =~ s/1022/1024/; }

      elsif ($line =~ /LINE_SAMPLES/)
        { $line =~ s/1024/1022/; }

      elsif ($line =~ /^\s*END\s*$/)
        { $not_done = 0; }

      elsif ($line =~ /MISSING_CONSTANT/)
        { $line =~ s/=.*$/= 0.00/; }

#      elsif ($line =~ /_DISPLAY_DIRECTION/)
#        { $line =~ s/'/"/g; }                     #'

      elsif ($line =~ /START_TIME\s*=\s*(\S+)\s/)
        { $start = $1; }

      elsif ($line =~ /STOP_TIME\s*=\s*(\S+)\s/)
        { $stop = $1; }

      elsif ($line =~ /TIME_FROM_CLOSEST_APPROACH\s*= (.*)</)
        {$line =~ s/=.*/= $cat/;}

      elsif ($line =~ /TARGET_NAME\s+=\s*"(.*)"/)
        { $target_name = $1; 
          $line =~ s/N\/A/CALIBRATION FIELD/;
          $target_name =~ s/N\/A/CALIBRATION FIELD/;
        }
      elsif ($line =~ /EXPOSURE_DURATION/)
        {printf LBL "EXPOSURE_DURATION              = %8.4f \<MS\>\n",$exposure ;
         next;}

      elsif ($line =~ /MISSION_NAME/)
        { $junk = <IMG>; }  # Discard extra INSTRUMENT_HOST_NAME line 

      elsif ($line =~ /FRAME_SEQUENCE_NUMBER\s*=\s*(\S+)\s/)
        { $frame = $1; }

      elsif ($line =~ /SCAN_MIRROR_RATE/)
        { $line =~ s/SEC/S/; }

      elsif ($line =~ /PHASE_ANGLE/)
      {if ($phase < 0) {printf LBL "PHASE_ANGLE                    = \"N\/A\"\n";}
       else {printf LBL "PHASE_ANGLE                    = %7.3f \<DEG\>\n", $phase;}
         next;}

      elsif ($line =~ /SC_TARGET_POSITION_VECTOR/)
         {if ($otx == -999999999.99) 
	   {printf LBL "SC_TARGET_POSITION_VECTOR = \(\"N\/A\", \"N\/A\", \"N\/A\"\)\n";
            next;}
         else
	   {printf LBL "SC_TARGET_POSITION_VECTOR = \( %13.1f, %13.1f, %13.1f \)\n", $otx,$oty,$otz;
	    next;}}

      elsif ($line =~ /SC_TARGET_VELOCITY_VECTOR/)
         {if ($otvx == -999.9999) 
	   {printf LBL "SC_TARGET_VELOCITY_VECTOR = \(\"N\/A\", \"N\/A\", \"N\/A\"\)\n";
            next;}
         else
	   {printf LBL "SC_TARGET_VELOCITY_VECTOR = \( %13.6f, %13.6f, %13.6f \)\n", $otvx,$otvy,$otvz;
	    next;}}

      elsif ($line =~ /SC_SUN_POSITION_VECTOR/)
	 {printf LBL "SC_SUN_POSITION_VECTOR = \( %13.1f, %13.1f, %13.1f \)\n", $osx,$osy,$osz;
	  next;}

      elsif ($line =~ /SC_SUN_VELOCITY_VECTOR/)
	 {printf LBL "SC_SUN_VELOCITY_VECTOR = \( %13.6f, %13.6f, %13.6f \)\n", $osvx,$osvy,$osvz;
	  next;}

      elsif ($line =~ /TARGET_CENTER_DISTANCE/)
      {   if ($otdist > 0) {printf LBL "TARGET_CENTER_DISTANCE    = %13.1f \<KM\>\n",$otdist ;} 
	  else {printf LBL "TARGET_CENTER_DISTANCE    = \"N\/A\"\n" ;} 
          if ($tsdist > 0) {printf LBL "SOLAR_DISTANCE            = %13.1f \<KM\>\n",$tsdist ;}
	  else {printf LBL "SOLAR_DISTANCE            = \"N\/A\"\n" ;} 
          printf LBL "SPACECRAFT_SOLAR_DISTANCE = %13.1f \<KM\>\n",$osdist ;
          next;}

      elsif ($line =~ /HORIZONTAL_PIXEL_SCALE/)
      {if ($pixscl < 0) {printf LBL "HORIZONTAL_PIXEL_SCALE         = \"N\/A\"\n";}
       else {printf LBL "HORIZONTAL_PIXEL_SCALE         = %7.3f <M/PIXEL>\n", $pixscl;}
         next;}

      elsif ($line =~ /VERTICAL_PIXEL_SCALE/)
        {if ($pixscl < 0) {printf LBL "VERTICAL_PIXEL_SCALE           = \"N\/A\"\n";}
         else {printf LBL "VERTICAL_PIXEL_SCALE           = %7.3f <M/PIXEL>\n", $pixscl;}
         next;}

      elsif ($line =~ /TWIST_ANGLE/)
        { printf LBL "TWIST_ANGLE                    = %8.4f \<DEG\>\n",$twist ;
	  printf LBL "TWIST_ANGLE_TYPE               = DEFAULT\n";
          next;}

     elsif ($line =~ /SMEAR_AZIMUTH\s*=\s*(.*)/)
       { $smear = 90.+ $1; 
         if ($smear < 0.) {$smear = $smear + 360.};
         if ($smear > 180.) {$smear = $smear - 180.};
         printf LBL "SMEAR_AZIMUTH                  = %8.4f \<DEG\>\n",$smear ;
         next;}

      elsif ($line =~ /RETICLE_POINT_RA/)
	 {printf LBL "RETICLE_POINT_RA     = \( %8.4f, %8.4f, %8.4f, %8.4f \)\n", $ra1,$ra2,$ra3,$ra4 ;
	  next;}

      elsif ($line =~ /RETICLE_POINT_DECLINATION/)
	 {printf LBL "RETICLE_POINT_DECLINATION = \( %8.4f, %8.4f, %8.4f, %8.4f \)\n", $dec1,$dec2,$dec3,$dec4 ;
	  next;}
  
      elsif ($line =~ /RIGHT_ASCENSION/)
        {printf LBL "RIGHT_ASCENSION                = %8.4f \<DEG\>\n", $ra;
         next;}

      elsif ($line =~ /DECLINATION/)
        {printf LBL "DECLINATION                    = %8.4f \<DEG\>\n", $dec;
         next;}


### Adding the WINDOW information 
      elsif ($line =~ /correspond to the transposed image geometry\.\"/)
        {if ($nwin > 0) 
          { printf LBL "    correspond to the transposed image geometry.\n";
            printf LBL "    \n";
            printf LBL "    Windows are located at pixels: \n";
            $iw = 1;
            while ($iw < $nwin+1) 
               {printf LBL "    WINDOW$iw   = $win[$iw]\n";
            $iw = $iw + 1;}
            $line = "    \"\n";
        }}

      elsif ($line =~ /LINE_DISPLAY_DIRECTION/)
        { printf LBL "  LINE_DISPLAY_DIRECTION       = \"UP\"\n";
          next;}

      elsif ($line =~ /SAMPLE_DISPLAY_DIRECTION/)
        { printf LBL "  SAMPLE_DISPLAY_DIRECTION     = \"RIGHT\"\n";
          next;}

      elsif ($line =~ /^\s*OBJECT\s+=/)
        { printf LBL "OBJECT                         = FITS_HEADER\n";
          printf LBL "  HEADER_TYPE                  = \"FITS\"\n";
          printf LBL "  BYTES                        = 2880\n";
          printf LBL "  RECORDS                      = 1\n";
          printf LBL "END_OBJECT                     = FITS_HEADER\n";
          printf LBL "\n";
        }

#      elsif ($line =~ /\*/)
#        { $line =~ s/\//\n\//; }

      if ($line =~ /=\s\s+/) 
        { $line =~ s/=\s+/= /;}
          
      $line =~ s/</ </;
      $line =~ s/\s+</ </;


      printf LBL $line;

      }
  close(IMG);
  close(LBL);


      @start_time = split(/T/,$start);
      @start_time1 = split(/Z/,$start_time[1]);
      @stop_time = split(/T/,$stop);
      @stop_time1 = split(/Z/,$stop_time[1]);

  # Writing the FITS header:

  printf HDR "SIMPLE  =                    T\n";
  printf HDR "BITPIX  =                  -32\n";
  printf HDR "NAXIS   =                    2\n";
  printf HDR "NAXIS1  =                 1022\n";
  printf HDR "NAXIS2  =                 1024\n";
  printf HDR "FRAMENUM= %20.20s / Frame sequence number\n", $frame;  
  printf HDR "MISSION = \'STARDUST\'\n";
  printf HDR "INSTRUM = \'NAVCAM\'\n";
  printf HDR "FILTER  = \'OPNAV\'\n";
  printf HDR "WAVELEN =                698.8 / Central wavelength of filter (nm)\n";
  printf HDR "DATE-OBS= %-20.20s / Date of observation\n", $start_time[0];
  printf HDR "STARTTIM= %-20.20s / Start time of observation\n", $start_time1[0];
  printf HDR "STOPTIME= %-20.20s / Start time of observation\n", $stop_time1[0];
  printf HDR "TARGET  = \'$target_name\'\n";
  printf HDR "EXPOSURE= %20.2f / Exposure time (ms)\n", $exposure;
  printf HDR "EXP_CMND= %20.2f / Commanded Exposure time (ms)\n", $expcmnd;
  printf HDR "RA      = %20.20s / Right Ascension (deg)\n",$ra;
  printf HDR "DEC     = %20.20s / Declination (deg)\n",$dec;
  if ($phase > 0)
  {printf HDR "PHASE   = %20.3f / Solar phase angle (deg)\n",$phase;}
  if ($otdist > 0)
  {printf HDR "RANGE   = %20.3f / Spacecraft to target distance (km)\n",$otdist;}
  printf HDR "CATIME  = %20.20s / Time from closest approach (days)\n", $cat;
  if ($pixscl > 0.0)
  {printf HDR "PIXSCALE= %20.3f / Pixel scale (m/pix)\n",$pixscl;}
  printf HDR "TWISTANG= %20.3f / Twist angle (deg)\n", $twist;
  printf HDR "CELESTN = %20.3f / Celestial north clock angle (deg)\n", $cnorth;
  printf HDR "ECLIPTN = %20.3f / Ecliptic north clock angle (deg)\n", $enorth;
  if ($sunvec > 0)
  {printf HDR "SUNVEC  = %20.3f / Sunward vector clock angle (deg)\n", $sunvec;}
  if ($xpos !~ -999)
  {printf HDR "X_POS   = %20.3f / Approximate X position of the target\n",$xpos};
  if ($xpos !~ -999)
  {printf HDR "Y_POS   = %20.3f / Approximate Y position of the target\n",$ypos};
  if ($nwin > 0)  {$iw = 1;
     while ($iw < $nwin+1) {
       printf HDR "WINDOW$iw = %-20.20s / Pixels defining window #$iw\n", $win[$iw] ;
		   $iw = $iw +1;}};
   printf HDR "END\n";
  close(HDR);

  # That's it for the label.

  return;
  }

#=============================================================================