CN102521811A - 基于各向异性扩散和互信息匀质性测度的sar图像降斑方法 - Google Patents

Abstract

本发明公开了一种基于各向异性扩散和互信息匀质性测度的SAR图像相干斑抑制算法，主要解决现有多数相干斑抑制算法在去除匀质区域中的噪声的同时无法有效保留图像中的结构特征的缺陷。其实现过程是：(1)对输入图像中的每个像素点，计算该点的互信息匀质性测度；(2)利用互信息匀质性测度构造各向异性扩散方程的扩散系数；(3)利用新建的扩散系数对均匀区域进行各向同性扩散，对边缘区域进行各向异性扩散；(4)对图像进行一定次数的迭代扩散后，输出最终的降噪图像。与现有的各向异性扩散降斑方法、空域及小波域降斑方法相比，本发明在增强均匀区域的噪声平滑效果的同时，还能有效保留图像中结构和细节特征，各方面的降噪性能提升显著，可用于SAR图像的分割、分类与目标识别。

Description

基于各向异性扩散和互信息匀质性测度的SAR图像降斑方法

技术领域

本发明属于图像处理技术领域，涉及一种SAR图像相干斑抑制方法，可用于对SAR图像分割、分类和目标识别。

背景技术

针对SAR图像的抑斑算法的基本目标是在抑制图像均匀区域斑点噪声的前提下，保留图像的边缘、纹理和强反射点目标。针对SAR图像中的斑点噪声为乘性的特点，人们发展了许多的降斑算法。这类算法大致可以分为空域滤波算法和变换域滤波算法。典型的空域滤波算法有Lee滤波器，Kuan滤波器和Frost滤波器；典型的变换域滤波算法包括小波域、Contourlet域、curvelet域滤波算法。然而，不管是空域滤波算法，还是小波域滤波算法，它们在抑制相干斑时是各向同性的，即对边缘附近各个方向平滑程度相同，因此必然会破坏图像中的边缘和细节信息。基于非线性偏微分方程的图像滤波方法能够在降噪的同时增强图像的边缘。YuYongjian等提出了一种各种异性扩散相干斑抑制算法——SRAD(speckle reduction anisotropicdiffusion)。SRAD能够在均匀区域进行各向同性扩散，而在边缘附近进行各向异性扩散，因此可以在去除均匀区域相干斑噪声的同时增强和保护图像的边缘，但经过对SRAD的深入研究，发现其在均匀区域的降斑效果却不佳。主要是因为SRAD将变差系数作为匀质性测度来辨别图像中均匀区域、异质区域以及边缘，而变差系数极易受到相干斑噪声的影响，因此在均匀区域进行扩散时，有可能因为受到噪声影响而停止扩散，这样就使得均匀区域内噪声平滑不够充分。虽然可以通过增加迭代次数来解决，但随着迭代次数的增加图像中的边缘和纹理信息会最终而变得模糊。

发明内容

本发明的目的在于克服SRAD的不足，提出一种基于互信息匀质性测度的各向异性扩散降斑方法，以提升各向异性扩散降斑方法在均匀区域的平滑效果并保留图像中的边缘和细节信息。

为实现上述目的，本发明的实现步骤如下：

步骤1：根据图像的格式(强度与幅度)求解相干斑变差系数τ₀，当输入图像为强度格式时，

当为幅度图像时，

步骤2：计算图像中每一个像素点的互信息匀质性测度τ_i，j：

步骤3：利用互信息匀质性测度构建各向异性扩散偏微分方程的扩散系数，即对像素点(i，j)，其对应的扩散系数表达式为

步骤4：基于步骤(3)得到的扩散系数，利用各向异性扩散偏微分方程

对图像中每一个像素点的灰度值进行迭代更新；

步骤5：更新迭代次数n：n＝n+1，并判断n是否超过设定的最大迭代次数，若是输出处理后的降噪图像；否则，返回步骤2。

本发明与现有技术相比的优点在于：本发明相对于经典的空域滤波器如Lee滤波器、Frost滤波器及变换域滤波方法，其优势在于不仅可以有效去除SAR图像中均匀区域的相干斑噪声，还能够保留图像中的重要纹理信息并增强图像的边缘；相对于以往的各向异性扩散抑斑方法，本发明的优势在于可以增强均匀区域内对噪声的平滑效果。本发明的特点是由于使用了更加精确且对噪声不敏感的互信息匀质性测度作为SAR图像匀质区域、异质区域及边缘区域的检测算子，因此可以在匀质区域内进行充分的各向同性扩散从而尽可能地平滑掉噪声，在边缘与纹理区域进行各向异性扩散从而尽可能地保留原始图像信息。

附图说明

