[求助]求取频响函数并分析的问题!!

simomo

简单的说就是把采集到的信号经过FFT后(得到Xi,Fi),
再进行自谱和互谱计算(GFX=Xi*Fi',GFF=Fi*Fi')
频响函数=GFX/GFF

最后的步骤是截取频响函数的前半部分

请问怎样编程才能实现截取前半部分?下面是上面说的简单的编城,有不少错误哈,高手忍一下,下面要怎样才能实现----截取前半部分????

1. clc
   
2. N=1000;
   
3. n=0:N-1;
   
4. Xk=fft(xn); %xn是采集到的响应信号
   
5. Fk=fft(fn); %fn是采集到的激励信号
   
6. GFX=(Fk').'*Xk;
   
7. FK=(Fk').'
   
8. GFX=Xk.*FK;
   
9. GFF=Fk.*FK;
   
10. TXF=GFX./GFF;
    
11. plot(,); %画这个幅频曲线要怎么画?

复制代码

谢谢各位了

[ 本帖最后由 happy 于 2006-11-6 18:10 编辑 ]

skynew

1. clc
   
2. close
   
3. 
   
4. Fs = 1000;          % Sampling frequency
   
5. T = 1/Fs;             % Sampling time
   
6. x;                         % x is excitation signal
   
7. y;                         % y is response signal
   
8. L = Fs*t_time;     % Length of signal
   
9. NFFT = 2^nextpow2(L);     % Next power of 2 from length of signal
   
10. X = fft(x,NFFT)/FsL;           % fft translation of x
    
11. Y = fft(y,NFFT)/FsL;          % fft translation of y
    
12. 
    
13. H = (Y.*conj(X))./(X.*conj(X));     % FRF
    
14. f = Fs/2*linspace(0,1,NFFT/2);     % x axis
    
15. plot(f,abs(H(1:NFFT/2)));              % plot FRF

复制代码

如有错误，或更好的方法，欢迎讨论
skynew2005@gmail.com

[ 本帖最后由 happy 于 2006-11-6 18:10 编辑 ]

skynew

% 如果要处理很多数据，就要把数据分段，分别作变换后，平均。
%
% 比如，要做400Hz内的频响，分辨率为1600线，即需要800Hz的采样频率，4分钟的测量数据，共3200个测量点。
% 用函数来实现作频响。

clc
close
clear

t_time = 10;
Fs = 800;

[t,x_state,y_out] = sim('tt_1',[0,t_time]);

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)));
    [f,H_s(kk,:)] = 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 [f,H] = 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

skynew

老版本的命令是tfe
新版本的命令是tfestimate
直接得出频响
计算效果很好