【求助】位移转加速度(谱转换法)

那啥,小弟初学matlab,比较菜

有这样一组数据(样本见附件)

0.3987 1.3757 -5.0000 -0.8385

0.3984 1.3759 -5.0000 -0.8220

0.3989 1.3765 -5.0000 -0.7935

0.3986 1.3760 -5.0000 -0.7515

0.3984 1.3759 -5.0000 -0.7213

0.3976 1.3756 -5.0000 -0.7024

0.3973 1.3753 -5.0000 -0.6691

0.3978 1.3756 -5.0000 -0.6296

0.3969 1.3754 -5.0000 -0.5919

... .... ... ...

说明:

1.采样速率20K，k=1000

2.文件长度22K

3.第1，2 列激光位移，第3列是干扰,(忽略)第4列加速度

现在要求用频谱转换法,将加速度转为位移

算法是这个样子的:

如果在时间T内采集N个数据,

x(n)经离散傅里叶变换后得到的X(k)是一个长度为N的复数序列(频谱),

X(k)=DFT[x(n)]=[(a0，jb0)，(a1，jb1 )，...,(an-1 ，jbn-1 )]

x(n)中各谐波分量的幅值、圆频率及初相角可由式(1)~(3)求出

Ak=sqrt(ak^2+bk^2) (1)

ωk=2*pi*k/T (2)

φ=atan(bk/ak) (3)

根据信号叠加原理，任何周期信号是若干简谐信号的叠加。加速度信号可用式(4)表示，其相应的位移用式(5)表示，两者幅值和相位之间的关系如式(6)

a=Aa0cos(ω0t+φa0t)+Aa0cos(ω0t+φa0t)+...+Aan-1cos(ωn-1t+φan-1t) (4)

d=Ad0cos(ω0t+φd0t)+Aa0cos(ω0t+φd0t)+...+Aan-1cos(ωn-1t+φdn-1t) (5)

Adi=Aai/ωi^2 , φdi=φai-pi

下面是我自己用matlab写该算法的实现代码,比较臭

%读取数据

a='d:\1.txt';

data=load(a);

[m,n]=size(data);

%提取加速度,

acce=data(:,4);

XK=fft(acce);

re=real(XK);

im=imag(XK);

%AK=zeros(m,1);

T=1;

for i=1:m

%计算幅值

AK_A(i,1)=sqrt(re(i,1)^2+im(i,1)^2);

%计算初相

PHY_A(i,1)=atan(im(i,1)/re(i,1));

%计算角频率

OMEGA(i,1)=2*pi*i/T;

end

for i=1:m

%计算幅值

AK_D(i,1)=AK_A(i,1)/(OMEGA(i,1)^2);

%计算初相

PHY_D(i,1)=PHY_A(i,1)-pi;

end

%计算位移

%%%%%%%%%%%%%%%%%%

t=1/20000;

%%%%%%%%%%%%%%%%%%

D=zeros(m,1);

D(i,1)=AK_D(i,1)*cos(OMEGA(i,1)*t+PHY_D(i,1));

for i=2:m

D(i,1)=D(i-1,1)+AK_D(i,1)*cos(OMEGA(i,1)*t+PHY_D(i,1));

end

plot(D)

现在的问题是:

1 通过该算法得到的位移D与数据里给出的两列位移都不匹配 (说实话我不太清楚两列位移各表示什么,但我算出来的D却和这两列数据毫无关系)

2.我的算法实现是不是存在问题,如果是,那我应该如何得到位移呢?

希望达人能帮帮忙,谢谢啦

上传的附件 1.txt (15.6 KB, 7 次查看)