Toast data and metadata improvements: (#1522)

- Support loading raw data as int32 DAQ units.

- When loading metadata, ensure that byte arrays are converted
  to strings as the starting point.  These get automatically
  converted to / from bytestrings during HDF5 save / load.

- The HKManager now excludes several large fields by default
  (WG / HWP encoders and ACU udp stream).

- HKManager converts all byte strings to unicode on load, and
  then separately converts to / from bytestrings during HDF5
  save / load.

- Small fix to always use "so3g compatibility mode" when
  converting boresight Az/El to RA/DEC on load.

Author: tskisner

sotodlib/toast/hkmanager.py

@@ -2,13 +2,17 @@
 # Full license can be found in the top level "LICENSE" file.
 
 import copy
+import re
 from collections.abc import MutableMapping
 
 import numpy as np
 from scipy.interpolate import CubicSpline, PchipInterpolator
 
-from toast.utils import Logger
-from toast.io.hdf_utils import load_meta_object, save_meta_object
+from toast.utils import Logger, replace_byte_arrays, array_equal
+from toast.io.hdf_utils import (
+    load_meta_object,
+    save_meta_object,
+)
 
 from ..io import hkdb
 
@@ -49,6 +53,17 @@ class HKManager(MutableMapping):
 
     """
 
+    # Exclude the massive HWP and WG housekeeping data by default, since
+    # this information is only used to build the HWP angle and WG calibration
+    # metadata.  Loading these fields dramatically increases the load time
+    # for normal use.  These can still be forcibly loaded by specifying them
+    # in the `plat_fields` list.
+    exclude_by_default = [
+        r"^hwp-.*",
+        r"^wg-.*",
+        r"^acu.acu_udp_stream.*",
+    ]
+
     def __init__(
         self,
         comm=None,
@@ -83,7 +98,7 @@ def __init__(
                 # Load data on just one process per group
                 self._load_hk_data(site_root, site_db, site_fields, site_aliases)
 
-        # Load the plot HK data
+        # Load the platform HK data
         if plat_root is not None:
             if plat_db is None:
                 raise RuntimeError("If plat_root is specified, plat_db is required")
@@ -99,6 +114,7 @@ def __init__(
     def _load_hk_data(self, root, db, fields, aliases):
         """Helper function to do the actual loading."""
         log = Logger.get()
+        exclude_pat = [re.compile(x) for x in self.exclude_by_default]
         conf = hkdb.HkConfig.from_dict(
             {
                 "hk_root": root,
@@ -120,7 +136,17 @@ def _load_hk_data(self, root, db, fields, aliases):
         )
         all_fields = hkdb.get_feed_list(test_spec)
         if fields is None or len(fields) == 0:
-            selected = list(all_fields)
+            # The user is requesting "all" fields.  We will prune this list to exclude
+            # fields related to the HWP and wire grid encoders, which are excessively
+            # large and rarely used for downstream analysis.
+            selected = list()
+            for fld in all_fields:
+                keep = True
+                for pat in exclude_pat:
+                    if pat.match(fld) is not None:
+                        keep = False
+                if keep:
+                    selected.append(fld)
         else:
             selected = fields
         # Load the data
@@ -133,7 +159,10 @@ def _load_hk_data(self, root, db, fields, aliases):
             hkdb=None,
         )
         result = hkdb.load_hk(lspec, show_pb=False)
-        self._internal.update(result.data)
+
+        converted = replace_byte_arrays(result.data)
+
+        self._internal.update(converted)
         if aliases is not None:
             for k, v in aliases.items():
                 if v not in self._internal:
@@ -159,12 +188,24 @@ def interp_cubic(self, field):
         itrp = CubicSpline(times, vals, extrapolate=True)
         return itrp(self._stamps)
 
-    def memory_use(self):
+    def memory_use(self, threshold=0):
+        """Return the total size of all fields and a string with the breakdown.
+
+        Args:
+            threshold (int):  Only include in the string objects larger than this.
+
+        Returns:
+            (tuple):  The (total size, info string)
+
+        """
         bytes = 0
+        fstr = ""
         for k, v in self._internal.items():
-            bytes += v[0].nbytes
-            bytes += v[1].nbytes
-        return bytes
+            fbytes = v[0].nbytes + v[1].nbytes
+            bytes += fbytes
+            if fbytes > threshold:
+                fstr += f"{k}: {fbytes} bytes\n"
+        return bytes, fstr
 
     # Mapping methods
 
@@ -204,14 +245,14 @@ def __eq__(self, other):
             log.verbose(f"  keys {self._internal} != {other._internal}")
             return False
         for k, v in self._internal.items():
-            try:
-                if not np.allclose(v, other._internal[k]):
-                    log.verbose(f"  data array {v} != {other._internal[k]}")
-                    return False
-            except Exception:
-                if self._internal[k] != other._internal[k]:
-                    log.verbose(f"  data value {v} != {other._internal[k]}")
-                    return False
+            times = v[0]
+            vals = v[1]
+            other_times = other._internal[k][0]
+            other_vals = other._internal[k][1]
+            if not array_equal(times, other_times, log_prefix=f"HKManager {k} times"):
+                return False
+            if not array_equal(vals, other_vals, log_prefix=f"HKManager {k} values"):
+                return False
         if self._aliases != other._aliases:
             log.verbose(f"  aliases {self._aliases} != {other._aliases}")
             return False

sotodlib/toast/ops/load_context.py

@@ -24,8 +24,7 @@
     Float,
 )
 from toast.ops.operator import Operator
-from toast.ops.pipeline import Pipeline
-from toast.utils import Logger
+from toast.utils import Logger, replace_byte_arrays
 from toast.dist import distribute_discrete
 from toast.observation import default_values as defaults
 
@@ -135,13 +134,9 @@ class LoadContext(Operator):
         help="Path to DB for site housekeeping",
     )
 
-    hk_site_fields = List(
-        list(), help="Restrict loading to only these site fields"
-    )
+    hk_site_fields = List(list(), help="Restrict loading to only these site fields")
 
-    hk_site_aliases = Dict(
-        dict(), help="Optional convenience aliases for site fields"
-    )
+    hk_site_aliases = Dict(dict(), help="Optional convenience aliases for site fields")
 
     hk_platform_root = Unicode(
         None,
@@ -330,6 +325,10 @@ class LoadContext(Operator):
 
     bandwidth = Float(0.2, help="Fractional bandwith used in analytic bandpass")
 
+    daq_units = Bool(
+        False, help="If True, convert raw data to original int32 DAQ units"
+    )
+
     def __init__(self, **kwargs):
         super().__init__(**kwargs)
 
@@ -355,6 +354,10 @@ def _exec(self, data, detectors=None, **kwargs):
                 msg = "Only one of the context or context_file should be specified"
                 raise RuntimeError(msg)
 
+        if self.daq_units and self.preprocess_config is not None:
+            msg = "Cannot convert raw signal to DAQ units if passing data through preprocessing"
+            raise RuntimeError(msg)
+
         # Build our detector selection dictionary.  Merge our explicit traits
         # with any pre-existing detector selection.
         dets_select = None
@@ -635,8 +638,14 @@ def _exec(self, data, detectors=None, **kwargs):
             # Read and communicate data
             self._load_data(ob, have_pointing, preproc_conf)
 
+            # Now that all metadata has been loaded, ensure that all byte strings
+            # are converted to unicode arrays.
+            ob._internal = replace_byte_arrays(ob._internal)
+
             # Optionally load housekeeping data
             if self.hk_site_root is not None or self.hk_platform_root is not None:
+                hk_timer = Timer()
+                hk_timer.start()
                 ob.hk = HKManager(
                     ob.comm.comm_group,
                     ob.shared[self.times].data,
@@ -649,6 +658,11 @@ def _exec(self, data, detectors=None, **kwargs):
                     plat_fields=self.hk_platform_fields,
                     plat_aliases=self.hk_platform_aliases,
                 )
+                log.debug_rank(
+                    f"LoadContext {obs_name} load House Keeping in",
+                    comm=comm.comm_group,
+                    timer=hk_timer,
+                )
 
             # Compute the boresight pointing and observatory position
             if have_pointing:
@@ -1052,9 +1066,14 @@ def _create_observation(
                 dtype=np.float64,
             )
         if ax_det_signal is not None:
-            ob.detdata.create(
-                self.det_data, dtype=np.float64, units=self.det_data_units
-            )
+            if self.daq_units:
+                ob.detdata.create(
+                    self.det_data, dtype=np.int32, units=u.dimensionless_unscaled
+                )
+            else:
+                ob.detdata.create(
+                    self.det_data, dtype=np.float64, units=self.det_data_units
+                )
             ob.detdata.create(self.det_flags, dtype=np.uint8)
 
         if meta is not None:
@@ -1118,6 +1137,7 @@ def _load_data(self, ob, have_pointing, pconf):
             ignore_preprocess_archive=self.ignore_preprocess_archive,
             context=self.context,
             context_file=self.context_file,
+            daq_units=self.daq_units,
         )
 
         log.debug_rank(
@@ -1315,6 +1335,7 @@ def _load_data(self, ob, have_pointing, pconf):
                     is_flag=(mask is not None),
                     flag_invert=do_invert,
                     flag_mask=mask,
+                    daq_units=self.daq_units,
                 )
 
         # Original wafer data no longer needed.  AxisManager does not seem to