CN110310310B - 一种针对航空影像配准的改进方法 - Google Patents

Abstract

本发明提出一种针对航空影像配准的改进方法，提高航空影像配准的时间效率，同时提高配准结果的精度，针对航空影像变化尺度较大，序列重合度低等特点，采用该改进方法能够有效地完成航空影像的配准任务。基于SURF特征提取方法生成航空影像的匹配点对；对提取的特征点对进行预处理并生成样本优先集；建立模型适应度函数代替RANSAC算法中人为设置阈值的刚性判别；通过改进RANSAC算法计算航空影像间的变换矩阵，最终获得较优的配准结果。

Description

一种针对航空影像配准的改进方法

技术领域

本发明属于计算机视觉中面向图像处理的技术领域，尤其涉及一种针对航空影像配准的改进方法。

背景技术

图像配准(Image Registration)过程包含设定一幅图片作为参考图片，对于其他的图片采用空间转换来匹配此参考图片，此过程亦称作图像匹配或者图像相关(ImageMatching or Image Correlation)。图像的空间转换就是寻找一种坐标变换模型，建立从一个图像坐标到另一个图像坐标之间的映射关系。常用变换模型包括四种变换：刚体变换、仿射变换、投影变换和非线性变换。由于空间转换为从一幅图像的点转换为另一幅图像中的点，因此，正确的空间转换参数是图像正确配准的关键。特征点对的误匹配、照相机的捕捉视角、镜头及传感器的畸变、摄像机的差别等均会导致图像的误匹配。介于后面影响因素均为硬件条件，不宜调整，故相关研究均致力于降低特征误匹配对图像配准产生的影响。

由于误匹配特征的存在,模型参数估计将产生较大的偏差,甚至出现估计失败。通常,模型参数估计可以分为线性法、迭代法和鲁棒法。线性法速度较快,但对于存在错误匹配点以及由噪声引起坏点的情况精确性较低。迭代法精度比线性算法高,但计算时间长,且误匹配点剔除效果不是很好。鲁棒法是通过制定错误信息排除策略,并利用正确信息进行模型参数估计的方法。目前,研究人员已经提出了众多鲁棒性参数估计算法,其中应用最为广泛的有极大似然估计(M-estimation)、最小中值(LMedS)和随机抽样一致性(RANSAC)3种算法。RANSAC算法与前两种相比具有算法结构简单、鲁棒性强、能够有效处理较多外点存在的情况等优势，故广泛应用于图像配准过程中的模型参数估计问题中。

RANSAC算法在外点比率较大时算法效率低下，且通过人为设置阈值刚性判断内外点，过大或过小均会对估计结果产生不利影响。本发明针对航空影像的变化尺度大，受光照影响明显且对配准过程实时性要求高等特点，对RANSAC算法进行了一定改进，在快速准确完成影像间变换模型参数估计的前提下实现航空影像配准。

发明内容

发明目的：针对以上现有技术存在的问题，本发明提出一种航空影像间的配准方法，该方法公开了一种改进方法来计算变换矩阵的优化方案，并以该优化方案为基础实现了航空影像间的精确配准，解决了航空影像的特殊性所面临的技术问题。

技术方案：为实现本发明的目的，本发明所采用的技术方案是：本发明公开了一种航空影像配准的改进方法，将改进的图像配准方法应用于变化尺度大，受光照影响明显且对配准过程实时性要求高等特点的航空影像配准。本发明的整体操作流程见图1所示。

(1)对航空影像进行SURF特征提取，包含四个步骤：建立尺度空间，精确定位及选择特征点，决定特征点的主方向，特征点描述。操作流程见图2所示。详细操作步骤如下：

(1.1)选择二阶标准高斯函数作为滤波器，高斯函数的标准差为尺度因子，将不同尺度的高斯核与图像进行卷积以计算二阶偏导数，构建不同尺度σ下的Hessian矩阵

其中，L_xx、L_yy、L_xy分为高斯滤波二阶导

同图像I＝(x,y)的卷积结果。

用盒子滤波近似代替二阶高斯滤波，并借助积分图像加速卷积，简化Hessian矩阵的求解过程。

