CN108520494A - 基于结构条件互信息的sar图像与可见光图像配准方法 - Google Patents

Abstract

本发明公开了一种基于结构条件互信息的SAR图像与可见光图像配准方法，主要解决现有技术不稳定且配准精度较低的问题，其实现如下：1)输入参考图像与待配准图像；2)分别计算参考图像与待配准图像的经过非局部均值滤波算法处理的相位一致性信息；3)用参考图像与待配准图像的相位一致性信息分别计算参考图像与待配准图像的相位一致性边缘强度图和方向图；4)根据相位一致性边缘强度图和方向图计算参考图像与待配准图像的结构条件互信息；5)在确定搜索空间下，记录结构条件互信息最大时所对应的变换参数；6)利用变换参数对待配准图像进行变换，得到配准结果。本发明配准稳定，且配准精度较高，可用于遥感图像融合与变化检测。

Description

基于结构条件互信息的SAR图像与可见光图像配准方法

技术领域

本发明属于遥感图像处理技术领域，更进一步涉及一种SAR图像与可见光图像的配准方法，可用于遥感图像融合与变化检测。

背景技术

遥感图像配准是对不同时间、不同视场、不同成像模式的两幅或多幅图像进行空间几何变换，使得各个图像在几何上能够匹配对应起来。其主要目的是减少或消除参考图像与待配准图像之间的几何形变，在遥感图像融合等领域有着广泛的应用。

S.Suri等人在其发表的论文“Mutual-information-based registration ofTerraSAR-X and Ikonos imagery in urban areas”(《IEEE TRANSACTIONS ONGEOSCIENCE AND REMOTE SENSING》文章编号：0196-2892(2010)48-0939-10)利用直方图分割的方式提取图像中的部分信息用来计算互信息，在不影响配准精度的情况下，提高了大幅SAR图像与可见光图像配准速度，但是该方法对噪声较大，灰度差异较大的图像无法得到稳定的配准结果。

J.P.W.Pluim等人在其发表的论文“Image registration by maximization ofcombined mutual information and gradient information”(《IEEE TRANSACTIONS ONMEDICAL IMAGING》文章编号：0278-0062(2000)19-0809-05)提出了一种基于梯度互信息的图像配准方法，该方法将图像的梯度信息作为空间信息与互信息相结合，提高了图像配准的稳定性，但是由于梯度对图像噪声的敏感性，使得该方法主要针对质量较好的图像进行配准，对于涉及噪声较大的SAR图像配准则配准结果不稳定且精度较低，影响后续遥感图像融合或变化检测的效果。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于结构条件互信息的SAR图像与可见光图像配准方法，以提高配准的稳定性和精度，便于后续使用。

本发明的技术方案是：首先对参考图像与待配准图像计算相位一致性信息，并通过非局部均值滤波对其处理，利用滤波后的相位一致性信息计算得到参考图像与待配准图像的相位一致性边缘强度图与方向图；然后计算参考图像与待配准图像的结构条件互信息，选择搜索空间中最大条件互信息值所对应的变换参数，对待配准图像进行变换，得到配准结果。其实现步骤包括如下：

(1)输入可见光图像与SAR图像，并将可见光图像作为参考图像V，SAR图像作为待配准图像S；

(2)利用非局部均值滤波算法分别对参考图像V和待配准图像S的相位一致性信息进行滤波，得到滤波后的参考图像相位一致性信息C_V和待配准图像的相位一致性信息C_S；

(3)用参考图像V的相位一致性信息C_V计算参考图像的相位一致性边缘强度图E_V和方向图Φ_V；用待配准图像S的相位一致性信息C_S计算待配准图像的相位一致性边缘强度图E_S和方向图Φ_S；

(4)计算参考图像V与待配准图像S的结构条件互信息：

(4a)分别对待配准图像S以及对应的相位一致性边缘强度图E_S及方向图Φ_S进行空间变换；

(4b)计算参考图像的相位一致性方向图Φ_V与待配准图像的相位一致性方向图Φ_S的角度差ΔΦ_V,S：

(4c)计算参考图像V与待配准图像S在角度差ΔΦ_V,S下的条件互信息CMI(V,S|ΔΦ_V,S)：

