CN101950413A - 基于非下采样Contourlet域MRF模型的SAR图像降斑方法 - Google Patents

Abstract

本发明公开了一种基于非下采样Contourlet域MRF模型的SAR图像降斑方法，主要解决传统多尺度几何分析应用于SAR图像降斑易产生匀质区域不够光滑，细节信息不够完整，点目标不够清晰的问题。其降斑过程为：对待降斑SAR图像作对数运算，得到图像数据；对图像数据进行非下采样Contourlet变换，得到变换系数C_j，l；估计图像数据噪声在非下采样Contourlet域的标准差σ_j，l；根据C_j，l和σ_j，l对掩码X_j，l进行初始估计并迭代更新，得到X_j，l的最终估计；利用X_j，l的最终估计，对C_j，l进行收缩，得到估计系数Y_j，l；对Y_j，l作指数运算，并对其进行非下采样Contourlet逆变换，得到最终的降斑图像。本发明具有降斑结果同质区域光滑，纹理信息和边缘信息的完整性保持较好，以及点目标清晰可见，可用于SAR图像目标识别。

Description

基于非下采样Contourlet域MRF模型的SAR图像降斑方法

技术领域

本发明属于图像处理领域，特别是一种涉及SAR图像降斑的方法，可应用于目标识别。

背景技术

合成孔径雷达SAR是一种高分辨率雷达体制，具有全天候、全天时、穿透能力强、可侧视成像等优点，具有广泛的应用。SAR可应用于农业，军事，导航，地理监视等诸多领域。它与其它遥感成像系统、光学成像系统是不同的。SAR成像过程中，由于相干照射而产生的衰落效应的影响，使图像内原本具有相同后向散射系数的匀质区域表现出颗粒状噪声，这种噪声称之为相干斑。相干斑的存在增加了图像解译的复杂性，大大降低了图像分割、目标分类以及其他信息提取的有效性。

对相干斑的抑制技术，即降斑，大体分为成像前的多视平滑预处理和成像后的滤波处理两类。早期的SAR成像处理中，大多采用多视处理技术，但随着SAR图像应用的不断拓展，对其空间分辨率的要求不断提高，多视处理技术已不能满足要求。成像后的滤波技术又可分为：空域滤波技术和频域滤波技术。

基于空域滤波技术的降斑方法，从20世纪80年代开始，出现了许多以空域滤波技术为基础的高分辨率SAR图像降斑算法，并获得了广泛的研究和应用。空域滤波技术都是假定噪声模型，然后在图像上取一个滑动窗，以窗内所有像素作为滤波器的输入值，估计无相干斑噪声的图像，它是基于局部统计特性进行滤波处理的。空域滤波方法有Enhanced Lee滤波，Kuan滤波，Frost滤波，以及Gamma Map滤波等。这些空域的滤波方法，难以保持图像的细节特征，其滤波性能的好坏，很大程度上依赖于所选滤波窗口的大小。

