diff --git a/dascore/proc/basic.py b/dascore/proc/basic.py index a802ce9aa..610511d8e 100644 --- a/dascore/proc/basic.py +++ b/dascore/proc/basic.py @@ -22,7 +22,16 @@ from dascore.exceptions import ParameterError from dascore.models import ArrayLike from dascore.utils.array import _apply_binary_ufunc -from dascore.utils.array_api import array_namespace, nan_reduce +from dascore.utils.array_api import ( + array_namespace, + asarray_like, + backend_name, + device, + is_numpy, + nan_reduce, + to_numpy, + warn_numpy_fallback, +) from dascore.utils.misc import _get_nullish from dascore.utils.patch import ( align_patch_coords, @@ -306,7 +315,7 @@ def abs(patch: PatchType) -> PatchType: >>> pa = dascore.get_example_patch() # generate example patch >>> out = pa.abs() # take absolute value of generated example patch data """ - return Abs()._apply(patch) + return Abs._call(patch) class Abs(PatchProcessor): @@ -334,7 +343,7 @@ def conj(patch: PatchType) -> PatchType: >>> dft = pa.dft(None) # multi-dim dft >>> conj = dft.conj() """ - return Conj()._apply(patch) + return Conj._call(patch) class Conj(PatchProcessor): @@ -366,7 +375,7 @@ def real(patch: PatchType) -> PatchType: >>> pa = dascore.get_example_patch() >>> out = pa.real() """ - return Real()._apply(patch) + return Real._call(patch) class Real(PatchProcessor): @@ -393,7 +402,7 @@ def imag(patch: PatchType) -> PatchType: >>> pa = dascore.get_example_patch() >>> out = pa.imag() """ - return Imag()._apply(patch) + return Imag._call(patch) class Imag(PatchProcessor): @@ -470,7 +479,7 @@ def normalize( >>> # Bit normalization (sign only) >>> bit_norm = patch.normalize(dim="time", norm="bit") """ - return Normalize(dim=dim, norm=norm)._apply(self) + return Normalize._call(self, dim=dim, norm=norm) class Normalize(PatchProcessor): @@ -552,7 +561,7 @@ def standardize( standardized_distance = patch.standardize('distance') ``` """ - return Standardize(dim=dim)._apply(self) + return Standardize._call(self, dim=dim) class Standardize(PatchProcessor): @@ -674,16 +683,83 @@ def fillna(patch: PatchType, value, include_inf=True) -> PatchType: >>> # Replace all occurrences of NaN with 5 >>> out = patch.fillna(5) """ - if include_inf: - to_replace = ~np.isfinite(patch.data) - else: - to_replace = pd.isnull(patch.data) - if not np.any(to_replace): # nothing nullish to fill - return patch - new_data = patch.data.copy() - new_data[to_replace] = value + return FillNa._call(patch, value=value, include_inf=include_inf) - return patch.new(data=new_data) + +class FillNa(PatchProcessor): + """Put a value where the data has none.""" + + value: Any + include_inf: bool = True + + @property + def needs_numpy(self) -> bool: + """ + Whether the value has a shape, or nothing means more than non-finite. + + A value with a shape is spent positionally, one element per null, + which `where` cannot say -- it would broadcast. And + `include_inf=False` asks `pandas.isnull` what counts as nothing, + which is not a question a backend answers. + + A value `np.ndim` refuses -- a ragged nested list -- is answered + for rather than raised on, so that a patch with nothing to fill + stays the no-op it has always been and the complaint comes from + numpy, at the fill, exactly as it used to. + """ + with suppress(ValueError): + return not self.include_inf or bool(np.ndim(self.value)) + return True + + def kernel(self, data, meta, out_meta): + """ + Return the data with its nulls filled, or the data unchanged. + + `where` rather than assigning into a copy: writing into an array + by boolean mask is not something the standard asks a backend for. + Handing the data straight back where there is nothing to fill is + what says the operation did nothing. + """ + xp = array_namespace(data) + replace = self._nulls(data) + if not xp.any(replace): + return data + value = xp.asarray(self.value, dtype=data.dtype, device=device(data)) + return xp.where(replace, value, data) + + def _nulls(self, data): + """ + Return where the data has nothing, in the backend's own terms. + + A boolean array is not asked: `isfinite` is undefined for one on + a strict backend, and where it is defined the answer is that + every value is finite, which is also true. + """ + xp = array_namespace(data) + if xp.isdtype(data.dtype, "bool"): + return xp.zeros(data.shape, dtype=xp.bool, device=device(data)) + return ~xp.isfinite(data) + + def numpy_kernel(self, data, meta, out_meta): + """ + Fill with numpy, for the two things the standard cannot say. + + A value with a shape is spent positionally, one element per null, + where `where` would broadcast it across the whole array. And + `pandas.isnull` is what `include_inf=False` means by nothing. + """ + if not is_numpy(data): + warn_numpy_fallback("fillna", backend_name(data)) + array = to_numpy(data) + replace = ~np.isfinite(array) if self.include_inf else pd.isnull(array) + if not np.any(replace): + return data + filled = np.array(array) + filled[replace] = self.value + return asarray_like(filled, data) + + +register_implementation("fillna", FillNa) @patch_function() @@ -941,12 +1017,44 @@ def flip(patch, *dims, flip_coords=True): >>> # Flip patch over all dimensions. >>> out = patch.flip(*patch.dims) """ - if not dims: - return patch # no-op - axes = tuple(patch.get_axis(name) for name in dims) - data = np.flip(patch.data, axis=axes) if dims else patch.data - coords = patch.coords.flip(*dims) if flip_coords else patch.coords - return patch.new(data=data, coords=coords) + return Flip._call(patch, dims=tuple(dims), flip_coords=flip_coords) + + +class Flip(PatchProcessor): + """Reverse a patch along one or more of its dimensions.""" + + dims: tuple[Any, ...] = () + flip_coords: bool = True + + def derive_meta(self, meta): + """ + Return the coordinates reversed along the same dimensions. + + Named no dimensions, the operation has nothing to reverse: the + metadata comes back untouched and so does the data, and it takes + both for `_apply` to hand the patch back. `flip_coords=False` + leaves the metadata alone too, but the data is still reversed, so + that one is a new patch wearing its old coordinates. + """ + if not self.dims: + return meta + # Resolved even when the coordinates stay put, so that a name the + # patch does not have is refused in the same terms either way. + for name in self.dims: + meta.get_axis(name) + if not self.flip_coords: + return meta + return meta.update(coords=meta.coords.flip(*self.dims)) + + def kernel(self, data, meta, out_meta): + """Return the data reversed along the axes the dimensions name.""" + if not self.dims: + return data + axes = tuple(meta.get_axis(name) for name in self.dims) + return array_namespace(data).flip(data, axis=axes) + + +register_implementation("flip", Flip) @patch_function(data_type="") @@ -973,8 +1081,54 @@ def full(patch, fill_value): >>> # Same thing, except for 0s. >>> zero_patch = patch.full(0.0) """ - array = np.full(patch.data.shape, fill_value) - return patch.update(data=array) + return Full._call(patch, fill_value=fill_value) + + +class Full(PatchProcessor): + """Replace every sample with one value.""" + + fill_value: Any + + @property + def needs_numpy(self) -> bool: + """ + Whether the fill value is one the standard does not promise. + + The standard names which python scalars a namespace must accept. + A numpy scalar is not one of them -- numpy fills with it and + keeps its dtype where a strict backend refuses it outright -- and + neither is an integer too large for a signed 64-bit, which numpy + widens to unsigned or to object and a backend overflows on. + """ + value = self.fill_value + if type(value) not in (int, float, bool, complex): + return True + return type(value) is int and not -(2**63) <= value < 2**63 + + def kernel(self, data, meta, out_meta): + """ + Return an array of one value, the shape the patch is. + + The data is here for its namespace and its device and nothing + else: what comes out depends on the shape and the value alone. + """ + xp = array_namespace(data) + return xp.full(meta.shape, self.fill_value, device=device(data)) + + def numpy_kernel(self, data, meta, out_meta): + """ + Fill with a value only numpy will take, then hand it back. + + Converted back to the data's own backend afterwards: which of + the two kernels ran is decided by the fill value, and that must + not decide what the patch is made of. + """ + if not is_numpy(data): + warn_numpy_fallback("full", backend_name(data)) + return asarray_like(np.full(meta.shape, self.fill_value), data) + + +register_implementation("full", Full) @patch_function() @@ -1028,16 +1182,42 @@ def demedian(patch, dim: str = "time"): >>> plt.show() # doctest: +SKIP >>> plt.close(fig) """ - axis = patch.get_axis(dim) - data = patch.data + return Demedian._call(patch, dim=dim) - # Compute median along axis, keep dims for broadcasting - med = np.nanmedian(data, axis=axis, keepdims=True) - new_data = data - med +class Demedian(PatchProcessor): + """Remove the median along a dimension.""" - # Return a new patch with updated data - return patch.new(data=new_data) + dim: str = "time" + + @property + def needs_numpy(self) -> bool: + """ + Always: the standard has no median which skips nulls. + + `nan_reduce` says so by not offering one. This is about the + kernel written here, not about the operation: a package which + registers a median of its own for its backend is chosen ahead of + this and is not held to it. + """ + return True + + def numpy_kernel(self, data, meta, out_meta): + """ + Return the data with the median of each slice taken out. + + By `np.nanmedian`; `needs_numpy` says why. The only kernel this + class has, so data from another backend make the trip to numpy + and back rather than the operation refusing them. + """ + if not is_numpy(data): + warn_numpy_fallback("demedian", backend_name(data)) + array = to_numpy(data) + median = np.nanmedian(array, axis=meta.get_axis(self.dim), keepdims=True) + return asarray_like(array - median, data) + + +register_implementation("demedian", Demedian) @patch_function() @@ -1091,7 +1271,7 @@ def demean(patch, dim: str = "time"): >>> plt.show() # doctest: +SKIP >>> plt.close(fig) """ - return Demean(dim=dim)._apply(patch) + return Demean._call(patch, dim=dim) class Demean(PatchProcessor): diff --git a/dascore/proc/coords.py b/dascore/proc/coords.py index 2005f7038..8fe5b0a44 100644 --- a/dascore/proc/coords.py +++ b/dascore/proc/coords.py @@ -237,7 +237,7 @@ def rename_coords(self: PatchType, **kwargs) -> PatchType: >>> pa2 = pa.rename_coords(distance='fragrance') >>> assert 'fragrance' in pa2.dims """ - return RenameCoords(**kwargs)._apply(self) + return RenameCoords._call(self, **kwargs) class RenameCoords(PatchProcessor): @@ -289,8 +289,26 @@ def update_coords(self: PatchType, **kwargs) -> PatchType: >>> pa2 = pa.update_coords(distance=new_dist) >>> assert np.allclose(pa2.coords.get_array('distance'), new_dist) """ - new_coord = self.coords.update(**kwargs) - return self.new(coords=new_coord, dims=new_coord.dims) + return UpdateCoords._call(self, **kwargs) + + +class UpdateCoords(PatchProcessor): + """ + Give a patch other coordinates. + + No kernel: which values a coordinate holds is not what the data are. + """ + + # The coordinates arrive under whatever names the caller used, so the + # fields cannot be known in advance; see `RenameCoords`. + model_config = ConfigDict(extra="allow", frozen=True) + + def derive_meta(self, meta): + """Return the coordinates with the given ones changed or added.""" + return meta.update(coords=meta.coords.update(**self._params())) + + +register_implementation("update_coords", UpdateCoords) @patch_function() @@ -771,7 +789,7 @@ def transpose(self: PatchType, *dims: str) -> PatchType: >>> # Set distance as the first dimension. >>> out = pa.transpose("distance", ...) """ - return Transpose(dims=tuple(dims))._apply(self) + return Transpose._call(self, dims=tuple(dims)) class Transpose(PatchProcessor): diff --git a/dascore/workflow/meta.py b/dascore/workflow/meta.py index 3e579fa41..382eb8580 100644 --- a/dascore/workflow/meta.py +++ b/dascore/workflow/meta.py @@ -8,9 +8,9 @@ told apart -- the metadata step never sees an array, and the kernel never sees a coordinate. -That separation is what makes an operation fusible: something which wants -to compile a chain of them can ask each one what it does to the metadata -without touching, or even holding, the data. +That separation is what lets a chain of operations be compiled into one +pass: something which wants to can ask each one what it does to the +metadata without touching, or even holding, the data. """ from __future__ import annotations diff --git a/dascore/workflow/processor.py b/dascore/workflow/processor.py index d0abf9aa8..d0a0bf1cf 100644 --- a/dascore/workflow/processor.py +++ b/dascore/workflow/processor.py @@ -49,7 +49,12 @@ from dascore.workflow.checks import attr_type, check_patch_attrs, check_patch_coords from dascore.workflow.meta import PatchMeta from dascore.workflow.serialize import digest -from dascore.workflow.task import _VERSION_KEY, Task, _resolve_default +from dascore.workflow.task import ( + _VERSION_KEY, + Task, + _resolve_default, + _take_ownership, +) # Stands in for the patch while a call is bound to a signature. The bind # only needs something to put in that slot; nothing ever looks at it. @@ -263,26 +268,80 @@ def plan_kernel(self, meta: PatchMeta, out_meta: PatchMeta): difference between them -- `transpose` wants the permutation which takes the old dimension order to the new. - The default finds a kernel registered for the data's backend, and - failing that the class's own `kernel`, which is written to the - array API standard and so runs on any of them. A class with no - kernel at all is a metadata-only operation and gets None. + A kernel registered for the data's backend wins, being someone + saying they took this operation on there, arguments and all. + Failing that, the class's own `kernel`, written to the array API + standard and so able to run on any backend -- unless + `needs_numpy` says these arguments are outside what the standard + promises, in which case the class's `numpy_kernel` answers for + them instead. A class with none of the three is a metadata-only + operation and gets None. + + Which kernel runs is settled here rather than inside a kernel so + that something reading a chain of operations can see what each + one got without running any of them. """ - if (found := _resolve_kernel(type(self), meta.backend)) is None: + fallback = self.needs_numpy + if (found := _resolve_kernel(type(self), meta.backend, fallback)) is None: return None return functools.partial(found, self, meta=meta, out_meta=out_meta) + @property + def needs_numpy(self) -> bool: + """ + Whether these arguments are outside what the standard promises. + + Some operations are portable for only part of what they accept: + the standard names which python scalars a namespace must take, + and `full` given a numpy scalar is asking for something outside + that. A class which says yes here supplies a `numpy_kernel` for + those arguments; the default is no, since most operations have + only the one kernel. + + Answered from the operation's own parameters and never from the + data, so the choice is made before anything is read. + + This says nothing about whether the operation can be fused -- + that is a property of the kernel which ends up running, not of + the operation, and it is not this class's to answer. A package + which registers its own kernel for a backend may well lower + these same arguments happily, and its kernel is chosen ahead of + the numpy one precisely so that it can. + """ + return False + def reconcile(self, data, meta: PatchMeta) -> PatchMeta: """ Return the metadata the data actually turned out to have. - Defining this says the operation cannot be fused: it is the one - step which has to see both halves at once. The default only - carries the data's dtype back, since a kernel may promote. + The one step which sees both halves at once, which is why an + operation which needs it cannot be described by metadata alone. + The default only carries the data's dtype back, since a kernel + may promote. """ dtype = getattr(data, "dtype", meta.dtype) return meta if dtype == meta.dtype else meta.update(dtype=dtype) + @classmethod + def _call(cls, patch: PatchType, /, **kwargs) -> PatchType: + """ + Build the operation from a patch function's arguments and run it. + + This is what a patch function's body calls. Built here rather + than by the body so that the arrays among the arguments are not + taken over: a task freezes what it is handed so its fingerprint + cannot come to describe values it no longer holds, but an + operation built inside a patch function is run once and thrown + away, and freezing would reach back and lock the caller's own + array for the rest of its life. + """ + token = _take_ownership.set(False) + try: + operation = cls(**kwargs) + finally: + _take_ownership.reset(token) + return operation._apply(patch) + def _apply(self, patch: PatchType) -> PatchType: """ Run the operation, metadata first and then the data. @@ -471,21 +530,28 @@ def decorate(func): return decorate -def _resolve_kernel(cls: type[PatchProcessor], backend: str): +def _resolve_kernel(cls: type[PatchProcessor], backend: str, fallback: bool = False): """ Return the kernel a class runs for one backend, or None if it has none. A kernel registered for the backend wins; failing that the class's own `kernel`, which is written to the array API standard and so runs on any of them. A class which defines neither is metadata-only. + + `fallback` says the arguments are outside what the standard promises, + so the class's `numpy_kernel` stands in for the generic one. A + registered kernel still wins over it: whoever registered it took this + backend on and gets to say what it does with these arguments. """ - # One class at a time, both questions asked of it before moving up: + # One class at a time, every question asked of it before moving up: # a subclass which wrote its own `kernel` means it, and a backend # kernel registered against its parent must not answer for it. for klass in cls.__mro__: contents = klass.__dict__ if (found := contents.get("_kernels", {}).get(backend)) is not None: return found + if fallback and (numpy_kernel := contents.get("numpy_kernel")) is not None: + return numpy_kernel if (generic := contents.get("kernel")) is not None: return generic return None diff --git a/dascore/workflow/task.py b/dascore/workflow/task.py index e312e7559..a034c76fb 100644 --- a/dascore/workflow/task.py +++ b/dascore/workflow/task.py @@ -31,6 +31,7 @@ import warnings import weakref from collections.abc import Callable, Mapping +from contextvars import ContextVar from functools import cached_property from pathlib import Path from typing import Any, ClassVar, Self @@ -65,6 +66,15 @@ _PARAMS_KEY = "params" +# Whether a task being built takes over the arrays it is handed. A patch +# function builds an operation, runs it and throws it away, so there is no +# fingerprint left to go stale and no reason to lock the caller's buffer for +# the rest of its life; `PatchProcessor._call` turns this off for that one +# case. A ContextVar rather than a module flag so two threads building tasks +# at once cannot see each other's answer. +_take_ownership: ContextVar[bool] = ContextVar("_take_ownership", default=True) + + class Task(DascoreBaseModel): """ Base class for a fingerprintable, serializable operation. @@ -90,7 +100,9 @@ class Task(DascoreBaseModel): @classmethod def _own_the_arrays(cls, data: Any) -> Any: """Take ownership of any array a task was handed.""" - return own_arrays(data) if isinstance(data, dict) else data + if not isinstance(data, dict) or not _take_ownership.get(): + return data + return own_arrays(data) # Bumped by a subclass whenever the same parameters should mean a # different answer, so that old fingerprints do not name the new diff --git a/docs/contributing/extending_dascore.qmd b/docs/contributing/extending_dascore.qmd index c9f3dc9a5..98e225c0f 100644 --- a/docs/contributing/extending_dascore.qmd +++ b/docs/contributing/extending_dascore.qmd @@ -194,6 +194,10 @@ def _abs_mybackend(processor, data, meta, out_meta): The kernel is chosen by the backend the data report, so nothing else about the operation changes — `patch.abs()` reaches it, and the name, arguments and fingerprint are what they always were. An operation with no registered kernel for a backend falls back to the class's own `kernel`, which is written to the standard. +Some operations are portable for only some of what they accept — `patch.full(1.5)` is something every backend can do, `patch.full(np.float32(1))` is not. Those declare it through `needs_numpy`, a property answered from the operation's arguments and never from the data, and supply a `numpy_kernel` for those arguments, which converts to NumPy and back and issues a `NumpyFallbackWarning`. Which of the two runs is settled before any data is read. + +A kernel you register for your backend wins over the NumPy one, including for the arguments `needs_numpy` names. That is deliberate: registering a kernel says you took that operation on, arguments and all, and DASCore has no way to know what your backend can express. `needs_numpy` describes the kernel written here, never the operation — an operation DASCore has to run through NumPy may be one your backend lowers happily. + Only a few operations have a processor today; the rest are plain functions, and adding one is a deliberate act rather than something every patch function needs. ## Related guides diff --git a/scripts/differential_check.py b/scripts/differential_check.py index 4bb71854d..71afab4ae 100644 --- a/scripts/differential_check.py +++ b/scripts/differential_check.py @@ -130,7 +130,9 @@ def make_arrays() -> dict: "norm_bit": lambda patch: patch.normalize("time", norm="bit"), "norm_l2_distance": lambda patch: patch.normalize("distance", norm="l2"), "demean_distance": lambda patch: patch.demean("distance"), + "demedian_distance": lambda patch: patch.demedian("distance"), "rename": lambda patch: patch.rename_coords(time="t"), + "full_bool": lambda patch: patch.full(True), "transpose_noop": lambda patch: patch.transpose(*patch.dims), "transpose_ell": lambda patch: patch.transpose(..., "distance"), "norm_l1": lambda patch: patch.normalize("time", norm="l1"), @@ -301,6 +303,47 @@ def get_calls() -> dict: "transpose_ell_last": lambda: patch.transpose(..., "distance"), "transpose_ell_first": lambda: patch.transpose("distance", ...), "rename_coords": lambda: patch.rename_coords(distance="depth"), + # The five operations lowered here -- flip, fillna, full, demedian + # and update_coords -- with the branches only they reach. + "flip_noop": lambda: patch.flip(), + "flip_complex": lambda: dft_patch.flip("ft_time"), + "flip_bad_dim": lambda: patch.flip("nope"), + # The same name refused with the coordinates left alone: the two + # take different routes through the class and must still agree. + "flip_bad_dim_no_coords": lambda: patch.flip("nope", flip_coords=False), + "fillna_nothing_to_do": lambda: patch.fillna(0), + "fillna_complex": lambda: dft_patch.fillna(0), + "fillna_null_noinf": lambda: null_patch.fillna(-1, include_inf=False), + # A value `np.ndim` refuses, on a patch with nothing to fill: a + # no-op, and not the error measuring the value would raise. + "fillna_ragged_noop": lambda: patch.fillna([1, [2, 3]]), + "full_complex": lambda: dft_patch.full(1 + 1j), + "full_bool": lambda: patch.full(True), + "full_on_int": lambda: int_patch.full(3), + # Values which are not plain python scalars: numpy takes these + # and keeps their dtype, and the standard refuses them. + "full_np_scalar": lambda: patch.full(np.float32(1)), + "full_np_int8": lambda: patch.full(np.int8(3)), + # No entry for an integer too large for a dtype: numpy answers + # with an object array, and hashing one hashes the pointers in + # it, so the same call never agrees with itself. `needs_numpy` is + # what pins that case, in tests/test_workflow/test_processor_seam.py. + # A value with a shape is spent on the nulls one element each, + # which is not what broadcasting it would do. + "fillna_array_value": lambda: null_patch.fillna( + np.arange(int((~np.isfinite(np.asarray(null_patch.data))).sum())).astype( + "float64" + ) + ), + "demedian_distance": lambda: patch.demedian("distance"), + "demedian_null": lambda: null_patch.demedian("time"), + "demedian_int": lambda: int_patch.demedian("time"), + "update_coords_new": lambda: patch.update_coords( + quality=("distance", np.arange(patch.shape[0], dtype="float64")) + ), + "update_coords_replace": lambda: patch.update_coords( + distance=patch.get_array("distance") * 2 + ), "rename_nondim": lambda: with_nondim.rename_coords(quality="grade"), "transpose_named": lambda: patch.transpose("time", "distance"), "squeeze": lambda: patch.select(distance=0, samples=True).squeeze(), diff --git a/tests/test_proc/test_basic.py b/tests/test_proc/test_basic.py index 2014eea61..d56768ed4 100644 --- a/tests/test_proc/test_basic.py +++ b/tests/test_proc/test_basic.py @@ -465,6 +465,16 @@ def test_inf_not_dropped(self, patch_with_inf): class TestFillNa: """Tests for replacing nullish values in a patch.""" + def test_a_boolean_patch_has_nothing_to_fill(self, random_patch): + """ + A boolean holds nothing which is not a value. + + Answered from the dtype rather than asked of the data, since + `isfinite` is undefined for a boolean array on a strict backend. + """ + patch = random_patch.new(data=np.asarray(random_patch.data) > 0) + assert patch.fillna(True) is patch + def test_fillna(self, patch_with_null): """Ensure we can fillna and keep the other values the same.""" patch = patch_with_null.fillna(0) diff --git a/tests/test_utils/test_array_api.py b/tests/test_utils/test_array_api.py index d2d0ae391..ddc7b1daa 100644 --- a/tests/test_utils/test_array_api.py +++ b/tests/test_utils/test_array_api.py @@ -20,6 +20,7 @@ to_numpy, ) from dascore.utils.misc import suppress_warnings +from dascore.warnings import NumpyFallbackWarning @pytest.fixture(scope="module") @@ -205,6 +206,13 @@ def _identity(patch): return patch +def _with_a_null(patch): + """Return the patch with a null in it, so filling one does something.""" + data = np.asarray(patch.data).copy() + data[0, 0] = np.nan + return patch.new(data=data) + + def _make_complex(patch): """Return the patch with complex data, so a conjugate means something.""" data = np.asarray(patch.data) @@ -228,6 +236,16 @@ class _Case(NamedTuple): call=lambda patch: patch.rename_coords(time="t") ), "dascore.proc.basic.abs": _Case(call=lambda patch: patch.abs()), + "dascore.proc.basic.flip": _Case(call=lambda patch: patch.flip("time")), + "dascore.proc.basic.full": _Case(call=lambda patch: patch.full(1.5)), + # fillna is given something to fill; with nothing null it hands the + # patch straight back, which says nothing about the backend. + "dascore.proc.basic.fillna": _Case( + call=lambda patch: patch.fillna(0.0), setup=_with_a_null + ), + "dascore.proc.coords.update_coords": _Case( + call=lambda patch: patch.update_coords(time=patch.get_array("time")) + ), # conj and real hand a real patch straight back, which says nothing # about the backend, so they are given something to actually do. "dascore.proc.basic.conj": _Case( @@ -286,6 +304,51 @@ def test_backend_preserved(self, name, random_patch, to_backend, backend): assert np.allclose(array, np.asarray(expected.data), equal_nan=True) +# Operations which are portable for only part of what they accept, with an +# argument which takes them off the portable path. Unlike ARRAY_API_CASES +# these are expected to convert to numpy and back, and to say so. +NUMPY_FALLBACK_CASES = { + "full_numpy_scalar": _Case(call=lambda patch: patch.full(np.int8(3))), + "fillna_pandas_nulls": _Case( + call=lambda patch: patch.fillna(0.0, include_inf=False), setup=_with_a_null + ), + "demedian": _Case(call=lambda patch: patch.demedian("time")), +} + + +class TestTheNumpyFallbacks: + """ + What an operation does with the half of its arguments it cannot lower. + + The data make the trip to numpy and back rather than the operation + refusing another backend, and the trip is announced: a caller who + handed over a dask array has just had the whole of it materialised. + """ + + @pytest.mark.parametrize("name", sorted(NUMPY_FALLBACK_CASES)) + def test_it_warns_and_stays_on_the_backend( + self, name, random_patch, to_backend, backend + ): + """The patch comes back as it went in, and the detour is announced.""" + case = NUMPY_FALLBACK_CASES[name] + numpy_patch = case.setup(random_patch) + patch = to_backend(numpy_patch) + with pytest.warns(NumpyFallbackWarning): + out = case.call(patch) + assert backend_name(out.data) == backend + expected = case.call(numpy_patch) + array = np.asarray(out.data) + assert array.dtype == expected.data.dtype + assert np.allclose(array, np.asarray(expected.data), equal_nan=True) + + @pytest.mark.parametrize("name", sorted(NUMPY_FALLBACK_CASES)) + def test_numpy_data_is_not_a_fallback(self, name, random_patch): + """Nothing was converted, so nothing is said about converting.""" + case = NUMPY_FALLBACK_CASES[name] + with warnings_as_errors(): + case.call(case.setup(random_patch)) + + class TestNanReduce: """Tests for reductions which ignore nan values.""" diff --git a/tests/test_workflow/test_patch_op.py b/tests/test_workflow/test_patch_op.py index cd7bd93e5..b37f6ab1d 100644 --- a/tests/test_workflow/test_patch_op.py +++ b/tests/test_workflow/test_patch_op.py @@ -353,8 +353,8 @@ def test_a_name_with_an_implementation(self): def test_a_star_args_group_by_hand(self): """Including one whose arguments cannot be passed by name.""" assert PatchOp( - name="flip", kwargs={"dims": ("time",), "flip_coords": True} - ) == dc.proc.flip.op("time") + name="sort_coords", kwargs={"coords": ("time",), "reverse": False} + ) == dc.proc.sort_coords.op("time") class TestAPositionalBeforeAStarArgs: @@ -616,6 +616,11 @@ def test_the_registry_gains_one_tag(self): "RenameCoords", "Standardize", "Transpose", + "Demedian", + "FillNa", + "Flip", + "Full", + "UpdateCoords", } def test_the_document_names_the_operation(self): diff --git a/tests/test_workflow/test_processor_seam.py b/tests/test_workflow/test_processor_seam.py index 274e65d4e..4a3ed31bd 100644 --- a/tests/test_workflow/test_processor_seam.py +++ b/tests/test_workflow/test_processor_seam.py @@ -10,13 +10,21 @@ from __future__ import annotations import pickle +from typing import Any import numpy as np import pytest import dascore as dc from dascore.exceptions import CoordDataError, ParameterError -from dascore.proc.basic import Abs, Normalize, _known_real +from dascore.proc.basic import ( + Abs, + Demedian, + FillNa, + Full, + Normalize, + _known_real, +) from dascore.workflow import PatchMeta, PatchProcessor, Task, register_kernel from dascore.workflow.processor import ( _resolve_kernel, @@ -168,6 +176,167 @@ class Nothing(PatchProcessor): assert Nothing()._apply(patch) is patch +class TestWhichKernelIsPlanned: + """ + Which of a class's kernels a call gets, and why. + + Settled from the operation's parameters before any array is read, so + a chain can be inspected without being run. + """ + + def test_one_kernel_means_no_question(self): + """Most operations are portable for everything they accept.""" + assert not Abs().needs_numpy + assert not Normalize(dim="time", norm="l2").needs_numpy + + def test_a_kernel_which_is_numpy_says_so(self): + """The standard has no median which skips nulls.""" + assert Demedian().needs_numpy + + def test_it_can_depend_on_the_arguments(self): + """ + Some operations are portable for only some of what they accept. + + `full` takes any value numpy would; the standard promises only + the plain python scalars, and only those which fit a dtype. + `fillna` given a value with a shape spends it positionally, + which `where` cannot say, and `include_inf=False` asks pandas + what counts as nothing, which no backend answers. + """ + assert not Full(fill_value=1.5).needs_numpy + assert Full(fill_value=np.float64(1.5)).needs_numpy + assert Full(fill_value=2**70).needs_numpy + assert not FillNa(value=0).needs_numpy + assert FillNa(value=[1, 2]).needs_numpy + assert FillNa(value=0, include_inf=False).needs_numpy + + def test_a_value_numpy_cannot_measure(self): + """ + A ragged value is answered for rather than raised on. + + `np.ndim` refuses it, and a property which raised would turn a + patch with nothing to fill from a no-op into an error. + """ + assert FillNa(value=[1, [2, 3]]).needs_numpy + + def test_the_answer_needs_no_data(self): + """ + Reached with no array anywhere, which is the whole point. + + A property which read the data would raise here rather than + answer, since the operation is never given a patch at all. + """ + assert not Full(fill_value=1.5).needs_numpy + assert Demedian(dim="time").needs_numpy + + def test_a_registered_kernel_is_not_held_to_it(self, patch): + """ + Someone else's backend may express what ours cannot. + + `Demedian` says `needs_numpy` because the median written here is + numpy's. A package which registers a median for its own backend + is answering a different question, and is chosen ahead of the + numpy one so that it can -- otherwise registering a kernel for + the operations DASCore finds hardest would buy nothing. + """ + + class Middling(Demedian): + """A `demedian` whose backend someone else claimed.""" + + @register_kernel(Middling, "numpy") + def _theirs(self, data, meta, out_meta): + """Answer with something no other kernel would.""" + return np.zeros(meta.shape) + + operation = Middling(dim="time") + meta = PatchMeta.from_patch(patch) + assert operation.needs_numpy + assert operation.plan_kernel(meta, meta).func is _theirs + + +class TestTheNumpyFallbacks: + """ + What the two half-portable operations do with the other half. + + The parity check covers these, but it is not what the coverage gate + runs, and an untested fallback is how a rewrite quietly narrows what + an operation accepts. + """ + + def test_a_value_with_a_shape_is_spent_positionally(self): + """One element per null, in order -- not broadcast.""" + patch = dc.get_example_patch("patch_with_null") + data = np.asarray(patch.data) + nulls = ~np.isfinite(data) + values = np.arange(int(nulls.sum()), dtype="float64") + expected = data.copy() + expected[nulls] = values + assert np.array_equal(np.asarray(patch.fillna(values).data), expected) + + def test_nothing_to_fill_hands_the_patch_back(self, patch): + """ + Whichever of the two fills was planned, an empty mask is a no-op. + + The identity is what the decorator reads as nothing having + happened, so no history is written and no id advances. + """ + assert patch.fillna(np.arange(3.0)) is patch + assert patch.fillna(0) is patch + + @pytest.mark.parametrize("value", [np.int8(3), 2**70]) + def test_a_value_the_standard_will_not_take_is_planned_onto_numpy( + self, patch, value + ): + """ + The plan says numpy, and it says so before any data is read. + + Asserted on the plan rather than on the dtype: on numpy the + portable fill answers a numpy scalar identically, so a result + cannot tell which kernel produced it. + """ + meta = PatchMeta.from_patch(patch) + planned = Full(fill_value=value).plan_kernel(meta, meta) + assert planned.func is Full.numpy_kernel + # And the dtype numpy keeps for it is what comes out. + assert patch.full(value).data.dtype == np.full((1,), value).dtype + + +class TestTheCallersArguments: + """ + What a patch function does to the arguments it was handed. + + A task freezes the arrays it is given so its fingerprint cannot come + to describe values it no longer holds. An operation built inside a + patch function is run once and thrown away, so there is no such + fingerprint -- and freezing would reach back and lock a buffer the + caller means to keep writing to. + """ + + def test_a_fill_array_stays_writable(self, patch): + """The values are read, not taken over.""" + values = np.arange(3.0) + patch.fillna(values) + assert values.flags.writeable + + def test_a_coordinate_array_stays_writable(self, patch): + """`update_coords` copies what it is given, and always has.""" + values = np.arange(patch.shape[patch.dims.index("time")], dtype="float64") + patch.update_coords(time=values) + assert values.flags.writeable + + def test_a_task_still_takes_ownership(self): + """The policy is off for the one case, not repealed.""" + values = np.arange(3.0) + + class Holding(Task): + """A task which holds whatever it is handed.""" + + value: Any = None + + Holding(value=values) + assert not values.flags.writeable + + class TestRegistrationRefuses: """What a class is turned away for, at import rather than at use."""