NSST域三维块匹配的SAR图像相干斑抑制算法

余 翔，李 娅，王诗言

(重庆邮电大学 通信与信息工程学院，重庆 400065)

0 引 言

星载合成孔径雷达(synthetic aperture radar，SAR)可实现高分辨率、全天时、全天候成像，在海洋观测、资源勘探和军事侦察领域具有很大的优势[1-3]。然而，SAR的相干成像机制会导致SAR图像中产生乘性的相干斑噪声，为SAR的后续应用如图像分割、目标检测及识别增加困难。因此，多年来许多科研工作者们致力于对SAR图像进行相干斑抑制研究，尽管已经出现了相当多的有效算法，但在相干斑去除及边缘纹理保持之间仍然难以兼顾。根据相干斑抑制发生在成像前后可分为两大类：成像前的多视处理和成像后的滤波处理。其中，多视处理方法较为单一，会导致SAR图像的几何分辨率降低，因此，人们常常将精力聚焦于滤波处理方法。

最早的滤波方式为空域滤波，主要有Lee滤波[4]、Frost滤波[5]、Gamma MAP滤波[6]及相关的改进增强方法[7-8]，该类算法复杂度较低，但无法保留SAR图像的细节信息。随着小波理论的发展和完善，基于小波变换的滤波被广泛应用于SAR图像相干斑抑制[9]，然而小波变换仅能表示图像3个方向的特征，它对于图像的边缘处理并不是最优的。鉴于上述缺点，各种多尺度变换方法应运而生，其中曲线波变换、轮廓波变换和剪切波(Shearlet)变换是使用较为广泛的方法。而Shearlet变换结合了曲线波变换和轮廓波变换的优点，不仅能够对图像进行最优的稀疏表示，而且在方向选择上十分灵活且易于实现，因此，得到了广泛的研究及应用[10-11]。NSST是Shearlet变换的扩展，由于没有进行下采样操作，NSST除了具有Shearlet变换所具有的优点以外，还具有平移不变性及更好的方向选择性[12-14]。

由于图像具有自相似性，文献[15]在空域滤波的基础上提出了非局部均值(non-local means，NLM)滤波方法，之后，李敬曼等[3,16]将该方法用于相干斑抑制并进行了改进。文献[17]将NLM与变换域滤波进行完美结合提出了三维块匹配(block-matching 3D，BM3D)算法，该算法在去除自然图像中的高斯噪声时十分有效。文献[18]基于BM3D框架提出了SAR-BM3D算法，该算法在相似块匹配阶段使用概率相似性度量[19](probabilistic patch-based，PPB)方法，在变换域利用局部线性最小均方误差准则获得更准确的收缩。目前SAR-BM3D算法被认为是相干斑抑制效果最好的算法之一，但同时它具有NLM所共有的缺点：相似块匹配效率过低，相似性度量容易受到搜索窗口大小和块大小的影响。在文献[20]中，作者基于灰度理论的灰色关联分析法提出了一种灰度理论块匹配方法，可用于定量表明图像块之间的相似性，该方法可以很好地克服之前已有块匹配方法的缺陷。文献[21]将灰度理论与加权核范数最小化结合应用到相干斑抑制中，取得了良好的效果。

受以上研究的启发，本文针对SAR图像相干斑抑制和细节信息保持之间难以平衡的问题，基于BM3D框架提出一种新的相干斑抑制算法。该算法相对于原有算法在块匹配阶段进行了改进，因为不管是BM3D中基于二范数的方法还是SAR-BM3D中的PPB方法都有如下缺陷：不仅匹配效率低，且相似性大小易受搜索框及块大小的影响。因此，本文采用更为先进的灰度理论进行相似性度量，实验的第1部分即是对相似块匹配方法的对比讨论。另外，为了更好地描述图像的细节，在变换阶段使用最优稀疏性的非下采样剪切波变换(non-subsampled Shearlet transform，NSST)将相似块组变换到NSST域，通过对NSST系数进行收缩从而滤除相干斑噪声。同BM3D思想一样，本文利用基本估计图像和原噪声图像进行联合滤波，从而得到具有良好相干斑抑制效果并成功保留图像细节信息的最终估计图像。通过实验分析，本文的相干斑抑制算法比现有的GT-WNNM、增强Lee滤波、Gamma MAP滤波、SAR-BM3D、SAR-NSST算法具有更好的细节保持能力，同时相干斑抑制效果相当。

1 NSST域三维块匹配的相干斑抑制算法

本节将介绍如何基于BM3D框架提出本文的相干斑抑制算法。假设相干斑完全发育，观察到的反向散射信号可以表示为

I(x,y)=R(x,y)·u(x,y)

(1)

(1)式中：(x,y)表示SAR图像像素点的二维坐标；I(x,y)为雷达回波信号的实际强度；R(x,y)为地物目标后向散射系数，表示雷达回波信号的真实强度；u(x,y)为相干斑噪声。通常认为R(x,y)与u(x,y)相互独立，对(1)式进行对数变换可将乘性噪声转化为加性噪声进行处理，但这会导致图像出现失真，可以通过改写(1)式的形式将乘性噪声转化为加性噪声。

I(x,y)=R(x,y)+[u(x,y)-1]R(x,y)=R(x,y)+u′(x,y)R(x,y)=R(x,y)+v(x,y)

(2)

由于R(x,y)和u(x,y)相互独立，并且u′(x,y)为零均值，因此，v(x,y)是方差取决于R(x,y)且均值为零的加性噪声，并且与R(x,y)不相关。得到这样的加性噪声模型之后，代表相干斑可以使用本文算法进行去噪。

SAR图像受相干斑的影响，每个像素点将产生随机乘性噪声。传统的NLM方法相似性度量易受块匹配参数影响，导致去噪图像边缘模糊且相似块匹配效率低。而基于灰度理论的块匹配方法能有效利用图像结构信息，可在样本数量贫乏的情况下探究图像块之间的关联性，因此，本文将灰度理论块匹配方法引入到BM3D算法中对SAR图像进行相干斑抑制，提高BM3D算法块匹配的鲁棒性。另外，小波变换仅可表示3个方向的细节，不能对图像进行最优的稀疏表示。因此，为了表示SAR图像更多方向的细节并使频域具有平移不变性，本文在变换阶段引入了NSST，将相似块组变换到最优稀疏的NSST域进行收缩从而滤除相干斑噪声。NSST域三维块匹配的相干斑抑制算法利用了SAR图像的自相似性和稀疏性，并对自相似性和稀疏性的利用做了改进，以达到去除相干斑的同时尽量保存SAR图像细节的目的。

1.1 基于灰度理论的相似性度量

传统的NLM中的相似性度量方法会降低相干斑抑制算法边缘的保持能力，而在SAR-BM3D算法中，虽然概率相似性度量使得边缘保持能力得到了一定的提升，但它的计算复杂度很高，且会导致图像过度平滑，边缘保持能力依旧不够理想。

在客观世界中，内部特征结构完全未知的系统通常表示为“黑色”，内部特征已得到充分研究的系统则表示为“白色”，而大部分系统属于“灰色”，即该系统的一部分信息已知，但另一部分有待探究。含噪图像包含原始信息和噪声，可被视为灰色系统，从而利用灰度理论进行噪声滤除操作。当某些像素受到污染或者干扰时，灰色理论所提出的权重函数能更健壮地表示2个图像块的相似性。灰色关联分析是用于分析灰色系统的重要方法，它能够有效利用图像的结构信息来衡量2个图像块的相关性。

设搜索框内的图像块组被定义为B=[B0,B1,…,Bi,…,Bj,…，Bs]，参考块(待恢复图像块)Bi=[Bi(1),Bi(2),…,Bi(M2-1),Bi(M2)]和匹配块Bj=[Bj(1),Bj(2),…,Bj(M2-1),Bj(M2)]是需要进行相似性度量的2个图像块，其中，M表示正方形图像块的长度。可利用(3)式计算它们的关联系数。

(3)

2个图像块之间的相似性权重定义为

(4)

(4)式中：i=(0,1,…，s-1,s)；j=(0,1,…，s-1,s)；k∈[1,M2],k∈Z，Z为整数。

1.2 NSST域相干斑抑制算法

Shearlet变换是将几何学与多分辨分析结合，利用经典仿射理论得到的新型多尺度几何变换方法。对于二维空间，具有复合扩张的仿射系统为

j,l∈Z,k∈Z2}

(5)

(5)式中：A=[4,0;0,2]，B=[1,1;0,1]分别为各向异性膨胀矩阵和剪切矩阵，且|detB|=1；函数ψj,k,l在频域中的支撑范围为

(6)

图1 Shearlet函数的支撑区间Fig.1 Support interval of Shearlet function

图2 Shearlet变换频域剖分图Fig.2 Tiling of frequency plane of Shearlet transform

(7)

(7)式中：j≥0；-2j≤l≤2j-1；k∈Z2；d=0,1。

NSST的离散化需要利用剪切滤波器(Shearlet filter，SF)和非下采样金字塔(non-subsampled pyramid，NSP)滤波器组来实现，其中，SF和NSP分别用于实现多方向分解和多尺度分解。首先经过k次NSP分解产生一个低频子带和k个高频子带，然后使用SF对各个子带进行多方向分解，最终得到多子带多方向的NSST系数。NSST将标准的SF从伪极化网络坐标映射到笛卡尔系统来实现平移不变性。NSST的3层分解示意图如图3所示。

图3 NSST的3层分解示意图Fig.3 Three level decomposition of NSST

本算法中的NSST域相干斑抑制具体步骤如下。

步骤1进行NSST变换得到各尺度各方向上的非下采样剪切波系数；

步骤2对变换系数进行处理(基础估计时使用硬阈值滤波，最终估计时使用联合维纳滤波)得到收缩的NSST系数；

步骤3NSST逆变换；

步骤4按照权重得到去噪的参考块。

1.3 本文算法流程

综上所述，基于NSST域三维块匹配的相干斑抑制算法流程如图4所示。

图4 基于NSST域三维块匹配的相干斑抑制算法流程图Fig.4 Flowchart of suppression algorithm based on block matching 3D in NSST domain

相干斑抑制算法总体上分为2个步骤：基础估计和最终估计。2个估计阶段的处理步骤基本一致，都经过了图4虚线框中的灰度相似块匹配、NSST域系数收缩、图像块聚集3个步骤。但2个阶段采用的滤波方法不同，基础估计阶段采用硬阈值滤波将信号和噪声分离，从而实现图像更稀疏的NSST表示；最终估计阶段输入了2个相似块组，它相对于第1阶段可以更好地利用相对干净的基础估计图像设计联合维纳滤波器，然后输入原始SAR图像执行实际的去噪过程。

1.3.1 基础估计

基础估计的步骤如下。

步骤1在N×N的搜索框内进行块分组，对于每个M×M含噪参考块利用(4)式计算所有匹配块与之的相似度，选择最相似的前Q-1个图像块，与参考块组成相似块组；

步骤2对相似块组使用NSST变换获得多尺度多方向的NSST系数，使用硬阈值滤波得到相应的收缩系数；

步骤3对处理后的变换系数进行NSST反变换并根据步骤1中的相似权重计算得到相对干净的去噪参考块；

步骤4遍历完整幅SAR图像后将所有去噪参考块按照位置进行聚集，输出基础估计SAR图像。

1.3.2 最终估计

最终估计的步骤如下。

步骤1对基础估计获得的SAR图像重复基础估计的步骤1得到相似块组。此时有2个相似块组，一个来自基础估计SAR图像，另一个来自原SAR图像，但位置由基础估计SAR图像相似块组的位置所决定；

步骤2在NSST域对2个相似块组进行联合维纳滤波并利用NSST反变换得到收缩系数；

步骤3按位置聚集得到最终估计SAR图像。

2 实验结果与分析

在本节中，将实验分为2部分：①对灰度相似性度量的优越性进行实验验证；②将本算法与其他相干斑抑制算法进行对比。实验中用到的图片如图5所示，分别为图5a Bedfordshire(256×256)，图5b Horsetrack(256×256)。实验采用Matlab语言在64位Windows10系统下进行，计算机硬件为Intel Core i5-7400处理器、8.00 GByte内存、3.00 GHz主频。2部分的实验都将从主观和客观进行评价，本文使用的4个客观评价指标分别为PSNR[21]、等效视数(equivalent numbers of looks，ENL)[21]、边缘保持指数(edge preservation index，EPI)[21、结构相似性指数(structural similarity，SSIM)[21]。其中，PSNR数值越大，表明算法具有更强的相干斑抑制能力；ENL越大，表明SAR图像的相干斑越弱，即该算法处理后的图像具有更好的平滑度；EPI越大，表明算法具有更强的边缘保持能力；SSIM越大，表明相干斑抑制后的图像与原图像越相似，失真越小。

图5 真实SAR图像Fig.5 Real SAR images

2.1 相似块匹配方法对比

为验证灰度相似性度量在本算法中的有效性，采取数学中的控制变量法，即只改变相似性度量方法，其他步骤均采用本文提出的方法，观察使用灰度相似性度量的情况下，本算法是否表现最好。实验中进行对比的相似性度量方法为：二范数[17]、PPB[19]、灰度理论。下面从主观和客观2方面评价。Bedfordshire和Horsetrack上各相似块匹配方法的对比分别如图6和图7所示，可以直观地看出，利用灰度理论进行相似性度量时，图片更好地抑制了相干斑，它能更准确地搜索相似块组并将它们变换到NSST域进行相干斑抑制。

图6 Bedfordshire上各相似块匹配方法的对比Fig.6 Comparison of various similar block matching methods on Bedfordshire

图7 Horsetrack上各相似块匹配方法的对比Fig.7 Comparison of various similar block matching methods on Horsetrack

为了更客观地体现灰度相似块匹配的优越性，表1列出了几种方法处理后的评价指标。从表1可以看出，除Horsetrack图上SSIM指标稍逊于PPB以外，GT方法的其他指标均表现更好。

表1 各相似块匹配方法评价指标对比Tab.1 Comparison of evaluation indexes of similar block matching methods

2.2 相干斑抑制算法对比

在本实验中，将本文算法(NSST-BM3D)与经典的增强Lee滤波[7]、Gamma MAP滤波[6]以及去噪效果良好的GT-WNNM[21]、SAR-BM3D[18]、SAR-NSST[12]进行对比。同时，为了充分说明本算法相对于SAR-NSST算法的优越性，本实验增加了真实SAR图像的 Kilauea(512×512)和Pentagon(590×490)原斑点图，如图8所示。下面分别从主观视觉和评价指标对本算法进行对比分析。

图8 真实SAR图像Fig.8 Real SAR images

2.2.1 主观评价

各相干斑抑制算法的对比如图9—图12所示。主观视觉上，可以看出，使用增强Lee滤波的效果最差；Gamma MAP滤波会残余大量噪声；使用SAR-BM3D会使图像过度平滑，不利于保留细节信息；GT-WNNM整体上优于SAR-BM3D，但是依旧会有边缘过度平滑的弊端；SAR-NSST滤波后的图像会产生许多不属于图像本身的纹理；而使用本文算法处理后的图像避免了上述缺点，较好地保持了图像边缘并抑制了相干斑。

图9 Bedfordshire上各相干斑抑制算法的对比Fig.9 Comparison of various speckle suppression algorithms on Bedfordshire

图11 Kilauea上各相干斑抑制算法的对比Fig.11 Comparison of various speckle suppression algorithms on Kilauea

图12 Pentagon上各相干斑抑制算法的对比Fig.12 Comparison of various speckle suppression algorithms on Pentagon

2.2.2 客观评价

各相干斑抑制算法评价指标如表2所示。

表2 各相干斑抑制算法评价指标Tab.2 Evaluation index of each speckle reduction algorithm

从表2的各相干斑抑制算法参考指标可以看出，本方法在整体上比之前最具有竞争力的SAR-BM3D、GT-WNNM算法表现更好。在PSNR指标上，除在Bedfordshire上本方法略低于SAR-NSST外，在剩余3个场景中本方法均表现最好；在ENL指标上，本方法除在Pentagon上稍低于GT-WNNM外，其他场景均优于GT-WNNM和SAR-NSST算法，而增强Lee滤波和SAR-BM3D在ENL指标上更高的原因是图像产生了很大的模糊，因而平滑度更好，但也使得自身的EPI指标降低。由于前3个场景经过SAR-NSST算法处理后图像产生了很多不属于图像本身的纹理，因此，在EPI指标上，本方法低于SAR-NSST，但均高于其他4种算法；而对于SSIM指标，本方法在4个场景中都是表现最好的算法。因此，综合4个评价指标进行考量，本算法达到了改进相干斑抑制效果的目的。

3 结束语

本文提出了一种新的基于BM3D的SAR图像相干斑抑制算法。为了适应于SAR图像特性，减少计算量并增强相似块匹配的鲁棒性，在块匹配阶段使用灰度理论进行相似性度量；为了获得更多方向的细节，本文在收缩阶段采用NSST变换，在NSST域采取硬阈值滤除相干斑噪声；为了更精确地估计，在第2阶段使用维纳滤波对基础估计图像和原噪声图像进行联合滤波。在主观视觉及PSNR、ENL、EPI、SSIM这4个指标上进行了2个基于真实SAR图像的测试实验，分别说明了灰度块匹配方法在本算法中的有效性以及本算法相对于GT-WNNM、增强Lee滤波、Gamma MAP滤波、SAR-BM3D、SAR-NSST这4个算法在整体上的优越性。