#!/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.

# Creates a PDS label for the FITS file and a header for the FITS file.

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

$SRCDIR = "../IMGFILES";
$NEWDIR = "../data";
$FITHDIR = "../IMGFILES/fithdir";

$EXPFILE = "../IMGFILES/exposure_corrections.tab";
$SPICEFILE1 = "../IMGFILES/af_geometry1.dat";
$SPICEFILE2 = "../IMGFILES/af_geometry2.dat";

$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 = "$FITHDIR/$file";
#    $fits =~ tr/A-Z/a-z/;
    $fits =~ s/\.img$/.fith/;

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


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

sub run_label

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

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

# define the arrays that hold the window corner positions
  @l1 = (0,0,0,0,0,0,0,0,0,0);
  @l2 = (0,0,0,0,0,0,0,0,0,0);
  @s1 = (0,0,0,0,0,0,0,0,0,0);
  @s2 = (0,0,0,0,0,0,0,0,0,0);

  $img = "$SRCDIR/$file";
  $label = "$NEWDIR/$file";
  $label =~ s/\.img/.lbl/;
  $fitsfile = $file;
  $fitsfile =~ s/img/fit/;

## get the filename for comparing to the geometry values being added
  $sfile = "$file";
  $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];
         $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);


  # Digging through the label file:

  $nwin = 0;
# go through the label and get the window information
### note that the positions of the Windows are computed to include the 
### overscan regions 
  $not_done = 1;
  while (($line=<IMG>) && $not_done)
      {if ($line =~ /OBJECT\s.*= WINDOW/)
        { $nwin = $nwin+1 ;
          {until ($line =~ /END_OBJECT/) {$line = <IMG>;
             if ($line =~ /FIRST_LINE\s*=\s*(.*)/) {$l1[$nwin] = $1+0}
             elsif ($line =~ /LINES\s*=\s*(.*)/) {$ln = $1};
             if ($line =~ /FIRST_LINE_SAMPLE\s*=\s*(.*)/) {$s1[$nwin] = $1+8}
             elsif ($line =~ /LINE_SAMPLES\s*=\s*(.*)/) {$sn = $1};}
	   next;};
          $l2[$nwin] = $l1[$nwin]+($ln-1);
          $s2[$nwin] = $s1[$nwin]+($sn-1);
	next;
        }

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

  close(IMG);

  open(IMG,$img) || die "Could not open source image pipe ($!)";

  $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 =~ /STANDARD_DEVIATION/);
      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/;}

      if ($line =~ /NUM_WIN\s*=\s*(.*)\s/)
         {$numwin = $1;
          next;}

      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/;
          next;} 
 
      elsif ($line =~ /LINE_SAMPLES/)
        {$line =~ s/=.*/= 1044/;}

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

      elsif ($line =~ /SAMPLE_BITS\s*=\s*(.*)\s/)
      {  if ($1 == 8) {printf LBL "  SAMPLE_TYPE        = MSB_INTEGER  \n";
	               printf LBL "  SAMPLE_BITS        = 16  \n";
                       printf LBL "  COMPRESSION_TYPE   = \"8_BIT\" \n";
                       $comp = "'8_BIT'";
	               $bpix = 16;}
	 else {printf LBL "  SAMPLE_TYPE        = MSB_INTEGER  \n";
               printf LBL "  SAMPLE_BITS        = 16 \n";
	       printf LBL "  COMPRESSION_TYPE   = \"NONE\" \n";
               $comp = "'NONE'";
               $bpix = 16;}
      next;}

#      elsif ($line =~ /LINE_PREFIX_BYTES/)
#        {$line =~ s/=.*/= $pbytes/;}
#      elsif ($line =~ /LINE_SUFFIX_BYTES/)
#        {$line =~ s/=.*/= 0/;}

      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*= (.*)</)
        { $cat_time = $1; }

      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\s*=\s*(.*)</)
        {printf LBL "EXPOSURE_DURATION    = %8.2f \<MS\>\n",$exposure ;
         next;}

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

#      elsif ($line =~ /MISSION_PHASE_NAME/)
#         { $line =~ s/CRUISE\ 2/CRUISE/;}

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

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

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

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

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

      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\s*= (.*)</)
        { $line =~ s/=.*/= $pixscl/;}

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

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

     elsif ($line =~ /SMEAR_AZIMUTH\s*=\s*(.*)</)
       { $smear = 180.- $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\s*=\s*\( (.*),(.*),(.*),(.*) \)/)
	 {printf LBL "RETICLE_POINT_RA     = \( %8.4f, %8.4f, %8.4f, %8.4f \)\n", $ra1,$ra2,$ra3,$ra4 ;
	  next;}

      elsif ($line =~ /RETICLE_POINT_DECLINATION\s*=\s*\( (.*),(.*),(.*),(.*) \)/)
	 {printf LBL "RETICLE_POINT_DECLINATION = \( %8.4f, %8.4f, %8.4f, %8.4f \)\n", $dec1,$dec2,$dec3,$dec4 ;
	  next;}
  
      elsif ($line =~ /SATURATED_PIXEL_COUNT/)
        { printf LBL "  LINE_DISPLAY_DIRECTION   = \"UP\"\n";
          printf LBL "  SAMPLE_DISPLAY_DIRECTION = \"RIGHT\"\n";
          printf LBL "  MISSING_CONSTANT    = 0\n";
      }

      elsif ($line =~ /IMAGE_HISTOGRAM/)
        { until ($line =~ /END_OBJECT/) {$line = <IMG>} ;
	   next;
       } 

### Adding the WINDOW information 
      elsif ($line =~ /non-null data\.\"/)
	{if ($nwin > 0) 
	  { printf LBL "    containing non-null data. \n";
            printf LBL "    For these images, the windows have been inserted into full-frame \n";
            printf LBL "    images, with zeros in the regions outside the windows\.\n";
            printf LBL "    The windows are located at pixels: \n";
	    $iw = 1;
            while ($iw < $nwin+1) 
               {printf LBL "    WINDOW$iw   = [$s1[$iw]:$s2[$iw],$l1[$iw]:$l2[$iw]]\n";
	    $iw = $iw + 1;}
	    $line = "    \"\n";
	}}

       elsif ($line =~ /OBJECT\s.*= WINDOW/)
         {until ($line =~ /END_OBJECT/) {$line = <IMG>;
					 next;} 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  = %20.0f\n", $bpix;
  printf HDR "NAXIS   =                    2\n";
  printf HDR "NAXIS1  =                 1044\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 "COMPRESS= %-20.20s / Data compression type used\n", $comp ;
  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.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) {
	 $windat = "\'\[$s1[$iw]\:$s2[$iw]\,$l1[$iw]\:$l2[$iw]\]\'";
       printf HDR "WINDOW$iw = %-20.20s / Pixels defining window #$iw\n", $windat ;
		   $iw = $iw +1;}};
   printf HDR "END\n";
  close(HDR);

  # That's it for the label.

  return;
  }

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