DawnVIR Driver

.gitignore

@@ -13,6 +13,7 @@ __pycache__/
 # Distribution / packaging
 .Python
 build/
+docs/
 develop-eggs/
 dist/
 downloads/
@@ -222,3 +223,4 @@ dmypy.json
 *.CUB
 *.img
 *.IMG
+print.prt

CHANGELOG.md

@@ -40,6 +40,7 @@ release.
 ### Added
 - Mariner10 IsisLabelNaifSpice driver, tests, and test data [#547](https://github.com/DOI-USGS/ale/pull/547)
 - DawnFC IsisLabelNaifSpice driver, tests, and test data [#567](https://github.com/DOI-USGS/ale/pull/567)
+- DawnVIR IsisLabelNaifSpice driver, tests, and test data [#566](https://github.com/DOI-USGS/ale/pull/566)
 - Updated the spiceypy version used in environment.yml [#552]
 
 ### Fixed

ale/drivers/dawn_drivers.py

@@ -1,22 +1,31 @@
-import pvl
+import re
 import spiceypy as spice
 import os
+import math
+import numpy as np
 
-from glob import glob
+from scipy.interpolate import CubicSpline
 
 import ale
 from ale.base import Driver
 from ale.base.label_isis import IsisLabel
 from ale.base.data_naif import NaifSpice
+from ale.base.data_isis import IsisSpice
 from ale.base.label_pds3 import Pds3Label
 from ale.base.type_distortion import NoDistortion
-from ale.base.type_sensor import Framer
+from ale.base.type_sensor import Framer, LineScanner
+from ale.base.data_isis import read_table_data
+from ale.base.data_isis import parse_table
+from ale.transformation import TimeDependentRotation
+from ale.transformation import ConstantRotation
 
 ID_LOOKUP = {
     "FC1" : "DAWN_FC1",
     "FC2" : "DAWN_FC2"
 }
 
+degs_per_rad = 57.2957795
+
 class DawnFcPds3NaifSpiceDriver(Framer, Pds3Label, NaifSpice, Driver):
     """
     Dawn driver for generating an ISD from a Dawn PDS3 image.
@@ -393,3 +402,306 @@ def ephemeris_center_time(self):
           center ephemeris time
         """
         return self.ephemeris_start_time + (self.exposure_duration_ms / 2.0)
+
+
+class DawnVirIsisLabelNaifSpiceDriver(LineScanner, IsisLabel, NaifSpice, NoDistortion, Driver):
+    @property
+    def instrument_id(self):
+        """
+        Returns the ID of the instrument
+
+        Returns
+        -------
+        : str
+          Name of the instrument
+        """
+        lookup_table = {'VIR': 'Visual and Infrared Spectrometer'}
+        return lookup_table[super().instrument_id]
+
+    @property
+    def sensor_name(self):
+        """
+        Returns the name of the instrument
+
+        Returns
+        -------
+        : str
+          Name of the sensor
+        """
+        return self.label["IsisCube"]["Instrument"]["InstrumentId"]
+
+    @property
+    def line_exposure_duration(self):
+      """
+      The exposure duration of the image, in seconds
+
+      Returns
+      -------
+      : float
+        Exposure duration in seconds
+      """
+      return self.label["IsisCube"]["Instrument"]["FrameParameter"][0]
+
+    @property
+    def focal_length(self):
+      return float(spice.gdpool('INS{}_FOCAL_LENGTH'.format(self.ikid), 0, 1)[0])
+
+    @property
+    def detector_center_sample(self):
+        """
+        Returns the center detector sample. Expects ikid to be defined. This should
+        be an integer containing the Naif Id code of the instrument.
+
+        Returns
+        -------
+        : float
+          Detector sample of the principal point
+        """
+        return super().detector_center_sample - 0.5
+
+    @property
+    def ikid(self):
+        """
+        Returns the Naif ID code for the instrument
+        Expects the instrument_id to be defined. This must be a string containing
+        the short name of the instrument.
+
+        Returns
+        -------
+        : int
+          Naif ID used to for identifying the instrument in Spice kernels
+        """
+        lookup_table = {
+          'VIS': -203211,
+          "IR": -203213
+        }
+        return lookup_table[self.label["IsisCube"]["Instrument"]["ChannelId"]]
+
+    @property
+    def sensor_frame_id(self):
+        """
+        Returns the FRAME_DAWN_VIR_{VIS/IR}_ZERO Naif ID code if there are no associated articulation kernels. 
+        Otherwise the original Naif ID code is returned.
+        Returns
+        -------
+        : int
+          Naif ID code for the sensor frame
+        """
+        if self.has_articulation_kernel:
+          lookup_table = {
+            'VIS': -203211,
+            "IR": -203213
+         }
+        else:
+          lookup_table = {
+            'VIS': -203221,
+            "IR": -203223
+          }
+        return lookup_table[self.label["IsisCube"]["Instrument"]["ChannelId"]]
+
+    @property
+    def housekeeping_table(self):
+        """
+        This table named, "VIRHouseKeeping", contains four fields: ScetTimeClock, ShutterStatus,
+        MirrorSin, and MirrorCos.  These fields contain the scan line time in SCLK, status of 
+        shutter - open, closed (dark), sine and cosine of the scan mirror, respectively.
+
+        Returns
+        -------
+        : dict
+          Dictionary with ScetTimeClock, ShutterStatus, MirrorSin, and MirrorCos
+        """
+        isis_bytes = read_table_data(self.label['Table'], self._file)
+        return parse_table(self.label['Table'], isis_bytes)
+
+    @property
+    def line_scan_rate(self):
+        """
+        Returns
+        -------
+        : list
+          Start lines
+        : list
+          Line times
+        : list
+          Exposure durations
+        """
+        line_times = []
+        start_lines = []
+        exposure_durations = []
+
+        line_no = 0.5
+
+        for line_midtime in self.hk_ephemeris_time:
+          if not self.is_calibrated:
+            line_times.append(line_midtime - (self.line_exposure_duration / 2.0) - self.center_ephemeris_time)
+            start_lines.append(line_no)
+            exposure_durations.append(self.label["IsisCube"]["Instrument"]["FrameParameter"][2])
+            line_no += 1
+
+        line_times.append(line_times[-1] + self.label["IsisCube"]["Instrument"]["FrameParameter"][2])
+        start_lines.append(line_no)
+        exposure_durations.append(self.label["IsisCube"]["Instrument"]["FrameParameter"][2])
+
+        return start_lines, line_times, exposure_durations
+
+    @property
+    def sensor_model_version(self):
+        """
+        Returns
+        -------
+        : int
+          ISIS sensor model version
+        """
+        return 1
+
+    @property
+    def optical_angles(self):
+        hk_dict = self.housekeeping_table
+
+        opt_angles = []
+        x = np.array([])
+        y = np.array([])
+        for index, mirror_sin in enumerate(hk_dict["MirrorSin"]):
+            shutter_status = hk_dict["ShutterStatus"][index].lower().replace(" ", "")
+            is_dark = (shutter_status == "closed")   
+
+            mirror_cos = hk_dict["MirrorCos"][index]
+
+            scan_elec_deg = math.atan(mirror_sin/mirror_cos) * degs_per_rad
+            opt_ang = ((scan_elec_deg - 3.7996979) * 0.25/0.257812) / 1000
+
+            if not is_dark:
+                x = np.append(x, index + 1)
+                y = np.append(y, opt_ang)
+
+            if not self.is_calibrated:
+                opt_angles.append(opt_ang)
+
+        cs = CubicSpline(x, y, extrapolate="periodic")
+
+        for i, opt_ang in enumerate(opt_angles):
+          shutter_status = hk_dict["ShutterStatus"][i].lower().replace(" ", "")
+          is_dark = (shutter_status == "closed")
+
+          if (is_dark):
+            if (i == 0):
+              opt_angles[i] = opt_angles[i+1]
+            elif (i == len(opt_angles) - 1):
+              opt_angles[i] = opt_angles[i-1]
+            else:
+              opt_angles[i] = cs(i+1)
+
+        return opt_angles
+
+    @property
+    def hk_ephemeris_time(self):
+
+        line_times = []
+        scet_times = self.housekeeping_table["ScetTimeClock"]
+        for scet in scet_times:
+          line_midtime = spice.scs2e(self.spacecraft_id, scet)
+          line_times.append(line_midtime)
+
+        return line_times
+
+    @property
+    def ephemeris_start_time(self):
+        """
+        Returns the starting ephemeris time of the image using the first et time from
+        the housekeeping table minus the line exposure duration / 2. 
+
+        Returns
+        -------
+        : double
+          Starting ephemeris time of the image
+        """
+        return self.hk_ephemeris_time[0] - (self.line_exposure_duration / 2)
+
+
+    @property
+    def ephemeris_stop_time(self):
+        """
+        Returns the ephemeris stop time of the image using the last et time from
+        the housekeeping table plus the line exposure duration / 2. 
+
+        Returns
+        -------
+        : double
+          Ephemeris stop time of the image
+        """
+        return self.hk_ephemeris_time[-1] + (self.line_exposure_duration / 2)
+
+    @property
+    def center_ephemeris_time(self):
+      return self.hk_ephemeris_time[int(len(self.hk_ephemeris_time)/2)]
+
+    @property
+    def is_calibrated(self):
+        """
+        Checks if image is calibrated.
+
+        Returns
+        -------
+        : bool
+        """
+        return self.label['IsisCube']['Archive']['ProcessingLevelId'] > 2
+
+    @property 
+    def has_articulation_kernel(self):
+        """
+        Checks if articulation kernel exists in pool of kernels.
+
+        Returns
+        -------
+        : bool
+        """
+        try:
+          regex = re.compile('.*dawn_vir_[0-9]{9}_[0-9]{1}.BC')
+          return any([re.match(regex, i) for i in self.kernels])
+        except ValueError:
+          return False
+
+    @property
+    def frame_chain(self):
+      frame_chain = super().frame_chain
+      frame_chain.add_edge(rotation=self.inst_pointing_rotation)
+      return frame_chain
+
+    @property
+    def inst_pointing_rotation(self):
+        """
+        Returns a time dependent instrument pointing rotation for -203221 and -203223 sensor frame ids.
+
+        Returns
+        -------
+        : TimeDependentRotation
+          Instrument pointing rotation
+        """
+        time_dep_quats = np.zeros((len(self.hk_ephemeris_time), 4))
+        avs = []
+
+        for i, time in enumerate(self.hk_ephemeris_time):
+          try:
+            state_matrix = spice.sxform("J2000", spice.frmnam(self.sensor_frame_id), time)
+          except:
+            rotation_matrix = spice.pxform("J2000", spice.frmnam(self.sensor_frame_id), time)
+            state_matrix = spice.rav2xf(rotation_matrix, [0, 0, 0])
+
+          opt_angle = self.optical_angles[i]
+
+          xform = spice.eul2xf([0, -opt_angle, 0, 0, 0, 0], 1, 2, 3)
+          xform2 = spice.mxmg(xform, state_matrix)
+
+          rot_mat, av = spice.xf2rav(xform2)
+          avs.append(av)
+
+          quat_from_rotation = spice.m2q(rot_mat)
+          time_dep_quats[i,:3] = quat_from_rotation[1:]
+          time_dep_quats[i, 3] = quat_from_rotation[0]
+
+        time_dep_rot = TimeDependentRotation(time_dep_quats, self.hk_ephemeris_time, 1, self.sensor_frame_id, av=avs)
+
+        return time_dep_rot
+
+

ale/formatters/formatter.py

@@ -2,6 +2,8 @@
 import numpy as np
 from scipy.interpolate import interp1d, BPoly
 
+import spiceypy as spice
+
 from networkx.algorithms.shortest_paths.generic import shortest_path
 
 from ale.transformation import FrameChain
@@ -140,6 +142,7 @@ def to_isd(driver):
     instrument_pointing['constant_frames'] = shortest_path(frame_chain, sensor_frame, destination_frame)
     constant_rotation = frame_chain.compute_rotation(destination_frame, sensor_frame)
     instrument_pointing['constant_rotation'] = constant_rotation.rotation_matrix().flatten()
+
     meta_data['instrument_pointing'] = instrument_pointing
 
     # interior orientation