数字信号处理教程-MATLAB释义与实现第07章.ppt

《数字信号处理教程-MATLAB释义与实现第07章.ppt》由会员分享，可在线阅读，更多相关《数字信号处理教程-MATLAB释义与实现第07章.ppt（144页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第七章第七章 FIRFIR滤波器设计滤波器设计 1FIRFIR滤波器设计滤波器设计本章的主要内容为：7.1 数字滤波器设计的基本要求 7.2 线性相位滤波器的特性 7.3 设计FIR滤波器的窗函数法 7.4 设计FIR滤波器频率样本法7.5 FIR滤波器的最优设计法 7.6 FIR滤波器设计的一些深入问题 2滤波器设计概述滤波器设计概述从本章起将着重研究信号分析的逆问题，那就是给定了输入输出指标，设计系统。这种可以人为调整参数来满足给定要求的硬件或软件系统，称为滤波器。这种问题没有唯一的解，因为可能有多种类型的滤波器都能满足给定的要求。要人为地进行选择，就要考虑更多的因素,所以设计比之分析更加

2、接近与工程实际。有两种类型的数字滤波系统。本章将讨论FIR数字滤波器的几种基本设计算法；下一章讨论IIR数字滤波器的设计问题。3滤波器设计概述滤波器设计概述本书中的设计大多数是针对选频型的滤波器。从幅频特性的要求而言，有低通、高通、带通和带阻滤波器；从相频特性来说，往往要求线性相位。在FIR滤波器设计中，还会涉及到像微分器或希尔伯特变换器之类的系统，尽管它们不是选频滤波器，但是也遵循选频滤波器的设计技术。在计算技术高度发展的今天，几乎所有的电子设计软件中都把滤波器设计作为其模块之一。那种用查表、手工计算的方法在工程实践中已经濒于淘汰。4滤波器设计概述滤波器设计概述 数字滤波器的设计要经过三个步

3、骤：确定指标确定指标：在设计一个滤波器前，必须有一些指标。这些指标要根据应用确定。模型逼近模型逼近：一旦确定了技术指标，就可利用前面学习过的基本原理和关系式，提出一个滤波器模型来逼近给定的指标体系。这是滤波器设计所要研究的主要问题。实现实现：上面的两步的结果得到的滤波器，通常是以差分方程、系统函数或脉冲响应来描述的。根据这个描述，可以如第六章所讨论的那样，用硬件或计算机软件来实现它。5滤波器设计概述滤波器设计概述数字滤波器的主要指标：数字滤波器的主要指标：在很多实际应用，象语音或音频信号处理中，数字滤波器用来实现选频功能。因此，指标的形式通常为频域中的幅度和相位响应。幅度指标幅度指标可以两种方

4、式给出。第一种，叫做绝对指标，它提出了对幅度响应函数|H(j)|的要求。这些指标一般可直接用于FIR滤波器。第二种方法叫做相对指标，它以分贝(dB)值的形式提出对幅频特性的要求。其值定义为：6FIRFIR滤波器设计概述滤波器设计概述相对幅度指标：由于定义中所包含的归一化，滤波器的相对幅频特性最高处的值为0dB，又由于定义中的负号，幅频特性小的地方，其dB值反而是正的。为了说明这些指标，以低通滤波器为例进行讨论。首先看下面两种指标图，它们分别用绝对指标和相对指标表示。7FIRFIR滤波器设计概述滤波器设计概述8FIRFIR滤波器设计概述滤波器设计概述绝对指标：绝对指标：图a是典型的低通幅特性绝对

5、指标：0,p段叫做通带，1是在理想通带中能接受的振幅波动（或容限）,011。s,段叫做阻带，2是阻带中能接受的振幅波动（或容限），021。p,s叫做过渡带，在此段上对幅度响应通常没有限制，也可以给些弱限制。相对指标（相对指标（dB）：）：低通的典型相对技术指标如图b所示，三种频带的定义不变，只是1，2 换成了Rp，As。其中 Rp是通带波动的dB数;As是阻带衰减的dB数；9FIRFIR滤波器设计概述滤波器设计概述绝对指标和相对指标间的变换关系：可以用MATLAB语句fgp711准确地求出绝对指标和相对指标的对应数值。并画出如图中所示曲线。10FIRFIR滤波器设计概述滤波器设计概述11FIR

