Lo L1B - Histogram Rate with Spin Cycle and ESA Mode #2482

sdhoyt · 2025-12-04T21:54:12Z

Change Summary

Overview

Added spin cycle and ESA mode fields to the L1B Histogram Rate product. This PR includes a re-write of the spin cycle field that should be at least close to working for DEs as well, but it has not been integrated with the DE code. The DE code is still using an older function and this should be updated in a later PR.
Code was added to drop any histograms tied to ASCs that did not complete 28 spins
The get_spin_number function was also updated to accept a list of times and to return a list of spin numbers

greglucas · 2025-12-05T15:13:49Z

imap_processing/lo/l1b/lo_l1b.py

+    # Find the closest start_acq for each direct event
+    # computes the index of the closest spin_met_per_asc for each science_met_per_asc
+    # so the resulting array will be of length len(science_met_per_asc), one index per
+    # ASC, but the value of each index will be the index of the closest spin data.
+    closest_start_acq_indices = np.abs(
+        science_met_per_asc[:, None] - spin_met_per_asc
+    ).argmin(axis=1)


I think we want this to be one-sided. If an event has a time stamp of 1000 and the next spin ASC has a timestamp of 1002, we actually want to place it in the previous ASC because that ASC hasn't started yet.
https://numpy.org/doc/stable/reference/generated/numpy.searchsorted.html

I had issues with the way searchsorted inserts instead of matches and certain edge cases that come with that. I think my change addressed the main issue of the spin time needing to be before the science time.

