对简单仿真信号进行EMD分解和后续Hilbert谱分析时的一些问题

Snikas

本人这个学期的毕业设计与HHT在轨道检测中的应用有关，刚接触这个方面的东西感觉还蛮有兴趣的，本人不久前将于德介老师的《机械故障诊断的Hibert-Huang变换方法》大致看了下，对EMD、IMF和HHT等概念有了一些基本的了解，这几天天天在泡论坛，想通过Matlab将上书中的一些简单仿真信号的处理自己实现一遍，例如书中P29的仿真信号：

x(t)=2*sin(30*pi*t)+4*sin(20*pi*t).*sin(2*pi*1/10*t)+sin(2*pi*5*t)，t∈[0,1]

我的matlab程序如下（很多参照了论坛中达人提供的程序）：

T=0:0.01:1;
x=2*sin(30*pi*T)+4*sin(20*pi*T).*sin(2*pi*1/10*T)+sin(2*pi*5*T);
plot(T,x)  %绘出仿真波形图，如图1
OPTIONS.DISPLAY=2;
OPTIONS.T=0:0.01:1;
imf=emd(x,OPTIONS);
emd_visu(x,0:0.01:1,imf)  %绘出各个IMF波形，如图2
[A,f,tt]=hhspectrum(imf);
[im,tt]=toimage(A,f);
disp_hhs(im,[],100);
colormap(flipud(gray))  %绘出Hilbert谱，如图3


     从仿真信号的函数可以看到理论上对信号EMD分解出来应该是一个15hz的正弦信号、一个5hz的正弦信号和一个10hz的调幅信号，可是在图2中一共得到了4个IMF和一个residue，而且前3个IMF的效果很差，比书上给的EMD分解图效果差很多，第2个IMF几乎没有调幅的特征了，我不知道是我做的EMD分解过程有误还是这就是EMD的端点效应问题？要是是端点效应的问题话，不知道用Matlab该如何实现抑制端点效应？
     还有我用disp_hhs(im,[],100)语句倒是把hilbert谱的纵坐标转换成了实际频率（采样频率100hz），或者可以单独把hilbert谱的横坐标转换成实际时间，但不知道怎么同时把两个坐标都转换成实际的，我想通过disp_hhs(im,tt,[],100)（tt=0:0.01:1）同时将两坐标转换过来，但会有错误，提示说输入参数数量太多了。。。搞了半天不知道该怎么解决这个问题。。。

[ 本帖最后由 Snikas 于 2008-4-12 19:50 编辑 ]

若菱

我也很奇怪这个问题，不过你可以修改disp_hhs.m这个函数的代码，将其中
if fs == 0
  imagesc(t,[0,0.5],im,[inf,0]);
  ylabel('normalized frequency')
else
imagesc(t,[0,0.5*fs],im,[inf,0]);改为imagesc(t/fs,[0,0.5*fs],im,[inf,0]);
  ylabel('frequency')
end
然后通过disp_hhs(im,[],fs);调用，则可以同时在图中表现对应的实际时间和频率了

Snikas

确实可行，谢谢哈:victory: 原程序确实有点问题哈，自己说明里都说了有四个输入参量，可是实际又不行。。。

green

disp_hhs(im,tt,[],100)中的tt好像不是0：0.01：1，因为程序中对时间进行了截取，可以这样处理(t-1)/fs

gb9813

为什么我在调试程序时总是出现下面的错误呢？
??? Undefined command/function 'instfreq'.

Error in ==> hhspectrum at 79
  f(i,:)=instfreq(an(i,:)',tt,l)';

Error in ==> Untitled5 at 8
[A,f,tt]=hhspectrum(imf);
是没有instfreq函数吗？

zhangnan3509

没有instfreq这个函数

jinnian

将error(nargchk(1,3,nargin));改为error(nargchk(1,4,nargin));
调用时用disp_hhs(im,t/fs,[],fs),也可以解决问题的！

nkdtxf

很感谢楼上的几位，让我学到了好多东西。知道了，调用的函数能修改。谢谢！

huangqiang

这个与设置的迭代次数及采样点数多少有关。
我分解的是下面的图。

dzkt

请问楼上参数是怎么设置的？

xiangyu537

参看相关文章，MSSP杂志上 ，有专门的文章讲述，参数的选择问题，楼上可以看一下

南宫轩宇

回复 1 # Snikas 的帖子

请问楼主如何才能让每一个imf图上的横纵坐标的值显示出来呢，谢谢哈，急！！！

南宫轩宇

请问楼主如何才能让每一个imf图上的横纵坐标的值显示出来呢，谢谢哈，急！！！

saige2

谢谢，有收获

yghit08

回复 7 # jinnian 的帖子

确实好用，谢谢！