基于matlab的fir數(shù)字濾波器

1、基于matlab的fir數(shù)字濾波器設(shè)計(jì)基于matlab的fir數(shù)字濾波器sea摘 要：介紹了應(yīng)用Matlab語言設(shè)計(jì)FIR數(shù)字濾波器時(shí)采用直接程序設(shè)計(jì)法。同時(shí)介紹了FIR數(shù)字濾波器幾種設(shè)計(jì)方法的函數(shù)調(diào)用格式；通過實(shí)例，給出了程序設(shè)計(jì)法的詳細(xì)步驟，并在Matlab的Simulink環(huán)境下，對(duì)所設(shè)計(jì)的濾波器進(jìn)行了仿真。關(guān)鍵詞：數(shù)字濾波器；Matlab引言1. MATLAB的特點(diǎn)MATLAB(Matrix Laboratory)機(jī)矩陣實(shí)驗(yàn)室是由美國MathWorks公司推出的一款集數(shù)值分析、矩陣計(jì)算、信號(hào)處理和圖形顯示于一體的工具軟件。它的命令語句功能非常強(qiáng)大,包含了大量高度集成的可直接調(diào)用的函數(shù)，

2、高效簡潔；另一方面，它又是一個(gè)開放系統(tǒng)，針對(duì)不同的學(xué)科，推出了不同的工具箱。自1984年推向市場以來，經(jīng)過十幾年的發(fā)展和競爭，現(xiàn)已成為國際認(rèn)可(IEEE)的最優(yōu)化的科技應(yīng)用軟件一。正是由于MATLAB具有良好的擴(kuò)展性以及強(qiáng)大的數(shù)據(jù)分析和處理能力，現(xiàn)已廣泛應(yīng)用于矩陣代數(shù)、數(shù)值計(jì)算、數(shù)字信號(hào)處理、振動(dòng)理論、神經(jīng)網(wǎng)絡(luò)控制、動(dòng)態(tài)仿真等領(lǐng)域。MATLAB是用M語言編程，盡管不能在M文件中直接調(diào)用C語言程序，但可以通過MATLAB提供的應(yīng)用編程接口(API)來與外部接口，在MATLAB環(huán)境中調(diào)用C語言或Fortran程序、輸入或輸出數(shù)據(jù)以及與其他軟件程序建立客戶/服務(wù)器關(guān)系。在MATLAB中調(diào)用C語言程序

3、，必須通過MEX文件來實(shí)現(xiàn)。2. FIR數(shù)字濾波器濾波就是有選擇性地提取或去掉(或削弱)某一段或某幾段頻率范圍內(nèi)的信號(hào),數(shù)字濾波器是一種用來過濾時(shí)間離散信號(hào)的數(shù)字系統(tǒng),它是通過對(duì)抽樣數(shù)據(jù)進(jìn)行數(shù)學(xué)處理來達(dá)到選頻目的。數(shù)字濾波器可分為IIR (無限沖激響應(yīng))和FIR (有限沖激響應(yīng))兩種結(jié)構(gòu)。隨著計(jì)算機(jī)技術(shù)和集成電路技術(shù)的發(fā)展，數(shù)字信號(hào)處理以其方便、靈活的特點(diǎn)，越來越引起人們的重視。數(shù)字濾波器是數(shù)字信號(hào)處理的重要內(nèi)容，數(shù)字濾波器的設(shè)計(jì)已成為數(shù)字信號(hào)處理研究中的熱點(diǎn)一。應(yīng)用Matlab語言可以快捷地設(shè)計(jì)出由軟件組成的數(shù)字濾波器,很容易通過對(duì)參數(shù)的修改進(jìn)行性能的優(yōu)化。FIR濾波器(即有限長沖擊響應(yīng)濾波

4、器, Finite Imp ulse Response Digital Filter)最大的優(yōu)點(diǎn)就是在滿足幅頻特性的同時(shí)，還可以獲得嚴(yán)格的線性相位特性，這使它在語音處理、圖像處理等要求高保真的數(shù)字信號(hào)處理中顯得十分重要。 FIR數(shù)字濾波器的設(shè)計(jì)原理1. FIR濾波器設(shè)計(jì)原理由于FIR 濾波器沖擊響應(yīng)hn是有限長序列，因此這種結(jié)構(gòu)可用非遞歸結(jié)構(gòu)來實(shí)現(xiàn)1。FIR 數(shù)字濾波器系統(tǒng)函數(shù)一般形式為：。FIR 濾波器數(shù)學(xué)表達(dá)式可用差分方程來表示：。式中:y（n）輸出序列；h（k）濾波器系數(shù)；n濾波器階數(shù)；x（k）輸入序列。應(yīng)用Matlab 設(shè)計(jì)FIR濾波器的主要任務(wù)就是根據(jù)給定的性能指標(biāo)，設(shè)計(jì)一個(gè)H（z）