CMI(V,S|ΔΦ_V,S)＝H(V,ΔΦ_V,S)+H(S,ΔΦ_V,S)-H(V,S,ΔΦ_V,S)-H(ΔΦ_V,S)其中，H(V,ΔΦ_V,S)表示参考图像V与角度差ΔΦ_V,S的联合信息熵，H(S,ΔΦ_V,S)表示待配准图像S与角度差ΔΦ_V,S的联合信息熵，H(V,S,ΔΦ_V,S)表示参考图像V、待配准图像S与角度差ΔΦ_V,S三者的联合信息熵，H(ΔΦ_V,S)表示角度差ΔΦ_V,S的信息熵；

(4d)利用参考图像V的相位一致性边缘强度图E_V和方向图Φ_V，待配准图像S的相位一致性边缘强度图E_S和方向图Φ_S，计算权重W(V,S)：

其中，x为矢量，表示图像中每一点的坐标，V∩S表示参考图像V与待配准图像S的重合区域；

(4e)计算参考图像V与待配准图像S的结构条件互信息SCMI(V,S)：

SCMI(V,S)＝W(V,S)·CMI(V,S|ΔΦ_V,S)

(5)在确定的搜索空间下，重复执行步骤(4)，记录整个搜索空间下结构条件互信息最大时所对应的变换参数；

(6)利用(5)中的变换参数对待配准图像S进行空间变换，得到配准结果。

本发明与现有技术相比具有如下优点：

第一，由于本发明利用了经过非局部均值滤波处理的相位一致性信息，克服了基于梯度方法的特征提取算子对SAR图像的噪声较为敏感的问题，保证提取的图像结构特征能够反映图像的真实信息。

第二，由于本发明在计算互信息中引入图像相位一致性边缘强度信息与方向信息，有效的结合了图像的空间结构信息，提高了图像配准的稳定性。

第三，由于本发明采用图像相位一致性信息与图像的灰度信息计算得到的结构条件互信息作为相似性度量，减少了陷入局部极值的可能性，提高了图像配准的精度。

附图说明

图1为本发明的实现流程图；

图2为本发明的实验仿真图。

图3为本发明的实验结果图与现有的常用两种方法的配准结果图。

具体实施方式

下面结合附图和实例对本发明的技术方案和效果做进一步描述：

参照图1，本发明的实施步骤如下。

步骤1，设定参考图像和待配准图像。

输入需要进行配准的两幅实测图像，即可见光图像与SAR图像，并将可见光图像作为参考图像V，SAR图像作为待配准图像S。

步骤2，分别对参考图像V和待配准图像S的相位一致性信息进行滤波，得到滤波后的参考图像相位一致性信息C_V和待配准图像的相位一致性信息C_S。

常用的滤波方法有均值滤波、中值滤波、非局部均值滤波等，本实例使用但不限于用非局部均值滤波算法对图像的相位一致性信息进行处理，其实现步骤如下：

2.1)利用非局部均值滤波算法对参考图像V相位一致性中的能量En_V进行滤波处理，即对能量En_V中的任一点x的值En_V(x)，通过该点邻域内的每一点与该点的权重关系按下式重新计算该点的值DEn_V(x)：

其中x为矢量，表示En_V中的任一点，I_V表示以点x为中心的邻域，y为矢量，表示邻域I_V中的任一点，w_V(x,y)表示点y在加权平均时的权重系数，v_x,k表示以点x为中心的图像块Δ_x中的第k个像素点的灰度值，v_y,k表示以y为中心的图像块Δ_y的第k个像素点的灰度值，α_k表示中心对称的高斯核，h表示控制平滑程度的滤波参数；

2.2)遍历参考图像V相位一致性中的能量En_V中的所有点，得到非局部均值滤波后的能量DEn_V；

2.3)按照下式，计算得到滤波后的参考图像任一点的相位一致性信息C_V(x)：

其中，x为矢量，表示图像中的某一像素点的位置，W_V(x)表示参考图像频率传播的权重系数，DEn_V(x)表示参考图像滤波后的复数形式能量，表示参考图像在尺度n下的幅度，T表示噪声阈值，ε表示一个很小的常数,ε＞0；

2.4)遍历参考图像中所有点，得到滤波后的参考图像相位一致性信息C_V；

