Lecture 04 — Image Pyramids
Lecture 4
Naïve Image Downsampling
Keep every other row and column, then repeat.
import matplotlib.pyplot as plt
y2 = I[::2, ::2]
y4 = y2[::2, ::2]
y8 = y4[::2, ::2]
fig, ax = plt.subplots(1, 3, figsize=(8, 4),
gridspec_kw={"width_ratios": [4, 2, 1]})
for a, y, s in zip(ax, [y2, y4, y8], [2, 4, 8]):
a.imshow(y, cmap="gray", vmin=0, vmax=255,
interpolation="nearest")
a.set_title(f"1/{s}")
a.axis("off")
plt.tight_layout()
plt.show()
Reminder
Camera / sensor hardware
CMOS sensor array
Images are a discrete, or sampled, representation of a continuous world.
Sampling
A simple example: a sine wave
How would you discretize this signal?
Sampling
Sample the signal at regularly spaced locations.
Sampling
The sampled signal consists of discrete values.
How many samples should I take?
Can I take as many samples as I want?
Sampling
The sampled signal consists of discrete values.
How many samples should I take?
Can I take as few samples as I want?
Undersampling
What happens when we take fewer samples?
Unsurprising effect: information is lost.
Undersampling
A different signal can agree with the same samples.
Unsurprising effect: information is lost.
Surprising effect: we can confuse the signal with one of lower frequency.
Undersampling
The same samples can also match a signal of higher frequency.
The samples alone do not uniquely identify the original signal.
We need an assumption about its highest frequency.
Aliasing
Aliasing: Undersampling disguises a signal as one of a lower frequency.
Aliasing in Textures
Aliasing in Photographs
These patterns are also known as moiré patterns.
Striped fabric
Building façades
AC grille
Temporal Aliasing
- The camera records separate snapshots; we do not see the motion between frames.
- Follow the dot: the wheel is rotating clockwise
- Without the dot, the four spokes are indistinguishable.
- A large clockwise step produces the same spoke positions as a small counterclockwise step.
The wheel can therefore appear to rotate backward—an example of temporal aliasing.
Anti-aliasing: Oversampling
How would you deal with aliasing?
Approach 1: Use more samples to represent fine detail.
Anti-aliasing: Smoothing
Approach 2: Smooth the image before downsampling.
Some fine detail is lost, but aliasing artifacts are reduced.
How would you smooth a signal?
Better Image Downsampling
Apply a Gaussian filter, then keep every other row and column. Repeat.
import cv2
import matplotlib.pyplot as plt
y2g = cv2.GaussianBlur(I, (3,3), 2)[::2, ::2]
y4g = cv2.GaussianBlur(y2g, (3,3), 2)[::2, ::2]
y8g = cv2.GaussianBlur(y4g, (3,3), 2)[::2, ::2]
fig, ax = plt.subplots(1, 3, figsize=(8, 4),
gridspec_kw={"width_ratios": [4, 2, 1]})
for a, y, s in zip(ax, [y2g, y4g, y8g], [2, 4, 8]):
a.imshow(y, cmap="gray", vmin=0, vmax=255,
interpolation="nearest")
a.set_title(f"1/{s}")
a.axis("off")
plt.tight_layout()
plt.show()
plt.close(fig)
Comparison of Downsampling Methods
Naïve Image Downsampling:
Smoothed Image Downsampling:
Image Pyramids
Why Use Image Pyramids?
Image compression
Texture mapping
Image blending
Denoising
Focal stack compositing
Multiscale detection
Multiscale registration
Also: Optical flow and feature tracking estimate motion from a coarse to fine paradigm.
Constructing a Gaussian Pyramid
cv2.pyrDown() combines Gaussian smoothing and downsampling.
import cv2
import matplotlib.pyplot as plt
G = [I]
for _ in range(3):
y = cv2.pyrDown(G[-1], borderType=cv2.BORDER_REFLECT)
G.append(y)
heights = [g.shape[0] for g in G[::-1]]
fig, ax = plt.subplots(4, 1, figsize=(10, 7.5),
gridspec_kw={"height_ratios": heights})
W = I.shape[1]
for a, i in zip(ax, range(3, -1, -1)):
h, w = G[i].shape
a.imshow(G[i], cmap="gray", vmin=0, vmax=255,
interpolation="nearest",
extent=((W-w)/2, (W+w)/2, h, 0))
a.set_xlim(0, W)
a.text(1.02, .5, rf"$G_{i}$", fontsize=28, transform=a.transAxes, va="center")
a.axis("off")
fig.subplots_adjust(left=.03, right=.87, bottom=.02, top=.98, hspace=.12)
plt.show()
What Happens at Coarser Levels?
What happens to the details of the image?
- Fine details get smoothed out as we move to higher, coarser levels.
- Small structures and texture become less distinguishable.
What is preserved at the higher levels?
- Mostly large uniform regions and broad intensity variations.
- The overall shape and arrangement of large structures.
Can We Recover the Original Image?
How would you reconstruct the original image from an upper level image?
- Exact reconstruction is not possible!
- Smoothing suppresses fine details.
- Downsampling retains fewer samples.
- Enlarging the coarse image cannot recover the missing information!
Blurring Removes Fine Detail
Gaussian smoothing suppresses detail before we even downsample.
What does the residual look like?
The residual is the original image minus the smoothed image.
What Did Smoothing Remove?
- The residual contains edges and fine texture suppressed by smoothing.
- Mid-gray represents zero; brighter and darker values represent positive and negative differences.
Can we make a pyramid that allows exact reconstruction?
Yes! Retain the residual information along with the coarse image. This is the idea behind a Laplacian pyramid.
Start with One Level
The original image equals an expanded coarse image plus a residual.
Do we need to store every blurred image as well as every residual?
No. We can recover each finer image from its residual and the next coarser image.
Laplacian Pyramid in OpenCV
Use floating-point arrays to preserve negative residual values.
G = [I.astype("float32")]
L = []
for _ in range(3):
y = cv2.pyrDown(G[-1], borderType=cv2.BORDER_REFLECT)
u = cv2.pyrUp(y, dstsize=G[-1].shape[::-1])
L.append(G[-1] - u)
G.append(y)
L.append(G[-1]) # Keep the smallest image too
heights = [im.shape[0] for im in L[::-1]]
fig, ax = plt.subplots(4, 1, figsize=(10, 7.5),
gridspec_kw={"height_ratios": heights})
W = I.shape[1]
for a, i in zip(ax, range(3, -1, -1)):
h, w = L[i].shape
s = max(float(abs(L[i]).max()), 1)
lo, hi = (0, 255) if i == 3 else (-s, s)
a.imshow(L[i], cmap="gray", vmin=lo, vmax=hi,
interpolation="nearest",
extent=((W-w)/2, (W+w)/2, h, 0))
a.set_xlim(0, W)
label = r"$L_3=G_3$" if i == 3 else rf"$L_{i}$"
a.text(1.02, .5, label, fontsize=28,
transform=a.transAxes, va="center")
a.axis("off")
fig.subplots_adjust(left=.03, right=.76, bottom=.02,
top=.98, hspace=.12)
plt.show()
What Must We Store?
What do we need to reconstruct the original image?
1. The residuals
\[ L_0,\quad L_1,\quad L_2 \]
These retain the detail lost between consecutive Gaussian levels.
2. The smallest Gaussian image
\[ L_3 = G_3 \]
This provides the starting point for reconstruction.
The intermediate Gaussian images can also be recovered as needed.
Residual contrast is scaled for display; mid-gray means zero.
Reconstruction uses the original signed arrays.
Reconstructing the Original Image
Start with the smallest image.
At each step:
- Expand to the next finer resolution.
- Add the residual at that level.
- Repeat until reaching the original resolution.
\[ G_i = \operatorname{expand}(G_{i+1}) + L_i \]
The residual restores the detail missing from the expanded image.
Keeping the residuals and smallest image allows exact reconstruction, apart from numerical roundoff.
Gaussian vs. Laplacian Pyramid
Which pyramid takes more space to store?
Their raw storage is equal when using the same data type.
A Laplacian pyramid can require more bits per value to represent signed differences and preserve precision.
They take about 33% more than the original image:
\[ N + \frac{N}{4} + \frac{N}{16} + \cdots \approx \frac{4}{3}N \]
Gaussian pyramid
Laplacian pyramid
Laplacian Pyramids Useful for Compression
Raw array size and compressed file size are different.
- Laplacian levels store differences between scales.
- In smooth regions, these differences are often zero or close to zero.
Frequently occurring values can be encoded using fewer bits
- The residuals can be inexpensive to store after compression, even though the raw pyramid contains more values than the original image.
Lossy compression
- We can round small residuals to zero or store them with less precision.
- This reduces storage further, at the cost of some image detail.
- Actual savings depend on the image and the encoding method.
Application: Selective Detail Enhancement
- Amplify foreground residuals, attenuate background residuals, then reconstruct:
\[ \widetilde{L}_i = \left[\alpha M_i+\beta(1-M_i)\right]L_i, \qquad \alpha=1.5,\quad \beta=0.25 \]
- \(M_i\) is the foreground mask at level \(i\): white selects the foreground, black selects the background, and gray gives a gradual transition.
Steerable Pyramids
- A multi-scale, multi-orientation image representation.
- Each pyramid level contains responses to structures at several orientations.
- A filter is steerable when the response at an arbitrary orientation can be synthesized from a small set of basis filters.
For first-order derivative basis responses:
\[ h_\theta = \cos(\theta)\,h_x + \sin(\theta)\,h_y \]
Applications
- Edge and contour analysis
- Texture and orientation analysis
- Image denoising and enhancement
- Image fusion
- Computational models of visual perception
A simplified first-order steerable pyramid: rows correspond to Gaussian pyramid levels \(G_0\), \(G_1\), and \(G_2\), while columns show orientation-selective responses synthesized from the basis responses \(h_x\) and \(h_y\).
Wavelets
A wavelet transform separates an image into spatial-frequency subbands:
- LL: coarse approximation
- LH: horizontal detail
- HL: vertical detail
- HH: diagonal detail
The LL image can be decomposed again, producing a multi-scale representation.
Applications
- Image compression — e.g., JPEG 2000
- Image denoising
- Feature extraction
- Image fusion
- Medical and scientific imaging
Haar wavelet decomposition: at each level, the input is separated into a low-frequency approximation LL and three directional detail bands. The LL band becomes the input to the next level.