EMD分解后如何求周期
本帖最后由 hmjerry 于 2011-3-16 15:04 编辑大家好，看论坛里都是做信号处理的，我是跨专业，没接触过信号处理，现在想用EMD方法做周期分析，因为用的是气象资料，无法给出采样频率，现在已经得到了IMF图，想求出各个IMF中的周期，请问该如何编程？
请大家帮帮忙！
谢谢各位！
数据、代码和结果如下
clear
z=
t=1982:1:2009;
x=
c=emd(z);%EMD 分解
=size(c);
for i=1:m-1
%--------------------------------------------------------------------
%计算各IMF的方差贡献率
%定义：方差为平方的均值减去均值的平方
%均值的平方
%imfp2=mean(c(i,:),2).^2
%平方的均值
%imf2p=mean(c(i,:).^2,2)
%各个IMF的方差
mse(i)=mean(c(i,:).^2,2)-mean(c(i,:),2).^2;
end;
for i=1:m-1
mse(i)=mean(c(i,:).^2,2)-mean(c(i,:),2).^2;
%方差百分比，也就是方差贡献率
mseb(i)=mse(i)/sum(mse)*100;
%显示各个IMF的方差和贡献率
end;
figure(1)
for i=1:m-1
disp(['imf',int2str(i)]) ;disp();
end;
subplot(m+1,1,1)
plot(t,z)
set(gca,'fontname','times New Roman')
set(gca,'fontsize',12.0)
set(gca,'xtick',x);
ylabel(['signal'])
for i=1:m-1
subplot(m+1,1,i+1);
set(gcf,'color','w')
plot(t,c(i,:),'k')
set(gca,'fontname','times New Roman')
set(gca,'fontsize',12.0)
set(gca,'xtick',x);
ylabel(['imf',int2str(i)])
end
subplot(m+1,1,m+1);
set(gcf,'color','w')
plot(t,c(m,:),'k')
set(gca,'fontname','times New Roman')
set(gca,'fontsize',12.0)
set(gca,'xtick',x);
ylabel(['r',int2str(m-1)])
xlabel('Time')
这么久了怎么都没人回复呢...大家帮帮忙啊！ 学习中，帮你顶。 回复 1 # hmjerry 的帖子
水平有限但好奇, 若如LZ所说无法给出采样频率, 可以求出各个对应的周期吗?
敢问LZ的数据是否为一年一笔数据, 若是, 那一年不就是采样週期(采样频率的倒数) 回复 4 # ChaChing 的帖子
谢谢你的问题
我用的是气温数据，原数据为一天一个值，处理成年平均值，现在所用为一年一个值。
想要计算的周期是很多年的变化周期，若给出一天或一年为采样频率，则周期为一天或一年，但此周期并不能代表长年的变化周期。
有人求算出每个IMF所对应的中心频率，然后求中心频率的倒数为周期，但个人觉得还不是很准确，所为问题一直没有得到解决，若有兴趣，可以一起讨论。 我用的是轴承故障数据，滚动体有9个，出现9个冲击波形即为一个周期。 困扰啊？如何解决？
页:
[1]