Honor the input mask in cosmicray_median - #979
Conversation
The `if hasattr(ccd, "mask"): data = ccd.data` branch was dead code, immediately overwritten by `xp.asarray(ccd)`, which silently drops the mask of a numpy masked array. Masked pixels are now never flagged as cosmic rays and are replaced by the local median for the purposes of detection and replacement, so they do not contaminate the filters. The existing masked test only passed because the mask was ignored; it now masks a subset of the injected cosmic rays and checks they are not flagged. Adds a nomask case and a CCDData-with-mask case. Fixes astropy#932 Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #979 +/- ##
==========================================
+ Coverage 95.90% 96.01% +0.10%
==========================================
Files 8 8
Lines 1637 1655 +18
==========================================
+ Hits 1570 1589 +19
+ Misses 67 66 -1 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
pre-commit.ci has been failing on main since the black 26.5.1 bump in astropy#919 because black now requires a blank line after the module docstring in ccdproc/__init__.py. This is the only formatting change it wants. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
The first implementation replaced masked pixels with the local median computed from the raw data, so a masked region wider than the median box still leaked its values into the median of neighbouring pixels. Fill the masked pixels with the mean of the unmasked data first, compute the local median, then refine once with that median, so the detection median never depends on the masked values. Add a regression test with a wide masked bright region and a marginal cosmic ray next to it, which the one-pass version misses. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
|
Pushed 8ac4a23: masked pixels are now excluded from the median filter exactly, rather than approximately. The first version replaced masked pixels with the local median of the raw data, which still leaked masked values into the median of neighbouring pixels whenever a masked region was wider than the median box. Now the masked pixels are first filled with the mean of the unmasked data (no dependence on the masked values), the local median is computed, and then refined once with that median. Cost: one extra New regression test |
There was a problem hiding this comment.
Pull request overview
Updates cosmicray_median to preserve and apply input masking semantics (especially for numpy.ma.MaskedArray inputs), ensuring masked pixels neither influence detection nor get flagged/replaced, and extends tests/documentation accordingly.
Changes:
- Implement explicit mask extraction/handling in
cosmicray_median’s array branch so masked pixels are excluded from detection and median filtering effects. - Expand/strengthen regression tests for masked-array and CCDData-with-mask behavior (including unioning masks).
- Document the behavior change in the docstring and release notes.
Reviewed changes
Copilot reviewed 4 out of 4 changed files in this pull request and generated 1 comment.
| File | Description |
|---|---|
CHANGES.rst |
Adds a bug-fix entry describing the corrected masked-input handling in cosmicray_median. |
ccdproc/core.py |
Implements real mask semantics in the array branch of cosmicray_median (preserve mask, exclude masked pixels from detection/filtering). |
ccdproc/tests/test_cosmicray.py |
Strengthens regression coverage for masked inputs, including partial masking, nomask, and CCDData mask union behavior. |
ccdproc/__init__.py |
Minor formatting change (blank line). |
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
Address review on astropy#979: - Forward the CCDData mask to the detection instead of only OR-ing it into the output mask, by factoring the array implementation into a helper shared by the array and CCDData branches. - Return early when every pixel is masked so the fill value is not a 0/0 that emits a RuntimeWarning. - Use American spelling (honor, neighbor) and make the docstring lead with the CCDData case. - Add tests for the CCDData neighbor-bias case and the all-masked case. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
The scipy.ndimage boundary escape moved from cosmicray_median into the new helper, which the baseline ratchet keys by function name. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
There was a problem hiding this comment.
Pull request overview
Copilot reviewed 5 out of 5 changed files in this pull request and generated no new comments.
Suppressed comments (2)
Previously missed (2) — in code that hasn't changed since the last review.
ccdproc/core.py:2147
- This two-pass fill does not actually exclude masked samples from the median. If masked pixels are a majority of a local window, the first median remains the global
fill, and the second pass preserves that value; an unmasked source beside or inside such a region can consequently be compared against the global mean rather than the median of the available local pixels. Please use a genuinely mask-aware rolling median (with an explicit fallback for windows containing no valid pixels), or narrow the documented semantics if this approximation is intentional.
fill = xp.sum(xp.where(keep, data, 0)) / n_keep
filt_data = xp.where(in_mask, fill, data)
marr = xp.asarray(ndimage.median_filter(filt_data, size=(mbox, mbox)))
filt_data = xp.where(in_mask, marr, data)
marr = xp.asarray(ndimage.median_filter(filt_data, size=(mbox, mbox)))
ccdproc/core.py:2122
- When
error_imageis omitted for a masked-array input, this still computes the standard deviation from every raw value, including masked pixels. A bright masked defect can therefore inflate the detection threshold enough to hide unmasked cosmic rays, contradicting the new guarantee that masked pixels do not affect detection. Compute the standard deviation from~in_maskonly (while retaining the all-masked early return) and add a regression case that omitserror_image.
if error_image is None:
error_image = xp.std(data)
mwcraig
left a comment
There was a problem hiding this comment.
Critical review of the diff (algorithm-semantics, numpy-regression, multi-backend, and test-quality passes, every number below re-run locally on the PR head 49cf28e).
What holds. No regression for unmasked inputs: 209 input combinations (dtypes × thresh/mbox/gbox/rbox × error_image forms, CCDData with/without uncertainty, memmap, non-contiguous, nomask, Quantity) are bit-identical to main except four intentional mask-honouring differences and a harmless C→K memory-order change. The headline fix for #932 works: masked bright values (bad column at 65535, bright blocks, bright edge rows) no longer leak into the median on a flat background (≤ 0.4σ from an exact masked median, zero mis-flags), masked pixels are never flagged and come back unchanged even with rbox>0, all-masked input is clean, NaN under the mask is harmless, jax/dask/numpy give identical CR counts, and the dask escape ratchet is green.
What doesn't, inline. Three real algorithm problems: (1) the two-pass fill is biased toward the global mean next to masked regions wider than ~mbox/2 on a non-flat background — measured false positives at ~15σ and a missed 6σ CR; (2) masked pixels seed gbox growth into clean pixels because the mask is applied after maximum_filter; (3) error_image=None still takes the std over masked pixels, so one masked saturated column turns detection off entirely. Then test coverage: 6 of 7 mutants of the new function survive the whole file, and none of the strict xfails fails for the stated reason (two already XPASS). Plus device/dask notes and a docstring that over-promises.
Housekeeping: the branch is based on a0b6610 and needs a rebase (the strict test_background_deviation_box failure is #977's, already fixed on main).
— Written by Claude at @mwcraig's direction.
Drop the approximate masked median (global-mean fill plus a second median pass) introduced earlier on this branch: it was biased toward the global mean next to wide masked regions on non-flat backgrounds and was beyond the scope of astropy#932. The median filter now runs on the data as-is, exactly as before. Masked pixels are still never flagged as cosmic rays -- the mask is now applied before the gbox growth step as well as after it, so masked pixels can neither seed growth nor be grown into -- are returned unchanged, and are excluded from the noise estimate when error_image is None so that a masked saturated column cannot inflate the threshold. Tests replaced accordingly; the two XPASSing array-api-strict markers are removed. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
|
Scoped back down in dc27fbb, per @mwcraig: this PR was meant to be the array-API fix for #932, and the approximate masked median had turned it into an algorithm change. The global-mean fill and second median pass are gone; the median filter runs on the data exactly as on — Written by Claude at @mwcraig's direction. |
…change - Remove the carried-over "scipy.ndimage ignores the mask" comment from _cosmicray_median_array, which now handles the mask. - Raise ValueError when the mask shape does not match the data shape. - Mention in CHANGES that masked pixels no longer seed gbox growth. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
…bool) Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA
cosmicray_medianhad dead code in its array branch:if hasattr(ccd, "mask"): data = ccd.datawas immediately overwritten bydata = xp.asarray(ccd), which silently drops the mask of anumpy.ma.MaskedArray(backend-dependently), so an input mask was ignored.This PR makes the mask handling deliberate and backend-independent, with intentionally narrow semantics:
numpy.ma.nomask), and theCCDDatabranch passesccd.maskstraight through to a shared_cosmicray_median_arrayhelper instead of round-tripping through the array branch;gboxgrowth step (so masked pixels cannot seed growth) and after it (so growth cannot extend into them);rbox > 0;error_imageisNone, the noise is estimated from the unmasked pixels only, so a masked saturated column cannot inflate the threshold and switch detection off;CCDDataoutput mask is built with|rather than+.The median filter itself still runs on the data as-is:
scipy.ndimageknows nothing about masks, so bright masked pixels can still influence the local median of their neighbors. That is documented in the Notes, and making the median mask-aware is #984 (an earlier revision of this branch attempted an approximate version and it was not good enough).Tests:
test_cosmicray_median_maskednow masks a subset of the injected cosmic rays and asserts those are not flagged while the rest are (fails onmain); new tests covergbox/rboxwith a mask (masked array andCCDData),error_image=Nonewith a masked saturated column, all-masked input,nomask, and theCCDDataunion mask. Removing the mask guard or the masked noise estimate is caught by the tests.Test matrix (branch merged with current
main): numpy 392 passed; dask 381 passed, escape ratchet clean; jax (x64)test_cosmicray.py46 passed / 2 xfailed; array-api-stricttest_cosmicray.py12 passed / 36 xfailed, no XPASS.Fixes #932
🤖 Generated with Claude Code
https://claude.ai/code/session_01DTN9DnPnLKK2u7knnMJ2gA