连续小波变换CWT以及MATALB例程

连续小波变换（CWT）以及MATALB例程2.1连续小波变换及其性质2.1.1连续小波基函数小波，即小区域的波，是一种特殊的长度有限、平均值为零的波形。小波的可容许条件：RC|||)(|2^小波特点：（一）“小”。即在时域都具有紧支集或近似紧支集。（二）正负交替的“波动性”。即直流分量为零。信号可分解为一系列由同一个母小波函数经平移与尺度伸缩得到的小波函数的叠加。将小波母函数进行伸缩和平移，就可以得到函数：小波函数基，它们是由同一母函数经伸缩和平移后得到的一组函数序列。)(t0;,),(1)(,aRaatata)(t伸缩和平移的含义1.尺度伸缩2.时间平移由于小波基函数在时间、频率域都具有有限或近似有限的定义域，显然，经过伸缩平移后的函数在时、频域仍是局部性的。小波基函数的窗口随尺度因子的不同而伸缩，当a逐渐增大时，基函数的时间窗口也逐渐增大，而其对应的频域窗口逐渐减小；反之亦然。Haar小波101/2()11/210ttt其它/224ˆ()sin/4iie1112()t01定量分析-时域假定小波母函数窗口宽度为△t，窗口中心为t0，则相应可求出连续小波的窗口中心为at0+τ，窗口宽度为a·△t。即信号限制在时间窗内：[at0+τ-△t·a/2,at0+τ+△t·a/2])(1)(,atata定量分析-频域同样，对于小波母函数的频域变换，其频域窗口中心为ω0，窗口宽度为△ω，则相应的连续小波的傅立叶变换为：其频域窗口中心为：窗口宽度为：信号在频域窗内：)()(21,aeaja0,1aaa1]211,211[00aaaa从上面的时频域的讨论可见，连续小波的时频域窗口中心及其宽度都随a的变化而伸缩，如果我们称△t·△ω为窗口函数的窗口面积，则：可见：连续小波基函数的窗口面积不随参数的变化而变化。atataa1,,几点结论：（1）尺度的倒数1/a在一定意义上对应于频率ω。即尺度越小，对应的频率越高。如果我们将尺度理解为时间窗口的话，则小尺度信号为短时间信号，大尺度信号为长时间信号。（2）在任何τ值上，小波的时频窗口大小△t和△ω都随频率ω（或a）的变化而变化。与短时傅立叶变换中的基不同。tjetgtg)()(,（3）在任何尺度a，时间点τ上，窗口面积保持不变，也可以说时间、尺度分辨率是相互制约的，不可能同时得到提高。（4）品质因素不随尺度变化而变化。0Q“恒Q性质”：假设(t)的中心为t0，有效宽度为Dt；()的中心为0，有效宽度为D；则a,b(t)提取的是f(t)在窗口[b+at0-aDt/2,b+at0+aDt/2]|中的性质，相应地从频域上说a,b()提取地是F()在窗口[0/a-D/(2a),0/a+D/(2a)]中的性质，因此对于小波来说时域窗口宽度和频域窗口宽度的乘积始终为DtD。2.1.2连续小波变换的定义和性质1.连续小波变换的定义将任意L2（R）空间中的函数f(t)在小波基下展开，称这种展开为函数f(t)的连续小波变换（CWT）。其表达式为：其中：dtattfattfaWTRaf)()(1)(),(),(,}0{,),(||)(21,RaRbatata从定义可以看出：小波变换和傅立叶变换一样，也是一种变换，为小波变换系数。也可见其与傅立叶变换的区别。),(aWTf逆变换若小波满足容许条件，则连续小波变换存在着逆变换。容许条件：逆变换公式：RdC|||)(|2dataaWTadaCdtaWTadaCtffaf)(1),(1)(),(1)(02,02说明：（1）必须满足“容许条件”，反变换才存在。（2）在实际应用中，对基本小波的要求往往不局限于满足容许条件，对还要施加所谓“正则性条件”，使在频域上表现出较好的局域性能。为了在频域上有较好的局域性，要求随a的减小而迅速减小，所以这就要求的前n阶原点距为0，且n值越高越好。即：)(t)(|),(|aWTf)(t值越大越好。且nnpdtttp,~1,0)(3.连续小波变换的再生核尺度和位移的连续变化的连续小波基函数构成了一组非正交的过渡完全基，小波展开系数之间有相关关系，采用如下描述RaadtttCaaK)()(1),;,(,,1.CWT系数具有很大的冗余，计算量比较大2.利用冗余性可以实现去噪和数据恢复的目的。重建核方程daaKaWadaaWff),;,(),(),(0002004.连续小波变换具有以下重要性质：(1)线性性：一个多分量信号的小波变换等于各个分量的小波变换之和。(2)平移不变性：若f(t)的小波变换为Wf(a，b)，则f(t－)的小波变换为Wf(a，b－)(3)伸缩共变性：若f(t)的小波变换为Wf(a，b)，则f(ct)的小波变换为),(1cbcaWcf0c(4)自相似性：对应不同尺度参数a和不同平移参数b的连续小波变换之间是自相似的。(5)冗余性：连续小波变换把一维信号变换到二维空间，因此在连续小波变换中存在信息表述的冗余度(redundancy)。小波变换的逆变换公式不是唯一的。小波变换的冗余性事实上也是自相似性的直接反映，它主要表现在以下两个面：①由连续小波变换恢复原信号的重构分式不是唯一的。也就是说，信号f(t)的小波变换与小波重构不存在一一对应关系，而傅里叶变换与傅里叶反变换是一一对应的。②小波变换的核函数即小波函数a，b(t)存在许多可能的选择(例如，它们可以是非正交小波、正交小波、双正交小波，甚至允许是彼此线性相关的)。小波变换在不同的(a，b)之间的相关性增加了分析和解释小波变换结果的困难，因此，小波变换的冗余度应尽可能减小，它是小波分析中的主要问题之一。在MATLAB中，可以用cwt函数实现对信号的连续小波变换。2.2几种常用的小波小波分类的标准（1）的支撑长度，即当时间或频率趋向于无穷大时，它们从一个有限值收敛到0。（2）对称性。它在图像处理中可以有效的避免移相。)()()()(、、、tt（3）的消失矩阵数。这对于压缩非常有用。（4）正则性。它在对信号或图像的重构获得较好的平滑效果作用上是非常有用的。具有对称性的小波不易产生相位畸变；具有好的正则性的小波，易于获得光滑的重构曲线和图像，从而减小误差。)()(tt和常用的小波1.Haar小波。其他,0121,1210,1)(ttt2.Daubechies（dbN）小波12002220110121)()2(sin)2(cos|)(|)(NkjkkNkNNkkkNehmPmCyCyPkk式中，则有：为二项式的系数，，其中，令3.MexicanHat(mexh)小波其函数为Gauss函数的二阶导数：2222222)()1()(eettt4.Morlet小波它是高斯包络下的单频率复正弦函数。是重构时的归一化常数CxCett)5cos()(222.3连续小波变换的步骤（1）选择小波函数及其尺度a值。（2）从信号的起始位置开始，将小波函数和信号进行比较，即计算小波系数。（3）沿时间轴移动小波函数，即改变参数b，在新的位置计算小波系数，直至信号的终点。（4）改变尺度a值，重复（2）、（3）步。2.4尺度和频率之间的关系a为尺度；△为采样间隔；Fc为小波的中心频率；Fa为伪频率。aFFca2.5应用实例例已知一信号f(t)＝3sin(100t)＋2sin(68t)＋5cos(72t)，且该信号混有白噪声，对该信号进行连续小波变换。小波函数取db3，尺度为1、1.2、1.4、1.6、…、3。其MATLAB程序如下：t＝0：0.01：1；f＝3*sin(100*pi*t)＋2*sin(68*pi*t)＋5*cos(72*pi*t)＋randn(1，length(t))；coefs＝cwt(f，［1：0.2：3］，db3，plot)；title(对不同的尺度小波变换系数值)；Ylabel(尺度)；Xlabel(时间)；程序输出结果如图所示。图1.11小波变换的系数用图所示的灰度值图表征，横坐标表示变换系数的系号，纵坐标表示尺度，灰度颜色越深，表示系数的值越大。绘图原理1.需要用到的小波工具箱中的三个函数cwt（），centfrq（），scal2frq（）COEFS=cwt(S,SCALES,‘wname’)说明：该函数能实现连续小波变换，其中S为输入信号，SCALES为尺度，wname为小波名称。FREQ=centfrq(‘wname’)说明：该函数能求出以wname命名的母小波的中心频率。F=scal2frq(A,‘wname’,DELTA)说明：该函数能将尺度转换为实际频率，其中A为尺度，wname为小波名称，DELTA为采样周期。注：这三个函数还有其它格式，具体可参阅matlab的帮助文档。2.尺度与频率之间的关系设a为尺度，fs为采样频率，Fc为小波中心频率，则a对应的实际频率Fa为Fa＝Fc×fs/a(1)显然，根据采样定理，为使小波尺度图的频率范围为(0,fs/2)，尺度范围应为(2*Fc,inf),其中inf表示为无穷大。在实际应用中，只需取尺度足够大即可。3.尺度序列的确定由式(1)可以看出，为使转换后的频率序列是一等差序列，尺度序列必须取为以下形式：c/totalscal,...,c/(totalscal-1),c(2)其中，totalscal是对信号进行小波变换时所用尺度序列的长度(通常需要预先设定好)，c为一常数。下面讲讲c的求法。根据式(1)容易看出，尺度c/totalscal所对应的实际频率应为fs/2，于是可得c=2×Fc×totalscal(3)将式(3)代入式(2)便得到了所需的尺度序列。4.时频图的绘制确定了小波基和尺度后，就可以用cwt求小波系数coefs（系数是复数时要取模），然后用scal2frq将尺度序列转换为实际频率序列f，最后结合时间序列t，用imagesc(t,f,abs(coefs))便能画出小波时频图。注意：直接将尺度序列取为等差序列，例如1:1:64，将只能得到正确的尺度－时间－小波系数图，而无法将其转换为频率－时间－小波系数图。这是因为此时的频率间隔不为常数。。下面给出一实际例子来说明小波时频图的绘制。所取仿真信号是由频率分别为100Hz和200Hz的两个正弦分量所合成的信号。clear;clc;fs=1024;%采样频率f1=100;f2=200;t=0:1/fs:1;s=sin(2*pi*f1*t)+sin(2*pi*f2*t);%两个不同频率正弦信号合成的仿真信号%%%%%%%%%%%%%%%%%小波时频图绘制wavename='cmor3-3';totalscal=256;%尺度序列的长度，即scal的长度wcf=centfrq(wavename);%小波的中心频率cparam=2*wcf*totalscal;%为得到合适的尺度所求出的参数a=totalscal:-1:1;scal=cparam./a;%得到各个尺度，以使转换得到频率序列为等差序列coefs=cwt(s,scal,wavename);%得到小波系数f=scal2frq(scal,wavename,1/fs);%将尺度转换为频率imagesc(t,f,abs(coefs));%绘制色谱图colorbar;xlabel('时间t/s');ylabel('频率f/Hz');title('小波时频图');