數(shù)字濾波器研究課程論文

1、研 究 生 課 程 論 文(2010-2011學(xué)年第一學(xué)期)數(shù)字濾波器的研究研究生： 提交日期： 研究生簽名：學(xué) 號學(xué) 院電子與信息學(xué)院課程編號S0809074課程名稱現(xiàn)代數(shù)字信號處理學(xué)位類別碩士任課教師教師評語：成績評定： 分 任課教師簽名： 年 月 日中文題目：數(shù)字濾波器的研究研究生姓名 摘要：無限沖激響應(yīng)（IIR）數(shù)字濾波器和有限沖激響應(yīng)（FIR）數(shù)字濾波器廣泛應(yīng)用于數(shù)字信號處理系統(tǒng)中。FIR 濾波器具有很好的線性相位特性，因此線性相位特性在工程實際中如數(shù)據(jù)通信、圖像處理等應(yīng)用領(lǐng)域具有比較重要的意義。IIR 數(shù)字濾波器具結(jié)構(gòu)簡單和實現(xiàn)容易的特點, 這對于相位要求不是特別嚴(yán)格的場合有很重要

2、的意義。利用MATLAB 實現(xiàn)IIR 低通數(shù)字濾波器的設(shè)計,并對被干擾的信號進(jìn)行濾波,達(dá)到了預(yù)期結(jié)果。關(guān)鍵詞: IIR；FIR；線性相位；數(shù)字濾波器；MATLAB1引言數(shù)字濾波器是現(xiàn)代數(shù)字信號處理的不可或缺的一部分，數(shù)字濾波器按其響應(yīng)形式可分為無限脈沖響應(yīng)（IIR）濾波器和有限脈沖響應(yīng)（FIR）濾波器和兩種。當(dāng)系統(tǒng)無嚴(yán)格相位要求時，IIR 數(shù)字濾波器與FIR 數(shù)字濾波器相比，可用較低的階數(shù)獲得較高的選擇性1。對于相同的設(shè)計指標(biāo)，F(xiàn)IR 濾波器所要求的階數(shù)比IIR濾波器高510 倍，而且信號的延遲也較大2。MATLAB 作為一種矩陣運(yùn)算為基礎(chǔ)的交互式程序語言，是進(jìn)行科學(xué)研究常用且必不可少的工具。

3、MATLAB著重針對科學(xué)計算、工程計算和繪圖的需求。它用解釋方式工作，鍵入程序立即得出結(jié)果，人機(jī)交互性能好，為科學(xué)人員所樂于接受。MATLAB 提供了豐富的函數(shù)，其中fir1 函數(shù)實現(xiàn)了加窗線性相位FIR 濾波器設(shè)計的經(jīng)典方法，fir1 主要用于常用的標(biāo)準(zhǔn)通帶濾波器設(shè)計，包括低通、帶通、高通和帶阻數(shù)字濾波器3。一般來講，數(shù)字濾波器的設(shè)計要經(jīng)過三步：確定設(shè)計指標(biāo)、模擬逼近和數(shù)字轉(zhuǎn)換。通常在設(shè)計濾波器之前，先根據(jù)具體的應(yīng)用確定一些技術(shù)指標(biāo)，然后就可以根據(jù)數(shù)學(xué)知識和濾波器的基本原理提出一個濾波器的模型來逼近給定的指標(biāo)，逼近的結(jié)果通常是得到以差分方程或脈沖響應(yīng)描述的濾波器，最后可以根據(jù)這個描述用硬件或

