python计算机视觉编程第二章 局部图像描述子

1.Harris角点检测

Harris角点检测算法（也称Harris&Stephens角点检测器）是一个极为简单的角点检测算法。该算法的主要思想是，如果像素周围显示存在多于一个方向的边，我们认为该点为兴趣点。该点就称为角点。

我们把图像域中点x上的对称半正定矩阵M_I=M_I(x)定义为：

其中∇I为包含导数I_x和I_y的图像梯度。由于该定义，M_I的秩为1，特征值为λ₁=|∇I|²和λ₂=0。现在对于图像的每个像素，我们可以计算出该矩阵。

选择权重矩阵W（通常为高斯滤波器G_σ），我们可以得到卷积：

该卷积的目的是得到M_I在周围像素上的局部平均。计算出的矩阵M_I^-有称为Harris矩阵。W的宽度决定了在像素x周围的感兴趣区域。像这样在区域附近对矩阵M_I^-取平均的原因是，特征值会依赖于局部图像特性而变化。如果图像梯度在该区域变化，那么M_I^-的第二个特征值将不再为0.如果图像的梯度没有变换，M_I^-的特征值也不会变化。

取决于该区域∇I的值，Harris矩阵M_I^-的特征值有三种情况： 1）如果λ₁和λ₂都是很大的正数，则该x点为角点； 2）如果λ₁很大，λ₂≈0，则该区域内存在一个边，该区域内的平均M_I^-的特征值不会变化太大； 3）如果λ₁≈λ₂≈0，则该区域为空。

在不需要实际计算特征值的情况下，为了把重要的情况和其他情况分开，Harris和Stephens在文献中引入了指示函数：

为了去除加权常数κ，我们通常使用商数：

作为指示器。

首先，我们先完成计算Harris角点检测器的响应函数：

def compute_harris_response(im, sigma=3):
    """在一幅灰度图像中，对每个像素计算Harris角点检测器响应函数"""
    
    # 计算导数
    imx = zeros(im.shape)
    filters.gaussian_filter(im, (sigma, sigma), (0,1), imx)    
    imy = zeros(im.shape)
    filters.gaussian_filter(im, (sigma, sigma), (1,0), imy)
    
    # 计算Harris矩阵的分量
    Wxx = filters.gaussian_filter(imx*imx, sigma)
    Wxy = filters.gaussian_filter(imx*imy, sigma)
    Wyy = filters.gaussian_filter(imy*imy, sigma)
    
    # 计算特征值和迹
    Wdet = Wxx * Wyy - Wxy**2
    Wtr = Wxx + Wyy
    
    return Wdet / Wtr

上面的函数会返回像素值为Harris响应函数值的一幅图像。现在，我们需要从这幅图像中挑选需要的信息。然后，选取像素值高于阈值的所有图像点；在加上额外的限制，即角点之间的间隔必须大于设定的最小距离。这种方法会产生很好的角点检测结果。为了实现该算法，我们获取所有候选像素点，以角点响应值递减的顺序排序，然后将距离已标记为角点位置过近的区域从候选像素点中删除：

def get_harris_points(harrisim, min_dist=10, threshold=0.1):
    """从一幅Harris响应图像中返回角点。min_dist为分割角点和图像边界的最小像素数目"""

    # 寻找高于阈值的候选角点
    corner_threshold = harrisim.max() * threshold
    harrisim_t = (harrisim > corner_threshold) * 1

    # 得到候选点的坐标
    coords = array(harrisim_t.nonzero()).T

    # 以及它们的Harris响应值
    candidate_values = [harrisim[c[0], c[1]] for c in coords]

    # 对候选点按照Harris响应值进行排序
    index = argsort(candidate_values)
    
    # 将可行点的位置保存到数组中
    allowed_locations = zeros(harrisim.shape)
    allowed_locations[min_dist:-min_dist,min_dist:-min_dist] = 1
    
    # 按照min_distance原则，选择最佳Harris点
    filtered_coords = []
    for i in index:
        if allowed_locations[coords[i,0], coords[i,1]] == 1:
            filtered_coords.append(coords[i])
            allowed_locations[(coords[i,0]-min_dist):(coords[i,0]+min_dist),
            (coords[i,1]-min_dist):(coords[i,1]+min_dist)] = 0
            
    return filtered_coords

def plot_harris_points(image,filtered_coords):
    """绘制图像中检测到的角点"""
    figure()
    gray()
    imshow(image)
    plot([p[1] for p in filtered_coords],[p[0] for p in filtered_coords],'*')
    axis('off')
    show()

if __name__=='__main__':
    im = array(Image.open('D:\\123\图像处理\Image Processing\Image Processing\jimei1.jpg').convert('L'))
    harrisim = compute_harris_response(im)
    filtered_coords = get_harris_points(harrisim,6)
    plot_harris_points(im, filtered_coords)

1.1在图像间寻找对应点

Harris角点检测器仅仅能够检测出图像中的兴趣点，但没有给出通过比较图像间的兴趣点来寻找匹配角点的方法。我们需要在每个点上加入描述子信息，并给出一个比较这些描述子的方法。

兴趣点描述子是分配给兴趣点的一个向量，描述该点附近的图像的表现信息。描述子越好，寻找到的对应点越好。我们用对应点或者点的对应来描述相同物体和场景点在不同图像上形成的像素点。

Harris角点的描述子通常是由周围图像像素块的灰度值，以及用于比较的归一化互相关矩阵构成的。图像的像素块由以该像素点为中心的周围矩形部分图像构成。

通常，两个（相同大小）像素块I₁(x)和I₂(x)的相关矩阵定义为：

其中，函数f随着相关方法的变化而变化。上式取像素块中所有像素位置x的和。对于互相关矩阵，函数f(I₁,I₂)=I₁I₂，因此，c(I₁,I₂)=I₁·I₂，其中·表示向量乘法。c(I₁,I₂)的值越大，像素块I₁和I₂的相似度越高。

归一化的互相关矩阵是互相关矩阵的一种变形，可以定义为：

其中，n为像素块中像素的数目，μ₁和μ₂表示每个像素块中平均像素值强度，σ₁和σ₂分别表示每个像素块中的标准差。通过减去均值和除以标准差，该方法对图像亮度变化具有稳健性。

def get_descriptors(image, filtered_coords, wid=5):
    "