Array-API-native k×k window filters for cosmicray_median, median_filter, background_deviation_filter and ccdmask

[mwcraig]

Depends on #1009 (rebased on top of it; merge that first).

Part of #971, section 3. The masked-median fix for #984 is a planned follow-up (PR-2) that builds on the NaN handling added here.

What changed and why

scipy.ndimage's window filters are numpy-only, so the four call sites that used them silently copied a jax or dask array to the host, and failed outright for an array on a non-default device. This adds ccdproc/_windowfilters.py, an implementation of the four window reductions ccdproc.core needs, written purely in terms of the array API standard:

helper	replaces
window_rank	ndimage.percentile_filter
window_median	ndimage.median_filter
window_any	ndimage.maximum_filter (on a boolean mask)
window_reduce	ndimage.generic_filter

The approach is deliberately literal: the array is padded once — the standard has no pad, so reflect and nearest are built from flip/concat of edge slices — and each window offset is taken as a shifted slice, with the offsets stacked into a trailing axis, which turns a window reduction into an ordinary reduction over that axis. The rank filters reduce with a new _nanfuncs._nanrank, factored out of nanmedian's sentinel-sort machinery; that is why they reproduce ndimage's rank exactly, including the upper-middle element of an even window.

The stack holds prod(size) copies of the input — 3.8 GiB for an 11x11 window over a 2048x2048 float64 image — so _windowed processes the output in bands of rows sized to a 256 MiB module-constant budget. Bands are cut from the padded array with the window overhang included, so each output pixel sees exactly the window it would have seen unbanded; band_rows stays on the private helpers (a test pins banded == unbanded) and is not exposed on any public function.

Dispatch policy: the numpy path is unchanged

Four wrappers in core.py — _dispatch_median_filter, _dispatch_percentile_filter, _dispatch_maximum_filter, _dispatch_generic_filter — send numpy input to scipy.ndimage exactly as before and everything else to _windowfilters. Concentrating the choice in one place per filter also leaves the seam #984 needs, where numpy input carrying a mask will start taking the native path too.

Call sites converted: _cosmicray_median_array's mbox median, gbox growth and rbox replacement; the public median_filter; background_deviation_filter (which gains an xp argument, matching background_deviation_box); and ccdmask's median and its two percentile filters. ccdmask's byblocks branch is untouched.

median_filter on a non-numpy array now accepts only size and mode. The argument list is bound against ndimage's own signature rather than duplicated, so a repeated or unknown argument still gets ndimage's error message, and footprint/origin/output/cval/axes raise TypeError naming the argument. The numpy path stays a verbatim passthrough.

cosmicray_median's CCDData branch now assigns through nccd._mask, as ccd_process already does: NDDataArray's mask setter runs the value through np.asarray, which was the real source of the "cosmicray_median numpy.asarray" escape the baseline attributed to ndimage.

Documented divergences from ndimage

• Integer input is promoted to the namespace's default real floating dtype; ndimage keeps an integer dtype.
• Only ndimage's reflect and nearest boundary modes exist; anything else raises rather than quietly behaving like reflect.
• A window needing more padding than its axis holds raises, where ndimage re-reflects.
• Cost is O(k2 log k2) per pixel for a k-by-k window, from a sort, against ndimage's O(k**2) selection.
• NaN handling. The rank filters exclude NaNs from a window and take the rank among the values that remain; ndimage sorts NaNs in with the values, above every real number. The one caller this reaches is ccdmask, whose input is a flat ratio that may well contain NaN — it opens by masking the non-finite pixels — so a ratio with NaN in it can give a slightly different mask off numpy. Every other caller filters finite data. The exclusion is deliberate: it is what will let cosmicray_median keep masked pixels out of its median in PR-2. test_windowfilters.py::test_ccdmask_window_filters_exclude_nan_off_numpy pins it in both directions, against an explicit sliding_window_view rank reference off numpy and against ndimage on numpy, asserting first that the two references really disagree.

All of these are in the docs/array_api.rst limitations list.

Verified backends

Dev environment: numpy 2.1.0, scipy 1.14.0, astropy 7.1.0, jax 0.5.0 (JAX_ENABLE_X64=1), dask 2025.7.0, array-api-strict 2.5 (on the non-default device1).

Affected files = test_nanfuncs.py test_windowfilters.py test_blocks.py test_ccdmask.py test_cosmicray.py test_ccdproc.py.

