马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
精华帖<Hilbert边界谱>中,版友"破凰"给出了一个较为完整的边际谱计算程序。鉴于该贴楼层过高、内容繁杂、版友查找信息不便,特此新开一帖对其中的讨论做一些阶段性总结。对实际信号的边际谱处理,还有待广大版友补充,谢谢!
一、EMD程序问题 有关程序的使用问题请参考:
二、原帖37楼版友“破凰"的程序 原帖由 破凰 于 2007-4-27 21:51 发表
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=[y1;y2;zeros(size(y1))]; %IMF集
%%%%%%%%%%%%%求边际谱
[A,fa,tt]=hhspectrum(y);
[E,tt1]=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 发表
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=[y1;y2;zeros(size(y1))]; %IMF集
%%%%%%%%%%%%%求边际谱
[A,fa,tt]=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=[y1;y2;zeros(size(y1))]; %IMF集
%%%%%%%%%%%%%求边际谱
[A,fa,tt]=hhspectrum(y); NN=500;%设置频率点数
[E,tt1]=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 发表
是的,求频率坐标系列时应该用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 发表
谢谢破凰,总算有人给以肯定的答案了,我有些疑问:
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分解的准确与否。(当然还有整周期采样问题)
六、本版一些边际谱相关的帖子
大家可以输入关键字“边际谱”搜索更多相关内容。 欢迎大家补充、指正。
[ 本帖最后由 zhlong 于 2007-10-22 17:23 编辑 ] |