声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 900|回复: 0

[HHT] 大家帮帮忙吧,我检查不出错误啊?

[复制链接]
发表于 2010-4-13 16:21 | 显示全部楼层 |阅读模式

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

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

x
我写了个程序关于插入与移除masking 信号的,检查了两天了还是没有检查出错误,请大家帮帮忙吧。load fivemeter.dat
t=fivemeter(:,1);
x=fivemeter(:,2);
plot(t,x)
title('signal')
xlabel('time(s)')
ylabel('wavesurface')
imf=emd(x);
k=size(imf,1);
for i=1:k
subplot(k+1,1,i)
plot(t,imf(i,:))
end
subplot(k+1,1,k+1)
plot(t,x,'r')
title('signal')
n=4096;
h=hilbert(imf(1,:));
%for j=1:k-1
%for kk=1:n
%h(kk,j)=hilbert((imf(j,:))',t)
%end
%end
fs=1/0.0247;
for kk=1:n
    a=((real(h)).^2+(imag(h)).^2).^0.5;
    thta=atan(imag(h)./real(h));
    f(kk)=(real(h(kk)).*(imag(kk+1)-imag(kk-1))-imag(kk).*(real(kk+1)-real(kk-1)))./2.*0.0247.*((real(kk)).^2+(imag(kk)).^2);
    fx(kk)=fs.*asin(f(kk)./fs);
end

f1=0;
for kk=1:n
f1=f1+a(kk).*((fx(kk)).^2)
kk=kk+1
end
f1end=f1(1,4097)
f2=0;
for kk=1:n
f2=f2+a(kk).*fx(kk)
kk=kk+1
end
f2end=f2(1,4097);
f=f1end./f2end;
fs=1./0.0247;
for kk=1:n
s(kk)=1.6.*sin(2.*pi.*f.*kk/fs);
x1(kk)=fivemeter(2,:)+s(kk)
x2(kk)=fivemeter(2,:)-s(kk)
end
%kk=1:n;
%kk=0:1/fs:n/fs;
%s(kk)=s(kk);
%x1(kk)=x1(kk);
%x2(kk)=x2(kk);
imf1=emd(x1);
imf2=emd(x2);
m=size(imf1,1);
n=size(imf2,1);
mink=min(m,n);
for i=1:mink
imf(i,:)=(imf1(i,:)+imf2(i,:))./2
end
for j=1:mink
subplot(mink+3,1,j)
plot(t,imf(j,:))
end
subplot(mink+3,1,mink+1)
plot(t,x,'r');
title('original signal')
subplot(mink+3,1,mink+2)
plot(t,x1,'g')
title('masking signal added')
subplot(mink+3,1,mink+3)
plot(t,x2,'m')
title('masking signal removed')
%以上是利用masking信号进行emd分解
一运行就显示??? Index exceeds matrix dimensions.并且很奇怪的是,所有的f1都为0,这是不可能的啊。所以请高手帮帮忙吧
回复
分享到:

使用道具 举报

您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-11-11 15:22 , Processed in 0.060571 second(s), 18 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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