《课程设计(论文)基于巴特沃斯法数字IIR带阻滤波器的设计》由会员分享,可在线阅读,更多相关《课程设计(论文)基于巴特沃斯法数字IIR带阻滤波器的设计(18页珍藏版)》请在金锄头文库上搜索。
1、基于巴特沃斯法数字IIR带阻滤波器的设计 1、 数字带阻IIR滤波器设计IIR数字滤波器在很多领域中有着广阔的应用。与FIR数字滤波器相比,它可以用较低的阶数获得高选择性,所用存储单元少,经济而效率高,在相同门级规模和相同时钟速度下可以提供更好的带外衰减特性。下面介绍一种设计实现IIR数字滤波器的方法。设计步骤如下:步骤1:将数字滤波器H(z)的技术指标p和s,通过tan(/2)转变为模拟滤波器G(s)的技术指标p和s,作归一化处理后,得到p1,ss/p;步骤2:化解为模拟原型滤波器G(s)的技术指标;步骤3:设计模拟原型滤波器G(p);步骤4:将G(p)转换为模拟滤波器的转移函数G(s);步
2、骤5:将G(s)转换成数字滤波器的转移函数H(z),s=(z1)(z+1)。所谓原型滤波器是指归一化的低通滤波器。本节主要讨论通过IIR数字滤波器的原型转换设计法和IIR数字滤波器的直接设计方法来设计数字高通、带通及带阻滤波器其转换方法主要有3种:一是直接由模拟低通滤波器转换成数字高通、带通或带阻滤波器;二是先由模拟低通滤波器转换成模拟高通、带通或带阻滤波器,然后再把它转换成相应的数字滤波器;三是将模拟低通滤波器先转换成数字低通滤波器,再通过变量代换变换成高通、带通或带阻滤波器。数字IIR滤波器的设计结构如图:模拟低通原型频率变换模拟离散化IIR数字滤波器图1、IIR滤波器的设计步骤2、带阻滤
3、波器的设计规范: 本设计中分别用h(t)、H(s)、Ha(j)表示模拟滤波器的单位脉冲相应、系统函数、频率响应函数,三者的关系如下:Ha(s) = LTha(t) =edtHa(j)=FTha(t) =edt可以用h(t)、H(s)、Ha(j)中任一个描述模拟滤波器,也可以用线性常系数微分方程描述模拟滤波器。但是设计模拟滤波器时,设计指标一般由幅频相应函数|Ha(j)|给出,而模拟滤波器设计就是根据设计指标,求系统函数H(s)。工程实际中通常用所谓的损耗函数即衰减函数A()来描述滤波器的幅频响应特性,对归一化幅频响应函数,A()定义如下:A()=-20lg|Ha(j)|=-10lg|Ha(j)
4、| dB模拟带阻滤波器的设计指标参数有、。和分别表示带阻滤波器的通带下边界频率和通带上边界频率;和分别表示带阻滤波器的阻带下边界频率和阻带上边界频率。1和2分别称为通带和阻带波纹幅度 图2、带阻滤波器带阻滤波器的技术指标要求必须是几何对称的。如果带阻滤波器不是几何对称的,为了使对称的带阻滤波器满足最小阻带衰减要求,在由非对称的带阻滤波器变换为对称的带阻滤波器时,保留非对称带阻滤波器的阻带截止频率调整两个通带截止频率中的一个,调整过程步骤如下:1、 计算=2、 计算=,如果,用代替;3、 如果,计算=,并用代替;4、 如果A A,选择A=min A, A。所有类型的变换,包括各自与低通原型的频率
5、变换如表1所示。用频率变换法设计一个模拟滤波器的步骤如下:1、确定低通、高通、带通和带阻模拟滤波器的技术要求。当要求的带通或带阻滤波器是非几何对称时,根据上面描述的过程使这些技术要求成几何对称的;2、根据表1中的变换关系,确定归一化低通滤波器的技术要求:通带截止频率为,阻带截止频率为,阻带衰减为A(dB),阻带衰减为A (dB);3、根据、A 和A,用巴特沃斯设计归一化低通滤波器;4、根据表1最右栏的变换关系,得到要求的非归一化模拟滤波器。 表1 模拟滤波器的频率变换滤波器类型归一化低通滤波器H(p)的技术指标要求要求的滤波器H(s), H(s)= H(p)低通H(p)低通H(s)P=低通H(
6、p)高通H(s)P=低通H(p)带通H(s)P=低通H(p)带阻H(s)P=下面我们用巴特沃斯滤波法来设计一个数字IIR带阻型滤波器:滤波器规范:通带下边界频率: 50HZ 通带上边界频率: 450HZ阻带下边界频率: 200HZ阻带上边界频率: 300HZ通带最大衰减: 3dB阻带最小衰减: 20dB 抽样频率: 1KHZ3、带阻滤波器的系数计算: 首先确定所要求的带阻滤波器是否为几何对称的:=2*50*2*450=4*2.25*10=2*200*2*300=4*6*10因为,而且设计的是带阻滤波器,所以需要调整这个带阻滤波器的通带上下边界频率。因为: =,=,则= 用值代替值,即令=2*1
7、33.3HZ。因此,所要求的带阻滤波器的指标调整为:=2*133.3/1000 rad/s=2*450/1000 rad/s=2*200 /1000rad/s=2*300 /1000rad/sB=2*266.7 /1000rad/s=2*245 /1000rad/sA= A=3 dBA= A=20 dB根据表1中间栏的变换关系式,将上述给定的带阻滤波器指标要求转化为相应的归一化低通技术要求,有 =1=2.6根据上面的技术要求,可以采用查表法或计算法来设计归一化的巴特沃斯低通滤波器,这里采用计算法。 由式 确定最小的的滤波器阶数,得N=1.9取N=2,计算巴特沃斯滤波器的归一化极点 p=ee k
8、=0,1 p=e=-0.5000+0.8660j p= e=-1.0000+0.0000j因此,低通滤波器的归一化传输函数H(p)为 H(p)= = 根据表1最右边栏的变换关系式,把归一化低通滤波器变成所要求的带阻滤波器的传输函数H(s) H(s) =H(p) =由双线性变换法求得数字巴特沃斯带阻滤波器的传递函数为: H(Z)= H(s)=软件设计部分Fs=1000;Ap=3;As=20;Wp=50/500,450/500;Ws=200/500,300/500;b,a=butter(N,Ws,stop)系统指标:传输函数系数A:系数B:4、带阻滤波器的实现结构: 令H(Z)=H(Z)H(Z):
9、 其中 H(Z)= H(Z)= 对应的定义滤波器的差分方程如下:滤波器部分1:w(n)=(1/s)x(n)-aw(n-2) y(n)= w(n)s/s+ w(n-2) s/s滤波器部分2:w(n)= y(n)- w(n-2)y(n)= w(n) s+ w(n-2) s系数和的精确值依赖于我们对H(Z)的多项式的分子和分母如何组对,以及实现多项式的二阶滤波器部分是如何排序的。最好的组对和排序只能通过有限字长分析来确定。实现的结构图如下:x(n) y(n) -a s/s s Z Z 图3、滤波器的实现结构图5、有效字长效应分析:由给定的性能规范,我们将假设使用的算法是定点的2的补码算术,每一个系数
10、被舍入量化到16位字长。这里我们主要的是评估不同量化误差对滤波器性能的影响,以及根据信噪比确定最好的滤波器结构用于实现,所要考虑的误差源主要是:1、 溢出误差2、 舍入误差3、 系数量化误差但由于DSP本身有限字长和精度的因素,加上IIR滤波器在结构上存在反馈回路,是递归型的,再者高阶滤波器参数的动态范围很大。这样一来造成两个后果:结果溢出和误差增大,从而导致算法无法在DSP上实现。解决此问题的有效方法是把高阶IIR数字滤波器简化成几个2阶滤波器来设计,即采用级联结构。因此系统函数为:H(z)=H1(z)H2(z)Hn(z)(3)上述4阶IIR数字滤波器采用级联的结构后变成2个滤波器级联的形式
11、,每个滤波器的系数表如下所示。量化前滤波器系数为了设计可实现的IIR滤波器:需要对上一步分解获得的二阶基本节的滤波器系数进行量化,即用一个固定的字长加以表示。量化过程中由于存在不同程度的量化误差,由此会导致滤波器的频率响应出现偏差,严重时会使IIR滤波器的极点移到单位圆之外,系统因而失去稳定性。为了获得最优的滤波器系数量化后滤波器的系数:除了系数存在量化误差,数字滤波器运算过程中有限字长效应也会造成误差,因此对滤波器中乘法器、加法器及寄存器的数据宽度要也进行合理的设计,以防止产生极限环现象和溢出振荡。量化过程中由于存在不同程度的量化误差,由此会导致滤波器的频率响应出现偏差,严重时会使IIR滤波
12、器的极点移到单位圆之外,系统因而失去稳定性图4、量化前后系统的零极点图 由图4可看出系统零极点都在单位圆内,系统稳定可实现。这两个可能的滤波器结构的每一个都有不同的比例因子,以及不同的信号舍入误差。这个步骤的目标是根据信噪比性能方面来确定最好的组对和排序。溢出和舍入误差是紧密相连的,所以伸缩变换和舍入分析应该同时进行。利用有限字长的分析程序,可以得到对于上面两个可能的滤波器的比例因子,它们是基于L 、L 和L的范数。在这个例子里,我们利用了L的范数,对于一个以两个标准部分串联实现的四阶滤波器,在伸缩变换以后,输出端的舍入噪声为:=3s| H(Z)H(Z)|+5s|H(z)|+3其中q是量化步长
13、或舍入,|.|是L范数的平方。H(Z)是第一级滤波器的传递函数,H(Z)是第二级滤波器的传递函数,s是第一个滤波器阶段的伸缩比例因子,s是第二个滤波器的比例因子。系统的量化误差的效应:我们检查那些给出的系数字长是否满足稳定和频率响应的规范。当极点不是非常靠近单位圆时,16位的系数字长对保持稳定是足够的。例如,对于第一个滤波器部分,实验程序表明2位就足以满足稳定要求,把系数量化成16位仅改变极点半径。采用16位系数字长时,滤波器的响应实际上和未量化的滤波器的响应是一样的。图5描绘了量化前后的滤波器的频率响应。图5、带阻滤波器量化前后的幅频响应-为量化前 量化后6、带阻滤波器的软件实现:由系统指标规范可编程:N为滤波器的阶数;zz和pz为系统的零点和极点;kz为系统的增益:Fs=1000;Ap=3;As=20;Wp=50/500,450/500;Ws=200/500,300/500;N,Wc=buttord(Wp,Ws,Ap,As);zz,pz,kz=butter(N,Ws,stop);b,a=butter(N,Ws,stop);Subplot(2,1,1)H,f=freqz(b,a,512,Fs);Plot(f,abs(H)Xlabel(频率/(HZ)Ylabel(幅度/(dB)Subplot(2,1,2)Zplane(b,a)