基于特征匹配和迭代优化的航拍图像拼接

Posted 2023-03-22

参考技术A 本文研究了无人机（UAV）遥感图像拼接过程中重叠区域的不匹配问题。为了解决这个问题，首先通过将双重匹配与随机抽样共识（RANSAC）方法相结合来过滤特征点。其次，为了保证每幅图像与全景照片的投影关系的一致性，我们提出了一种局部拼接的方法。为了避免随着图像数量的增加透视变化累积而导致图像倾斜的错误，我们建立了图像旋转坐标系，并将图像之间的关系限制为平移和旋转。用坐标原点的相对位置来表示平移距离，通过迭代求解最优旋转角度。最后，图像的重叠部分通过线性加权融合。通过实验结果验证，本文提出的方法在大量图像的情况下能够保证更快的处理速度和更高的处理精度，从而达到理想的拼接效果。

近年来，随着计算机视觉的不断进步，图像拼接技术在海洋和矿产勘探、遥感勘探、医学成像、效果生成、虚拟现实等方面得到了广泛的应用。许多航拍遥感图像可以通过配备摄像头的无人机在地面拍摄得到。通常，由于无人机飞行高度、相机焦距等因素，单幅图像存在信息量少、全局分辨率低等问题。因此，要获得广角高分辨率的照片，就需要研究全景图像拼接技术。Brown 在 2003 年引入了著名的 AutoSitich 算法，很快就被用于商业产品，如 Photoshop。但是，该算法假定图像的重叠区域没有深度变化。2013 年，萨拉戈萨 J 等人。将图像拆分为密集的网格，并为每个网格使用单个更改，称为网格变形。该方法在一定程度上解决了图像变形、尺寸缩放、重定向等问题。

图像拼接技术一般分为图像几何校正、图像预处理、图像对齐、图像融合四个步骤。由于相机镜头的畸变，需要对无人机的图像进行校正，使得到的图像满足个别地图的投影关系。图像预处理是几乎所有图像处理技术的重要组成部分，包括去噪、灰度变化等。这个过程可以降低匹配难度，提高匹配精度。然而，对于无人机遥感图像的拼接，图像匹配和图像融合是成功的关键。

图像匹配技术是图像拼接的基础。1975年米尔格拉姆提出了计算机拼接技术。于是，在重叠区域寻找最优接缝线就成为一个重要的研究方向。同年，Kuslin 提出了一种相位相关方法，通过傅里叶变化将图像转换到频域，并利用功率谱计算平移。1987 年，Reddy 提出了一种扩展的相位相关方法，该方法可以计算图像的平移和旋转关系并解决图像缩放问题。图像拼接的另一个分支是基于图像特征。1988年Harris提出经典的Harris点检测算法，它使用特定的旋转不变性哈里斯点进行特征匹配。2004 年，Lowe 提出了一种完美的尺度不变特征变换算法（SIFT），对平移、旋转、尺度缩放、不均匀光照等图像领域应用最广泛的技术具有良好的适应性。C Aguerrebere 根据输入图像的 SNR 条件给出的问题难度级别显示不同的行为区域。Wu通过建立模型，将深度学习和进化算法应用于遥感图像的拼接，实现概率意义上的全局优化。

图像融合技术是遥感图像拼接技术中的另一项核心技术，分为像素级融合、特征级融合、决策级融合。像素级融合仍然是现阶段最常用的图像融合方法之一。

对于无人机的遥感图像，存在图像数量多、光照条件多变等问题。每次拼接过程中的小错误都难以避免。随着图像数量的增加，误差不断累积，图像拼接后期会出现图像失真和重影。S Bang 创建高质量全景图，过滤掉视频的模糊帧，选择关键帧，并校正相机镜头失真。Zhang 提出了基于 STIF 的 GA-SIFT 并给出了一种自适应阈值方法来解决计算量大和拼接时间长的问题。李明基于动态规划解决无人机侧视问题寻找最佳接缝线。然而，当图像数量逐渐增加时，现有的拼接算法存在误差累积。

也有一些基于网格变形的图像拼接算法，但计算量太大。在本文中，图像被匹配两次以过滤特征点以提高准确度。拼接问题对应于通过坐标系转换的旋转角度解，应用高斯-牛顿迭代计算最优旋转角度。此外，我们练习局部匹配方法以减少错误并使用加权融合来实现过度平滑。

