脉冲响应不变法设计IIR数字滤波器

1、第一章 摘要 本设计采用脉冲响应不变法设计 IIR数字低通滤波器，所以在根据脉冲响应不变法设计 出来的IIR数字滤波器会很好地重现原模拟滤波器的频率特性；数字滤波器的单位脉冲响 应完全模仿模拟滤波器的单位脉冲响应，时域特性逼近好，但容易产生频谱混叠现象， 只适合低通和带通滤波器的设计，不适合高通和带阻滤波器的设计。 关键词：数字滤波器；脉冲响应不变法；频率混叠 第二章 引言 数字滤波器可以满足滤波器对幅度和相位特性的严格要求，精确度高，有高度的可编 程性，灵活性好，并且它没有电压飘移、温度漂移及噪声等问题，基本不受环境影响， 稳定性好等。正是由于数字滤波器的以上优点，使得数字滤波器广泛应用于语

2、音处理、 图像处理、模式识别、频谱分析、医学仪器等领域。 第三章 设计原理 3.1 数字滤波器 数字滤波器是一种用来过滤时间离散信号的数字系统，通过对抽样数据进行数学处理 来达到频域滤波的目的。可以设计系统的频率响应，让它满足一定的要求，从而对通过 该系统的信号的某些特定的频率成分进行过滤，这就是滤波器的基本原理。如果系统是 一个连续系统，则滤波器称为模拟滤波器。如果系统是一个离散系统，则滤波器称为数 字滤波器。 数字滤波实质上是一种运算过程， 实现对信号的运算处理。 输入数字信号 （数字序列） 通过特定的运算转变为输出的数字序列，因此，数字滤波器本质上是一个完成特定运算 的数字计算过程，也可

3、以理解为是一台计算机。描述离散系统输出与输入关系的卷积和 差分方程只是给数字信号滤波器提供运算规则，使其按照这个规则完成对输入数据的处 理。时域离散系统的频域特性 ： Y ejX ej H ej 其中Y ej 、X ej分别是数字滤波器的输出序列和输入序列的频域特性(或称为 频谱特性)，H ej是数字滤波器的单位取样响应的频谱，又称为数字滤波器的频域响 应。输入序列的频谱 X ej经过滤波后X ej H ej ，因此，只要按照输入信号频谱 的特点和处理信号的目的，适当选择H ej ，使得滤波后的X ej H ej满足设计 的要求，这就是数字滤波器的滤波原理。 3.2脉冲响应不变法 脉冲响应不变

4、法是实现模拟滤波器数字化的一种直观而常用的方法。它特别适合于对 滤波器的时域特性有一定要求的场合。 脉冲响应不变法的核心原理是通过对连续函数ha(t)等间隔采样得到离散序列 ha(nT)。令h(n)= ha(nT)，功采样间隔。它是一种时域上的转换方法。 一个模拟滤波器的传递函数可以用有理分式表达式表示为： M (s) a。as a” s(3-1) H a(S)N b。bs bzs 通过反拉普拉斯变换我们就可以得到它的冲激相应： 1 ha(U L【Ha(S)l( 3-2) 脉冲响应不变法就是要保证脉冲响应不变，即： (3-3) (3-4) h(n) ha( nT) 对上式的冲激相应序列 h(n

5、)作Z变换，就可以得到数字滤波器的传递函数: H(z) ZTha (nT) 般来说，H a(S)的分母多项式阶次总是大于分子多项式的阶次。假定H a(S)的没 有多重极点，则式(1 )就可分解为: Ha(S) (3-5) 式中Si，Ai均为复数，S是H a(s)的极点。其拉氏变换为 N Si.T n ha(t)Ae u(n) i 1 般来说，Ha(S)的分母多项式阶次总是大于分子多项式的阶次。数字滤波器的传 递函数H(Z)经过合并简化，成为一般形式的有理分式传递函数 i CoCiZ dodiZ L ClZ dNZ N (3-6) 精选范本，供参考! 在讨论采样序列z变换与模拟信号拉氏变换之间关

6、系的有关章节中，我们已经知道 H(z)st Ha(s jy K) (3-7) ST 按照z=e的关系，每一个s平面上宽度为2 /T的水平条带将重迭映射到z平面 上。因此脉冲响应不变法将s平面映射到z平面，不是一个简单的对应的关系。对 于高采样频率(T小)的情况，数字滤波器在频域可能有极高的增益。为此我们采用 H(z) TAi SiT1 e z (3-8) N h(n) Tha( nT)TAeSTnu( n)(3-9) i 1 在脉冲响应不变法设计中，模拟频率与数字频率之间的转换关系是线性的 (*t)。同时，它可以保持脉冲响应不变，h(t) = ha(nT)。因此，这一方法往往用 于低通时域数字

7、滤波器设计及相应的模拟系统数字仿真设计。 3.3巴特沃斯低通滤波器 没有起伏，而在组频带 巴特沃斯滤波器的特点是同频带内的频率响应曲线最为平坦, 则逐渐下降为零。在振幅的对数对角频率的波特图上，从某一边界见频率开始，振幅随 着角频率的增加而逐渐减少，趋向于负无穷大。 一阶巴特沃斯滤波器的衰减率为每倍频 20分贝，二阶巴特沃斯滤波器的衰减率为 每倍频12分贝，三阶的衰减率为每分贝 18分贝，如此类推，巴特沃斯滤波器的振幅对 角频率单调下降，并且滤波器的结束越高，在组频带振幅衰减速度越快，其他滤波器高 阶的振幅对角频率图和低阶数的振幅对角频率有不同的形状。 2 S )N (3-10) 上述函数的特