backend	affected files	full pytest ccdproc
numpy	591 passed	1203 passed, 5 skipped
jax	586 passed, 5 skipped	1196 passed, 36 skipped
dask	586 passed, 5 skipped	1197 passed, 16 skipped (run with docs, under the escape ratchet)
array-api-strict	560 passed, 31 xfailed	1166 passed, 10 skipped, 32 xfailed

No failures and no XPASS on any backend. Against the pre-PR branch point, array-api-strict loses 11 xfails, all of them tests that only failed because they went through ndimage.

The escape-baseline ratchet was run the way CI's py312-alldeps-dask-enforce job runs it (CCDPROC_ARRAY_LIBRARY=dask CCDPROC_LOG_ARRAY_ESCAPES=1 CCDPROC_ENFORCE_ESCAPE_BASELINE=1 pytest ccdproc docs). Before the baseline edit it reported exactly three entries as "not hit this run" — _cosmicray_median_array, background_deviation_filter and cosmicray_median; after deleting them it reports OK: no library escapes outside the baseline.

Deviations from the agreed plan

1. Banding is not skipped on dask. The plan said to skip it and let dask's chunking bound memory. That made the dask suite unusable: test_cosmicray_median_rbox took 202 s and test_background_deviation_filter 230 s. The cause is that slicing a chunked band into prod(size) differently-offset windows makes dask realign every one of them — a 21x21 window over a 100x100 image builds a graph of 585,149 tasks, against 912 when the band is first collapsed to a single chunk. So each band is now rechunk(-1)ed (guarded by is_dask_namespace), and banding stays on for dask, since it is the band budget that makes collapsing safe. The dask affected-file run went from 496 s to ~15 s.
2. dask's PerformanceWarning about the chunk-count multiplication is suppressed inside _stack_from_padded, matched by message so dask need not be imported. That multiplication is the algorithm and no caller can act on it, and ccdproc's pytest configuration turns warnings into errors, so leaving it would fail every dask test touching a window filter. Documented in the module docstring.
3. Test-side numpy-isms had to be fixed for the marker removals to mean anything. Dropping the ndimage backend_xfails exposed that add_cosmicrays built its rays by copying the image to numpy — impossible for an array on a non-default device — and that several assertions used .sum(), .mean(), .std() or np.array() on backend arrays. add_cosmicrays now blends a host-side overlay in with a single where; the assertions use xp.* or _to_numpy.
4. One extra core.py fix: nccd.mask | crarr coerced the namespace array to numpy once crarr stopped being a numpy array, so the existing mask is brought into the data's namespace first.
5. One extra marker pruned: test_cosmicray_lacosmic_detects_inconsistent_units began XPASSing on array-api-strict once add_cosmicrays worked (it raises before ever reaching astroscrappy), so its marker is gone, per the "prune on XPASS" note in docs/array_api.rst.
6. The helpers are N-D rather than 2-D-only. The padding was already axis-generic, so generality cost nothing.

Not verified

• The Sphinx docs build was not run locally. sphinx-astropy is not installed in the dev environment. The first CI run caught two nitpicky cross-references to the private _windowfilters module in public docstrings; those are now spelled as literals. docs/array_api.rst was checked with docutils (clean apart from Sphinx-only roles) and the .. _scipy.ndimage: link target the new bullets use was added.
• CuPy. Not installed, so untested, as usual for this repo.
• The XPASS-driven marker removals were confirmed only on macOS, not in CI — docs/array_api.rst asks for CI confirmation before deleting a marker, so that one is worth a second look on the CI logs.

🤖 Generated with Claude Code

https://claude.ai/code/session_01F6f9L1GyMHrKnvLrrgWbWN

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 98.20%. Comparing base (9d18599) to head (aa71121).
⚠️ Report is 5 commits behind head on main.

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #1010      +/-   ##
==========================================
+ Coverage   97.97%   98.20%   +0.23%     
==========================================
  Files           9       11       +2     
  Lines        1927     2177     +250     
==========================================
+ Hits         1888     2138     +250     
  Misses         39       39              

Flag	Coverage Δ
dask	97.19% <97.52%> (+<0.01%)	⬆️
jax	97.37% <97.52%> (+0.02%)	⬆️
numpy	97.38% <94.34%> (-0.50%)	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.