双线性变换法课件.ppt

1、1数字信号处理数字信号处理by Zaiyue YangCSE, ZJU, 20162第第5章章 数字滤波器的设计数字滤波器的设计5.1 数字滤波器的基本概念数字滤波器的基本概念5.2 模拟滤波器的设计模拟滤波器的设计5.3 用脉冲响应不变法设计用脉冲响应不变法设计IIR数字低通滤波器数字低通滤波器5.4 用双线性变换法设计用双线性变换法设计IIR数字低通滤波器数字低通滤波器5.5 数字高通、带通和带阻滤波器的设计数字高通、带通和带阻滤波器的设计3SpecificationsDesired IIRButterworth, Chebyshev, Ellipse, Bessel, etc.脉冲响应不

2、变法阶跃响应不变法双线性变换法上一堂课学了些什么？LP to HP, BP, BSHow to design analog lowpass filtersButterworth：Chebyshev：Elliptic：它们的特点是什么？如何确定它们的参数？4221()1()aNcHj2221()1()aNpHjC2221()1()aNpHjUUsing MatlabButterworth： buttap, buttord, butterChebyshev： cheb1ap, cheb1ord, cheby1, cheb2ap, cheb2ord, cheby2Elliptic： ellipap,

3、 ellipord, ellip5How to transform LP to HP/BP/BSFrequency band transformations: LP to HP: LP to BP: LP to BS:6ph 220wB 220wB Using Matlab 仍然调用设计模拟低通滤波器的函数 注意通带和阻带的截止频率 LP: p s HP: s p BP: sl pl pu su BS: pl sl su pu 注意加high、stop7Analog to Digital/ S to ZGeneral Principle:Preserve certain “aspects” o

4、f analog and digital filters8AspectsTransformationsShape of the impulse responseImpulse invariance *Shape of step responseStep invarianceSystem function representationBilinear *Goal:Design digital filter H(z) for a given analog filter Ha(s)95.3 脉冲响应不变法脉冲响应不变法Preserve the shape of the impulse respons

5、eThe digital filter impulse response looks “similar” to that of an analog filter( )( )|()at nTah nh th nTh(n)在采样点上等于ha(t)脉冲响应不变法是一种时域逼近方法数字滤波器H(z)的单位脉冲响应模拟滤波器Ha(s)的单位冲击响应采样时间Goal: 已知Ha(s)，寻找H(z)，使得h(n) = ha(nT)How to solve it?1. 从从Ha(s)到到ha(t)2. 从从ha(t)到到h(n)3. 从从h(n)到到H(z)101( )NiaiiAHsss 111( )( )

6、( ) 1iNs taaaiih tLSHs LStLSHsAe1( )()iNs nTaiih nh nTAe0011101( )( )1iiiNs nTnninniNNs nTniisTiniH zh n zAezAAezez结论: 已知 求得111( )NiaiiAHsss11( )1iNisTiAH zez条件: 模拟滤波器Ha(s)只有单阶极点，且分母多项式的阶次高于分子多项式的阶次注意: 对于同样Ha(s)，采样时间T不同，H(z)不同Ha(s)极点与H(z)极点的映射关系12 例5.3.1 ：已知模拟滤波器的传输函数Ha(s)为 用脉冲响应不变法将Ha(s)转换成数字滤波器的系统

7、函数H(z)。 解 首先将Ha(s)写成部分分式：20.5012( )0.64490.7079aHsss0.32240.3224( )0.32240.77720.32240.7772ajjHssjsj极点为12(0.32240.772),(0.32240.7772)sjsj 那么H(z)的极点为131212,s Ts Tzeze经过整理，得到： 设T=1s时用H1(z)表示，T=0.1s时用H2(z)表示，则111212120.3276( )1 1.03280.2470.0485( )1 1.93070.9375zH zzzzHzzz脉冲响应不变法的sz映射关系：脉冲响应不变法对应的脉冲响应不

