CN104637037A - 一种基于非本地分类稀疏表示的sar图像降噪方法 - Google Patents

Abstract

本发明公开了一种基于非本地分类稀疏表示的SAR图像降噪方法。属于数字图像处理技术领域。它是一种基于分类系数估计与稀疏表示相结合的SAR图像降噪方法。首先按照图像的同质区和异质区进行分类，而后分别采用非本地均值方法和3D变换域硬阈值收缩方法实现对同质区和异质区的降噪，并由降噪结果获得非本地系数。为进一步提高降噪性能，利用所估计出的非本地系数对稀疏表示模型中的系数进行约束，以更好的逼近真实图像的系数。最后,利用最终估计的系数重构出降噪后的图像。本发明对图像的同质区和异质区采用不同方法进行处理，不仅可有效的实现同质区抑斑，同时能保留异质区中存在的细节，可用于SAR图像降噪。

Description

一种基于非本地分类稀疏表示的SAR图像降噪方法

技术领域

本发明属于数字图像处理技术领域，它特别涉及基于图像非本地相似性与稀疏表示相结合的降噪方法，用于SAR图像降噪处理。

背景技术

合成孔径雷达SAR被广泛的应用于多个领域，并在国民经济发展与国防安全保卫等领域中发挥着重要的作用。但SAR图像的相干成像机理使其成像后的图像中存在严重的相干斑噪声，对SAR图像的后处理，如图像分割、边缘检测和特征提取等带来了一定影响，不利于SAR图像的解译。因此，对相干斑噪声的有效抑制便成为SAR图像处理的一个重要环节。

信号的稀疏表示能有效提取信号的本质特征，基于稀疏表示的图像降噪因其良好的降噪效果而得到国内外学者广泛关注。其中一些方法通过将稀疏表示与其他降噪方法相结合来提高降噪的性能，如将稀疏降噪与均值滤波，全微分相结合等。尤其是将稀疏降噪与非本地冗余特性相结合的降噪方法取得了较好的结果，如BM3D及其改进方法在降噪的整体效果方面达到了领先水平。

发明内容

本发明的目在于针对现有SAR图像降噪中存在的不足，提出一种基于非本地相似性与稀疏表示相结合的SAR降噪方法。该方法充分考虑了SAR图像中不同区域内的特征与噪声分布的情况，首先对图像分类，然后不同的方法对不同的图像块进行降噪以估计出每一图像块对应的非本地系数，最后利用非本地系数调节本地稀疏系数使最终的稀疏解更接近于真实图像的系数。因此，该方法可以提高稀疏降噪的性能。包括以下步骤：

步骤一、图像内容的分类

对于图像中的每个像素点，计算以该点为中心的图像块的方差，然后将其与预先设置好的阈值相比较，若大于阈值，则该中心像素点被标记为异质区，若小于阈值，则标记为同质区。

步骤二、分类估计非本地系数

首先，根据目标块分类后的属性设定相应的块大小与搜索区域。然后，由式(2)计算出目标块与搜索范围内候选图像块的相似度：

$\mathrm{\ω}(y_i,y_j)=\exp [-\underset{k\∈\mathrm{\Ω}}{\mathrm{\Σ}}((2L-1)\log (\frac{z_i\left(k\right)}{z_j\left(k\right)}+\frac{z_j\left(k\right)}{z_i\left(k\right)})+\mathrm{\ζL}\frac{{|{\hat{x}}_i\left(k\right)-{\hat{x}}_j\left(k\right)|}^2}{{\hat{x}}_i\left(k\right){\hat{x}}_j\left(k\right)})]$     式(2)

其中y_i和y_j分别表示目标块和候选块，ω(y_i,y_j)为y_i和y_j之间的相似度，k为块内像素，L为SAR图像视数，ζ为权重调节参数，表示上一次降噪结果，初始降噪时为0。在求得候选图像块与目标图像块之间的相似度后，从中选择相似度最大的K个图像块作为目标块的相似块。对于异质区目标块，将其本身及相似块联合构造为3D矩阵，并对其进行小波变换与硬阈值操作，再经过小波逆变换得到目标块的降噪结果；对于同质区的目标块，将其对应的相似块进行加权平均以实现目标块的降噪处理。

为获得稀疏表示中更适合该图像的字典，分别对同质区和异质区的目标块进行Kmeans聚类，并对每一类图像块使用主成分分析法训练相应的自适应学习字典D，最后将降噪后的目标块在该字典上投影，获得同质区目标块的非本地稀疏系数η_i和异质区稀疏系数γ_i。

步骤三、结合非本地系数与稀疏表示的SAR图像降噪

