Matlab和雙線性變換的濾波器設計

1、武漢理工大學MATLAB課程設計報告書 題 目: MATLAB課程設計基于Matlab和雙線性變換的濾波器設計 初始條件：MATLAB仿真軟件數字信號處理與圖像處理基礎知識雙線性變換法的原理和算法巴特沃斯高通濾波器的性能指標要求完成的主要任務: 利用MATLAB仿真軟件系統(tǒng)結合雙線性變換法設計一個數字巴特沃斯高通、低通和帶通IIR濾波器。說明：自己確定濾波器的不同參數得到不同的結果，最后自己分析原因。課程設計題目由指導教師提供，每人隸屬一組完成任務，每組完成的內容不能雷同（按學號分組）時間安排：序號階段內容所需時間1方案設計1天2軟件設計2天3系統(tǒng)調試1天4答辯1天合 計5天指導教師簽名： 年

2、 月 日系主任（或責任教師）簽名： 年 月 日目 錄摘 要數字濾波器是數字信號處理的基礎，用來對信號進行過濾、檢測與參數估計等處理，在通信、圖像、語音、雷達等許多領域都有著十分廣泛的應用。尤其在圖像處理、數據壓縮等方面取得了令人矚目的進展和成就。鑒于此，數字濾波器的設計就顯得尤為重要。 此報告重點介紹了用雙線性不變法設計IIR數字濾波器的基本流程，比較了各種設計方法的優(yōu)缺點，總結了模擬濾波器的性能特征。最后以雙線性不變法設計了一個高通巴特沃斯FIR數字濾波器，介紹了設計步驟，然后在Matlab環(huán)境下進行了仿真與調試，實現(xiàn)了設計目標。關鍵字：數字濾波器 巴特沃斯 MATLAB 雙線性變法Abst

3、ract Matlab for research and application of digital filter provides an intuitive, efficient and convenient weapon. He matrix operations based on the calculation, visualization, program design integrated into an interactive working environment. In particular, Matlab toolbox to enable researchers in v

4、arious fields can be easy and intuitive application of scientific research and engineering. This report introduced with emphasis of the basic flow of designing the IIR digit filter by the bilinear political reform, compared with each kind of design method's good and bad points, summarized analog

5、 filter's performance characteristic. Finally design one by the bilinear political reform to pass Butterworth high IIR digit filter, introduced the design procedure, then has carried on the simulation and the debugging under the Matlab environment, has achieved the project objective. Keyword: di

6、gital filter Butterworth MATLAB1 數字濾波器1.1 數字濾波器的基本概念數字濾波器是數字信號處理的重要基礎，是對信號都是過濾檢測與參數估計等處理過程中，它是使用最為廣泛的一種線性系統(tǒng)。數字濾波器處理的對象是經由采樣期間將模擬信號轉換而得到的數字信號。數字濾波器是指完成信號濾波處理功能的，用有限精度算法實現(xiàn)的離散時間線性非時變系統(tǒng)。數字濾波器的輸入是一組數字量。它本身既可以是用數字硬件裝配而成的一臺用于完成給定運算的專用數字計算機，也可以是將所需的運算編寫的程序通過計算機來執(zhí)行。數字濾波器具有穩(wěn)定性高、精度高、靈活性大等突出優(yōu)點。隨著數字技術的發(fā)展，用數字技術實現(xiàn)

7、濾波器的功能愈來愈受到人們的重視，并得到了廣泛的應用。數字濾波器的數學運算通常有兩種實現(xiàn)方式。一種是頻域法，即利用FFT快速運算方法對輸入信號進行離散傅里葉變換，分析其頻譜，然后再根據所希望的頻率特性進行濾波，再利用傅里葉反變換得到時域信號。這種方法具有較好的頻域選擇特性和靈活性，并且由于信號頻率與所希望的頻率特性是簡單的相乘關系，所以它比計算等價的時域卷積要快得多。另一種方法是時域法，這種方法是通過離散的抽樣數據做差分數學運算來達到濾波目的的。1.2 數字濾波器的分類數字濾波器按照不同的分類方法，可分為許多種，但總體來講可以分成兩大類。一類稱為經典濾波器，即一般濾波器，特點是輸入信號中有用的