8、变法对应的s平面到平面到z平面的映射关系为：平面的映射关系为：141( )NiaiiiAHssss的极点11( )1iiNsTisTiAH zeez的极点Wild guess!sTze脉冲响应不变法的性质：1. 因果稳定的因果稳定的Ha(s)得到因果稳定的得到因果稳定的H(z)15设那么，,jsjzre jTj TreeeTreT sTze 即： =0, r=1 0, r0, r1s平面z平面虚轴单位圆左半平面单位圆内右半平面单位圆外数字频率模拟频率模拟频率和数字频率转换关系16 例：用脉冲响应不变法设计数字低通滤波器，指标参数如下， 解：（1）将数字频率转换为模拟频率0.2 ,1,0.35

9、,40ppssadBadB（2）设计模拟滤波器Ha(s)（3）将Ha(s) 转换成H(z)0.2,10.35,40pppsssadBTTadBTT Matlab实现实现T=1;wp=0.2*pi/T; ws=0.35*pi/T; ap=1; as=40;N,wc=buttord(wp,ws,ap,as,s);B,A=butter(N,wc,s);Bz,Az=impinvar(B,A);Freqz(Bz,Az);00.10.20.30.40.50.60.70.80.91-1000-800-600-400-2000Normalized Frequency ( rad/sample)Phase (d

10、egrees)00.10.20.30.40.50.60.70.80.91-150-100-50050Normalized Frequency ( rad/sample)Magnitude (dB)脉冲响应不变法的性质：2. 映射的周期性：变化2/T的整数倍，映射值不变182,jM TsTTj TTTzeeeeeM为任意整数物理意义： 将s平面沿着虚轴分割成一条条宽为2/T的水平带，每条水平带都对应着整个z平面。线性关系： 模拟频率从-/T变化到/T ，数字频率线性地从-变到19图5.3.1 z=esT, s平面与z平面之间的映射关系脉冲响应不变法的缺点：若ha(t)的频带不在/T、-/T之间，

11、那么在水平带边界（即/T的奇数倍）附近产生频谱混叠20为什么？（从采样定理考虑）采样频率2sT 22ssT 采样频率小于2倍信号频率，不满足采样定理，因此在频域会有混叠采样周期T越小，混叠越少，但T不可能无限逼近于021图5.3.2 脉冲响应不变法的频率混叠现象22图5.3.3 例5.3.1的幅度特性H1(z), T=1sH2(z),T=0.1s从图中能观察到什么现象？总结：脉冲响应不变法的特点23频率变换关系是线性的，即=T，如果不存在频谱混叠，该方法可以很好地重现模拟滤波器的频响特性；数字滤波器的单位脉冲响应完全模仿模拟滤波器的单位冲击响应波形，时域特性逼近好会产生不同程度的频谱混叠失真，

12、适合于带限（低通、带通）滤波器，不适合于非带限（高通，带阻）滤波器245.4 双线性变换法双线性变换法目的：为了克服模拟频率在/T奇数倍附近的频谱混叠办法：先将s平面上整个虚轴压缩到/T之间，再用zesT转换到z平面注意：1. 压缩一定是非线性的2. 可以有多种压缩方式25图5.4.1 双线性变换法的映射关系s1和1为中间变量非线性压缩脉冲响应不变法26一种常用的非线性压缩方法 正切变换： 121tan()2TT 式中T仍是采样间隔，当由经过0变化到+时，1从/T经过0变化到/T，实现了s平面上整个虚轴完全压缩到s1平面上虚轴的/T之间的转换。结论：双线性变换关系1122 1,21szTszT

13、zsT271s Tze112 11s Ts TesTe如何得到？111111112211112211sin()21222 12tan()121cos()2jTjTjTjTjTjTjTeeejjTTTTTeTee 令s=j, s1=j1112 11zsTz22sTzsT28双线性变换法的频率变换关系121tan2TT 图5.4.2 双线性变换法的频率变换关系21tan2T 1T 脉冲不变法的特性双线性变换法的性质：1. 使s左半平面映射到z平面单位圆的内部，因此Ha(s)因果稳定，那么H(z)也因果稳定2. 非线性频率变换关系 ，消除了频谱混叠现象3. 然而，非线性关系使得数字滤波器频响曲线不能

