Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
55 changes: 44 additions & 11 deletions imap_processing/hit/l1a/hit_l1a.py
Original file line number Diff line number Diff line change
Expand Up @@ -347,6 +347,10 @@ def find_complete_mod10_sets(mod_vals: np.ndarray) -> np.ndarray:
"""
Find start indices where mod values match [0,1,2,3,4,5,6,7,8,9] pattern.

To be a valid set, the previous set needs to be complete as well. Since the
first set will only be used for its livetime values, it is valid if it
matches the pattern.

Parameters
----------
mod_vals : np.ndarray
Expand All @@ -369,7 +373,19 @@ def find_complete_mod10_sets(mod_vals: np.ndarray) -> np.ndarray:
# Use sliding windows to find pattern matches
sw_view = np.lib.stride_tricks.sliding_window_view(mod_vals, window_size)
matches = np.all(sw_view == np.arange(window_size), axis=1)
return np.where(matches)[0]
starts = np.where(matches)[0]
# If no matches, return an empty array
if len(starts) == 0:
return np.array([], dtype=int)
# Find where the previous set is complete as well.
valid = np.isin(starts - 10, starts)
# If the first complete set is the first set,
# that set is considered complete since it's only used for livetime
# values.
valid[0] = starts[0] == 0
# Return the start inds of complete sets where the previous set is complete
# as well.
return starts[valid]


def subset_sectored_counts(
Expand Down Expand Up @@ -423,6 +439,14 @@ def subset_sectored_counts(
else:
start_indices = start_indices[start_indices >= bin_size]

# Sets can also be dropped when the previous set is incomplete (e.g.
# due to a dropped frame), so check again after filtering.
if start_indices.size == 0:
raise ValueError(
"No data to process - no complete sets of sectored counts with "
"a complete set of livetime values from the previous 10 minutes."
)

# Subset data for complete sets of sectored counts.
# Each set of sectored counts is 10 minutes long, so we take the indices
# starting from the start indices and extending to the bin size of 10.
Expand Down Expand Up @@ -519,21 +543,30 @@ def subset_livetime(dataset: xr.Dataset) -> xr.Dataset:
"Epoch values are empty. Cannot proceed with livetime subsetting."
)

# Get index positions of epoch[0] and epoch[-1] in epoch_livetime
start_idx = np.where(epoch_livetime_vals == epoch_vals[0])[0][0]
end_idx = np.where(epoch_livetime_vals == epoch_vals[-1])[0][0]
# Get the inds of the epoch_livetime where it matches the filtered epoch vals.
match_inds = np.where(np.isin(epoch_livetime_vals, epoch_vals))[0]

if start_idx < 10:
# Each epoch must match exactly one epoch_livetime value, otherwise
# livetime will be misaligned with the sectored counts.
if match_inds.size != epoch_vals.size:
raise ValueError(
"Start index for livetime is less than 10. This indicates that the "
"dataset is too small to shift livetime correctly."
f"Found {match_inds.size} epoch_livetime matches for "
f"{epoch_vals.size} epoch values. Livetime cannot be aligned with "
"sectored counts."
)

# Compute shifted indices by 10 minutes
start_trimmed = max(start_idx - 10, 0)
end_trimmed = max(end_idx - 10, 0)
# Shift 10 frames (10 minutes) earlier to get the livetime indices
livetime_inds = match_inds - 10

if np.any(livetime_inds < 0):
raise ValueError(
"Livetime indices are negative. This indicates that the "
"dataset is too small to shift livetime correctly."
)

return dataset.isel(epoch_livetime=slice(start_trimmed, end_trimmed + 1))
# Then return the dataset where there are only complete sets of
# livetime values corresponding to the sectored counts.
return dataset.isel(epoch_livetime=livetime_inds)


def filter_dataset_to_processing_day(
Expand Down
186 changes: 114 additions & 72 deletions imap_processing/tests/hit/test_hit_l1a.py
Original file line number Diff line number Diff line change
Expand Up @@ -202,17 +202,50 @@ def test_subset_livetime():
)
with pytest.raises(
ValueError,
match="Start index for livetime is less than 10. This indicates that the "
match="Livetime indices are negative. This indicates that the "
"dataset is too small to shift livetime correctly.",
):
subset_livetime(dataset)

# Test case 4: Non-contiguous epochs (sets dropped in the middle). Each
# epoch should map to the livetime value 10 frames earlier.
epoch_livetime = np.arange(60)
epoch = np.concatenate([epoch_livetime[10:20], epoch_livetime[40:50]])
dataset = xr.Dataset(
{
"livetime_counter": ("epoch_livetime", np.arange(60) * 2),
},
coords={"epoch": epoch, "epoch_livetime": epoch_livetime},
)
trimmed = subset_livetime(dataset)
expected_epoch_livetime = epoch - 10
np.testing.assert_array_equal(
trimmed["epoch_livetime"].values, expected_epoch_livetime
)
np.testing.assert_array_equal(
trimmed["livetime_counter"].values, expected_epoch_livetime * 2
)

# Test case 5: An epoch value is missing from epoch_livetime, so livetime
# can't be aligned with the sectored counts.
dataset = xr.Dataset(
{
"livetime_counter": ("epoch_livetime", np.arange(60)),
},
coords={
"epoch": np.array([10, 11, 100]),
"epoch_livetime": np.arange(60),
},
)
with pytest.raises(ValueError, match="Livetime cannot be aligned"):
subset_livetime(dataset)


def test_find_complete_mod10_sets():
"""Test the find_complete_mod10_sets function."""

complete_set = np.arange(10)
# Test case 1: Valid pattern exists throughout
mod_vals = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
mod_vals = np.tile(complete_set, 2)
expected_indices = np.array([0, 10])
result = find_complete_mod10_sets(mod_vals)
assert np.array_equal(result, expected_indices), (
Expand All @@ -228,8 +261,23 @@ def test_find_complete_mod10_sets():
)

# Test case 3: Pattern in the middle
mod_vals = np.array([5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4])
expected_indices = np.array([5])
mod_vals = np.concat(
[
np.array([5, 6, 7, 8, 9]),
complete_set,
complete_set,
np.array(
[
0,
1,
2,
3,
4,
]
),
]
)
expected_indices = np.array([15])
result = find_complete_mod10_sets(mod_vals)
assert np.array_equal(result, expected_indices), (
f"Expected {expected_indices}, got {result}"
Expand All @@ -243,6 +291,35 @@ def test_find_complete_mod10_sets():
f"Expected {expected_indices}, got {result}"
)

# Test case 5.
# incomplete, complete, complete, incomplete, complete
# A valid set has to have complete set before so only the set
# starting at 19 is valid.
mod_vals = np.concat(
[
np.array([0, 1, 2, 3, 4, 5, 6, 8, 9]),
complete_set,
complete_set,
np.array([0, 1, 2, 3, 4, 6, 7, 8, 9]),
complete_set,
]
)
expected_indices = np.array([19])
result = find_complete_mod10_sets(mod_vals)
assert np.array_equal(result, expected_indices), (
f"Expected {expected_indices}, got {result}"
)

# Test case 6: Data starts with 10 frames that are not a complete set
# (leftover frame + set missing a frame). The first complete set at index
# 10 is invalid because the previous 10 frames aren't a complete set.
mod_vals = np.array([9, 0, 1, 2, 3, 4, 5, 6, 8, 9] + list(range(10)) * 2)
expected_indices = np.array([20])
result = find_complete_mod10_sets(mod_vals)
assert np.array_equal(result, expected_indices), (
f"Expected {expected_indices}, got {result}"
)


def test_calculate_uncertainties():
"""Test the calculate_uncertainties function.
Expand Down Expand Up @@ -415,77 +492,42 @@ def validate_subset(l1a_counts_dataset):
l1a_counts_dataset = create_l1a_counts_dataset(np.arange(105, 135))
validate_subset(l1a_counts_dataset)

