Sampling and Aliasing

A camera sensor, a renderer, or a resizing routine records an image only at the points of a grid. Sampling theory says when those samples capture the image completely and what goes wrong when they do not. The failure, aliasing, makes detail finer than the grid reappear as false, coarser patterns: moiré on fabric and screens, jagged edges in renderings, and outputs of convolutional networks that change when the input moves by one pixel.

Definition

Sampling turns a continuous image $f(x, y)$ into a discrete one by recording its values on a regular grid, $f[m, n] = f(m\,\Delta x,\, n\,\Delta y)$, with sampling rates $1/\Delta x$ and $1/\Delta y$ samples per unit length. Aliasing is the effect by which a spatial frequency above half the sampling rate produces exactly the same samples as a lower frequency, its alias, so that after sampling the two cannot be told apart. Half the sampling rate is called the Nyquist frequency.

Intuition

Sample a sinusoid at a few points per cycle and the samples trace the wave faithfully. Sample it at fewer than two points per cycle and a slower wave passes through the very same points; nothing in the samples says which one was there. The fine pattern has not vanished. It has been relabeled as a coarse one that was never in the scene.

The same happens in time. A spoked wheel filmed at 24 frames per second can appear to turn slowly or even backwards, the wagon-wheel effect, because between frames each spoke moves nearly to the position of its neighbor. In images, a striped shirt photographed at the wrong distance shows broad, wavy color bands, and a high-contrast diagonal edge rendered with one sample per pixel becomes a staircase.

The cure follows from the picture: remove the detail the grid cannot represent before sampling, so that nothing is left to masquerade as something else.

Formal Definition

Shannon (1949) stated the sampling theorem for signals of one variable: a function containing no frequencies higher than $W$ is completely determined by its values at points spaced $1/(2W)$ apart. He reconstructed it by placing a $\sin(x)/x$ pulse at every sample,

$$f(t) = \sum_{n=-\infty}^{\infty} f\!\left(\frac{n}{2W}\right) \operatorname{sinc}(2Wt - n), \qquad \operatorname{sinc}(x) = \frac{\sin \pi x}{\pi x},$$

noted that mathematicians, citing Whittaker, had given the theorem earlier in other forms, and called $1/(2W)$ the Nyquist interval after Nyquist’s work on telegraphy.

The theorem extends directly to images on a rectangular grid. If the spectrum $F(u, v)$ of $f$ vanishes for $|u| \ge 1/(2\Delta x)$ or $|v| \ge 1/(2\Delta y)$, the image is recovered exactly by a product of sinc functions:

$$f(x, y) = \sum_{m} \sum_{n} f[m, n]\; \operatorname{sinc}\!\left(\frac{x}{\Delta x} - m\right) \operatorname{sinc}\!\left(\frac{y}{\Delta y} - n\right).$$

In units of cycles per sample, the band limit is $\pm\tfrac{1}{2}$ along each axis.

Aliasing is visible in the frequency domain. Sampling makes the spectrum periodic: the spectrum of the sampled image is the sum of copies of $F$ shifted by every multiple of the sampling rate along each axis. When $F$ is band-limited, the copies do not overlap and an ideal low-pass filter recovers the original. When it is not, the copies overlap, and a frequency $u$ is indistinguishable from $u - k/\Delta x$ for every integer $k$. In cycles per sample, a component at frequency $u$ appears at

$$u_{\text{alias}} = u - \operatorname{round}(u),$$

and likewise for $v$. A 2D frequency can therefore alias to a different orientation as well as a different scale, which is why moiré bands often run at an angle to the true pattern. The Fourier transform article in this series develops the frequency-domain view in full.

Anti-Aliasing and Reconstruction

Prefiltering

