MAG CDF updates (#349)

* Updating mag CDF generation to generate two files

* PR updates, naming changes

imap_processing/cli.py

@@ -11,6 +11,7 @@
 """
 
 import argparse
+import logging
 import sys
 from abc import ABC, abstractmethod
 from json import loads
@@ -35,6 +36,8 @@
 from imap_processing.swe.l1a.swe_l1a import swe_l1a
 from imap_processing.swe.l1b.swe_l1b import swe_l1b
 
+logger = logging.getLogger(__name__)
+
 
 def _parse_args():
     """Parse the command line arguments.
@@ -313,17 +316,23 @@ def process(self):
                     f"Unexpected dependencies found for MAG L1A:"
                     f"{file_paths}. Expected only one dependency."
                 )
-            filename = imap_data_access.ScienceFilePath.generate_from_inputs(
-                "mag", "l1a", "raw", self.start_date, self.end_date, self.version
-            )
-            mag_l1a(file_paths[0], filename.construct_path())
-            print(f"Generated file: {filename.construct_path()}")
+            filename_norm = imap_data_access.ScienceFilePath.generate_from_inputs(
+                "mag", "l1a", "raw-norm", self.start_date, self.end_date, self.version
+            ).construct_path()
+            filename_burst = imap_data_access.ScienceFilePath.generate_from_inputs(
+                "mag", "l1a", "raw-burst", self.start_date, self.end_date, self.version
+            ).construct_path()
+            mag_l1a(file_paths[0], filename_norm, filename_burst)
 
             if self.upload_to_sdc:
-                print(f"Uploading file from: {filename.construct_path()}")
                 # TODO: figure out data_dir, because now this fails.
                 #  Should switch to using IMAP_DATA_DIR env var.
-                imap_data_access.upload(filename.construct_path())
+                if filename_norm.exists():
+                    logging.info(f"Uploading file: {filename_norm}")
+                    imap_data_access.upload(filename_norm)
+                if filename_burst.exists():
+                    logging.info(f"Uploading file: {filename_burst}")
+                    imap_data_access.upload(filename_burst)
 
 
 class Swapi(ProcessInstrument):

imap_processing/mag/l0/decom_mag.py

@@ -71,58 +71,117 @@ def export_to_xarray(l0_data: list[MagL0]):
     ----------
     l0_data: list[MagL0]
         A list of MagL0 datapoints
