基于cwt的时频分析
本例介绍了如何利用连续小波变换(CWT)对信号进行时频联合分析。该实例讨论了瞬态的定位，其中CWT优于短时傅里叶变换(STFT)。实例还展示了如何使用逆CWT合成时频局部化信号逼近。
在一些例子中，将CWT与STFT进行了比较。必须有信号处理工具箱™才能运行光谱图的例子。
调制信号的时频分析
加载一个二次啁啾信号，并显示其频谱图。信号的频率在t= 0时开始于约500hz，在t=2时下降到100hz，在t=4时增加到500hz。采样频率为1khz。
负载quadchirp；Fs = 1000;[S,F,T] =频谱图(四啁啾，100,98,128,fs);helperCWTTimeFreqPlot (S T F,“冲浪”，二次啁啾的STFT，“秒”，“赫兹”）
使用CWT获得该信号的时频图。
[cfs,f] = cwt(quadchirp,fs，“WaveletParameters”[14200]);helperCWTTimeFreqPlot (cfs tquad f,“冲浪”，二次啁啾的CWT，“秒”，“赫兹”）
带有凹凸小波的CWT产生的时频分析非常类似于STFT。
频率和幅度调制在自然信号中经常出现。使用CWT获得一个回声定位脉冲由一个大棕色蝙蝠(Eptesicus Fuscus)发射的时频分析。采样间隔为7微秒。使用每八度32个声音的bump小波。感谢伊利诺伊大学贝克曼中心的Curtis Condon、Ken White和Al Feng提供的蝙蝠数据，并允许在这个例子中使用它。
负载batsignalt = 0:DT:(数字(蝙蝠信号)*DT)-DT;[cfs,f] = cwt(蝙蝠信号，“撞”1 / DT,“VoicesPerOctave”、32);helperCWTTimeFreqPlot (cfs t * 1 e6, f / 1000“冲浪”，“蝙蝠回声定位的CWT”，.．.微秒的，“赫兹”）
获取并绘制蝙蝠数据的STFT。
[S,F,T] =光谱图(蝙蝠信号，50,48,128,1/DT);helperCWTTimeFreqPlot (S, t . * 1 e6, f / 1 e3,“冲浪”，蝙蝠回声定位(STFT)，.．.微秒的，“赫兹”）
对于模拟信号和自然调制信号，CWT提供了类似于STFT的结果。
利用CWT检测振荡中的瞬态
在时频分析中，在某些情况下，CWT可以提供比短时间傅里叶变换更有信息的时频变换。当信号被瞬态信号损坏时，就会出现这样的情况。这些瞬变现象的出现和消失往往具有物理意义。因此，除了描述信号中的振荡分量外，能够定位这些瞬态是很重要的。为了模拟这种情况，创建一个由频率为150和200 Hz的两个正弦波组成的信号。采样频率为1khz。正弦波有不相交的时间支撑。万博1manbetx150赫兹的正弦波发生在100到300毫秒之间。200赫兹的正弦波发生在700毫秒到1秒之间。此外，有两个瞬态在222和800毫秒。 The signal is corrupted by noise.
rng默认的；Dt = 0.001;T = 0:dt:1-dt;addNoise = 0.025*randn(size(t));x = cos(2 *π* 150 * t) * (t > = 0.1 & t < 0.3) +罪(2 *π* 200 * t) * (t > 0.7);x = x+addNoise;X ([222 800]) = X ([222 800]) +[-2 2];图;情节(t。* 1000 x);包含(的毫秒）;ylabel (“振幅”）;
放大这两个瞬态，看看它们代表了150和200 Hz振荡中的扰动。
次要情节(2,1,1)情节(t(184:264)。* 1000 x (184:264));ylabel (“振幅”)标题(瞬态的)轴紧；次要情节(2,1,2)情节(t(760:840)。* 1000 x (760:840));ylabel (“振幅”)轴紧；包含(的毫秒）
利用解析Morlet小波得到并绘制CWT。
图;[cfs,f] = cwt(x,1/dt，“爱”）;轮廓(t * 1000 f, abs (cfs))网格在C = colorbar;包含(的毫秒) ylabel (“频率”) c.Label.String =“级”；
分析型Morlet小波的频率局部化较凹凸小波差，但时间局部化较好。这使得Morlet小波成为一种较好的瞬态定位方法。绘制平方级精细比例系数，以证明瞬态的局部化。
图;Wt = cwt(x,1/dt，“爱”）;情节(t。* 1000、abs (wt (1:)) ^ 2);包含(的毫秒）;ylabel (“级”）;网格在；标题(解析Morlet小波——精细尺度系数）;持有从
小波收缩以实现瞬态的时间定位，同时拉伸以实现150和200 Hz振荡的频率定位。
STFT只能将瞬变定位到窗口的宽度。您必须以分贝(dB)为单位绘制STFT，以便能够可视化瞬态。
[S,F,T] =谱图(x,50,48,128,1000);冲浪(t . * 1000 F, 20 * log10 (abs (S)),“edgecolor”，“没有”）;视图(0,90);轴紧；阴影插值函数；colorbar包含(“时间”), ylabel (“赫兹”）;标题(“STFT”）
瞬态只在宽带功率增加时出现在STFT中。比较在第一次瞬态出现之前(居中为183 msec)和之后(居中为223 msec)从STFT获得的短时功率估计。
图;情节(20 * log10 (abs (S (:, 80))));持有在；情节(20 * log10 (abs (S (:, 100))),“r”）;传奇(“T = 183毫秒”，'T = 223毫秒')包含(“赫兹”）;ylabel (“数据库”）;持有从；
STFT不提供与CWT相同程度的瞬态本地化能力。
利用逆CWT去除时间局部频率分量
创建一个由指数加权正弦波组成的信号。有两个25赫兹的分量——一个以0.2秒为中心，一个以0.5秒为中心。有两个70hz分量——一个以0.2秒为中心，一个以0.8秒为中心。前25-Hz和70-Hz分量同时出现。
T = 0:1/2000:1-1/2000;Dt = 1/2000;X1 = sin(50* *t).*exp(-50* *(t-0.2).^2);X2 = sin(50*pi*t).*exp(-100*pi*(t-0.5).^2);X3 = 2*cos(140*pi*t).*exp(-50*pi*(t-0.2).^2);X4 = 2*sin(140*pi*t).*exp(-80*pi*(t-0.8).^2);X = x1+x2+x3+x4;图;情节(t, x)标题(叠加信号的）
获取并显示CWT。
类(x, 2000);标题(用默认莫尔斯小波分析CWT）;
通过将CWT系数归零，去除发生在大约0.07秒到0.3秒之间的25hz分量。使用逆CWT (icwt)重建信号的近似值。
[cfs,f] = cwt(x,2000);T1 = .07;T2 = .33;F1 = 19;F2 = 34;cfs(f > F1 & f < F2, t> T1 & t < T2) = 0;Xrec = icwt(cfs);
显示重构信号的CWT。初始的25hz组件被移除。
类(xrec, 2000)
绘制原始信号和重建信号。
次要情节(2,1,1);情节(t, x);标题(原始信号的）;次要情节(2,1,2);情节(t, xrec)标题(“去除第一个25赫兹分量的信号”）;
最后，将重建信号与没有以0.2秒为中心的25赫兹分量的原始信号进行比较。
Y = x2+x3+x4;图;xrec情节(t)在情节(t y“r——”)传说(逆CWT近似，“原始信号无25赫兹”）;持有从
通过解析CWT确定准确的频率
当你用解析小波得到正弦波的小波变换时，解析CWT系数实际上编码了频率。
为了说明这一点，考虑从人耳获得的耳声发射。耳声发射(OAEs)是由耳蜗(内耳)发出的，它们的存在表明听力正常。加载并绘制OAE数据。数据以20千赫采样。
负载dpoae情节(t。* 1000,dpoaets)包含(的毫秒) ylabel (“振幅”）
发射是由一个开始于25毫秒，结束于175毫秒的刺激引起的。根据实验参数，发射频率应为1230 Hz。获得并绘制CWT。
类(dpoaets“撞”, 20000);包含(的毫秒）;
你可以通过找到频率最接近1230 Hz的CWT系数并检查它们的大小作为时间的函数来研究OAE的时间演化。将震级与时间标记一起画出来，标明诱发刺激的开始和结束。
[dpoaeCWT,f] = cwt(dpoaets，“撞”, 20000);[~，idx1230] = min(abs(f-1230));cfsOAE = dpoaeCWT(idx1230，:);情节(t。* 1000、abs (cfsOAE));持有在AX = gca;plot([25 25]，[AX.YLim(1) AX.YLim(2)]，“r”) plot([175 175]，[AX.YLim(1) AX.YLim(2)]，“r”)包含(的毫秒), ylabel (“级”）;标题(“CWT系数幅度”）
在诱发刺激的开始和OAE之间有一定的延迟。一旦唤醒刺激终止，OAE的强度立即开始衰减。
另一种隔离发射的方法是使用逆CWT在时域重构频率局部化近似。
通过在频率范围[1150 1350]Hz内反转CWT来重建频率局域发射近似。将原始数据与重建图和指示诱发刺激开始和结束的标记一起绘制出来。
xrec = icwt(dpoaeCWT，“撞”f (1150 1350));图(t.*1000,dpoaets)保持在；xrec情节(t。* 1000年,“r”) AX = gca;ylim = AX.YLim;ylim情节(25 [25],“k”) plot([175 175]，ylim，“k”)包含(的毫秒) ylabel (“振幅”)标题(“发射的频域重构”）
在时域数据中，你可以清楚地看到在触发刺激的应用和终止时，发射是如何启动和关闭的。
值得注意的是，即使为重建选择了一个频率范围，分析小波变换实际上编码了发射的确切频率。为了证明这一点，采用从解析CWT重建的发射近似的傅里叶变换。
XDFT = fft(xrec);Freq = 0:2e4/numel(xrec):1e4;XDFT = XDFT (1:numel(xrec)/2+1);图绘制(频率、abs (xdft));包含(“赫兹”）;ylabel (“级”)标题(基于cw的信号逼近的傅里叶变换）;[~，maxidx] = max(abs(xdft));流(频率为%4.2f Hz\n频率(maxidx));
频率为1230.00 Hz
结论及进一步阅读
在这个例子中，您学习了如何使用CWT来获得1-D信号的时频分析，使用解析小波类．您看到了CWT提供与STFT类似结果的信号示例，以及CWT可以提供比STFT更多可解释结果的示例。最后，您学习了如何使用重建信号的时间尺度(频率)局部近似icwt．
有关CWT的更多信息，请参阅小波工具箱文档。
参考资料:马拉特，S。《信号处理的小波漫游:稀疏方式》，学术出版社，2009年。
附录
本例中使用了以下helper函数。
您也可以从以下列表中选择一个网站: