I-ALiRT - Swapi Rolling Average (#2608)

Author: laspsandoval

imap_processing/ialirt/l0/process_swapi.py

@@ -115,6 +115,60 @@ def optimize_pseudo_parameters(
     return sol[0]
 
 
+def geometric_mean(
+    swapi_met_list: list,
+    pseudo_speed_list: list,
+    pseudo_proton_density_list: list,
+    pseudo_proton_temperature_list: list,
+) -> tuple:
+    """
+    Find moving geometric mean of SWAPI data.
+
+    Parameters
+    ----------
+    swapi_met_list : list
+        Mission elapsed time for first time measurement of SWAPI sweep.
+    pseudo_speed_list : list
+        Pseudo speed for SWAPI sweep.
+    pseudo_proton_density_list : list
+        Pseudo proton density for SWAPI sweep.
+    pseudo_proton_temperature_list : list
+        Pseudo proton temperature for SWAPI sweep.
+
+    Returns
+    -------
+    avg_swapi_met : float
+        Average swapi met value.
+    avg_proton_density : float
+        Average proton density value.
+    avg_pseudo_speed : float
+        Average pseudo speed value.
+    avg_proton_temperature : float
+        Average proton temperature value.
+    """
+    met_arr = np.asarray(swapi_met_list)
+
+    # If any of the values are equal to nan then do not include that index.
+    valid = (
+        ~np.isnan(pseudo_speed_list)
+        & ~np.isnan(pseudo_proton_density_list)
+        & ~np.isnan(pseudo_proton_temperature_list)
+    )
+
+    pseudo_speed_arr = np.asarray(pseudo_speed_list)[valid]
+    avg_pseudo_speed = np.exp(np.mean(np.log(pseudo_speed_arr)))
+
+    density_arr = np.asarray(pseudo_proton_density_list)[valid]
+    avg_proton_density = np.exp(np.mean(np.log(density_arr)))
+
+    temperature_arr = np.asarray(pseudo_proton_temperature_list)[valid]
+    avg_proton_temperature = np.exp(np.mean(np.log(temperature_arr)))
+
+    avg_swapi_met = np.mean(met_arr[valid])
+
+    return avg_swapi_met, avg_proton_density, avg_pseudo_speed, avg_proton_temperature
+
+
 def process_swapi_ialirt(
     unpacked_data: xr.Dataset, calibration_lut_table: pd.DataFrame
 ) -> list[dict]:
@@ -145,6 +199,10 @@ def process_swapi_ialirt(
     sci_dataset["met"] = met
     incomplete_groups = []
     swapi_data = []
+    pseudo_proton_speed_list = []
+    pseudo_proton_density_list = []
+    pseudo_proton_temperature_list = []
+    swapi_met_list = []
 
     # Extract energy values from the calibration lookup table file
     calibration_lut_table["timestamp"] = pd.to_datetime(
@@ -210,16 +268,44 @@ def process_swapi_ialirt(
             raw_coin_rate.squeeze(), count_rate_error.squeeze(), energy_passbands
         )
 
-        swapi_data.append(
-            _populate_instrument_header_items(met)
-            | {
-                "instrument": "swapi",
-                "swapi_epoch": int(met_to_ttj2000ns(mid_measurement)),
-                "swapi_pseudo_proton_speed": Decimal(f"{pseudo_speed:.3f}"),
-                "swapi_pseudo_proton_density": Decimal(f"{pseudo_density:.3f}"),
-                "swapi_pseudo_proton_temperature": Decimal(f"{pseudo_temperature:.3f}"),
-            }
-        )
+        pseudo_proton_speed_list.append(pseudo_speed)
+        pseudo_proton_density_list.append(pseudo_density)
+        pseudo_proton_temperature_list.append(pseudo_temperature)
+        swapi_met_list.append(mid_measurement)
+
+        # Begin averaging after 1 minute has passed (5 sweeps) and make certain that
+        # the data is sequential (~12 s cadence).
+        if len(swapi_met_list) >= 5 and np.all(
+            np.isclose(np.diff(swapi_met_list[-5:]), 12.0, atol=0.05)
+        ):
+            (
+                avg_swapi_met,
+                avg_pseudo_proton_density,
+                avg_pseudo_proton_speed,
+                avg_pseudo_proton_temperature,
+            ) = geometric_mean(
+                swapi_met_list[-5:],
+                pseudo_proton_speed_list[-5:],
+                pseudo_proton_density_list[-5:],
+                pseudo_proton_temperature_list[-5:],
+            )
+
+            swapi_data.append(
+                _populate_instrument_header_items(met)
+                | {
+                    "instrument": "swapi",
+                    "swapi_epoch": int(met_to_ttj2000ns(avg_swapi_met)),
+                    "swapi_pseudo_proton_speed": Decimal(
+                        f"{avg_pseudo_proton_speed:.3f}"
+                    ),
+                    "swapi_pseudo_proton_density": Decimal(
+                        f"{avg_pseudo_proton_density:.3f}"
+                    ),
+                    "swapi_pseudo_proton_temperature": Decimal(
+                        f"{avg_pseudo_proton_temperature:.3f}"
+                    ),
+                }
+            )
     if incomplete_groups:
         logger.info(
             f"The following swapi groups were skipped due to "

imap_processing/tests/external_test_data_config.py

@@ -126,6 +126,7 @@
     ("iois_1_packets_2025_284_05_53_38", "ialirt/data/l0/"),
     ("iois_1_packets_2025_284_05_54_39", "ialirt/data/l0/"),
     ("iois_1_packets_2025_344_05_57_56", "ialirt/data/l0/"),
+    ("iois_1_packets_2025_344_05_59_58", "ialirt/data/l0/"),
     ("imap_recon_od005_20250925_20251014_v01.bsp", "spice/test_data/"),
     ("imap_2025_283_2025_284_001.ah.bc", "spice/test_data/"),
 

imap_processing/tests/ialirt/unit/conftest.py

@@ -22,19 +22,16 @@ def sc_packet_path():
 @pytest.fixture
 def swapi_postlaunch_sc_packet_path():
     """Returns the spacecraft packet directory."""
-    packet_path = (
-        imap_module_directory
-        / "tests"
-        / "ialirt"
-        / "data"
-        / "l0"
-        / "iois_1_packets_2025_344_05_57_56"
-    )
     xtce_ialirt_path = (
         imap_module_directory / "ialirt" / "packet_definitions" / "ialirt.xml"
     )
 
-    return packet_path, xtce_ialirt_path
+    directory = imap_module_directory / "tests" / "ialirt" / "data" / "l0"
+    filenames = [
+        "iois_1_packets_2025_344_05_57_56",
+        "iois_1_packets_2025_344_05_59_58",
+    ]
+    return tuple(directory / fname for fname in filenames), xtce_ialirt_path
 
 
 @pytest.fixture

imap_processing/tests/ialirt/unit/test_process_swapi.py

@@ -3,10 +3,12 @@
 import numpy as np
 import pandas as pd
 import pytest
+import xarray as xr
 
 from imap_processing import imap_module_directory
 from imap_processing.ialirt.l0.process_swapi import (
     count_rate,
+    geometric_mean,
     optimize_pseudo_parameters,
     process_swapi_ialirt,
 )
@@ -257,16 +259,114 @@ def test_optimize_parameters():
             )
 
 
+def test_geometric_mean():
+    """Test geometric_mean function."""
+
+    swapi_met_list = [12, 24, 36, 48, 60]
+
+    pseudo_proton_speed_list = [400, 420, 440, 460, 480]
+    pseudo_proton_density_list = [5.0, 6.0, 7.0, 8.0, 9.0]
+    pseudo_proton_temperature_list = [60000, 62000, 64000, 66000, 68000]
+
+    avg_swapi_met, avg_density, avg_speed, avg_temperature = geometric_mean(
+        swapi_met_list,
+        pseudo_proton_speed_list,
+        pseudo_proton_density_list,
+        pseudo_proton_temperature_list,
+    )
+
+    expected_density = np.exp(np.mean(np.log(pseudo_proton_density_list)))
+    expected_speed = np.exp(np.mean(np.log(pseudo_proton_speed_list)))
+    expected_temperature = np.exp(np.mean(np.log(pseudo_proton_temperature_list)))
+    expected_met = np.mean(swapi_met_list)
+
+    assert np.isclose(avg_density, expected_density)
+    assert np.isclose(avg_speed, expected_speed)
+    assert np.isclose(avg_temperature, expected_temperature)
+    assert np.isclose(avg_swapi_met, expected_met)
+
+
+def test_geometric_mean_nan():
+    """Test geometric_mean function."""
+
+    swapi_met_list = [12, 24, 36, 48, 60]
+
+    pseudo_proton_speed_list = [400, 420, 440, 460, np.nan]
+    pseudo_proton_density_list = [5.0, 6.0, 7.0, 8.0, np.nan]
+    pseudo_proton_temperature_list = [60000, 62000, 64000, 66000, np.nan]
+
+    avg_swapi_met, avg_density, avg_speed, avg_temperature = geometric_mean(
+        swapi_met_list,
+        pseudo_proton_speed_list,
+        pseudo_proton_density_list,
+        pseudo_proton_temperature_list,
+    )
+
+    expected_density = np.exp(np.mean(np.log(pseudo_proton_density_list[0:4])))
+    expected_speed = np.exp(np.mean(np.log(pseudo_proton_speed_list[0:4])))
+    expected_temperature = np.exp(np.mean(np.log(pseudo_proton_temperature_list[0:4])))
+    expected_met = np.mean(swapi_met_list[0:4])
+
+    assert np.isclose(avg_density, expected_density)
+    assert np.isclose(avg_speed, expected_speed)
+    assert np.isclose(avg_temperature, expected_temperature)
+    assert np.isclose(avg_swapi_met, expected_met)
+
+
+def test_geometric_gaps():
+    """Test geometric_mean function."""
+
+    swapi_met_list = [0, 12, 24, 36, 240, 252, 264, 272, 284]
+
+    bool_check = len(swapi_met_list) >= 5 and np.all(
+        np.isclose(np.diff(swapi_met_list[-5:]), 12.0, atol=0.05)
+    )
+    assert not bool_check
+
+    pseudo_proton_speed_list = [400, 420, 440, 460, 480, 500, 520, 540, 560]
+    pseudo_proton_density_list = [5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0]
+    pseudo_proton_temperature_list = [
+        60000,
+        62000,
+        64000,
+        66000,
+        68000,
+        70000,
+        72000,
+        74000,
+        76000,
+    ]
+
+    avg_swapi_met, avg_density, avg_speed, avg_temperature = geometric_mean(
+        swapi_met_list[4::],
+        pseudo_proton_speed_list[4::],
+        pseudo_proton_density_list[4::],
+        pseudo_proton_temperature_list[4::],
+    )
+
+    expected_density = np.exp(np.mean(np.log(pseudo_proton_density_list[4::])))
+    expected_speed = np.exp(np.mean(np.log(pseudo_proton_speed_list[4::])))
+    expected_temperature = np.exp(np.mean(np.log(pseudo_proton_temperature_list[4::])))
+    expected_met = np.mean(swapi_met_list[4::])
+
+    assert np.isclose(avg_density, expected_density)
+    assert np.isclose(avg_speed, expected_speed)
+    assert np.isclose(avg_temperature, expected_temperature)
+    assert np.isclose(avg_swapi_met, expected_met)
+
+
 @pytest.mark.external_test_data
 def test_process_spacecraft_packet(
     esa_unit_conversion_table, swapi_postlaunch_sc_packet_path
 ):
     """Tests spacecraft packet processing."""
 
     packet_path, xtce_ialirt_path = swapi_postlaunch_sc_packet_path
-    postlaunch_sc_xarray_data = packet_file_to_datasets(
-        packet_path, xtce_ialirt_path, use_derived_value=False
-    )[478]
+    xarray_data = tuple(
+        packet_file_to_datasets(packet, xtce_ialirt_path, use_derived_value=False)[478]
+        for packet in packet_path
+    )
+    postlaunch_sc_xarray_data = xr.concat(xarray_data, dim="epoch")
 
     postlaunch_sc_xarray_data["swapi_version"].data = np.full_like(
         postlaunch_sc_xarray_data["swapi_version"].data, 2
@@ -275,19 +375,4 @@ def test_process_spacecraft_packet(
         postlaunch_sc_xarray_data, esa_unit_conversion_table
     )
 
-    assert len(swapi_product) == 4
-
-    key_names = [
-        "apid",
-        "met",
-        "met_in_utc",
-        "ttj2000ns",
-        "swapi_pseudo_proton_density",
-        "swapi_pseudo_proton_speed",
-        "swapi_pseudo_proton_temperature",
-    ]
-
-    for key in key_names:
-        assert swapi_product[0][key] is not None, (
-            f"The expected attribute {key} was not filled in the result dict."
-        )
+    assert len(swapi_product) == 0