首先将目标块在字典D上投影得到本地稀疏系数α_i，即：

α_i＝D^Ty_i    式(3)

然后利用式(4)和式(5)求解式(1)中的权重参数λ₁和λ₂或λ₁和λ₃：

${\mathrm{\λ}}_1=\frac{{\left|\right|y_i\left|\right|}_2}{{\mathrm{L\σ}}_i}$     式(4)

${\mathrm{\λ}}_2=\frac{{\left|\right|y_i\left|\right|}_2}{{\mathrm{L\δ}}_i},{\mathrm{\λ}}_3=\frac{{\left|\right|y_i\left|\right|}_2}{{\mathrm{L\χ}}_i}$     式(5)

其中σ_i表示α_i的标准差，δ_i和χ_i分别表示α_i-η_i和α_i-γ_i的标准差，||y||₂表示图像块的2范数。在D、α_i、λ₁、λ₂、λ₃、η_i、γ_i确定后，利用迭代收缩算法求解式(1),得到最终的稀疏解最后，使用式(6)将重构得到降噪后的图像：

$\hat{x}=D\widehat{\mathrm{\α}}={\left(\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{R}_i^T{R}_i\right)}^{-1}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left(R_i^T{D}_i{\widehat{\mathrm{\α}}}_i\right)$     式(6)

将降噪后图像与初始含噪图像进行凸组合后的结果作为新的含噪图像，然后重复步骤二和步骤三的降噪处理，当相邻两次降噪所获图像间的平均绝对偏差MAD小于设置的阈值时结束循环。

本发明的创新点是利用图像的非本地冗余信息提高图像稀疏降噪的性能；利用不同的方法实现对不同类型图像块的降噪以获得非本地稀疏系数；并利用估计出的非本地系数对稀疏表示模型中的系数进行约束，以进一步的提高降噪性能，并将该方法用于SAR图像降噪。

本发明的有益效果：稀疏降噪的过程中结合图像的非本地冗余特性，一方面可以训练出更适合图像表示的字典，另一方面可调节本地稀疏系数使其更接近于干净图像的稀疏系数，从本质上提高稀疏降噪的性能；在求解非本地系数时，根据图像块内容的不同而选择不同的参数和方法以提高非本地系数的估计精度，从而更好的调节本地稀疏系数；从稀疏表示降噪原理出发，获得乘性噪声下的稀疏降噪模型及参数，以更好的保持SAR图像的固有特征，有利于SAR图像的后续处理。

本发明主要采用仿真实验的方法进行验证，所有步骤、结论都在MATLAB8.0上验证正确。

附图说明

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

图2是本发明仿真中使用的真实含噪的田园SAR图像；

其中白色矩形区域为选择的同质区；

图3是PPB方法对图2的降噪结果图；

图4是BM3D方法对图2的降噪结果图；

图5是本发明方法对图2的降噪结果图。

具体实施方式

参照图1，本发明基于非本地分类稀疏表示的SAR图像降噪方法，具体步骤包括如下：

步骤一、图像内容的分类

首先，对SAR图像进行简单的高斯空域滤波，方差为1.2；然后计算8×8大小图像块的方差c＝σ/μ。设置阈值其中σ_n为相干斑的标准差，n为块内像素点个数。若c大于T，则该中心像素点属于异质区，并设k_i＝1；反之，属于同质区，并设k_i＝0。

步骤二、分类估计非本地系数

对属于同质区的目标块，设置块大小为9×9，搜索范围的大小为15×15，当目标块属于异质区时设置块的大小为3×3，搜索范围的大小为39×39；首先通过式(2)计算出搜索范围内K个候选块的之间相似度，其中ζ＝1，K＝16，并将这些候选块设置为初始相似块的集合，对相似度大小进行排序，然后逐个计算搜索范围内剩余候选块的相似度，若小于当前相似块集合中候选块的相似度则丢弃，否则将该块插入相似块集合并更新相似块集合。

对异质区相似块组成的3D矩阵使用式(7)估计目标块的非本地系数：

${\hat{X}}_i={\mathrm{\Γ}}_{3D}^{-1}\left(\mathrm{\τ}\left({\mathrm{\Γ}}_{3D}\left(Y_i\right)\right)\right)$     式(7)

其中Y_i表示相似块组成的3D矩阵，Γ_3D和分别表示小波变换和小波逆变换，τ表示阈值处理，且阈值的大小由局部噪声方差决定。当目标块属于同质区时，采用对相似块进行平均加权的方法估计非本地系数，如式(8)所示：

