6.7 Image formation as a linear system🔗⧉

In the previous chapters, we have focused on linear operators from an input image to an output image. We now turn to a theme we foreshadowed with the corner camera example in Chapter 6.1, Introduction to forward and inverse linear methods. In general, as the Chapter 2.4, Image measurements as integrals over the plenoptic function established, pixel values are integrals of the incoming light. For fixed optics and geometry, incoherent light integration gives a linear map from incident radiance to expected, unclipped sensor values. The choice of unknown matters: the same measurements are generally nonlinear in scene depth or material parameters. We can model and discretize the scene into a list of numbers (a vector $\mathbf{x}$) which represent, for example, the color (radiance) of discretized light rays, or discrete pixels on a flat billboard at the focusing distance. We then discretize the image the same way into a vector $\mathbf{y}$, and the camera can be modeled as a matrix $A$:

The expected measurement is $A\mathbf{x}$; $\mathbf n$ is the recorded value minus that expectation. Writing noise additively does not imply that it is independent of the signal: photon shot-noise variance depends on the expected measurement (Chapter 2.8, Noise, signal-to-noise ratio and dynamic range). A calibrated offset must be subtracted or included explicitly. Each row of $A$ contains one measurement’s weights; each column shows how one scene element contributes across measurements.

The content of this chapter is the realization that an imaging system can be described by which weights you put in $A$. A pinhole, a lens, a coded aperture, a stereo pair: they differ only in their matrix. In turn, the properties of an imaging matrix inform us about the characteristics of an imaging system.

This provides a linear algebra perspective on big lessons. We said image formation is selecting and combining rays: a pinhole selects a single ray for each pixel, a lens combines a whole cone of them. The matrix $A$ is exactly that selection-and-combination, written down. A row that is zero everywhere except one entry selects one scene element; a row full of small weights combines many. Read a camera's matrix and you have read what the camera does to light.

To keep $A$ small enough to look at, we can work in flatland, a two-dimensional world with a one-dimensional scene and a one-dimensional sensor, and we discretize both, so that $A$ is a modest matrix we can draw as a picture (Figure 6.7.1⛶ fullscreen). Everything generalizes to the real 2-D-scene, 2-D-sensor case; the matrices just get larger (a 2-D image flattened into a long vector, $A$ correspondingly huge), and they are too big to draw but conceptually similar.

How we feed $A$ depends on how we choose to write down "the scene." There are two standard parameterizations, and the difference between them is the spine of this chapter. The first treats the scene as a flat billboard at a fixed depth: a single row of colored points. It is the obvious choice, it makes the matrices concrete, and it has one sharp limitation we will hit head-on. The second treats the scene as ray space (the set of all rays, the light field), and it is the parameterization that removes the limitation. We take them in turn.

6.7.1 Linear systems from a billboard to an image🔗⧉

Start with the simplest possible scene: a flat billboard standing at some chosen depth, facing the camera, painted with a 1-D pattern of color. Discretize it into a row of points, stack their values into $\mathbf{x}$, and ask what matrix each flatland camera applies. Because the billboard is flat and at a known depth, each system's $A$ is a small matrix we can render as a heat-map (bright where the weight is large, dark where it is zero) and the shape of that heat-map is the character of the camera (Figure 6.7.2⛶ fullscreen).

