图像映射与全景拼接

图像拼接(image mosaic)是一个日益流行的研究领域，他已经成为照相绘图学、计算机视觉、图像处理和计算机图形学研究中的热点。图像拼接解决的问题一般式，通过对齐一系列空间重叠的图像，构成一个无缝的、高清晰的图像，它具有比单个图像更高的分辨率和更大的视野。

图像拼接基本流程

1.基础流程

①针对某个场景中拍摄多张/序列图像
②计算第二张图像与第一张图像之间的变换关系
③将第二张图像叠加到第一张图像的坐标系中
④变换后的融合/合成
⑤在多图场景中，重复上述过程

2.计算图像间的变换结构

①针对两张及多张图片提取特征
②特征匹配
找出输入的两张图片（或者多张）的特征点（SIFT特征点）然后再匹配上这些特征点。
SIFT相关原理见（链接: link.）

③根据图像变换特点，选取合适的变换结构
2D图像变换，详细变换原理见下方

④根据DLT等方法计算变换结构

DLT(Direct Linear Transform)算法
直接线性变换解法是建立像点的“坐标仪坐标”和相应物点的物方空间坐标直接的线性关系的解法。
直接线性变换解法的特点：
1、不归心、不定项
2.不需要内外方位元素的起始值
3.物方空间需布置一组控制点
4.特别适合于处理非量测相机所摄影像
5.本质是一种空间后交-欠交解法

图像拼接相关原理

1.2D图像变换原理

根据图像的变换特点，我们要选择不同的变换结构，通常，2D图像的变换原理有以下的几种类型：

2.特征点错误匹配干扰

在进行图像特征点匹配过程中，难免会出现特征点匹配错误的问题，如果使用错误的特征点进行图像拼接，可能会导致比较大的误差。

为此，我们可以采用曲线拟合的方法，在匹配的特征点中找到一条合适的曲线尽量消除噪声点以减小误差。

方式：
1.直线拟合
在空间中给定若干点（已匹配的特征点），任取两点构成一条直线，给定一个阈值，计算inliers的数量（在线上的点的数量），选取最多的那条直线为拟合直线。

2.圆拟合
两点确定一条直线，三点则能确定一个圆，原理同直线拟合。

3.曲线拟合
给定曲线多项式，求解多项式的系数。

3.APAP算法

产生“鬼影”的原因在于图像映射是全局的单应性变换，然后图像场景中各个物体往往具有不同的深度，如果采用处于不同深度物体的特征点进行全局单应性变换，由于此时图像中的物体无法满足近似于同一平面的条件，计算得到的单应性矩阵会有较大的误差，仅仅由一个全局的单应性变换无法完全描述两幅图像之间的变换关系。
在图像拼接融合的过程中，受客观因素的影响，拼接融合后的图像可能会存在上图中的“鬼影现象”以及图像间过度不连续等问题。2013年，Julio Zaragoza等人发表了一种新的图像配准算法APAP（As-Projective-As-Possible Image Stitching with Moving DLT）解决鬼影现象可以采用APAP算法。

算法流程：
1.提取两张图片的SIFT特征点
2.对两张图片的特征点进行匹配
3.匹配后，使用RANSAC算法进行特征点对的筛选，排除错误点。筛选后的特征点基本能够一一对应。
4.使用DLT算法，将剩下的特征点对进行透视变换矩阵的估计。
5.因为得到的透视变换矩阵是基于全局特征点对进行的，即一个刚性的单应性矩阵完成配准。为提高配准的精度，APAP将图像切割成无数多个小方块，对每个小方块进行单应性矩阵变换。

APAP虽然能够较好地完成配准，但非常依赖于特征点对。若图像高频信息较少，特征点对过少，配准将完全失效，并且对大尺度的图像进行配准，其效果也不是很好，一切都决定于特征点对的数量。

4.RANSAC方法图像拼接

RANSAC是“RANdom SAmple Consensus（随机抽样一致）”的缩写。它可以从一组包含“局外点”的观测数据集中，通过迭代方式估计数学模型的参数。它是一种不确定的算法——它有一定的概率得出一个合理的结果；为了提高概率必须提高迭代次数。

RANSAC的基本假设是：
（1）数据由“局内点”组成，例如：数据的分布可以用一些模型参数来解释；
（2）“局外点”是不能适应该模型的数据；
（3）除此之外的数据属于噪声。
局外点产生的原因有：噪声的极值；错误的测量方法；对数据的错误假设。
RANSAC也做了以下假设：给定一组（通常很小的）局内点，存在一个可以估计模型参数的过程；而该模型能够解释或者适用于局内点。

算法基本步骤

1.选择4对匹配特征点（选择4对特征点因为单应性矩阵有8个自由度，需要4个特征点，每个特征点都能构造2个方程，然后求解单应性矩阵，但是任意3个点或4个点不能在同一条直线上）
2.根据直接线性变换解法DLT计算单应性矩阵H
3.对所匹配点，计算映射误差
4.根据误差阈值，确定inliers数量
5.针对最大的inliers集合，重新计算单应性矩阵H

