CN103514606A - 一种异源遥感影像配准方法 - Google Patents

Abstract

本发明公开了一种异源遥感影像配准方法，核心思想是基于多尺度匹配以直线交点为基元、采用联合Voronoi图与谱图的点匹配方法，集成迭代的特征提取与匹配策略，克服了现有方法严重依赖特征提取、可靠性差、精度不高等问题。实现过程为：对原始影像进行多尺度分析，在最粗尺度上提取直线并获取交点；对交点集用联合Voronoi图与谱图的点匹配方法获取同名点对；检测匹配结果是否合格:若合格，转入下一步，否则自适应调整参数，重新进行直线提取与点集匹配；对待配准影像进行初变换，且分别提取直线特征，寻找同名直线段，获取候选同名点对；利用KNN图得到精确的匹配点对，求解变换参数。本发明主要用于可见光、红外和合成孔径雷达（SAR）等异源遥感影像的配准。

Description

一种异源遥感影像配准方法

技术领域

本发明属于遥感影像处理技术领域，涉及一种异源遥感影像配准方法，尤其是涉及一种基于迭代的线特征和联合Voronoi图与谱图点匹配的异源遥感影像配准方法。

背景技术

图像配准技术是将相同地区，在不同时刻、不同视角、不同传感器或不同光照条件下拍摄的图像进行空间对准的过程。其中，异源遥感图像的自动配准一直是遥感图像处理领域需要解决的难题，因为不同传感器的影像能够反映地物的不同特征，如光学图像反映了地物在可见光和近红外波长范围内的反射光谱特性，SAR图像反映了地物对电磁波的后向散射特性，红外影像反映了地物的热辐射特性，使得异源遥感影像所反映的地物信息往往具有一定的互补性，因而对这些异源数据根据应用需求进行数据融合，可以为决策者们提供更加准确、全面和丰富的信息，其中高精度的影像配准是数据融合的重要前提。然而，由于异源遥感影像的成像机理不同，影像间往往存在着较大的辐射和几何差异，甚至遭受严重的噪声干扰，同一场景在异源影像上可能呈现完全不同的图像，这些差异导致传统的图像配准方法一般无法直接应用在异源图像的配准。

目前异源遥感影像的配准方法大致可分为两类：基于区域的配准方法和基于特征的配准方法。

（1）基于区域的配准方法：

通常选择图像的某一区域或者整幅图像作为一种特殊的图像特征，其核心是匹配相似度的度量。常用的匹配相似度量包括：归一化互相关方法（normalized cross correlation coefficientmethods）、傅立叶方法（Fourier methods）、互信息方法（mutual information methods）、交叉累积剩余熵（cross-cumulative residual entropy）。基于区域的配准方法大多利用影像的灰度信息，然而异源遥感影像的不同成像机理，使得同一地物灰度可能完全不同，因此大多基于区域的配准方法难以获得较好的配准结果，甚至有学者断言，基于区域配准方法不可能直接应用于异源遥感影像的配准（Zitova,2003;Inglada,2004;Hel-Or,2011）,且一般寻求相似度最大的过程是一个庞大的搜索过程，计算效率一般较慢。

（2）基于特征的配准方法：

通常先分别从两幅图像中提取一些几何特征作为配准基元，然后对配准基元进行匹配。由于基于特征的匹配方法利用特征对象进行匹配，并不直接对图像的灰度信息进行操作，受光照、噪声、几何形变等影响较小，算法稳健性更强，因而更适用于异源遥感影像间的配准。

常用的三种特征为：点、线和面。

点特征一般包括角点、高曲率点等，通常采用各种角点检测算法提取。然而点特征不易精确定位，且点特征含有的信息有限，匹配困难。线特征和面特征含有更多信息，但是图像中并不一定能够提取足够的面特征，如城区影像，且SAR影像受噪声影响、红外影像对比度不明显等导致难以找到稳定的面特征，同时由于分割原因即使存在面特征也很难获得满意的特征，因而基于面特征的方法具有一定的局限性；而影像中的线特征通常比较丰富，采用线特征作为配准基元将是一个不错的选择。

