对简单仿真信号进行EMD分解和后续Hilbert谱分析时的一些问题
本人这个学期的毕业设计与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∈
我的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
=hhspectrum(imf);
=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,,im,);
ylabel('normalized frequency')
else
imagesc(t,,im,);改为imagesc(t/fs,,im,);
ylabel('frequency')
end
然后通过disp_hhs(im,[],fs);调用，则可以同时在图中表现对应的实际时间和频率了
回复 4楼 的帖子
确实可行，谢谢哈:victory: 原程序确实有点问题哈，自己说明里都说了有四个输入参量，可是实际又不行。。。回复 楼主 的帖子
disp_hhs(im,tt,[],100)中的tt好像不是0：0.01：1，因为程序中对时间进行了截取，可以这样处理(t-1)/fs 为什么我在调试程序时总是出现下面的错误呢？??? Undefined command/function 'instfreq'.
Error in ==> hhspectrum at 79
f(i,:)=instfreq(an(i,:)',tt,l)';
Error in ==> Untitled5 at 8
=hhspectrum(imf);
是没有instfreq函数吗？
回复 5楼 的帖子
没有instfreq这个函数 将error(nargchk(1,3,nargin));改为error(nargchk(1,4,nargin));调用时用disp_hhs(im,t/fs,[],fs),也可以解决问题的！ 很感谢楼上的几位，让我学到了好多东西。知道了，调用的函数能修改。谢谢！ 这个与设置的迭代次数及采样点数多少有关。
我分解的是下面的图。
对简单仿真信号进行EMD分解和后续Hilbert谱分析时的一些问题
请问楼上参数是怎么设置的？ 参看相关文章，MSSP杂志上 ，有专门的文章讲述，参数的选择问题，楼上可以看一下 回复 1 # Snikas 的帖子请问楼主如何才能让每一个imf图上的横纵坐标的值显示出来呢，谢谢哈，急！！！ 请问楼主如何才能让每一个imf图上的横纵坐标的值显示出来呢，谢谢哈，急！！！{:{39}:} 谢谢，有收获 回复 7 # jinnian 的帖子
确实好用，谢谢！
页:
[1]
2