CN101980286B - 结合双树复小波与双变量模型的sar图像降斑方法 - Google Patents

Abstract

本发明公开了一种结合双树复小波与双变量模型的的SAR图像降斑方法，主要解决现有SAR图像降斑方法不能很好的抑制斑点噪声和丢失部分边缘和细节信息的问题。其实现过程为：对原始SAR图像进行双树复小波分解，得到在各个尺度上分解系数的实部和虚部；利用非对数加性噪声模型，求解得到噪声系数的方差；根据利用局部邻域窗口，求得复小波系数实部边缘方差和虚部边缘方差；利用最大后验估计求得阈值收缩函数，对双树复小波分解的系数进行阈值收缩；对缩减后的系数运行双树复小波重构，得到降斑后图像。本发明具有有效去除SAR图像斑点噪声和边缘保持性好的优点，可用于对边缘信息丰富、细节丰富的SAR图像，特别是含有机场、跑道、道路的SAR图像降斑。

Description

结合双树复小波与双变量模型的SAR图像降斑方法

技术领域

本发明属于图像处理技术领域，涉及图像噪声抑制，具体地说是一种复小波域的SAR图像降斑方法，可用于合成孔径雷达图像的斑点噪声的抑制。

背景技术

合成孔径雷达(SAR)是一种高分辨率成像雷达。它具有全天候、多极化、多视角、多俯角数据获取能力及对一些地物的穿透能力，不仅广泛地被应用在军事上，在农业、气象、地形地貌、灾情监控等民用上也有大量的应用。但由于SAR发射的是相干波，这些相干波经过与地物的相干作用，特别是地物的后向散射作用，使目标回波信号产生了衰减，这种衰减表现在图像上就是相干斑点噪声。因此如何抑制图像中的相干斑噪声，提高图像的解读能力以及获得更多的信息成为一个重要的问题。

SAR图像降斑的首要目标是在滤除斑点噪声的同时，尽可能的保持图像的细节信息。斑点噪声是一种复杂的乘性噪声模型的信号。对于斑点噪声这种特殊的性质，在过去的近二十年中，人们提出了很多经典的SAR降斑方法，如Lee滤波，增强Lee滤波，Kuan滤波等等。这些方法是用一个已定义好的滤波器窗来估计局域斑点噪声的方差，并进行滤波处理，其结果通常过分的平滑了边缘细节信息，这些方法均在一定程度上收到了较好的效果。1995年，美国学者Donoho把小波理论引入到图像去噪中，提出了小波软阈值方法。小波软阈值方法是一种非线性的算法，依然存在破坏图像细节信息的问题，对图像的辐射特性保持也不好。

离散小波变换由于具有缺少平移不变性与较差的方向选择性的缺点，最近英国学者Kingsbury等人提出了双树复小波变换，在图像去噪中的应用已初步显示出其明显的优势：与离散小波变换相比，双树复小波变换由于其具有近似平移不变性和更多的方向选择性，有效的解决离散小波变换中出现的振铃效应。但是这种双树复小波变换降斑方法并没有充分考虑图像的几何特征以及SAR图像在复小波域的统计特性和系数之间的局部相关性，降斑后的SAR图像平滑区域的斑点噪声滤除不充分，同时图像的细节和边缘信息部分丢失。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种结合双树复小波与双变量模型的SAR图像降斑方法，以充分滤除降斑后的SAR图像平滑区域的斑点噪声，完整的保留图像的细节和边缘信息。

实现本发明目的的技术思路是结合双树复小波变换的平移不变性和多方向选择性，利用复小波系数实虚部双变量模型对双树复小波变换分解高频系数进行自适应萎缩，得到滤除斑点噪声、保留细节信息的SAR图像。其具体实现步骤包括如下：

(1)对原始SAR图像I进行双树复小波分解，得到在尺度j上的分解复系数y_j，其实部和虚部分别为y_r，j，y_i，j；

(2)利用非对数加性噪声模型，求解各尺度上的噪声方差

(3)利用局部邻域窗口，分别求解在尺度j上复小波系数的实部边缘标准差σ_r，j和虚部边缘标准差σ_i，j；

(4)分别对在尺度j上的复系数的实部和虚部进行阈值收缩，求得估计的无噪小波系数

(5)对缩减后的系数

运行双树复小波重构，得到降斑后的图像

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

