Skip to content

FIX: Lo update spin time calculation #2635

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
greglucas merged 1 commit into from
Jan 29, 2026
+108 −130
Merged

FIX: Lo update spin time calculation #2635

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
81 changes: 24 additions & 57 deletions imap_processing/lo/l1b/lo_l1b.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,11 @@
get_pointing_times,
interpolate_repoint_data,
)
from imap_processing.spice.spin import get_spin_data, get_spin_number
from imap_processing.spice.spin import (
get_spin_data,
get_spin_number,
interpolate_spin_data,
)
from imap_processing.spice.time import (
epoch_to_fractional_doy,
et_to_utc,
Expand Down Expand Up @@ -166,6 +170,7 @@
"spin_cycle",
]
# -------------------------------------------------------------------
DE_CLOCK_TICK_S = 4.096e-3 # seconds per DE clock tick


def lo_l1b(
Expand Down Expand Up @@ -277,13 +282,9 @@ def l1b_de(
avg_spin_durations_per_cycle = get_avg_spin_durations_per_cycle(spin_data)
# set the spin cycle for each direct event
l1b_de = set_spin_cycle(pointing_start_met, l1a_de, l1b_de)
# get spin start times for each event
spin_start_time = get_spin_start_times(l1a_de, l1b_de, spin_data)

# get the absolute met for each event
l1b_de = set_event_met(
l1a_de, l1b_de, spin_start_time, avg_spin_durations_per_cycle
)
l1b_de = set_event_met(l1a_de, l1b_de)
# set the epoch for each event
l1b_de = set_each_event_epoch(l1b_de)
# Set the ESA mode for each direct event
Expand Down Expand Up @@ -794,7 +795,7 @@ def _check_sufficient_spins(spin_data: xr.Dataset) -> np.ndarray:


def get_spin_start_times(
l1a_de: xr.Dataset, l1b_de: xr.Dataset, spin_data: xr.Dataset
l1a_de: xr.Dataset,
) -> xr.DataArray:
"""
Get the start time for the spin that each direct event is in.
Expand All @@ -807,59 +808,31 @@ def get_spin_start_times(
----------
l1a_de : xr.Dataset
The L1A DE dataset.
l1b_de : xr.Dataset
The L1B DE dataset.
spin_data : xr.Dataset
The L1A Spin dataset.

Returns
-------
spin_start_time : xr.DataArray
spin_start_time : np.ndarray
The start time for the spin that each direct event is in.
"""
# align l1a_de packets with spin_data packets
de_to_spin_indices = match_science_to_spin_asc(
l1a_de["epoch"].values, spin_data["epoch"].values
)
# Repeat this for each direct event based on the de_count
de_to_spin_indices = np.repeat(de_to_spin_indices, l1a_de["de_count"])

# There are 28 spins per epoch (1 aggregated science cycle)
# Set the spin_cycle_num to the spin number relative to the
# start of the ASC
spin_cycle_num = l1b_de["spin_cycle"] % 28
asc_starts = spin_data["acq_start_sec"] + spin_data["acq_start_subsec"] * 1e-6
avg_spin_durations = get_avg_spin_durations_per_cycle(spin_data)
# Get the actual spin start times from the spin data
# Use the individual spin start times rather than calculating from ASC averages
spin_start_times = interpolate_spin_data(l1a_de["shcoarse"].values)[
"spin_start_met"
].values

# Calculate the time based off of the start of the acquisition period
# then using an average spin duration to calculate the offset within the ASC
# NOTE: We don't want to use the spin start times directly from the ASC spin packet
# because we are using an average spin_cycle for the ASC and there are
# times when only half the spins were completed in an ASC. This allows us
# to still get a valid spin_cycle start time for each direct event, even
# if the average spin_cycle was after the end of the acquisition period.
# TODO: Can we do even better by knowing how many esa_steps and spins were complete?
# i.e. change the spin_cycle calculation
spin_start_time = (
asc_starts[de_to_spin_indices]
+ spin_cycle_num * avg_spin_durations[de_to_spin_indices]
)
return spin_start_time
return spin_start_times


def set_event_met(
l1a_de: xr.Dataset,
l1b_de: xr.Dataset,
spin_start_time: xr.DataArray,
avg_spin_durations: xr.DataArray,
) -> xr.Dataset:
"""
Get the event MET for each direct event.

Each direct event is converted from a data number to engineering unit in seconds.
de_eu_time de_dn_time / 4096 * avg_spin_duration
where de_time is the direct event time Data Number (DN) and avg_spin_duration
is the average spin duration for the ASC that the event was measured in.
time_from_start_of_spin = de_time * DE_CLOCK_TICK_S
Copy link
Collaborator Author

@greglucas greglucas Jan 29, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we be adding a half-tick here?
(de_time + 0.5) * DE_CLOCK_TICK_S

Copy link
Contributor

@subagonsouth subagonsouth Jan 29, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For Hi, each DE is tagged with the exact clock time (not adding 1/2 of a tick). When computing spin-phase or geometry, the 1/2 tick is added.

where de_time is the direct event time Data Number (DN).

The direct event time is the time of direct event relative to the start of the spin.
The event MET is the sum of the start time of the spin and the
Expand All @@ -871,26 +844,18 @@ def set_event_met(
The L1A DE dataset.
l1b_de : xr.Dataset
The L1B DE dataset.
spin_start_time : np.ndarray
The start time for the spin that each direct event is in.
avg_spin_durations : xr.DataArray
The average spin duration for each epoch.

Returns
-------
l1b_de : xr.Dataset
The L1B DE dataset with the event MET.
"""
counts = l1a_de["de_count"].values
de_time_asc_groups = np.split(l1a_de["de_time"].values, np.cumsum(counts)[:-1])
de_times_eu = []
for i, de_time_asc in enumerate(de_time_asc_groups):
# DE Time is 12 bit DN. The max possible value is 4095
# divide by 4096 to get fraction of a spin duration
de_times_eu.extend(de_time_asc / 4096 * avg_spin_durations[i].values)
# get spin start times for each event
spin_start_times = get_spin_start_times(l1a_de)

# spin start + offset based on de_time ticks
l1b_de["event_met"] = xr.DataArray(
spin_start_time + de_times_eu,
spin_start_times + l1a_de["de_time"].values * DE_CLOCK_TICK_S,
dims=["epoch"],
# attrs=attr_mgr.get_variable_attributes("epoch")
)
Expand Down Expand Up @@ -1344,12 +1309,14 @@ def set_pointing_bin(l1b_de: xr.Dataset) -> xr.Dataset:
# first column: radius (Not needed)
# second column: longitude
lons = direction[:, 1]
# shift to 0-360 range (spin-phase 0 should be in bin 0)
lons = (lons + 360) % 360
# third column: latitude
lats = direction[:, 2]

# Define bin edges
# 3600 bins, 0.1° each
lon_bins = np.linspace(-180, 180, 3601)
lon_bins = np.linspace(0, 360, 3601)
# 40 bins, 0.1° each
lat_bins = np.linspace(-2, 2, 41)

Expand Down
150 changes: 77 additions & 73 deletions imap_processing/tests/lo/test_lo_l1b.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
from imap_processing.cdf.imap_cdf_manager import ImapCdfAttributes
from imap_processing.cdf.utils import load_cdf
from imap_processing.lo.l1b.lo_l1b import (
DE_CLOCK_TICK_S,
calculate_de_rates,
calculate_histogram_rates,
calculate_star_sensor_profile_for_group,
Expand Down Expand Up @@ -219,7 +220,9 @@ def l1a_hist():
"imap_processing.lo.l1b.lo_l1b.cartesian_to_latitudinal",
return_value=np.zeros((2000, 3)),
)
@patch("imap_processing.lo.l1b.lo_l1b.interpolate_spin_data")
def test_lo_l1b_de(
mock_interpolate_spin_data,
mock_frame_transform,
mock_lo_instrument_pointing,
mocked_get_pointing_times,
Expand All @@ -229,6 +232,15 @@ def test_lo_l1b_de(
anc_dependencies,
):
# Arrange
# Mock the spin data to provide spin start times
# Create a DataFrame covering the time range of the test data
num_events = 2000
mock_spin_df = pd.DataFrame(
{
"spin_start_met": np.ones(num_events),
}
)
mock_interpolate_spin_data.return_value = mock_spin_df

# Add l1b_nhk dependency with pivot angle information
l1b_nhk = xr.Dataset(
Expand Down Expand Up @@ -490,94 +502,82 @@ def test_spin_cycle(mock_get_spin_number):
np.testing.assert_array_equal(spin_cycle_data["spin_cycle"], spin_cycle_expected)


def test_get_spin_start_times():
@patch("imap_processing.lo.l1b.lo_l1b.interpolate_spin_data")
def test_get_spin_start_times(mock_interpolate_spin_data):
# Arrange
l1b_de = xr.Dataset(
# Mock the spin data to return specific spin start times
mock_spin_df = pd.DataFrame(
{
"spin_cycle": ("epoch", [0, 1, 2, 3, 4]),
},
coords={
"epoch": [
0,
1,
2,
3,
4,
]
},
"spin_start_met": [10.5, 10.5, 30.1, 30.1, 30.1],
}
)
mock_interpolate_spin_data.return_value = mock_spin_df

l1a_de = xr.Dataset(
{
"shcoarse": ("epoch", [0, 1]),
"shcoarse": ("epoch", [15, 35]),
"de_count": ("epoch", [2, 3]),
"met": ("epoch", [0, 1]), # MET per time epoch, not per direct event
"de_time": ("direct_event", [0000, 1000, 2000, 3000, 4000]),
"de_time": ("direct_event", [0, 1000, 2000, 3000, 4000]),
},
coords={"epoch": [0, 1], "direct_event": [0, 1, 2, 3, 4]},
)
spin = xr.Dataset(
{
"shcoarse": ("epoch", [0, 1]),
"acq_start_sec": (
"epoch",
[20, 25],
),
"acq_start_subsec": (
"epoch",
[0, 0],
),
"acq_end_sec": (
"epoch",
[25, 30],
),
"acq_end_subsec": (
"epoch",
[0, 0],
),
"num_completed": (
"epoch",
[28, 14],
),
}
)

# Expected: shcoarse 15 should match spin at index 0 (10 < 15 < 20)
# shcoarse 35 should match spin at index 2 (30 < 35 < 40)
# Repeated by de_count: [2, 3] -> [index0, index0, index2, index2, index2]
spin_start_times_expected = np.array(
[20, 20 + 5 / 28, 25 + 5 * 2 / 14, 25 + 5 * 3 / 14, 25 + 5 * 4 / 14]
[10.5, 10.5, 30.1, 30.1, 30.1] # 10 + 0.5e6*1e-6 # 30 + 0.1e6*1e-6
)
spin_start_times = get_spin_start_times(l1a_de, l1b_de, spin)

# Act
spin_start_times = get_spin_start_times(l1a_de)

# Assert
np.testing.assert_allclose(
spin_start_times,
spin_start_times_expected,
atol=1e-4,
)


def test_set_event_met():
@patch("imap_processing.lo.l1b.lo_l1b.interpolate_spin_data")
def test_set_event_met(mock_interpolate_spin_data):
# Arrange
# Mock the spin data
mock_spin_df = pd.DataFrame(
{
"spin_start_met": [10, 10, 30, 30, 30],
}
)
mock_interpolate_spin_data.return_value = mock_spin_df

l1b_de = xr.Dataset()
l1a_de = xr.Dataset(
{
"shcoarse": ("epoch", [15, 35]),
"de_count": ("epoch", [2, 3]),
"de_time": ("direct_event", [0000, 1000, 2000, 3000, 4000]),
"de_time": ("direct_event", [0, 1000, 2000, 3000, 4000]),
},
coords={
"epoch": [0, 1],
"direct_event": [
0,
1,
2,
3,
4,
],
"direct_event": [0, 1, 2, 3, 4],
},
)
avg_spin_durations = xr.DataArray([5, 10])
spin_start_times = xr.DataArray([10, 20, 30, 40, 50])
expected_event_met = np.array([10, 21.2207, 34.8828, 47.3242, 59.7656])

# shcoarse 15 -> spin_start 10, shcoarse 35 -> spin_start 30
# event_met = spin_start + de_time * DE_CLOCK_TICK_S
expected_event_met = np.array(
[
10 + 0 * DE_CLOCK_TICK_S, # 10.0
10 + 1000 * DE_CLOCK_TICK_S, # 14.096
30 + 2000 * DE_CLOCK_TICK_S, # 38.192
30 + 3000 * DE_CLOCK_TICK_S, # 42.288
30 + 4000 * DE_CLOCK_TICK_S, # 46.384
]
)

# Act
l1b_de = set_event_met(l1a_de, l1b_de, spin_start_times, avg_spin_durations)
l1b_de = set_event_met(l1a_de, l1b_de)

# Assert
np.testing.assert_allclose(
Expand All @@ -586,24 +586,25 @@ def test_set_event_met():
atol=1e-4,
)

def test_set_each_event_epoch():
l1b_de = xr.Dataset(
{
"event_met": ("epoch", [10, 20, 30, 40, 50]),
},
coords={
"epoch": [0, 1, 2, 3, 4],
},
)
epoch_expected = met_to_ttj2000ns(np.array([10, 20, 30, 40, 50]))

l1b_de = set_each_event_epoch(l1b_de)
def test_set_each_event_epoch():
l1b_de = xr.Dataset(
{
"event_met": ("epoch", [10, 20, 30, 40, 50]),
},
coords={
"epoch": [0, 1, 2, 3, 4],
},
)
epoch_expected = met_to_ttj2000ns(np.array([10, 20, 30, 40, 50]))

np.testing.assert_allclose(
l1b_de["epoch"].values,
epoch_expected,
atol=1e-4,
)
l1b_de = set_each_event_epoch(l1b_de)

np.testing.assert_allclose(
l1b_de["epoch"].values,
epoch_expected,
atol=1e-4,
)


def test_set_avg_spin_durations_per_event():
Expand Down Expand Up @@ -838,6 +839,8 @@ def test_set_direction(mock_lo_instrument_pointing, imap_ena_sim_metakernel):
)
@patch(
"imap_processing.lo.l1b.lo_l1b.cartesian_to_latitudinal",
# Longitudes: -180 -> 180, 0 -> 0, 90 -> 90, 180 -> 180
# After shift to 0-360: 180, 0, 90, 180
return_value=np.array([[0, -180, -2], [0, 0, 0], [0, 90, 1], [0, 180, 2]]),
)
def test_pointing_bins(mock_cartesian_to_latitudinal, mock_frame_transform):
Expand All @@ -859,7 +862,8 @@ def test_pointing_bins(mock_cartesian_to_latitudinal, mock_frame_transform):
)

expected_pointing_lats = np.array([0, 20, 30, 40])
expected_pointing_lons = np.array([0, 1800, 2700, 3600])
# Longitude bins are now in 0-360 range after the shift
expected_pointing_lons = np.array([1800, 0, 900, 1800])

# Act
l1b_de = set_pointing_bin(l1b_de)
Expand Down
Loading