python计算机视觉-图像处理基础章节之第二章 图像局部描述符

---

前言

本文大量的参考了<<python计算机视觉编程>>一书的第二章，以及大量的优质博客，温馨提示，本文内容过多，可依靠博客目录学习。

---

一、Harris角点检测算法

1.角点(corner points)是什么?

角点具有的特征:
局部窗口沿各方向移动，均产生明显变化的点；图像局部曲率突变的点；轮廓之间的交点；
对于同一场景，即使视角发生变化，通常具备稳定性质的特征；
该点附近区域的像素点无论在梯度方向上还是其梯度幅值上有着较大变化；

2.角点检测算法基本思想是什么?

该算法基本思想是:
使用一个固定小窗口在图像上进行任意方向上的滑动，比较滑动前与滑动后两种情况，窗口中的像素灰度变化程度，如果存在任意方向上的滑动，都有着较大灰度变化，那么我们可以认为该窗口中存在角点。

3.什么是好的角点检测算法?

●检测出图像中“真实的”角点
●准确的定位性能
●很高的稳定性
●具有对噪声的鲁棒性
●具有较高的计算效率

4.角点的特征如何用数学方法来刻画呢？

当窗口发生[u，v]移动的时候，滑动前和滑动后对应的窗口中的像素点灰度变化描述如下:

所以得到，当在平坦区域时，窗口向任意方向移动，无灰度变化；在边缘区域时，沿着边缘方向移动，无灰度变化；在角点区域时，沿着任意方向移动，有明显灰度变化

经过一系列E(u，v)表达式的演化， E(u,v)表达式可以更新为：

其中矩阵M是2x2矩阵，可由图像的导数求得:

5.如何用度量角点响应?

角点响应函数R


取决于该区域vI的值，Harris矩阵Mi的特征值有三种情况:
1)如果 λ1和 λ 2都是很大的正数，则该x点为角点;
2)如果λ1很大,λ 2 很小(接近于0)，则该区域内存在一个边，该区域内的平均Mi的特征值不会变化太大;
3)如果λ1λ2 都是很小的数 (λ1≈λ2≈0), 则该区域为空。

6.Harris角点检测器的响应函数是什么?

Harris角点检测器的响应函数会返回像素值为 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.py 文件中：
代码如下：

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

7. Harris角点检测算法优点：

1. 旋转不变性，椭圆转过一定角度但是其形状保持不变（特征值保持不变）。
2. 对于图像灰度的仿射变化具有部分的不变性，由于仅仅使用了图像的一介导数，对于图像灰度平移变化不变；
3. 对于图像灰度尺度变化不变。

8.Harris角点检测算法缺点：

1. 它对尺度很敏感，不具备几何尺度不变性。
2. 提取的角点是像素级的

9.Harris角点检测实例

代码:

# -*- codeing =utf-8 -*-
# @Time : 2021/3/22 19:03
# @Author : ArLin
# @File : demo1.py
# @Software: PyCharm
# -*- coding: utf-8 -*-
from pylab import *
from PIL import Image
from PCV.localdescriptors import harris


# 读入图像
im = array(Image.open('../JMU/14.jpg').convert('L'))


# 检测harris角点
harrisim = harris.compute_harris_response(im)


# Harris响应函数
harrisim1 = 255 - harrisim


figure()
gray()


#画出Harris响应图
subplot(141)
imshow(harrisim1)
print (harrisim1.shape)
axis('off')
axis('equal')


threshold = [0.01,0.05,0.1]
for i, thres in enumerate(threshold):


    filtered_coords = harris.get_harris_points(harrisim, 6, thres)
    subplot(1, 4, i+2)
    imshow(im)
    print(im.shape)
    plot([p[1] for p in filtered_coords], [p[0] for p in filtered_coords], '*')
    axis('off')



show()

结果展示:

使用 Harris 角点检测器检测角点：（a）为 Harris 响应函数；（b-d）分别为使用阈值 0.01、0.05 和 0.1 检测出的角点

结果分析:

增大α的值，将减小角点响应值R，降低角点检测的灵性，减少被检测角点的数量；减小α值，将增大角点响应值R，增加角点检测的灵敏性，增加被检测角点的数量。使用阈值 0.01、0.05 和 0.1 检测出的角点依次减少。

10.图像中寻找对应点：

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

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

代码：

# -*- codeing =utf-8 -*-
# @Time : 2021/3/22 19:23
# @Author : ArLin
# @File : demo2.py
# @Software: PyCharm
# -*- coding: utf-8 -*-
from pylab import *
from PIL import Image
from PCV.localdescriptors import harris
from PCV.tools.imtools import imresize

im1 = array(Image.open("../JMU/9.jpg").convert("L"))
im2 = array(Image.open("../JMU/12.jpg").convert("L"))

# 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 = harris.compute_harris_response(im1, 5)
filtered_coords1 = harris.get_harris_points(harrisim, wid+1)
d1 = harris.get_descriptors(im1, filtered_coords1, wid)

harrisim = harris.compute_harris_response(im2, 5)
filtered_coords2 = harris.get_harris_points(harrisim, wid+1)
d2 = harris.get_descriptors(im2, filtered_coords2, wid)

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

figure()
gray()
harris.plot_matches(im1, im2, filtered_coords1, filtered_coords2, matches)
show()

结果展示:

结果分析：

虽然大部分得到了对应点的匹配，
但是我们很明显也能看出该算法的结果中也存在一些不正确的匹配，这是因为，与现代的一些方法相比，图像像素块的互相关矩阵具有较弱的描述性。实际运用中，我们通常使用更稳健的方法来处理这些对应匹配。这些描述符还有一个问题，它们不具有尺度不变性和旋转不变性，而算法中像素块的大小也会影响对应匹配的结果。
而且，Harris的算法效率也不高。

---

二、SIFT（尺度不变特征变换）

1.为什么要提出SIFT(目的和意义)?:

1999年David G. Lowe教授总结了基于特征不变技术的检测方法，在图.
像尺度空间基础上，提出了对图像缩放、旋转保持不变性的图像局部
特征描述算子-SIFT (尺度不变特征变换)，该算法在2004年完善。

2.SIFT算法可以解决什么问题?

●目标的旋转、缩放、平移(RST)
●图像仿射/投影变换(视点viewpoint )
●弱光照影响( ilumination)
●部分目标遮挡(occlusion)
●杂物场景( clutter)
●噪声

3.哪些点是SIFT中要查找的关键点(特征点) ?

1. 尺度
2. 方向
3. 位移
4. 光照

4.SIFT特性：

1. 特性独特性，也就是特征点可分辨性高，类似指纹，适合在海量数据中匹配。
2. 多量性，提供的特征多。
3. 高速性，就是速度快。
4. 可扩展，能与其他特征向量联合使用。

5.什么是尺度空间?

尺度空间理论最早于1962年提出，其主要思想是通过对原始图像进行尺度变换，获得图像多尺度下的空间表示。从而实现边缘、角点检测和不同分辨率上的特征提取，以满足特征点的尺度不变性。
尺度空间中各尺度图像的模糊程度逐渐变大，能够模拟人在距离目标由近到远时目标在视网膜上的形成过程。
尺度越大图像越模糊。

6.SIFT算法如何实现特征匹配?

在不同尺度空间上查找特征点（关键点）的问题)
1、提取关键点;
2、对关键点附加详细的信息(局部特征)，即描述符;
3、通过特征点(附带上特征向量的关键点)的两两比较找出相互匹配的若干对特征点，建立景物间的对应关系。

7.Sift特征的生成一般包括以下几个步骤：

（1）构建尺度空间，检测极值点，获得尺度不变性；

（2）特征点过滤并进行精确定位；

（3）为特征点分配方向值；

（4）计算变换参数
当两幅图像的Sift特征向量生成以后，下一步就可以采用关键点特征向量的欧式距离来作为两幅图像中关键点的相似性判定度量。取图1的某个关键点，通过遍 历找到图像2中的距离最近的两个关键点。在这两个关键点中，如果次近距离除以最近距离小于某个阙值，则判定为一对匹配点。

8.检测兴趣点

