數(shù)字濾波器設(shè)計

1、 1.1 第一章 緒論引言 隨著信息時代的到來，數(shù)字信號處理的理論和技術(shù)更加成熟，它應(yīng)用廣泛，滲透到各個重要的科學(xué)領(lǐng)域。作為數(shù)字信號處理的一項重大突破就是數(shù)字濾波技術(shù)。數(shù)字濾波在通信，語音編碼，雷達等許多領(lǐng)域中有著十分廣泛的應(yīng)用，目前，數(shù)字濾波器的設(shè)計圖像處理，數(shù)據(jù)壓縮等方面的應(yīng)用取得了令人矚目的進展和成就，近年來迅速發(fā)展起來的小波理論，由于其局部分析性能的友誼，在圖像處理中的應(yīng)用研究，尤其是在圖像壓縮，圖像去噪等方面的應(yīng)用研究，收到了越來越多的關(guān)注。由于數(shù)字濾波器的概念比較抽象，加上其數(shù)值計算又比較繁瑣，所以借助好的計算機軟件來進行輔助設(shè)計，是數(shù)字濾波器研究領(lǐng)域的一個發(fā)展趨勢。這樣的軟件有很

2、多種，其中最具代表性的就是MATLAB。MATLAB是美國MathWorks公司自20世紀80年代中期推出的數(shù)學(xué)軟件，優(yōu)秀的數(shù)值計算能力和卓越的數(shù)據(jù)可視化能力使其很快在數(shù)學(xué)軟件中脫穎而出。MATLAB是“矩陣實驗室”（MATrix LABoratoy）的縮寫，它是由美國Mathworks公司于1984年正式推出的，是一種以矩陣運算為基礎(chǔ)的交互式程序語言，專門針對科學(xué)、工程計算及繪圖的需求。MATLAB是功能強大的科學(xué)及工程計算軟件，它不但具有以矩陣計算為基礎(chǔ)的強大數(shù)學(xué)計算和分析功能，而且還具有豐富的可視化圖形表現(xiàn)功能和方便的程序計算能力。MATLAB的應(yīng)用領(lǐng)域極為廣泛，除數(shù)學(xué)計算和分析外，還被

3、廣泛地應(yīng)用于自動控制、系統(tǒng)仿真、數(shù)字信號處理、圖形圖像分析、數(shù)理統(tǒng)計、人工智能、虛擬現(xiàn)實技術(shù)、通信工程、金融系統(tǒng)等領(lǐng)域，因此，MATLAB是面向21世紀的計算機程序設(shè)計及科學(xué)計算語言。 Simulink是MATLAB最重要的組件之一，它提供一個動態(tài)系統(tǒng)建模、仿真和綜合分析的集成環(huán)境。在該環(huán)境中，無需大量書寫程序，而只需要通過簡單直觀的鼠標操作，就可構(gòu)造出復(fù)雜的系統(tǒng)。Simulink具有適應(yīng)面廣、結(jié)構(gòu)和流程清晰及仿真精細、貼近實際、效率高、靈活等優(yōu)點，并基于以上優(yōu)點Simulink已被廣泛應(yīng)用于控制理論和數(shù)字信號處理的復(fù)雜仿真和設(shè)計。同時有大量的第三方軟件和硬件可應(yīng)用于或被要求應(yīng)用于Simuli

4、nk。MATLAB/Simulink實現(xiàn)數(shù)字濾波器的好處是效率高，直觀而且容易修改。有高性能數(shù)值計算的高級算法，特別適合矩陣代數(shù)領(lǐng)域；有大量事先定義的數(shù)學(xué)函數(shù)，并且有很強的用戶自定義函數(shù)的能力；有強大的繪圖功能以及具有教育、科學(xué)和藝術(shù)學(xué)的圖解和可視化的二維、三維圖；基于HTML的完整的幫助功能；適合個人應(yīng)用的強有力的面向矩陣(向量)的高級程序設(shè)計語言；與其它語言編寫的程序結(jié)合和輸入輸出格式化數(shù)據(jù)的能力； 數(shù)字濾波器是指完成信號濾波處理功能的，用有限精度算法實現(xiàn)的離散線性非時變系統(tǒng)，其輸入是一組（由模擬信號取樣和量化的)數(shù)字量，其輸出是經(jīng)過變換或說處理的另一組數(shù)字量。因此，它本身既可以是用數(shù)字硬

5、件裝配成德一臺完成給定運算的專用數(shù)字計算機，也可將所需的運算編成程序，讓通用計算機來執(zhí)行。第二章 數(shù)字濾波器數(shù)字濾波器是指輸入輸出均為數(shù)字信號，通過一定運算關(guān)系改變輸入信號所含頻率成分的相對比例或者濾除某些頻率成分的器件。它工作在數(shù)字信號域，它處理的對象是經(jīng)由采樣器件將模擬信號轉(zhuǎn)換而得到的數(shù)字信號。數(shù)字濾波器一般由寄存器、延時器、加法器和乘法器等基本數(shù)字電路實現(xiàn)。數(shù)字濾波器一詞出現(xiàn)在60年代中期。由于電子計算機技術(shù)和大規(guī)模集成電路的發(fā)展，數(shù)字濾波器已可用計算機軟件實現(xiàn)，也可用大規(guī)模集成數(shù)字硬件實時實現(xiàn)。數(shù)字濾波器具有穩(wěn)定性高、精度高、靈活性大等突出優(yōu)點。隨著數(shù)字技術(shù)的發(fā)展，用數(shù)字技術(shù)實現(xiàn)濾波器