Let me know if you think this addresses your comment (or if you think the searchsorted should be a lot simpler than I'm making it sound) and also your thoughts on handling the case were the times match exactly. I was having decimal issues for exact matches and added a tolerance to check against. You might have a better way to handle that though.

I'm not sure why we need to worry about floating point issues, I think that'd all just fall out. That seems to add quite some complexity / nuance IMO.

My thinking for how this would work is with the following example, where you'd want to use side="right" I think to get what you're after if I'm understanding correctly. You might also have to subtract by 1 to get the index you want since it is where to insert it into the array.

>>> de_times = [0, 0.5, 0.9, 1, 1.0001, 1.5] >>> spin_times = [0, 1, 2, 3] >>> np.searchsorted(spin_times, de_times) array([0, 1, 1, 1, 2, 2]) >>> np.searchsorted(spin_times, de_times, side="right") array([1, 1, 1, 2, 2, 2])

I did try that, but was having a more difficult time getting around the floating point issue I mentioned. Maybe it's something wrong with how I'm setting up my test. This is what I'm getting:

spin_met_per_asc Out[2]: array([100.]) science_met_per_asc Out[3]: array([100.]) np.searchsorted(spin_met_per_asc, science_met_per_asc, side="right") - 1 Out[4]: array([-1])

the 100 for science is getting interpreted as before the spin time because of the floating point. My code handled the case where a science met is before the first spin met by dropping times when the searchsorted is < 0, but in this case, they should be equal so I wouldn't want it to get dropped.

Search sorted needs to know where to put it. In a 1element array do you want it before or after? This does bring up the question of whether we want to put a right edge on the furthest out spins if say we got DEs way after ASCs because those ASCs weren't available or valid...

greglucas · 2025-12-05T15:17:36Z

imap_processing/lo/l1b/lo_l1b.py

+    valid_asc_mask = []
+    for spin_idx in closest_start_acq_indices:
+        spin_count = spin_data["num_completed"].isel(epoch=spin_idx)
+        valid_asc_mask.append(spin_count >= 28)


I think you can avoid the loop here with numpy indexing.

Suggested change

valid_asc_mask = []

for spin_idx in closest_start_acq_indices:

spin_count = spin_data["num_completed"].isel(epoch=spin_idx)

valid_asc_mask.append(spin_count >= 28)

# It isn't a valid cycle if we didn't complete all 28 spins in an ASC

valid_cycles = spin_data["num_completed"] == 28

valid_asc_mask = valid_cycles.values[closest_start_acq_indices]

greglucas · 2025-12-05T15:19:27Z

imap_processing/lo/l1b/lo_l1b.py

+        closest_start_acq_indices = closest_start_acq_indices[valid_asc_mask]
+
+    closest_start_acq_per_asc = acq_start.isel(
+        epoch=xr.DataArray(closest_start_acq_indices)


Suggested change

epoch=xr.DataArray(closest_start_acq_indices)

epoch=closest_start_acq_indices

greglucas · 2025-12-05T15:24:47Z

imap_processing/lo/l1b/lo_l1b.py

+    # Get the spin cycle number from the spin data for each direct event
+    spin_start_num_per_asc = np.atleast_1d(get_spin_number(closest_start_acq_per_asc))


Could you put in a comment about "spin_cycle" being the average spin number for that ESA step during that ASC? It was not clear to me after our discussion on Monday that the cycle was an average value because we repeat the ESA step pattern twice with 0-13 containing the first sequence and 14-27 containing the second sequence so we are averaging over those two sequences.

…asc functions

greglucas · 2025-12-11T03:00:58Z

imap_processing/lo/l1b/lo_l1b.py

+    # computes the index of the closest spin_met_per_asc for each science_met_per_asc
+    # so the resulting array will be of length len(science_met_per_asc), one index per
+    # ASC, but the value of each index will be the index of the closest spin data.
+    # Find the index of the spin data that starts just before or at each science epoch


I thought we weren't sure which packet got timestamped first, so we wanted to do this as the absolute value closest like you originally had?

greglucas · 2025-12-11T03:02:15Z

imap_processing/lo/l1b/lo_l1b.py

+    closest_start_acq_indices = time_diff_masked.argmin(axis=1)
+
+    return closest_start_acq_indices, time_diff


Do you need to return time_diff here, or can you set all time_diff < 0 to an index of -1 or something that indicates it is a bad index for use later instead?

greglucas · 2025-12-11T03:13:35Z

imap_processing/lo/l1b/lo_l1b.py

+        temp_mask = np.copy(valid_mask)
+        temp_mask[valid_mask] &= ~insufficient_spins
+        valid_mask = temp_mask


Rather than this subsetting of valid indices, what if each of these if functions was its own function / handler so you are operating on the full arrays each time? Then at the end do something like combine all the masks.

valid_mask = valid_cycles & valid_times & valid_prior_spins

greglucas · 2025-12-11T03:17:33Z

imap_processing/lo/l1b/lo_l1b.py

+    spin_met_per_asc = spin_data["shcoarse"].values.astype(np.float64)
+    science_met_per_asc = ttj2000ns_to_met(l1a_science["epoch"]).astype(np.float64)
+
+    closest_start_acq_indices, time_diff = match_science_to_spin_asc(


I don't think this is actually the start_acq indices, because you are matching on the SHCOARSE values which were timestamped after the full acquisition of each of the packets right?

Suggested change

closest_start_acq_indices, time_diff = match_science_to_spin_asc(

science_to_spin_indices, time_diff = match_science_to_spin_asc(

yeah, you're right, I'll update that name

greglucas · 2025-12-11T03:31:47Z

imap_processing/tests/lo/test_lo_l1b.py

+    # science_met[2]=300 -> searchsorted returns 3, -1 = 2 (spin_met[2]=250)
+    assert l1b_hist["spin_cycle"][0, 0] == 7  # 0 + 7 + 0
+    assert l1b_hist["spin_cycle"][1, 0] == 35  # 28 + 7 + 0
+    assert l1b_hist["spin_cycle"][2, 0] == 63  # 56 + 7 + 0


Can you test a different energy step here too?

Suggested change

assert l1b_hist["spin_cycle"][2, 0] == 63 # 56 + 7 + 0

assert l1b_hist["spin_cycle"][2, 2] == 67 # 56 + 7 + 2*2

greglucas · 2025-12-11T03:31:59Z

imap_processing/tests/lo/test_lo_l1b.py

+    # science_met[0]=100 -> searchsorted returns 1, -1 = 0 (spin_met[0]=50)
+    # science_met[1]=200 -> searchsorted returns 2, -1 = 1 (spin_met[1]=150)
+    # science_met[2]=300 -> searchsorted returns 3, -1 = 2 (spin_met[2]=250)


No more searchsorted anymore.

sdhoyt added 2 commits December 3, 2025 16:40

added spin cycle code and setting esa_mode

4255475

updated docstrings

6de406a

sdhoyt added this to the December 2025 milestone Dec 4, 2025

sdhoyt requested review from greglucas, lacoak21 and subagonsouth December 4, 2025 21:54

sdhoyt self-assigned this Dec 4, 2025

sdhoyt added Ins: Lo Related to the IMAP-Lo instrument Level: L1 Level 1 processing labels Dec 4, 2025

sdhoyt added this to IMAP Dec 4, 2025

sdhoyt requested review from KenWFairchild, frahmanif, hi1011 and nschwadron December 4, 2025 21:54

greglucas reviewed Dec 5, 2025

View reviewed changes

sdhoyt added 2 commits December 5, 2025 15:44

fixed spin, histogram ASC matching

0f8c981

refactored new spin cycle function, created asc match and find valid …

64494a9

…asc functions

greglucas reviewed Dec 11, 2025

View reviewed changes

back to using abs value

cd58a7a

	epoch=xr.DataArray(closest_start_acq_indices)
	epoch=closest_start_acq_indices

		# Get the spin cycle number from the spin data for each direct event
		spin_start_num_per_asc = np.atleast_1d(get_spin_number(closest_start_acq_per_asc))

		closest_start_acq_indices = time_diff_masked.argmin(axis=1)

		return closest_start_acq_indices, time_diff

	closest_start_acq_indices, time_diff = match_science_to_spin_asc(
	science_to_spin_indices, time_diff = match_science_to_spin_asc(

	assert l1b_hist["spin_cycle"][2, 0] == 63 # 56 + 7 + 0
	assert l1b_hist["spin_cycle"][2, 2] == 67 # 56 + 7 + 2*2

Lo L1B - Histogram Rate with Spin Cycle and ESA Mode #2482

Are you sure you want to change the base?

Lo L1B - Histogram Rate with Spin Cycle and ESA Mode #2482

Uh oh!

Conversation

sdhoyt commented Dec 4, 2025

Change Summary

Overview

Uh oh!

Choose a reason for hiding this comment

Uh oh!

sdhoyt Dec 5, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

sdhoyt Dec 5, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

sdhoyt Dec 5, 2025 •

edited

Loading

sdhoyt Dec 5, 2025 •

edited

Loading