使用开源工具包 VLFeat 提供的二进制文件来计算图像的 SIFT 特征 。VLFeat 工具包可以从 http://www.vlfeat.org/ 下载，二进制文件可以在所有主要的平台上运行。

代码:

def process_image(imagename, resultname, params="--edge-thresh 10 --peak-thresh 5"):
    """处理一幅图像，然后将结果保存在文件中"""
    
    if imagename[-3:] != 'pgm':
        #创建一个pgm文件
        im = Image.open(imagename).convert('L')
        im.save('tmp.pgm')
        imagename = 'tmp.pgm'
        
    cmmd= str("sift"+imagename+"--output="+resultname+" "+params)
    os.system(cmmd)
    print('processed',imagename,'to',resultname)

由于该二进制文件需要的图像格式为灰度.pgm，所以如果图像为其他各是，我们需要首先将其转换成.pgm格式文件。其中数据的每一行前4个数值依次表示兴趣点的坐标、尺度和方向角度，后面紧跟着的是对应描述符的128维向量。

下面函数用于从上面输出文件中，将特征读取到Numpy数组中的函数。

def read_features_from_file(filename):
    """读取特征值属性值，然后将其以矩阵形式返回"""
    
    f = loadtxt(filename)
    return f[:,:4],f[:,4:] # 特征位置，描述子

读取特征后，通过在图像上绘制出它们的位置，可以将其可视化。下面的函数可以实现这个功能：

def plot_features(im, locs, circle=False):
    """显示带有特征的图像
        输入：im（数组图像），locs（每个特征的行、列、尺度和方向角度）"""


    def draw_circle(c,r):
        t = arange(0,1.01,.01)*2*pi
        x = r*cos(t) + c[0]
        y = r*sin(t) + c[1]
        plot(x,y,'b',linewidth=2)


    imshow(im)
    if circle:
        for p in locs:
            draw_circle(p[:2],p[2])
    else:
        plot(locs[:,0],locs[:,1],'ob')
    axis('off')
    return

代码演示:

# -*- codeing =utf-8 -*-
# @Time : 2021/3/22 19:34
# @Author : ArLin
# @File : demo3.py
# @Software: PyCharm
# -*- coding: utf-8 -*-
from PIL import Image
from pylab import *
from PCV.localdescriptors import sift
from PCV.localdescriptors import harris
import os






# 添加中文字体支持2
from matplotlib.font_manager import FontProperties
font = FontProperties(fname=r"d:\JsVison\font\SimSun.ttc", size=14)


imname = '../JMU/7.jpg'
im = array(Image.open(imname).convert('L'))
sift.process_image(imname, 'empire.sift')
l1, d1 = sift.read_features_from_file('empire.sift')




figure()
gray()
subplot(131)
sift.plot_features(im, l1, circle=False)
title(u'SIFT特征',fontproperties=font)
subplot(132)
sift.plot_features(im, l1, circle=True)
title(u'用圆圈表示SIFT特征尺度',fontproperties=font)


# 检测harris角点
harrisim = harris.compute_harris_response(im)


subplot(133)
filtered_coords = harris.get_harris_points(harrisim, 6, 0.1)
imshow(im)
plot([p[1] for p in filtered_coords], [p[0] for p in filtered_coords], '*')
axis('off')
title(u'Harris角点',fontproperties=font)

show()

结果展示:

结果分析:

sift和Harris角点的两种算子选择了不同的坐标
SIFT算法的实质是在不同的尺度空间上查找关键点(特征点)，并计算出关键点的方向。SIFT所查找到的关键点是一些十分突出，不会因光照，仿射变换和噪音等因素而变化的点，如角点、边缘点、暗区的亮点及亮区的暗点等。

9.关键点检测一_DOG

DoG ( Difference of Gaussian)函数

DoG在计算上只需相邻高斯平滑后图像相减，因此简化了计算!.

10.描述子

上面讨论的兴趣点（关键点）位置描述子给出了兴趣点的位置和尺度信息。为了实现旋转不变性，基于每个点周围图像梯度的方向和大小，SIFT 描述子又引入了参考方向。SIFT 描述子使用主方向描述参考方向。主方向使用方向直方图（以大小为权重）来度量。