• Pinhole. Each sensor pixel sees, through the tiny hole, exactly one point of the billboard. The matrix is essentially a permutation (one bright entry per row, zero elsewhere), a re-indexing of the scene onto the sensor. In this geometrical toy model it introduces no blur; a physical pinhole has finite size and diffraction, and it throws away nearly all the light (every ray that does not pass through the hole is discarded). Sharp because each row selects; dim because it selects only one ray.
• Pinspeck. Leave an otherwise open view and insert a single small occluder instead of a single small hole, and you get the pinhole's photographic complement: an "anti-pinhole." Its matrix is the unoccluded transport matrix minus the blocked pattern (all-ones minus a pinhole pattern only in the equal-weight toy model): each pixel sums all but one of the scene's contributions. This is the principle behind accidental and occluder-based cameras, where an everyday object's shadow inadvertently encodes an image of the room.
• Lens. Each pixel now integrates a range of billboard points (the cone of rays the aperture admits), so the matrix is a banded blur matrix: a diagonal smeared into a band. On the focal plane the band collapses to the diagonal (sharp, $A\approx$ identity) and off it the band widens (the defocus disc, the circle of confusion from the lens chapter). The lens combines where the pinhole selects: far more light, at the cost of a single plane of focus.
• Corner camera. Let an occluding edge (the corner of a wall) cut across the scene. As you scan past the edge, each successive pixel sees one more sliver of the hidden region than the last, so the measurement is a cumulative sum of the scene and $A$ is essentially lower-triangular. Inverting the linear system, i.e. differentiating the measurement, undoes the cumulative sum and recovers the scene's angular structure, which is how the faint penumbra at a corner can be read as an image of what is around it.
• Stereo pair. Two pinholes side by side give two images; stack them and $A$ simply gets taller: more rows, two viewpoints of the same scene. At a known billboard depth the two blocks have known shifts. Estimating unknown depth requires correspondence or a search over depth-dependent operators; it is not a linear solve for billboard colors (Chapter 36.1, Multiple view geometry). The lesson here: more rows means more measurements, and more measurements can carry extra information: here, depth.
• Coded aperture. Replace the lens's clear aperture with a patterned mask. Now $A$ is a coded blur: still a banded matrix, but with a deliberately chosen pattern of weights rather than a smooth disc. The point of the design is to make that blur invertible: a clear disc blurs away certain spatial frequencies irretrievably, whereas a code can improve conditioning over a selected band, at a cost in light throughput, so some detail is more stably recoverable and (as a bonus) the scene's depth read off from which code best deblurs it. We are only naming it here; the design and the deconvolution are the subject of Part 33, Computational optics and coded imaging.

Six cameras, six matrices, one picture each. That gallery is the content of half this book in miniature: design the matrix, and you have designed the camera.

Compare the defocused lens with the corner camera. The lens attenuates some scene variations strongly, so inversion amplifies their noise. The corner camera’s cumulative-sum toy model has a simple difference inverse, but simplicity does not make it noise-free: differencing independent measurements doubles their noise variance, and physical angular-derivative scaling can amplify it further. Coded apertures change these tradeoffs; none guarantees recovery of all frequencies under arbitrary noise.

In practice, most restoration problems let us drop the billboard-and-depth machinery entirely and take a simpler view of $\mathbf{x}$. For denoising, demosaicking, deblurring, and super-resolution, the unknown we actually want is just an idealized image: the sharp, noise-free, fully-sampled picture a perfect pinhole at the sensor's own resolution would have recorded. There is no billboard, no depth, no ray space, only that clean image $\mathbf{x}$ and an operator $A$ that degraded it, adding noise, dropping color samples, blurring, or subsampling. This is the setting for the bulk of the book, and the one to hold in mind for the rest of this section: $\mathbf{x}$ is the image we wish we had, $\mathbf{y}$ is the one the camera gave us, and $A$ is the damage in between.

6.7.2 How solvable is an imaging problem?🔗⧉

Those two examples, a defocused lens you can barely invert and a corner camera you invert in one line, are a preview of the question this section addresses: given $A$ and a measurement $\mathbf{y}=A\mathbf{x}+\mathbf{n}$, how hard is it to recover $\mathbf{x}$?

Before trying to invert anything, inspect $A$ to predict how the inverse will behave. Two properties tell you almost everything: the operator's shape and its conditioning (Figure 6.7.4⛶ fullscreen).

A square, full-rank matrix determines one exact solution, though it can be sensitive to noise. A tall matrix has more measurements than unknowns, but only full column rank guarantees a unique least-squares solution. Repeated independent measurements can reduce noise; duplicate geometry does not fill an existing null space. A wide matrix necessarily has a null space. If $A\mathbf z=0$, then $\mathbf x$ and $\mathbf x+\mathbf z$ give the same noiseless data. Demosaicking illustrates this: one measured channel per pixel does not determine three unrestricted color values without additional assumptions.

Before the general case, look at the easy one, because it organizes everything that follows: $A$ square and diagonal. Then the equation decouples completely (each measurement is just one unknown times one diagonal entry, $y_i = a_{ii}\,x_i + n_i$) and inverting is nothing but per-element division, $x_i = y_i / a_{ii}$. You can read the trouble straight off the diagonal: wherever an entry $a_{ii}$ is large the inversion is gentle, and wherever it is near zero the division blows the noise up. No subspaces, no entanglement, just a list of independent gains. This case matters because many operators can be made diagonal by a change of basis. For the convolutions that dominate imaging (blur, defocus, motion) the Fourier basis does exactly this (each frequency is rescaled on its own; → Chapter 6.3, Fourier), turning deblurring into per-frequency division and ill-conditioning into simply the frequencies the blur zeroed out. Diagonalizable operators are the ones we truly understand.

