數(shù)字信號處理-數(shù)字濾波器設計方法論

1、1第第5章章 數(shù)字濾波器的設計數(shù)字濾波器的設計 5.1 數(shù)字濾波器的基本概念數(shù)字濾波器的基本概念5.2 模擬濾波器的設計模擬濾波器的設計5.3 用脈沖響應不變法設計用脈沖響應不變法設計IIR數(shù)字低通濾波器數(shù)字低通濾波器5.4 用雙線性變換法設計用雙線性變換法設計IIR數(shù)字低通濾波器數(shù)字低通濾波器5.5 數(shù)字高通、帶通和帶阻濾波器的設計數(shù)字高通、帶通和帶阻濾波器的設計2數(shù)字濾波器：是指輸入輸出均為數(shù)字信號，通過一定運算關系改變輸入信號所含頻率成分的相對比例或者濾除某些頻率成分的器件。 高精度、穩(wěn)定、體積小、重量輕、靈活，不要求阻抗匹配，可實現(xiàn)特殊濾波功能優(yōu)點：5.1 數(shù)字濾波器的基本概念3 1.

2、數(shù)字濾波器的分類(1)總體分類： 經典濾波器，信號與干擾分占不同頻帶：選頻濾波器 現(xiàn)代濾波器，信號與干擾頻帶重疊：維納濾波器、卡爾曼濾波器、自適應濾波器等(2)濾波特性：低通、高通、帶通、帶阻等(3)網絡結構：IIR、FIR1001( )( )( )1MrrNnrNknkkb zH zH zh n za z4理想低通、高通、帶通、帶阻數(shù)字濾波器幅度特性 )(ejH)(ejH)(ejH)(ejH0低通0高通0帶通0帶阻222222225 2.數(shù)字濾波器的技術指標 我們通常用的數(shù)字濾波器一般屬于選頻濾波器。假設數(shù)字濾波器的傳輸函數(shù)H(ej)用下式表示： 數(shù)字低通濾波器的技術要求()()()jjjH

3、 eH ee 通帶紋波幅度阻帶紋波幅度通帶截止頻率3dB通帶截止頻率阻帶截止頻率6 通帶內和阻帶內允許的衰減一般用dB數(shù)表示，通帶內允許的最大衰減用ap表示，阻帶內允許的最小衰減用as表示，ap和as分別定義為：00()20lgdB()()20lgdB()psjpjjsjH eaH eH eaH e(5.1.3) (5.1.4) 如將|H(ej0)|歸一化為1，(5.1.3)和(5.1.4)式則表示成：20lg() dB20lg() dBpsjpjsaH eaH e (5.1.5) (5.1.6) 7 3. 數(shù)字濾波器設計方法概述 SpecificationsDesired IIRButter

4、worth, Chebyshev, elliptic, Bessel, etc.脈沖響應不變法階躍響應不變法雙線性變換法 IIR濾波器和FIR濾波器的設計方法是不同，IIR濾波器設計需要借助于模擬濾波器來完成 回想控制器設計問題85.2 模擬濾波器的設計模擬濾波器的設計 模擬濾波器的理論和設計方法已發(fā)展得相當成熟，且有若干典型的模擬濾波器供我們選擇，如巴特沃斯(Butterworth)濾波器、切比雪夫(Chebyshev)濾波器、橢圓(Elliptic)濾波器、貝塞爾(Bessel)濾波器等，這些濾波器都有嚴格的設計公式、現(xiàn)成的曲線和圖表供設計人員使用。 9 圖5.2.1 各種理想模擬濾波器的

5、幅頻特性 )(jaH低通帶通帶阻高通)(jaH)(jaH)(jaH000c10 1.模擬低通濾波器的設計指標及逼近方法模擬低通濾波器的設計指標及逼近方法 模擬低通濾波器的設計指標有ap, p,as和s。 p：通帶截止頻率s：阻帶截止頻率 ap：通帶最大衰減系數(shù)as：阻帶最小衰減系數(shù)22()10lg()apapHjaHj(5.2.1) (5.2.2) 22()10lg()asasHjaHj 將=0處幅度已歸一化到1，即|Ha(0)|=1,得到210lg()papaHj 210lg()sasaHj (5.2.3) (5.2.4) 11 技術指標給定后，設計一個傳輸函數(shù)Ha(s)，希望其幅度平方函數(shù)

6、滿足給定的指標ap和as。 一般濾波器的單位沖激響應為實數(shù)，因此2()( )()()()aaasjaaHjHs HsHjHj (5.2.5) 注意：1. 從數(shù)學上講，有無限多種Ha(s)能滿足指標2. 實際中，針對結構給定的Ha(s) （如：Butterworth, Chebyshev, Elliptic等），選取合適的參數(shù)12 幅度平方函數(shù)：221()1()aNcHj(5.2.6) 兩個參數(shù)：N, c2.巴特沃斯低通濾波器的設計方法巴特沃斯低通濾波器的設計方法13 將|Ha(j)|2寫成s的函數(shù)： 21( )()1()aaNcHs Hssj(5.2.7) N階Butterworth濾波器，|

7、Ha(j)|2=Ha(s)Ha(-s)有2N個極點，極點sk為：12121()22( 1)()NkckjNcsje (5.2.8) 三階巴特沃斯濾波器極點分布|Ha(j)|2的2N個極點均勻分布在半徑為c的圓上14 為形成穩(wěn)定的濾波器，2N個極點中只取s平面左半平面的N個極點構成Ha(s)，LHP( )()NcakHsss例：設N=3，極點有6個，它們分別為(5.2.10) 22330121133345jjcccjjcccsessesesse 取s平面左半平面的極點s0, s1, s2組成Ha(s)： 32233( )()()()cajjcccHssss15 歸一化：由于各濾波器的幅頻特性不同

8、，為使設計統(tǒng)一，需要將所有的頻率歸一化這里采用對3dB截止頻率c歸一化，歸一化后的Ha(s)表示為 令歸一化復變量p=s/c，pk=sk/c，得到歸一化巴特沃斯的傳輸函數(shù)101( )()aNkkccHsss(5.2.11) 101( )()aNkkHppp(5.2.12) 16 式中，pk為歸一化極點，用下式表示： 將極點表示式(5.2.13)代入(5.2.12)式，得到的Ga(p)的分母是p的N階多項式，用下式表示： 121()22,0,1,1kjkNkcspekN(5.2.13)(5.2.15) 11101( )aNNNGppbpb pb歸一化后Ga(p)的N個極點均勻分布在左半單位圓上表

9、6.2.1 巴特沃斯歸一化低通濾波器參數(shù) 20求參數(shù)N: 通帶指標: 阻帶指標:/10/10101()101psapNas/10/10101lglg101psapasN/1021 ()10papNc N可能有小數(shù)部分，應取大于等于N的最小整數(shù)/1021 ()10saNsc/10/10101lglg101psapasN(5.2.18)21求參數(shù) c: 3dB截止頻率c可由通帶指標/阻帶指標分別得到: 通帶指標10.12(101)paNcp 阻帶指標(5.2.19)(5.2.20) 10.12(101)saNcs 22低通巴特沃斯濾波器的設計步驟總結： (2) 根據(jù)N求出歸一化極點pk，以及歸一化

10、傳輸函數(shù)Ga(p):/10/10101lglg101psapasN (4) Ga(p)去歸一化:121()22kjNkpe101( )()aNkkGppp(3) 求c:110.10.122(101)or(101)psaaNNcpcs 將p=s/c, pk=sk/c代入Ga(p)，得到實際的傳輸函數(shù)Ha(s) (1) 根據(jù)技術指標p, ap, s和as，求濾波器的階數(shù)N:23 例5.2.1 已知通帶截止頻率fp=5kHz，通帶最大衰減p=2dB，阻帶截止頻率fs=12kHz，阻帶最小衰減s=30dB，按照以上技術指標設計巴特沃斯低通濾波器。 解： (1) 確定階數(shù)N 0.10.11010.024

11、210122.42lg0.02424.255lg2.4psaaspffN 24 (2) 其極點為3455016523754,jjjjjpepepepepe歸一化傳輸函數(shù)為401( )()akkGppp25 上式分母可以展開成為五階多項式，或者將共軛極點放在一起，形成因式分解形式。這里不如直接查表5.2.1；由N=5，直接查表得到： 極點：-0.3090j0.9511,-0.8090j0.5878; -1.00005432432101( )aHppa pa pa pa pa 其中，a0=1.0000, a1=3.2361, a2=5.2361, a3=5.2361, a4=3.236126 (3

12、) 求3dB截止頻率c 按照(5.2.19)式，得到：10.12(101)25.2755/paNcpkrad s (4) Ga(p)去歸一化，將p=s/c, pk=sk/c代入Ha(p)中得到：55423324543210( )accccccHssasasasasa27Matlab實現(xiàn)常用函數(shù)：buttap, buttord, butter%Butterworth_Lowwp=2*pi*5000;ws=2*pi*12000;Ap=2;As=30;N,wc=buttord(wp,ws,Ap,As,s);B,A=butter(N,wc,s);fk=0:14000/512:14000;wk=2*pi

13、*fk;figure(1); freqs(B,A,wk);Hk=freqs(B,A,wk);figure(2);plot(fk/1000,20*log10(abs(Hk);xlabel(Frequency(kHz);ylabel(Magnitude(dB);axis(0,14,-40,5);28103104-200-1000100200Frequency (rad/s)Phase (degrees)10310410-210-1100Frequency (rad/s)Magnitude02468101214-40-35-30-25-20-15-10-505Frequency(kHz)Magnit

14、ude(dB)293.切比雪夫低通濾波器的設計方法Motivation：Butterworth的頻率特性曲線在通帶和阻帶內都是單調遞減若通帶（阻帶）邊界處滿足指標，則通帶（阻帶）內會有較大余量更有效的方法在通帶（或阻帶、兩者）之內等波紋(equiripple)，可降低濾波器階數(shù)Chebyshev和Elliptic濾波器都是根據(jù)此原則設計的30切比雪夫低通濾波器：I型通帶等波紋，阻帶單調遞減II型通帶單調遞減，阻帶等波紋31 0 1，表示通帶內幅度波動的程度，愈大，波動幅度愈大；p為通帶截止頻率；CN(x)為N階切比雪夫多項式cos(arccos )1( )cosh(arcosh )1NNxxC

15、xNxx22221()1()aNpAHjC 切比雪夫I型低通濾波器幅度平方函數(shù)：(5.2.24) cosh2xxeex雙曲函數(shù)：p給定，兩個參數(shù)和N32當： N=0時，C0(x)=1 N=1時，C1(x)=x N=2時，C2(x)=2x2-1 N=3時，C3(x)=4x3-3x可歸納出高階切比雪夫多項式的遞推公式為 CN+1(x)=2xCN(x) - CN-1(x) (5.2.25) CN(x)特性： (1)切比雪夫多項式的過零點在|x|1的范圍內； (2)當|x|1時，|CN(x)|1,在|x|1時，CN(x)是雙曲線函數(shù)，隨x單調上升。33N=0,4,5切比雪夫多項式曲線Observati

16、ons:CN(x)有N個過零點N為偶數(shù)時，在x正半軸有N/2個過零點N為奇數(shù)時，在x正半軸有(N+1)/2個過零點CN(x)在x正半軸的過零點對應著|Ha(j)|=134切比雪夫型濾波器幅頻特性 N=3N=435|Ha(j)|的性質：x = 0(或 = 0)x = 1(或 = p)0 x 1(或0p) |Ha(j)|在1和 之間振蕩x 1 (或 p) |Ha(j)|單調遞減到02( 0)1for odd1( 0)for even1aaHjNHjN21( 1)for all 1aHjN21136 切比雪夫型與巴特沃斯低通的A2()曲線N相同時，切比雪夫I型性能更好37|Ha(s)|的極點分布：N

17、階|Ha(s)|2有2N個極點，令分布在左平面的N個穩(wěn)定極點為|Ha(s)|的極點。通過繁瑣的計算可知極點sj為，,0,.,1cos,sin2122iiiipiipiisjiNabiN 極點位于橢圓之上11112,22111NNNNab短軸為ap ，在實軸上長軸為bp，在虛軸上38三階切比雪夫濾波器的極點分布三個極點，i分別為:0 = 2/31 = 2 = 4/339歸一化與去歸一化11( )()aNiiGpcpp(5.2.38)c為待定系數(shù)，根據(jù)幅度平方函數(shù)(5.2.24)式可得：c = 2N-1111( )2()aNNiiGppp(5.2.40) 歸一化：令p = s/p， pi = si

18、/p ，1111( )2()2()NNppaNNNNiipiiHsssspcossin,0,.,1iiipajbiN去歸一化：(5.2.39)40指標與參數(shù)、N之間的關系決定通帶內最大衰減ap|Ha(j)|2最大值：|Ha(j)|2最小值：通帶內最大衰減ap ：N影響過渡帶寬度，以及通帶內波動的疏密（N等于通帶內極值的總個數(shù)）2max10lg10lg 1minpa2max()1aHj221min()1aHj0.12101pa(5.2.27)(5.2.26)求得41和N決定阻帶內最大衰減as222max10lg10lg 1()()ssNPasaCHj(5.2.44)2221()1()assNPH

19、jC0.122101()sasNPC 0.12101pa 0.10.1101arcosharcosh101spasaPN結合0.10.1101()cosharcosh101spassNaPPCN求得N42(1) 給定技術指標ap, p, as和s切比雪夫型濾波器設計步驟(2) 計算參數(shù)和N0.1101pa0.10.1101arcosharcosh101spasaPN(3) 求歸一化極點pi及系統(tǒng)Ga(p)21cossin,0,.,122iiiiipajbiNN111( )2()aNNiiGppp(4) 去歸一化得到實際系統(tǒng)Ha(s)11( )2()NpaNNipiHssp43 例5.2.2設計

20、低通切比雪夫濾波器，要求通帶截止頻率fp=3kHz，通帶最大衰減ap=0.1dB，阻帶截止頻率fs=12kHz，阻帶最小衰減as=60dB。 解：(1) 濾波器的技術要求： 0.1,2,60,2pppsssadBfadBf (2) 求階數(shù)N和：0.10.10.10.011016553,4101arcosh(6553)9.474.65arcosh(4)2.061011010.1526sppasapaN 44 (3) 求Ga(p):5(5 1)11( )0.1526 2()aiiGppp可求出N=5時的極點pi，得到：22( )12.442(0.5389)(0.33311.1949)(0.8720

21、0.6359)aGpppppp (4) 將Ga(p)去歸一化，得到：7261427141( )(1.0158 10 )(6.2788 104.2459 10 )11.6437 102.2595 10aHssssss45Matlab實現(xiàn)常用函數(shù)：cheb1ap, cheb1ord, cheby1, cheb2ap, cheb2ord, cheby2%Chebyshev1_Lowwp=2*pi*3000;ws=2*pi*12000;Ap=0.1;As=60;N,wp=cheb1ord(wp,ws,Ap,As,s);B,A= cheby1(N, Ap,wp,s);fk=0:14000/512:140

22、00;wk=2*pi*fk;figure(1);freqs(B,A,wk);Hk=freqs(B,A,wk);figure(2);plot(fk/1000,20*log10(abs(Hk);xlabel(Frequency(kHz);ylabel(Magnitude(dB);axis(0,14,-70,5);46103104-200-1000100200Frequency (rad/s)Phase (degrees)10310410-410-2100Frequency (rad/s)Magnitude02468101214-70-60-50-40-30-20-100Frequency(kHz)

23、Magnitude(dB)474.橢圓低通濾波器簡介Motivation：Butterworth的頻率特性曲線在通帶和阻帶內都是單調遞減Chebyshev-I在通帶等波紋， Chebyshev-II在阻帶等波紋Elliptic在通帶和阻帶內都是等波紋幾種濾波器比較：Butterworth濾波器是Chebyshev濾波器的特殊形式Chebyshev濾波器是Elliptic濾波器的特殊形式Elliptic濾波器性能最好48 UN(x)為N階雅可比橢圓函數(shù)。2221()1()aNpHjU橢圓低通濾波器幅度平方函數(shù)：p給定，兩個參數(shù)和N。的求法與切比雪夫濾波器相似，N的求法如下： 2121210.12

24、2011,1sin101spasK k KkNK kKkdkkK xx49Observations：橢圓濾波器通帶和阻帶波紋幅度固定時，階數(shù)越高，過渡帶越窄；橢圓濾波器階數(shù)固定時，通帶和阻帶波紋幅度越小，過渡帶就越寬General Conclusions：階數(shù)越高，性能越好，系統(tǒng)越復雜；波紋幅度和過渡帶寬度是一對矛盾50Matlab實現(xiàn)常用函數(shù)：ellipap, ellipord, ellip%Elliptic_Lowwp=2*pi*3000;ws=2*pi*12000;Ap=0.1;As=60;N,wp= ellipord(wp,ws,Ap,As,s);B,A= ellip(N, Ap,As

25、,wp,s);fk=0:14000/512:14000;wk=2*pi*fk;figure(1);freqs(B,A,wk);Hk=freqs(B,A,wk);figure(2);plot(fk/1000,20*log10(abs(Hk);xlabel(Frequency(kHz);ylabel(Magnitude(dB);axis(0,14,-80,5);51103104-200-1000100200Frequency (rad/s)Phase (degrees)10310410-610-410-2100Frequency (rad/s)Magnitude02468101214-80-70-

26、60-50-40-30-20-100Frequency(kHz)Magnitude(dB)52幾種濾波器比較：N, ap, as相同時的過渡帶寬度：Butterworth Chebyshev I = Chebyshev II EllipticButterworth：單調下降Chebyshev I/II：通帶/阻帶等波紋，過渡帶和阻帶/通帶單調下降Elliptic：通帶、阻帶都等波紋相位逼近情況： Butterworth和Chebyshev在大約3/4通帶接近線性相位特性；Elliptic在大約1/2通帶接近線性相位特性復雜性：在滿足相同幅頻響應指標下的階數(shù)N，Butterworth Cheby

27、shev I/II Elliptic535. 頻率變換與模擬高通、帶通和帶阻濾波器設計目標：設計指標為as、ap的模擬高通、帶通或帶阻濾波器Hd(s)54步驟：歸一化低通濾波器G (p)的截止頻率根據(jù)as,ap設計G (p)得到Hd(s)Hd(s)的截止頻率頻率變換公式Butterworth, Chebyshev, Elliptic頻率變換公式：p = F(s)1( )( )( )( )( )( )dp F sds FpHsG pG pHs55(1) 低通到高通的頻率變換公式phps(5.2.50) phphjj 復頻域： 頻域： (5.2.51) 高通濾波器通帶截止頻率 歸一化低通濾波器G(

28、p)頻率高通濾波器Hd(s)頻率通帶0,1ph, +阻帶s,+0,sh高通濾波器阻帶截止頻率 ph 由對稱性可得56低通與高通濾波器的幅度特性57例： 設計高通濾波器, ph =200Hz, sh=100Hz，幅度特性單調下降，ph處最大衰減為3dB，阻帶最小衰減as=15dB。 解 （1）高通技術要求： ph=200Hz, ap=3dB; sh=100Hz, as=15dB （2）頻率轉換： p=1, 2phssh58 （3）設計歸一化低通G(p)。采用巴特沃斯濾波器，故0.10.1321010.18,2101lg0.182.473lg21( )221psspNG pppp （4）求模擬高通

29、Hd(s)：33223( )( )22phdpphphphssHsG psss 59Matlab實現(xiàn)常用函數(shù)：lp2hp%Butt_LP2HPwp=1;ws=2;Ap=3;As=15;N,wc=buttord(wp,ws,Ap,As,s);B,A=butter(N,wc,s);wph=2*pi*200;BH,AH=lp2hp(B,A,wph);figure(1); freqs(B,A);figure(2);freqs(BH,AH);60102103-200-1000100200Frequency (rad/s)Phase (degrees)10210310-610-410-2100Freque

30、ncy (rad/s)Magnitude10-1100101-200-1000100200Frequency (rad/s)Phase (degrees)10-110010110-310-210-1100Frequency (rad/s)Magnitude61其實可以直接來%Butt_HPwp=2*pi*200;ws=2*pi*100;Ap=0.1;As=40;N,wc=buttord(wp,ws,Ap,As,s);B,A=butter(N,wc,high,s);freqs(B,A);102103-200-1000100200Frequency (rad/s)Phase (degrees)10

31、210310-610-410-2100Frequency (rad/s)Magnitude62(2) 低通到帶通的頻率變換公式220wspB s(5.2.53) 220wB 復頻域： 頻域： (5.2.54) 帶通濾波器通帶的高、低截止頻率 歸一化低通濾波器G(p)頻率高通濾波器Hd(s)頻率通帶0,1-1,0pl,00,pu阻帶s,+-,-s0,slsu,0,puplwpuplB 對數(shù)意義下的中心頻率63低通與帶通濾波器的幅度特性 64例：設計模擬帶通濾波器，通帶帶寬Bw=2200rad/s，中心頻率0=21000rad/s，通帶內最大衰減ap=3dB，阻帶sl=2830rad/s, su=

32、21200rad/s，阻帶最小衰減as=15dB。 解 (1) 模擬帶通的技術要求： 0=21000rad/s, ap=3dB sl=2830rad/s, su=21200rad/s as=15dB, Bw=2200rad/s; (2) 頻率轉換： p=1,2202201.8741.8331.833slwslssuwsuBB 65 （3）設計模擬歸一化低通濾波器G(p)，采用巴特沃斯型，有0.10.1321010.18,1.833101lg0.182.833lg1.8331( )221psspNG pppp （4）求模擬帶通Hd(s)： 2()( )( )plpupupldspsHsG p 66Matlab實現(xiàn)%Butt_BPwp=2*pi*905,1105;%905*1105 approximates 1000*1000w