From 00b70ad7448ac2f7024881669d1ef28203ba691a Mon Sep 17 00:00:00 2001 From: Tenzin Choedon <36522642+tech3371@users.noreply.github.com> Date: Wed, 31 Jul 2024 12:35:29 -0600 Subject: [PATCH] SWAPI l2 Implementation - counts to rate (#721) --- .../config/imap_swapi_global_cdf_attrs.yaml | 11 ++- .../cdf/config/imap_swapi_variable_attrs.yaml | 44 ++++++++- imap_processing/cli.py | 10 ++ imap_processing/swapi/l1/swapi_l1.py | 13 ++- imap_processing/swapi/l2/__init__.py | 0 imap_processing/swapi/l2/swapi_l2.py | 97 +++++++++++++++++++ imap_processing/tests/swapi/conftest.py | 16 +++ .../tests/swapi/test_swapi_decom.py | 24 ++--- imap_processing/tests/swapi/test_swapi_l1.py | 19 ++-- imap_processing/tests/swapi/test_swapi_l2.py | 21 ++++ 10 files changed, 221 insertions(+), 34 deletions(-) create mode 100644 imap_processing/swapi/l2/__init__.py create mode 100644 imap_processing/swapi/l2/swapi_l2.py create mode 100644 imap_processing/tests/swapi/conftest.py create mode 100644 imap_processing/tests/swapi/test_swapi_l2.py diff --git a/imap_processing/cdf/config/imap_swapi_global_cdf_attrs.yaml b/imap_processing/cdf/config/imap_swapi_global_cdf_attrs.yaml index c43fc7e78..ed8ce687e 100644 --- a/imap_processing/cdf/config/imap_swapi_global_cdf_attrs.yaml +++ b/imap_processing/cdf/config/imap_swapi_global_cdf_attrs.yaml @@ -13,14 +13,21 @@ instrument_base: &instrument_base imap_swapi_l1_sci: <<: *instrument_base # NOTE: Right now, this Data_level is required to produce valid CDF - Data_level: 1 + Data_level: "1" Data_type: L1_SCI>Level-1 Science data in 1 minute resolution Logical_source: imap_swapi_l1_sci-1min Logical_source_description: SWAPI Instrument Level-1 Science Data in 1 minute resolution imap_swapi_l1_hk: <<: *instrument_base - Data_level: 1 + Data_level: "1" Data_type: L1_HK>Level-1B Housekeeping data Logical_source: imap_swapi_l1_hk Logical_source_description: SWAPI Instrument Level-1 Housekeeping Data + +imap_swapi_l2_sci: + <<: *instrument_base + Data_level: "2" + Data_type: L2_SCI>Level-2 Science data in 1 minute resolution + Logical_source: imap_swapi_l2_sci-1min + Logical_source_description: SWAPI Instrument Level-1 Science Data in 1 minute resolution diff --git a/imap_processing/cdf/config/imap_swapi_variable_attrs.yaml b/imap_processing/cdf/config/imap_swapi_variable_attrs.yaml index 083050eb1..7746a5e86 100644 --- a/imap_processing/cdf/config/imap_swapi_variable_attrs.yaml +++ b/imap_processing/cdf/config/imap_swapi_variable_attrs.yaml @@ -24,7 +24,7 @@ counts_default: &counts_default DEPEND_1: energy DISPLAY_TYPE: spectrogram LABL_PTR_1: energy_label - FILLVAL: -9223372036854775808 + FILLVAL: 4294967295 FORMAT: I5 UNITS: counts VALIDMIN: 0 @@ -37,7 +37,7 @@ compression_default: &compression_default DEPEND_1: energy DISPLAY_TYPE: spectrogram LABL_PTR_1: energy_label - FILLVAL: -9223372036854775808 + FILLVAL: 4294967295 FORMAT: I1 UNITS: ' ' VALIDMIN: 0 @@ -58,47 +58,87 @@ uncertainty_default: &uncertainty_default VAR_TYPE: data SCALETYP: linear +rate_default: &rate_default + DEPEND_0: epoch + DEPEND_1: energy + DISPLAY_TYPE: spectrogram + LABL_PTR_1: energy_label + FILLVAL: -1.0000000E+31 + FORMAT: E19.5 + UNITS: counts + VALIDMIN: 0.0 + VALIDMAX: 1.7976931348623157e+308 # TODO: find actual value + VAR_TYPE: data + SCALETYP: linear + pcem_counts: <<: *counts_default CATDESC: Primary Channel Electron Multiplier (CEM) Counts FIELDNAM: Primary CEM Counts + LABLAXIS: pcem_cnts scem_counts: <<: *counts_default CATDESC: Secondary Channel Electron Multiplier (CEM) Counts FIELDNAM: Secondary CEM Counts + LABLAXIS: scem_cnts coin_counts: <<: *counts_default CATDESC: Coincidence Counts FIELDNAM: Coincidence CEM Counts + LABLAXIS: coin_cnts pcem_flags: <<: *compression_default CATDESC: Primary Channel Electron Multiplier (CEM) compression flag FIELDNAM: Primary CEM Flag + LABLAXIS: pcem_flag scem_flags: <<: *compression_default CATDESC: Secondary Channel Electron Multiplier (CEM) compression flag FIELDNAM: Secondary CEM Flag + LABLAXIS: scem_flag coin_flags: <<: *compression_default CATDESC: Coincidence compression flag FIELDNAM: Coincidence Flag + LABLAXIS: coin_flag pcem_uncertainty: <<: *uncertainty_default CATDESC: Primary Channel Electron Multiplier (CEM) uncertainty FIELDNAM: Primary CEM Uncertainty + LABLAXIS: pcem_unc scem_uncertainty: <<: *uncertainty_default CATDESC: Secondary Channel Electron Multiplier (CEM) uncertainty FIELDNAM: Secondary CEM Uncertainty + LABLAXIS: scem_unc coin_uncertainty: <<: *uncertainty_default CATDESC: Coincidence uncertainty FIELDNAM: Coincidence Uncertainty + LABLAXIS: coin_unc + +pcem_rate: + <<: *rate_default + CATDESC: Primary Channel Electron Multiplier (CEM) Rates + FIELDNAM: Primary CEM Rates + LABLAXIS: pcem_rate + +scem_rate: + <<: *rate_default + CATDESC: Secondary Channel Electron Multiplier (CEM) Rates + FIELDNnam: Secondary CEM Rates + LABLAXIS: scem_rate + +coin_rate: + <<: *rate_default + CATDESC: Coincidence Rates + FIELDNAM: Coincidence Rates + LABLAXIS: coin_rate diff --git a/imap_processing/cli.py b/imap_processing/cli.py index 30024fba4..ad7cb7763 100644 --- a/imap_processing/cli.py +++ b/imap_processing/cli.py @@ -51,6 +51,7 @@ from imap_processing.mag.l1b.mag_l1b import mag_l1b from imap_processing.mag.l1c.mag_l1c import mag_l1c from imap_processing.swapi.l1.swapi_l1 import swapi_l1 +from imap_processing.swapi.l2.swapi_l2 import swapi_l2 from imap_processing.swe.l1a.swe_l1a import swe_l1a from imap_processing.swe.l1b.swe_l1b import swe_l1b from imap_processing.ultra.l1a import ultra_l1a @@ -725,6 +726,15 @@ def do_processing(self, dependencies: list) -> list[xr.Dataset]: ) # process data datasets = [swapi_l1(dependencies[0], self.version)] + elif self.data_level == "l2": + if len(dependencies) > 1: + raise ValueError( + f"Unexpected dependencies found for SWAPI L2:" + f"{dependencies}. Expected only one dependency." + ) + # process data + l1_dataset = load_cdf(dependencies[0]) + datasets = [swapi_l2(l1_dataset, self.version)] return datasets diff --git a/imap_processing/swapi/l1/swapi_l1.py b/imap_processing/swapi/l1/swapi_l1.py index c43d91c26..cd3c29f62 100644 --- a/imap_processing/swapi/l1/swapi_l1.py +++ b/imap_processing/swapi/l1/swapi_l1.py @@ -488,17 +488,16 @@ def process_swapi_science(sci_dataset: xr.Dataset, data_version: str) -> xr.Data ) # Add other global attributes + # TODO: add others like below once add_global_attribute is fixed cdf_manager.add_global_attribute("Data_version", data_version) - cdf_manager.add_global_attribute( - "sweep_table", f"{sci_dataset['sweep_table'].data[0]}" - ) - cdf_manager.add_global_attribute( - "plan_id", f"{sci_dataset['plan_id_science'].data[0]}" - ) + l1_global_attrs = cdf_manager.get_global_attributes("imap_swapi_l1_sci") + l1_global_attrs["Sweep_table"] = f"{sci_dataset['sweep_table'].data[0]}" + l1_global_attrs["Plan_id"] = f"{sci_dataset['plan_id_science'].data[0]}" + l1_global_attrs["Apid"] = f"{sci_dataset['pkt_apid'].data[0]}" dataset = xr.Dataset( coords={"epoch": epoch_time, "energy": energy, "energy_label": energy_label}, - attrs=cdf_manager.get_global_attributes("imap_swapi_l1_sci"), + attrs=l1_global_attrs, ) dataset["swp_pcem_counts"] = xr.DataArray( diff --git a/imap_processing/swapi/l2/__init__.py b/imap_processing/swapi/l2/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/imap_processing/swapi/l2/swapi_l2.py b/imap_processing/swapi/l2/swapi_l2.py new file mode 100644 index 000000000..f59826659 --- /dev/null +++ b/imap_processing/swapi/l2/swapi_l2.py @@ -0,0 +1,97 @@ +"""SWAPI L2 processing module.""" + +import logging + +import xarray as xr + +from imap_processing.cdf.imap_cdf_manager import ImapCdfAttributes + +logger = logging.getLogger(__name__) + + +TIME_PER_BIN = 0.167 # seconds + + +def swapi_l2(l1_dataset: xr.Dataset, data_version: str) -> xr.Dataset: + """ + Produce science data to L2. + + To process science data to L2, we need to: + - convert counts to rates. This is done by dividing the counts by the + TIME_PER_BIN time. TIME_PER_BIN is the exposure time per energy bin which is + obtained by dividing the time for one complete sweep + (12 s, coarse + fine sweep) by the total energy steps (72), + i.e., TIME_PER_BIN = 12/72 = 0.167 s. This will be constant. + + - update uncertainty. Calculate new uncertainty value using + SWP_PCEM_ERR data from level one and divide by TIME_PER_BIN. Eg. + SWP_PCEM_UNC = SWP_PCEM_ERR / TIME_PER_BIN + Do the same for SCEM and COIN data. + + Parameters + ---------- + l1_dataset : xarray.Dataset + The L1 data input. + data_version : str + Version of the data product being created. + + Returns + ------- + data : xarray.Dataset + Processed data to L2. + """ + # Load the CDF attributes + cdf_manager = ImapCdfAttributes() + cdf_manager.add_instrument_global_attrs("swapi") + cdf_manager.load_variable_attributes("imap_swapi_variable_attrs.yaml") + + # Copy over only certain variables from L1 to L2 dataset + l1_data_keys = [ + "epoch", + "energy", + "energy_label", + "swp_pcem_flags", + "swp_scem_flags", + "swp_coin_flags", + ] + l2_dataset = l1_dataset[l1_data_keys] + + # Update L2 specific attributes + l2_dataset.attrs["Data_version"] = data_version + l2_global_attrs = cdf_manager.get_global_attributes("imap_swapi_l2_sci") + l2_dataset.attrs["Data_level"] = l2_global_attrs["Data_level"] + l2_dataset.attrs["Data_type"] = l2_global_attrs["Data_type"] + l2_dataset.attrs["Logical_source"] = l2_global_attrs["Logical_source"] + l2_dataset.attrs["Logical_source_description"] = l2_global_attrs[ + "Logical_source_description" + ] + + # convert counts to rate + l2_dataset["swp_pcem_rate"] = l1_dataset["swp_pcem_counts"] / TIME_PER_BIN + l2_dataset["swp_scem_rate"] = l1_dataset["swp_scem_counts"] / TIME_PER_BIN + l2_dataset["swp_coin_rate"] = l1_dataset["swp_coin_counts"] / TIME_PER_BIN + # update attrs + l2_dataset["swp_pcem_rate"].attrs = cdf_manager.get_variable_attributes("pcem_rate") + l2_dataset["swp_scem_rate"].attrs = cdf_manager.get_variable_attributes("scem_rate") + l2_dataset["swp_coin_rate"].attrs = cdf_manager.get_variable_attributes("coin_rate") + + # update uncertainty + l2_dataset["swp_pcem_unc"] = l1_dataset["swp_pcem_err"] / TIME_PER_BIN + l2_dataset["swp_scem_unc"] = l1_dataset["swp_scem_err"] / TIME_PER_BIN + l2_dataset["swp_coin_unc"] = l1_dataset["swp_coin_err"] / TIME_PER_BIN + # update attrs + l2_dataset["swp_pcem_unc"].attrs = cdf_manager.get_variable_attributes( + "pcem_uncertainty" + ) + l2_dataset["swp_scem_unc"].attrs = cdf_manager.get_variable_attributes( + "scem_uncertainty" + ) + l2_dataset["swp_coin_unc"].attrs = cdf_manager.get_variable_attributes( + "coin_uncertainty" + ) + + # TODO: add thruster firing flag + # TODO: add other flags + logger.info("SWAPI L2 processing complete") + + return l2_dataset diff --git a/imap_processing/tests/swapi/conftest.py b/imap_processing/tests/swapi/conftest.py new file mode 100644 index 000000000..6b91ba97c --- /dev/null +++ b/imap_processing/tests/swapi/conftest.py @@ -0,0 +1,16 @@ +import pytest + + +@pytest.fixture(scope="session") +def swapi_test_data_path(imap_tests_path): + return imap_tests_path / "swapi/" + + +@pytest.fixture(scope="session") +def swapi_l0_test_data_path(swapi_test_data_path): + return swapi_test_data_path / "l0_data/" + + +@pytest.fixture(scope="session") +def swapi_l0_validation_data_path(swapi_test_data_path): + return swapi_test_data_path / "l0_validation_data/" diff --git a/imap_processing/tests/swapi/test_swapi_decom.py b/imap_processing/tests/swapi/test_swapi_decom.py index f6c5e406e..f3aff0273 100644 --- a/imap_processing/tests/swapi/test_swapi_decom.py +++ b/imap_processing/tests/swapi/test_swapi_decom.py @@ -10,16 +10,15 @@ @pytest.fixture(scope="session") -def decom_test_data(): +def decom_test_data(swapi_l0_test_data_path): """Read test data from file""" - test_folder_path = "tests/swapi/l0_data" - packet_files = list(imap_module_directory.glob(f"{test_folder_path}/*.pkts")) + test_file = "imap_swapi_l0_raw_20231012_v001.pkts" + packet_file = imap_module_directory / swapi_l0_test_data_path / test_file packet_definition = ( f"{imap_module_directory}/swapi/packet_definitions/swapi_packet_definition.xml" ) data_list = [] - for packet_file in packet_files: - data_list.extend(decom_packets(packet_file, packet_definition)) + data_list.extend(decom_packets(packet_file, packet_definition)) return data_list @@ -39,12 +38,11 @@ def test_number_of_packets(decom_test_data): assert len(aut_packets) == expected_aut_packets -def test_swapi_sci_data(decom_test_data): +def test_swapi_sci_data(decom_test_data, swapi_l0_validation_data_path): """This test and validate raw data of SWAPI raw science data.""" # read validation data - test_data_path = imap_module_directory / "tests/swapi/l0_validation_data" raw_validation_data = pd.read_csv( - test_data_path / "idle_export_eu.SWP_SCI_20231012_125245.csv", + swapi_l0_validation_data_path / "idle_export_eu.SWP_SCI_20231012_125245.csv", index_col="SHCOARSE", ) @@ -77,12 +75,11 @@ def test_swapi_sci_data(decom_test_data): ) -def test_swapi_hk_data(decom_test_data): +def test_swapi_hk_data(decom_test_data, swapi_l0_validation_data_path): """This test and validate raw data of SWAPI raw housekeeping data.""" # read validation data - test_data_path = imap_module_directory / "tests/swapi/l0_validation_data" raw_validation_data = pd.read_csv( - test_data_path / "idle_export_raw.SWP_HK_20231012_125245.csv", + swapi_l0_validation_data_path / "idle_export_raw.SWP_HK_20231012_125245.csv", index_col="SHCOARSE", ) @@ -113,12 +110,11 @@ def test_swapi_hk_data(decom_test_data): assert value.raw_value == validation_data[key] -def test_swapi_aut_data(decom_test_data): +def test_swapi_aut_data(decom_test_data, swapi_l0_validation_data_path): """This test and validate raw data of SWAPI raw autonomy data.""" # read validation data - test_data_path = imap_module_directory / "tests/swapi/l0_validation_data" raw_validation_data = pd.read_csv( - test_data_path / "idle_export_raw.SWP_AUT_20231012_125245.csv", + swapi_l0_validation_data_path / "idle_export_raw.SWP_AUT_20231012_125245.csv", index_col="SHCOARSE", ) diff --git a/imap_processing/tests/swapi/test_swapi_l1.py b/imap_processing/tests/swapi/test_swapi_l1.py index 11b7b41d3..5c221eef1 100644 --- a/imap_processing/tests/swapi/test_swapi_l1.py +++ b/imap_processing/tests/swapi/test_swapi_l1.py @@ -19,10 +19,10 @@ @pytest.fixture(scope="session") -def decom_test_data(): +def decom_test_data(swapi_l0_test_data_path): """Read test data from file""" - test_file = "tests/swapi/l0_data/imap_swapi_l0_raw_20231012_v001.pkts" - packet_files = imap_module_directory / test_file + test_file = "imap_swapi_l0_raw_20231012_v001.pkts" + packet_files = imap_module_directory / swapi_l0_test_data_path / test_file packet_definition = ( f"{imap_module_directory}/swapi/packet_definitions/swapi_packet_definition.xml" ) @@ -316,13 +316,14 @@ def test_process_swapi_science(decom_test_data): assert cdf_path.name == cdf_filename -def test_swapi_l1_cdf(): +def test_swapi_l1_cdf(swapi_l0_test_data_path): """Test housekeeping processing and CDF file creation""" - l0_data_path = ( - f"{imap_module_directory}/tests/swapi/l0_data/" - "imap_swapi_l0_raw_20231012_v001.pkts" - ) - processed_data = swapi_l1(l0_data_path, data_version="v001") + test_packet_file = swapi_l0_test_data_path / "imap_swapi_l0_raw_20231012_v001.pkts" + processed_data = swapi_l1(test_packet_file, data_version="v001") + + assert processed_data[0].attrs["Apid"] == f"{SWAPIAPID.SWP_SCI}" + assert processed_data[0].attrs["Plan_id"] == "0" + assert processed_data[0].attrs["Sweep_table"] == "0" # Test CDF File # sci cdf file diff --git a/imap_processing/tests/swapi/test_swapi_l2.py b/imap_processing/tests/swapi/test_swapi_l2.py new file mode 100644 index 000000000..850e3f964 --- /dev/null +++ b/imap_processing/tests/swapi/test_swapi_l2.py @@ -0,0 +1,21 @@ +import numpy as np + +from imap_processing.cdf.utils import write_cdf +from imap_processing.swapi.l1.swapi_l1 import swapi_l1 +from imap_processing.swapi.l2.swapi_l2 import TIME_PER_BIN, swapi_l2 + + +def test_swapi_l2_cdf(swapi_l0_test_data_path): + """Test housekeeping processing and CDF file creation""" + l0_data_path = swapi_l0_test_data_path / "imap_swapi_l0_raw_20231012_v001.pkts" + processed_data = swapi_l1(l0_data_path, data_version="v001") + l1_dataset = processed_data[0] + + l2_dataset = swapi_l2(l1_dataset, data_version="v001") + l2_cdf = write_cdf(l2_dataset) + assert l2_cdf.name == "imap_swapi_l2_sci-1min_20100101_v001.cdf" + + # Test uncertainty variables are as expected + np.testing.assert_array_equal( + l2_dataset["swp_pcem_unc"], l1_dataset["swp_pcem_err"] / TIME_PER_BIN + )