14、保真地模仿模拟滤波器的频响曲线形状2921tan2T 当 的刻度是均匀的，的刻度是均匀的，而是随着的增加越来越密集30图5.4.3 双线性变换法幅度和相位特性的非线性映射31 例5.4.1试分别用脉冲响应不变法和双线性不变法将图5.4.4所示的RC低通滤波器转换成数字滤波器。 解 首先按照图6.4.4写出该滤波器的传输函数Ha(s)为1( ),aHssRC 利用脉冲响应不变法转换，数字滤波器的系统函数H1(z)为 11( )1TH zez图5.4.4 RC低通滤波器32 利用双线性变换法转换，数字滤波器的系统函数H2(z)为 111121212112(1)( )( )12,22azsTzzHz

15、Hsa zTTTTH1(z)和H2(z)的网络结构分别如图6.4.5(a),(b)所示。图5.4.5 例5.4.1图H1(z)和H2(z)的网络结构 (a)H1(z); (b)H2(z) 33图5.4.6例5.4.1 图Ha(s)、H1(z)和H2(z)的幅频特性 能从图中观察到什么？34 下面我们总结利用模拟滤波器设计IIR数字低通滤波器的步骤。 (1)确定数字低通滤波器的技术指标：通带截止频率p、通带衰减ap、阻带截止频率s、阻带衰减as。 (2)将数字低通滤波器的技术指标转换成模拟低通滤波器的技术指标。如采用脉冲响应不变法，边界频率的转换关系为： ,pspsTT 如果采用双线性变换法，边

16、界频率的转换关系为2121tan(),tan()22ppssTT 35 (3)按照模拟低通滤波器的技术指标设计模拟低通滤波器。 (4)将模拟滤波器Ha(s)，从s平面转换到z平面，得到数字低通滤波器系统函数H(z)。 例5.4.2 设计低通数字滤波器，要求在通带内频率低于0.2rad时，容许幅度误差在1dB以内；在频率0.3到之间的阻带衰减大于15dB。指定模拟滤波器采用巴特沃斯低通滤波器。试分别用脉冲响应不变法和双线性变换法设计滤波器。36 解 (1) 用脉冲响应不变法设计数字低通滤波器。 数字低通的技术指标为 p=0.2rad, ap=1dB; s=0.3rad, as=15dB 模拟低通

17、的技术指标为 T=1s, p=0.2rad/s, ap=1dB; s=0.3rad/s, as=15dB37 设计巴特沃斯低通滤波器。先计算阶数N及3dB截止频率c。 0.10.10.31.50.21010.092101lg0.0925.8846lg1.5psspN 38 取N=6。为求3dB截止频率c,将p和ap代入(6.2.17)式，得到c=0.7032rad/s。 根据阶数N=6，查表6.2.1，得到归一化传输函数为23456( )11 3.86377.46419.14167.46413.8637aHppppppp 为去归一化，将p=s/c代入Ha(p)中，得到实际的传输函数Ha(s),

18、 3962652433425665432( )3.86377.46419.14167.46413.86370.12092.7163.6913.1791.8250.1210.1209accccccHsssssssssssss 用脉冲响应不变法将Ha(s)转换成H(z)。首先将Ha(s)进行分解，然后根据 得到：1112121120.28710.44662.14281.1454( )10.12970.69491 1.06910.36991.85580.630410.99720.2570zzH zzzzzzzzisTize40图5.4.7 例5.4.2图用脉冲响应不变法设计的数字低通滤波器的幅度特性

19、41 (2) 用双线性变换法设计数字低通滤波器。 数字低通技术指标仍为 p=0.2rad, ap=1dB; s=0.3rad, as=15dB 模拟低通的技术指标为21tan,122tan0.10.65/ ,12tan0.151.019/ ,15ppPpssTTrad sdBrad sdB 42 设计巴特沃斯低通滤波器。阶数N计算如下：0.10.11.0191.5680.651010.092101lg0.0925.3066lg1.568psspN 取N=6。为求c,将s和as代入(5.2.18)式中，得到c=0.7662rad/s。43 根据N=6，查表5.2.1得到的归一化传输函数Ha(p)