Now conditioning, which is subtler and where most of the real pain lives, because even a square, technically invertible $A$ can be useless. The plain way to see a general operator (one no single change of basis diagonalizes) is through its singular values: the set of gains $\sigma$ by which it stretches or shrinks the different directions of the scene, the diagonal case's "gains" generalized to any matrix. Inverting $A$ means multiplying each direction by $1/\sigma$. If the singular values run from large down to nearly zero (a large condition number $\kappa(A)=\sigma_{\max}/\sigma_{\min}$) then the directions with tiny $\sigma$ get multiplied by an enormous $1/\sigma$, and any noise that happened to land on those directions is amplified without mercy. A lens blur is the textbook villain: it attenuates high spatial frequencies almost to nothing (tiny singular values up there), so naive deblurring divides the noisy high frequencies by almost zero and the image dissolves into a blizzard. Figure 6.7.4⛶ fullscreen shows this in one dimension, and it is the picture to keep: the true signal, blurred and lightly noised; the operator's singular values decreasing toward zero; the naive inverse exploding into garbage; and a regularized inverse recovering the signal cleanly.

Singular values alone do not specify recoverability. Their scale depends on the units of the scene and measurement, and the noise level must be compared with the attenuated signal in each measured direction. Whitening accounts for unequal or correlated measurement noise before assessing the relative gains.

Hadamard calls a problem well-posed when a solution exists, is unique, and depends continuously on the data Hadamard 1902. A finite full-rank system can satisfy all three and still be badly conditioned: continuity does not guarantee small error amplification. Rank deficiency creates nonuniqueness, and noisy overdetermined equations can be inconsistent. These are distinct reasons to replace exact inversion with a statistical estimate.

A useful way to read the operator is to ask what happens to differences between candidate scenes. If two scenes produce exactly the same measurement, their difference lies in the null space; no noiseless inversion can tell them apart. If the measurements differ only slightly, noise may hide the difference even though the operator is technically invertible. Adding more equations helps only when they constrain directions that the earlier measurements left ambiguous. Repeating the same exposure can reduce independent noise, but it does not by itself add a missing spatial frequency or a new view.

6.7.3 Regularization🔗⧉

When a problem is ill-posed, no amount of solver cleverness recovers what the null space erased or the noise drowned; the missing information has to come from outside the data. Regularization is the linear-algebra way of supplying it, and it is the cure the regularized inverse of Figure 6.7.4⛶ fullscreen (panel d) already quietly used, the single most important idea in the linear-algebra view of recovery. Rather than chase the data exactly, which means chasing the noise, modify the problem to make it balanced, encouraging the solution to have desirable characteristics.

Least squares minimizes $\lVert A\mathbf{x}-\mathbf{y}\rVert^2$. Its objective is a convex quadratic, with a unique minimizer only when $A$ has full column rank. If the rank is deficient, it has a flat set of minimizers. Noise does not always prevent an exact fit: a square invertible matrix fits even noisy measurements exactly, which is part of the problem.

Regularization is then simply a second term added to that score, a penalty that grows when $\mathbf{x}$ looks implausible, so the winning scene has to balance fitting the data against staying believable:

The first term insists the answer explain the measurement; the second term, $R(\mathbf{x})$, penalizes implausible scenes (ones that are too wiggly, or not smooth) and the knob $\lambda$ sets how much plausibility is allowed to override data fit. The classic choice just penalizes scenes that are large or that vary too fast, $R(\mathbf{x})=\lVert\mathbf{x}\rVert^2$ or $\lVert\nabla\mathbf{x}\rVert^2$, which gently reins in the runaway inversion and is the standard first resort Tikhonov 1963; a smarter penalty charges for total wiggle while sparing genuine edges, $R(\mathbf{x})=\lVert\nabla\mathbf{x}\rVert_1$, so recovered images stay sharp instead of smearing. We return to these choices, and to why they work, once we have the tools to appreciate them; for now the shape is the whole point. This single line (a data term plus a regularizer) is the skeleton of nearly every restoration method in the book, and inverting imaging operators this way is a field of its own Bertero & Boccacci 1998.

