Skip to content

v1.2.10

Choose a tag to compare

@github-actions github-actions released this 14 Apr 16:58
· 24 commits to main since this release

Highlights

This release is a correctness overhaul of every analytical reconstruction path
(parallel / fan / cone) and of the autograd adjoint path shared by the
iterative examples. It also ships a unified architecture where each geometry
has its own dedicated voxel-driven FBP/FDK gather kernel, separate from the
pure Siddon adjoint used by autograd.

Bug fixes

Analytical FBP / FDK amplitude bugs

Prior to this release, all three analytical reconstruction examples produced
results that were off by large constant factors:

Example Pre-fix reco range Post-fix reco range MSE improvement
fbp_parallel.py [0, 6.52] [−0.02, 1.01] ~580×
fbp_fan.py [0, 6.33] [−0.08, 1.01] ~900×
fdk_cone.py [0, 10.08] [−0.08, 1.00] ~600×

Root causes:

  • Cone and fan used (sdd/U)^2 with U = sdd + x·sin − y·cos instead of
    the correct FDK weight (sid/U)^2 with U = sid + x·sin − y·cos.
  • All three geometries were missing the 1/(2π) Fourier-convention
    constant from the FBP/FDK reconstruction formula.
  • Cone and fan used a Siddon ray-driven scatter for the analytical path,
    which is neither a true adjoint nor a classical FBP/FDK gather.

Autograd adjoint bug (cone + fan)

ConeProjectorFunction.backward, ConeBackprojectorFunction.forward,
FanProjectorFunction.backward and FanBackprojectorFunction.forward
were all passing distance_weight=1.0 to the shared Siddon backward
kernel. This meant the autograd backward was not the true adjoint P^T
of the forward projector — it had a (sdd/U)^2 per-voxel factor baked in,
biasing the gradient by ~2–3× depending on voxel position. Any iterative
reconstruction that relied on autograd (including iterative_reco_cone.py
and iterative_reco_fan.py) was running on a biased gradient flow.

Parallel beam was unaffected (its backward kernel had no distance_weight
parameter to begin with).

This release removes the distance_weight parameter from both
_cone_3d_backward_kernel and _fan_2d_backward_kernel entirely, so the
dead code path cannot be accidentally reintroduced, and fixes the four
autograd call sites to use the pure adjoint.

New features

Voxel-driven gather kernels for every geometry

Three new CUDA kernels, all under a dedicated fastmath=False decorator
for FDK-grade accuracy:

  • _parallel_2d_fbp_backproject_kernel — no distance weighting (no source).
  • _fan_2d_fbp_backproject_kernel — (sid/U)^2 + linear detector interp.
  • _cone_3d_fdk_backproject_kernel — (sid/U)^2 + bilinear detector interp.

Each kernel is voxel-driven: one thread per output pixel/voxel, loops over
views, computes the projected detector coordinate, interpolates the filtered
sinogram, weights and accumulates.

parallel_weighted_backproject (new public helper)

Mirrors fan_weighted_backproject and cone_weighted_backproject.
Applies the 1/(2π) Fourier-convention constant so a unit-density disk
reconstructs to amplitude 1. Exported from diffct.

ramp_filter_1d — new backward-compatible kwargs

  • sample_spacing (default 1.0): physical detector cell pitch. Output
    is rescaled by 1/sample_spacing for physical-unit correctness.
  • pad_factor (default 1): zero-pad to pad_factor * N before the
    FFT to suppress circular-convolution wrap-around. 2 is recommended
    for FBP/FDK.
  • window (default None): frequency-domain apodization. Options:
    None / "ram-lak", "hann", "hamming", "cosine",
    "shepp-logan".
  • use_rfft (default True): faster real-FFT path for real-valued
    inputs.

Existing ramp_filter_1d(sino, dim=1) calls keep working unchanged.

Analytical FBP / FDK scale factors

Each analytical helper now applies the correct analytical constant
automatically so reconstructions are already amplitude-calibrated:

  • parallel_weighted_backproject: 1 / (2π).
  • fan_weighted_backproject: sdd / (2π · sid).
  • cone_weighted_backproject: sdd / (2π · sid).

Testing

Test suite expanded from 6 to 27 tests:

  • test_adjoint_inner_product.py (new): the definitive check.
    For random x and y, asserts ⟨A x, y⟩ = ⟨x, A^T y⟩ for parallel,
    fan and cone autograd pairs. Permanently guards against any regression
    that would reintroduce the distance_weight=1.0 bug.
  • test_fdk_cone_accuracy.py (new): 128³ Shepp-Logan RMSE / amplitude
    bounds.
  • test_fdk_cone_offsets.py (new): detector, center and combined
    offsets with amplitude assertions.
  • test_cone_projector_autograd.py (new): gradient finiteness and
    non-zero sanity.
  • test_fbp_fan_accuracy.py (new): 256×256 Shepp-Logan RMSE.
  • test_fbp_fan_offsets.py (new): three offset configurations.
  • test_fbp_parallel_accuracy.py (new): parallel FBP RMSE + offset.

Examples (unified structure)

examples/fdk_cone.py, examples/fbp_fan.py and
examples/fbp_parallel.py have been rewritten with a consistent 8-step
layout and detailed inline comments documenting every geometry variable
(what it is, units, typical values, available options). Each example
prints raw MSE, clamped MSE, and the reconstruction / phantom data
ranges.

The iterative reconstruction examples (iterative_reco_parallel.py,
iterative_reco_fan.py, iterative_reco_cone.py) are unchanged — they
silently benefit from the fixed autograd adjoint.

Documentation

  • docs/source/api.rst: adds autofunction entries for
    parallel_weighted_backproject, documents the new ramp_filter_1d
    options, and adds an "Analytical FBP / FDK architecture" section
    explaining the scale-factor derivation.
  • docs/source/fdk_cone_example.rst: formula updated to reflect the
    voxel-gather kernel and (sid/U)^2 weight.
  • docs/source/fbp_fan_example.rst and fbp_parallel_example.rst:
    completely rewritten to reference the new analytical helpers and fix
    a sign-convention inconsistency in the math section.

Compatibility

This release is source-compatible with 1.2.9 for every public entry point.
The distance_weight parameter was internal and is removed from the
private backward kernels; no public API touched that argument.

PyPI

Install with pip install diffct==1.2.10 (PyPI page).