20、与脉冲响应不变法得到的相同。为去归一化，将p=s/c代入Ha(p)，得实际的Ha(s)， 用双线性变换法将Ha(s)转换成数字滤波器H(z)：222( )0.2024(0.3960.5871)(1.0830.5871)(1.4800.5871)aHsssssss111211 6121212( )( )0.0007378(1)1(1 1.2680.7051)(1 1.0100.358) 1 0.90440.2155azszH zHszzzzzzzMatlab实现实现T=1;wpz=0.2*pi; wsz=0.3*pi;wp=2/T*tan(wpz/2); ws=2/T*tan(wsz/2); a

21、p=1; as=15;N,wc=buttord(wp,ws,ap,as,s);B,A=butter(N,wc,s);Bz,Az=bilinear(B,A,1/T);45图5.4.8 例5.4.2图用双线性变换法设计的数字低通滤波器的幅度特性465.5 数字高通、带通和带阻滤波器的设计数字高通、带通和带阻滤波器的设计 例如高通数字滤波器等。具体设计步骤如下： (1) 确定所需类型数字滤波器的技术指标。 (2) 将所需类型数字滤波器的技术指标转换成所需类型模拟滤波器的技术指标，转换公式为 21tan2T 47 (3)将所需类型模拟滤波器技术指标转换成模拟低通滤波器技术指标(具体转换公式参考本章6.

22、2节)。 (4)设计模拟低通滤波器。 (5)将模拟低通通过频率变换，转换成所需类型的模拟滤波器。 (6)采用双线性变换法，将所需类型的模拟滤波器转换成所需类型的数字滤波器。48 例5.5.1 设计一个数字高通滤波器，要求通带截止频率p=0.8rad，通带衰减不大于3dB，阻带截止频率s=0.44rad,阻带衰减不小于15dB。希望采用巴特沃斯型滤波器。 解 (1)数字高通的技术指标为 p=0.8rad,p=3dB; s=0.44rad,s=15dB49 (2) 模拟高通的技术指标计算如下： 令T=1，则有1tan6.155/ ,321tan1.655/ ,32pppsssrad s adBra

23、d s adB (3)模拟低通滤波器的技术指标计算如下：1/ ,36.1553.719/ ,151.655ppssrad sdBrad sdB50 (4)设计归一化模拟低通滤波器G(p)。模拟低通滤波器的阶数N计算如下：0.10.11010.18031013.7191.3042psspN51 查表5.2.1，得到归一化模拟低通传输函数G(p)为21( )21G ppp (5) 将模拟低通转换成模拟高通。将上式中G(p)的变量换成p/s，得到模拟高通Ha(s)： 222( )( )|2papppssHsG pss (6)用双线性变换法将模拟高通Ha(s)转换成数字高通H(z)：111212121

