CN113222820B - 一种位姿信息辅助的航空遥感图像拼接方法 - Google Patents
一种位姿信息辅助的航空遥感图像拼接方法 Download PDFInfo
- Publication number
- CN113222820B CN113222820B CN202110553778.5A CN202110553778A CN113222820B CN 113222820 B CN113222820 B CN 113222820B CN 202110553778 A CN202110553778 A CN 202110553778A CN 113222820 B CN113222820 B CN 113222820B
- Authority
- CN
- China
- Prior art keywords
- image
- coordinate system
- navigation
- remote sensing
- images
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 46
- 239000011159 matrix material Substances 0.000 claims abstract description 77
- 230000032669 eclosion Effects 0.000 claims abstract description 11
- 238000007500 overflow downdraw method Methods 0.000 claims abstract description 9
- 238000006073 displacement reaction Methods 0.000 claims abstract description 8
- 238000007781 pre-processing Methods 0.000 claims abstract description 6
- 230000009466 transformation Effects 0.000 claims description 26
- 239000013598 vector Substances 0.000 claims description 25
- 238000013519 translation Methods 0.000 claims description 17
- 230000004927 fusion Effects 0.000 claims description 14
- 238000006243 chemical reaction Methods 0.000 claims description 12
- 230000008859 change Effects 0.000 claims description 9
- 230000003044 adaptive effect Effects 0.000 claims description 8
- 230000003287 optical effect Effects 0.000 claims description 7
- 238000004364 calculation method Methods 0.000 claims description 5
- 238000012545 processing Methods 0.000 claims description 5
- 238000012216 screening Methods 0.000 claims description 5
- 238000007499 fusion processing Methods 0.000 claims description 4
- 230000001131 transforming effect Effects 0.000 claims description 4
- 238000005259 measurement Methods 0.000 claims description 3
- 230000004044 response Effects 0.000 claims description 3
- 230000007704 transition Effects 0.000 claims description 3
- 238000007562 laser obscuration time method Methods 0.000 claims 1
- 238000002474 experimental method Methods 0.000 description 8
- 238000004458 analytical method Methods 0.000 description 5
- 238000001514 detection method Methods 0.000 description 3
- 238000010276 construction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000012876 topography Methods 0.000 description 2
- 230000001419 dependent effect Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000007429 general method Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000035800 maturation Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T3/00—Geometric image transformations in the plane of the image
- G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
- G06T3/4038—Image mosaicing, e.g. composing plane images from plane sub-images
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
- G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Image Processing (AREA)
Abstract
本发明公开了一种位姿信息辅助的航空遥感图像拼接方法，涉及航空摄影遥感和图像拼接领域。首先获取航空遥感图像数据和导航位姿数据并进行预处理；然后根据位姿信息求解航带内图像位移量，获得相邻图像重叠区域；在图像重叠区域内提取图像特征进行匹配，求解图像单应性矩阵进行图像配准；分析图像配准质量，若判定配准失效，根据导航数据转换成的相机位姿直接求解单应性矩阵，进行图像配准；最后采用缝合带羽化融合方法进行融合，根据位姿信息进行全景图地理定位。本发明利用机载导航系统的位姿信息辅助航空遥感图像拼接，提高了遥感图像拼接速度，而且能实现特征稀少图像的拼接，最终获得具有地理信息的遥感全景图，具有较大的实用价值。
Description
技术领域
本发明涉及航空摄影遥感和图像拼接领域，具体涉及一种位姿信息辅助的航空遥感图像拼接方法。本发明的主要应用是多场景的航空遥感摄影图像拼接，包括地貌单一、特征稀少的遥感图像的快速拼接。
背景技术
随着航空遥感技术的成熟和相机传感器性能的不断提高，航空摄影遥感图像在生态环境监测、城市规划建设、土地综合利用、地理信息系统建设等领域中广泛应用。受限于高分辨率和大视角的矛盾，遥感图像必须进行拼接才能用于后续研究，其拼接质量直接影响了后续研究及应用的正确性和效果。因此，高性能的遥感图像拼接技术具有广阔的应用场景。现有的遥感图像拼接方法主要有基于特征的图像拼接方法、基于灰度模板匹配的拼接方法、基于频域变换的图像拼接方法等，这些方法计算量都较大，导致遥感图像拼接速度慢；相应的一些改进方法也对硬件依赖大，不适用于一些低成本的应用。另外，现有的图像拼接方法均依赖于图像信息，例如图像特征、图像模板等，导致在一些特殊场景下出现失效。比如在沙漠、水面、草原等场景获取的遥感图像，使用现有的图像拼接算法无法实现拼接。因此，对于复杂的应用场景，现有遥感图像拼接方法的拼接速度受限且适应性较差，引入机载相机位姿信息进行辅助，可以解决传统图像拼接方法依赖单一图像信息的局限性，提高遥感图像拼接速度和适应性。
发明内容
本发明提供一种位姿信息辅助的航空遥感图像拼接方法，以提高航空遥感图像拼接速度，并解决地貌单一、特征稀少的遥感图像拼接易失效的问题。
本发明采取的技术方案为：
一种位姿信息辅助的航空遥感图像拼接方法，具体步骤如下：
S1：遥感图像数据与导航位姿数据的作业获取与预处理；
S2：根据位姿信息计算航带相邻图像的位移量，求解图像重叠区域；
S3：基于SURF算法提取重叠区域图像特征，定位特征点并生成特征描述子后进行特征匹配，求解图像变换单应性矩阵，实现图像配准：
基于SURF算法提取特征点，
根据自适应阈值求取特征匹配，
使用RANSAC算法进行误配筛选，最终求解图像配准单应性矩阵；
S4：分析基于重叠区域特征的图像配准结果的配准质量，若配准失效，则将导航系统位姿信息转换为相机的旋转矩阵和平移向量，直接求解图像间单应性矩阵，实现图像配准；
S5：基于单应性矩阵进行图像变换，使用基于缝合带羽化的图像融合方法进行融合处理；
S6：根据像点和对面点共线关系，将导航位姿信息变换为地面地理坐标，获得包含图像地理信息的遥感图像全景图。
步骤S1中的数据的获取与预处理，包括以下步骤：
(1)对相机传感器和导航系统进行标定实验，包括标定时间差，进行图像、位姿数据频率同步，并标定获得相机的内参数，用于后续拼接；
(2)航空遥感载体搭载相机和导航系统，通过作业获取图像、位姿数据，选择图像数据类型进行存储，并根据导航算法解算获得作业中的经度、纬度、高度和三个姿态角；
(3)使用自适应直方图均衡方法对图像进行增强处理，校正曝光，增加细节。
(4)对由传感器质量、大气、载体姿态、地形起伏、地球自转等因素引起的图像畸变进行校正处理，恢复图像的地貌特征，为遥感图像拼接做准备。
步骤S2中的相邻图像重叠区域求解，包括以下步骤：
(1)根据航带内图像I1和图像I2拍摄时刻的经纬度信息，求解图像I1和图像I2中像主点的地理坐标(X1,Y1)、(X2,Y2)，表明了图像在东向、北向的移动量
(2)已知航向角ψ，根据旋转变换关系，可以求得图像I1和图像I2在图像右上的x、y方向上的移动量。其中图像坐标系
(3)根据已知的飞行航高H、相机焦距f，像元尺寸dxy，可求得图像中一个像素对应地面的尺寸
(4)计算得相邻图像的像元移动量，最后即可根据图像尺寸计算得到相邻图像重叠区域
(5)同一航带内运行速度变化较小，图像x、y方向上的位移偏差不大，因此根据(1)～(4)步骤计算前四组图像得重叠区域，求取平均值，作为整条航线上图片的重叠区域大小。其中，通过航向角的改变幅度，判断航线的变化。
步骤S3中的基于重叠区域特征的图像配准，包括以下步骤：
(1)基于SURF算法提取特征点。包括构建Hessian矩阵，通过判别式检测特征点，舍去对比度低的点和边缘响应点；生成图像尺度空间；定位特征点，构建64维的特征描述算子。
(2)根据自适应阈值求取特征匹配对假设两张图像分别有n1、n2个特征点，图像特征点的距离用欧式距离度量特征向量相似性来表示，求取最近邻距离dmx与次近邻距离dmy。
其中，所有特征点m＝1,…,n1，(xm1,xm2,…,xmn)T为特征点的特征描述向量，n为特征描述向量的维数。
在计算最近邻特征距离和次近邻距离后，定义距离比值
计算图像I1到图像I2的匹配特征，再计算图像I2到图像I1的匹配特征，得到两个匹配集合A、B。根据集合A、B中的距离比值，求解公共子集C。此时公共子集C的最小距离比值即为自适应阈值。此时根据特征点的特征距离比值和阈值可以判断特征点对匹配是否成功。
(3)使用RANSAC算法进行误配筛选，最终求解图像配准的单应性矩阵。随机取出m对(至少4对)不共线的匹配特征点，计算变换矩阵H，记为模型M；计算初始匹配集中的内点，根据全部内点计算模型M的统计误差error，迭代更新模型M和error，最小error的模型为最优模型，此时的H为图像配准单应性矩阵。
步骤S4中的图像配准结果的配准质量分析。首先，设定特征匹配对数量有效阈值，若特征匹配对数量过少，则认定配准失败。其次，设定配准误差阈值，假如特征匹配对数量足够且求解出单应性矩阵，则求解图像配准的配准误差，若配准误差过高，则同样认定为图像配准失败。
步骤S4中的基于位姿信息直接求解单应性矩阵的方法，包括以下步骤：
(1)由姿态数据计算图像I1与图像I2相机坐标系的旋转矩阵R。
遥感平台搭载的导航系统直接输出的数据，是载体的导航信息，包括在WGS-84大地坐标系下的纬度、经度、高度，载体坐标系相对导航坐标系的三个旋转角度：航向角、横滚角、俯仰角。需要将导航位姿数据转换为相机拍摄两张照片的位姿关系，包括旋转运动R矩阵和平移运动t向量。假设相机在两个时刻拍摄两张图像I1、I2，相机光心分别位于点OC1和点OC2。
旋转矩阵R的转换关系为导航坐标系(n)→载体坐标系(b)→相机坐标系(C)。定义导航坐标系为东北天导航坐标系，载体坐标系为右上载体坐标系，x轴指向地理东向、y轴指向地理北向，z轴指向天。相机坐标系为x轴、y轴和图像坐标系保持一致，z轴指向物点。图像I1相机坐标系和图像I2相机坐标系之间的旋转矩阵R为
其中，定义航向角ψ，北偏东为正；俯仰角θ，抬头为正；横滚角γ，右倾为正，且
(2)由姿态数据计算图像I1与图像I2相机坐标系的平移向量t。平移向量t的转换关系为：WGS-84大地坐标系→WGS-84地心空间直角坐标系(e)→导航坐标系(n)→载体坐标系(b)→世界坐标系(C1)→相机坐标系(C2)。WGS-84坐标系原点为地球质心，X轴指向零度子午面和赤道的交点，Z轴指向协议地球极，Y轴由右手定则得到。根据WGS-84大地坐标和空间直角坐标转换关系计算拍摄时刻空间直角坐标X、Y、Z。
其中，卯酉圈半径a为参考椭球体长半轴，b为参考椭球体短半轴；/>
图像I1相机坐标系和图像I2相机坐标系之间的平移向量t为
其中，R为公式中旋转矩阵，
(3)根据位姿信息求解两幅图像的单应性矩阵。以图像I1的相机坐标系作为世界坐标系，即作为图像配准的基准坐标系，图像I1和图像I2之间的单应性矩阵为
其中，由于航拍图像的拍摄高度较高，在地面起伏不大的情况下，h由光心O1的高度表示。
(4)求解航带内多图像间单应性矩阵。将图像I2变换到图像I1平面，得到I1和I2组成的大图Ip1，I1坐标系与拼接结果大图Ip1坐标系仅为平移关系。假设输入图片尺寸为宽w(pixel)，高h(pixel)，则图像I2的四个角点在I2像素坐标系中坐标为A2(0,0)，B2(w,0)，C2(w,h)，D2(0,h)，根据单应性变换到图I1像素坐标系中，假设I1、I2的四个角点在I1像素坐标系中坐标分别为A2,1(uA2,vA2)，B2,1(uB2,vB2)，C2,1(uC2,vC2)，D2,1(uD2,vD2)；A1(uA1,vA1)，B1(uB1,vB1)，C1(uC1,vC1)D1(uD1,vD1)。拼接结果大图Ip1的尺寸宽W、高H为
I1像素坐标系与Ip1像素坐标系的变换关系为
其中，
图像I2与Ip1的变换关系为
图像I3与Ip1的变换关系为
同理图像间单应性矩阵皆可求得，可以求解多图像的单应性矩阵。
步骤S5中的基于缝合带羽化的图像融合方法，包括以下步骤：
(1)根据单应性矩阵计算变换后图像四个端点的坐标，求解缝合线坐标。设置缓冲宽度，得到缝合带区域。
(2)在缝合带区域内运用羽化方法分别对上下边界、左右边界图像像素进行处理。设置过渡步长：依次遍历融合：第n行/列的像素值为n·sI1(x,y)+(1-n·s)I2(x,y)。
其中，D1-2表示缝合线区域即融合区域，其他区域皆保留图像I12。α、β表示融合时的权重，β＝1-α。其中Whalf表示缝合线的半宽度，d∈(0,2Whalf)表示融合时距缝合线区域的初始边界的距离。上下边融合时，在缝合线区域由上向下过渡时d从0变化到2Whalf，缝合带中接缝线上侧的区域为I1、I2的羽化融合，缝合线下侧区域为I1。
步骤S6中的地理坐标求解方法，具体方法为：
以图I1相机坐标系为世界坐标系，变换后的全景图在图I1坐标系平面中，图像中像点像素坐标和地理坐标的关系为
其中，图I1时刻相机光心在地理坐标系的坐标为(XO,YO,0)，求的像点对应坐标(X,Y,Z)。
图I1像素坐标(u0,v0,1)与全景图像素坐标的关系为
经纬度与地理坐标的转换关系，由高斯-克吕格投影得到。根据以上公式，在获得全景图之后，根据像点在全景图的像素坐标，即可求得对应地面点在中国WGS-84公里网的地理坐标，实现对全景图的地理粗定位，实现地面物体的大致定位，对遥感图像的后续分析与检测等提供帮助。
本发明的实现效果和优点包括：
(1)引入机载导航系统的位姿信息，辅助实现了航空遥感图像无缝拼接，解决传统图像拼接方法依赖单一图像信息的局限性；
(2)通过使用位姿信息缩小特征提取空间，提高了遥感图像拼接的速度，降低了拼接误差，遥感图像拼接性能总体得到提升；
(3)在一般方法失效时，本方法根据导航位姿信息得到相机旋转、平移量，直接求解图像单应性矩阵进行拼接，解决了传统图像拼接方法无法拼接特征稀少图像的问题；
(4)拼接得到具有地理信息的全景图，丰富了遥感图像拼接结果的信息，对航空遥感图像的后续分析与检测提供了帮助，有助于应急救援等应用。
附图说明
图1是本发明的位姿信息辅助遥感图像拼接方法的总体流程图；
图2是本发明的基于位姿信息求解图像重叠区域示意图；
图3是本发明的位姿辅助特征匹配的遥感图像拼接方法流程图；
图4是本发明的基于位姿信息图像单应性矩阵的空间模型；
图5是本发明的位姿信息直接求解单应性矩阵的拼接方法流程图；
图6是本发明的基于缝合带羽化的图像融合方法示意图。
具体实施方式
以下结合附图和具体实施案例对本发明作进一步说明。
实施例1
一种位姿信息辅助的航空遥感图像拼接方法，具体流程图如图1所示，步骤内容如下：
S1：遥感图像数据与导航位姿数据的作业获取与预处理；
S2：根据位姿信息计算航带相邻图像的位移量，求解图像重叠区域；
S3：提取重叠区域图像特征，定位特征点并生成特征描述子后进行特征匹配，求解图像变换单应性矩阵，实现图像配准；
S4：分析基于重叠区域特征的图像配准结果的配准质量，若配准失效，则将导航系统位姿信息转换为相机的旋转矩阵和平移向量，直接求解图像间单应性矩阵，实现图像配准；
S5：基于单应性矩阵进行图像变换，使用基于缝合带羽化的图像融合方法进行融合处理；
S6：根据像点和对面点共线关系，将导航位姿信息变换为地面地理坐标，获得包含图像地理信息的遥感图像全景图。
步骤S1中的数据的获取与预处理，包括以下步骤：
(1)对相机传感器和导航系统进行标定实验，包括标定时间差，进行图像、位姿数据频率同步，并标定获得相机的内参数，用于后续拼接；
(2)航空遥感载体搭载相机和导航系统，通过作业获取图像、位姿数据，选择图像数据类型进行存储，并根据导航算法解算获得作业中的经度、纬度、高度和三个姿态角；
(3)使用自适应直方图均衡方法对图像进行增强处理，校正曝光，增加细节。
(4)对由传感器质量、大气、载体姿态、地形起伏、地球自转等因素引起的图像畸变进行校正处理，恢复图像的地貌特征，为遥感图像拼接做准备。
步骤S2中的相邻图像重叠区域求解示意图如图2所示，包括以下步骤：
(1)根据航带内图像I1与图像I2拍摄时刻的经纬度信息，求解图像I1与图像I2中像主点的地理坐标(X1,Y1)、(X2,Y2)，表明了图像在东向、北向的移动量
(2)已知航向角ψ，根据旋转变换关系，可以求得图像I1与图像I2在图像右上的x、y方向上的移动量。其中图像坐标系
(3)根据已知的飞行航高H、相机焦距f，像元尺寸dxy，可求得图像中一个像素对应地面的尺寸
(4)计算得相邻图像的像元移动量，最后即可根据图像尺寸计算得到相邻图像重叠区域
(5)同一航带内运行速度变化较小，图像x、y方向上的位移偏差不大，因此根据(1)～(4)步骤计算前四组图像得重叠区域，求取平均值，作为整条航线上图片的重叠区域大小。其中，通过航向角的改变幅度，判断航线的变化。
步骤S3中的基于重叠区域特征的图像配准，包括以下步骤：
(1)基于SURF算法提取特征点。包括构建Hessian矩阵，通过判别式检测特征点，舍去对比度低的点和边缘响应点；生成图像尺度空间；定位特征点，构建64维的特征描述算子。
(2)根据自适应阈值求取特征匹配对假设两张图像分别有n1、n2个特征点，图像特征点的距离用欧式距离度量特征向量相似性来表示，求取最近邻距离dmx与次近邻距离dmy。
其中，所有特征点m＝1,…,n1，(xm1,xm2,…,xmn)T为特征点的特征描述向量，n为特征描述向量的维数。
在计算最近邻特征距离和次近邻距离后，定义距离比值
计算图像I1到图像I2的匹配特征，再计算图像I2到图像I1的匹配特征，得到两个匹配集合A、B。根据集合A、B中的距离比值，求解公共子集C。此时公共子集C的最小距离比值即为自适应阈值。此时根据特征点的特征距离比值和阈值可以判断特征点对匹配是否成功。
(3)使用RANSAC算法进行误配筛选，最终求解图像配准的单应性矩阵。随机取出m对(至少4对)不共线的匹配特征点，计算变换矩阵H，记为模型M；计算初始匹配集中的内点，根据全部内点计算模型M的统计误差error，迭代更新模型M和error，最小error的模型为最优模型，此时的H为图像配准单应性矩阵。基于重叠区域特征匹配的图像拼接方法具体流程图如图3所示。
步骤S4中的图像配准结果的配准质量分析。首先，设定特征匹配对数量有效阈值，若特征匹配对数量过少，则认定配准失败。其次，设定配准误差阈值，假如特征匹配对数量足够且求解出单应性矩阵，则求解图像配准的配准误差，若配准误差过高，则同样认定为图像配准失败。若配准失败，根据位姿直接求解单应性矩阵实现图像拼接，流程图如图5所示。
步骤S4中的基于位姿信息直接求解单应性矩阵的方法，包括以下步骤：
(1)由姿态数据计算图像I1与图像I2相机坐标系的旋转矩阵R。
遥感平台搭载的导航系统直接输出的数据，是载体的导航信息，包括在WGS-84大地坐标系下的纬度、经度、高度，载体坐标系相对导航坐标系的三个旋转角度：航向角、横滚角、俯仰角。需要将导航位姿数据转换为相机拍摄两张照片的位姿关系，包括旋转运动R矩阵和平移运动t向量。假设相机在两个时刻拍摄两张图像I1、I2，相机光心分别位于点OC1和点OC2，物点P在连续两幅图像中成像分别为像点p1、p2，在各自相机坐标系中的坐标为(u1,v1,1)T、(u2,v2,1)T，如图4所示。
旋转矩阵R的转换关系为导航坐标系(n)→载体坐标系(b)→相机坐标系(C)。定义导航坐标系为东北天导航坐标系，载体坐标系为右上载体坐标系，x轴指向地理东向、y轴指向地理北向，z轴指向天。相机坐标系为x轴、y轴和图像坐标系保持一致，z轴指向物点。图像I1相机坐标系和图像I2相机坐标系之间的旋转矩阵R为
其中，定义航向角ψ，北偏东为正；俯仰角θ，抬头为正；横滚角γ，右倾为正，且
(2)由姿态数据计算图像I1与图像I2相机坐标系的平移向量t。平移向量t的转换关系为：WGS-84大地坐标系→WGS-84地心空间直角坐标系(e)→导航坐标系(n)→载体坐标系(b)→世界坐标系(C1)→相机坐标系(C2)。WGS-84坐标系原点为地球质心，X轴指向零度子午面和赤道的交点，Z轴指向协议地球极，Y轴由右手定则得到。根据WGS-84大地坐标和空间直角坐标转换关系计算拍摄时刻空间直角坐标X、Y、Z。
其中，卯酉圈半径a为参考椭球体长半轴，b为参考椭球体短半轴；/>
图像I1相机坐标系和图像I2相机坐标系之间的平移向量t为
其中，R为公式中旋转矩阵，
(3)根据位姿信息求解两幅图像的单应性矩阵。以图像I1的相机坐标系作为世界坐标系，即作为图像配准的基准坐标系，图像I1和图像I2之间的单应性矩阵为
其中，由于航拍图像的拍摄高度较高，在地面起伏不大的情况下，h由光心O1的高度表示。
(4)求解航带内多图像间单应性矩阵。将图像I2变换到图像I1平面，得到I1和I2组成的大图Ip1，I1坐标系与拼接结果大图Ip1坐标系仅为平移关系。假设输入图片尺寸为宽w(pixel)，高h(pixel)，则图像I2的四个角点在I2像素坐标系中坐标为A2(0,0)，B2(w,0)，C2(w,h)，D2(0,h)，根据单应性变换到图I1像素坐标系中，假设I1、I2的四个角点在I1像素坐标系中坐标分别为A2,1(uA2,vA2)，B2,1(uB2,vB2)，C2,1(uC2,vC2)，D2,1(uD2,vD2)；A1(uA1,vA1)，B1(uB1,vB1)，C1(uC1,vC1)D1(uD1,vD1)。拼接结果大图Ip1的尺寸宽W、高H为
I1像素坐标系与Ip1像素坐标系的变换关系为
其中，
图像I2与Ip1的变换关系为
图像I3与Ip1的变换关系为
同理图像间单应性矩阵皆可求得，可以求解多图像的单应性矩阵。
步骤S5中的基于缝合带羽化的图像融合方法，包括以下步骤：
(1)根据单应性矩阵计算变换后图像四个端点的坐标，求解缝合线坐标。设置缓冲宽度，得到缝合带区域，如图6所示。
(2)在缝合带区域内运用羽化方法分别对上下边界、左右边界图像像素进行处理。设置过渡步长：依次遍历融合：第n行/列的像素值为n·sI1(x,y)+(1-n·s)I2(x,y)。
其中，D1-2表示缝合线区域即融合区域，其他区域皆保留图像I12。α、β表示融合时的权重，β＝1-α。其中Whalf表示缝合线的半宽度，d∈(0,2Whalf)表示融合时距缝合线区域的初始边界的距离。上下边融合时，在缝合线区域由上向下过渡时d从0变化到2Whalf，缝合带中接缝线上侧的区域为I1、I2的羽化融合，缝合线下侧区域为I1。
步骤S6中的地理坐标求解方法，具体方法为：
以图I1相机坐标系为世界坐标系，变换后的全景图在图I1坐标系平面中，图像中像点像素坐标和地理坐标的关系为
其中，图I1时刻相机光心在地理坐标系的坐标为(XO,YO,0)，求的像点对应坐标(X,Y,Z)。
图I1像素坐标(u0,v0,1)与全景图像素坐标的关系为
经纬度与地理坐标的转换关系，由高斯-克吕格投影得到。根据以上公式，在获得全景图之后，根据像点在全景图的像素坐标，即可求得对应地面点在中国WGS-84公里网的地理坐标，实现对全景图的地理粗定位，实现地面物体的大致定位，对遥感图像的后续分析与检测等提供帮助。
实施例2
本发明提供一种位姿信息辅助的航空遥感图像拼接方法，包括以下步骤：
首先，设计作业实验的拍摄时间、航线等；提前做好相机传感器和导航系统的标定实验，，包括标定相机和导航系统时间差，标定获得相机的内参数；进行作业实验，获取遥感图像数据和导航解算位姿数据；对获取图像进行增强、畸变校正处理。
其次，根据位姿信息求解航带图像位移量，进而求解相邻图像重叠区域；在重叠区域内通过SURF进行特征提取，并使用自适应阈值和RANSAC进行特征匹配，求解图像变换的单应性矩阵；
再次，通过特征匹配数量和配准误差，判断图像配准质量；若配准失效，将导航系统位姿信息转换为相机的旋转矩阵和平移向量，直接求解图像间单应性矩阵，并推广求解多图像单应性矩阵；
最后，根据求解的单应性矩阵进行图像变换，使用基于缝合带羽化的图像融合方法进行融合处理；根据像点和对面点共线关系，将导航位姿信息变换为地面地理坐标，获得包含图像地理信息的遥感图像全景图。
本发明利用航空遥感搭载导航系统的位姿信息获得机载相机位姿信息，辅助遥感图像配准、拼接，最后通过融合获得具有地理信息的全景图。与现有技术相比，本发明实现了复杂场景下的遥感图像快速拼接，为航空遥感应用提供了一种高性能、高适应性的图像拼接方法。
Claims (7)
1.一种位姿信息辅助的航空遥感图像拼接方法，其特征在于：包括以下步骤：
S1：遥感图像数据与导航位姿数据的作业获取与预处理；
S2：根据位姿信息计算航带相邻图像的位移量，求解图像重叠区域；
S3：基于SURF算法提取重叠区域图像特征，定位特征点并生成特征描述子后进行特征匹配，求解图像变换单应性矩阵，实现图像配准：
基于SURF算法提取特征点，
根据自适应阈值求取特征匹配，
使用RANSAC算法进行误配筛选，最终求解图像配准单应性矩阵；
S4：分析基于重叠区域特征的图像配准结果的配准质量，若配准失效，则将导航系统位姿信息转换为相机的旋转矩阵和平移向量，直接求解图像间单应性矩阵，实现图像配准；
S5：基于单应性矩阵进行图像变换，使用基于缝合带羽化的图像融合方法进行融合处理；
S6：根据像点和对面点共线关系，将导航位姿信息变换为地面地理坐标，获得包含图像地理信息的遥感图像全景图。
2.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S2中的相邻图像重叠区域简化求解方法，具体实现步骤如下：
(1)根据航带内图像I1和图像I2拍摄时刻的经纬度信息，求解图像I1和图像I2中像主点的地理坐标(X1,Y1)、(X2,Y2)，表明了图像在东向、北向的移动量
(2)已知航向角ψ，根据旋转变换关系，可以求得图像I1和图像I2在图像右上的x、y方向上的移动量，其中图像坐标系
(3)根据已知的飞行航高H、相机焦距f，像元尺寸dxy，可求得图像中一个像素对应地面的尺寸
(4)计算得相邻图像的像元移动量，最后即可根据图像尺寸计算得到相邻图像重叠区域
(5)同一航带内运行速度变化较小，图像x、y方向上的位移偏差不大，因此根据(1)～(4)步骤计算前四组图像得重叠区域，求取平均值，作为整条航线上图片的重叠区域大小，其中，通过航向角的改变幅度，判断航线的变化。
3.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S3中的基于重叠区域特征匹配的图像配准方法，具体实现步骤如下：
(1)基于SURF算法提取特征点，包括构建Hessian矩阵，通过判别式检测特征点，舍去对比度低的点和边缘响应点；生成图像尺度空间；定位特征点，构建64维的特征描述算子；
(2)根据自适应阈值求取特征匹配对假设两张图像分别有n1、n2个特征点，图像特征点的距离用欧式距离度量特征向量相似性来表示，求取最近邻距离dmx与次近邻距离dmy
其中，所有特征点m＝1,…,n1，(xm1,xm2,…,xmn)T为特征点的特征描述向量，n为特征描述向量的维数；
在计算最近邻特征距离和次近邻距离后，定义距离比值
计算图像I1到图像I2的匹配特征，再计算图像I2到图像I1的匹配特征，得到两个匹配集合A、B，根据集合A、B中的距离比值，求解公共子集C，此时公共子集C的最小距离比值即为自适应阈值，此时根据特征点的特征距离比值和阈值可以判断特征点对匹配是否成功；
(3)使用RANSAC算法进行误配筛选，最终求解图像配准的单应性矩阵，随机取出m对，至少4对，不共线的匹配特征点，计算变换矩阵H，记为模型M；计算初始匹配集中的内点，根据全部内点计算模型M的统计误差error，迭代更新模型M和error，最小error的模型为最优模型，此时的H为图像配准单应性矩阵。
4.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S4中的基于位姿信息直接求解单应性矩阵的方法，具体实现步骤如下：
(1)由姿态数据计算图像I1与图像I2相机坐标系的旋转矩阵R：
遥感平台搭载的导航系统直接输出的数据，是载体的导航信息，包括在WGS-84大地坐标系下的纬度、经度、高度，载体坐标系相对导航坐标系的三个旋转角度：航向角、横滚角、俯仰角，需要将导航位姿数据转换为相机拍摄两张照片的位姿关系，包括旋转运动R矩阵和平移运动t向量，假设相机在两个时刻拍摄两张图像I1、I2，相机光心分别位于点OC1和点OC2；
旋转矩阵R的转换关系为导航坐标系(n)→载体坐标系(b)→相机坐标系(C)，定义导航坐标系为东北天导航坐标系，载体坐标系为右上载体坐标系，x轴指向地理东向、y轴指向地理北向，z轴指向天，相机坐标系为x轴、y轴和图像坐标系保持一致，z轴指向物点，图像I1相机坐标系和图像I2相机坐标系之间的旋转矩阵R为
其中，定义航向角ψ，北偏东为正；俯仰角θ，抬头为正；横滚角γ，右倾为正，且
(2)由姿态数据计算图像I1与图像I2相机坐标系的平移向量t，平移向量t的转换关系为：WGS-84大地坐标系→WGS-84地心空间直角坐标系(e)→导航坐标系(n)→载体坐标系(b)→世界坐标系(C1)→相机坐标系(C2)，WGS-84坐标系原点为地球质心，X轴指向零度子午面和赤道的交点，Z轴指向协议地球极，Y轴由右手定则得到，根据WGS-84大地坐标和空间直角坐标转换关系计算拍摄时刻空间直角坐标X、Y、Z
其中，卯酉圈半径a为参考椭球体长半轴，b为参考椭球体短半轴；/>
图像I1相机坐标系和图像I2相机坐标系之间的平移向量t为
其中，R为公式中旋转矩阵，
(3)根据位姿信息求解两幅图像的单应性矩阵，以图像I1的相机坐标系作为世界坐标系，即作为图像配准的基准坐标系；
(4)求解航带内多图像间单应性矩阵，将图像I2变换到图像I1平面，得到I1和I2组成的大图Ip1，I1坐标系与拼接结果大图Ip1坐标系仅为平移关系，假设输入图片尺寸为宽w(pixel)，高h(pixel)，则图像I2的四个角点在I2像素坐标系中坐标为A2(0,0)，B2(w,0)，C2(w,h)，D2(0,h)，根据单应性变换到图I1像素坐标系中，假设I1、I2的四个角点在I1像素坐标系中坐标分别为A2,1(uA2,vA2)，B2,1(uB2,vB2)，C2,1(uC2,vC2)，D2,1(uD2,vD2)；A1(uA1,vA1)，B1(uB1,vB1)，C1(uC1,vC1)D1(uD1,vD1)，拼接结果大图Ip1的尺寸宽W、高H为
I1像素坐标系与Ip1像素坐标系的变换关系为
其中，
图像I2与Ip1的变换关系为
图像I3与Ip1的变换关系为
同理图像间单应性矩阵皆可求得，可以求解多图像的单应性矩阵。
5.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S5中的基于缝合带羽化的图像融合方法，具体实现步骤如下：
(1)根据单应性矩阵计算变换后图像四个端点的坐标，求解缝合线坐标，设置缓冲宽度，得到缝合带区域；
(2)在缝合带区域内运用羽化方法分别对上下边界、左右边界图像像素进行处理，设置过渡步长：s，依次遍历融合：第n行/列的像素值为n·sI1(x,y)+(1-n·s)I2(x,y)
其中，D1-2表示缝合线区域即融合区域，其他区域皆保留图像I12，α、β表示融合时的权重，β＝1-α，其中Whalf表示缝合线的半宽度，d∈(0,2Whalf)表示融合时距缝合线区域的初始边界的距离，上下边融合时，在缝合线区域由上向下过渡时d从0变化到2Whalf，缝合带中接缝线上侧的区域为I1、I2的羽化融合，缝合线下侧区域为I1。
6.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S3和步骤S4中所述的单应性矩阵的解算是两种独立方案，当S3图像配准失效后进行S4的单应性矩阵求解与图像配准，能适用于复杂场景的遥感图像拼接。
7.根据权利要求1所述的一种位姿信息辅助的航空遥感图像拼接方法，其特征在于步骤S4中所述的单应性矩阵的解算不需要图像特征、模板信息，通过位姿信息直接求解，可以拼接地貌单一、特征稀少的遥感图像。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110553778.5A CN113222820B (zh) | 2021-05-20 | 2021-05-20 | 一种位姿信息辅助的航空遥感图像拼接方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110553778.5A CN113222820B (zh) | 2021-05-20 | 2021-05-20 | 一种位姿信息辅助的航空遥感图像拼接方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113222820A CN113222820A (zh) | 2021-08-06 |
| CN113222820B true CN113222820B (zh) | 2024-05-07 |
Family
ID=77093352
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110553778.5A Active CN113222820B (zh) | 2021-05-20 | 2021-05-20 | 一种位姿信息辅助的航空遥感图像拼接方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113222820B (zh) |
Families Citing this family (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113781444B (zh) * | 2021-09-13 | 2024-01-16 | 北京理工大学重庆创新中心 | 基于多层感知机校正的快速拼接航拍图像的方法和系统 |
| CN113705734B (zh) * | 2021-09-30 | 2022-12-09 | 中国电子科技集团公司第五十四研究所 | 一种基于多传感器和地心的遥感图像特征点高程获取方法 |
| CN114519671B (zh) * | 2022-02-16 | 2022-11-01 | 天津中科无人机应用研究院 | 无人机遥感影像动态快速拼接方法 |
| CN116363185B (zh) * | 2023-06-01 | 2023-08-01 | 成都纵横自动化技术股份有限公司 | 地理配准方法、装置、电子设备和可读存储介质 |
| CN117036666B (zh) * | 2023-06-14 | 2024-05-07 | 北京自动化控制设备研究所 | 基于帧间图像拼接的无人机低空定位方法 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2014193670A2 (en) * | 2013-05-29 | 2014-12-04 | Capso Vision, Inc. | Reconstruction of images from an in vivo multi-camera capsule |
| CN109325913A (zh) * | 2018-09-05 | 2019-02-12 | 北京悦图遥感科技发展有限公司 | 无人机图像拼接方法及装置 |
| WO2020037615A1 (zh) * | 2018-08-23 | 2020-02-27 | 深圳市大疆创新科技有限公司 | 云台系统及其图像处理方法、无人机 |
| CN111798373A (zh) * | 2020-06-11 | 2020-10-20 | 西安视野慧图智能科技有限公司 | 一种基于局部平面假设及六自由度位姿优化的快速无人机图像拼接方法 |
-
2021
- 2021-05-20 CN CN202110553778.5A patent/CN113222820B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2014193670A2 (en) * | 2013-05-29 | 2014-12-04 | Capso Vision, Inc. | Reconstruction of images from an in vivo multi-camera capsule |
| WO2020037615A1 (zh) * | 2018-08-23 | 2020-02-27 | 深圳市大疆创新科技有限公司 | 云台系统及其图像处理方法、无人机 |
| CN109325913A (zh) * | 2018-09-05 | 2019-02-12 | 北京悦图遥感科技发展有限公司 | 无人机图像拼接方法及装置 |
| CN111798373A (zh) * | 2020-06-11 | 2020-10-20 | 西安视野慧图智能科技有限公司 | 一种基于局部平面假设及六自由度位姿优化的快速无人机图像拼接方法 |
Non-Patent Citations (3)
| Title |
|---|
| Extended Ellipsoidal Outer-Bounding Set-Membership Estimation for Nonlinear Discrete-Time Systems with Unknown-but-Bounded Disturbances;Liu Y., Zhao Y., Wu F.;Discrete Dynamics in Nature and Society;第11卷(第1期);全文 * |
| 一种基于无人机位姿信息的航拍图像拼接方法;程争刚;张利;;测绘学报(06);全文 * |
| 视觉辅助惯性定位定姿技术研究;史俊;杨功流;陈雅洁;万振塬;;航空计算技术(01);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113222820A (zh) | 2021-08-06 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN113222820B (zh) | 一种位姿信息辅助的航空遥感图像拼接方法 | |
| CN108648240B (zh) | 基于点云特征地图配准的无重叠视场相机姿态标定方法 | |
| US7233691B2 (en) | Any aspect passive volumetric image processing method | |
| CN110675450B (zh) | 基于slam技术的正射影像实时生成方法及系统 | |
| CN106373088B (zh) | 大倾斜低重叠率航空图像的快速拼接方法 | |
| CN108765298A (zh) | 基于三维重建的无人机图像拼接方法和系统 | |
| US9098229B2 (en) | Single image pose estimation of image capture devices | |
| CN107808362A (zh) | 一种基于无人机pos信息与图像surf特征结合的图像拼接方法 | |
| CN111882612A (zh) | 一种基于三维激光检测车道线的车辆多尺度定位方法 | |
| CN109446973B (zh) | 一种基于深度神经网络图像识别的车辆定位方法 | |
| CN109596121B (zh) | 一种机动站自动目标检测与空间定位方法 | |
| US8547375B2 (en) | Methods for transferring points of interest between images with non-parallel viewing directions | |
| CN111507901A (zh) | 基于航带gps及尺度不变约束的航拍图像拼接定位方法 | |
| CN110555813B (zh) | 一种无人机遥感影像的快速几何校正方法及系统 | |
| CN111161154A (zh) | 一种无人机视频实时快速正射拼接系统及其拼接方法 | |
| CN113313659A (zh) | 一种多机协同约束下高精度图像拼接方法 | |
| CN114372992A (zh) | 一种基于动平台的边缘角点检测四目视觉算法 | |
| CN115451964A (zh) | 基于多模态混合特征的船舶场景同时建图与定位方法 | |
| CN112750075A (zh) | 一种低空遥感影像拼接方法以及装置 | |
| CN117036666B (zh) | 基于帧间图像拼接的无人机低空定位方法 | |
| Zhao et al. | Digital Elevation Model‐Assisted Aerial Triangulation Method On An Unmanned Aerial Vehicle Sweeping Camera System | |
| Li et al. | Learning dense consistent features for aerial-to-ground structure-from-motion | |
| CN114199250A (zh) | 一种基于卷积神经网络的景象匹配导航方法及装置 | |
| CN108986025B (zh) | 基于不完备姿轨信息的高精度异时图像拼接及校正方法 | |
| CN116907511B (zh) | 一种将管道坐标转换为图像坐标的方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |