FIR等波纹高通滤波器报告

课程设计（报告）题目高通等波纹FIR滤波器的设计与分析学院名称指导教师班 级学 号学生姓名2011年6月14日摘要：FIR数字滤波器的设计方法有窗函数法，频率采样法以及等波纹逼近法等，其中等波纹逼近法为最优化设计，在同样的技术指标下，用这种方法设计得到的滤波器要比窗函数法和频率采样法得到的滤波器的长度均要小，而且设计过程简单易行。等波纹逼近法使用Parks-McClellan算法，其原理是应用最大误差最小化准则，使设计的FIR数字滤波器的频域特性与理想滤波器的频域特性的最大误差为最小，并且它的幅度在通带和阻带内具有等波动性，误差在整个频带内均匀分布。在本课程设计中，借助MATLAB，设计出高通等波纹FIR滤波器，再分析不同字长对系统性能的影响及不同网络结构对其性能指标的影响。关键字：Parks-McClellan算法等波纹字长网络结构Abstract:FIRdigitalfilterdesignmethodsarewindowmethod,frequencysamplingmethod,rippleapproximationmethodandsoon,inwhichsuchrippleapproximationmethodistheoptimizeddesign.Withthesametechnicalindicators,thelengthofthefilteristhesmallestbyusingthismethod,andthedesignprocessissimple.RippleapproximationmethodisbasedontheParks-McClellanalgorithm,theprincipleisappliedtominimizeerrorcriteriontothemaximum,andmakethemaximumerrorofthefrequencycharacteristicsofFIRdigitalfiltercomparedwithidealfilter,andmaketheamplitude'srippleapproximateinit'spassbandandstopband,andtheerrorisevenlydistributedthroughouttheband.Inthispaper,byUsingMATLAB,first,workingoutrippleapproximationFIRhigh-passfilter;second,analyzingsystemperformanceinfluencedbydifferentword-lengtheffectanddifferentnetworkstructures.Keywords:Parks-McClellanalgorithmRippleapproximationWord-lengthNetworkstructures目录TOC\o"1-5"\h\z目录 3\o"CurrentDocument"一、 FIR滤波器的基本认识 4\o"CurrentDocument"二、 线性相位FIR数字滤波器及其特点 4\o"CurrentDocument"线性相位FIRDF 41.2线性相位FIR数字滤波器的零点分布特点 5\o"CurrentDocument"三、 等波纹最佳逼近法概述 5四、 Remez交换算法 6\o"CurrentDocument"五、 基于MATLAB的数字高通FIR最佳等波纹的实现 7\o"CurrentDocument"remezord和remez函数介绍 75.2高通滤波器设计指标 75.3MATLAB编程（程序见附录一FIR高通等波纹滤波器设计程序 7）设计滤波器 7\o"CurrentDocument"六、 使用FDAtool对滤波器进行设计与分析 9\o"CurrentDocument"七、 不同量化字长对系统的影响 10\o"CurrentDocument"八、 不同网络结构对滤波器性能指标的影响 148.1直接型结构 148.2线性相位结构 158.3级联型结构 168.4编写量化误差对直接型和级联型结构系统性能指标的影响程序程序 178.4.1直接型 17\o"CurrentDocument"8.4.2级联型 18\o"CurrentDocument"九、 课程设计心得体会 20\o"CurrentDocument"参考文献 21附录MATLAB程序清单 21附录一FIR高通等波纹滤波器设计程序 21附录二量化误差对直接型和级联型结构系统性能指标的影响程序 22附录三量化子函数quant2 23弓I言：数字滤波器(digitalfilter)是由数字乘法器、加法器和延时单元组成的一种装置，在通信、图像、语音、雷达等许多领域都有着十分广泛的应用。在数字信号处理中，数字滤波占有极其重要的地位。目前对数字滤波器的设计有多种方法。其中Matlab软件已成为设计数字滤波器的强有力工具。利用Matlab信号处理工具箱可以快速有效地实现由软件组成的常规数字滤波器的设计、分析和仿真，极大地减轻了工作量，有利于滤波器设计的最优化。一、 FIR滤波器的基本认识滤波器的单位脉冲响应h(n)长度为N,那么H(z)=£^h(n)Z-nn=0H(z)是z-i的N-l次多项式(FIR滤波器的阶数是N-1)该系统在z平面上有N-1个零点，在z=0处有一个N-1重极点。系统永远稳定。与IIR滤波器比较，设计方法差别很大，由于h(n)长度为有限的，因此FIR滤波器设计通常是找出这有限个值二、 线性相位FIR数字滤波器及其特点1.1线性相位FIRDF设FIRDF的单位脉冲响应h(n)长度为N,则其频率响应函数为HS)=H(3)eje(w)g式中，H(3)称为幅度特性函数，是w的实函数(可以取负值)，与|H(ej3)gI不同；e(3)称为相位特性函数。线性相位 0㈣=-®T1-11 为常数如果e(3)满足下式：0@)0@)=0-TW，0001-12 是起始相位严格地说，此时e(3)不具有线性相位特性，但以上两种情况都满足群时延d0（①）_是一个常数，即—页一一T也称这种情况为线性相位。一般称满足式1-11是第一类线性相位；满足式1-12为第二类线性相位。表1线性相位条件对FIRDF时域约束呦）对称性W)长度N幅度函数馬（2）相位特性可以设计的滤波器类型第一类线性相位偶对称奇数3=0*兀*2兀三点偶对称N-1-€0 2四种滤波器都可设计第一类线性相位偶对称偶数关于3=71奇对称，关于3=0、1F偶对称N-1-co 2不能实现髙通和带阻滤波特性「第二类线性相位奇对称奇数3=0*兀*2兀三点奇对称N-1一群/2_& 2只能设计带通滤波器第二类线性相位奇对称偶数关于°=71偶对称，关于3=0*1f奇对称N-1-^：1-3 1可设计高通、带通滤波器，不能设计低通和带阻1・2线性相位FIR数字滤波器的零点分布特点零点倒数成对出现；零点共轭成对出现。三、等波纹最佳逼近法概述等波纹最佳逼近法是一种优化设计法，它克服了窗函数设计法和频率采样法的缺点，使最大误差（即波纹的峰值）最小化，并在整个逼近频段上均匀分布。用等波纹最佳逼近法设计的FIR数字滤波器的幅频响应在通带和阻带都是等波纹的，而且可以分别控制通带和阻带波纹幅度。这就是等波纹的含义。最佳逼近是指在滤波器长度给定的条件下，使加权误差波纹幅度最小化。与窗函数设计法和频率采样法比较，由于这种设计法使滤波器的最大逼近误差均匀分布，所以设计的滤波器性能价格比最高。阶数相同时，这种设计法使滤波器的最大逼近误差最小，即通带最大衰减最小，阻带最小衰减最大；指标相同时，这种设计法使滤波器阶数最低。实现FIR数字滤波器的等波纹最佳逼近法的MATLAB信号处理工具函数为remez和remezord。Remez函数采用数值分析中的remez多重交换迭代算法求解等波纹最佳逼近问题，求的满足等波纹最佳逼近准则的FIR数字滤波器的单位脉冲响应h(n)。由于切比雪夫和雷米兹对解决该问题做出了贡献，所以又称之为切比雪夫逼近法和雷米兹逼近法。Remez交换算法是利用交替定理来求解的一种最优化算法，这种方法的框图如下图所示，分三步求逼近问题的解。(1)对于给定的一组频率(3，0W3W3或3W3W刀｝k=k k p p k0,1,…,r,求解满足下式的。值，即W(ejWk)[H(ejWk)—P(ejWk)]=(―1)k5 k=0,1,,rd利用已求出的C和给定的r个频率点，求解P(ej3)。由于已求出P(ej3)，则可计算E(ejw)=W>(ejw)[HH(ej)—P(ejw)]d上述算法的结果仅得到P(ej3)，还要经过傅里叶反变换，才能得到冲激响应h(n)。因此，根据交替定理，最优化滤波器设计程序有以下几个步骤：输入部分，包括滤波器的性能要求及滤波器类型。前者指的是所需的频率响应的幅度函数为D(w)，加权函数W(w)和滤波器单位抽样响应长度N,后者是要指出所需设计的是带通(包括低通、带通、高通、带阻等)滤波器或是微分器或是希尔伯特变换器。用公式表示逼近问题，也就是表示加权逼近误差E(e).用Remez多次交换算法，求逼近问题的解。计算滤波器的单位抽样响应。五、基于MATLAB的数字高通FIR最佳等波纹的实现5.1remezord和remez函数介绍[M,fo,m0,w]=remezord[f,m,rip]功能：利用remezord函数可以通过估算得到滤波器的近似阶数M,归一化频率带边界fo,频带内幅值m0及各个频带内的加权系w。输入参数f为频带边缘频率，m为各个频带所期望的幅度值，rip是各个频带允许的最大波动。h=remez[M,fo,m0,w,‘ftype']功能：利用remez函数可以得到最优化设计的数字FIR滤波器的h(n)系数，输入参数M是滤波器的阶数。ftype是所设计的滤波器类型，它除了可以设计普通的滤波器外，它还可以设计数字希尔伯特变换器以及数字微分器。实际设计中，由于remezord函数可跑高估或低估滤波器的阶数n,因此在得到滤波器的系数后，必须检查其阻带最小衰减是否满足设计要求。如果此时的技术指标不能满足设计要求，则必须提高滤波器的阶数到n+1,n+2等。5・2高通滤波器设计指标逼近通带[wp,n],通带最大衰减：apdB;逼近阻带：[0,ws],阻带最小衰减：asdB。在确定设计指标时，出于以下考虑：阶数越高，逼近精度越高；为增强可分析性，先用remezord函数估算出滤波器阶数M使阻带最小衰减小于要求，再根据线性相位特性情况增加阶数；为了后续的系统结构分析方便，阶数不能设置过高。最终选定如下指标：通带截止频率=0.75 通带最大衰减ap=1dB阻带截止频率二0.4 阻带最小衰减as=40dB5・3MATLAB编程(程序见附录一FIR高通等波纹滤波器设计程序)设计滤波器

M=14,N=M+1=15时，.得出了如下的衰减函数，此时阻带最小衰减＜40dB,不满足要求，说明用remez估算的M略小，需要增加M。由hn曲线可以看出，hn关于(N+1)/2偶对称，该FIR滤波器属于第一类线性相位特性。因而要使N为偶数，M需增加2。实际脉冲响应h(n)_E-0.2■°%----------r~~~~~~~~~~1-_E-0.2■°%----------r~~~~~~~~~~1--_—-—--n1 11 11 11 11 11 1----------n ___q1 1 11 1 11 1 11 1 11 1 11 1 111111111111111■ 111111111111111111111 1 T 11 1 11 1 11 1 19 ? [ 「更:1 11 11 11 11 1111111111111111 1 1 1 11 1: : d1 11 1>< 1 T 11 1 11 1 1b : : :1 1 11 1 1101214令M=M+2=16此时观察衰减函数曲线可看出，阻带最小衰减＞40dB，满足了要求。另外给出其幅频相应与相频相应。实际脉冲响应h(n).2J

o.2J-Oo.4

a101214160.3 0.4 0.5 0.6NormalizedFrequency(xnrad/sample)0.70.8 0.9mp)apn-FBewmseLIdo-0DO0.30.4 0.5 0.6NormalizedFrequency(x“rad/sample)0.70.8 0.9-2D00000实际脉冲响应h(n).2J

o.2J-Oo.4

a101214160.3 0.4 0.5 0.6NormalizedFrequency(xnrad/sample)0.70.8 0.9mp)apn-FBewmseLIdo-0DO0.30.4 0.5 0.6NormalizedFrequency(x“rad/sample)0.70.8 0.9-2D00000-6六、使用FDAtool对滤波器进行设计与分析首先在命令窗口键人FDAtool命令，调出FDAtool界面，在FilterType选项中选择Highpass（高通）；在DesignMethod中选择FIR滤波器；接着在FIR中选择Equiripple（等波纹）法。然后在FilterOrder中选择阶数（使用最小阶数时,阻带最小衰减同样出现＜40dB的情况，需增加阶数）。在FrenquencySpecification（频率设定）选项中，Unit设为Normalized即归一化频率，wstop:0.4,wpass:0.75。再对MagnitudeSpecification（幅度设定）进行设定，使Apass（通带波纹）1dB；Astop（阻带衰减）40dB，最后点击下方的DesignFilter（滤波器设计）即可得到所设计的FIR滤波器。设计完成后，可以通过菜单选项Analysis来分析滤波器的幅频响应、相频响应特性、单位脉冲响应、零极点分布图以及舍入功率密度谱。20-12-44-76-108-140andPhase20-12-44-76-108-140andPhaseNormalizedFrequency(xtiradfeample)■wmofeoEoCTBLId2000Structure: Direct-FormFIROrder: 16Sections: 1Stable: YesSource: DesignedStoreFilter...FilterManager...七、不同量化字长对系统的影响数字计算机要用二进制编码信号进行运算，二进制编码的位数（或字长）有限，带来了各种量化误差，形成有限精度的运算。如果信号值用b+1位二进制数表示，其中b为表示小数部分，能表示的最小单位称为量阶，用q表示，q=2人（-b）。对于超过b位的部分进行尾数处理。舍入法将尾数第b+1位逢1进位，逢0不进位，b+1位以后的数略去。舍入法统计平均值为0,方差为qA2/12oMATLAB按二进制双精度格式表示数，表示一个数用8字节（64位二进制数），所以MATLAB的量化误差可以忽略不计。实际工程中要把用MATLAB设

计的滤波器付诸实现，须米用字长一般为8、16、32位的DSP芯片，因而用MATLAB完成滤波器设计后，必须考虑实际系统的有效字长。用FDAtool分析不同字长对系统性能的影响设置字长为8：幅频和相频特性零极点图tIKdA-I4SUD14SE-11T-I111111111 1 1| | |1 1 1| | |tIKdA-I4SUD14SE-11T-I111111111 1 1| | |1 1 1| | |i 丨⑪厂"「吨］ ■ :「口 0 : ;11.■: : f□■暇I 丨 丨::■□ : J F…v ； V 1 :埠 1 : ::何 乂 L J 「亠」_Ej Zj 1 1 .一誌 i j L■ 1■ 1> 1■ 1111111111111 11口 J 1 I-2265O5-1.5-40RealPart单位脉冲响应QJpntEn-QJpntEn-舍入噪声功率谱msgod-8400.20.3 0.4 0.50.6 0.70.8 0.9-8I:-82msgod-8400.20.3 0.4 0.50.6 0.70.8 0.9-8I:-82由上面的幅频相频特性和零极点图可以明显看出，字长为8时，量化误差较大,对系统性能的影响较大。设置字长为16幅频和相频特性00.10.2 0.30.4 0.50.6 0.7NormalizedFrequency(xnrad/sample)0.80.900.10.2 0.30.4 0.50.6 0.7NormalizedFrequency(xnrad/sample)0.80.9零极点图-3 -2 -1 0 12 3 4RealPart单位脉冲响应epn七-dE-=J4.3.2Jepn七-dE-=J4.3.2JO.11O.O..311-舍入噪声功率谱0.4 0.5 0.6NormalizedFrequen匚竽(xurad/sample)0.70.8 0.956-U设置字长为32幅频和相频特性零极点图12 3 4URealPart单位脉冲响应-S5七-一L1LE=]-0.4 0.5 0.6NormalizedFrequen匚竽(xurad/sample)0.70.8 0.956-U设置字长为32幅频和相频特性零极点图12 3 4URealPart单位脉冲响应-S5七-一L1LE=]-舍入噪声功率谱〔〔LDP)god小结：显然，字长越长，量化误差越小。这从幅频相频特性和零极点图可直观看出；单位脉冲响应不太明显；而从舍入噪声功率谱可定性地看出，字长越长，波动越剧烈，说明量化噪声越大。八、不同网络结构对滤波器性能指标的影响运算量化误差不仅和字长有关，还与网络结构有关。8.1直接型结构由前面运算出的hn二［-0.003206790257856,0.004922580543198470.0157426771800879,-0.0199854762313814,-0.0325727662487369,0.0715881229387441,0.0496124695097384,-0.304969178702688,0.443290020816201,-0.304969178702688,0.0496124695097384,0.0715881229387441,-0.0325727662487369,-0.0199854762313814,0.0157426771800879,0.00492258054319847,-0.003206790257856］，N=17。可得出直接型结构为：H(z)=^+1hnz-n=—0.003206790257856+0.00492258054z-1+0.0157426771800879z-2n=00.00492258054319847z-i5+0.003206790257856z-i6表示系统输入输出关系的差分方程可写作:y(n)=£1h(m)x(n-m)m=0直接由差分方程得出的实现结构如图所示:

x(n)-1-1-1x(n)-1-1-1h(D户h(1ph(驴h(N-2)^ h(N-1)11——=■_——=■ H_——HH——■血)图FIR直接型网络结构假设定点乘法运算按b位进行量化，量化误差用e(n)表示。系统中所有乘法支1路都要引入一个噪声源，每个噪声源的方差均为b2=q2，q=2-b。e12e二e+e+e+ e+ef0 1 2 15 16q2q2丄q2丄 十 十1212q2_17q2十12—~L2~8.2线性相位结构线性相位结构是FIR系统直接型结构的简化网络结构，特点是网络具有线性相位特性，比直接型结构节约了一半的乘法器。若系统具有线性相位，应满足：h(n)=土h(N-n-1)〃+〃代表第一类线性相位滤波器；“-”代表第二类线性相位滤波器。N为偶数时，H(z)=肾-h(n)[z-n土z-(N-n-1)]n=0N_12-1 N-1N-1N为奇数时，H(z)=£h(n)[z-n土z-(N-n-1)]+h(—)z-22n=0N为奇数N为偶数由hn,得出线性相位结构为由于乘法器减小到了8个，因而量化误差也减小了将近一半。e二e+ee+eTOC\o"1-5"\h\zf0 1 6 78q2 2q2b2= = —f12 38.3级联型结构将H(z)进行因式分解，并将共轭成对的零点放在一起，形成一个系数为视数的二阶形式，这样级联型网络就是由一阶或二阶因子构成的级联结构，其中每一个因式都用直接型实现。借助MATLAB中的tf2sos将直接型转换为级联型结构，得出H(z)二乞2h(n)[-0.003206790257856+0.00492258054319847z-i+ +0.0496124695097384z-6-0.304969178702688z-7]n=0+0.443290020816201z-8H(z)=-0.003206790257856(1-0.197030972073081z-1-5.07657538969303z-2)(1+0.038811788843542(1-1.94432422612996z-1+0.999999999999994z-2)(1+0.77334303393771z-1+0.248584289912379z-2)因生成的S矩阵共8行，则H(z)=H(z)H(z)……H(z)H(z)TOC\o"1-5"\h\z1 2 7 8e (n)=[e (n)+e (n)+e (n)]*[h (n)*h (n) * h (n)*h (n)]f 0 1 2 2 3 7 8+[e(n)+e(n)+e(n)]*[h(n) *h(n)*h(n)] 3 4 5 3 7 8+[e(n)+e(n)+e(n)]*h(n)+[e(n)+e(n)+e(n)]18 19 20 8 21 22 23dzdzH(z-1)H(z-1)]z7 8z2 1 fG2= q2x f[H(z)H(z)f12 2丐 2 3cH7(z)H8(z川H2(z-1)H3(z-1)2 1fq2 1fq2xf12 2njc2q212"zH(z)H(z-1)上+8 8z小结：量化误差直接型网络结构最大；线性相位结构的量化误差约为直接型的一半；级联型的比直接型的小。8.4编写量化误差对直接型和级联型结构系统性能指标的影响程序程序（见附录二量化误差对直接型和级联型结构系统性能指标的影响程序,而其中要调用编写的量化子函数quant2（见附录三量化子函数quant2））做出8位量化对直接型和级联型结构零极点和幅频特性的影响。为了明显地观察量化效应对各种结构的影响，都进行8位量化。该程序的设计思想是，分别对每种网络结构的系数进行量化：对于直接型，直接对hn进行量化。对于级联型，要转换成级联型，再对S矩阵量化。因为需要作出量化前后衰减函数及相频特性曲线，而使用直接型系数即hn作图最方便，因而量化后在转为直接型。由于MATLAB用双精度格式表示数，量化误差比起8位量化可忽略不计，所以可以采用上述方法。8・4・1直接型量化前后零极点图量化前零极点图量化零极点图-2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5RealPart-2 -1.5 -1 -0.5 0 0.5 1 1.5 2RealPart量化前后衰减函数量化前后相频特性8.4.2级联型量化前后零极点量化前零极点图量化零极点图411111111-1—111411111111-1—1113111-1—1113111-1—1112111-1_112111_1_11o io i11■ 严弋沁 -/ : Q:°：16◎0■tj T t/- -tJ -； E□ T - -^3 T 七：i\O: ?\:/-1- 、一存 __-1. 、---J一费 _-2O :1111_ 1 _1 1 11111-2o :1111_ 1 _1 1 —1111-3-4111-1-1111111-1-1111111-3-4111-1-1111111-1-1111111-2-1012-2-1012RealPart RealPart量化前后衰减函数量化前后相频特性600小结：比较两种结构量化前后的性能指标，可以明显地观察出，级联型量化误差比直接型小。九、课程设计心得体会在本学期的理论课上，主要进行了前三章数字信号处理基础知识的学习。在滤波器的设计与分析时，还需要掌握第五章网络结构的知识；第六章数字滤波器的基本知识；第七章FIR滤波器的知识、等波纹FIR滤波器的设计方法、FDAtool的使用简介；第九章字长效应和不同网络结构的量化效应。因而，对滤波器设计分析前投入了很大的时间精力看教材，打好理论基础。只有对基本知识有了系统的把握，在设计与分析滤波器时才能理清思路。当然，仅靠课本是远远不够的。比如，等波纹滤波器实现的理论基础课本仅做了基本思想的介绍，还需要到网上查阅资料。更重要的是，由于需要借助MATLAB进行滤波器的设计与分析，许多思路是以如何能用MATLAB程序实现为出发点的。以本次课程设计为契机，巩固和提高了MATLAB的使用技能，尤其是大量的特性分析图的运用使本人熟练掌握了子图、双变量及各种特性图的作图方法。本课程设计编写许多程序，程序的编写会经常出现问题，有时是因为符号输入错误等低级失误，而大部分问题还是由于对所调用函数的理解不到位。其中量化子函数quant原先是参考课本的，但对级联型网络结构S矩阵量化是总是出错，于是修改为quant2，使用中再没出现问题。通过对滤波器不同网络结构的分析，进一步巩固了系统函数及系统性能的相关问题，也体会到了量化效应对系统性能指标的影响，这些都是对本学期理论知识的巩固与提升。本次课程设计获益匪浅，感谢陈老师及对本人提出过帮助的同学！参考文献高西全丁玉美•数字信号处理（第三版）•西安电子科技大学出版社JohnG.ProakisDimitrisG.Manolakis.（张晓林肖创柏译）.数字信号处理原理、算法与应用.电子工业出版社附录MATLAB程序清单附录一FIR高通等波纹滤波器设计程序%FIR高通等波纹滤波器ap=l；as=40;%通带最大衰减ldB,阻带最小衰减40dBdl=(10人(ap/20)-l)/(10人(as/20)+l);%通带振荡波纹d2=10人(-as/20);%阻带振荡波纹f=[0.4,0.75];m=[0,l];%与f对应的幅度向量rip=[d2,dl];[M,fO,mO,W]=remezord(f,m,rip);%求出阶数、频率、幅值、误差加权向量M=M+2hn=remez(M,fO,mO,W);n=0:M;subplot(221);stem(n,hn);gridon;title('实际脉冲响应h(n)');xlabel('n');ylabel('