通过改变盒子滤波的大小，将不同大小的盒子滤波同原始图像做卷积，即可生成不同尺度的金字塔图像。生成的金字塔图像分为4个组，每一组包含4张尺度不同的图像。对于L×L大小的盒子滤波，其对应的尺度为

尺度空间金字塔的最底层由9×9的盒子滤波输出得到，则对应尺度为σ₀＝1.2。为保证盒子滤波的结构不变，后续盒子滤波器的大小至少以6个像素值步长增加。在生成金字塔图像的过程中，图像大小保持不变，不同组之间的图像通过改变盒子滤波尺寸而获得，每一组的第一层盒子滤波尺寸是上一组的第二层盒子滤波的尺寸。第1组中，相邻盒子滤波尺寸L以6个像素的步长逐步增加；第2组中，相邻盒子滤波尺寸L以12个像素的步长逐步增加；第3组中，相邻盒子滤波尺寸L以24个像素的步长逐步增加；第4组中，相邻盒子滤波尺寸L以48个像素的步长逐步增加，最终，生成不同尺度的金字塔图像。

(1.2)选择上述每组图像的第2层和第3层卷积图像中的像素点的Hessian响应，即Hessian矩阵行列式的近似值，分别与上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应进行比较，若该点的Hessian响应大于或小于上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应值，该点为最大或最小极值点，则保留为初步特征点。图2中标记为‘x’的像素点的Hessian响应若大于3×3×3周围像素的Hessian响应，则可确定该点为该区域的初步特征点。然后，在尺度空间采用三维线性插值法对初步特征点进行插值运算得到亚像素级别的极值点，该极值点为稳定的特征点，同时可获得该特征点的位置及所在的尺度。

(1.3)以稳定的特征点为中心，计算半径为6s的邻域内的点在水平x和垂直y方向的Haar小波响应，即Haar小波特征值，s为特征点所属尺度空间的尺度值，并对所求的响应值赋高斯权重系数σ₁＝2s，并加权；圆心角为60°的扇形以特征点为中心环绕一周，旋转步长为10°，计算该扇形处于每个角度时所包含图像点的Haar小波响应之和，由于每一点都有x、y两个方向的响应，因此扇形区域中所有点的响应之和构成一个新的矢量，遍历整个圆形区域，选择最长矢量的方向为该特征点的主方向。图2中标记出了中心特征点的主方向。

(1.4)以特征点为中心，将坐标轴旋转到特征点的主方向，即特征点的主方向为y轴。按照主方向选取边长为20s的正方形区域，将区域分为4×4个子区域。以特征点为中心，采样步长取该特征点所属尺度空间的尺度值大小是，在每个子区域内统计5s×5s范围内的Haar小波响应，相对于主方向的水平、垂直方向的Haar小波响应分别记作d_x、d_y，同样赋予高斯权重系数σ₂＝3.3s，并加权。将每个子区域中的Haar小波响应及响应的绝对值相加得到∑d_x、∑d_y、∑|d_x|、∑|d_y|，即在每个子区域可获得四维矢量V_sub＝(∑d_x,∑|d_x|,∑d_y,∑|d_y|)，最后，形成4×4×4＝64维的描述向量。

(2)由于航空序列影像间的重叠度较低，采用最近邻算法容易产生较多误匹配点，故在提取SURF特征点之后，对特征点匹配进行预处理以消除误匹配点。预处理之后，设定阈值生成样本优先集作为RANSAC算法的初始采样集。

(2.1)预处理过程具体为通过建立评测标准η(x)＝1-d₁(x)/d₂(x)，来剔除明显误匹配点完成预处理步骤，其中：x为特征描述向量，η(x)为特征点匹配调整系数，d₁(x)为左影像中特征点通过遍历找出的右影像中与之最近的特征点之间的欧式距离，d₂(x)为左影像中该特征点通过遍历找出的右影像中与之次近的特征点之间的欧式距离，其中，左影像为参考影像，右影像为待配准影像。若η(x)＞0.3，则将左右影像中的特征点视为一对特征点对。

(2.2)预处理后，得到特征点对样本集X＝{x_i|i＝1,2,…,N}，其中，N为总特征点对数，设第i对特征点对之间的距离为d(i)，i＝1,2,…,N。根据

计算标准差估计值σ_d，将满足d(i)≤σ_d条件的所有特征点对加入样本优先集X'＝{x'_j|j＝1,2,…,N'}，其中，N'为符合d(i)≤σ_d条件的特征点对数。