${\hat{x}}_i=\frac{\underset{j\∈C_i}{\mathrm{\Σ}}\mathrm{\ω}(y_i,y_j)\·y_j}{\underset{j\∈C_i}{\mathrm{\Σ}}\mathrm{\ω}(y_i,y_j)}$     式(8)

其中表示图像块y_i的非本地估计。

利用Kmeans聚类算法，分别将同质区和异质区的目标块分为15类，针对每类图像块，利用主成分分析法训练自适应字典D，将降噪后的图像块在字典D上投影得到非本地稀疏系数η_i或γ_i。

步骤三、结合估计系数与稀疏表示的SAR图像降噪

首先利用式(3)求解目标块在字典D上投影的本地稀疏系数α_i，采用式(4)和式(5)求解式(1)中的权重参数λ₁，λ₂和λ₃，并结合步骤二中估计出的非本地稀疏系数η_i或γ_i，并由式(9)计算出式(1)的解

${\widehat{\mathrm{\&alpha;}}}_i=\left\{\begin{array}{cc}{S}_{{\mathrm{\&lambda;}}_1,{\mathrm{\&lambda;}}_i,{\mathrm{\&beta;}}_i}\left({\mathrm{\&alpha;}}_i\right)& {\mathrm{\&beta;}}_i\&GreaterEqual;0\\ {-S}_{{\mathrm{\&lambda;}}_1,{\mathrm{\&lambda;}}_i,{\mathrm{\&beta;}}_i}\left({\mathrm{\&alpha;}}_i\right)& {\mathrm{\&beta;}}_i<0\end{array}\right.$     式(9)

其中λ_i表示k_i确定后，式(1)右端非本地约束项的权重参数，β_i为对应的非本地系数，的表达式如下：

$S_{{\mathrm{\&lambda;}}_1,{\mathrm{\&lambda;}}_i,{\mathrm{\&beta;}}_i}\left({\mathrm{\&alpha;}}_i\right)=\left\{\begin{array}{cc}{\mathrm{\&alpha;}}_i+{\mathrm{\&lambda;}}_1+{\mathrm{\&lambda;}}_i& {\mathrm{\&alpha;}}_i<{-\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i\\ 0& {-\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i\&le;{\mathrm{\&alpha;}}_i\&le;{\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i\\ {\mathrm{\&alpha;}}_i-{\mathrm{\&lambda;}}_1+{\mathrm{\&lambda;}}_i& {\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i<{\mathrm{\&alpha;}}_i<{\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i+{\mathrm{\&beta;}}_i\\ {\mathrm{\&beta;}}_i& {\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i+{\mathrm{\&beta;}}_i\&le;{\mathrm{\&alpha;}}_i\&le;{\mathrm{\&lambda;}}_1+{\mathrm{\&lambda;}}_i+{\mathrm{\&beta;}}_i\\ {\mathrm{\&alpha;}}_i-{\mathrm{\&lambda;}}_1-{\mathrm{\&lambda;}}_i& {\mathrm{\&alpha;}}_i>{\mathrm{\&lambda;}}_1+{\mathrm{\&lambda;}}_i+{\mathrm{\&beta;}}_i\end{array}\right.$     式(10)

在得到稀疏解后，利用式(6)重构稀疏解得到降噪后图像将与含噪图像y进行凸组合：

$y^{\&prime;}=\mathrm{\&rho;y}+(1-\mathrm{\&rho;})\hat{x}$     式(11)

其中ρ为0.02，对结合后图像y′重复执行步骤二、步骤三，并在每次循环之前通过式(12)计算图像的平均绝对偏差：

$\mathrm{MAD}=\frac{1}{N}{\left|\right|{\hat{x}}^k-{\hat{x}}^{k-1}\left|\right|}_1$     式(12)

其中N为图像像素值个数，和分别为第k-1和k次的降噪结果；当平均绝对偏差MAD小于设置门限ε时停止循环，得到最终降噪图像，其中门限ε设置为2.5e-04。

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

一、实验条件和内容

实验条件：实验使用的输入图像为图2，像素大小为256×256。实验中各降噪算法都使用MATLAB语言编程实现。

实验内容：在上述实验条件下，使用PPB方法、SAR-BM3D方法与本发明方法进行对比。降噪能力的客观评价结果用同质区等效视数ENL衡量。

实验1：用本发明方法和现有的PPB方法及SAR-BM3D方法分别对图2进行降噪，其中PPB方法的图像块大小为7×7，搜索窗大小为21×21，迭代次数为25，降噪结果如图3所示；SAR-BM3D方法的图像块大小为8×8，搜索窗大小为39×39，降噪结果如图4所示；本发明中同质区的图像块大小为9×9，搜索范围的大小为15×15，异质区的块大小为3×3，搜索范围的大小为39×39，降噪结果如图5所示。通过图3，图5之间的对比可以看出，本发明方法的同质区降噪效果与PPB方法相接近，但本发明方法在细节纹理保持方面有明显优势；将图3，图5与图4对比，能够发现SAR-BM3D的抑斑能力最弱，但在纹理细节保持方面最好，并且本发明方法的细节保持能力介于PPB方法和SAR-BM3D方法之间。

表1 同质区1，2的ENL指标

	noise	PPB方法	SAR-BM3D方法	本发明方法
					同质区1	9.99	182.17	72.83	190.49
同质区2	10.21	81.02	54.73	80.39

表1给出了图2中的两个矩形同质区降噪前后的ENL情况，其中ENL值提高越多表示降噪效果越强。可以看出本发明方法与PPB方法相近，且均远高于SAR-BM3D方法，此结果与降噪效果图相吻合。

上述实验表明，本发明降噪方法在降噪后图像的视觉效果及客观评价指标上都较好，由此可见本发明对SAR图像降噪是有效的。

Claims (1)

1.一种基于非本地分类稀疏表示的SAR图像降噪方法，其特征在于具体步骤如下：

步骤一、图像内容的分类

SAR图像中不同区域内噪声分布不同，并且不同区域间图像本身的特性也有所不同，为实现对不同图像内容采用不同方法降噪，将图像块的方差与预设的阈值进行比较，当方差大于阈值时视为异质区，小于阈值时视为同质区；

步骤二、分类估计非本地系数

结合分类后SAR图像中同质区和异质区的特点，选择不同的方法对其进行降噪处理；首先以图像的每个像素点为中心，依次抽取图像块，若图像块属于异质区，则根据异质区图像块包含奇异结构和细节边缘而难以在像素点附近范围内找到真正相似图像块的特点，选择较小的块尺寸并在较大的邻域范围内搜索与该目标块相似的图像块，并在得到相似图像块后使用3D变换域硬阈值收缩方法实现对目标图像块的降噪；对于同质区目标块，由于同质区内的像素点所包含的真实信号服从相同的分布，因此设置较大的块尺寸和较小的搜索范围以寻找相似块，并利用加权平均的方法实现对同质区目标块的降噪；然后利用主成分分析法训练出自适应学习字典，最后将降噪后的图像块在对应的自适应学习字典上投影得到图像块的非本地稀疏系数；

步骤三、结合非本地系数与稀疏表示的SAR图像降噪

为进一步提高图像降噪的性能，利用步骤二得到的非本地系数调节稀疏解使其更接近于真实解，即将非本地系数代入式(1)求解：

$\widehat{\mathrm{\&alpha;}}=\arg \min {\left|\right|y-\mathrm{D\&alpha;}\left|\right|}_2^2+{\mathrm{\&lambda;}}_1\stackrel{N}{\mathrm{\&Sigma;}}{\left|\right|{\mathrm{\&alpha;}}_i\left|\right|}_1+{\mathrm{\&lambda;}}_2\stackrel{N}{\mathrm{\&Sigma;}}{k}_i{\left|\right|{\mathrm{\&alpha;}}_i-{\mathrm{\&gamma;}}_i\left|\right|}_1+{\mathrm{\&lambda;}}_3\stackrel{N}{\mathrm{\&Sigma;}}(1-k_i){\left|\right|{\mathrm{\&alpha;}}_i{-\mathrm{\&eta;}}_i\left|\right|}_1$   式(1)

其中y表示含噪图像，D为字典，N为图像所包含的像素点个数，k_i为分类后像素所属类型的标记，α_i为本地稀疏系数，γ_i和η_i为非本地稀疏系数，λ₁，λ₂和λ₃为权重调节参数；式(1)右端前两项为图像的本地稀疏表示，后两项为使本地稀疏系数α_i更接近于真实解的约束项；由于图像内容仅分为同质区和异质区，k_i仅为1或0，因此针对每个图像块，求解式(1)时均会从最后两项中选择出一项对进行α_i约束，此时的优化问题转化为二范数求解问题；在确定了参数λ₁，λ₂和λ₃后，先将目标图像块在字典上投影得到本地稀疏系数α_i，然后利用非本地系数γ_i或η_i调节α_i，得到最终稀疏解并将重构得到降噪后图像；将重构后的降噪图像与原含噪图像非凸组合后继续执行步骤二和步骤三，当前后两次降噪结果的平均绝对偏差MAD小于设定的阈值时停止循环，得到最终降噪图像。