计算机视觉 （二）

这次的博客是关于图像局部特征描述子。

特征点检测、匹配都是一些图像处理的重要步骤，特征点必须具有不变性：位移、旋转、尺度、光照、曝光等。特征包括：角点、梯度特征点、边缘特征、纹理特征。

Harris角点检测器

角点：局部窗口沿各方向移动，均发生明显变化的点，图像局部曲率突变的点。

Harris角点检测基本思想，从图像局部的小窗口观察图像特征，窗口向任意方向移动都导致图像灰度的明显变化即为角点，如下图红色窗口。

 

根据泰勒展开式

 

以上的过程的数学表达：α=0.04~0.06

|R|为小数值就是平坦区，角点为大数值正数（要设置一个数值来判断大小），边缘为大数值负数

ps：为了消除参数α的影响，也可以采用商来计算。

Q：那角点检测是否就是直接选择响应值R大的点呢？

A：不是，因为响应值R大的点可能会集中在一个区域内——一个真正角点周围的临近点R值也会很大。所以我们还需要一步非极大值抑制来提取R的局部极大值。

from pylab import *
from PIL import Image
from pylab import *
from numpy import *
from scipy.ndimage import filters
import cv2 


def compute_harris_response(im,sigma=3):
    #计算响应值

    
    # 计算导数
    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)
    
    
    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

以上就是计算响应值R 的函数。

根据响应值R对候选角点进行排序，其中还包括非极大值抑制：

def get_harris_points(harrisim,min_dist=10,threshold=0.1):

    
    # 找到R高于阈值的点
    corner_threshold = harrisim.max() * threshold
    harrisim_t = (harrisim > corner_threshold) * 1
    
    # 坐标和R
    coords = array(harrisim_t.nonzero()).T
    candidate_values = [harrisim[c[0],c[1]] for c in coords]
    
    # 对R排序
    index = argsort(candidate_values)[::-1]
    allowed_locations = zeros(harrisim.shape)
    allowed_locations[min_dist:-min_dist,min_dist:-min_dist] = 1
    
    # 非极大值抑制
    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()
    imshow(image)
    plot([p[1] for p in filtered_coords],
                [p[0] for p in filtered_coords],'*')
    axis('off')
    

 结果：

在图像中寻找对应点对

检测出图像中的角点后，还需要对不同图片检测出的角点进行匹配：给检测出的角点打上标签，根据标签对角点进行匹配，这些标签叫做角点的描述子 ——是分配给角点的一个向量，这个向量通常是由周围图像像素块的灰度值，以及用于比较的归一化互相关矩阵构成的，通常，两个相同大小像素块的相关矩阵定义为：

 c值越大，两个像素块的相关度越高，下列为归一化后的互相关矩阵，归一化可以保证不在一个规模的数据之间变得具有可比性：

 

def get_descriptors(image,filtered_coords,wid=5):
    """返回角点周围2*wid +1个像素的值 ，计算每个角点的描述子"""
    
    desc = []
    for coords in filtered_coords:
        patch = image[coords[0]-wid:coords[0]+wid+1,
                            coords[1]-wid:coords[1]+wid+1].flatten()
        desc.append(patch)
    
    return desc
def match(desc1,desc2,threshold=0.5):
    """ 对角点描述子，进行归一化互相关，匹配两幅图片的角点. """
    
    n = len(desc1[0])
    
    # pair-wise distances
    d = -ones((len(desc1),len(desc2)))
    for i in range(len(desc1)):
        for j in range(len(desc2)):
            d1 = (desc1[i] - mean(desc1[i])) / std(desc1[i])
            d2 = (desc2[j] - mean(desc2[j])) / std(desc2[j])
            ncc_value = sum(d1 * d2) / (n-1) 
            if ncc_value > threshold:
                d[i,j] = ncc_value
            
    ndx = argsort(-d)
    matchscores = ndx[:,0]
    
    return matchscores

以上是单项地从第一幅图像向第二幅图像进行匹配，若也从第二幅图像向第一幅图像匹配一次，就可以得到更稳定的结果：

def match_twosided(desc1,desc2,threshold=0.5):
    """ 进行双向匹配 """
    
    matches_12 = match(desc1,desc2,threshold)
    matches_21 = match(desc2,desc1,threshold)
    
    ndx_12 = where(matches_12 >= 0)[0]
    
    # 去除不对称的匹配点对
    for n in ndx_12:
        if matches_21[matches_12[n]] != n:
            matches_12[n] = -1
    
    return matches_12

连接匹配角点对，实现匹配可视化：

def appendimages(im1,im2):
    
    
    rows1 = im1.shape[0]    
    rows2 = im2.shape[0]
    
    if rows1 < rows2:
        im1 = concatenate((im1,zeros((rows2-rows1,im1.shape[1]))),axis=0)
    elif rows1 > rows2:
        im2 = concatenate((im2,zeros((rows1-rows2,im2.shape[1]))),axis=0)
    
    return concatenate((im1,im2), axis=1)
    
    
def plot_matches(im1,im2,locs1,locs2,matchscores,show_below=True):
    """ 显示带有匹配对的图像们 """
    
    im3 = appendimages(im1,im2)
    if show_below:
        im3 = vstack((im3,im3))
    
    imshow(im3)
    
    cols1 = im1.shape[1]
    for i,m in enumerate(matchscores):
        if m>0:
            plot([locs1[i][1],locs2[m][1]+cols1],[locs1[i][0],locs2[m][0]],'c',linewidth = 0.3,color = 'r')
    axis('off')
    plt.savefig('C:\Python\imgs\p3.jpg',dpi=1080)

 测试：

def imresize(im,sz):
    pil_im = Image.fromarray(uint8(im))
    
    return array(pil_im.resize(sz))
# Figure 2-2下面的图
img1 = array(Image.open("C:\Python\imgs\p1.jpg").convert('RGB'))
img2 = array(Image.open("C:\Python\imgs\p2.jpg").convert('RGB'))

im1 = array(Image.open("C:\Python\imgs\p1.jpg").convert('L'))
im2 = array(Image.open("C:\Python\imgs\p2.jpg").convert('L'))

img1 = imresize(img1, (img1.shape[1]//2, img1.shape[0]//2))
img2 = imresize(img2, (img2.shape[1]//2, img2.shape[0]//2))
# resize让匹配更快点
im1 = imresize(im1, (im1.shape[1]//2, im1.shape[0]//2))
im2 = imresize(im2, (im2.shape[1]//2, im2.shape[0]//2))


wid = 5
harrisim = compute_harris_response(im1, 5)
filtered_coords1 = get_harris_points(harrisim, wid+80)
d1 = get_descriptors(im1, filtered_coords1, wid)

harrisim = compute_harris_response(im2, 5)
filtered_coords2 = get_harris_points(harrisim, wid+80)
d2 = get_descriptors(im2, filtered_coords2, wid)

print ('starting matching')
matches = match_twosided(d1, d2)

figure()
gray() 
plot_matches(img1, img2, filtered_coords1, filtered_coords2, matches)
show()

结果：

 但是对于旋转和尺度不同的特征点Harris角点特征检测仍然有问题。