(3)图3为RANSAC算法的原理图。以图3为例，随机选择图中两个点生成一条直线，位于距离该直线一定值范围内的点称为支撑点，并统计支撑点数目。重复选择两个点生成一系列直线，统计不同直线所产生的支撑点数目，并选择包含支撑点最多的直线为最终计算结果。属于该直线的支撑点称为内点，剩余点则为外点。图3显示出RANSAC算法优于最小二乘法。类似于两点生成一条直线的RANSAC理论，将RANSAC算法应用到四对特征点生成变换矩阵的情况。

在本发明中，通过建立模型适应度函数对原始RANSAC算法进行改进，建立模型适应度函数

其中/>

代替根据余差R_i大小所产生的内点比率来对变换矩阵H进行选择，其中，x_i＝(x_i,y_i)与x'_i＝(x′_i,y′_i)为一对特征点对，/>

为样本集X中样本对由优先样本集X'采样所求得的变换矩阵H的余差，d()为两点之间的欧式影像距离，v＝R_max-R_min，R_max和R_min为余差最大值和最小值，/>

为标准差估计，取值为/>

其中，m为影像间变换矩阵参数估计所需最小样本数，i＝1,2,…,N。

将改进的RANSAC算法用于计算影像间的变换矩阵，该算法可以有效地剔除外点，即错误匹配点，并利用内点，即正确匹配点，实现变换矩阵参数的最优估计。改进RANSAC算法具体步骤可描述为：

(3.1)输入匹配点样本集X＝{x_i|i＝1,2,…,N}，样本优先集X'＝{x'_j|j＝1,2,…,N'}，样本中内点比率ε＝0.8，用作计算采样次数，置信率η₀，内外点判定阈值τ，影像间变换矩阵参数估计所需最小样本数m，迭代次数l＝0；

(3.2)计算采样次数k＝log(1-η₀)/log(1-ε^m)，并从样本优先集X'中随机抽取四个点对(x_a,x'_a)建立方程组

经过k次采样，即k次随机抽取四个点对，求得影像间得变换矩阵集Ω_H＝(H₁,H₂,…H_k)，变换矩阵H的求解步骤入下：