基于频域滤波技术的降斑方法，从20世纪90年代以后，是以频域处理技术为基础的小波方法给信号处理领域带来了崭新的思想，并被有效地用于SAR图像相干斑抑制中。小波变换具有时频局部化特性，为SAR图像的降斑提供了另一条思路。它可以从不同的分辨率空间来描述图像的局部特征，使得信号和噪声在小波变换域中表现出不同的特征，从而较易区分信号和噪声。基于小波变换的图像去噪方法中，以Donoho提出的硬阈值和软阈值去噪方法为起源，得到了广泛的应用。基于Donoho的阈值思想，逐渐又出现了类似的系数收缩方法，经典的收缩方法有VisuShrink，SureShrink，OracleShrink，BayesShrink等。但是，这些方法主要是针对加性噪声模型的。针对SAR图像中相干斑乘性噪声的算法，通常都是通过对数运算将其转换为加性噪声，然后利用上述的收缩方法对系数进行处理，再通过指数运算重构图像。上述方法都是在假设图像信息独立性的条件下进行的一系列处理，忽视了图像内部相关性，于是对于一致性较好的区域，噪声抑制效果仍然不够理想，目标的保持也不够清晰。为了进一步挖掘图像内部相关性，人们将马尔可夫随机场模型与小波变换相结合，建立了小波马尔可夫模型，并用于SAR图像抑斑。Dusan Gleich等人将MRF引入二进小波变换域，对小波域中的无噪图像用广义高斯马尔可夫随机场来建模，并用于图像降斑，参见Wavelet-Based Despeckling of SAR Images Using Gauss-Markov Random Fields.IEEE Trans.on Geoscience and Remote Sensing，2007，45(12)：4127-4143。经过比较分析，小波变换的方法对边缘有着较好的保留效果，但是由于小波变换有限的方向性，它仅能捕获图像的水平、垂直及对角三个方向的信息，而对于其他方向却无能为力，这就导致处理后的图像中细节信息丢失严重，纹理信息和边缘信息不完整。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种基于非下采样Contourlet域MRF模型的SAR图像降斑方法，利用非下采样Contourlet变换完备的方向性信息提高边缘分辨精度和区域一致性，保证降斑后图像纹理信息和边缘信息的完整性以及目标的清晰性，提高SAR图像降斑的质量。

实现本发明目的的技术方案是：使用高效的多尺度几何分析工具非下采样Contourlet变换，充分挖掘图像中的方向信息，构建出具有严密结构的非下采样Contourlet域MRF模型。其具体实现步骤包括如下：

(1)对待降斑SAR图像作对数运算，得到对数运算后的图像数据；

(2)对图像数据进行非下采样Contourlet变换，得到变换系数C_j，l，

其中，j表示尺度，范围为1～J，J为非下采样Contourlet变换的分解层数；l表示方向子带，范围为1～L，L表示每个尺度中的方向数；C_j，l＝{c₁，...，c_h，...，c_H}，c_h表示C_j，l里的第h个元素，1≤h≤H，H表示待降斑SAR图像像素点的总个数；

(3)利用Monte Carlo方法，估计图像数据噪声在非下采样Contourlet域的标准差σ_j，l；

(4)利用迭代条件模型ICM算法，根据变换系数C_j，l和标准差σ_j，l对掩码X_j，l进行初始估计并迭代更新，得到掩码的最终估计：

(4a)估计MRF模型掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值，

其中，当c_h≥σ_j，l时，x_h＝1；当c_h＜σ_j，l时，x_h＝-1，1≤h≤H；

(4b)计算重要性度量M_j，l＝{m₁，...，m_h，...，m_H}，

其中，m_j＝abs(c_h)，abs表示取绝对值，1≤h≤H；

(4c)计算似然比Lhood_j，l＝{ξ₁，...，ξ_h，...，ξ_H}和先验比Prior_j，l＝{μ₁，...，μ_h，...，μ_H}，以及似然比和先验比的乘积R_j，l＝{r₁，...，r_h，...，r_H}，

其中：

${\mathrm{\ξ}}_h=\frac{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=1)}{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=-1)}$

