<H1 style="text-align: center">ECMM426 - Computer Vision / ECMM441 - Machine Vision (Professional)</H1>
<H1 style="text-align: center"></H1>
<H2 style="text-align: center">Workshop 2</H2>
<H2 style="text-align: center">Corner Detection</H2>

## Imports


```
import urllib
import matplotlib, cv2
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = (20.0, 10.0)
plt.rcParams['figure.dpi'] = 72

print('OpenCV version: {}'.format(cv2.__version__))
```

## Download Images
Download chessboard image


```
import os
if not os.path.exists('chessboard.zip'):
  !wget --no-check-certificate https://empslocal.ex.ac.uk/people/staff/ad735/ECMM426/chessboard.zip
  !unzip -q chessboard.zip
```

## Read and show the image


```
# read an image
img = cv2.imread('chessboard.png')

# convert image to grayscale
gray_img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY).astype(np.float32)

# show the image
plt.imshow(gray_img, cmap='gray')
```

## Add border to detect all the corners


```
gray_img = cv2.copyMakeBorder(gray_img.copy(), 10, 10, 10, 10, cv2.BORDER_CONSTANT, value=[0, 0, 0])
# show the image
plt.imshow(gray_img, cmap='gray')
```

## Harris corner detection
We know that corners are regions in the image with large variation in intensity in all the directions. One early attempt to find these corners was done by **Chris Harris & Mike Stephens** in their paper **A Combined Corner and Edge Detector** in 1988, so now it is called the **Harris Corner Detector**. He took this simple idea to a mathematical form. It basically finds the difference in intensity for a displacement of (u,v) in all directions. This is expressed as below:

$$E(u,v) = \sum_{x,y} \underbrace{w(x,y)}_\text{window function} \, [\underbrace{I(x+u,y+v)}_\text{shifted intensity}-\underbrace{I(x,y)}_\text{intensity}]^2$$

The window function is either a rectangular window or a Gaussian window which gives weights to pixels underneath. For sake of simplicity we can consider:

$$E(u,v) = \sum_{x,y} [\underbrace{I(x+u,y+v)}_\text{shifted intensity}-\underbrace{I(x,y)}_\text{intensity}]^2$$

We have to maximize this function $E(u,v)$ for corner detection. 

Taylor's series expansion

$$f(x)=f(a+x-a)=\sum_{n=0}^{\infty} \frac{f^{(n)}(a)}{n!}(x-a)^n=f(a)+\frac{f^{(1)}(a)}{1!}(x-a)+\frac{f^{(2)}(a)}{2!}(x-a)^2+\ldots$$

Applying Taylor's expansion to the above equation

$$\begin{align}
I(x+u,y+v)&=I(x,y)+\frac{\partial I}{\partial x}u+\frac{\partial I}{\partial y}v+\frac{\partial^2 I}{\partial x^2}u^2+\frac{\partial^2 I}{\partial y^2}v^2+\ldots\\
&=I(x,y)+I_xu+I_yv+I_{xx}u^2+I_{yy}v^2+\ldots
\end{align}$$

Therefore,

$$E(u,v) \approx \sum_{x,y}(I_xu+I_yv)^2 \approx \sum_{x,y}I_xI_xu^2+2\sum_{x,y}I_xI_yuv+\sum_{x,y}I_yI_yv^2$$

 and using some mathematical steps (please refer to any standard text books you like for full derivation), we get the final equation as:

$$E(u,v) \approx \begin{bmatrix} u & v \end{bmatrix} M \begin{bmatrix} u \\ v \end{bmatrix}$$

where

$$M = \sum_{x,y} \begin{bmatrix}I_x I_x & I_x I_y \\ I_x I_y & I_y I_y \end{bmatrix}$$

Here, $I_x$ and $I_y$ are image derivatives in x and y directions respectively.

## Image gradient


```
# construct the Sobel x-axis kernel
sobelX = np.array(([-1, 0, 1], [-2, 0, 2], [-1, 0, 1]), dtype=np.float)
# construct the Sobel y-axis kernel
sobelY = np.array(([-1, -2, -1], [0, 0, 0], [1, 2, 1]), dtype=np.float)
# convolution
I_x = cv2.filter2D(gray_img, -1, sobelX)
I_y = cv2.filter2D(gray_img, -1, sobelY)
# plot
plt.subplot(1,3,1),plt.imshow(gray_img, cmap = 'gray')
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1,3,2),plt.imshow(I_x, cmap = 'gray')
plt.title('$I_x$'), plt.xticks([]), plt.yticks([])
plt.subplot(1,3,3),plt.imshow(I_y, cmap = 'gray')
plt.title('$I_y$'), plt.xticks([]), plt.yticks([])
plt.show()
```

