CN105160686A - 一种基于改进sift算子的低空多视角遥感影像匹配方法 - Google Patents

Abstract

本发明公开了一种基于改进SIFT算子的低空多视角遥感影像匹配方法，首先采用优化的DoG算子对第一影像和第二影像进行特征点检测，然后利用基于局部区域采样模拟的SIFT描述子对特征点进行描述形成特征向量，最后采用NNDR策略形成初始匹配特征点，并利用基于极几何约束的RANSAC算法对匹配点进行提纯。本发明能够有效地解决低空遥感影像匹配中存在的多视角以及弱纹理等问题。

Description

一种基于改进SIFT算子的低空多视角遥感影像匹配方法

技术领域

本发明属于遥感影像处理技术领域，涉及一种遥感影像匹配方法，涉及一种基于改进SIFT算子的低空多视角遥感影像匹配方法。

背景技术

影像匹配一直是摄影测量和遥感领域中的研究重点和热点，是影像定向，正射影像制作，影像自动配准以及拼接、三维重建的关键步骤之一，直接影响到后续产品制作的精度和效果。与传统遥感影像相比，低空遥感影像的获取环境更加复杂，特别是从不同视角获取的影像，存在着严重的几何畸变、相似纹理、纹理断裂、阴影、遮挡等问题，大大减小了同一目标在具有重叠区域的不同影像上的相似性，给低空影像的匹配带来极大的困难。因此，如何从具有严重几何畸变的低空多视角影像中获得稳定可靠的匹配点具有重大意义。

影像匹配算法主要分为两类：基于影像灰度的匹配和基于影像特征的匹配。基于影像灰度的匹配方法直接对影像上一定窗口范围内的像素进行操作，通过窗口范围内像素的灰度相关性来进行同名点匹配，该类方法虽然在影像辐射和几何形变小的情况下能获得较高的匹配精度，然而其对影像的灰度变化比较敏感，不适应于低空多视角影像中的几何畸变、相似纹理、纹理断裂、阴影、遮挡等问题。基于影像特征的匹配方法通过提取影像中的点、线、面等几何特征来进行匹配，能较好地补充基于灰度的匹配方法的不足。

目前，基于特征的匹配方法中，应用最广泛的是Lowe提出的SIFT(ScaleInvariantFeatureTransformation)算法，能较好地处理影像几何变形中的旋转、尺度缩放和平移变化的问题，但是在大视角影像对的问题上，仍然无法取得令人满意的效果。为了改善针对不同视角影像的匹配效果，学者们对SIFT算子进行多种研究和改进，提出了ASIFT(Affine-SITT)，ISIFT(Iterative-SITT)以及MM-SIFT(Multi-resolutionMSERbasedSIFT)等多视角影像匹配方法。ASIFT根据相机仿射模型，对经度、纬度和尺度三个参数进行采样，以模拟影像视角变换，从而建立仿射空间影像，然后在所有模拟影像的上进行SIFT特征点检测和匹配，因而该方法具有完全仿射不变性，即使在较大的视角变化情况下，也能获得较好的匹配结果，与此同时，其时间复杂度随之增加，限制了其应用范围。为了平衡匹配结果和时间效率，ISIFT方法利用迭代估计的几何变换模型模拟其中一幅影像，然后利用仅有的一幅模拟影像与另一幅影像进行SIFT匹配，大大提高了匹配效率。但是当影像间几何畸变较大时，该算法会因几何模型求解失败而失效。MM-SIFT根据影像中对应局部区域满足同一个仿射模型的理论，对两幅影像中的局部椭圆区域进行圆形归一化，然后在归一化的椭圆区域内进行SIFT匹配。该方法大大提高了时间效率，而且在影像间存在大视角变化的情况下，仍然能够得到一些正确的匹配点，然而其对提取的局部区域依赖性大，其采用MSER算子提取局部区域时，并没有最佳阈值，且现有的区域分割算法用于多视角影像尚不成熟。

发明内容