6、的功能愈來愈受到人們的注意和廣泛使用。這里所說的數(shù)字濾波器是指理想帶通，低通等的頻率選擇數(shù)字濾波器。數(shù)字濾波器設(shè)計的一個重要步驟是確定一個可實現(xiàn)的傳輸函數(shù)H(z)，這個確定傳輸函數(shù)H(z)的過程稱為數(shù)字濾波器設(shè)計。數(shù)字濾波器的一般設(shè)計過程為：（1）按照實際需要，確定濾波器的性能要求（通常在頻域內(nèi)給定數(shù)字濾波的性能要求）。（2）尋找一滿足預(yù)定性能要求的離散時間線性系統(tǒng)。（3）用有限精度的運算實現(xiàn)所設(shè)計的系統(tǒng)。（4）通過模擬，驗證所設(shè)計的系統(tǒng)是否符合給定性能要求。數(shù)字濾波器是對數(shù)字信號實現(xiàn)濾波的線性時不變系統(tǒng)。數(shù)字濾波實質(zhì)上是一種運算過程，實現(xiàn)對信號的運算處理。輸入數(shù)字信號（數(shù)字序列）通過特定的運

7、算轉(zhuǎn)變?yōu)檩敵龅臄?shù)字序列，因此，數(shù)字濾波器本質(zhì)上是一個完成特定運算的數(shù)字計算過程，也可以理解為是一臺計算機。描述離散系統(tǒng)輸出與輸入關(guān)系的卷積和差分方程只是給數(shù)字信號濾波器提供運算規(guī)則，使其按照這個規(guī)則完成對輸入數(shù)據(jù)的處理。時域離散系統(tǒng)的頻域特性:,數(shù)字濾波器可以理解為是一個計算程序或算法，將代表輸入信號的數(shù)字時間序列轉(zhuǎn)化為代表輸出信號的數(shù)字時間序列，并在轉(zhuǎn)化的過程中，使信號按預(yù)定的形式變化。數(shù)字濾波器有多種分類，從數(shù)字濾波器功能上分可分為低通、高通、帶阻、帶通濾波器，數(shù)字濾波器有低通、高通、帶通、帶阻和全通等類型。它可以是時不變的或時變的、因果的或非因果的、線性的或非線性的。根據(jù)數(shù)字濾波器沖激響

8、應(yīng)的時域特征，可將數(shù)字濾波器分為兩種，即無限長沖激響應(yīng)濾波器（IIR）和有限長沖激響應(yīng)濾波器（FIR）。數(shù)字濾波器指標：一般來說，濾波器的幅頻特性是分段常數(shù)的，以低通為例，在通帶內(nèi)逼近于1，阻帶內(nèi)逼近與0，實際設(shè)計的濾波器并非是銳截止的通帶和阻帶兩個范圍，兩者之間總有一個過渡帶。在設(shè)計濾波器時事先給定幅頻特性允許誤差，在通帶范圍內(nèi)幅度響應(yīng)以誤差逼近于1，在阻帶內(nèi)幅度響應(yīng)以誤差逼近于0。 （1）式中wc和wr分別為通帶邊界頻率和阻帶邊界頻率，wr-wc為過渡帶。在具體的技術(shù)指標中往往用通帶波動來表示，用最小阻帶衰減At來表示，其具體的對應(yīng)公式這里就不詳述了。IIR數(shù)字濾波器設(shè)計IIR DF的沖激

9、響應(yīng)h(n)是無限長的，其輸入輸出的關(guān)系為： (2)系統(tǒng)函數(shù)為= (3)設(shè)計無限長單位脈沖響應(yīng)（IIR）數(shù)字濾波器一般可有三種方法。第一種方法，先設(shè)計一個合適的模擬濾波器，然后將其數(shù)字話，即將S平面映射到Z平面得到所需的數(shù)字濾波器。模擬濾波器的設(shè)計技巧非常成熟，不僅得到的是閉合形式的公式，而且設(shè)計系數(shù)已經(jīng)表格化了。因此，由模擬濾波器設(shè)計數(shù)字濾波器的方法準確，簡便，得到普遍采用。對于這種方法，工程上有兩種常見得變換法脈沖響應(yīng)不變法及雙線性變換法。第二種方法，在Z平面直接設(shè)計IIR數(shù)字濾波器，給出閉合形式的公式，或者以所希望的濾波器響應(yīng)作為依據(jù)，直接在Z平面上通過多次選定極點和零點的位置，以逼近該

10、響應(yīng)。第三種方法，利用最優(yōu)化技術(shù)設(shè)計參數(shù)，選定極點和零點在Z平面上的合適位置，在某種最優(yōu)化準則意義上逼近所希望的響應(yīng)。但一般不能得到濾波器的系數(shù)（即零，極點的位置）作為給定響應(yīng)的閉合形式函數(shù)表達式。優(yōu)化設(shè)計需要完成大量的迭代運算，這種設(shè)計法實際上也是IIR濾波器的直接設(shè)計。本文著重介紹由模擬濾波器設(shè)計相應(yīng)的IIR數(shù)字濾波器的方法。（1）脈沖響應(yīng)不變法 脈沖響應(yīng)不變法是使數(shù)字濾波器的單位脈沖響應(yīng)序列h(n)逼近模擬濾波器的沖激響應(yīng)，讓h(n)正好等于的采樣值。設(shè)已有滿足要求的模擬濾波器, 則可®®® 因為:,公式導(dǎo)出:具體轉(zhuǎn)換如下:設(shè)(以一階極點為例) (4)作拉氏

11、反變換,得采樣得作Z變換,得 (5)與極點關(guān)系為: (6) 一般對應(yīng)關(guān)系 (7)圖1 S平面到Z平面變換示意圖所以, 模擬系統(tǒng)穩(wěn)定因果®數(shù)字系統(tǒng)穩(wěn)定因果。按照脈沖響應(yīng)不變法，從S平面到Z平面的映射不是單值關(guān)系，而是先將在S平面沿虛軸作周期嚴拓，再按照映射關(guān)系將映射到Z平面，得到，因此，脈沖響應(yīng)不變法只適用于帶限的濾波器（如低通、帶通）。 在Matlab中利用M文件impinvar可以對模擬傳輸函數(shù)實行脈沖響應(yīng)不變法。（2）雙線性變換法 脈沖響應(yīng)不變法不適帶阻和高通濾波器的設(shè)計，因為高頻帶為通帶，前述方法易引起混頻。故希望:平面虛軸Ûz平面單位圓一周, 且應(yīng)有,因為, 所以選

12、變換 (8)其中可取任意正常數(shù), 后面將導(dǎo)出.設(shè)計思路: .®設(shè)計出模擬濾波器®轉(zhuǎn)化成數(shù)字濾波器. 圖2 數(shù)字域頻率與模擬域頻率的對應(yīng)關(guān)系轉(zhuǎn)化公式推導(dǎo)如下:因只關(guān)心頻率轉(zhuǎn)換, 故可設(shè), 則有, (稱為雙線性變換) (9)所以模擬濾波器轉(zhuǎn)換成數(shù)字濾波器的公式為 (10)由雙線性變換公式, 可得, 視為兩復(fù)平面變換, 再由.可得 (11)從上式可得:時, 上虛軸Û上單位圓周。時, 上左半平面Û上單位圓內(nèi)。時, 上右半平面Û上單位圓外。故若模擬濾波器穩(wěn)定，則雙線性變換后數(shù)字濾波器也穩(wěn)定。由于雙線性變換法是一種單值映射，因此消除了頻率混疊的現(xiàn)象。雙線性

13、變換法的缺點是模擬頻率與數(shù)字頻率間的非線性，這種非線性關(guān)系要求被變換的連續(xù)系統(tǒng)的幅度響應(yīng)是分段常數(shù)型的（某一頻率范圍內(nèi)幅度響應(yīng)近似于某一常數(shù)），不然所映射出的數(shù)字頻率響應(yīng)相對于原來的模擬頻率響應(yīng)會產(chǎn)生變形。為解決雙線性變換中的頻率非線性關(guān)系，我們采用預(yù)畸的方法，即，其中K=2/T。在Matlab中利用M文件bilinear可以對模擬傳輸函數(shù)實行雙線性變換法。 MATLAB中IIR數(shù)字濾波器的設(shè)計過程包括兩步。第一步，根據(jù)給定指標，確定濾波器的階數(shù)N和頻率縮放因子Wn。第二步，利用這些參數(shù)和給定的波紋參數(shù)，確定傳輸函數(shù)的關(guān)系。階數(shù)估計：利用雙線性變換法設(shè)計數(shù)字濾波器時，首先要對IIR數(shù)字濾波器的

14、階數(shù)進行估計，相應(yīng)的M文件為：buttord用于巴特沃斯濾波器，cheb1ord用于切比雪夫1型濾波器，cheb2ord用于切比雪夫2型濾波器，ellipord用于橢圓濾波器。濾波器的設(shè)計：對于基于雙線性變換法的IIR濾波器設(shè)計，對應(yīng)于四種逼近技術(shù)（即巴特沃斯、切比雪夫1型和2型及橢圓逼近），MATLAB工具箱中有相應(yīng)的函數(shù)。特別地可以用到下面的M文件：butter用于巴特沃斯濾波器的設(shè)計，cheby1用于切比雪夫1型濾波器的設(shè)計，cheby2用于切比雪夫2型濾波器的設(shè)計，ellip用于橢圓濾波器的設(shè)計。這些函數(shù)的輸出可以是濾波器傳輸函數(shù)分子和分母的系數(shù)向量，也可以是濾波器的零極點向量和標量增

15、益因子。同時，利用zp2tf可以由濾波器的零極點向量和標量增益因子得到傳輸函數(shù)分子和分母的系數(shù)向量。相應(yīng)地，利用函數(shù)zp2sos可以得到傳輸函數(shù)分子和分母系數(shù)向量的二次項因子。在計算出傳輸函數(shù)的系數(shù)之后，可以利用M文件freqz來計算頻率響應(yīng)。2.2 FIR數(shù)字濾波器設(shè)計FIR DF的沖激響應(yīng)h(n)是有限長的，M階FIR DF可以表示為： (12)其系統(tǒng)函數(shù)為: (13)與IIR數(shù)字濾波器的設(shè)計不同，F(xiàn)IR濾波器的設(shè)計與模擬濾波器的設(shè)計沒有任何聯(lián)系。因此，F(xiàn)IR濾波器的設(shè)計基于對指定幅度響應(yīng)的直接逼近，并通常要求其具有線性相位響應(yīng)。為了保證濾波器具有線性相位特性，濾波器系數(shù)必須滿足條件：h(

16、n)=h(M-1-n)。目前關(guān)于FIR濾波器的設(shè)計方法主要有三種，即窗函數(shù)法，頻率取樣法和切比雪夫等波紋逼近的最優(yōu)化設(shè)計方法。一般應(yīng)用較多的是第一種和第三種方法。這是因為窗函數(shù)法比較簡單，可應(yīng)用現(xiàn)成的窗函數(shù)公式，在技術(shù)指標要求不嚴格的情況下市比較靈活的。最優(yōu)化設(shè)計法必須借助計算機計算，但是它能得到最佳的等波紋的線性相位FIR濾波器。目前切比雪夫等波紋的線性相位FIR濾波器的計算機機助設(shè)計程序已經(jīng)比較完善，由于采用了REMEZ迭代算法，所以設(shè)計效率也很高，在應(yīng)用中越來越占優(yōu)勢。（1） 窗函數(shù)法 一般設(shè)計過程總是先給定一理想的濾波器頻率響應(yīng)，然后設(shè)計一個FIR濾波器，用它的頻率響應(yīng)來逼近理想的。這

