声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 18975|回复: 55

[HHT] HHT边际谱的定义以及理解

  [复制链接]
发表于 2005-7-22 09:39 | 显示全部楼层 |阅读模式

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

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

x
(1)我对HHT边际谱的定义理解为就是将时频分布对时间求积分,具体地说就是将频率相同的时刻所对应的幅值进行累加。这样理解正确吗?
(2) 如果我的理解是正确的,那么写程序时关键的一步就是判断在各种IMF下哪些时刻的频率是相同的。在实际计算中,计算机的精度导致了各种IMF下各时刻的频率值几乎不可能是完全一致的,那么我要用什么样的算法来做边际谱呢?哪位高手能提供一下具体的算法或者有现成的code可以下载?

[ 本帖最后由 yejet 于 2006-9-17 12:57 编辑 ]

本帖被以下淘专辑推荐:

回复
分享到:

使用道具 举报

发表于 2005-7-22 21:54 | 显示全部楼层
Marginal spectrum是对Hilbert spectrum即H(A,w,t)进行时间积分得到的。Hilbert spectrum是指H(A,w,t)时、频、幅三维图。
发表于 2005-7-27 14:02 | 显示全部楼层
Fourier变换是将任何信号分解为正弦信号的加权和,而每一个正弦信号对应着一个固定的频率(Fourier频率)和固定的幅值,因此,用Fourier 变换分析频率不随时间变化的平稳信号是十分有效的。但对于频率随时间变化的非平稳信号,Fourier 变换就无能为力了。

HHT是历史上首次对Fourier变换的基本信号和频率定义作的创造性的改进。他们不再认为组成信号的基本信号是正弦信号,而是一种称为固有模态函数的信号,也就是满足以下两个条件的信号: (1) 整个信号中,零点数与极点数相等或至多相差1 ; (2) 信号上任意一点,由局部极大值点确定的包络线和由局部极小值点确定的包络线的均值均为零,即信号关于时间轴局部对称。

无论Hilbert谱中的频率还是边际谱中的频率(即瞬时频率) ,其意义都与Fourier分析中的频率(即Fourier 频率) 完全不同,但在Fourier分析中,某一频率处能量的存在,代表一个正弦或余弦波在整个时间轴上的存在,而边际谱h中某一频率处能量的存在仅代表在整个时间轴上可能有这样一个频率的振动波在局部出现过,h越大,代表该频率出现的可能性越大。

[ 本帖最后由 yejet 于 2006-9-17 12:57 编辑 ]
发表于 2005-8-23 16:56 | 显示全部楼层
目前最好的贴
发表于 2005-9-19 16:58 | 显示全部楼层

感谢 meil32的回贴

我正急需emd分解的matlab的程序呢,meil32的回贴借了我的燃眉之急。
发表于 2005-9-20 11:13 | 显示全部楼层
HHT问题很多,下面是别人提出来的
(1)因为分解的不稳定,分解所得的固有模态函数并非是系统的固有(一定)模态。
(2)且正因为分解的不稳定,分解结果一般也就很难用确切的物理过程来解释了。
(3)此外,再对一般情况下不稳定分解所得的本来就无明确物理意义的分量进行诸如波内调频(intra-wave frequency modulation)的细致分析,无疑亦是在做多此一举的无用功了。

[ 本帖最后由 yejet 于 2006-9-17 12:57 编辑 ]
发表于 2005-11-13 15:28 | 显示全部楼层
to:meil32,你的源程序里好像缺一个函数吧,extr(m),缺了这个不能正常运行啊,麻烦再传一下了,谢谢!
发表于 2005-12-4 20:45 | 显示全部楼层

to mail32

对你的程序里却两个函数,是不能正常运行的,你要是有的话请传上来,或发到邮箱wugeren@126.com

[ 本帖最后由 yejet 于 2006-9-17 12:57 编辑 ]
发表于 2005-12-22 09:25 | 显示全部楼层

HHT著作已经出版,

HHT著作已经出版,那位高人去搞个电子版的 !

http://www.worldscibooks.com/mathematics/5862.html
http://www.amazon.com/gp/product ... s&v=glance&n=283155

