巴特沃斯高通数字滤波器

巴特沃斯高通数字滤波器 ------------------------------------------------------------------------------------------------ 巴特沃斯高通数字滤波器 数字信号处理课程 老 师 陈 忠 泽 老 师 学 院 电气学院 班 级 电子信息工程0 81班 学 号 20084470110 姓 名 何 依 阳 二0一一年五月 目录: 一、IIR数字高通滤波器的设计 1、数字滤波器的概述 2、数字滤波器的设计步骤 3、设计方法 4、IIR巴特沃斯数字高通滤波器的实例计算 二、软件仿真工具及实现环境简介 1、计算机辅助设计方法 2、 MATLAB直接设计IIR巴特沃斯数字高通滤波器 三、滤波器结构对数字滤波器性能指标的影响 1、 IIR系统的基本网络结构 (1)直接型 (2)级联型 —————————————————————————————————————— ------------------------------------------------------------------------------------------------ (3)并联型 四、有限字长运算在网络结构中对数字滤波器的影响 1 、 运算量化效应对数字滤波器的影响 2 、 参数的字长对数字滤波器性能指标的影响 2.1 、系数量化对数字滤波器的影响 五 、运用MATLAB的辅助工具FDATOOL画出系统函数图像 六 、设计心得 IIR数字高通滤波器的设计 一、IIR数字高通滤波器的设计 1、数字滤波器的概述 所谓数字滤波器，是指输入、输出均为数字信号，通过数值运算处理改变输入信号所含频率成分的相对比例，或者滤除某些频率成分的数字器件或程序。 2、 数字滤波器的设计步骤 设计一个IIR数字滤波器主要包括下面5个步骤: (1) 确定滤波器要求的规范指标。 (2) 选择合适的滤波器系数的计算(如图一图所示)。 (3) 用一个适当的结构来表示滤波器(实现结构)。 (4) 有限字长效应对滤波器性能的影响分析。 (5) 用软件或硬件来实现滤波器。 本次设计的IIR数字滤波器系数的计算是根据已知的模拟滤波器的特性转换到等价的数字滤波器。本次设计用双线性变换法得到数字—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 滤波器。而且，双线性变换法得到的数字滤波器保留了模拟滤波器的幅度响应特性。 3 、设计方法 频率变换法设计思想: 1 、从归一化模拟低通原型出发，先在模拟域内经频率变换成为所需类型的模拟滤波器;然后进行双线性变换，由 S 域变换到 Z 域，而得到所需类型的数字滤波器。 2 、先进行双线性变换，将模拟低通原型滤波器变换成数字低通滤波器;然 后在 Z 域内经数字频率变换为所需类型的数字滤波器。 4 、IIR巴特沃斯数字高通滤波器的实例计算 (1)设计要求 设计一个巴特沃斯数字高通滤波器,要求通带截止频率 减不大于1dB,阻带截止频率?p?0.32? rad ,通带衰?p?0.16?rad，阻带衰减不小于15dB (2)确定数字高通的技术指标: 通带截止频率 ?p?0.32? ap?1 rad， 通带最小衰减 阻带截止频率 dB rad， ?p?0.16? 阻带最大衰减 as?15 dB ?ph?tan ?sh1?p?0.5494rad/s,?p?1dB2; 1?tan?s?0.2566rad/s,?s?15dB2。 (4)模拟低通滤波器的技术指标计算如下:对通带边界频率归—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 一化，由于本设计的为1dB截止频率，所以 ?p??c?1,?p?1dB 把?p?1?s?sh带入和-???p?ph?sh求得归一化巴特沃斯低通滤波器的阻带截止频率为: ?s???p?ph?sh?2.141 ， ?s?15 dB (5)设计归一化巴特沃斯模拟滤波器G(p)。 ksp???10.8757: ?sp? N??s?2.141?p: ?3.1349 : lgksplg?sp 所以取N=4，根据巴特沃斯归一化低通滤波器参数表(见附录)可得归一化模拟低通原型系统函数G(p)为: G(p)?1 p4?2.613p3?3.4142p2?2.6131p?1 1?22(p?0.7654p?1)(p?1.8478p?1) (6)利用频率变换公式p??p?phs ?p?ph s将G(p)转换成模拟高通HHP(s): Ha(s)?G(p)|p? 把 S4?4?ph+2.613?3phS+3.4142?2phS2+2.6131?phS3+S4?ph?0.5494代入此式可得: 4 s?2340.0911?0.4333s?1.0305s?1.4356s?s 用双线性变换法将模拟高通Ha(s)转换成数字高通H(z): —————————————————————————————————————— ------------------------------------------------------------------------------------------------ 令 21?S?T1?z?1?1 , 0.3647?1.4587z?1?2.1881z?2?1.4587z?3?0.3647z?4 H(z)?Ha(s)|1?z?1??1?2?3?4s?2?11?2.0578z?1.8545z?0.7895z?0.1331z1 ?z 二 、 软件仿真工具及实现环境简介 1、 计算机辅助设计方法 在优秀科技应用软件MATLAB的信号处理工具箱中提供了一整套模拟，数字滤波器的设计命令和运算函数，方便准确，简单容行使得设计人员除了可按上述传统设计步骤快速的进行较复杂高阶选频滤波器的计算、分析外，还可通过原型变换直接进行各种典型数字滤波器设计，即应用MATLAB设计工具从模拟原型直接变换成满足原定频域指标要求有数字滤波器。 2、 MATLAB直接设计IIR巴特沃斯数字高通滤波器 MATLAB编程如下: fs=5000; wp=800*2/fs; ws=400*2/fs; rp=1; rs=15; Nn=128; —————————————————————————————————————— ------------------------------------------------------------------------------------------------ [N,wn]=buttord(wp,ws,rp,rs) [b,a]=butter(N,wn,'high') freqz(b,a,Nn,fs) MATLAB运算结果如下: N = 4 wn = 0.2388 b = 0.3647 -1.4587 2.1881 -1.4587 0.3647 a = 1.0000 -2.0578 1.8545 -0.7895 0.1331 三、 滤波器结构对数字滤波器性能指标的影响分析 1 、 IIR系统的基本网络结构 IIR系统的基本网络结构有三种，即直接型、级联型和并联型。 (1)直接型 N阶差分方程如下: MN y(n)??bix(n?i)??aiy(n?i) i?0i?1 对应的系统函数为 H(z)??ibz?iM 1??aiz?i i?1i?0N 设M=N=2，按照差分方程可以直接画出网络结构如下图(a)所示。图中第一部分系统函数用H1(z)表示，第二部分用H2(z)表示，那—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 么H(z)?H1(z)H2(z)，当然也可以写成H(z)?H2(z)H1(z),按照该式，相当 于将下图(a)中两部分流图交换位置，如下图(b)所示。该图中节点 变量w1=w2，因此前后两部分的延时支路可以合并，形成如下图 (c) 所示的网络结构流图，我们将下图 (c)所示的这类流图称为IIR直接型 网络结构。 M=N=2时的系统函数为 b0?b1z?1?b2z?2 H(z)?1?a1z?1?a2z?2 接按照H(z)画出直接型结构流图。 对照下图 (c)的各支路的增益系数与H(z)分母分子多项式的系数 可见，可以直 IIR网络直接型结构 由bz和az写出数字滤波器系统函数: 0.3647?1.4587z?2.1881z?1.4587z?0.3647zH(z)?1?2.0578z?1?1.8545z?2? 0.7895z?3?0.1331z?4 有H(z)写出差分方程如下: ?1?2?3?4 y(n)?2.0578y(n?1)?1.8545y(n?2)?0.7895y(n?3)?0.1331y(n?4) ?0.3647x(n)?1.4587x(n?1)?2.1881x(n?2)?1.4587x(n?3)?0.3647x(n?4) 直接型网络结构如下图: 直接型网络结构图 在后面的分析中我们将发现，直接型系统对滤波器的性能控制作 —————————————————————————————————————— ------------------------------------------------------------------------------------------------ 用不明显，极点对系数的变化不灵敏，易出现不稳定或较大误差，而且运算的累积误差较大。因此，在设计时一般不选用直接型。 (2)级联型 级联型结构是将系统传递函数H(z)写成具有实系数的二阶节的乘积。将分子和分母多项式分解为各自的根，然后将一对复数共轭根(或者任意两个实数根)组合成二阶多项式。 在直接型表示的系统函数H(z)中，分子、分母均为多项式，且多项式的系 数一般为实数。现将分子、分母多项式分别进行因式分解，得到: H(z)?A?(1?C ?(1?d r?1r?1NMrz?1)z?1) r 上式中, A是常数; Cr和dr分别表示H(z)的零点和极点。由于多项式的系数是实数，Cr和dr是实数或者是共轭成对的复数，将共轭成对的零点(极点)放在一起，形成一个二阶多项式，其系数仍为实数;再将分子、分母均为实系数的二阶多项式放在一起，形成一个二阶网络Hj(z)。Hj(z)如下式: Hj(z)? 上式中Hj(z)?0j??1jz?1??2jz?21?a1jz?1?a2jz?2 Hj(z)表示一个一阶或二阶的数字网络的子系统函数，每个的网络结构均采用前面介绍的直接型网络结构，如下图所示，H(z)则由k个子系统级联构成。 一阶和二阶直接型网络结构 —————————————————————————————————————— ------------------------------------------------------------------------------------------------ 级联型的MATLAB的表示与实现: 将数字滤波器系统函数H(z): 0.3647?1.4587z?1?2.1881z?2?1.4587z?3?0.3647z?4 H(z)? ?1?2?3?41?2.0578z?1.8545z?0.7895z?0.1331z 用MATLAB转换成级联型的程序如下: bz=[0.3647 -1.4587 2.1881 -1.4587 0.3647]; az=[1.0000 -2.0578 1.8545 -0.7895 0.1331]; [S,G]=tf2sos(bz,az) freqz(bz,az) 说明:[S,G]=tf2sos(bz,az):实现直接型到级联型的变换。B和A分 别为直接型系统函数的分子和分母多项式系数向量。返回L级二阶级 联型结构的系数矩阵S和增益常数G。 MATLAB运算结果如下: S = 1.0000 -2.1807 1.1988 1.0000 -0.8976 0.2272 1.0000 -1.8190 0.8341 1.0000 -1.1602 0.5859 G = 0.3647 由S和G写出数字滤波器级联型系统函数: 1?2.1807z?1?1.1988z?21?1.819z?1?0.8341z?2 H(z)?0.3647()()?1?2?1?2 1?0.8976z?0.2272z1?1.1602z?0.5859z 级联型网络结构如下图: 在级联型结构中，每一级分子的系数确定一对零点，分母的系数 确定一对极点，因为子网络的零极点也即整体网络的零极点，所以整 —————————————————————————————————————— ------------------------------------------------------------------------------------------------ 个系统的零极点都可以准确的由每一级的系数来调整和控制，这样便于调整滤波器的频率响应性能，其灵敏度特性优于直接型和正准型结构。其次，级联结构具有最少的存储器。并联级联型网络结构图 支路的极点也是整个网络的极点，而并联支路的零点却不是整个网络的零点，因此并联网络能独立的调整系统的极点位置，但不能控制零点。并联结构的灵敏度由于直接型和正准型，运算累积误差比级联型小。 (3)并联型 如果将级联形式的H(z)展成部分分式形式，则得到: H(z)?H1(z)?H2(z)???Hk(z) 对应的网络结构为这k个子系统并联。上式中，Hi(z)通常为一阶网络或二阶网络，网络系统均为实数。二阶网络的系统函数一般为 ?0i??1iz?1 Hi(z)?1??1iz?1??2iz?2 式中，?0i、?1i、?1i和?2ii都是实数。如果?1i=?2i=0，则构成一阶网络。由上式，其输出Y(z)表示为 Y(z)?H1(z)X(z)?H2(z)X(z)???Hk(z)X(z) 上式表明将x(n)送入每个二阶(包括一阶)网络后，将所有输出加起来得到输出y(n)。 在并联型结构中，每一个一阶网络决定一个实数极点，每一个二阶网络决定一对共轭极点，因此调整极点位置方便，但调整零点位置不如级联型方便。另外，各个基本网络是并联的，产生的运算误差互—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 不影响，不像直接型和级联型那样有误差积累，因此，并联形式运算误差最小。由于基本网络并联，可同时对输入信号进行运算，因此并联型结构与直接型和级联型比较，其运算速度最高。 四、有限字长运算在网络结构中对数字滤波器的影响 1 、 运算量化效应对数字滤波器的影响 在实现数字滤波器时，将遇到相乘与求和运算。在定点制运算中，每一次乘法运算之后都要作一次舍入(或截留)处理，研究定点实现相乘运算的流图如下图所示。图(a)表示无限精度乘积y(n);图(b)表示有限精度乘积y(n),[?]表示? 舍入运算。采用统计分析方法时，可以将舍入误差作为独立噪声e(n)叠加在信号上，如图(c)所示。 x(n)y(n)x(n)y(n)?(n)yx (n)?(n)y a理想相乘 b实际相乘的非线性流图 c统计模型的线性流图 定点相乘运算的流图表示 显然，采用统计分析方法后，实际的输出可以表示为: 对于舍入处理，e(n)的均值为零，方差为: 2?e?(n),y(n)?e(n)yq2?12 现在以一个一阶IIR巴特沃斯数字高通滤波器为例来讨论分析方法。表示其输入与输出关系的差分方程为 y(n)?ay(n?1)?x(n)， n?0 式中|a|<1。它含有乘积项a?y(n?1) —————————————————————————————————————— ------------------------------------------------------------------------------------------------ ，这将引入一个舍入噪声，其统计分析流图示于下图。 一阶IIR滤波器的幅频特性 整个系统可以当作线性系统来处理。输出噪声ef(n)是由噪声源e(n)造成的输出误差，可由量化误差通过线性系统的方法求得输出噪声。由于e(n)叠加在输入端，因此 式中h(n)=anu(n)是一阶系统的单位脉冲响应,由式 ef(n)?e(n)?h(n)?e(n)?anu(n) 2 f????h(m)h(l)?(m?l)???2f2em?0l???2em?0?h(m)2? (?2f为输出噪声方差)和 1?1dz????H(z)H(z)c2?jz 2e 可求得输出噪声的方差 22m?2??afe? m?0? 式中H(z)为一阶IIR系统的传递函数，即 zH(z)?z?a 有此可以求得 1q22?2b ??????b222q?21?a12(1?a)3(1?a) ，(,) 2f2e 由此可见，字长b越大，数字滤波器输出端的噪声越小。 2 、参数的字长对数字滤波器性能指标的影响 2.1 、系数量化对数字滤波器的影响 系统对输入信号进行处理时需要若干参数或者称为系数，这些系数都要存储在有限位数的寄存器中，因此存在系数的量化效应。系数—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 的量化效应误差直接影响系统函数的零、极点位置，如果发生了偏移，会使系统的频率响应偏移理论设计的频率响应，不满足实际需要。 下面分析系数量化误差对极零点位置的影响。如果极零点位置改变了，严重时不仅IIR系统的频率响应会发生变化，还会影响系统的稳定性。因此研究极点位置的改变更加重要。为了表示系数量化对极点位置的影响，引入极点位置灵敏度的概念，所谓极点灵敏度, 是指每个极点对系数偏差的敏感程度。相应的还有零点位置灵敏度，分析方法相同。下面讨论系数量化对极点位置的影响。 现分析一个N阶直接型结构的IIR滤波器的传递函数 H(z)??kbz?kM 上式表示了一个N阶直接型结构的IIR数字滤波器的系统函数，该滤波器的极点k?11??akz?kk?0NB(z)?A(z)都在单位圆内聚集在z=1附近。系数ak和bk必须用有限位二进制数进行量化，存储在有限长的寄存器中，经过量化后的系数误差为?ar和?br，量化后的系数用ar和br，即 ?? ?r?ar??ar a ??b??bbrrr 则实际的系统函数可表示为: ?(z)?H?z?kb?k ?kz?k1??a k?1k?0NM 显然，系数量化后的频率响应已不同于原来设计的频率响应。用—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 直接型结构来实现该滤波器时，系数ak和bk都将直接出现在信号流程图中，其中ak影响着极点的位置。当由于系数量化误差使一个极点从单位圆内移动到单位圆上或单位圆外时，滤波器的稳定性即受到破坏。所以，只要有一个系数由于量化产生很微小的误差，就有可能使系统失去稳定。反馈支路的阶次N越高，使滤波器失去稳定的系数量化误差的绝对值就越小，则越容易使滤波器变得不稳定。 ?r?ar??ar和 设滤波器的传输函数H(z)，系数ak 和bk经舍入量化后为a ??b??bbrrr，这里?ar和?br是量化误差。 分母多项式有N个极点，用 的传输函数为: pi(i=1,2,… N)表示。这样，实际的滤波器 H(z)??z?bk k?0 NM?k? ?z?kb?kk?0?1[1?(pi??pi)z]?i?1NM?kz?k1??ak?1 上式中，?pi是第i个极点的偏差，称为极点误差，它应该和各个系数偏差都有关，它和各系数偏差的关系用下式表示: ?pi?pi??ak, k?1?akNi?1,?,N ?pi?pi a?a?p上式中，?ak的大小决定着系数k的偏差k对极点偏差i的影响程度。?ak ?pi —————————————————————————————————————— ------------------------------------------------------------------------------------------------ p?a?p越大，k对i影响也越大。称?ak为极点i对系数变化的灵敏度。 五 、运用MATLAB的辅助工具FDATOOL画出系统函数图像 系统函数系数向量经过1位二进制舍入量化前后幅频特性和相位响应 系统函数系数向量经过4位二进制舍入量化前后幅频特性和相位响应 系统函数系数向量经过1位二进制舍入量化前后零极点特性图 系统函数系数向量经过4位二进制舍入量化前后零极点特性图 由MATLAB运算1位和4位二进制舍入量化后结果可以看出，因为系数的量化，使极点位置发生变化，系数量化的相对误差(p-pq)不到10%，极点位置的相对误差ap达到了70%。问题不但是数量的变化，算一下极点的模，可以发现所有根的模都变大了，说明量化后的极点离单位圆稍远一些，如系统函数系数向量经过1和4位二进制舍入量化前后的零极点特性图如图(c)、(d)所示。这致使数字高通滤波器的幅频特性降低，运用MATLAB的辅助工具FDATOOL画出1位和4位量化墙后的幅频特性曲线分别如图(a)、(b)中的实线和虚线所示，这说明由于系数量化效应，使极点位置发生了变化，从而改变了原来设计的频率响应特性。 另外，我们还注意到4位二进制舍入量化后极点都在单位圆内部，但在系统函数系数向量经过1位二进制舍入量化后，原来较小的系数相对误差变化较大，使滤波器性能偏离原设计指标要求，使本来稳定的系统变成了不稳定滤波器。 —————————————————————————————————————— ------------------------------------------------------------------------------------------------ 从以上分析可以看出，系数的量化效应误差直接影响系统函数的零、极点位置，如果发生了偏移，会使系统的频率响应偏移理论设计的频率响应，不满足实际需要。因此，在设计滤波器的时候应尽量选择合理的系数量化，以减小量化误差，以便设计出较为理想的实用滤波器。 六 、设计心得: 在课设之前， 我对 MATLAB 软件，特别是滤波器设计中的函数基本上处于一种模糊状态。但是通过与同学不断的交流，最后完成了这次课设，对滤波器的设计有了比较清楚的了解。 这次课设的完成，真的不太容易。我的理论基础不是很扎实，所以完成起来要比其他同学费劲。在课程设计的过程中，我学到了很多东西，比如设计滤波器的一些基本函数的用法， 各种模拟滤波器的特性，设计滤波器的一些基本方法。但更为重要的是，我对于解决一个问题的思路更加清晰，找到了属于自己的方法。 当然，在设计的过程中，不可能避免的遇到了很多问题，如刚开始思路比较混乱，没有明确的方向。也由于一些基本的概念了解的不够清楚，比如模拟频带变换，我理解为是做模拟滤波器，而老师的课程设计题上有一句话似乎要求是要做数字滤波器。所以，在这个问题上，走了很多弯路。但通过与同学交流终于弄清楚数字滤波器而不是模拟滤波器。模拟频带变换只不过是一种变化方式或者说一种变化过程。是由模拟低通变为模拟高通，再由模拟高通变为数字高通的过程。 总的来说，这次课程设计让我对 MATLAB 有了更深刻的了解，对数字滤波器的—————————————————————————————————————— ------------------------------------------------------------------------------------------------ 设计流程有了大致的了解，掌握了一些设计滤波器的基本方法，提高 了理论用于实践的能力，掌握了更多专业相关的使用知识与技能。同 时，也暴露了我很多的不足，在以后的学习中， 将进一步发现并克 服缺点。 ——————————————————————————————————————