There are two kinds of "hard," and they are constantly confused. Some difficulty is merely numerical: the system is enormous and sparse, slow or unstable to solve by brute force. That kind yields to better tools (conjugate gradient, multigrid, FFT-based diagonalization, preconditioning) all of which live in Chapter 6.9, Efficient linear solvers, and none of which change what is recoverable; they just reach it faster and more stably. But conditioning and rank are properties of $A$ itself. No solver, however clever, recovers a scene direction the null space erased or the noise drowned. The distinction is the whole ethic of the recovery problems ahead: "hard to compute" is an engineering problem; "impossible to recover" is a fact about the operator and the noise, and the only things that change it are more measurements or stronger assumptions about what the scene can look like, full stop.

Consider a measurement that is just the average of two unknown pixels. An observed value of one half is consistent with (0,1), (1,0), and (0.5,0.5), among many other pairs. Least squares gives all these exact fits the same data cost. A penalty favoring neighboring values to be similar selects (0.5,0.5) in this example. That selection comes from the assumption, not from a second measurement.

Now suppose both pixels are measured individually with noise. The data term can have a unique minimum, yet matching the noisy values exactly may preserve unwanted fluctuations. Regularization trades some measurement fit for the chosen image property. Thus missing information and noisy information are different motivations for regularization; neither is captured by saying that a quadratic always has one useful minimum.

If a square matrix is invertible, unregularized least squares agrees with direct inversion: the data cost reaches zero at the exact solution. The minimization formulation therefore includes the familiar solve while also accommodating inconsistent measurements and explicit preferences among ambiguous solutions. Adding a regularizer changes the problem deliberately; it is not simply a faster way to compute the same inverse.

6.7.4 Linear systems from ray space (plenoptic function or light field) to an image🔗⧉

So far we have looked at a simple imaging model where our billboard parameterization is tied to a single depth. The matrix $A$ for the lens was a band whose width came from how far the billboard sat from the focal plane; the matrix for the stereo pair encoded a disparity that came from the billboard's distance. Move the billboard nearer or farther and every one of these matrices changes: the lens band changes width, and the stereo shift changes. A billboard matrix describes the imaging of one flat scene at one depth, and a real scene has many depths at once. We need a way to write down the scene that does not bake in a depth. That is exactly what the light field does.

The fix is to stop describing the scene as a surface and start describing it as rays, the plenoptic-function view introduced in Chapter 2.4, Image measurements as integrals over the plenoptic function. Instead of "the color of the billboard at position $u$," take as the unknown "the radiance of the ray arriving at the camera from direction so-and-so, through aperture position so-and-so." The full set of those rays is ray space, essentially the plenoptic function from Chapter 2.4, Image measurements as integrals over the plenoptic function, or its close cousin the light field (Part 31, Light fields and plenoptic cameras), and stacking all of it into $\mathbf{x}$ gives a parameterization that knows nothing about depth: depth is no longer a property of the operator; it is structure inside $\mathbf{x}$ (a depth edge is a particular pattern of which rays carry which color).

Once the scene is ray space, the imaging system becomes a fixed linear projection of that ray space down to the image, and the projection no longer depends on scene depth. Within this radiometric model, one fixed camera configuration describes different incident light fields, because the scene's depth has moved out of $A$ and into $\mathbf{x}$. The same cameras from the billboard gallery reappear, now as different ways of slicing the one light field (Figure 6.7.5⛶ fullscreen):

• a pinhole keeps only a thin slice of ray space: the rays through one point;
• a lens integrates ray space along a sheared line: it sums all the rays in the aperture that converge to one image point, and the shear of that line is the focus setting (refocusing is re-shearing);
• a stereo pair or light-field camera keeps more of ray space (more rows), which is exactly why it can recover depth and refocus;
• a coded aperture keeps a coded projection of the aperture's rays;
• a focus sweep averages projections over time. Wavefront coding changes wave phase and requires a wave-optics model; an incoherent intensity PSF can still describe its image formation under appropriate assumptions.

Every camera is a different projection of the same underlying light field. That single picture (drawn here after Levin, Freeman and Durand) ties the whole "cameras differ by how they slice the light field" theme together, and it is the reason the operator view is worth installing now.