4、軟件實現(xiàn)，至此完成一個濾波器設(shè)計的全過程3。2.數(shù)字濾波器的設(shè)計2.1模擬橢圓濾波器對加噪語音信號的濾波：設(shè)計模擬濾波器的一般方法是：首先根據(jù)技術(shù)指標(biāo)確定濾波器的傳輸函數(shù)H(s）,然后綜合電路網(wǎng)絡(luò)實現(xiàn)該傳遞函數(shù)。而設(shè)計濾波器傳輸函數(shù)的關(guān)鍵是找到逼近函數(shù)，在各種濾波器響應(yīng)中，橢圓函數(shù)響應(yīng)最優(yōu)越4-6。首先，打開的是一個摩托車引擎信號的波形文件，它是一個頻率范圍大致分布在5000-6000Hz的一個低頻語音信號。對于該語音信號，施加高斯白噪聲，設(shè)計一個1dB截止頻率6000Hz 最小阻帶衰減為100dB的橢圓低通濾波器。其原始信號波形、加噪聲后信號波形與濾波后信號波形分別圖1和圖2：圖1.加噪聲后

5、語音信號波形與頻譜圖2.原始語音信號與濾波后語音信號波形與頻譜對比 通過matlab仿真可以發(fā)現(xiàn)，原始語音信號基本可以從高斯白噪聲中恢復(fù)，但是濾波后的語音信號比原始語音信號有明顯的波形延遲，這是模擬橢圓低通濾波器不足的地方。2.2 IIR濾波器的結(jié)構(gòu)與設(shè)計這里采用的是直接II型結(jié)構(gòu)實現(xiàn)IIR 濾波器，如式（1）所示： （1）此結(jié)構(gòu)便于準(zhǔn)確地實現(xiàn)濾波器的零點、極點，也便于調(diào)整濾波器的頻率響應(yīng)性能7；另一個優(yōu)點是所需的存儲單元較少，在硬件實現(xiàn)時甚至還可以用一個二階節(jié)進(jìn)行時分復(fù)用，進(jìn)一步降低了對現(xiàn)場可編程邏輯陣列硬件資源的要求8。IIR濾波器的系數(shù)在做強(qiáng)噪音下的語音增強(qiáng)，需要通過共振峰對語音進(jìn)行端點

6、識別。前兩個共振峰對區(qū)別不同語音有非常重要的作用9；1 kHz以下基本為第一共振峰的范圍,對語音感知、語意識別作用比較重要的第二共振峰基本在1 kHz之外9。因此對帶噪語音信號加入預(yù)處理環(huán)節(jié),即進(jìn)行數(shù)字濾波(通帶下限1 kHz,Butterworth濾波器高通,阻帶衰減3 dB，語音信號采樣頻率為8000 Hz)。利用MATLAB提供的FDATool便可直接得到直接II型各個子系由于硬件當(dāng)中對小數(shù)進(jìn)行運(yùn)算耗費(fèi)資源比較大，需要在計算精度與速度之間優(yōu)化選擇10。下面用雙線性變換設(shè)計一個0.4dB截止頻率為10KHz且在30KHz處有最小阻帶衰減為50dB的數(shù)字巴特沃茲低通濾波器，其抽樣率為100K

7、Hz；設(shè)計的具體步驟：1： 先設(shè)計一個符合以上指標(biāo)的模擬巴特沃茲低通濾波器，作為與要設(shè)計的數(shù)字濾波器的參照濾波器。先用手算計算出濾波器所需的階數(shù)是Ns=7, 與Ns,wn=buttord(wp,ws,Rp,Rs,s);得出的結(jié)果一樣。2： 雙線性變換公式： （2）其中，是模擬的頻率，是數(shù)字的歸一化頻率，T是積分步長。用公式（2）預(yù)畸所求數(shù)字濾波器的的數(shù)字頻率指標(biāo)（數(shù)字濾波器通帶截止頻率0.2，阻帶截止頻率0.6），得到一個等價的模擬低通濾波器的頻率指標(biāo)。3： 得到原型模擬低通濾波器的指標(biāo)后作頻率響應(yīng)，從而得模擬的傳輸函數(shù)H(s)，然后對傳輸函數(shù)進(jìn)行雙線性變換，從而得到數(shù)字濾波器的傳輸函數(shù)H(z

8、).結(jié)果分析：圖3.模擬原型與數(shù)字濾波器幅度相應(yīng)比較圖4. 模擬原型與數(shù)字濾波器相位相應(yīng)比較由圖1可以看出模擬原型低通濾波器的通帶截止頻率為10KHz對應(yīng)著數(shù)字低通濾波器的截止頻率0.2 rad與雙線性變換結(jié)果一致。由圖2可以看出，相位響應(yīng)就有很大差別了，這是由雙線性變換的頻率畸變所致。2.3 FIR濾波器的結(jié)構(gòu)與設(shè)計對于IIR濾波器必須明確要求推導(dǎo)出來的傳輸函數(shù)是穩(wěn)定的，而另一方面，對于FIR濾波器的設(shè)計，由于濾波器的傳輸函數(shù)是以1 /z的多項式表示的，所以FIR數(shù)字濾波器都是穩(wěn)定的11。利用傅里葉加窗級數(shù)法設(shè)計線性相位低通FIR濾波器（沒有使用到fir1函數(shù)）設(shè)計指標(biāo)：通帶截止頻率4rad

9、/s，阻帶截止頻率6rad/s，最大通帶衰減為0.2dB,最小阻帶衰減為42dB，抽樣率為18rad/s.圖5.由截短得到的理想FIR低通濾波器的沖擊響應(yīng)圖6.各種窗函數(shù)設(shè)計的低通濾波器結(jié)果分析：Hamming窗函數(shù)的長度是30通帶波紋是Rp = 0.0394dB最小阻帶衰減時As = 52dBHanning窗函數(shù)的長度是29通帶波紋是Rp= 0.0711dB最小阻帶衰減時As = 44dBBlackman窗函數(shù)的長度是52通帶波紋是Rp= 0.0027dB最小阻帶衰減時As = 75dB均滿足濾波器的性能指標(biāo)。結(jié)論在電子系統(tǒng)中，具有嚴(yán)格的線性相位特性的FIR 數(shù)字濾波器被廣泛使用，是信號處理

10、的基本構(gòu)成之一。利用matlab的強(qiáng)大運(yùn)算功能，基于matlab信號處理工具箱（signal processing toolbox）的數(shù)字濾波器設(shè)計法可以快速有效的設(shè)計由軟件組成的常規(guī)數(shù)字濾波器，設(shè)計方便、快捷，極大的減輕了工作量。在設(shè)計過程中可以對比濾波器特性，隨時更改程序參數(shù)，以達(dá)到濾波器設(shè)計的最優(yōu)化。利用matlab設(shè)計數(shù)字濾波器在信號處理軟件和微機(jī)保護(hù)中，有著廣泛的應(yīng)用前景。本文通過調(diào)用simulink中的功能模塊構(gòu)成數(shù)字濾波器的仿真框圖，編寫matlab語言進(jìn)行仿真，從結(jié)果可以看出各種濾波器的性能并反映實際的情況。參考文獻(xiàn) 1 李惠瓊，孟令軍，王宏濤.基于Nios 的IIR數(shù)字濾波器

11、的設(shè)計J.通信技術(shù),2009,42(08)2 王秀敏,汪毓鐸,張洋,楊世華. 通信系統(tǒng)中FIR數(shù)字濾波器的設(shè)計研究J. 通信技術(shù), 2009,(09) . 3梁辰 基于MATLAB的FIR數(shù)字濾波器的設(shè)計J. 機(jī)械設(shè)計與制造，2010，（12）4 史燕,楊小雪. 基于Matlab和Multisim的綜合性實驗橢圓濾波器設(shè)計與仿真J. 北華航天工業(yè)學(xué)院學(xué)報, 2008,(S1) . 5 王靖,李永全. 數(shù)字橢圓濾波器的Matlab設(shè)計與實現(xiàn)J. 現(xiàn)代電子技術(shù), 2007,(06) . 6 劉抒珍,童子權(quán),任麗軍,劉小紅. DDS波形合成技術(shù)中低通橢圓濾波器的設(shè)計J. 哈爾濱理工大學(xué)學(xué)報, 200

12、4,(05) .7 楊曉琳. 循環(huán)碼理論及其譯碼算法研究設(shè)計距離為n 的二元BCH 碼構(gòu)造及其B-M 迭代譯碼算法實現(xiàn)D.成都：成都理工大學(xué),2008.8 王宇. BCH 編碼在GPS 探空儀中的應(yīng)用J.信息安全與通信保密,2010，(07)9 趙亞梅, 杜紅棉, 張志杰. 基于MATLAB一種IIR數(shù)字帶通濾波器的設(shè)計與仿真J. 微計算機(jī)信息, 2007,(13) 10楊世華,王秀敏,陳豪威;基于DSP Builder和FPGA的IIR濾波器設(shè)計J,2010,(12)11 LAWRENCE.R.R FIR Digital Filter Design Techniques Using Weig

13、hted Chebyshev Approximation PROCEEDINGS OF THE IEEE, VOL. 63, NO. 4, APRIL 1975.附仿真源程序：程序1：y1,Fs,bits=wavread(engine.wav,2048);%Fs是音頻采樣率，bits是每個采樣點包含的bit數(shù)目sound(y1,Fs,bits);hold onY1=fft(y1,4096);f=Fs*(0:4095)/4096;y2=y1+0.2*randn(2048,1);figure(1)subplot(211)plot(y2);xlabel (時間序號n); ylabel (振幅);le

14、gend(疊加噪聲后信號波形,1)axis (0 ,2100, -1 , 1);Y2=fft(y2,4096);subplot(212)plot(f,fftshift(abs(Y2);xlabel (信號頻率Hz); ylabel (振幅);legend(疊加噪聲后信號頻譜,1)axis (4000,7000, 0 , 150);Ap=1;%通帶波紋As=100;%最小阻帶衰減fp=200;fs=300;Ft=20000;wp=2*pi*fp/Ft;%歸一化通帶截止角頻率ws=2*pi*fs/Ft;%歸一化阻帶截止角頻率n,wn=ellipord(wp,ws,Ap,As);num,den=el

15、lip(n,Ap,As,wn);y3=filter(num,den,y2);Y3=fft(y3,4096);figure(2)subplot(211)plot(y1,r-);hold onplot(y3,g-);xlabel (時間序號n); ylabel (振幅);legend(原始信號波形,濾波后信號波形,1)axis (0 ,2100, -1 , 1);subplot(212)plot(f,fftshift(abs(Y1),r-);hold onplot(f,fftshift(abs(Y3),g-);xlabel (信號頻率Hz); ylabel (振幅);legend(原始信號頻譜,濾

16、波后信號頻譜,1)axis (4000,7000, 0 , 150);程序2：%用雙線性變換設(shè)計一個0.4dB截止頻率為10KHz且在30KHz處有最小阻帶衰減為50dB的數(shù)字巴特沃茲低通濾波器，其抽樣率為100KHz；clear all;fp=10000;fs=30000;Fs=100000;Rp=0.4;Rs=50;wp=2*pi*fp;ws=2*pi*fs;Ns,wn=buttord(wp,ws,Rp,Rs,s);bs,as=butter(Ns,wn,s);%計算模擬原型濾波器的頻率響應(yīng)w = 0: 200: 80000*pi;h = freqs(bs,as,w);%預(yù)畸所求數(shù)字濾波器的

17、數(shù)字頻率指標(biāo)，得到一個等價的的模擬濾波器的指標(biāo)Wp=tan(wp/(2*Fs)Ws=tan(ws/(2*Fs)Nz,Wn=buttord(Wp,Ws,Rp,Rs,s);b,a=butter(Nz,Wn,s);%進(jìn)行雙線性變換，得到所求的數(shù)字IIR傳輸函數(shù)Bz,Az=bilinear(b,a,0.5);H,omega = freqz(Bz,Az,1024);figure(1)subplot (211);plot (w/(2*pi*1000),20*log10(abs(h);xlabel(Frequency, KHz); ylabel(Gain, dB);title(模擬原型濾波器幅度響應(yīng))sub

18、plot (212);plot (omega/pi,20*log10(abs(H);xlabel(omega/pi); ylabel(Gain, dB);title(數(shù)字濾波器幅度響應(yīng))figure(2)subplot (211);phase=180*angle(h)/pi;plot (w/(2*pi*1000),phase);xlabel(Frequency, KHz); ylabel(phase,);title(模擬原型濾波器相位響應(yīng))subplot (212);Phase=180*angle(H)/pi;plot(omega/pi,Phase);xlabel(omega/pi); yla

19、bel(Phase,);title(數(shù)字濾波器相位響應(yīng))程序3：%利用窗函數(shù)發(fā)設(shè)計線性相位低通FIR濾波器wp=2*pi*4/18; %歸一化通帶截止頻率 ws = 2*pi*6/18; %歸一化阻帶截止頻率 tr_width =ws-wp; %確定過度帶寬N1=ceil (6.64*pi/tr_width)+1;% 確定濾波器階數(shù) ，ceil是向上取整函數(shù)N2=ceil (6.22*pi/tr_width)+1;% 確定濾波器階數(shù) ，ceil是向上取整函數(shù)N3=ceil (11.12*pi/tr_width)+1;% 確定濾波器階數(shù) ，ceil是向上取整函數(shù)n1 =0:1:N1-1;n2 =

20、0:1:N2-1;n3 =0:1:N3-1;wc=(ws+wp)/2; %理想低通的截止頻率hd1=ideal_lp(wc, N1);%調(diào)用理想低通濾波器的單位沖擊響應(yīng)函數(shù)hd2=ideal_lp(wc, N2);%調(diào)用理想低通濾波器的單位沖擊響應(yīng)函數(shù)hd3=ideal_lp(wc, N3);%調(diào)用理想低通濾波器的單位沖擊響應(yīng)函數(shù)W1=(hamming(N1);W2=hanning(N2);W3=blackman(N3);h1= hd1.*W1;h2= hd2.*W2;h3= hd3.*W3;db1,mag1,pha1,w1=freqz_m(h1,1);%計算hamming窗低通濾波器的幅度響應(yīng)db2,mag2,pha2,w2=freqz_m(h2,1);%計算hanning窗低通濾波器的幅度響應(yīng)db3,mag3,pha3,w3=freqz_m(h3,1);%計算blackman窗低通濾波器的幅度響應(yīng)delta_w=2*pi/18;Rp1= - (min(db1(1:1: wp/delta_w +1) %hamming窗實際通帶波動As1= - round(max(db1(ws/delta_w +1: 1: 501 )% hamming窗最小阻帶衰減Rp2= - (min(db2(1:1: wp/delta_w +1) %hanni