SIFT特征点不仅在空间尺度和旋转上保持不变，而且在光照和视角变化的条件下，还具有优异的抗干扰能力和良好的稳定性。为了实现空间尺度的不变性，SIFT特征点可以根据物体远看小而模糊，反之大而清晰的特点，建立高斯金字塔模型。差分金字塔 (DoG) 是通过计算金字塔中相邻两层图像之间的差异来获得的。使用函数拟合在 DOG 空间中测试极值。通过对确定场中基于SIFT特征点的梯度信息进行统计，选择加权幅度最大的梯度方向作为主梯度方向。通过将特征点与其主梯度方向相关联，可以解决图像特征点的旋转不变性问题。最后，利用特征点周围像素的信息建立一个128维的向量作为特征点的描述符。

提取特征点后，需要对两幅图像的特征点进行匹配。通过特征点成对匹配，可以计算出两个特征点对应的描述符之间的欧氏距离，选择欧氏距离最小的点作为匹配点对。为了减少不匹配的发生率， 被用作正确匹配的阈值。具有大于 的描述符欧几里得距离的匹配点对被消除。

RANSAC 是特征点匹配中最常用的方法之一。它首先从匹配结果中随机选择四对特征点并计算单应矩阵。其次，根据上一步得到的单应矩阵，计算第一幅图像在第二幅图像中的重投影坐标，并计算该坐标与第二幅图像中匹配点对坐标的距离。通过设置距离阈值记录所有匹配点对中正确匹配特征点对的个数。重复上面的过程，最终留下与最多点对数的正确匹配。

高斯-牛顿迭代是求解非线性最小二乘优化问题的算法之一，可以描述为：

我们选择一个初始值，然后不断更新当前优化变量以减小目标函数值。高斯-牛顿迭代的主要思想是对函数 进行一阶泰勒展开，计算 及其雅可比矩阵 对应的函数值。使用 和 计算 的增量，直到 足够小。

加权平均法是图像融合中简单有效的方法之一。第一幅图像和第二幅图像重叠区域的像素值由两幅图像像素的加权求和得到，表示为：

其中：越接近 img1， 的值越大。 的值从1逐渐变为0，重叠区域从第一幅图像逐渐过渡到第二幅图像，从而实现画面的平滑过渡。照片的加权平均融合因其直观的简单性和快速的运行速度而被广泛使用和图像拼接。

对于两幅图像的拼接，由于无人机的遥感相机通常安装在一个稳定的平台上，通过选择合适的坐标系，将图像对齐问题转化为单幅图像旋转问题，如图1所示。

此外，大多数具有相关高光频的常用相机通常在连续帧之间具有较大的重叠区域。因此，在图像拼接过程中，第 幅图像在全局位置上的投影关系，不仅受第 幅图像的影响，还与 图像相关。为了保证图像变形的一致性，首先将 张图像拼接在一起，然后将结果整合到整幅图像中。大量的实验测试证明，当i设置为3时效果最佳。整个过程如图图2。

图像中的特征点有很多种，本文使用最常见的SIFT特征点。我们提取并匹配两张输入图像的特征点，结果如下所示。

特征点的匹配精度直接影响旋转角度的计算，因此使用前必须对特征点对进行过滤。鉴于过滤特征点的方法很多，本文先将左图与右图进行匹配，再将右图与左图进行匹配。两次相同结果配对的匹配点将被保留。在此基础上，使用RANSAC方法对结果进行优化，成功匹配了上图中的121个特征点。

从无人机拍摄的两张照片之间通常存在旋转和平移。为了独立优化旋转角度，我们首先建立如图 5所示的坐标系。

以图像匹配成功的特征点坐标值的平均值作为该坐标的原点，坐标轴与像素坐标系的两个坐标轴平行。根据公式（3），特征点从图像坐标系转换为图像旋转坐标系：

其中 为滤波后的特征对的总数， 为特征点在原始图像坐标系中的坐标值，并且 是新的值。

在计算图像的旋转角度之前，我们首先需要分析图像的缩放比例。由飞行高度引起的尺寸变化将在轴上具有相同的缩放比例。因此，根据所有特征点与图像旋转坐标系原点的欧氏距离比，可以计算出两幅图像之间的缩放比例，对图像进行缩放和改变。

图像缩放后，计算图像旋转的角度。高斯牛顿迭代的方式计算旋转角度的最优解。首先设置目标函数：

通过迭代选择最优的 使得：

使用误差函数 的泰勒展开进行迭代。

其中

根据

我们可以发现增量值 每次迭代。最终，当我们计算出的 满足条件时，停止迭代过程。可以使用最佳旋转角度和旋转中心来求解图像的变换矩阵。

