唐 璐,吴亚娟,李 天
(西华师范大学计算机学院,四川南充 637001)
合成孔径雷达(SAR)图像是通过从卫星或飞机向目标发射微波并使用目标反射的信号获得图像,由于它是通过主动传感器获取雷达图像,故含有不同于常规光学传感器获得的一般图像特征的噪声,即散斑噪声.散斑噪声具有与乘性噪声相似的特性,并呈瑞利分布[1],导致图像边缘信息和关键特征丢失,图像质量降低.SAR图像去噪的常用方法有非线性滤波、变换域滤波和空间域滤波等.自适应平均滤波器和中值滤波器等非线性滤波器虽能有效降低散斑噪声,但会使图像的重要区域和边缘区域模糊而丢失信息,因此不适用于散斑噪声的去除.为了解决该问题,Tomasi等提出了双边滤波器[2],通过考虑相邻像素的空间特性及相邻像素的阴影值来确定滤波器的权重从而保持局部边界,但这也增加了计算复杂度.为了降低计算复杂度,加快进程执行,Zhang提出的自适应双边滤波器[3],不仅允许高斯范围核的中心和宽度随像素变化,计算复杂度得以降低,而且在图像增强和去噪方面也不逊于传统的双边滤波器.Gavaskar提出了一种快速自适应双边滤波器[4],它用一个多项式替换直方图,用一个积分替换有限范围空间和,不仅使进程执行快20倍,且在图像去噪方面也具有有效性.
在变换域滤波的方法中,离散小波变换对于SAR图像去噪十分有效.小波是一种适用于非平稳连续信号展开的函数,它们构成小波变换的核,使信号从时域映射到频域,主要特点是在不同的频率下可以获得不同的分辨率.Pesquet提出了平稳小波变换,使小波分解时不变[5],提高了小波在信号去噪和增强中的功率.对于平稳或冗余变换,在每个尺度[6]上进行滤波卷积之前,先进行上采样,而不是下采样.因此,它可以产生更精确的频率定位信息,这种变换的冗余性有助于识别信号中的显著特征,特别是识别噪声.张玉叶[7]提出的一种平稳小波方向能量阈值滤波的SAR图像去噪方法,去噪效果较好.
在空间域滤波的方法中,非局部均值算法是SAR图像去噪的热门算法.Buades在2005年提出了非局部均值算法(也称为NLM)[8],采用图像的块相似性代替传统的单像素相似性来构造权重,充分利用图像中的冗余信息,在去噪的同时最大限度地保持图像的细节特征.针对SAR图像的乘性噪声,文献[9]引入均值比和变异系数来构造权函数.文献[10]将图像对数转化为加性噪声,并结合变异系数进行NLM去噪.文献[11][12]结合双边滤波的思想,且Zhang等人在考虑块相似性时考虑了中心像素的相似性,也取得了良好的滤波效果.
在上述研究的基础上,本文融合平稳小波变换(SWT)、快速自适应双边滤波(FABF)、小波阈值去噪和优化贝叶斯非局部均值算法(OBNLM),提出改进快速自适应双边滤波(SFOBNLM)算法,以提高SAR图像质量.
SFOBNLM算法的流程如图1所示.首先,针对SAR图像固有的散斑噪声,利用平稳小波变换将图像分解为近似子带和细节子带,对近似子带采用快速自适应双边滤波(其执行效率高,运算速度快,去噪效果佳);对细节子带采用软阈值滤波,保留图像特征信息.然后,利用小波逆变换对去噪后的图像进行重构,使重构后的图像噪声减少且特征信息明显.最后结合OBNLM算法进一步去除噪声,强化边缘细节信息,输出滤波图像.
图1 改进快速自适应双边滤波(SFOBNLM)算法流程
应用平稳小波变换,将图像实际分解为如图2(a)所示的四个子带.这四个子带由小波变换分析滤波器组的低通和高通滤波器的可分离应用而产生;分解提供了对应于不同分辨率的方向的子带.LH1、HL1和HH1子带代表最细尺度的小波系数(即细节图像),LL1子带代表粗尺度的小波系数(即近似图像).利用相似滤波器组对子带LL1进行进一步分解,得到下一个粗级小波系数.这些结果在二级小波分解中如图2(b)所示.本文将SWT应用于二级分解图像中,对噪声图像进行SWT处理,将图像分解为近似子带和细节子带,对近似子带进行快速自适应双边滤波处理.对于第二次分解,将SWT应用于一级分解图像的近似子带上,然后对细节子带进行小波阈值处理,并利用小波逆变换对恢复后的图像进行重构.
图2 图像的SWT分解
Gavaskar提出的一种快速的自适应双边滤波算法[4],其复杂度不随空间滤波器宽度的变化而变化,从而加快图像处理执行速度,且能有效进行纹理滤波.SAR图像经平稳小波变换后产生的近似子带进行快速自适应双边滤波,可有效快速地去除子带中大量噪声.
快速自适应双边滤波是通过用一个多项式替换直方图,用一个积分替换有限范围空间和,使用解析函数近似滤波器对图像进行处理.特别地,多项式是通过匹配其矩阵与目标直方图的矩阵进行拟合的,解析函数是通过分部积分递归计算的.本文使用文献[2]中的经典双边过滤器的定义:
设f:I→R为输入图像,I⊂Z2是图像域,则输出图像g:I→R为:
这里Ω是一个以原点为中心的窗口,窗口通常设 置 为Ω=[-3ρ,3ρ]2,η(i)表 示 归 一 化 指 数,φi:R→R是局部高斯值域核,式(1)中的中心θ(i)和式(3)中的宽度σ(i)是空间变化函数,w(j)表示空间域权值.
从式(1)(2)可以看出,经典双边滤波器的蛮力计算要求每像素进行O(ρ2)运算,其运算代价大,使得实时实现具有挑战性.本文采用的快速自适应双边滤波器模型[6]计算复杂度为O(1),执行效率高,运算速度快,图像去噪效果佳,模型如下:
对于i∈I,设:
对于t∈Λi,定义加权直方图:
其中δ是Kronecker增量,即δ(t)=1,t=0或δ(t)=0,t≠0.
定义边界Λi,设:
定义在区间[αi,βi]上的函数pi近似hi,根据式(6)并用积分替换式(1)(2)中的有限和,即将式(1)(2)近似如下:
若pi是一个多项式,设pi(t)=c0+c1t+…+cNtN,通过多项式直方图,得到滤波器输出:
本文将平稳小波变换分离出的近似子带通过快速自适应双边滤波器进行处理,首先将近似子带图像用一个多项式替换该图像的直方图进行表示,通过快速卷积进行拟合,最后用一个积分替换有限范围空间和,将分部积分递归计算输出结果.计算过程快速有效,输出的结果图像在去噪的同时细节信息也得以保留.
小波阈值化的动机[13]是由于小波变换擅长能量压缩,信号通过小波变换后,系数较小的为噪声组成部分,系数较大的为信号的重要特征.这些系数可以阈值化而不影响图像的重要特征.阈值法是一种简单的非线性方法,每次只对一个小波系数进行处理.在最基础的形式中,每个系数都是通过与一个阈值进行比较来设定阈值的.如果系数小于阈值,则将其设为零;否则,它将保持原样或进行修改.对阈值结果进行小波零变换和小波逆变换,可以得到噪声较小的信号基本特征.在这项工作中,分别使用了硬阈值和软阈值的方法进行去噪,实验表明软阈值去噪的效果更好.
在硬阈值法中,将每个小波系数的绝对值与阈值h进行比较.如果小于或等于h,则将其设为零,因为它们大多对应噪声.边缘相关系数通常在阈值以上,因此被保留.阈值的选择方法有很多种,一种常用的方法是计算每个子频带的标准差,将阈值设置为相应子频带中标准差(o)的倍数,因此h的值为2o.硬阈值表示为:
其中,dik表示小波系数,h为阈值.执行硬阈值化的函数如图3(a)所示.
软阈值通常表示为:
其中,sign(dik)表示signum函数.执行软阈值化的函数如图3(b)所示.
图3 阈值函数
基于平稳小波变换后处理的图像存在图像模糊的特征,为了提高图像质量,改进图像分析方法的性能,文献[10]提出了一种用于超声图像去斑点的非局部均值(NLM)滤波器,文献[11]将贝叶斯框架应用到NLM滤波中,实现了图像斑点噪声的消除.因此,本文采用OBNLM算法对图像进行进一步滤波.与传统公式相比,优化贝叶斯非局部均值(OBNLM)滤波器采用了分块方法和Pearson距离计算权重.
NLM滤波的基本公式[14]表示为:
其中,NL(u)(xi)是像素xi的恢复强度和像素u(xi)在“搜索量”Vi的加权平均值,w(xi,xj)是像素恢复xi指定的强度值u(xj)的权重.
根据文献[15]中NLM过滤器的贝叶斯解释,将图像划分为大小为P=(2α+1)d的块Bik(d是图像维度),考虑到块Bik之间的交叉点不为空(即2α≥n),分块NLM定义为:
其中:p(u(Bik)|u(Bj)表示u(Bik)|u(Bj)的分布,p(u(Bj))表示斑块分布.
基于贝叶斯公式,文献[11]提出利用噪声分布计算图像中斑块间距离的新方法:
这种距离允许平滑明亮区域,而不是黑暗区域.在搜索位置选择像素有助于提高过滤器的速度,并更好地保持对比度.
为了验证本文提出去噪算法效果,采用来自数据集Set12的两幅灰度图像(如图4的(a)、(b)所示,Set12的两幅灰度图像house和lena,大小分别为256×256和512×512)和三幅真实的SAR图像(如图4的(c)(d)(e)所示三幅真实的SAR图像,大小分别为256×256、256×256和500×500)进行比较.在两幅灰度图像中添加乘性噪声,去噪图像采用全参考度量(PSNR、FSIM、SSIM、EPF)对图像的去噪效果进行评估;对三幅真实的SAR图像采用非参考度量(ENL)对本文算法进行评价.所有的实验都在Intel Core i5处理器、16 GB RAM和64位操作系统Windows 7电脑的MATLAB 2020a中实现.
图4 实验采用的2幅灰度图像和3幅SAR图像
用于验证的评估指标有:
(1)峰值信噪比:
其中d为大小为M×N[16]的原始图像与重建图像的变化量.峰值信噪比值越高,表明图像中的噪声水平越低.
(2)结构相似性指数:
其中:L、C、S分别为亮度、对比度、结构变化量[17].该指数越接近1,相似性越高.
(3)特征相似度指数:
其中:S(x)为x位置的局部相似度,PCm(x)为x位置f1和f2之间的最大相位同余值[18].该指数越接近1,相似性越高.
(4)边缘保持因子:
其中:Δx和Δy分别为x和y的边缘滤波图像和为它们的均值[19].接近1的值表示边缘保持最佳.
(5)等效视数:
其中:μf为均值,σf为均匀区域f[20]的标准差.等效视数的值越高,说明区域越平滑,即散斑的出现越少.
在house图像和lena图像中添加噪声等级L=4的乘性噪声,其后用SFOBNLM去噪算法与经典的去噪算法[21-22](如Lee滤波、Frost滤波、Mean滤波、FABF、NLM、FNLM和ANLM等),对加噪的图像进行去噪处理,图5和图6分别展示了经典方法与本文方法去噪后的效果.明显可见,采用FABF算法、Lee滤波、Frost滤波和Mean滤波去噪后的图像,虽然散斑噪声较加噪图像有所减少,但整体图像受到噪声的影响仍然较大,图像的纹理结构、局部轮廓和背景信息均较为模糊;采用NLM、FNLM和ANLM算法去噪后的图像,散斑噪声明显减少,图像的纹理结构和局部轮廓较为清晰,但背景信息仍受噪声影响较为粗糙;采用SFOBNLM算法去噪后的图像,无明显的散斑噪声,具有良好的纹理结构、局部轮廓保持能力强,背景信息清晰,说明本文的去噪效果和细节保持都明显优于其他算法.采用表1、表2的全参考度量评价指标,进一步对去噪结果进行评估,对比表明本文算法的PSNR指数最高,SSIM、FSIM、EPF更接近1,去噪效果更好.
表1 经典算法与SFOBNLM算法对house图像去噪效果的全参考度量评价
表2 经典算法与SFOBNLM算法对lena图像去噪效果的全参考度量评价
图5 house加噪图像和去噪结果图像
图6 lena加噪图像和去噪结果图像
针对三幅SAR图像生成特征,由于没有不带散斑图像的真值图,因此用等效视数ENL评价图像的平滑效果,其值越高,区域越平滑,去噪效果越好.如图7的所示,蓝色方框中的像素被用来估计ENL.表3为各算法滤波后图像的ENL结果,从表3看出,不同的SAR图像通过本文算法滤波,都可以得到较大的ENL值,表明本文算法在平滑SAR图像上有较好的结果.
图7 三幅SAR图像(蓝色方框内的像素被用以估计ENL)
表3 三幅SAR图像不同区域不同算法的ENL估计
本文基于平稳小波变换,对噪声图像进行自适应双边滤波和小波阈值去噪,并结合优化贝叶斯非局部均值算法对图像进一步去噪,提出了一种基于快速自适应双边滤波和OBNLM的SAR图像去噪算法,并将本文算法与经典的乘性噪声去噪算法进行了比较.通过视觉感受和评价指标的数据呈现,清楚表明不论是加噪的灰度图像,还是真实的SAR图像,采用本文算法都可以达到较好的去噪效果.因此本文提出的算法具有良好的散斑抑制能力和细节保留能力,是一种有效的SAR图像散斑噪声抑制算法.