17、種逼近中最直接的方法，是在時域中用FIR濾波器的單位脈沖響應(yīng)h(n)去逼近理想的單位脈沖響應(yīng)。因而，先由的IDTFT導(dǎo)出 (14)由于是矩形頻率特性，故一定是無限長的序列，且是非因果的。然而FIR濾波器是有限長的，所以用有限長的h(n)來逼近無限長的，最簡單的方法是截取中最重要的一段，將無限長的截取成長度為M的有限長序列，等效于再上施加了一個長度為M的矩形窗口，更為一般的，可以用一個長度為M的窗口函數(shù)w(n)來截取，即 (15)這一方法通常稱為窗函數(shù)法，窗口函數(shù)的形狀及長度M的選擇是窗函數(shù)法的關(guān)鍵。下面我們一低通為例，了解一下窗函數(shù)法的運用:提出希望頻率響應(yīng)函數(shù)（低通） 圖3 理想低通濾波器的

18、頻響線性相位, 具有片斷特點, 即算出(無限長) 圖4 理想低通的單位脈沖響應(yīng)（無限長的一部分）加窗,長, 得 (*)要線性相位, 就要關(guān)于偶對稱,而關(guān)于偶對稱, 故要求所以要求關(guān)于偶對稱. 圖5 窗函數(shù) 圖6 加窗后的單位脈沖響應(yīng)再回過來檢驗是否滿足精度要求. 圖7 圖4的脈沖響應(yīng)的頻響 圖8 理想頻響與實際頻響的對比若基本滿足, 則依截取的, 制硬件, 編軟件.為便于選擇使用, 將5種常見的窗函數(shù)基本參數(shù)如表1所示。表1 5種常見的窗函數(shù)基本參數(shù)類型窗函數(shù)的旁瓣峰過渡帶寬度加窗后濾波器的阻帶最小衰減rectwin-134p/N-21bartlet三角-258p/N-25hanning-31

19、8p/N-44hamming-418p/N-53blackman-5712p/N-74（2） 頻率取樣法 窗口設(shè)計法事從時域出發(fā)，把理想的用一定形狀的窗口函數(shù)截取成有限長的h(n)，以此h(n)來近似理想的，從而頻率響應(yīng)也近似于理想的頻率響應(yīng)。我們知道一個有限長序列可以通過其頻譜的相同長度的等間隔采樣值準確地恢復(fù)原有的序列。頻率采樣法便是從頻域出發(fā)，對理想的頻率響應(yīng)加以等間隔采樣 (16)然后以此作為實際FIR濾波器的頻率特性的離散樣本H(k)，即 (17)由H(k)通過IDFT可求出有限長序列h(n)為 (18)利用M個頻率的離散樣本H(k)同樣可求出FIR濾波器的系統(tǒng)函數(shù)H(z)及頻率響應(yīng)

20、。 (19)令可得到濾波器的頻率響應(yīng)。如果設(shè)計的是線性相位的FIR數(shù)字濾波器，其采樣值H(k)的相位的幅度一定要滿足特定的約束條件，這個設(shè)計時一定要注意。（3） 最優(yōu)化設(shè)計法 最優(yōu)化設(shè)計法事以最佳一致逼近（最大誤差最小化）理論為基礎(chǔ)，利用雷米茲算法設(shè)計的具有等波紋特性的設(shè)計方法。具體設(shè)計步驟如下： 對設(shè)計指標進行歸一化處理。 確定remezord函數(shù)所需要的參數(shù)。包括歸一化邊界頻率、各頻帶的幅度要求和波紋要求等。歸一化邊界頻率總是從0開始到1結(jié)束，故只需遞增列出中間的邊界頻率；頻帶幅度要求不含過渡區(qū)，個數(shù)是邊界頻率個數(shù)的一半加1；波紋要求是頻帶內(nèi)幅度允許的波動要求，與分貝間的關(guān)系是： (20)