6、FIR滤波器设计概述滤波器设计概述相位指标相位指标主要是线性相位条件。理想的滤波器相位应该和频率成正比，也就是()=-。因此满足此式的系统对信号中所有的频率分量都具有相同的时间延迟(无量纲)，因此无波形失真。略低一点的要求是相位和频率成线性关系，即()=0-。IIR滤波器无法做到线性相位，通常用群迟延g()来评价。其定义为相位对频率的导数 ，它愈接近常数愈好。12FIRFIR滤波器设计概述滤波器设计概述上述指标是针对低通滤波器的，其他类型的选频滤波器，像高通和带通滤波器，它们的技术指标可作类似的规定，并且可以由低通滤波器演化而得。本章讨论FIR数字滤波器的设计和逼近。FIR滤波器在设计和实现上

7、具有如下的优越性：(a).相位响应可为严格的线性；(b).不存在稳定性问题；(c).只包含实数算法，不需要递推运算。137.2 7.2 线性相位滤波器的特性线性相位滤波器的特性 线性相位条件要求滤波器分子系数bn满足对称性条件。设滤波器的系数长度为N，则这些系数应关于中心点=n0=(N-1)/2对称。偶对称时，h(n)=h(N-n-1)；而奇对称时，h(n)=-h(N-n-1)。再考虑到N可以为奇数或偶数，总共有四种类型的线性相位FIR滤波器。在讨论线性相位滤波器频率响应时需要引进幅特性正负号的概念。以往常设幅特性(Magnitude Response)为正数，因为幅特性的反号可以用相特性加减

8、来补偿。当相位特性要求线性，不得随便增减时，幅特性就必须分出正负,称为符幅特性(Amplitude Response)。14线性相位滤波器的特性线性相位滤波器的特性为了计算符幅特性，需要直接作解析计算，求出频率特性的解析式，并将它分解为幅度和相位特性，即 ，如果()满足线性相位条件，这时的A()就是符幅特性。对于长度为N的序列h(n)，先考虑两个对称的系数h(n)和h(N-1-n)，它们的傅立叶变换为（注意运用样本值相等的条件使两项合一）：15I I型线性相位滤波器的特性型线性相位滤波器的特性下面分别四种类型进行推导：下面分别四种类型进行推导：1 1类类型型I I，系数，系数对对称，称，N N

9、奇数奇数序列h的下标n从0起算，到N-1为止，中点位置是L=(N-1)/2。h(L)是一个孤项，其他可配成L对系数。将它们的符幅特性加起来，得到：16I I型线性相位滤波器的特性型线性相位滤波器的特性为了使四种类型滤波器的公式简明统一，引进两个参数。一个是令 ，它就是(7.1.5)式中的值，用以反映延迟的样本数或群延迟，它可以是分数。另一个是 。即把值向下取整，因为要把L用作下标，它必须是整数。当N为奇数时，=L。在推导中，还注意前一个求和号是从0到L-1,后一个求和号是从0到L,那是把孤项也当成成对项求和，然后再减去一个孤项h(L)。这样做的目的也是使四种类型滤波器的公式统一。17I I型线

10、性相位滤波器的特性型线性相位滤波器的特性显然，为相角项，相角-,它与成严格的线性关系；A()就是符幅项。可见，其符幅特性由L+1个余弦项叠加组成。n=0时的符幅项为2h(0)cos(L)，在=02频段间波动L个周期；余弦函数在=0，和2处都不等于零，因此类型I线性相位滤波器既可以用作低通滤波器(在=0处，幅特性不为零)；也可以用作高通滤波器(在=处，幅特性不为零)；且可用作带通和带阻滤波器。18IIII型线性相位滤波器的特性型线性相位滤波器的特性类型类型IIII，系数对称，系数对称，N N偶数偶数如果N为偶数，那末全部系数都可以配对，不会出现中心点的单项，一共有N/2组对称出现的系数，频率特性

11、就成为 此时=(N-1)/2将不是整数,也就是说，对称中心将在两个样本点的中间。其相角特性仍为()=-,它仍与成严格的线性关系。19IIII型线性相位滤波器的特性型线性相位滤波器的特性用整数 来表示求和号的上限。得出其符幅特性的表达式可见，它的符幅特性也由L+1个余弦分量构成。这样，不管N为奇数或偶数时，L都是整数，求和的项数都是L+1。由于=(N-1)/2是分数，所有余弦分量cos(-n)中都含一个0.5。在=处，就会出现符幅特性必定为零的情况。所以类型II不能用作高通和带阻滤波器。20IIIIII型线性相位滤波器的特性型线性相位滤波器的特性3 3类型类型IIIIII，系数反对称，系数反对称

