关于HHT谱和边际谱的理解和疑问,请大家帮忙分析一下啊
本帖最后由 cwb 于 2014-9-3 23:14 编辑我有三个版本的程序:National Taiwan Central University, Matlab File Exchange Center,和G Rilling(2007版)的。我试过了后两个版本:Matlab File Exchange Center是分解出IMF,然后画出 时间-频率图,而幅度好像是用灰度(我不确定)表示的,没有给出画幅度-频率-时间三维图的画法,也没有给出边际谱的画法;而G rilling的程序也是一样,我在论坛上看到的画A-f-t谱图和边际谱图的的方法是:
>> x=xlsread('matlab');
>> y=x(:,2);
>> imf=emd(y);%y是行向量或者列向量没关系吧?
>>=hhspectrum(imf);
>> =toimage(A,fa,tt,length(tt));
>> disp_hhs(E,tt1) %画Hilbert-Huang spectrum
>>fs=500000;%采样频率500kHz
>> for k=1:size(E,1)
bjp(k)=sum(E(k,:))*1/fs;
end
>>f=(0:(length(y)-3))/(length(y))*(fs/2);%频率
>>plot(f,bjp); %画边际谱
其实我不理解为什么按照上面求得的bjp(k)=sum(E(k,:))*1/fs;就是边际谱?难道函数toimage的功能就是排序吗?在这里,我的数据y是9000个点,得到的E是8998*8998的矩阵,难道E的每一行就代表同一个瞬时频率下的8998个瞬时幅值吗?
还有一个问题,有解释说边际谱的纵轴并不是实际幅值,而是代表某个频率出现的概率的大小,如果是这样的话,sum(bjp)是否应该等于1呀?可是我求得的值是sum(bjp)=0.0034,完全不等于1嘛,晕掉了。。。
还有,我用下面的程序求瞬时幅度-时间谱和瞬时频率-时间谱可以吗?
igure;
for k=1:length(imf) %幅度-时间谱
plot(tt,A(k,:));
hold on;
end
figure;
for k=1:length(imf) %频率-时间谱
plot(tt,fa(k,:));
hold on;
end
实际上就是把所有的瞬时幅度叠加在一张图上,把所有的瞬时频率叠加在另一张图上,这涉及到我对Hilbert-Huang Spectrum的理解,我不知道对不对。
另外,我试图用disp_hhs(E,tt1)画图时,出现错误了,??? Out of memory. Type HELP MEMORY for your options.
Error in ==> disp_hhs at 67
im = 10*log10(im/M);这是什么意思?很奇怪,好像和调用toimage时的参数length(tt)有关,如果我去掉参数length(tt),运行=toimage(A,fa,tt,);之后再运行disp_hhs(E,tt1)就可以出图,但这个时候的E变成400*8998的矩阵了。。。如果要完整的画出Hilbert-Huang Spectrum,一定要用disp_hhs吗,这画出来的图背景是蓝色的(如下图所示),横坐标是时间,纵坐标是normalized frequency,那幅度怎么表示呢?纵轴怎么转化成实际的频率呢?
file:///C:/Users/CAOWEN~1/AppData/Local/Temp/msohtml1/01/clip_image002.jpg
唉,不懂的太多,请大家不吝赐呀!谢谢了!!!
怎么么有大侠相助呀,。。 试试National Taiwan Central University中的Fast EEMD,时频图像更清晰。 shuihai707 发表于 2014-9-6 11:22
试试National Taiwan Central University中的Fast EEMD,时频图像更清晰。
嗯嗯嗯!好的,我下好了,就是还没有试,我来试试! 为何大家都不给解答呀,,郁闷,,,, 在实际信号的处理时,我也遇到了Out of memory的问题,以及画边际谱的时候f与bjp维数不匹配的问题,但在仿真信号且采样率较低时就没有这个问题,你是怎么解决的呢
页:
[1]