8、頻率成分和希望濾除的頻率成分各占有不同的頻帶，通過一個合適的選頻濾波器達到濾波的目的。但對于一般濾波器，如果信號和干擾的頻帶互不重疊，則不能完成對干擾的有效濾除，這是需要采用另一類所謂的現(xiàn)代濾波器，例如維納濾波器、卡爾曼濾波器、自適應濾波器等最佳濾波器，這些濾波器可按照隨機信號內部的一些統(tǒng)計分布規(guī)律，從干擾中最佳地提取信號。從功能上，一般數字濾波器可以分為低通、高通、帶通、帶阻和全通等，此種分類方法是和模擬濾波器一樣的。1.3 數字濾波器的MATLAB實現(xiàn)MATLAB工具箱未濾波器的設計應用提供了豐富而簡便的方法，如函數方法和圖形工具方法等，使原來非常繁瑣復雜的程序設計變成了簡單的函數調用，為

9、濾波器的設計和實現(xiàn)開辟了廣闊的天地。數字濾波器的一般設計過程為：按照實際需要，確定濾波器的性能要求；用一個因果的、穩(wěn)定的離散線性時不變系統(tǒng)，去逼近這一性能指標；用有限精度的運算實現(xiàn)所設計的系統(tǒng)；通過模擬，驗證所設計的系統(tǒng)是否符合給定性能要求。Dd2 雙線性變換法2.1 雙線性變換法知識簡介脈沖響應不變法的主要缺點是產生頻率響應的混疊失真。這是因為從S平面到平面是多值的映射關系所造成的。為了克服這一缺點，可以采用非線性頻率壓縮方法，將整個頻率軸上的頻率范圍壓縮到-/T/T之間，再用z=esT轉換到Z平面上。也就是說，第一步先將整個S平面壓縮映射到S1平面的-/T/T一條橫帶里；第二步再通過標準變

10、換關系z=es1T將此橫帶變換到整個Z平面上去。這樣就使S平面與Z平面建立了一一對應的單值關系，消除了多值變換性，也就消除了頻譜混疊現(xiàn)象。由圖2-1看出，在零頻率附近，模擬角頻率與數字頻率之間的變換關系接近于線性關系；但當進一步增加時，增長得越來越慢，最后當時，終止在折疊頻率=處，因而雙線性變換就不會出現(xiàn)由于高頻部分超過折疊頻率而混淆到低頻部分去的現(xiàn)象，從而消除了頻率混疊現(xiàn)象。圖2-12.2 雙線性變換法設計數字濾波器原理將S平面j軸壓縮變換到s1平面j軸上的-/T到/T一段，可以采用以下變換關系： 這樣變換到，變換到，可將上式寫成令，解析延拓到整個s平面和s1平面，可得再將平面通過以下變換關

11、系映射到z平面，即從而得到s平面和z平面的單值映射關系為一般來說，為了使模擬濾波器的某一頻率與數字濾波器的任一頻率有對應關系，可引入待定常數c，將代入到上式，可得在MATLAB中，雙線性Z變換可以通過bilinear函數實現(xiàn)，其調用格式為：Bz，Azbilinear(B，A，F(xiàn)s)；其中B，A為模擬濾波器傳遞函數G(s)的分子分母多項式的系數向量，而Bz，Az為數字濾波器的傳遞函數H(z)的分子分母多項式的系數向量。3 設計任務及方案選擇3.1 設計任務與要求利用MATLAB仿真軟件系統(tǒng)結合雙線性變換法設計一個數字巴特沃斯高通、低通和帶通FIR濾波器。說明：自己確定濾波器的不同參數得到不同的結

