Computer Vision

Before detecting anything, we need a precise mental model of what an image is. A grayscale image is a discrete 2D function $I: \mathbb{Z}^2 \to [0, 255]$ (or $[0,1]$ after normalization). A color image has three such functions stacked (R, G, B channels). Pixels are samples of a continuous scene, quantized in both space (the grid) and amplitude (integer intensity).

Convolution is the central operation. Convolving image $I$ with kernel $K$ produces output $O$: $$O(x,y) = \sum_{u}\sum_{v} I(x-u, y-v) \cdot K(u,v)$$ Think of sliding $K$ over every pixel location, taking a weighted sum of the neighborhood. This one operation implements blurring, sharpening, edge detection, and (stacked in layers) entire deep networks.

The Gaussian filter is the canonical smoother. Its kernel is sampled from $G_{\sigma}(x,y) = \frac{1}{2\pi\sigma^2}e^{-(x^2+y^2)/2\sigma^2}$. Applying it blurs the image, attenuating high-frequency noise. The parameter $\sigma$ controls the spatial scale: small $\sigma$ preserves fine detail; large $\sigma$ blurs broadly. This is not arbitrary — the Gaussian is the only kernel that satisfies all the axioms of a linear scale-space (causality, isotropy, linearity), meaning it's the theoretically correct way to look at an image at multiple scales.

import numpy as np
from scipy.ndimage import convolve

def gaussian_kernel(sigma, size=None):
    """Build a 2D Gaussian kernel."""
    if size is None:
        size = int(6 * sigma + 1) | 1  # odd size, ~3-sigma radius
    k = size // 2
    y, x = np.mgrid[-k:k+1, -k:k+1]
    g = np.exp(-(x**2 + y**2) / (2 * sigma**2))
    return g / g.sum()  # normalize so weights sum to 1

kernel = gaussian_kernel(sigma=2.0)
smoothed = convolve(image.astype(float), kernel)

Edges are locations where image intensity changes rapidly — they correspond to object boundaries, surface creases, and lighting discontinuities. Detecting them is the first step toward understanding scene structure.

Gradient-based detection. The image gradient $\nabla I = (I_x, I_y)$ measures how fast intensity changes in each direction. We approximate derivatives with finite differences. The Sobel operator uses a slightly smoothed version: $$K_x = \begin{bmatrix}-1&0&1\-2&0&2\-1&0&1\end{bmatrix}, \quad K_y = K_x^T$$ Edge magnitude is $|\nabla I| = \sqrt{I_x^2 + I_y^2}$ and direction is $\theta = \text{atan2}(I_y, I_x)$.

The Canny edge detector (1986) is the gold standard. It does four steps: (1) Gaussian smoothing to suppress noise; (2) Sobel gradient computation; (3) non-maximum suppression — thin edges to single-pixel width by keeping only local maxima along the gradient direction; (4) hysteresis thresholding — use two thresholds $T_{high}$ and $T_{low}$: pixels above $T_{high}$ are definite edges; pixels between thresholds are edges only if connected to a definite edge. This produces clean, thin, connected edge chains.

The Laplacian of Gaussian (LoG) is an alternative: smooth with Gaussian then apply the Laplacian $\nabla^2 I = I_{xx} + I_{yy}$. Edges correspond to zero crossings of the LoG response. The Difference of Gaussians (DoG) — subtracting two Gaussian-blurred images at different $\sigma$ — approximates LoG cheaply and is central to SIFT's blob detection.

import cv2
import numpy as np

def canny_demo(img_gray):
    # Step 1: Gaussian smoothing (cv2.Canny does this internally)
    blurred = cv2.GaussianBlur(img_gray, (5, 5), sigma=1.4)
    # Step 2-4: Canny (non-max suppression + hysteresis)
    edges = cv2.Canny(blurred, threshold1=50, threshold2=150)
    return edges

# Manual Sobel (illuminates what's happening)
Ix = cv2.Sobel(img_gray, cv2.CV_64F, 1, 0, ksize=3)
Iy = cv2.Sobel(img_gray, cv2.CV_64F, 0, 1, ksize=3)
magnitude = np.sqrt(Ix**2 + Iy**2)
direction = np.arctan2(Iy, Ix)

