{ "cells": [ { "cell_type": "markdown", "id": "2e2a9a77", "metadata": {}, "source": [ "# Lesson 12: Differentiation and Edge Detection\n", "\n", "An edge is, informally, a place where intensity changes quickly. That's a statement about a *derivative*. This lesson builds up discrete image derivatives — Prewitt, Sobel, and Scharr operators — as convolution kernels (Lesson 10), uses them inside the classic **Canny** edge detector, and finishes by turning the resulting edge maps into actual shapes with the **Hough transform**." ] }, { "cell_type": "code", "execution_count": null, "id": "7ab986e0", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import cv2\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "id": "91a8e89f", "metadata": {}, "source": [ "## From derivatives to gradients\n", "\n", "Treat the image as a function $I(x,y)$. Its gradient $\\nabla I = (I_x, I_y)$ points in the direction of steepest intensity increase, where $I_x = \\partial I/\\partial x$ and $I_y = \\partial I/\\partial y$ are the **partial derivatives** of image intensity with respect to the horizontal and vertical directions — how fast the intensity changes as you move one pixel right, or one pixel down. Two numbers summarize the gradient at each pixel:\n", "\n", "- **magnitude** $\\|\\nabla I\\| = \\sqrt{I_x^2 + I_y^2}$ — how sharp the change is (large at edges, near zero on flat regions)\n", "- **direction** $\\theta = \\mathrm{atan2}(I_y, I_x)$ — which way intensity is increasing fastest (perpendicular to the edge)\n", "\n", "$I_x$ and $I_y$ are themselves computed by convolving the image with small derivative kernels." ] }, { "cell_type": "code", "execution_count": null, "id": "ae46f2f3", "metadata": {}, "outputs": [], "source": [ "def make_test_image(size=200):\n", " img = np.zeros((size, size), dtype=np.uint8)\n", " cv2.rectangle(img, (30, 30), (100, 100), 200, -1)\n", " cv2.circle(img, (140, 140), 40, 150, -1)\n", " pts = np.array([[150, 30], [190, 90], [110, 90]], dtype=np.int32)\n", " cv2.fillPoly(img, [pts], 25) # low-contrast triangle, for the threshold demo later\n", " return img\n", "\n", "img = make_test_image()\n", "plt.imshow(img, cmap='gray')\n", "plt.title('Test image (edges at several orientations)')\n", "plt.axis('off')\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "829560b3", "metadata": {}, "source": [ "## The naive derivative: just a difference\n", "\n", "The simplest $I_x$ estimate is a 1D central-difference kernel $[-1,\\ 0,\\ 1]$. It works, but with no smoothing in the perpendicular direction, it's maximally sensitive to noise." ] }, { "cell_type": "code", "execution_count": null, "id": "a49de564", "metadata": {}, "outputs": [], "source": [ "simple_x = np.array([[-1, 0, 1]], dtype=np.float64)\n", "simple_y = simple_x.T\n", "\n", "gx = cv2.filter2D(img.astype(np.float64), cv2.CV_64F, simple_x)\n", "gy = cv2.filter2D(img.astype(np.float64), cv2.CV_64F, simple_y)\n", "magnitude = np.sqrt(gx**2 + gy**2)\n", "\n", "plt.imshow(magnitude, cmap='gray')\n", "plt.title('Gradient magnitude, simple [-1, 0, 1] kernel')\n", "plt.axis('off')\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "db062970", "metadata": {}, "source": [ "## Prewitt, Sobel, and Scharr: derivative + smoothing\n", "\n", "Real operators combine a difference in one direction with a *smoothing* average in the perpendicular direction — this is what makes them robust to noise. They differ only in how they weight that smoothing:\n", "\n", "| Operator | $x$-kernel | Perpendicular weighting |\n", "|---|---|---|\n", "| **Prewitt** | $\\begin{bmatrix}-1&0&1\\\\-1&0&1\\\\-1&0&1\\end{bmatrix}$ | uniform: 1, 1, 1 |\n", "| **Sobel** | $\\begin{bmatrix}-1&0&1\\\\-2&0&2\\\\-1&0&1\\end{bmatrix}$ | binomial: 1, 2, 1 |\n", "| **Scharr** | $\\begin{bmatrix}-3&0&3\\\\-10&0&10\\\\-3&0&3\\end{bmatrix}$ | 3, 10, 3 |\n", "\n", "Each $y$-kernel is just the transpose of its $x$-kernel. Sobel's binomial weights approximate a small Gaussian, which is why it's the most commonly used default. Scharr's weights were chosen (numerically optimized) specifically to make the *estimated gradient direction* as rotationally accurate as possible — i.e. as close to a true continuous derivative as an integer $3\\times3$ kernel can get — which matters most in applications like optical flow that depend on precise gradient angles rather than just edge location." ] }, { "cell_type": "code", "execution_count": null, "id": "baf2ceb3", "metadata": {}, "outputs": [], "source": [ "prewitt_x = np.array([[-1, 0, 1]] * 3, dtype=np.float64)\n", "prewitt_y = prewitt_x.T\n", "\n", "operators = {\n", " 'Simple': (simple_x, simple_y),\n", " 'Prewitt': (prewitt_x, prewitt_y),\n", " 'Sobel': None, # use cv2.Sobel directly below\n", " 'Scharr': None, # use cv2.Scharr directly below\n", "}\n", "\n", "img_f = img.astype(np.float64)\n", "magnitudes = {}\n", "magnitudes['Simple'] = magnitude\n", "magnitudes['Prewitt'] = np.sqrt(cv2.filter2D(img_f, cv2.CV_64F, prewitt_x)**2 +\n", " cv2.filter2D(img_f, cv2.CV_64F, prewitt_y)**2)\n", "magnitudes['Sobel'] = np.sqrt(cv2.Sobel(img_f, cv2.CV_64F, 1, 0, ksize=3)**2 +\n", " cv2.Sobel(img_f, cv2.CV_64F, 0, 1, ksize=3)**2)\n", "magnitudes['Scharr'] = np.sqrt(cv2.Scharr(img_f, cv2.CV_64F, 1, 0)**2 +\n", " cv2.Scharr(img_f, cv2.CV_64F, 0, 1)**2)\n", "\n", "fig, axes = plt.subplots(1, 4, figsize=(13, 3.5))\n", "for ax, (name, mag) in zip(axes, magnitudes.items()):\n", " ax.imshow(mag, cmap='gray')\n", " ax.set_title(name, fontsize=10)\n", " ax.axis('off')\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "adc4fe66", "metadata": {}, "source": [ "On a clean, low-noise image like this, the four results look nearly identical — the practical differences between them are subtle and show up mainly under noise or when precise sub-pixel angle matters, not as a visibly different edge map." ] }, { "cell_type": "markdown", "id": "62eaf1f5", "metadata": {}, "source": [ "### Noise robustness: why the perpendicular smoothing matters\n", "\n", "To make the benefit of that perpendicular smoothing concrete, we measure each kernel's response to pure noise on a flat (edge-free) region. To compare fairly, we first rescale each kernel so it has the *same gain* on an ideal step edge (i.e. divide by the sum of its positive weights) — otherwise a kernel with larger raw coefficients would look noisier just from its overall scale, not its shape." ] }, { "cell_type": "code", "execution_count": null, "id": "f6fc8f29", "metadata": {}, "outputs": [], "source": [ "rng = np.random.default_rng(0)\n", "flat_noisy = np.full((150, 150), 128.0) + rng.normal(0, 15, (150, 150))\n", "\n", "kernels = {\n", " 'Simple [-1,0,1]': simple_x,\n", " 'Prewitt': prewitt_x,\n", " 'Sobel': np.array([[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]], dtype=np.float64),\n", " 'Scharr': np.array([[-3, 0, 3], [-10, 0, 10], [-3, 0, 3]], dtype=np.float64),\n", "}\n", "\n", "print(f'{\"kernel\":>18} {\"noise std (gain-normalized)\":>30}')\n", "for name, k in kernels.items():\n", " gain = k[k > 0].sum()\n", " response = cv2.filter2D(flat_noisy, cv2.CV_64F, k / gain)\n", " print(f'{name:>18} {response.std():>30.2f}')" ] }, { "cell_type": "markdown", "id": "fbc4627b", "metadata": {}, "source": [ "The plain difference kernel is noticeably noisier than any of the $3\\times3$ operators — averaging over 3 rows (or columns) while differencing cuts down the noise response substantially, which is exactly the point of building smoothing into the derivative kernel." ] }, { "cell_type": "markdown", "id": "b8293efa", "metadata": {}, "source": [ "## From gradients to edges: the Canny detector\n", "\n", "Simply thresholding gradient magnitude gives thick, noisy edge blobs. The **Canny** edge detector (Canny, 1986★) refines this with a multi-stage pipeline:\n", "\n", "1. **Smooth** with a Gaussian, to suppress noise before differentiating.\n", "2. **Compute gradients** (Sobel, internally) to get magnitude and direction at every pixel.\n", "3. **Non-maximum suppression**: at each pixel, keep the gradient magnitude only if it's a local maximum *along the gradient direction* — this thins wide gradient ridges down to single-pixel-wide lines.\n", "4. **Double thresholding + hysteresis**: pixels above a high threshold are definite edges; pixels below a low threshold are discarded; pixels in between are kept only if they connect to a definite edge. This links up weak-but-real edge segments while suppressing isolated noise responses." ] }, { "cell_type": "code", "execution_count": null, "id": "1237d94f", "metadata": {}, "outputs": [], "source": [ "calvin = cv2.imread('../img/calvin.png', cv2.IMREAD_GRAYSCALE)\n", "\n", "calvin_sobel_mag = np.sqrt(cv2.Sobel(calvin.astype(np.float64), cv2.CV_64F, 1, 0, ksize=3)**2 +\n", " cv2.Sobel(calvin.astype(np.float64), cv2.CV_64F, 0, 1, ksize=3)**2)\n", "naive_edges = (calvin_sobel_mag > 150).astype(np.uint8) * 255 # crude: just threshold |gradient|\n", "canny_edges = cv2.Canny(calvin, threshold1=80, threshold2=160)\n", "\n", "fig, axes = plt.subplots(1, 3, figsize=(10, 3.5))\n", "for ax, im, title in zip(axes, [calvin, naive_edges, canny_edges],\n", " ['Original', 'Naive: threshold |gradient|\\n(thick, blobby)', 'Canny\\n(thin, connected)']):\n", " ax.imshow(im, cmap='gray')\n", " ax.set_title(title, fontsize=9)\n", " ax.axis('off')\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "bdb9c90d", "metadata": {}, "source": "

Image source: National Gallery of Art

" }, { "cell_type": "markdown", "id": "b3f03aa5", "metadata": {}, "source": [ "### Threshold sensitivity\n", "\n", "Canny's two thresholds trade off completeness against noise. Too low, and noise gets picked up as spurious edges; too high, and real (but low-contrast) edges are missed." ] }, { "cell_type": "code", "execution_count": null, "id": "62999371", "metadata": {}, "outputs": [], "source": [ "noisy_calvin = np.clip(calvin.astype(np.float64) + rng.normal(0, 8, calvin.shape), 0, 255).astype(np.uint8)\n", "\n", "threshold_pairs = [(20, 60), (80, 160), (150, 220)]\n", "\n", "fig, axes = plt.subplots(1, len(threshold_pairs), figsize=(10, 3.5))\n", "for ax, (lo, hi) in zip(axes, threshold_pairs):\n", " edges = cv2.Canny(noisy_calvin, lo, hi)\n", " ax.imshow(edges, cmap='gray')\n", " ax.set_title(f'thresholds=({lo}, {hi})', fontsize=9)\n", " ax.axis('off')\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "6ce7639d", "metadata": {}, "source": [ "With low thresholds on the noisy image, speckle noise survives as spurious edge fragments scattered everywhere, especially in the flat background. At the highest thresholds, faint detail (the window panes, the fine hatching that shades the curtain and robe) disappears completely, while the figure's strong dark outline and the books' edges survive — a real example of Canny's threshold trading off noise rejection against sensitivity to faint-but-genuine edges." ] }, { "cell_type": "markdown", "id": "4db89191", "metadata": {}, "source": [ "## From edges to lines: the Hough transform\n", "\n", "Edge detection finds *where* intensity changes sharply, but it doesn't know that a scattered set of edge pixels forms a straight line. The **Hough transform** (Hough, 1962; Duda & Hart, 1972) answers a different question: given a set of edge points, which ones are consistent with lying on a common line?\n", "\n", "Parameterize a line not as $y=mx+b$ (which blows up for vertical lines) but by its distance from the origin and the angle of its normal:\n", "\n", "$$\\rho = x\\cos\\theta + y\\sin\\theta$$\n", "\n", "For a *fixed* edge point $(x,y)$, this equation traces out a sinusoidal curve in $(\\rho,\\theta)$ space as $\\theta$ varies — every $(\\rho,\\theta)$ pair on that curve describes a line through $(x,y)$. The key idea: if several edge points are **colinear**, their sinusoids all cross at the *same* $(\\rho,\\theta)$ — the parameters of the line they share." ] }, { "cell_type": "code", "execution_count": null, "id": "eb34e8fc", "metadata": {}, "outputs": [], "source": [ "line_points = [(30, 270), (150, 150), (270, 30)] # three points on the same line, x + y = 300\n", "\n", "# the line's true (rho, theta): its normal direction is (1,1)/sqrt(2), i.e. theta=45 degrees,\n", "# and rho is the line's distance from the origin along that normal\n", "theta_true = np.radians(45)\n", "rho_true = 300 * np.cos(theta_true)\n", "\n", "thetas = np.linspace(0, np.pi, 500)\n", "plt.figure(figsize=(6, 3.5))\n", "for x, y in line_points:\n", " rho = x * np.cos(thetas) + y * np.sin(thetas)\n", " plt.plot(np.degrees(thetas), rho, label=f'({x},{y})')\n", "plt.axvline(45, color='gray', linestyle='--', linewidth=1)\n", "plt.scatter([np.degrees(theta_true)], [rho_true], color='black', zorder=5, marker='x', s=80,\n", " label=f'true line: rho={rho_true:.1f}, theta=45°')\n", "plt.xlabel('theta (degrees)')\n", "plt.ylabel('rho')\n", "plt.legend(fontsize=8)\n", "plt.title('Each point traces a sinusoid; colinear points cross at one (rho, theta)')\n", "plt.show()\n", "\n", "# confirm all three curves genuinely pass through that point, not just visually\n", "for x, y in line_points:\n", " rho_at_true_theta = x * np.cos(theta_true) + y * np.sin(theta_true)\n", " print(f'point ({x:>3},{y:>3}): rho at theta=45 deg = {rho_at_true_theta:.2f} (true rho = {rho_true:.2f})')" ] }, { "cell_type": "markdown", "id": "3035bb40", "metadata": {}, "source": [ "### The accumulator: voting for lines\n", "\n", "In practice, every edge pixel votes for its entire sinusoid of possible $(\\rho,\\theta)$ lines, into a discretized 2D **accumulator** array. The bin with the most votes is the line most consistent with the data. We build this from scratch for a single synthetic line, and check it against `cv2.HoughLines`." ] }, { "cell_type": "code", "execution_count": null, "id": "acd5e9fd", "metadata": {}, "outputs": [], "source": [ "line_img = np.zeros((300, 300), dtype=np.uint8)\n", "cv2.line(line_img, (30, 270), (270, 30), 255, 3) # the same x+y=300 line as above\n", "\n", "ys, xs = np.nonzero(line_img)\n", "theta_bins = np.linspace(0, np.pi, 180, endpoint=False)\n", "max_rho = int(np.hypot(*line_img.shape))\n", "accumulator = np.zeros((2 * max_rho, len(theta_bins)), dtype=np.int32)\n", "\n", "for x, y in zip(xs, ys):\n", " rho_vals = (x * np.cos(theta_bins) + y * np.sin(theta_bins)).astype(int) + max_rho\n", " accumulator[rho_vals, np.arange(len(theta_bins))] += 1\n", "\n", "peak_rho_idx, peak_theta_idx = np.unravel_index(np.argmax(accumulator), accumulator.shape)\n", "manual_rho, manual_theta = peak_rho_idx - max_rho, theta_bins[peak_theta_idx]\n", "\n", "cv_lines = cv2.HoughLines(line_img, 1, np.pi / 180, threshold=100)\n", "cv_rho, cv_theta = cv_lines[0, 0]\n", "\n", "print(f'true line parameters: rho={rho_true:.1f}, theta={np.degrees(theta_true):.1f} deg')\n", "print(f'manual accumulator peak: rho={manual_rho}, theta={np.degrees(manual_theta):.1f} deg')\n", "print(f'cv2.HoughLines result: rho={cv_rho:.0f}, theta={np.degrees(cv_theta):.1f} deg')\n", "\n", "# recover two far-apart points on the detected line, to draw it back on the image\n", "a, b = np.cos(cv_theta), np.sin(cv_theta)\n", "x0, y0 = a * cv_rho, b * cv_rho\n", "pt1 = (int(x0 + 1000 * (-b)), int(y0 + 1000 * a))\n", "pt2 = (int(x0 - 1000 * (-b)), int(y0 - 1000 * a))\n", "overlay = cv2.cvtColor(line_img, cv2.COLOR_GRAY2BGR)\n", "cv2.line(overlay, pt1, pt2, (0, 0, 255), 3)\n", "\n", "fig, axes = plt.subplots(1, 3, figsize=(12, 3.5))\n", "axes[0].imshow(line_img, cmap='gray')\n", "axes[0].set_title('Input: single line', fontsize=10)\n", "axes[0].axis('off')\n", "\n", "# most bins get only 1-2 votes from individual points' sinusoids, while the peak bin gets one\n", "# vote per point on the line (~240 here) -- on a linear color scale that peak washes out\n", "# everything else, so we log-compress the counts to make the fainter sinusoid trails visible too\n", "axes[1].imshow(np.log1p(accumulator), cmap='hot', aspect='auto',\n", " extent=[0, 180, -max_rho, max_rho], origin='lower')\n", "axes[1].scatter([np.degrees(theta_true)], [rho_true], edgecolor='cyan', facecolor='none', s=150,\n", " linewidth=1.5, label='true (rho, theta)')\n", "axes[1].set_xlabel('theta (degrees)')\n", "axes[1].set_ylabel('rho')\n", "axes[1].legend(fontsize=7)\n", "axes[1].set_title('Accumulator (log-scaled)', fontsize=10)\n", "\n", "axes[2].imshow(cv2.cvtColor(overlay, cv2.COLOR_BGR2RGB))\n", "axes[2].set_title('Recovered line (red),\\noverlaid on input', fontsize=10)\n", "axes[2].axis('off')\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "9f72dd09", "metadata": {}, "source": [ "As you can see, the basic Hough algorithm returns an infinite line (red output) rather than a line segment (white input)." ] }, { "cell_type": "markdown", "id": "8958832e", "metadata": {}, "source": [ "### `cv2.HoughLinesP`: the practical version\n", "\n", "The probabilistic Hough transform additionally returns line *segment* endpoints (not just infinite lines), and is what's normally used in practice. We test it on a scene with two lines of known orientation, run through an actual Canny edge detector first (Hough transforms operate on edge maps, not raw images)." ] }, { "cell_type": "code", "execution_count": null, "id": "1239643f", "metadata": {}, "outputs": [], "source": [ "shapes_with_lines = np.zeros((300, 300), dtype=np.uint8)\n", "cv2.line(shapes_with_lines, (30, 270), (270, 30), 200, 2) # true angle: -45 degrees\n", "cv2.line(shapes_with_lines, (50, 50), (250, 50), 200, 2) # true angle: 0 degrees\n", "cv2.circle(shapes_with_lines, (220, 220), 30, 150, -1) # a distractor shape, no straight edges\n", "\n", "edges = cv2.Canny(shapes_with_lines, 80, 160)\n", "segments = cv2.HoughLinesP(edges, 1, np.pi / 180, threshold=60, minLineLength=80, maxLineGap=5)\n", "\n", "# color each segment by which of the two true lines (-45 or 0 degrees) it's closer to\n", "# (colors given in BGR, since `vis` gets converted to RGB before display below)\n", "color_diagonal, color_horizontal = (255, 0, 0), (80, 80, 255) # red, blue\n", "\n", "vis = cv2.cvtColor(shapes_with_lines, cv2.COLOR_GRAY2BGR)\n", "print(f'{\"segment\":>20} {\"angle (deg)\":>12}')\n", "for x1, y1, x2, y2 in segments.reshape(-1, 4):\n", " angle = np.degrees(np.arctan2(y2 - y1, x2 - x1))\n", " print(f'({x1},{y1})-({x2},{y2})'.rjust(20), f'{angle:>12.1f}')\n", " color = color_diagonal if abs(angle - (-45)) < abs(angle - 0) else color_horizontal\n", " cv2.line(vis, (x1, y1), (x2, y2), color, 2)\n", "\n", "fig, axes = plt.subplots(1, 3, figsize=(9, 3.5))\n", "for ax, im, title in zip(axes, [shapes_with_lines, edges, vis],\n", " ['Original', 'Canny edges', 'Detected segments\\n(red = -45°, blue = 0°)']):\n", " ax.imshow(im, cmap='gray' if im.ndim == 2 else None)\n", " ax.set_title(title, fontsize=9)\n", " ax.axis('off')\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "3492a104", "metadata": {}, "source": [ "Every detected segment lands almost exactly on one of the two true line angles (-45 and 0 degrees) — each line typically yields two or three near-duplicate detections rather than one, because a drawn line has finite width, so Canny finds an edge on *each* side of the stroke, a few pixels apart. Real Hough-based pipelines usually merge nearby detections in a post-processing step for exactly this reason. The circle — which has no straight edges at all — correctly produces no line detections whatsoever.\n", "\n", "### Beyond lines: Hough circles\n", "\n", "The same voting idea generalizes to any shape describable by a small number of parameters. A circle needs 3 (center $x$, center $y$, radius $r$), so `cv2.HoughCircles` accumulates votes in a 3D parameter space instead of 2D." ] }, { "cell_type": "code", "execution_count": null, "id": "bc3f6c64", "metadata": {}, "outputs": [], "source": [ "circle_img = np.zeros((200, 200), dtype=np.uint8)\n", "true_center, true_radius = (100, 100), 50\n", "cv2.circle(circle_img, true_center, true_radius, 255, 1)\n", "\n", "detected = cv2.HoughCircles(circle_img, cv2.HOUGH_GRADIENT, dp=1, minDist=50,\n", " param1=50, param2=20, minRadius=30, maxRadius=70)\n", "cx, cy, r = detected[0, 0]\n", "\n", "print(f'true circle: center={true_center}, radius={true_radius}')\n", "print(f'detected circle: center=({cx:.1f}, {cy:.1f}), radius={r:.1f}')\n", "\n", "vis = cv2.cvtColor(circle_img, cv2.COLOR_GRAY2BGR)\n", "cv2.circle(vis, (int(cx), int(cy)), int(r), (0, 0, 255), 1)\n", "cv2.drawMarker(vis, (int(cx), int(cy)), (0, 0, 255), cv2.MARKER_CROSS, 10)\n", "plt.imshow(cv2.cvtColor(vis, cv2.COLOR_BGR2RGB))\n", "plt.title('Detected circle (red) overlaid on the true one')\n", "plt.axis('off')\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "f21078d2", "metadata": {}, "source": [ "### Peeking inside: the circle accumulator\n", "\n", "`cv2.HoughCircles` doesn't expose its internal accumulator, but building one from scratch is the same voting idea as the line accumulator above, with one more dimension: a circle needs a center *and* a radius, so each edge pixel votes for every $(c_x, c_y, r)$ triple it could be consistent with, for every candidate radius $r$." ] }, { "cell_type": "code", "execution_count": null, "id": "abe9451a", "metadata": {}, "outputs": [], "source": [ "circle_edges = cv2.Canny(circle_img, 50, 100)\n", "edge_ys, edge_xs = np.nonzero(circle_edges)\n", "\n", "H, W = circle_img.shape\n", "radii = np.arange(30, 71)\n", "vote_thetas = np.linspace(0, 2 * np.pi, 72, endpoint=False) # angle around each candidate circle\n", "cos_t, sin_t = np.cos(vote_thetas), np.sin(vote_thetas)\n", "\n", "circle_accumulator = np.zeros((len(radii), H, W), dtype=np.int32)\n", "for ri, r in enumerate(radii):\n", " # for this radius, every edge point votes for every center that's r away from it\n", " cand_cx = np.round(edge_xs[:, None] - r * cos_t[None, :]).astype(int)\n", " cand_cy = np.round(edge_ys[:, None] - r * sin_t[None, :]).astype(int)\n", " valid = (cand_cx >= 0) & (cand_cx < W) & (cand_cy >= 0) & (cand_cy < H)\n", " np.add.at(circle_accumulator[ri], (cand_cy[valid], cand_cx[valid]), 1)\n", "\n", "peak_r_idx, peak_cy, peak_cx = np.unravel_index(np.argmax(circle_accumulator), circle_accumulator.shape)\n", "print(f'true circle: center={true_center}, radius={true_radius}')\n", "print(f'accumulator peak: center=({peak_cx}, {peak_cy}), radius={radii[peak_r_idx]}')\n", "\n", "radius_vs_cx = circle_accumulator.max(axis=1) # max over cy\n", "radius_vs_cy = circle_accumulator.max(axis=2) # max over cx\n", "\n", "fig, axes = plt.subplots(1, 2, figsize=(9, 4))\n", "axes[0].imshow(radius_vs_cx, aspect='auto', cmap='hot', extent=[0, W, radii[-1], radii[0]])\n", "axes[0].set_xlabel('cx'); axes[0].set_ylabel('radius')\n", "axes[0].set_title('radius vs. cx (max over cy)', fontsize=10)\n", "axes[1].imshow(radius_vs_cy, aspect='auto', cmap='hot', extent=[0, H, radii[-1], radii[0]])\n", "axes[1].set_xlabel('cy'); axes[1].set_ylabel('radius')\n", "axes[1].set_title('radius vs. cy (max over cx)', fontsize=10)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "d03ee0cc", "metadata": {}, "source": [ "To avoid the difficulty of visualizing a 3D array, we instead take 2D slices: the maximum vote count over $c_y$ (giving radius vs. $c_x$) and over $c_x$ (giving radius vs. $c_y$). Each slice shows a bright X: every edge point's votes trace two diagonal lines as the candidate radius grows (one for centers to its left, one to its right), and all of them cross at the true center and radius, which is exactly the brightest point in each plot." ] }, { "cell_type": "markdown", "id": "05e078d4", "metadata": {}, "source": [ "### Exercise\n", "\n", "1. Increase the noise standard deviation in `flat_noisy` and re-run the noise-robustness comparison. Does the relative ordering of the four kernels change?\n", "2. `cv2.Canny` accepts an `apertureSize` parameter for its internal Sobel step (default 3). Try `apertureSize=5` and compare the result to the default on the noisy image.\n", "3. Canny's hysteresis step needs a *connected* path of above-low-threshold pixels between a weak edge and a strong one. Construct a small binary example (by hand, as a NumPy array) where a real edge is broken into two segments with a 1-pixel gap, and confirm that hysteresis fails to link them even though a human viewer would clearly see one edge.\n", "4. Lower `cv2.HoughLinesP`'s `threshold` argument well below 60 on `shapes_with_lines`. Do you start getting spurious short segments from noise in the circle's Canny boundary, even though it has no straight edges? What does raising `minLineLength` do to fix this?\n", "5. Add a third line to `shapes_with_lines` at a shallow angle (e.g. from `(30, 200)` to `(270, 150)`), and predict its angle before checking `HoughLinesP`'s output against your prediction.\n", "6. `cv2.HoughCircles`'s `param2` argument is roughly an accumulator vote threshold (like `HoughLines`'s `threshold`) — lower values report more (and more spurious) circles. Sweep `param2` from 30 down to 10 on a noisy version of `circle_img` (add Gaussian noise first) and observe when false-positive circles start appearing." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "name": "python", "version": "3.x" } }, "nbformat": 4, "nbformat_minor": 5 }