针对现有影像匹配技术存在的不足，本发明的目的是提供一种时间效率高，且对影像几何畸变较稳定的基于SIFT改进的低空多视角影像匹配方法。

本发明所采用的技术方案是：一种基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于，包括以下步骤：

步骤1：使用自适应阈值的DoG算子检测第一影像和第二影像中的特征点；

步骤2：基于改进的SIFT描述子对特征点进描述；

对于第一影像和第二影像上提取的特征点，在其邻域内分别采用基于圆形区域采样和基于椭圆区域采样的灰度值对特征点进行描述，最终形成128维特征向量；

步骤3：基于特征点对应的特征向量之间的最近邻距离与次近邻距离的比值(NNDR，NearestNeighborDistanceRatio)这一策略获取初始匹配点集；

步骤4：基于极几何约束的RANSAC算法对初始匹配点集进行提纯，得到最终匹配点。

作为优选，步骤1中所述的使用自适应阈值的DoG算子检测第一影像和第二影像中的特征点，是根据以特征点为中心的一定窗口范围内的灰度信息自适应确定阈值T_contrast，其表达式为：

$T\_contrast=\frac{1}{4\·r^2}\underset{i=m-r}{\overset{i=m+r}{\Σ}}\underset{j=n-r}{\overset{j=n+r}{\Σ}}\sqrt{\frac{1}{N^2}\underset{k=i-\frac{N-1}{2}}{\overset{k=i+\frac{N-1}{2}}{\Σ}}\underset{l=j-\frac{N-1}{2}}{\overset{l=j+\frac{N-1}{2}}{\Σ}}\left(x_{kl}-{\stackrel{\‾}{x}}_N\right)};$

其中，(m,n)是特征点在影像中的坐标，N是窗口大小，x_kl是像素在DoG尺度空间的灰度值，是窗口范围内像素的灰度平均值，r由影像的尺度来决定，可按如下公式计算：

$r=4\·{\mathrm{\σ}}_0\·2^{\frac{2S_i}{M_I}};$

其中，S_i为当前特征点所在的尺度，M_I表示DoG尺度空间中影像的数量，σ₀表示影像平滑的程度。

作为优选，步骤2中所述的采用基于圆形区域采样的灰度值对特征点进行描述，对于第一影像的任意特征点(x_i,y_i)，其邻域内的圆形区域采样点的表达式为：

$x_i^{jk}=x_i+r_{ij}sin\left({\mathrm{\θ}}_k\right);$

$y_i^{jk}=y_i+r_{ij}cos\left({\mathrm{\θ}}_k\right);$

其中，r_ij为采样半径，θ_k为采样角度。

作为优选，步骤2中所述的采用基于椭圆区域采样的灰度值对特征点进行描述，对于第二影像的任意特征点(x′_i,y′_i)，其邻域内的椭圆区域采样点(x′_i ^jk,y′_i ^jk)的表达式为：

x′_i ^jk＝x′_i+a_ijcos(θ_k)cos(θ₀)+b_ijcos(θ_k)sin(θ₀)；

y′_i ^jk＝y′_i+a_ijsin(θ_k)sin(θ₀)-b_ijcos(θ_k)cos(θ₀)；

其中，a_ij为椭圆采样长轴，b_ij为椭圆采样短轴，b_ij＝a_ij/γ，γ≥1,θ₀为椭圆方向，θ_k为采样角度。

作为优选，采样点的梯度计算表达式如下：