## Gradient covariances


```
# Gradient covariance and light smoothing
I_x_I_x = cv2.GaussianBlur(I_x*I_x, (3, 3), sigmaX=1, sigmaY=1)
I_y_I_y = cv2.GaussianBlur(I_y*I_y, (3, 3), sigmaX=1, sigmaY=1)
I_x_I_y = cv2.GaussianBlur(I_x*I_y, (3, 3), sigmaX=1, sigmaY=1)

# plot
plt.subplot(1,4,1),plt.imshow(gray_img, cmap = 'gray')
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1,4,2),plt.imshow(I_x_I_x, cmap = 'gray')
plt.title('$I_xI_x$'), plt.xticks([]), plt.yticks([])
plt.subplot(1,4,3),plt.imshow(I_y_I_y, cmap = 'gray')
plt.title('$I_yI_y$'), plt.xticks([]), plt.yticks([])
plt.subplot(1,4,4),plt.imshow(I_x_I_y, cmap = 'gray')
plt.title('$I_xI_y$'), plt.xticks([]), plt.yticks([])
plt.show()
```

## Harris response calculation

Then comes the main part. After this, they created a score, basically an equation, which determines if a window can contain a corner or not.

$$R = \text{det}(M) - k(\text{trace}(M))^2$$

where
* $\text{det}(M) = \lambda_1 \lambda_2$
* $\text{trace}(M) = \lambda_1 + \lambda_2$
* $\lambda_1$ and $\lambda_2$ are the eigenvalues of $M$


```
k = 0.05
# determinant
detA = I_x_I_x * I_y_I_y - I_x_I_y ** 2
# trace
traceA = I_x_I_x + I_y_I_y
R = detA - k * traceA ** 2
```

## Find corners using R > 0


```
y, x = np.where(R > 0)
plt.figure()
plt.imshow(gray_img, cmap='gray')
plt.scatter(x, y, c='r', s=40)
plt.show()
```

##Find edges using R < 0


```
y, x = np.where(R < 0)
plt.figure()
plt.imshow(gray_img, cmap='gray')
plt.scatter(x, y, c='r', s=10)
plt.show()
```

## Find flat regions using R = 0


```
y, x = np.where(R == 0)
plt.figure()
plt.imshow(gray_img, cmap='gray')
plt.scatter(x, y, c='r', s=1)
plt.show()
```

## Non-maximum suppression


```
from skimage.feature import corner_peaks
# non-maximal suppression
corners = corner_peaks(R)
x, y = corners[:, 1], corners[:, 0]
plt.figure()
plt.imshow(gray_img, cmap='gray')
plt.scatter(x, y, c='r', s=40)
plt.show()
```

## `HarrisCorner` function


```
from skimage.feature import corner_peaks
def HarrisCorner(img, ksize=5, sigmaX=2, sigmaY=2, k=0.05, non_max_suppr=True):
  # convert image to grayscale
  gray_img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY).astype(np.float32)
  # light Gaussian smoothing
  gray_img = cv2.GaussianBlur(gray_img, (ksize, ksize), sigmaX=sigmaX, sigmaY=sigmaY)
  # construct the Sobel x-axis kernel
  sobelX = np.array(([-1, 0, 1], [-2, 0, 2], [-1, 0, 1]), dtype=np.float)
  # construct the Sobel y-axis kernel
  sobelY = np.array(([-1, -2, -1], [0, 0, 0], [1, 2, 1]), dtype=np.float)
  # convolution
  I_x = cv2.filter2D(gray_img, -1, sobelX)
  I_y = cv2.filter2D(gray_img, -1, sobelY)
  # gradient covariances and light Gaussian smoothing
  I_x_I_x = cv2.GaussianBlur(I_x*I_x, (ksize, ksize), sigmaX=sigmaX, sigmaY=sigmaY)
  I_y_I_y = cv2.GaussianBlur(I_y*I_y, (ksize, ksize), sigmaX=sigmaX, sigmaY=sigmaY)
  I_x_I_y = cv2.GaussianBlur(I_x*I_y, (ksize, ksize), sigmaX=sigmaX, sigmaY=sigmaY)
  # determinant
  detA = I_x_I_x * I_y_I_y - I_x_I_y ** 2
  # trace
  traceA = I_x_I_x + I_y_I_y
  # Harris criteria
  R = detA - k * traceA ** 2
  if non_max_suppr:
    corners = corner_peaks(R)
    x, y = corners[:, 1], corners[:, 0]
  else:
    y, x = np.where(R > 0)
  return (x, y)
```

