Week 1 · Day 1 Lab — An image is an array

AIMS Cameroon · Introduction to Computer Vision

Today you will make friends with the single data structure that every method in this course operates on: a NumPy array of numbers.

Learning goals

1. Load images in three different ways and inspect shape, dtype and value range.
2. Index, crop, flip and rotate images using array slicing only.
3. Quantify how "nuisance" transformations change the raw pixels.
4. Build a centre–surround (difference-of-Gaussians) filter and reproduce the Mach-band illusion.

Cells marked ✏️ Exercise contain TODOs for you. Work in pairs; ask lots of questions.

0. Environment check

# 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. Loading and inspecting an image

skimage.data ships a few classic test images, so this notebook works offline.

img = data.chelsea()          # a cat
print("shape :", img.shape)   # (rows, cols, channels) = (H, W, C)
print("dtype :", img.dtype)
print("range :", img.min(), "to", img.max())
show(img, titles=["data.chelsea()"])

✏️ Exercise 1.1 — Read single pixels

NumPy indexes images as img[row, col, channel], with row 0 at the top.

1. Store the RGB value of the pixel at row 100, column 200 in px.
2. Store the red channel of the 10×10 patch whose top-left corner is (140, 170) in patch.
3. Store the number of pixels in the image (not counting channels) in n_pixels.

px = None        # TODO
patch = None     # TODO
n_pixels = None  # TODO
print(px, patch.shape, n_pixels)

assert patch.shape == (10, 10) and n_pixels == 300 * 451
print("✓ Exercise 1.1 passed")

2. Three ways to load an image — and the BGR trap

We save the cat to disk, then read it back with scikit-image, OpenCV and Pillow.

from PIL import Image
skimage.io.imsave("cat.png", img)

a = skimage.io.imread("cat.png")          # RGB
b = cv2.imread("cat.png")                 # BGR (!)
c = np.array(Image.open("cat.png"))       # RGB
print("skimage == PIL :", np.array_equal(a, c))
print("skimage == cv2 :", np.array_equal(a, b))
print("skimage == cv2 reversed channels:", np.array_equal(a, b[..., ::-1]))
show(a, b, titles=["skimage (RGB)", "cv2 shown as if RGB"])

Practitioner's note. OpenCV stores colour images in BGR order for historical reasons. Mixing libraries without cv2.cvtColor(img, cv2.COLOR_BGR2RGB) is one of the most common silent bugs in vision pipelines: nothing crashes, but a model sees scrambled colours.

Also note: cv2.imread returns None (no error!) if the path is wrong. Always assert img is not None.

3. Geometry with slicing only

✏️ Exercise 3.1

Using only slicing and transpose (no np.flip, np.rot90, cv2):

• mirror: flip left–right
• upside: flip top–bottom
• rot90: rotate 90° counter-clockwise
• face: crop the cat's face, roughly rows 30–260 and columns 120–340

mirror = None  # TODO
upside = None  # TODO
rot90 = None   # TODO  hint: transpose swaps rows and columns
face = None    # TODO

assert np.array_equal(mirror, np.fliplr(img))
assert np.array_equal(upside, np.flipud(img))
assert np.array_equal(rot90, np.rot90(img))
show(mirror, upside, rot90, face, titles=["mirror", "upside down", "rot90", "face"], size=3)
print("✓ Exercise 3.1 passed")

4. Grayscale, the perceptual way

✏️ Exercise 4.1

Implement to_gray(rgb) using the ITU-R BT.601 luminance weights $Y = 0.299R + 0.587G + 0.114B$. Return a float64 array with the same value range as the input.

def to_gray(rgb):
    # TODO: one line with a matrix product or a weighted sum of channels
    raise NotImplementedError

g = to_gray(img)
ref = color.rgb2gray(img) * 255      # skimage uses slightly different (Rec. 709) weights
print("max |difference| vs skimage:", np.abs(g - ref).max().round(2), "grey levels")
print("pixel (220, 40, 30) ->", to_gray(np.array([[[220, 40, 30]]]))[0, 0].round(1))
show(g, titles=["BT.601 luminance"])

5. Same cat, different pixels

On Day 1 we argued that recognition is hard because nuisance factors change the pixels a lot while leaving the meaning untouched. Let's measure this.

✏️ Exercise 5.1

Write pixel_corr(a, b) returning the Pearson correlation between the grayscale versions of two equally sized images (flatten them, then use np.corrcoef).

def pixel_corr(a, b):
    # TODO
    raise NotImplementedError

h, w = img.shape[:2]
variants = {
    "illumination x0.35": np.clip(img * 0.35, 0, 255).astype(np.uint8),
    "rotate 35°": (transform.rotate(img, 35, mode="edge") * 255).astype(np.uint8),
    "shift 20 px": np.roll(img, 20, axis=1),
    "occlusion": img.copy(),
    "noise": util.img_as_ubyte(util.random_noise(img, var=0.04, rng=1)),
    "blur": cv2.GaussianBlur(img, (0, 0), 6),
}
variants["occlusion"][60:260, 120:330] = [60, 110, 60]
for name, v in variants.items():
    print(f"{name:20s} corr = {pixel_corr(img, v):+.2f}")
show(*variants.values(), titles=list(variants), size=2.6)

Discuss with your neighbour

1. Why is the correlation for the darkened image exactly 1.0? (Hint: what is correlation invariant to?) Would a non-linear change (e.g. gamma) also give 1.0? Try img ** 0.5.
2. A 20-pixel shift is invisible to you. What does it do to the correlation? What does this suggest about comparing images pixel by pixel?

6. A retina in five lines: centre–surround filtering

Retinal ganglion cells respond like a difference of Gaussians (DoG): an excitatory centre minus an inhibitory surround.

✏️ Exercise 6.1

1. Write gaussian_kernel(size, sigma) returning a normalised 2D Gaussian (size odd, weights sum to 1).
2. Build dog = gaussian_kernel(15, 1.0) - gaussian_kernel(15, 3.0).

def gaussian_kernel(size, sigma):
    # TODO: build a grid of (x, y) offsets with np.meshgrid, apply the Gaussian formula, normalise
    raise NotImplementedError

dog = None  # TODO

assert abs(gaussian_kernel(9, 2).sum() - 1) < 1e-9
fig, ax = plt.subplots(1, 2, figsize=(9, 3.5))
ax[0].imshow(dog, cmap="RdBu_r"); ax[0].set_title("DoG kernel"); ax[0].axis("off")
ax[1].plot(dog[7], "o-"); ax[1].axhline(0, color="k", lw=0.5); ax[1].set_title("centre row")
plt.show()

✏️ Exercise 6.2 — Mach bands

Build a staircase image (8 uniform vertical strips getting brighter), filter it with your DoG using ndi.convolve(..., mode="nearest"), and plot one row of the response. Where do the overshoots occur, and how does this explain the light/dark bands you see at each step?

stairs = None  # TODO: shape (60, 320), values in [0, 1]
resp = None    # TODO
# TODO: plot stairs and one row of resp

7. Reflection (write 2–3 sentences each)

1. Name three nuisance factors that matter for your sticky-note application.
2. Which part of today's lab is most like a simple cell? What operation would you add to mimic a complex cell?
3. You found one image library that silently returned the wrong channel order. What defensive check will you add to your future pipelines?