变换矩阵H为大小为3×3的矩阵，h₃₃＝1。剩余变换矩阵H的8个未知数需要8个方程来求解，每对特征点对(x_a,x'_a)可提供两个方程，因此，四对特征点对即可正确求解出变换矩阵H。

(3.3)对变换矩阵集Ω_H＝(H₁,H₂,…H_k)计算不同变换矩阵H_h，h＝1,2,…,k，所对应的余差R_ih，并利用函数F(H)来检验当前变换矩阵H的模型适应度，遍历k个变换矩阵H之后，选择函数F(H)最大值所对应的变换矩阵H，记为H_max；

(3.4)H_max对应的余差R_i，i＝1,2,…,N，与设定的阈值τ作比较，若H_max对应的余差项R_i小于τ的特征点对为内点，反之为外点，更新内点比率

n_inlier为采样集X中的内点数目；

(3.5)迭代次数l＝l+1，与采样次数k进行比较，若l＞k，停止循环，并输出内点比率最高的变换矩阵作为最终结果，否则重复步骤(3.2)-(3.5)。

有益效果：与现有技术相比，本发明的技术方案具有以下有益技术效果：

1、采用SURF特征算子提取算法用于航空影像配准，能够快速生成相对可靠的特征用于后续处理步骤；

2、对匹配点对进行预处理之后，设定阈值生成样本优先集作为RANSAC算法的初始采样集，提高运行效率；

3、建立模型适用度函数用于RANSAC算法中变换矩阵的可靠性分析，提高算法的鲁棒性。

附图说明

图1是本发明的整体流程图；

图2是提取SURF特征算子的过程；

图3是RANSAC原理的示意图；

图4是不同计算方法所获得的变换矩阵；

图5是基于SIFT和SURF特征算子进行点对匹配的对比图；

图6是利用本发明对航空影像配准的结果示意图。

具体实施方式

下面结合附图对发明的技术方案进行详细说明。下面通过参考附图描述的实施方式是示例性的，仅用于解释本发明，而不能解释为对本发明的限制。

本领域的技术人员可以理解，除非另外定义，这里使用的所有术语，包括技术术语和科学术语，具有与本发明所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样定义，不会用理想化或过于正式的含义来解释。

本发明公开了一种航空影像配准的改进方法，将改进的图像配准方法应用于变化尺度大，受光照影响明显且对配准过程实时性要求高等特点的航空影像配准。本发明的整体操作流程见图1所示。

(1)对航空影像进行SURF特征提取，包含四个步骤：建立尺度空间，精确定位及选择特征点，决定特征点的主方向，特征点描述。操作流程见图2所示。详细操作步骤如下：

(1.1)选择二阶标准高斯函数作为滤波器，高斯函数的标准差为尺度因子，将不同尺度的高斯核与图像进行卷积以计算二阶偏导数，构建不同尺度σ下的Hessian矩阵

其中，L_xx、L_yy、L_xy分为高斯滤波二阶导

同图像I＝(x,y)的卷积结果。

用盒子滤波近似代替二阶高斯滤波，并借助积分图像加速卷积，简化Hessian矩阵的求解过程。

通过改变盒子滤波的大小，将不同大小的盒子滤波同原始图像做卷积，即可生成不同尺度的金字塔图像。生成的金字塔图像分为4个组，每一组包含4张尺度不同的图像。对于L×L大小的盒子滤波，其对应的尺度为

尺度空间金字塔的最底层由9×9的盒子滤波输出得到，则对应尺度为σ₀＝1.2。为保证盒子滤波的结构不变，后续盒子滤波器的大小至少以6个像素值步长增加。在生成金字塔图像的过程中，图像大小保持不变，不同组之间的图像通过改变盒子滤波尺寸而获得，每一组的第一层盒子滤波尺寸是上一组的第二层盒子滤波的尺寸。第1组中，相邻盒子滤波尺寸L以6个像素的步长逐步增加；第2组中，相邻盒子滤波尺寸L以12个像素的步长逐步增加；第3组中，相邻盒子滤波尺寸L以24个像素的步长逐步增加；第4组中，相邻盒子滤波尺寸L以48个像素的步长逐步增加，最终，生成不同尺度的金字塔图像。

(1.2)选择上述每组图像的第2层和第3层卷积图像中的像素点的Hessian响应，即Hessian矩阵行列式的近似值，分别与上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应进行比较，弱该点的Hessian响应大于或小于上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应值，该点为最大或最小极值点，则保留为初步特征点。图2中标记为‘x’的像素点的Hessian响应若大于3×3×3周围像素的Hessian响应，则可确定该点为该区域的初步特征点。然后，在尺度空间采用三维线性插值法对初步特征点进行插值运算得到亚像素级别的极值点，该极值点为稳定的特征点，同时可获得该特征点的位置及所在的尺度。

(1.3)以稳定的特征点为中心，计算半径为6s的邻域内的点在水平x和垂直y方向的Haar小波响应，即Haar小波特征值，s为特征点所属尺度空间的尺度值，并对所求的响应值赋高斯权重系数σ₁＝2s，并加权；圆心角为60°的扇形以特征点为中心环绕一周，旋转步长为10°，计算该扇形处于每个角度时所包含图像点的Haar小波响应之和，由于每一点都有x、y两个方向的响应，因此扇形区域中所有点的响应之和构成一个新的矢量，遍历整个圆形区域，选择最长矢量的方向为该特征点的主方向。图2中标记出了中心特征点的主方向。

(1.4)以特征点为中心，将坐标轴旋转到特征点的主方向，即特征点的主方向为y轴。按照主方向选取边长为20s的正方形区域，将区域分为4×4个子区域。以特征点为中心，采样步长取该特征点所属尺度空间的尺度值大小是，在每个子区域内统计5s×5s范围内的Haar小波响应，相对于主方向的水平、垂直方向的Haar小波响应分别记作d_x、d_y，同样赋予高斯权重系数σ₂＝3.3s，并加权。将每个子区域中的Haar小波响应及响应的绝对值相加得到∑d_x、∑d_y、∑|d_x|、∑|d_y|，即在每个子区域可获得四维矢量V_sub＝(∑d_x,∑|d_x|,∑d_y,∑|d_y|)，最后，形成4×4×4＝64维的描述向量。

(2)由于航空序列影像间的重叠度较低，采用最近邻算法容易产生较多误匹配点，故在提取SURF特征点之后，对特征点匹配进行预处理以消除误匹配点。预处理之后，设定阈值生成样本优先集作为RANSAC算法的初始采样集。

(2.1)预处理过程具体为通过建立评测标准η(x)＝1-d1(x)/d2(x)，来剔除明显误匹配点完成预处理步骤，其中：x为特征描述向量，η(x)为特征点匹配调整系数，d₁(x)为左影像中特征点通过遍历找出的右影像中与之最近的特征点之间的欧式距离，d₂(x)为左影像中该特征点通过遍历找出的右影像中与之次近的特征点之间的欧式距离，其中，左影像为参考影像，右影像为待配准影像。若η(x)＞0.3，则将左右影像中的特征点视为一对特征点对。

(2.2)预处理后，得到特征点对样本集X＝{x_i|i＝1,2,…,N}，其中，N为总特征点对数，设第i对特征点对之间的距离为d(i)，i＝1,2,…,N。根据

计算标准差估计值σ_d，将满足d(i)≤σ_d条件的所有特征点对加入样本优先集X'＝{x'_j|j＝1,2,…,N'}，其中，N'为符合d(i)≤σ_d条件的特征点对数。

(3)图3为RANSAC算法的原理图。以图3为例，随机选择图中两个点生成一条直线，位于距离该直线一定值范围内的点称为支撑点，并统计支撑点数目。重复选择两个点生成一系列直线，统计不同直线所产生的支撑点数目，并选择包含支撑点最多的直线为最终计算结果。属于该直线的支撑点称为内点，剩余点则为外点。图3显示出RANSAC算法优于最小二乘法。类似于两点生成一条直线的RANSAC理论，将RANSAC算法应用到四对特征点生成变换矩阵的情况。

在本发明中，通过建立模型适应度函数对原始RANSAC算法进行改进，建立模型适应度函数

其中/>

代替根据余差R_i大小所产生的内点比率来对变换矩阵H进行选择，其中，x_i＝(x_i,y_i)与x'_i＝(x′_i,y′_i)为一对特征点对，/>

为样本集X中样本对由优先样本集X'采样所求得的变换矩阵H的余差，d()为两点之间的欧式影像距离，v＝R_max-R_min，R_max和R_min为余差最大值和最小值，/>

为标准差估计，取值为/>

其中，m为影像间变换矩阵参数估计所需最小样本数，i＝1,2,…,N。

将改进的RANSAC算法用于计算影像间的变换矩阵，该算法可以有效地剔除外点，即错误匹配点，并利用内点，即正确匹配点，实现变换矩阵参数的最优估计。改进RANSAC算法具体步骤可描述为：

(3.1)输入匹配点样本集X＝{x_i|i＝1,2,…,N}，样本优先集X'＝{x'_j|j＝1,2,…,N'}，样本中内点比率ε＝0.8，用作计算采样次数，置信率η₀，内外点判定阈值τ，影像间变换矩阵参数估计所需最小样本数m，迭代次数l＝0；

(3.2)计算采样次数k＝log(1-η₀)/log(1-ε^m)，并从样本优先集X'中随机抽取四个点对(x_a,x'_a)建立方程组

