找到栅栏的尖锐点之角点检测邻居家新建了一排栅栏，上面尖尖的角非常可怕，你能写一个程序检测出这些尖尖角的位置，以免下次不会

一、问题背景

邻居家有一排栅栏，上面尖尖的角非常恐怖，你能写一个程序检测出这些尖尖角的位置，以免下次不会撞上吗？

二、Harris 角点检测算法思想

使用 Harris角点检测子可以检测出图像中的角点，其基本思想为：

利用滑动窗口的思想，让以当前检测点为中心的小窗口在附近滑动。如下图所示，当滑动窗口在所有方向移动时，窗口内的像素灰度出现了较大的变化，就可能是角点。

注：如果是平坦点，则在附件滑动时，像素基本没有变化，如下最左侧图，如果是一条边，则在沿着边的方向滑动不会产生变化，如果时候角点，则无论如何移动，像素值都会发生变化。

根据上面的思想，我们计算出在点(x, y)处移动窗口到(x+u, y+v)点，产生的像素值变化 $E(u, v)$ 记为：

E(u, v) = \sum_{x, y}{w(x, y){[I(x+u, y+v) - I(x, y)]}^2}

其中w(x,y)是一种加权函数，可以设为高斯函数。

根据泰勒公式，将其在(0, 0)展开到二阶，得到：

E(u, v) = E(0, 0) + (u, v) \begin{bmatrix}E_1(0, 0) \\ E_2(0, 0) \end{bmatrix} + {1 \over 2} (u, v) \begin{bmatrix}E_{1, 1}(0, 0) & E_{1, 2}(0, 0) \\ E_{2, 1}(0, 0) & E_{2, 2}(0, 0) \end{bmatrix} \begin{bmatrix}u \\ v \end{bmatrix}

我们可以计算得到:

E(0, 0) = 0 \\ E_1(0, 0) = E_2(0, 0) = 0 \\ E_{1, 1} = 2\sum_{x, y} {w(x, y)I_1(x, y)*I_1(x, y)} \\ ...

代入 $E(U, V)$ 可得到:

E(u, v) = {1 \over 2} (u, v) \begin{bmatrix}E_{1, 1}(0, 0) & E_{1, 2}(0, 0) \\ E_{2, 1}(0, 0) & E_{2, 2}(0, 0) \end{bmatrix} \begin{bmatrix}u \\ v \end{bmatrix} = (u, v)M\begin{bmatrix}u \\ v \end{bmatrix} \propto M

其中矩阵M定义为：

M = \sum_{x, y} \begin{bmatrix}I_1^2 & I_1I_2 \\ I_2I_1 & I_2^2 \end{bmatrix} \\

注： $I_1$ 表示 $I(x, y)$ 对其第一个参数 $x$ 求导，亦可表示为 $I_x$ 。

假设角点是互相垂直的，即 $I_1 * I_2 = I_x * I_y = 0$ ，则有

M = \begin{bmatrix}\lambda_1 & 0 \\ 0 & \lambda_2 \end{bmatrix} \\

我们使用 $\lambda_1, \lambda_2$ 分别表示 $x, y$ 方向上的像素变化强度，只有当两者都非常大时，我们就认为检测到了角点。

所以我们来构造一个函数，来判断 $\lambda_1, \lambda_2$ 是否同时很大:

R(\lambda_1, \lambda_2) = |M| - k(tr(M))^2 = \lambda_1 × \lambda_2 - k(\lambda_1 + \lambda_2)^2

其中 $k$ 为一超参数。

三、Harris 角点检测算法步骤

计算图片在x, y方向上的导数dx, dy
计算每一个滑动窗口内的二阶矩矩阵M
计算响应函数 $R(\lambda_1, \lambda_2)$
根据阈值判断当前点是否为角点

我们使用一、中的思想，根据上面步骤来检测下图中栅栏上的角点

四、Harris 角点检测算法实现

'''
img: 灰度图
win_size: 滑动窗口大小
k: 常数
thresh: 判断为角点的阈值
'''
def detect_corners(img:np.ndarray, win_size:int, k:float, thresh:int):
    # 1. 计算图片在`x, y`方向上的导数`dx, dy`
    dy, dx = np.gradient(img)

    Ixx = dx ** 2, Ixy = dy * dx, Iyy = dy ** 2

    offset = win_size//2
    height = img.shape[0]
    width = img.shape[1]

    corners = [] # 记录角点坐标

    for y in range(offset, height-offset):
        for x in range(offset, width-offset):
            # 2. 计算每一个滑动窗口内的二阶矩矩阵`M`
            win_Ixx = Ixx[y-offset:y+offset+1, x-offset:x+offset+1].sum()
            win_Ixy = Ixy[y-offset:y+offset+1, x-offset:x+offset+1].sum()
            win_Iyy = Iyy[y-offset:y+offset+1, x-offset:x+offset+1].sum()

            det = win_Ixx * win_Iyy - win_Ixy * win_Ixy
            tr = win_Ixx + win_Iyy
            
            # 3. 计算响应函数 $R(\lambda_1, \lambda_2)$
            R = det - k * tr
            
            # 4. 根据阈值判断当前点是否为角点
            if R > thresh:
                corners.append((x, y))

    return corners

使用上面的方法我们上面侧栅栏图片，得到检测的效果如下所示。

以上完整代码见 detect_corners