Analysis of reconstruction from noisy discrete generalized Radon data

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} + \eta_{j,k}$. We show that the reconstruction error $N_\epsilon^{\text{rec}}=\mathcal A \eta_{j,k}$ satisfies $N^{\text{rec}}(\check x;x_0)=\lim_{\epsilon\to0}N_\epsilon^{\text{rec}}(x_0+\epsilon\check x)$, $\check x\in D$. Here $x_0$ is a fixed point, $D\subset\mathbb{R}^n$ is a bounded domain, and $\eta_{j,k}$ are independent, but not necessarily identically distributed, random variables. $N^{\text{rec}}$ and $N_\epsilon^{\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 $\eta_{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.