FIR数字滤波器设计.docVIP

FIR数字滤波器设计摘要 数字滤波器的输入输出均为数字信号，信号通过数字滤波器后，可以改变频率成分的相对比例或滤除某些频率成分。数字滤波器可以分为IIR数字滤波器和FIR数字滤波器。文中只介绍FIR数字滤波器的设计，可以根据所给定的频率特性直接设计FIR数字滤波器。FIR数字滤波器在保证幅度特性满足要求的同时，能够做到严格的线性特性。本文采用了窗函数法、频率采样法以及基于firls函数和remez函数的最优化方法设计FIR滤波器。在文中对FIR滤波器进行了详细的理论分析，并且对应于每种方法都给出了设计实例。通过编写MATLAB语言程序，运行程序，得到幅频和相频特性图。对于窗函数和firls函数设计的滤波器，还通过建立Simulink系统模块进行仿真，观察滤波器滤波情况。关键词 窗函数法；频率采样法；最优化方法；仿真Design of FIR Digital FilterAbstract Firstly, the paper gives a brief description of digital filter and MATLAB language.The digital filter can be divided into IIR digital filter and FIR digital filter.The paper only give the digital filter of FIR, as for the FIR filter, we can design it directly based on the giving frequency,in the paper it uses the design of the window function，frequency sampling method and optimization method.As for all filter design method above, the paper all describes in detail. And for each method, the paper also gives some examples. At last, through the preparation of the relevant procedures in MATLAB simulation software, operating procedures and records, observation of frequency and phase of plans to FIR digital filter in the MATLAB software in the implementation process.Keywords 目 次引言11数字滤波器21.1数字滤波器的概述21.2数字滤波器的分类21.3数字滤波器设计指标32 FIR数字滤波器设计基础42.1 FIR数字滤波器的特点42.2 FIR数字滤波器的线性相位条件：42.3 FIR数字滤波器的基本结构53 FIR数字滤波器设计93.1 窗函数法设计FIR数字滤波器93.2频率采样法设计FIR数字滤波器153.3最优化方法设计FIR数字滤波器16结 论21参考文献22致 谢2323引言随着信息时代与数字技术的发展，数字信号处理己逐渐发展成为当今极其重要的学科与技术领域之一。数字信号处理在通信、语音、图像、自动控制雷达、军事、航空航天、医疗和家用电器等众多领域得到了广泛的应用。在数字信号处理的基本方法中，通常会涉及到变换、滤波、频谱分析、调制解调和编码解码等处理。其中，滤波是应用非常广泛的一个环节，数字滤波器的相关理论也一直都是人们研究的重点之一。数字滤波器是数字信号处理的重要基础，在对信号的滤波、检测及参数的估计等信号应用中，数字滤波器是使用最为广泛的一种线性系统。数字滤波器根据其单位冲击响应函数的时域特性可分为两类：无限冲击响应(IIR)数字滤波器和有限冲击响应(FIR)数字滤波器。与IIR数字滤波器相比,FIR数字滤波器的实现是非递归的，稳定性好，精度高；更重要的是FIR数字滤波器在满足幅度响应要求的同时，可以获得严格的线性相位。因此，它在高保真的信号处理中，如数字音频、图像处理、数据传输和生物医学等领域得到广泛应用。1数字滤波器1.1数字滤波器的概述所谓数字滤波器，是指输入输出均为数字信号，通过一定的运算关系，改变输入信号中所含频率成分的相对比例，或则滤除某些频率成分的器件3。数字滤波器具有稳定性高、精度高、灵活性大等突出优点。对于数字滤波器而言，若系统函数为H(z),其脉冲响应为h(n),输入时间序列为x(n),则它们在时域内的关系式如下： (1-1)在Z域内，输入和输出存在如下关系： (1-2)式中， X(z)、Y(z)分别为x(n)和y(n)的Z变换。在频域内，输入和输出则存在如下关系： (1-3)式中，是数字滤波器的频率特性；、分别为x(n)和y(n)的频谱，而为数字角频率。1.2数字滤波器的分类数字滤波器可以有很多种分类方法，但总体上可分为两大类。一类称为经典滤波器，即一般的滤波器，其特点是输入信号中的有用成分和希望滤除的成分占用不同的频带，通过合适的选频滤波器可以实现滤波4。例如，若输入信号中有干扰，信号和干扰的频带互不重叠，则可滤出信号中的干扰得到纯信号。但是，如果输入信号中信号和干扰的频带相互重叠，则干扰就不能被有效的滤除。另一类称为现代滤波器，如维纳滤波器、卡尔曼滤波器等，其输入信号中有用信号和希望滤除的频带成分重叠。对于经典滤波器，从频域上也可以分为低通、高通、带通和带阻滤波器。从时域特性上看，数字滤波器还可以分为有限脉冲响应（FIR，finite impulse response）数字滤波器和无限脉冲响应（IIR, infinite impulse response）数字滤波器5。对于有限脉冲响应（FIR）数字滤波器，其输出y(n)只取决于有限个过去和现在的输入，x(n),x(n-1),x(n-m),滤波器的输入输出关系可表示为 (1-4)对于无限脉冲响应（IIR）数字滤波器，它的输出不仅取决于过去和现在的输入，而且还取决于过去的输出，其差分方程为 (1-5)该差分方程的单位冲激响应是无限延续的。1.3数字滤波器设计指标设数字滤波器的传输函数用下式表示： (1-6)式中，H(e)为幅频特性，为相频特性6。幅频特性表示信号通过滤波器后各频率成分的衰减情况，相频特性则反映各频率成分通过滤波器后在时间上的延时情况。通常，选频滤波器的指标要求都以幅频特性给出，对相频特性不作要求，如果需要对输出波形有严格要求，如语音合成、波形传输等，则要求设计线性相位数字滤波器7。数字滤波器的参数指标是、和。和分别称为通带截止频率和阻带截止频率。通带和阻带内允许的衰减一般用分贝数表示，通带内允许的最大衰减用表示，阻带内允许的最小衰减用表示，和分别定义为 dB (1-7） dB (1-8)式中均假定已被归一化为12 FIR数字滤波器设计基础2.1 FIR数字滤波器的特点FIR滤波器在保证幅度特性的同时，很容易做到严格的线性相位特性。在数字滤波器中，FIR滤波器的最主要特点是没有反馈回路，故不存在不稳定的问题；同时，在幅度特性可以任意设置的同时，保证了精确的线性相位。稳定和线性相位是FIR滤波器的突出优点。另外还有以下特点：设计方式是线性的；硬件容易实现；滤波器过渡过程具有有限区间；相对IIR滤波器而言，阶次较高，其延迟也要比同样性能的IIR滤波器大得多8。2.2 FIR数字滤波器的线性相位条件：设滤波器单位脉冲响应的长度为N，系统函数为 (2-1)由此式可见，H(z)是的(N-1)次多项式,它在Z平面上有(N-1)个零点，原点z=0是(N-1)阶重极点，位于r =1的单位圆内，系统永远稳定。稳定性和线性相位特性是FIR滤波器的突出优点。 FIR滤波器的设计任务是选择有线长度的h(n)，使传输函数H(e)满足要求。 线性相位条件：对于长度为N的h(n)，传输函数为 (2-2) (2-3)式中,称为幅度特性，称为相位特性。线性相位是指相位函数满足如下特性：或, 是起始相位，为常数，一般称第一种情况为第一类线性相位，称第二种情况为第二类线性相位。满足第一类线性相位的充要条件是：h(n)为实序列，并且对(N-1)/2偶对称，即；满足第二类线性相位的充要条件是：h(n)为实序列，并且对(N-1)/2奇对称，即。2.3 FIR数字滤波器的基本结构FIR滤波器的基本结构有以下几种：直接型、级联型、线性相位型、频率采样型。1. 直接型设FIR滤波器的单位冲击响应h(n)为一个长度为N的序列，则滤波器系统函数为： （2-4）表示这一系统输入输出关系的差分方程为 (2-5)直接由差分方程可得出对应的网络结构如图2-1所示： 图2-1 FIR滤波器的直接型结构直接型结构的优点：简单直观，乘法运算量较少。缺点：调整零点较难。2级联型当需要控制滤波器的传输零点时，可将H(z)分解为实系数二阶因子的乘积形式： (2-6)式中，为的变换，为实数。级联型结构如图2-2所示：图2-2 FIR滤波器的级联型结构该结构的优点：调整零点比直接型方便。缺点：中的系数比直接型多，因而需要的乘法器多。当的阶次高时，也不易分解。3.线性相位型结构FIR滤波器的线性相位结构有偶对称和奇对称，不论为偶对称还是奇对称都有：当N为偶数时，系统函数为 (2-7)当N为奇数时，系统函数为 (2-8)对这两种情况，都可以用FIR直接型实现，其信号流图如图2-3所示。(a)N为偶数(b)N为奇数图2-3 线性相位型结构 这种结构在本质上是直接型，但乘法次数比直接型省了一半。4频率采样型频率采样型结构是一种用系数将滤波器参数化的一种实现结构。一个有限长序列可以由相同长度频域采样值惟一确定。系统函数在单位圆上作N等分取样就是单位取样相应h(n)的离散傅里叶变换。与系统函数之间的关系可用内插公式表示： (2-9)式中 这样，是由梳状滤波器和N个一阶网络的并联结构进行级联而成的，其网络结构（信号流图）如图2-3所示。是一个梳妆网络，其零点为, k= 0, 1，2，N-1刚好和极点一样，等间隔地分布在单位圆上。理论上，极点和零点相互抵消，保证了网络的稳定性。图2-5 FIR滤波器的频率采样结构频率采样结构的优点：1)在频率采样点，只要调整就可以有效地调整频响特性。2）只要长度N相同，对于任何频响，其梳状滤波器部分和N个一阶网络部分完全相同，只是各支路增益不同。相同部分便于标准化、模块化。缺点：1）寄存器长度都是有限的，零、级点可能不能正好抵消，造成系统不稳定。2）当N很大时，其结构很复杂，需要的乘法器和延时单元很多。3 FIR数字滤波器设计 FIR数字滤波器的设计方法有窗函数法、频率采样法和基于firls函数和remez函数的最优化方法。MATLAB语言中的数字信号处理工具箱，提供了一些滤波器的函数，使FIR滤波器的运算更加方便和快捷。在MATLAB中提供的滤波函数有fir1(),此函数以经典的方法实现加窗线性相位FIR数字滤波器设计，可以设计出低通、高通、带通和带阻滤波器；fir2函数设计的FIR滤波器，其滤波的频率特性由矢量f和m决定，f和m分别为滤波器的期望幅频响应的频率相量和幅值相量。Firls()和remez()的基本格式用于设计I型和II型线性相位FIR滤波器，I型和II型的区别是偶函数还是奇函数。freqz()用于求数字滤波器的频率响应。并且提供了各种窗函数的函数，比如，hamming()是海明窗函数，hanning()是汉宁窗函数，kaiser()是凯泽窗函数，使在设计的过程中，不用自己重新设计窗函数。Simulink是MATLAB众多工具包中的一员，对于建模，Simulink提供了一个图形化的用户界面（GUI）。Simulink包括一个复杂的由接收器、信号源、线性和非线性组件以及连接件组成的模块库。定义完一个模型后，就可以通过Simulink的菜单或者在MATLAB的命令窗口输入命令对它进行仿真。使用Scopes或者其他的显示模块，可以在运行仿真时观察到仿真的结果。另外，还可以在仿真时改变参数，并且立即就可看到变化。 3.1 窗函数法设计FIR数字滤波器设我们所要设计的FIR滤波器的传输函数是(e),(n)是与其对应的单位脉冲响应，因此 （3-1） （3-2)如果我们能够在已知的情况下，求出，经过Z变换可得到滤波器的系统函数。通常情况下理想数字滤波器的单位脉冲相应是无限长的，且是非因果序列。获得有限脉冲响应滤波器的一种可能方法是对截取一段来近似代替，可是这样会改变原来的滤波器指标，出现吉布斯效应误差。窗函数法就是用被称为窗函数的有限加权序列w(n)来修正式(1)的傅里叶基数以求得要求的有限脉冲响应序列，即 （3-3）w(n)是有限长序列，当nN-1时，w(n)=0。 这种方法的重点在于选择某种合适的窗函数。要求窗函数主瓣宽度尽可能窄，以获得最小的过渡带；旁瓣相对值尽可能小，以使得通带波纹小，并且阻带衰减大。下面介绍几种常用的窗函数：1.矩形窗(Rectangle Window) （3-4）其频率函数为： （3-5）2.三角形窗(Bartlett Window) （3-6）其频率函数为： （3-7)3.汉宁(Hanning)窗，又称升余弦窗 （3-8）利用傅里叶变换得到频率函数为： （3-9）当时，所以窗函数的幅度函数为 (3-10)4.汉明(Hamming)窗，又称改进的升余弦窗 （3-11）其幅度函数为： （3-12）5.布莱克曼(Blankman)窗，又称二阶升余弦窗 （3-13)其幅度函数为： (3-14）6.凯泽(Kaiser)窗 （3-15）其中：是一个可自由选择的参数，I0( x)是第一类修正零阶贝塞尔函数10.上述窗函数的基本参数如下表窗函数旁瓣峰值幅度/db过渡带宽阻带最小衰减/db矩形窗-134/N-21三角形窗-268/N-25汉宁窗-318/N-44汉明窗-408/N-53布莱克曼窗-5712/N-74凯泽窗-5710/N-80窗函数法设计滤波器的步骤：1）根据技术要求确定待求滤波器的单位取样响应。2）根据对过渡带和阻带衰减的要求，选择窗函数的形式，并估计窗口长度N。3）计算滤波器的单位取样响应h(n)： （3-16)式中,是前面所选择好的窗函数。4）检验技术指标是否满足要求。根据下式计算： (3-17)如果不满足要求，根据具体情况重复步骤(2)(3)(4)步，直到满足要求为止。本文以一个FIR滤波器的设计为例说明如何使用MATLAB设计数字滤波器 设计实例：用窗函数法设计线性相位FIR低通数字滤波器，要求通带截止频率Wp=0.4*,阻带截止频率Ws=0.5*, 通带衰减不大于3db，阻带衰减不小于40db。1.程序如下：Wp=0.4*pi;Ws=0.5*pi;Wdel=Ws-Wp;N=ceil(8*pi/Wdel);Wn=(0.4+0.5)*pi/2;window=hanning(N+1);b=fir1(N,Wn/pi,window);freqz(b,1,512)程序执行后得幅频和相频如下图所示：2用Simulink仿真对窗函数设计的低通滤波器进行滤波：FIR低通滤波器的Simulink仿真模块图如下所示：输入信号：s(t)=sin(30t)+sin(56t)+sin(70t)数字滤波器参数设置：已知滤波器的阶数n=80。本例中，首先在Filter Type中选择lowpass；在Design Method选项中选择FIR Window，接着在Window选项中选取hann；指定Filter Order项中的Specify order为80；采样频率Fs=100Hz,截止频率Fc=22.5HZ。设置完以后点击窗口下方的Design Filter，在窗口上方就会看到所设计滤波器的幅频响应和相频响应曲线。设计界面如下图所示：信号源参数设置：三个信号源的频率分别是15Hz、28Hz、35Hz ： 运行仿真，得波形如下图所示： 15HZ信号波形 28HZ信号波形 35HZ信号波形 滤波前的叠加波形 滤波后波形由上图可知，15Hz、28Hz、35Hz叠加后通过截止频率为22.5HZ的FIR低通滤波器，输出与15Hz信号波形相似，其它的高频被阻隔。3.2频率采样法设计FIR数字滤波器1.对理想滤波器的系统函数Hd(z)进行频率采样以得到系统的理想频响Hd(ejw)的等间隔采样值H(k)。H(k)实际上是所要求的滤波器的单位采样响应(h(n)的离散傅里叶变换（DFT），如下试： (3-18) (3-19)为了减小H（k）的通带边缘由于抽样点的变化而引起的起伏振荡，可以增加过渡点，加宽过渡带以减小通带的起伏。每一个抽样值产生一个与sin()/sin()成正比，并位移（2k）/N的频率响应，而H（k）与内插函数的线性组合就是FIR滤波器的频率响应，增加一点过渡可以使阻带衰减提高到-4454dB,二点过渡衰减-6575,三点过渡衰减-8595dB.如果不能使过渡带太宽，同时要求增大阻带衰减，可以增加取样点数N，但这样会增加计算量、延时和误差。频率取样型FIR滤波器设计步骤：(1)给定理想滤波器频率响应Hd(ejw)。(2)根据过渡带宽和阻带衰减确定过渡点数和h(n)的长度N。 (3-20)(3)由IFFT计算IDFT得到： (3-21)设计实例：率采样法设计一个带通滤波器，满足：低阻带边缘：w1s=0.2*;低通带边缘：w1p=0.35*;高通带边缘：w2p=0.65*;高阻带边缘：w2s=0.8*。设计过渡带中的频率样本值为t1和t2，取t1=0.109021,t2=0.59417456。 设计程序如下：M=40;al=(M-1)/2;l=0:M-1;t1=0.109021;t2=0.59417456;Hrs=zeros(1,5),t1,t2,ones(1,7),t2,t1,zeros(1,9),t1,t2,ones(1,7),t2,t1,zeros(1,4);k1=0:floor(M-1)/2);k2=floor(M-1)/2)+1:M-1;angh=-al*(2*pi)/M*k1,al*(2*pi)/M*(M-k2);H=Hrs.*exp(j*angh);h=real(ifft(H,M);freqz(h,1,512,1000)实验得幅频相频特性如下图所示：3.3最优化方法设计FIR数字滤波器MATLAB7.0信号处理工具箱提供了通用的工具函数firls和remez，采用不同的优化方法可以设计最优滤波器和多频带滤波器。 firls是fir1和fir2函数的扩展，其基本设计准则是利用最小二乘法使期望的频率响应和实际的频率响应之间的整体误差最小。Remez函数则利用Parks-McClellan算法设计线性相位FIR滤波器。该算法能使期望频率响应和实际频率响应之间的最大误差最小。Firls函数设计的滤波器在频率响应上出现等波纹，因此也叫等波纹滤波器。函数firls、remez调用的语法规则相同。(1)firls。其调用格式为. b=firls(n,f,a). b=firls(n,f,a,w). b=firls(n,f,a,ftype). b=firls(n,f,a,w,ftype)该函数用于设计n阶FIR滤波器，其幅频特性由f和a向量确定，f是频率点向量，其范围为01，是递增向量；a是指定频率点的幅度响应，与f必须同长度；b为函数返回的滤波器系数，长度为n+1,且具有偶对称关系： b(k)=-b(n+2-k), k=1,2,n+1 b=firls(n,f,a,w)则使用权系数w给误差加权。 b=firls(n,f,a,w,ftype)调用形式中，参数ftype用于指定所设计的滤波器类型，ftype=Hilbert,为奇对称的线性相位滤波器，返回的滤波器系数满足b(k)=-b(n+2-k), k=1,2,n+1；ftype=differentiator,则采用特殊加权技术，生成奇对称的线性相位滤波器，使低频段误差大大小于高频段误差。(2)函数remez的基本调用格式为 b=remez(n,f,a)其中各项意义同函数firls所述。设计实例：分别使用函数firls和remez函数设计20阶带通滤波器，阻带频率为ws1=0.3*,ws2=0.7*,通带频率为wp1=0.4*,wp2=0.6*,并画出幅频特性曲线。1.程序如下所示：f=0 0.3 0.4 0.6 0.7 1;m=0 0 1 1 0 0;n=20;b=firls(n,f,m); %firls函数设计FIR数字滤波器h,w=freqz(b); bb=remez(n,f,m); %remez函数设计FIR数字滤波器hh,w=freqz(bb);axes(position,0.2 0.2 0.5 0.5); plot(w/pi,abs(h),b:,w/pi,abs(hh),r-);xlabel(frequency);ylabel(magnitude);grid on;实验的幅频特性曲线如下：2.用Simulink仿真对用firls函数设计的带通滤波器进行滤波：仿真模块如下： 输入信号：s(t)=sin(20t)+sin(50t)+sin(90t) 数字滤波器参数设置：已知滤波器的阶数n=20。本例中，首先在Filter Type中选择bandpass；在Design Method选项中选择FIR least-squares;指定Filter Order项中的Specify order为20；采样频率F