Merge master into dev

benchmarks/readme.md

@@ -17,6 +17,7 @@ pytest benchmarks/ --codspeed
 pytest benchmarks/test_patch_benchmarks.py --codspeed
 pytest benchmarks/test_io_benchmarks.py --codspeed
 pytest benchmarks/test_spool_benchmarks.py --codspeed
+pytest benchmarks/test_import_benchmarks.py --codspeed
 ```
 
 ## Benchmark Structure
@@ -27,6 +28,7 @@ Benchmarks are now organized as pytest tests in the `benchmarks/` directory:
 - `test_io_benchmarks.py` - File I/O operations benchmarks
 - `test_spool_benchmarks.py` - Spool chunking and selection benchmarks
 - `test_lookup_benchmarks.py` - In-memory lookups on hot paths (format resolution, remote-cache and IO handle resolution, repeat spool access). These are deliberately small: a change of a few microseconds per lookup is invisible in the end-to-end benchmarks above, because one file read costs far more than the lookups it makes.
+- `test_import_benchmarks.py` - Import benchmarks (dascore's own modules; third party dependencies stay warm)
 
 Each benchmark uses the `@pytest.mark.benchmark` decorator to automatically measure performance.
 

benchmarks/test_import_benchmarks.py

@@ -0,0 +1,53 @@
+"""Benchmarks for importing dascore using pytest-codspeed."""
+
+from __future__ import annotations
+
+import importlib
+import sys
+import warnings
+
+import pint
+import pytest
+
+
+def _dascore_module_names():
+    """Get the names of all currently imported dascore modules."""
+    return [x for x in sys.modules if x == "dascore" or x.startswith("dascore.")]
+
+
+class TestImportBenchmarks:
+    """
+    Benchmarks for re-executing dascore's modules.
+
+    Only dascore's own modules are removed from the module cache, so these
+    measure the cost of executing dascore's module bodies, not the one-time
+    cost of importing third party dependencies. A fresh interpreter (eg a
+    CLI call) also pays the latter, but it cannot be measured here; a
+    subprocess falls outside the region CodSpeed instruments. Instead, see
+    tests/test_imports.py for the guards which keep slow dependencies out of
+    the import chain entirely.
+    """
+
+    @pytest.fixture()
+    def restore_dascore_modules(self):
+        """Put the original dascore modules back after re-importing them."""
+        # Import here so third party dependencies are warm before timing,
+        # even when this file is the only one collected.
+        importlib.import_module("dascore")
+        saved = {x: sys.modules[x] for x in _dascore_module_names()}
+        # Re-importing dascore makes (and installs) a new pint registry.
+        registry = pint.get_application_registry().get()
+        # dascore adds a warning filter on import; catch_warnings undoes that.
+        with warnings.catch_warnings():
+            yield
+        for name in _dascore_module_names():
+            del sys.modules[name]
+        sys.modules.update(saved)
+        pint.set_application_registry(registry)
+
+    @pytest.mark.benchmark
+    def test_reimport_dascore(self, restore_dascore_modules):
+        """Time re-importing the top-level dascore module."""
+        for name in _dascore_module_names():
+            del sys.modules[name]
+        importlib.import_module("dascore")

benchmarks/test_patch_benchmarks.py

@@ -190,6 +190,48 @@ def test_slope_mute(self, example_patch):
         patch.slope_mute(slopes=(1000, 3000))
 
 
+class TestPatchConstructionBenchmarks:
+    """Benchmarks for the patch construction paths."""
+
+    @pytest.fixture(scope="module")
+    def new_data(self, example_patch):
+        """Data for constructing new patches; built outside the timed call."""
+        return np.asarray(example_patch.data) * 2
+
+    @pytest.fixture(scope="module")
+    def decimated(self, example_patch):
+        """A coord manager with one dimension shortened."""
+        coord = example_patch.get_coord("time")[::2]
+        return example_patch.coords.update(time=coord)
+
+    @pytest.mark.benchmark
+    def test_new_data_only(self, example_patch, new_data):
+        """Time new when only data changes; coords and attrs are reused."""
+        example_patch.new(data=new_data)
+
+    @pytest.mark.benchmark
+    def test_new_with_coords(self, example_patch, decimated):
+        """Time new when the coords change, so attrs must be rebuilt."""
+        example_patch.new(data=example_patch.data[:, ::2], coords=decimated)
+
+    @pytest.mark.benchmark
+    def test_new_with_coords_and_attrs(self, example_patch, new_data):
+        """Time new when both coords and attrs are passed."""
+        patch = example_patch
+        patch.new(data=new_data, coords=patch.coords, attrs=patch.attrs)
+
+    @pytest.mark.benchmark
+    def test_patch_init(self, example_patch, new_data):
+        """
+        Time the normal constructor.
+
+        This is a control; it should not move, since the strict path is
+        deliberately left alone.
+        """
+        patch = example_patch
+        dc.Patch(data=new_data, coords=patch.coords, dims=patch.dims, attrs=patch.attrs)
+
+
 class TestTransformBenchmarks:
     """Benchmarks for patch transform operations."""
 
@@ -353,6 +395,11 @@ def test_rolling_large_roller_mean(self, big_roller):
         """Time rolling mean calculation."""
         big_roller.mean()
 
+    @pytest.mark.benchmark
+    def test_rolling_mean_full_call(self, example_patch):
+        """Time a complete rolling mean, including roller construction."""
+        example_patch.rolling(time=5, samples=True, center=True).mean()
+
 
 class TestAlignBenchmarks:
     """Benchmarks for align_to_coord operation."""

dascore/core/coords.py

@@ -163,6 +163,18 @@ def _reduce_time_like(func, data):
     return np.atleast_1d(out)
 
 
+def _validate_new_length(length) -> int:
+    """Ensure a requested coordinate length is a non-negative integer."""
+    # bool is an int subclass; True/False are never a sensible length.
+    if isinstance(length, bool) or not isinstance(length, int | np.integer):
+        msg = f"change_length requires an integer length, not {length!r}."
+        raise ParameterError(msg)
+    if length < 0:
+        msg = f"change_length requires a non-negative length, not {length}."
+        raise ParameterError(msg)
+    return int(length)
+
+
 def _get_dtype(value, dtype):
     """Get the data type based on the first argument."""
     if dtype is not None and dtype != "":
@@ -399,7 +411,7 @@ def _select_by_sample_array(self, array):
             msg = "Using an array input for select with samples requires integer dtype."
             raise CoordError(msg)
         # Filter out bad indices
-        if self.ndim > 1:
+        if self.ndim != 1:
             msg = "Select only works on 1D coords."
             raise CoordError(msg)
         inds = np.arange(len(self))
@@ -1148,7 +1160,12 @@ def change_length(self, length: int) -> Self:
         Parameters
         ----------
         length
-            The output length.
+            The output length. Must be a non-negative integer.
+
+        Raises
+        ------
+        ParameterError
+            If length is not a non-negative integer.
         """
         msg = f"Coordinate type {self.__class__} does not implement change_length"
         raise NotImplementedError(msg)
@@ -1308,7 +1325,7 @@ def change_length(self, length: int) -> Self:
         if self.ndim != 1:
             msg = "change_length only works on 1D coords."
             raise CoordError(msg)
-        return get_coord(shape=(length,))
+        return get_coord(shape=(_validate_new_length(length),))
 
     def to_summary(self, dims=()) -> CoordSummary:
         """Get the summary info about the coord."""
@@ -1668,6 +1685,7 @@ def change_length(self, length: int) -> Self:
         """
         # CoordRange is always 1D by construction; keep as an internal invariant.
         assert self.ndim == 1, "Can only change length for 1D coords."
+        length = _validate_new_length(length)
         if len(self) == length:
             return self
         # Only the sample count changes; start/step are already valid.

dascore/io/dasvader/core.py

@@ -25,8 +25,9 @@ class DASVaderV1(FiberIO):
     Notes
     -----
     Legacy DASVader files may contain anonymous JLD2 object references. DASCore
-    detects those files and raises `DASVaderCompatibilityError` with compatibility
-    instructions instead of failing inside `h5py`. A known working stack for
+    reads these references when supported by HDF5 and raises
+    `DASVaderCompatibilityError` with compatibility instructions when
+    dereferencing fails. A known working stack for
     such legacy files is `h5py<3.16` with `HDF5 1.14.x`.
     """
 

dascore/io/dasvader/utils.py

@@ -4,7 +4,6 @@
 
 import h5py
 import numpy as np
-from h5py import h5r
 from h5py.h5r import Reference, dereference
 
 import dascore as dc
@@ -77,15 +76,18 @@ def _raise_legacy_ref_error(h5, field_name: str) -> None:
 
 
 def _dereference(h5, value, field_name: str):
-    """Resolve an HDF5 reference, rejecting legacy anonymous DASVader refs."""
+    """Resolve an HDF5 reference or raise a clear compatibility error."""
     if not isinstance(value, Reference):
         return value
-    if h5r.get_name(value, h5.id) is None:
-        _raise_legacy_ref_error(h5, field_name)
     try:
         return h5[value]
     except KeyError:
-        return h5py.Dataset(dereference(value, h5.id))
+        # The high-level lookup fails for some references HDF5 can still
+        # resolve directly, so try that before giving up.
+        try:
+            return h5py.Dataset(dereference(value, h5.id))
+        except Exception:
+            _raise_legacy_ref_error(h5, field_name)
 
 
 # --- Metadata parsing

dascore/proc/basic.py

@@ -318,6 +318,10 @@ def normalize(
     """
     Normalize a patch along a specified dimension.
 
+    NaN values are ignored when computing the norm. They remain NaN in the
+    output but do not affect any other sample. Slices with a norm of zero,
+    meaning they contain nothing but zeros and NaN, are returned unscaled.
+
     Parameters
     ----------
     dim
@@ -348,26 +352,24 @@ def normalize(
     data = self.data
     if norm in {"l1", "l2"}:
         order = int(norm[-1])
-        norm_values = np.linalg.norm(self.data, axis=axis, ord=order)
+        # Equivalent to np.linalg.norm, but skips NaN rather than letting a
+        # single null blank every sample sharing its slice. The float exponent
+        # promotes ints so the powers cannot overflow a narrow dtype.
+        norm_values = np.nansum(np.abs(data) ** float(order), axis=axis) ** (1 / order)
+        divisor = np.expand_dims(norm_values, axis=axis)
     elif norm == "max":
-        norm_values = np.max(np.abs(data), axis=axis)
+        divisor = np.expand_dims(np.nanmax(np.abs(data), axis=axis), axis=axis)
     elif norm == "bit":
-        pass
+        divisor = np.abs(data)
     else:
         msg = (
             f"Norm value of {norm} is not supported. "
             f"Supported values are {('l1', 'l2', 'max', 'bit')}"
         )
         raise ValueError(msg)
-    if norm == "bit":
-        new_data = np.divide(
-            data, np.abs(data), out=np.zeros_like(data), where=np.abs(data) != 0
-        )
-    else:
-        expanded_norm = np.expand_dims(norm_values, axis=axis)
-        new_data = np.divide(
-            data, expanded_norm, out=np.zeros_like(data), where=expanded_norm != 0
-        )
+    # A zero divisor means there is nothing but zeros and nulls to scale, so
+    # divide those by one; the zeros stay zero and the nulls stay null.
+    new_data = data / np.where(divisor == 0, 1, divisor)
     return self.new(data=new_data)
 
 
@@ -385,6 +387,9 @@ def standardize(
     where u is the mean of the training samples or zero if with_mean=False,
     and s is the standard deviation of the training samples or one if with_std=False.
 
+    NaN values are ignored when computing the mean and standard deviation. They
+    remain NaN in the output but do not affect any other sample.
+
     Parameters
     ----------
     dim
@@ -406,8 +411,8 @@ def standardize(
     """
     axis = self.get_axis(dim)
     data = self.data
-    mean = np.mean(data, axis=axis, keepdims=True)
-    std = np.std(data, axis=axis, keepdims=True)
+    mean = np.nanmean(data, axis=axis, keepdims=True)
+    std = np.nanstd(data, axis=axis, keepdims=True)
     new_data = (data - mean) / std
     return self.new(data=new_data)
 
@@ -810,6 +815,10 @@ def demedian(patch, dim: str = "time"):
     """
     Remove the median along a given dimension of a DASCore patch.
 
+    NaN values are ignored when computing the median, consistent with
+    [Patch.median](`dascore.proc.aggregate.median`). They remain NaN in the
+    output but do not affect any other sample.
+
     Parameters
     ----------
     patch :
@@ -856,7 +865,7 @@ def demedian(patch, dim: str = "time"):
     data = patch.data
 
     # Compute median along axis, keep dims for broadcasting
-    med = np.median(data, axis=axis, keepdims=True)
+    med = np.nanmedian(data, axis=axis, keepdims=True)
 
     new_data = data - med
 
@@ -869,6 +878,10 @@ def demean(patch, dim: str = "time"):
     """
     Remove the mean along a given dimension of a DASCore patch.
 
+    NaN values are ignored when computing the mean, consistent with
+    [Patch.mean](`dascore.proc.aggregate.mean`). They remain NaN in the output
+    but do not affect any other sample.
+
     Parameters
     ----------
     patch :
@@ -915,7 +928,7 @@ def demean(patch, dim: str = "time"):
     data = patch.data
 
     # Compute mean along axis, keep dims for broadcasting
-    mea = np.mean(data, axis=axis, keepdims=True)
+    mea = np.nanmean(data, axis=axis, keepdims=True)
 
     new_data = data - mea