声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1337|回复: 4

[编程技巧] 请问各位大侠如何绘制小波的频谱(例如morlet小波)

[复制链接]
发表于 2010-9-22 20:50 | 显示全部楼层 |阅读模式

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

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

x
请给出源程序,不胜感激
回复
分享到:

使用道具 举报

发表于 2010-10-13 19:51 | 显示全部楼层
关注,借贴同问
发表于 2010-10-19 16:45 | 显示全部楼层
本帖最后由 happy 于 2010-10-19 16:48 编辑

本程序来自网络,版权归原作者所有,由于原作者无法考究,故此处没有标明
如原作者看到本贴,请和管理员联系添加相关说明或做其他处理

  1. %%======================%%
  2. %%小波包分析程序
  3. %%=======================%%
  4. clc;clear all;close all;
  5. % %-----------------输入实际信号---------------
  6. N=8192;
  7. fs=20480;
  8. [DATAfile DATApath]=uigetfile('.txt','输入信号');%显示一个取文件的窗口
  9. FILENAME=[DATApath,DATAfile]
  10. DATA=load(FILENAME);

  11. Data=DATA(1:N);
  12. DataLen=length(Data);%数据长度
  13. t=[0:(DataLen-1)]'/fs;%时间数组
  14. % -------------------小波包分解---------------
  15. n=3;%分解层数
  16. rr=wpdec(Data,n,'db15');%小波包分解

  17. % ---------对特征频带波形进行分析-----
  18. for i=1:2^n
  19.      Datar(i,:)=wprcoef(rr,[n,i-1]);%求取每个频带的小波包系数
  20. end
  21. Datar=[Datar(1,:);Datar(2,:);Datar(4,:);Datar(3,:);Datar(7,:);Datar(8,:);Datar(6,:);Datar(5,:)];%数据调整

  22. f=(0:N/2-1)*fs/N;
  23. for i=1:2^n
  24.     DataFFT(i,:)=abs(fft(Datar(i,:)))/(N/2);  %%fft
  25. end
  26. % %-----------小波包分解计算能量-------
  27. for i=0:(2^n-1)
  28.     rcfs=wprcoef(rr,[n i]);
  29.     r(i+1,:)=rcfs;
  30. end
  31. r_check=[r(1,:);r(2,:);r(4,:);r(3,:);r(7,:);r(8,:);r(6,:);r(5,:)];
  32. for i=1:2^n
  33.     xx(i)=sum(r_check(i,:).^2)/N;%求能量
  34. end
  35. temp=sum(xx);
  36. for i=1:2^n
  37.     xx(i)=xx(i)/temp;%%能量归一化
  38. end
  39. %---------画小波分频带图---------
  40. % figure()
  41. % for i=1:2^n/2
  42. %     subplot(4,1,i);
  43. %     plot(t,Datar(i,:));ylabel(i);axis([min(t) max(t) min(Datar(i,:)) max(Datar(i,:))]);
  44. %     if i==1
  45. %        title('小波包分解')
  46. %     end
  47. % end
  48. % xlabel('time');
  49. % grid on
  50. % figure()
  51. % for i=2^n/2+1:2^n
  52. %     subplot(4,1,i-2^n/2);
  53. %     plot(t,Datar(i,:));ylabel(i);axis([min(t) max(t) min(Datar(i,:)) max(Datar(i,:))]);
  54. %         if i== 2^n/2+1
  55. %            title('小波包分解')
  56. %         end
  57. % end
  58. % xlabel('time');
  59. % grid on
  60. %%-----------对每个频带求取频谱--------

  61. figure()
  62. set(gcf,'Name','ff1')
  63. for i=1:2^n/2
  64.     subplot(4,1,i);
  65.     plot(f,DataFFT(i,1:N/2));ylabel(i);axis([0 10240 0 max(DataFFT(i,1:N/2))]);grid on
  66.     if i==1
  67.        title('小波包各频带频谱图')
  68.     end
  69. end
  70. xlabel('frequency')

  71. figure()
  72. set(gcf,'Name','ff2')
  73. for i=2^n/2+1:2^n
  74.     subplot(4,1,i-2^n/2);
  75.     plot(f,DataFFT(i,1:N/2));ylabel(i);axis([0 10240 0 max(DataFFT(i,1:N/2))]);grid on
  76.     if i==2^n/2+1
  77.        title('小波包各频带频谱图')
  78.     end
  79. end
  80. xlabel('frequency')

  81. %-------------画出能量谱-----------
  82. figure()
  83. set(gcf,'Name','能量谱')
  84. width=0.5;
  85. bar(xx,width);
  86. colormap(cool);
  87. title('小波包能量谱图');
  88. xlabel('频 带');
  89. ylabel('能量归一化');
  90. axis([0 9 0 1]);
  91. grid on
复制代码

评分

1

查看全部评分

 楼主| 发表于 2010-10-19 22:01 | 显示全部楼层
回复 happy 的帖子

万分感谢
 楼主| 发表于 2010-10-25 19:28 | 显示全部楼层
回复 happy 的帖子

你好能否再为我解答几个问题?

小波多分辨分析方面的问题
http://forum.vibunion.com/forum- ... fromuid-151591.html


您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-5-6 10:57 , Processed in 0.262300 second(s), 18 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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