12、，N N奇数奇数当系数序列反对称时，在(7.2.1)式的推导过程中，h(N-1-n)项将反号而有h(N-1-n)=-h(n)，两个反对称项合成的结果是：所以只要把类型I 结果中的cos换成sin，得到：21IIIIII型线性相位滤波器的特性型线性相位滤波器的特性当N为奇数时，=L=整数。相特性为()=/2-(N-1)/2=/2-,它虽然仍与成线性关系，但多了一个常数项/2。注意反对称序列的对称中心n=L处，序列的值h(L)=0。所以可以放入求和号中，得到它的符幅特性为：正弦函数在=0，和2处都等于零，因此类型III线性相位滤波器的符幅特性在=0，和2处都等于零。它不能用作低通，也不能用作高通滤

13、波器。故不宜用于选频滤波器。22IVIV型线性相位滤波器的特性型线性相位滤波器的特性4 4类型类型IVIV，系数反对称，系数反对称，N N偶数偶数与类型II滤波器相仿，若N为偶数，则共有N/2组数值相同反对称出现的系数，不出现单项，频率特性式就成为此处仍有=(N-1)/2，因此将不是整数,带有小数部分0.5。故对称中心将在两个样本点的中间。它的相角特性为()=/2-,虽然仍与成线性关系，但多了常数项/2，故不通过原点。23IVIV型线性相位滤波器的特性型线性相位滤波器的特性类型IV的滤波器符幅特性为：因为L是整数，它的符幅特性也由L+1=N/2个正弦分量构成。因为sin(0)=0，在=0处，符

14、幅特性必定为A(0)=0。这样的特性显然不能用作低通和带阻滤波器。实际上这种类型的滤波器也不适用于选频滤波。把这四种类型的线性相位滤波器的相角特性和符幅特性用对比的方法画在一张图上如表。242526符幅特性的计算举例符幅特性的计算举例例7.2.1 设系统的系数向量为b=-1,2,4,2,-1，a=1，求其幅频特性及符幅特性，并作比较。解：先用手工计算，因为a=1,说明这是一个FIR滤波器，观察系数向量b，其长度N=5，阶次为M=N-1=4，传递函数为H(z)=-1+2z-1+4z-2+2 z-3 z 4注意h(n)=bn，其频率响应为：调用abs(H)函数只能求出常为正的幅特性。27幅特性和符

15、幅特性的算法对比幅特性和符幅特性的算法对比要求出符幅特性，必须把对称系数配对计算。将n=0,4组成一对，n=1,3组成一对，n=2处是对称中心，h(2)独立计算，于是有：由于对称配置，的高次和低次项搭配成为余弦，合成后的相位项都是 ，故其符幅特性是这三个幅度项的代数和：把本题的幅特性和符幅特性画成曲线如下图28幅特性和符幅特性的结果对比幅特性和符幅特性的结果对比左上是幅特性，右上是符幅特性；可见符幅特性取绝对值，就得幅特性。下方是相特性，右图是线性相位，左图的第二个突跳与幅特性反号对应。29计算符幅特性的计算符幅特性的MATLABMATLAB程序程序本书编了一个通用的程序，它能自动判别滤波器的

16、类型并能在(7.2.2a)(7.2.5a)四个公式中，自动选择正确的符幅特性公式进行计算。把这个函数程序命名为amplres.m，其调用方法：A,w,type,tao=amplres(h)其中：h=FIR滤波器的脉冲响应或分子系数向量 A滤波器的符幅特性 w取的频率向量，在0到pi之间取501个点 type线性相位滤波器的类型 tao符幅特性的群迟延30计算符幅特性的通用子程序子程序amplres.m的核心语句为：N=length(h);tao=(N-1)/2;L=floor(n-1)/2);n=1:L+1;w=0:500*pi/500;%滤波器频率向量if all(abs(h(n)-h(N-

17、n+1)1e-10)%若系数对称 A=2*h(n)*cos(N+1)/2-n)*w)-mod(N,2)*h(L+1);%对称条件下计算A（两种类型）elseif all(abs(h(n)+h(N-n+1)1e-10)&(h(L+1)*mod(N,2)=0)%系数为反对称条件下计算A（两种类型）A=2*h(n)*sin(N+1)/2-n)*w);end 31线性相位滤波器的特性在阅读这个程序时，注意以下几点：（1）。首先要注意程序中的n比公式中的n大1，这是因为MATLAB下标从1开始，而公式中的下标从零开始；公式中的求和指标从n=0到L在程序中是n=1到L+1，依此类推。（2）。关于判断序列对