1)本发明由于采用非对数加性噪声模型，利用这种模型可以避免在乘性模型转化为加性模型时，因取对数操作带来的对原图像的辐射特性保持的不足。因此能够更加充分的保留原图像的辐射特性。

2)本发明由于利用复小波系数实虚部双变量模型，充分考虑到了SAR图像本身的方向特性与局部特性，保留了更加丰富的图像边缘和细节信息，充分滤除了SAR图像平滑区域的斑点噪声。

3)仿真结果表明，本发明方法较其它几种现有的经典的SAR图像降斑方法，在平滑区域的平滑效果以及边缘保持能力方面都有显著的提高。

附图说明

图1是本发明的流程框图；

图2是本发明与现有两种降斑方法应用于X波段幅度SAR图像的降斑仿真结果对比图；

图3是本发明与现有两种降斑方法应用于Ku波段强度SAR图像的降斑仿真结果对比图。

具体实施方式

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

步骤1，对输入SAR原始图像进行双树复小波分解。

输入原始SAR图像记为I，该SAR图像本身就是被斑点噪声污染的图像，因此不需要像研究自然图像去噪一样，给原图加上一个随机噪声或者某种特性的噪声，可直接对该图像进行降斑处理，对输入的原始SAR图像I进行双树复小波分解，得到一个与一个低频图像和J个尺度，每个尺度有6个高频图像，在尺度j上高频图像的复小波系数记为y_j

y_j＝y_r，j+i·y_i，j                 (1)

其中y_r，j为复小波系数实部，y_i，j为复小波系数虚部。

步骤2，利用非对数加性噪声模型，求解各尺度上的噪声方差

(2a)利用非对数加性噪声模型，将输入的原始SAR图像表示为：

I＝RX＝X+(R-1)X＝X+N               (2)

其中R代表相干斑，其均值为1，方差为

X代表地物的真实后向散射强度，N是将要滤除的加性噪声；从输入的原始SAR图像中的每一点取方形窗I(k)，窗口大小为k×k，计算每个局部方形窗的噪声方差：

${\mathrm{\σ}}_{N\left(k\right)}^2=({\mathrm{\σ}}_{I\left(k\right)}^2+m_{I\left(k\right)}^2)\·{\mathrm{\σ}}_R^2/({\mathrm{\σ}}_R^2+1)---\left(3\right)$

其中m_I(k)，

分别表示局部窗口I(k)的均值和方差，

为相干斑R方差，k取值为3，5，7；

(2b)对每个局部方形窗的噪声方差

求平均值，得到各尺度上噪声方差为：

${\mathrm{\σ}}_n^2=\mathrm{mean}\left({\mathrm{\σ}}_{N\left(k\right)}^2\right).---\left(4\right)$

步骤3，利用局部邻域窗口，分别求解在尺度j上复小波系数的实部边缘标准差σ_r，j和虚部边缘标准差σ_i，j。

分别对第j尺度上复小波系数的实部y_r，j和虚部y_i，j中的每一个点取方形窗N(l)，窗口大小为l×l，利用下式分别计算在尺度j上复小波系数的实部边缘标准差σ_r，j和虚部边缘标准差σ_i，j：

${\mathrm{\σ}}_{r,j}=\sqrt{\max (1/M\·{\mathrm{\Σ}}_{m\∈N\left(l\right)}{y}_{r,j}^2\left(m\right)-{\mathrm{\σ}}_n^2,0)}---\left(5\right)$

${\mathrm{\σ}}_{i,j}=\sqrt{\max (1/M\·{\mathrm{\Σ}}_{m\∈N\left(l\right)}{y}_{i,j}^2\left(m\right)-{\mathrm{\σ}}_n^2,0)}---\left(6\right)$

其中M是方形窗N(l)中系数的个数，l取值为3，5，7。

步骤4，分别对第j尺度上的复系数的实部和虚部进行阈值收缩，求得估计的无噪小波系数

(4a)设任一尺度的复小波系数的实部和虚部近似满足如下分布：

$p_y\left(y\right)=\frac{3}{2{\mathrm{\π\σ}}^2}\·\exp (-\frac{\sqrt{3}}{\mathrm{\σ}}\sqrt{y_{r,j}^2+y_{i,j}^2})---\left(7\right)$

其中σ为复小波系数边缘标准差，y_r，j与y_i，j分别为在第j尺度上复小波系数的实部和虚部；