图1是本发明对SAR图像进行降噪的流程图；

图2是实验用的真实SAR图像以及本发明对原始图像的降噪图；

图3是实验用的真实SAR图像以及本发明对原始图像的降噪图.

具体实施方式

如图1所示，本发明的具体实施步骤如下：

步骤1根据图像的格式(强度与幅度)求解相干斑变差系数τ₀，当输入图像为强度格式时，

当为幅度图像时，

步骤2计算图像中每一个像素点的互信息匀质性测度τ_i，j。依据Bruno Aiazzi等人于2004年发表在IEEE Transaction on Geoscience and Remote Sensing上的文章《Information-TheoreticHeterogeneity Measurement for SAR Imagery》，τ_i，j是按照如下步骤进行计算的：

2a)在一个3×3的滑动窗口内计算出像素点(i，j)的局部均值

与标准差

${\widehat{\mathrm{\μ}}}_g=\frac{1}{9}\underset{k=1}{\overset{9}{\mathrm{\Σ}}}{I}_k$ ${\widehat{\mathrm{\σ}}}_g=\sqrt{\frac{1}{8}\underset{k=1}{\overset{9}{\mathrm{\Σ}}}{(I_k-{\widehat{\mathrm{\μ}}}_g)}^2};$

2b)计算出所有像素点的

和

后，以

和

为坐标轴，作出该图像所有像素点对应的

和的二维分布图；

2c)对和

的二维分布图进行256×256的网格划分，得到二维直方图h(i，j)，并对h(i，j)进行规范化；

2d)对二维的离散直方图h(i，j)进行高斯平滑得到

与

的离散联合概率密度 $p(i,j)=p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g);$

2e)根据

计算

与

的二维条件概率密度函数

$p\left({\widehat{\mathrm{\σ}}}_g\right|{\widehat{\mathrm{\μ}}}_g)=\frac{p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}{p\left({\widehat{\mathrm{\μ}}}_g\right)}=\frac{p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}{\underset{{\widehat{\mathrm{\μ}}}_g}{\mathrm{\Σ}}p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}$

2f)根据

计算像素点(i，j)的互信息匀质性测度τ(i，j)：

$\mathrm{\τ}(i,j)=-\mathrm{\α}\log \left\{p\left(\widehat{\mathrm{\σ}}(i,j)\right|\widehat{\mathrm{\μ}}(i,j))\right\}$

其中，α为一常数，对强度图像

对幅度图像

步骤3 利用互信息匀质性测度构建各向异性扩散偏微分方程的扩散系数，即对像素点(i，j)，其对应的扩散系数表达式为

步骤4基于步骤(3)得到的扩散系数，利用各向异性扩散偏微分方程

对图像中每一个像素点的灰度值进行迭代更新，按照以下步骤进行的：

4a)首先对各向异性扩散偏微分方程

进行离散化，得到如下离散表达式：

$I_{i,j}^{n+1}=I_{i,j}^n+\frac{\mathrm{\Δt}}{\left|{\mathrm{\η}}_{i,j}\right|}\underset{p\∈{\mathrm{\η}}_{i,j}}{\mathrm{\Σ}}C\left(\mathrm{\τ}\right){\▿I}_{i,j,p}^n$

其中η_s表示像素点s的邻域，|η_s|表示邻域内点的个数，p表示η_s中点，Δt为时间步长。

4b)利用前向差分代替步骤3a)中的离散偏微分方程中的梯度算子，得到各个像素点灰度值的更新方程：

$\begin{array}{c}{I}_s^{n+1}=I_s^n+\frac{\mathrm{\Δt}}{\left|{\mathrm{\η}}_s\right|}[c_{i+1,j}(I_{i+1,j}-I_{i,j})+c_{i,j}(I_{i-1,j}-I_{i,j})\\ +c_{i,j+1}(I_{i,j+1}-I_{i,j})+c_{i,j}(I_{i,j-1}-I_{i,j})]\end{array}$

本发明的有效性可以通过以下实验得到验证：

仿真实验的环境为：硬件环境为AMD双核CUP 2.91GHz，2G RAM内存；软件环境为MATLAB R2007a。

实验内容包括：

分别采用本发明给出的方法、SRAD、Lee滤波器对真实SAR图像进行降噪。实验结果如图2，图3所示，其中图2(a)、3(a)是原始SAR图像，图2(b)、3(b)是Lee滤波器的降斑结果，图2(c)、3(c)是SRAD的降斑结果，图2(d)、3(d)是本发明的降斑结果。为了定量地比较各种降斑方法的性能，本发明选取了均值MEAN，标准差STD，等效视数ENL作为性能评价指标。表1和表2分别给出了各滤波器对图2(a)和图3(a)的降噪性能指标参数计算结果。从降噪后的图像以及各表中可以看出，本发明给出的降噪方法相比Lee滤波器和SRAD，更能够有效平滑同质区域的相干斑噪声并降低图像的标准差；另一方面，虽然在均值保持方面要稍次于SRAD，但滤波后的均值还是非常接近于原始图像的均值。

