SIFT算法应用：全景图像拼接

什么是全景图像拼接

全景图像拼接顾名思义就是将多张存在重叠部分的图像拼成一幅全景图像。个人能理解为：两幅图像之间存在匹配的点的情况下，可以通过特征匹配得到对应点，将这些对应点坐标重合而保留两张图像其他的部分，就可以得到两幅图像的拼接结果，以这种方式循环，通过两两拼接的方式最终拼接成一幅全景图。

---

实现全景图像拼接需要哪些步骤？

可以简单分为以下几步：

1.根据给定图像/集，实现特征匹配。
2.通过匹配特征计算图像之间的变换结构。
3.利用图像变换结构，实现图像映射。
4.针对叠加后的图像，采用APAP之类的算法，对齐特征点。(图像配准)
5.通过图割方法，自动选取拼接缝。
6.根据multi-band blending策略实现融合。

---

相关介绍

1.RANSAC

RANSAC 是“RANdom SAmple Consensus”（随机一致性采样）的缩写。该方法是用来找到正确模型来拟合带有噪声数据的迭代方法。给定一个模型，例如点集之间的单应性矩阵，RANSAC 基本的思想是，数据中包含正确的点和噪声点，合理的模型应该能够在描述正确数据点的同时摒弃噪声点。

RANSAC 是个非常有用的算法，我们将在估计单应性矩阵中使用它。

2.稳健的单应性矩阵估计

任何模型中都可以使用 RANSAC 模块，在使用 RANSAC 模块时，我们只需要在相应 Python 类中实现 fit() 和 get_error() 方法，剩下就是正确地使用 ransac.py。我们这里使用可能的对应点集来自动找到用于全景图像的单应性矩阵。
可以使用 SIFT 特征自动找到匹配对应
代码如下：

featname = ['D:/SHIJUE/SY5/jmu/p' + str(i + 1) + '.sift' for i in range(5)]  # 图片路径记得修改
imname = ['D:/SHIJUE/SY5/jmu/p' + str(i + 1) + '.jpg' for i in range(5)]

# extract features and match
l = {}
d = {}
for i in range(5):
    sift.process_image(imname[i], featname[i])
    l[i], d[i] = sift.read_features_from_file(featname[i])

matches = {}
for i in range(4):
    matches[i] = sift.match(d[i + 1], d[i])

可以得到如下一组结果：




显然，并不是所有图像中的对应点对都是正确的。实际上，SIFT 是具有很强稳健性的描述子，能够比其他描述子，例如图像块相关的 Harris 角点，产生更少的错误的匹配。但是该方法仍然远非完美

可以将 RANSAC 算法应用于对应点对上:

def convert_points(j): 
 ndx = matches[j].nonzero()[0]
 fp = homography.make_homog(l[j+1][ndx,:2].T)
 ndx2 = [int(matches[j][i]) for i in ndx]
 tp = homography.make_homog(l[j][ndx2,:2].T)
 return fp,tp
# 估计单应性矩阵
model = homography.RansacModel()
fp,tp = convert_points(1)
H_12 = homography.H_from_ransac(fp,tp,model)[0] # im1 到 im2 的单应性矩阵
fp,tp = convert_points(0)
H_01 = homography.H_from_ransac(fp,tp,model)[0] # im0 到 im1 的单应性矩阵
tp,fp = convert_points(2) # 注意：点是反序的
H_32 = homography.H_from_ransac(fp,tp,model)[0] # im3 到 im2 的单应性矩阵
tp,fp = convert_points(3) # 注意：点是反序的
H_43 = homography.H_from_ransac(fp,tp,model)[0] # im4 到 im3 的单应性矩阵

在该例子中，图像 2 是中心图像，也是我们希望将其他图像变成的图像。图像 0 和图像 1 应该从右边扭曲，图像 3 和图像 4 从左边扭曲。在每个图像对中，由于匹配是从最右边的图像计算出来的，所以我们将对应的顺序进行了颠倒，使其从左边图像开始扭曲。因为我们不关心该扭曲例子中的正确点对，所以仅需要该函数的第一个输出（单应性矩阵）。