$=\left\{\begin{array}{cc}\exp (-\mathrm{\&alpha;}),& m_h/{\mathrm{\&sigma;}}_{j,l}<(1-\mathrm{\&delta;});\\ \exp (\mathrm{\&alpha;}\&CenterDot;(m_h/{\mathrm{\&sigma;}}_{j,l}-1)/\mathrm{\&delta;}),& (1-\mathrm{\&delta;})\&le;m_h/{\mathrm{\&sigma;}}_{j,l}\&le;(1+\mathrm{\&delta;});\\ \exp \left(\mathrm{\&alpha;}\right),& m_h/{\mathrm{\&sigma;}}_{j,l}>(1+\mathrm{\&delta;}).\end{array}\right.$

${\mathrm{\&mu;}}_h=\frac{p(x_h=1|x_{\&PartialD;h})}{p(x_h=-1|x_{\&PartialD;h})}=\exp (\mathrm{\&gamma;}\&CenterDot;\underset{i\&Element;\&PartialD;h}{\mathrm{\&Sigma;}}{x}_i)$

r_H＝ξ_h·μ_h

其中，(m_h|x_h＝1)表示在X_j，l条件下，当x_h＝1时出现M_j，l中的m_h的概率，

(m_h|x_h＝-1)表示在X_j，l条件下，当x_h＝-1时出现M_j，l中的m_h的概率，exp表示指数运算，α＝5，δ＝0.9，

表示以h为中心的方形区域，i为该方形区域中的位置，

表示在

区域中掩码存在的条件下出现x_h＝1的概率，

表示在

区域中掩码存在的条件下出现x_h＝-1的概率，γ为邻域影响因子，γ＝0.2；

(4d)对掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值进行更新，如果似然比和先验比的乘积r_h大于1，则把掩码的初始值中的x_h更新为1，如果似然比和先验比的乘积r_h小于等于1，则把掩码的初始值中的x_h更新为-1，1≤h≤H；

(4e)统计掩码的更新数目，判断掩码更新的终止条件，如果掩码的更新数目大于H/2，则认为不满足收敛条件，转至步骤4c)；如果掩码的更新数目小于等于H/2，则认为满足收敛条件，其对x_h的更新值为最终估计，转至步骤5)，H表示待降斑SAR图像像素点的总个数；

(5)利用掩码的最终估计，对变换系数C_j，l＝{c₁，...，c_h，...，c_H}进行收缩，得到估计系数Y_j，l＝{y₁，...，y_h，...y_H}，其中

y_h示Y_j，l里的第h个元素， 和

是根据掩码的最终估计计算得到的，计算公式与步骤4c)的公式相同，1≤h≤H，H表示待降斑SAR图像像素点的总个数；

(6)对估计系数Y_j，l作指数运算，并对指数运算后的结果系数进行非下采样Contourlet逆变换，得到降斑图像。

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

1、本发明由于使用多尺度几何分析工具非下采样Contourlet变换的冗余特性，能够对SAR图像纹理区域的方向性细节信息进行提取，保证了图像纹理信息和边缘信息的完整性以及目标的清晰性；

2、本发明将非下采样Contourlet变换的多方向性和MRF模型相结合，尺度内对应关系严密，对同质区域有效地进行了降斑；

3、仿真结果表明，本发明方法较现有的Gamma Map滤波和Enhanced Lee滤波方法能更有效的进行图像降斑。

附图说明

图1是本发明的流程示意图；

图2是本发明使用的一幅2视幅度待降斑SAR图像；

图3是本发明对图2的仿真结果对比图；

图4是本发明使用的一幅4视强度待降斑SAR图像；

图5是本发明对图4的仿真结果对比图。

具体实施方式

参照图1，本发明的具体实现步骤如下：

步骤一、输入待降斑SAR图像，并对它作对数运算，得到对数运算后的图像数据。

这一步骤的目标是将待降斑SAR图像中的乘性噪声转换为加性噪声，输入待降斑SAR图像如图2和图4所示。

步骤二、对图像数据进行非下采样Contourlet变换，得到变换系数C_j，l。

该非下采样Contourlet变换采用“maxflat”非下采样塔形分解和“dmaxflat7”非下采样方向滤波器组分解，通过对图像数据进行上述非下采样Contourlet变换，得到变换系数C_j，l，这里采用三层非下采样Contourlet变换，从第一个尺度到第三个尺度，每个尺度的方向数依次为8，4，8；其中，j表示尺度，范围为1～J，J为非下采样Contourlet变换的分解层数；l表示方向子带，范围为1～L，L表示每个尺度中的方向数；C_j，l＝{c₁，...，c_h，...，c_H}，c_h表示C_jl里的第h个元素，1≤h≤H，H表示待降斑SAR图像像素点的总个数，实验中采用的图像大小为256×256，像素点的总个数H＝65536。

步骤三、利用Monte Carlo方法，估计图像数据噪声在非下采样Contourlet域的标准差σ_j，l。