由于拍摄图像时光线不均匀，连续两张图像之间可能存在一些颜色差异。此外，图像旋转不可避免地存在小误差，因此我们练习线性加权融合以消除两幅图像之间的拼接线和色度变化。图像的重叠是按距离加权的，这样拼接结果自然是从img1到img2过度了。

我们利用OpenCV的功能从遥感图像中提取SIFT特征点并进行匹配。从Stitch拼接功能、基于透视变化的图像拼接结果以及本文的拼接速度的对比可以看出，本文采用的方法具有一定的优越性。

从表1数据可以看出，在拼接少量图像时，三种算法的拼接结果相似，没有出现明显的拼接误差。但是，Stitcher 算法比其他两种拼接方法花费的时间要多得多。

图 11很明显，随着图像数量的增加，基于透视变换的图像拼接算法出现了严重的失配。然而，本文采用的方法取得了比较满意的结果，因为在无人机拍摄的图像中，地面上的所有特征都可以近似地视为在同一平面上。根据透视变换，无人机的远近抖动会引入图像拼接导致错误。图像数量的不断增加会导致错误的积累，从而导致严重的失配。另外，这使得程序中断，从而无法完成所有60幅图像的拼接。假设同一平面上图片的仿射变化会更符合无人机遥感图像的实际情况。最后，可以通过线性加权融合来解决误差问题，以提高拼接效果。考虑到stitch算法耗时过长，本文不会对两者进行比较。

在上面的图 12 中，使用 100 张图像来测试本文中的方法。图像的仿射变换是通过计算围绕图像特征点中心的旋转角度来进行的。变换后的图像采用线性加权融合后，可以得到大量图像数据处理后的结果。拼接自然，符合人类视觉体验。

我们在网络上跑了一组数据，结果如下。

鉴于以上实验结果，该方法具有一定的抗干扰能力，可以高速运行。与高度集成的Stitcher和基于透视变换的图像拼接结果相比，我们可以发现，基于透视变化的图像拼接结果随着图像数量的增加而逐渐变差。然而，尽管拼接效果很好，但 Stitcher 需要更长的处理时间。

在本文中，我们研究了无人机遥感图像的拼接技术，主要贡献可以总结如下：

通过实验结果可以看出，本文提出的方法比现有方法具有更好的实时性，对于相机平面与成像平面平行的情况具有更好的拼接效果。

计算机视觉项目实战-基于特征点匹配的图像拼接

😊😊😊欢迎来到本博客😊😊😊
本次博客内容将继续讲解关于OpenCV的相关知识
🎉作者简介：⭐️⭐️⭐️目前计算机研究生在读。主要研究方向是人工智能和群智能算法方向。目前熟悉深度学习（keras、pytorch、yolo），python网页爬虫、机器学习、计算机视觉（OpenCV）、群智能算法。然后正在学习深度学习的相关内容。以后可能会涉及到网络安全相关领域，毕竟这是每一个学习计算机的梦想嘛！
📝目前更新：🌟🌟🌟目前已经更新了关于网络爬虫的相关知识、机器学习的相关知识、目前正在更新计算机视觉-OpenCV的相关内容。
💛💛💛本文摘要💛💛💛

本文我们将继续讲解计算机视觉领域项目-基于特征点匹配的图像拼接。 文章目录 🌟项目前言 🌟项目讲解前期准备 ⭐️图像特征检测Harris原理 ⭐️图像特征检测SIFT原理 ⭐️图像特征匹配实战 🌟项目详解

🌟项目前言

之前我们介绍过基于OpenCv的特征匹配操作，我们通过特征匹配可以精确的找到目标。本节我们继续探索基于特征匹配还可以做哪些事情。我们都在拍一个集体的过程中使用过苹果手机的全图效果进行拍照留念。那么苹果手机这个效果它是基于什么技术来做的呢？没错其实就是特征匹配。他是实时拍取多个照片，然后使用特征匹配操作继续两个图像之间特征点的匹配，然后生成转换矩阵，最后转换成效果图，我们本次博客就是要介绍一下这个操作如何使用OpenCv进行实现。

🌟项目讲解前期准备

⭐️图像特征检测Harris原理

• 角点：在图像的角度来看，无论是沿着水平方向还是竖直方向进行移动时候，灰度级会发生变化，而且这个变化是非常迅速的，我们称这个图像就是一个角点。
• 边界：在图像的角度上来看，指把图像按照水平或者垂直方向移动的时候，只有一个一个方向变化的比较明显，另一个方向变化就比较微弱，这样的我们就称之为边界。
• 平面：平面就是说无论向垂直方向还是水平方向移动，图像的灰度值都不会发生迅速的变化，这个就是平面。