根据特征描述的不同，线特征可以分为以下三类：曲线特征、边缘特征和直线段特征。然而，线特征在提取过程中容易出现断裂、提取不完整等问题，直接使用线特征不利于后续匹配。例如，轮廓分明、未断裂的曲线或边缘很难在SAR影像和红外影像上提取，而在光学影像上通常可以提取完整的曲线和边缘；另外，由于受光照影响、视角和成像方式的不同，SAR影像、红外影像和光学影像上提取的直线段特征通常在位置、长度、数量等方面也存在着较大差异，在这些情况下很难定义一种相似性测度来度量两个特征集间的相似性。然而，尽管提取的直线段特征是断裂的，位置不同，长度也不尽相同，但同名直线段是在同一条直线上的。因此，可以将直线段特征转化为直线交点作为匹配基元，它不需要产生交点的两对同名直线在两幅图像上处于同一位置。

基于直线交点匹配的关键是如何找到正确的匹配点对。基于谱图的点集匹配方法被认为是最有效的方法之一。然而，传统的谱图点集匹配方法存在两个问题：（1）对噪声敏感：噪声通常来源于影像获取和特征提取过程中，它会给点的位置带来偏差；（2）对局外点敏感：局外点的存在容易干扰正确匹配关系的判断。传统谱图方法是从整体出发，对点集结构和空间分布进行分析，加入点集的局部约束可以提高点集匹配精度。

同时，几乎所有的基于特征的配准方法很大程度上依赖于特征提取算法：特征匹配结果受特征提取的影响。虽然存在成功的异源影像匹配算法，但它们通常是对特定的影像有效，很难成为通用的方法。因此，多次特征提取与匹配的策略是需要的。

发明内容

为了解决上述的技术问题，本发明提出了一种基于迭代的线特征和联合Voronoi图与谱图点匹配方法的异源遥感影像配准方法；本发明是基于直线交点的配准方法，为了提高方法的健壮性和准确性，本发明同时提出了迭代的特征提取与匹配策略和多层次的匹配策略。本发明主要将可见光、红外和合成孔径雷达（SAR）等异源遥感影像的配准作为研究目标。

本发明所采用的技术方案是：一种异源遥感影像配准方法，其特征在于，包括以下步骤：

步骤1：对基准影像和待配准影像利用小波变换进行多尺度分析，得到最粗尺度基准影像和最粗尺度待配准影像，所述的基准影像和待配准影像为光学影像、红外影像或SAR影像；

步骤2：在最粗尺度基准影像和最粗尺度待配准影像上进行直线提取并获取交点集；

步骤3：对步骤2所得交点集采用联合Voronoi图与谱图的点匹配方法获取同名点对；

步骤4.判断，所述的同名点对的数量是否小于预定阀值m？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则进入步骤5；

步骤5.对所述的同名点对进行均方根误差计算，并判断，所述的均方根误差是否大于预定阀值n？

如果是，则排除具有最大误差的同名点对，对剩余的同名点对重新计算均方根误差，然后进入步骤6；

如果否，则进入步骤7；

步骤6：判断，所述的剩余的同名点对的数量是否小于预定阀值m？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则回转执行所述的步骤5；

步骤7：将所述的同名点对映射到原始影像上，计算初始变换参数并得到初变换影像；

步骤8：对基准影像和待配准影像进行预处理，包括对光学影像和红外影像进行高斯滤波处理，对SAR影像进行Frost滤波处理；

步骤9：在对步骤8处理后的影像上分别提取直线特征，根据初始变换参数计算寻找同名直线段，进而获得候选同名点对；

步骤10：利用KNN图从结构上得到精确的匹配点对，对匹配点对采用多项式变换模型求解变换参数，得到最终配准结果。

作为优选，步骤2中所述的在最粗尺度基准影像和最粗尺度待配准影像上进行直线提取并获取交点集，其具体实现包括以下子步骤：

步骤2.1：低层特征提取：在最粗尺度基准影像和最粗尺度待配准影像上提取线特征，其中光学影像和红外影像采用LSD算法，SAR影像先采用基于Gaussian-Gamma-Shaped Bi窗口算法提取边缘特征，再采用Hough变换提取直线特征；

步骤2.2：低层直线段预处理：对所有直线按长度从大到小排序，从最长的直线起统计斜率相近的线段，并计算线段间的距离，若距离小于阈值d，刚认为是属于同一条直线段，删除除长度最长的直线外的所有直线，只保留同一方向上最长的直线；

步骤2.3：低层交点求取：对预处理后的直线段求交点，为了防止斜率相近的线段相交得到极点，对线集中夹角在10°内的线段不进行交点求取，将最终获得的交点集作为控制点集。