## Download some other images


```
import os
if not os.path.exists('images.zip'):
  !wget --no-check-certificate https://empslocal.ex.ac.uk/people/staff/ad735/ECMM426/images.zip
  !unzip -q images.zip
```

## Apply on the other images


```
img = cv2.imread('images/books.jpg')
x, y = HarrisCorner(img, sigmaX=3, sigmaY=3, k=0.1)
x_, y_ = HarrisCorner(img, sigmaX=3, sigmaY=3, k=0.1, non_max_suppr=False)
# plot
plt.subplot(1, 3, 1), plt.imshow(img)
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 2), plt.imshow(img)
plt.scatter(x, y, c='r', s=40)
plt.title('Corners with nonmax suppression'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 3), plt.imshow(img)
plt.scatter(x_, y_, c='r', s=40)
plt.title('Corners without nonmax suppression'), plt.xticks([]), plt.yticks([])
plt.show()
```


```
img = cv2.imread('images/piet.png')
x, y = HarrisCorner(img, sigmaX=3, sigmaY=3, k=0.1)
x_, y_ = HarrisCorner(img, sigmaX=3, sigmaY=3, k=0.1, non_max_suppr=False)
# plot
plt.subplot(1, 3, 1), plt.imshow(img)
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 2), plt.imshow(img)
plt.scatter(x, y, c='r', s=40)
plt.title('Corners with nonmax suppression'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 3), plt.imshow(img)
plt.scatter(x_, y_, c='r', s=40)
plt.title('Corners without nonmax suppression'), plt.xticks([]), plt.yticks([])
plt.show()
```

## Qualitative comparison with OpenCV implementation
Note the difference between the plotting of our implementation and OpenCV implementation


```
img = cv2.imread('images/books.jpg')
img_cv = img.copy()
x, y = HarrisCorner(img, ksize=3, sigmaX=3, sigmaY=3, k=0.01)
dst = cv2.cornerHarris(np.float32(cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)), 3, 3, 0.1)
# Results are marked through the dilated corners 
dst = cv2.dilate(dst, None, iterations = 3) 
# Reverting back to the original image, with optimal threshold value 
img_cv[dst > 0.01 * dst.max()] = [255, 0, 0]
# plot
plt.subplot(1, 3, 1), plt.imshow(img)
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 2), plt.imshow(img)
plt.scatter(x, y, c='r', s=40)
plt.title('Corners (our implementation)'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 3), plt.imshow(img_cv)
plt.title('Corners (opencv implementation)'), plt.xticks([]), plt.yticks([])
plt.show()
```


```
img = cv2.imread('images/piet.png')
img_cv = img.copy()
x, y = HarrisCorner(img, ksize=3, sigmaX=3, sigmaY=3, k=0.1)
dst = cv2.cornerHarris(np.float32(cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)), 2, 3, 0.1)
# Results are marked through the dilated corners 
dst = cv2.dilate(dst, None, iterations = 10) 
# Reverting back to the original image, with optimal threshold value 
img_cv[dst > 0.01 * dst.max()] = [255, 0, 0]
# plot
plt.subplot(1, 3, 1), plt.imshow(img)
plt.title('Original'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 2), plt.imshow(img)
plt.scatter(x, y, c='r', s=40)
plt.title('Corners (our implementation)'), plt.xticks([]), plt.yticks([])
plt.subplot(1, 3, 3), plt.imshow(img_cv)
plt.title('Corners (opencv implementation)'), plt.xticks([]), plt.yticks([])
plt.show()
```

bf05f805740684f100c899e5816a7285a1de5118
ipynb
Jupyter Notebook
Notebooks/ECMM426_ECMM441_Corner_Detection (1).ipynb
Aarif1430/monthly-challenges
23990bb3e8930eeb9dc4308c73b43b77a6a3ab49
MIT
2021-02-03T19:19:27.000Z
2021-02-03T19:19:27.000Z
Notebooks/ECMM426_ECMM441_Corner_Detection (1).ipynb
Aarif1430/monthly-challenges
23990bb3e8930eeb9dc4308c73b43b77a6a3ab49
MIT
Notebooks/ECMM426_ECMM441_Corner_Detection (1).ipynb
Aarif1430/monthly-challenges
23990bb3e8930eeb9dc4308c73b43b77a6a3ab49
MIT
Qwen/Qwen-72B
1. YES
2. YES