5、 ,使其逼近這一指標(biāo)，進(jìn)而計(jì)算并確定濾波器的系數(shù)b(n)，再將所設(shè)計(jì)濾波器的幅頻響應(yīng)、相頻響應(yīng)曲線作為輸出,與設(shè)計(jì)要求進(jìn)行比較，對(duì)設(shè)計(jì)的濾波器進(jìn)行優(yōu)化。FIR濾波器的主要設(shè)計(jì)方法有窗函數(shù)法、最優(yōu)化設(shè)計(jì)法及約束最小二乘逼近法2。在濾波器傳統(tǒng)設(shè)計(jì)中，要得到其幅頻和相頻響應(yīng)特性，需要根據(jù)這些方法進(jìn)行大量的計(jì)算，這使得濾波器的設(shè)計(jì)緩慢，周期變長，不利于設(shè)計(jì)的優(yōu)化。Matlab信號(hào)處理工具箱中提供了基于濾波器設(shè)計(jì)方法的工具函數(shù)，編程中可根據(jù)設(shè)計(jì)要求直接調(diào)用相應(yīng)的函數(shù)，方便快捷，見表1 。表1 FIR數(shù)字濾波器設(shè)計(jì)方法及函數(shù)調(diào)用格式設(shè)計(jì)方法調(diào)用方法調(diào)用格式說明窗函數(shù)法Fir1Fir2b=fir1(n,Wn

6、)b=fir2（n，f，m）n為階數(shù)；Wn為截止頻率f，m為期望幅頻向量和幅值向量最優(yōu)化設(shè)計(jì)法FirlsRemezRemezordb=firls(n,f,a)bremez(n,f,a,w) n ,fo ,ao ,w =remezord (f ,a ,dev)兩者僅算法不同，f為頻率點(diǎn)向量，n為指定頻率點(diǎn)幅度響應(yīng)，w為權(quán)系數(shù)Fo為歸一化頻率邊界，ao為頻帶內(nèi)幅值，w為權(quán)向量最小二乘法FirclsFircls1b=fircls(n,f,a,up,lo)b=fircls1(n,wo,dp,ds)up,lo為每個(gè)頻帶上邊界和下邊界頻率，f，a為期望幅頻特性的頻率向量和幅值向量。Wo為截至頻率，dp為離

7、幅值1的最大偏差，ds為阻帶離幅值0的最大偏差。2. 用窗函數(shù)法設(shè)計(jì)FIR數(shù)字濾波器窗函數(shù)設(shè)計(jì)3技術(shù)是FIR濾波器設(shè)計(jì)的主要方法之一，由于其運(yùn)算簡單、精度高，已成為工程中應(yīng)用最廣泛的設(shè)計(jì)方法。設(shè)計(jì)濾波器的基本思想是讓待設(shè)計(jì)的實(shí)際濾波器逼近理想特性。理想低通濾波器的頻率特性應(yīng)為：振幅特性在通帶內(nèi)為1，阻帶內(nèi)為0；在通帶內(nèi)的相位特性與w成線性關(guān)系，即它對(duì)應(yīng)的理想單位沖激相應(yīng)hd(n)為可見，理想低通濾波器的單位沖激響應(yīng)hd (n)是無限長的非因果序列。而我們要設(shè)計(jì)的是FIR濾波器，其h( n)必然是有限長的，為了構(gòu)造物理上可以實(shí)現(xiàn)的長度為N 的因果線性相位濾波器，必須將hd( n) 截取長度為N