5.最佳拼接缝

由于图像重叠会导致鬼影问题，所以可以采取找最佳拼接缝的方法，找到两张图片之间最佳的一条拼接线，左边是一张图，右边是另一张图，就很好的解决了鬼影问题。
方法：最大流最小割
最大流
给定指定的一个有向图，其中有两个特殊的点源S(Sources)和汇T(Sinks)，每条边有指定的容量(Capacity)，求满足条件的从S到T的最大流(MaxFlow)。
最小割
割是网络中定点的一个划分，它把网络中的所有顶点划分成两个顶点集合S和T，其中源点s∈S,汇点t∈T。记为CUT(S,T)，满足条件的从S到T的最小割（Min cut）。

不断的找从S到E的一条路径，每找到一条路径，就让该路径上的所有边的权重减去这条路径上权重值最小边的权重，同时在这条路径上构造一个方向的有向边，权重值设置为减去的值，另外如果减完后权重为0，则删除这条边。算法一直循环知道无法找到S到E的路径借结束，最后切割那些为空或是流已满的边

参考：link

6.代码部分

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

"""
This is the panorama example from section 3.3.
"""


featname = ['D:/pythonCode/test/data/testimages/' + str(i + 1) + '.sift' for i in range(5)]
imname = ['D:/pythonCode/test/data/testimages/' + str(i + 1) + '.jpg' for i in range(5)]
# extract features and m
# 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

# 扭曲图像
delta = 100  # 用于填充和平移 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_42, "uint8"))
axis('off')
show()

图像拼接代码测试

1.测试一
测试图片

运行结果：

2.测试二

测试图片

运行结果：

结果分析

1.进行了7组图片的测试，上方结果只展示了两种，其他五组效果较差，分析原因可能是因为图片拍摄时角度和位置变化较大，导致特征匹配时匹配点较少，使得拼接结果较差，可以间接说明该代码在进行角度和位置同时改变时的适应能力较差。
2.在除了图片的问题以外，一开始运行代码时，出现了以下的运行结果，原因是代码中delta（用于平移和填充的数据）较大，导致图片结果较长，可以调节小一些。

3.在进行测试时，有一组图片拼接结果最差（下图）

从结果来看，基本就是乱拼，看了一下几张特征匹配的图，特征匹配的效果还不错（错误的匹配点较少，但总的匹配数量较少），分析原因可能是图片的建筑物差异点（特征点）不够明显，加上数量较少，多图直接拼接容易弄错。

4.运行过程中的错误：
（1）：

原因：图片较大
修改方式：①更改图片大小
②在系统环境变量中调节虚拟内存大小
（2）：

原因：图片大小不同
修改方式：图片大小规格改为相同

5.对于第二组测试的结果：
第二组图片是用别人拍的照片进行测试的（自己拍的效果都比较差QAQ）
自己的运行结果

但是在别人的电脑上运行的结果

可以看出，拼接出的分割线基本相同，但前一个的分割线有明显的的图线划分，原因目前只能解释不同电脑也有所差异（ε=(´ο｀*)）之后如果清楚原因会再补充。

总结，该算法大体上较为成熟，拼接时得到的效果较好，只是对拍摄图片的角度和位置不能有太大的变动（类似于手机拍长图，只能转动方位），但在拼接线方面还是有明显的痕迹，对于特征点数量要求较高，需要匹配度很高（但该要求在进行图像拼接时比较合理），即在图片较好的情况下能进行不错的拼接，没有出现“鬼影”。

# coding=utf-8
from PIL import Image, ImageDraw, ImageFont
import cv2
import numpy as np


def jigsaw(imgs, direction="horizontal", gap=0):
    imgs = [Image.fromarray(img) for img in imgs]
    w, h = imgs[0].size
    if direction == "horizontal":
        result = Image.new(imgs[0].mode, ((w+gap)*len(imgs)-gap, h))
        for i, img in enumerate(imgs):
            result.paste(img, box=((w+gap)*i, 0))
    elif direction == "vertical":
        result = Image.new(imgs[0].mode, (w, (h+gap)*len(imgs)-gap))
        for i, img in enumerate(imgs):
            result.paste(img, box=(0, (h+gap)*i))
    else:
        raise ValueError("The direction parameter has only two options: horizontal and vertical")
    return np.array(result)


if __name__ == '__main__':
    img1 = cv2.imread("output/frankfurt_000000_000294_leftImg8bit.png")
    img2 = cv2.imread("output/frankfurt_000000_000576_leftImg8bit.png")
    img3 = cv2.imread("output/frankfurt_000000_001016_leftImg8bit.png")
    img = jigsaw([img1, img2, img3])
    cv2.imwrite("test.png", img)