+
+    Returns
+    -------
+    norm_data : xr.Dataset
+        xarray dataset for generating burst data CDFs
+    burst_data : xr.Dataset
+        xarray dataset for generating burst data CDFs
     """
     # TODO split by mago and magi using primary sensor
-    # TODO split by norm and burst
-    norm_data = defaultdict(list)
-    burst_data = norm_data.copy()
-
-    for datapoint in l0_data:
-        if datapoint.ccsds_header.PKT_APID == Mode.NORMAL:
-            for key, value in dataclasses.asdict(datapoint).items():
-                if key != "ccsds_header":
-                    norm_data[key].append(value)
-        if datapoint.ccsds_header.PKT_APID == Mode.BURST:
-            burst_data["SHCOARSE"].append(datapoint.SHCOARSE)
-            burst_data["raw_vectors"].append(datapoint.VECTORS)
+    norm_data = []
+    burst_data = []
+
+    for packet in l0_data:
+        if packet.ccsds_header.PKT_APID == Mode.NORMAL:
+            norm_data.append(packet)
+        if packet.ccsds_header.PKT_APID == Mode.BURST:
+            burst_data.append(packet)
+
+    norm_dataset = None
+    burst_dataset = None
+
+    if len(norm_data) > 0:
+        norm_dataset = generate_dataset(norm_data)
+    if len(burst_data) > 0:
+        burst_dataset = generate_dataset(burst_data)
+
+    return norm_dataset, burst_dataset
+
+
+def generate_dataset(l0_data: list[MagL0]):
+    """
+    Generate a CDF dataset from the sorted L0 MAG data.
+
+    Used to create 2 similar datasets, for norm and burst data.
+
+    Parameters
+    ----------
+    l0_data : list[MagL0]
+        List of sorted L0 MAG data.
+
+    Returns
+    -------
+    dataset : xr.Dataset
+        xarray dataset with proper CDF attributes and shape.
+    """
+    vector_data = np.zeros((len(l0_data), len(l0_data[0].VECTORS)))
+    shcoarse_data = np.zeros(len(l0_data))
+
+    support_data = defaultdict(list)
+
+    for index, datapoint in enumerate(l0_data):
+        vector_len = len(datapoint.VECTORS)
+        if vector_len > vector_data.shape[1]:
+            # If the new vector is longer than the existing shape, first reshape
+            # vector_data and pad the existing vectors with zeros.
+            vector_data = np.pad(
+                vector_data,
+                (
+                    (
+                        0,
+                        0,
+                    ),
+                    (0, vector_len - vector_data.shape[1]),
+                ),
+                "constant",
+                constant_values=(0,),
+            )
+        vector_data[index, :vector_len] = datapoint.VECTORS
+
+        shcoarse_data[index] = calc_start_time(datapoint.SHCOARSE)
+
+        # Add remaining pieces to arrays
+        for key, value in dataclasses.asdict(datapoint).items():
+            if key not in ("ccsds_header", "VECTORS", "SHCOARSE"):
+                support_data[key].append(value)
 
     # Used in L1A vectors
-    direction_norm = xr.DataArray(
-        np.arange(len(norm_data["VECTORS"][0])),
+    direction = xr.DataArray(
+        np.arange(vector_data.shape[1]),
         name="Direction",
         dims=["Direction"],
         attrs=mag_cdf_attrs.direction_attrs.output(),
     )
 
-    norm_epoch_time = xr.DataArray(
-        [calc_start_time(shcoarse) for shcoarse in norm_data["SHCOARSE"]],
+    epoch_time = xr.DataArray(
+        shcoarse_data,
         name="Epoch",
         dims=["Epoch"],
         attrs=ConstantCoordinates.EPOCH,
     )
 
     # TODO: raw vectors units
-    norm_raw_vectors = xr.DataArray(
-        norm_data["VECTORS"],
+    raw_vectors = xr.DataArray(
+        vector_data,
         name="Raw_Vectors",
         dims=["Epoch", "Direction"],
         attrs=mag_cdf_attrs.mag_vector_attrs.output(),
     )
 
     # TODO add norm to attrs somehow
-    norm_dataset = xr.Dataset(
-        coords={"Epoch": norm_epoch_time, "Direction": direction_norm},
+    output = xr.Dataset(
+        coords={"Epoch": epoch_time, "Direction": direction},
         attrs=mag_cdf_attrs.mag_l1a_attrs.output(),
     )
 
-    norm_dataset["RAW_VECTORS"] = norm_raw_vectors
-
-    # TODO: retrieve the doc for the CDF description (etattr(MagL0, "__doc__", {}))
+    output["RAW_VECTORS"] = raw_vectors
 
-    for key, value in norm_data.items():
+    for key, value in support_data.items():
         # Time varying values
         if key not in ["SHCOARSE", "VECTORS"]:
-            norm_datarray = xr.DataArray(
+            output[key] = xr.DataArray(
                 value,
                 name=key,
                 dims=["Epoch"],
@@ -135,6 +194,5 @@ def export_to_xarray(l0_data: list[MagL0]):
                     display_type="no_plot",
                 ).output(),
             )
-            norm_dataset[key] = norm_datarray
 
-    return norm_dataset
+    return output

imap_processing/mag/l1a/mag_l1a.py

@@ -1,23 +1,34 @@
 """Methods for decomming packets, processing to level 1A, and writing CDFs for MAG."""
+import logging
 from pathlib import Path
 
 from imap_processing.cdf.utils import write_cdf
 from imap_processing.mag.l0 import decom_mag
 
+logger = logging.getLogger(__name__)
 
-def mag_l1a(packet_filepath, output_filepath):
+
+def mag_l1a(packet_filepath, output_filepath_norm, ouptput_filepath_burst):
     """
     Process MAG L0 data into L1A CDF files at cdf_filepath.
 
     Parameters
     ----------
-    packet_filepath:
+    packet_filepath :
         Packet files for processing
