harris角點檢測算法實作

算法流程：

        1、将圖像轉換為灰階圖像；

        2、利用Sobel濾波器求出 海森矩陣 (Hessian matrix) :

        3、将高斯濾波器分别作用于Ix²、Iy²、IxIy；

        4、計算每個像素的 R= det(H) - k(trace(H))²。det(H)表示矩陣H的行列式，trace表示矩陣H的迹。通常k的取值範圍為[0.04,0.16]；

        5、滿足 R>=max(R) * th 的像素點即為角點。th常取0.1。

python代碼實作：

#!/usr/bin/env python3
# coding: utf-8

import cv2
import copy
import numpy as np

def bgr2gray(img):
    gray = 0.2126 * img[..., 2] + 0.7152 * img[..., 1] + 0.0722 * img[..., 0]
    gray = gray.astype(np.uint8)
    return gray

def sobel_filtering(gray):
    # get shape
    img_h, img_w = gray.shape

    # sobel kernel
    sobel_y = np.array(((1, 2, 1),
                        (0, 0, 0),
                        (-1, -2, -1)), dtype=np.float32)

    sobel_x = np.array(((1, 0, -1),
                        (2, 0, -2),
                        (1, 0, -1)), dtype=np.float32)

    # padding
    tmp = np.pad(gray, (1, 1), 'edge')

    # prepare
    ix = np.zeros_like(gray, dtype=np.float32)
    iy = np.zeros_like(gray, dtype=np.float32)

    # get differential
    for y in range(img_h):
        for x in range(img_w):
            ix[y, x] = np.mean(tmp[y: y + 3, x: x + 3] * sobel_x)
            iy[y, x] = np.mean(tmp[y: y + 3, x: x + 3] * sobel_y)

    ix2 = ix ** 2
    iy2 = iy ** 2
    ixy = ix * iy
    return ix2, iy2, ixy

def gaussian_filtering(I, k_size=3, sigma=3):
    # get shape
    img_h, img_w = I.shape

    # gaussian
    i_t = np.pad(I, (k_size // 2, k_size // 2), 'edge')

    # gaussian kernel
    K = np.zeros((k_size, k_size), dtype=np.float32)
    for x in range(k_size):
        for y in range(k_size):
            _x = x - k_size // 2
            _y = y - k_size // 2
            K[y, x] = np.exp(-(_x ** 2 + _y ** 2) / (2 * (sigma ** 2)))
    K /= (sigma * np.sqrt(2 * np.pi))
    K /= K.sum()

    # filtering
    for y in range(img_h):
        for x in range(img_w):
            I[y, x] = np.sum(i_t[y: y + k_size, x: x + k_size] * K)
    return I

def corner_detect(img, ix2, iy2, ixy, k=0.04, th=0.1):
    # prepare output image
    out = copy.deepcopy(img)

    # get R
    R = (ix2 * iy2 - ixy ** 2) - k * ((ix2 + iy2) ** 2)

    # detect corner
    out[R >= np.max(R) * th] = [255, 0, 0]
    out = out.astype(np.uint8)
    return out

def harris_corner(img):
    # 1. grayscale
    gray = bgr2gray(img)

    # 2. get difference image
    ix2, iy2, ixy = sobel_filtering(gray)

    # 3. gaussian filtering
    ix2 = gaussian_filtering(ix2, k_size=3, sigma=3)
    iy2 = gaussian_filtering(iy2, k_size=3, sigma=3)
    ixy = gaussian_filtering(ixy, k_size=3, sigma=3)

    # 4. corner detect
    out = corner_detect(img, ix2, iy2, ixy)
    return out

def main():
    # Read image
    img = cv2.imread("test.jpg")
    img = cv2.resize(img, (512, 512))
    img = img.astype(np.float32)

    # Harris corner detection
    out = harris_corner(img)
    cv2.imwrite("out.jpg", out)
    print("proc ok.")

if __name__ == "__main__":
    main()
           

結果：

 opencv python版本實作：

        對于每一個像素(x,y)，在(blockSize*blockSize)鄰域内，計算梯度圖的協方差矩陣M(x,y) ，然後通過上面第二步中的角點響應函數得到結果圖。圖像中的角點可以為該結果圖的局部最大值。即可以得到輸出圖中的局部最大值，這些值就對應圖像中的角點。

        Harris 角點檢測函數：

cv2.cornerHarris(src, blockSize, ksize, k, dst, borderType)
           

        src：資料類型為float32的輸入圖像

        blockSize：角點檢測中要考慮的領域大小

        ksize：Sobel求導中使用的視窗大小

        k：Harris 角點檢測方程中的自由參數，取值參數為 [0,04，0.06].

        dst：目标圖像

        borderType：邊界類型

import cv2
import numpy as np

def test():
    img = cv2.imread("test.jpg")
    img = cv2.resize(img, (512, 512))
    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    gray = np.float32(gray)
    # 輸入圖像必須是 float32，最後一個參數在 0.04 到 0.05 之間
    dst = cv2.cornerHarris(gray, 3, 3, 0.04)
    # result is dilated for marking the corners, not important
    # dst = cv2.dilate(dst, None)
    # Threshold for an optimal value, it may vary depending on the image.
    img[dst > 0.1 * dst.max()] = [255, 0, 0]
    cv2.imwrite("out2.jpg", img)
    print("proc ok.")
           

結果：

相關連結：

1、harris角點檢測原理

2、SHI-TOMASI角點檢測

參考文章：

1、python實作Harris角點檢測算法（本文代碼來源）

https://www.jb51.net/article/201957.htm

2、opencv實作版本python

https://blog.csdn.net/yukinoai/article/details/88759615

3、opencv c++版本實作

https://blog.csdn.net/zhu_hongji/article/details/81235643