6.7.5 Space–time: imaging across time🔗⧉

So far $\mathbf{x}$ has held a static scene. Following our plenoptic perspective, add time as another axis of $\mathbf{x}$, a moving scene captured over an exposure, and the picture is unchanged: it is still $\mathbf{y}=A\mathbf{x}$, now over space and time. The previous chapter already noted that the exposure is an integral over time; here that integral is just more rows and columns of $A$. A moving point smears along its trajectory, and motion blur is a linear operator whose kernel is the path the point traced during the exposure.

A flutter shutter modulates exposure in time to improve the response for an assumed motion family. A parabolic camera sweep can make blur approximately invariant to object speed along the sweep direction over a designed range. Neither makes arbitrary motion known or invertible: finite exposure, off-axis motion, occlusions, and light loss still matter. See Code in time and Chapter 30.1, Motion blur, temporal sampling, and resampling.

6.7.6 Other forms of linearity, and multilinearity🔗⧉

We have treated $A$ as a linear operator on the scene radiance $\mathbf{x}$, taking the illumination as given. But radiance is not the only quantity imaging is linear in. Light adds and scales, so the image is equally linear in the illumination: double every source and every pixel doubles; light a scene by the sum of two illuminants and you get the sum of the two images. That superposition is what makes light stages, photometric stereo, and relighting work.

The image is linear in still other quantities: for a fixed scene a convolutional blur is linear in its kernel exactly as it is linear in the image (the fact that makes blind deblurring so slippery), and a pixel is linear in the sensor's spectral response. Each is a different axis along which the same measurement is a linear map.

The complication is that these linearities usually multiply. A pixel is the product of what the surface reflects and what the light delivers; a blurred image is the product (convolution) of a sharp image and a kernel. Such a measurement is bilinear: linear in each factor with the other held fixed, but not jointly, and more generally multilinear when several unknowns combine. Joint least-squares objectives for bilinear models are generally non-convex even though each factor alone is linear, so when both factors are unknown we cannot simply invert a matrix; we alternate between them, or factor the product apart with priors. Everything earlier in this chapter (conditioning, invertibility, the reliance on priors) carries straight over, now doubled.

Color constancy separates illuminant spectrum and surface reflectance inside the sensor’s wavelength integral. That integral is bilinear in the two spectra when the other transport factors are fixed. A channel-wise RGB product is only an approximation, as discussed in Chapter 4.2, White balance and white point. The ambiguity is visible even in a scalar product: multiplying one factor by a constant and dividing the other by it leaves the measurement unchanged. Priors or calibration must constrain that scale and other ambiguities.

Brainard and Freeman used a Bayesian formulation to reason about color constancy, while Tenenbaum and Freeman used bilinear models to separate style from content. The recurring issue is how observations constrain two interacting factors. Gray-world assumptions, reflectance statistics and calibration constrain different aspects of the ambiguity; merely naming both factors does not make them separately observable. Chapter 7.1, Imaging as Bayesian inference with probabilistic priors develops how such assumptions enter an estimate.

6.7.7 Designing measurements for recovery🔗⧉

Assess a camera by the quantities it must recover under a specified photon budget, resolution, and noise model. A better condition number alone is insufficient if the design discards most of the light. Coded optics and additional exposures alter both the measurement operator and the noise statistics.

There are two complementary ways to make $A$ carry more recoverable information, and they correspond to the book's two later parts on the subject:

• Coded imaging and computational optics (Part 33, Computational optics and coded imaging) keep the number of measurements fixed but redesign $A$ to be better-conditioned (a coded aperture, wavefront coding, the lattice-focal lens Levin et al. 2009 (4D frequency / lattice-focal)) so that the same amount of data inverts more cleanly (and, often, reveals depth as a side effect).
• Light fields and multiple-exposure photography (Part 31, Light fields and plenoptic cameras, HDR merging) instead make $A$ taller (take more measurements, more views or more exposures) so the system is more determined and simply records more about the scene.

For shift-invariant systems with matching boundary conditions, Chapter 6.3, Fourier makes these tradeoffs visible frequency by frequency. More general cameras can require singular-value or probabilistic analysis, including the light-field comparisons of Levin and colleagues Levin et al. 2009 (4D frequency / lattice-focal). The relevant question is which scene variations become distinguishable in the measurements, and at what noise level.