12、果，最后自己分析原因。假設一個信號：其中f1=100Hz，f2=550Hz。試將該信號分別通過低通、高通、帶通濾波器，對輸出信號進行分析比較。Step1:Step2:Step3:Step4:3.2 方案設計論證與選擇本次濾波器設計的主要方法是先設計低通模擬濾波器，然后轉換為高通、帶通或帶阻數字濾波器。對于其他如高通，帶通，則通過頻率變換轉換為設計相應的高通，帶通等。在設計的全過程的各個步驟，matlab都提供相應的工具箱函數，使得IIR數字濾波器設計變得非常簡單?？偟膩碚f，我的設計思路主要有以下兩種：方案一：從歸一化模擬低通原型出發(fā)，先在模擬域內經頻率變換成為所需類型的模擬濾波器；然后進行雙線

13、性變換，由S域變換到Z域，而得到所需類型的數字濾波器。歸一化模擬低通原型數字高,帶通或帶阻模擬高,帶通或帶阻 模擬域 雙線性變換法頻率變換 圖3-1 先頻率變換再離散方案二：先進行雙線性變換，將模擬低通原型濾波器變換成數字低通濾波器；然后在Z域內經數字頻率變換為所需類型的數字濾波器。 歸一化模擬低通原型數字高,帶通或帶阻數字原型低通 數字域 雙線性變換法 頻率變換圖3-2 先離散再頻率變換以上兩種方案都可以，我最后選擇了第一種方案進行設計，即先在模擬域內經頻率變換成為所需類型的模擬濾波器；然后進行雙線性變換，由S域變換到Z域，而得到所需類型的數字濾波器。因為高通，帶通濾波器的設計方法是先將要設

14、計的濾波器的技術特性指標通過頻率轉換成模擬低通濾波器的技術指標，再根據這些性能指標設計出低通濾波器的傳遞函數，為了設計程序的簡潔，故選擇了方案一。4 程序設計與調試4.1 設計要求若實驗設計要求得巴特沃斯高通濾波器的性能指標為：通帶邊界頻率fp：fp=400Hz；通帶波紋 Rp：Rp=1dB；阻帶邊界頻率fs：fs=300Hz；阻帶衰減 Rs： Rs=30dB；采樣頻率 Fs： Fs=2000Hz；要求巴特沃斯低通濾波器通帶邊界頻率與阻帶邊界頻率與高通的正好相反。模擬帶通濾波器中心頻率為f0=200Hz，帶寬為Bw=400Hz，4.2 程序設計與說明4.2.1 模擬低通濾波首先應用雙線性變換作

15、為數字低通濾波器到模擬高通濾波器的頻率預畸變處理： 利用上式確定的技術指標和已知的參數Rp、Rs，就可以進行模擬原型低通濾波器設計了。該部分的程序設計與注釋如下： Fs=2000; Ts=1/Fs; %采樣頻率和采樣周期wp=400*2*pi;ws=300*2*pi; %低通、高通數字濾波器截止頻率Rp=1;Rs=30; %通帶波紋和阻帶衰減wpl=Fs*2*tan(ws*Ts/2); %按頻率轉換公式進行轉換wsl=Fs*2*tan(wp*Ts/2);N,Wn=buttord(wpl,wsl,Rp,Rs,'s'); %計算模擬濾波器的最小階數z,p,k=buttap(N);

