fftfortimedata.m

来自「数字信号处理用的快速傅里叶变化」· M 代码 · 共 36 行

fprintf('START\n');
clear;
%****changable parameters when run each time******
  fin1=fopen('8_1-5.txt','r');%txt file: time  data
  period=1;  %time period number of one Time Series Datas,each period is FFTed independently. 
  start_at=2;%the first data to be ffted 
  t_point=2048;%number of data at time domain within each time period.
%*************************************************
figure(1);
for i=1:1:start_at-1
  line=fgetl(fin1);
end
for k=1:1:period
  for i=1:1:t_point
    t(i)=fscanf(fin1,'%g',1);
    data(i)=fscanf(fin1,'%g', 1);
  end
  delta_t=t(2)-t(1);
  smp_frq=1/delta_t;
  data=data-mean(data);

  Y = fft(data,t_point);
  Pyy = Y.* conj(Y) / t_point;
  f = smp_frq*(0:t_point/2)/t_point;
  subplot(period,1,k);
  %subplot(2,1,1);
  pmax=max(Pyy(1:t_point/2+1));
  Pyy(1:t_point/2+1)=Pyy(1:t_point/2+1)/pmax;
  plot(f,Pyy(1:t_point/2+1))
  axis([0 50 0 1.2]);
  set(gca,'fontname','times new roman','fontsize',10);%for axis
  xlabel('frequency (Hz)','fontname','times new roman','fontsize',12)%,'fontweight','bold');
  ylabel('Spectrum','fontname','times new roman','fontsize',12)%,'fontweight','bold');
end%k
fprintf('OK!\n');