作为优选，步骤3中所述的对步骤2所得交点集采用联合Voronoi图与谱图的点匹配方法获取同名点对，其具体实现包括以下子步骤：

步骤3.1：建立Voronoi图：设基准影像交点集为P=[p₁,p₂,…,p_m],待配准影像交点集为Q=[q₁,q₂,…,q_n]，首先采用构建Delaunay三角网，再利用对偶的方法由三角网生成Voronoi图；

步骤3.2：生成邻近矩阵：假设生成的Voronoi多边形分别为：VP=[vp₁,vp₂,…,vp_m]，VQ=[vq₁,vq₂,…,vq_n]，计算每个交点集对应的Voronoi多边形间的Laplace矩阵，得到基准影像对应的邻近矩阵HP：

$\mathrm{HP}=H_{\mathrm{ij}}=\left[h_{\mathrm{ij}}\right]=\left\{\begin{array}{cc}{-r}_{\mathrm{ij}}^2& \mathrm{if}(i\≠j)\\ -\underset{k\≠i}{\mathrm{\Σ}}{h}_{\mathrm{ik}}& \mathrm{if}(i=j)\end{array}\right.$

其中，

$r_{\mathrm{ij}}=e^{-\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)/2{\mathrm{\σ}}^2}$

$\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)=\underset{a\∈{\mathrm{vp}}_i}{\max }\left\{\underset{b\∈{\mathrm{vp}}_j}{\min }\right\{d(a,b)\left\}\right\}$

a,b为任意多边形vp_i,vp_j上的点，d(a,b)为点a,b之间的欧式距离；hausdorff(vp_i,vp_j)为多边形vp_i,vp_j间的Hausdorff距离，用来度量两个多边形间的相似性；σ为控制参数；

采用上述的方法，得到待配准影像对应的邻近矩阵HQ；

步骤3.3：特征值分解：

HP=UΔ₁U^T

HQ=VΔ₂V^T

其中，Δ₁=diag(λ₁,λ₂,…,λ_m)，Δ₂=diag(γ₁,γ₂,…,γ_n)，其对角元素分别为HP，HQ的特征值，U=(u₁,u₂,…,u_m)，u_i为HP的对应于特征值λ_i的特征向量，V=(v₁,v₂,…,v_n)，v_i为HQ的对应于特征值γ_i的特征向量；若交点集P中含有m个特征点，Q中含有n个特征点，且m＜n，则删除V中的n-m列，以达到V与U的维数一致，反之亦然；

步骤3.4：特征干扰与符号校正：

对初始邻近矩阵添加特征干扰：

K₁=Udiag(δ,δ²,…,δ^m)U^Tδ=rand(0,1)

K₂=Vdiag(ζ,ζ²,…,ζ^m)V^Tζ=rand(0,1)

对特征向量进行符号校正：

