{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center>\n",
    "<img src=\"../ontario-tech-univ-logo.png\" width=\"22%\">\n",
    "</center>\n",
    "\n",
    "# CSCI 3240U: Computer Vision I\n",
    "\n",
    "## Lab 7: Interest Points and Local Descriptors\n",
    "\n",
    "Faisal Z. Qureshi  \n",
    "Faculty of Science, Ontario Tech University  \n",
    "Oshawa ON Canada  \n",
    "<http://vclab.science.ontariotechu.ca>\n",
    "\n",
    "*Fall 2026*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Copyright information\n",
    "\n",
    "&copy; Faisal Qureshi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## License\n",
    "\n",
    "<a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\"><img alt=\"Creative Commons Licence\" style=\"border-width:0\" src=\"https://i.creativecommons.org/l/by-nc/4.0/88x31.png\" /></a><br />This work is licensed under a <a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\">Creative Commons Attribution-NonCommercial 4.0 International License</a>."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## About this notebook\n",
    "\n",
    "This notebook is the starter for **Part 1 of Lab 7, corner detection**.  Read\n",
    "the handout first.  The sections below follow its numbered tasks.\n",
    "\n",
    "| Section here | Handout |\n",
    "|---|---|\n",
    "| Task 1: Harris corner detector | Part 1, Task 1 |\n",
    "| Non-maximum suppression | Deliverables (no numbered task asks for it, but the Deliverables do) |\n",
    "| Task 2: Shi-Tomasi detector | Part 1, Task 2 |\n",
    "| Task 3: compare the two detectors | Part 1, Task 3 |\n",
    "| Task 4: rotation repeatability | Part 1, Task 4 |\n",
    "| Optional extra: `img_match` | Part 4, optional and **not marked** |\n",
    "\n",
    "The other parts of the lab live elsewhere.  Part 2 (optical flow and tracking)\n",
    "has no starter notebook; write it yourself against\n",
    "`../lab1-setup/traffic-short.mp4`.  Part 3 (the SIFT-like descriptor) has its\n",
    "own starter, `sift-like-descriptor.ipynb`.\n",
    "\n",
    "Cells marked **TODO** are yours to write.  Everything else runs as given and is\n",
    "there to save you plumbing.\n",
    "\n",
    "Submit a single notebook, executed top to bottom, **via Canvas**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "import os\n",
    "\n",
    "import cv2\n",
    "import numpy as np\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "print('OpenCV', cv2.__version__)\n",
    "print('NumPy', np.__version__)\n",
    "print('Matplotlib', matplotlib.__version__)"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Helper code\n",
    "\n",
    "Nothing below is graded.  These helpers load images, draw points, rotate an\n",
    "image by a known angle, move point coordinates between the original and the\n",
    "rotated frame, and plot a repeatability curve.  Read them before you use them.\n",
    "\n",
    "Two conventions matter and they are easy to mix up.\n",
    "\n",
    "- OpenCV point coordinates are `(x, y)`, that is `(column, row)`.\n",
    "- NumPy array indices are `[row, column]`.\n",
    "\n",
    "Every point array in this notebook is `(N, 2)` holding `(x, y)`, in float."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "DATA_FOLDER = '.'\n",
    "\n",
    "\n",
    "def load_rgb(filename):\n",
    "    \"\"\"Read an image and return it as RGB, ready for plt.imshow.\n",
    "\n",
    "    cv2.imread gives you BGR. Passing that straight to plt.imshow swaps red and\n",
    "    blue, which is the classic first bug of every OpenCV lab.\n",
    "    \"\"\"\n",
    "    path = os.path.join(DATA_FOLDER, filename)\n",
    "    bgr = cv2.imread(path)\n",
    "    if bgr is None:\n",
    "        raise FileNotFoundError(path)\n",
    "    return cv2.cvtColor(bgr, cv2.COLOR_BGR2RGB)\n",
    "\n",
    "\n",
    "def load_gray(filename):\n",
    "    \"\"\"Read an image as a single channel uint8 array.\"\"\"\n",
    "    path = os.path.join(DATA_FOLDER, filename)\n",
    "    gray = cv2.imread(path, cv2.IMREAD_GRAYSCALE)\n",
    "    if gray is None:\n",
    "        raise FileNotFoundError(path)\n",
    "    return gray\n",
    "\n",
    "\n",
    "def show(img, title=None, size=(10, 10), cmap=None):\n",
    "    \"\"\"Display one image.  Pass an RGB array, not a BGR one.\"\"\"\n",
    "    plt.figure(figsize=size)\n",
    "    if title:\n",
    "        plt.title(title)\n",
    "    if img.ndim == 2 and cmap is None:\n",
    "        cmap = 'gray'\n",
    "    plt.imshow(img, cmap=cmap)\n",
    "    plt.axis('off')\n",
    "\n",
    "\n",
    "def draw_points(img, pts, colour=(255, 0, 0), radius=4, thickness=-1):\n",
    "    \"\"\"Draw (x, y) points on a copy of an RGB image and return the copy.\"\"\"\n",
    "    out = img.copy()\n",
    "    if pts is None:\n",
    "        return out\n",
    "    for x, y in np.asarray(pts, dtype=float).reshape(-1, 2):\n",
    "        cv2.circle(out, (int(round(x)), int(round(y))), radius, colour, thickness)\n",
    "    return out"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def rotate_image(img, angle_deg):\n",
    "    \"\"\"Rotate about the image centre onto a canvas large enough to hold it all.\n",
    "\n",
    "    Returns the rotated image and the 2x3 affine matrix M that maps a point\n",
    "    in the original image to its location in the rotated image.\n",
    "    \"\"\"\n",
    "    h, w = img.shape[:2]\n",
    "    centre = (w / 2.0, h / 2.0)\n",
    "    M = cv2.getRotationMatrix2D(centre, angle_deg, 1.0)\n",
    "\n",
    "    cos, sin = abs(M[0, 0]), abs(M[0, 1])\n",
    "    new_w = int(round(h * sin + w * cos))\n",
    "    new_h = int(round(h * cos + w * sin))\n",
    "\n",
    "    # Shift so the rotated content lands inside the new canvas.\n",
    "    M[0, 2] += new_w / 2.0 - centre[0]\n",
    "    M[1, 2] += new_h / 2.0 - centre[1]\n",
    "\n",
    "    rotated = cv2.warpAffine(img, M, (new_w, new_h), flags=cv2.INTER_LINEAR)\n",
    "    return rotated, M\n",
    "\n",
    "\n",
    "def transform_points(pts, M):\n",
    "    \"\"\"Apply a 2x3 affine matrix to an (N, 2) array of (x, y) points.\"\"\"\n",
    "    pts = np.asarray(pts, dtype=np.float64).reshape(-1, 2)\n",
    "    if len(pts) == 0:\n",
    "        return pts\n",
    "    homogeneous = np.hstack([pts, np.ones((len(pts), 1))])\n",
    "    return homogeneous @ M.T\n",
    "\n",
    "\n",
    "def inverse_affine(M):\n",
    "    \"\"\"Invert a 2x3 affine matrix.  Use this to map rotated-frame points back.\"\"\"\n",
    "    return cv2.invertAffineTransform(M)\n",
    "\n",
    "\n",
    "def inside(pts, shape, margin=0):\n",
    "    \"\"\"Boolean mask: which (x, y) points fall inside an image of this shape.\"\"\"\n",
    "    pts = np.asarray(pts, dtype=float).reshape(-1, 2)\n",
    "    h, w = shape[:2]\n",
    "    return ((pts[:, 0] >= margin) & (pts[:, 0] < w - margin) &\n",
    "            (pts[:, 1] >= margin) & (pts[:, 1] < h - margin))\n",
    "\n",
    "\n",
    "def plot_repeatability(angles, curves, tol):\n",
    "    \"\"\"Plot repeatability against rotation angle.\n",
    "\n",
    "    curves is a dict mapping a label to a list of rates, one per angle.\n",
    "    \"\"\"\n",
    "    plt.figure(figsize=(7, 4))\n",
    "    for label, rates in curves.items():\n",
    "        plt.plot(angles, rates, marker='o', label=label)\n",
    "    plt.xlabel('rotation (degrees)')\n",
    "    plt.ylabel('repeatability')\n",
    "    plt.title('Repeatability within %g px' % tol)\n",
    "    plt.ylim(0.0, 1.05)\n",
    "    plt.grid(alpha=0.3)\n",
    "    plt.legend()"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Test image\n",
    "\n",
    "`cn-tower-1.jpg` is a reasonable default.  `cb.png` is a synthetic checkerboard\n",
    "and is useful while you are debugging, because you know where its corners are."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "img_file = 'cn-tower-1.jpg'\n",
    "# img_file = 'cb.png'\n",
    "\n",
    "img = load_rgb(img_file)\n",
    "gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)\n",
    "\n",
    "print('image shape', img.shape)\n",
    "show(img, 'Test image')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 1: Harris corner detector\n",
    "\n",
    "Both detectors in this lab are built from the structure tensor\n",
    "\n",
    "$$\n",
    "M = \\sum_{(x,y) \\in W} w(x,y)\n",
    "\\left[\n",
    "\\begin{array}{cc}\n",
    "I_x^2 & I_x I_y \\\\\n",
    "I_x I_y & I_y^2\n",
    "\\end{array}\n",
    "\\right],\n",
    "$$\n",
    "\n",
    "and they differ only in how they turn its eigenvalues into a score.  Harris uses\n",
    "$R = \\det(M) - k\\,\\mathrm{tr}(M)^2$ with $k$ around 0.04 to 0.06.\n",
    "\n",
    "The cell below runs OpenCV's `cornerHarris` and shows the raw response.  Use it\n",
    "as a reference.  You may also build $M$ yourself from the $I_x$ and $I_y$\n",
    "machinery of Lab 4; the two responses should agree up to a scale factor."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "src = np.float32(gray)\n",
    "\n",
    "blocksize = 2   # neighbourhood used to accumulate the structure tensor\n",
    "ksize = 3       # Sobel aperture used for the derivatives\n",
    "k = 0.04        # Harris parameter in R = det(M) - k (trace M)^2\n",
    "\n",
    "harris_response = cv2.cornerHarris(src, blocksize, ksize, k)\n",
    "\n",
    "# A crude threshold on the raw response. Note the blobs: each real corner fires\n",
    "# over a small neighbourhood, so one corner produces many detections. That is\n",
    "# the problem non-maximum suppression solves.\n",
    "crude = img.copy()\n",
    "crude[harris_response > 0.01 * harris_response.max()] = [255, 0, 0]\n",
    "\n",
    "plt.figure(figsize=(20, 15))\n",
    "plt.subplot(2, 1, 1)\n",
    "plt.title('Harris response R')\n",
    "plt.imshow(harris_response, cmap='gray')\n",
    "plt.axis('off')\n",
    "plt.subplot(2, 1, 2)\n",
    "plt.title('Thresholded response, no suppression')\n",
    "plt.imshow(crude)\n",
    "plt.axis('off')\n",
    "\n",
    "print('pixels above threshold:', int((harris_response > 0.01 * harris_response.max()).sum()))"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO (Task 1)\n",
    "\n",
    "Turn the response map into a list of interest points.\n",
    "\n",
    "Write `harris_keypoints`.  It must return an `(N, 2)` float array of `(x, y)`\n",
    "locations, sorted by response so that `max_points` keeps the strongest ones.\n",
    "Threshold the response relative to its maximum, then suppress non-maxima with\n",
    "the `nms` function you write in the next section."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def harris_keypoints(gray_img, max_points=300, rel_threshold=0.01,\n",
    "                     blocksize=2, ksize=3, k=0.04, nms_window=7):\n",
    "    \"\"\"Return an (N, 2) array of Harris interest points as (x, y).\"\"\"\n",
    "    # TODO\n",
    "    #   1. compute the response (cv2.cornerHarris, or your own structure tensor)\n",
    "    #   2. threshold it at rel_threshold * response.max()\n",
    "    #   3. apply non-maximum suppression over an nms_window neighbourhood\n",
    "    #   4. sort the survivors by response, strongest first, and keep max_points\n",
    "    raise NotImplementedError('Task 1: harris_keypoints')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Non-maximum suppression\n",
    "\n",
    "The Deliverables require your detections to be non-maximum suppressed.  No\n",
    "numbered task says so in as many words, so it is easy to miss.  Without it you\n",
    "report a blob of detections per corner, and every count you make afterwards,\n",
    "including the Task 4 repeatability, is meaningless.\n",
    "\n",
    "The idea is one line: keep a pixel only if it is the largest response in its own\n",
    "neighbourhood.  A max filter over the response makes this easy.  `cv2.dilate`\n",
    "with a rectangular structuring element *is* a max filter, so\n",
    "`response == cv2.dilate(response, kernel)` marks the local maxima.  Combine that\n",
    "with your threshold.\n",
    "\n",
    "Note that `cv2.dilate` used on its own, as many online tutorials do, does the\n",
    "opposite of suppression: it spreads each response outwards so the marks are\n",
    "easier to see.  Do not confuse the two uses."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def nms(response, threshold, window=7):\n",
    "    \"\"\"Non-maximum suppression on a response map.\n",
    "\n",
    "    Keep a pixel if its response exceeds `threshold` and is the maximum inside\n",
    "    a window x window neighbourhood centred on it.\n",
    "\n",
    "    Returns an (N, 2) float array of (x, y) locations.\n",
    "    \"\"\"\n",
    "    # TODO\n",
    "    raise NotImplementedError('non-maximum suppression')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "# Check your suppression once nms is written. Both counts print, so you can see\n",
    "# how much the suppression removed. Expect a large reduction.\n",
    "thr = 0.01 * harris_response.max()\n",
    "raw_count = int((harris_response > thr).sum())\n",
    "pts_nms = nms(harris_response, thr, window=7)\n",
    "\n",
    "print('above threshold      :', raw_count)\n",
    "print('after suppression    :', len(pts_nms))\n",
    "\n",
    "plt.figure(figsize=(20, 8))\n",
    "plt.subplot(1, 2, 1)\n",
    "plt.title('No suppression')\n",
    "plt.imshow(crude)\n",
    "plt.axis('off')\n",
    "plt.subplot(1, 2, 2)\n",
    "plt.title('After suppression')\n",
    "plt.imshow(draw_points(img, pts_nms, colour=(255, 0, 0), radius=4))\n",
    "plt.axis('off')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 2: Shi-Tomasi detector\n",
    "\n",
    "Shi and Tomasi score a window by $R = \\min(\\lambda_1, \\lambda_2)$, the smaller\n",
    "eigenvalue of the same $M$.  Part 2 of this lab explains why that particular\n",
    "score earns the title *Good Features to Track*.\n",
    "\n",
    "`cv2.goodFeaturesToTrack` implements it.  Its `minDistance` argument performs\n",
    "the suppression for you, so its output is already sparse.  Keep that in mind\n",
    "when you compare it against Harris: compare like with like, by giving both\n",
    "detectors the same suppression window and the same point budget."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "corners = cv2.goodFeaturesToTrack(gray,\n",
    "                                  maxCorners=25,\n",
    "                                  qualityLevel=0.01,\n",
    "                                  minDistance=10)\n",
    "\n",
    "# np.int0 was removed in NumPy 2.0. Use np.intp, or just keep the floats and\n",
    "# round when you draw. Here we keep floats, which is what the helpers expect.\n",
    "pts_shi_demo = corners.reshape(-1, 2).astype(np.float64)\n",
    "\n",
    "print('Shi-Tomasi points:', len(pts_shi_demo))\n",
    "show(draw_points(img, pts_shi_demo, colour=(0, 255, 0), radius=5),\n",
    "     'Shi-Tomasi corner detector')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO (Task 2)\n",
    "\n",
    "Write `shi_tomasi_keypoints` with the same signature style as\n",
    "`harris_keypoints`, so that Task 3 can call both the same way and compare them\n",
    "fairly.  Return an `(N, 2)` float array of `(x, y)` locations.\n",
    "\n",
    "You may call `cv2.goodFeaturesToTrack`, or compute $\\min(\\lambda_1,\\lambda_2)$\n",
    "yourself with `cv2.cornerMinEigenVal` and reuse your `nms`.  If you use\n",
    "`goodFeaturesToTrack`, set `minDistance` to match your Harris suppression\n",
    "window, otherwise the comparison in Task 3 measures your parameters rather than\n",
    "the detectors."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def shi_tomasi_keypoints(gray_img, max_points=300, rel_threshold=0.01,\n",
    "                         nms_window=7):\n",
    "    \"\"\"Return an (N, 2) array of Shi-Tomasi interest points as (x, y).\"\"\"\n",
    "    # TODO\n",
    "    raise NotImplementedError('Task 2: shi_tomasi_keypoints')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 3: compare the two detectors\n",
    "\n",
    "Run both detectors on the same image with matched settings, draw both point\n",
    "sets on one figure, and say what differs."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "MAX_POINTS = 300\n",
    "NMS_WINDOW = 7\n",
    "REL_THRESHOLD = 0.01\n",
    "\n",
    "pts_harris = harris_keypoints(gray, max_points=MAX_POINTS,\n",
    "                              rel_threshold=REL_THRESHOLD,\n",
    "                              nms_window=NMS_WINDOW)\n",
    "pts_shi = shi_tomasi_keypoints(gray, max_points=MAX_POINTS,\n",
    "                               rel_threshold=REL_THRESHOLD,\n",
    "                               nms_window=NMS_WINDOW)\n",
    "\n",
    "print('Harris     :', len(pts_harris), 'points')\n",
    "print('Shi-Tomasi :', len(pts_shi), 'points')\n",
    "\n",
    "overlay = draw_points(img, pts_harris, colour=(255, 0, 0), radius=6, thickness=1)\n",
    "overlay = draw_points(overlay, pts_shi, colour=(0, 255, 0), radius=3, thickness=-1)\n",
    "show(overlay, 'Harris (red circles) and Shi-Tomasi (green dots)', size=(14, 14))"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO (Task 3)\n",
    "\n",
    "\"They look about the same\" is a fair first impression, and it is not an answer.\n",
    "Quantify it.\n",
    "\n",
    "Write a function that, given the two point sets and a tolerance in pixels,\n",
    "reports what fraction of Harris points have a Shi-Tomasi point within that\n",
    "tolerance, and the fraction the other way round.  Then write a short paragraph\n",
    "in the markdown cell below: where do the two agree, where do they disagree, and\n",
    "what kind of image structure sits at the points only one of them found?"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def overlap_fraction(pts_a, pts_b, tol=3.0):\n",
    "    \"\"\"Fraction of points in pts_a with a point of pts_b within tol pixels.\"\"\"\n",
    "    # TODO\n",
    "    raise NotImplementedError('Task 3: overlap_fraction')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Your Task 3 write-up goes here.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 4: rotation repeatability\n",
    "\n",
    "Rotate the image by a known angle, detect again, and ask whether the *same\n",
    "physical points* come back.\n",
    "\n",
    "The measurement has to be made on **locations**, not on counts.  A detector that\n",
    "finds 300 points before the rotation and 300 after has told you nothing: they\n",
    "could be 300 completely different points.  So:\n",
    "\n",
    "1. detect in the original image, giving `pts_ref`;\n",
    "2. rotate the image by $\\theta$ with `rotate_image`, which also hands you the\n",
    "   affine matrix $M$;\n",
    "3. detect in the rotated image, giving `pts_rot`;\n",
    "4. map `pts_rot` back into the original frame with the inverse of $M$;\n",
    "5. count how many of `pts_ref` have a mapped-back point within a few pixels.\n",
    "\n",
    "Step 4 is pure geometry and is provided.  Step 5 is the measurement and is\n",
    "yours.\n",
    "\n",
    "One detail decides whether your number means anything.  `rotate_image` enlarges\n",
    "the canvas, so the rotated image carries blank corners and a long straight\n",
    "border that was never in the scene.  Detections there are artefacts.  Keep only\n",
    "the mapped-back points that land inside the original image, and be ready to\n",
    "discard a border margin as well.  `inside()` does both."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "ANGLE = 30.0\n",
    "\n",
    "rot_img, M = rotate_image(img, ANGLE)\n",
    "rot_gray = cv2.cvtColor(rot_img, cv2.COLOR_RGB2GRAY)\n",
    "M_inv = inverse_affine(M)\n",
    "\n",
    "print('original shape', img.shape[:2], ' rotated shape', rot_img.shape[:2])\n",
    "\n",
    "# Sanity check on the geometry: send a point out and bring it back.\n",
    "probe = np.array([[100.0, 200.0], [400.0, 50.0]])\n",
    "round_trip = transform_points(transform_points(probe, M), M_inv)\n",
    "print('round trip error (px):', np.abs(round_trip - probe).max())\n",
    "\n",
    "show(rot_img, 'Rotated by %g degrees' % ANGLE, size=(12, 12))"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO (Task 4)\n",
    "\n",
    "Write `repeatability`.  Given the reference points and the mapped-back points,\n",
    "return the fraction of reference points that were recovered.\n",
    "\n",
    "Decide, and state in your write-up, how you handle these:\n",
    "\n",
    "- a mapped-back point should be used at most once, otherwise one lucky\n",
    "  detection can claim several reference points;\n",
    "- a reference point sitting on the image border loses part of its neighbourhood\n",
    "  once the image is rotated, so decide whether to exclude a margin;\n",
    "- what tolerance you use, and why. A few pixels is the usual choice.\n",
    "\n",
    "The driver cell below already drops mapped-back points that fall outside the\n",
    "original image, so `repeatability` receives two point sets in the same frame."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def repeatability(pts_ref, pts_mapped, tol=3.0):\n",
    "    \"\"\"Fraction of pts_ref matched by a distinct point of pts_mapped within tol.\"\"\"\n",
    "    # TODO\n",
    "    raise NotImplementedError('Task 4: repeatability')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "# Driver. This plumbing is given; the number it reports comes from your\n",
    "# repeatability function, so it is only as good as that.\n",
    "TOL = 3.0\n",
    "angles = [0, 10, 20, 30, 45, 60, 90]\n",
    "curves = {'Harris': [], 'Shi-Tomasi': []}\n",
    "\n",
    "detectors = {'Harris': harris_keypoints, 'Shi-Tomasi': shi_tomasi_keypoints}\n",
    "\n",
    "for name, detect in detectors.items():\n",
    "    ref = detect(gray, max_points=MAX_POINTS, rel_threshold=REL_THRESHOLD,\n",
    "                 nms_window=NMS_WINDOW)\n",
    "    for angle in angles:\n",
    "        rotated, Ma = rotate_image(img, angle)\n",
    "        rotated_gray = cv2.cvtColor(rotated, cv2.COLOR_RGB2GRAY)\n",
    "        found = detect(rotated_gray, max_points=MAX_POINTS,\n",
    "                       rel_threshold=REL_THRESHOLD, nms_window=NMS_WINDOW)\n",
    "        back = transform_points(found, inverse_affine(Ma))\n",
    "        back = back[inside(back, img.shape)]   # drop points off the original canvas\n",
    "        curves[name].append(repeatability(ref, back, tol=TOL))\n",
    "    print(name, ['%.2f' % r for r in curves[name]])\n",
    "\n",
    "plot_repeatability(angles, curves, TOL)"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "# Look at the 30 degree case. Reference points in red, points detected in the\n",
    "# rotated image and mapped back in green. Where the two coincide the detector\n",
    "# was repeatable.\n",
    "ref_h = harris_keypoints(gray, max_points=MAX_POINTS,\n",
    "                         rel_threshold=REL_THRESHOLD, nms_window=NMS_WINDOW)\n",
    "rot_h = harris_keypoints(rot_gray, max_points=MAX_POINTS,\n",
    "                         rel_threshold=REL_THRESHOLD, nms_window=NMS_WINDOW)\n",
    "back_h = transform_points(rot_h, M_inv)\n",
    "back_h = back_h[inside(back_h, img.shape)]\n",
    "\n",
    "vis = draw_points(img, ref_h, colour=(255, 0, 0), radius=6, thickness=1)\n",
    "vis = draw_points(vis, back_h, colour=(0, 255, 0), radius=3, thickness=-1)\n",
    "show(vis, 'Original detections (red) and rotated detections mapped back (green)',\n",
    "     size=(14, 14))"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Your Task 4 write-up goes here.**  Report the repeatability at 30 degrees for\n",
    "both detectors, with the tolerance you used.  Then account for the points that\n",
    "were lost.  Interpolation, the discrete detection window and the image boundary\n",
    "all cost you something.  If your number is close to 100 percent, check that you\n",
    "compared locations rather than counts."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Part 4 of the handout: deciding whether two images match\n",
    "\n",
    "**This part is optional and is not marked.**  Everything above is.  Come back to\n",
    "it if you have time, or after Lab 8, which does descriptor matching properly\n",
    "inside RANSAC."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### SIFT keypoints\n",
    "\n",
    "The code below finds locations suitable for SIFT, each with a scale and an\n",
    "orientation.\n",
    "\n",
    "SIFT moved into the main OpenCV module in version 4.4.  Call\n",
    "`cv2.SIFT_create()`.  The old `cv2.xfeatures2d.SIFT_create()` does not exist in\n",
    "`opencv-python` and raises an `AttributeError`."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "sift = cv2.SIFT_create()\n",
    "kp = sift.detect(gray, None)\n",
    "\n",
    "print('SIFT keypoints:', len(kp))\n",
    "\n",
    "# Draw each keypoint with its scale and orientation.\n",
    "sift_vis = cv2.drawKeypoints(gray, kp, None,\n",
    "                             flags=cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)\n",
    "show(sift_vis, 'SIFT key points')"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO (optional, not marked)\n",
    "\n",
    "Use local feature descriptors from two images to compute a match score for the\n",
    "pair.  Here is a recipe.  **You are encouraged not to use the built-in matcher.\n",
    "Experiment with different matching schemes instead.**\n",
    "\n",
    "1. Compute $n$ SIFT descriptors from image 1.\n",
    "2. Compute $m$ SIFT descriptors from image 2.\n",
    "3. Each descriptor is a 128-dimensional vector, so pick a distance between two\n",
    "   vectors, say Euclidean, and a threshold that decides whether two vectors\n",
    "   match.  Using that distance, count the matched vectors between the two sets.\n",
    "4. Report the matched count as a fraction of the smaller of $n$ and $m$.\n",
    "\n",
    "__Complete the following function__, which returns a number between 0.0 and 1.0."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "def img_match(filename1, filename2):\n",
    "\n",
    "    # TO DO\n",
    "\n",
    "    return 0.0"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will use this function as follows."
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "source": [
    "filename1 = 'cn-tower-1.jpg'\n",
    "filename2 = 'cn-tower-2.jpg'\n",
    "\n",
    "print ('Match score is', img_match(filename1, filename2))"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Submission\n",
    "\n",
    "Submit a single Jupyter notebook **via Canvas**.  Run it top to bottom before\n",
    "you submit, so every figure and every number is visible in the file.  Code that\n",
    "has not been executed cannot be marked."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center>\n",
    "    <tr>\n",
    "    <td><img src=\"../ontario-tech-univ-logo.png\" width=\"25%\"></img></td>\n",
    "    </tr>\n",
    "</center>"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