---

实现

图像集

同一场景下不同拍摄点拍摄的图像集合

源码

# -*- coding: utf-8 -*-
from pylab import *
from numpy import *
from PIL import Image

# If you have PCV installed, these imports should work
from PCV.geometry import homography, warp
from PCV.localdescriptors import sift

np.seterr(invalid='ignore')
"""
This is the panorama example from section 3.3.
"""

# set paths to data folder
featname = ['D:/SHIJUE/SY5/jmu/p' + str(i + 1) + '.sift' for i in range(5)]  # 图片路径记得修改
imname = ['D:/SHIJUE/SY5/jmu/p' + str(i + 1) + '.jpg' for i in range(5)]

# extract features and match
l = {}
d = {}
for i in range(5):
    sift.process_image(imname[i], featname[i])
    l[i], d[i] = sift.read_features_from_file(featname[i])

matches = {}
for i in range(4):
    matches[i] = sift.match(d[i + 1], d[i])

# visualize the matches (Figure 3-11 in the book)
for i in range(4):
    im1 = array(Image.open(imname[i]))
    im2 = array(Image.open(imname[i + 1]))
    figure()
    sift.plot_matches(im2, im1, l[i + 1], l[i], matches[i], show_below=True)


# function to convert the matches to hom. points
def convert_points(j):
    ndx = matches[j].nonzero()[0]
    fp = homography.make_homog(l[j + 1][ndx, :2].T)
    ndx2 = [int(matches[j][i]) for i in ndx]
    tp = homography.make_homog(l[j][ndx2, :2].T)

    # switch x and y - TODO this should move elsewhere
    fp = vstack([fp[1], fp[0], fp[2]])
    tp = vstack([tp[1], tp[0], tp[2]])
    return fp, tp


# estimate the homographies
model = homography.RansacModel()

fp, tp = convert_points(1)
H_12 = homography.H_from_ransac(fp, tp, model)[0]  # im 1 to 2

fp, tp = convert_points(0)
H_01 = homography.H_from_ransac(fp, tp, model)[0]  # im 0 to 1

tp, fp = convert_points(2)  # NB: reverse order
H_32 = homography.H_from_ransac(fp, tp, model)[0]  # im 3 to 2

tp, fp = convert_points(3)  # NB: reverse order
H_43 = homography.H_from_ransac(fp, tp, model)[0]  # im 4 to 3

# warp the images
delta = 1600  # for padding and translation

im1 = array(Image.open(imname[1]), "uint8")
im2 = array(Image.open(imname[2]), "uint8")
im_12 = warp.panorama(H_12, im1, im2, delta, delta)

im1 = array(Image.open(imname[0]), "f")
im_02 = warp.panorama(dot(H_12, H_01), im1, im_12, delta, delta)

im1 = array(Image.open(imname[3]), "f")
im_32 = warp.panorama(H_32, im1, im_02, delta, delta)

im1 = array(Image.open(imname[4]), "f")
im_42 = warp.panorama(dot(H_32, H_43), im1, im_32, delta, 2 * delta)

figure()
imshow(array(im_02, "uint8"))
axis('off')
savefig("example1.png", dpi=300)
show()

拼接后的全景图

可以看出拼接后的全景图存在图像边缘过渡不自然的情况，可能由于原始图像存在亮度差异，导致拼接缝明显。对于图片拼接不自然问题，应进行图片融合。理论上可使用Graph cuts方法解决。

其他测试

突发奇想，试试将一张图像分成几部分进行重新拼接，比较与原始图像的差异

原图：

分成5个图像：

拼接结果：

该拼接的图像与原始图像相较，除部分图像丢失外无明显差别； 与上一拼接图像相较，不存在明显的拼接缝

---

总结

全景图像拼接通过SIFT得到特征匹配，通过RANSAC筛选正确的特征匹配得到透视矩阵，利用透视变换完成两幅图像的拼接，循环这个流程迭代所有图像就可以完成多张图像的全景图的简单拼接。