基于MATLAB函数的滤波器设计技术

采用MATLAB滤波器函数快速设计数字滤波器的方法IIR滤波器设计•基于IIR滤波器的特性，必须首先理解滤波器的3种描述形式并建立其关系：•冲击响应：线性系统时域基本的输入输出描述•频率响应：包括幅值响应和相位响应，是线性时不变系统频域基本的输入输出描述•零、极点图：线性系统频域基本的传输函数描述常用四类IIR滤波器•巴特沃斯（设计函数butter/buttap/buttord）•切比雪夫I型（设计函数cheb1ap/cheby1/cheby1ord）•切比雪夫II型（设计函数cheb2ap/cheby2/cheby2ord）•椭圆（设计函数ellipap/ellip/ellipord）butter：巴特沃思模拟和数字滤波器设计函数•[b,a]=butter(n,Wn)•[b,a]=butter(n,Wn,’ftype’)•[b,a]=butter(n,Wn,’ftype’,’s’)•[z,p,k]=butter(…)•[A,B,C,D]=butter(…)•参数含义：n－滤波器阶数；Wn－归一化截止频率（若为二元向量则给出2n阶数字带通滤波器）ftype=high－高通数字滤波器ftype=stop,Wn为二元向量－带阻数字滤波器•返回变量b,a－滤波器分子分母多项式系数z,p,k－滤波器零极点及增益buttap:巴特沃斯模拟低通滤波器原型设计函数•[z,p,k]=buttap(n)•参数含义：n－滤波器阶数•返回变量z,p,k－滤波器零点、极点及增益buttord:巴特沃斯滤波器阶数和截止频率计算函数•[n,Wn]=buttord(Wp,Ws,Rp,Rs)•[n,Wn]=buttord(Wp,Ws,Rp,Rs,’s’)•参数含义：Wp－通带转折频率Ws－阻带转折频率Rp－通带波纹系数Rs－阻带衰减系数cheb1ap:切比雪夫Ⅰ型模拟低通滤波器原型•[z,p,k]=cheb1ap(n,Rp)•参数含义：n－阶数Rp－通带波纹系数(dB)•返回变量：z－空矩阵（没有零点）k－增益cheby1:切比雪夫Ⅰ型低通滤波器设计函数•[b,a]=cheby1(n,Rp,Wn)•[b,a]=cheby1(n,Rp,Wn,’ftype’)•[b,a]=cheby1(n,Rp,Wn,’s’)•[b,a]=cheby1(n,Rp,Wn,’ftype’,’s’)•[z,p,k]=cheby1(…)•[A,B,C,D]=cheby1(…)•参数含义：n－滤波器阶数Rp－通带波纹系数Wn－归一化截止频率（若为二元向量则给出2n阶数字带通滤波器）ftype=high－高通数字滤波器ftype=stop,Wn为二元向量－带阻数字滤波器•返回变量b,a－滤波器分子分母多项式系数z,p,k－滤波器零极点及增益A,B,C,D－滤波器的状态空间型描述参数cheb1ord:切比雪夫Ⅰ型滤波器的阶数计算函数•[n,Wn]=cheb1ord(Wp,Ws,Rp,Rs)•[n,Wn]=cheb1ord(Wp,Ws,Rp,Rs,’s’)•参数含义：Wp－通带转折频率Ws－阻带转折频率Rp－通带波纹系数Rs－阻带衰减系数cheb2ap:切比雪夫Ⅱ型模拟低通滤波器原型•[z,p,k]=cheb2ap(n,Rp)•参数含义：n－阶数Rs－阻带波纹系数(dB)•返回变量：z,p,k－零点、极点、增益cheby2:切比雪夫Ⅱ型低通滤波器设计函数•[b,a]=cheby2(n,Rs,Wn)•[b,a]=cheby2(n,Rs,Wn,’ftype’)•[b,a]=cheby2(n,Rs,Wn,’s’)•[b,a]=cheby2(n,Rs,Wn,’ftype’,’s’)•[z,p,k]=cheby2(…)•[A,B,C,D]=cheby2(…)•参数含义：n－滤波器阶数Rs－阻带波纹系数Wn－归一化截止频率（若为二元向量则给出2n阶数字带通滤波器）ftype=high－高通数字滤波器ftype=stop,Wn为二元向量－带阻数字滤波器•返回变量b,a－滤波器分子分母多项式系数z,p,k－滤波器零极点及增益A,B,C,D－滤波器的状态空间型描述参数cheb2ord:切比雪夫Ⅱ型滤波器的阶数计算函数•[n,Wn]=cheb2ord(Wp,Ws,Rp,Rs)•[n,Wn]=cheb2ord(Wp,Ws,Rp,Rs,’s’)•参数含义：Wp－通带转折频率Ws－阻带转折频率Rp－通带波纹系数Rs－阻带衰减系数ellipap:椭圆模拟低通滤波器原型•[z,p,k]=ellipap(n,Rp,Rs)•参数含义：n－阶数Rp－通带波纹系数(dB)Rs－阻带波纹系数(dB)•返回变量：z,p,k－零点、极点、增益ellip:椭圆滤波器设计函数•[b,a]=ellip(n,Rp,Rs,Wn)•[b,a]=ellip(n,Rp,Rs,Wn,’ftype’)•[b,a]=ellip(n,Rp,Rs,Wn,’s’)•[b,a]=ellip(n,Rp,Rs,Wn,’ftype’,’s’)•[z,p,k]=ellip(…)•[A,B,C,D]=ellip(…)•参数含义：n－滤波器阶数Rp－通带波纹系数(dB)Rs－阻带波纹系数(dB)Wn－归一化截止频率（若为二元向量则给出2n阶数字带通滤波器）ftype=high－高通数字滤波器ftype=stop,Wn为二元向量－带阻数字滤波器•返回变量b,a－滤波器分子分母多项式系数z,p,k－滤波器零极点及增益A,B,C,D－滤波器的状态空间型描述参数ellipord:切比雪夫I型滤波器的阶数计算函数•[n,Wn]=ellipord(Wp,Ws,Rp,Rs)•[n,Wn]=ellipord(Wp,Ws,Rp,Rs,’s’)•参数含义：Wp－通带转折频率Ws－阻带转折频率Rp－通带波纹系数Rs－阻带衰减系数说明•MATLAB设计函数的规范化频率设置不是通常的归一化概念！•当设计的滤波器转折频率为Wp(rad)和Ws(rad)时，它们必须被π除。(*fs/2π=*fNyq/π)•例如：滤波器技术指标要求当抽样频率为5kHz时，系统截止频率为1kHz时，在MATLAB设计函数中必须以0.4赋值：ωc=2π(1000/5000)=2π(0.2)=0.4π典型IIR滤波器的比较类型特点|H(jΩ)|(Ω=1)butterΩ=0和Ω=∞处提供最好的Taylor级数逼近，整体平滑单调cheby1理想与实际频响绝对误差最小，阻带响应最平坦10-Rp/20cheby2优点：理想与实际频响绝对误差最小，通带响应无波纹最平坦缺点：阻带逼近零比cheby1慢，且阶次n为偶数时不逼近零10-Rs/20ellip满足相同要求所需阶次最低，相同指标（n,Rp,Rs）过渡带最窄10-Rp/20yulewalk:基于最小均方意义的递归滤波器设计函数•[b,a]=yulewalk(n,f,m)•designsrecursiveIIRdigitalfiltersusingaleast-squaresfittoaspecifiedfrequencyresponse.•fisavectoroffrequencypoints,specifiedintherangebetween0and1,where1correspondstohalfthesamplefrequency(theNyquistfrequency).Thefirstpointoffmustbe0andthelastpoint1,withallintermediatepointsinincreasingorder.Duplicatefrequencypointsareallowed,correspondingtostepsinthefrequencyresponse.•misavectorcontainingthedesiredmagnituderesponseatthepointsspecifiedinf.(给定频率点的理想幅值响应)•fandmmustbethesamelength.•plot(f,m)displaysthefiltershape.FIR滤波器设计设计步骤：确定指标→算法逼近→系统实现FIR滤波器设计指标•FIR滤波器均可满足线性相位要求，故只需考虑幅度指标：幅度指标的定义：•绝对指标：由幅度响应函数|H(ejω)|给出约束[0,ωp]定义为通带，ωp为通带截止频率δ1为通带波纹/理想通带响应最大容限[ωs,π]定义为阻带，ωs为阻带截止频率δ2为阻带波纹/理想阻带响应最大容限[ωp,ωs]定义为过渡带，Δω=ωs-ωp表示•相对指标：以分贝(dB)的形式给出约束其中：Rp为通带波纹，As为阻带衰减1maxmax1211()dB20lg0()(1)()120lg0(0)20lg0111jjjpsHeHeHeRA基于窗函数的设计•设一理想频率选择性滤波器的频响为Hd(ejω)|，其冲击响应序列为hd(n)，对hd(n)对称截断得到h(n)。•基于窗函数的设计思想就是从冲击响应序列h(n)入手，选择合适的窗函数和理想滤波器。•加窗运算：对hd(n)进行对称（双边）截断得到一个长度为N的因果线性相位FIR滤波器的冲击响应序列h(n)。•决定窗函数设计法设计的滤波器频响与理想滤波器逼进程度的两个因素：1、窗函数频响的主瓣宽度Δ2、窗函数频响的旁瓣峰值A•窗函数的性质：1、长度为N时主瓣宽度正比于1/N2、通过N可调节过渡带宽度：NΔω=C3、旁瓣峰值与N无关4、不同窗函数相比，旁瓣幅值减小，主瓣宽度增加5、旁瓣产生波动，在通带和阻带中形状相似1,()0,1()()()2sin[()]11sin[()]2()ccjcjdcjjjndddjjnccceHehnIDTFTHeHeedneednn()(),01()012()()()()()1()()()()()2ddjjjjjddhnnNhnNhnhnwnwnwNnHeHeWeHeWed常用窗函数性能比较窗函数名过渡带宽度最小阻带衰减/dB近似值(2π/N)准确值(2π/N)矩形窗20.921汉宁窗(hanning)43.144海明窗(hamming)43.353布莱克曼窗(blackman)65.574从应用的角度考虑，海明窗具有较好的综合性能。采用hamming窗的FIR滤波器使99.963%的能量集中在主瓣内，旁瓣峰值小于主瓣峰值的1%，有效地消除频谱泄漏，并能在一定程度上减少栅栏效应是各种窗函数中性能较好的一种，因此被广泛应用。kaiser窗•定义：•I0为第一类修正零阶贝塞尔(Bessel)函数•Kaiser窗函数是近似于给定旁瓣电平，是主板具有最大能量意义下的最佳窗函数。•β是取决于N的参数，一般根据阻带要求得到。对于期望阻带衰减As(dB)，As50dB时有经验公式：β=0.1102As-0.95872002012111(),01/21!mmnINwnnNIxIxm窗函数和设计参数的选择•滤波器长度估计N、C取决于窗函数的选择通带波纹和阻带波纹越匹配，窗长度越小该公式结果较为保守。•截止频率的选择•凯泽窗函数归一化过渡带滤波器阶次形状参数spCN12cpsps7.9612.286sAN0.40.11028.7500.5842210.7886212150021ssssssAAdBAAdBAdBAdB凯泽窗函数形状参数β与所设计FIR滤波器性能指标的关系β预计阻带衰减/dB实际阻带衰减/dB5.654566.063647.072728.081819.090901