Aliasing cannot be undone after sampling, because the alias and the original produce identical samples. It must be prevented by low-pass filtering before or during sampling, attenuating everything above the new Nyquist frequency. The ideal prefilter is a perfect low-pass filter, whose kernel is a sinc. It is unusable in practice: the sinc has infinite extent, decays slowly, and its sharp cutoff causes ringing around edges. Practical prefilters, such as Gaussians, binomial kernels, and windowed sincs, trade a little attenuation of legitimate detail below the Nyquist frequency against incomplete suppression of frequencies just above it.

Prefiltering takes several forms:

• In cameras, the pixel’s light-collecting area already averages the image over a small region. Many cameras add an optical low-pass filter in front of the sensor that blurs the image slightly. Some high-resolution cameras omit it for sharper images, at the risk of moiré on fine repetitive patterns.
• When downsampling a digital image, the image is smoothed by convolution with a low-pass kernel before discarding samples. This is the smoothing step of every image pyramid; for subsampling by a factor $s$, the kernel must remove frequencies above $1/(2s)$ cycles per original pixel.
• When warping, a transformation that shrinks part of the image must be prefiltered by the local scale factor, or sample from a prefiltered pyramid, as described in image warping.
• In rendering, where the image is defined by an algorithm that computes point samples and cannot be filtered in advance, renderers take several samples per pixel (supersampling), often at jittered positions, and average them, which turns structured aliasing into less objectionable noise.

Reconstruction

The reverse step, turning samples into a continuous image (for display, resizing, or warping), is interpolation, and it is a convolution of the samples with a reconstruction kernel. Mitchell and Netravali (1988) distinguish aliasing caused by undersampling, which they call prealiasing, from postaliasing, caused by a reconstruction filter that fails to remove the spectral copies. They showed that a near-ideal windowed sinc greatly reduces postaliasing but rings at edges, and studied a family of piecewise cubic filters to balance blurring against ringing. Nearest-neighbor, bilinear, and cubic interpolation, compared in image warping, are increasingly close approximations to the ideal reconstruction filter.

Examples

The code samples two oblique cosine gratings at every fourth pixel, once directly and once after Gaussian smoothing, and reads off the strongest frequency and its amplitude in the subsampled image.

import numpy as np
from scipy import ndimage

N, factor = 500, 4
y, x = np.mgrid[0:N, 0:N].astype(float)


def grating(fx, fy):
    """Oblique cosine grating; fx, fy in cycles per pixel."""
    return np.cos(2 * np.pi * (fx * x + fy * y))


def dominant_frequency(img):
    """Frequency (cycles/sample, x then y) and amplitude of the strongest non-DC component."""
    F = np.fft.fft2(img) / img.size
    F[0, 0] = 0
    iy, ix = np.unravel_index(np.abs(F).argmax(), F.shape)
    fy, fx = np.fft.fftfreq(img.shape[0])[iy], np.fft.fftfreq(img.shape[1])[ix]
    return (fx, fy), 2 * np.abs(F[iy, ix])


sigma = factor / 2   # Gaussian prefilter width, in original pixels
for fx, fy in [(0.30, 0.20), (0.05, 0.03)]:
    img = grating(fx, fy)
    naive = img[::factor, ::factor]
    pre = ndimage.gaussian_filter(img, sigma, mode="wrap")[::factor, ::factor]
    (ax, ay), a_naive = dominant_frequency(naive)
    _, a_pre = dominant_frequency(pre)
    print(f"grating ({fx}, {fy}) cycles/px = ({fx * factor:.2f}, {fy * factor:.2f}) cycles/sample after subsampling")
    print(f"  naive:       strongest frequency ({ax:+.2f}, {ay:+.2f}) cycles/sample, amplitude {a_naive:.3f}")
    print(f"  prefiltered: amplitude {a_pre:.5f} "
          f"(Gaussian response {np.exp(-2 * np.pi**2 * sigma**2 * (fx**2 + fy**2)):.5f})")

Output:

