【文章內(nèi)容簡介】 就可以把它們變成數(shù)字濾波器了。這些變換均是復(fù)值映射，許多文獻(xiàn)對此都有研究，根據(jù)數(shù)字濾波器所保持的模擬濾波器的不同特性，研究出不同的變換技術(shù)。其中，最重要的有兩種：脈沖響應(yīng)不變法（保持脈沖響應(yīng)不變，又叫沖激響應(yīng)不變法）和雙線性Z變換法（保持系統(tǒng)函數(shù)不變）。 脈沖響應(yīng)不變法脈沖響應(yīng)不變法的設(shè)計(jì)原理是使數(shù)字濾波器的單位抽樣響應(yīng)序列h(n)，模仿模擬濾波器的脈沖響應(yīng)g(t)。設(shè)系統(tǒng)傳遞函數(shù)為G(s)的模擬濾波器的單位脈沖響應(yīng)為g(t),并將脈沖響應(yīng)g(t)進(jìn)行等間隔采樣，使得數(shù)字濾波器的單位抽樣響應(yīng)h(n)剛好等于g(t)的采樣值，即： (26)其中的為采樣周期。G(s)是模擬濾波器的系統(tǒng)傳遞函數(shù)，又令H(z)是數(shù)字濾波器的系統(tǒng)傳遞函數(shù)。采樣信號的拉式變換與相應(yīng)的采樣序列Z變換的映射關(guān)系為： (27)所以系統(tǒng)函數(shù)G(s)和H(z)的關(guān)系為： (28)式(28)的物理意義為首先將模擬濾波器的系統(tǒng)函數(shù)G(s)作周期的延拓，在經(jīng)過式(27)的映射變換，映射到Z平面上，從而得到數(shù)字濾波器的系統(tǒng)函數(shù)H(z)。且模擬和數(shù)字頻率滿足下列關(guān)系：ω=ΩT。經(jīng)過式(27)的映射，s平面的左半平面映射為Z平面的單位圓內(nèi)，因此，一個因果的和穩(wěn)定的模擬濾波器映射成因果的和穩(wěn)定的數(shù)字濾波器。經(jīng)過以上分析，按照脈沖響應(yīng)不變法，通過模擬濾波器的系統(tǒng)傳遞函數(shù)G(s)，可直接求得數(shù)字濾波器系統(tǒng)函數(shù)H(Z),其設(shè)計(jì)具體步驟歸納如下：(1)利用ω=ΩT（可由關(guān)系式推出），將數(shù)字濾波器指標(biāo)，轉(zhuǎn)換為模擬濾波器指標(biāo)，(2)根據(jù)指標(biāo)，來設(shè)計(jì)模擬濾波器G(s)(3)利用部分分式展開法，把G(s)展成 (29)(4)最后把模擬極點(diǎn)轉(zhuǎn)換為數(shù)字極點(diǎn)，得到數(shù)字濾波器： (210)根據(jù)上述理論，將舉例在MATLAB環(huán)境下用函數(shù)實(shí)現(xiàn)脈沖響應(yīng)不變法設(shè)計(jì)一數(shù)字低通濾波器。其函數(shù)為[b,a]=impinvar(c,d,T)，其中，b表示數(shù)字濾波器自變量為的分子多項(xiàng)式，a表示數(shù)字濾波器自變量為的分母多項(xiàng)式，c表示模擬濾波器自變量為s的分子多項(xiàng)式,d表示模擬濾波器自變量為s的分母多項(xiàng)式，T表示采樣變換參數(shù)。%設(shè)模擬低通巴特沃斯濾波器，通帶紋波為Rp=1dB,通帶上限角頻率ωp=,阻帶下限角頻率ωs=，帶阻最小衰減αs=15dB，根據(jù)該低通模擬濾波器，利用沖激響應(yīng)不變法設(shè)計(jì)相應(yīng)的數(shù)字低通濾波器wp=*pi。ws=*pi。Rp=1。As=15。T=1。 %采樣周期%性能指標(biāo)Rip=10^(Rp/20)。Atn=10^(As/20)。OmgP=wp*T。OmgS=ws*T。[N,OmgC]=buttord(OmgP,OmgS,Rp,As,39。s39。)。 %選取模擬濾波器的階數(shù)[cs,ds]=butter(N,OmgC,39。s39。)。 %設(shè)計(jì)出所需的模擬低通濾波器[b,a]=impinvar(cs,ds,T)。 %應(yīng)用脈沖響應(yīng)不變法進(jìn)行轉(zhuǎn)換 [db,mag,pha,grd,w]=freqz_m(b,a)。 %求得相對，絕對頻響及相位、群延遲響應(yīng)%繪制各條曲線subplot(2,2,1)。plot(w/pi,mag)。title(39。幅頻特性39。)。xlabel(39。w(/pi)39。)。ylabel(39。|H(jw)|39。)。axis([0,1,0,])。set(gca,39。XTickMode39。,39。manual39。,39。XTick39。,[0 1])。set(gca,39。YTickMode39。,39。manual39。,39。YTick39。,[0 Atn Rip 1])。gridsubplot(2,2,2)。plot(w/pi,db)。title(39。幅頻特性(dB)39。)。xlabel(39。w(/pi)39。)。ylabel(39。dB39。)。axis([0,1,40,5])。set(gca,39。XTickMode39。,39。manual39。,39。XTick39。,[0 1])。set(gca,39。YTickMode39。,39。manual39。,39。YTick39。,[40 As Rp 0])。gridsubplot(2,2,3)。plot(w/pi,pha/pi)。title(39。相頻特性39。)。xlabel(39。w(/pi)39。)。ylabel(39。pha(/pi)39。)。axis([0,1,1,1])。set(gca,39。XTickMode39。,39。manual39。,39。XTick39。,[0 1])。gridsubplot(2,2,4)。plot(w/pi,grd)。title(39。群延遲39。)。xlabel(39。w(/pi)39。)。ylabel(39。Sample39。)。axis([0,1,0,12])。set(gca,39。XTickMode39。,39。manual39。,39。XTick39。,[0 1])。grid%本設(shè)計(jì)用到的M文件函數(shù)function[db,mag,pha,grd,w]=freqz_m(b,a)%濾波器幅值響應(yīng)（絕對、相對）、相位響應(yīng)及群延遲%Usage:[db,mag,pha,grd,w]=freqz_m(b,a) %500點(diǎn)對應(yīng)[0,pi]%db 相對幅值響應(yīng)； mag 絕對幅值響應(yīng)； pha 相位響應(yīng)； grd 群延遲響應(yīng)；%w 采樣頻率； b 系統(tǒng)函數(shù)H(z)的分子項(xiàng)%a 系統(tǒng)函數(shù)H(z)的分母項(xiàng)[H,w]=freqz(b,a,500)。 %500點(diǎn)的幅頻響應(yīng)mag=abs(H)。db=20*log10(mag/max(mag))。pha=angle(H)。grd=grpdelay(b,a,w)。圖27 脈沖響應(yīng)不變法設(shè)計(jì)數(shù)字低通濾波器的特性(基于MATLAB實(shí)現(xiàn))總結(jié)以上，脈沖響應(yīng)不變法的優(yōu)點(diǎn)是頻率坐標(biāo)變換是線性的，即ω=ΩT，如不考慮頻率混疊現(xiàn)象，用這種方法設(shè)計(jì)數(shù)字濾波器會很好的重現(xiàn)原模擬濾波器的頻率響應(yīng)。另外一個優(yōu)點(diǎn)是數(shù)字濾波器的單位脈沖響應(yīng)完全模仿模擬濾波器的單位沖激響應(yīng)，時(shí)域逼近好。但其也具有很大的缺點(diǎn)，若抽樣頻率不高或其它原因?qū)a(chǎn)生混疊失真，不能重現(xiàn)原模擬濾波器頻率響應(yīng)。所以，脈沖響應(yīng)不變法適合低通、帶通濾波器設(shè)計(jì)，不適合高通、帶阻濾波器的設(shè)計(jì)。 雙線性Z變換法利用脈沖響應(yīng)不變法設(shè)計(jì)數(shù)字濾波器時(shí)，由于ω=ΩT的頻率關(guān)系是根據(jù)推導(dǎo)的，所以是軸每隔2π/T便映射到單位圓上一周，引起了頻域混疊。為克服這一現(xiàn)象，人們找到了另一種映射關(guān)系： (211)此關(guān)系稱為雙線性Z變換法。雙線性Z變換法的基本思路是：首先將整個s平面壓縮到s1平面的一條帶寬為2π/T（叢π/T到π/T）的橫帶里，然后通過標(biāo)準(zhǔn)的變化關(guān)系將橫帶變換成整個Z平面上去，這樣就得到s平面與Z平面間的一一對應(yīng)的單值關(guān)系，整個過程如圖28所示： jΩ jΩ1 jIm(Z) π/T 0 б 0 б 0 1 бπ/T s平面 s1平面 Z平面圖28 雙線性Z變換法的映射關(guān)系由式(211)得 (212)及 (213) (214)式(211)及式(212)給出了s和z之間的映射關(guān)系，而式(213) 和式(214)給出了Ω和ω之間的映射關(guān)系，但這是一種非映射關(guān)系，雙線性Z變換法正是利用了正切函數(shù)的非線性特點(diǎn)，把整個jΩ軸壓縮到了單位圓的一周上。在MATLAB中，雙線性Z變換可通過bilinear函數(shù)實(shí)現(xiàn),其調(diào)用格式為：[Bz,Az]=bilinear(B,A,Fs)，其中B,A為模擬濾波器的傳遞函數(shù)G(s)的分子分母多項(xiàng)式系數(shù)分量，而Bz,Az為數(shù)字濾波器的傳遞函數(shù)H(Z)的分子分母多項(xiàng)式的系數(shù)分量。%使用雙線性Z變換設(shè)計(jì)一低通數(shù)字濾波器，fp=100HZ,fs=300HZ，αp=3dB,αs=20dB,抽樣頻率為Fs=1000HZfp=100。fs=300。Fs=1000。rp=3。rs=20。 wp=2*pi*fp/Fs。ws=2*pi*fs/Fs。Fs=Fs/Fs %使Fs為一wap=tan(wp/2)。was=tan(ws/2)。[n,wn]=buttord(wap,was,rp,rs,39。s39。)[z,p,k]=buttap(n)。[bp,ap]=zp2tf(z,p,k)[bs,as]=lp2lp(bp,ap,wap)[bz,az]=bilinear(bs,as,Fs/2)[H,w]=freqz(bz,az,256,Fs*1000)。圖29 雙線性Z變換設(shè)計(jì)數(shù)字低通濾波器的特性(基于MATLAB實(shí)現(xiàn)) 頻帶變換及其MATLAB實(shí)現(xiàn)除了低通數(shù)字濾波器之外，實(shí)際中還常常需要高通、帶通、及帶阻數(shù)字濾波器。設(shè)計(jì)這三種濾波器的方法有很多，例如基于模擬濾波器轉(zhuǎn)換法的IIR數(shù)字濾波器設(shè)計(jì)，基于直接數(shù)字域法的IIR數(shù)字濾波器設(shè)計(jì)等等。下面將具體介紹基于濾波器轉(zhuǎn)換法的IIR數(shù)字濾波器設(shè)計(jì)—直接由模擬低通濾波器轉(zhuǎn)換成數(shù)字高通、帶通或帶阻濾波器。 模擬低通濾波器轉(zhuǎn)換成數(shù)字高通濾波器若已知模擬低通濾波器的系統(tǒng)傳遞函數(shù)為G(s)，則模擬低通濾波器s平面到數(shù)字高通濾波器的Z平面的變換公式為： (215)頻率變換公式為： (216)從而得到數(shù)字高通濾波器的傳遞函數(shù)表達(dá)式為： (217)接下來，基于MATLAB利用巴特沃斯模擬濾波器舉例，設(shè)計(jì)一數(shù)字高通濾波器，,通帶內(nèi)衰減不大于1dB,，阻帶內(nèi)衰減不小于15dB,采樣周期為Ts=1。clear all。Wp=*pi。Ws=*pi。Ap=1。As=15。[N,wn]=buttord(Wp/pi,Ws/pi,Ap,As) %計(jì)算巴特沃斯濾波器階次和截至頻率[b,a]=butter(N,wn,39。high39。)。 %頻率變換法設(shè)計(jì)巴特沃斯高通濾波器[b0,B,A]=dir2cas(b,a) %數(shù)字高通濾波器級聯(lián)型[db,mag,pha,grd,w]=freqz_m(b,a) %數(shù)字濾波器響應(yīng)subplot(211)。plot(w/pi,mag)。title(39。高通數(shù)字濾波器幅頻響|H(ejOmega)|39。)subplot(212)。plot(w/pi,db)。title(39。高通數(shù)字濾波器幅頻響(db)39。)其MATLAB運(yùn)行結(jié)果為：b0 =B = A = 所以由巴特沃斯模擬濾波器設(shè)計(jì)的數(shù)字高通濾波器的傳遞函數(shù)表達(dá)式為：圖210 基于巴特沃斯模擬濾波器設(shè)計(jì)的數(shù)字高通濾波器的幅頻響應(yīng)曲線 模擬低通濾波器轉(zhuǎn)換成數(shù)字帶通濾波器若已知模擬低通濾波器的系統(tǒng)傳遞函數(shù)為G(s)，則模擬低通濾波器s平面到數(shù)字帶通濾波器的Z平面的變換公式為： (218)頻率變換公式為： (219)其中： (220)式中為下截止頻率，為上截止頻率。從而可得數(shù)字帶通濾波器的傳遞函數(shù)表達(dá)式為： (221)接下來，基于MATLAB利用巴特沃斯模擬濾波器舉例，設(shè)計(jì)一數(shù)字帶通濾波器，,帶通內(nèi)衰減不大于3dB，，阻帶內(nèi)衰減不小于18dB。clear all。Wp=[*pi,*pi]。Ws=[*pi,*pi]。Ap=3。As=18。[N,wn]=buttord(Wp/pi,Ws/pi,Ap,As)。 %計(jì)算巴特沃斯濾波器階次和截至頻率[b,a]=butter(N,wn,39。bandpass39。)。 %頻率變換法設(shè)計(jì)巴特沃斯帶通濾波器[b0,B,A]=dir2cas(b,a) %數(shù)字帶通濾波器級聯(lián)型[db,mag,pha,grd,w]=freqz_m(b,a)。 %數(shù)字濾波器響應(yīng)subplot(211)。plot(w/pi,mag)。title(39。帶通數(shù)字濾波器幅頻響|H(ejOmega)|39。)subplot(212)。plot(w/pi,db)。title(39。帶通數(shù)字濾波器幅頻響(db)39。)其MATLAB運(yùn)行結(jié)果為：b0 =B = A = 所以由巴