1、經(jīng)典數(shù)字濾波器及其設(shè)計4.1 引 言數(shù)字濾波器(Digital Filter, DF在信號的過濾、檢測和參數(shù)估計等方面起著重要的作用。信號往往夾雜著噪聲及無用信號成分，必須將這些干擾成分濾除。數(shù)字濾波器對信號進(jìn)行篩選，可通過特定頻段的信號。一般來說，噪聲信號往往是高頻信號，而經(jīng)典濾波器正是假定有用信號與噪聲信號具有不同的頻段，所以利用經(jīng)典濾波器可以去除噪聲。但如果有用信號和無用信號，或有用信號和噪聲的頻譜相互重疊，那么經(jīng)典濾波器則不能實現(xiàn)理想的濾波性能。現(xiàn)代濾波器的作用是從含有噪聲的信號中估計出信號的某些特征或信號本身，一旦信號被估計出，那么估計出來的信號與原信號相比會有更高的信噪比。這類濾波

2、器主要有維納濾波器、卡爾曼濾波器、線性預(yù)測濾波器以及自適應(yīng)濾波器等。4.2 經(jīng)典數(shù)字濾波器原理數(shù)字濾波是數(shù)字信號分析中最重要的組成部分之一，與模擬濾波相比，它具有精度和穩(wěn)定性高、系統(tǒng)函數(shù)容易改變、靈活性強(qiáng)、便于大規(guī)模集成和可實現(xiàn)多維濾波等優(yōu)點(diǎn)。在信號的過濾、檢測和參數(shù)的估計等方面，經(jīng)典數(shù)字濾波器是使用最廣泛的一種線性系統(tǒng)。數(shù)字濾波器的作用是利用離散時間系統(tǒng)的特性對輸入信號波形(或頻譜 進(jìn)行加工處理，或者說利用數(shù)字方法按預(yù)定的要求對信號進(jìn)行變換。若濾波器的輸入、輸出都是離散時間信號，那么該濾波器的單位沖激響應(yīng)h(n也必然是離散的，這種濾波器稱為數(shù)字濾波器。當(dāng)用硬件實現(xiàn)一個DF 時，所需的元件是乘

3、法器、延時器和相加器；而用MA TLAB 軟件實現(xiàn)時，它僅僅需要線性卷積程序就可以實現(xiàn)。眾所周知，模擬濾波器(Analog Filter ，AF 只能用硬件來實現(xiàn)，其元件有電阻R ，電感L ，電容C 及運(yùn)算放大器等。因此，DF 的實現(xiàn)要比AF 容易得多，并且更容易獲得較理想的濾波性能。 數(shù)字濾波器的作用是對輸入信號進(jìn)行濾波，就如同信號通過系統(tǒng)一樣。對于線性時不變系統(tǒng)，其時域輸入輸出關(guān)系是：(4-1若y(n、x(n的傅里葉變化存在，則輸入輸出的頻域關(guān)系是：(4-2當(dāng)輸入信號x(n通過濾波器h(n后，其輸出y(n中不再含有的頻率成分，僅使的信號成分通過，其中是濾波器的轉(zhuǎn)折頻率。經(jīng)典數(shù)字濾波器按照單

4、位取樣響應(yīng)h(n的時域特性可分為無限沖激響應(yīng)(IIR，Infinite Impulse Response 系統(tǒng)和有限沖激響應(yīng)(FIR，F(xiàn)inite Impulse Response 系統(tǒng)。如果單位取樣響應(yīng)是時寬無限的h(n，則稱之為IIR 系統(tǒng)；而如果單位取樣響應(yīng)是時寬有限的h(n，則稱之為FIR 系統(tǒng)。數(shù)字濾波器按照實現(xiàn)的方法和結(jié)構(gòu)形式分為遞歸型或非遞歸型兩類。遞歸型數(shù)字濾波器的當(dāng)前輸出y(n是輸入x(n的當(dāng)前值和以前各輸入值x(n，x(n1 ，. ，及以前各輸出值y(n，y(n1 ，. 的函數(shù)。一個N 階遞歸型數(shù)字濾波器(IIR濾波器 的差分方程為：(4-3其中，式(4-3中的系數(shù)至少有一

5、項不為零。說明必須將延時的輸出序列進(jìn)行反饋。 遞歸系統(tǒng)的傳統(tǒng)函數(shù)定義為：(4-4遞歸系統(tǒng)的傳遞函數(shù)H(z在Z 平面上不僅有零點(diǎn)，而且有極點(diǎn)。非遞歸型數(shù)字濾波器當(dāng)前的輸出值y(n僅為當(dāng)前及以前的輸入序列的函數(shù)，而與以前的各個輸出值無關(guān)，因此從結(jié)構(gòu)上看非遞歸系統(tǒng)沒有反饋環(huán)路。一個N 階的非遞歸型數(shù)字濾波器(FIR濾波器 的差分方程為：(4-5差分方程式(4-5中的系數(shù)等于單位取樣響應(yīng)的序列值h(n，其系統(tǒng)函數(shù)H(z可以表示為以下形式：(4-6H(z是的多項式，因此它的極點(diǎn)只能在Z 平面的原點(diǎn)上。這兩類濾波器無論是在性能上還是在設(shè)計方法上都有著很大的區(qū)別。FIR 濾波器可以對給定的頻率特性直接進(jìn)行設(shè)

6、計，而IIR 濾波器目前最通用的方法是利用已經(jīng)很成熟的模擬濾波器的設(shè)計方法來進(jìn)行設(shè)計。1. 特征頻率濾波器的頻率參數(shù)主要有：通帶截頻為通帶與過渡帶的邊界點(diǎn)，在該點(diǎn)信號增益下降到規(guī)定的下限。阻帶截頻為阻帶與過渡帶的邊界點(diǎn)，在該點(diǎn)信號衰耗下降到規(guī)定的下限。轉(zhuǎn)折頻率為信號功率衰減到1/2(約3dB 時的頻率，在很多情況下，也常以fc 作為通帶或阻帶截頻。當(dāng)電路沒有損耗時，固有頻率，就是其諧振頻率，復(fù)雜電路往往有多個固有頻率。2. 增益與衰耗濾波器在通帶內(nèi)的增益并非常數(shù)。對低通濾波器通帶增益一般指=0時的增益；高通指時的增益；帶通則指中心頻率處的增益。對帶阻濾波器，應(yīng)給出阻帶衰耗，衰耗定義為增益的倒數(shù)

7、。通帶增益變化量指通帶內(nèi)各點(diǎn)增益的最大變化量，如果以dB 為單位，則指增益dB 值的變化量。3. 阻尼系數(shù)與品質(zhì)因數(shù)阻尼系數(shù)是表征濾波器對角頻率為0信號的阻尼作用，是濾波器中表示能量衰耗的一項指標(biāo)，它是與傳遞函數(shù)的極點(diǎn)實部大小相關(guān)的一項系數(shù)。它可由傳遞函數(shù)的分母多項式系數(shù)求得：(4-7式中，表示傳遞函數(shù)的分母多項式系數(shù)。的倒數(shù)Q=1/稱為品質(zhì)因數(shù)，是評價帶通與帶阻濾波器頻率選擇特性的一個重要指標(biāo)，Q 為：(4-8式中的為帶通或帶阻濾波器的3dB 帶寬，為中心頻率，在很多情況下中心頻率與固有頻率相等。4. 靈敏度濾波電路由許多元件構(gòu)成，每個元件參數(shù)值的變化都會影響濾波器的性能。濾波器某一性能指標(biāo)

8、對某一元件參數(shù)變化的靈敏度記作，定義為：(4-9靈敏度是濾波電路設(shè)計中的一個重要參數(shù)，可以用來分析元件實際值偏離設(shè)計值時，電路實際性能與設(shè)計性能的偏差程度；也可以用來估計在使用過程中元件參數(shù)值變化時，電路性能變化情況。該靈敏度與測量儀器或電路系統(tǒng)靈敏度概念不同，該靈敏度越小，標(biāo)志著電路容錯能力越強(qiáng)，穩(wěn)定性也越高。5. 群時延函數(shù)當(dāng)濾波器幅頻特性滿足設(shè)計要求時，為保證輸出信號失真度不超過允許范圍，對其相頻特性也應(yīng)提出一定要求。在濾波器設(shè)計中，常用群時延函數(shù)評價信號經(jīng)濾波后相位失真程度。越接近常數(shù)，信號相位失真越小。6. 濾波器設(shè)計的步驟不論是IIR 濾波器還是FIR 濾波器的設(shè)計都包括三個步驟：

9、(1 按照實際任務(wù)的要求，確定濾波器的性能指標(biāo)。(2 用一個因果、穩(wěn)定的離散線性時不變系統(tǒng)的系統(tǒng)函數(shù)去逼近這一性能指標(biāo)。根據(jù)不同的要求可以用IIR 系統(tǒng)函數(shù)，也可以用FIR 系統(tǒng)函數(shù)去逼近。(3 利用有限精度算法實現(xiàn)系統(tǒng)函數(shù)，包括結(jié)構(gòu)選擇、字長選擇等。IIR 數(shù)字濾波器設(shè)計的最通用的方法是借助于模擬濾波器的設(shè)計方法。模擬濾波器設(shè)計已經(jīng)有一套相當(dāng)成熟的方法，它不但有完整的設(shè)計公式，而且還有較為完整的圖表以供查詢，因此充分利用這些已有的資源將會給數(shù)字濾波器的設(shè)計帶來很大的方便。IIR 數(shù)字濾波器的設(shè)計步驟是：(1 按一定的規(guī)則將給出的數(shù)字濾波器的技術(shù)指標(biāo)轉(zhuǎn)換為模擬低通濾波器的技術(shù)指標(biāo)。(2 根據(jù)轉(zhuǎn)

10、換后的技術(shù)指標(biāo)設(shè)計模擬低通濾波器G(s。(3 再按一定的規(guī)則將G(s轉(zhuǎn)換成H(s。若所設(shè)計的數(shù)字濾波器是低通的，那么上述設(shè)計工作可以結(jié)束，若所設(shè)計是高通、帶通或帶阻濾波器，那么還需進(jìn)行以下步驟。(4 將高通、帶通或帶阻數(shù)字濾波器的技術(shù)指標(biāo)先轉(zhuǎn)換為低通模擬濾波器的技術(shù)指標(biāo)，然后按照上述步驟設(shè)計出低通G(s，再將G(s轉(zhuǎn)換為所需的H(z。4.3 IIR 數(shù)字濾波器的結(jié)構(gòu)與設(shè)計IIR 濾波器系統(tǒng)函數(shù)的極點(diǎn)可以在單位圓內(nèi)的任何位置，實現(xiàn)IIR 濾波器的階次較低，所用的存儲單元較少，效率高，又由于IIR 數(shù)字濾波器能夠保留一些模擬濾波器的優(yōu)良特性，因此應(yīng)用很廣。設(shè)計IIR 數(shù)字濾波器的方法主要有基于沖激

11、響應(yīng)不變法的IIR 數(shù)字濾波器設(shè)計，基于雙線性Z 變換法的IIR 數(shù)字濾波器設(shè)計，數(shù)字高通、帶通及帶阻IIR 濾波器設(shè)計，基于MATLAB 函數(shù)直接設(shè)計IIR 數(shù)字濾波器。沖激響應(yīng)不變法的設(shè)計原理是利用數(shù)字濾波器的單位抽樣響應(yīng)序列H(z來逼近模擬濾波器的沖激響應(yīng)g(t。按照沖激響應(yīng)不變法的原理，通過模擬濾波器的系統(tǒng)傳遞函數(shù)G(s，可以直接求得數(shù)字濾波器的系統(tǒng)函數(shù)H(z，其轉(zhuǎn)換步驟如下：1 利用=T(可由關(guān)系式推導(dǎo)出 ，將，轉(zhuǎn)換成，而，不變；2 求解低通模擬濾波器的傳遞函數(shù)G(s；3 將模擬濾波器的傳遞函數(shù)G(s轉(zhuǎn)換為數(shù)字濾波器的傳遞函數(shù)H(z。盡管通過沖激響應(yīng)不變法求取數(shù)字濾波器的系統(tǒng)傳遞函數(shù)

12、比較方便，并具有良好的時域逼近特性，但若G(s不是帶限的，或是抽樣頻率不高，那么在中將發(fā)生混疊失真，數(shù)字濾波器的頻率響應(yīng)不能重現(xiàn)模擬濾波器的頻率響應(yīng)。只有當(dāng)模擬濾波器的頻率響應(yīng)在超過折疊頻率后的衰減很大時，混疊失真才很小，此時采樣脈沖響應(yīng)不變法設(shè)計的數(shù)字濾波器才能滿足設(shè)計的要求，這是沖激響應(yīng)不變法的一個嚴(yán)重的缺點(diǎn)。【例4-1】設(shè)模擬低通巴特沃斯濾波器，通帶紋波為1dB ，通帶上限角頻率=0.2，阻帶下限角頻率=0.3，阻帶最小衰減=15dB，根據(jù)該低通模擬濾波器，利用沖激響應(yīng)不變法設(shè)計響應(yīng)的數(shù)字低通濾波器，并繪出設(shè)計后的數(shù)字濾波器的特性曲線。實現(xiàn)例4-1的MA TLAB 程序如例程4-1所示，

13、程序的運(yùn)行結(jié)果如圖4-1所示。例程4-1 利用沖激響應(yīng)不變法設(shè)計數(shù)字低通濾波器%利用模擬巴特沃斯濾波器設(shè)計數(shù)字低通濾波器%沖激響應(yīng)不變法wp=0.2*pi;ws=0.3*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,'s' %選取模擬濾波器的階數(shù)cs,ds=butter(N,OmgC,'s' %設(shè)計出所需的模擬低通濾波器b,a=impinvar(cs,ds,T; %應(yīng)用脈沖響應(yīng)不變法進(jìn)行轉(zhuǎn)換%求得相對、絕對

14、頻響及相位、群遲延響應(yīng)db,mag,pha,grd,w=freqz_m(b,a;%下面繪出各條曲線subplot(2,2,1;plot(w/pi,mag;title('幅頻特性'xlabel('w(/pi'ylabel('|H(jw|'axis(0,1,0,1.1;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;set(gca,'YTickMode','manual','YTick',0 Atn

15、 Rip 1;gridsubplot(2,2,2;plot(w/pi,db;title('幅頻特性(db'xlabel('w(/pi'ylabel('dB'axis(0,1,-40,5;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;set(gca,'YTickMode','manual','YTick',-40 -As -Rp 0;gridsubplot(2,2,3;plot(w/pi,pha

16、/pi;title('相頻特性'xlabel('w(/pi'ylabel('pha(/pi'axis(0,1,-1,1;set(gca,'XTickMode','manual','XTick',0 0.2 0.3 0.5 1;gridsubplot(2,2,4;plot(w/pi,grd;title('群延時'xlabel('w(/pi'ylabel('Sample'axis(0,1,0,12;set(gca,'XTickMode',&

17、#39;manual','XTick',0 0.2 0.3 0.5 1;grid%所使用的M 文件函數(shù)functiondb,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 %w 采樣頻率； b 系統(tǒng)函數(shù)H(z的分子項(對FIR ，b=h %a 系統(tǒng)函數(shù)H(z的分母項(對FIR ，a=1H,w=freqz(b,a,500; %500點(diǎn)的復(fù)頻響應(yīng) mag=abs(H; %絕對幅值響應(yīng) db=20*log10(mag/max(ma

18、g; %相對幅值響應(yīng) pha=angle(H; %相位響應(yīng)grd=grpdelay(b,a,w; %群延遲響應(yīng)圖4-1 利用沖激響應(yīng)不變法設(shè)計的數(shù)字低通濾波器的特性曲線由于的頻率映射關(guān)系是根據(jù)推導(dǎo)的，所以使j軸每隔2/Ts便映射到單位圓上一周，利用沖激響應(yīng)不變法設(shè)計數(shù)字濾波器時可能會導(dǎo)致上述的頻域混疊現(xiàn)象。為了克服這一問題，需要找到由s 平面到z 平面的另外的映射關(guān)系，這種關(guān)系應(yīng)保證：1 s平面的整個j軸僅映射為z 平面單位圓上的一周；2 若G(s是穩(wěn)定的，由G(s映射得到的H(z也應(yīng)該是穩(wěn)定的；3 這種映射是可逆的，既能由G(s得到H(z，也能由H(z得到G(s；4 如果G(j0=1，那么。

19、雙線性Z 變換滿足以上4個條件的映射關(guān)系，其變換公式為(4-10雙線性Z 變換的基本思路是：首先將整個s 平面壓縮到一條從/Ts到/Ts的帶寬為2/Ts的橫帶里，然后通過標(biāo)準(zhǔn)的變換關(guān)系將橫帶變換成整個z 平面上去，這樣就得到s 平面與z平面間的一一對應(yīng)的單值關(guān)系。在MA TLAB 中，雙線性Z 變換可以通過bilinear 函數(shù)實現(xiàn)，其調(diào)用格式為：Bz，Azbilinear(B，A ，F(xiàn)s ；其中B ，A 為模擬濾波器傳遞函數(shù)G(s的分子分母多項式的系數(shù)向量，而Bz ，Az 為數(shù)字濾波器的傳遞函數(shù)H(z的分子分母多項式的系數(shù)向量。【例4-2】使用雙線性Z 變換法設(shè)計一個低通數(shù)字濾波器，給定的數(shù)

20、字濾波器的技術(shù)指標(biāo)為fp 100Hz ，fs=300Hz，p=3dB，p=20dB，抽樣頻率為Fs=1000Hz。實現(xiàn)例4-2的MATLAB程序如例程4-2所示，程序的運(yùn)行結(jié)果如圖4-2所示。例程4-2 基于雙線性Z 變換法設(shè)計數(shù)字濾波器clear all;fp=100;fs=300;Fs=1000;rp=3;rs=20;%wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Fs=Fs/Fs;%Firstly to finish frequency prewarping;wap=tan(wp/2;was=tan(ws/2;n,wn=buttord(wap,was,rp,rs,'s&

21、#39;%Note:'s'!z,p,k=buttap(n;bp,ap=zp2tf(z,p,kbs,as=lp2lp(bp,ap,wap%Note:s=(2/Ts(z-1(z+1;Ts=1,that is 2fs=1,fs=0.5;bz,az=bilinear(bs,as,Fs/2h,w=freqz(bz,az,256,Fs*1000;plot(w,abs(h;grid on;【運(yùn)行結(jié)果】n =2wn =0.4363bp =0 0 1ap =1.0000 1.4142 1.0000bs =0.1056as =1.0000 0.4595 0.1056bz =0.0675 0.134

22、9 0.0675az =1.0000 -1.1430 0.4128圖4-2 基于雙線性Z 變換法設(shè)計的數(shù)字濾波器的特性曲線除了低通數(shù)字濾波器之外，實際中還常常需要高通、帶通及帶阻數(shù)字濾波器，這三種數(shù)字濾波器的設(shè)計步驟如下：步驟1：將數(shù)字濾波器H(z的技術(shù)指標(biāo)p和s，通過tan(/2轉(zhuǎn)變?yōu)槟M濾波器G(s的技術(shù)指標(biāo)p和s，作歸一化處理后，得到p1，ss/p；步驟2：化解為模擬原型濾波器G(s的技術(shù)指標(biāo)；步驟3：設(shè)計模擬原型濾波器G(p；步驟4：將G(p轉(zhuǎn)換為模擬濾波器的轉(zhuǎn)移函數(shù)G(s；步驟5：將G(s轉(zhuǎn)換成數(shù)字濾波器的轉(zhuǎn)移函數(shù)H(z，s=(z1(z+1。所謂原型濾波器是指歸一化的低通濾波器。本節(jié)

23、主要討論通過IIR 數(shù)字濾波器的原型轉(zhuǎn)換設(shè)計法和IIR 數(shù)字濾波器的直接設(shè)計方法來設(shè)計數(shù)字高通、帶通及帶阻濾波器其轉(zhuǎn)換方法主要有3種：一是直接由模擬低通濾波器轉(zhuǎn)換成數(shù)字高通、帶通或帶阻濾波器；二是先由模擬低通濾波器轉(zhuǎn)換成模擬高通、帶通或帶阻濾波器，然后再把它轉(zhuǎn)換成相應(yīng)的數(shù)字濾波器；三是將模擬低通濾波器先轉(zhuǎn)換成數(shù)字低通濾波器，再通過變量代換變換成高通、帶通或帶阻濾波器。1. 直接由模擬低通濾波器轉(zhuǎn)換成數(shù)字高通、帶通或帶阻濾波器(1 模擬低通濾波器轉(zhuǎn)換成數(shù)字高通濾波器【例4-3】 利用巴特沃斯模擬濾波器設(shè)計數(shù)字高通濾波器，要求通帶截止頻率為0.6，通帶內(nèi)衰減不大于1dB ，阻帶起始頻率為0.4，阻

24、帶內(nèi)衰減不小于15dB ，采樣周期Ts 1。該例的MA TLAB 實現(xiàn)如例程4-3所示，該數(shù)字濾波器的幅頻特性曲線如圖4-3所示。例程4-3 基于巴特沃斯模擬濾波器設(shè)計數(shù)字高通濾波器Wp=0.6*pi;Ws=0.4*pi;Ap=1;As=15;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %計算巴特沃斯濾波器階次和截止頻率b,a=butter(N,wn,'high' %頻率變換法設(shè)計巴特沃斯高通濾波器db,mag,pha,grd,w=freqz_m(b,a; %數(shù)字濾波器響應(yīng)Plot(w/pi,mag;Title('數(shù)字濾波器幅頻響應(yīng)|H(ejOmeg

25、a|'圖4-3 基于巴特沃斯模擬濾波器設(shè)計的數(shù)字高通濾波器幅頻響應(yīng)曲線(2 模擬低通濾波器轉(zhuǎn)換成數(shù)字帶通濾波器【例4-4】利用巴特沃斯模擬濾波器設(shè)計數(shù)字帶通濾波器，要求通帶上下截止頻率為0.4、0.3，通帶內(nèi)衰減不大于3dB ，阻帶上下起始頻率為0.5、0.2，阻帶內(nèi)衰減不小于18dB 。該例的MA TLAB 實現(xiàn)如例程4-4所示，該數(shù)字濾波器的幅頻特性曲線如圖4-4 所示。例程4-4 基于巴特沃斯模擬濾波器設(shè)計數(shù)字帶通濾波器Wp=0.3*pi,0.4*pi;Ws=0.2*pi,0.5*pi;Ap=3;As=18;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %計算巴

26、特沃斯濾波器階次和截止頻率b,a=butter(N,wn,'bandpass' %頻率變換法設(shè)計巴特沃斯帶通濾波器db,mag,pha,grd,w=freqz_m(b,a; %數(shù)字濾波器響應(yīng)Plot(w/pi,mag;Title('數(shù)字濾波器幅頻響應(yīng)|H(ejOmega|'圖4-4 基于巴特沃斯模擬濾波器設(shè)計的數(shù)字帶通濾波器的幅頻響應(yīng)曲線(3 模擬低通濾波器轉(zhuǎn)換成數(shù)字帶阻濾波器【例4-5】利用巴特沃斯模擬濾波器設(shè)計數(shù)字帶阻濾波器，要求通帶上下截止頻率為0.8和該例的MA TLAB 實現(xiàn)如例程4-5所示，該數(shù)字濾波器的幅頻特性曲線如圖4-5所示。例程4-5 基于巴

27、特沃斯模擬濾波器設(shè)計數(shù)字帶阻濾波器Wp=0.4*pi,0.7*pi;Ws=0.2*pi,0.8*pi;Ap=1;As=30;N,wn=buttord(Wp/pi,Ws/pi,Ap,As; %計算巴特沃斯濾波器階次和截止頻率b,a=butter(N,wn,'stop' %頻率變換法設(shè)計巴特沃斯帶通濾波器db,mag,pha,grd,w=freqz_m(b,a; %數(shù)字濾波器響應(yīng)plot(w/pi,mag;Title('數(shù)字濾波器幅頻響|H(ejOmega| '圖4-5 基于巴特沃斯模擬濾波器設(shè)計數(shù)字帶阻濾波器2. 先由模擬低通濾波器轉(zhuǎn)換成模擬高通、帶通或帶阻濾波器

28、，然后再把它轉(zhuǎn)換成相應(yīng)的數(shù)字濾波器首先由頻率變換將模擬低通原型濾波器變換成模擬高通、帶通或帶阻濾波器，然后根據(jù)雙線性Z 變換將模擬高通濾波器變換成相應(yīng)的高通、帶通或帶阻數(shù)字濾波器。其模擬低通原型濾波器變換成模擬高通、帶通或帶阻濾波器的關(guān)系式如表4-1所示。表4-1 模擬低通原型濾波器變換成模擬高通、帶通或帶阻濾波器頻率變換模擬濾波器類型頻率變換公式高通帶通帶阻3. 將模擬低通濾波器先轉(zhuǎn)換成數(shù)字低通濾波器，再通過變量代換變換成高通、帶通或帶阻濾波器該方法中的模擬低通濾波器先轉(zhuǎn)換成數(shù)字低通濾波器，可以利用沖激響應(yīng)不變法和雙線性變換法來實現(xiàn)，再通過變量代換變換成數(shù)字高通、帶通或帶阻濾波器。【例4-6

29、】試用雙線性變換法設(shè)計一個帶通橢圓濾波器使其幅頻特性逼近于一個具有如下技術(shù)指標(biāo)的模擬帶通濾波器：Wp1=10Hz,Wp2=20Hz,在通帶內(nèi)的最大衰減為0.5dB, 在阻帶內(nèi)的最大衰減為50dB ，抽樣頻率為100Hz 。MATLAB 實現(xiàn)如例程4-6所示，該數(shù)字濾波器的幅頻特性曲線如圖4-6所示。wp1=10;wp2=20;Fs=100;rp=0.5;rs=50;wp1=2*pi*wp1;wp2=2*pi*wp2;Bw=wp2-wp1;Wo=sqrt(wp2*wp1;z,p,k=ellipap(7,rp,rs;A,B,C,D=zp2ss(z,p,k;At,Bt,Ct,Dt=lp2bp(A,B

30、,C,D,Wo,Bw;At1,Bt1,Ct1,Dt1=bilinear(At,Bt,Ct,Dt,Fs;num,den=ss2tf(At1,Bt1,Ct1,Dt1;H,W=freqz(num,den;plot(W*Fs/(2*pi,abs(H;grid;xlabel('頻率/Hz'ylabel('幅值'圖4-6 橢圓濾波器的頻率響應(yīng)4.4 FIR 數(shù)字濾波器的結(jié)構(gòu)與設(shè)計由于IIR 數(shù)字濾波器能夠保留一些模擬濾波器的優(yōu)良特性，因此應(yīng)用很廣。但是這些特性是以犧牲線性相位頻率特性為代價的，即用Butterworth 、chelbchev 和橢圓法設(shè)計的數(shù)字濾波器逼近理想

31、的濾波器的幅度頻率特性，得到的濾波器往往是非線性的。在許多電子系統(tǒng)中，對幅度頻率特性和線性相位特性都有較高的要求，所以IIR 濾波器在這些系統(tǒng)中往往難以勝任。有限長單位沖激響應(yīng)(FIR數(shù)字濾波器具有以下優(yōu)良的特點(diǎn)： 可在設(shè)計任意幅度頻率特性濾波器的同時，保證精確、嚴(yán)格的線性相位特性。 FIR數(shù)字濾波器的單位沖激響應(yīng)h(n是有限長的，可以用一個固定的系統(tǒng)來實現(xiàn)，因而FIR數(shù)字濾波器可以做成因果穩(wěn)定系統(tǒng)。 允許設(shè)計多通帶(多阻帶 系統(tǒng)。根據(jù)FIR 濾波器的幅頻特性，可以將其分為以下4種情況，其中h(n為FIR 濾波器的沖激響應(yīng)：1. h(n為偶對稱，且N 為奇數(shù)(型濾波器型濾波器的幅頻函數(shù)可以通過

32、例程4-7的MA TLAB 語言來實現(xiàn)。例程4-7 型濾波器的幅頻函數(shù)function Hr,w,a,L=hr_type1(h;%計算所設(shè)計的1型濾波器的振幅響應(yīng)%Hr=振幅響應(yīng)%a=1型濾波器的系數(shù)%L=Hr的階次%h=1型濾波器的單位沖擊響應(yīng)M=length(h;L=(M-1/2;a=h(L+1 2*h(L:-1:1;n=0:1:L;w=0:1:500'*2*pi/500;Hr=cos(w*n*a'2. h(n為偶對稱，且N 為偶數(shù)(型濾波器型濾波器的幅頻函數(shù)可以通過例程4-8 的MA TLAB 語言來實現(xiàn)。例程4-8 型濾波器的幅頻函數(shù)function Hr,w,b,L=

33、hr_type2(h;%計算所設(shè)計的2型濾波器的振幅響應(yīng)%Hr=振幅響應(yīng)%b=2型濾波器的系數(shù)%L=Hr的階次%h=2型濾波器的單位沖擊響應(yīng)M=length(h;L=M/2;b= 2*h(L:-1:1;n=1:1:L;n=n-0.5;w=0:1:500'*2*pi/500;Hr=cos(w*n*b'3. h(n為奇對稱，且N 為奇數(shù)(型濾波器型濾波器的幅頻函數(shù)可以通過例程4-9的MA TLAB 語言來實現(xiàn)。例程4-9 型濾波器的幅頻函數(shù)function Hr,w,c,L=hr_type3(h;%計算所設(shè)計的3型濾波器的振幅響應(yīng)%Hr=振幅響應(yīng)%b=3型濾波器的系數(shù)%L=Hr的階

34、次%h=3型濾波器的單位沖擊響應(yīng)M=length(h;L=(M-1/2;c= 2*h(L+1:-1:1;n=0:1:L;w=0:1:500'*2*pi/500;Hr=sin(w*n*c'4. h(n為奇對稱，且N 為偶數(shù)(型濾波器型濾波器的幅頻函數(shù)可以通過例程4-10 的MATLAB 語言來實現(xiàn)。例程4-10 型濾波器的幅頻函數(shù)function Hr,w,d,L=hr_type4(h;%計算所設(shè)計的4型濾波器的振幅響應(yīng)%Hr=振幅響應(yīng)%b=4型濾波器的系數(shù)%L=Hr的階次%h=4型濾波器的單位沖擊響應(yīng)M=length(h;L=M/2;d= 2*h(L:-1:1;n=1:1:L;

35、n=n-0.5;w=0:1:500'*2*pi/500;Hr=sin(w*n*d'針對FIR 數(shù)字濾波器的結(jié)構(gòu)特點(diǎn)，目前主要采用窗函數(shù)法、頻率采樣法和最優(yōu)化法等三種方法設(shè)計FIR 數(shù)字濾波器。窗函數(shù)法就是設(shè)計FIR 數(shù)字濾波器的最簡單的方法。它在設(shè)計FIR 數(shù)字濾波器中有很重要的作用，正確地選擇窗函數(shù)可以提高設(shè)計數(shù)字濾波器的性能，或者在滿足設(shè)計要求的情況下，減小FIR 數(shù)字濾波器的階次。常用的窗函數(shù)有以下幾種：矩形窗(Rectangular window、三角窗(Triangular window 、漢寧窗(Hanning window 、海明窗(Hamming window

36、、布拉克曼窗(Blackman window、切比雪夫窗(Chebyshev window、巴特里特窗(Bartlett window及凱塞窗(Kaiser window。在MA TLAB 中，實現(xiàn)矩形窗的函數(shù)為boxcar 和rectwin ，其調(diào)用格式如下：w=boxcar(Nw=rectwin(N其中N 是窗函數(shù)的長度，返回值w 是一個N 階的向量，它的元素由窗函數(shù)的值組成。實際上，w=boxcar(N等價于w=ones(N,1。在MA TLAB 中，實現(xiàn)三角窗的函數(shù)為triang ，調(diào)用格式為：w=triang(N在MA TLAB 中，實現(xiàn)漢寧窗的函數(shù)為hann ，調(diào)用格式如下：w=h

37、ann(Nw=hann(N,'sflag'Hann 函數(shù)中的參數(shù)sflag 為采樣方式，其值可取symmetric(默認(rèn)值 或periodic 。當(dāng)sflag symmetric 時，為對稱采樣；當(dāng)sflag periodic 時，為周期采樣，此時hann 函數(shù)計算N+1個點(diǎn)的窗，但是僅返回前N 個點(diǎn)。在MA TLAB 中，實現(xiàn)海明窗的函數(shù)為hamming ，調(diào)用格式分別如下：w=hamming (Nw=hamming (N,'sflag'其中sflag 的用法同上。在MA TLAB 中，實現(xiàn)布拉克曼窗的函數(shù)為blackman ，調(diào)用格式如下：w=blackma

38、n (Nw=blackman (N,'sflag'在MA TLAB 中，實現(xiàn)切比雪夫窗的函數(shù)為chebwin ，調(diào)用格式為：w=chebwin (N,r其中r 表示切比雪夫窗函數(shù)的傅里葉變換旁瓣幅度比主瓣低rdB(其默認(rèn)值為100dB ，且旁瓣是等紋波的。在MA TLAB 中，實現(xiàn)巴特里特窗的函數(shù)為bartlett ，調(diào)用格式為：w=bartlett (N在MA TLAB 中，實現(xiàn)凱塞窗的函數(shù)為kaiser ，調(diào)用格式為：w=kaiser (N,beta其中beta 為窗函數(shù)的參數(shù)。各種窗函數(shù)的性能比較可見表4-2，在設(shè)計FIR 濾波器的過程中可以根據(jù)要求選擇合適的窗函數(shù)：表4

39、-2 各種窗函數(shù)的性能比較窗 函 數(shù)第一旁瓣相對于主瓣衰減/dB主 瓣 寬阻帶最小衰減/dB矩形窗134/N21三角窗258/N25漢寧窗318/N44海明窗418/N53布拉克曼窗5712/N74凱塞窗可調(diào)可調(diào)可調(diào)切比雪夫窗可調(diào)可調(diào)可調(diào)1. 數(shù)字低通濾波器的窗函數(shù)設(shè)計計算、型理想低通濾波器的單位沖激響應(yīng)hd(n的MATLAB 實現(xiàn)如例程4-11所示。 例程4-11 、型理想低通濾波器的單位沖激響應(yīng)計算function hd=ideal_lp1(Wc,N%compute the ideal lowpass fiter unit pulse respondence hd(n%wc: cutoff

40、 frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(Wc*m./(pi*m;計算、型理想低通濾波器的單位沖激響應(yīng)hd(n的MATLAB 實現(xiàn)如例程4-12所示。 例程4-12 、型理想低通濾波器的單位沖激響應(yīng)計算function hd=ideal_lp2(Wc,N%compute the ideal lowpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window leng

41、th%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(Wc*m./(pi*m;【例4-7】根據(jù)下列技術(shù)指標(biāo)，設(shè)計一個FIR 數(shù)字低通濾波器：=0.2，=0.4，=0.25bB，=50bB選擇一個適當(dāng)?shù)拇昂瘮?shù)，確定單位沖激響應(yīng)，繪出所設(shè)計的濾波器的幅度響應(yīng)。根據(jù)窗函數(shù)最小阻帶衰減的特性表4-2，可采用海明窗和布拉克曼窗可提供大于50dB 的衰減。本例選擇海明窗，其過渡帶為6.6/N，因此具有較小的階次。該例的MA TLAB 實現(xiàn)如例程4-13 所示，其設(shè)計濾波器的沖激響應(yīng)與幅度響應(yīng)如圖4-7所示。 例程4

42、-13 利用海明窗設(shè)計型數(shù)字低通濾波器clear all;Wp=0.2*pi;Ws=0.4*pi;tr_width=Ws-Wp; %過渡帶寬度N=ceil(6.6*pi/tr_width+1 %濾波器長度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通濾波器的截止頻率hd=ideal_lp1(Wc,N; %理想低通濾波器的單位沖激響應(yīng)w_ham=(hamming(N' %海明窗h=hd.*w_ham; %截取得到實際的單位脈沖響應(yīng)db,mag,pha,w=freqz_m2(h,1; %計算實際濾波器的幅度響應(yīng)delta_w=2*pi/1000;Ap=-(min(db(1:1:W

43、p/delta_w+1 %實際通帶紋波As=-round(max(db(Ws/delta_w+1:1:501 %實際阻帶紋波subplot(221stem(n,hdtitle('理想單位脈沖響應(yīng)hd(n'subplot(222stem(n,w_hamtitle('海明窗w(n'subplot(223stem(n,htitle('實際單位脈沖響應(yīng)hd(n'subplot(224plot(w/pi,dbtitle('幅度響應(yīng)(dB'axis(0,1,-100,10%-functiondb,mag,pha,w=freqz_m2(b,a%

44、濾波器的幅值響應(yīng)(相對、絕對 、相位響應(yīng)%db:相對幅值響應(yīng)%mag：絕對幅值響應(yīng)%pha: 相位響應(yīng)%w 采樣頻率；%b 系統(tǒng)函數(shù)H(z的分子項(對FIR ，b=h%a 系統(tǒng)函數(shù)H(z的分母項(對FIR ，a=1H,w=freqz(b,a,1000,'whole'H=(H(1:1:501'w=(w(1:1:501'mag=abs(H; %絕對幅值響應(yīng)db=20*log10(mag+eps/max(mag; %相對幅值響應(yīng)pha=angle(H; %相位響應(yīng)【運(yùn)行結(jié)果】N =34Ap =0.0477As =52由此結(jié)果可知，所設(shè)計的低通濾波器為型濾波器，它的通帶

45、紋波和阻帶紋波均滿足設(shè)計要求。圖4-7 型FIR 數(shù)字低通濾波器沖激響應(yīng)與幅度響應(yīng)2. 數(shù)字高通濾波器的窗函數(shù)設(shè)計數(shù)字高通濾波器的設(shè)計分兩種情況，即無相移和相移為±/2兩種情況。當(dāng)無相移時，N 為奇數(shù)時，所設(shè)計的FIR 數(shù)字高通濾波器為型濾波器；當(dāng)N 為偶數(shù)時，為型濾波器。當(dāng)相移為±/2時，N 為奇數(shù)時，所設(shè)計的FIR 數(shù)字高通濾波器為型濾波器；當(dāng)N 為偶數(shù)時，為型濾波器、型理想高通濾波器的單位沖激響應(yīng)hd(n計算的MATLAB 實現(xiàn)如例程4-14所示。 例程4-14 、型理想高通濾波器的單位沖激響應(yīng)計算function hd=ideal_hp1(Wc,N%compute

46、the ideal highpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(pi*m-sin(Wc*m./(pi*m;、型理想高通濾波器的單位沖激響應(yīng)hd(n計算的MATLAB 實現(xiàn)如例程4-15所示。 例程4-15 、型理想高通濾波器的單位沖激響應(yīng)hd(n計算function hd=ideal_hp2(Wc,N%compute the ideal h

47、ighpass fiter unit pulse respondence hd(n%wc: cutoff frequency%N: window length%hd: unit pulse respondencealpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(pi*m-cos(Wc*m./(pi*m;【例4-8】根據(jù)下列技術(shù)指標(biāo)，設(shè)計一個FIR 數(shù)字高通濾波器：=0.6，=0.4，=0.25bB，=50bB選擇一個適當(dāng)?shù)拇昂瘮?shù)，確定單位沖激響應(yīng)，繪出所設(shè)計濾波器的幅度響應(yīng)。根據(jù)窗函數(shù)最小阻帶衰減的特性表4-2所示，選擇漢寧窗可達(dá)到 44dB 的最小阻帶衰

48、減，其過渡帶為6.2/N，因此具有較小的階次。該例的MA TLAB 實現(xiàn)如例程4-16所示，其設(shè)計的濾波器的沖激響應(yīng)與幅度響應(yīng)如圖4-8所示。 例程4-16 利用漢寧窗設(shè)計型數(shù)字高通濾波器clear all;Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws; %過渡帶寬度N=ceil(6.2*pi/tr_width %濾波器長度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通濾波器的截止頻率hd=ideal_hp1(Wc,N; %理想低通濾波器的單位沖激響應(yīng)w_han=(hanning(N' %漢寧窗h=hd.*w_han; %截取得到實際的單位脈沖響應(yīng)db

49、,mag,pha,w=freqz_m2(h,1; %計算實際濾波器的幅度響應(yīng)delta_w=2*pi/1000;Ap=-(min(db(Wp/delta_w+1:1:501 %實際通帶紋波As=-round(max(db(1:1:Ws/delta_w+1 %實際阻帶紋波subplot(221stem(n,hdtitle('理想單位脈沖響應(yīng)hd(n'subplot(222stem(n,w_hantitle('漢寧窗w(n'subplot(223stem(n,htitle('實際單位脈沖響應(yīng)hd(n'subplot(224plot(w/pi,dbti

50、tle('幅度響應(yīng)(dB'axis(0,1,-100,10【運(yùn)行結(jié)果】N = 31Ap = 0.0887As = 44由此結(jié)果可知，所設(shè)計的高通濾波器為型濾波器，它的通帶紋波和阻帶紋波均滿足設(shè)計要求。圖4-8 型FIR 數(shù)字高通濾波器沖激響應(yīng)與幅度響應(yīng)【例4-9】根據(jù)下列技術(shù)指標(biāo)，設(shè)計一個具有+/2相移的FIR 數(shù)字高通濾波器：=0.6，=0.4，=0.25bB，=50bB選擇一個適當(dāng)?shù)拇昂瘮?shù)，確定單位沖激響應(yīng)，繪出所設(shè)計的濾波器的幅度響應(yīng)。根據(jù)窗函數(shù)最小阻帶衰減的特性表4-2，選擇漢寧窗可達(dá)到 44dB 的最小阻帶衰減，其過渡帶為6.2/N，因此具有較小的階次。該例的MA T

51、LAB 實現(xiàn)如例程4-17所示，其設(shè)計的濾波器的沖激響應(yīng)與幅度響應(yīng)如圖4-9所示。 例程4-17 利用漢寧窗設(shè)計型數(shù)字高通濾波器clear all;Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws; %過渡帶寬度N=ceil(6.2*pi/tr_width+1 %濾波器長度n=0:1:N-1;Wc=(Ws+Wp/2; %理想低通濾波器的截止頻率hd=ideal_hp2(Wc,N; %理想低通濾波器的單位沖激響應(yīng)w_han=(hanning(N' %漢寧窗h=hd.*w_han; %截取得到實際的單位脈沖響應(yīng)db,mag,pha,w=freqz_m2(h,1; %計算

52、實際濾波器的幅度響應(yīng)delta_w=2*pi/1000;Ap=-(min(db(Wp/delta_w+1:1:501 %實際通帶紋波As=-round(max(db(1:1:Ws/delta_w+1 %實際阻帶紋波subplot(221stem(n,hdtitle('理想單位脈沖響應(yīng)hd(n'subplot(222stem(n,w_hantitle('漢寧窗w(n'subplot(223stem(n,htitle('實際單位脈沖響應(yīng)hd(n'subplot(224plot(w/pi,dbtitle('幅度響應(yīng)(dB'axis(0,

53、1,-100,10【運(yùn)行結(jié)果】N = 32Ap = 0.0711As = 44由此結(jié)果可知，所設(shè)計的高通濾波器為型濾波器，它的通帶紋波和阻帶紋波均滿足設(shè)計要求。圖4-9 型FIR 數(shù)字高通濾波器沖激響應(yīng)與幅度響應(yīng)3. 數(shù)字帶通濾波器的窗函數(shù)設(shè)計數(shù)字帶通濾波器的設(shè)計也分兩種情況討論，即無相移和相移為±/2兩種情況。當(dāng)無相移時，N 為奇數(shù)時，所設(shè)計的FIR 數(shù)字帶通濾波器為型濾波器；當(dāng)N 為偶數(shù)時，為型濾波器。當(dāng)相移為±/2時，N 為奇數(shù)時，所設(shè)計的FIR 數(shù)字帶通濾波器為型濾波器；當(dāng)N 為偶數(shù)時，為型濾波器、型理想帶通濾波器的單位沖激響應(yīng)hd(n計算的MATLAB 實現(xiàn)如例程

54、4-18所示。 例程4-18 、型理想帶通濾波器的單位沖激響應(yīng)計算function hd=ideal_bp1(Wcl,Wch,N%compute the ideal bandpass fiter unit pulse respondence hd(n%wcl: low cutoff frequency%wch: high cutoff frequency%N: window length%hd: unit pulse respondencealpha= (N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(Wch*m-sin(Wcl*m./(pi*m;、型理想帶通濾波器的

55、單位沖激響應(yīng)hd(n計算的MATLAB 實現(xiàn)如例程4-19所示。 例程4-19 、型理想帶通濾波器的單位沖激響應(yīng)hd(n計算function hd=ideal_bp2(Wcl,Wch,N%compute the ideal bandpass fiter unit pulse respondence hd(n%wcl: low cutoff frequency%wch: high cutoff frequency%N: window length%hd: unit pulse respondencealpha= (N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=cos(Wch*

56、m-cos(Wcl*m./(pi*m;【例4-10】根據(jù)下列技術(shù)指標(biāo)，設(shè)計一個FIR 數(shù)字帶通濾波器：=0.4，=0.6，=0.2，=0.8.=1bB，=1bB，=60bB ，=60bB.選擇一個適當(dāng)?shù)拇昂瘮?shù)，確定單位沖激響應(yīng)，繪出所設(shè)計的濾波器的幅度響應(yīng)。根據(jù)窗函數(shù)最小阻帶衰減的特性表4-2，選擇布拉克曼窗可達(dá)到 75dB 的最小阻帶衰減，其過渡帶為11/N。該例的MA TLAB 實現(xiàn)如例程4-20所示，其設(shè)計的濾波器的沖激響應(yīng)與幅度響應(yīng)如圖4-10所示。例程4-20 利用布拉克曼窗設(shè)計型數(shù)字帶通濾波器clear all;Wpl=0.4*pi;Wph=0.6*pi;Wsl=0.2*pi;Wsh=0.8*pi;tr_width=min(Wpl-Wsl,(Wsh-Wph; %過渡帶寬度N=ceil(11*pi/tr_width+1 %濾波器長度n=0:1:N-1;