CN101957991A - 一种遥感图像配准方法 - Google Patents

Abstract

本发明涉及图像处理技术领域，具体是指一种遥感图像配准算法，其包括以下操作步骤：(1)从存储介质中读取原始遥感图像数据；(2)对各原始图像进行预处理；(3)进行二值化处理实现图像分割；(4)对二值化图像进行腐蚀和膨胀处理；(5)对腐蚀和膨胀处理后的二值图像运用Roberts交叉梯度算子，提取单值边缘；(6)对边缘图像计算匹配度，寻找最佳匹配点；(7)进行图像间的变换；(8)输出配准信息和配准图像。本发明可应用于多光谱数据融合和遥感图像拼接中。

Description

一种遥感图像配准方法

技术领域

本发明涉及图像处理技术，具体指一种遥感图像配准方法，应用于多光谱数据融合和遥感图像拼接。

背景技术

图像配准技术是根据应用需求发展起来的图像分析技术，在医疗影像、遥感图像等领域有广泛的应用。图像配准算法可概括为两类，即基于区域的配准和基于特征的配准。如文献《Infrared Image Registration and High-Resolution Reconstruction Using Multiple Translationally shifted Aliased Video Frames》所述，传统图像配准算法复杂，计算量大；又如文献《图像配准及其在目标精确定位中的应用》所述，两类算法均对场景变化敏感，使算法不具鲁棒性。

发明内容

本发明的目的是提供一种计算量较小，对缓慢变化场景具有鲁棒性的遥感图像配准算法。

本发明在基于区域的配准算法基础上，对待配准图像进行变换处理，对匹配度计算准则进行优化，降低计算量，提高在变化场景情况下算法的鲁棒性。算法步骤如下：(1)从存储介质中读取原始图像数据；(2)对各原始图像进行预处理；(3)进行二值化处理实现图像分割；(4)对二值化图像进行腐蚀和膨胀处理；(5)对腐蚀和膨胀处理后的二值图像运用Roberts交叉梯度算子，提取单值边缘；(6)对边缘图像计算匹配度，寻找最佳匹配点；(7)对任一图像进行双线性插值和平移变换；(8)输出配准图像。

其中步骤(2)中预处理算法包括基于统计参数的灰度信息调整和基于常量统计量法的非均匀性校正，表述如下：

A.有两幅灰度级不一致的灰度图像f _a 和f _b ，坐标(x，y)处的灰度值分别为f_a(x，y)和f_b(x，y)，

和

分别为f _a 和f _b 的均值，var(f_a)和var(f_b)分别为f _a 和f _b 的方差，f_a′(x，y)和f_b′(x，y)分别为校正后坐标(x，y)处的灰度值，则基于图像统计值的灰度校正方法为

${f_a}^{\′}(x,y)=\frac{f_a(x,y)-\stackrel{\‾}{f_a}}{\mathrm{var}\left(f_a\right)}\·\mathrm{var}\left(f_a\right)+\stackrel{\‾}{f_a}---\left(1\right)$

f_b′(x，y)＝f_b(x，y)         (2)

B.对于具有N个像元的线列探测器，假设只存在增益非均匀性。对每一像元i，计算其响应均值

${\stackrel{\‾}{V}}_i=\underset{j}{\mathrm{\Σ}}{V}_{i,j}---\left(3\right)$

其中，j＝1，2，...，M，M为采样次数，V_i，j为第i像元第j个采样点的响应值。然后获得所有像元的平均响应值

$\stackrel{\‾}{V}=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{V}}_i---\left(4\right)$

则增益校正因子为

$k_i={\stackrel{\‾}{V}}_i/\stackrel{\‾}{V}---\left(5\right)$

校正后第i像元第j个采样点响应为

V_i，j′＝k_i×V_i，j       (6)

对预处理后的结果图像用灰度值矩阵分别表示为I₁(i，j)和I₂(i，j)，(i，j)为坐标。

步骤(3)中的基于全局阈值二值化算法步骤如下：

A.灰度图像宽和高分别为W和H，用其灰度值矩阵I(i，j)表示，(i，j)为矩阵坐标。计算图像的均值作为初始门限T_old，即：

设置叠代次数为Q，门限之差的阈值为T₀，当前为第q次叠代，初始值为q＝0。

B.用T_old分割图像，产生两个集合：G₁＝{I(i，j)：I(i，j)＞T_odd}，G₂＝{I(i，j)：I(i，j)≤T_odd}。

C.对区域集合G₁和G₂，计算其所有像素的平均灰度值，分别记为μ₁和μ₂。

D.计算新的门限值： $T_{\mathrm{new}}=\frac{{\mathrm{\μ}}_1+{\mathrm{\μ}}_2}{2}.$

E.当满足q＜Q，|T_old-T_new|≥T₀时，T_old＝T_new，q＝q+1，重复2)～4)步，直至q≥Q或|T_old-T_new|＜T₀，此时T_new即为二值化处理的全局门限。