-    output_filepath:
-        Directory for output
+    output_filepath_norm :
+        Full directory and filename for raw-norm CDF file
+    ouptput_filepath_burst :
+        Full directory and filename for raw-burst CDF file
     """
     mag_l0 = decom_mag.decom_packets(packet_filepath)
 
-    mag_datasets = decom_mag.export_to_xarray(mag_l0)
+    mag_norm, mag_burst = decom_mag.export_to_xarray(mag_l0)
+
+    if mag_norm is not None:
+        write_cdf(mag_norm, Path(output_filepath_norm))
+        logging.info(f"Created CDF file at {output_filepath_norm}")
 
-    write_cdf(mag_datasets, Path(output_filepath))
+    if mag_burst is not None:
+        write_cdf(mag_burst, Path(ouptput_filepath_burst))
+        logging.info(f"Created CDF file at {output_filepath_norm}")

imap_processing/mag/mag_cdf_attrs.py

@@ -26,7 +26,9 @@
 )
 
 mag_l1a_attrs = GlobalDataLevelAttrs(
+    # TODO: data type should include "norm" and "burst" L1A-norm>Level-1A-normal-rate
     "L1A>Level-1A",
+    # Should also include data type
     logical_source="imap_mag_l1a",
     logical_source_desc="IMAP Mission MAG Instrument Level-1A Data.",
     instrument_base=mag_base,
@@ -50,6 +52,8 @@
 
 # TODO: display type, catdesc, units, format, label_axis
 
+# TODO: update descriptor to be more accurate for L1A raw
+# TODO: does raw value need "counts"
 mag_vector_attrs = ScienceAttrs(
     validmin=GlobalConstants.INT_FILLVAL,
     validmax=GlobalConstants.INT_MAXVAL,

imap_processing/tests/mag/test_mag_decom.py

@@ -1,8 +1,10 @@
 from pathlib import Path
 
 import pandas as pd
+from cdflib.xarray import cdf_to_xarray
 
 from imap_processing.cdf import global_attrs
+from imap_processing.cdf.utils import write_cdf
 from imap_processing.mag.l0.decom_mag import decom_packets, export_to_xarray
 
 
@@ -12,7 +14,6 @@ def test_mag_decom():
     l0 = decom_packets(burst_test_file)
 
     expected_output = pd.read_csv(current_directory / "mag_l0_test_output.csv")
-
     for index, test in enumerate(l0):
         assert test.ccsds_header.PKT_APID == expected_output["PHAPID"][index]
         assert test.ccsds_header.SRC_SEQ_CTR == expected_output["PHSEQCNT"][index]
@@ -35,15 +36,50 @@ def test_mag_decom():
     assert len(l0) == len(expected_output.index)
 
 
-def test_mag_raw_cdf():
+def test_mag_raw_xarray():
     current_directory = Path(__file__).parent
     burst_test_file = current_directory / "mag_l0_test_data.pkts"
     l0 = decom_packets(str(burst_test_file))
 
-    output_data = export_to_xarray(l0)
+    norm_data, burst_data = export_to_xarray(l0)
     required_attrs = list(
         global_attrs.GlobalInstrumentAttrs("", "", "").output().keys()
     )
 
-    assert all([item in list(output_data.attrs.keys()) for item in required_attrs])
-    assert all([item is not None for _, item in output_data.attrs.items()])
+    assert all([item in list(norm_data.attrs.keys()) for item in required_attrs])
+    assert all([item is not None for _, item in norm_data.attrs.items()])
+
+    assert all([item in list(burst_data.attrs.keys()) for item in required_attrs])
+    assert all([item is not None for _, item in burst_data.attrs.items()])
+
+    expected_norm_len = 17
+    assert norm_data.dims["Epoch"] == expected_norm_len
+
+    expected_burst_len = 19
+    assert burst_data.dims["Epoch"] == expected_burst_len
+
+
+def test_mag_raw_cdf_generation(tmp_path):
+    current_directory = Path(__file__).parent
+    test_file = current_directory / "mag_l0_test_data.pkts"
+    l0 = decom_packets(str(test_file))
+
+    norm_data, burst_data = export_to_xarray(l0)
+
+    test_data_path_norm = tmp_path / "mag_l1a_raw-normal_20210101_20210102_v01-01.cdf"
+
+    assert not test_data_path_norm.exists()
+    output = write_cdf(norm_data, test_data_path_norm)
+    assert test_data_path_norm.exists()
+
+    input_xarray = cdf_to_xarray(output)
+    assert input_xarray.attrs.keys() == norm_data.attrs.keys()
+
+    test_data_path_burst = tmp_path / "mag_l1a_raw-burst_20210101_20210102_v01-01.cdf"
+
+    assert not test_data_path_burst.exists()
+    output = write_cdf(burst_data, test_data_path_burst)
+    assert test_data_path_burst.exists()
+
+    input_xarray = cdf_to_xarray(output)
+    assert input_xarray.attrs.keys() == burst_data.attrs.keys()