The Central Idea
Let
\[ \mathbf{x} \]represent the unknown image or field.
A measurement system produces data
\[ \mathbf{y} \]through a forward operator:
\[ \boxed{ \mathbf{y} = A\mathbf{x} + \boldsymbol{\varepsilon} }. \]Here:
- \(A\) describes the imaging physics;
- \(\mathbf{x}\) is the unknown object;
- \(\mathbf{y}\) is the observed data;
- \(\boldsymbol{\varepsilon}\) represents noise and measurement error.
What Reconstruction Means
Reconstruction seeks an estimate
\[ \hat{\mathbf{x}} \]that explains the measurements while respecting physical or structural constraints.
In ideal arithmetic, one might hope to write
\[ \mathbf{x} = A^{-1}\mathbf{y}. \]In real imaging problems, the operator may be rectangular, singular, ill-conditioned, incomplete, or too large to invert directly.
Images as Vectors
A digital image with \(m\) rows and \(n\) columns can be rearranged into a vector:
\[ \mathbf{x} \in \mathbb{R}^{mn}. \]A three-dimensional volume can likewise be stacked into a long vector.
This allows imaging problems to be studied using linear algebra, even though the unknown has a spatial interpretation.
Forward Operators
Different imaging systems use different forward operators.
Examples include:
- convolution for blur;
- line integrals for tomography;
- Fourier sampling in magnetic resonance imaging;
- wave propagation in diffraction imaging;
- projection and attenuation models in transmission imaging.
The Observation Operator Matters
Reconstruction quality depends on what the instrument actually measures.
Two imaging systems can observe the same object but produce very different inverse problems because their operators
\[ A_1 \quad\text{and}\quad A_2 \]preserve different spatial features.
Blur as Convolution
A blurred one-dimensional signal can be modeled as
\[ y(t) = \int k(t-s) x(s) \,ds + \varepsilon(t). \]In shorthand,
\[ \boxed{ y = k*x + \varepsilon }. \]The function \(k\) is the blur kernel or point-spread function.
Point-Spread Function
The point-spread function describes how the imaging system records an ideal point source.
A narrow point-spread function preserves fine detail.
A broad point-spread function spreads a point over many pixels and produces stronger blur.
Blur in the Frequency Domain
Convolution becomes multiplication after a Fourier transform:
\[ \widehat{y}(\omega) = \widehat{k}(\omega) \widehat{x}(\omega) + \widehat{\varepsilon}(\omega). \]Formally,
\[ \widehat{x}(\omega) = \frac{ \widehat{y}(\omega) }{ \widehat{k}(\omega) }. \]Why Deblurring Is Ill Posed
Blur kernels often suppress high-frequency image content.
Thus
\[ |\widehat{k}(\omega)| \]may become very small for large \(|\omega|\).
Dividing by these small values greatly amplifies high-frequency noise.
Deconvolution
Recovering an image from convolution data is called deconvolution.
A practical deconvolution method must balance:
\[ \boxed{ \text{sharpness} \leftrightarrow \text{noise amplification} }. \]Tomographic Imaging
In tomography, measurements often represent integrals of an unknown field along paths.
A simplified continuous model is
\[ y(\ell) = \int_{\ell} x(\mathbf{r}) \,ds. \]Each measurement combines information from many spatial locations.
The Radon Transform Idea
In two-dimensional tomography, line integrals through an object can be viewed as samples of the Radon transform.
Measurements taken at many angles encode information about the interior structure.
Reconstruction approximately inverts this transform.
Projection Data
A set of line-integral measurements at one viewing angle forms a projection.
Multiple projections at different angles provide complementary information about the object.
Limited angular coverage can leave some structures poorly determined.
Discretized Tomography
After discretization, tomography often takes the form
\[ \boxed{ \mathbf{y} = A\mathbf{x} + \boldsymbol{\varepsilon} }. \]Each row of \(A\) describes how one measurement depends on the image pixels or voxels.
System Matrix
The matrix
\[ A \]is often called the system matrix.
Its entries may represent:
- ray-pixel intersection lengths;
- blur weights;
- Fourier sampling coefficients;
- sensor sensitivities;
- wave-propagation effects.
Direct Reconstruction
Some imaging problems admit direct analytic or transform-based inversion formulas.
Direct methods can be fast, but they may be sensitive to incomplete data or noise.
Filtered Backprojection Idea
In classical tomography, a direct reconstruction can be built by filtering each projection and then spreading, or backprojecting, it through the image domain.
Filtering restores spatial frequencies suppressed by projection geometry.
Backprojection accumulates the information from multiple viewing directions.
Backprojection Alone Blurs
Simply smearing each projection back across the image produces a blurred reconstruction.
The filtering step compensates for that blur.
Iterative Reconstruction
Instead of applying a direct inverse formula, iterative methods repeatedly improve an image estimate.
A generic iteration has the form
\[ \mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \Delta\mathbf{x}^{(k)}. \]The update is chosen to reduce mismatch between measured and predicted data.
Data-Fidelity Objective
A basic least-squares reconstruction minimizes
\[ \boxed{ \Phi(\mathbf{x}) = \|\mathbf{y}-A\mathbf{x}\|_2^2 }. \]If the problem is ill posed, however, minimizing this term alone may amplify noise.
Regularized Reconstruction
A common reconstruction objective is
\[ \boxed{ \Phi_{\lambda}(\mathbf{x}) = \|\mathbf{y}-A\mathbf{x}\|_2^2 + \lambda R(\mathbf{x}) }. \]The regularizer
\[ R(\mathbf{x}) \]encodes plausible structure in the recovered image.
Quadratic Image Regularization
One simple choice is
\[ R(\mathbf{x}) = \|L\mathbf{x}\|_2^2. \]If \(L\) computes first or second differences, the penalty discourages rapid spatial variation and promotes smoothness.
Total Variation Reconstruction
Total variation regularization uses a penalty based on image gradients:
\[ R(\mathbf{x}) \approx \sum_j \|\nabla x_j\|. \]It suppresses small oscillatory noise while allowing relatively sharp edges.
Sparse Reconstruction
Some images are not sparse in pixel space but are sparse in another basis.
Let
\[ \mathbf{x} = \Psi\boldsymbol{\alpha}, \]where \(\boldsymbol{\alpha}\) contains transform coefficients.
If most coefficients are near zero, one may minimize
\[ \|\mathbf{y}-A\Psi\boldsymbol{\alpha}\|_2^2 + \lambda \|\boldsymbol{\alpha}\|_1. \]Compressed-Sensing Idea
If an image has a sparse representation and the measurement system satisfies suitable conditions, accurate recovery may be possible from fewer measurements than a naive dimension count would suggest.
The essential ingredients are:
- sparsity or compressibility;
- informative measurements;
- appropriate reconstruction constraints.
Nonnegativity
Many imaging quantities must satisfy
\[ x_j\ge0. \]Examples include:
- intensity;
- density;
- concentration;
- emission strength.
Enforcing nonnegativity can improve reconstruction stability and physical plausibility.
Bounded Reconstructions
Some image values have known physical limits:
\[ x_{\min} \le x_j \le x_{\max}. \]Including these bounds reduces the admissible reconstruction space.
Resolution
Resolution describes the ability of an imaging system to distinguish nearby or fine-scale features.
Resolution is limited by:
- sensor spacing;
- sampling density;
- blur;
- viewing geometry;
- noise;
- regularization.
Resolution Is Not the Same as Pixel Size
An image may be stored on a very fine grid while the data contain much less fine-scale information.
Small pixels do not automatically imply high physical resolution.
Point Resolution
One way to study resolution is to reconstruct an idealized point source.
If the reconstructed point spreads over a broad neighborhood, the imaging system has limited local resolution.
Resolution Matrix
In a linear reconstruction,
\[ \hat{\mathbf{x}} = B\mathbf{y}. \]If noiseless data satisfy
\[ \mathbf{y}=A\mathbf{x}_{\text{true}}, \]then
\[ \hat{\mathbf{x}} = BA\mathbf{x}_{\text{true}}. \]Define
\[ \boxed{ R = BA }. \]This is a resolution matrix.
Ideal Resolution
If
\[ R=I, \]every true component is reconstructed exactly in the noiseless linear model.
In real regularized inverse problems,
\[ R \ne I. \]The departure from the identity reflects smoothing, cross-talk, or unresolved structure.
Artifacts
A reconstruction artifact is a feature that appears in the recovered image but does not correspond faithfully to the underlying object.
Artifacts can arise from:
- limited data;
- incorrect forward models;
- strong noise;
- motion;
- undersampling;
- poor regularization choices.
Limited-Angle Artifacts
In tomography, missing viewing angles create weakly observed spatial directions.
Structures aligned with those directions may become blurred, elongated, or distorted.
Undersampling Artifacts
If the measurement grid is too sparse, high-frequency structure may be misrepresented as lower-frequency structure.
This is related to aliasing.
Ring and Streak Artifacts
Systematic detector errors or incomplete projection geometry can produce structured artifacts such as rings or streaks.
Their geometry often reveals the source of the reconstruction error.
Model Mismatch
Suppose the true imaging process is
\[ \mathbf{y} = A_{\text{true}} \mathbf{x} + \boldsymbol{\varepsilon}, \]but reconstruction uses
\[ A_{\text{model}}. \]Even small operator mismatch can produce systematic reconstruction bias in an ill-conditioned problem.
Motion Error
If the object moves while measurements are being collected, the data may no longer correspond to one fixed image
\[ \mathbf{x}. \]Motion correction may need to be incorporated into the forward model.
Blind Deconvolution
Sometimes both the image and blur kernel are unknown:
\[ \mathbf{y} = k*\mathbf{x} + \boldsymbol{\varepsilon}. \]Estimating both \(k\) and \(\mathbf{x}\) is called blind deconvolution.
This is a more difficult inverse problem because image and blur effects can trade off against one another.
Identifiability in Imaging
Imaging systems can contain exact or approximate null spaces.
If
\[ A\mathbf{v} = \mathbf{0}, \]then adding
\[ \mathbf{v} \]to the image does not change the measurements.
That image direction is unidentifiable from the available data.
Near-Null Image Features
If
\[ \|A\mathbf{v}\| \ll \|\mathbf{v}\|, \]then feature direction \(\mathbf{v}\) has only a weak effect on the data.
Reconstructing such features is highly sensitive to noise.
Singular Values in Imaging
The singular values of \(A\) quantify how strongly different image directions influence the observations.
Small singular values correspond to weakly observed spatial patterns.
Regularization suppresses these unstable directions.
Noise Propagation
If a linear reconstruction is
\[ \hat{\mathbf{x}} = B\mathbf{y} \]and
\[ \operatorname{Cov}(\mathbf{y}) = \Sigma_y, \]then
\[ \boxed{ \operatorname{Cov}( \hat{\mathbf{x}} ) = B\Sigma_yB^T }. \]Reconstruction algorithms therefore transform measurement noise into spatial uncertainty patterns.
Spatially Varying Uncertainty
Not every reconstructed pixel is equally reliable.
Geometry, sensor coverage, and regularization can make uncertainty much larger in some spatial regions than others.
Uncertainty Maps
A reconstruction may be accompanied by an uncertainty map showing the estimated variance or credible interval width at each location.
This can be more informative than displaying the reconstructed image alone.
Bayesian Imaging
Bayesian reconstruction combines a likelihood
\[ p( \mathbf{y} \mid \mathbf{x} ) \]with an image prior
\[ p(\mathbf{x}). \]The posterior is
\[ \boxed{ p( \mathbf{x} \mid \mathbf{y} ) \propto p( \mathbf{y} \mid \mathbf{x} ) p( \mathbf{x} ) }. \]MAP Reconstruction
The maximum a posteriori image is
\[ \boxed{ \hat{\mathbf{x}}_{\text{MAP}} = \operatorname*{arg\,max}_{\mathbf{x}} p( \mathbf{x} \mid \mathbf{y} ) }. \]Under Gaussian assumptions, this often reduces to regularized least squares.
Posterior Mean Imaging
The posterior mean is
\[ \mathbb{E} [ \mathbf{x} \mid \mathbf{y} ]. \]It averages over reconstruction uncertainty rather than selecting only the most probable image.
Image Priors
Priors can encode assumptions such as:
- smoothness;
- edge preservation;
- sparsity;
- nonnegativity;
- known texture statistics;
- known anatomical or geometric structure.
Prior Influence
In weakly observed image directions, the prior can have substantial influence on the reconstruction.
A visually sharp image may therefore reflect strong prior assumptions rather than strong measurement support.
Data Fidelity Versus Prior Structure
A reconstruction balances two questions:
- Does the image reproduce the measured data?
- Does the image satisfy plausible structural assumptions?
Reliable reconstruction requires both.
Iterative Gradient Reconstruction
For
\[ \Phi(\mathbf{x}) = \|\mathbf{y}-A\mathbf{x}\|_2^2, \]the gradient is
\[ \nabla\Phi = -2A^T ( \mathbf{y}-A\mathbf{x} ). \]A gradient step has the form
\[ \mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \alpha A^T ( \mathbf{y}-A\mathbf{x}^{(k)} ). \]Backprojection as an Adjoint Operation
The matrix
\[ A^T \]maps measurement residuals back into image space.
In many imaging systems, this is physically interpreted as a backprojection or adjoint operation.
Landweber Iteration
The update
\[ \mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \alpha A^T ( \mathbf{y}-A\mathbf{x}^{(k)} ) \]is a classical Landweber-type iteration for linear inverse problems.
Early stopping can act as regularization.
Early Stopping in Imaging
Early iterations tend to reconstruct strongly observed image components.
Later iterations may increasingly fit weak noisy components.
Stopping before full convergence can therefore reduce noise amplification.
Conjugate-Gradient Ideas
Quadratic reconstruction problems can often be solved more efficiently using conjugate-gradient methods rather than simple gradient descent.
These methods exploit the geometry of the linear system to accelerate convergence.
Proximal Reconstruction
Nonsmooth penalties such as L1 and total variation are often handled by proximal methods.
These methods alternate between a data-fit step and a regularization step.
Operator Splitting
Large imaging objectives can be decomposed into simpler subproblems.
Operator-splitting methods are useful when the objective contains:
- multiple penalties;
- constraints;
- nonsmooth terms;
- large forward operators.
Nonlinear Imaging
Some imaging systems obey a nonlinear forward map:
\[ \mathbf{y} = F(\mathbf{x}) + \boldsymbol{\varepsilon}. \]Examples include strong multiple scattering, nonlinear sensor responses, and coupled wave-interaction problems.
Reconstruction then requires nonlinear optimization or probabilistic inference.
Jacobian in Nonlinear Imaging
Local sensitivity is described by
\[ J = \frac{ \partial F }{ \partial\mathbf{x} }. \]The Jacobian maps small image perturbations to small measurement perturbations.
Adjoint-Based Gradients
In large nonlinear imaging problems, explicitly storing the full Jacobian may be impossible.
Adjoint methods can compute the gradient of a scalar objective without forming every Jacobian entry.
Image Registration
Multiple images may need to be spatially aligned before they can be compared or combined.
Estimating the transformation that aligns images is itself an inverse problem.
Multi-Modal Reconstruction
Different imaging modalities may observe complementary properties of the same object.
Joint reconstruction can combine measurements
\[ \mathbf{y}^{(1)} = A_1\mathbf{x} + \boldsymbol{\varepsilon}^{(1)} \]and
\[ \mathbf{y}^{(2)} = A_2\mathbf{x} + \boldsymbol{\varepsilon}^{(2)}. \]Complementary operators can improve identifiability.
Data Consistency
A reconstructed image should reproduce the measurements within their expected uncertainty:
\[ A\hat{\mathbf{x}} \approx \mathbf{y}. \]Large unexplained residuals indicate that the image, forward model, or uncertainty model may be inadequate.
Residual Images and Sinograms
In imaging, residuals can be displayed in measurement-space geometry.
Structured residual patterns can reveal:
- miscalibration;
- motion;
- incorrect blur kernels;
- missing physics;
- detector defects.
Visual Quality Is Not Enough
A reconstruction may look plausible while being quantitatively wrong.
Visual inspection should be accompanied by:
- data-consistency checks;
- resolution analysis;
- uncertainty analysis;
- synthetic recovery tests;
- validation against known structures.
Synthetic Phantom Tests
A known test object, often called a phantom, can be used to generate synthetic measurements:
\[ \mathbf{y}_{\text{syn}} = A\mathbf{x}_{\text{true}} + \boldsymbol{\varepsilon}. \]Reconstruction quality can then be evaluated against the known truth.
Resolution Versus Noise Experiments
Reconstruct the same phantom under several noise levels and regularization strengths.
This reveals how spatial detail degrades as uncertainty increases.
Common Error: Treating Pixel Size as Resolution
A fine reconstruction grid can create an illusion of detail.
True resolution is determined by the information content of the measurements and the reconstruction method.
Common Error: Inverting Small Singular Values Directly
Small singular values correspond to weakly observed image directions.
Direct inversion can amplify measurement noise catastrophically.
Regularization or truncation is often necessary.
Common Error: Over-Regularizing
Excessive regularization can erase real edges, small objects, or fine textures.
A smooth image is not automatically an accurate image.
Common Error: Under-Regularizing
Too little regularization allows weak inverse directions to fit noise.
The result may contain grain, oscillation, ringing, or unstable detail.
Common Error: Using the Wrong Blur Kernel
Deconvolution depends directly on the assumed point-spread function.
A mismatched kernel can create systematic sharpening errors or artifacts.
Common Error: Ignoring Boundary Conditions
Convolution and finite-difference regularizers require assumptions about image boundaries.
Periodic, zero, reflective, and replicated boundary rules can produce different edge behavior.
Common Error: Trusting One Reconstruction Algorithm
If several reasonable methods produce very different images, the data may not strongly determine the disputed features.
Comparing algorithms can expose weakly constrained structure.
Common Error: Ignoring Uncertainty
A reconstruction is an estimate, not an exact picture of the hidden object.
Weakly observed regions and fine features should be accompanied by appropriate uncertainty or resolution statements.
A Reliable Imaging Workflow
- Define the unknown image or field precisely.
- Write the forward imaging operator.
- Model the observation noise and calibration process.
- Inspect rank, singular values, or local sensitivity.
- Choose an image grid consistent with actual information content.
- Select physical constraints such as nonnegativity or bounds.
- Choose a regularizer that matches plausible image structure.
- Select the regularization strength using a principled criterion.
- Compute the reconstruction with a suitable direct or iterative method.
- Check data consistency and residual structure.
- Assess resolution and uncertainty.
- Run synthetic phantom recovery tests.
- Check sensitivity to modeling and regularization assumptions.
- Report which image features are strongly supported and which are uncertain.
Practice 1: Forward Imaging Model
Write the standard linear imaging model relating unknown image \(\mathbf{x}\), system matrix \(A\), data \(\mathbf{y}\), and noise \(\boldsymbol{\varepsilon}\).
Practice 2: Blur
What mathematical operation commonly models spatially invariant blur?
Practice 3: Deconvolution Instability
Why does dividing by a very small Fourier-domain blur factor amplify noise?
Practice 4: Tomography
What kind of measurements are commonly used in classical tomographic reconstruction?
Practice 5: Regularized Objective
Write a standard regularized least-squares imaging objective.
Practice 6: Nonnegativity
Give one reason a nonnegativity constraint can be useful in imaging.
Practice 7: Resolution
Does decreasing pixel size automatically improve physical image resolution?
Practice 8: Null Space
If
\[ A\mathbf{v} = \mathbf{0}, \]what does this imply about image component \(\mathbf{v}\)?
Practice 9: Noise Covariance
If
\[ \hat{\mathbf{x}} = B\mathbf{y} \]and
\[ \operatorname{Cov}(\mathbf{y}) = \Sigma_y, \]what is the covariance of the reconstructed image?
Practice 10: Over-Regularization
What visual or quantitative effect can excessive smoothness regularization produce?
Practice Answers
Practice 1
\[ \boxed{ \mathbf{y} = A\mathbf{x} + \boldsymbol{\varepsilon} }. \]Practice 2
Spatially invariant blur is commonly modeled by convolution:
\[ y=k*x. \]Practice 3
If
\[ |\widehat{k}(\omega)| \]is very small, then
\[ \frac{ \widehat{\varepsilon}(\omega) }{ \widehat{k}(\omega) } \]becomes very large. Weakly transmitted frequencies therefore amplify measurement noise during inversion.
Practice 4
Classical tomography commonly uses line-integral or projection measurements taken from multiple viewing angles.
Practice 5
\[ \boxed{ \Phi_{\lambda}(\mathbf{x}) = \|\mathbf{y}-A\mathbf{x}\|_2^2 + \lambda R(\mathbf{x}) }. \]Practice 6
Many physical image quantities such as intensity, density, and concentration cannot be negative. The constraint removes physically impossible solutions and can improve stability.
Practice 7
No. Smaller pixels increase numerical sampling density, but physical resolution is limited by the measurement operator, noise, geometry, and regularization.
Practice 8
The image direction \(\mathbf{v}\) is unobservable. Adding it to a candidate image does not change the measurements, so that component cannot be uniquely determined from the data.
Practice 9
\[ \boxed{ \operatorname{Cov}( \hat{\mathbf{x}} ) = B\Sigma_yB^T }. \]Practice 10
Excessive smoothness regularization can blur edges, suppress small structures, and remove real fine-scale detail.
What to Remember
Imaging and reconstruction are inverse problems in which the unknown is a spatial field, image, or volume inferred from indirect measurements.
- The forward operator determines what image information reaches the measurements.
- Blur is commonly modeled by convolution.
- Tomography reconstructs internal structure from projection or line-integral data.
- Small singular values correspond to weakly observed image directions.
- Direct inversion can amplify noise severely.
- Regularization balances data fidelity against plausible image structure.
- Quadratic, total-variation, sparse, and constraint-based regularizers serve different purposes.
- Pixel size is not the same as physical resolution.
- Artifacts can arise from limited data, model mismatch, motion, undersampling, and poor regularization.
- Iterative methods use forward and adjoint operations to reduce measurement mismatch.
- Uncertainty and resolution should accompany reconstructed images.
- Synthetic phantom tests are essential validation tools.
The central principle is:
\[ \boxed{ \text{reconstruct only the spatial detail that the measurements can actually support} }. \]Where We Go Next
Imaging and reconstruction bring together the major themes of inverse modeling: identifiability, least squares, likelihoods, Bayesian estimation, regularization, ill-posedness, measurement error, and computational optimization.
With the core 5.19 topics now assembled, the next step is to connect the individual lessons through the section index and ensure the sequence moves cleanly from forward modeling through inference, stabilization, sequential estimation, and reconstruction.