21、 利用remezord函數(shù)確定remez所需參數(shù)。 調(diào)用remez函數(shù)進行設(shè)計。 利用freqz函數(shù)驗算技術(shù)指標是否滿足要求。3 數(shù)字濾波器的MATLAB設(shè)計3.1 直接程序設(shè)計法（1）IIR的直接程序設(shè)計法 例如欲設(shè)計一數(shù)字（IIR）帶阻濾波器，其數(shù)字域指標為：數(shù)字阻帶邊緣頻率分別為0.4和0.7，數(shù)字通帶邊緣頻率為0.25和0.8，通帶波動為1db最小阻帶衰減為40db。 此題的MATLAB程序為：ws=0.4*pi 0.7*pi; %數(shù)字阻帶邊緣頻率wp=0.25*pi 0.8*pi; %數(shù)字通帶邊緣頻率rp=1 ; %通帶波動（db）as=40; %阻帶衰減n,wn=cheb2ord(

22、wp/pi,ws/pi,rp,as);%根據(jù)給定指標，確定濾波器的階數(shù)N和頻率縮放因子Wnb,a=cheby2(n,as,ws/pi,'stop');%返回的b,a分別為H（z)的分子、分母。h,w=freqz(b,a,512);%返回的h,w分別為濾波器的頻率響應(yīng)及其頻率plot(w/pi,abs(h);%畫出頻率響應(yīng)（以w/pi為橫軸）grid;xlabel('w/pi');ylabel('幅值');title('頻率響應(yīng)');程序運行結(jié)果為：圖9 所設(shè)計的帶阻濾波器的頻率響應(yīng)在設(shè)計中如果該濾波器的特性不滿足要求，原有的參數(shù)必

23、須做相應(yīng)的調(diào)整，在程序中只需對參數(shù)做新的設(shè)定就可以得到所需的濾波器。接下來我們來看看此題所設(shè)計的濾波器的濾波效果：S為含有3個頻率成分的信號（歸一化頻率（w/2）分別為0.1、0.3、0.45），用所設(shè)計的濾波器濾除歸一化頻率為0.3的成分。n=0:100;s1=sin(pi*0.2*n);s2=sin(pi*0.6*n);s3=sin(pi*0.9*n);s4=s1+s3;s=s1+s2+s3;sf=filter(b,a,s);subplot(311)stem(n,s);title('濾波前的信號');subplot(312);stem(n,sf);title('濾波

24、后的信號');subplot(313);stem(n,s4);title('想要保留的信號');程序運行的結(jié)果為：圖10 采用filter函數(shù)進行數(shù)字濾波前后信號比較示意圖由圖可以看出，濾波后的信號與想要保留的信號基本一致（相位有些許偏差，但基本一致），所以我們可以說該濾波器基本滿足了以上所提出的濾波要求。（2）FIR的直接程序設(shè)計法例如欲設(shè)計一個線性相位數(shù)字（FIR）帶通濾波器，其數(shù)字域指標為：數(shù)字通帶邊界頻率為0.35和0.65，數(shù)字阻帶邊界頻率為0.2和0.8，通帶波動為1db，最小阻帶衰減為60db。 FIR數(shù)字濾波器的窗函數(shù)法 此題的MATLAB程序為：ws1

25、=0.2*pi;wp1=0.35*pi;wp2=0.65*pi;ws2=0.8*pi;as=60;tr=min(wp1-ws1),(ws2-wp2);M=ceil(11*pi/tr)+1; %濾波器的階數(shù)，程序運行后M=75n=0:1:M-1;r=(M-1)/2;%r為群時延wc1=(ws1+wp1)/2;wc2=(wp2+ws2)/2;hd=sin(wc2*(n-r)+eps)./(pi*(n-r)+eps)-sin(wc1*(n-r)+eps)./(pi*(n-r)+eps);%hd為理想濾波器的脈沖響應(yīng)w_bla=(blackman(M)'%長度為M的blackman窗h=hd.

26、*w_bla;%h為濾波器的實際脈沖響應(yīng)stem(n,h);title('濾波器的實際單位脈沖響應(yīng)');freqz(h,1,512);title('幅度響應(yīng)和相位響應(yīng)');圖11 所設(shè)計的濾波器的實際單位脈沖響應(yīng)由上圖可知濾波器的實際脈沖響應(yīng)h是偶對稱的，即h(n)=h(M-1-n),故該濾波器滿足FIR線性相位的條件，該濾波器是線性相位的FIR濾波器。圖12 所設(shè)計的帶通濾波器的幅度和相位響應(yīng)由濾波器的相位特性也可以看出該濾波器是線性相位的FIR濾波器。接下來我們來看看此題所設(shè)計的濾波器的濾波效果：S為含有3個頻率成分的信號（歸一化頻率（w/2）分別為0.05

27、、0.2、0.45），用所設(shè)計的濾波器濾除歸一化頻率為0.05和0.45的成分。l=0:100;s1=sin(0.1*pi*l);s2=sin(0.4*pi*l);s3=sin(pi*0.9*l);s=s1+s2+s3;sf=filter(h,1,s);subplot(311)stem(l,s);title('濾波前的信號');subplot(312);stem(l,sf);title('濾波后的信號');subplot(313);stem(l,s2);title('想要保留的信號');圖13 采用filter函數(shù)進行數(shù)字濾波前后信號比較示意圖由

28、上圖可知濾波后的信號和想要保留的信號的幅度和頻率基本一致（濾波后的信號相對于想要保的信號有一個相位延遲，這是線性相位FIR濾波器的群延遲引起的，此濾波器留的群延遲r=(M-1)/2=37），所以我們可以說該濾波器基本滿足了以上所提出的濾波要求。 FIR數(shù)字濾波器的頻率采樣法 此題的MATLAB程序為：M=40;%取濾波器的階數(shù)為40al=(M-1)/2;%群時延n=0:M-1;T2=0.59417456;T1=0.109021;Hrs=zeros(1,5),T1,T2,ones(1,7),T2,T1,zeros(1,9),T1,T2,ones(1,7),T2,T1,zeros(1,4);%采樣

29、值的幅值k1=0:floor(M-1)/2);k2=floor(M-1)/2)+1:M-1;angH=-al*(2*pi)/M*k1,al*(2*pi)/M*(M-k2);%采樣值的相位H=Hrs.*exp(j*angH);h=real(ifft(H,M);%長度為M的單位脈沖響應(yīng)stem(n,h);title('濾波器的實際單位脈沖響應(yīng)');freqz(h,1,512);title('幅度響應(yīng)和相位響應(yīng)');圖14 所設(shè)計的濾波器的實際單位脈沖響應(yīng)由圖14可知濾波器的實際脈沖響應(yīng)h是偶對稱的，即h(n)=h(M-1-n),故該濾波器滿足FIR線性相位的條件，該

30、濾波器是線性相位的FIR濾波器。圖15 所設(shè)計的帶通濾波器的幅度和相位響應(yīng)由濾波器的相位特性也可以看出該濾波器是線性相位的FIR濾波器。此濾波器的群延時為al=(M-1)/2=19.5。接下來我們來看看此題所設(shè)計的濾波器的濾波效果：S為含有3個頻率成分的信號（歸一化頻率（w/2）分別為0.05、0.2、0.45），用所設(shè)計的濾波器濾除歸一化頻率為0.05和0.45的成分。l=0:100;s1=sin(0.1*pi*l);s2=sin(0.4*pi*l);s3=sin(pi*0.9*l);s=s1+s2+s3;sf=filter(h,1,s);subplot(311);stem(l,s);tit

31、le('濾波前的信號');subplot(312);stem(l,sf);title('濾波后的信號');subplot(313);stem(l,s2);title('想要保留的信號');圖16 采用filter函數(shù)進行數(shù)字濾波前后信號比較示意圖同上面分析相似，濾波后的信號和想要保留的信號的幅度和頻率基本一致（濾波后的信號相對于想要保的信號有一個相位延遲，這是線性相位FIR濾波器的群延遲引起的，此濾波器留的群延遲（r=(M-1)/2=19.5），所以我們可以說該濾波器基本滿足了以上所提出的濾波要求。 FIR數(shù)字濾波器的最優(yōu)設(shè)計法此題的MATLAB

32、程序為：%設(shè)計指標ws1=0.2*pi;wp1=0.35*pi;wp2=0.65*pi;ws2=0.8*pi;rp=1;as=60;%設(shè)置邊界頻率和幅度要求F=ws1/pi,wp1/pi,wp1/pi,ws2/pi;A=0,1,0;%設(shè)置各頻帶的波紋要求devp=(10(rp/20)-1)/(10(rp/20)+1);devs=10(-as/20);dev=devs,devp,devs;%確定remez參數(shù)，其中濾波器的階數(shù)為（N+1），程序運行后得到N=26N,Fo,Ao,W=remezord(F,A,dev);%調(diào)用remez函數(shù)進行設(shè)計h=remez(N,Fo,Ao,W);n=0:N;s

33、tem(n,h);title('濾波器的單位沖激響應(yīng)');freqz(h,1,512);title('幅度響應(yīng)和相位響應(yīng)');圖17 所設(shè)計的濾波器的實際單位脈沖響應(yīng)由圖17可知濾波器的實際脈沖響應(yīng)h是偶對稱的，即h(n)=h(N-n),故該濾波器滿足FIR線性相位的條件，該濾波器是線性相位的FIR濾波器。圖18 所設(shè)計的帶通濾波器的幅度和相位響應(yīng)由濾波器的相位特性也可以看出該濾波器是線性相位的FIR濾波器。此濾波器的群延時為al=(N)/2=13。接下來我們來看看此題所設(shè)計的濾波器的濾波效果：S為含有3個頻率成分的信號（歸一化頻率（w/2）分別為0.05、0.

34、2、0.45），用所設(shè)計的濾波器濾除歸一化頻率為0.05和0.45的成分。l=0:100;s1=sin(0.1*pi*l);s2=sin(0.4*pi*l);s3=sin(pi*0.9*l);s=s1+s2+s3;sf=filter(h,1,s);subplot(311)stem(l,s);title('濾波前的信號');subplot(312);stem(l,sf);title('濾波后的信號');subplot(313);stem(l,s2);title('想要保留的信號');圖19采用filter函數(shù)進行數(shù)字濾波前后信號比較示意圖同上面分析

35、相似，濾波后的信號和想要保留的信號的幅度和頻率基本一致（濾波后的信號相對于想要保的信號有一個相位延遲，這是線性相位FIR濾波器的群延遲引起的，此濾波器留的群延遲r=(N)/2=13），所以我們可以說該濾波器基本滿足了以上所提出的濾波要求。3.2 利用信號處理工具箱SPTOOL設(shè)計法SPTOOL是信號處理工具箱中一個具有交互式圖形用戶界面的信號處理工具，專門用于完成常用的數(shù)字信號處理任務(wù)。通過這個工具，只需要鼠標簡單的操縱鼠標（點擊或拖動），就可以完成載入、觀察、分析、實現(xiàn)和設(shè)計數(shù)字濾波器并進行譜分析等甚至十分復(fù)雜的數(shù)字信號處理任務(wù)，而不需要用戶對數(shù)字濾波器的設(shè)計原理非常熟悉。下面是用SPTOO

36、L工具設(shè)計的IIR濾波器對上述信號進行濾波。SPTOOL使用步驟：導(dǎo)入信號、濾波器、頻譜，設(shè)計濾波器，對信號濾波，分析輸入輸出信號的譜。1 IIR濾波器的設(shè)計創(chuàng)建并導(dǎo)入信號源， 在MATLAB命令窗口輸入命令：fs=4000;（這里取采樣頻率為4000Hz）;n=0:100;s1=sin(pi*0.2*n);s2=sin(pi*0.6*n);s3=sin(pi*0.9*n);s4=s1+s3;s=s1+s2+s3;此時，變量fs、n、s、s1、s2、s3、s4將顯示在MATLAB的workspace列表中。在MTLAB命令窗口中輸入sptool，即可彈出SPTOOL的主界面,如圖所示： 圖20

37、 SPTOOL的主界面這里我們首先來設(shè)計前面的IIR帶阻濾波器，前面的數(shù)字指標轉(zhuǎn)化為模擬指標為：通帶頻率分別為500Hz和1600Hz，阻帶頻率為800Hz和1400Hz，通帶波動為1db，最小阻帶衰減為40db。以上三個信號轉(zhuǎn)換為模擬域的頻率分別為400Hz，1200Hz，1800Hz。1）導(dǎo)入信號。使用【file/import】可以導(dǎo)入信號、濾波器和譜。信號的來源可以是MATLAB工作空間變量或MAT數(shù)據(jù)文件，也可以是在DATA前面的文本框中直接輸入信號的數(shù)據(jù)。這里信號源仍然利用前面程序設(shè)計法中的混頻信號，將信號s和采樣頻率fs=4kHz導(dǎo)入并命名為sig1，將信號s2導(dǎo)入命名為sig2,

38、將信號s4導(dǎo)入命名為sig4。2）濾波器的設(shè)計。單擊Filter欄中的New按鈕打開濾波器設(shè)計工具。在界面的最上面制定濾波器的名字、采樣頻率和設(shè)計方法。這里命名濾波器名稱為filt1，采樣頻率fs=4kHz和設(shè)計方法為chebyshev typeII IIR；左側(cè)的Specification面板用來指定設(shè)計指標；minimum order,type為bandstop，通帶頻率分別為500Hz和1600Hz，阻帶頻率為800Hz和1400Hz，通帶波動為1db，最小阻帶衰減為40db。完成以上輸入之后單擊Apply按鈕，則完成了濾波器的設(shè)計；中間的Frequency Response面板則根據(jù)設(shè)

39、計指標自動繪出濾波器的幅頻響應(yīng)，而右側(cè)的Measurements面板顯示了濾波器設(shè)計完成后的實測參數(shù)。此時在Filters欄中，選中filt1按下View按鈕即彈出界面，我們還可以觀察到濾波器的相位響應(yīng)、幅頻相位響應(yīng)、脈沖響應(yīng)、階躍響應(yīng)等等（只須點擊界面上的相應(yīng)快捷按鈕）。圖21 所設(shè)計的濾波器的幅值響應(yīng)3）將濾波器應(yīng)用到sig1信號序列。分別在Signals、Filters欄中選擇sig2、filt1單擊Filters欄列表下的Apply按鈕，在彈出的Apply Filter對話框中將輸出信號命名為Sig3（濾波后信號）。同時按下Alt與Shift鍵選中Sig3與Sig2點擊Signals欄

40、下的View可觀察到它們的時域波形(其中紅色是信號sig3，藍色是信號sig2)，如圖所示。由圖可觀察出信號s2通過該帶阻濾波器后得到很大的衰減，幾乎衰減為0。圖22 信號s2濾波前后的比較分別在Signals、Filters欄中選擇sig1、filt1單擊Filters欄列表下的Apply按鈕，在彈出的Apply Filter對話框中將輸出信號命名為Sig5（濾波后信號）。同時按下Alt與Shift鍵選中Sig4與Sig5點擊Signals欄下的View可觀察到它們的時域波形(其中紅色是信號sig5，藍色是信號sig4)，如圖所示。由圖可觀察出信號s通過濾波器后所得的信號與s4基本吻合，只是

41、相位上有點偏差，即該帶阻濾波器基本濾除了信號s2，濾波效果良好。圖23 濾波后的信號與s4的比較4）進行頻譜分析。在Signals中選擇sig1，單擊Spectra欄下的Create按鈕，在彈出的Spectra Viewer界面中選擇method為FFT，nfft=512，單擊Apply按鈕生成s的頻譜spect1同理得到sig4的頻譜spect4，sig5的頻譜spect5,同時按下Alt與Shift鍵選中spect1與spect5，點擊Spectra欄下的View可觀察到它們的頻譜,下圖是sig1和sig5的頻譜，其中紅色代表濾波后信號的頻譜（sig5），藍色代表濾波前信號的頻譜（sig1

42、）。圖24 濾波前后信號的頻譜的比較由頻譜可以看出，混頻信號s中頻率為1200Hz的部分基本被濾除了，頻率為400Hz和1800Hz的部分基本被保留下來，未受濾波的影響。同時按下Alt與Shift鍵選中spect4與spect5，點擊Spectra欄下的View可觀察到它們的頻譜,下圖是sig4和sig5的頻譜，其中紅色代表濾波后信號的頻譜（sig5），藍色代表濾波前信號的頻譜（sig4），如下圖。圖25 濾波后的信號與信號s4頻譜的比較由上圖可以看出，s4與濾波后的信號頻譜基本一致，所以我們可以得出結(jié)論：該帶阻濾波器的濾波效果較好。2 FIR濾波器的設(shè)計 與1的設(shè)計過程類似，這里就簡述過程。

43、 創(chuàng)建并導(dǎo)入信號源， 在MATLAB命令窗口輸入命令：fs=4000Hz; l=0:100; s1=sin(0.1*pi*l); s2=sin(0.4*pi*l); s3=sin(pi*0.9*l); s=s1+s2+s3; s4=s1+s3;前面的數(shù)字指標轉(zhuǎn)化為模擬指標為：通帶頻率分別為700Hz和1300Hz，阻帶頻率為400Hz和1600Hz，通帶波動為1db，最小阻帶衰減為60db，以上三個信號轉(zhuǎn)換為模擬域的頻率分別為200Hz，800Hz，1800Hz。1）導(dǎo)入信號。將信號s和采樣頻率fs=4kHz導(dǎo)入并命名為sig1，將信號s2導(dǎo)入命名為sig2,將信號s4導(dǎo)入命名為sig4。2）

44、 濾波器的設(shè)計。圖26 所設(shè)計的濾波器的幅值響應(yīng)3） 將濾波器應(yīng)用到sig1信號、sig4信號，生成的信號分別命名為sig5和sig3。圖27 信號s4濾波前后的比較上圖中sig4信號與sig3信號比較（其中紅色代表sig3，藍色代表sig4）可知信號s4通過該帶通濾波器后得到大大衰減，幾乎衰減為0。圖28 濾波后的信號與信號s2的比較上圖中sig2和sig5比較（其中紅色代表sig5，藍色代表sig2），可知信號s通過該帶通濾波器后的信號與s2基本具有相同的頻率和幅值（只是相位上有一點誤差），所以我們認為該濾波器較好的濾掉了s4，濾波性能優(yōu)良。4）進行頻譜分析。分別生成sig1的頻譜spec

45、t1，sig2的頻譜spect2，sig5的頻譜spect5。下圖是sig1和sig5的頻譜，其中紅色代表濾波后信號的頻譜（sig5），藍色代表濾波前信號的頻譜（sig1）。由圖可看出該濾波器基本濾除了頻率為200Hz和1800Hz的信號，保留了頻率為800Hz的信號。圖29 濾波前后信號的頻譜的比較下圖是sig2和sig5的頻譜，其中紅色代表濾波后信號的頻譜（sig5），藍色代表濾波前信號的頻譜（sig2）。由圖可看出濾波后的信號與s2具有基本一致的頻譜。圖30 濾波后的信號與信號s2頻譜的比較3.3 濾波器分析設(shè)計工具箱FDATOOL法FDATool（Filter Design &

46、 Analysis Tool）是MATLAB信號處理工具箱里專用的濾波器設(shè)計分析工具，MATLAB6.0以上的版本還專門增加了濾波器設(shè)計工具箱（Filter Design Toolbox）。FDATool可以設(shè)計幾乎所有的基本的常規(guī)濾波器，包括FIR和IIR的各種設(shè)計方法。它操作簡單，方便靈活。FDATool界面總共分兩大部分，一部分是Design Filter，在界面的下半部，用來設(shè)置濾波器的設(shè)計參數(shù)，另一部分則是特性區(qū)，在界面的上半部分，用來顯示濾波器的各種特性。Design Filter部分主要分為：Filter Type（濾波器類型）選項，包括Lowpass（低通）、Highpass（

47、高通）、Bandpass（帶通）、Bandstop（帶阻）和特殊的FIR濾波器。Design Method（設(shè)計方法）選項，包括IIR濾波器的Butterworth（巴特沃思）法、Chebyshev Type I（切比雪夫I型）法、 Chebyshev Type II（切比雪夫II型） 法、Elliptic（橢圓濾波器）法和FIR濾波器的Equiripple法、Least-Squares（最小乘方）法、Window（窗函數(shù)）法。Filter Order（濾波器階數(shù)）選項，定義濾波器的階數(shù)，包括Specify Order（指定階數(shù)）和Minimum Order（最小階數(shù)）。在Specify Or

48、der中填入所要設(shè)計的濾波器的階數(shù)（N階濾波器，Specify OrderN-1），如果選擇Minimum Order則MATLAB根據(jù)所選擇的濾波器類型自動使用最小階數(shù)。Frenquency Specifications選項，可以詳細定義頻帶的各參數(shù)，包括采樣頻率Fs和頻帶的截止頻率。它的具體選項由Filter Type選項和Design Method選項決定，例如Bandpass（帶通）濾波器需要定義Fstop1（下阻帶截止頻率）、Fpass1（通帶下限截止頻率）、Fpass2（通帶上限截止頻率）、Fstop2（上阻帶截止頻率），而Lowpass（低通）濾波器只需要定義Fstop1、Fpa

49、ss1。采用窗函數(shù)設(shè)計濾波器時，由于過渡帶是由窗函數(shù)的類型和階數(shù)所決定的，所以只需要定義通帶截止頻率，而不必定義阻帶參數(shù)。Magnitude Specifications選項，可以定義幅值衰減的情況。例如設(shè)計帶通濾波器時，可以定義Wstop1（頻率Fstop1處的幅值衰減）、Wpass（通帶范圍內(nèi)的幅值衰減）、Wstop2（頻率Fstop2處的幅值衰減）。當采用窗函數(shù)設(shè)計時，通帶截止頻率處的幅值衰減固定為6db，所以不必定義。Window Specifications選項，當選取采用窗函數(shù)設(shè)計時，該選項可定義，它包含了各種窗函數(shù)?，F(xiàn)用FDATOOL設(shè)計一個滿足上述要求的IIR濾波器：在命令窗口

50、里輸入fdatool，即可彈出fdatool的主界面，如圖31 fdatool的主界首先在Filter Type中選擇bandstop，接著在Design Method選項中選擇chebyshev typeII IIR,再在指定Filter Order項中的minimum order，然后令units為normalized(0to1),數(shù)字阻帶邊緣頻率分別為0.4和0.7，數(shù)字通帶邊緣頻率為0.25和0.8，通帶波動為1db最小阻帶衰減為40db。設(shè)置完以后點擊Design Filter即可得到所設(shè)計的FIR濾波器。通過菜單選項Analysis可以在特性區(qū)看到所設(shè)計濾波器的幅頻響應(yīng)、相頻響應(yīng)、

51、零極點配置和濾波器系數(shù)等各種特性。圖32 所設(shè)計的帶阻濾波器的幅值響應(yīng) 在設(shè)計過程中，可以對比濾波器幅頻相頻特性和設(shè)計要求，隨時調(diào)整參數(shù)和濾波器類型，以便得到最佳效果。其它類型的FIR濾波器和IIR濾波器也都可以使用FDATool來設(shè)計。 同樣地，我們也可以設(shè)計滿足前面指標的FIR濾波器，其幅值響應(yīng)為：圖33所設(shè)計的帶通濾波器的幅值響應(yīng)在設(shè)計過程中，可以對比濾波器幅頻相頻特性和設(shè)計要求，隨時調(diào)整參數(shù)和濾波器類型，以便得到最佳效果。其它類型的FIR濾波器和IIR濾波器也都可以使用FDATool來設(shè)計。4 Simulink仿真 本文通過調(diào)用Simuink中的功能模塊構(gòu)成數(shù)字濾波器的仿真框圖，在仿真

52、過程中，可以雙擊各功能模塊，隨時改變參數(shù)，獲得不同狀態(tài)下的仿真結(jié)果。在Simulink環(huán)境下,通過設(shè)置濾波器功能模塊的參數(shù)（digital filter design模塊的參數(shù)）來直接設(shè)計濾波器。 我們先設(shè)計滿足上述帶阻IIR濾波器參數(shù)的濾波器，再進行仿真。找到Simuink功能模塊，communications blockset/comm filters/filter designs link/digital filter design,找到digital filter design模塊。分別在sources和sinks中找到正玄信號示波器，再在math operation中找到相加器。信號n=0:100;s1=sin(pi*0.2*n);s2=sin(pi*0.6*n);s3=sin(pi*0.9*n);s4=s1+s3;s=s1+s2+s3;修改正玄信號的參數(shù)sine type/time ba