8、的一段，或者說用一個(gè)有限長度的窗口函數(shù)序列w ( n) 來截取hd ( n)，即h( n) = w ( n) hd ( n)，因而窗函數(shù)序列的形狀及長度的選擇就是關(guān)鍵。工程實(shí)際中常用的窗函數(shù)有6種，即矩形窗、三角形窗、漢寧(Hanning)窗、哈明(Hamming)窗、布萊克曼(Blackman) 窗和凱澤( Kaiser) 窗。它們之間的性能比較如表2所示。實(shí)際的FIR數(shù)字濾波器的頻率響應(yīng)H() 逼近理想濾波器頻率響應(yīng)Hd()的好壞，取決于窗函數(shù)的頻率特性W( w)。表2 6種窗函數(shù)的特性窗函數(shù)旁瓣峰值/dB近似過渡帶寬精確過渡帶寬阻帶最小衰減/dB矩形窗-134/N1.8/N21三角形窗-

9、258/N6.1/N25漢寧窗-318/N6.2/N44哈明窗-418/N6.6/N53布萊克曼窗-5712/N11/N74凱澤窗（）-5710/N803. 窗函數(shù)法設(shè)計(jì)步驟3.1用窗函數(shù)設(shè)計(jì)FIR濾波器的步驟為: (1)選擇窗函數(shù)類型能夠和長度，寫出窗函數(shù)w（n）表達(dá)式。根據(jù)阻帶最小衰減選擇窗函數(shù)w(n)的類型，再根據(jù)過渡帶寬確定所選窗函數(shù)的長度N。(2)構(gòu)造希望逼近的頻率響應(yīng)函數(shù)。根據(jù)設(shè)計(jì)需要，一般選擇相應(yīng)的線性相位理想濾波器（理想低通，理想高通，理想帶通，理想帶阻）。應(yīng)當(dāng)注意，理想濾波器的截止頻率wc點(diǎn)（對(duì)低通濾波器）近似為最終設(shè)計(jì)的FIRDFD的-6dB頻率。(3)計(jì)算：(4)加窗得到

10、設(shè)計(jì)結(jié)果：3.2實(shí)際設(shè)計(jì)時(shí)，一般采用MATLAB工具箱函數(shù)實(shí)現(xiàn)。步驟（1）由設(shè)計(jì)者根據(jù)設(shè)計(jì)指標(biāo)完成，步驟（2）-（4）的解題過程可調(diào)用設(shè)計(jì)函數(shù)fir1來實(shí)現(xiàn)。Fir1時(shí)用窗函數(shù)法設(shè)計(jì)線性相位FIRDF的工具箱函數(shù)，調(diào)用格式如下：hn=fir1(N,wc,ftype,window) fir實(shí)現(xiàn)線性相位FIR濾波器的標(biāo)準(zhǔn)窗函數(shù)法設(shè)計(jì)。“標(biāo)準(zhǔn)”時(shí)指再設(shè)計(jì)低通、高通、帶通、帶阻FIR濾波器時(shí)，分別取相應(yīng)的理想低通、高通、帶通、帶阻濾波器，故而設(shè)計(jì)的濾波器的頻率響應(yīng)稱為標(biāo)注頻率響應(yīng)。hn=fir1(N,wc)可得到6dB截至頻率為wc的N階（單位脈沖響應(yīng)h(n)長度為N+1）FIR低通濾波器，默認(rèn)（缺省

