Moravec corner detection

Moravec corner detection is one of the earliest corner detection algorithm.

We will cover the intuition, mathematics formulation, and implementation of the Moravec corner detection.

Let's first import the necessary libraries.

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
from scipy import signal

We will consider the following image

img = np.array([
              [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 100, 100, 100, 100, 100, 100, 100, 100],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 100, 100, 100, 100, 100, 100, 100, 100, 100],
              [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
              [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], 
])
plt.figure(figsize=(7, 6))
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Let's create the patch image of size 3 * 3

plt.figure(figsize=(7, 6))
rect = Rectangle((6,5),2,2,linewidth=2,edgecolor='r',facecolor='none')
plt.gca().add_patch(rect)
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

We can shift the current patch image by 1 pixel into 8 directions.

When we shift the patch image, we can calculate the sum of the square distance between the patch image and shifted patch image as a difference between the two images.

plt.figure(figsize=(7, 6))
rect = Rectangle((6,5),2,2,linewidth=2,edgecolor='r',facecolor='none')
plt.gca().add_patch(rect)
rect = Rectangle((7,6),2,2,linewidth=2,edgecolor='b',facecolor='none')
plt.gca().add_patch(rect)
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Three cases need to be considered

Case I: Patch Image in flat region
All shifts will result in only a small change.

plt.figure(figsize=(7, 6))
rect = Rectangle((6,5),2,2,linewidth=2,edgecolor='r',facecolor='none')
plt.gca().add_patch(rect)
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Case II: Patch image in Edge region
Shift along the edge will result in a small change, but a shift perpendicular to the edge will result in a large change.

plt.figure(figsize=(7, 6))
rect = Rectangle((2,5),2,2,linewidth=2,edgecolor='r',facecolor='none')
plt.gca().add_patch(rect)
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Case III: Patch image in Corner region
All shifts will result in a large change

plt.figure(figsize=(7, 6))
rect = Rectangle((2,1),2,2,linewidth=2,edgecolor='r',facecolor='none')
plt.gca().add_patch(rect)
plt.imshow(img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Moravec corner detection algorithm

• A corner is detected when the minimum change produced by any of the shifts is larger than a certain threshold.
• Here the shifts comprise {(1, 0), (1, 1), (0, 1), (-1, 1)}

For every pixel location

In the above equation, p = 2k + 1 and q = 2l + 1

If the V is above a certain threshold, the corner is detected in that pixel location.

window_size = 3
offset = window_size // 2
res_img = np.zeros((img.shape))
for img_y in range(offset, img.shape[0]-offset):
    for img_x in range(offset,img.shape[1]-offset):
val = [0, 0, 0, 0]
for j in range(-offset, offset -1 + 1):
            for i in range(-offset, offset + 1):
                val[0] += (img[img_y + j, img_x + i] - img[img_y + j + 1, img_x + i])**2
        
        p = 2 * (offset + 1) + 1
        q = 2 * (offset -1 + 1) + 1
        val[0] = val[0] / (p * (q - 1))
for j in range(-offset, offset+ 1):
            for i in range(-offset, offset -1 + 1):
                val[1] += (img[img_y + j, img_x + i] - img[img_y + j, img_x + i + 1])**2
p = 2 * (offset -1 + 1) + 1
        q = 2 * (offset + 1) + 1
        val[1] = val[1] / ((p -1)* q)
for j in range(-offset, offset -1 + 1):
            for i in range(-offset, offset -1 + 1):
                val[2] += (img[img_y + j, img_x + i] - img[img_y + j + 1, img_x + i + 1])**2
p = 2 * (offset -1 + 1) + 1
        q = 2 * (offset -1 + 1) + 1
        val[2] = val[2] / ((p -1) * (q - 1))
for j in range(-offset, offset -1 + 1):
            for i in range(-offset, offset -1 + 1):
                val[3] += (img[img_y + j + 1, img_x + i] - img[img_y + j, img_x + i + 1])**2
p = 2 * (offset -1 + 1) + 1
        q = 2 * (offset -1 + 1) + 1
        val[3] = val[3] / ((p -1) * (q - 1))
min_val = min(val)
if min_val > 2000:
            res_img[img_y, img_x] = min_val
plt.figure(figsize=(7, 6))
plt.imshow(res_img, cmap='gray')
plt.xticks(list(range(0, img.shape[1])))
plt.yticks(list(range(0, img.shape[0])))
plt.grid()
plt.show()

Corner point is detected at (3, 2) and (3, 12)

Also, we can represent the change E produced by a shift (u, v) as

where w is the image window whose value is 1 within a specified rectangular region, and 0 elsewhere

References:

1. https://www.ece.lsu.edu/gunturk/EE7700/A_Combined_Corner_and_Edge_Detector.pdf