Claims (3)

1.利用各向异性扩散和互信息匀质性测度的SAR图像降斑方法，其特征在于如下步骤：

步骤1：根据图像的格式(强度与幅度)求解相干斑变差系数τ₀，当输入图像为强度格式时，

当为幅度图像时，

步骤2：计算图像中每一个像素点的互信息匀质性测度τ_i，j；

步骤3：利用互信息匀质性测度构建各向异性扩散偏微分方程的扩散系数，即对像素点(i，j)，其对应的扩散系数表达式为

步骤4：基于步骤(3)得到的扩散系数，利用各向异性扩散偏微分方程

对图像中每一个像素点的灰度值进行迭代更新；

步骤5：更新迭代次数n：n＝n+1，并判断n是否超过设定的最大迭代次数，若是，输出处理后的降噪图像；否则，返回步骤(2)。

2.根据权利要求1所述的方法，其中步骤(2)所述的计算图像的互信息匀质性测度τ(i，j)，依据Bruno Aiazzi等人于2004年发表在IEEE Transaction on Geoscience and Remote Sensing上的文章《Information-Theoretic Heterogeneity Measurement for SAR Imagery》，是按照如下步骤进行计算的：

2a)在一个3×3的滑动窗口内计算出像素点(i，j)的局部均值与标准差

${\widehat{\mathrm{\μ}}}_g=\frac{1}{9}\underset{k=1}{\overset{9}{\mathrm{\Σ}}}{I}_k$ ${\widehat{\mathrm{\σ}}}_g=\sqrt{\frac{1}{8}\underset{k=1}{\overset{9}{\mathrm{\Σ}}}{(I_k-{\widehat{\mathrm{\μ}}}_g)}^2};$

2b)计算出所有像素点的和

后，以

和

为坐标轴，绘出该图像所有像素点对应的

和

的二维分布图；

2c)对

和

的二维分布图进行256×256的网格划分，得到二维直方图h(i，j)，并对h(i，j)进行规范化；

2d)对二维的离散直方图h(i，j)进行高斯平滑得到与

的离散联合概率密度 $p(i,j)=p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g);$

2e)根据

计算

与

的二维条件概率密度函数

$p\left({\widehat{\mathrm{\σ}}}_g\right|{\widehat{\mathrm{\μ}}}_g)=\frac{p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}{p\left({\widehat{\mathrm{\μ}}}_g\right)}=\frac{p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}{\underset{{\widehat{\mathrm{\μ}}}_g}{\mathrm{\Σ}}p({\widehat{\mathrm{\μ}}}_g,{\widehat{\mathrm{\σ}}}_g)}$

2f)根据

计算像素点(i，j)的互信息匀质性测度τ(i，j)：

$\mathrm{\τ}(i,j)=-\mathrm{\α}\log \left\{p\left(\widehat{\mathrm{\σ}}(i,j)\right|\widehat{\mathrm{\μ}}(i,j))\right\}$

其中，α为一常数，对强度图像

对幅度图像

3.根据权利要求1所述的方法，其中步骤(4)所述的利用各向异性扩散偏微分方程对图像中每一个像素点进行灰度值更新，是按照以下步骤进行的：

3a)首先对各向异性扩散偏微分方程

进行离散化，得到如下离散表达式：

$I_{i,j}^{n+1}=I_{i,j}^n+\frac{\mathrm{\Δt}}{\left|{\mathrm{\η}}_{i,j}\right|}\underset{p\∈{\mathrm{\η}}_{i,j}}{\mathrm{\Σ}}C\left(\mathrm{\τ}\right){\▿I}_{i,j,p}^n$

其中η_s表示像素点s的邻域，|η_s|表示邻域内点的个数，p表示η_s中点，Δt为时间步长。

3b)利用前向差分代替步骤3a)中的离散偏微分方程中的梯度算子，得到各个像素点灰度值的更新方程：

$\begin{array}{c}{I}_s^{n+1}=I_s^n+\frac{\mathrm{\Δt}}{\left|{\mathrm{\η}}_s\right|}[c_{i+1,j}(I_{i+1,j}-I_{i,j})+c_{i,j}(I_{i-1,j}-I_{i,j})\\ +c_{i,j+1}(I_{i,j+1}-I_{i,j})+c_{i,j}(I_{i,j-1}-I_{i,j})]\end{array}$

Images