[ 本帖最后由 yejet 于 2006-9-17 12:57 编辑 ]
发表于 2005-12-29 22:38 | 显示全部楼层
这两本书那位高手能介绍买的地方,谢谢了!!!
发表于 2006-2-26 09:46 | 显示全部楼层
非常感谢你的程序

[ 本帖最后由 zhlong 于 2007-10-15 22:29 编辑 ]
发表于 2006-2-26 16:54 | 显示全部楼层
谢谢,我正在求这个程序,谢谢!
发表于 2006-2-27 21:15 | 显示全部楼层

回复:(qlisikefed)[转帖]HHT边际谱

其中缺少的extr试试下面的代码

  1. %extr
  2. function [indmin,indmax,indzer]=extr(x)
  3. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  4. watchon;

  5. FILENAME=[DATApath,DATAfile];
  6. x=load(FILENAME);

  7. numzer=1;
  8. nummin=1;
  9. nummax=1;
  10. for g=2:(length(x)-1)

  11. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  12. indzer(numzer)=g-1;
  13. numzer=numzer+1;
  14. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  15. indzer(numzer)=g;
  16. numzer=numzer+1;
  17. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  18. indzer(numzer)=(2*g-1)/2;
  19. numzer=numzer+1;
  20. end

  21. if (x(g-1)>x(g))&(x(g)>x(g+1))
  22. continue;
  23. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  24. nummin=nummin+1;
  25. elseif (x(g-1) continue;
  26. elseif (x(g-1)x(g+1))
  27. indmax(nummax)=g;
  28. nummax=nummax+1;
  29. end
  30. end
  31. %subplot(2,3,1)
  32. plot(indzer,x(indzer),'r-+')
  33. title('indzer')
  34. pause
  35. % subplot(2,3,2)
  36. plot(indmin,x(indmin),'b-o')
  37. title('indmin')
  38. pause
  39. %subplot(2,3,3)
  40. plot(indmax,x(indmax),'g-*')
  41. title('indmax')
  42. pause
  43. %subplot(2,3,4)
  44. plot(1:length(x),x,'k-s')
  45. title('x')
  46. pause

  47. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  48. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  49. %subplot(2,3,5)
  50. plot(1:length(x),envmin,'y-d')
  51. title('envmin')
  52. pause
  53. %subplot(2,3,6)
  54. plot(1:length(x),envmax,'c-h')
  55. title('envmax')
  56. pause
  57. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  58. title('mean((envmax+envmin)/2)')%extr
  59. function [indmin,indmax,indzer]=extr(x)
  60. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  61. watchon;

  62. FILENAME=[DATApath,DATAfile];
  63. x=load(FILENAME);

  64. numzer=1;
  65. nummin=1;
  66. nummax=1;
  67. for g=2:(length(x)-1)

  68. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  69. indzer(numzer)=g-1;
  70. numzer=numzer+1;
  71. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  72. indzer(numzer)=g;
  73. numzer=numzer+1;
  74. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  75. indzer(numzer)=(2*g-1)/2;
  76. numzer=numzer+1;
  77. end

  78. if (x(g-1)>x(g))&(x(g)>x(g+1))
  79. continue;
  80. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  81. nummin=nummin+1;
  82. elseif (x(g-1) continue;
  83. elseif (x(g-1)x(g+1))
  84. indmax(nummax)=g;
  85. nummax=nummax+1;
  86. end
  87. end
  88. %subplot(2,3,1)
  89. plot(indzer,x(indzer),'r-+')
  90. title('indzer')
  91. pause
  92. % subplot(2,3,2)
  93. plot(indmin,x(indmin),'b-o')
  94. title('indmin')
  95. pause
  96. %subplot(2,3,3)
  97. plot(indmax,x(indmax),'g-*')
  98. title('indmax')
  99. pause
  100. %subplot(2,3,4)
  101. plot(1:length(x),x,'k-s')
  102. title('x')
  103. pause

  104. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  105. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  106. %subplot(2,3,5)
  107. plot(1:length(x),envmin,'y-d')
  108. title('envmin')
  109. pause
  110. %subplot(2,3,6)
  111. plot(1:length(x),envmax,'c-h')
  112. title('envmax')
  113. pause
  114. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  115. title('mean((envmax+envmin)/2)')%extr
  116. function [indmin,indmax,indzer]=extr(x)
  117. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  118. watchon;

  119. FILENAME=[DATApath,DATAfile];
  120. x=load(FILENAME);

  121. numzer=1;
  122. nummin=1;
  123. nummax=1;
  124. for g=2:(length(x)-1)

  125. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  126. indzer(numzer)=g-1;
  127. numzer=numzer+1;
  128. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  129. indzer(numzer)=g;
  130. numzer=numzer+1;
  131. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  132. indzer(numzer)=(2*g-1)/2;
  133. numzer=numzer+1;
  134. end

  135. if (x(g-1)>x(g))&(x(g)>x(g+1))
  136. continue;
  137. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  138. nummin=nummin+1;
  139. elseif (x(g-1) continue;
  140. elseif (x(g-1)x(g+1))
  141. indmax(nummax)=g;
  142. nummax=nummax+1;
  143. end
  144. end
  145. %subplot(2,3,1)
  146. plot(indzer,x(indzer),'r-+')
  147. title('indzer')
  148. pause
  149. % subplot(2,3,2)
  150. plot(indmin,x(indmin),'b-o')
  151. title('indmin')
  152. pause
  153. %subplot(2,3,3)
  154. plot(indmax,x(indmax),'g-*')
  155. title('indmax')
  156. pause
  157. %subplot(2,3,4)
  158. plot(1:length(x),x,'k-s')
  159. title('x')
  160. pause

  161. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  162. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  163. %subplot(2,3,5)
  164. plot(1:length(x),envmin,'y-d')
  165. title('envmin')
  166. pause
  167. %subplot(2,3,6)
  168. plot(1:length(x),envmax,'c-h')
  169. title('envmax')
  170. pause
  171. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  172. title('mean((envmax+envmin)/2)')%extr
  173. function [indmin,indmax,indzer]=extr(x)
  174. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  175. watchon;

  176. FILENAME=[DATApath,DATAfile];
  177. x=load(FILENAME);

  178. numzer=1;
  179. nummin=1;
  180. nummax=1;
  181. for g=2:(length(x)-1)

  182. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  183. indzer(numzer)=g-1;
  184. numzer=numzer+1;
  185. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  186. indzer(numzer)=g;
  187. numzer=numzer+1;
  188. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  189. indzer(numzer)=(2*g-1)/2;
  190. numzer=numzer+1;
  191. end

  192. if (x(g-1)>x(g))&(x(g)>x(g+1))
  193. continue;
  194. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  195. nummin=nummin+1;
  196. elseif (x(g-1) continue;
  197. elseif (x(g-1)x(g+1))
  198. indmax(nummax)=g;
  199. nummax=nummax+1;
  200. end
  201. end
  202. %subplot(2,3,1)
  203. plot(indzer,x(indzer),'r-+')
  204. title('indzer')
  205. pause
  206. % subplot(2,3,2)
  207. plot(indmin,x(indmin),'b-o')
  208. title('indmin')
  209. pause
  210. %subplot(2,3,3)
  211. plot(indmax,x(indmax),'g-*')
  212. title('indmax')
  213. pause
  214. %subplot(2,3,4)
  215. plot(1:length(x),x,'k-s')
  216. title('x')
  217. pause

  218. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  219. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  220. %subplot(2,3,5)
  221. plot(1:length(x),envmin,'y-d')
  222. title('envmin')
  223. pause
  224. %subplot(2,3,6)
  225. plot(1:length(x),envmax,'c-h')
  226. title('envmax')
  227. pause
  228. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  229. title('mean((envmax+envmin)/2)')%extr
  230. function [indmin,indmax,indzer]=extr(x)
  231. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  232. watchon;

  233. FILENAME=[DATApath,DATAfile];
  234. x=load(FILENAME);

  235. numzer=1;
  236. nummin=1;
  237. nummax=1;
  238. for g=2:(length(x)-1)

  239. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  240. indzer(numzer)=g-1;
  241. numzer=numzer+1;
  242. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  243. indzer(numzer)=g;
  244. numzer=numzer+1;
  245. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  246. indzer(numzer)=(2*g-1)/2;
  247. numzer=numzer+1;
  248. end

  249. if (x(g-1)>x(g))&(x(g)>x(g+1))
  250. continue;
  251. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  252. nummin=nummin+1;
  253. elseif (x(g-1) continue;
  254. elseif (x(g-1)x(g+1))
  255. indmax(nummax)=g;
  256. nummax=nummax+1;
  257. end
  258. end
  259. %subplot(2,3,1)
  260. plot(indzer,x(indzer),'r-+')
  261. title('indzer')
  262. pause
  263. % subplot(2,3,2)
  264. plot(indmin,x(indmin),'b-o')
  265. title('indmin')
  266. pause
  267. %subplot(2,3,3)
  268. plot(indmax,x(indmax),'g-*')
  269. title('indmax')
  270. pause
  271. %subplot(2,3,4)
  272. plot(1:length(x),x,'k-s')
  273. title('x')
  274. pause

  275. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  276. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  277. %subplot(2,3,5)
  278. plot(1:length(x),envmin,'y-d')
  279. title('envmin')
  280. pause
  281. %subplot(2,3,6)
  282. plot(1:length(x),envmax,'c-h')
  283. title('envmax')
  284. pause
  285. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  286. title('mean((envmax+envmin)/2)')%extr
  287. function [indmin,indmax,indzer]=extr(x)
  288. [DATAfile DATApath]=uigetfile('*.txt','输入信号');
  289. watchon;

  290. FILENAME=[DATApath,DATAfile];
  291. x=load(FILENAME);

  292. numzer=1;
  293. nummin=1;
  294. nummax=1;
  295. for g=2:(length(x)-1)

  296. if (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)==0)
  297. indzer(numzer)=g-1;
  298. numzer=numzer+1;
  299. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g)==0)
  300. indzer(numzer)=g;
  301. numzer=numzer+1;
  302. elseif (x(g-1)*x(g)<0)|(x(g-1)*x(g)==0)&(x(g-1)~=0)&(x(g)~=0)
  303. indzer(numzer)=(2*g-1)/2;
  304. numzer=numzer+1;
  305. end

  306. if (x(g-1)>x(g))&(x(g)>x(g+1))
  307. continue;
  308. elseif (x(g-1)>x(g))&(x(g) indmin(nummin)=g;
  309. nummin=nummin+1;
  310. elseif (x(g-1) continue;
  311. elseif (x(g-1)x(g+1))
  312. indmax(nummax)=g;
  313. nummax=nummax+1;
  314. end
  315. end
  316. %subplot(2,3,1)
  317. plot(indzer,x(indzer),'r-+')
  318. title('indzer')
  319. pause
  320. % subplot(2,3,2)
  321. plot(indmin,x(indmin),'b-o')
  322. title('indmin')
  323. pause
  324. %subplot(2,3,3)
  325. plot(indmax,x(indmax),'g-*')
  326. title('indmax')
  327. pause
  328. %subplot(2,3,4)
  329. plot(1:length(x),x,'k-s')
  330. title('x')
  331. pause

  332. envmax = interp1(indmax,x(indmax),1:length(x),'spline');
  333. envmin = interp1(indmin,x(indmin),1:length(x),'spline');
  334. %subplot(2,3,5)
  335. plot(1:length(x),envmin,'y-d')
  336. title('envmin')
  337. pause
  338. %subplot(2,3,6)
  339. plot(1:length(x),envmax,'c-h')
  340. title('envmax')
  341. pause
  342. plot(1:length(x),mean((envmax+envmin)/2),'m-p')
  343. title('mean((envmax+envmin)/2)')
复制代码

[ 本帖最后由 yejet 于 2006-9-17 12:58 编辑 ]
发表于 2006-3-17 13:32 | 显示全部楼层

我也正准备研究一下

本人东南大学无线电系大三学生,正准备做个SRTP,各位大侠多多指教.MEIL32的程序我收下了,缺的两个函数有时间发到syyu514@126.com.谢谢各位。

[ 本帖最后由 yejet 于 2006-9-17 12:58 编辑 ]
发表于 2006-3-21 18:47 | 显示全部楼层
我刚刚下载
多谢了!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-11-25 01:00 , Processed in 0.060310 second(s), 20 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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