FIR数字滤波器的Matlab实现.pdf
第第7 FIR章章数字滤波器的数字滤波器的Matlab实现实现7.1 实验目的学习用窗函数法设计 FIR 数字滤波器的原理及其设计步骤；学习编写数字滤波器的设计程序的方法，并能进行正确编程；根据给定的滤波器指标，给出设计步骤。7.2 实验原理及实例分析FIR 低通数字滤波器的设计原理如果系统的冲激响应hd(n)为，那么系统的输入/输出关系为：y(n)x(n)hd(n)对于低通滤波器，只要设计出低通滤波器的冲激响应函数，就可以由上式得到系统的输出了。假设所希望的数字滤波器的频率响应为Hd(ejw)，它是频域的周期函数，周期为 2，那么它与Hd(ejw)相对应的傅立叶系数为hd(n)-1以hd(n)为冲激响应的数字滤波器将具有频域响但是将hd(n)作为滤波器脉冲响应有两个问题：(1)(2)Hd(ejw)ejnw dwHd(ejw)。它是无限长的，与 FIP 滤波器脉冲响应有限长这一前提不一致它是非因果的，hd(n)0,n 0对此，要采取以下的措施，(1)将hd(n)截短2将其往右平移，N 1由此得到h2(n)的实际频域响应Hd(ejw)h2(n)ejnw，与理想频域响应Hd(ejw)相近，n 0但不完全一致。理论证明上述现象是对hd(n)进行简单截短处理的必然结果，一般称为吉布斯现象，为尽可能的减少吉布斯现象，应对hd(n)进行加窗截取，即以h(n)hd(n)WN(n)作为 FIR滤波器的系数。常用的窗函数有矩形窗、海明窗和布莱克曼窗等。用窗函数法设计 FIR 滤波器Matlab 设计 FIR 滤波器有多种方法和对应的函数，见表方法窗方法多带方法最小一乘法任意响应法描述使用窗函数和逆傅立叶变换实现包含子带频率域使用最小二乘法将整个频率域上的错误 几率压缩到最小使用任意响应，包括非线性相位以及复滤 波器firrcos 等余弦法使用三角函数的低通响应窗函数方法不仅在数字滤波器的设计中占有重要的地位，同时可以用于功率谱的估计，从根本上讲，使用窗函数的目的就是消除由无限序列的截短而引起的响。窗函数设计线性相位 FIR 滤波器步骤如下：(1)(2)确定数字滤波器的性能要求，临界频率根据性能要求，合理选择单位脉冲响应响应Hd(ejw)的幅频特性和相频特性(3)求理想单位脉冲响应hd(n)，在实际计算中，可对Hd(ejw)采样，并对其求 IDFT 的hMGibbs 现象所带来的影cremez 等7-1。表 7 1 matlab 设计 FIR 滤波器的方法和函数函数fir1,fir2,kaiserord 等firls,remez 等fircls,fircls1 等wk，滤波器单位脉冲响应长度Nh(n)的奇偶对称性，从而确定理想频率(n),用hM(n)代替(n)(4)选择适当的窗函数 wn，根据h(n)hd(n)WN(n)求所需设计的 FIR 滤波 器单位脉冲响应(5)求Hd(ejw)，分析其幅频特性，假设不满足要求，可适当改变窗函数形式或长度 N，重复上述设计过程，以得到满意的结果调用格式：a)b)b=fir1(N,w n)b=fir1(N,wn,hign N 为滤波器的阶次，wn 是通带截止频率，其值在 01 之间，1 对应抽样频率的一半；b 是设计好的滤波器系数h(n)。wn 是 1xL 的向量，那么对于格式 a),假设 wn 是标量，那么用来设计低通滤波器，假设用来设计带通(带阻)滤波器。其格式为：b=fir1(N,wn,DC-1 或 b=fir1(N,wn,DC-0 其中，前者保证第一个带为通带，后者保证 第一个带为阻带。b用来设计高通滤波器。用于产生窗函数的 matlab 文件有如下 8 个:bartlett.m(巴特利特窗)blackman.m(布莱克曼窗)rectangle.m(矩形窗)hamming.m汉明窗 hanning.m汉宁窗 triang.m(三角窗)chebwin.m切比雪夫窗 Kaiser.m 凯泽窗【实例 7-1，使用具有 30dB 波纹的 Chebyshev 窗。解：源程序如下：b=fir1(34,0.48,high,chebw in(35,30);freqz(b,1,512)其响应波形如图 7-1 所示。0|山-BOD：盲100D三-15002D0D001020 30 40 506070 90 91NEirnahredl FrECjLiEntv i：E sample)图 7 1 带通 FIR 滤波器【实例 7-2】设计一个 24 阶 FIR 带通滤波器，通带频率0.35 w 0.65解：源程序如下：wn=0.35 0.65;N=24;b=fir1(2*N,w n,DC-0);%DC-0表示第一个带为阻带，这样设计出来的滤波器为带通滤波器。freqz(b,1,512)EQH.P-PU?Z号JOSCI80 10.20.21313D.4D.4：.5.50.50.5Q.7Q.7 O.E 0.9 O.E 0.9NoiiiiAln*d riiequenc/|xjtrad-rtanipl*)Tom X TO 1 1 d%d%I I1 11 11iI I1 11 11 11 1i ii i 111 11 1i i ta ta1 11 1|m Tm T4 4 i ii i1 1I II I F F_L L _!_ _ _-_ _二-斗-r厂T Ti i1 1 i i i ii iK K1 1E E1 11 11 11 1a a|-1 1 r-r-:0)5I0 1 OJ2 0303D.405050.6Q.7Q.7 O.S O.S OiS1otJcim-iliz-ifLl riq jeiicy iQLrmd-lauipl m)图 7 2 带通 FIR 滤波器【实例 7-3】设计具有下面指标的低通FIR 滤波器Wp0.2,Rpws0.3,AS小的过渡带。解：MATLAB 源程序为%数字滤波器指标wp=0.2*pi;ws=0.3*pi;tr_width=ws-wp;%通带宽度N=ceil(6.6*pi/tr_width)+1;%滤波器长度if rem(N,2)=0N=N+1;endNw=N;%滤波器长度为奇数n=0:1:Nw-1;wc=(ws+wp)/2;%截止频率alpha=(Nw-1)/2;0.25dB,50dB由于其最小阻带衰减为50dB，因此可以选择 hamming 窗来实现这个滤波器，因为它具有较m=n-alpha+;%如果不加 0.00001，当 n=alpha 时，m 为 0，下式分母为 0，不能执行hd=sin(wc*m)./(pi*m);%理想滤波器的冲激响应%生成 hamming 窗w_ham=(hamming(Nw);%频域图像的绘制h=hd.*w_ham;%滤波器实际冲激响应freqz(h,1)figure(2);subplot(2,2,1),stem(n,hd);title(理想脉冲响应)axis(O Nw-1-0.3 0.3);xlabel(n);ylabel(hd(n)subplot(2,2,2),stem(n,w_ham);title(hamming 窗)axis(0 Nw-1-0.3 1.2);xlabel(n);ylabel(w(n)subplot(2,2,3),stem(n,h);title(实际脉冲响应)axis(0 Nw-1-0.3 0.3);xlabel(n);ylabel(h(n)其响应波形如图 7-2 所示。j一FD-50pOJ-100-15001 11 1i i1 1.H H4 4I Ia aI II I 1 1I I*ii iiV V1 1 1 11 1|1|11 1i i1 1II IIit it1 11 1 4 41 1G Gi iP P1 111ii ii1*1*：_-1 1-1 1-w w-*7-*7-:-|i i厂甘怕亍卩丫厂甘怕亍卩丫的丫协的丫协7 7丫汁丫汁辛辛0 0了了V V1 11 11 11 11 1A AH H1 11 1 I I1 1I-I-H Hh hH H B B 1 1I I I I 1 1A AH H Ii IiI Ifafa 1 11 1F Fh h1 1r r 1 10 1 D 20 3 C 40 0 50 0 510 0 70 80 9TIIorm orm I I ire d F requonire d F requon匚y(xp rodacmple|y(xp rodacmple|h h1 1 r iir iit tii ii r ri ii i Pi hi h-1000-1000-1500-1500-2000-2000.:1 11 1II II|l|lI I1 1h hi iii iii iibib1 1i ii i 1 1I IF Fi iVW0.20.20.30.30.40.40.50.50.60.60.70.70.60.60 90 9f Jcrrmal ipd F rpqiLJpnCA-一rad smplp:0.10.1理牌脉7中响时0 00-0-0-n n比际眛?中响丿竝图 7 3 hamming 窗函数设计 FIR 滤波器实例 7-4】设带通滤波器的指标为wp10.35,wp20.8,ws10.2,ws20.65,As60dB选择 Blackman 窗来实现这个滤波器。解：MATLAB 源程序为%数字滤波器指标ws1=0.2*pi;wp1=0.35*pi;ws2=0.65*pi;wp2=0.8*pi;As=60;tr_width=min(wp1-ws1),(wp2-ws2);%过渡带带宽N=ceil(11*pi/tr_width)+1;if rem(N,2)=0N=N+1;endNw=N;n=0:1:Nw-1;wc1=(ws1+wp1)/2;wc2=(ws2+wp2)/2;alpha=(Nw-1)/2;m=n-alpha+0.00001;%如果不加 0.00001，当 n=alpha 时，m 为 0，下式分母为 0，不能执行hd=sin(wc2*m)./(pi*m)-sin(wc1*m)./(pi*m);%理想滤波器的冲激响应%生成 blackman 窗 w_bla=(blackman(M);h=hd.*w_bla;%频域图像的绘制freqz(h,1)figure(2);subplot(2,2,1),stem(n,hd);title(idael impulse response)axis(0 Nw-1-0.4 0.5);xlabel(n);ylabel(hd(n)subplot(2,2,2),stem(n,w_bla);title(blackman window)axis(0 Nw-1 0 1.1);xlabel(n);ylabel(w(n)subplot(2,2,3),stem(n,h);title(actual impulse response)axis(0 Nw-1-0.4 0.5);xlabel(n);ylabel(h(n)其响应波形如图 7-4 所示。0-1O0-15000.10 20.30 40 0 60G G 0 7 O.g 0 G1Nomnalied Frequency(XTT rd/smple)20000-2000-4000-600%0L10.2F0.30.40.5C.60 70.60.91 Jcrmlied FreqLiencyx 一rad/saiTiple|idel impulse responsehiackman window图 7-4 blackman 窗函数设计 FIR 滤波器7.3编程练习1.设计一个 48 阶 FIR 带阻滤波器，通带频率2.用窗函数设计一个线性相位0 w 0.35,0.65 wFIR 高通滤波器，并满足性能指标：通带边界频率wp 0.66阻带边界频率ws 0.5，阻带衰减不小于 40db。提示：参考例 7 3，4;首先根据阻带衰减，选择适宜的窗函数，见课本根据 P345,式 7-74 写出理想滤波器的冲激响应hd。P342，表 7-3。再
