DawnVIR Driver (#566)

* Draft DawnVIR driver

* Adds interpolate and framechain updates

* fixes line scan rate

* Fixed FrameChain

* updated gitignore

* Fixed Body Rotation

* fixed rebase deletion

* added tests

* failing tests - no alespiceroot

* fix tests

* Initial fixes for dawn virs

* Fixed tests

* addressed PR feedback

* Fixed dawn vir isd

---------

Co-authored-by: acpaquette <[email protected]>

.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