力锤激励试验中怎么用matlab实现集中平均法？

zyy2011

本人证在用matlab处理力锤激励试验的数据，要求识别出振型和频率。
采用的是多点激励单点响应法，每个激励点锤击4次，请问：
1. 每个点得到的频响函数怎么用matlab实现集中平均，能否提供集中平均法的程序？
2. 得到的多个点频响后，是用每个点频响数据进行识别固有频率和振型？还是先用什么方法合成一个综合的频响，再进行频率和振型的识别？
欢迎赐教！

zyy2011

另外，还有为啥得到的谱图只有到25Hz？
采用频率为100.16Hz，FFT长度1024

程序如下：

%随机信号谱分析
%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
close all hidden
format long
%%%%%%%%%%%%%%%%%%%%%%
fni=input('随机信号谱分析-输入数据文件名:','s');
fid=fopen(fni,'r');
fun = fscanf(fid,'%d',1);     %分析内容(1自谱,2互谱,3频响函数,4相干函数)
sf = fscanf(fid,'%f',1);      %采样频率
nfft = fscanf(fid,'%d',1);    %FFT长度
win = fscanf(fid,'%d',1);     %窗函数(1矩形,2汉宁,3海明,4布莱克曼,5三角)
fg = fscanf(fid,'%d',1);      %图形方式(1实部、虚部图,2幅值、相位图)
fno = fscanf(fid,'%s',1);     %输出数据文件名
a = fscanf(fid,'%f',[2,inf]); %按列读入数据,第1列激励,第2列响应
x=a(1,:);                     
y=a(2,:)-a(1,:);
status=fclose(fid);
f=0:sf/nfft:sf/2-sf/nfft;
%加窗
switch win
  case 1  %矩形窗
    w=boxcar(nfft);
  case 2  %汉宁窗
    w=hanning(nfft);
  case 3  %海明窗
    w=hamming(nfft);
  case 4  %布莱克曼窗
    w=blackman(nfft);
  case 5  %三角窗
    w=triang(nfft);
  otherwise   
    w=boxcar(nfft);
end
switch fun
  case 1  %计算自谱
    z=psd(y,nfft,sf,w,nfft/2);
  case 2  %计算互谱
    z=csd(x,y,nfft,sf,w,nfft/2);
  case 3  %计算频响函数
    z=tfe(x,y,nfft,sf,w,nfft/2);
  case 4  %计算相干函数
    z=cohere(x,y,nfft,sf,w,nfft/2);
  otherwise   
   ;
end
%绘制曲线图
nn=1:nfft/4;
  subplot(2,1,1);
if fg==1
  plot(f(nn),real(z(nn)));
  ylabel('实部');             %添加纵向坐标轴的标注   
else
  plot(f(nn),abs(z(nn)));
  ylabel('幅值');
end  
xlabel('频率 (Hz)');          %添加横向坐标轴的标注
grid on;                      %在图幅上添加坐标网格
if fun>1&fun<4
  subplot(2,1,2);
  if fg==1
    plot(f(nn),imag(z(nn)));
    ylabel('虚部');           %添加纵向坐标轴的标注   
  else
    plot(f(nn),angle(z(nn)));
    ylabel('相位');
  end  
  xlabel('频率 (Hz)');        %添加横向坐标轴的标注
  grid on;                    %在图幅上添加坐标网格
end
fid=fopen(fno,'w');           %以写的方式打开文件或建立一个新文件
if fun==1|fun==4              %自谱和相干函数
  for k=1:nfft/2
    fprintf(fid,'%f        %f\r\n',f(k),abs(z(k))); %每行输出2个实型数据，f为频率，z为幅值
  end
else                          %互谱和频响函数
  for k=1:nfft/2
    fprintf(fid,'%f        %f        %f\r\n',f(k),real(z(k)),imag(z(k))); %每行输出3个实型数据
  end
end
status=fclose(fid);