%Samp3_18
dt=0.02;N=512;
n=0:N-1; t=n*dt; f=n/(N*dt);   %時(shí)間序列及頻率序列
f1=3; f2=10;   %信號(hào)的頻率成分
x=.5*sin(2*pi*f1*t)+cos(2*pi*f2*t);
subplot(2,2,1), plot(t,x);    %繪制原來(lái)的信號(hào)
title('原始信號(hào)的時(shí)間域');xlabel('時(shí)間/s');
y=fft(x);      %對(duì)原信號(hào)作FFT變換
subplot(2,2,2), plot(f, abs(y)*2/N)   %繪制原信號(hào)的振幅譜
xlabel('頻率/Hz'), ylabel('振幅')
xlim([0 50]);title('原始振幅譜')
f1=8;f2=15;    %要濾去頻率的上限和下限
yy=zeros(1,length(y));  %設(shè)置與y相同元素的數(shù)組
for m=0:N-1   %將頻率落在該頻率范圍及其大于Nyquist頻率的波濾去
   if(m/(N*dt)>f1&m/(N*dt)<f2)...   %小于Nyquist頻率的濾波范圍
|(m/(N*dt)>(1/dt-f2)&m/(N*dt)<(1/dt-f1))  %大于Nyquist頻率的濾波范圍
%1/dt為一個(gè)頻率周期
      yy(m+1)=0.;       %置在此頻率范圍內(nèi)的振動(dòng)振幅為零
   else
      yy(m+1)=y(m+1);     %其余頻率范圍的振動(dòng)振幅不變
   end
end
subplot(2,2,4),plot(f,abs(yy)*2/N)    %繪制濾波后的振幅譜
xlim([0 50]);xlabel('頻率/Hz');ylabel('振幅')
gstext=sprintf('自 %4.1f -%4.1f Hz 的頻率被濾除',f1,f2); %將濾波范圍顯示作為標(biāo)題
title(gstext)
subplot(2,2,3),plot(t,real(ifft(yy)))   %繪制濾波后的數(shù)據(jù)運(yùn)用ifft變換回時(shí)間域并繪圖。
title('通過(guò)IFFT回到時(shí)間域');
xlabel('時(shí)間/s');