一种基于剪切波系数处理的SAR图像降噪方法

摘要

本发明公开了一种基于剪切波系数处理的SAR图像降噪方法。属于数字图像处理技术领域。它是利用图像剪切波变换后系数具有的稀疏特性，首先建立基于图像剪切波系数的稀疏表示模型，而后通过分段正交匹配追踪StOMP算法实现统计均值意义上稀疏表示系数的无偏估计，并将稀疏表示后的剪切波系数重构为降噪后图像；为弥补稀疏表示中丢失部分系数对图像细节的损失，并利用这部分系数对应的剪切波函数具有提取图像边缘细节的能力，针对图像在丢失系数对应的剪切波函数空间中投影重构的结果，结合基于能量泛函的总变分TV方法进一步迭代去噪，最终得到细节丰富的去噪图像，既抑制了SAR图像斑点噪声又保持了图像的细节纹理，可用于SAR图像降噪。

说明书

技术领域

本发明属于数字图像处理技术领域，它特别涉及SAR图像降噪方法，用于 对SAR图像进行降噪处理。

背景技术

合成孔径雷达在成像中具有全天时，全天候的探测与侦察跟踪能力，能有 效的识别伪装和穿透掩盖物，因此SAR图像被广泛应用于航空摄影测量与遥感、 卫星海洋观测、航天侦察、图像匹配制导、深空探测等方面。但由SAR相干成 像机理给图像带来的斑点噪声，对目标识别及图像压缩等后期处理带来不利影 响。能否有效的滤除斑点噪声已成为图像后续解译的重要前提。

由于SAR图像拥有丰富的纹理和边缘，因此在有效滤除斑点噪声的同时充 分保留图像的纹理和边缘，是SAR图像降噪处理的重点。近年来，由于稀疏表 示在自然图像降噪中具有的良好性能，正逐渐成为图像降噪的有效方法之一。 但此类方法在稀疏表示的同时，难以避免的存在图像细节损失的问题。为提高 图像细节保持能力，已有方法大多局限于对稀疏降噪后的图像进行处理，但由 于稀疏表示中损失的信息无法恢复，其改善程度有限，不能明显的修复降噪图 像损失的边缘细节。

发明内容

本发明为克服现有SAR图像降噪方法存在图像边缘细节信息损失的不足， 以获得细节清晰的降噪SAR图像，提供了一种基于剪切波系数处理的SAR图像 降噪方法。该方法充分考虑到剪切波系数的特性，首先建立剪切波系数稀疏表 示降噪模型，然后使用TV方法对图像进一步修复，不仅可以很好的抑制SAR 图像的斑点噪声，而且解决了降斑图像的细节纹理保持问题，因此该方法可以 有效的实现SAR图像降噪。包括以下步骤：

步骤一、图像噪声模型转换

使用非对数加性噪声模型，将输入SAR图像中均值为1的乘性噪声转化为 均值为0的加性噪声。

步骤二、剪切波域稀疏降噪

首先对噪声图像进行剪切波变换得到图像的剪切波系数w。为实现系数的 稀疏表示，用测量矩阵Φ对系数w进行变换，即y＝Φw。进一步假设y在测量矩 阵Φ下的稀疏逼近表示为z。得到如下最优化稀疏表示模型：

$g\left(a\right)={\left|\right|y-\mathrm{\Phi z}\left|\right|}_2^2+\gamma {\left|\right|z\left|\right|}_0$  式(1)

当式(1)取最小值时对应的最优解z即为原系数w的稀疏表示且稀疏 表示的均值均为干净图像剪切波系数均值的无偏估计。使用StOMP算法求解 式(1)，首先将y投影到测量矩阵Φ的每个原子上，然后将大于阈值的投影值对 应原子放入索引集I₁，通过索引集中原子对y进行逼近后得到残差R¹y。为采用 更多的原子进一步逼近y，可对残差R¹y采用相同的方式进一步分解。若第m次 迭代的残差为R^my，则R^m+1y可表示为:

R^my＝<R^my,φ_rm>φ_rm+R^m+1y  式(2)

其中φ_rm表示第m次迭代时，从测量矩阵Φ中选出的原子，将这些所选原子并入 I_m-1得到更新的I_m。如果满足终止条件||R^my||²≤ε或|I_m|>L则迭代结束，否则继续 迭代。当迭代结束时可以通过式(3)得到最优解z

${\left(z_m\right)}_{I_m}={\left({\Phi }_{I_m}^T{\Phi }_{I_m}\right)}^{-1}{\Phi }_{I_m}^{T}y$  式(3)

优化解z即为原系数w的稀疏表示将及w稀疏表示时丢弃的小系数进行剪 切波反变换可得到降噪后图像u_s及残差图像u_d，并将稀疏表示中丢弃的小系数对 应的剪切波空间记为M。

步骤三、TV降噪与细节修复

根据剪切波系数特性，结合步骤二处理结果，建立基于能量泛函的TV模型：

$F\left(u\right)=\underset{\Omega }{\int }\phi \left(\right||▿P_S\left(u\right)|\left|\right)\mathrm{dxdy}+\frac{\lambda }{2}\underset{\Omega }{\int }{(u-u_0)}^2\mathrm{dxdy}$  式(4)

其中P_S(u)如式(5)所示：

$P_S\left(u\right)=\underset{j,l,k\in M}{\Sigma }<u,{\psi }_{j,l,k}>{\psi }_{j,l,k}$  式(5)

其中M为j,k,l所属集合，该集合由步骤二稀疏表示中丢失系数对应的剪切波函 数确定，P_S(u)表示图像u在这部分剪切波函数空间上投影重构的结果。当 时式(3)的能量泛函取到极小值，此时的u即为TV处理后的去噪图像。 由式(3)可得如下所示：