F.采用上述全局门限对图像进行二值化处理，得到坐标(i，j)处的二值信息为：

$E(i,j)=\left\{\begin{array}{cc}0& I(i,j)>T_{\mathrm{new}}\\ 1& I(i,j)\≤T_{\mathrm{new}}\end{array}\right.---\left(7\right)$

步骤(4)中的腐蚀和膨胀算子分别定义为：

A被B膨胀：

B对A的腐蚀： $\mathrm{A\ΘB}=\left\{z\right|{\left(B\right)}_z\⊆A\}---\left(9\right)$

其中，A和B是二维整数空间Z²中的集合，

为B的反射，(A)_z＝{c|c＝a+z，a∈A}表示为集合A平移到点z＝(z₁，z₂)，z∈Z²。定义模版B为

$B=\left[\begin{array}{ccc}1& *& *\\ *& 1& *\\ *& *& 1\end{array}\right],$ *表示不关心该元素的值(10)

依次对二值图像进行P次腐蚀和P次膨胀(P＞1，由用户指定)，有效消除闭合区域内的微小空洞，即图像中的细节。

步骤(5)中的Roberts交叉梯度算子分别为：

$G_1=\left[\begin{array}{cc}-1& 0\\ 0& 1\end{array}\right],$ $G_2=\left[\begin{array}{cc}0& -1\\ 1& 0\end{array}\right]---\left(11\right)$

边缘图像信息为：E′(i，j)＝|E(i，j)-E(i-1，j-1)|+|E(i，j-1)-E(i-1，j)|

步骤(6)中的匹配度计算步骤如下：

A.有两幅二值边缘图像E₁′(i，j)和E₂′(i，j)，其中E ₁ ′(i，j)的宽和高分别为W₁、H₁，E ₂ ′(i，j)的宽和高分别为W₂、H₂。对每一图像分别定义行向量R和列向量L为：

$R_1\left(i\right)=\underset{j=1}{\overset{W_1}{\mathrm{\Σ}}}{E_1}^{\′}(i,j),$ $L_1\left(j\right)=\underset{i=1}{\overset{H_1}{\mathrm{\Σ}}}{E_1}^{\′}(i,j),i=\mathrm{1,2},...,H_1,j=\mathrm{1,2},...,W_1---\left(12\right)$

$R_2\left(i\right)=\underset{j=1}{\overset{W_2}{\mathrm{\Σ}}}{E_2}^{\′}(i,j),$ $L_2\left(j\right)=\underset{i=1}{\overset{H_2}{\mathrm{\Σ}}}{E_2}^{\′}(i,j),i=\mathrm{1,2},...,H_2,j=\mathrm{1,2},...,W_2---\left(13\right)$

B.计算图像间行向量间和列向量间的相关系数CorrR和CorrL：

$\mathrm{CorrR}\left(m\right)=\underset{i}{\mathrm{\Σ}}{R}_2\left(i\right)R_1(i+m),$ $\mathrm{CorrL}\left(n\right)=\underset{j}{\mathrm{\Σ}}{L}_2\left(j\right)L_1(j+n)---\left(14\right)$

其中，m，n为实数，通过对R₁和L₁的线性插值可以得到R₁(i+m)和L₁(j+n)

分别表示为m，n的整数部分。

C.分别求CorrR和CorrL中的最大值，计CorrR(M)＝max{CorrR(m)}，CorrL(N)＝max{CorrL(n)}，则(M，N)为最佳匹配点。

步骤(7)中的图像变换包括对图像的双线性插值和平移变换：

A.双线性插值

分别记g，h为M，N的整数部分，即

则有

I₁(i+M，j+N)＝(1-M+g)×(1-N+h)×I₁(i+g，j+h)

+(1-M+g)×(N-h)×I₁(i+g，j+h+1)        (16)

+(M-g)×(1-N+h)×I₁(i+g+1，j+h)

+(M-g)×(N-h)×I₁(i+g+1，j+h+1)

B.平移变换I₂(i，j)＝I₁(i+M，j+N)       (17)

本发明由于采用以上技术方案，具有以下优点：(1)本发明运用腐蚀和膨胀两种形态学方法对二值化图像进行处理，可以消除二值图像中的微小区域或闭合区域内的小空洞，通常这些微小区域或闭合区域内的小空洞是由图像细节或光谱特征引起的，消除其影响可避免误配准；(2)基于二值图像的边缘提取可以保证获得单像素宽度的边界，并且不需要采用复杂的算子进行提取，大大减少了计算量；(3)匹配度(相似度)的计算从二维空间降到一维空间进行，避免了大量计算；(4)降维处理只与目标边界数有关，而与目标位置无关，在处理运动场景时具有一定的鲁棒性；(5)算法可以达到亚像素级的配准精度。

说明书附图

图1是本发明的流程图；

图2是预处理后的待配准的两幅图像；

图3是图像配准结果；

具体实施方式

下面结合附图和实施例，对本发明进行详细描述。

本发明是一种遥感图像配准算法，其具体操作流程如图1所示：

首先，运行图像预处理程序，读入原始图像数据，进行图像预处理(包括灰度校正和非均匀性校正)，处理结果分别如图2中两幅图所示。两幅图像中云层目标场景存在着缓慢的变化。

然后运行图像配准过程程序，对预处理后的图像进行配准。由于待配准的图像一般为两幅以上(包含两幅)，在本发明的实施例中以两幅图像为例。

下一步是进行图像二值化操作。二值化是基于全局阈值进行图像分割得到的。

然后对二值图像进行形态学处理。首先运用模版B对二值图像进行N＝3次腐蚀处理，然后再运用模版B对上述结果图进行N＝3次膨胀操作。

接下来对上述结果图像运用Roberts交叉梯度算子提取单像素边缘。Roberts交叉梯度算子为两个2*2模版，梯度幅值用绝对和表示，相对于其他边缘提取算子，其运算量得到大幅降低。由于在本发明中只需要利用大区域的边界信息，此单边界满足需求，且节省了运算时间。

最后是进行匹配度计算，得到匹配参数。匹配度的计算基于一维向量的相关系数，定义为CorrR(m)＝∑R₂(i)R₁(i+m)，一维向量通过对边缘图像的各行或各列累加得到。相关系数最大的位置即为匹配参数，设为(M，N)。利用匹配参数对图像进行双线性插值及图像平移变换：I₂(i，j)＝I₁(i+M，j+N)，获得匹配结果。匹配参数计算结果为：M＝-8.7，N＝-4.6，图像匹配结果如图3所示。

Claims (3)

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

(1)从存储介质中读取原始图像数据；

(2)原始图像预处理：对原始图像数据进行基于统计参数的灰度信息调整和基于常量统计量法的非均匀性校正；

(3)进行二值化处理实现图像分割；

(4)对二值化图像进行腐蚀和膨胀处理；

(5)对腐蚀和膨胀处理后的二值图像运用Roberts交叉梯度算子，提取单值边缘；

(6)对单值边缘图像计算匹配度，寻找最佳匹配点；

(7)对任一图像进行双线性插值和平移变换；

(8)输出配准图像。

2.如权利要求1所述的一种遥感图像配准方法，其特征在于：步骤4中所述的对二值化图像进行腐蚀和膨胀处理的方法是：

A被B膨胀：

B对A的腐蚀： $\mathrm{A\ΘB}=\left\{z\right|{\left(B\right)}_z\⊆A\}---\left(\text{2}\right)$

其中，A和B是二维整数空间Z²中的集合，

为B的反射，(A)_z＝{c|c＝a+z，a∈A}表示为集合A平移到点z＝(z₁，z₂)，z∈Z²，定义模版B为

$B=\left[\begin{array}{ccc}1& *& *\\ *& 1& *\\ *& *& 1\end{array}\right],$ *表示不关心该元素的值(3)

依次对二值图像进行P次腐蚀和P次膨胀，P＞1，由用户指定，消除闭合区域内的微小空洞，即图像中的细节。

3.如权利要求1所述的一种遥感图像配准方法，其特征在于：步骤6中所述的匹配度计算步骤如下：

对于两幅二值边缘图像E₁′(i，j)和E₂′(i，j)，其中E ₁ ′(i，j)的宽和高分别为W₁、H₁，E ₂ ′(i，j)的宽和高分别为W₂、H₂，对每一图像分别定义行向量R和列向量L为：

$R_1\left(i\right)=\underset{j=1}{\overset{W_1}{\mathrm{\Σ}}}{E_1}^{\′}(i,j),$ $L_1\left(j\right)=\underset{i=1}{\overset{H_1}{\mathrm{\Σ}}}{E_1}^{\′}(i,j),i=\mathrm{1,2},...,H_1,j=\mathrm{1,2},...,W_1---\left(4\right)$

$R_2\left(i\right)=\underset{j=1}{\overset{W_2}{\mathrm{\Σ}}}{E_2}^{\′}(i,j),$ $L_2\left(j\right)=\underset{i=1}{\overset{H_2}{\mathrm{\Σ}}}{E_2}^{\′}(i,j),i=\mathrm{1,2},...,H_2,j=\mathrm{1,2},...,W_2---\left(5\right)$

A.计算图像间行向量间和列向量间的相关系数CorrR和CorrL：

$\mathrm{CorrR}\left(m\right)=\underset{i}{\mathrm{\Σ}}{R}_2\left(i\right)R_1(i+m),$ $\mathrm{CorrL}\left(n\right)=\underset{j}{\mathrm{\Σ}}{L}_2\left(j\right)L_1(j+n)---\left(6\right)$

其中，m，n为实数，通过对R₁和L₁的线性插值可以得到R₁(i+m)和L₁(j+n)

分别表示为m，n的整数部分。

B.分别求CorrR和CorrL中的最大值，计CorrR(M)＝max{CorrR(m)}，

CorrL(N)＝max{CorrL(n)}，则(M，N)为最佳匹配点。

Images