SIFT特征点提取与匹配

SIFT特征点提取与匹配

一、 SIFT算法

1、算法简介

  尺度不变特征转换即SIFT (Scale-invariant feature transform)是一种计算机视觉的算法。它用来侦测与描述影像中的局部性特征，它在空间尺度中寻找极值点，并提取出其位置、尺度、旋转不变量，此算法由 David Lowe在1999年所发表，2004年完善总结。

   其应用范围包含物体辨识、机器人地图感知与导航、影像缝合、3D模型建立、手势辨识、影像追踪和动作比对。

  局部影像特征的描述与侦测可以帮助辨识物体，SIFT特征是基于物体上的一些局部外观的兴趣点而与影像的大小和旋转无关。对于光线、噪声、些微视角改变的容忍度也相当高。基于这些特性，它们是高度显著而且相对容易撷取，在母数庞大的特征数据库中，很容易辨识物体而且鲜有误认。使用 SIFT特征描述对于部分物体遮蔽的侦测率也相当高，甚至只需要3个以上的SIFT物体特征就足以计算出位置与方位。在现今的电脑硬件速度下和小型的特征数据库条件下，辨识速度可接近即时运算。SIFT特征的信息量大，适合在海量数据库中快速准确匹配。

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

二、图像特征提取与匹配

图像特征提取是图像分析与图像识别的前提，它是将高维的图像数据进行简化表达最有效的方式，从一幅图像的M×N×3的数据矩阵中，我们看不出任何信息，所以我们必须根据这些数据提取出图像中的关键信息，一些基本元件以及它们的关系。
特征点主要指的就是能够在其他含有相同场景或目标的相似图像中以一种相同的或至少非常相似的不变形式表示图像或目标。因此图像特征点具有两点基本要求：

差异性：特征点应该呈现出区别于非特征点的明显特征
重复性：对应同一三维点的特征点应该在不同视角中被重复检测到

Harris Corner Detector

harris 基本思想是在角点小的窗口内，如果往任何方向移动窗口都会引起比较大的灰度变化，那么这个点就可以认为是角点。

λ1≈λ2≈0 像素点在光滑的区域(flat)
λ1>0,λ2≈0 像素点在边缘区域(edge)
λ1,λ2>0 像素在角点区域(corner)

三.运行结果与代码

from pylab import *
from PIL import Image
from numpy import *
import os

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)

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

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')

if name == ‘main’:
# imname = ‘raccoon.jpg’
# im1 = array(Image.open(imname).convert(‘L’))
# process_image(imname, ‘raccoon.sift’)
# l1, d1 = read_features_from_file(‘raccoon.sift’)
#
# figure()
# gray()
# plot_features(im1, l1, circle=True)
# show()

im1f = r'zhongshan21.jpg'
im2f = r'zhongshan11.jpg'
im1 = array(Image.open(im1f))
im2 = array(Image.open(im2f))

process_image(im1f, 'out_sift_1.txt')
l1, d1 = read_features_from_file('out_sift_1.txt')
figure()
gray()
subplot(121)
plot_features(im1, l1, circle=False)

process_image(im2f, 'out_sift_2.txt')
l2, d2 = read_features_from_file('out_sift_2.txt')
subplot(122)
plot_features(im2, l2, circle=False)

matches = match_twosided(d1, d2)
print('{} matches'.format(len(matches.nonzero()[0])))
figure()
gray()
plot_matches(im1, im2, l1, l2, matches, show_below=True)
show()

从这两个结果可以看出，该算法对特征点的检测效果不错，但是在图像匹配方向，效果一般，甚至出现几乎相同的角度都无法有效匹配的情况。

地理匹配
为了创建显示可能图像组的图，如果匹配的数目高于一个阙值，我们使用边来连接相应的图像节点。为了使图像看起来漂亮，我们需要将每幅图像尺度化为缩略图形式，下面是具体实现代码：

具体代码：

-- coding: utf-8 --

import json
import os
import urllib
from pylab import *
from PIL import Image
from PCV.localdescriptors import sift
from PCV.tools import imtools
import pydot

if name == ‘main’:
download_path = “zhongshan1.jpg”
path = “zhongshan1”

imlist = imtools.get_imlist('zhongshan11.jpg')
nbr_images = len(imlist)

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

# 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) + '.jpg'
            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((100, 100))
            filename = path + str(j) + '.jpg'
            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_jpg('22.jpg')