$v_i^{\′}=\left\{\begin{array}{cc}{v}_i& \mathrm{if}\left(\right||u_i+v_i||>{\left|\right|u}_i-v_i|\left|\right)\\ -v_i& \mathrm{otherwise}\end{array}\right.;$

步骤3.5：角度一致性约束：假设基准影像中的两条直线段

生成控制点P_o，待配准影像中的两条直线段

生成控制点P_s，如果直线段

间的夹角与直线段

间的夹角之差过大，则点P_o与P_s不可能成为同名点，因而，将角度一致性约束定义为：

$A(P_o,P_s)=\left\{\begin{array}{cc}1& \mathrm{if}\left(\right|\mathrm{\&theta;}(L_i^O,L_j^O)-\mathrm{\&theta;}(L_m^S,L_n^S)|<t)\\ 0& \mathrm{otherwise}\end{array}\right.$

其中，t为控制变形敏感性的阈值，t取值为5°；

步骤3.6：匹配判断：

Z=A*(U*V'^T)

其中V’为V经符号校正后的矩阵；如果Zi_j是它所有行与列的最大元素，则认为点集P中的第i个点与点集Q中的第j个点匹配。

作为优选，步骤4中所述的预定阀值m=3。

作为优选，步骤5中所述的预定阀值n=5。

作为优选，步骤4中所述的自适应调整直线提取参数，所述的直线提取参数包括：线支持区搜索过程中的角度容忍值τ、检测窗口的宽度α及两个检测窗口间的间隔β，其初始值设定为τ=22.5，α=3，β=1.5，在迭代参数调整过程中，若基准影像中检测到更多直线段，则固定基准影像参数不变，将待配准影像参数变量α与β均增加2，以检测到更多直线段；否则，固定待配准影像参数不变，将基准影像参数变量τ增加7.5。

作为优选，步骤7中所述的将所述的同名点对映射到原始影像上，计算初始变换参数并得到初变换影像；其具体实现过程为：设控制点在最粗尺度影像上的坐标为(x,y),最粗尺度影像相对于原始影像经过了L次降分，则该点映射到原始影像上的坐标为（2^L-2^L-1）*(x,y)；初始变换模型采用仿射变换模型，利用最粗尺度上的同名点对映射后的相应点对作为同名点对计算变换参数，并对待配准影像进行变换。

作为优选，步骤9中所述的在对步骤8处理后的影像上分别提取直线特征，根据初始变换参数计算寻找同名直线段，进而获得候选同名点对；其具体实现包括以下子步骤：

步骤9.1：对基准影像与待配准影像采用与低层同样的直线提取方法与直线段预处理方法；

步骤9.2：通过定义两条直线段匹配代价函数判断基准影像与待配准影像上的直线段是否为同名线段，代价函数为：

$C_{\mathrm{r\&theta;}}(s_i,o_j)=\frac{1}{\sqrt{2}}{\{{\left[\frac{d_r(s_i,o_j)}{d_{r\max }}\right]}^2+{\left[\frac{d_{\mathrm{\&theta;}}(s_i,o_j)}{d_{\mathrm{\&theta;}\max }}\right]}^2\}}^{1/2}$

其中，d_r(s_i,o_j)为直线段o_j到直线段s_i的距离，d_θ(s_i,o_j)为直线段o_j与直线段s_i之间的角度差异，d_rmax、d_θmax分别为d_r、d_θ的最大容许值，用于归一化作用，值设定为d_rmax=5像素，d_θmax=3弧度；

如果两条直线段为同名线段，则由上式得到的值趋于0；如果值大于1，则认为两条直线不可能同名；对于待配准影像上的每条直线段s(i),基准影像上与s(i)斜率相近的直线段均用来计算与s(i)的匹配代价，如果基准影像上的所有直线与s(i)的匹配代价均大于1，则认为s(i)为孤立线段，否则，基准影像上与s(i)具有最小匹配代价的线段被认为与s(i)是同名线段；进而将任意两对同名线段的交点作为候选同名点对。

本发明提出的一种基于直线交点的配准方法，配合迭代的特征提取与匹配策略和多层次的匹配策略，提高了异源遥感影像配准方法的健壮性和准确性；本发明主要用于可见光、红外和合成孔径雷达（SAR）等异源遥感影像的配准。

附图说明

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

具体实施方式

以下将结合附图和具体实施例对本发明做进一步的阐述。

请见图1，本发明所采用的技术方案是：一种异源遥感影像配准方法，其特征在于，包括以下步骤：

步骤1：对基准影像和待配准影像利用小波变换进行多尺度分析，得到最粗尺度基准影像和最粗尺度待配准影像，使得分解后的最粗尺度影像大小不超过256×256像素；所述的基准影像和待配准影像为光学影像、红外影像或SAR影像。

步骤2：在最粗尺度基准影像和最粗尺度待配准影像上进行直线提取并获取交点集，其具体实现包括以下子步骤：

步骤2.1：低层特征提取：在最粗尺度基准影像和最粗尺度待配准影像上提取线特征，其中光学影像和红外影像采用LSD算法，SAR影像先采用基于Gaussian-Gamma-Shaped Bi窗口算法提取边缘特征，再采用Hough变换提取直线特征；

步骤2.2：低层直线段预处理：对所有直线按长度从大到小排序，从最长的直线起统计斜率相近的线段，并计算线段间的距离，若距离小于阈值d，刚认为是属于同一条直线段，删除除长度最长的直线外的所有直线，只保留同一方向上最长的直线；

步骤2.3：低层交点求取：对预处理后的直线段求交点，为了防止斜率相近的线段相交得到极点，对线集中夹角在10°内的线段不进行交点求取，将最终获得的交点集作为控制点集。步骤3：对步骤2所得交点集采用联合Voronoi图与谱图的点匹配方法获取同名点对，其具体实现包括以下子步骤：

步骤3.1：建立Voronoi图：设基准影像交点集为P=[p₁,p₂,…,p_m],待配准影像交点集为Q=[q₁,q₂,…,q_n]，首先采用构建Delaunay三角网，再利用对偶的方法由三角网生成Voronoi图；

步骤3.2：生成邻近矩阵：假设生成的Voronoi多边形分别为：VP=[vp₁,vp₂,…,vp_m]，VQ=[vq₁,vq₂,…,vq_n]，计算每个交点集对应的Voronoi多边形间的Laplace矩阵，得到基准影像对应的邻近矩阵HP：

$\mathrm{HP}=H_{\mathrm{ij}}=\left[h_{\mathrm{ij}}\right]=\left\{\begin{array}{cc}{-r}_{\mathrm{ij}}^2& \mathrm{if}(i\&NotEqual;j)\\ -\underset{k\&NotEqual;i}{\mathrm{\&Sigma;}}{h}_{\mathrm{ik}}& \mathrm{if}(i=j)\end{array}\right.$

其中，

$r_{\mathrm{ij}}=e^{-\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)/2{\mathrm{\&sigma;}}^2}$

$\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)=\underset{a\&Element;{\mathrm{vp}}_i}{\max }\left\{\underset{b\&Element;{\mathrm{vp}}_j}{\min }\right\{d(a,b)\left\}\right\}$

a,b为任意多边形vp_i,vp_j上的点，d(a,b)为点a,b之间的欧式距离；hausdorff(vp_i,vp_j)为多边形vp_i,vp_j间的Hausdorff距离，用来度量两个多边形间的相似性；σ为控制参数；

采用上述的方法，得到待配准影像对应的邻近矩阵HQ；

步骤3.3：特征值分解：

HP=UΔ₁U^T

HQ=VΔ₂V^T

其中，Δ₁=diag(λ₁,λ₂,…,λ_m)，△₂=diag(γ₁,γ₂,…,γ_n)，其对角元素分别为HP，HQ的特征值，U=(u₁,u₂,…,u_m)，u_i为HP的对应于特征值λ_i的特征向量，V=(v₁,v₂,…,v_n)，v_i为HQ的对应于特征值γ_i的特征向量；若交点集P中含有m个特征点，Q中含有n个特征点，且m＜n，则删除V中的n-m列，以达到V与U的维数一致，反之亦然；

步骤3.4：特征干扰与符号校正：由于在进行特征值分解时可能带来特征值多重性问题，因而可能带来特征向量顺序变化的问题；同时，在计算特征向量时存在符号模糊的问题，即Hu_i=λ_iu_i同时H(-u_i)=λ_i(-u_i)；

对初始邻近矩阵添加特征干扰：

K₁=Udiag(δ,δ²,…,δ^m)U^Tδ=rand(0,1)

K₂=Vdiag(ζ,ζ²,…,ζ^m)V^Tζ=rand(0,1)

由于δ与ζ的取值都较小，轻微的扰动操作会使得特征值互不相同且不会对特征向量有明显影响；

对特征向量进行符号校正：

$v_i^{\&prime;}=\left\{\begin{array}{cc}{v}_i& \mathrm{if}\left(\right||u_i+v_i||>{\left|\right|u}_i-v_i|\left|\right)\\ -v_i& \mathrm{otherwise}\end{array}\right.;$

步骤3.5：角度一致性约束：假设基准影像中的两条直线段

生成控制点P_o，待配准影像中的两条直线段

生成控制点P_s，如果直线段

间的夹角与直线段

间的夹角之差过大，则点P_o与P_s不可能成为同名点，因而，将角度一致性约束定

$A(P_o,P_s)=\left\{\begin{array}{cc}1& \mathrm{if}\left(\right|\mathrm{\&theta;}(L_i^O,L_j^O)-\mathrm{\&theta;}(L_m^S,L_n^S)|<t)\\ 0& \mathrm{otherwise}\end{array}\right.$

其中，t为控制变形敏感性的阈值，t值越大，A对图像间变形的兼容性越好，但是精度会越低，t取值为5°；

步骤3.6：匹配判断：

Z=A*(U*V^'T)

其中V’为V经符号校正后的矩阵；如果Zi_j是它所有行与列的最大元素，则认为点集P中的第i个点与点集Q中的第j个点匹配。

步骤4.判断，所述的同名点对的数量是否小于预定阀值m=3？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则进入步骤5；

其中，自适应调整直线提取参数，所述的直线提取参数包括：线支持区搜索过程中的角度容忍值τ、检测窗口的宽度α及两个检测窗口间的间隔β，其中边缘检测结果直接影响着Hough变换获取直线结果，而在所采用的基于比率的边缘检测方法中，较为关键的参数为：检测窗口的宽度α及两个检测窗口间的间隔β，因而在参数调整过程中将这三个参数作为变量，其初始值设定为τ=22.5，α=3，β=1.5；在迭代参数调整过程中，若基准影像中检测到更多直线段，则固定基准影像参数不变，将待配准影像参数变量α与β均增加2，以检测到更多直线段；否则，固定待配准影像参数不变，将基准影像参数变量τ增加7.5。

步骤5：对所述的同名点对进行均方根误差（RMSE）计算，并判断，所述的均方根误差是否大于预定阀值n=5？

如果是，则排除具有最大误差的同名点对，对剩余的同名点对重新计算RMSE，然后进入步骤6；

如果否，则进入步骤7。

步骤6：判断，所述的剩余的同名点对的数量是否小于预定阀值m=3？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则回转执行所述的步骤5。

步骤7：将所述的同名点对映射到原始影像上，计算初始变换参数并得到初变换影像；其具体实现过程为：设控制点在最粗尺度影像上的坐标为(x,y),最粗尺度影像相对于原始影像经过了L次降分，则该点映射到原始影像上的坐标为（2^L-2^L-1）*(x,y)；初始变换模型采用仿射变换模型，利用最粗尺度上的同名点对映射后的相应点对作为同名点对计算变换参数，并对待配准影像进行变换。

步骤8：对基准影像和待配准影像进行预处理，包括对光学影像和红外影像进行高斯滤波处理，对SAR影像进行Frost滤波处理。

步骤9：在对步骤8处理后的影像上分别提取直线特征，根据初始变换参数计算寻找同名直线段，进而获得候选同名点对；其具体实现包括以下子步骤：

步骤9.1：对基准影像与待配准影像采用与低层同样的直线提取方法与直线段预处理方法；

步骤9.2：通过定义两条直线段匹配代价函数判断基准影像与待配准影像上的直线段是否为同名线段，代价函数为：

$C_{\mathrm{r\&theta;}}(s_i,o_j)=\frac{1}{\sqrt{2}}{\{{\left[\frac{d_r(s_i,o_j)}{d_{r\max }}\right]}^2+{\left[\frac{d_{\mathrm{\&theta;}}(s_i,o_j)}{d_{\mathrm{\&theta;}\max }}\right]}^2\}}^{1/2}$

其中，d_r(s_i,o_j)为直线段o_j到直线段s_i的距离，d_θ(s_i,o_j)为直线段o_j与直线段s_i之间的角度差异，d_rmax、d_θmax分别为d_r、d_θ的最大容许值，用于归一化作用，值设定为d_rmax=5像素，d_θmax=3弧度；

如果两条直线段为同名线段，则由上式得到的值趋于0；如果值大于1，则认为两条直线不可能同名；对于待配准影像上的每条直线段s(i),基准影像上与s(i)斜率相近的直线段均用来计算与s(i)的匹配代价，如果基准影像上的所有直线与s(i)的匹配代价均大于1，则认为s(i)为孤立线段，否则，基准影像上与s(i)具有最小匹配代价的线段被认为与s(i)是同名线段；进而将任意两对同名线段的交点作为候选同名点对。

步骤10：利用KNN图从结构上得到精确的匹配点对，对匹配点对采用多项式变换模型求解变换参数，得到最终配准结果。

Claims (8)

1.一种异源遥感影像配准方法，其特征在于，包括以下步骤：

步骤1：对基准影像和待配准影像利用小波变换进行多尺度分析，得到最粗尺度基准影像和最粗尺度待配准影像，所述的基准影像和待配准影像为光学影像、红外影像或SAR影像；

步骤2：在最粗尺度基准影像和最粗尺度待配准影像上进行直线提取并获取交点集；

步骤3：对步骤2所得交点集采用联合Voronoi图与谱图的点匹配方法获取同名点对；

步骤4.判断，所述的同名点对的数量是否小于预定阀值m？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则进入步骤5；

步骤5.对所述的同名点对进行均方根误差计算，并判断，所述的均方根误差是否大于预定阀值n？

如果是，则排除具有最大误差的同名点对，对剩余的同名点对重新计算均方根误差，然后进入步骤6；

如果否，则进入步骤7；

步骤6：判断，所述的剩余的同名点对的数量是否小于预定阀值m？

如果是，则进行自适应调整直线提取参数，回转执行所述的步骤2，重新进行直线提取与交点集匹配，直到迭代次数达到预定阈值；

如果否，则回转执行所述的步骤5；

步骤7：将所述的同名点对映射到原始影像上，计算初始变换参数并得到初变换影像；

步骤8：对基准影像和待配准影像进行预处理，包括对光学影像和红外影像进行高斯滤波处理，对SAR影像进行Frost滤波处理；

步骤9：在对步骤8处理后的影像上分别提取直线特征，根据初始变换参数计算寻找同名直线段，进而获得候选同名点对；

步骤10：利用KNN图从结构上得到精确的匹配点对，对匹配点对采用多项式变换模型求解变换参数，得到最终配准结果。

2.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤2中所述的在最粗尺度基准影像和最粗尺度待配准影像上进行直线提取并获取交点集，其具体实现包括以下子步骤：

步骤2.1：低层特征提取：在最粗尺度基准影像和最粗尺度待配准影像上提取线特征，其中光学影像和红外影像采用LSD算法，SAR影像先采用基于Gaussian-Gamma-Shaped Bi窗口算法提取边缘特征，再采用Hough变换提取直线特征；

步骤2.2：低层直线段预处理：对所有直线按长度从大到小排序，从最长的直线起统计斜率相近的线段，并计算线段间的距离，若距离小于阈值d，刚认为是属于同一条直线段，删除除长度最长的直线外的所有直线，只保留同一方向上最长的直线；

步骤2.3：低层交点求取：对预处理后的直线段求交点，为了防止斜率相近的线段相交得到极点，对线集中夹角在10°内的线段不进行交点求取，将最终获得的交点集作为控制点集。

3.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤3中所述的对步骤2所得交点集采用联合Voronoi图与谱图的点匹配方法获取同名点对，其具体实现包括以下子步骤：

步骤3.1：建立Voronoi图：设基准影像交点集为P=[p₁,p₂,…,p_m],待配准影像交点集为Q=[q₁,q₂,…,q_n]，首先采用构建Delaunay三角网，再利用对偶的方法由三角网生成Voronoi图；

步骤3.2：生成邻近矩阵：假设生成的Voronoi多边形分别为：VP=[vp₁,vp₂,…,vp_m]，VQ=[vq₁,vq₂,…,vq_n]，计算每个交点集对应的Voronoi多边形间的Laplace矩阵，得到基准影像对应的邻近矩阵HP：

$\mathrm{HP}=H_{\mathrm{ij}}=\left[h_{\mathrm{ij}}\right]=\left\{\begin{array}{cc}{-r}_{\mathrm{ij}}^2& \mathrm{if}(i\&NotEqual;j)\\ -\underset{k\&NotEqual;i}{\mathrm{\&Sigma;}}{h}_{\mathrm{ik}}& \mathrm{if}(i=j)\end{array}\right.$

其中，

$r_{\mathrm{ij}}=e^{-\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)/2{\mathrm{\&sigma;}}^2}$

$\mathrm{hausdorff}({\mathrm{vp}}_i,{\mathrm{vp}}_j)=\underset{a\&Element;{\mathrm{vp}}_i}{\max }\left\{\underset{b\&Element;{\mathrm{vp}}_j}{\min }\right\{d(a,b)\left\}\right\}$

a,b为任意多边形vp_i,vp_j上的点，d(a,b)为点a,b之间的欧式距离；hausdorff(vp_i,vp_j)为多边形vp_i,vp_j间的Hausdorff距离，用来度量两个多边形间的相似性；σ为控制参数；

采用上述的方法，得到待配准影像对应的邻近矩阵HQ；

步骤3.3：特征值分解：

HP=UΔ₁U^T

HQ=VΔ₂V^T

其中，Δ₁=diag(λ₁,λ₂,…,λ_m)，Δ₂=diag(γ₁,γ₂,…,γ_n)，其对角元素分别为HP，HQ的特征值，U=(u₁,u₂,…,u_m)，u_i为HP的对应于特征值λ_i的特征向量，V=(v₁,v₂,…,v_n)，v_i为HQ的对应于特征值γ_i的特征向量；若交点集P中含有m个特征点，Q中含有n个特征点，且m＜n，则删除V中的n-m列，以达到V与U的维数一致，反之亦然；

步骤3.4：特征干扰与符号校正：

对初始邻近矩阵添加特征干扰：

K₁=Udiag(δ,δ²,…,δ^m)U^Tδ=rand(0,1)

K₂=Vdiag(ζ,ζ²,…,ζ^m)V^Tζ=rand(0,1)

对特征向量进行符号校正：

$v_i^{\&prime;}=\left\{\begin{array}{cc}{v}_i& \mathrm{if}\left(\right||u_i+v_i||>{\left|\right|u}_i-v_i|\left|\right)\\ -v_i& \mathrm{otherwise}\end{array}\right.;$

步骤3.5：角度一致性约束：假设基准影像中的两条直线段

生成控制点P_o，待配准影像中的两条直线段

生成控制点P_s，如果直线段

间的夹角与直线段

间的夹角之差过大，则点P_o与P_s不可能成为同名点，因而，将角度一致性约束定义为：

$A(P_o,P_s)=\left\{\begin{array}{cc}1& \mathrm{if}\left(\right|\mathrm{\&theta;}(L_i^O,L_j^O)-\mathrm{\&theta;}(L_m^S,L_n^S)|<t)\\ 0& \mathrm{otherwise}\end{array}\right.$

其中，t为控制变形敏感性的阈值，t取值为5°；

步骤3.6：匹配判断：

Z=A*(U*V'^T)

其中V’为V经符号校正后的矩阵；如果Z_ij是它所有行与列的最大元素，则认为点集P中的第i个点与点集Q中的第j个点匹配。

4.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤4中所述的预定阀值m=3。

5.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤5中所述的预定阀值n=5。

6.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤4中所述的自适应调整直线提取参数，所述的直线提取参数包括：线支持区搜索过程中的角度容忍值τ、检测窗口的宽度α及两个检测窗口间的间隔β，其初始值设定为τ=22.5，α=3，β=1.5，在迭代参数调整过程中，若基准影像中检测到更多直线段，则固定基准影像参数不变，将待配准影像参数变量α与β均增加2，以检测到更多直线段；否则，固定待配准影像参数不变，将基准影像参数变量τ增加7.5。

7.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤7中所述的将所述的同名点对映射到原始影像上，计算初始变换参数并得到初变换影像；其具体实现过程为：设控制点在最粗尺度影像上的坐标为(x,y),最粗尺度影像相对于原始影像经过了L次降分，则该点映射到原始影像上的坐标为（2^L-2^L-1）*(x,y)；初始变换模型采用仿射变换模型，利用最粗尺度上的同名点对映射后的相应点对作为同名点对计算变换参数，并对待配准影像进行变换。

8.根据权利要求1所述的异源遥感影像配准方法，其特征在于：步骤9中所述的在对步骤8处理后的影像上分别提取直线特征，根据初始变换参数计算寻找同名直线段，进而获得候选同名点对；其具体实现包括以下子步骤：

步骤9.1：对基准影像与待配准影像采用与低层同样的直线提取方法与直线段预处理方法；

步骤9.2：通过定义两条直线段匹配代价函数判断基准影像与待配准影像上的直线段是否为同名线段，代价函数为：

$C_{\mathrm{r\&theta;}}(s_i,o_j)=\frac{1}{\sqrt{2}}{\{{\left[\frac{d_r(s_i,o_j)}{d_{r\max }}\right]}^2+{\left[\frac{d_{\mathrm{\&theta;}}(s_i,o_j)}{d_{\mathrm{\&theta;}\max }}\right]}^2\}}^{1/2}$

其中，d_r(s_i,o_j)为直线段o_j到直线段s_i的距离，d_θ(s_i,o_j)为直线段o_j与直线段s_i之间的角度差异，d_rmax、d_θmax分别为d_r、d_θ的最大容许值，用于归一化作用，值设定为d_rmax=5像素，d_θmax=3弧度；

如果两条直线段为同名线段，则由上式得到的值趋于0；如果值大于1，则认为两条直线不可能同名；对于待配准影像上的每条直线段s(i),基准影像上与s(i)斜率相近的直线段均用来计算与s(i)的匹配代价，如果基准影像上的所有直线与s(i)的匹配代价均大于1，则认为s(i)为孤立线段，否则，基准影像上与s(i)具有最小匹配代价的线段被认为与s(i)是同名线段；进而将任意两对同名线段的交点作为候选同名点对。

Images