Wiener-Butterworth unmatched back projector for fast RL deconvolution#52
Merged
Conversation
…ectors Port the reference Guo et al. (2020) BackProjector.m to a device-agnostic create_backprojector(). It builds a spatial-domain back-projector PSF from a forward PSF, supporting the traditional, gaussian, butterworth, wiener, and wiener-butterworth (WB) types. The WB projector flattens the spectral product |FT(f)·FT(b)| across the passband so Richardson-Lucy converges in ~1-2 iterations instead of 10-50. Works on NumPy (CPU) or CuPy (GPU) arrays; FWHM is computed on tiny host line profiles, grids are created on the input device. Exported from cubic.preprocessing. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…RL deconvolution Add a `backprojector` argument to richardson_lucy_xp (and thread it through decon_xpy, richardson_lucy_iter, and deconv_iter_num_finder). When given, a faithful unmatched Richardson-Lucy update is used (ported from the reference DeconSingleView.m): both projectors normalized to sum=1, no H^T1 term, no correction<0 clip, estimate floored at 1e-3 each iteration. This lets an unmatched back projector (e.g. Wiener-Butterworth) reach resolution-limited quality in ~1-2 iterations. The default (backprojector=None) path is byte-for-byte the existing matched RL, verified identical on CPU and GPU. The unmatched path is circulant only: it raises NotImplementedError for noncirc/mask and ValueError for the skimage implementation or a back projector whose shape differs from the PSF. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Add a 3D example notebook (and its paired script) demonstrating the Wiener-Butterworth unmatched back projector on the astrocyte-nuclei stack: building the back projector, the flattened combined transfer function, and matched-vs-WB convergence tracked by PSNR/SSIM (vs a matched-converged reference) and absolute FSC resolution. WB reaches matched-converged quality in ~1-2 iterations vs ~10 for the matched projector. Linked from the README example table. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Construct the flipped-PSF OTF, beta_fp cutoff gain, and Wiener filter only in the branches that consume them, so traditional/gaussian/butterworth no longer allocate unused full-size complex arrays (matters for large 3D GPU volumes). Factor the shared Butterworth mask into a helper and collapse the cutoff-gain reduction to a single multi-axis max. Outputs are bit-identical for every back-projector type on CPU and GPU. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Collapse the duplicated device check into one call, and reject pad_psf=True with a back projector up front (it must stay PSF-shaped) with a clear error instead of a misleading downstream shape-mismatch. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Contributor
There was a problem hiding this comment.
Sorry @alxndrkalinin, your pull request is larger than the review limit of 150000 diff characters
The trailing-edge guard in _fwhm_1d checked `i == n - 1`, but the search loop exits with `i == n` when a profile never returns below half-maximum before the array boundary, so it fell through to an out-of-bounds `y[i]` (IndexError) instead of returning NaN. Mirror the leading-edge guard (`i >= n`) so an unresolved profile yields NaN, and raise a clear ValueError in create_backprojector when the FWHM cannot be estimated (PSF too wide for its array) rather than crashing downstream. Add regression tests. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
There was a problem hiding this comment.
Pull request overview
This PR adds support for unmatched back projectors (including the Wiener-Butterworth design) to accelerate Richardson–Lucy (RL) deconvolution in cubic while preserving the existing matched-path behavior when the feature is not used.
Changes:
- Introduces
create_backprojector(cubic.preprocessing.create_backprojector) to generate several back-projector PSFs (traditional/gaussian/butterworth/wiener/WB) in a device-agnostic way. - Threads an optional
backprojectorthrough the RL implementation (richardson_lucy_xp) and the higher-level deconvolution wrappers (decon_xpy,richardson_lucy_iter,deconv_iter_num_finder) with appropriate validation/guardrails. - Adds unit tests for back projector generation and unmatched RL behavior, plus a new example script and README entry.
Reviewed changes
Copilot reviewed 9 out of 10 changed files in this pull request and generated 3 comments.
Show a summary per file
| File | Description |
|---|---|
cubic/preprocessing/backprojector.py |
New back-projector PSF generator (including Wiener-Butterworth). |
cubic/preprocessing/richardson_lucy_xp.py |
Adds an unmatched RL branch driven by a user-supplied back projector. |
cubic/preprocessing/deconvolution.py |
Plumbs backprojector through deconvolution entry points with validation/errors. |
cubic/preprocessing/__init__.py |
Exports create_backprojector from the preprocessing package. |
tests/preprocessing/test_backprojector.py |
New tests covering back-projector types, normalization, and CPU/GPU parity. |
tests/preprocessing/test_deconvolution.py |
New tests for unmatched RL execution, parity, and expected error cases. |
examples/scripts/deconvolution_wb_backprojector_3d.py |
New example script demonstrating matched vs WB convergence behavior. |
README.md |
Adds a README entry linking the new WB deconvolution notebook. |
cubic/__init__.py |
Bumps package version to 0.8.0a2. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
…izing A PSF that sums to zero or contains NaNs would otherwise propagate NaNs silently through the FFTs and return an invalid back projector (the downstream FWHM check is skipped entirely for res_flag=2). Validate the sum up front and raise a clear ValueError instead. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…o xpy The signature defaults to implementation="xpy" but the docstring listed "skimage" as the default. Align the docstring with the actual default. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
Adds an unmatched Wiener-Butterworth (WB) back projector to Richardson-Lucy (RL) deconvolution, following Guo et al. 2020 (Nat. Biotechnol. 38:1337; reference code
eguomin/regDeconProject). Replacing the matched back projector (flipped PSF / conjugate OTF) with one that flattens the spectral product|FT(f)·FT(b)|across the passband makes RL converge in ~1–2 iterations instead of 10–50. The speedup is iteration reduction (each iteration still costs the same 4 FFTs).What's included
create_backprojector(cubic/preprocessing/backprojector.py) — device-agnostic (NumPy/CuPy) port ofBackProjector.m, supportingtraditional,gaussian,butterworth,wiener, andwiener-butterworthtypes. Exported fromcubic.preprocessing.backprojectortoggle onrichardson_lucy_xp, threaded throughdecon_xpy,richardson_lucy_iter, anddeconv_iter_num_finder. WhenNone(default) the existing matched path is byte-for-byte unchanged; when given, a faithful unmatched RL update is used (noH^T1term, nocorrection<0clip, estimate floored at 1e-3, both projectors sum-normalized — perDeconSingleView.m). The unmatched path is circulant only (raisesNotImplementedErrorfornoncirc/mask,ValueErrorfor the skimage implementation or a shape-mismatched back projector).examples/notebooks/deconvolution_wb_backprojector_3d.ipynb(+ paired script + README row) demonstrating matched-vs-WB convergence by PSNR/SSIM (vs a matched-converged reference) and absolute FSC resolution.0.8.0a2.Validation
traditional == flip/sum, FWHM ≈ 2.3548·σ, combined-transfer coefficient of variationCV_WB ≪ CV_matched, CPU/GPU parity ≤ 4e-8.alpha=0.05(default); WB degrades past ~3 iterations (noise amplification).ruff,ruff format,mypy, and the fullpytestsuite (397 passed, GPU active) all green.Notes
fix/deconv-data-zenodois behindmainand touches the same files; it would clobber this work if it lands — reconcile before merging that branch.🤖 Generated with Claude Code