已知输入输出 求传递函数 四种方法比较
看了论坛上很多关于测量传递函数的帖子:具有代表性的:
http://forum.vibunion.com/thread-110984-1-1.html
里面有H1,H2,Hv的程序
我自己找了个例子实现了一下,发现这三个方法均差与互谱除以自谱的方法。
互谱和自谱分别用的是matlab中的psd,csd来做。
有兴趣的可以运行一下:
%% 传递函数
clear
close all
num=1;
den=;
sys=tf(num,den);
%% 输入
nfft=1024;
x=randn(1,nfft);
fs=1000;
t=(1:nfft)/fs;
f1=(0:nfft/2-1)/(nfft/2)*fs/2;
%% 求输出
y=lsim(sys,x,t);
y=y';%%GXF./GXX
X=fft(x);
Y=fft(y);
X=X(1:end/2);
Y=Y(1:end/2);
Sxx=X.*conj(X);
Syx=Y.*conj(X);
Txy=Syx./Sxx;
figure;plot(f1,abs(Txy)/max(abs(Txy)));
ang=angle(Txy);
for i=1:length(ang)
if ang(i)<0
ang(i)=ang(i)+2*pi;
end
end
figure;plot(f1,ang);%% GXF./GXX
Syy=Y.*conj(Y);
Sxy=X.*conj(Y);
Txy=Syy./Sxy;
figure;plot(f1,abs(Txy)/max(abs(Txy)));
ang=angle(Txy);
for i=1:length(ang)
if ang(i)<0
ang(i)=ang(i)+2*pi;
end
end
figure;plot(f1,ang);
%%Hv
for ii=1:nfft/2;
G=;
=eig(G);
orig_lambda=diag(d);
=sort(real(orig_lambda));
lambda=orig_lambda(I);
psi=xx(:,I);
Hv(1,ii)=-psi(2,1)/psi(1,1);
end;
Txy=Hv;
figure;plot(f1,abs(Txy)/max(Txy));
ang=angle(Txy);
for i=1:length(ang)
if ang(i)<0
ang(i)=ang(i)+2*pi;
end
end
figure;plot(f1,ang);%%PSD 互谱/自谱
window=gausswin(nfft/4)';
=psd(x,nfft/4,fs,window);
=csd(x,y,nfft/4,fs,window);
Txy=Sxy./Sxx;
figure;plot(f1,abs(Txy)/max(abs(Txy)));
ang=angle(Txy);
for i=1:length(ang)
if ang(i)<0
ang(i)=ang(i)+2*pi;
end
end
figure;plot(f1,ang);
自己是第一次涉足,不对的地方还请提宝贵意见 本帖最后由 westrongmc 于 2015-1-20 16:46 编辑
sunyuxinhe 发表于 2015-1-19 21:12
自己是第一次涉足,不对的地方还请提宝贵意见
所谓的利用matlab的 PSD 互谱/自谱中,进行了四次平均,加了窗函数;
而前面的H1,H2,Hv算子中的互谱、自谱没有进行平均,没有加窗。
呵呵,不在一个条件下,可以对比吗? westrongmc 发表于 2015-1-20 13:54
所谓的利用matlab的 PSD 互谱/自谱中,进行了四次平均,加了窗函数;
而前面的H1,H2,Hv算子的中的互谱 ...
多谢所提建议!我也是刚开始接触,看来平均和窗函数确实很重要! 正好在做,只是刚做到滤波截断。过来学习了
页:
[1]