原理图，这里第一个图表示的就是平面灰度值没有明显变化，第二个图就是要给边界灰度值水平方向变化明显垂直方向灰度值变化并不明显，第三个图表示的就是一个角点，无论水平还是垂直方向都很明显。主要看灰度级的变化结果：

边界：一个特征值大，一个特征值小，自相关函数在某一个方向上大，在其他方向上小。
平面：两个特征都小，且近似相等。
角点：两个特征都大，且近似相等，自相关函数在所有方向都大。

在OpenCV当中我们使用，cv2.cornerHarris()来进行角点检测。
其中参数都有：

• img： 数据类型为 float32 的入图像
• blockSize： 角点检测中指定区域的大小
• ksize： Sobel求导中使用的窗口大小
• k： 取值参数为 [0,04,0.06]
  我们使用几个小图像去做一下角点检测：
  

import cv2 
import numpy as np

img = cv2.imread('white-black.webp')
print ('img.shape:',img.shape)
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
# gray = np.float32(gray)
dst = cv2.cornerHarris(gray, 2, 3, 0.04)
print ('dst.shape:',dst.shape)
img[dst>0.01*dst.max()]=[0,0,255]
cv2.imshow('dst',img) 
cv2.waitKey(0) 
cv2.destroyAllWindows()

这里有一些点检测的不是特别好，然后我们用黑白棋盘来看一下。

这个效果堪称完美！！！

⭐️图像特征检测SIFT原理

SIFT是指尺度空间：是指在一定的范围内，无论物体是大是小，人眼都可以进行一个识别，然后计算机要去识别却很难，所以要让计算机能够对物体进行一个在不同尺度下都存在一个统一的认知，就要考虑图像在不同的尺度在都存在的特点，尺度空间的获取一般使用高斯模糊来实现-高斯滤波。
我们再说说SIFT有什么优点：
1、具有较好的稳定性和不变性，能够适应旋转、尺度缩放、亮度的变化，能在一定程度上不受视角变化、仿射变换、噪声的干扰。
2、区分性好，能够在海量特征数据库中进行快速准确的区分信息进行匹配
3、多量性，就算只有单个物体，也能产生大量特征向量
4、高速性，能够快速的进行特征向量匹配
5、可扩展性，能够与其它形式的特征向量进行联合
我们都知道如果图片经过高斯滤波操作之后呢，他会变模糊，那么为什么要这么做呢？因为当我们从很近的看一个人的时候，他是清晰的，那么从很远看的时候他就是模糊的。所以我们为了模拟这个过程，就用高斯滤波来进行相同了一个模拟。
不同σ的高斯函数决定了对图像的平滑程度，越大的σ值对应的图像越模糊。
因此我们要介绍一个金字塔，高斯差分金字塔。
我们需要做一个多分辨率的金字塔，对于金字塔的每一层都要做高斯滤波。

根据意思就是5个输入的高斯图像，相邻的进行像素值相减，得到4张差分后的结果。那么我们想要找什么呢？找SIFT，就是特征点，那么什么样的点被认为是特征点呢？通常关键点数值较大、差分结果较大的、极值里面较大的。是不是有点像之前讲的图像金字塔。

⭐️图像特征匹配实战

🌟项目详解

前面我们简单的将特征匹配介绍了一下，然后我们开始讲一下我们本次博客的内容。
主函数：

from Stitcher import Stitcher
import cv2
# 读取拼接图片
imageA = cv2.imread("left_02.jpg")
imageB = cv2.imread("right_02.jpg")
# 把图片拼接成全景图
stitcher = Stitcher()
(result, vis) = stitcher.stitch([imageA, imageB], showMatches=True)
# 显示所有图片
cv2.imshow("Image A", imageA)
cv2.imshow("Image B", imageB)
cv2.imshow("Keypoint Matches", vis)
cv2.imshow("Result", result)
cv2.waitKey(0)
cv2.destroyAllWindows()

在主函数这里我们把两张图像导入进来，注意这里的图像宽度必须一致，因为后续要进行拼接。长度可以保持不一致。然后我们进入拼接全景图的操作当中。
Stitcher()函数部分：

import numpy as np
import cv2