经过k次采样，即k次随机抽取四个点对，求得影像间得变换矩阵集Ω_H＝(H₁,H₂,…H_k)，变换矩阵H的求解步骤入下：

变换矩阵H为大小为3×3的矩阵，h₃₃＝1。剩余变换矩阵H的8个未知数需要8个方程来求解，每对特征点对(x_a,x'_a)可提供两个方程，因此，四对特征点对即可正确求解出变换矩阵H。

(3.3)对变换矩阵集Ω_H＝(H₁,H₂,…H_k)计算不同变换矩阵H_h，h＝1,2,…,k，所对应的余差R_ih，并利用函数F(H)来检验当前变换矩阵H的模型适应度，遍历k个变换矩阵H之后，选择函数F(H)最大值所对应的变换矩阵H，记为H_max；

(3.4)H_max对应的余差R_i，i＝1,2,…,N，与设定的阈值τ作比较，若H_max对应的余差项R_i小于τ的特征点对为内点，反之为外点，更新内点比率

n_inlier为采样集X中的内点数目；

(3.5)迭代次数l＝l+1，与采样次数k进行比较，若l＞k，停止循环，并输出内点比率最高的变换矩阵作为最终结果，否则重复步骤(3.2)-(3.5)。

利用本发明中的方法，采用一组无人机航拍的影像进行配准测试。附图4显示了不同变换矩阵的计算方法对影像变化产生的影像，其中(b)为本发明中的改进变换矩阵方法生成的结果图。对比基于SIFT和SURF特征算子的提取过程，结果证明为了满足航空影像配准的实时性及可靠性，SURF更适用于航空影像配准的特征提取。附图5为基于SIFT和SURF特征算子进行点对匹配的对比图。结合本发明的技术方案最终实现的航空影像配准结果如附图6所示。

其中，本发明的技术方案所提取的SURF特征算子可抵抗尺度、光照等变化带来的影响，经过预处理后可获得可靠性较高的匹配点对，并用于改进的RANSAC算法计算航空影像间的变换矩阵，提供了较优的航空影像配准结果。

上述实施例仅为本发明技术方案的一种实现方式，不构成对本发明实施例的限定，本领域的技术人员在本发明公开的航空影像配准方案的基础上，能够将其应用到其它的配准领域中。

Claims (1)

1.一种针对航空影像配准的改进方法，其特征在于，该方法包括如下步骤：

(1)对航空影像进行SURF特征提取，包含四个步骤：建立尺度空间，选择特征点，决定特征点的主方向，并进行特征点描述；

(2)在参考影像上寻找与上述特征点匹配的点形成点对集合，并对点对集合进行预处理以消除误匹配的点生成样本优先集；

(3)构建适应度函数，通过样本优先集计算航空影像和参考影像之间的最佳转换矩阵实现影像上特征点对的映射；

步骤(1)中，对航空影像进行SURF特征提取，包含四个步骤：建立尺度空间，精确定位及选择特征点，决定特征点的主方向，特征点描述的方法如下：

(1.1)选择二阶标准高斯函数作为滤波器，高斯函数的标准差为尺度因子，将不同尺度的高斯核与图像进行卷积以计算二阶偏导数，构建不同尺度σ下的Hessian矩阵

其中，L_xx、L_yy、L_xy分为高斯滤波二阶导/>

同图像I＝(x,y)的卷积结果，通过二阶高斯滤波器设置盒子滤波器，借助积分图像加速卷积，简化Hessian矩阵的求解过程；

(1.2)通过改变盒子滤波的大小，将不同大小的盒子滤波同原始图像做卷积，即可生成不同尺度的金字塔图像，生成的金字塔图像分为4个组，每一组包含4张尺度不同的图像；

(1.3)选择上述每组图像的第2层和第3层卷积图像中的像素点的Hessian响应，分别与上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应进行比较，若该点的Hessian响应大于或小于上一尺度、下一尺度及本尺度周围的26个邻域像素点的Hessian响应值，该点为最大或最小极值点，则保留为初步特征点，在尺度空间采用三维线性插值法对初步特征点进行插值运算得到亚像素级别的极值点，该极值点为稳定的特征点，同时可获得该特征点的位置及所在的尺度；

(1.3)以稳定的特征点为中心，计算半径为6s的邻域内的点在水平x和垂直y方向的Haar小波响应，即Haar小波特征值，s为特征点所属尺度空间的尺度值，并对所求的响应值赋高斯权重系数σ₁＝2s，并加权；圆心角为60°的扇形以特征点为中心环绕一周，旋转步长为10°，计算该扇形处于每个角度时所包含图像点的Haar小波响应之和，由于每一点都有x、y两个方向的响应，因此扇形区域中所有点的响应之和构成一个新的矢量，遍历整个圆形区域，选择最长矢量的方向为该特征点的主方向；

(1.4)以稳定特征点为中心，将坐标轴旋转到特征点的主方向，即特征点的主方向为y轴，按照主方向选取边长为20s的正方形区域，将区域分为4×4个子区域，以特征点为中心，采样步长取该特征点所属尺度空间的尺度值大小s，在每个子区域内统计5s×5s范围内的Haar小波响应，相对于主方向的水平、垂直方向的Haar小波响应分别记作d_x、d_y，同样赋予高斯权重系数σ₂＝3.3s，并加权，将每个子区域中的Haar小波响应及响应的绝对值相加得到∑d_x、∑d_y、∑|d_x|、∑|d_y|，即在每个子区域可获得四维矢量V_sub＝(∑d_x,∑|d_x|,∑d_y,∑|d_y|)；

所述步骤(2)，在参考影像上寻找与上述特征点匹配的点形成点对集合，并对点对集合进行预处理以消除误匹配的点生成样本优先集的方法如下：

(2.1)通过建立评测标准η(x)＝1-d₁(x)/d₂(x)，来剔除明显误匹配点完成预处理步骤，其中：x为特征描述向量，η(x)为特征点匹配调整系数，d₁(x)为左影像中特征点通过遍历找出的右影像中与之最近的特征点之间的欧式距离，d₂(x)为左影像中该特征点通过遍历找出的右影像中与之次近的特征点之间的欧式距离，其中，左影像为参考影像，右影像为待配准影像，若η(x)＞M，则将左右影像中的特征点视为一对特征点对，M为预设的阈值；

(2.2)记所获得的特征点对样本集X＝{x_i|i＝1,2,…,N}，其中，N为总特征点对数，设第i对特征点对之间的距离为d(i)，i＝1,2,…,N，根据

计算标准差估计值σ_d，将满足d(i)≤σ_d条件的所有特征点对加入样本优先集X'＝{x'_j|j＝1,2,…,N'}，其中，N'为符合d(i)≤σ_d条件的特征点对数；

步骤(3)中，构建适应度函数，通过样本优先集计算航空影像和参考影像之间的最佳转换矩阵实现特征点对的映射的方法如下：

(3.1)建立模型适应度函数

其中

代替根据余差R_i大小所产生的内点比率来对变换矩阵H进行选择，其中，x_i＝(x_i,y_i)与x'_i＝(x_i',y_i')为一对特征点对，

为样本集X中样本对由优先样本集X'采样所求得的变换矩阵H的余差，d()为两点之间的欧式影像距离，v＝R_max-R_min，R_max和R_min为余差最大值和最小值，/>

为标准差估计，取值为/>

其中，m为影像间变换矩阵参数估计所需最小样本数，i＝1,2,…,N；

(3.2)输入匹配点对样本集X＝{x_i|i＝1,2,…,N}，样本优先集对X'＝{x'_j|j＝1,2,…,N'}，设置样本中内点比率ε，用作计算采样次数，置信率η₀，内外点判定阈值τ，影像间变换矩阵参数估计所需最小样本数m，迭代次数l＝0；

(3.3)计算采样次数k＝log(1-η₀)/log(1-ε^m)，并从样本优先集X'中随机抽取四个点对(x_a,x'_a)建立方程组

