小波变换的matlab实现.ppt
1,第4章 小波变换的matlab实现,2,1.Matlab中小波种类,15种经典类小波:Harr小波、Morlet小波、Mexican hat小波、Gaussian小波正交小波:db小波、对称小波、Coiflets小波、Meyer小波双正交小波查看命令 wavemngr(read,1),3,小波分析示例,一维连续小波 1.coefs=cwt(s,scale,wname)2.coefs=cwt(s,scale,wname,plot)c=cwt(noissin,1:48,db4,plot);,4,C=cwt(noissin,2:2:128,db4,plot),5,图形接口方式(GUI)命令:wavemenu,6,7,8,一维离散小波分解,命令:dwt格式:cA1,cD1=dwt(X,wname)cA1,cD1=dwt(X,Lo_D,Hi_D),举例:load leleccum;s=leleccum(1:3920);ls=length(s);cA1,cD1=dwt(s,db1);,9,原始信号,低频系数,高频系数,10,系数重构,命令:upcoef格式:1.Y=upcoef(O,X,wname,N)2.Y=upcoef(O,X,wname,N,L)3.Y=upcoef(O,X,Lo_R,Hi_R,N)4.Y=upcoef(O,X,Lo_R,Hi_R,N,L)5.Y=upcoef(O,X,wname)6.Y=upcoef(O,X,Lo_R,Hi_R)O=a 低频,O=d 高频,11,举例:A1=upcoef(a,cA1,db1,1,ls);D1=upcoef(d,cD1,db1,1,ls);subplot(1,2,1);plot(A1);title(Approximation A1)subplot(1,2,2);plot(D1);title(Detail D1),12,逆变换恢复信号,命令:idwt格式:1.X=idwt(cA,cD,wname)2.X=idwt(cA,cD,Lo_R,Hi_R)3.X=idwt(cA,cD,wname,L)4.X=idwt(cA,cD,Lo_R,Hi_R,L),13,举例:A0=idwt(cA1,cD1,db1,ls);,14,多尺度一维分解,命令:wavedec格式:C,L=wavedec(X,N,wname)C,L=wavedec(X,N,Lo_D,Hi_D),15,C,L=wavedec(s,3,db1);,16,低频系数提取,命令:appcoef格式:1.A=appcoef(C,L,wname,N)2.A=appcoef(C,L,wname)3.A=appcoef(C,L,Lo_R,Hi_R,N)4.A=appcoef(C,L,Lo_R,Hi_R),17,高频系数提取,命令:detcoef格式:1.A=detcoef(C,L,N)2.A=detcoef(C,L),18,举例 cA3=appcoef(C,L,db1,3);cD3=detcoef(C,L,3);cD2=detcoef(C,L,2);cD1=detcoef(C,L,1);,19,重构系数,命令:wrcoef格式:1.X=wrcoef(type,C,L,wname,N)2.X=wrcoef(type,C,L,Lo_R,Hi_R,N)3.X=wrcoef(type,C,L,wname)4.X=wrcoef(type,C,L,Lo_R,Hi_R)type=a 低频,type=d 高频,20,A3=wrcoef(a,C,L,db1,3);D1=wrcoef(d,C,L,db1,1);D2=wrcoef(d,C,L,db1,2);D3=wrcoef(d,C,L,db1,3);,21,重构原始信号,命令:waverec格式:1.X=waverec(C,L,wname)2.X=waverec(C,L,Lo_R,Hi_R)例子:A0=waverec(C,L,db1);重构最大误差:Err=max(abs(s-A0),22,23,图形接口方式(GUI),24,25,26,27,28,29,2.二维离散小波,单尺度分解dwt2格式:1.cA1,cH1,cV1,cD1=dwt2(X,wname)2.cA1,cH1,cV1,cD1=dwt2(X,Lo_D,Hi_D)cA1,cH1水平;cV1垂直;cD1对角应用:load wbarb;figure(1);image(X);colormap(map);colorbar;cA1,cH1,cV1,cD1=dwt2(X,bior3.7),30,重构系数,命令:upcoef2格式:1.Y=upcoef2(O,X,wname,N,S)2.Y=upcoef2(O,X,Lo_R,Hi_R,N,S)3.Y=upcoef2(O,X,wname,N)4.Y=upcoef2(O,X,Lo_R,Hi_R,N)5.Y=upcoef2(O,X,wname)6.Y=upcoef2(O,X,Lo_R,Hi_R)O:a低频;h水平;v垂直;d对角,31,A1=upcoef2(a,cA1,bior3.7,1);H1=upcoef2(h,cH1,bior3.7,1);V1=upcoef2(v,cV1,bior3.7,1);D1=upcoef2(d,cD1,bior3.7,1);figure(2);colormap(map);subplot(2,2,1);image(wcodemat(A1,192);title(Approximation A1)subplot(2,2,2);image(wcodemat(H1,192);title(Horizontal Detail H1)subplot(2,2,3);image(wcodemat(V1,192);title(Vertical Detail V1)subplot(2,2,4);image(wcodemat(D1,192);title(Diagonal Detail D1),32,33,二维逆变换,命令:idwt2格式:1.X=idwt2(cA1,cH1,cV1,cD1,bior3.7);2.X=idwt2(cA1,cH1,cV1,cD1,bior3.7);3.X=idwt2(cA1,cH1,cV1,cD1,bior3.7);4.X=idwt2(cA1,cH1,cV1,cD1,bior3.7);应用:Xsyn=idwt2(cA1,cH1,cV1,cD1,bior3.7);,34,多尺度二维小波,命令:wavedec2格式:1.C,S=wavedec2(X,N,wname)2.C,S=wavedec2(X,N,Lo_D,Hi_D),35,C,S=wavedec2(X,2,bior3.7);%图像的多尺度二维小波分解,36,提取低频系数,命令:appcoef2格式:1.A=appcoef2(C,S,wname,N)2.A=appcoef2(C,S,wname)3.A=appcoef2(C,S,Lo_R,Hi_R)4.A=appcoef2(C,S,Lo_R,Hi_R,N)cA2=appcoef2(C,S,bior3.7,2);%从上面的C中提取第二层的低频系数,37,提取高频系数,命令:detcoef2格式:A=detcoef2(type,C,S,wname,N)说明:Type:h 水平;v垂直;d对角 cH2=detcoef2(h,C,S,2);cV2=detcoef2(v,C,S,2);cD2=detcoef2(d,C,S,2);cH1=detcoef2(h,C,S,1);cV1=detcoef2(v,C,S,1);cD1=detcoef2(d,C,S,1);,38,重构系数,命令:wrcoef2格式:1.X=wrcoef2(type,C,S,wname,N)2.X=wrcoef2(type,C,S,Lo_R,Hi_R,N)3.X=wrcoef2(type,C,S,wname)4.X=wrcoef2(type,C,S,Lo_R,Hi_R,N)说明:Type:a低频;h 水平;v垂直;d对角,39,A2=wrcoef2(a,C,S,bior3.7,2);H1=wrcoef2(h,C,S,bior3.7,1);%重构第1、2层的高频信号 V1=wrcoef2(v,C,S,bior3.7,1);D1=wrcoef2(d,C,S,bior3.7,1);H2=wrcoef2(h,C,S,bior3.7,2);V2=wrcoef2(v,C,S,bior3.7,2);D2=wrcoef2(d,C,S,bior3.7,2);,40,41,重构原始信号,命令:waverec2格式:X=waverec2(C,S,wname)X=waverec2(C,S,Lo_R,Hi_R)应用:X0=waverec2(C,S,bior3.7);,42,2D图形接口,43,显示,44,小波分析用于信号处理,常用信号的小波分析信号的特征提取信号处理GUI进行信号处理,45,正弦波的线性组合,S(t)=sin(2t)+sin(20t)+sin(200t),46,间断点检测波形未来预测各分信号的频率识别信号从近似到细节的迁移,47,分段信号,S(t)=sin(0.03t)t=1:500 或sin(0.3t)t=500:1000,信号抑制信号未来预测,48,信号的特征提取,信号的突变点往往是它的重要特征信号的频率谱和它的幅值等表征了信号的许多信息。信号的连续性(即信号的奇异性)分析、信号的频率谱分析和幅值谱分析不可或缺。小波分析进行特征提取时,两种处理方法,即边界的处理和滤波。利用小波分析得到低频和高频部分。,49,检测信号的突变点,提取了信号的近似特征a和细节特征d。在原始信号图像上,无法得知原始信号导数的不连续性。,50,信号的奇异点检测,【定义】在某一尺度 下,如果存在一点 使得,则称点 是局部极值点,且 在 上有一个模极大值(过零)点。如果对 的某一领域内的任意点,则称 为小波变换模极大值(过零)点。尺度空间中所有的模极大值点的连线称为模极大值线。,51,定理:设 为一严格的整数,为具有 阶消失矩、次连续可微和紧支集的小波,(为某一实数区间),若存在尺度,使得,没有局部极大值点,则在区间 上是一致Lipschitz(为任一小的正数)。一般来讲,函数在某一点的Lipschitz指数 表征了该点的奇异性大小,越大,该点的光滑度越高,越小,该点的奇异性越大。,52,当小波函数可看做某一平滑函数的一阶导数时,信号小波变换模的局部极值点对应于信号的突变点(或边缘);当小波函数可看做某一平滑函数的二阶导数时,信号小波变换模的过零点,也对应于信号的突变点(或边缘)。因此,采用检测小波变换系数模的过零点和局部极值点的方法可以检测信号的边缘位置。比较而言,采用局部边缘进行检测更具有优越性。,53,信号的奇异性通常可以分为两种情况:第一种类型的间断点:信号在某一时刻,其幅值发生突变,引起信号的不连续,信号的突变处是间断点;第二种类型的间断点:信号外观上很光滑,其幅值没有突变,但在信号的一阶微分上有突变产生,且一阶微分是不连续的。,54,55,56,信号自相似性的检测,小波系数与自相似性:小波分解可通过计算信号和小波的“自相似指数”得到。自相似指数也就是小波系数,如果自相似指数大,则信号的自相似程度就高,反之亦然。如果一个信号在不同的尺度上都相似于它本身,那么,其“自相似指数”,或者小波系数在不同的尺度上也是相似的。,57,58,信号发展趋势的识别,59,在某一频率区间上信号的识别,60,信号抑制与衰减,消失矩:如果(其中 是小波函数)()的平均值为0,那么该小波有n+1个消失矩,并且可以利用该小波对n次多项式信号进行抑制。,61,62,信号消噪与提取弱信号,其中,为含噪信号,为有用信号,为噪声信号。,消噪的三个步骤:1.一维信号的小波分解2.小波分解高频系数的阈值量化3.一维小波重构,63,信号消噪处理,命令:wden格式:1.sd=wden(s,tptr,sorh,scal,n,wavename)说明:2.tptr指定阈值选取规则;3.sorh指定选取软阈值(sorh=s)或硬阈值(sorh=h)4.scal=one 基本模式 scal=sln 未知尺度的基本模式,且仅根据第一层的小波分解系数来估计噪声的层次,并只进行一次估计,以此来变换阈值的尺度。scal=mln 非白噪声的基本模式,且在每个不同的小波分解层次上都估计噪声的层次,以此来变换阈值的尺度。,64,命令:wdencmp格式:xd=wdencmp(opt,x,wavename,n,thr,sorh,keepapp)其中:(1)opt=gbl,thr0,则阈值为全局阈值 opt=lvd,thr是向量,则阈值是在各层上大小不同的数值。(2)keepapp=1,不对小波分解后的低频系数做处理 keepapp=0,对小波分解后的低频系数也进行阈值化处理(3)x是待处理信号。(4)xd是处理后信号。(5)其他同wden参数,65,小波分析进行消噪处理3种方法,默认阈值消噪处理。该法利用函数ddencmp生成信号的默认阈值,然后利用函数wdencmp进行消噪处理。给定阈值消噪处理。在实际的消噪处理中,阈值往往通过经验公式获得,且这种阈值比默认阈值的可信度高。在进行阈值量化处理时可用函数wthresh。强制消噪处理。该法是将小波分解结构中的高频系数全置0,即滤掉所有高频部分,然后进行小波重构。这种方法比较简单,且消噪后的信号比较平滑,但容易丢失信号中的有用成分。,66,67,68,阈值选取规则,命令:wthresh格式:yt=wthresh(y,sorh,thr)说明:该函数根据参数sorh的取值返回输入分解系数的软阈值或硬阈值。硬阈值对应于最简单的处理方法,而软阈值具有很好的数学特性,并且所得到的理论结果是可用的。,69,70,71,信号压缩,步骤1.信号的小波分解2.对高频系数进行阈值量化处理。对第一到第N层的高频系数,均可选择不同的阈值,并用硬阈值进行系数的量化。对量化的系数进行小波重构。压缩与消噪主要区别:第2步。有效的信号压缩方法:1.对信号进行小波尺度的扩展,并保留绝对值最大的系数;2.根据分解后各层的效果来确定某一层的阈值,且这些阈值是互不相同。,72,73,小波分析图像处理,图像的小波分解与重构,H,L,图像小波分解示意图,74,第1级 L1斜线细节,第1级 L1水平细节,第1级 L1垂直细节,第2级 L2细节,近似图象,第3级 L3,75,小波分解数据流示意图,76,小波重构数据流示意图,77,78,79,80,图像边缘失真的处理,补零:假设在原始支撑以外的信号以零补足。其缺点是人为地在边界处制造了不连续。边界对称化:假设通过对称的边界值复制可以恢复信号和图像的原始支撑以外的信号和图像。缺点是在边界处,人为地制造了一阶导数的不连续性,但该方法通常对图像处理非常有效。一阶平滑填补:假设通过简单的一阶导数外插(填补时,对前两个值和后两个值使用线性拟合)能够从信号或图像的原始支撑之外恢复信号或图像。该方法通常对光滑信号较为有效。,81,零阶平滑填补:假设通过简单的常数外插能够从信号或图像的原始支撑之外恢复信号或图像。对于信号延拓来说,该方法是位于左边的第一个值和右边的最后一个值的重复。周期性填补1:这里假设通过周期延拓恢复信号或图像原始支撑以外的信号或图像。其缺点是在边界处人为地制造了不连续性。周期性填补2:如果信号的长度是奇数,首先给信号加一个采样点,其值等于最后一个值,接下来对信号进行周期性填补(1),即在两端对信号进行最小周期延拓。对于图像,采取同样的方法。这种模式可生成最小长度的小波分解,但为了确保完全重构,在逆变换中也应采用同样的延拓模式。总结:前5种,会存在一定的冗余,因此在任意一种模式中,在逆变换中都能确保对信号和图像的完全重构。,82,补零,对称,平滑填补,83,图像消噪的步骤,二维图像信号的小波分解 选择合适的小波和恰当的分解层次(记为N),然后对待分析的二维图像信号X进行N层分解计算。对分解后的高频系数进行阈值量化。对于分解的每一层,选择一个恰当的阈值,并对该层高频系数进行软阈值量化处理。二维小波的重构图像信号根据小波分解后的第N层近似(低频系数)和经过阈值量化处理后的各层细节(高频系数),来计算二维信号的小波重构。,84,85,86,图像压缩,87,小波变换图像压缩,88,对图像小波分解后,可得到一系列不同分辨率的子图像(它们所对应的频率不相同)。而对于图像来说,表征它的最主要的部分是低频部分,而高频部分大部分点的数值均接近于0,而且频率越高,这种现象越明显。因此,利用小波分解去掉图像的高频部分而仅保留图像的低频部分是一种最简单的图像压缩方法。,89,90,压缩前图像的大小:Name Size Bytes Class Attributes X 256x256 524288 double 第一次压缩图像的大小:Name Size Bytes Class Attributes ca1 135x135 145800 double 第二次压缩图像大小:Name Size Bytes Class Attributes ca2 75x75 45000 double,91,%装载并显示原始图像load flujet;subplot(1,2,1);image(X);colormap(map);title(原始图像);axis square;%=%首先利用小波db3对图像X进行2层分解c,l=wavedec2(X,2,db3);%=%全局阈值thr,sorh,keepapp=ddencmp(cmp,wv,X);%=%进行压缩处理:对所有高频系数进行同样的阈值量化处理Xcmp,cxc,lxc,perf0,perfl2=wdencmp(gbl,c,l,db3,2,thr,sorh,keepapp);%=%图示压缩结果subplot(1,2,2);image(Xcmp);colormap(map);title(压缩后的图像);axis square;disp(小波分解系数中置0的系数个数百分比:);perf0disp(压缩后保留能量百分比:);perfl2,92,小波分解系数中置0的系数个数百分比:perf0=89.5226压缩后保留能量百分比:perfl2=99.9831,93,Wdencmp函数,%装载并显示原始图像load wbarb;subplot(1,2,1);image(X);colormap(map);title(原始图像);%=%采用默认的全局阈值对图像进行压缩thr,sorh,keepapp,crit=ddencmp(cmp,wp,X);Xc=wpdencmp(X,sorh,3,bior3.1,crit,thr,keepapp);subplot(1,2,2);image(Xc);colormap(map);title(全局阈值压缩图像);,