16、%設計模擬原型低通濾波器b,a=zp2tf(z,p,k); %變零極點增益形式為系統(tǒng)傳遞函數形式H,w=freqs(b,a); %求模擬低通濾波器的頻率響應figure;plot(w,abs(H);grid; %繪制低通頻率特性title('模擬低通原型');xlabel('w');ylabel('幅度');4.2.2 模擬到數字的轉換 （1）從低通向低通的轉換b,a=lp2lp(Bap,Aap,Wn);功能是把模擬濾波器原型轉換成截至頻率為 Wn 的低通濾波器。再由得到的模擬低通濾波器通過函數Bz，Azbilinear(B，A，F(xiàn)s)得到數字低

17、通濾波器。該部分的程序設計如下：e1,f1=lp2lp(b,a,300*2*pi) ; %原型低通轉換為低通模擬濾波器h1,w1=freqs(e1,f1); %求模擬低通濾波器的頻率響應figure;plot(w1/(2*pi),abs(h1); %繪制高通頻率特性grid;title('模擬低通濾波器');xlabel('頻率Hz');ylabel('幅度');r1,q1=bilinear(e1,f1,1000); %模擬低通濾波器轉換為數字濾波器s1,w2=freqz(r1,q1); %求數字低通濾波器的頻率響應figure;plot(w2/

18、pi, 20*log10(abs(s1);grid; %繪制數字低通濾波器的頻率響應title('數字低通濾波器');xlabel('w/pi');ylabel('幅度(dB)');（2）從低通向高通的轉換b,a=lp2hp(Bap,Aap,Wn);功能是把模擬濾波器原型轉換成截至頻率為 Wn 的模擬高通濾波器，再進一步轉換成數字高通濾波器。程序設計如下：e2,f2=lp2hp(b,a,300*2*pi) ; %原型低通模擬濾波器轉換為高通模擬濾波器h2,w3=freqs(e2,f2); %求模擬高通濾波器的頻率響應figure;plot(w3/

19、(2*pi),abs(h2); %繪制高通頻率特性grid;title('模擬高通濾波器');xlabel('頻率Hz');ylabel('幅度');r2,q2=bilinear(e2,f2,1000); %模擬高通濾波器轉換為數字濾波器s2,w4=freqz(r2,q2); %求數字高通濾波器的頻率響應figure;plot(w4/pi, 20*log10(abs(s2);grid; %繪制數字高通濾波器的頻率響應title('數字高通濾波器');xlabel('w/pi');ylabel('幅度(dB)

20、');（3）從低通向帶通的轉換b,a=lp2bp(Bap,Aap,W0，Bw);功能是把模擬濾波器原型轉換成中心頻率為 W0、帶寬為Bw的模擬帶通濾波器，在進一步轉換成數字帶通濾波器。程序設計如下：W0=200*2*pi;Bw=200;e3,f3=lp2bp(b,a,W0,Bw) ; %原型低通轉換為帶通模擬濾波器h3,w5=freqs(e3,f3); %求模擬帶通濾波器的頻率響應figure;plot(w5/(2*pi),abs(h3); %繪制帶通頻率特性grid;title('模擬帶通濾波器');xlabel('頻率Hz');ylabel('

21、;幅度');r3,q3=bilinear(e3,f3,1000); %模擬濾波器轉換為數字濾波器s3,w6=freqz(r3,q3); %求數字帶通濾波器的頻率響應figure;plot(w6/pi, 20*log10(abs(s3);grid; %繪制數字帶通濾波器的頻率響應title('數字帶通濾波器');xlabel('w/pi');ylabel('幅度(dB)');4.1.3 對輸入信號進行濾波讀取輸入信號程序：n=0:180;dt=1/Fs;t=n*dt; %時間序列f1=100;f2=550; %輸入信號頻率x=sin(2*p

22、i*f1*t)+cos(2*pi*f2*t)+sin(2*pi*1*t); %輸入信號figure;plot(t,x);grid; %繪制輸入信號title('輸入信號');xlabel('時間t');ylabel('幅度');將輸入信號分別通過得到的低通、高通、帶通濾波器，調用函數y=filter(bz,az,x)對輸入信號進行濾波。y=filter(r1,q1,x); %對輸入信號進行低通濾波 figure;plot(t,y);grid; %繪制輸出信號title('低通濾波輸出信號');xlabel('時間t'

23、;);ylabel('幅度')再進行高通、帶通濾波只需將調用函數filter中參量改為r2、q2和r3、q3即可。4.3 仿真結果與分析將設計好的程序在matlab軟件下進行仿真，得到仿真波形。4.3.1 仿真結果記錄圖4-1 模擬低通原型幅頻響應（a）低通模擬濾波器的頻率特性曲線 （b）巴特沃斯低通數字濾波器的頻率特性曲線圖4-2 （a）高通模擬濾波器的頻率特性曲線 （b）巴特沃斯高通數字濾波器的頻率特性曲線圖4-3（a）帶通模擬濾波器的頻率特性曲線 （b）巴特沃斯帶通數字濾波器的頻率特性曲線圖4-4圖4-5 輸入信號波形圖4-6 經過巴特沃斯低通濾波器濾波后的輸出波形圖4-

24、7 經過巴特沃斯高通濾波器濾波后的輸出波形圖4-8 經過巴特沃斯高通濾波器濾波后的輸出波形4.3.2 仿真結果分析由仿真結果可知，設計出的數字低通，高通，帶通的幅頻響應均符合各自的特性，從模擬濾波器特性曲線可以看出，所設計的高通濾波器在大于400Hz為通帶，其衰減均小于1dB；小于300Hz為阻帶，其衰減大于30dB，而設計的低通濾波器情況正好相反，符合本次設計的要求。當濾波器輸入100Hz和550Hz兩種信號后，低通濾波器可以濾除頻率為550Hz的高頻信號，而低頻信號100Hz的信號可以通過濾波器，驗證了該低通濾波器的功能。而高通濾波器可以濾除頻率為100Hz的信號，而高頻信號550Hz的信

25、號可以通過濾波器。4.3.3 程序調試中出現(xiàn)的問題帶通濾波器的仿真結果雖然符合其帶通特性，但要求得到的性能指標存在較大誤差，由于時間緊迫，將在課后進一步尋找原因及解決方案。一個線性是不變系統(tǒng)的典型表達式有狀態(tài)空間型、傳遞函數型、零極增益型、極點留數型共四種，在開始的設計過程中我發(fā)現(xiàn)在把零極點增益型向傳遞函數型的轉換過程中多用了一步，即先把零極點增益型轉換為狀態(tài)空間型，進行相應的頻率變換和雙線性變換后，再由狀態(tài)空間型轉換為傳遞函數型，進而繪制出其頻率響應圖，但是經仿真分析后發(fā)現(xiàn)，性能并不理想，所以并沒有采用。5心得體會本次課程設計是基于數字濾波器的MATLAB軟件應用，在本學期中，我在數字信號處

26、理課程中學習了關于本次課程設計的一些內容，對于MATLAB軟件在數字濾波器的設計中的一些函數也有所認識，通過著一周的時間更是讓自己在這個方面的學習更進一步，使自己在這方面的學習深度有所加深。本次用雙線性變換法對巴特沃斯低通，高通及帶通濾波器的設計，要求我能夠用所學理論知識對巴特沃斯高通濾波器進行設計，并能夠對設計的濾波器的頻率特性進行分析。在設計過程中，首先，我通過查找資料等途徑，對關于所要設計的題目進行更加全面深入的學習，認真分析其中的原理，弄明白設計的思路；其次，根據自己選擇的性能指標要求，結合自己通過學習掌握的知識和資料進行程序設計，根據自己提出的方案利用軟件實現(xiàn)；最后，運行程序，分析所

27、得結果是否符合自己的設計要求，對結果進行分析，實現(xiàn)由理論到實踐，再由實踐結果驗證理論的過程。總的來說，這次課程設計讓我對MATLAB有了更深刻的了解，對數字濾波器的設計流程有了大致的了解，掌握了一些設計濾波器的基本方法，提高了理論用于實踐的能力，掌握了更多專業(yè)相關的使用知識與技能。同時，也暴露了我很多的不足，在以后的學習中，將進一步發(fā)揚有點，克服缺點。 參考文獻1袁潤章，自蔓延高溫合成技術研究進展.武漢:武漢工業(yè)大學出版社,19942Schacht E.Industrial polysaccharides. Amsterdam:Elsevier Science,19873郭景坤.陶瓷材料的強化

28、與增韌新途徑的探索.無機材料學報,1998,13910;23-264Spriggs G E.A history of fine grained hardmetal. Int J of Refractory Metal and Hard Material,1995,13:241-2555王為民.二硼化鈦陶瓷的自蔓延高溫合成與加工:博士論文.武漢工業(yè)大學,19886吳中慶,張基田.SiC顆粒與液體鋁的界石層.廣州:92秋季中國材料科學研討會,C-MRS,1992.603-609.附錄本次課設總程序如下：Fs=2000; Ts=1/Fs; %采樣頻率和采樣周期wp=400*2*pi;ws=300*2

29、*pi; %低通、高通數字濾波器截止頻率Rp=1;Rs=30; %通帶波紋和阻帶衰減wpl=Fs*2*tan(ws*Ts/2); %按頻率轉換公式進行轉換wsl=Fs*2*tan(wp*Ts/2);N,Wn=buttord(wpl,wsl,Rp,Rs,'s'); %計算模擬濾波器的最小階數z,p,k=buttap(N); %設計模擬原型低通濾波器b,a=zp2tf(z,p,k); %變零極點增益形式為系統(tǒng)傳遞函數形式H,w=freqs(b,a); %求模擬低通濾波器的頻率響應figure;plot(w,abs(H);grid; %繪制低通頻率特性title('模擬低通原

30、型');xlabel('w');ylabel('幅度');e1,f1=lp2lp(b,a,300*2*pi) ; %原型低通轉換為低通模擬濾波器h1,w1=freqs(e1,f1); %求模擬低通濾波器的頻率響應figure;plot(w1/(2*pi),abs(h1); %繪制高通頻率特性grid;title('模擬低通濾波器');xlabel('頻率Hz');ylabel('幅度');r1,q1=bilinear(e1,f1,1000); %模擬低通濾波器轉換為數字濾波器s1,w2=freqz(r1,q

31、1); %求數字低通濾波器的頻率響應figure;plot(w2/pi, 20*log10(abs(s1);grid; %繪制數字低通濾波器的頻率響應title('數字低通濾波器');xlabel('w/pi');ylabel('幅度(dB)');e2,f2=lp2hp(b,a,300*2*pi) ; %原型低通模擬濾波器轉換為高通模擬濾波器h2,w3=freqs(e2,f2); %求模擬高通濾波器的頻率響應figure;plot(w3/(2*pi),abs(h2); %繪制高通頻率特性grid;title('模擬高通濾波器');

32、xlabel('頻率Hz');ylabel('幅度');r2,q2=bilinear(e2,f2,1000); %模擬高通濾波器轉換為數字濾波器s2,w4=freqz(r2,q2); %求數字高通濾波器的頻率響應figure;plot(w4/pi, 20*log10(abs(s2);grid; %繪制數字高通濾波器的頻率響應title('數字高通濾波器');xlabel('w/pi');ylabel('幅度(dB)');W0=200*2*pi;Bw=200;e3,f3=lp2bp(b,a,W0,Bw) ; %原型低

33、通轉換為帶通模擬濾波器h3,w5=freqs(e3,f3); %求模擬帶通濾波器的頻率響應figure;plot(w5/(2*pi),abs(h3); %繪制帶通頻率特性grid;title('模擬帶通濾波器');xlabel('頻率Hz');ylabel('幅度');r3,q3=bilinear(e3,f3,1000); %模擬濾波器轉換為數字濾波器s3,w6=freqz(r3,q3); %求數字帶通濾波器的頻率響應figure;plot(w6/pi, 20*log10(abs(s3);grid; %繪制數字帶通濾波器的頻率響應n=0:180;

34、dt=1/Fs;t=n*dt; %時間序列f1=100;f2=550; %輸入信號頻率x=sin(2*pi*f1*t)+cos(2*pi*f2*t)+sin(2*pi*1*t); %輸入信號figure;plot(t,x);grid; %繪制輸入信號title('輸入信號');xlabel('時間t');ylabel('幅度');y=filter(r1,q1,x); %對輸入信號進行低通濾波 figure;plot(t,y);grid; %繪制輸出信號title('低通濾波輸出信號');xlabel('時間t');ylabel('幅度')y=filter(r2,q2,x); %對輸入信號進行高通濾波 figure;plot(t,y);grid; %繪制輸出信號title('高通通濾波輸出信號'