2.5)利用非局部均值滤波算法对待配准图像S相位一致性中的能量En_S进行滤波处理，即对能量En_S中的任一点x的值En_S(x)，通过该点邻域内的每一点与该点的权重关系按下式重新计算该点的值DEn_S(x)：

其中x为矢量，表示En_S中的任一点，I_S表示以点x为中心的邻域，y为矢量，表示邻域I_S中的任一点，w_S(x,y)表示点y在加权平均时的权重系数，s_x,k表示以点x为中心的图像块Δ_x中的第k个像素点的灰度值，s_y,k表示以y为中心的图像块Δ_y的第k个像素点的灰度值，α_k表示中心对称的高斯核，h表示控制平滑程度的滤波参数；

2.6)遍历待配准图像S相位一致性中的能量En_S中的所有点，得到非局部均值滤波后的能量DEn_S；

2.7)按照下式，计算得到滤波后的待配准图像任一点的相位一致性信息C_S(x)：

其中，x为矢量,表示图像中的某一像素点的位置，W_S(x)表示待配准图像频率传播的权重系数，DEn_S(x)表示待配准图像滤波后的复数形式能量，表示待配准图像在尺度n下的幅度，T表示噪声阈值，ε表示一个很小的常数,ε＞0；

2.8)遍历待配准图像中所有点，得到滤波后的待配准图像相位一致性信息C_S。

步骤3，分别计算参考图像V和待配准图像S的相位一致性边缘强度图E_V、E_S和方向图Φ_V、Φ_S。

3.1)用参考图像V的相位一致性信息C_V计算参考图像的相位一致性边缘强度图E_V和方向图Φ_V：

其中θ表示6个不同的方向，分别是0°、30°、60°、90°、120°、150°，C_V(θ)表示参考图像V在θ方向的滤波后相位一致性信息，a_V、b_V、c_V表示中间变量；

3.2)用待配准图像S的相位一致性信息C_S计算待配准图像的相位一致性边缘强度图E_S和方向图Φ_S：

其中θ表示6个不同的方向，分别是0°、30°、60°、90°、120°、150°，C_S(θ)表示待配准图像S在θ方向的滤波后相位一致性信息，a_S、b_S、c_S表示中间变量。

步骤4，计算参考图像V与待配准图像S的结构条件互信息。

4.1)分别对待配准图像S以及对应的相位一致性边缘强度图E_S及方向图Φ_S进行空间变换；

4.2)计算参考图像的相位一致性方向图Φ_V与待配准图像的相位一致性方向图Φ_S的角度差ΔΦ_V,S：

4.3)计算参考图像V与待配准图像S在角度差ΔΦ_V,S下的条件互信息CMI(V,S|ΔΦ_V,S)：

CMI(V,S|ΔΦ_V,S)＝H(V,ΔΦ_V,S)+H(S,ΔΦ_V,S)-H(V,S,ΔΦ_V,S)-H(ΔΦ_V,S)其中，H(V,ΔΦ_V,S)表示参考图像V与角度差ΔΦ_V,S的联合信息熵，H(S,ΔΦ_V,S)表示待配准图像S与角度差ΔΦ_V,S的联合信息熵，H(V,S,ΔΦ_V,S)表示参考图像V、待配准图像S与角度差ΔΦ_V,S三者的联合信息熵，H(ΔΦ_V,S)表示角度差ΔΦ_V,S的信息熵；

4.4)利用参考图像V的相位一致性边缘强度图E_V和方向图Φ_V，待配准图像S的相位一致性边缘强度图E_S和方向图Φ_S，计算权重W(V,S)：

其中，x为矢量，表示图像中每一点的坐标，V∩S表示参考图像V与待配准图像S的重合区域；

4.5)计算参考图像V与待配准图像S的结构条件互信息SCMI(V,S)：

SCMI(V,S)＝W(V,S)CMI(V,S|ΔΦ_V,S)。

步骤5，在确定的搜索空间下，重复执行步骤(4)，记录整个搜索空间下结构条件互信息最大时所对应的变换参数。

步骤6，利用(5)中的变换参数对待配准图像S进行空间变换，即根据变换参数对待配准图像S进行平移、旋转操作，得到最终的配准结果。

下面结合仿真结果对本发明的效果做进一步的描述。

1.仿真条件：

本发明的实验平台为Intel(R)Pentium(R)CPU G3240@3.10GHz内存4.0GB，运行Windows7的PC机，编程语言为Matlab R2016b。

