pr response

imap_processing/pointing_frame_handler.py → imap_processing/spice/pointing_frame_handler.py

@@ -7,9 +7,8 @@
 1. imap_science_0001.tf - pointing frame kernel
 2. imap_sclk_0000.tsc - spacecraft clock kernel
 3. imap_wkcp.tf - spacecraft frame kernel
-4. de430.bsp - standard SPICE planetary ephemeris kernel
-5. naif0012.tls - standard NAIF leapsecond kernel
-6. imap_spin.bc - test attitude kernel available at:
+4. naif0012.tls - standard NAIF leapsecond kernel
+5. imap_spin.bc - test attitude kernel available at:
    https://lasp.colorado.edu/galaxy/display/IMAP/Data
 These need to be placed in tests/pointing_frame/test_data.
 
@@ -31,7 +30,7 @@
 # TODO : Add multiple pointings to the pointing frame.
 
 
-def get_coverage(ck_kernel: str) -> tuple[float, float, np.ndarray]:
+def get_et_times(ck_kernel: str) -> tuple[float, float, np.ndarray]:
     """
     Create the pointing frame.
 
@@ -87,6 +86,13 @@ def average_quaternions(et_times: np.ndarray) -> tuple[np.ndarray, list[np.ndarr
     aggregate = np.zeros((4, 4))
 
     for tdb in et_times:
+        # we use a quick and dirty method here for grabbing the quaternions
+        # from the attitude kernel.  Depending on how well the kernel input
+        # data is built and sampled, there may or may not be aliasing with this
+        # approach.  If it turns out that we need to pull the quaternions
+        # directly from the CK there are several routines that exist to do this
+        # but it's not straight forward.  We'll revisit this if needed.
+
         # Rotation matrix from IMAP spacecraft frame to ECLIPJ2000.
         # https://spiceypy.readthedocs.io/en/main/documentation.html#spiceypy.spiceypy.pxform
         body_rots = spice.pxform("IMAP_SPACECRAFT", "ECLIPJ2000", tdb)
@@ -175,7 +181,7 @@ def create_pointing_frame() -> Path:
     # Furnish the kernels.
     with spice.KernelPool(kernels):
         # Get timerange for the pointing frame kernel.
-        et_start, et_end, et_times = get_coverage(str(ck_kernel[0]))
+        et_start, et_end, et_times = get_et_times(str(ck_kernel[0]))
         # Create a rotation matrix
         rotation_matrix, _ = create_rotation_matrix(et_times)
 

imap_processing/tests/pointing_frame/test_pointing_frame_handler.py → imap_processing/tests/spice/test_pointing_frame_handler.py

@@ -9,11 +9,11 @@
 import pytest
 import spiceypy as spice
 
-from imap_processing.pointing_frame_handler import (
+from imap_processing.spice.pointing_frame_handler import (
     average_quaternions,
     create_pointing_frame,
     create_rotation_matrix,
-    get_coverage,
+    get_et_times,
 )
 
 
@@ -24,12 +24,11 @@ def kernel_path(tmp_path):
     test_dir = (
         Path(sys.modules[__name__.split(".")[0]].__file__).parent
         / "tests"
-        / "pointing_frame"
+        / "spice"
         / "test_data"
     )
 
     kernels = [
-        "de430.bsp",
         "naif0012.tls",
         "imap_science_0001.tf",
         "imap_sclk_0000.tsc",
@@ -57,25 +56,26 @@ def create_kernel_list(kernel_path):
 
 @pytest.fixture()
 def et_times(create_kernel_list):
-    """Tests get_coverage function."""
+    """Tests get_et_times function."""
     kernels, ck_kernel = create_kernel_list
 
     with spice.KernelPool(kernels):
-        et_start, et_end, et_times = get_coverage(str(ck_kernel[0]))
+        et_start, et_end, et_times = get_et_times(str(ck_kernel[0]))
 
     return et_times
 
 
 @pytest.mark.xfail(reason="Will fail unless kernels in pointing_frame/test_data.")
-def test_get_coverage(create_kernel_list):
-    """Tests get_coverage function."""
+def test_get_et_times(create_kernel_list):
+    """Tests get_et_times function."""
     kernels, ck_kernel = create_kernel_list
 
     with spice.KernelPool(kernels):
-        et_start, et_end, et_times = get_coverage(str(ck_kernel[0]))
+        et_start, et_end, et_times = get_et_times(str(ck_kernel[0]))
 
     assert et_start == 802008069.184905
     assert et_end == 802094467.184905
+    assert len(et_times) == 57599
 
 
 @pytest.mark.xfail(reason="Will fail unless kernels in pointing_frame/test_data.")
@@ -89,6 +89,7 @@ def test_average_quaternions(et_times, create_kernel_list):
     # Generated from MATLAB code results
     q_avg_expected = np.array([-0.6838, 0.5480, -0.4469, -0.1802])
     np.testing.assert_allclose(q_avg, q_avg_expected, atol=1e-4)
+    assert len(z_eclip_time) == 57599
 
 
 @pytest.mark.xfail(reason="Will fail unless kernels in pointing_frame/test_data.")
@@ -120,7 +121,7 @@ def test_create_pointing_frame(monkeypatch, kernel_path, create_kernel_list):
     kernels = [str(file) for file in kernel_path.iterdir()]
 
     with spice.KernelPool(kernels):
-        et_start, et_end, et_times = get_coverage(str(ck_kernel[0]))
+        et_start, et_end, et_times = get_et_times(str(ck_kernel[0]))
 
         rotation_matrix_1 = spice.pxform("ECLIPJ2000", "IMAP_DPS", et_start + 100)
         rotation_matrix_2 = spice.pxform("ECLIPJ2000", "IMAP_DPS", et_start + 1000)
@@ -134,22 +135,20 @@ def test_create_pointing_frame(monkeypatch, kernel_path, create_kernel_list):
     )
     np.testing.assert_allclose(rotation_matrix_1, rotation_matrix_expected, atol=1e-4)
 
+    # Verify imap_dps.bc has been created.
+    assert (kernel_path / "imap_dps.bc").exists()
 
-@pytest.mark.xfail(reason="Will fail unless kernels in pointing_frame/test_data.")
-def test_z_axis(create_kernel_list):
+
+@pytest.mark.skip(reason="Plotting only")
+def test_declination_plot(create_kernel_list, et_times):
     """Tests Inertial z axis and provides visualization."""
     kernels, ck_kernel = create_kernel_list
 
     with spice.KernelPool(kernels):
-        et_start, et_end, et_times = get_coverage(str(ck_kernel[0]))
-
         # Converts rectangular coordinates to spherical coordinates.
         q_avg, z_eclip_time = average_quaternions(et_times)
-        z_avg_expected = spice.q2m(list(q_avg))[:, 2]
         _, z_avg = create_rotation_matrix(et_times)
 
-        assert np.array_equal(z_avg, z_avg_expected)
-
         # Create visualization
         declination_list = []
         for time in z_eclip_time:
@@ -169,7 +168,7 @@ def test_z_axis(create_kernel_list):
         np.full(len(et_times), avg_declination * 180 / np.pi),
         "-r",
         linewidth=2,
-        label="mean z-axis for pointing frame",
+        label="mean declination for pointing frame",
     )
     plt.xlabel("Ephemeris Time")
     plt.ylabel("Spacecraft Spin Axis Declination")