声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 2228|回复: 3

[HHT] g。rilling 的emd_separation运行错误,求指导。

[复制链接]
发表于 2017-3-2 15:37 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
新人一枚。。最近学习g.rilling的emd程序在emd_separation,比例设在4和1时都可以运行出论文中的结果。如下:
4.png 振幅比为4 1.png 振幅比为1
可是当比例设置为0.25时,运行不出论文中的结果,matlab显示:
Warning: Divide by zero.
> In emd>stop_sifting at 350
  In emd at 183
  In emd_separation at 42
在实验比例为0.3时,可以运行出结果,如图: 0.3.png
0.25就运行不出来了。。求问大神。。为啥呢

本帖被以下淘专辑推荐:

回复
分享到:

使用道具 举报

 楼主| 发表于 2017-3-2 16:07 | 显示全部楼层
程序如下
% EMD_SEPARATION.M
%
% P. Flandrin, Mar. 13, 2003 - modified Mar. 2, 2006
%
% computes and displays an error measure in the EMD
% separation of two tones
%
% produces Figure 4 in
%
% G. Rilling, P. Flandrin and P. Gon峚lv弒
% "On Empirical Mode Decomposition and its algorithms"
% IEEE-EURASIP Workshop on Nonlinear Signal and Image Processing
% NSIP-03, Grado (I), June 2003

N = 256;% # of data samples
t = 1:N;
tt = N/4:3*N/4;

Nf = 129;% # of tested fequencies
f = linspace(0,.5,Nf);

rapp = 0.29;% amplitude ratio between modes
a1 = sqrt(rapp);
a2 = 1/a1;

x1 = a1*cos(2*pi*f'*t);
x2 = a2*cos(2*pi*f'*t);

se = zeros(Nf);

for k1 = 5:Nf-1
       
        for k2 = 2:k1-1;
               
                y1 = x1(k1,:);
                y2 = x2(k2,:);
               
                sy1 = sum((y1(tt)).^2);
                sy2 = sum((y2(tt)).^2);
                sy = sum((y1(tt)).^2+(y2(tt)).^2);
       
                imf = emd(y1+y2);
        if size(imf,1) == 1
            imf(2,:) = 0;
        end
                se(k1,k2) = sqrt((sy1*sum((imf(1,tt)-y1(tt)).^2) + sy2*sum((imf(2,tt)-y2(tt)).^2))/(sy1+sy2)/sy);       

                [k1 k2 size(imf)]
                               
        end

end

imagesc(f,f,flipud(se'))
axis('square')
axis([0 .5 0 .5])
hold on
plot([f(2) f(Nf-1)],[f(Nf-1) f(Nf-1)])
plot([f(Nf-1) f(Nf-1)],[f(2) f(Nf-1)])
plot([f(2) f(Nf-1)],[f(Nf-1) f(2)])
set(gca,'YTick',[]);set(gca,'XTick',[])
xlabel('f_1 > f_2')
ylabel('f_2')
colormap(jet)
%colormap(flipud(gray))
colorbar('vert')
hold off
发表于 2017-3-3 08:55 | 显示全部楼层
你这几张图的结果都正确吗? 你这个比例我没理解?
 楼主| 发表于 2017-3-3 10:01 | 显示全部楼层
sovereign 发表于 2017-3-3 08:55
你这几张图的结果都正确吗? 你这个比例我没理解?

我是对照着g.riling 03年的论文中的结果运行的。。。论文中设置了三个比例4,1,0.25。。我运行之后,只有在4和1的情况下才能算出结果。0.25就提示错误。。
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-10 17:00 , Processed in 0.082344 second(s), 23 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表