2.仿真内容与结果分析：

仿真1，用本发明和现有两种方法分别对一组大小为231×314的SAR图像与可见光图像进行配准，结果如图2，其中图2(a)是可见光图像，图2(b)是对应的SAR图像，图2(c)是利用归一化互信息进行配准的结果的轮廓重叠图,图2(d)是利用梯度互信息进行配准的结果的轮廓重叠图,图2(e)是利用本发明进行配准的结果的轮廓重叠图。

从图2中可以看出，利用现有的归一化互信息与梯度互信息的方法配准，其结果存在误差，而本发明的配准误差很小，表明本发明能够得到正确的配准结果。

对图2中三种方法所对应的配准结果计算均方根误差，结果如表1，

表1本发明和归一化互信息与梯度互信息配准结果的均方根误差对比

仿真类型	本发明	归一化互信息	梯度互信息
				均方根误差	0.97	2.56	2.81

从表1可见本发明的配准误差远低于现有的两种方法。

仿真2，用本发明和现有两种方法分别对一组大小为410×419的SAR图像与可见光图像进行配准，结果如图3，其中图3(a)是SAR图像，图3(b)是对应的可见光图像，图3(c)是利用归一化互信息进行配准的结果的轮廓重叠图,图3(d)是利用梯度互信息进行配准的结果的轮廓重叠图,图3(e)是利用本发明进行配准的结果的轮廓重叠图。

从图3中可以看出，梯度互信息的配准结果明显错误，归一化互信息结果存在微小误差，本发明的配准结果正确，说明本发明能够得到正确的配准结果。

对图3中三种方法所对应的配准结果计算均方根误差，结果如表2，表2本发明和归一化互信息与梯度互信息配准结果的均方根误差对比

仿真类型	本发明	归一化互信息	梯度互信息
				均方根误差	0.92	1.13	错误

从表2可见，本发明的配准结果的误差最低。

以上两组实验仿真表明，本发明能够对SAR图像与可见光图像进行稳定配准，且配准结果的精度高。

Claims (5)

1.一种基于结构条件互信息的SAR图像与可见光图像配准方法，包括：

(1)输入可见光图像与SAR图像，并将可见光图像作为参考图像V，SAR图像作为待配准图像S；

(2)利用非局部均值滤波算法分别对参考图像V和待配准图像S的相位一致性信息进行滤波，得到滤波后的参考图像相位一致性信息C_V和待配准图像的相位一致性信息C_S；

(3)用参考图像V的相位一致性信息C_V计算参考图像的相位一致性边缘强度图E_V和方向图Φ_V；用待配准图像S的相位一致性信息C_S计算待配准图像的相位一致性边缘强度图E_S和方向图Φ_S；

(4)计算参考图像V与待配准图像S的结构条件互信息：

(4a)分别对待配准图像S以及对应的相位一致性边缘强度图E_S及方向图Φ_S进行空间变换；

(4b)计算参考图像的相位一致性方向图Φ_V与待配准图像的相位一致性方向图Φ_S的角度差ΔΦ_V,S：

(4c)计算参考图像V与待配准图像S在角度差ΔΦ_V,S下的条件互信息CMI(V,S|ΔΦ_V,S)：

CMI(V,S|ΔΦ_V,S)＝H(V,ΔΦ_V,S)+H(S,ΔΦ_V,S)-H(V,S,ΔΦ_V,S)-H(ΔΦ_V,S)

其中，H(V,ΔΦ_V,S)表示参考图像V与角度差ΔΦ_V,S的联合信息熵，H(S,ΔΦ_V,S)表示待配准图像S与角度差ΔΦ_V,S的联合信息熵，H(V,S,ΔΦ_V,S)表示参考图像V、待配准图像S与角度差ΔΦ_V,S三者的联合信息熵，H(ΔΦ_V,S)表示角度差ΔΦ_V,S的信息熵；

(4d)利用参考图像V的相位一致性边缘强度图E_V和方向图Φ_V，待配准图像S的相位一致性边缘强度图E_S和方向图Φ_S，计算权重W(V,S)：

其中，x为矢量，表示图像中每一点的坐标，V∩S表示参考图像V与待配准图像S的重合区域；