(3a)对图像数据进行一层正交小波变换，提取正交小波变换后的对角方向子带系数D，用Donoho的鲁棒中值估计方法计算出图像数据的噪声标准差

median表示取中值，abs表示取绝对值；

(3b)产生一幅大小与待降斑SAR图像相同，均值为0，标准差为σ的高斯白噪声图像；

(3c)对高斯白噪声图像进行非下采样Contourlet变换，计算噪声变换系数N_j，l的标准差，该非下采样Contourlet变换方法采用“maxflat”非下采样塔形分解和“dmaxflat7”非下采样方向滤波器组分解，通过对高斯白噪声图像进行非下采样Contourlet变换，得到噪声变换系数N_j，l，再计算噪声变换系数N_j，l在非下采样Contourlet域的标准差：σ_j，l＝median(abs(N_j，l-median(N_j，l)))/0.6745。

步骤四、利用迭代条件模型ICM算法，根据变换系数C_j，l和标准差σ_j，l对掩码X_j，l进行初始估计并迭代更新，得到掩码的最终估计。

(4a)估计MRF模型掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值，

其中，当c_h≥σ_j，l时，x_h＝1；当c_h＜σ_j，l时，x_h＝-1，x_h表示X_j，l里的第h个元素，1≤h≤H；

(4b)计算重要性度量M_j，l＝{m₁，...，m_h，...，m_H}，

实验中，取变换系数的绝对值作为重要性度量，即M_j，l＝abs(C_j，l)，其中m_h＝abs(c_h)，m_h表示M_j，l里的第h个元素，1≤h≤H，此外，重要性度量还可以用Lipschitz指数来表示；

(4c)计算似然比Lhood_j，l＝{ξ₁，...，ξ_h，...，ξ_H}和先验比Prior_j，l＝{μ₁，...，μ_h，...，μ_H}，以及似然比和先验比的乘积R_j，l＝{r₁，...，r_h，...，r_H}，

其中：

${\mathrm{\&xi;}}_h=\frac{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=1)}{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=-1)}$

$=\left\{\begin{array}{cc}\exp (-\mathrm{\&alpha;}),& m_h/{\mathrm{\&sigma;}}_{j,l}<(1-\mathrm{\&delta;});\\ \exp (\mathrm{\&alpha;}\&CenterDot;(m_h/{\mathrm{\&sigma;}}_{j,l}-1)/\mathrm{\&delta;}),& (1-\mathrm{\&delta;})\&le;m_h/{\mathrm{\&sigma;}}_{j,l}\&le;(1+\mathrm{\&delta;});\\ \exp \left(\mathrm{\&alpha;}\right),& m_h/{\mathrm{\&sigma;}}_{j,l}>(1+\mathrm{\&delta;}).\end{array}\right.$

${\mathrm{\&mu;}}_h=\frac{p(x_h=1|x_{\&PartialD;h})}{p(x_h=-1|x_{\&PartialD;h})}=\exp (\mathrm{\&gamma;}\&CenterDot;\underset{i\&Element;\&PartialD;h}{\mathrm{\&Sigma;}}{x}_i)$

r_h＝ξ_h·μ_h

其中，(m_h|x_h＝1)表示在X_j，l条件下当x_h＝1时，出现M_j，l中的m_h的概率，

(m_h|x_h＝-1)表示在X_j，l条件下当x_h＝-1时，出现M_j，l中的m_h的概率，exp表示指数运算，α＝5，δ＝0.9，

表示以h为中心的方形区域，i为该方形区域中的位置，

表示在

区域中掩码存在的条件下出现x_h＝1的概率，

表示在

区域中掩码存在的条件下出现x_h＝-1的概率，γ为邻域影响因子，γ＝0.2；

(4d)对掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值进行更新，如果似然比和先验比的乘积r_h大于1，则把掩码的初始值中的x_h更新为1，如果似然比和先验比的乘积r_h小于等于1，则把掩码的初始值中的x_h更新为-1，1≤h≤H；

(4e)统计掩码的更新数目，判断掩码更新的终止条件，如果掩码的更新数目大于H/2，则认为不满足收敛条件，转至步骤4c)；如果掩码的更新数目小于等于H/2，则认为满足收敛条件，其对x_h的更新值为最终估计，转至步骤5)，H表示待降斑SAR图像像素点的总个数。

