fir_kai.m

利用凯泽窗产生FIR滤波器

function [kai_hd,kai_h,w_kai,kai_db,kai_ap,kai_as]=FIR_kai(fp,fs,ap,as,Fs);

   clear; 
   clc;
   Fs=10000;
   fp=2000;
   fs=2500;
   ap=0.005;
   as=0.005;
   Wp=2*pi*fp/Fs;
   Ws=2*pi*fs/Fs;
   As=-20*log10(as);
   tr_width=Ws-Wp;   %计算过度带宽 
   M=ceil((As-7.95)*2*pi/(2.286*tr_width))+1;%按凯泽窗计算滤波器长度 
   disp(['滤波器的长度',num2str(M)]); 
   beta=0.1102*(As-8.7); 
   n=[0:1:M-1]; 
   disp(['线性相位斜率为',num2str(beta)]); 
   Wc=(Ws+Wp)/2;    %求理想低通滤波器的截止频率
   kai_hd =ideal_lp(Wc,M); 
   w_kai=(kaiser(M,beta))'; %求凯泽窗窗函数
   kai_h= kai_hd.*w_kai; 
   [kai_db,kai_mag,kai_w]=freqz_m(kai_h,[1]); 
   delta_w=2*pi/1000; 
    
   kai_ap=-min(kai_db(1:1:Wp/delta_w+1)); 
   disp(['实际通带波动为kai_ap=',num2str(1-10^(-kai_ap/20))]); 
   kai_as=-round(max(kai_db(Ws/delta_w+1:1:501))); 
   disp(['最小阻带衰减为kai_as=',num2str(10^(-kai_as/20))]); 
   
   figure(1);
   subplot(2,1,1); 
   stem(n,kai_hd); 
   axis([0 M-1 -0.1 0.2]); 
   title('理想脉冲响应'); 
   ylabel('hd（n）'); 
   subplot(2,1,2); 
   stem(n,kai_h); 
   title('实际脉冲响应'); 
   axis([0 M-1 -0.1 0.2]); 
   xlabel('n'); 
   ylabel('h(n)'); 
   
   figure(2);
   subplot(2,1,1); 
   stem(n,w_kai); 
   axis([0 M-1 0 1.1]); 
   title('凯泽窗'); 
   ylabel('wd(n)'); 
   subplot(2,1,2); 
   plot(kai_w/pi,kai_db); 
   title('幅度响应/db'); 
   grid; 
   axis([0 1 -110 10]);  
   xlabel('以pi为单位的频率'); 
   ylabel('分贝数/db') 
   
   function[db,mag,w]=freqz_m(b,a);  
 %db=[0到pi]区间内的相对振幅；  
 %mag=[0:pi]区间内的绝对振幅；  
 %pha=[0:pi]区间内的相位响应；  
 %grd=[0:pi]区间内的群延迟；  
 %w=[0;pi]区间内的501个频率样本向量;  
 %b=直接型滤波器分子多项式系数(对FIR:b=h);  
 %a=直接型滤波器分母多项式系数(对FIR:a=[1]);  
  
 [H,w]=freqz(b,a,1000,'whole');  
 H=(H(1:1:501))';  
 w=(w(1:1:501))';  
 mag=abs(H);  
 db=20*log10((mag+eps)/max(mag));  

  
function hd=ideal_lp(Wc,M);%[hd]=ideal_lp(Wc,M)；  
 %hd=0 to M-1之间的理想脉冲响应；  
 %Wc=截止频率(弧度);  
 %M=理想滤波器的长度;  
alpha=(M-1)/2;  
n=[0:1:(M-1)];  
m=n-alpha+eps; %加一个小数以避免零作除数;  
hd=sin(Wc*m)./(pi*m);