# Test with partial data in the middle of the dataset
# Test with a gap in the middle of the dataset. The set after the gap
# (130-139) is complete but the 10 frames before it are not a complete
# set, so there is no valid livetime and no data to process.
l1a_counts_dataset = create_l1a_counts_dataset(
[
100,
101,
102,
103,
104,
105,
106,
107,
108,
109,
110,
111,
112,
113,
114,
120,
121,
122,
123,
124,
130,
131,
132,
133,
134,
135,
136,
137,
138,
139,
]
list(range(100, 115)) + list(range(120, 125)) + list(range(130, 140))
)
validate_subset(l1a_counts_dataset)
with pytest.raises(
ValueError,
match="No data to process - no complete sets of sectored counts with "
"a complete set of livetime values from the previous 10 minutes.",
):
subset_sectored_counts(l1a_counts_dataset, packet_date="20100106")

# Test with partial data at the start, middle, and end of the dataset
# Test with partial data at the start, middle, and end of the dataset.
# The complete set (130-139) follows partial data, so it is dropped.
# No complete set is preceded by another complete set, so no valid start
# indices are found.
l1a_counts_dataset = create_l1a_counts_dataset(
[
105,
106,
107,
108,
109,
110,
111,
112,
113,
114,
115,
116,
117,
118,
119,
120,
121,
122,
130,
131,
132,
133,
134,
135,
136,
137,
138,
139,
140,
141,
]
list(range(105, 123)) + list(range(130, 142))
)
with pytest.raises(
ValueError,
match="No data to process - valid start indices not found for "
"complete sectored counts.",
):
subset_sectored_counts(l1a_counts_dataset, packet_date="20100106")

# Test with a dropped frame in the middle of the dataset. The set with the
# dropped frame (110-119, missing 119) and the set after it (120-129) are
# dropped. Only the following set (130-139) is kept.
l1a_counts_dataset = create_l1a_counts_dataset(
list(range(100, 119)) + list(range(120, 140))
)
subset_dataset = subset_sectored_counts(l1a_counts_dataset, packet_date="20100106")
np.testing.assert_array_equal(
subset_dataset["hdr_minute_cnt"].values, np.arange(130, 140)
)
validate_subset(l1a_counts_dataset)

Expand Down
Loading