步骤五、利用掩码的最终估计，对变换系数C_j，l＝{c₁，...，c_h，...，c_H}进行收缩，得到估计系数Y_j，l＝{y₁，...，y_h，...y_H}。

根据掩码的最终估计，重新计算似然比

和先验比

，以及似然比和先验比的乘积

，计算公式与步骤4c)的公式相同，再对变换系数进行收缩，收缩公式为

y_h表示Y_j，l里的第h个元素，1≤h≤H，H表示待降斑SAR图像像素点的总个数。

步骤六、对估计系数Y_j，l作指数运算，并对指数运算后的结果系数进行非下采样Contourlet逆变换，通过非下采样Contourlet逆变换重构出降斑图像。

本发明的效果可通过以下仿真进一步说明：

1仿真内容：应用本发明方法和现有的Gamma Map滤波方法，Enhanced Lee滤波方法分别对两幅SAR图像进行降斑实验，并从同质区域的内部一致性、纹理信息和边缘信息的完整性以及目标的清晰性评价这些方法的各自性能。

2仿真实验结果

A)2视幅度SAR图像的实验结果

用本发明方法以及Gamma Map滤波方法，Enhanced Lee滤波方法对图2的2视幅度SAR图像进行降斑，其效果比较如图3所示。其中图3(a)为Gamma Map滤波方法对图2进行降斑得到的结果；图3(b)为Enhanced Lee滤波方法对图2进行降斑得到的结果；图3(c)为本发明方法对图2进行降斑得到的结果。由图3可见，本发明中，边缘清晰细节信息完整，点目标保持较好。相比之下，Gamma Map滤波方法在消除平滑区域相干斑上比较有效，但是图像出现比较严重的模糊现象，边缘信息丢失严重。Enhanced Lee滤波方法处理后的图像中，斑点噪声没有得到很好的抑制。

B)4视强度SAR图像的实验结果

用本发明方法以及Gamma Map滤波方法，Enhanced Lee滤波方法对图4的4视强度SAR图像进行降斑，其效果比较如图5所示。其中图5(a)为Gamma Map滤波方法对图4进行降斑得到的结果；图5(b)为Enhanced Lee滤波方法对图4进行降斑得到的结果；图5(c)为本发明方法对图4进行降斑得到的结果。从图5可见，本发明更好地抑制了斑点噪声，匀质区域更加光滑，边缘信息保存良好；相比之下，Gamma Map滤波方法处理后，图像边缘信息丢失严重。Enhanced Lee滤波方法处理后，斑点噪声的抑制不够明显。

Claims (2)

1.一种基于非下采样Contourlet域MRF模型的SAR图像降斑方法，包括如下步骤：

1)对待降斑SAR图像作对数运算，得到对数运算后的图像数据；

2)对图像数据进行非下采样Contourlet变换，得到变换系数C_j，l，

其中，j表示尺度，范围为1～J，J为非下采样Contourlet变换的分解层数；l表示方向子带，范围为1～L，L表示每个尺度中的方向数；C_j，l＝{c₁，...，c_h，...，c_H}，c_h表示C_j，l里的第h个元素，1≤h≤H，H表示待降斑SAR图像像素点的总个数；

3)利用Monte Carlo方法，估计图像数据噪声在非下采样Contourlet域的标准差σ_j，l；

4)利用迭代条件模型ICM算法，根据变换系数C_j，l和标准差σ_j，l对掩码X_j，l进行初始估计并迭代更新，得到掩码的最终估计：

4a)估计MRF模型掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值，

其中，当c_h≥σ_j，l时，x_h＝1；当c_h＜σ_j，l时，x_h＝-1，1≤h≤H；

4b)计算重要性度量M_j，l＝{m₁，...，m_h，...，m_H}，

其中，m_h＝abs(c_h)，abs表示取绝对值，1≤h≤H；