$d_{\mathrm{\θ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{j(k-1)},y_i^{j(k-1)});$

$d_{\mathrm{\γ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{(j-1)k},y_i^{(j-1)k});$

$m(x_i^{jk},y_i^{jk})=\sqrt{d_{\mathrm{\θ}}^2(x_i^{jk},y_i^{jk})+d_{\mathrm{\γ}}^2(x_i^{jk},y_i^{jk})}.$

其中，d_.(.,.)表示该点θ和γ方向的梯度值，m(.,.)表示该点的梯度方向，采样点的灰度值则采用双线性内插的方式求取，记为

作为优选，采用迭代求解的策略确定椭圆参数γ和θ₀，具体步骤如下：

步骤2.1：从第一影像中随机选取M个特征点，其特征向量记为

步骤2.2：在第二影像中进行椭圆采样时，设γ＝γ^s，θ₀＝0，获得所有特征点特征向量，记为使用最邻近最小距离方法匹配第一影像和第二影像，记正确匹配个数为n；

步骤2.3：改变γ和θ₀值，其中γ∈[γ^s,γ^e]，θ₀∈[0°,90°]，获得多种模拟情况下的D′和相应的正确匹配点n_i；

步骤2.4：比较多种情况下的匹配结果，当获得max(n_i)时，对应的γ^*和θ^*即为最优解；

步骤2.5：在第一影像中挑选max(n_i)时的特征点作为初始点，重新进行步骤2.1-步骤2.5，直到γ^*和θ^*的值趋于稳定。

作为优选，步骤2中完成邻近点采样后，采样点的灰度值由双线性内插方法求取。

与现有技术相比，本发明具有如下特点和有益效果，

1、在特征点提取阶段采用自适应阈值，使得DoG算子能够适应不同纹理的影像，从而保证了后续匹配结果；

2、在特征描述阶段，采用椭圆采样的方式去模拟视角变化，增加了同一目标在不同视角影像上的相似性，大大改善了多视角及弱纹理区域的匹配效果；

3、采用对局部区域进行模拟而不是对整幅影像进行处理，在提高匹配精度的同时，兼顾了时间效率。

附图说明

图1：为本发明实施例的流程图；

图2：为本发明实施例中的圆形区域采样示意图；

图3：为本发明实施例中的椭圆区域采样示意图。

具体实施方式

为了便于本领域普通技术人员理解和实施本发明，下面结合附图及实施例对本发明作进一步的详细描述，应当理解，此处所描述的实施示例仅用于说明和解释本发明，并不用于限定本发明。

本发明提出的一种基于改进SIFT算子的低空多视角遥感影像匹配方法，首先采用优化的DoG算子在参考影像和待匹配影像上分别提取特征点；其次采用基于局部区域模拟的SIFT描述子进行特征描述；然后基于最近距离与次邻近距离的比值策略对特征点进行匹配，并采用基于极几何的RANSAC对初始匹配结果进行提纯，得到最终匹配点。

具体实施时，本发明可采用计算机软件技术实现自动运行流程。为详细说明本发明技术方案，请见图1，本发明提供的一种基于改进SIFT算子的低空多视角遥感影像匹配方法，包括以下步骤：

步骤1：使用自适应阈值的DoG算子检测第一影像和第二影像中的特征点；

为了适应不同细节的影像，本发明不采用原始DoG中的固定阈值，而是根据以特征点为中心的一定窗口范围内的灰度信息自适应确定阈值T_contrast，其表达式为：

$T\_contrast=\frac{1}{4\·r^2}\underset{i=m-r}{\overset{i=m+r}{\Σ}}\underset{j=n-r}{\overset{j=n+r}{\Σ}}\sqrt{\frac{1}{N^2}\underset{k=i-\frac{N-1}{2}}{\overset{k=i+\frac{N-1}{2}}{\Σ}}\underset{l=j-\frac{N-1}{2}}{\overset{l=j+\frac{N-1}{2}}{\Σ}}\left(x_{kl}-{\stackrel{\‾}{x}}_N\right)};$

其中，(m,n)是特征点在影像中的坐标，N是窗口大小，x_kl是像素在DoG尺度空间的灰度值，是窗口范围内像素的灰度平均值，r由影像的尺度来决定，可按如下公式计算：

$r=4\·{\mathrm{\σ}}_0\·2^{\frac{2S_i}{M_I}};$

其中，S_i为当前特征点所在的尺度，M_I表示DoG尺度空间中影像的数量，σ₀表示影像平滑的程度。

步骤2：基于改进的SIFT描述子对特征点进描述；

对于第一影像和第二影像上提取的特征点，在其邻域内分别采用基于圆形区域采样和基于椭圆区域采样的灰度值对特征点进行描述，最终形成128维特征向量；

本发明的研究对象是低空多视角遥感影像，由于成像条件的复杂性，不同视角影像之间存在着严重的几何畸变,使得整体影像不满足同一个仿射变换,然而,具有重叠区域的多视角影像对应的局部平面区域可以用同一个仿射变换表示。在影像几何中,圆和椭圆之间的关系可以用来刻画一个仿射变换,因此,本发明在特征描述阶段，一幅影像中的特征点描述范围采用圆形区域，另一幅影像中的描述范围采用模拟的椭圆区域,目的是模拟影像的视角变换，增加特征点邻域内的灰度信息的相似性，从而改善多视角影像的匹配效果。

当对其中一幅影像的特征点用圆形区域来描述时，对每一个特征点(x_i,y_i)，其局部圆形描述区域的大小可用半径表示，其表达式如下。

$w_i=(3{\mathrm{\σ}}_i\×\sqrt{2}\×5+1)/2;$

其中，σ_i为特征点所在尺度。

与原始SIFT直接使用相邻整数坐标内的灰度值不同，本发明使用采样和插值的方式来获取邻域内的灰度信息，采样方式如图2所示。对局部圆形描述区域内的半径和角度的采样表达式分别为：

r_ij＝jw_i/s；

θ_k＝kπ/180；

其中，s为总的半径采样次数，j＝1,...,s，k为采样角度，k＝1,...,360°/n,...,360°。

根据采样半径和采样角度，可计算采样点的坐标：

$x_i^{jk}=x_i+r_{ij}sin\left({\mathrm{\θ}}_k\right);$

$y_i^{jk}=y_i+r_{ij}cos\left({\mathrm{\θ}}_k\right);$

其中，r_ij为采样半径，θ_k为采样角度。

采样点的灰度值则采用双线性内插的方式求取，记为双线性插值是一种常见的方法，在这里不再赘述。

当对另外一幅影像采用椭圆区域描述特征点时，特征点椭圆区域包含的内容应该尽量与同名点圆形区域包含的内容相同以保证匹配效果。对于每一个特征点(x′_i,y′_i)，根据长轴a_ij，短轴b_ij以及椭圆方向θ₀三个椭圆参数，即可确定椭圆区域的采样范围。a_ij的计算方式与r_ij相同，设b_ij＝a_ij/γ，则控制椭圆区域变化的参数为γ和θ₀。对椭圆进行采样时，采样方式如图3所示，每一个采样点(x′_i ^jk,y′_i ^jk)的坐标可根据如下公式计算：

x′_i ^jk＝x′_i+a_ijcos(θ_k)cos(θ₀)+b_ijcos(θ_k)sin(θ₀)；

y′_i ^jk＝y′_i+a_ijsin(θ_k)sin(θ₀)-b_ijcos(θ_k)cos(θ₀)；

其中j＝1,...,s,θ₀∈[0°,180°]，θ_k＝kπ/180，k＝1,...,360°/n,...,360°。

完成邻近点采样后，椭圆区域采样点的灰度值同样采用双线性内插的方式获取。

在计算梯度时，与原始SIFT不同，不使用相邻特征点，如(x_i,y_i)，(x_i,y_i-1)，(x_i-1,y_i)等，而是采用上述相邻采样点，如 等，以圆形描述区域为例，梯度具体计算表达式如下：

$d_{\mathrm{\θ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{j(k-1)},y_i^{j(k-1)});$

$d_{\mathrm{\γ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{(j-1)k},y_i^{(j-1)k});$

$m(x_i^{jk},y_i^{jk})=\sqrt{d_{\mathrm{\θ}}^2(x_i^{jk},y_i^{jk})+d_{\mathrm{\γ}}^2(x_i^{jk},y_i^{jk})}.$

其中，d_.(.,.)表示该点θ和γ方向的梯度值，m(.,.)表示该点的梯度方向。

根据以上公式获得的灰度梯度，可按照原始SIFT的计算方式计算梯度方向直方图，进而获得各特征点的特征向量。

上述基于局部区域模拟采样的SIFT描述子能否有效解决大视角下的几何形变问题，关键在于椭圆区域对圆形区域的模拟程度，即控制椭圆区域变化的参数为γ和θ₀。本发明采用迭代策略求解γ和θ₀的值，以保证匹配效果，具体步骤如下：

步骤2.1：从第一影像中随机选取M个特征点，其特征向量记为

步骤2.2：在第二影像中进行椭圆采样时，设γ＝γ^s，θ₀＝0，获得所有特征点特征向量，记为使用最邻近最小距离方法匹配第一影像和第二影像，记正确匹配个数为n；

步骤2.3：改变γ和θ₀值，其中γ∈[γ^s,γ^e]，θ₀∈[0°,90°]，获得多种模拟情况下的D′和相应的正确匹配点n_i；

步骤2.4：比较多种情况下的匹配结果，当获得max(n_i)时，对应的γ^*和θ^*即为最优解；

步骤2.5：在第一影像中挑选max(n_i)时的特征点作为初始点，重新进行步骤2.1-步骤2.5，直到γ^*和θ^*的值趋于稳定。

步骤3：采用NNDR匹配策略匹配特征点

采用最近邻与次近邻比值(NNDR，NearestNeighborDistanceRatio)的策略进行匹配，具体实施时，计算各特征描述符之间的欧式距离，对于每一个特征，如果与其最相近的两个特征的距离满足最小距离与次小距离的比值小于一定的阈值△，则认为该匹配对是一对正确的匹配对；

步骤4：利用基于极几何约束的RANSAC剔除错误匹配点

完成特征点匹配后，在初始匹配结果中，会存在一些错误匹配。传统方法是利用RANSAC(RandomSampleConsensus)估计影像间的仿射变换矩阵来剔除错误匹配，但是对于大视角影像来说，整幅影像并不适用于同一个仿射变换，因此，在具体实施时，采用极几何约束的方法，通过估计影像像对之间的基本矩阵来剔除错误匹配。

应当理解的是，本说明书未详细阐述的部分均属于现有技术。

应当理解的是，上述针对较佳实施例的描述较为详细，并不能因此而认为是对本发明专利保护范围的限制，本领域的普通技术人员在本发明的启示下，在不脱离本发明权利要求所保护的范围情况下，还可以做出替换或变形，均落入本发明的保护范围之内，本发明的请求保护范围应以所附权利要求为准。

Claims (7)

1.一种基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于，包括以下步骤：

步骤1：使用自适应阈值的DoG算子检测第一影像和第二影像中的特征点；

步骤2：基于改进的SIFT描述子对特征点进描述；

对于第一影像和第二影像上提取的特征点，在其邻域内分别采用基于圆形区域采样和基于椭圆区域采样的灰度值对特征点进行描述，最终形成128维特征向量；

步骤3：基于特征点对应的特征向量之间的最近邻距离与次近邻距离的比值这一策略获取初始匹配点集；

步骤4：基于极几何约束的RANSAC算法对初始匹配点集进行提纯，得到最终匹配点。

2.根据权利要求1所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：步骤1中所述的使用自适应阈值的DoG算子检测第一影像和第二影像中的特征点，是根据以特征点为中心的一定窗口范围内的灰度信息自适应确定阈值T_contrast，其表达式为：

$T\_contrast=\frac{1}{4\·r^2}\underset{i=m-r}{\overset{i=m+r}{\Σ}}\underset{j=n-r}{\overset{j=n+r}{\Σ}}\sqrt{\frac{1}{N^2}\underset{k=i-\frac{N-1}{2}}{\overset{k=i+\frac{N-1}{2}}{\Σ}}\underset{l=j-\frac{N-1}{2}}{\overset{l=j+\frac{N-1}{2}}{\Σ}}\left(x_{kl}-{\stackrel{\‾}{x}}_N\right)};$

其中，(m,n)是特征点在影像中的坐标，N是窗口大小，x_kl是像素在DoG尺度空间的灰度值，是窗口范围内像素的灰度平均值，r由影像的尺度来决定，可按如下公式计算：

$r=4\·{\mathrm{\σ}}_0\·2^{\frac{2S_i}{M_I}};$

其中，S_i为当前特征点所在的尺度，M_I表示DoG尺度空间中影像的数量，σ₀表示影像平滑的程度。

3.根据权利要求1所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：步骤2中所述的采用基于圆形区域采样的灰度值对特征点进行描述，对于第一影像的任意特征点(x_i,y_i)，其邻域内的圆形区域采样点的表达式为：

$x_i^{jk}=x_i+r_{ij}sin\left({\mathrm{\θ}}_k\right);$

$y_i^{jk}=y_i+r_{ij}cos\left({\mathrm{\θ}}_k\right);$

其中，r_ij为采样半径，θ_k为采样角度。

4.根据权利要求1所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：步骤2中所述的采用基于椭圆区域采样的灰度值对特征点进行描述，对于第二影像的任意特征点(x′_i,y′_i)，其邻域内的椭圆区域采样点(x′_i ^jk,y′_i ^jk)的表达式为：

x′_i ^jk＝x′_i+a_ijcos(θ_k)cos(θ₀)+b_ijcos(θ_k)sin(θ₀)；

y′_i ^jk＝y′_i+a_ijsin(θ_k)sin(θ₀)-b_ijcos(θ_k)cos(θ₀)；

其中，a_ij为椭圆采样长轴，b_ij为椭圆采样短轴，b_ij＝a_ij/γ，γ≥1,θ₀为椭圆方向，θ_k为采样角度。

5.根据权利要求1所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：步骤2中所述的基于区域采样的特征点梯度计算表达式如下：

$d_{\mathrm{\θ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{j(k-1)},y_i^{j(k-1)});$

$d_{\mathrm{\γ}}(x_i^{jk},y_i^{jk})=I(x_i^{jk},y_i^{jk})-I(x_i^{(j-1)k},y_i^{(j-1)k});$

$m(x_i^{jk},y_i^{jk})=\sqrt{d_{\mathrm{\θ}}^2(x_i^{jk},y_i^{jk})+d_{\mathrm{\γ}}^2(x_i^{jk},y_i^{jk})};$

其中，d(.,.)表示该点θ和γ方向的梯度值，m(.,.)表示该点的梯度方向，采样点的灰度值则采用双线性内插的方式求取，记为其中*为通用符。

6.根据权利要求4所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：采用迭代求解的策略确定椭圆参数γ和θ₀，具体步骤如下：

步骤2.1：从第一影像中随机选取M个特征点，其特征向量记为

步骤2.2：在第二影像中进行椭圆采样时，设γ＝γ^s，θ₀＝0，获得所有特征点特征向量，记为使用最邻近距离方法匹配第一影像和第二影像，记正确匹配个数为n；

步骤2.3：改变γ和θ₀值，其中γ∈[γ^s,γ^e]，θ₀∈[0°,90°]，获得多种模拟情况下的D′和相应的正确匹配点n_i；

步骤2.4：比较多种情况下的匹配结果，当获得max(n_i)时，对应的γ^*和θ^*即为最优解；

步骤2.5：在第一影像中挑选max(n_i)时的特征点作为初始点，重新进行步骤2.1-步骤2.5，直到γ^*和θ^*的值趋于稳定。

7.根据权利要求1所述的基于改进SIFT算子的低空多视角遥感影像匹配方法，其特征在于：步骤2中完成邻近点采样后，采样点的灰度值由双线性内插方法求取。