Week 1 · Day 2 Lab — Representation and colour¶

AIMS Cameroon · Introduction to Computer Vision

Learning goals

  1. Simulate sampling and quantisation and measure their cost with PSNR.
  2. Produce aliasing, then prevent it.
  3. Implement grayscale and RGB→HSV conversion from the formulas.
  4. Build a colour-based segmenter that survives a shadow.
  5. Reproduce (and fix) the classic dtype and channel-order bugs.
In [ ]:
# If you are on Google Colab or Kaggle, everything below is pre-installed.
# Locally:  pip install numpy matplotlib scipy scikit-image scikit-learn opencv-python
import numpy as np
import matplotlib.pyplot as plt
import cv2, skimage, scipy
from skimage import data, color, util, exposure, filters, feature, measure, transform
from scipy import ndimage as ndi

plt.rcParams["figure.dpi"] = 100
rng = np.random.default_rng(0)

def show(*images, titles=None, cmap="gray", size=4, **kw):
    """Display images side by side."""
    n = len(images)
    fig, axes = plt.subplots(1, n, figsize=(size * n, size))
    axes = np.atleast_1d(axes)
    for i, (ax, im) in enumerate(zip(axes, images)):
        ax.imshow(im, cmap=cmap if np.ndim(im) == 2 else None, **kw)
        ax.axis("off")
        if titles:
            ax.set_title(titles[i])
    plt.tight_layout(); plt.show()

print("numpy", np.__version__, "| opencv", cv2.__version__, "| scikit-image", skimage.__version__)

1. Quantisation¶

✏️ Exercise 1.1¶

Implement quantise(img, bits) mapping a uint8 image to $L=2^b$ levels with $q(v)=\lfloor vL/256\rfloor$, then rescale back to $[0,255]$ for display: $q\cdot 255/(L-1)$. Also implement psnr(ref, img).

In [ ]:
def quantise(img, bits):
    # TODO
    raise NotImplementedError

def psnr(ref, img):
    # TODO: 10 log10(255^2 / MSE)
    raise NotImplementedError
In [ ]:
cam = data.camera()
bits = range(1, 9)
scores = [psnr(cam, quantise(cam, b)) for b in bits]
plt.plot(list(bits), scores, "o-"); plt.xlabel("bits"); plt.ylabel("PSNR (dB)"); plt.grid(alpha=.3)
plt.title("Each extra bit buys roughly 6 dB"); plt.show()
show(*[quantise(cam, b) for b in (8, 4, 2, 1)], titles=["8 bits", "4 bits", "2 bits", "1 bit"], size=3)

Question. The slope of the curve is ≈ 6 dB/bit. Can you derive this? (Hint: the quantisation error is roughly uniform on an interval of width $\Delta=256/L$, so its variance is $\Delta^2/12$.)

2. Sampling and aliasing¶

✏️ Exercise 2.1¶

Downsample the zone plate by a factor of 6 in two ways:

  • naive: keep every 6th pixel
  • safe: Gaussian-blur first (σ ≈ 3), then keep every 6th pixel

Then also try cv2.resize(..., interpolation=cv2.INTER_AREA).

In [ ]:
n = 512
yy, xx = np.mgrid[-1:1:n * 1j, -1:1:n * 1j]
zone = 0.5 + 0.5 * np.cos(np.pi * 140 * (xx ** 2 + yy ** 2))
In [ ]:
naive = None  # TODO
safe = None   # TODO
area = None   # TODO
show(zone, naive, safe, area, titles=["original", "naive", "blur + subsample", "INTER_AREA"], size=3)

Question. The naive result shows rings that are not in the original. Using Nyquist ($f_s > 2f_{\max}$), explain where they come from. Which σ would you choose for a downsampling factor of 2? Of 10?

3. Colour spaces from the formulas¶

✏️ Exercise 3.1 — RGB → HSV¶

Implement rgb_to_hsv(rgb) for an RGB image with values in $[0,1]$, returning H, S, V all in $[0,1]$ (H = degrees/360, like scikit-image):

$V=\max(R,G,B)$, $\Delta=V-\min(R,G,B)$, $S=\Delta/V$ (0 if $V=0$),

$H = 60° \times \begin{cases}((G-B)/\Delta) \bmod 6 & V=R\\ (B-R)/\Delta+2 & V=G\\ (R-G)/\Delta+4 & V=B\end{cases}$, and $H=0$ if $\Delta=0$.

Vectorise it: no Python loops over pixels. Use np.where or boolean masks, and be careful about division by zero.

In [ ]:
def rgb_to_hsv(rgb):
    # TODO
    raise NotImplementedError