11.匹配描述子

对于将一幅图像中的特征匹配到另一幅图像的特征，一种稳健的准则（同样是由Lowe提出的）是使用者两个特征距离和两个最匹配特征距离的比率。相比于图像中的其他特征，该准则保证能够找到足够相似的唯一特征。使用该方法可以使错误的匹配数降低。下面的代码实现了匹配函数。

代码:

def match(desc1, desc2):
    """对于第一幅图像的每个描述子，选取其在第二幅图像中的匹配
        输入：desc1（第一幅图像中的描述子），desc2（第二幅图像中的描述子）"""


    desc1 = array([d/linalg.norm(d) for d in desc1])
    desc2 = array([d/linalg.norm(d) for d in desc2])


    dist_ratio = 0.6
    desc1_size = desc1.shape


    matchscores = zeros((desc1_size[0],1), 'int')
    desc2t = desc2.T    #预先计算矩阵转置
    for i in range(desc1_size[0]):
        dotprods = dot(desc1[i,:], desc2t) #向量点乘
        dotprods = 0.9999*dotprods
        # 反余弦和反排序，返回第二幅图像中特征的索引
        index = argsort(arccos(dotprods))


        # 检查最近邻的角度是否小于dist_ratio乘以第二近邻的角度
        if arccos(dotprods)[index[0]] < dist_ratio * arccos(dotprods)[index[1]]:
            matchscores[i] = int(index[0])


    return matchscores

def match_twosided(desc1,decs2):
    """双向对称版本的match"""


    matches_12 = match(desc1, decs2)
    matches_21 = match(decs2, decs2)


    ndx_12 = matches_12.nonzero()[0]


    # 去除不对称匹配
    for n in ndx_12:


        if matches_21[int(matches_12[n])] != n:
            matches_12[n] = 0


    return matches_12

def appendimages(im1, im2):
    """返回将两幅图像并排拼接成的一幅新图像"""


    # 选取具有最少行数的图像，然后填充足够的空行
    row1 = im1.shape[0]
    row2 = im2.shape[0]


    if row1 < row2:
        im1 = concatenate((im1,zeros((row2-row1,im1.shape[1]))), axis=0)
    elif row1 > row2:
        im2 = concatenate((im2,zeros((row1-row2,im2.shape[1]))), axis=0)


    # 如果这些情况都没有，那么他们的行数相同，不需要进行填充


    return concatenate((im1,im2), axis=1)

def plot_matches(im1, im2, locs1, locs2, matchscores, show_below=True):
    """显示一幅带有连接匹配之间连线的图片
        输入：im1，im2（数组图像），locs1，locs2（特征位置），matchscores（match的输出），
        show_below（如果图像应该显示再匹配下方）"""


    im3 = appendimages(im1,im2)
    if show_below:
        im3 = vstack((im3,im3))


    imshow(im3)


    cols1 = im1.shape[1]
    for i in range(len(matchscores)):
        if matchscores[i] > 0:
            plot([locs1[i, 0], locs2[matchscores[i, 0], 0] + cols1], [locs1[i, 1], locs2[matchscores[i, 0], 1]], 'c')
    axis('off')

结果展示:

相同景物不同角度


相同图片不同分辨率:

Harris角点检测结果:

对比分析：

从匹配的准确率看来，相较harris算法。sift的算法具有更高的匹配性，显而易见地，sift算法具有较高的准确度及稳健性（对图片亮度不敏感）；
从时间上来说，sift算法的效率远高与harris算法，所以sift算法具有更高的计算效率；
总的来说，两个算法由于选择特征点的方法不同，而显示出了差异bo
sift算法综合了许多前人已经整理出的算法的优点，
最后才得到了现在的算法（更优）。

以及在用SIFT检测时，注意要用相同大小的图片，否则在程序中无法运行。

12.SIFT的不足:

SIFT采用henssian矩阵获取图像局部最值还是十分稳定的，
但是在求主方向阶段太过于依赖局部区域像素的梯度方向，
有可能使得找到的主方向不准确，后面的特征向量提取以及匹配都严重依赖于主方向，即使不大偏差角度也可以造成后面特征匹配的放大误差，从而匹配不成功；
另外图像金字塔的层取得不足够紧密也会使得尺度有误差，后面的特征向量提取同样依赖相应的度，发明者在这个问题上的折中解决方法是取适量的层然后进行插值。
SIFT是一种只利用到灰度性质的算法，忽略了色彩信息.

三、匹配地理标记图像

代码:

# -*- codeing =utf-8 -*-
# @Time : 2021/3/23 10:45
# @Author : ArLin
# @File : demo8.py
# @Software: PyCharm
# -*- coding: utf-8 -*-
from pylab import *
from PIL import Image
from PCV.localdescriptors import sift
from PCV.tools import imtools
import pydot


""" This is the example graph illustration of matching images from Figure 2-10.
To download the images, see ch2_download_panoramio.py."""


download_path = "D:\python\pytharm\demo\pythonProject\JsVision\jmujpg"  # set this to the path where you downloaded the panoramio images
path = "D:\python\pytharm\demo\pythonProject\JsVision\jmujpg"  # path to save thumbnails (pydot needs the full system path)
      
# list of downloaded filenames
imlist = imtools.get_imlist(download_path)
nbr_images = len(imlist)


# extract features
featlist = [imname[:-3] + 'sift' for imname in imlist]
for i, imname in enumerate(imlist):
    sift.process_image(imname, featlist[i])


matchscores = zeros((nbr_images, nbr_images))


for i in range(nbr_images):
    for j in range(i, nbr_images):  # only compute upper triangle
        print('comparing ', imlist[i], imlist[j])
        l1, d1 = sift.read_features_from_file(featlist[i])
        l2, d2 = sift.read_features_from_file(featlist[j])
        matches = sift.match_twosided(d1, d2)
        nbr_matches = sum(matches > 0)
        print ('number of matches = ', nbr_matches)
        matchscores[i, j] = nbr_matches
print ("The match scores is: \n", matchscores)


# copy values
for i in range(nbr_images):
    for j in range(i + 1, nbr_images):  # no need to copy diagonal
        matchscores[j, i] = matchscores[i, j]


#可视化


threshold = 2  # min number of matches needed to create link


g = pydot.Dot(graph_type='graph')  # don't want the default directed graph


for i in range(nbr_images):
    for j in range(i + 1, nbr_images):
        if matchscores[i, j] > threshold:
            # first image in pair
            im = Image.open(imlist[i])
            im.thumbnail((100, 100))
            filename = path + str(i) + '.png'
            im.save(filename)  # need temporary files of the right size
            g.add_node(pydot.Node(str(i), fontcolor='transparent', shape='rectangle', image=filename))


            # second image in pair
            im = Image.open(imlist[j])
            im.thumbnail((200, 200))
            filename = path + str(j) + '.png'
            im.save(filename)  # need temporary files of the right size
            g.add_node(pydot.Node(str(j), fontcolor='transparent', shape='rectangle', image=filename))
            g.add_edge(pydot.Edge(str(i), str(j)))
g.write_png('seeingJMU3.png')

结果展示:

原照片集:

测试1:

测试2:

结果分析:

测试1用了七张图，大小没有经过任何的调整，所以存在5M以上大小，图片也基本都是4032*2687这个大小左右，匹配时间较长，但是准确度非常高，在角度，光线，树叶等影响下，仍然具有非常高的匹配性
测试2的21张图，经过设定500&500大小的图片后，匹配效率迅猛提升，但是同样可以看到因为经过压缩后，失去了一定的准确性

---

四、总结

总的来说，这篇博客的内容较多，其实在学习的过程中也耗费了大量的时间，参考了大量的优质博客，在实现的过程中也遇到了很多的问题，最后在网络搜索，同学，助教帮助，老师讲解的情况加还是顺利的解决了。
可以提升的地方:
会再捋一遍这两个算法的逻辑，争取更有条理的讲解说明
接下来可能也会把图片的相关资源以及遇到的问题以及解决过程发在网上。