声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 35938|回复: 75

[小波] 小波时频(尺度)图的绘制原理与实现

  [复制链接]
发表于 2007-12-17 17:03 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
应广大版友的需要,下面将介绍小波时频(尺度)图的绘制原理,并举例加以说明。

-、绘制原理

    1.需要用到的小波工具箱中的三个函数

    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/4,c/2,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('小波时频图');

     程序运行结果如下:
wavsc.jpg
     说明:(1)应用时只须改变wavename和totalscal两个参数即可。
              (2)在这个例子中,最好选用复的morlet小波,其它小波的分析效果不好,而且morlet小波的带宽参数和中心频率取得越大,时频图上反映的时频聚集性越好。




[ 本帖最后由 破凰 于 2007-12-17 17:29 编辑 ]

点评

赞成: 5.0
赞成: 5
  发表于 2016-9-15 18:47

评分

3

查看全部评分

回复
分享到:

使用道具 举报

发表于 2007-12-17 17:13 | 显示全部楼层

回复 #1 破凰 的帖子

防止灌水和保护帖子的需要,现在将帖子关闭如有其他问题和建议,请专门发帖说明。:handshake
发表于 2009-3-26 20:13 | 显示全部楼层

回复 楼主 破凰 的帖子

请教你一个一直困扰我的问题'cmor3-3'这里面的贷款系数和中心频率如何确定呢?
发表于 2009-3-30 22:21 | 显示全部楼层

楼主 破凰 的帖子

楼主,我按照你所给的上述方法用db小波对一个时变信号做时频法分析,但结果和原文出入很大!能不能给解释呢!
初来乍到,传文件是给弄错了,结果给弄到我的空间了,烦劳你看看!先谢啦!
发表于 2009-4-26 18:04 | 显示全部楼层
我觉得这一段写的不是很简单明了,绕晕了
改为这样理解,大家看是不是好一些
Fa=Fc*fs/a;
根据实际信号的频率分布范围确定要选择的频率范围fmin,fmax;
但是对于fs采样频率,只有小于或等于fs/2的频率成分,才满足采样定理,所以设置的频率范围fmax<=fs/2
如果想得到频率为等距间隔,那么直接:
f=linspace(fmin,fmax,n);
a=Fc*fs./f;
这样就可以先确定fmin,fmax 来确定a的范围,而且又保证f是等距间隔

我只是觉得第3段写的不易理解,其它都很好
楼主很多帖子我都看过,很佩服啊,

3.尺度序列的确定

由式(1)可以看出,为使转换后的频率序列是一等差序列,尺度序列必须取为以下形式:

c/totalscal,...,c/(totalscal-1),c/4,c/2,c (2)

其中,totalscal是对信号进行小波变换时所用尺度序列的长度(通常需要预先设定好),c为一常数。

下面讲讲c的求法。

根据式(1)容易看出,尺度c/totalscal所对应的实际频率应为fs/2,于是可得

c=2×Fc/totalscal (3)

将式(3)代入式(2)便得到了所需的尺度序列。
发表于 2009-9-14 13:27 | 显示全部楼层

回复 楼主 破凰 的帖子

楼主,你的例程怎么运行不了,我是复制过去的。
??? Attempt to execute SCRIPT cwt as a function.

Error in ==> example_cwt at 15
coefs=cwt(s,scal,wavename);       %得到小波系数
发表于 2010-8-27 11:56 | 显示全部楼层
回复 cuixuepeng 的帖子 你没有小波工具箱


   
发表于 2010-8-30 13:44 | 显示全部楼层
本帖最后由 xinglong-liu 于 2010-8-30 13:47 编辑

回复 chihewanle
  1. 如果想得到频率为等距间隔,那么直接:
  2. f=linspace(fmin,fmax,n);
  3. a=Fc*fs./f;
  4. 这样就可以先确定fmin,fmax 来确定a的范围,而且又保证f是等距间隔
复制代码
小波的频宽是随中心分析频率变化的,也就是,中心分析频率越高,频宽就越宽。因此,当采用等间隔的分析频率设置时,会出现分析频率的严重重叠,分析效率降低。
发表于 2010-8-30 20:37 | 显示全部楼层
  

              
发表于 2010-8-30 20:38 | 显示全部楼层
  c=2×Fc/totalscal       应该是        c=2×Fc×totalscal      

点评

赞成: 5.0
赞成: 5
对的!  发表于 2012-5-28 14:55
发表于 2010-9-8 20:44 | 显示全部楼层
很受用,支持,学习
发表于 2010-10-11 15:14 | 显示全部楼层
学习了,新手,所以对于此图的实际意义不太了解,觉得没有那种频率-幅值图给人的视觉效果好,可能我还是不太了解小波的缘故吧.得努力了.
发表于 2010-10-18 21:38 | 显示全部楼层
能否给出不使用cwt函数的小波变换,s变换、广义s变换的MATLAB程序,非常感谢
发表于 2010-10-19 09:37 | 显示全部楼层
回复 破凰 的帖子

弱弱的问一下,为什么最后信号S的时频图上一个时间对应两个频率,两个正弦信号合成后的信号S的频率要怎么求……
发表于 2011-1-9 10:51 | 显示全部楼层
新手,有个问题想请教一下,你的程序画处理的时频图,为什么频率的坐标是倒着的?谢谢了。
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-13 16:09 , Processed in 0.081076 second(s), 24 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表