经典数字滤波器及其设计

1、经典数字滤波器及其设计4.1 引 言数字滤波器(Digital Filter, DF在信号的过滤、检测和参数估计等方面起着重要的作用。信号往往夹杂着噪声及无用信号成分，必须将这些干扰成分滤除。数字滤波器对信号进行筛选，可通过特定频段的信号。一般来说，噪声信号往往是高频信号，而经典滤波器正是假定有用信号与噪声信号具有不同的频段，所以利用经典滤波器可以去除噪声。但如果有用信号和无用信号，或有用信号和噪声的频谱相互重叠，那么经典滤波器则不能实现理想的滤波性能。现代滤波器的作用是从含有噪声的信号中估计出信号的某些特征或信号本身，一旦信号被估计出，那么估计出来的信号与原信号相比会有更高的信噪比。这类滤波

2、器主要有维纳滤波器、卡尔曼滤波器、线性预测滤波器以及自适应滤波器等。4.2 经典数字滤波器原理数字滤波是数字信号分析中最重要的组成部分之一，与模拟滤波相比，它具有精度和稳定性高、系统函数容易改变、灵活性强、便于大规模集成和可实现多维滤波等优点。在信号的过滤、检测和参数的估计等方面，经典数字滤波器是使用最广泛的一种线性系统。数字滤波器的作用是利用离散时间系统的特性对输入信号波形(或频谱 进行加工处理，或者说利用数字方法按预定的要求对信号进行变换。若滤波器的输入、输出都是离散时间信号，那么该滤波器的单位冲激响应h(n也必然是离散的，这种滤波器称为数字滤波器。当用硬件实现一个DF 时，所需的元件是乘

3、法器、延时器和相加器；而用MA TLAB 软件实现时，它仅仅需要线性卷积程序就可以实现。众所周知，模拟滤波器(Analog Filter ，AF 只能用硬件来实现，其元件有电阻R ，电感L ，电容C 及运算放大器等。因此，DF 的实现要比AF 容易得多，并且更容易获得较理想的滤波性能。 数字滤波器的作用是对输入信号进行滤波，就如同信号通过系统一样。对于线性时不变系统，其时域输入输出关系是：(4-1若y(n、x(n的傅里叶变化存在，则输入输出的频域关系是：(4-2当输入信号x(n通过滤波器h(n后，其输出y(n中不再含有的频率成分，仅使的信号成分通过，其中是滤波器的转折频率。经典数字滤波器按照单

4、位取样响应h(n的时域特性可分为无限冲激响应(IIR，Infinite Impulse Response 系统和有限冲激响应(FIR，Finite Impulse Response 系统。如果单位取样响应是时宽无限的h(n，则称之为IIR 系统；而如果单位取样响应是时宽有限的h(n，则称之为FIR 系统。数字滤波器按照实现的方法和结构形式分为递归型或非递归型两类。递归型数字滤波器的当前输出y(n是输入x(n的当前值和以前各输入值x(n，x(n1 ，. ，及以前各输出值y(n，y(n1 ，. 的函数。一个N 阶递归型数字滤波器(IIR滤波器 的差分方程为：(4-3其中，式(4-3中的系数至少有一

5、项不为零。说明必须将延时的输出序列进行反馈。 递归系统的传统函数定义为：(4-4递归系统的传递函数H(z在Z 平面上不仅有零点，而且有极点。非递归型数字滤波器当前的输出值y(n仅为当前及以前的输入序列的函数，而与以前的各个输出值无关，因此从结构上看非递归系统没有反馈环路。一个N 阶的非递归型数字滤波器(FIR滤波器 的差分方程为：(4-5差分方程式(4-5中的系数等于单位取样响应的序列值h(n，其系统函数H(z可以表示为以下形式：(4-6H(z是的多项式，因此它的极点只能在Z 平面的原点上。这两类滤波器无论是在性能上还是在设计方法上都有着很大的区别。FIR 滤波器可以对给定的频率特性直接进行设

6、计，而IIR 滤波器目前最通用的方法是利用已经很成熟的模拟滤波器的设计方法来进行设计。1. 特征频率滤波器的频率参数主要有：通带截频为通带与过渡带的边界点，在该点信号增益下降到规定的下限。阻带截频为阻带与过渡带的边界点，在该点信号衰耗下降到规定的下限。转折频率为信号功率衰减到1/2(约3dB 时的频率，在很多情况下，也常以fc 作为通带或阻带截频。当电路没有损耗时，固有频率，就是其谐振频率，复杂电路往往有多个固有频率。2. 增益与衰耗滤波器在通带内的增益并非常数。对低通滤波器通带增益一般指=0时的增益；高通指时的增益；带通则指中心频率处的增益。对带阻滤波器，应给出阻带衰耗，衰耗定义为增益的倒数

7、。通带增益变化量指通带内各点增益的最大变化量，如果以dB 为单位，则指增益dB 值的变化量。3. 阻尼系数与品质因数阻尼系数是表征滤波器对角频率为0信号的阻尼作用，是滤波器中表示能量衰耗的一项指标，它是与传递函数的极点实部大小相关的一项系数。它可由传递函数的分母多项式系数求得：(4-7式中，表示传递函数的分母多项式系数。的倒数Q=1/称为品质因数，是评价带通与带阻滤波器频率选择特性的一个重要指标，Q 为：(4-8式中的为带通或带阻滤波器的3dB 带宽，为中心频率，在很多情况下中心频率与固有频率相等。4. 灵敏度滤波电路由许多元件构成，每个元件参数值的变化都会影响滤波器的性能。滤波器某一性能指标

8、对某一元件参数变化的灵敏度记作，定义为：(4-9灵敏度是滤波电路设计中的一个重要参数，可以用来分析元件实际值偏离设计值时，电路实际性能与设计性能的偏差程度；也可以用来估计在使用过程中元件参数值变化时，电路性能变化情况。该灵敏度与测量仪器或电路系统灵敏度概念不同，该灵敏度越小，标志着电路容错能力越强，稳定性也越高。5. 群时延函数当滤波器幅频特性满足设计要求时，为保证输出信号失真度不超过允许范围，对其相频特性也应提出一定要求。在滤波器设计中，常用群时延函数评价信号经滤波后相位失真程度。越接近常数，信号相位失真越小。6. 滤波器设计的步骤不论是IIR 滤波器还是FIR 滤波器的设计都包括三个步骤：

9、(1 按照实际任务的要求，确定滤波器的性能指标。(2 用一个因果、稳定的离散线性时不变系统的系统函数去逼近这一性能指标。根据不同的要求可以用IIR 系统函数，也可以用FIR 系统函数去逼近。(3 利用有限精度算法实现系统函数，包括结构选择、字长选择等。IIR 数字滤波器设计的最通用的方法是借助于模拟滤波器的设计方法。模拟滤波器设计已经有一套相当成熟的方法，它不但有完整的设计公式，而且还有较为完整的图表以供查询，因此充分利用这些已有的资源将会给数字滤波器的设计带来很大的方便。IIR 数字滤波器的设计步骤是：(1 按一定的规则将给出的数字滤波器的技术指标转换为模拟低通滤波器的技术指标。(2 根据转

10、换后的技术指标设计模拟低通滤波器G(s。(3 再按一定的规则将G(s转换成H(s。若所设计的数字滤波器是低通的，那么上述设计工作可以结束，若所设计是高通、带通或带阻滤波器，那么还需进行以下步骤。(4 将高通、带通或带阻数字滤波器的技术指标先转换为低通模拟滤波器的技术指标，然后按照上述步骤设计出低通G(s，再将G(s转换为所需的H(z。4.3 IIR 数字滤波器的结构与设计IIR 滤波器系统函数的极点可以在单位圆内的任何位置，实现IIR 滤波器的阶次较低，所用的存储单元较少，效率高，又由于IIR 数字滤波器能够保留一些模拟滤波器的优良特性，因此应用很广。设计IIR 数字滤波器的方法主要有基于冲激

11、响应不变法的IIR 数字滤波器设计，基于双线性Z 变换法的IIR 数字滤波器设计，数字高通、带通及带阻IIR 滤波器设计，基于MATLAB 函数直接设计IIR 数字滤波器。冲激响应不变法的设计原理是利用数字滤波器的单位抽样响应序列H(z来逼近模拟滤波器的冲激响应g(t。按照冲激响应不变法的原理，通过模拟滤波器的系统传递函数G(s，可以直接求得数字滤波器的系统函数H(z，其转换步骤如下：1 利用=T(可由关系式推导出 ，将，转换成，而，不变；2 求解低通模拟滤波器的传递函数G(s；3 将模拟滤波器的传递函数G(s转换为数字滤波器的传递函数H(z。尽管通过冲激响应不变法求取数字滤波器的系统传递函数

12、比较方便，并具有良好的时域逼近特性，但若G(s不是带限的，或是抽样频率不高，那么在中将发生混叠失真，数字滤波器的频率响应不能重现模拟滤波器的频率响应。只有当模拟滤波器的频率响应在超过折叠频率后的衰减很大时，混叠失真才很小，此时采样脉冲响应不变法设计的数字滤波器才能满足设计的要求，这是冲激响应不变法的一个严重的缺点。【例4-1】设模拟低通巴特沃斯滤波器，通带纹波为1dB ，通带上限角频率=0.2，阻带下限角频率=0.3，阻带最小衰减=15dB，根据该低通模拟滤波器，利用冲激响应不变法设计响应的数字低通滤波器，并绘出设计后的数字滤波器的特性曲线。实现例4-1的MA TLAB 程序如例程4-1所示，

13、程序的运行结果如图4-1所示。例程4-1 利用冲激响应不变法设计数字低通滤波器%利用模拟巴特沃斯滤波器设计数字低通滤波器%冲激响应不变法wp=0.2*pi;ws=0.3*pi;Rp=1;As=15;T=1;%性能指标Rip=10(-Rp/20;Atn=10(-As/20;OmgP=wp*T;OmgS=ws*T;N,OmgC=buttord(OmgP,OmgS,Rp,As,'s' %选取模拟滤波器的阶数cs,ds=butter(N,OmgC,'s' %设计出所需的模拟低通滤波器b,a=impinvar(cs,ds,T; %应用脉冲响应不变法进行转换%求得相对、绝对

14、频响及相位、群迟延响应db,mag,pha,grd,w=freqz_m(b,a;%下面绘出各条曲线subplot(2,2,1;plot(w/pi,mag;title('幅频特性'xlabel('w(/pi'ylabel('|H(jw|'axis(0,1,0,1.1;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;set(gca,'YTickMode','manual','YTick',0 Atn

15、 Rip 1;gridsubplot(2,2,2;plot(w/pi,db;title('幅频特性(db'xlabel('w(/pi'ylabel('dB'axis(0,1,-40,5;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;set(gca,'YTickMode','manual','YTick',-40 -As -Rp 0;gridsubplot(2,2,3;plot(w/pi,pha

16、/pi;title('相频特性'xlabel('w(/pi'ylabel('pha(/pi'axis(0,1,-1,1;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;gridsubplot(2,2,4;plot(w/pi,grd;title('群延时'xlabel('w(/pi'ylabel('Sample'axis(0,1,0,12;set(gca,'XTickMode',&

17、#39;manual','XTick',0 0.2 0.3 0.5 1;grid%所使用的M 文件函数functiondb,mag,pha,grd,w=freqz_m(b,a%滤波器的幅值响应(相对、绝对 、相位响应及群延迟%Usage:db,mag,pha,grd,w=freqz_m(b,a %500点对应0,pi %w 采样频率； b 系统函数H(z的分子项(对FIR ，b=h %a 系统函数H(z的分母项(对FIR ，a=1H,w=freqz(b,a,500; %500点的复频响应 mag=abs(H; %绝对幅值响应 db=20*log10(mag/max(ma

18、g; %相对幅值响应 pha=angle(H; %相位响应grd=grpdelay(b,a,w; %群延迟响应图4-1 利用冲激响应不变法设计的数字低通滤波器的特性曲线由于的频率映射关系是根据推导的，所以使j轴每隔2/Ts便映射到单位圆上一周，利用冲激响应不变法设计数字滤波器时可能会导致上述的频域混叠现象。为了克服这一问题，需要找到由s 平面到z 平面的另外的映射关系，这种关系应保证：1 s平面的整个j轴仅映射为z 平面单位圆上的一周；2 若G(s是稳定的，由G(s映射得到的H(z也应该是稳定的；3 这种映射是可逆的，既能由G(s得到H(z，也能由H(z得到G(s；4 如果G(j0=1，那么。

19、双线性Z 变换满足以上4个条件的映射关系，其变换公式为(4-10双线性Z 变换的基本思路是：首先将整个s 平面压缩到一条从/Ts到/Ts的带宽为2/Ts的横带里，然后通过标准的变换关系将横带变换成整个z 平面上去，这样就得到s 平面与z平面间的一一对应的单值关系。在MA TLAB 中，双线性Z 变换可以通过bilinear 函数实现，其调用格式为：Bz，Azbilinear(B，A ，Fs ；其中B ，A 为模拟滤波器传递函数G(s的分子分母多项式的系数向量，而Bz ，Az 为数字滤波器的传递函数H(z的分子分母多项式的系数向量。【例4-2】使用双线性Z 变换法设计一个低通数字滤波器，给定的数

20、字滤波器的技术指标为fp 100Hz ，fs=300Hz，p=3dB，p=20dB，抽样频率为Fs=1000Hz。实现例4-2的MATLAB程序如例程4-2所示，程序的运行结果如图4-2所示。例程4-2 基于双线性Z 变换法设计数字滤波器clear all;fp=100;fs=300;Fs=1000;rp=3;rs=20;%wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Fs=Fs/Fs;%Firstly to finish frequency prewarping;wap=tan(wp/2;was=tan(ws/2;n,wn=buttord(wap,was,rp,rs,'s&

21、#39;%Note:'s'!z,p,k=buttap(n;bp,ap=zp2tf(z,p,kbs,as=lp2lp(bp,ap,wap%Note:s=(2/Ts(z-1(z+1;Ts=1,that is 2fs=1,fs=0.5;bz,az=bilinear(bs,as,Fs/2h,w=freqz(bz,az,256,Fs*1000;plot(w,abs(h;grid on;【运行结果】n =2wn =0.4363bp =0 0 1ap =1.0000 1.4142 1.0000bs =0.1056as =1.0000 0.4595 0.1056bz =0.0675 0.134

22、9 0.0675az =1.0000 -1.1430 0.4128图4-2 基于双线性Z 变换法设计的数字滤波器的特性曲线除了低通数字滤波器之外，实际中还常常需要高通、带通及带阻数字滤波器，这三种数字滤波器的设计步骤如下：步骤1：将数字滤波器H(z的技术指标p和s，通过tan(/2转变为模拟滤波器G(s的技术指标p和s，作归一化处理后，得到p1，ss/p；步骤2：化解为模拟原型滤波器G(s的技术指标；步骤3：设计模拟原型滤波器G(p；步骤4：将G(p转换为模拟滤波器的转移函数G(s；步骤5：将G(s转换成数字滤波器的转移函数H(z，s=(z1(z+1。所谓原型滤波器是指归一化的低通滤波器。本节

23、主要讨论通过IIR 数字滤波器的原型转换设计法和IIR 数字滤波器的直接设计方法来设计数字高通、带通及带阻滤波器其转换方法主要有3种：一是直接由模拟低通滤波器转换成数字高通、带通或带阻滤波器；二是先由模拟低通滤波器转换成模拟高通、带通或带阻滤波器，然后再把它转换成相应的数字滤波器；三是将模拟低通滤波器先转换成数字低通滤波器，再通过变量代换变换成高通、带通或带阻滤波器。1. 直接由模拟低通滤波器转换成数字高通、带通或带阻滤波器(1 模拟低通滤波器转换成数字高通滤波器【例4-3】 利用巴特沃斯模拟滤波器设计数字高通滤波器，要求通带截止频率为0.6，通带内衰减不大于1dB ，阻带起始频率为0.4，阻

24、带内衰减不小于15dB ，采样周期Ts 1。该例的MA TLAB 实现如例程4-3所示，该数字滤波器的幅频特性曲线如图4-3所示。例程4-3 基于巴特沃斯模拟滤波器设计数字高通滤波器Wp=0.6*pi;Ws=0.4*pi;Ap=1;As=15;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %计算巴特沃斯滤波器阶次和截止频率b,a=butter(N,wn,'high' %频率变换法设计巴特沃斯高通滤波器db,mag,pha,grd,w=freqz_m(b,a; %数字滤波器响应Plot(w/pi,mag;Title('数字滤波器幅频响应|H(ejOmeg

25、a|'图4-3 基于巴特沃斯模拟滤波器设计的数字高通滤波器幅频响应曲线(2 模拟低通滤波器转换成数字带通滤波器【例4-4】利用巴特沃斯模拟滤波器设计数字带通滤波器，要求通带上下截止频率为0.4、0.3，通带内衰减不大于3dB ，阻带上下起始频率为0.5、0.2，阻带内衰减不小于18dB 。该例的MA TLAB 实现如例程4-4所示，该数字滤波器的幅频特性曲线如图4-4 所示。例程4-4 基于巴特沃斯模拟滤波器设计数字带通滤波器Wp=0.3*pi,0.4*pi;Ws=0.2*pi,0.5*pi;Ap=3;As=18;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %计算巴

26、特沃斯滤波器阶次和截止频率b,a=butter(N,wn,'bandpass' %频率变换法设计巴特沃斯带通滤波器db,mag,pha,grd,w=freqz_m(b,a; %数字滤波器响应Plot(w/pi,mag;Title('数字滤波器幅频响应|H(ejOmega|'图4-4 基于巴特沃斯模拟滤波器设计的数字带通滤波器的幅频响应曲线(3 模拟低通滤波器转换成数字带阻滤波器【例4-5】利用巴特沃斯模拟滤波器设计数字带阻滤波器，要求通带上下截止频率为0.8和该例的MA TLAB 实现如例程4-5所示，该数字滤波器的幅频特性曲线如图4-5所示。例程4-5 基于巴

27、特沃斯模拟滤波器设计数字带阻滤波器Wp=0.4*pi,0.7*pi;Ws=0.2*pi,0.8*pi;Ap=1;As=30;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %计算巴特沃斯滤波器阶次和截止频率b,a=butter(N,wn,'stop' %频率变换法设计巴特沃斯带通滤波器db,mag,pha,grd,w=freqz_m(b,a; %数字滤波器响应plot(w/pi,mag;Title('数字滤波器幅频响|H(ejOmega| '图4-5 基于巴特沃斯模拟滤波器设计数字带阻滤波器2. 先由模拟低通滤波器转换成模拟高通、带通或带阻滤波器

28、，然后再把它转换成相应的数字滤波器首先由频率变换将模拟低通原型滤波器变换成模拟高通、带通或带阻滤波器，然后根据双线性Z 变换将模拟高通滤波器变换成相应的高通、带通或带阻数字滤波器。其模拟低通原型滤波器变换成模拟高通、带通或带阻滤波器的关系式如表4-1所示。表4-1 模拟低通原型滤波器变换成模拟高通、带通或带阻滤波器频率变换模拟滤波器类型频率变换公式高通带通带阻3. 将模拟低通滤波器先转换成数字低通滤波器，再通过变量代换变换成高通、带通或带阻滤波器该方法中的模拟低通滤波器先转换成数字低通滤波器，可以利用冲激响应不变法和双线性变换法来实现，再通过变量代换变换成数字高通、带通或带阻滤波器。【例4-6

29、】试用双线性变换法设计一个带通椭圆滤波器使其幅频特性逼近于一个具有如下技术指标的模拟带通滤波器：Wp1=10Hz,Wp2=20Hz,在通带内的最大衰减为0.5dB, 在阻带内的最大衰减为50dB ，抽样频率为100Hz 。MATLAB 实现如例程4-6所示，该数字滤波器的幅频特性曲线如图4-6所示。wp1=10;wp2=20;Fs=100;rp=0.5;rs=50;wp1=2*pi*wp1;wp2=2*pi*wp2;Bw=wp2-wp1;Wo=sqrt(wp2*wp1;z,p,k=ellipap(7,rp,rs;A,B,C,D=zp2ss(z,p,k;At,Bt,Ct,Dt=lp2bp(A,B

30、,C,D,Wo,Bw;At1,Bt1,Ct1,Dt1=bilinear(At,Bt,Ct,Dt,Fs;num,den=ss2tf(At1,Bt1,Ct1,Dt1;H,W=freqz(num,den;plot(W*Fs/(2*pi,abs(H;grid;xlabel('频率/Hz'ylabel('幅值'图4-6 椭圆滤波器的频率响应4.4 FIR 数字滤波器的结构与设计由于IIR 数字滤波器能够保留一些模拟滤波器的优良特性，因此应用很广。但是这些特性是以牺牲线性相位频率特性为代价的，即用Butterworth 、chelbchev 和椭圆法设计的数字滤波器逼近理想

31、的滤波器的幅度频率特性，得到的滤波器往往是非线性的。在许多电子系统中，对幅度频率特性和线性相位特性都有较高的要求，所以IIR 滤波器在这些系统中往往难以胜任。有限长单位冲激响应(FIR数字滤波器具有以下优良的特点： 可在设计任意幅度频率特性滤波器的同时，保证精确、严格的线性相位特性。 FIR数字滤波器的单位冲激响应h(n是有限长的，可以用一个固定的系统来实现，因而FIR数字滤波器可以做成因果稳定系统。 允许设计多通带(多阻带 系统。根据FIR 滤波器的幅频特性，可以将其分为以下4种情况，其中h(n为FIR 滤波器的冲激响应：1. h(n为偶对称，且N 为奇数(型滤波器型滤波器的幅频函数可以通过

32、例程4-7的MA TLAB 语言来实现。例程4-7 型滤波器的幅频函数function Hr,w,a,L=hr_type1(h;%计算所设计的1型滤波器的振幅响应%Hr=振幅响应%a=1型滤波器的系数%L=Hr的阶次%h=1型滤波器的单位冲击响应M=length(h;L=(M-1/2;a=h(L+1 2*h(L:-1:1;n=0:1:L;w=0:1:500'*2*pi/500;Hr=cos(w*n*a'2. h(n为偶对称，且N 为偶数(型滤波器型滤波器的幅频函数可以通过例程4-8 的MA TLAB 语言来实现。例程4-8 型滤波器的幅频函数function Hr,w,b,L=

33、hr_type2(h;%计算所设计的2型滤波器的振幅响应%Hr=振幅响应%b=2型滤波器的系数%L=Hr的阶次%h=2型滤波器的单位冲击响应M=length(h;L=M/2;b= 2*h(L:-1:1;n=1:1:L;n=n-0.5;w=0:1:500'*2*pi/500;Hr=cos(w*n*b'3. h(n为奇对称，且N 为奇数(型滤波器型滤波器的幅频函数可以通过例程4-9的MA TLAB 语言来实现。例程4-9 型滤波器的幅频函数function Hr,w,c,L=hr_type3(h;%计算所设计的3型滤波器的振幅响应%Hr=振幅响应%b=3型滤波器的系数%L=Hr的阶

34、次%h=3型滤波器的单位冲击响应M=length(h;L=(M-1/2;c= 2*h(L+1:-1:1;n=0:1:L;w=0:1:500'*2*pi/500;Hr=sin(w*n*c'4. h(n为奇对称，且N 为偶数(型滤波器型滤波器的幅频函数可以通过例程4-10 的MATLAB 语言来实现。例程4-10 型滤波器的幅频函数function Hr,w,d,L=hr_type4(h;%计算所设计的4型滤波器的振幅响应%Hr=振幅响应%b=4型滤波器的系数%L=Hr的阶次%h=4型滤波器的单位冲击响应M=length(h;L=M/2;d= 2*h(L:-1:1;n=1:1:L;

35、n=n-0.5;w=0:1:500'*2*pi/500;Hr=sin(w*n*d'针对FIR 数字滤波器的结构特点，目前主要采用窗函数法、频率采样法和最优化法等三种方法设计FIR 数字滤波器。窗函数法就是设计FIR 数字滤波器的最简单的方法。它在设计FIR 数字滤波器中有很重要的作用，正确地选择窗函数可以提高设计数字滤波器的性能，或者在满足设计要求的情况下，减小FIR 数字滤波器的阶次。常用的窗函数有以下几种：矩形窗(Rectangular window、三角窗(Triangular window 、汉宁窗(Hanning window 、海明窗(Hamming window

36、、布拉克曼窗(Blackman window、切比雪夫窗(Chebyshev window、巴特里特窗(Bartlett window及凯塞窗(Kaiser window。在MA TLAB 中，实现矩形窗的函数为boxcar 和rectwin ，其调用格式如下：w=boxcar(Nw=rectwin(N其中N 是窗函数的长度，返回值w 是一个N 阶的向量，它的元素由窗函数的值组成。实际上，w=boxcar(N等价于w=ones(N,1。在MA TLAB 中，实现三角窗的函数为triang ，调用格式为：w=triang(N在MA TLAB 中，实现汉宁窗的函数为hann ，调用格式如下：w=h

37、ann(Nw=hann(N,'sflag'Hann 函数中的参数sflag 为采样方式，其值可取symmetric(默认值 或periodic 。当sflag symmetric 时，为对称采样；当sflag periodic 时，为周期采样，此时hann 函数计算N+1个点的窗，但是仅返回前N 个点。在MA TLAB 中，实现海明窗的函数为hamming ，调用格式分别如下：w=hamming (Nw=hamming (N,'sflag'其中sflag 的用法同上。在MA TLAB 中，实现布拉克曼窗的函数为blackman ，调用格式如下：w=blackma

38、n (Nw=blackman (N,'sflag'在MA TLAB 中，实现切比雪夫窗的函数为chebwin ，调用格式为：w=chebwin (N,r其中r 表示切比雪夫窗函数的傅里叶变换旁瓣幅度比主瓣低rdB(其默认值为100dB ，且旁瓣是等纹波的。在MA TLAB 中，实现巴特里特窗的函数为bartlett ，调用格式为：w=bartlett (N在MA TLAB 中，实现凯塞窗的函数为kaiser ，调用格式为：w=kaiser (N,beta其中beta 为窗函数的参数。各种窗函数的性能比较可见表4-2，在设计FIR 滤波器的过程中可以根据要求选择合适的窗函数：表4

39、-2 各种窗函数的性能比较窗 函 数第一旁瓣相对于主瓣衰减/dB主 瓣 宽阻带最小衰减/dB矩形窗134/N21三角窗258/N25汉宁窗318/N44海明窗418/N53布拉克曼窗5712/N74凯塞窗可调可调可调切比雪夫窗可调可调可调1. 数字低通滤波器的窗函数设计计算、型理想低通滤波器的单位冲激响应hd(n的MATLAB 实现如例程4-11所示。 例程4-11 、型理想低通滤波器的单位冲激响应计算function hd=ideal_lp1(Wc,N%compute the ideal lowpass fiter unit pulse respondence hd(n%wc: cutoff

40、 frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(Wc*m./(pi*m;计算、型理想低通滤波器的单位冲激响应hd(n的MATLAB 实现如例程4-12所示。 例程4-12 、型理想低通滤波器的单位冲激响应计算function hd=ideal_lp2(Wc,N%compute the ideal lowpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window leng

41、th%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(Wc*m./(pi*m;【例4-7】根据下列技术指标，设计一个FIR 数字低通滤波器：=0.2，=0.4，=0.25bB，=50bB选择一个适当的窗函数，确定单位冲激响应，绘出所设计的滤波器的幅度响应。根据窗函数最小阻带衰减的特性表4-2，可采用海明窗和布拉克曼窗可提供大于50dB 的衰减。本例选择海明窗，其过渡带为6.6/N，因此具有较小的阶次。该例的MA TLAB 实现如例程4-13 所示，其设计滤波器的冲激响应与幅度响应如图4-7所示。 例程4

42、-13 利用海明窗设计型数字低通滤波器clear all;Wp=0.2*pi;Ws=0.4*pi;tr_width=Ws-Wp; %过渡带宽度N=ceil(6.6*pi/tr_width+1 %滤波器长度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通滤波器的截止频率hd=ideal_lp1(Wc,N; %理想低通滤波器的单位冲激响应w_ham=(hamming(N' %海明窗h=hd.*w_ham; %截取得到实际的单位脉冲响应db,mag,pha,w=freqz_m2(h,1; %计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(1:1:W

43、p/delta_w+1 %实际通带纹波As=-round(max(db(Ws/delta_w+1:1:501 %实际阻带纹波subplot(221stem(n,hdtitle('理想单位脉冲响应hd(n'subplot(222stem(n,w_hamtitle('海明窗w(n'subplot(223stem(n,htitle('实际单位脉冲响应hd(n'subplot(224plot(w/pi,dbtitle('幅度响应(dB'axis(0,1,-100,10%-functiondb,mag,pha,w=freqz_m2(b,a%

44、滤波器的幅值响应(相对、绝对 、相位响应%db:相对幅值响应%mag：绝对幅值响应%pha: 相位响应%w 采样频率；%b 系统函数H(z的分子项(对FIR ，b=h%a 系统函数H(z的分母项(对FIR ，a=1H,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; %相对幅值响应pha=angle(H; %相位响应【运行结果】N =34Ap =0.0477As =52由此结果可知，所设计的低通滤波器为型滤波器，它的通带

45、纹波和阻带纹波均满足设计要求。图4-7 型FIR 数字低通滤波器冲激响应与幅度响应2. 数字高通滤波器的窗函数设计数字高通滤波器的设计分两种情况，即无相移和相移为±/2两种情况。当无相移时，N 为奇数时，所设计的FIR 数字高通滤波器为型滤波器；当N 为偶数时，为型滤波器。当相移为±/2时，N 为奇数时，所设计的FIR 数字高通滤波器为型滤波器；当N 为偶数时，为型滤波器、型理想高通滤波器的单位冲激响应hd(n计算的MATLAB 实现如例程4-14所示。 例程4-14 、型理想高通滤波器的单位冲激响应计算function hd=ideal_hp1(Wc,N%compute

46、the ideal highpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(pi*m-sin(Wc*m./(pi*m;、型理想高通滤波器的单位冲激响应hd(n计算的MATLAB 实现如例程4-15所示。 例程4-15 、型理想高通滤波器的单位冲激响应hd(n计算function hd=ideal_hp2(Wc,N%compute the ideal h

47、ighpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(pi*m-cos(Wc*m./(pi*m;【例4-8】根据下列技术指标，设计一个FIR 数字高通滤波器：=0.6，=0.4，=0.25bB，=50bB选择一个适当的窗函数，确定单位冲激响应，绘出所设计滤波器的幅度响应。根据窗函数最小阻带衰减的特性表4-2所示，选择汉宁窗可达到 44dB 的最小阻带衰

48、减，其过渡带为6.2/N，因此具有较小的阶次。该例的MA TLAB 实现如例程4-16所示，其设计的滤波器的冲激响应与幅度响应如图4-8所示。 例程4-16 利用汉宁窗设计型数字高通滤波器clear all;Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws; %过渡带宽度N=ceil(6.2*pi/tr_width %滤波器长度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通滤波器的截止频率hd=ideal_hp1(Wc,N; %理想低通滤波器的单位冲激响应w_han=(hanning(N' %汉宁窗h=hd.*w_han; %截取得到实际的单位脉冲响应db

49、,mag,pha,w=freqz_m2(h,1; %计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(Wp/delta_w+1:1:501 %实际通带纹波As=-round(max(db(1:1:Ws/delta_w+1 %实际阻带纹波subplot(221stem(n,hdtitle('理想单位脉冲响应hd(n'subplot(222stem(n,w_hantitle('汉宁窗w(n'subplot(223stem(n,htitle('实际单位脉冲响应hd(n'subplot(224plot(w/pi,dbti

50、tle('幅度响应(dB'axis(0,1,-100,10【运行结果】N = 31Ap = 0.0887As = 44由此结果可知，所设计的高通滤波器为型滤波器，它的通带纹波和阻带纹波均满足设计要求。图4-8 型FIR 数字高通滤波器冲激响应与幅度响应【例4-9】根据下列技术指标，设计一个具有+/2相移的FIR 数字高通滤波器：=0.6，=0.4，=0.25bB，=50bB选择一个适当的窗函数，确定单位冲激响应，绘出所设计的滤波器的幅度响应。根据窗函数最小阻带衰减的特性表4-2，选择汉宁窗可达到 44dB 的最小阻带衰减，其过渡带为6.2/N，因此具有较小的阶次。该例的MA T

51、LAB 实现如例程4-17所示，其设计的滤波器的冲激响应与幅度响应如图4-9所示。 例程4-17 利用汉宁窗设计型数字高通滤波器clear all;Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws; %过渡带宽度N=ceil(6.2*pi/tr_width+1 %滤波器长度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通滤波器的截止频率hd=ideal_hp2(Wc,N; %理想低通滤波器的单位冲激响应w_han=(hanning(N' %汉宁窗h=hd.*w_han; %截取得到实际的单位脉冲响应db,mag,pha,w=freqz_m2(h,1; %计算

52、实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(Wp/delta_w+1:1:501 %实际通带纹波As=-round(max(db(1:1:Ws/delta_w+1 %实际阻带纹波subplot(221stem(n,hdtitle('理想单位脉冲响应hd(n'subplot(222stem(n,w_hantitle('汉宁窗w(n'subplot(223stem(n,htitle('实际单位脉冲响应hd(n'subplot(224plot(w/pi,dbtitle('幅度响应(dB'axis(0,

53、1,-100,10【运行结果】N = 32Ap = 0.0711As = 44由此结果可知，所设计的高通滤波器为型滤波器，它的通带纹波和阻带纹波均满足设计要求。图4-9 型FIR 数字高通滤波器冲激响应与幅度响应3. 数字带通滤波器的窗函数设计数字带通滤波器的设计也分两种情况讨论，即无相移和相移为±/2两种情况。当无相移时，N 为奇数时，所设计的FIR 数字带通滤波器为型滤波器；当N 为偶数时，为型滤波器。当相移为±/2时，N 为奇数时，所设计的FIR 数字带通滤波器为型滤波器；当N 为偶数时，为型滤波器、型理想带通滤波器的单位冲激响应hd(n计算的MATLAB 实现如例程

54、4-18所示。 例程4-18 、型理想带通滤波器的单位冲激响应计算function hd=ideal_bp1(Wcl,Wch,N%compute the ideal bandpass fiter unit pulse respondence hd(n%wcl: low cutoff frequency%wch: high cutoff frequency%N: window length%hd: unit pulse respondencealpha= (N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(Wch*m-sin(Wcl*m./(pi*m;、型理想带通滤波器的

55、单位冲激响应hd(n计算的MATLAB 实现如例程4-19所示。 例程4-19 、型理想带通滤波器的单位冲激响应hd(n计算function hd=ideal_bp2(Wcl,Wch,N%compute the ideal bandpass fiter unit pulse respondence hd(n%wcl: low cutoff frequency%wch: high cutoff frequency%N: window length%hd: unit pulse respondencealpha= (N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(Wch*

56、m-cos(Wcl*m./(pi*m;【例4-10】根据下列技术指标，设计一个FIR 数字带通滤波器：=0.4，=0.6，=0.2，=0.8.=1bB，=1bB，=60bB ，=60bB.选择一个适当的窗函数，确定单位冲激响应，绘出所设计的滤波器的幅度响应。根据窗函数最小阻带衰减的特性表4-2，选择布拉克曼窗可达到 75dB 的最小阻带衰减，其过渡带为11/N。该例的MA TLAB 实现如例程4-20所示，其设计的滤波器的冲激响应与幅度响应如图4-10所示。例程4-20 利用布拉克曼窗设计型数字带通滤波器clear all;Wpl=0.4*pi;Wph=0.6*pi;Wsl=0.2*pi;Wsh=0.8*pi;tr_width=min(Wpl-Wsl,(Wsh-Wph; %过渡带宽度N=ceil(11*pi/tr_width+1 %滤波器长度n=0:1:N-1;