class Stitcher:
    #拼接函数
    def stitch(self, images, ratio=0.75, reprojThresh=4.0,showMatches=False):
        #获取输入图片
        (imageB, imageA) = images
        (kpsA, featuresA) = self.detectAndDescribe(imageA)
        (kpsB, featuresB) = self.detectAndDescribe(imageB)
        M = self.matchKeypoints(kpsA, kpsB, featuresA, featuresB, ratio, reprojThresh)
        # 如果返回结果为空，没有匹配成功的特征点，退出算法
        if M is None:
            return None
        # 否则，提取匹配结果
        # H是3x3视角变换矩阵      
        (matches, H, status) = M
        # 将图片A进行视角变换，result是变换后图片
        result = cv2.warpPerspective(imageA, H, (imageA.shape[1] + imageB.shape[1], imageA.shape[0]))
        #self.cv_show('result', result)
        # 将图片B传入result图片最左端
        result[0:imageB.shape[0], 0:imageB.shape[1]] = imageB
        #self.cv_show('result', result)
        # 检测是否需要显示图片匹配
        if showMatches:
            # 生成匹配图片
            vis = self.drawMatches(imageA, imageB, kpsA, kpsB, matches, status)
            # 返回结果
            return (result, vis)
        # 返回匹配结果
        return result

首先我们将程序写成一个类，然后在类中使用self.detectAndDescribe()函数检测A、B图片的SIFT关键特征点，并计算特征描述子。

    def detectAndDescribe(self, image):
        # 将彩色图片转换成灰度图
        gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)
        descriptor = cv2.SIFT_create()
        (kps, features) = descriptor.detectAndCompute(image, None)
        kps = np.float32([kp.pt for kp in kps])
        return (kps, features)

首先我们将图像转换成了灰度图，然后我们建立了使用SIFT检测特征点的检测器。然后对两张图片进行特征点检测，然后将结果转化成numpy数组并且返回。
然后在使用self.matchKeypoints()对特征点进行匹配。

    def matchKeypoints(self, kpsA, kpsB, featuresA, featuresB, ratio, reprojThresh):
        # 建立暴力匹配器
        matcher = cv2.BFMatcher()
        # 使用KNN检测来自A、B图的SIFT特征匹配对，K=2
        rawMatches = matcher.knnMatch(featuresA, featuresB, 2)
        matches = []
        for m in rawMatches:
            # 当最近距离跟次近距离的比值小于ratio值时，保留此匹配对
            if len(m) == 2 and m[0].distance < m[1].distance * ratio:
            # 存储两个点在featuresA, featuresB中的索引值
                matches.append((m[0].trainIdx, m[0].queryIdx))
        # 当筛选后的匹配对大于4时，计算视角变换矩阵
        if len(matches) > 4:
            # 获取匹配对的点坐标
            ptsA = np.float32([kpsA[i] for (_, i) in matches])
            ptsB = np.float32([kpsB[i] for (i, _) in matches])
            # 计算视角变换矩阵
            (H, status) = cv2.findHomography(ptsA, ptsB, cv2.RANSAC, reprojThresh)
            # 返回结果
            return (matches, H, status)
        # 如果匹配对小于4时，返回None
        return None

通过上述代码我们可以计算出来两张图象的视角变换矩阵，然后返回结果。最后进行可视化操作。

    def drawMatches(self, imageA, imageB, kpsA, kpsB, matches, status):
        # 初始化可视化图片，将A、B图左右连接到一起
        (hA, wA) = imageA.shape[:2]
        (hB, wB) = imageB.shape[:2]
        vis = np.zeros((max(hA, hB), wA + wB, 3), dtype="uint8")
        vis[0:hA, 0:wA] = imageA
        vis[0:hB, wA:] = imageB

        # 联合遍历，画出匹配对
        for ((trainIdx, queryIdx), s) in zip(matches, status):
            # 当点对匹配成功时，画到可视化图上
            if s == 1:
                # 画出匹配对
                ptA = (int(kpsA[queryIdx][0]), int(kpsA[queryIdx][1]))
                ptB = (int(kpsB[trainIdx][0]) + wA, int(kpsB[trainIdx][1]))
                cv2.line(vis, ptA, ptB, (0, 255, 0), 1)

        # 返回可视化结果
        return vis

我们来看一下效果图：
左一图：

右一图：

最终效果图：

在用一个经典案例：

🔎支持：🎁🎁🎁如果觉得博主的文章还不错或者您用得到的话，可以免费的关注一下博主，如果三连收藏支持就更好啦！这就是给予我最大的支持！