小波阈值图像降噪及MATLAB仿真

小波阈值图像降噪及MATLAB仿真图像信号在生成和传输过程中常常因受到各种噪声的干扰和影响而使图像降质，这对后续图像的处理(如分割、压缩和图像理解等)将产生不利影响。为了抑制噪声，改善图像质量，便于更高层次的处理，必须对图像进行降噪预处理。小波降噪的方法有多种，如利用小波分解与重构的方法滤波降噪、利用小波变换模极大值的方法去噪、利用信号小波变换后空域相关性进行信噪别离、非线性小波阈值方法去噪、平移不变量小波降噪法，以及多小波降噪等等。归结起来主要有三类:模极大值检测法、阈值降噪法和屏蔽(相关)降噪法。其中最常用的就是阈值法去噪，其根本思想就是利用图像小波分解后，各个子带图像的不同特性选取不同的阈值，从而到达较好的降噪目的。而且，小波变换本身是一种线形变换，因而对于类似于高斯噪声的效果较好。线性运算往往还会造成边缘模糊，小波分析技术正因其独特的时频局部化特性在图像信号和噪声信号的区分以及有效去除噪声并保存有用信息等方面较之传统的降噪具有明显的优势，且在降噪的同时实现了图像一定程度的压缩和边缘特征的提取。所以小波降噪具有无可比较的优越性。小波降噪主要优点有：低熵性，小波系数的稀疏分布，使得图象变换后的熵降低；多分辨率，由于采用了多分辨率的方法，所以可以非常好地刻画信号的非平稳特征，如边缘、尖峰、断点等；去相关性，因为小波变换可以对信号进行去相关，且噪声在变换后有白化趋势，所以小波频域比时域更利于降噪；选基灵活性，由于小波变换可以灵活选择变换基，从而对不同应用场合、不同的研究对象，可以选用不同的小波函数，以获得最正确的效果。一、阈值去噪法简述1992年，斯坦福大学的DonohoDL和Johnstone教授提出一种具有良好的统计优化特性的降噪方法，称作“WaveletShrinkage〞〔即阈值收缩法〕。该方法的主要思想是：基于图像和噪声在经小波变换后具有不同的统计特性：图像本身的能量对应着幅值较大的小波系数,主要集中在高频〔〕；噪声能量那么对应着幅值较小的小波系数，并分散在小波变换后的所有系数中。根据该特征，设置一个阈值门限，认为大于该阈值的小波系数的主要成份为有用的信号，给予收缩后保存；小于该阈值的小波系数，主要成份为噪声，予以剔除，这样就可以到达降噪的目的。降噪时，通常认为低通系数含有大量的图像能量，一般不作处理，只对剩余三个高通局部进行处理。因此，一次阈值降噪并不能完全去除噪声，还需要对未作处理的低频局部〔〕再次进行小波分解和阈值去噪，直到实际图像与估计图像的偏差到达最小值。但是，随着分解和降噪次数的增加，小波系数中的噪声能量越来越少，并且趋于分散，降噪的效果将逐渐降低。一般来说，进行3-4层小波分解和降噪就可以到达满意的降噪效果。小波阈值降噪的根本思路是：〔1〕先对含噪信号做小波变换，得到一组小波系数；〔2〕通过对进行阈值处理，得到估计系数，使得与两者的差值尽可能小；〔3〕利用进行小波重构，得到估计信号即为降噪后的信号。Donoho提出了一种非常简洁的方法对小波系数进行估计。对连续做几次小波分解后，有空间分布不均匀信号各尺度上小波系数在某些特定位置有较大的值，这些点对应于原始信号的奇变位置和重要信息，而其他大局部位置的较小；对于白噪声，它对应的小波系数在每个尺度上的分不都是均匀的，并随尺度的增加，系数的幅值减小。因此，通常的降噪方法是寻找一个适宜的数作为阈值〔门限〕，把低于λ的小波函数〔主要由信号引起〕，设为零，而对于高于的小波函数〔主要由信号引起〕，那么予以保存或进行收缩，从而得到估计小波系数，它可理解为根本由信号引起的，然后对进行重构，就可以重构原始信号。估计小波系数的方法如下，取：〔4-1〕定义：〔4-2〕称之为硬阈值估计方法。一般软阈值估计定义为〔4-3〕两种阈值方法各有差异,前者具有连续性,在数学上易于处理,而后者更接近实际情况。阈值化处理的关键是阈值的选择,如阈值太小,降噪后仍留有噪声,但阈值如果太大,重要的信号与图像特征会被滤掉,引起偏差。二、基于MATLAB的小波降噪函数简介常用的图像降噪方式是小波阈值降噪方法。这是一种实现简单而效果较好的降噪方法，阈值降噪方法的思想很简单，就是对小波分解后的各层系数模大于和小于某阈值的系数分别进行处理，然后利用处理后的小波系数重构出降噪后的图像。在阈值降噪中，阈值函数表达了对小波分解系数的不同处理策略和不同的估计方法。常用的阈值函数有硬阈值函数和软阈值函数。硬阈值函数可以很好地保存图像边缘等局部特征，但图像会出现伪吉布斯效应等视觉失真等现象；而软阈值处理相对较光滑，但可能会造成边缘模糊等失真现象，为此人们提出了半软阈值函数。小波阈值降噪方法处理阈值的选取，另一个关键因素是阈值的具体估计。如果阈值太小，降噪后的图像仍然存在噪声；相反如果阈值太大，重要图像特征有被滤掉，引起偏差。从直观上讲,对于给定的小波系数，噪声越大，阈值就越大。MATLAB中实现图像的降噪，主要是阈值获取和图像降噪实现两个方面。1阈值获取MATLAB中实现阈值获取的函数有ddencmp、select、wbmpen和wdcbm2。这里主要介绍函数ddencmp。函数ddencmp的功能是获取降噪或压缩的默认值。该函数是降噪和压缩的导向函数，它给出一维或二维信号使用小波或小波包进行降噪和压缩一般过程的所有默认值。其语法格式为：[THR,SORH,KEEPAPP,CRIT]=ddencmp〔IN1，IN2，X〕[THR,SORH,KEEPAPP]=ddencmp〔IN1，‘wv’，X〕[THR,SORH,KEEPAPP,CRIT]=ddencmp〔IN1，‘wp’，X〕2阈值降噪MATLAB中实现阈值降噪的函数有wden、wdencmp、wpdencmp、wthresh、wpthcoef和wthcoef2。这里主要介绍函数wdencmp。其语法格式为：[XC,CXC,LXC,PERF0,PERFL2]=wdencmp('gbl',X,'wname',N,THR,SORH,KEEPAPP)[XC,CXC,LXC,PERF0,PERFL2]=wdencmp('lvd',X,'wname',N,THR,SORH)[XC,CXC,LXC,PERF0,PERFL2]=wdencmp('lvd',C,L,'wname',N,THR,SORH)函数wdencmp的功能是使用小波进行降噪。该函数是二维小波降噪的导向函数。它使用小波对信号或图像执行降噪过程。wname是所用的小波函数。gbl〔global〕表示每层都采用同一个阈值进行处理。lvd表示每层用不同的阈值进行处理。N表示小波分解的层数。THR为阈值向量，长度为N。SORH表示选择软阈值或硬阈值〔分别取值为‘s’和‘h’〕。参数KEEPAPP取值为1是，那么低频系数不进行阈值量化，反之，那么低频系数要进行阈值量化。XC是降噪后的信号，[CXC,LXC]是XC的小波分解结构，PHRF0和PERFL2是恢复和压缩L2的范数百分比。如果[C，L]是x的小波分解结构，那么PERFL2=100〔CXC向量的范数/C向量的范数〕2；如果X是一维信号，小波wname是一个正交小波，那么PERFL2=。三、用小波对信号进行降噪使用小波进行一维信号的降噪对含噪一维信号降噪，实现信噪别离。含噪的正弦信号：对含噪信号进行分解，从分解中得出高频与低频系数。高频系数对应细节信号，低频系数对应逼近信号。由高频系数中检测噪声，低频系数中识别各分量信号的不同频率。然后，再从去噪的高频信号与逼近的低频信号重构目标信号。降噪的程序如下：%生成含噪正弦信号N=1024;t=1:N;sig=sin(0.03*t);figure(1);subplot(211);plot(t,sig);title('正弦信号');%叠加噪声x=sig+randn(1,N);subplot(212);plot(t,x);title('含噪正弦信号');%一维小波分解，使用'haar'进行4层分解[c,l]=wavedec(x,4,'haar');%重构第1-4层逼近信号a4=wrcoef('a',c,l,'haar',4);a3=wrcoef('a',c,l,'haar',3);a2=wrcoef('a',c,l,'haar',2);a1=wrcoef('a',c,l,'haar',1);%显示各层逼近信号figure(2);subplot(411);plot(a4);ylabel('a4');subplot(412);plot(a3);ylabel('a3');subplot(413);plot(a2);ylabel('a2');subplot(414);plot(a1);ylabel('a1');%重构第1-4层细节信号d4=wrcoef('d',c,l,'haar',4);d3=wrcoef('d',c,l,'haar',3);d2=wrcoef('d',c,l,'haar',2);d1=wrcoef('d',c,l,'haar',1);%显示各层细节信号figure(3);subplot(411);plot(d4);ylabel('d4');subplot(412);plot(d3);ylabel('d3');subplot(413);plot(d2);ylabel('d2');subplot(414);plot(d1);ylabel('d1');2、使用小波对图像的降噪处理对含有高斯白噪声的图像，使用sym4小波进行分解，再用阈值法进行降噪处理，最后重构图像，实现降噪目的。在下面的例子中，分别采用了软阈值和硬阈值两种方法，观察降噪后的图像并比较两种方法的优缺点。程序如下：%装载原始图像I=imread('D:\MATLAB7\toolbox\images\imdemos\cameraman.tif');nbc=size(I,1);X=im2double(I);%产生噪声图像init=2055415866;randn('seed',init);x=X+randn(size(X))/10;%使用sym4执行图像的2层小波分解wname='sym4';lev=2;[c,l]=wavedec2(x,lev,wname);sigma_s=0.054779;%图像降噪时，使用wbmpen函数选择阈值alpha=2;thr_s=wbmpen(c,l,sigma_s,alpha);%使用软阈值和保存的低频信号，进行图像降噪keepapp=1;xds=wdencmp('gbl',x,wname,lev,thr_s,'s',keepapp);sigma_h=0.062818;thr_h=wbmpen(c,l,sigma_h,alpha);%使用硬阈值和保存的低频信号，进行图像降噪xdh=wdencmp('gbl',x,wname,lev,thr_h,'h',keepapp);%画出原始图像和降噪后的图像figure(1);subplot(221);imshow(I,nbc)