18、称和反对称的条件语句，似乎该用all(h(n)-h(N-n+1)=0，程序中则采用了all(abs(h(n)-h(N-n+1)1e-10)。这是因为生成h(n)和h(N-n+1)中会有计算误差，为了防止造成误判，要允许两者有一点差别；32线性相位滤波器的特性（3）。计算符幅特性时，用行向量2h(n)与列向量cos(-n)相乘，实现求和号下的相乘累加运算。（4）。程序中采用了mod(N,2)函数，来区分奇偶的不同效果，避免用多个条件转移语句。在N为奇数时它等于0，而在N为偶数时它等于1。用mod(N,2)去乘一个变量，相当于在N为偶数时把这项变量取消。所以使适用于四种情况的计算程序程序变得如此简

19、练。33符幅特性算例给出由序列h0=3,-1,2,-3 构成的四种类型的系数向量：(a),h1=3,-1,2,-3,5,-3,2,-1,3;(b),h2=3,-1,2,-3,-3,2,-1,3;(c),h3=3,-1,2,-3,0,3,-2,1,-3;(d),h4=3,-1,2,-3,03,-2,1,-3;分别求出它们的符幅特性曲线，进行比较。解：对(a)的语句如下，其余类推。h1=3,-1,2,-3,5,-3,2,-1,3;A1,w1,typea,tao1=amplres(h1);typea 34符幅特性算例的结果运行程序得到typea=1,typeb=2,typec=3 和typed=4，

20、相应的符幅特性曲线见图。可验证表7.2.1:类型I幅特性关于=对称，在=0和=处可以取任何值；类型II幅特性关于=反对称，在=0可以取任何值，在=处必定为零；类型III幅特性关于=反对称，在=0和=处都必定为零；类型IV幅特性关于=对称，在=0处必定为零，而=处可以取任何值；35线性相位滤波器的特性线性相位滤波器的特性36线性相位滤波器的零极点分布线性相位滤波器的零极点分布零极点应该从FIR滤波器正幂传递函数来判断。所以它具有与零点数目相同的M个极点，它们集中在z平面的原点处，成为M重极点。可以证明，对于具有对称性的系数向量，如果z1是(7.2.13)的分子多项式的根，那末它的倒数q1=1/z

21、1也是这个多项式的根。37线性相位滤波器的零极点分布线性相位滤波器的零极点分布如果z1是复数零点，根据实系数多项式的性质，它的共軛数z1*也是复数零点；现在又证明了它的倒数1/z1也是零点，因而倒数的共軛数z1*也是复数零点。四个复数零点必成组出现，如图。387.3 7.3 设计设计FIRFIR滤波器的窗函数法滤波器的窗函数法第五章中研究了理想低通滤波器的频率特性应为：振幅(符幅)特性在通带内为1，阻带内为0；在通带内的相位特性与成线性关系。即 对应的理想单位脉冲响应(滤波器系数向量)为理想滤波器之所以不能实现，因为(1)其长度是无限的;(2)它是非因果序列。39设计设计FIRFIR滤波器的窗

22、函数法滤波器的窗函数法为了构造特性接近理想低通滤波器的实际滤波器，应当使它具备如下特性：(1).因果序列，在n0处无样本，起点在n=0处；(2).有限长度，长度为N(阶次为N-1),系数下标从0到N-1;(3).系数对称，以保证线性相位。对称中心应在n0=(N-1)/2处；也即把理想脉冲响应右移n0;用长N对称中心在n0=(N-1)/2的矩形窗RN(n)与hd(n)相乘，截取出的实际系数序列h(n)，如下图40矩形窗截取实际滤波器系数矩形窗截取实际滤波器系数41矩形窗截取的频率特性计算矩形窗截取的频率特性计算上图为理想滤波器的脉冲响应hd(n)和其频谱；中图为矩形窗函数RN(n)和其频谱；下图

23、时域为hd(n)和RN(n)的乘积，其频谱为上两频谱的卷积。出现波动。42矩形窗截取的频率特性计算矩形窗截取的频率特性计算求理想脉冲函数加窗后的频谱特性(图7.3.1(f)可以利用频域卷积定理。时域相乘h(n)=hd(n)RN(n),对应于频域卷积，做卷积必须考虑函数的正负号，所以其频谱需要用符幅特性来表示：其中Hd(j)为hd(n)的频谱WRN(j)为RN(n)的频谱，用符幅特性表示。43矩形窗截取的频率特性计算矩形窗截取的频率特性计算分析其幅度，可得截断后滤波器的频率特性这一积分的几何意义见图。实际的H()和理想的Hd()有两点差别：（1）通带内增加了波动，最大的峰值在=c-2/N处。阻带

24、内产生了余振，最大的负峰在=c+2/N处；（2）在理想特性不连续点=c附近形成过渡带。过渡带的宽度近似等于ARN()主瓣宽度。44(a).理想滤波器频谱；(b).矩形窗函数符幅特性；此时相当于=0，它与(a)卷积的结果为矩形窗截断区的正面积减负面积，设为1。(c).=c时卷积输出应为阴影面积的积分，此处应为0.5。(d).c时少一块负面积，卷积输出大于1；(e)L+1，即方程数比未知数的数目多，形成了所谓超定方程组。那就不可能找到精确满足这些方程组的系数d(n)，只能找到最近似地满足这些方程的最小二乘解。人们总想用最少的系数，在尽量多的频点上近似逼近预期幅特性。这所以种情况是最常见的。如果K0

25、的频段引入90度相移，而对0的频段引入90度相移，所以也称它为90度相移器。想把信号带宽减小一半，可以利用希尔波特变换器的这个特性。137希尔波特变换器设计希尔波特变换器设计设x(n)的频谱为 ，让它通过希尔波特变换器，则输出信号y(n)的频谱为：然后把x(n)和y(n)组成复数信号z(n)=x(n)+jy(n)，则z(n)的频谱函数为:138希尔波特变换器设计希尔波特变换器设计它在0的频段内的频谱.，所以其频带宽度减小了一半。由于频谱的共軛对称性，它也保留了原来频谱中的所有信息，可以在接收端完全恢复原来的信号内容。同样由于非因果和非绝对可加的问题，理想希尔波特变换器是无法实现的。先按因果性的

26、要求，把脉冲序列移后-0.5(N-1)拍，0.5(N-1)必须是整数，以避免x(n)+jy(n)时间匹配的困难，所以实际上只取N为奇数。139希尔波特变换器设计希尔波特变换器设计实际使用希尔波特变换器的框图如下：因为在变换器中有迟延，y(n)延迟了L拍，需要把x(n)也迟延L拍，才能把两者正交叠加为z(n)。140希尔波特变换器设计希尔波特变换器设计因此考虑因果性后的变换器预期特性为：它的IDFT应该就是此滤波器的系数向量：很容易证明，其系数序列是反对称的。而序列的长度N为奇数，故它是类型III的FIR线性相位滤波器。141希尔波特变换器设计算例希尔波特变换器设计算例例7.6.2 用矩形窗和汉

27、宁窗设计一个10阶的希尔波特变换器。解：设计希尔波特变换器的程序语句很简单。只是检验语句有些不同。因为它在正负频率区间的相位特性不同，要检验这一点，就必须计算它在全频域的频率特性，因此在freqz函数中加变元whole；另外，要看出希尔波特变换器构成的相位突变，就不希望出现因MATLAB取 相 角 主 值 带 来 的 突 变。因 此 用 了unwrap(angle(H)代替通常的angle(H)。本例的程序hc762如下。142希尔波特变换器设计算例希尔波特变换器设计算例N=input(滤波器长度(输入奇数);tau=(N-1)/2;n=0:2*tau+1e-10;hd=2*(sin(n-tau)*pi/2).2./(n-tau)/pi;%矩形窗截断的系数向量hhn=hd.*hann(N);%汉宁窗截断的系数向量H,w=freqz(hd,1,whole);%求频率响应Hhn,w=freqz(hhn,1,whole);程序运行的结果如下图。它说明了III型滤波器的特征恰好符合希尔波特变换器的需要。143希尔波特变换器设计算例希尔波特变换器设计算例144