(4b)求解在尺度j上无噪信号的最大后验估计MAP估计，得到第j尺度上复小波系数实部的收缩函数

和虚部的收缩函数

分别为：

${\hat{w}}_{r,j}=\frac{\mathrm{soft}(\sqrt{y_{r,j}^2+y_{r,j}^2}-\frac{\sqrt{3}{\mathrm{\σ}}_n^2}{\mathrm{\σ}})}{\sqrt{y_{r,j}^2+y_{r,j}^2}}\·y_{r,j}---\left(8\right)$

${\hat{w}}_{i,j}=\frac{\mathrm{soft}(\sqrt{y_{i,j}^2+y_{i,j}^2}-\frac{\sqrt{3}{\mathrm{\σ}}_n^2}{\mathrm{\σ}})}{\sqrt{y_{i,j}^2+y_{i,j}^2}}\·y_{i,j}---\left(9\right)$

其中soft(g)定义为：

$\mathrm{soft}\left(g\right)=\left\{\begin{array}{cc}0& g<0\\ g& g\&GreaterEqual;0\end{array}\right.;---\left(10\right)$

(4c)求解在尺度j上复小波实部阈值T_r，j和虚部阈值T_i，j中的较大值T_j

T_j＝max(T_r，j，T_i，j)                            (11)

其中

σ_r，j，σ_i，j分别为第j尺度上复小波系数的实部边缘标准差和虚部边缘标准差；

(4d)在第j尺度上利用下式进行阈值收缩，计算得到缩减后的第j尺度上复小波系数：

${\hat{w}}_j=\mathrm{soft}\left(\right|y_j|-T_j)\&CenterDot;e^{i\&CenterDot;\mathrm{\&theta;}\left(y_j\right)}---\left(12\right)$

其中θ(y_j)表示y_j方向的弧度值。

步骤5，将缩减后的J个尺度的6个方向的高频系数与低频图像输入双树复小波重构滤波器组进行重构，最终得到重构的降斑图像

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

1.实验内容

实验1，将现有的Gamma-MAP降斑方法、经典双变量模型的降斑方法和本发明降斑方法应用于X波段的幅度SAR图像进行降斑。

实验2，将Gamma-MAP降斑方法、经典双变量模型的降斑方法和本发明降斑方法应用于Ku波段的强度SAR图像进行降斑。

本发明利用等效视数ENL，图像均值M和标准差V作为评价SAR图像降斑性能的客观标准。ENL越高，说明平滑区域的降斑性能越好，降斑后的均值越接近原图像均值越好，降斑后标准差越低，说明平滑效果越好。对比结果进一步说明了本发明在降噪方面的优越性。

2.实验结果

实验1的结果如图2所示，其中图2(a)为原SAR图像，图2(b)为Gamma-MAP方法降斑结果，图2(c)为经典双变量模型的方法降斑图，图2(d)为本发明降斑结果图。图2(a)中所示矩形区域1、2为表1中计算ENL所需的同质区域，表1列出了仿真内容1的仿真结果所得的均值、方差和等效视数对比结果。

表1：不同降斑方法客观指标评价：等效视数ENL、均值M、标准差V

从表1中可以看出，本发明获得了比现有其他方法更高的等效视数，降斑后的图像均值非常接近原图像均值，降斑后的标准差最低。因此本发明取得了比其他降斑方法更优的平滑降斑效果。

实验2的结果如图3所示，其中图3(a)为原SAR图像，图3(b)为经典Gamma-MAP方法降斑结果，图3(c)为经典双变量模型的方法降斑图，图3(d)为本发明降斑结果图。图3(a)中所示矩形区域3、4为表2中计算ENL所需的同质区域。表2列出了仿真内容(2)的仿真结果所得的均值、方差和等效视数对比结果。

表2：不同降斑方法客观指标评价：等效视数ENL、均值M、标准差V

从表2中可以看出，本发明获得了比其他方法更高的等效视数，降斑后的图像均值非常接近原图像均值，降斑后的标准差最低。因此本发明取得了比其他降斑方法更优的平滑降斑效果。

Claims (1)

1.一种结合双树复小波与双变量模型的的SAR图像降斑方法，包括如下步骤：

(1)对原始SAR图像I进行双树复小波分解，得到在尺度j上的分解复系数y_j，其实部和虚部分别为y_r，j，y_i，j；

(2)利用非对数加性噪声模型，求解各尺度上的噪声方差

(2a)根据非对数加性噪声模型，对原图像中I的每一点取方形窗I(k)，窗口大小为k×k，求每个局部方形窗的噪声方差：

${\mathrm{\&sigma;}}_{N\left(k\right)}^2=({\mathrm{\&sigma;}}_{I\left(k\right)}^2+m_{I\left(k\right)}^2)\&CenterDot;{\mathrm{\&sigma;}}_R^2/({\mathrm{\&sigma;}}_R^2+1)$

其中m_I(k)，

分别表示局部窗口I(k)的均值和方差，

为原图I的相干斑噪声方差，k取3，5，7；

(2b)对每个局部方形窗的噪声方差

求平均值，得到各尺度上的噪声方差：

${\mathrm{\&sigma;}}_n^2=\mathrm{mean}\left({\mathrm{\&sigma;}}_{N\left(k\right)}^2\right);$

(3)利用局部邻域窗口，分别求解在尺度j上复小波系数的实部边缘标准差σ_r，j和虚部边缘标准差σ_i，j：

分别对第j尺度上复小波系数的实部y_r，j和虚部y_i，j中的每一个点取方形窗N(l)，窗口大小为l×l，分别求解在尺度j上复小波系数的实部边缘标准差σ_r，j和虚部边缘标准差σ_i，j

${\mathrm{\&sigma;}}_{r,j}=\sqrt{\max (1/M\&CenterDot;{\mathrm{\&Sigma;}}_{m\&Element;N\left(l\right)}{y}_{r,j}^2\left(m\right)-{\mathrm{\&sigma;}}_n^2,0)}$

${\mathrm{\&sigma;}}_{i,j}=\sqrt{\max (1/M\&CenterDot;{\mathrm{\&Sigma;}}_{m\&Element;N\left(l\right)}{y}_{i,j}^2\left(m\right)-{\mathrm{\&sigma;}}_n^2,0)}$

其中M是方形窗N(l)中系数的个数，l取3，5，7；

(4)分别对在尺度j上的复系数的实部和虚部进行阈值收缩，求得估计的无噪小波系数

(4a)设任一尺度的复小波系数的实部和虚部近似满足如下分布：

$p_y\left(y\right)=\frac{3}{2\mathrm{\&pi;}{\mathrm{\&sigma;}}^2}\&CenterDot;\exp (-\frac{\sqrt{3}}{\mathrm{\&sigma;}}\sqrt{y_{r,j}^2+y_{i,j}^2})$

其中σ为复小波系数边缘标准差，y_r，j与y_i，j分别为在第j尺度上复小波系数的实部和虚部；

(4b)求解尺度j上无噪信号的最大后验估计，得到第j尺度上复小波系数实部收缩函数和虚部的收缩函数

分别为：

${\hat{w}}_{r,j}=\frac{\mathrm{soft}(\sqrt{y_{r,j}^2+y_{r,j}^2}-\frac{\sqrt{3}{\mathrm{\&sigma;}}_n^2}{\mathrm{\&sigma;}})}{\sqrt{y_{r,j}^2+y_{r,j}^2}}\&CenterDot;y_{r,j}$

${\hat{w}}_{i,j}=\frac{\mathrm{soft}(\sqrt{y_{i,j}^2+y_{i,j}^2}-\frac{\sqrt{3}{\mathrm{\&sigma;}}_n^2}{\mathrm{\&sigma;}})}{\sqrt{y_{i,j}^2+y_{i,j}^2}}\&CenterDot;y_{i,j}$

其中soft(g)定义为：

$\mathrm{soft}\left(g\right)=\left\{\begin{array}{cc}0& g<0\\ g& g\&GreaterEqual;0\end{array}\right.$

(4c)求解在尺度j上复小波实部阈值T_r，j和虚部阈值T_i，j中的较大值T_j

T_j＝max(T_r，j，T_i，j)

其中

σ_i，j分别为第j尺度上复小波系数的实部边缘标准差和虚部边缘标准差；

(4d)在第j尺度上利用下式进行阈值收缩，计算得到缩减后的复小波系数：

${\hat{w}}_j=\mathrm{soft}\left(\right|y_j|-T_j)\&CenterDot;e^{i\&CenterDot;\mathrm{\&theta;}\left(y_j\right)}$

其中θ(y_j)表示y_j方向的弧度值；

(5)对缩减后的系数

运行双树复小波重构，得到降斑后图像

Images