实验四 窗函数法设计FIR 数字滤波器一、实验目的1、掌握窗函数法设计FIR 数字滤波器的原理及具体方法。

2、掌握频率取样法设计FIR 数字滤波器的原理和基本方法。

3、学习利用窗函数法和频率取样法设计低通、带通、高通、带阻数字滤波器。

二、实验环境计算机、MATLAB 软件 三、实验基础理论窗函数设计FIR 滤波器 1.基本原理窗函数设计法的基本思想为，首先选择一个适当的理想的滤波器()j d H e ω，然后用窗函数截取它的单位脉冲响应(n)d h ，得到线性相位和因果的FIR 滤波器。

这种方法的重点是选择一个合适的窗函数和理想滤波器，使设计的滤波器的单位脉冲响应逼近理想滤波器的单位脉冲响应。

2.设计步骤（1）给定理想滤波器的频率响应()j d H e ω，在通带上具有单位增益和线性相位，在阻带上具有零响应。

一个带宽为()c c ωωπ<的低通滤波器由下式给定：πωωωωωωω≤<=≤=-||,0)(,||,)(c j d c ja j d e H e e H其中α为采样延迟，其作用是为了得到一个因果系统。

（2）确定这个滤波器的单位脉冲响应)())(sin()(a n a n n h c d --=πω为了得到一个(n)h 长度为N 的因果的线性相位FIR 滤波器，我们令21-=N a （3）用窗函数截取(n)d h 得到所设计FIR 数字滤波器：)()()(n R n h n h N d = 3.窗函数的选择常用的窗函数有矩形（Rectangular ）窗，汉宁（Hanning ）窗，海明（Hamming ）窗、布莱克曼（Blackman ）窗、凯瑟（Kaiser ）窗等表4-1 MATLAB 中产生窗函数的命令表4-2 常用窗函数的特性00()[]I n I ωβ⎡⎢⎣⎦=其中[]0I x 是修正的零阶贝塞尔函数，参数β控制最小阻带衰减，这种窗函数对于相同的N 可以提供不同的过渡带宽。

由于贝塞尔函数比较复杂，这种窗函数的设计方程很难推导，然而幸运的是，有一些经验设计方程可以直接使用。

已知给定的指标,,R p st p s A ωω和 ，滤波器长度N 和凯瑟窗参数β可以按如下凯瑟窗方程给出过渡带带宽：st p ωωω∆=-7.9512.285s A N ω-≈+∆0.40.1102(8.7),500.5842(21)0.07886(21),2150s s s s s A A A A A β-≥⎧⎪=⎨-+-<<⎪⎩频率取样设计FIR 滤波器 1.基本原理频率取样法从频域出发，把理想的滤波器()j d H e ω等间隔采样得到()d H k ,将()d H k 作为实际设计滤波器的()H k ：2()()()|0,1,,1j d k NH k H k H e k N ωπω====-L得到()H k 以后可以由()H k 来确定唯一确定滤波器的单位脉冲响应()h n ,()j H e ω可以由()H k 求得：10()[()]2()()()N j k h n IDFT H k H e H k k Nωπφω-===-∑其中()x φ为内插函数：12sin(/2))sin(/2)N j N e N ωωφωω--=⋅（有()H k 求得的频率响应()j H e ω将逼近()j d H e ω。

如果我们设计的是线性相位FIR 滤波器，则()H k 的幅度和相位满足线性相位滤波器的约束条件。

我们将()H k 表示为如下形式()()()=()()j k j k r H k H k e H k e θθ=当()h n 为实数，则*()()H k H N k =-由此得到()()r r H k H N k =-即()/2r H k k N =以为中心偶对称。

在利用线性相位条件可知，对于1型和2型线性相位滤波器：121()0,,22()121()()1,122N kN k N k N N N k k N N πθπ⎧--⎢⎥-=⎪⎢⎥⎪⎣⎦=⎨--⎢⎥⎪-=+-⎢⎥⎪⎣⎦⎩L L对于3型和4型线性相位滤波器121()0,,222()121()()1,1222N kN k N k N N N k k N N ππθππ⎧--⎢⎥±-=⎪⎢⎥⎪⎣⎦=⎨--⎢⎥⎪-=+-⎢⎥⎪⎣⎦⎩L m L2.设计步骤（1）由给定的理想滤波器给出()r H k 和()k θ。

（2）由()()()=()()j k j k r H k H k eH k e θθ=求得()H k（3）根据()H k 求得()h n 或()j H e ω四、实验内容1、设计一个数字低通FIR 滤波器，其技术指标如下：0.2,0.25p p R dB ωπ==0.3,50st s A dB ωπ==分别采用矩形窗、汉宁窗、海明窗、布莱克曼窗、凯瑟窗设计该滤波器。

结合实验结果，分别讨论采用上述方法设计的数字滤波器是否都能满足给定指标要求。

（1）矩形窗 程序代码：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp; N=ceil(1.8*pi/tr_width) n=0:N-1;wc=(wst+wp)/2; alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha)); w_boxcar=boxcar(N)'; h=hd.*w_boxcar; subplot(221);stem(n,hd,'filled');axis tight ;xlabel('n');ylabel('hd(n)'); [Hr,w1]=zerophase(h); subplot(222); plot(w1/pi,Hr);axis;xlabel('\omega/\pi');ylabel('H(\omega)'); subplot(223);stem(n,h,'filled');axis tight ;xlabel('n');ylabel('h(n)'); [H,w]=freqz(h,1); subplot(224);plot(w/pi,20*log10(abs(H)/max(H)));axis tight ;xlabel('\omega/\pi');ylabel('dB'); grid on ;MATLAB 图形：（2）汉宁窗 程序代码：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp; N=ceil(6.2*pi/tr_width) n=0:N-1;wc=(wst+wp)/2; alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha)); w_boxcar=hanning(N)'; h=hd.*w_boxcar; subplot(221);stem(n,hd,'filled');axis tight ;xlabel('n');ylabel('hd(n)'); [Hr,w1]=zerophase(h); subplot(222); plot(w1/pi,Hr);axis;xlabel('\omega/\pi');ylabel('H(\omega)'); subplot(223);stem(n,h,'filled');axis tight ;xlabel('n');ylabel('h(n)'); [H,w]=freqz(h,1); subplot(224);plot(w/pi,20*log10(abs(H)/max(H)));n h d (n )0.51-0.500.511.5ω/πH (ω)n h (n)0.20.40.60.8-80-60-40-20ω/πd Baxis tight ;xlabel('\omega/\pi');ylabel('dB'); grid on ; MATLAB 图形：（3）海明窗 程序代码：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp; N=ceil(6.6*pi/tr_width) n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha)); w_boxcar=hamming(N)'; h=hd.*w_boxcar; subplot(221);stem(n,hd,'filled');axis tight ;xlabel('n');ylabel('hd(n)'); [Hr,w1]=zerophase(h); subplot(222); plot(w1/pi,Hr);axis;xlabel('\omega/\pi');ylabel('H(\omega)'); subplot(223);stem(n,h,'filled'); axis tight ; xlabel('n');nh d (n )0.51-0.500.511.5ω/πH (ω)nh (n )0.20.40.60.8-120-100-80-60-40-200ω/πd Bylabel('h(n)'); [H,w]=freqz(h,1); subplot(224);plot(w/pi,20*log10(abs(H)/max(H)));axis tight ;xlabel('\omega/\pi');ylabel('dB'); grid on ; MATLAB 图形：（4）布莱克曼窗 程序代码：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp; N=ceil(11*pi/tr_width) n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha)); w_boxcar=blackman(N)'; h=hd.*w_boxcar; subplot(221);stem(n,hd,'filled');axis tight ;xlabel('n');ylabel('hd(n)'); [Hr,w1]=zerophase(h); subplot(222); plot(w1/pi,Hr);axis;xlabel('\omega/\pi');ylabel('H(\omega)'); subplot(223);stem(n,h,'filled');axis tight ;xlabel('n');ylabel('h(n)');nh d (n )-0.500.511.5ω/πH (ω)nh (n )0.20.40.60.8-100-50ω/πd B[H,w]=freqz(h,1); subplot(224);plot(w/pi,20*log10(abs(H)/max(H)));axis tight ;xlabel('\omega/\pi');ylabel('dB'); grid on ;MATLAB 图形为：（5）凯瑟窗 程序代码：wp=0.2*pi;wst=0.3*pi;tr_width=wst-wp;As=50; N=ceil((As-7.95)/(2.285*tr_width))+1; beta=0.1102*(As-8.7); n=0:N-1;wc=(wst+wp)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha)); w_boxcar=kaiser(N,beta)'; h=hd.*w_boxcar; subplot(221);stem(n,hd,'filled');axis tight ;xlabel('n');ylabel('hd(n)'); [Hr,w1]=zerophase(h); subplot(222); plot(w1/pi,Hr);axis;xlabel('\omega/\pi');ylabel('H(\omega)'); subplot(223);stem(n,h,'filled');nh d (n )-0.500.511.5ω/πH (ω)nh (n )-150-100-50ω/πd Baxis tight ;xlabel('n');ylabel('h(n)'); [H,w]=freqz(h,1); subplot(224);plot(w/pi,20*log10(abs(H)/max(H)));axis tight ;xlabel('\omega/\pi');ylabel('dB'); grid on ; MATLAB 图形：2、设计一个数字带通FIR 滤波器，其技术指标如下： 下阻带边缘：10.2,60st s A dB ωπ== 下通带边缘：10.35,1p p R dB ωπ== 上通带边缘：20.65,1p p R dB ωπ== 上阻带边缘：20.8,60st s A dBωπ==程序代码：wp1=0.2*pi;Rp1=1; wst1=0.35*pi;A1=60; width1=wst1-wp1;N1=ceil(11*pi/width1)+1; n1=0:(N1-1);wc1=(wp1+wst1)/2; alpha=(N1-1)/2;wp2=0.65*pi;Rp2=1;wst2=0.8*pi;A2=60;nh d (n )-0.50.511.5ω/πH (ω)nh (n )-100-50ω/πd Bwidth2=wst2-wp2;N2=ceil(11*pi/width2)+1; n2=0:(N2-1);wc2=(wp2+wst2)/2; alpha=(N2-1)/2;hd=(wc2/pi)*sinc((wc2/pi)*(n2-alpha))-(wc1/pi)*sinc((wc1/pi)*(n1-alpha)); w_w=blackman(N1)'; h=hd.*w_w; subplot(221);stem(n1,h,'filled'); subplot(222); [H,w]=freqz(h,1);plot(w/pi,20*log10(abs(H)/max(abs(H)))); subplot(223);[Hr,w1]=zerophase(h); plot(w1/pi,Hr); subplot(224);stem(n1,hd,'filled'); [Hr,wl]=zerophase(h); grid on ; MATLAB 图形：3.采用频率取样法设计FIR 数字低通滤波器，满足以下指标2040608000.51-150-100-500.51-0.500.511.50204060800.2,0.250.3,50p p st s R dB A dBωπωπ====（1）取N=20，过渡带没有样本。