11、參數(shù)windows）選用hammiing窗。其單位脈沖響應(yīng)h(n)為：h(n)=h(n+1),n=0,1,2,N而且滿足線性相位條件：h(n)=h(N-1-n)其中wc為對(duì)歸一化的數(shù)字頻率，0wc1。hn=fir1(N,wc,ftype)可設(shè)計(jì)高通和帶阻濾波器。·當(dāng)ftype=high時(shí)，設(shè)計(jì)高通FIR濾波器；·當(dāng)ftype=stop時(shí)，設(shè)計(jì)帶阻FIR濾波器。hn=fir1(N,wc,ftype,window)通過選擇wc、ftype和window參數(shù)，可以設(shè)計(jì)各種加窗濾波器。Fir2可以指定任何形狀的Hd（），用help命令查閱其調(diào)用格式。4. 用Matlab對(duì)FIR數(shù)字

12、濾波器進(jìn)行仿真設(shè)計(jì)4.1 編寫流程4我們以低通數(shù)字濾波器設(shè)計(jì)為例，來介紹用MATLAB對(duì)FIR數(shù)字濾波器進(jìn)行仿真過程：設(shè)計(jì)55階低通濾波器，截止頻率在0.3，并用濾波器對(duì)信號(hào)sin (2*pi*15*t)+0.5*sin(2*pi*90*t)+0.2*sin(2*pi*300*t)濾波(信號(hào)采樣頻率為600Hz)。1）用RLS設(shè)計(jì)數(shù)字低通濾波器B=fircls1(n,wo,dp,ds)函數(shù)返回長度為n+1 線性相位低通FIR濾波器，截止頻率為wo，取值0-1之間。通帶幅度偏離1.0最大值為dp，阻帶偏離0最大值為ds，n為濾波器階數(shù)。2）建立新的M文件選擇“File”菜單項(xiàng)，單擊“New M-

13、file”按鈕，便建立一個(gè)M文件5。3）編寫程序,調(diào)試程序,在MATLAB 平臺(tái)上得出實(shí)驗(yàn)結(jié)果編寫完程序，單擊“Debug”按鈕，即可觀察到所設(shè)計(jì)濾波器幅頻特性和濾波前后的波形。4.2 實(shí)例下面給出兩個(gè)實(shí)例說明用MATLAB設(shè)計(jì)FIR濾波器的優(yōu)點(diǎn)，及濾波器的濾波特性。先用最優(yōu)化設(shè)計(jì)法設(shè)計(jì)一個(gè)濾波器：例1： 設(shè)計(jì)一個(gè)最小階數(shù)的低通濾波器，采樣頻率fs=2000 Hz，通帶截止頻率為500Hz，阻帶的截至頻率為600Hz，阻帶最小衰減為40dB，通帶的最大衰減為3dB。在設(shè)計(jì)之前應(yīng)先確定用說明方法設(shè)計(jì)，本例可選擇等波紋的最優(yōu)化設(shè)計(jì)法。程序如下：fs=2000; %采樣頻率rp=3; %通帶波紋rs

14、=40; %阻帶波紋f=500 600; %截止頻率a=1 0; %期望的幅度dev=(10(rp/20)-1)/(10(rp/20)+1) 10(-rs/20);n,fo,ao,w=remezord(f,a,dev,fs);b=remez(n,fo,ao,w); %調(diào)用最優(yōu)設(shè)計(jì)法中remez 函數(shù)freqz(b,1,1024,fs)程序運(yùn)行后，計(jì)算機(jī)輸出該濾波器的幅頻及相頻響應(yīng)特性，如圖1所示。圖1 濾波器輸出的幅頻及相頻響應(yīng)特性上圖中幅頻特性曲線從500Hz開始向下折，這是因?yàn)轭}設(shè)中通帶截止頻率為500Hz，同樣圖中幅頻曲線在600Hz降為最低，是因?yàn)轭}設(shè)中阻帶的截至頻率為600Hz。用窗

