fre_sample_fir.m

经典《信号与系统》教程的matlab例程,对深入理解信号与系统相关概念有很大帮助

% Name:fre_sample_fir.m
clear,close all;
N=input('Type in the length N=');
w=-pi:2*pi/N:pi-2*pi/N;
%Hd=u(w+(2*pi/N)*(N/4))-u(w-(2*pi/N)*(N/4));
Hd=u(w+(2*pi/N)*(N/4))-u(w-(2*pi/N)*(N/4));
w1=0:2*pi/N:2*pi-2*pi/N;
H=u(w1)-u(w1-(2*pi/N)*(N/4));
k=0:N-1;
phai=-pi*k*(N-1)/N;
H(10)=0.5;H(11)=0.1;
H=H.*exp(i*phai);
kc=fix(N/2);
for i=0:kc
    H(N-(i))=conj(H(i+2));
end
%Hd(10)=0.5;Hd(11)=0.1;

subplot(221)
stem(w1/pi,abs(H),'.');grid on
axis([0,1,0,1])
title('The sampled frequency response of Hd')
h=ifft(H,N);
subplot(222)
n=0:N-1;
stem(n,h,'.');grid on
title('The impulse response h(n) of an FIR filter')
axis([0,N-1,min(h),max(h)])
k=0:N-1;
M=exp(-j*w*(N-1)/2);
%L=length(w)
%
subplot(223)
W=-pi:0.01:pi;
fa1=0;
for i=0:N-1
    if kc==N/2
        m=-(N+2)/2+i;
    else
        m=-(N+2)/2+i+1;
    end
    fa=abs(H(i+1))*sin(N*((W-2*pi*m/N)/2+eps))./sin((W-2*pi*m/N)/2+eps);
    fa1=fa1+abs(H(i+1))*sin(N*((W-2*pi*m/N)/2+eps))./sin((W-2*pi*m/N)/2+eps);
end
fa1=(1/N)*fa1;M=length(fa1);
fa2=zeros(1,M);
for i=0:ceil(M/2);
    fa2(i+ceil(M/2)-1)=fa1(M-i);
    %fa2(M-i)=fa1(i+1);
end
plot(W/pi,(fa2));grid on
axis([0,1,min(fa1),max(fa1)]);hold on
H=fft(h,N);%./exp(-i*w*(N-1)/2);
plot(k*2/N,abs(H),'r.');grid on
axis([0,1,min(fa1),max(fa1)]);hold on

subplot(224)
n=0:N-1;
stem(n,abs(H),'.');grid on
title('The Hg(k) of an FIR filter')
axis([0,N-1,min(abs(H)),max(abs(H))])
%==========================================
figure(2)
plot(W/pi,20*log10(abs(fa2)/max(abs(fa2))));grid on
axis([0,1,-80,0]);