边际谱相关讨论的阶段性总结
精华帖<Hilbert边界谱>中，版友"破凰"给出了一个较为完整的边际谱计算程序。鉴于该贴楼层过高、内容繁杂、版友查找信息不便，特此新开一帖对其中的讨论做一些阶段性总结。对实际信号的边际谱处理，还有待广大版友补充，谢谢！一、EMD程序问题 目前多数版友手头的程序都是网上下载的法国人的程序包，下面的边际谱程序也是基于该程序包的。该程序包可以到法国人的主页http://perso.ens-lyon.fr/patrick.flandrin/emd.html下载，也可以到我的E盘http://zhlong.ys168.com下载。 有关程序的使用问题请参考：送给搞EMD或者HHT但没有下载到完整程序的朋友求助EMD程序如何使用
二、原帖37楼版友“破凰"的程序原帖由 破凰 于 2007-4-27 21:51 发表 http://www.chinavib.com/forum/images/common/back.gif
clear;
fs=1000; %fs为采样频率;tspan=2;
t=1/fs:1/fs:tspan;N=length(t); %采样点数y1=5*sin(2*pi*241*t);
y2=3*sin(2*pi*73*t);
y=; %IMF集
%%%%%%%%%%%%%求边际谱
=hhspectrum(y);
=toimage(A,fa,tt,length(tt));
E=flipud(E);% 这句要或不要与你下载的程序包的版本有关，它的作用是使得你的边际谱左右翻转
% 如果觉得结果不对试着添加或删除这句for k=1:size(E,1)
bjp(k)=sum(E(k,:))*1/fs*1/tspan;
end
f=(0:N-3)/N*(fs/2);
plot(f,bjp);
xlabel('频率 / Hz');
ylabel('幅值'); 注意程序中红色的注释部分。%程序已按本帖七楼“破凰”版友的建议做了适当修改
三、关于频率点数的设置问题
原帖由 chiu 于 2007-7-27 13:21 发表 http://www.chinavib.com/forum/images/common/back.gif
clear;
fs=1000; %fs为采样频率;
N=1000; %采样点数
t=1/fs:1/fs:1;
y1=2*sin(60*pi*t);
y2=5*sin(90*pi*t);
y=; %IMF集
%%%%%%%%%%%%%求边际谱
=hhspectrum(y ...
<Hilbert边界谱>一帖中113楼chiu版友提出了设定边际谱频率点数的问题。因为在“破凰”的程序中边际谱的频率点数是N-3，所以当信号很长时，势必会使得程序中的时频距阵E过大，影响计算速度。为此将程序稍做修改如下，其中NN用于设定频率点数。NN越大，频率分辨率越高，占用的内存也越多。
clear;
fs=1000; %fs为采样频率;tspan=2;
t=1/fs:1/fs:tspan;N=length(t); %采样点数y1=5*sin(2*pi*241*t);
y2=3*sin(2*pi*73*t);
y=; %IMF集
%%%%%%%%%%%%%求边际谱
=hhspectrum(y);NN=500;%设置频率点数
=toimage(A,fa,tt,NN);E=flipud(E);% 注释同上
for k=1:NN
bjp(k)=sum(E(k,:))*1/fs*1/tspan;
end
f=(0:NN-1)/NN*(fs/2);
plot(f,bjp);
xlabel('频率 / Hz');
ylabel('幅值');
四、整周期采样对边际谱的影响
原帖由 破凰 于 2007-6-12 16:59 发表 http://www.chinavib.com/forum/images/common/back.gif
是的，求频率坐标系列时应该用fs/2而不是fs。
f=(0:N-3)/N*fs/2;
当不是对信号进行整周期采样时，效果就没那么好了，尤其是幅值很不准！
版友“破凰”提到整周期采样对边际谱幅值的影响问题。大家可以将上面程序中的信号改为如下试试得到的边际谱幅值如何。y1=5*sin(2*pi*241.5*t);
y2=3*sin(2*pi*73.5*t); 画出图来大家会发现频率点差不多对应，但幅值不对了。这个主要是因为hhspectrum函数计算瞬时频率和瞬时幅值时是基于信号希尔伯特变换的，用到了hilbert函数，而hilbert中又用到了fft函数。非整周期采样会给fft带来泄漏，最后传递到边际谱。详细讨论参考hilbert变换怪现象！ 9楼和10楼。
五、边际谱与EMD分解结果之间的关系原帖由 form 于 2007-5-22 10:47 发表 http://www.chinavib.com/forum/images/common/back.gif
谢谢破凰，总算有人给以肯定的答案了，我有些疑问：
1）按照上面程序运行，出来的谱线反而在右侧，去掉E=flipud(E);才对应了那副图，是否因为E的频率是反的，但size（E，1）算得是大小，因此这句for k=1:siz ...
版友form指出EMD虚假分量会影响边际谱，应该是这样的。边际谱是由HHT时频图沿时间轴积分得到的，所以得到准确幅值谱的前提是时频图准确。而EMD分解过程通常存在这样一些问题：高频分解不彻底、许多无意义的低频干扰，这些都影响HHT时频图，进而影响边际谱。所以要得到准确边际谱的前提是EMD分解准确。事实上大家可以试验一下，比如y1=5*sin(2*pi*241*t);y2=3*sin(2*pi*73*t);y=y1+y2;，然后对y做EMD分解，即使只利用前两阶IMF做边际谱，同样会发现边际谱幅值不准确。除非分解得到的这两阶IMF刚好和y1、y2完全一样。我想这也是为什么大家看到有关边际谱的论文中，只给出相对幅值，而没有给出准确幅值的原因之一。 综上，大家在讨论边际谱的时候，如果觉得得到的边际谱不尽人意，问题的本质个人觉得应该不在边际谱程序上，而在于信号EMD分解的准确与否。（当然还有整周期采样问题）
六、本版一些边际谱相关的帖子
Hilbert边界谱[转帖]HHT边际谱我边际谱作的图怎么不对啊
大家可以输入关键字“边际谱”搜索更多相关内容。欢迎大家补充、指正。
[ 本帖最后由 zhlong 于 2007-10-22 17:23 编辑 ] 占位编辑
回复 #2 zhlong 的帖子
支持一下，顶一下！感谢zhlong的辛苦，在这炎炎夏日给我们送来了凉爽的风。但是有一个问题，解决实际信号怎么办？好像一直没有成型的办法。总也没有讨论出来。我建议把题目改成 阶段总结，比较好，最后怎么也要有个结果啊
尤其是对于非平稳信号，很难说整周期，这些程序得到的效果恐怕就要商榷了
[ 本帖最后由 zhangnan3509 于 2007-8-10 21:03 编辑 ]
回复 #3 zhangnan3509 的帖子
我的观点还是：对实际信号EMD分解结果不对的话，无法保证边际谱的效果，所以根本问题应该不在边际谱，而在于EMD分解得到的IMFs。当然如果能从边际谱着手保证其准确性，而不用管EMD分解效果，那当然最好了。ps：已按版主意见更改标题:@) 。
回复 #4 zhlong 的帖子
你说的按道理是这样，实际上能让EMD有漂亮的分解结果 的并不多。我们不可能总得到y=; %IMF集 这样好的IMF集啊回复 #5 zhangnan3509 的帖子
这样好的IMF集只是为了说明边际谱程序的正确性。对于实际信号的边际谱，确实有待探讨。
原帖由 MVH 于 2005-9-20 11:13 发表 http://www.chinavib.com/forum/images/common/back.gif
HHT问题很多，下面是别人提出来的
(1)因为分解的不稳定，分解所得的固有模态函数并非是系统的固有(一定)模态。
(2)且正因为分解的不稳定，分解结果一般也就很难用确切的物理过程来解释了。
(3)此外，再对一般情况下不稳定分解所得的本来就无明确物理意义的分量进行诸如波内调频(intra-wave frequency modulation)的细致分析，无疑亦是在做多此一举的无用功了。
对于上面的论断，我们可以保留自己的观点。但也说明了HHT方法的重点研究方向或是弱点在哪里。
[ 本帖最后由 zhlong 于 2007-8-10 21:14 编辑 ] 版主太不厚道了，沙发都抢了:lol
总结得很好！
ps：黄在原文中指出，边际谱的幅值代表某个频率成分的能量，能量越大，该频率的振动波在整个时间轴出现的可能性也越大。时间取得越长，振动能量越大，边际谱的幅值也就越大。因此，我认为边际谱的幅值和信号的幅值没有直接联系。在我的例子中，信号的幅值与边际谱的幅值相等只是一种偶然，那是取样时间为1s的特殊情况。取样时间为2s，幅值将加倍（大家可以试试）。
为了正确表征信号幅值，可以在边际谱定义的基础上除以取样时间。
============================
已按破凰版友意见在1楼程序中进行了适当修改。
再次对破凰版友表示感谢！
by zhlong
============================
[ 本帖最后由 zhlong 于 2007-8-10 22:09 编辑 ]
回复 #7 破凰 的帖子
站位是应该的，也许以后有什么补充呢！所以也不算什么不厚道啊应该说这也不能怨大伙，已经很详细了，问题是EMD分解的根本问题就是IMF定义上有问题啊
这个解决了就好很多啊
因此现阶段也只能这样了，那我们也可以尝试在EMD分解之前做一些预处理把一些能解除的干扰降低一些，或许也是一个思路吧
回复 #7 破凰 的帖子
（边际谱的幅值和信号的幅值没有直接联系。）这是正确的，http://forum.vibunion.com/forum/thread-45278-1-1.html 15楼解释过了。你说的问题最好用程序来说明，这样更清楚
回复 #7 破凰 的帖子
^_^，破凰确实高，对HHT很有研究。实际信号得到边际谱效果不好、很多时候无法反应信号频率成份，破凰有何高见？不要吝啬啊:lol
[ 本帖最后由 zhlong 于 2007-8-10 21:55 编辑 ] http://forum.vibunion.com/forum/thread-43809-1-3.html 13，14楼有我现在用的HHT谱，边际谱，时频图的程序，并非rilling的，大家相互印证吧。我事先声明，原因大家都知道，有问题几乎都是因为EMD分解IMF的问题，这个办法也不是解决所有问题的良方妙药，其实跟破凰作的程序是原理相同的，只是方式不同罢了。
要是应用过程效果不尽如人意，大家不要怪我啊:lol
[ 本帖最后由 zhangnan3509 于 2007-8-10 21:42 编辑 ]
回复 #9 zhangnan3509 的帖子
呵呵，在例子中将点数改为2000就可以了。刚才我只是随便写点，就偷懒了哈！
回复 #10 zhlong 的帖子
不敢当哦！在形成时频图时（也就是用toimage时），频率点取少些，可以减少能量泄漏，效果可能会好些。
还有就是zhangnan3509 版主说的先做一些预处理。
[ 本帖最后由 破凰 于 2007-8-10 21:58 编辑 ]
回复 #11 zhangnan3509 的帖子
你在求边际谱时，也是除了取样时间的吧。回复 #13 破凰 的帖子
你说的取少些点，这有时候不现实，多少算是少？有时候要处理的数据2万算一般多的，这时候怎么办？[ 本帖最后由 zhlong 于 2007-8-10 22:10 编辑 ]