8、点是等距离分布在半径为的圆上。 因此，极点用下式表示为 j j(2k 1) Sk ce 2 e k 0,1,2, N 1 (3-11) H a(S)的表示式: Ha(s) n c (s Sk) k 0 (3-12) 为了使设计公式和图表统 将频率归一化。巴特沃斯滤波器采用 3dB截止频率 归一化，归一化后的系统函数为 Ga( (3-13) c s c ,- c ,称为归一化频率，p称为归一化复变量，这 样巴特沃斯滤波器的归一化低通原型系统函数为 (3-14) 1 N! (P Pk) k 0 式中，pk S c，为归一化极点，用下式表示： 1 2k 1 pk e 2 2Nk 0,1,2, N 1

9、(3-15) 第四章设计过程 利用模拟滤波器来设计IIR低通数字滤波器是常用的方法，称之为模拟一数字转换法。 4.1设计步骤 利用在MATLAB设计IIR数字滤波器可分以下几步来实现 (1) 按一定规则将数字滤波器的技术指标转换为模拟低通滤波器的技术指标； (2) 根据转换后的技术指标使用滤波器阶数函数，确定滤波器的最小阶数N和截止 频率Wc ； (3) 利用最小阶数 N产生模拟低通滤波原型； (4) 利用截止频率 Wc把模拟低通滤波器原型转换成模拟带通原型； (5) 利用冲激响应不变法或双线性不变法把模拟滤波器转换成数字滤波器。 4.2频率混叠现象 数字滤波器的频率响应是模拟滤波器频率响应的

10、周期延拓。只有当模拟滤波器的 频率响应是限带的，且带限于折叠频率以内时，才能使数字滤波器的频率响应在折叠频 率以内重现模拟滤波器的频率响应，而不产生混叠失真。即 1 H(ej )Ha jT|w|n 11(4-1) 但是，任何一个实际的模拟滤波器频率响应都不是严格限带的(非理想)，变换 后就会产生周期延拓分量的频谱交叠，即产生频率响应的混叠失真。这时数字滤波器的 频响就不同于原模拟滤波器的频响，而是有一定的失真。当模拟滤波器的频率响应在折 叠频率以上处衰减越大、越快时，变换后频率响应混叠失真就越小。这时，采用脉冲响 应不变法设计的数字滤波器才能得到良好的效果 精选范本，供参考! 图4.1脉冲响应

11、不变法的频率混叠现象 总结以上，脉冲响应不变法的优点是频率坐标变换是线性的，即W=，如果不考虑频 率混叠现象，用这种方法设计的数字滤波器会很好地重现原模拟滤波器的频率特性。另 一个优点是数字滤波器的单位脉冲响应完全模仿模拟滤波器的单位脉冲响应，时域特性 逼近好。其缺点是会产生频谱混叠现象，使数字滤波器的频响偏离模拟滤波器的频响。 脉冲响应不变法适合低通、带通滤波器的设计，不适合高通、带阻滤波器的设计。 第五章程序设计 5.1设计巴特沃斯模拟滤波器 5.1.1模拟滤波器设计程序 clear %采样频率为周期倒数 Fs=1000; Wp=200*pi; Ws=600*pi; %设置归一化通带和阻带

12、截止频率 Ap=3; As=20; %设置通带最大和最小衰减 N,Wc=buttord(Wp,Ws,Ap,As, s）;%调用 butter 函数确定巴特沃斯滤波器阶数 B,A=butter(N,Wc, s); %调用 butter 函数设计巴特沃斯滤波器 W=(0:pi:1000*pi); %指定一段频率值 hs=freqs(B,A,W); %计算模拟滤波器的幅频响应 plot(W/pi/2,abs(hs)/abs(hs(1); %绘出巴特沃斯模拟滤波器的扶贫特性曲线 grid on ; title（ 巴特沃斯模拟滤波器 xlabel（ 频率 /Hz ）; ylabel（ 归一化幅值 ）;

13、5.1.2 用分贝显示幅值 clear Fs=1000; %采样频率为周期倒数 Wp=200*pi; Ws=600*pi; %设置归一化通带和阻带截止频率 Ap=3; As=20; %设置通带最大和最小衰减 N,Wc=buttord(Wp,Ws,Ap,As, s）;%调用 butter 函数确定巴特沃斯滤波器阶数 B,A=butter(N,Wc, s); %调用 butter 函数设计巴特沃斯滤波器 W=(0:pi:1000*pi); %指定一段频率值 hs=freqs(B,A,W); %计算模拟滤波器的幅频响应 hs0=abs(hs)/abs(hs(1) hs1=20*log10(hs0)

14、plot(W/pi/2,hs1);%绘出巴特沃斯模拟滤波器的幅频特性曲线 grid on ; title( 巴特沃斯模拟滤波器 ); xlabel( 频率 /Hz ); ylabel( 分贝 /db ); 5.2 设计巴特沃斯数字滤波器 5.2.1 数字滤波器设计程序 clear %采样频率为周期倒数 Fs=1000; Wp=200*pi; Ws=600*pi; Ap=3; As=20; N,Wc=buttord(Wp,Ws,Ap,As, B,A=butter(N,Wc, s); W=(0:pi:1000*pi); D,C=impinvar(B,A,Fs); W=(0:0.001*pi:pi)

15、 Hz=freqz(D,C,W); plot(W/pi,abs(Hz)/abs(Hz(1); grid on ; title( 巴特沃斯模拟滤波器 ); xlabel( 频率 /Hz ); ylabel( 归一化幅值 ); %设置归一化通带和阻带截止频率 %设置通带最大和最小衰减 s);%调用 butter 函数确定巴特沃斯滤波器阶数 %调用 butter 函数设计巴特沃斯滤波器 %指定一段频率值 %调用脉冲响应不变法 %返回频率响应 %绘出巴特沃斯数字低通滤波器的扶贫特性曲线 精选范本 ,供参考！ Fs=1000; Wp=200*pi; Ws=600*pi; Ap=3; As=20; N,W

16、c=buttord(Wp,Ws,Ap,As, B,A=butter(N,Wc, s); W=(0:pi:1000*pi); hs=freqs(B,A,W); D,C=impinvar(B,A,Fs); W=(0:0.001*pi:pi); Hz=freqz(D,C,W); Hz0=abs(Hz)/abs(Hz(1); %绘出巴特沃斯数字低通滤波器的幅频特性曲线 5.2.2 用分贝显示幅值 clear %采样频率为周期倒数 %设置归一化通带和阻带截止频率 %设置通带最大和最小衰减 s);%调用 butter 函数确定巴特沃斯滤波器阶数 %调用 butter 函数设计巴特沃斯滤波器 %指定一段频率

17、值 %计算出模拟滤波器的幅频响应 %调用脉冲响应不变法 %返回频率响应 Hz1=20*log10(Hz0); plot(W/pi,Hz1); grid on ; title( 巴特沃斯数字滤波器 ); xlabel( 频率 /Hz ); ylabel( 分贝 /db ); 5.3 程序中涉及到的函数介绍 Buttord 功能：用于计算巴特沃斯数字滤波器的阶数 N 和 3dB 截止频率 wc 。 N,wc=buttord(wp , ws, a p , a s) 调用参数 wp ， ws 分别为数字滤波器的通带、阻带截止频率的归一化值，要求： 0 w wp w 1 , 0 w ws w 11表示数

18、字频率 pi。文论坛 a p , a分别为通带最大衰减和组带最小衰减(dB)。 当 ws w wp 时，为高通滤波器 ; 当 wp 和 ws 为二元矢量时，为带通或带阻滤波器，这时 wc 也是二元向量。 N ， wc 作为 butter 函数的调用参数。 N, Q c=buttord( Q p, Q s, a p, a s, s) 用于计算巴特沃斯模拟滤波器的阶数N和3dB截止频率Q c。 Q p , Q s均为实际模拟角频率。 Butter 功能：利用 butter 函数直接设计各种形式的数字滤波器 (也可以设计模拟滤波器 ) B,A = BUTTER(N,Wn)，设计一个阶数为 n,频率为

19、 Wn的低通滤波器； B,A = BUTTER(N,Wn,ftype) 可以设计高通，带阻滤波器，其中 ftype 参数的形式可以 指定何种滤波器， ftype为high 时，设计一个阶数为n,频率为 Wn的高通滤波器；ftype为stop 时，得到滤波器阶数为 2*n ，频率范围为 Wn = W1,W2 的带阻滤波器； 中文论坛 Z,P,K = BUTTER(.) ，若返回值是三个，则分别是滤波器的零点，极点和增益，其他函 数的参数同上面 。 Plot(X,Y) 功能：绘制线性二维图形 当 X,Y 均为实数向量时，并且为维数相同， X=X(i),Y=Y(i) ，则 plot(X,Y) 先描述

20、点 (X(i),Y(i) ，然后依次画线 ； 当 X,Y 均为复数向量时，则不考虑虚数部分 ； 当X,Y均为实数矩阵时，并且为维数相同，plot依次按照对应的列画出线，矩阵有几列 就有几条线 当 X,Y 一个为向量，一个为矩阵时，并且向量的维数等于矩阵的行数或者列数，则把矩 阵按照向量的方向分解为几个向量，在与向量配对分别画图，矩阵分解几个向量就有几 条线。 Grid on 2 维 /3 维绘图设置网格线 grid 功能：来设置当前坐标系的网格线的开 / 闭。 grid on 功能：对当前坐标系添加主要网格线。 grid off 功能：对当前坐标系移除主要 / 细节网格线。 grid(axes

21、_handle,.) 输入： axes_handle 坐标系句柄 功能：对指定的坐标系设置网格线是否显示。 与信号处理 grid minor 功能：对当前坐标系设置添加细节网格 Impinvar 功能：用“脉冲响应不变法”将模拟滤波器离散化。 BZ,AZ =impinvar(B,A,Fs) 把具有B,A模拟滤波器传递函数模型转换成采样频率为Fs(Hz)的数字滤波器的传递函 数模型BZ,AZ。采样频率Fs的默认值为Fs=1o Freqz 、Freqs Freqz 是计算数字滤波器的频率响应的函数 Freqs 是计算模拟滤波器的频率响应的函数 H,F = freqz(B,A,N,Fs) 其中 B/

22、A提供滤波器系数B为分子 A为分母(bO + b1ZA-1 +.)/(a0 + a1ZA-1 +.) N 表示选取单位圆的上半圆等间距的 N 个点作为频响输出 ; Fs 为采样频率，该参数可以省略 Matlab 中文论坛 H 为 N 个点处的频率响应复值输出向量，其模即为频响幅值曲线幅值 2Olog1O(abs(H)DB ，其幅角 angle(H) 即为频响相位曲线相位值。 F 为与第 N 点处对应的频率值 f(Hz) ，如果 Fs 参数省略时， 则频率值 w 为 rad/sample ， w = 2*pi*f/Fs book.iLoveM H,F=freqz(B,A,N,whole) 中文论

23、坛 其参数意义如上，只是调用时角频率的范围由 O-pi 扩展到了 O-2*pi ，F 返回的是离散系 统频率响应 H(exp(j*w) 在 O-2*pi 范围内 N 个频率等分点的对应 样值的列向量。 freqz(B,A,N) ，调用时并不返回系统频率响应的样值，而是直接绘制系统的幅频和相频 特性曲线。但要注意，该调用绘制的是对数幅频特性曲线 freqz(B,A,N,whole) ，调用时绘制 O-2*pi 频率范围内系统对数幅频和相频特性曲线 freqs 的功能与 freqz 类似。 Filter 功能：FILTER是一维数字滤波器 Y = FILTER(B,A,X) ，输入 X 为滤波前序

24、列， Y 为滤波结果序列， B/A 提供滤波器系数， B 为分子， A 为分母 整个滤波过程是通过下面差分方程实现的： a(1)*y(n)= b(1)*x(n)+ b(2)*x(n-1)+ . + b(nb+1)*x(n-nb)-a(2)*y(n-1)- . - a(na+1)*y(n-na) 中文论坛 Y,Zf = FILTER(B,A,X,Zi),输入X为滤波前序列，Y为滤波结果序列，B/A提供滤波器 系数，B为分子，A为分母，并输入 Zi指定X的初始状态，Zf为最终状态矢量 第六章运行波形 6.1巴特沃斯模拟滤波器幅频响应 巴特沃斯模拟滤波器运行波形图如下 a !： 领宝dr 精选范本，供参考! 图6.1巴特沃斯模拟低通滤波器归一化幅频响应 图6.2巴特沃斯模拟低通滤波器幅频响应（分贝） 如图6.1、6.2所示，模拟低通滤波器的幅频响应与所给参数基本相符。300Hz处的幅 值衰减为20db ，100Hz处的幅值衰减为1db。通带的指标与要求略有出入，这是因为在 程序中计算的时候调用的公式保证了阻带衰减，而通带的指标有富余。 6.2巴特沃斯数字滤波器幅频响应 数字滤波器的幅频响应图横轴为数字域频率，数字域频率为模拟角频率对采样频率的 归一化频率。程序设计中，横轴为数字域频率除以n后的值，所以横轴为0-1。 0.2X0 &O.& .S