Opencv Hessian角点检测

对thorino.jpg进行Hessian角点检测吧！

角点检测是检测边缘上的角点。

角点是曲率变大的点，下式定义了高斯曲率：
K=\frac{\det(H)}{(1+{I_x}^2+{I_y}^2)^2}

其中：

• $$\det(H)=I_{xx}\ I_{yy}-{I_{xy}}^2$$；
• $$H$$表示Hessian矩阵。图像的二次微分（通过将Sobel滤波器应用于灰度图像计算得来）。对于图像上的一点，按照下式定义：
  • $$I_x$$：应用$$x$$方向上的Sobel滤波器；
  • $$I_y$$：应用$$y$$方向上的Sobel滤波器；
  • $$H=\left[\begin{matrix}I_{xx}&I_{xy}\\I_{xy}&I_{yy}\end{matrix}\right]$$

在Hessian角点检测中，\det{H}将极大点视为j角点。

如果中心像素与其8-近邻像素相比值最大，则中心像素为极大点。

解答中，角点是\det(H)为极大值，并且大于\max(\det(H))\cdot 0.1的点。

输入 (thorino.jpg)	输出

python实现：

import cv2
import numpy as np
import matplotlib.pyplot as plt

# Hessian corner detection
def Hessian_corner(img):

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

    ## Sobel
    def Sobel_filtering(gray):
        # get shape
        H, W = gray.shape

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

        sobelx = 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(H):
            for x in range(W):
                Ix[y, x] = np.mean(tmp[y : y  + 3, x : x + 3] * sobelx)
                Iy[y, x] = np.mean(tmp[y : y + 3, x : x + 3] * sobely)

        Ix2 = Ix ** 2
        Iy2 = Iy ** 2
        Ixy = Ix * Iy

        return Ix2, Iy2, Ixy



    ## Hessian
    def corner_detect(gray, Ix2, Iy2, Ixy):
        # get shape
        H, W = gray.shape

        # prepare for show detection
        out = np.array((gray, gray, gray))
        out = np.transpose(out, (1,2,0))

        # get Hessian value
        Hes = np.zeros((H, W))

        for y in range(H):
            for x in range(W):
                Hes[y,x] = Ix2[y,x] * Iy2[y,x] - Ixy[y,x] ** 2

        ## Detect Corner and show
        for y in range(H):
            for x in range(W):
                if Hes[y,x] == np.max(Hes[max(y-1, 0) : min(y+2, H), max(x-1, 0) : min(x+2, W)]) and Hes[y, x] > np.max(Hes) * 0.1:
                    out[y, x] = [0, 0, 255]

        out = out.astype(np.uint8)

        return out


    # 1. grayscale
    gray = BGR2GRAY(img)

    # 2. get difference image
    Ix2, Iy2, Ixy = Sobel_filtering(gray)

    # 3. corner detection
    out = corner_detect(gray, Ix2, Iy2, Ixy)

    return out


# Read image
img = cv2.imread("thorino.jpg").astype(np.float32)

# Hessian corner detection
out = Hessian_corner(img)

cv2.imwrite("out.jpg", out)
cv2.imshow("result", out)
cv2.waitKey(0)