[求助]求取频响函数并分析的问题!!
简单的说就是把采集到的信号经过FFT后(得到Xi,Fi),再进行自谱和互谱计算(GFX=Xi*Fi',GFF=Fi*Fi')
频响函数=GFX/GFF
最后的步骤是截取频响函数的前半部分
请问怎样编程才能实现截取前半部分?下面是上面说的简单的编城,有不少错误哈,高手忍一下,下面要怎样才能实现----截取前半部分????
clc
N=1000;
n=0:N-1;
Xk=fft(xn); %xn是采集到的响应信号
Fk=fft(fn); %fn是采集到的激励信号
GFX=(Fk').'*Xk;
FK=(Fk').'
GFX=Xk.*FK;
GFF=Fk.*FK;
TXF=GFX./GFF;
plot(,); %画这个幅频曲线要怎么画?
谢谢各位了
[ 本帖最后由 happy 于 2006-11-6 18:10 编辑 ] clc
close
Fs = 1000; % Sampling frequency
T = 1/Fs; % Sampling time
x; % x is excitation signal
y; % y is response signal
L = Fs*t_time; % Length of signal
NFFT = 2^nextpow2(L); % Next power of 2 from length of signal
X = fft(x,NFFT)/FsL; % fft translation of x
Y = fft(y,NFFT)/FsL; % fft translation of y
H = (Y.*conj(X))./(X.*conj(X)); % FRF
f = Fs/2*linspace(0,1,NFFT/2); % x axis
plot(f,abs(H(1:NFFT/2))); % plot FRF
如有错误，或更好的方法，欢迎讨论
skynew2005@gmail.com
[ 本帖最后由 happy 于 2006-11-6 18:10 编辑 ]
如果有多组数据
% 如果要处理很多数据，就要把数据分段，分别作变换后，平均。%
% 比如，要做400Hz内的频响，分辨率为1600线，即需要800Hz的采样频率，4分钟的测量数据，共3200个测量点。
% 用函数来实现作频响。
clc
close
clear
t_time = 10;
Fs = 800;
= sim('tt_1',);
x = y_out(:,1); % excitation data
y = y_out(:,2); % response data
span = Fs/2; % span of fft
lines = 6400; % lines of fft
L_data = lines*2; % 每一段的数据长度
L_n = round(length(y)/L_data); % 总数据共分段数
for kk = 1:L_n
x_s = x(L_data*(kk-1)+1:min(L_data*kk,length(x)));
y_s = y(L_data*(kk-1)+1:min(L_data*kk,length(y)));
= cal_data_FRF(x_s,y_s,span,lines);
end
for kk = 1:length(H_s(1,:))
H(kk) = mean(H_s(:,kk));
end
plot(f,H)
function = cal_data_FRF(x,y,span,lines)
% calculate the frequency response function of the experiment data
% x, is the excitation signal
% y, is the response signal
% span, is the span of frequency of FFT ,for exzample 400,800,1600...
% lines, is the resolving of the FFT
Fs = span*2; % Sampling frequency
T = 1/Fs; % Sampling time
L = lines*2; % Length of signal
NFFT = 2^nextpow2(L); % Next power of 2 from length of signal
X = fft(x,NFFT)/L; % fft translation of x
Y = fft(y,NFFT)/L; % fft translation of y
HH = (Y.*conj(X))./(X.*conj(X)); % FRF
f = Fs/2*linspace(0,1,NFFT/2); % x axis
H = abs(HH(1:NFFT/2)); % y axis
有时候做出的曲线不光滑，大家一般怎么处理的？
欢迎讨论
skynew2005@gmail.com
现在找到一种更好的方法
老版本的命令是tfe新版本的命令是tfestimate
直接得出频响
计算效果很好
页:
[1]