In [ ]:
ast = util.img_as_float(data.astronaut())
mine, ref = rgb_to_hsv(ast), color.rgb2hsv(ast)
err = np.abs(mine - ref)
err[..., 0] = np.minimum(err[..., 0], 1 - err[..., 0])   # hue is circular
print("max error:", err.max())
assert err.max() < 1e-6, "your HSV does not match skimage yet"
print("✓ rgb_to_hsv matches skimage")
show(ast, mine[..., 0], mine[..., 1], mine[..., 2], titles=["RGB", "H", "S", "V"], size=3)

4. A segmenter that survives a shadow¶

Below is a synthetic scene: red fruit on green foliage, lit by light that fades from left to right.

✏️ Exercise 4.1¶

  1. rgb_mask: pixels with $R > 0.5$.
  2. hsv_mask: pixels whose hue is red ($H<0.05$ or $H>0.95$) and $S>0.5$.
  3. Compute the fraction of true fruit pixels recovered by each (truth is given).
In [ ]:
H, W = 240, 480
scene = np.zeros((H, W, 3)); scene[...] = [0.25, 0.55, 0.2]
yy, xx = np.mgrid[0:H, 0:W]
truth = np.zeros((H, W), bool)
for cx in [60, 180, 300, 420]:
    truth |= (xx - cx) ** 2 + (yy - 120) ** 2 < 45 ** 2
scene[truth] = [0.85, 0.12, 0.10]
lit = np.clip(scene * np.linspace(1.0, 0.25, W)[None, :, None] + rng.normal(0, 0.01, scene.shape), 0, 1)
show(lit, titles=["scene"], size=6)
In [ ]:
hsv = None       # TODO
rgb_mask = None  # TODO
hsv_mask = None  # TODO
# TODO: print recall = (mask & truth).sum() / truth.sum() for both masks

Prove it. Show on paper that $S$ and $H$ do not change under $(R,G,B)\mapsto(kR,kG,kB)$ for $k>0$. Which pixels does the argument fail for?

5. Bugs you will meet in real life¶

5.1 Integer overflow¶

In [ ]:
a = np.full((2, 2), 200, np.uint8); b = np.full((2, 2), 100, np.uint8)
print("uint8 a + b      :", (a + b)[0, 0])        # wraps modulo 256
print("cv2.add(a, b)    :", cv2.add(a, b)[0, 0])  # saturates
print("float a + b      :", a.astype(float) + b)

✏️ Exercise 5.1¶

brighten_bad(img, 80) below is buggy. Write brighten(img, delta) that adds delta and saturates at 0 and 255, returning uint8. Show the two results side by side on data.camera().

In [ ]:
def brighten_bad(img, delta):
    return img + np.uint8(delta)
In [ ]:
def brighten(img, delta):
    # TODO
    raise NotImplementedError

show(brighten_bad(cam, 80), brighten(cam, 80), titles=["buggy (wrap-around)", "saturating"])

5.2 Layout: HWC vs CHW¶

NumPy/OpenCV use (H, W, C); PyTorch uses (C, H, W) (with a batch dimension in front).

In [ ]:
x = data.astronaut()
chw = x.transpose(2, 0, 1)
print(x.shape, "->", chw.shape, "-> back", chw.transpose(1, 2, 0).shape)

6. Storage and normalisation¶

✏️ Exercise 6.1¶

  1. storage_mb(M, N, C, bits) returns raw size in megabytes ($10^6$ bytes).
  2. Use it for: a Full-HD RGB frame; a 16-bit 512×512 CT slice; a 13-band 16-bit 10980×10980 Sentinel-2 tile.
  3. standardise(img) applies the ImageNet statistics to an RGB uint8 image and returns float32.
In [ ]:
def storage_mb(M, N, C, bits):
    # TODO
    raise NotImplementedError

MEAN = np.array([0.485, 0.456, 0.406]); STD = np.array([0.229, 0.224, 0.225])
def standardise(img):
    # TODO
    raise NotImplementedError
In [ ]:
z = standardise(data.astronaut())
print("per-channel mean:", z.reshape(-1, 3).mean(0).round(2), " std:", z.reshape(-1, 3).std(0).round(2))
print("red pixel 180 ->", round(float((180 / 255 - 0.485) / 0.229), 3))

7. Stretch — counting objects with connectivity¶

  1. Smooth data.coins(), build a foreground mask (try Otsu: filters.threshold_otsu; does it work everywhere? If not, try the watershed recipe below).
  2. Count coins with measure.label(mask, connectivity=1) (4-connectivity) and connectivity=2 (8-connectivity). Do they differ? Why or why not?
In [ ]:
# TODO