Analysis of reconstruction from noisy discrete generalized Radon data

TL;DR

The paper develops a comprehensive local probabilistic analysis for tomographic reconstruction from noisy, discretized generalized Radon data in $\mathbb{R}^n$, where the data are modeled as $g=\mathcal{R}f+\eta$ and the reconstruction operator $\mathcal{A}$ is a Fourier integral operator with a phase linear in frequency. It proves that the reconstruction error, after applying $\mathcal{A}$ to the noise, converges (as the sampling step $\varepsilon\to0$) to a zero-mean Gaussian random field $N^{\text{rec}}$, with an explicit covariance given by $\text{Cov}(\check x,\check y)=C(\check x-\check y)$ computed from the principal symbol and the geometry of the transform via $G$ and $\Psi$. The results extend prior work on noisy 2D Radon data to a broad class of generalized Radon transforms (integrating over submanifolds of codimension $N$) and general FIO reconstruction operators, providing a rigorous, local probabilistic description of how discretization and noise affect reconstructions at native resolution. A numerical cone-beam CT example in $\mathbb{R}^3$ demonstrates excellent agreement between theory and simulation, highlighting the practical utility for assessing local detectability and resolution. Overall, the work delivers explicit local error statistics, enabling precise inference about small features near boundaries in tomographic images and informing robust reconstruction strategies under noise.

Abstract

We consider a wide class of generalized Radon transforms $\mathcal R$, which act in $\mathbb{R}^n$ for any $n\ge 2$ and integrate over submanifolds of any codimension $N$, $1\le N\le n-1$. Also, we allow for a fairly general reconstruction operator $\mathcal A$. The main requirement is that $\mathcal A$ be a Fourier integral operator with a phase function, which is linear in the phase variable. We consider the task of image reconstruction from discrete data $g_{j,k} = (\mathcal R f)_{j,k} + η_{j,k}$. We show that the reconstruction error $N_ε^{\text{rec}}=\mathcal A η_{j,k}$ satisfies $N^{\text{rec}}(\check x;x_0)=\lim_{ε\to0}N_ε^{\text{rec}}(x_0+ε\check x)$, $\check x\in D$. Here $x_0$ is a fixed point, $D\subset\mathbb{R}^n$ is a bounded domain, and $η_{j,k}$ are independent, but not necessarily identically distributed, random variables. $N^{\text{rec}}$ and $N_ε^{\text{rec}}$ are viewed as continuous random functions of the argument $\check x$ (random fields), and the limit is understood in the sense of probability distributions. Under some conditions on the first three moments of $η_{j,k}$ (and some other not very restrictive conditions on $x_0$ and $\mathcal A$), we prove that $N^{\text{rec}}$ is a zero mean Gaussian random field and explicitly compute its covariance. We also present a numerical experiment with a cone beam transform in $\mathbb{R}^3$, which shows an excellent match between theoretical predictions and simulated reconstructions.

Figures (5)

• Figure 1: Illustration of cone beam geometry.
• Figure 2: Illustration of the case when the set $Y_1(x,\check x)$ is an interval for a general source trajectory (left panel). The right panel illustrates the case of a circular source trajectory.
• Figure 3: Illustration of the case when the set $Y_2(x,\check x)$ has positive measure. The values $s\in Y_2(x,\check x)$ parametrize a line segment (shown in red), which is a subset of the source trajectory. The vector $\check x$ (shown in blue) is parallel to the line segment. The lines $L(x,s)$, $s\in Y_2(x,\check x)$, are all the same and contain the segment.
• Figure 4: Observed PDF and predicted PDF for the random variable $N_\epsilon^{\text{rec}}(\check x=0)$ (i.e., at $x_0$ in the original coordinates). The $x$-axis represents the values of the random variable.
• Figure 5: Observed PDF (left), predicted PDF (center), and their difference (right) for the random vector $(N_\epsilon^{\text{rec}}(\check x_1),N_\epsilon^{\text{rec}}(\check x_2))$. The display range is [-2E-5,1.7E-4] for all three figures. The horizontal and vertical axes represent the values of $N_\epsilon^{\text{rec}}(\check x_1)$ and $N_\epsilon^{\text{rec}}(\check x_2)$, respectively.

Theorems & Definitions (18)

• Definition 2.2
• Definition 2.7
• Theorem 3.1
• Theorem 3.2
• Definition 3.3
• Theorem 3.4
• Lemma 7.1
• Lemma 7.2
• Lemma 7.3
• Definition 8.1: Khoshnevisan2002