Revision 109e5a3f551db7f7aada4123c975a3aced2a0b40 (click the page title to view the current version)
Changes from 109e5a3f551db7f7aada4123c975a3aced2a0b40 to current
---
title: Corner Detection
categories: session
---
# Briefing
## Corners and Feature Points
**Reading** Ma 4.3 and 4.A and Szelisky Chapter 3.2
## Differentiation
**Warning** Ma starts Chapter 4 by discussing *tracking*,
which means that motion as a function of time is considered as well
as the image as a function of spatial co-ordinates.
This is a lot of concepts and quantities to process at the same time.
We will instead start by discussing features in a still image.
When we have a good idea of what features are and how they behave,
we shall introduce motion.
## Harris Feature Detector
**Briefing** [Corner Lecture]()
# Exercises
## Setup
## Learning Objectives
We will use scipy for a small part of the exercises,
if you haven't already, run 'pip install scipy'.
1. What makes a feature in visual terms?
2. What makes a feature in mathematical terms?
3. How do we differentiate a sampled signal?
4. How does the Harris corner detector work?
We need an image to work with, you can either load an image from disk or
capture a new image from the webcam with the below code.
We use the same setup as in [Image Filters]().
```python
# Imports
import cv2 as cv
import numpy as np
## Exercise 1
# Capture a single frame
vid = cv.VideoCapture(0)
image = None
Learning goal: 1D derivatives and 1D convolutions
# Capture frames until we click the space button, use that image
while True:
_, frame = vid.read()
**Part 1**<br>
Extract one row from the grayscale image, and visualize it (e.g. with matplotlib.pyplot).<br>
Does the values correspond to what you would expect from the row?
# Display the resulting frame
cv.imshow("frame", frame)
# Click space to capture frame
k = cv.waitKey(1)
if k % 256 == 32:
# SPACE pressed
image = frame
break
# After the loop release the cap object
vid.release()
# Destroy all the windows
cv.destroyAllWindows()
```
As we will be working with a lot of different images in this exercise,
it is recommended to save it to disk, we can do that with
```python
# Imports
from pathlib import Path
# First we create a path to an images directory
p = Path.cwd() / "images" # <--- current working directory + /images
if not p.is_dir():
p.mkdir()
# Then we save the image to the directory with name "frame.jpg"
cv.imwrite(str(p / "frame.jpg"), image)
```
Now we convert the image to gray-scale
<details>
<summary>Hint (Click to expand)</summary>
</br>
```python
# Convert frame to grayscale
image_gray = cv.cvtColor(image, cv.COLOR_BGR2GRAY)
# Load image im as usual and convert to greyscale.
row = im[50,:] # Row number 50
n = row.size
# Save gray to images
cv.imwrite(str((p / "gray.jpg")), image_gray)
import matplotlib.pyplot as plt
plt.plot(range(n),row,color="blue",label="Pixel Row")
plt.show()
```
## Exercise 1, Sobel/Derivative filter
The first exercise is to implement a Sobel-filter.
Recall from the theory that we need to implement two 3x3 kernels to convolve with the original image.
This can be done using scipy.signal.convolve2d ( https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.convolve2d.html )
(Note: For a larger challenge you can also try implementing your own algorithm for the convolve function using numpy)
</details>
Code answers from here on will be collapsed, we recommend that you try to
implement them yourself before reading an answer.
**Part 2**
Convolve the row with a $[1/2,-1/2]$ kernel and visualize.<br>
Do the values make sense?
<details>
<summary>Hint (Click to expand)</summary>
<summary>Manual with Numpy (Click to expand)</summary>
</br>
Use `scipy.signal.convolve2d(<image>, <filter>, boundary='symm', mode='same')`
</details>
To implement convolution manually using a loop over a numpy array, you can
create a new 1D array with `np.zeros(<shape>)` and fill it in by
iterating over the row and the kernel.
Remember that the resulting array should be smaller than the original.
</details>
<details>
<summary>Solution (Click to expand)</summary>
<summary>Using the API (Click to expand)</summary>
</br>
If you are not able to get a result using numpy, use scipy.signal
(can also be used to compare your result):
```python
from scipy import signal
sobel_x = np.array([[-1, 0, 1],
[-2, 0, 2],
[-1, 0, 1]])
sobel_y = np.array([[-1, -2, -1],
[0, 0, 0],
[1, 2, 1]])
grad_x = signal.convolve2d(image_gray, sobel_x, boundary='symm', mode='same')
grad_y = signal.convolve2d(image_gray, sobel_y, boundary='symm', mode='same')
row_d = signal.convolve(row, kernel)
# or
# row_d = signal.correlate(row, cv.flip(kernel, -1))
```
</details>
You should then show the images using cv.imshow or save using cv.imwrite, as we did earlier.
Discuss the results.
You can compare the results of your implementation with the built in function
`cv.Sobel`.
```python
grad_x_cv = cv.Sobel(image_gray, cv.CV_32F, 1, 0, ksize=3) # gradient along x axis,
grad_y_cv = cv.Sobel(image_gray, cv.CV_32F, 0, 1, ksize=3) # gradient along y axis,
```
Using `cv.Sobel`, adjust the kernel-size to see how that changes the result.
<details>
<summary>Cross-correlation (Click to expand)</summary>
</br>
If you use cross-correlation instead of convolution, flip the kernel.
</details>
Compute the magnitude and orientation of the derivatives using
```python
magnitude = np.sqrt((grad_x ** 2) + (grad_y ** 2))
orientation = np.arctan2(grad_y, grad_x) * (180 / np.pi) % 180
```
## Exercise 2
Show/Save the images, and discuss.
Learning goals: 2D Derivatives
## Exercise 2, Eigenvalues and eigenvectors
**Part 1**<br>
Apply the sobel operator
$$G_x = 1/8 \begin{bmatrix}
1 & 0 & -1 \\
2 & 0 & -2 \\
1 & 0 & -1 \\
\end{bmatrix}
$$
to the entire grayscale image
using `scipy.signal.convolve2d`, `scipy.signal.correlate2d` or `cv.filter2d`.
Either show the image or write it to file with `cv.imwrite`.
This should give you the derivative $I_x$ of the image $I$ with
respect to $x$.
TODO: ?
<details>
<summary>Hint</summary>
</br>
1. Look at the exercises from last Friday. We apply the Sobel filter just like a blurring
filter.
2. You can define the matrix $G_x$ as given above using the `numpy.array()` function.
</details>
As this is not a trivial step, we will calculate the EigenVals and Vecs using
`cv.cornerEigenValsAndVecs`, note that this function also applies a sobel-filter
so we will call this function using the grayscale image.
+ How does this relate to the 1D derivative you did in exercise 1.2?
+ What are the minimum and maximum values of the $I_x$ matrix?
```python
bsize = 3 # the size of neighbourhood considered for corner detection
ksize = 3 # Aperture parameter of the Sobel derivative used.
eigenv = cv.cornerEigenValsAndVecs(image_gray, bsize, ksize)
```
**Part 2**<br>
Show $I_x$ as an image.
You probably have negative pixel values, so you may have to
scale the image.
The cornerEigenValsAndVecs algorithm will
+ Try to take the absolute values of the luminence values.
+ Try to scale the luminences into the $0\ldots255$ range,
e.g. by adding $255$ and dividing by two.
+ What does the different visualisations tell you?
+ You may scale further to use the full $0\ldots255$ range and thus
increas contrast.
1. Calculate the derivatives $dI/dx$ and $dI/dy$ using the sobel filter (as we did in exercise 1)
2. For each pixel $p$, take a neighborhood $S(p)$ of blockSize `bsize*bsize`,
Calculate the covariation matrix $M$ of the derivatives over the neighborhood as:
**Part 3**<br>
Repeat Parts 1 and 2 with vertical derivation, i.e. use $G_y$ instead
of $G_x$.
$$M = \begin{bmatrix}
\sum_{S(p)}(dI/dx)^2 & \sum_{S(p)}(dI/dx)(dI/dy) \\
\sum_{S(p)}(dI/dx)(dI/dy) & \sum_{S(p)}(dI/dy)^2
$$G_y = 1/8 \begin{bmatrix}
1 & 2 & 1 \\
0 & 0 & 0 \\
-1 & -2 & -1 \\
\end{bmatrix}
$$
After this the eigenvectors and eigenvalues of $M$ are calculated and returned.
+ Compare the images. What differences can you make out?
This should result in a $h*w*6$ array with values ($\lambda_1,\lambda_2,x_1,y_1,x_2,y_2$) where
## Exercise 3
* $\lambda_1,\lambda_2$ are the non-sorted eigenvalues of $M$
* $x_1,y_1$ are the eigenvectors corresponding to $\lambda_1$
* $x_2,y_2$ are the eigenvectors corresponding to $\lambda_2$
Learning Goal: find rotation invariant heuristics for edges
Analyze? TODO: How can students learn from this exercise?
In Exercise 2, we calculated $I_x$ and $I_y$ which give a lot of
edge information. Now we want to aggregate this information over
a Window.
## Exercise 3, Harris-criterion
Note that $I_x$ and $I_y$ are matrices with the same dimensions as
the original image. The index of an entry in these matrices will
be denoted $\mathbf{x}$ below, and we are going to make more
*pseudo-images* with the same dimensions.
Recall from the theory that we can calculate the harris-criterion with:
$$\lambda_1 * \lambda_2 - k * (\lambda_1 + \lambda_2)^2$$
### 3.1
Create a new image of the same size as image_gray, and calculate the
harris criterion for all pixels.
For every point $\mathbf{x}$ we calculate the matrix
$$G(\mathbf{x}) = \begin{bmatrix}
\sum I_x^2 & \sum I_xI_y \\
\sum I_xI_y & \sum I_y^2
\end{bmatrix},$$
where the summations are made over a window, say a $5\times5$ window,
around $\mathbf{x}$.
<details>
<summary>Solution (Click to expand)</summary>
Note that this is not a *pseudo-image*. For each $\mathbf{x}$ we
have $2\times2$ matrix and not just a scalar.
A simple (naive) solution, is to create an empty image and
iterate over all pixels like this:
### 3.2
(Note that this is a very computationally expensive solution,
and takes a long time to complete)
For each pixel position $\mathbf{x}$ calculate the eigenvectors/-values
of $G(\mathbf{x})$.
We know that large eigenvalues indicate features, and we want
to visualise this information.
Several variants are possible, and you may need only one or two
to get the picture.
```python
# Create an empty image of size image_gray.shape
c_x = np.zeros(image_gray.shape, dtype=np.float32)
You can make matrices containing, for each $\mathbf{x}$
# Define k, e.g. 0.04
k = 0.04 # Harris free parameter
+ The maximum of the eigenvalues of $G(\mathbf{x})$
+ The minimum of the eigenvalues of $G(\mathbf{x})$
+ The sum of the eigenvalues of $G(\mathbf{x})$
+ The product of the eigenvalues of $G(\mathbf{x})$
# Iterate over image and set values
for i in range(c_x.shape[0]):
for j in range(c_x.shape[1]):
lambda_1 = eigenv[i, j, 0]
lambda_2 = eigenv[i, j, 1]
c_x[i, j] = lambda_1 * lambda_2 - k * np.power((lambda_1 + lambda_2), 2)
```
1. Scale these matrices so that they can be interpreted as
grey scale images, and visualise them.
2. Compare these visulisations to the edge plots from Exercise 2.
3. What do you see?
Another, faster solution:
## Exercise 4
```python
# Define k, e.g. 0.04
k = 0.04 # Harris free parameter
Learning goals: Introduction to Harris corner detector
# Extract eigen-values
lambdas_1 = eigenv[:,:,0]
lambdas_2 = eigenv[:,:,1]
If you have completed Exercise 3, you have done 90% of the
implementation of the Harris detector.
The missing part is the heuristics
$$C(G)=\lambda_{1}\lambda_{2} - k\cdot (\lambda _{1}+\lambda _{2})^{2}
=\det(G)-k\cdot \mathrm {tr} (G)^{2}$$
and the threshold used.
You can choose if you want to complete the implementation
of your very own Harris detector or use OpenCV's implementation.
# Compute harris-criterion
c_x = lambdas_1 * lambdas_2 - k * np.power((lambdas_1 + lambdas_2), 2)
```
</details>
Take a look at some of the values you calculated above and compare them to the
original image (specifically where you expect the algorithm to detect corners).
You can slice the image with the `c_x[:,:]` command, the syntax is
[starting_row(including):ending_row(not including),starting_column(including):ending_column(not including)]
, e.g. if you want to look at the top left corner around ~10% from the edge on a
480x480 image, you can use `c_x[40:60, 40:60]` or for the same but the bottom right
corner, use e.g. `c_x[-60:-40,-60:-40]`.
## Exercise 4, Simple corner-detector
**Part 1**<br>
For the last exercise for corner-detection, you should draw a circle around all the
corners you found on the original image/frame.
Consider the grayscale image we have been working with, and the gradient magnitude from 3.3.<br>
Where do you expect to find corners?
Start by picking a threshold T based on the information you found in the last exercise
(or by using `T = min_cx + (max_cx - min_cx) * 50 / 100)` ). Now draw a circle around
all pixels exceeding this value, using `cv.circle`.
Apply opencv's built in harris-detector.<br>
E.g. `cv.cornerHarris(img_gray, block_size, kernel_size, k)` with `block_size = 2`, `kernel_size = 5` and `k = 0.06`.<br>
Here, block_size is the size of neighbourhood considered for corner detection, kernel_size is the size of the sobel derivative kernel, while k is the harris free parameter.
Adjust the threshold T and try with different values to see how this changes the detection.
Make a copy of the original image (with colors) and make circles around any corners found by the harris-detector.<br>
Example code for drawing circles is added below.
<details>
<summary>Solution (Click to expand)</summary>
A simple solution, is to iterate over all pixels and check compare the value against
the threshold.
<summary>Code</summary>
</br>
```python
min_cx = np.min(c_x)
max_cx = np.max(c_x)
T = min_cx + (max_cx - min_cx) * 50 / 100)
f_circles = frame
cx = cv.cornerHarris(img_gray, bsize, ksize, k)
# Drawing a circle around corners (naive, slow)
for i in range(f_circles.shape[0]):
for j in range(f_circles.shape[1]):
if c_x[i, j] > T:
cv.circle(f_circles, (j, i), 5, (0, 0, 255), 2)
T = 0.1 # Threshold
c_image = img
for i in range(c_image.shape[0]):
for j in range(c_image.shape[1]):
if c_x[i, j] > T:
cv.circle(c_image, (j, i), 2, (0, 0, 255), 2)
```
</details>
Another solution is to create a list of all pixels above the threshold before
iterating over them:
Save/visualize the result, how does it compare with your expectations?
```python
min_cx = np.min(c_x)
max_cx = np.max(c_x)
T = min_cx + (max_cx - min_cx) * 50 / 100)
f_circles = frame
**Part 2**<br>
Adjust the threshold `T` when drawing circles, what does this do?
mask = c_x > T
corners = np.argwhere(mask).tolist()
**Part 3**<br>
Adjust the kernel_size (must be positive and odd),
block_size and/or $k$, and observe how they change the result.
for corner in corners:
cv.circle(f_circles, tuple(corner), 5, (0, 0, 255), 2)
```
## Optional Exercises
</details>
### Laplacian of Gaussian
Now do the same to an the original image based on the result from `cv.cornerHarris` and compare with your own result.
Test different variations of LoG filters on your test images.
You can construct Gaussians as we did in the last session on
[Image Filters]() and the Laplacian as suggested in the
[Corner Lecture]() today.
Using `cv.cornerHarris`, adjust block_size, aperture_size and the harris parameter (k),
to see how they affect the result.
How does LoG compare to the techniques you ghave tested above?
## Exercise 5, SIFT
### Building on Exercise 1
## Exercise 6, Feature matching
Repeat part 2 with an image column (instead of row) and the transpose of the kernel.
# Debrief
### Building on Exercise 2
Using the same method as in exercise 1.2 and 1.3, compute the gradient of all
rows and columns, and compute the magnitude (as in exercise 3.3), compare with the
magnitude from 3.3.