Edges are useful but ambiguous for matching: a point on a straight edge looks the same as its neighbors along the edge — the aperture problem. Corners (and more generally, interest points) are locatable in two dimensions and thus matchable across images.

The Harris corner detector (Harris & Stephens, 1988) measures how much the image changes when a small window is shifted in any direction. For a window $W$ centered at $(x,y)$, the structure tensor (or second moment matrix) is: $$M = \sum_{(u,v)\in W} \begin{bmatrix} I_x^2 & I_x I_y \ I_x I_y & I_y^2 \end{bmatrix}$$ The eigenvalues $\lambda_1, \lambda_2$ of $M$ tell the story:

• Both small: flat region
• One large, one small: edge
• Both large: corner

To avoid expensive eigenvalue computation, Harris uses the response: $$R = \det(M) - k \cdot (\text{trace}(M))^2 = \lambda_1\lambda_2 - k(\lambda_1+\lambda_2)^2$$ with $k \approx 0.04$–$0.06$. Points with $R > $ threshold are corners.

Shi-Tomasi (1994) simplifies: use $R = \min(\lambda_1, \lambda_2)$ directly. A point is a corner if the smaller eigenvalue exceeds a threshold. This is what OpenCV's goodFeaturesToTrack uses.

Key property: Harris is rotation-invariant (the eigenvalues don't change under rotation) but not scale-invariant. If you zoom in, a corner may look like an edge. Fixing this requires automatic scale selection — the foundation of SIFT.

import cv2
import numpy as np

def harris_corners(img_gray, k=0.04, threshold_ratio=0.01):
    img = img_gray.astype(np.float32)
    # Compute image gradients
    Ix = cv2.Sobel(img, cv2.CV_64F, 1, 0, ksize=3)
    Iy = cv2.Sobel(img, cv2.CV_64F, 0, 1, ksize=3)
    # Structure tensor elements (smoothed)
    sigma = 1.5
    Ixx = cv2.GaussianBlur(Ix*Ix, (0,0), sigma)
    Iyy = cv2.GaussianBlur(Iy*Iy, (0,0), sigma)
    Ixy = cv2.GaussianBlur(Ix*Iy, (0,0), sigma)
    # Harris response
    det_M   = Ixx*Iyy - Ixy**2
    trace_M = Ixx + Iyy
    R = det_M - k * trace_M**2
    # Threshold and non-maximum suppression
    corners = (R > threshold_ratio * R.max())
    return corners, R

SIFT (Lowe, 1999/2004) is one of the most influential algorithms in computer vision. It detects and describes local image features that are invariant to scale, rotation, and partially invariant to illumination change and viewpoint. Understanding SIFT teaches you the full pipeline from raw pixels to matchable descriptors.

Stage 1 — Scale-Space Extrema Detection. Build a Gaussian scale-space pyramid: repeatedly blur with Gaussians at doubling scales $\sigma, 2\sigma, 4\sigma\ldots$, grouped into octaves (each octave halves the image resolution). Compute the Difference of Gaussians (DoG) between consecutive Gaussian levels — this approximates the LoG (blob detector). Find 3D local extrema in $(x, y, \sigma)$ space: a candidate keypoint must be larger or smaller than all 26 neighbors (8 in the same scale level, 9 above, 9 below).

Stage 2 — Keypoint Localization. Fit a 3D quadratic to refine the extremum to sub-pixel/sub-scale accuracy (Taylor expansion of the DoG function). Reject low-contrast points ($|D(\hat{x})| < 0.03$) and edge-like points (using the ratio of principal curvatures from the Hessian of DoG).

Stage 3 — Orientation Assignment. For each keypoint, compute gradient magnitudes and directions in a neighborhood at the keypoint's scale. Build a 36-bin histogram of orientations (weighted by gradient magnitude and a Gaussian). The dominant peak becomes the canonical orientation. This makes the descriptor rotation-invariant.

Stage 4 — Descriptor Computation. Take a 16×16 pixel patch around the keypoint (at the keypoint's scale), divided into a 4×4 grid of cells. In each cell, build an 8-bin gradient orientation histogram. Concatenate all 4×4×8 = 128 values. Normalize to unit length to reduce illumination effects; clamp values at 0.2 and renormalize to reduce nonlinear illumination effects.

The result: a 128-D vector that describes the local appearance so robustly that SIFT matches can be found across significant viewpoint changes, scale differences, and lighting variations. Matching is done by nearest-neighbor search in descriptor space, using Lowe's ratio test: accept a match only if the distance to the nearest neighbor is $< 0.8$ times the distance to the second nearest.

import cv2

def sift_match_demo(img1, img2):
    """Detect SIFT keypoints and match between two images."""
    sift = cv2.SIFT_create(nfeatures=0, contrastThreshold=0.04)

    # Detect and compute 128-D descriptors
    kp1, des1 = sift.detectAndCompute(img1, None)
    kp2, des2 = sift.detectAndCompute(img2, None)

    # BFMatcher with L2 norm, then Lowe ratio test
    bf = cv2.BFMatcher(cv2.NORM_L2, crossCheck=False)
    raw_matches = bf.knnMatch(des1, des2, k=2)

    good = []
    for m, n in raw_matches:
        if m.distance < 0.75 * n.distance:  # Lowe ratio test
            good.append(m)

    # good[i].queryIdx, good[i].trainIdx index into kp1, kp2
    print(f"{len(good)} good matches from {len(kp1)} / {len(kp2)} keypoints")
    return good, kp1, kp2

Segmentation partitions an image into regions that correspond to meaningful parts of the scene. Unlike edge detection (which finds boundaries) or keypoints (which find sparse points), segmentation assigns every pixel a label.

K-Means clustering treats each pixel as a point in feature space (e.g., RGB or RGB + xy position) and clusters them. It's fast and simple but ignores spatial coherence and requires choosing $k$.

Graph-based segmentation (Shi & Malik, Normalized Cuts, 2000; Felzenszwalb & Huttenlocher, 2004) models the image as a graph $G = (V, E)$ where vertices are pixels and edge weights encode similarity (color, texture, gradient). Segmentation becomes a graph partitioning problem. Normalized cuts minimizes the cut cost normalized by the total edge weight within each partition — this is equivalent to a generalized eigenvalue problem on the graph Laplacian, giving globally coherent segments. The Felzenszwalb-Huttenlocher method is a greedy efficient alternative using minimum spanning trees.

Watershed segmentation treats the image gradient magnitude as a topographic surface. Imagine flooding the terrain from local minima: the ridges where floodwaters from different minima meet are segment boundaries. Implemented with priority-queue flooding from gradient minima. Tends to over-segment; markers (seeds) are often provided manually or from SIFT-style detectors.

GrabCut (Rother et al., 2004) combines graph cuts with Gaussian Mixture Models. The user draws a rough rectangle around the foreground; the algorithm alternates between fitting GMMs to foreground/background pixels and computing min-cut to update the segmentation. It's interactive and produces excellent results from minimal user input.

Modern approach: semantic segmentation. Fully Convolutional Networks (FCN, Long et al., 2015) replace the fully connected layers of a classification CNN with 1×1 convolutions and upsample with transposed convolutions, producing a per-pixel class label map. U-Net adds skip connections between encoder and decoder for fine-grained boundary recovery. DeepLab uses atrous (dilated) convolutions to increase the field of view without losing resolution.

# Normalized Cuts — conceptual skeleton
import numpy as np
from scipy.sparse.linalg import eigsh
from scipy.sparse import diags

def normalized_cuts(W):  # W: N×N affinity matrix (sparse)
    """Return Fiedler vector for 2-way Ncut segmentation."""
    d = np.array(W.sum(axis=1)).flatten()  # degree
    D = diags(d)
    D_inv_sqrt = diags(1.0 / np.sqrt(d))
    # Normalized graph Laplacian: L_sym = D^{-1/2}(D-W)D^{-1/2}
    L_sym = D_inv_sqrt @ (D - W) @ D_inv_sqrt
    # Fiedler vector = eigenvector for 2nd smallest eigenvalue
    vals, vecs = eigsh(L_sym, k=2, which='SM')
    fiedler = D_inv_sqrt @ vecs[:, 1]  # convert back
    # Threshold Fiedler vector → binary segment labels
    labels = (fiedler > 0).astype(int)
    return labels

Object detection must answer: where in the image is the object, and what class is it? Classical approaches slide a window over the image at multiple scales, classify each window, and merge detections.

Viola-Jones face detector (2001) was the first real-time face detector. Its key ideas:

1. Haar-like features: simple rectangular filter differences (sums of pixels in adjacent rectangles). Fast to compute via the integral image (each Haar feature computed in O(1) after O(N) preprocessing).
2. AdaBoost: train a cascade of weak classifiers (each a single Haar feature threshold) into a strong classifier. With ~6000 features, only ~20 are needed per stage.
3. Attentional cascade: arrange stages from fastest to most discriminative. Windows that fail early (definitely not a face) are rejected immediately. On a typical image, >99% of windows are rejected at stage 1–2, making the detector fast despite exhaustive sliding window search.

HOG (Histogram of Oriented Gradients, Dalal & Triggs, 2005) became the dominant pedestrian/person descriptor. For each detection window: compute gradient magnitude and direction at each pixel; divide into 8×8 cells; build a 9-bin orientation histogram per cell (unsigned, 0–180°); normalize over 2×2 blocks of cells (L2 norm). A 64×128 window gives a ~3780-D HOG vector. Train a linear SVM on HOG features. HOG+SVM on sliding window was state-of-the-art for pedestrian detection until deep learning.

Deformable Part Models (DPM, Felzenszwalb et al., 2010) extended HOG by modeling objects as a root HOG filter plus part filters with spatial deformation costs (spring model). Matched human performance on PASCAL VOC 2010. DPM was the best classical detector before deep convolutional nets.

# HOG feature extraction — conceptual
def hog_cell(patch_8x8):
    """Compute 9-bin unsigned gradient histogram for one 8x8 cell."""
    Ix = patch_8x8[:, 1:] - patch_8x8[:, :-1]  # x-derivative
    Iy = patch_8x8[1:, :] - patch_8x8[:-1, :]  # y-derivative
    # pad to same size
    Ix = np.pad(Ix, ((0,0),(0,1)))
    Iy = np.pad(Iy, ((0,1),(0,0)))
    mag = np.sqrt(Ix**2 + Iy**2)
    angle = np.degrees(np.arctan2(Iy, Ix)) % 180  # unsigned
    bins = np.arange(0, 180, 20)  # 9 bins × 20°
    hist, _ = np.histogram(angle, bins=9, range=(0, 180),
                           weights=mag)
    return hist

# Use cv2.HOGDescriptor for production use

The ImageNet Large Scale Visual Recognition Challenge (ILSVRC) tested classification of 1.2 million images into 1000 classes. In 2012, AlexNet (Krizhevsky, Sutskever, Hinton) achieved 15.3% top-5 error vs. 26.2% for the next-best method — a stunning 11-point gap. This moment marks the inflection point where handcrafted features were abandoned in favor of end-to-end learned representations.

Why CNNs work for vision:

• Local connectivity: each neuron sees only a small receptive field, matching the local structure of natural images.
• Weight sharing: the same filter is applied everywhere (translation equivariance). This dramatically reduces parameters and encodes the prior that useful features (edges, corners, textures) can appear anywhere.
• Pooling: max-pooling provides translation invariance within small windows and builds hierarchical representations.
• Depth: stacking layers composes simple features into complex ones. Layer 1 detects Gabor-like edges; layer 2 detects corners and junctions; layer 3 detects object parts; later layers detect object categories.

Key architectural milestones:

• AlexNet (2012): 5 conv layers + 3 FC layers, ReLU activations, dropout for regularization, trained on 2 GPUs.
• VGGNet (2014): All 3×3 filters, 16–19 layers. Showed depth matters more than filter size.
• GoogLeNet/Inception (2014): Inception modules compute multiple filter sizes in parallel, concatenate. Very deep (22 layers) but parameter-efficient via 1×1 bottleneck convolutions.
• ResNet (2015): The insight that was missing: residual connections $y = F(x) + x$. Instead of learning the full mapping, learn the residual (what to add to the input). This solves the vanishing-gradient problem, enabling 50, 101, 152 layer networks. ResNet-152 achieved 3.6% top-5 error — below human-level.
• Vision Transformer (ViT, 2020): Apply standard Transformer self-attention to non-overlapping 16×16 image patches. No convolutions at all. With enough data (or strong augmentation), ViT matches or beats ResNet.

Transfer learning is the practical payoff: a ResNet trained on ImageNet has learned universal visual features. Fine-tune its last layer on your 500-image medical dataset and get excellent results — you are leveraging millions of training examples you didn't collect.

import torch
import torch.nn as nn

# The residual block — the key idea of ResNet
class ResidualBlock(nn.Module):
    def __init__(self, channels):
        super().__init__()
        self.conv1 = nn.Conv2d(channels, channels, 3, padding=1, bias=False)
        self.bn1   = nn.BatchNorm2d(channels)
        self.conv2 = nn.Conv2d(channels, channels, 3, padding=1, bias=False)
        self.bn2   = nn.BatchNorm2d(channels)
        self.relu  = nn.ReLU(inplace=True)

    def forward(self, x):
        residual = x                       # identity shortcut
        out = self.relu(self.bn1(self.conv1(x)))
        out = self.bn2(self.conv2(out))
        out = out + residual               # add skip connection
        return self.relu(out)              # activations AFTER addition

# Transfer learning: fine-tune ResNet on a new task
import torchvision.models as models
model = models.resnet50(weights='IMAGENET1K_V2')  # pretrained
for p in model.parameters():                      # freeze backbone
    p.requires_grad = False
model.fc = nn.Linear(2048, num_classes)           # replace head
# Now train only model.fc on your dataset

Classification tells you what is in the image. Detection tells you what and where, producing bounding boxes and class labels for every object instance.

R-CNN (2014): Run Selective Search to generate ~2000 region proposals, warp each to a fixed size, run a CNN on each, classify with an SVM, and refine boxes with a linear regressor. Accurate but slow (47 seconds/image).

Fast R-CNN (2015): Run CNN once on the whole image to get a feature map. Project each proposal onto the feature map, extract a fixed-size feature with RoI Pooling, then classify and regress. 25× faster than R-CNN.

Faster R-CNN (2015): Eliminate Selective Search. Add a Region Proposal Network (RPN) that slides over the feature map and predicts objectness scores + box offsets for $k$ anchor boxes at each position (anchors are prior boxes of various scales/aspect ratios). RPN and detection head share the feature map — the whole network is end-to-end trainable. Near real-time at ~5–7 fps.

YOLO (You Only Look Once, 2016): Reformulate detection as a single regression problem. Divide the image into an $S \times S$ grid. Each cell predicts $B$ bounding boxes (x,y,w,h,confidence) and $C$ class probabilities — a single forward pass, no proposals. YOLO v1: 45 fps. YOLOv8 (2023): state-of-the-art accuracy and speed. The anchor-free YOLO variants (FCOS, YOLOv8) directly predict box centers and sizes without prior anchors.

Feature Pyramid Networks (FPN, 2017): Multi-scale detection is hard because small objects need high-res features (early layers) and large objects need semantically rich features (late layers). FPN builds a top-down pathway that merges semantic information from deep layers back into shallow high-res layers, giving a pyramid of feature maps each rich in both semantics and spatial detail.

Loss functions: detection combines classification loss (cross-entropy or focal loss) with localization loss (smooth L1 / IoU-based loss). IoU (Intersection over Union) is the core metric: $\text{IoU} = \frac{|A \cap B|}{|A \cup B|}$. Detections with IoU > 0.5 with a ground-truth box count as true positives; evaluation uses mean Average Precision (mAP) across IoU thresholds and classes.

# YOLO-style single-cell prediction (simplified)
# Grid cell at (i,j) predicts: [x, y, w, h, conf, cls_0, ..., cls_C]
import torch
import torch.nn as nn

class YOLOHead(nn.Module):
    """Single prediction head cell — simplified."""
    def __init__(self, in_ch, num_anchors, num_classes):
        super().__init__()
        self.pred = nn.Conv2d(
            in_ch,
            num_anchors * (5 + num_classes),  # 5 = x,y,w,h,conf
            kernel_size=1
        )

    def forward(self, feature_map):
        # feature_map: [B, in_ch, H, W]
        raw = self.pred(feature_map)   # [B, A*(5+C), H, W]
        return raw

# IoU between two boxes [x1,y1,x2,y2]
def iou(box_a, box_b):
    xa = max(box_a[0], box_b[0]); ya = max(box_a[1], box_b[1])
    xb = min(box_a[2], box_b[2]); yb = min(box_a[3], box_b[3])
    inter = max(0, xb-xa) * max(0, yb-ya)
    area_a = (box_a[2]-box_a[0]) * (box_a[3]-box_a[1])
    area_b = (box_b[2]-box_b[0]) * (box_b[3]-box_b[1])
    return inter / (area_a + area_b - inter + 1e-6)

Beyond bounding boxes, we often want pixel-level understanding.

Semantic segmentation assigns every pixel a class label (road, sky, person, car) but does not distinguish between separate instances of the same class. Instance segmentation further separates individual object instances (person #1, person #2).

Fully Convolutional Networks (FCN, 2015) were the first end-to-end trainable semantic segmentation networks. The trick: replace the FC layers of VGG/AlexNet with 1×1 convolutions (preserving spatial information), then upsample with transposed convolutions (also called deconvolutions or strided transpose convolutions) to restore the original resolution. Skip connections fuse fine spatial detail from early layers with semantic information from deep layers.

U-Net (2015) was developed for biomedical image segmentation and became ubiquitous. Its architecture: an encoder (contracting path) of conv+pool blocks halves spatial size while doubling channels, followed by a decoder (expanding path) of transposed convolutions. Crucially, skip connections concatenate encoder feature maps directly to the corresponding decoder level, letting the network recover fine boundaries. With data augmentation, U-Net achieves excellent results even on small datasets (30 training images in the original paper).

DeepLab series (2015–2018) introduces two key ideas: atrous (dilated) convolutions — insert zeros into the kernel to expand the receptive field without losing resolution or adding parameters. A dilation rate of $r$ gives the kernel a receptive field of $(2r+1)\times(2r+1)$ with the same number of weights as a $3\times3$. ASPP (Atrous Spatial Pyramid Pooling) applies atrous convolutions at multiple dilation rates in parallel, capturing multi-scale context. DeepLabv3+ achieved state-of-the-art on PASCAL VOC 2012 for years.

Mask R-CNN (2017) extends Faster R-CNN to instance segmentation. After RoI Pooling, an additional small FCN head predicts a binary mask for each proposed region, in parallel with the class and box heads. The key improvement over RoI Pooling: RoI Align uses bilinear interpolation instead of quantization, preserving spatial alignment critical for masks. Mask R-CNN is clean, modular, and was quickly extended to keypoint detection (human pose estimation).

# U-Net building block in PyTorch
import torch
import torch.nn as nn
import torch.nn.functional as F

class DoubleConv(nn.Module):
    def __init__(self, in_ch, out_ch):
        super().__init__()
        self.block = nn.Sequential(
            nn.Conv2d(in_ch, out_ch, 3, padding=1, bias=False),
            nn.BatchNorm2d(out_ch), nn.ReLU(inplace=True),
            nn.Conv2d(out_ch, out_ch, 3, padding=1, bias=False),
            nn.BatchNorm2d(out_ch), nn.ReLU(inplace=True),
        )
    def forward(self, x): return self.block(x)

class UNetDecoder(nn.Module):
    def __init__(self, in_ch, skip_ch, out_ch):
        super().__init__()
        self.up   = nn.ConvTranspose2d(in_ch, in_ch//2, 2, stride=2)
        self.conv = DoubleConv(in_ch//2 + skip_ch, out_ch)
    def forward(self, x, skip):
        x = self.up(x)
        # Align spatial dims (handles odd sizes)
        x = F.interpolate(x, size=skip.shape[-2:])
        x = torch.cat([skip, x], dim=1)  # skip connection
        return self.conv(x)

Having detected and described keypoints (SIFT, ORB, etc.), how do we use correspondences to estimate geometric relationships between images? This is the foundation of image stitching, 3D reconstruction, and augmented reality.

Homography. If two images are related by a pure rotation, or if the scene is planar, corresponding points are related by a homography (projective transformation): $p' \sim H p$ where $H$ is a $3\times3$ matrix (8 degrees of freedom, since scale is free). Given 4 point correspondences, $H$ can be computed exactly via Direct Linear Transform (DLT). With more correspondences, use least squares.

The problem: outliers. SIFT matching produces many outliers (incorrect correspondences) even after Lowe's ratio test. Solving DLT on outlier-contaminated data gives garbage.

RANSAC (Random Sample Consensus, Fischler & Bolles 1981) is the answer. Algorithm:

1. Randomly sample the minimum number of points needed to fit the model (4 for homography).
2. Fit the model to the sample.
3. Count inliers: all correspondences consistent with this model (reprojection error < threshold $\epsilon$).
4. Repeat $N$ times; return the model with most inliers.
5. Refit using all inliers.

How many iterations $N$? With inlier fraction $w$ and sample size $s$, the probability that at least one sample is all-inliers is $1 - (1-w^s)^N$. Set this to 0.99 and solve for $N$: $N = \log(1-0.99)/\log(1-w^s)$. With 50% inliers and $s=4$: $N \approx 72$ iterations.

Fundamental matrix / Essential matrix. For two uncalibrated cameras, corresponding points satisfy the epipolar constraint $p'^T F p = 0$ where $F$ is the $3\times3$ fundamental matrix (7 DOF). For calibrated cameras, the essential matrix $E = K'^T F K$ encodes relative rotation and translation. RANSAC + 8-point algorithm robustly estimates $F$ even with 50% outliers — this is the core of Structure from Motion (SfM).

import cv2
import numpy as np

def panorama_stitch(img1, img2):
    """Stitch two images using SIFT + RANSAC homography."""
    sift = cv2.SIFT_create()
    kp1, d1 = sift.detectAndCompute(img1, None)
    kp2, d2 = sift.detectAndCompute(img2, None)

    bf = cv2.BFMatcher(cv2.NORM_L2)
    matches = [m for m,n in bf.knnMatch(d1,d2,k=2)
               if m.distance < 0.75*n.distance]

    pts1 = np.float32([kp1[m.queryIdx].pt for m in matches])
    pts2 = np.float32([kp2[m.trainIdx].pt for m in matches])

    # RANSAC-robust homography estimation
    H, mask = cv2.findHomography(
        pts1, pts2,
        method=cv2.RANSAC,
        ransacReprojThreshold=4.0  # epsilon in pixels
    )
    n_inliers = mask.sum()
    print(f"{n_inliers}/{len(matches)} inliers")

    h, w = img2.shape[:2]
    warped = cv2.warpPerspective(img1, H, (w*2, h))
    warped[0:h, 0:w] = img2
    return warped, H

Static images give us snapshots; video gives us motion. Optical flow estimates the per-pixel velocity field $(u(x,y), v(x,y))$ that describes how each pixel moved between frames.

The optical flow constraint equation. Assuming brightness constancy ($I(x,y,t) = I(x+u, y+v, t+1)$) and linearizing: $$I_x u + I_y v + I_t = 0$$ where $I_x, I_y$ are spatial gradients and $I_t$ is the temporal gradient. This is one equation in two unknowns — the aperture problem again. We need an additional constraint.

Lucas-Kanade (1981) assumes the flow is constant within a small window $W$. Stacking equations over the window gives an overdetermined system $Ax = b$ where $A$ is the structure tensor $M$ from Harris! The solution $[u, v]^T = (A^T A)^{-1} A^T b$ is reliable only where $M$ has two large eigenvalues — i.e., at corners. Lucas-Kanade is a sparse optical flow method, tracking selected keypoints.

Horn-Schunck (1981) adds a global smoothness penalty $\int\int (|\nabla u|^2 + |\nabla v|^2) dx dy$, producing a dense flow field via iterative relaxation on every pixel.

Deep optical flow. FlowNet (2015) trained a CNN to predict dense flow end-to-end. PWC-Net (2018) uses a pyramid + warping + cost-volume architecture inspired by classical methods but fully learned, achieving real-time dense flow.

Video understanding. Two-stream CNNs (Simonyan & Zisserman, 2014) process RGB frames and optical flow in parallel, fusing for action recognition. 3D convolutions (C3D, I3D) extend conv kernels to the temporal dimension, learning spatio-temporal features. Modern video transformers (Video Swin, 2021) apply shifted-window attention across both space and time.

import cv2
import numpy as np

def lucas_kanade_tracking(prev_gray, curr_gray, prev_pts):
    """
    Track keypoints using Lucas-Kanade pyramidal optical flow.
    prev_pts: Nx1x2 float32 array of points to track.
    """
    lk_params = dict(
        winSize  = (21, 21),
        maxLevel = 3,           # pyramidal (3 scales)
        criteria = (cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_COUNT,
                    30, 0.01)
    )
    curr_pts, status, err = cv2.calcOpticalFlowPyrLK(
        prev_gray, curr_gray, prev_pts, None, **lk_params
    )
    # status[i]=1 means point i was found
    good_prev = prev_pts[status == 1]
    good_curr = curr_pts[status == 1]
    flow_vectors = good_curr - good_prev  # (dx, dy) per point
    return good_curr, good_prev, flow_vectors

Rigorous evaluation separates science from folklore. Every vision task has canonical metrics.

Classification: Top-1 accuracy (fraction of images where the predicted class is correct) and Top-5 accuracy (correct class is in the model's top 5 predictions). ImageNet benchmark uses both.

Detection: Mean Average Precision (mAP). For each class, sort detections by confidence, compute precision-recall curve, compute Average Precision (AP) as the area under the PR curve. mAP is the mean over classes (and, in COCO, over IoU thresholds from 0.5 to 0.95 in steps of 0.05).

Segmentation: Mean IoU (mIoU) per class: $\frac{TP}{TP+FP+FN}$ averaged over classes and images. Also pixel accuracy.

Optical flow: End-Point Error (EPE) — average Euclidean distance between predicted and ground-truth flow vectors.

The road ahead. Vision has moved from detecting edges to understanding scenes at a level competitive with humans:

• Foundation models: CLIP (Radford et al., 2021) trains on 400M image-text pairs with contrastive loss, learning visual concepts aligned with language. Zero-shot classification: embed an image and class names into the same space; pick the closest class. No task-specific fine-tuning needed.
• Segment Anything Model (SAM, 2023): a foundation model for segmentation trained on 1 billion masks. Accepts prompts (points, boxes, text) and segments arbitrary objects zero-shot.
• Diffusion models (DALL-E 2, Stable Diffusion) invert the perception problem: generate photorealistic images from text. They encode visual knowledge at a level that enables stunning controllable generation.
• Open challenges: 3D scene understanding, video reasoning, robustness to distribution shift, and making models that generalize the way humans do with a handful of examples.

def mean_iou(pred_mask, true_mask, num_classes):
    """
    Compute mean IoU over all classes.
    pred_mask, true_mask: H×W integer arrays with class labels.
    """
    iou_per_class = []
    for c in range(num_classes):
        pred_c = (pred_mask == c)
        true_c = (true_mask == c)
        intersection = (pred_c & true_c).sum()
        union        = (pred_c | true_c).sum()
        if union == 0:
            continue  # class absent from both — skip
        iou_per_class.append(intersection / union)
    return sum(iou_per_class) / len(iou_per_class)

# COCO-style AP (simplified AUC of precision-recall)
def average_precision(recalls, precisions):
    """Area under PR curve via trapezoidal rule."""
    import numpy as np
    r = np.concatenate([[0.0], recalls, [1.0]])
    p = np.concatenate([[1.0], precisions, [0.0]])
    # Monotonic decreasing envelope
    for i in range(len(p)-2, -1, -1):
        p[i] = max(p[i], p[i+1])
    return np.trapz(p, r)