15、函數(shù)設(shè)計(jì)一個(gè)濾波器：例2： 設(shè)計(jì)一個(gè)低通濾波器，性能指標(biāo)為通帶01000Hz，阻帶截止頻率為2000Hz，通帶波動(dòng)1%，阻帶波動(dòng)5%，采樣頻率為10000Hz，采用Kaiser窗。主要程序如下：fc=1000;fs=10000;n,Wn,beta,ftype=kaiserord(1000 2000,1 0,0.05 0.05 ,fs); %得出濾波器的階數(shù)w=2*fc/fs %模擬到數(shù)字濾波器的技術(shù)指標(biāo)的轉(zhuǎn)換window=kaiser(n+1,beta); %使用kaiser窗函數(shù)b=fir1(n,w,window); %使用標(biāo)準(zhǔn)頻率響應(yīng)的加窗設(shè)計(jì)函數(shù)fir1freqz(b,1,512); %

16、數(shù)字濾波器頻率響應(yīng);t=(1:100)/fs;s=sin(2*pi*t*800)+sin(2*pi*t*3000)+sin(2*pi*t*4000); %混和正弦波信號(hào)sf=filter(b,1,s); %對(duì)信號(hào)s進(jìn)行濾波濾波器長度為13，歸一化截止頻率為0.2，凱塞窗控制旁瓣的參數(shù)為1.5099，以及傳遞給函數(shù)firl用于指定濾波器類型的ftype為low，通過圖2-圖4可以看出它能滿足設(shè)計(jì)要求。 圖2 濾波器的幅頻與相頻特性圖3 濾波前混和正弦波信號(hào)的波形圖4 濾波后的波形 用最小二乘法設(shè)計(jì)一個(gè)fir濾波器例3：有一個(gè)采樣頻率為1000Hz的信號(hào)源，sin(30t)+0.5*sin(180

17、t)+0.2*sin(600t)，應(yīng)用約束最小二乘法設(shè)計(jì)一個(gè)帶通濾波器，通帶為0.2 0.4，對(duì)該信號(hào)進(jìn)行濾波。程序如下： t=0:1/1000:1; sig=sin(2*pi*15*t)+0.5*sin(2*pi*90*t)+0.2*sin(2*pi*300*t); plot(t,sig); n=50; f=0 0.2 0.4 1; a=0 1 0; up=0.02 1.02 0.01; lo=-0.02 0.98 -0.01; b=fircls(n,f,a,up,lo); newsig=fftfilt(b,sig); figure ft=t(301:350); ns=newsig(301:

18、350); zns=interp(ns,12); znt=interp(ft,12); plot(znt,zns);該題濾波器的幅頻特性曲線和相頻特性曲線及濾波前后的波形如圖5、圖6和圖7所示： 圖5 濾波器的幅頻特性曲線和相頻特性曲線。圖6 濾波前的波形圖7 濾波后的波形上文分別介紹了fir濾波器設(shè)計(jì)的三種方法：窗函數(shù)法、最優(yōu)化設(shè)計(jì)法和最小二乘法設(shè)計(jì)fir濾波器的程序步驟，并對(duì)程序進(jìn)行了仿真，通過波形圖，展示了濾波器的特性曲線和很好的仿真效果。5 結(jié)論：FIR濾波器的應(yīng)用十分廣泛，當(dāng)今的許多信號(hào)處理系統(tǒng)，如圖像處理等系統(tǒng)要求具有線性相位特性。在這方面FIR濾波器有其獨(dú)特的優(yōu)點(diǎn)，能很容易地設(shè)計(jì)出具有嚴(yán)格線性特性的FIR系統(tǒng)。此外，F(xiàn)IR濾波器的沖激響應(yīng)是有限長序列，其系統(tǒng)函數(shù)為一個(gè)多項(xiàng)式，它所含的極點(diǎn)多為原點(diǎn)，所以FIR濾波器永遠(yuǎn)是穩(wěn)定的。綜上所述，在Matlab中進(jìn)行數(shù)字濾