一、算法目的
在同一位置(即图像的相机位置相同)拍摄两张以上图片,这些图片是单应性相关的,即图片之间有相同的拍摄区域。基于此将图片进行缝补,拼成一个大的图像来创建全景图像。
二、基本原理
要实现两张图片的简单拼接,其实只需找出两张图片中相似的点 (至少四个,因为 homography 矩阵的计算需要至少四个点), 计算一张图片可以变换到另一张图片的变换矩阵 (homography 单应性矩阵),用这个矩阵把那张图片变换后放到另一张图片相应的位置 ,就是相当于把两张图片中定好的四个相似的点重合在一起。如此,就可以实现简单的全景拼接。
三、实现步骤
1、读入连续图片并使用SIFT特征查找匹配对应点对
import sift featname = ['D:/LearnSrc/PythonSrc/Homework/Lecture05/images'+str(i+1)+'.sift' for i in range(5)] imname = ['D:/LearnSrc/PythonSrc/Homework/Lecture05/images'+str(i+1)+'.jpg' for i in range(5)] 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算子存在错误的匹配点,因此需要用Ransac算法剔除错误匹配点。
2、利用RANSAC算法计算变换矩阵
2.1 RANSAC算法
RANSAC是"RANdom SAmple Consensus"(随机一致采样)的缩写。该方法是用来找到正确模型来拟合带有噪声数据的迭代方法。给定一个模型,例如点集之间的单应性矩阵。基本的思想是:数据中包含正确的点和噪声点,合理的模型应该能够在描述正确数据点的同时摒弃噪声点。
求解单应性矩阵:
class RansacModel(object): def __init__(self,debug=False): self.debug = debugdef fit(self, data): """ 计算选取四个对应的单应性矩阵 """ # 将其转置,来调用H_from_points()计算单应性矩阵 data = data.T #映射的起始点 fp = data[:3,:4] # 映射的目标点 tp = data[3:,:4] #计算单应性矩阵然后返回 return H_from_points(fp,tp) def get_error( self, data, H): """ 对所有的对应计算单应性矩阵,然后对每个变换后的点,返回相应的误差 """ data = data.T #映射的起始点 fp = data[:3] # 映射的目标点 tp = data[3:] #变换fp fp_transformed = dot(H,fp) #归一化齐次坐标 for i in range(3): fp_transformed[i] /= fp_transformed[2] return sqrt( sum((tp-fp_transformed)**2,axis=0) )
以上可看出,这个类包含fit()方法,仅接受由ransac算法选择的4个对应点对(data中的前4个点对),然后拟合一个单应性矩阵。get_error()方法对每个对应点对使用该单应性矩阵,然后返回相应的平方距离之和。因此ransac算法能够判定正确与错误的点。在实际中,还需在距离上使用一个阈值来决定合理的单应性矩阵是哪些。
3.拼接图像:
估计出图像间的单应性矩阵后(使用RANSAC算法),需要将所有的图像扭曲到一个公共的图像平面上。通常,这里的公共平面为中心图像平面(否则需要进行大量变形)。一种方法是创建一个很大的图像,比如图像中全部填充0,使其和中心图像平行,然后将所有的图像扭曲到上面。由于我所有的图像是由照相机水平旋转拍摄的,因此可使用一个较简单的步骤:将中心图像左边或右边的区域填充0,以便为扭曲的图像腾出空间。
代码:
from array import arrayfrom numpy import dotfrom numpy.ma import vstackfrom pylab import *from PIL import Image# If you have PCV installed, these imports should workfrom PCV.geometry import homography, warpfrom PCV.localdescriptors import sift"""This is the panorama example from section 3.3."""# set paths to data folderfeatname = ['D:/LearnSrc/PythonSrc/Homework/Lecture05/images/0' + str(i + 1) + '.sift' for i in range(5)] imname = ['D:/LearnSrc/PythonSrc/Homework/Lecture05/images/0' + str(i + 1) + '.jpg' for i in range(5)]# extract features and matchl = {} 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. pointsdef 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 homographiesmodel = homography.RansacModel() fp, tp = convert_points(1) H_12 = homography.H_from_ransac(fp, tp, model)[0] # im 1 to 2fp, tp = convert_points(0) H_01 = homography.H_from_ransac(fp, tp, model)[0] # im 0 to 1tp, fp = convert_points(2) # NB: reverse orderH_32 = homography.H_from_ransac(fp, tp, model)[0] # im 3 to 2tp, fp = convert_points(3) # NB: reverse orderH_43 = homography.H_from_ransac(fp, tp, model)[0] # im 4 to 3# warp the imagesdelta = 2000 # for padding and translationim1 = 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_42, "uint8")) axis('off') savefig("example5.png", dpi=300) show()
四、实验结果
注:图片不能太大,要拼接的图片的分辨率必须都相同。
第一次运行结果:
图片:
以第一张图为开始的结果:
由第三张图为中心拼接的结果:
从上面个结果可以看出,拼接的效果很不好。想过是不是因为拍摄的图片不好,周围建筑太相似而造成这种结果,所以重新拍摄了3张新照片测试:
运行结果:
可以发现右边的拼接效果比左边好很多,回头看前几次的运行结果,同样也是右边的效果比左边好一些。
看代码似乎与拼接的方向有关,但改成从右边往左拼接也没得到比较好的效果。
最后将图片的顺序改成倒序,并以第三张为中心拼接,得到的效果最好....同时也可以看出上面说的太多相似建筑而导致拼接效果不好的想法是不成立的。
从上面的结果可以看到,这种方法可以得到一个比较好的拼接效果。但仍然存在许多问题,比如图片有些拼接部位倾斜拉伸严重,有些地方的拼接缝比较明显。还有上面出现的断开的情况。这可能是由于拍摄时没有站在固定的一个点而造成的结果。
选的图片不是固定点拍摄、以及拼接的顺序不对的话,怎么做拼接效果都会很差。。