Harris Corner detection using OpenCV

HarrisCornerDetection.py

import numpy as np
import cv2
import glob

# Kernel parameters
sigma = 1.4
k_size = (3, 3)

# Corners parameters
r = 1
thickness = 1

alpha = 0.06
threshold = 600000
image_name = "<your_file_name>.png"

def getImageDerivative(I, H_p):
    '''
    Applies the H_p kernel over I and creates a new image I'
    '''
    w, h = I.shape
    I_d = np.zeros((w, h), np.float32)
    # For each pixel apply the H_p kernel
    for rows in range(w):
        for cols in range(h):
            if rows >= 1 or rows <= w-2 and cols >= 1 or cols <= h-2:
                for i in range(3):
                    for j in range(3):
                        I_d[rows][cols] += H_p[i][j] * I[rows - i-1][cols - j-1]
            else:
                I_d[rows][cols] = I[rows][cols]
    return I_d

def HarrisCornerDetection(I):
    '''
    Returns a list of the strongest corners found in the image I
    '''
    ## Prefilter (smooth) of the original image I' = I * H_p
    H_px = np.array([   [-1, 0, 1],
                        [-2, 0, 2],
                        [-1, 0, 1]
                    ])
    H_py = np.array([
                        [-1, -2, -1],
                        [ 0,  0,  0],
                        [ 1,  2,  1]
                    ])

    ### I_x <-- I' * H_px
    I_dX = getImageDerivative(I, H_px)
    ### I_y <-- I' * H_py
    I_dY = getImageDerivative(I, H_py)

    w, h = I.shape

    for i in range(w):
        for j in range(h):
            if  I_dY[i][j] < 0:
                I_dY[i][j] *= -1
                # ImgY[i][j] = 0
            if  I_dX[i][j] < 0:
                I_dX[i][j] *= -1
                # ImgX[i][j] = 0
    
    ## Compute local structure matrix M(u, v) = ( A, C 
    #                                             C, B )
    M_x  = np.square(I_dX)
    M_y  = np.square(I_dY)
    M_xy = np.multiply(I_dX, I_dY)
    M_yx = np.multiply(I_dY, I_dX)

    ## Blur each component of the structure matrix: M* = ( A*, C*
    #                                                      C*, B*)
    M_x  = cv2.GaussianBlur(M_x, k_size, sigma)
    M_y  = cv2.GaussianBlur(M_y, k_size, sigma)
    M_xy = cv2.GaussianBlur(M_xy, k_size, sigma)
    M_yx = cv2.GaussianBlur(M_yx, k_size, sigma)

    ## Find corners strength
    strengths = np.zeros((w, h), np.float32)
    for row in range(w):
        for col in range(h):
            M = np.array([  [M_x[row][col],  M_xy[row][col]],
                            [M_yx[row][col], M_y[row][col]]    ])
            strengths[row][col] = np.linalg.det(M) - (alpha * np.square(np.trace(M)))
    strengths

    corners = []
    ## Delete weak corners
    for row in range(w):
        for col in range(h):
            if strengths[row][col] > threshold:
                max = strengths[row][col]
                skip = False
                for nrow in range(5):
                    for ncol in range(5):
                        if row + nrow - 2 < w and col + ncol - 2 < h:
                            if strengths[row + nrow - 2][col + ncol - 2] > max:
                                skip = True
                                break
                if not skip: # Draw corner point in image
                    corners.append((row, col))
    
    return corners

##### MAIN #####
I_original = cv2.imread(image_name)
cv2.imshow("Original", I_original)
# Original image
I = cv2.cvtColor(I_original, cv2.COLOR_BGR2GRAY)
# Gray scale image
I_gray = cv2.cvtColor(I, cv2.COLOR_GRAY2RGB)
cv2.imshow("Gray image",I_gray)

# Corner detection
corners = HarrisCornerDetection(I)

# Draw corners
for corner in corners:
    cv2.circle(I_gray, (corner[1], corner[0]), r, (0,255,0), thickness)

# Display and save image with corners
cv2.imshow("Image corners", I_gray)
cv2.imwrite("Output_" + image_name + ".png", I_gray)

# End program when windows are closed
cv2.waitKey(0)
cv2.destroyAllWindows()