a＝1,2,3,4，经过k次采样，即k次随机抽取四个点对，求得影像间得变换矩阵集Ω_H＝(H₁,H₂,…H_k)，变换矩阵H的求解步骤入下：

变换矩阵H为大小为3×3的矩阵，h₃₃＝1，剩余变换矩阵H的8个未知数需要8个方程来求解，每对特征点对(x_a,x'_a)可提供两个方程，因此，四对特征点对即可正确求解出变换矩阵H；

(3.4)对变换矩阵集Ω_H＝(H₁,H₂,…H_k)计算不同变换矩阵H_h，h＝1,2,…,k，所对应的余差R_ih，并利用函数F(H)来检验当前变换矩阵H的模型适应度，遍历k个变换矩阵H之后，选择函数F(H)最大值所对应的变换矩阵H，记为H_max；

(3.5)使用H_max计算输入匹配点对样本集X＝{x_i|i＝1,2,…,N}对应的余差R_i，i＝1,2,…,N，与设定的阈值τ作比较，若计算的每个特征点对所对应的余差项R_i小于τ的特征点对为内点，反之为外点，n_inlier为采样集X中的内点数目，更新内点比率

(3.6)迭代次数l＝l+1，与采样次数k进行比较，若l＞k，停止循环，并输出内点比率最高的变换矩阵作为最终结果，否则重复步骤(3.3)-(3.6)。

Images