Filtered Backpropagation

From CT to RF -- Filtered Backpropagation

Plain backpropagation (matched filter) produces blurred images because the k-space samples are not uniformly distributed: samples are denser near the origin and sparser at high spatial frequencies. Filtered backpropagation corrects this non-uniformity by applying a frequency- domain weight before backprojecting.

The idea originated in X-ray computed tomography (CT), where the Fourier Slice Theorem shows that back-projection over-represents low spatial frequencies by a factor 1/∣f∣1/|f|. The classical Ram-Lak filter ∣f∣|f| compensates exactly. For RF imaging, we adapt this principle to the non-uniform k-space coverage of multi-static multi-frequency systems.

Definition:

Filtered Backpropagation (FBP)

The filtered backpropagation estimator applies a k-space weighting before the adjoint operation:

c^FBP(p)=∑m=1MWm ym ejκs,rm⋅p,\hat{c}_{\text{FBP}}(\mathbf{p}) = \sum_{m=1}^{M} W_m \, y_m \, e^{j {\kappa_{\mathbf{s},\mathbf{r}}}_{m} \cdot \mathbf{p}},

where WmW_m is a filter weight that compensates for the local density of k-space samples near κs,rm{\kappa_{\mathbf{s},\mathbf{r}}}_{m}.

In matrix form:

c^FBP=AHD−1y,\hat{\mathbf{c}}_{\text{FBP}} = \mathbf{A}^{H} \mathbf{D}^{-1} \mathbf{y},

where D=diag(d1,…,dM)\mathbf{D} = \text{diag}(d_1, \ldots, d_M) is a diagonal weighting matrix with dmd_m encoding the local k-space sampling density at measurement mm.

Interpretation: FBP is a weighted adjoint -- a pseudo-inverse restricted to the sampled k-space. Where the matched filter uses AH\mathbf{A}^{H}, FBP uses AHD−1\mathbf{A}^{H} \mathbf{D}^{-1}, partially inverting the non-uniform sampling.

Theorem: FBP as a Restricted Pseudo-Inverse

Let A=FΩΦ\mathbf{A} = \mathbf{F}_\Omega \boldsymbol{\Phi} where FΩ\mathbf{F}_\Omega selects rows of the 2D DFT matrix corresponding to the sampled k-space locations Ω\Omega, and Φ\boldsymbol{\Phi} is a discretization matrix.

For uniform k-space sampling (∣Ω∣|\Omega| points on a regular grid), the filtered backpropagation estimator satisfies:

c^FBP=Φ−1FΩHDΩ−1y=Φ−1F−1[FΩHDΩ−1y]zero-filled,\hat{\mathbf{c}}_{\text{FBP}} = \boldsymbol{\Phi}^{-1} \mathbf{F}_\Omega^H \mathbf{D}_\Omega^{-1} \mathbf{y} = \boldsymbol{\Phi}^{-1} \mathbf{F}^{-1} [\mathbf{F}_\Omega^H \mathbf{D}_\Omega^{-1} \mathbf{y}]_{\text{zero-filled}},

which is the zero-filled inverse FFT after density compensation -- the standard gridding reconstruction used in MRI and diffraction tomography.

For non-uniform k-space sampling, the density compensation factors dmd_m can be computed via Voronoi tessellation of the k-space locations or iteratively via the method of Pipe and Menon (1999).

Proof
Fourier structure

For a translation-invariant system in the far field, each measurement samples the Fourier transform of the reflectivity at wavenumber κs,rm{\kappa_{\mathbf{s},\mathbf{r}}}_{m}:

ym=c^(κs,rm)+wm.y_m = \hat{c}({\kappa_{\mathbf{s},\mathbf{r}}}_{m}) + w_m.

Back-projection as partial IDFT

The plain backpropagation computes:

c^BP(p)=∑m∈Ωc^(κs,rm)ejκs,rm⋅p,\hat{c}_{\text{BP}}(\mathbf{p}) = \sum_{m \in \Omega} \hat{c}({\kappa_{\mathbf{s},\mathbf{r}}}_{m}) e^{j {\kappa_{\mathbf{s},\mathbf{r}}}_{m} \cdot \mathbf{p}},

which is an IDFT restricted to Ω\Omega. If Ω\Omega is non-uniformly sampled, low-frequency components (near κs,r=0\kappa_{\mathbf{s},\mathbf{r}} = 0) are over-represented.

Density compensation

Weighting by Wm=1/dmW_m = 1/d_m where dmd_m is the local sampling density yields a uniformly-weighted IDFT:

c^FBP(p)=∑m∈Ω1dmc^(κs,rm)ejκs,rm⋅p,\hat{c}_{\text{FBP}}(\mathbf{p}) = \sum_{m \in \Omega} \frac{1}{d_m} \hat{c}({\kappa_{\mathbf{s},\mathbf{r}}}_{m}) e^{j {\kappa_{\mathbf{s},\mathbf{r}}}_{m} \cdot \mathbf{p}},

which correctly approximates c(p)c(\mathbf{p}) up to the bandwidth limit of Ω\Omega.

,

Definition:

Ram-Lak and Hamming Filters

The density compensation weight is often expressed as a 1D filter applied to each projection before backprojection:

Ram-Lak (ramp) filter:

HRL(f)=∣f∣,∣f∣≤fmax⁡.H_{\text{RL}}(f) = |f|, \quad |f| \leq f_{\max}.

This is the theoretically exact compensator for the 1/∣f∣1/|f| density of radial sampling.

Hamming-windowed ramp filter:

HHam(f)=∣f∣⋅[0.54+0.46cos⁡(πf/fmax⁡)].H_{\text{Ham}}(f) = |f| \cdot \bigl[0.54 + 0.46\cos(\pi f / f_{\max})\bigr].

The Hamming window suppresses the high-frequency amplification of the ramp filter, reducing noise and ringing at the cost of some resolution.

Filter Noise amplification Resolution Edge preservation
Ram-Lak High Best Good
Hamming Moderate ∼\sim1.3×\times wider Very good
Hann Low ∼\sim1.5×\times wider Moderate
Shepp-Logan Low-moderate ∼\sim1.1×\times wider Best

Example: Filtered Backpropagation for Circular Aperture

A monostatic radar rotates through 360 degrees around a scene (L=360L = 360 aperture positions, 1 degree spacing), transmitting Nf=128N_f = 128 frequencies over bandwidth B=1B = 1 GHz at fc=10f_c = 10 GHz.

(a) Show that the k-space coverage forms a filled disk.

(b) Derive the density compensation weight for this geometry.

(c) Compare the PSFs of plain backpropagation and filtered backpropagation.

Solution
k-space coverage

At angle ϕ\phi and frequency ff, the monostatic wavenumber is κs,r=(4πf/c)[cos⁡ϕ,sin⁡ϕ]T\kappa_{\mathbf{s},\mathbf{r}} = (4\pi f/c)[\cos\phi, \sin\phi]^T. As ϕ\phi sweeps 0--360 degrees and ff sweeps fc±B/2f_c \pm B/2, the coverage is an annular disk with inner radius 4π(fc−B/2)/c4\pi(f_c - B/2)/c and outer radius 4π(fc+B/2)/c4\pi(f_c + B/2)/c. For B≪fcB \ll f_c, this approximates a filled disk.

Density compensation

In polar k-space coordinates (k,ϕ)(k, \phi), the Jacobian is kk. The sampling density at radius kk is proportional to 1/k1/k (more samples near the origin due to the radial geometry). The density compensation weight is W(k)∝kW(k) \propto k, which is exactly the Ram-Lak filter ∣f∣|f| expressed in k-space.

PSF comparison
  • Plain BP: Low-frequency bias produces a blurred PSF with suppressed edges and enhanced smooth features.
  • FBP with Ram-Lak: The PSF becomes a 2D sinc-like function with well-defined mainlobe and controlled sidelobes. Resolution matches the diffraction limit c/(2B)c/(2B) isotropically.

Filtered vs. Unfiltered Backpropagation

Compares plain backpropagation with filtered backpropagation for a point-scatterer scene.

Left: Plain backpropagation (matched filter). Right: Filtered backpropagation with the selected filter.

Observe how the ramp filter sharpens the PSF mainlobe and reduces the low-frequency bias. The Hamming-windowed filter trades some resolution for reduced noise amplification.

Parameters
20
60

When Is FBP Sufficient?

Filtered backpropagation is a direct (non-iterative) method -- a single pass through the data produces the image. This makes it extremely fast (O(Nlog⁡N)O(N \log N) with FFT) but limits its capabilities:

  • FBP assumes uniform or smoothly-varying k-space density. For highly non-uniform coverage (e.g., limited-angle aperture), the density compensation factors become large, amplifying noise.
  • FBP cannot impose priors such as sparsity or positivity.
  • FBP struggles with missing k-space regions (angular gaps).

For these scenarios, iterative methods (LASSO, ADMM, OAMP) provide superior results by incorporating regularization, at the cost of higher computation. FBP remains the standard first-pass algorithm and the baseline for comparison.

Common Mistake: Density Compensation Can Amplify Noise

Mistake:

Applying the Ram-Lak filter ∣f∣|f| without considering the noise level. At high spatial frequencies, ∣f∣|f| amplifies both signal and noise, potentially degrading the image quality.

Correction:

Use a windowed ramp filter (Hamming, Shepp-Logan) or truncate the ramp at a frequency below the Nyquist limit. The optimal cutoff depends on the SNR: lower SNR requires more aggressive windowing. Alternatively, use iterative methods that naturally handle the noise through regularization.

🔧Engineering Note

Non-Uniform FFT for Efficient FBP

When k-space samples are non-uniformly distributed (the common case in multi-static RF imaging), the direct summation for FBP has complexity O(MN)O(MN). The non-uniform FFT (NUFFT) accelerates this to O(Mlog⁡M+Nlog⁡N)O(M \log M + N \log N) by:

  1. Gridding: Convolving the non-uniform samples onto a uniform grid using a Kaiser-Bessel or min-max interpolation kernel.
  2. FFT: Applying a standard inverse FFT on the gridded data.
  3. Deapodization: Dividing by the Fourier transform of the gridding kernel to remove the convolution blur.

NUFFT implementations (e.g., FINUFFT, torchkbnufft) achieve near-FFT speed with controllable approximation error, making FBP practical for large-scale RF imaging.

Practical Constraints
  • •

    NUFFT accuracy depends on the oversampling factor (typically 2×2\times) and kernel width (6--8 points)

  • •

    Gridding introduces interpolation error that limits dynamic range to approximately 60--80 dB

Filtered backpropagation (FBP)

An image reconstruction method that applies a frequency-domain filter (typically the ramp filter ∣f∣|f|) before backpropagation to compensate for non-uniform k-space sampling density. Produces sharper images than plain backpropagation.

Density compensation

Weighting k-space samples by the inverse of the local sampling density to achieve approximately uniform spectral representation before inverse Fourier transform.

Key Takeaway

Filtered backpropagation corrects the non-uniform k-space density inherent in RF imaging by applying a ramp-like filter before the adjoint operation. It sharpens the PSF compared to plain backpropagation and is the RF analogue of filtered back-projection in CT. However, it remains a linear, non-iterative method that cannot exploit sparsity or handle missing k-space regions -- limitations that motivate the iterative algorithms of Ch 14.

The Point-Spread Function (PSF)Adaptive Beamforming for Imaging