diff --git a/imap_processing/ena_maps/utils/corrections.py b/imap_processing/ena_maps/utils/corrections.py index 6346b8198..71afb5536 100644 --- a/imap_processing/ena_maps/utils/corrections.py +++ b/imap_processing/ena_maps/utils/corrections.py @@ -468,11 +468,23 @@ def add_spacecraft_velocity_to_pset( f"add_spacecraft_velocity_to_pset does not support PSETs with " f"Logical_source: {pset.attrs['Logical_source']}" ) - et = ttj2000ns_to_et(pset["epoch"].values[0] + pointing_duration_ns / 2) - # Get spacecraft state in HAE frame - sc_state = geometry.imap_state(et, ref_frame=geometry.SpiceFrame.IMAP_HAE) - sc_velocity_vector = sc_state[3:6] + # Handle case where pointing duration is zero or negative to avoid invalid + # ephemeris time (this is used, for example, for empty psets due to + # goodtimes filtering) + if pointing_duration_ns <= 0: + logger.warning( + "Pointing duration is zero or negative. " + "Setting spacecraft velocity to zero." + ) + sc_velocity_vector = np.zeros(3) # Zero velocity vector + else: + # Compute ephemeris time (J2000 seconds) of PSET midpoint + et = ttj2000ns_to_et(pset["epoch"].values[0] + pointing_duration_ns / 2) + + # Get spacecraft state in HAE frame + sc_state = geometry.imap_state(et, ref_frame=geometry.SpiceFrame.IMAP_HAE) + sc_velocity_vector = sc_state[3:6] # Store spacecraft velocity as DataArray pset["sc_velocity"] = xr.DataArray( diff --git a/imap_processing/lo/l1c/lo_l1c.py b/imap_processing/lo/l1c/lo_l1c.py index 435495715..1204fe51f 100644 --- a/imap_processing/lo/l1c/lo_l1c.py +++ b/imap_processing/lo/l1c/lo_l1c.py @@ -87,11 +87,22 @@ def lo_l1c(sci_dependencies: dict, anc_dependencies: list) -> list[xr.Dataset]: logical_source = "imap_lo_l1c_pset" l1b_de = sci_dependencies["imap_lo_l1b_de"] l1b_goodtimes_only = filter_goodtimes(l1b_de, anc_dependencies) - # TODO: Need to handle case where no good times are found - # Set the pointing start and end times based on the first epoch - pointing_start_met, pointing_end_met = get_pointing_times( - ttj2000ns_to_met(l1b_goodtimes_only["epoch"][0].item()) - ) + + # Handle case where no good times are found after filtering, + # which would lead to an empty dataset + if len(l1b_goodtimes_only["epoch"]) == 0: + logging.warning( + "No good times found in L1B DE dataset after filtering. " + "Creating empty PSET dataset with zero counts and exposure time." + ) + # Set dummy pointing start and end METs + pointing_start_met = 0.0 + pointing_end_met = 0.0 + else: + # Set the pointing start and end times based on the first epoch + pointing_start_met, pointing_end_met = get_pointing_times( + float(ttj2000ns_to_met(l1b_goodtimes_only["epoch"][0].item()).item()) + ) pset = xr.Dataset( coords={"epoch": np.array([met_to_ttj2000ns(pointing_start_met)])}, @@ -132,13 +143,19 @@ def lo_l1c(sci_dependencies: dict, anc_dependencies: list) -> list[xr.Dataset]: ) # Get the start and end spin numbers based on the pointing start and end MET + if 0 == pointing_start_met: + start_spin_number = 0 + end_spin_number = 0 + else: + start_spin_number = get_spin_number(pset["pointing_start_met"].item()) + end_spin_number = get_spin_number(pset["pointing_end_met"].item()) pset["start_spin_number"] = xr.DataArray( - [get_spin_number(pset["pointing_start_met"].item())], + [start_spin_number], dims="epoch", attrs=attr_mgr.get_variable_attributes("start_spin_number"), ) pset["end_spin_number"] = xr.DataArray( - [get_spin_number(pset["pointing_end_met"].item())], + [end_spin_number], dims="epoch", attrs=attr_mgr.get_variable_attributes("end_spin_number"), ) @@ -228,7 +245,7 @@ def filter_goodtimes(l1b_de: xr.Dataset, anc_dependencies: list) -> xr.Dataset: l1b_de : xarray.Dataset Filtered L1B Direct Event dataset. """ - # the goodtimes are currently the only ancillary file needed for L1C processing + # The goodtimes are one of several ancillary files needed for L1C processing goodtimes_table_df = lo_ancillary.read_ancillary_file( next(str(s) for s in anc_dependencies if "good-times" in str(s)) ) @@ -777,6 +794,17 @@ def calculate_exposure_times( # Calculate total pointing duration in seconds total_pointing_duration = pointing_end_met - pointing_start_met + if total_pointing_duration <= 0: + logging.warning( + "Pointing duration is zero or negative. Exposure times will be zero." + ) + # Return zero exposure times with correct shape and dimensions + zero_exposure = np.zeros(PSET_SHAPE, dtype=np.float32) + return xr.DataArray( + data=zero_exposure, + dims=PSET_DIMS, + ) + # Get representative spins from the pointing period representative_spins = get_representative_spin_times( pointing_start_met, pointing_end_met, n_representative_spins @@ -1063,8 +1091,8 @@ def set_background_rates( # TODO: This assumes that the backgrounds will never change mid-pointing. # Is that a safe assumption? pointing_bg_df = background_df[ - (background_df["GoodTime_start"] <= pointing_start_met) - & (background_df["GoodTime_end"] >= pointing_end_met) + (background_df["GoodTime_start"] < pointing_end_met) + & (background_df["GoodTime_end"] > pointing_start_met) ] # convert the bin start and end resolution from 6 degrees to 0.1 degrees @@ -1138,6 +1166,24 @@ def set_pointing_directions( hae_latitude : xr.DataArray The HAE latitude for each spin and off angle bin. """ + # Handle case where epoch is empty + if epoch == met_to_ttj2000ns(0): + return xr.DataArray( + data=np.zeros( + (1, len(SPIN_ANGLE_BIN_CENTERS), len(OFF_ANGLE_BIN_CENTERS)), + dtype=np.float64, + ), + dims=["epoch", "spin_angle", "off_angle"], + attrs=attr_mgr.get_variable_attributes("hae_longitude"), + ), xr.DataArray( + data=np.zeros( + (1, len(SPIN_ANGLE_BIN_CENTERS), len(OFF_ANGLE_BIN_CENTERS)), + dtype=np.float64, + ), + dims=["epoch", "spin_angle", "off_angle"], + attrs=attr_mgr.get_variable_attributes("hae_latitude"), + ) + et = ttj2000ns_to_et(epoch) # create a meshgrid of spin and off angles using the bin centers spin, off = np.meshgrid( diff --git a/imap_processing/tests/lo/test_lo_l1c.py b/imap_processing/tests/lo/test_lo_l1c.py index 71fdb971d..297eea389 100644 --- a/imap_processing/tests/lo/test_lo_l1c.py +++ b/imap_processing/tests/lo/test_lo_l1c.py @@ -295,6 +295,38 @@ def test_filter_goodtimes(l1b_de, anc_dependencies): xr.testing.assert_equal(l1b_goodtimes_only, l1b_goodtimes_onl_expected) +def test_lo_l1c_no_goodtimes( + l1b_de_spin, + anc_dependencies, + use_fake_repoint_data_for_time, + use_fake_spin_data_for_time, + repoint_met, +): + # Arrange + data = {"imap_lo_l1b_de": l1b_de_spin} + use_fake_spin_data_for_time(511000000) + use_fake_repoint_data_for_time(np.arange(511000000, 511000000 + 86400 * 5, 86400)) + expected_logical_source = "imap_lo_l1c_pset" + + # Act + output_dataset = lo_l1c(data, anc_dependencies)[0] + + # Assert + assert expected_logical_source == output_dataset.attrs["Logical_source"] + # Verify that pivot_angle is passed through from l1b_de + assert "pivot_angle" in output_dataset + assert output_dataset["pivot_angle"].values[0] == 45.0 + expected_counts = np.zeros((1, 7, 3600, 40)) + np.testing.assert_array_equal(output_dataset["h_counts"], expected_counts) + np.testing.assert_array_equal(output_dataset["o_counts"], expected_counts) + np.testing.assert_array_equal(output_dataset["doubles_counts"], expected_counts) + np.testing.assert_array_equal(output_dataset["triples_counts"], expected_counts) + np.testing.assert_array_equal(output_dataset["h_background_rates"], expected_counts) + np.testing.assert_array_equal(output_dataset["o_background_rates"], expected_counts) + expected = np.zeros((1, 3600, 40)) + np.testing.assert_array_equal(output_dataset["hae_latitude"], expected) + + def test_create_pset_counts(l1b_de): # Arrange expected_counts = np.zeros((1, 7, 3600, 40))