4c)计算似然比Lhood_j，l＝{ξ₁，...，ξ_h，...，ξ_H}和先验比Prior_j，l＝{μ₁，...，μ_h，...，μ_H}，以及似然比和先验比的乘积R_j，l＝{r₁，...，r_h，...，r_H}，

其中：

${\mathrm{\&xi;}}_h=\frac{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=1)}{p_{M_{j,l}|X_{j,l}}\left(m_h\right|x_h=-1)}$

$=\left\{\begin{array}{cc}\exp (-\mathrm{\&alpha;}),& m_h/{\mathrm{\&sigma;}}_{j,l}<(1-\mathrm{\&delta;});\\ \exp (\mathrm{\&alpha;}\&CenterDot;(m_h/{\mathrm{\&sigma;}}_{j,l}-1)/\mathrm{\&delta;}),& (1-\mathrm{\&delta;})\&le;m_h/{\mathrm{\&sigma;}}_{j,l}\&le;(1+\mathrm{\&delta;});\\ \exp \left(\mathrm{\&alpha;}\right),& m_h/{\mathrm{\&sigma;}}_{j,l}>(1+\mathrm{\&delta;}).\end{array}\right.$

${\mathrm{\&mu;}}_h=\frac{p(x_h=1|x_{\&PartialD;h})}{p(x_h=-1|x_{\&PartialD;h})}=\exp (\mathrm{\&gamma;}\&CenterDot;\underset{i\&Element;\&PartialD;h}{\mathrm{\&Sigma;}}{x}_i)$

r_h＝ζ_h·μ_h

其中，

(m_h|x_h＝1)表示在X_j，l条件下，当x_h＝1时出现M_j，l中的m_h的概率，

(m_h|x_h＝-1)表示在X_j，l条件下，当x_h＝-1时出现M_j，l中的m_h的概率lexp表示指数运算，α＝5，δ＝0.9，

表示以h为中心的方形区域，i为该方形区域中的位置，

表示在

区域中掩码存在的条件下出现x_h＝1的概率，

表示在

区域中掩码存在的条件下出现x_h＝-1的概率，γ为邻域影响因子，γ＝0.2；

4d)对掩码X_j，l＝{x₁，...，x_h，...，x_H}的初始值进行更新，如果似然比和先验比的乘积r_h大于1，则把掩码的初始值中的x_h更新为1，如果似然比和先验比的乘积r_h小于等于1，则把掩码的初始值中的x_h更新为-1，1≤h≤H；

4e)统计掩码的更新数目，判断掩码更新的终止条件，如果掩码的更新数目大于H/2，则认为不满足收敛条件，转至步骤4c)；如果掩码的更新数目小于等于H/2，则认为满足收敛条件，其对x_h的更新值为最终估计，转至步骤5)，H表示待降斑SAR图像像素点的总个数；

5)利用掩码的最终估计，对变换系数C_j，l＝{c₁，...，c_h，...，c_H}进行收缩，得到估计系数Y_j，l＝{y₁，...，y_h，...y_H}，其中

y_h表示Y_i，l里的第h个元素，

和是根据掩码的最终估计计算得到的，计算公式与步骤4c)的公式相同，1≤h≤H，H表示待降斑SAR图像像素点的总个数；

6)对估计系数Y_j，l作指数运算，并对指数运算后的结果系数进行非下采样Contourlet逆变换，得到降斑图像。

2.根据权利要求1所述的图像降斑方法，其中步骤3)按如下步骤进行：

3a)对图像数据进行正交小波变换，用Donoho的鲁棒中值估计方法估计出图像数据的噪声标准差σ；

3b)产生一幅大小与待降斑SAR图像相同，均值为0，标准差为σ的高斯白噪声图像；

3c)对高斯白噪声图像进行非下采样Contourlet变换，得到噪声变换系数N_j，l；计算噪声变换系数N_j，l在非下采样Contourlet域的标准差：σ_j，l＝median(abs(N_j，l-median(N_j，l)))/0.6745，其中median表示取中值。

Images