$\frac{\partial u}{\partial t}=▿(\frac{{\phi }^{\prime }\left(\right||▿\mathrm{Ps}\left(u\right)|\left|\right)}{\left|\right|▿\mathrm{Ps}\left(u\right)\left|\right|}▿\mathrm{Ps}\left(u\right)-\lambda (u-u_0)$  式(6)

式(6)第一项为扩散项，令ρ(x)＝φ'(x)/x，其中则ρ(x)可对扩散过 程中的平滑度进行控制。

使用最速下降法通过式(7)求解式(6)

u^k+1＝u^k+Δt[η(Ps(u^k))-λ(u^k-u₀)]  式(7)

其中Δt为迭代步长，k为第k次迭代，u₀为原始含噪 图像，迭代初始值为u_s+Δt[η(u_d)]，当前后两次迭代间的MAD小于某一门限时结 束迭代，获得最终的降噪图像。

本发明的创新点是利用图像剪切波系数的稀疏性与提取图像边缘的能力， 利用稀疏表示中保留系数对干净图像剪切波系数无偏估计的优势进行降噪，并 针对稀疏表示过程中图像损失的边缘纹理信息，利用剪切波小系数对图像噪声 和边缘细节的区分特性，针对图像在该空间的投影结果采用TV处理，进一步实 现了对图像的降噪及细节修复。

本发明的有益效果：使用基于统计最优的降噪模型来稀疏表示剪切波系数 而不是传统的软硬阈值处理方法，同时在求解模型最优化问题时，采用了StOMP 算法既保证了求解精度又提高了运算效率，因此整体降噪效果可达到较高水平； 并建立了基于能量泛函的TV处理，在降噪的基础上进一步对稀疏表示过程中损 失的信息进行修复，因此该方法既能在降噪的同时保持图像纹理细节，还具有 较高的运算效率。

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

附图说明

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

图2是本发明仿真使用的真实含噪的丘陵SAR图像；

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

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

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

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

具体实施方式

参照图1，本发明基于剪切波系数处理的SAR图像降噪方法，具体步骤包 括如下：

步骤一、图像噪声模型转换

采用式(8)非对数噪声转化模型将SAR图像乘性噪声模型转化为加性噪 声模型

I＝RX＝X+(R-1)X＝X+N  式(8)

其中I为被噪声污染的图像强度，R代表相干斑噪声，X代表地物真实的后向散 射强度，N＝(R-1)X为零均值的加性噪声。

步骤二、剪切波域稀疏降噪

将噪声图像进行剪切波变换，得到水平锥和垂直锥各5个尺度下的剪切波 系数，对各尺度下的系数进行按列稀疏。首先构造随机测量矩阵Φ＝{φ_r}_r∈Γ，然 后对其归一化。假设w表示系数的列，用测量矩阵Φ对w进行变换，即y＝Φw， 然后建立最优化稀疏表示模型并使用StOMP算法求解稀疏 模型。StOMP算法初始时稀疏逼近值z₀＝0，残差R⁰y＝y，所选原子下标索引集 I₀为空。当迭代次数m＝1时，首先求R⁰y与Φ各原子的内积，即R⁰y在各原子上 的投影，然后引入硬阈值t₁并通过(9)式选择与R⁰y最匹配的几个原子

I₁＝{i:|<φ_r,R⁰y>|≥t₁,r∈Γ}  式(9)

正交化所选原子，重新将R⁰y投影到正交化后的原子上，得到第一次稀疏逼近后 的残差R¹y，判断索引集中原子下标个数是否小于设置的稀疏度L，若小于继续 迭代。当循环结束时通过式(3)计算出最终的z_m。L的大小为w中大于阈值的 元素个数，而阈值是根据噪声大小的不同在0.3-0.8之间取值。通过剪切波反变 换重构出降噪后的图像u_d。根据稀疏后的系数找到稀疏表示过程中丢弃的小系 数及它们对应的剪切波所属集合M，并将这部分小系数重构后得到残差图像u_s。

步骤三、TV降噪与细节修复

首先建立基于能量泛函的TV模型，当F(u)的导数等于0时，此能量泛函取 得最小值。式(6)的离散化形式如式(7)，其中迭代步长Δt＝0.1，λ＝0.15。通过u⁰＝u_s+Δt[η(u_d)]求得迭代初始值u⁰， 将u⁰带入式(7)迭代循环，循环过程中的Ps(u^k)通过式(5)求得。在每次循环 之前通过式(10)计算图像的平均绝对偏差

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

其中N为图像像素值大小；当MAD小于设置门限ε时，停止迭代，得到最终修 复后图像u^k，否则继续修复u^k，其中门限ε设置为2.5e-04。

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

一、实验条件和内容

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

实验内容：在上述实验条件下，使用KSVD方法、KSVD和普通TV相结 合的KSVD_TV方法与本发明方法进行对比实验。去噪效果的客观评价结果用 同质区均值μ，方差σ²，等效视数ENL衡量。

实验1：用本发明方法和现有的KSVD及KSVD_TV方法分别对图2进行 降噪，其中KSVD方法重叠图像块大小为8×8，字典大小为64×256；KSVD_TV 方法TV迭代次数为10，步长为0.2，扩散项与保真项之间权重λ为0.15，得到 的降噪效果分别如图3和图4所示；本发明中β＝5，λ＝0.15，迭 代步长Δt＝0.1，迭代次数为20，降噪结果如图5所示。对比图3和图4可以看 出本发明降噪效果更强，同时拥有更丰富的边缘细节信息。

实验2：在图2中选择两块矩形同质区域，用本发明和KSVD及KSVD_TV 方法对它们分别去噪。用均值μ，方差σ²，等效视数ENL作为去噪效果的评价 指标，并列在表1中。

表1 同质区1，2不同降噪算法的性能参数

表1的结果表明，本发明方法降噪前后图像的均值改变量为几种他方法中 最小的，说明降噪时保持雷达辐射特性能力较强；方差的降低及NEL的大幅提 高表明了本发明方法具有很强的降噪能力。

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