24、0.13260.26530.1326( )( )1 0.73940.2699azszzzH zHszzMatlab实现实现%DIGI_BUTT_HPwpz=0.8; wsz=0.44; ap = 3; as = 15;N,wc=buttord(wpz,wsz,ap,as)Bz,Az=butter(N,wc,high)fk=0:1/200:1;Hk=freqz(Bz,Az,fk*pi);plot(fk,20*log10(abs(Hk);xlabel(Frequency(pi);ylabel(Magnitude(dB);Matlab实现实现Output: N = 2; wc = 0.6978; B

25、z = 0.1326 -0.2653 0.1326; Az = 1.0000 0.7394 0.269900.10.20.30.40.50.60.70.80.91-100-90-80-70-60-50-40-30-20-100Frequency(pi)Magnitude(dB)54 例5.5.2设计一个数字带通滤波器，通带范围为0.3rad到0.4rad，通带内最大衰减为3dB，0.2rad以下和0.5rad以上为阻带，阻带内最小衰减为18dB。采用巴特沃斯型模拟低通滤波器。 解 (1)数字带通滤波器技术指标为 通带上截止频率 u=0.4rad 通带下截止频率 l=0.3rad55 阻带上截止

26、频率 s2=0.5rad 阻带下截止频率 s1=0.2rad 通带内最大衰减p=3dB，阻带内最小衰减s=18dB。 56 (2) 模拟带通滤波器技术指标如下： 设T=1，则有2211012tan1.453/212tan1.019/212tan2/212tan0.650/21.217/0.434/uullssssulwulrad srad srad srad srad sBrad s (通带中心频率) (带宽) 57 (3) 模拟归一化低通滤波器技术指标： 归一化阻带截止频率222022.902sswsB归一化通带截止频率p=1p=3dB,s=18dB58 (4) 设计模拟低通滤波器：0.10

27、.11010.1271012.902lg0.1271.9402lg2.902psspN 查表5.2.1，得到归一化低通传输函数G(p),21( )21G ppp59 (5) 将归一化模拟低通转换成模拟带通： (6)通过双线性变换法将Ha(s)转换成数字带通滤波器H(z)。220()( )( )ulaspsHsG p 11121zsz60 将上式代入(5)中的转换公式，得11221 221 20021211212224(1)(1)()2(1)()5.484.57.4816.3135.18806190.868(1)1zsululzszzpszzzzzzz 将上面的p等式代入G(p)中，得 2412

28、340.021(12)( )1 1.4912.8481.681.273zzHzzzzzMatlab实现实现%Digi_Butt_BPwp=.3,.4;ws=.2,.5;Ap=3;As=18;N,wc=buttord(wp,ws,Ap,As);Bz,Az=butter(N,wc);fk=0:1/200:1;Hk=freqz(Bz,Az,fk*pi);plot(fk,20*log10(abs(Hk);xlabel(Frequency(pi);ylabel(Magnitude(dB);NBzAzMatlab实现实现Output: N = 2; Bz = 0.0213 0 -0.0426 0 0.02

29、13; Az = 1.0000 -1.6303 2.2183 -1.2919 0.632000.10.20.30.40.50.60.70.80.91-120-100-80-60-40-200Frequency(pi)Magnitude(dB)63 例5.5.3设计一个数字带阻滤波器，通带下限频率l=0.19,阻带下截止频率s1=0.198，阻带上截止频率s2=0.202，通带上限频率u=0.21，阻带最小衰减s=13dB，l和u处衰减p=3dB。采用巴特沃斯型。 解 (1) 数字带阻滤波器技术指标： l=0.19rad,u=0.21rad,p=3dB; s1=0.198rad,s2=0.202

30、rad,s=13dB64 (2) 模拟带阻滤波器的技术指标： 设T=1，则有1122112tan0.615/ ,2tan0.685/22112tan0.615/ ,2tan0.685/22lluussssrad srad srad srad s 阻带中心频率平方为 20=lu=0.421阻带带宽为 Bw=u-l=0.07rad/s65 (3) 模拟归一化低通滤波器的技术指标： 按照(5.2.48)式，有 p=1,p=3dB222204.434,13wssssBdB66 (4) 设计模拟低通滤波器：0.10.11010.2291014.434lg0.2290.991lg4.434psspN 22

31、0220( )( )wwasBpssBpsHsG p(5) 将G(p)转换成模拟阻带滤波器Ha(s)：67 (6) 将Ha(s)通过双线性变换，得到数字阻带滤波器H(z)。1121 221 202221 221 21200112122(1)4(1)(1)2(1)4(1)(1)0.969(1 1.619)( )( )1 1.5690.939zszzBpzzsBzBpszzzzH zG pzzMatlab实现实现%Digi_Butt_BSwp=.19,.21;ws=.198,.202;Ap=3;As=13;N,wc=buttord(wp,ws,Ap,As);Bz,Az=butter(N,wc,stop);fk=0:1/200:1;Hk=freqz(Bz,Az,fk*pi);plot(fk,20*log10(abs(Hk);xlabel(Frequency(pi);ylabel(Magnitude(dB);NBzAzMatlab实现实现Output: N = 1; Bz = 0.9732 -1.5748 0.9732; Az = 1.0000 -1.5748 0.946500.10.20.30.40.50.60.70.80.91-60-50-40-30-20-100Frequency(pi)Magnitude(dB)