(4e)计算参考图像V与待配准图像S的结构条件互信息SCMI(V,S)：

SCMI(V,S)＝W(V,S)·CMI(V,S|ΔΦ_V,S)

(5)在确定的搜索空间下，重复执行步骤(4)，记录整个搜索空间下结构条件互信息最大时所对应的变换参数；

(6)利用(5)中的变换参数对待配准图像S进行空间变换，得到配准结果。

2.根据权利要求1所述的方法，其特征在于，步骤(2)利用非局部均值滤波算法对参考图像V的相位一致性信息进行滤波，按如下步骤进行：

2a)利用非局部均值滤波算法对参考图像V相位一致性中的能量En_V进行滤波处理，即对能量En_V中的任一点x的值En_V(x)，通过该点邻域内的每一点与该点的权重关系按下式重新计算该点的值DEn_V(x)：

其中x为矢量，表示En_V中的任一点，I_V表示以点x为中心的邻域，y为矢量，表示邻域I_V中的任一点，w_V(x,y)表示点y在加权平均时的权重系数，v_x,k表示以点x为中心的图像块Δ_x中的第k个像素点的灰度值，v_y,k表示以y为中心的图像块Δ_y的第k个像素点的灰度值，α_k表示中心对称的高斯核，h表示控制平滑程度的滤波参数；

2b)遍历参考图像V相位一致性中的能量En_V中的所有点，得到非局部均值滤波后的能量DEn_V；

2c)按照下式，计算得到滤波后的参考图像任一点的相位一致性信息C_V(x)：

其中，x为矢量，表示图像中的某一像素点的位置，W_V(x)表示参考图像频率传播的权重系数，DEn_V(x)表示参考图像滤波后的复数形式能量，表示参考图像在尺度n下的幅度，T表示噪声阈值，ε表示一个很小的常数,ε＞0；

2d)遍历参考图像中所有点，得到滤波后的参考图像相位一致性信息C_V。

3.根据权利要求1所述的方法，其特征在于，步骤(2)利用非局部均值滤波算法对待配准图像S的相位一致性信息进行滤波，其实现如下：

3a)利用非局部均值滤波算法对待配准图像S相位一致性中的能量En_S进行滤波处理，即对能量En_S中的任一点x的值En_S(x)，通过该点邻域内的每一点与该点的权重关系按下式重新计算该点的值DEn_S(x)：

其中x为矢量，表示En_S中的任一点，I_S表示以点x为中心的邻域，y为矢量，表示邻域I_S中的任一点，w_S(x,y)表示点y在加权平均时的权重系数，s_x,k表示以点x为中心的图像块Δ_x中的第k个像素点的灰度值，s_y,k表示以y为中心的图像块Δ_y的第k个像素点的灰度值，α_k表示中心对称的高斯核，h表示控制平滑程度的滤波参数；

3b)遍历待配准图像S相位一致性中的能量En_S中的所有点，得到非局部均值滤波后的能量DEn_S；

3c)按照下式，计算得到滤波后的待配准图像任一点的相位一致性信息C_S(x)：

其中，x为矢量，表示图像中的某一像素点的位置，W_S(x)表示待配准图像频率传播的权重系数，DEn_S(x)表示待配准图像滤波后的复数形式能量，表示待配准图像在尺度n下的幅度，T表示噪声阈值，ε表示一个很小的常数,ε＞0；

3d)遍历待配准图像中所有点，得到滤波后的待配准图像相位一致性信息C_S。

4.根据权利要求1所述的方法，其特征在于，步骤(3)中用参考图像V的相位一致性信息C_V计算参考图像的相位一致性边缘强度图E_V和方向图Φ_V，按如下公式计算：

其中θ表示6个不同的方向，分别是0°、30°、60°、90°、120°、150°，C_V(θ)表示参考图像V在θ方向的滤波后相位一致性信息，a_V、b_V、c_V表示中间变量。

5.根据权利要求1所述的方法，其特征在于，步骤(3)中用待配准图像S的相位一致性信息C_S计算待配准图像的相位一致性边缘强度图E_S和方向图Φ_S，按如下公式计算：

其中θ表示6个不同的方向，分别是0°、30°、60°、90°、120°、150°，C_S(θ)表示待配准图像S在θ方向的滤波后相位一致性信息，a_S、b_S、c_S表示中间变量。