grating (0.3, 0.2) cycles/px = (1.20, 0.80) cycles/sample after subsampling
  naive:       strongest frequency (-0.20, +0.20) cycles/sample, amplitude 1.000
  prefiltered: amplitude 0.00003 (Gaussian response 0.00003)
grating (0.05, 0.03) cycles/px = (0.20, 0.12) cycles/sample after subsampling
  naive:       strongest frequency (+0.20, +0.12) cycles/sample, amplitude 1.000
  prefiltered: amplitude 0.76461 (Gaussian response 0.76456)

The first grating, at $(1.2, 0.8)$ cycles per new sample, is far beyond the Nyquist limit of $\tfrac{1}{2}$. Subsampled directly, it survives at full amplitude as a grating at $(0.2, -0.2)$, the value predicted by subtracting the nearest integer from each component (a cosine at $(-0.2, 0.2)$ is the same pattern). Its stripes are coarser and have swung from about 34 degrees on one side of vertical to 45 degrees on the other: a moiré pattern. Gaussian prefiltering with $\sigma = 2$ pixels reduces it to $3 \times 10^{-5}$ of its amplitude, matching the filter’s frequency response $e^{-2\pi^2\sigma^2(u^2 + v^2)}$. The second grating is below the new Nyquist frequency and is sampled correctly either way, but the prefilter also attenuates it to 76% of its amplitude. That is the price of a practical filter: a wider Gaussian suppresses aliases better and blurs legitimate detail more.

Common Misconceptions

• “Sampling at twice the highest frequency is always enough.” The theorem requires the spectrum to vanish at and above the Nyquist frequency; a sinusoid exactly at it can be sampled at its zero crossings and disappear. It also assumes ideal sinc reconstruction from infinitely many samples.
• “Real images are band-limited.” An image of finite extent with sharp edges cannot be strictly band-limited, as Shannon noted for signals limited in time. In practice the goal is to make the energy above the Nyquist frequency negligible, not zero.
• “Aliasing can be removed afterwards with a blur.” Blurring after sampling softens the alias but cannot separate it from genuine content at the same frequency. Only filtering before sampling prevents it.
• “Downsampling is taking every $s$-th pixel.” Without prefiltering, that is exactly what produces aliasing, as the example and the image pyramid article show.
• “Convolutional networks are shift-invariant by construction.” Convolution is, but the subsampling in strided layers and pooling is not.

Where It Is Used

• Image acquisition and demosaicing. A color sensor with a color filter array samples each color more sparsely than the full pixel grid, so each color channel has a lower Nyquist frequency than the sensor as a whole. Fine detail then aliases into false colors and color moiré, which demosaicing algorithms try to avoid by exploiting correlations between channels.
• Resizing and multi-scale processing. Image pyramids, thumbnails, and coarse-to-fine methods all subsample, and all depend on smoothing first.
• Geometric transformation. Warping, texture mapping, and rectification combine reconstruction (interpolation) with resampling, and minification requires prefiltering.
• Rendering. Mitchell and Netravali credit Crow (1977) with identifying aliasing as the cause of artifacts in synthetic images, and quote a 1934 paper on picture scanning by Mertz and Gray that already described moiré and the rule of thumb of cutting off frequencies above half the scanning rate. Supersampling, mipmapping, and careful reconstruction filters are the standard defenses.
• Convolutional networks. Zhang (2019) showed that max-pooling, average pooling, and strided convolution break shift-equivariance because they subsample without regard to the sampling theorem. Writing each as a dense operation followed by subsampling, he inserted a low-pass filter before the subsampling step (blur kernels from $[1, 1]$ to the binomial $[1, 4, 6, 4, 1]$) and reported both greater shift-invariance and higher ImageNet accuracy across several architectures. Azulay and Weiss (2019) showed that small translations or rescalings can drastically change the predictions of modern CNNs. They attributed this to architectures that ignore the sampling theorem and to data augmentation that helps only for images resembling the training set, and found antialiasing and more augmentation to be only partial remedies.