41FIR数字滤波器设计的窗函数法-数字信号处理

第4章有限冲激响应数字滤波器设计窗函数法设计FIR数字滤波器频率抽样法最佳一致逼近法FIR数字滤波器数字滤波器设计：由给定的系统频率特性,确定M和N及系数ai,bjLTI系统iiNijjMjzazbzH101)(若ai等于零，则系统为FIR数字滤波器。若ai至少有一个非零，则系统为IIR数字滤波器。FIR系统与IIR系统比较数字信号处理研究的内容主要有以下三方面1．FIR系统只有零点，要取得好的衰减特性，一般要求h（z）的阶次要高，也即M要大。优点：①系统总是稳定的。②易实现线性相位。③允许设计多通带（或多阻带）滤波器，目前，FIRDF的设计方法主要是建立在对理想滤波器的频率特性作某种近似的基础上的。FIR系统与IIR系统比较2．IIR系统零、极点全部包括，所以易取得比较好的通带与阻带衰减特性，对于IIRDF的设计方法，主要是借助与模拟滤波器的设计方法，这些面向极点的设计方法不适用于仅包含零点的FIR系统。对于FIRDF的设计方法包括窗函数法、频率抽样法及最佳一致逼近法.4．1FIR数字滤波器设计的窗函数法对于一个理想低通数字滤波器，其频率特性为Hd（ejω），现假定其幅频特性|Hd（ejω）|=1，相频特性Ф(ω)=0，那么，该滤波器的单位抽样响应hd（n）是以hd（0）为对称的sinc函数，hd（0）=ωc/π，由于n的取值范围为-∞~+∞，所以这样的系统是非因果的，因此是物理不可实现的。4．1FIR数字滤波器设计的窗函数法如果我们将hd（n）截短，例如仅取hd（-M/2）…hd（0）…hd（M/2），并将截短后的hd（n）移位，得那么hd（n）是因果的，且为有限长，长度为M+1，令M0nnz)n(h)z(HMnMnhnhd,1,0),2()(4．1FIR数字滤波器设计的窗函数法可得所设计的滤波器的转移函数，H(z)的频率响应将近似Hd(ejw)，且是线性相位的。若Mnnhnhd,1,0)()()2()2sin(21)(2MnMndeenhcnjMjdcc2)(M4．1FIR数字滤波器设计的窗函数法设计一FIR低通滤波器，所希望的频率响应Hd（ejω）在０≤｜ω｜≤0.25π之间为１，在0.25π＜｜ω｜≤π之间为０，分别取Ｍ＝10，20，40，观察其幅频响应的特点由上式可得当M=10时，求得h（0）＝h（10）＝-０.045,h（1）＝h（9）＝０,h（2）＝h（8）＝０.075,h（3）＝h（7）＝０.1592,h（4）＝h（6）＝０.2251,h（5）=0.25)2Mn(]25.0)2Mnsin[()n(hd4．1FIR数字滤波器设计的窗函数法吉布斯现象的产生是由于对hd（n）突然截短的结果，将无穷长的hd（n）仅取长为0~M，等于在hd（n）上施加了长为(M+1)的矩形窗口，加窗的结果，等于Hd（ejω）和矩形窗频谱的卷积，Hd（ejω），矩形窗频谱D（ejω）为sinc(ω)函数，它有着较大的边瓣，正是这些边瓣在和Hd（ejω）卷积时产生了吉布斯现象。实际上，在对hd（n）自然截短时，H（ejω）是对Hd（ejω）在最小平方意义上的逼近。4．1FIR数字滤波器设计的窗函数法为了减少吉布斯现象，我们应选取边瓣较小的窗函数，如用汉宁窗或汉明窗来代替矩形窗。由于本设计方法是选取某种好的窗函数对hd（n）予以截短来得到h（n）的，所以该设计方法称为“窗函数法”又由于Hd（ejω）是以２π为周期的周期函数，hd（n）可看作Hd（ejω）在频域展开为傅立叶级数时的系数，所以该方法又称“傅立叶级数法”由上述设计过程可以看出，FIRDF设计的窗函数法不象IIRDF的设计那样能精确指定通带、阻带的边缘频率ωp,ωs，也不能精确地给定在这两个带内的衰减αp，αs，而是仅给出通带截止频率ωc（即ωp），其他几个参数是靠h（n）的长度M及所使用的窗函数的性能来决定的。选定了窗函数后，可以不断地变化M，以检查通带、阻带是否达到了希望的技术要求，直到满意时为止。这些工作在计算机上是极易实现的。4．1FIR数字滤波器设计的窗函数法4．2窗函数在信号处理中不可避免地要遇到数据截短问题，例如，在应用DFT时，数据x（n）总是有限长，在FIR滤波器设计中遇到了对理想滤波器抽样响应hd（n）的截短问题，在功率谱估计中也要遇到对自相关函数的截短问题，总之，我们在实际工作中所能处理的离散序列总是有限长，把一个长序列变成有限长短序列不可避免地要用到窗函数，因此，窗函数本身的研究及应用是信号处理中的一个基本问题。设xn为一长序列，是长度为N的窗函数，则)()()(nnxnxNdeWeXeXjjjN)()(21)()()2/sin()2/sin(e)()e(2/)1(j10jNeNWNNnnjn吉伯斯现象窗函数的频谱)2/sin()2/sin(e)e(2/)1(jjNWN2NN04NW3N)(WW矩形窗的幅度函数W(W)1.W(e)的主瓣宽度4/N2.旁瓣相对衰减为常数)π3()0(lg20NWWA=13dB3dB带宽，最大边瓣峰值A（dB）,边瓣谱峰渐近衰减速度D（dB/oct）一个理想的窗函数，应该具有最小的B和A及最大的D除以上三个指标外，我们对窗函数还有一些共同的要求ω(n)应是非负的实偶函数，且ω(n)从对称中心开始，应是非递增的。几种常用的窗函数矩形窗汉宁窗（Hanning）汉明窗（Hamming）布莱克曼窗（Blackman）凯泽窗（Kaiser）常用窗函数矩形w=boxcar(M+1))其他001][MkkwAp=20lg(1dp)0.82dB,As=20lg(ds)21dB10.50.091.09WWcNπ8.1dp=ds=0.09常用窗函数Hann(汉宁)窗(w=hanning(M+1))其他00)/π2cos(5.05.0][MkMkkw01020301w[k]kAp0.056dB,As44dBWWcNπ2.60.006411.00641-0.0064dp=ds=0.0064常用窗函数Hamming(汉明)窗(w=hamming(M+1))其他00)/π2cos(46.054.0][MkMkkw01020301w[k]kAp0.019dB,As53dBdp=ds=0.0022WWcNπ0.71.00220.002211-0.0022常用窗函数Blackman窗(w=blackman(M+1))其他00)/π4cos(08.0)/π2cos(5.042.0][MkMkMkkw01020301w[k]kAp0.0017dB，As74dBdp=ds=0.0002WWcNπ4.111.00020.00021-0.0002常用窗函数性质窗的类型主瓣宽度近似过渡带宽度dp，dsAp(dB)As(dB)矩形4/N1.8/N0.090.8221Hann8/N6.2/N0.00640.05644Hamming8/N7/N0.00220.01953Blackman12/N11.4/N0.00020.001774][][][dnwnhnhN例：用非矩形窗设计截频为Wc=0.5的线性相位理想低通，并将设计结果与利用矩形窗设计结果比较。解：(1)确定线性相位FIR滤波器类型W其他01)(cjdeAd(W)=M/2(2)确定理想滤波器的幅度函数Ad()和相位d()(3)计算IDTFT得hd[n])]2/([sinπ][ccdMncnh(4)截断hd[n]可以是汉宁窗、哈明窗或Blackman窗矩形窗和汉宁窗设计的低通对比00.250.50.751-80-60-44-200SquareHanning例：用非矩形窗设计截频为Wc=0.5的线性相位理想低通，并将设计结果与利用矩形窗设计结果比较。矩形窗和哈明窗设计的低通对比00.250.50.751-80-60-52-200SquareHamming例：用非矩形窗设计截频为Wc=0.5的线性相位理想低通，并将设计结果与利用矩形窗设计结果比较。矩形窗和Blackman窗设计的低通对比00.250.50.751-100-75-60-40-200SquareBlackman例：用非矩形窗设计截频为Wc=0.5的线性相位理想低通，并将设计结果与利用矩形窗设计结果比较。Kaiser(凯泽)窗(w=kaiser(N,beta))MkIMkIkw0,)()]/21[1(][020是一可调参数，调节窗函数的形状。I0(x):themodifiedzeroth-orderBesselfunction.I0(x)可用幂级数表示为210!)2/(1)(nxxInn一般求20项就能达到所需精度。定义Kaiser(凯泽)窗——取不同值时窗的形状0246810121416182000.51=00246810121416182000.51=10246810121416182000.51=3Kaiser(凯泽)窗A=20lg(min{dp，ds})21,285.295.7spAAMWW21,05021),21(07886.0)21(5842.050),7.8(1102.04.0AAAAAAN与的确定Kaiser(凯泽)窗Kaiser窗设计FIR滤波器的MATLAB函数[M,Wc,beta,ftype]=kaiserord(f,a,dev)f:表示需设计的FIR滤波器的频带。a:B个元素的向量,表示FIR滤波器在B个频带中的幅度值。一般对通带取值为1，阻带取值为0。dev:B个元素的向量,表示FIR滤波器在B个频带中的波动值。返回参数M及beta，分别表示FIR滤波器阶数M及Kaiser窗的参数。返回参数Wc和ftype是函数fir1的调用参数。h=fir1(M,Wc,'ftype',window)M表示滤波器的阶数。Wc表示理想FIR滤波器的B个频带。ftype表示滤波器的类型，缺省值为空。window是一长度为N+1的向量,Hamming窗。Kaiser(凯泽)窗获得FIR滤波器的单位脉冲响应例：用Kaiser窗设计满足下列指标的I型线性相位FIR低通滤波器。Wp=0.3,Ws=0.5，Ap=0.1dB,As=40dB。解：(1)由给定指标确定待逼近理想低通的截频Wc由于理想低通滤波器的|H(ejW)|在截频Wc处收敛于0.5，因此常将截频Wc取在过渡带的中点Wc=(Wp+Ws)/2=0.4解：(2)由给定指标确定Kaiser窗的参数N和0114.0101p05.0pAd01.010s05.0sAdA=20lg(min{dp，ds})=As=40dB3.22285.295.7spWWAMI型线性相位滤波器阶数必须是偶数，取M=243953.3)21(07886.0)21(5842.04.0AA例：用Kaiser窗设计满足下列指标的I型线性相位FIR低通滤波器。Wp=0.3,Ws=0.5，Ap=0.1dB,As=40dB。解：(3)设计截频Wc=0.4的I型线性相位FIR低通滤波器其他01)(cd)]5.0([SaπccMkWWd(W)=0.5MWdee)(π21][j)(jdππddkAkh()][)12(π4.0Sa4.0][][][2525dkwkkwkhkh例：用Kaiser窗设计满足下列指标的I型线性相位FIR低通滤波器。Wp=0.3,Ws=0.5，Ap=0.1dB,As=40dB。解：0122400.20.