《數(shù)字信號處理-基于數(shù)值計算》課件-第5章

第5章有限長單位脈沖響應(yīng)數(shù)字濾波器的設(shè)計5.1線性相位FIR數(shù)字濾波器的特性5.2窗函數(shù)設(shè)計法5.3頻率采樣法5.4FIR數(shù)字濾波器的最優(yōu)化設(shè)計5.5IIR與FIR數(shù)字濾波器的比較5.6MATLAB中的濾波器設(shè)計工具FDATool

5.1線性相位FIR數(shù)字濾波器的特性

5.1.1線性相位的條件

線性相位意味著一個系統(tǒng)的相頻特性是頻率的線性函數(shù)，即

j(ω)=－αω

(5.1.1)

式中：α為常數(shù)，此時通過這一系統(tǒng)的各頻率分量的時延為一相同的常數(shù)，系統(tǒng)的群時延為(5.1.2)對于單位脈沖響應(yīng)h(n)為實數(shù)且長度為N的FIR濾波器，假設(shè)它具有線性相位,即j(ω)=－αω,則其DTFT可以寫為(5.1.3)式中：H(ω)是正或負的純實函數(shù)，稱為幅度函數(shù)；j(ω)=－αω是相位函數(shù)。由于復(fù)數(shù)等式相等，等式的第2部分和等式右邊的實部與虛部應(yīng)當(dāng)各自相等，因此實部與虛部的比值應(yīng)當(dāng)相等，即可以導(dǎo)出

(5.1.4)將式(5.1.4)兩邊交叉相乘，再將等式右邊各項移到左邊，應(yīng)用三角函數(shù)的恒等關(guān)系可得：滿足上式的條件是h(n)為實數(shù)，且滿足：(5.1.5)即：h(n)為實數(shù)且偶對稱，稱為第一類線性相位條件。

另外一種情況是，除了上述的線性相位外，還有一附加的相位，即

j(ω)=β－αω

(5.1.6)

利用類似的推導(dǎo)方法，可以導(dǎo)出新的線性相位條件為

(5.1.7)即：h(n)為實數(shù)且奇對稱，稱為第二類線性相位條件。兩類的線性相位差別只在于第二類線性相位具有一個附加相移，如圖5.1.1所示。圖5.1.1線性相位特性5.1.2線性相位FIR濾波器的幅度特性

當(dāng)FIR濾波器滿足線性相位時，由于h(n)為實數(shù)且具有偶對稱或者奇對稱，而h(n)的點數(shù)又可以分為奇數(shù)和偶數(shù)兩種情況，組合起來有4種情況，分別稱為線性相位1~4型。下面分別討論這4種情況下的幅度函數(shù)特性。

1.h(n)偶對稱，N為奇數(shù)——1型

線性相位1型濾波器的h(n)示意圖如圖5.1.2所示，條件是：h(n)=h(N－1－n)，代入式(5.1.3)，得圖5.1.2線性相位1型濾波器h(n)示意圖令，則令(5.1.10)則(5.1.11)由于cosnω對于頻率ω=0、π、2π呈偶對稱，因此H(ω)對于頻率ω=0、π、2π也呈偶對稱。

【例5.1.1】

N=5，h(0)=h(1)=h(3)=h(4)=－1/2，h(2)=2，求幅度函數(shù)H(ω)。

解

N為奇數(shù)且h(n)滿足偶對稱關(guān)系，則有a(0)=h(2)=2；a(1)=2h(3)=－1；a(2)=2h(4)=－1H(ω)=2－cosω－cos2ω=2－(cosω+cos2ω)其幅度函數(shù)曲線如圖5.1.3所示，由圖可見該曲線對于ω=π呈偶對稱。圖5.1.3幅度函數(shù)曲線

2.h(n)偶對稱，N為偶數(shù)——2型

線性相位2型濾波器的h(n)示意圖如圖5.1.4所示，條件是：h(n)=h(N－1－n)，代入式(5.1.3)，得圖5.1.4線性相位2型濾波器h(n)示意圖

(5.1.12)令m=N/2－n，代入(5.1.12)式，得(5.1.13)改寫為由于cos(ω(n－1/2))對ω=π呈奇對稱，所以H(ω)對ω=π也呈奇對稱；而且當(dāng)ω=π時，cos(π(n－1/2))=0，因此H(π)=0，即H(z)在z=－1處必有一個零點。因此，2型濾波器不能用來設(shè)計高通、帶阻濾波器。

3.h(n)奇對稱，N為奇數(shù)——3型

線性相位3型濾波器的h(n)示意圖如圖5.1.5所示，條件是：h(n)=－h(huán)(N－1－n)、，代入式(5.1.3)，得則有(5.1.14)圖5.1.5線性相位3型濾波器h(n)示意圖

(5.1.15)令

，得改寫為(5.1.16)由于sin(nω)在ω=0、π、2π處都為0，并對這些點呈奇對稱，則H(ω)在ω=0、π、2π處都為0，并對這些點呈奇對稱，H(z)在z=±1處必有零點。因此，其不能用作低通、高通和帶阻濾波器的設(shè)計。

4.h(n)奇對稱，N為偶數(shù)——4型

線性相位4型濾波器的h(n)示意圖如圖5.1.6所示，條件是：h(n)=－h(huán)(N－1－n)，代入式(5.1.3)，得則有圖5.1.6線性相位4型濾波器h(n)示意圖令，則有(5.1.17)由于在ω=0、2π處為零，所以H(ω)在ω=0、2π處為零，即H(z)在z=1上有零點，并對ω=0、2π呈奇對稱。因此，其不能用作低通、帶阻濾波器的設(shè)計。對這四種線性相位濾波器的幅度特性和相位特性歸納見表5.1.1。5.1.3線性相位FIR濾波器的零點特性

由于線性相位FIR數(shù)字濾波器的單位脈沖響應(yīng)具有對稱特性，即h(n)=±h(N－1－n)，則令N－1－n=m，代入上式，得所以有

H(z)=±z－(N－1)H(z－1)

(5.1.18)若z=zi是H(z)的零點，則z=zi－1也一定是H(z)的零點。由于h(n)是實數(shù)，H(z)的零點還必須共軛成對出現(xiàn)，所以z=z*i及z=1/z*也必是零點。

所以線性相位濾波器的零點必須是互為倒數(shù)的共軛對，這種共軛對共有四種可能的情況：

(1)零點既不在單位圓上，也不在實軸上，則有四個互為倒數(shù)共軛對：zi

、z*i、1/zi、1/z*i，如圖5.1.7(a)所示。

(2)零點在單位圓上，但不在實軸上，因倒數(shù)就是自己的共軛，所以有一對共軛零點：zi、z*i

，如圖5.1.7(b)所示。

(3)零點不在單位圓上，但在實軸上，是實數(shù)，共軛就是自己，所以有一對互為倒數(shù)的零點：zi、1/zi

，如圖5.1.7(c)所示。

(4)零點既在單位圓上，又在實軸上，共軛和倒數(shù)都合為一點，所以成單出現(xiàn)，只有兩種可能：zi=1或zi=－1，如圖5.1.7(d)所示。

在設(shè)計濾波器時，應(yīng)根據(jù)實際需要選擇合適類型，并在設(shè)計時遵循其約束條件。圖5.1.7線性相位FIR濾波器的零點特性

5.2窗函數(shù)設(shè)計法

FIR數(shù)字濾波器的設(shè)計同樣也采用逼近法。假設(shè)希望設(shè)計得到的濾波器(稱為目標(biāo)濾波器)的頻率響應(yīng)為Hd(ejω)，那么FIR濾波器的設(shè)計就是尋找一個傳遞函數(shù)去逼近它，對于長度為N的FIR數(shù)字濾波器，其單位脈沖響應(yīng)為h(n)，則其傳遞函數(shù)為

(5.2.1)常用的逼近方法有三種：窗函數(shù)設(shè)計法(時域逼近)、頻率采樣法(頻域逼近)、最優(yōu)化設(shè)計(等波紋逼近)。本節(jié)介紹窗函數(shù)設(shè)計法。窗函數(shù)設(shè)計法是從單位脈沖響應(yīng)序列著手，使設(shè)計的濾波器h(n)逼近目標(biāo)濾波器的單位脈沖響應(yīng)hd(n)。目標(biāo)濾波器通常都采用理想濾波器模型，可以根據(jù)所要設(shè)計的濾波器的類型和邊界頻率來確定，因此hd(n)可以通過傅氏反變換獲得，即(5.2.2)一般來說，理想頻響是分段恒定，在邊界頻率點處有突變點，所以，這樣得到的理想單位脈沖響應(yīng)hd(n)往往都是無限長序列，而且是非因果的。因此，濾波器的設(shè)計問題變?yōu)樵鯓佑靡粋€有限長的序列去近似無限長的hd(n)。最簡單的辦法是直接截取一段hd(n)代替h(n)，即

h(n)=w(n)hd(n)

(5.2.3)5.2.1矩形窗函數(shù)設(shè)計法

假設(shè)要設(shè)計一個截止頻率為ωc的具有線性相位的FIR數(shù)字低通濾波器，則該FIR數(shù)字濾波器的矩形窗函數(shù)法設(shè)計步驟如下：

(1)計算hd(n)。

以理想低通數(shù)字濾波器模型作為所要設(shè)計的目標(biāo)濾波器，則濾波器的傳遞函數(shù)為(5.2.4)式中：α為低通濾波器的群時延，則其單位脈沖響應(yīng)為(5.2.5)理想低通的單位脈沖響應(yīng)和幅度函數(shù)特性如圖5.2.1所示，圖(a)為理想低通的單位脈沖響應(yīng)，圖(b)為它的幅度函數(shù)。由圖(a)可以看出理想低通的單位脈沖響應(yīng)以n=α呈偶對稱。

(2)計算h(n)。由式(5.2.5)可見，它是一個中心點在α的偶對稱的無限長非因果序列。如果截取一段n=0～N－1的hd(n)作為h(n)，則為保證所得到的是線性相位FIR濾波器，應(yīng)必須滿足：圖5.2.1理想低通的單位脈沖響應(yīng)和幅度函數(shù)特性因此，可得式中：wR(n)=RN(n)。

(3)計算分析H(ejω)。對h(n)作傅立葉變換可得

H(ejω)=Hd(ejω)*WR(ejω)

(5.2.7)式中：WR(ejω)為矩形窗口函數(shù)的頻譜，因此有用幅度函數(shù)和相位函數(shù)來表示，則有

WR(ejω)=WR(ω)e－jωα

線性相位部分為e－jωα，α=(N－1)/2表示延時一半長度。幅度函數(shù)為

(5.2.8)矩形窗函數(shù)的時域波形及其幅度函數(shù)波形如圖5.2.2所示：圖(a)為矩形窗函數(shù)的時域波形；圖(b)為其幅度函數(shù)波形。由矩形窗函數(shù)的時域波形可以看出，長度為N的矩形窗是以n=(N－1)/2呈偶對稱的。圖5.2.2矩形窗函數(shù)及其幅度函數(shù)理想低通的頻率響應(yīng)也可以寫成幅度函數(shù)和相位函數(shù)的表示形式，即

Hd(ejω)=Hd(ω)ejωα

(5.2.9)

式(5.2.9)中幅度函數(shù)為則所以幅度函數(shù)為(5.2.10)正好是理想濾波器幅度函數(shù)與窗函數(shù)幅度函數(shù)的卷積。其卷積過程如圖5.2.3所示。圖5.2.3矩形窗的卷積過程下面討論4個特殊頻率點的卷積結(jié)果：

(1)當(dāng)ω=0時，幅度響應(yīng)H(0)等于Hd(θ)和WR(－θ)兩個函數(shù)相乘的積分，即圖5.2.3(a)中陰影部分的面積。

(2)當(dāng)ω=ωc時，Hd(θ)正好與WR(ω－θ)的一半重疊，如圖5.2.3(b)所示，因此H(ωc)/H(0)=0.5。

(3)當(dāng)ω=ωc－2π/N時，WR(ω－θ)的主瓣全部在Hd(θ)的通帶內(nèi)，而右邊具有負面積的第一旁瓣全部移出通帶，如圖5.2.3(c)所示，因此卷積結(jié)果出現(xiàn)最大值，頻響出現(xiàn)正肩峰。

(4)當(dāng)ω=ωc+2π/N

時，WR(ω－θ)的主瓣全部在Hd(θ)的通帶之外，Hd(θ)的通帶內(nèi)第一旁瓣(負數(shù))起著主導(dǎo)作用，如圖5.2.3(d)所示，卷積值為最負的值，出現(xiàn)負的肩峰。

整個卷積的結(jié)果如圖5.2.3(e)所示。由圖可以看出加窗后，窗函數(shù)對理想低通濾波器的特性具有以下影響：

(1)改變了理想頻響的邊沿特性，形成過渡帶寬近似等于WR(ω)的主瓣寬度4π/N。

(2)過渡帶兩旁產(chǎn)生肩峰和余振(帶內(nèi)、帶外起伏)，取決于WR(ω)的旁瓣，旁瓣多，余振多；旁瓣相對值大，肩峰強，與N無關(guān)(決定于窗口形狀)。

(3)N增加，過渡帶寬減小，肩峰值不變。因為(5.2.11)式中：x=Nω/2，所以N的改變不能改變主瓣與旁瓣的比例關(guān)系，只能改變WR(ω)的絕對值大小和起伏的密度。當(dāng)N增加時，幅值變大，頻率軸變密，而最大肩峰永遠為8.95%，這種現(xiàn)象稱為吉布斯(Gibbs)效應(yīng)。5.2.2幾種常用的窗函數(shù)

1.矩形窗

上面已講過，不再贅述。

在MATLAB中，函數(shù)Boxcar(N)用于生成長度為N的矩形窗。

2.漢寧窗

漢寧窗也叫升余弦窗，其函數(shù)表達式為(5.2.12)利用傅氏變換的移位特性，漢寧窗頻譜的幅度函數(shù)W(ω)可用矩形窗的幅度函數(shù)表示為幅度譜為其實質(zhì)是三部分矩形窗頻譜相加，使旁瓣互相抵消，能量集中在主瓣，旁瓣大大減小，但是其主瓣寬度增加了1倍，為8π/N。漢寧窗函數(shù)頻譜如圖5.2.4所示。圖5.2.4漢寧窗函數(shù)頻譜在MATLAB中，函數(shù)Hanning(N)用于生成長度為N的漢寧窗，例如：Window=hanning(5)，運行得到：

Window=[0.25000.75001.00000.75000.2500]

3.漢明窗(改進的升余弦窗)

漢明窗函數(shù)的表達式為(5.2.13)它是對漢寧窗的改進，在主瓣寬度(對應(yīng)第一零點的寬度)相同的情況下，旁瓣進一步減小，可使99.96%的能量集中在窗譜的主瓣內(nèi)。在MATLAB中，函數(shù)Hamming(N)用于生成長度為N的漢明窗，例如：Window=hamming(5)，運行得到：

Window=[0.08000.54001.00000.54000.0800]

4.布萊克曼窗(二階升余弦窗)

布萊克曼窗函數(shù)的表達式為(5.2.14)由式可見，它增加了一個二次諧波余弦分量，降低了旁瓣，但主瓣寬度進一步增加，為12π/N。頻譜的幅度函數(shù)為在MATLAB中，函數(shù)Blackman(N)用于生成長度為N的布萊克曼窗，例如：Window=Blackman(5)，運行得到：

Window=[0.00000.34001.00000.34000.0000]四種窗函數(shù)的時域波形如圖5.2.5所示，頻譜如圖5.2.6所示。圖5.2.5四種常用窗函數(shù)的時域波形圖5.2.6四種常用窗函數(shù)的頻譜在相同的技術(shù)指標(biāo)要求下，采用不同窗函數(shù)設(shè)計，當(dāng)N=51、

ωc=0.5π時數(shù)字低通濾波器的頻譜如圖5.2.7所示，它們的特性參數(shù)如表5.2.1所示。圖5.2.7相同指標(biāo)下四種窗函數(shù)設(shè)計的濾波器特性比較

5.凱塞窗

以上四種窗函數(shù)都是以增加主瓣寬度為代價來降低旁瓣，而且它們的主瓣寬度和旁瓣衰減都固定，難以實現(xiàn)精確設(shè)計。凱塞窗則可自由選擇主瓣寬度和旁瓣衰減，其函數(shù)表達式為(5.2.15)凱塞窗函數(shù)波形如圖5.2.8所示，通常取4<β<9，當(dāng)β=5.44時，接近漢明窗函數(shù)；當(dāng)β=8.5時，接近布萊克曼窗函數(shù)；當(dāng)β=0時，為矩形窗函數(shù)。表5.2.2給出了不同參數(shù)β時,用凱塞窗函數(shù)設(shè)計的濾波器的特性參數(shù),以便設(shè)計時查用。圖5.2.8凱塞窗函數(shù)波形若給定濾波器的過渡帶寬Δω和阻帶最小衰減At，則參數(shù)β、N可由以下經(jīng)驗公式計算：

At=－20lgδ2

(5.2.16)(5.2.17)(5.2.18)在MATLAB中，函數(shù)Kaiser(N，β)用于生成長度為N的布萊克曼窗，例如：Window=Kaiser(5，4.538)，運行得到：

Window=[0.05530.58861.00000.58860.0553]5.2.3窗函數(shù)法設(shè)計舉例

根據(jù)上述討論，如果給定FIR數(shù)字濾波器的技術(shù)指標(biāo)，則采用窗函數(shù)法設(shè)計具有線性相位的FIR數(shù)字濾波器時，可按以下設(shè)計步驟來完成:

(1)確定所要設(shè)計的濾波器的類型和技術(shù)指標(biāo)。

(2)確定目標(biāo)濾波器的頻率特性和hd(n)。以理想濾波器為模型并考慮線性相位確定目標(biāo)濾波器的頻率特性函數(shù)Hd(ejω)，通過求反傅立葉變換得到目標(biāo)濾波器的單位脈沖響應(yīng)hd(n)。當(dāng)目標(biāo)濾波器模型較為復(fù)雜時，直接用反傅立葉變換法計算hd(n)不容易求得。這時可用IDFT代替連續(xù)傅立葉反變換來求得，具體方法是：選擇一個足夠大的整數(shù)M，并對Hd(ejω)進行M點等間隔采樣，得

(5.2.19)再作M點的IDTFT，得(5.2.20)當(dāng)M>>N時，有hM(n)≈hd(n)。

(3)確定窗函數(shù)。查表5.2.1和表5.2.2，根據(jù)阻帶衰減指標(biāo)選擇窗函數(shù)，選擇原則是：在符合指標(biāo)要求的情況下，選擇最簡單的窗函數(shù)。

(4)確定濾波器的最小階數(shù)。查表5.2.1和表5.2.2，根據(jù)技術(shù)指標(biāo)要求的過渡帶和所選用的窗函數(shù)的過渡帶寬度，計算出窗的寬度N，也是濾波器的最小階數(shù)。

(5)確定濾波器階數(shù)。根據(jù)濾波器的類型和線性相位的約束條件，選擇線性相位類型，并以此確定濾波器的階數(shù)和群延時。

(6)計算h(n)=hd(n)wR(n)，即為所要設(shè)計濾波器的單位脈沖響應(yīng)。如有要求，還可進一步求出濾波器的系統(tǒng)函數(shù)。

【例5.2.1】

設(shè)計一線性相位FIR數(shù)字低通濾波器，通帶截止頻率為0.3π，阻帶邊界頻率為0.5π，通帶波動小于2dB，阻帶衰減不小于50dB。

解法一：通過計算設(shè)計。

(1)所需設(shè)計的是低通濾波器，技術(shù)指標(biāo)是：ωp=0.3π,ωr=0.5π，δ=2dB，At=50dB。

(2)以理想低通濾波器為目標(biāo)濾波器模型，確定Hd(ejω),經(jīng)IDTFT得

(3)因At=50dB，根據(jù)表5.2.1選擇漢明窗函數(shù)，其過渡帶為6.6π/N。因Δω=0.2π，故N=33。

(4)由于是低通，選擇情況1進行設(shè)計，N要求為奇數(shù)，則濾波器的階數(shù)確定為32，即N=33，α=(N－1)/2=16。(5)解法二：用MATLAB輔助設(shè)計。

wp=0.3*pi；wr=0.5*pi；

deltaw=wr-wp；

N0=ceil(6.6*pi/deltaw)；%選擇漢明窗，求濾波器長度N0

N=N0+mod(N0+1，2)；%選擇1型濾波器，確保N為奇數(shù)

alfa=(N-1)/2；

wc=(wp+wr)/2；

windows=hamming(N)；%漢明窗

n=0：N-1；

hd=sin(wc*(n-alfa+eps))./(pi*(n-alfa+eps))；%eps微小值，避免0作除數(shù)

h=hd.*windows′

freqz(h，1)；運行結(jié)果：

h=[

0.0015-0.0000-0.0025-0.00230.00330.0078-0.0000-0.0151

-0.01260.01680.0361-0.0000-0.0655-0.05750.09030.3001

0.40000.30010.0903-0.0575-0.0655-0.00000.03610.0168

-0.0126-0.0151-0.00000.00780.0033-0.0023-0.0025-0.00000.0015]解法三：用MATLAB的fir1函數(shù)設(shè)計。

MATLAB中基于窗函數(shù)的設(shè)計函數(shù)為fir1，它用于設(shè)計標(biāo)準的低通、高通、帶通、帶阻線性相位數(shù)字濾波器，其函數(shù)具有以下幾種調(diào)用形式：

b=fir1(n，Wn)

b=fir1(n，Wn，′ftype′)

b=fir1(n，Wn，window)

b=fir1(n，Wn，′ftype′，window)

函數(shù)的形式參數(shù)n為所要設(shè)計的濾波器的階數(shù)。

函數(shù)的形式參數(shù)Wn為所要設(shè)計的濾波器的歸一化邊界頻率，對于低通和高通，歸一化邊界頻率為一元向量，而對于帶通和帶阻，歸一化邊界頻率為二元向量；如果是多元向量，則為多帶濾波器。函數(shù)的形式參數(shù)ftype為所要設(shè)計的濾波器的類型。當(dāng)不指定ftype時，即采用第1種函數(shù)形式，歸一化邊界頻率為一元向量時默認為低通，為二元向量時默認為帶通，采用的是漢明窗函數(shù)設(shè)計。當(dāng)指定ftype時，即采用第2種函數(shù)形式。當(dāng)指定ftype為“high”時，表示設(shè)計高通濾波器；為“stop”時，表示設(shè)計帶阻濾波器；為“DC-1”時，表示設(shè)計多帶濾波器，并且第一個頻帶為通帶；為“DC-0”時，表示設(shè)計多帶濾波器，并且第一個頻帶為阻帶。

函數(shù)的形式參數(shù)window為所設(shè)計方法中使用的窗函數(shù)的類型。N點窗函數(shù)的調(diào)用是boxcar(N)、hanning(N)、hamming(N)、blackman(N)、kaiser(N，β)。如果不指定window選項，則默認的是采用漢明窗函數(shù)設(shè)計。

另外，對于高通和帶阻濾波器，fir1函數(shù)設(shè)計的結(jié)果都是偶數(shù)階的。利用fir1函數(shù)進行本例題的設(shè)計，根據(jù)方法一中的(1)~(4)步驟計算可知濾波器的階數(shù)為32，歸一化邊界頻率為0.4，采用漢明窗函數(shù)，則程序代碼如下：

n=32；

Wn=0.4；

b=fir1(n，Wn，hamming(n+1))

[b1，w1]=freqz(b，1)；

plot(w1/pi，20*log10(abs(b1)))；

axis([0，1，-80，10])；

grid；

xlabel(

′歸一化頻率/＼pi′)

ylabel(′幅度/dB′)設(shè)計結(jié)果：

b=[0.0015-0.0000-0.0025-0.00230.00330.0078-0.0000-0.0151

-0.01260.01680.0361-0.0000-0.0654-0.05750.09020.2998

0.39970.29980.0902-0.0575-0.0654-0.00000.03610.0168-0.0126

-0.0151-0.00000.00780.0033-0.0023-0.0025-0.00000.0015]

所得到的濾波器特性分析結(jié)果如圖5.2.9所示，由圖可見，其特性完全符合題目的技術(shù)指標(biāo)要求。圖5.2.9漢明窗函數(shù)設(shè)計的濾波器的幅頻特性

【例5.2.2】

用凱塞窗設(shè)計一線性相位FIR低通濾波器，通帶邊界頻率為0.3π，阻帶邊界頻率為0.5π，通帶波動小于2dB，阻帶衰減不小于50dB。

解根據(jù)題目可知所需設(shè)計的是低通濾波器，其技術(shù)指標(biāo)為：通帶邊界頻率ωp=0.3π，阻帶邊界頻率ωr=0.5π，δ=2dB，At=50dB。根據(jù)指標(biāo)要求，其理想低通的截止頻率應(yīng)為則理想低通濾波器的單位脈沖響應(yīng)為

濾波器的過渡帶為：Δω=ωr－ωp=0.2π。

β=0.1102(50－8.7)=4.55則設(shè)計的濾波器的單位脈沖響應(yīng)為與漢明窗函數(shù)設(shè)計相比，濾波器的階數(shù)變?yōu)?9，比較??；但是需要用到零階修正貝塞爾函數(shù)，計算復(fù)雜，需要計算機輔助計算。可采用MATLAB計算，程序代碼如下：

windows=kaiser(30，4.55)；

nn=[0：1：29]；

alfa=(30-1)/2；

hd=sin(0.4*pi*(nn-alfa))./(pi*(nn-alfa))；

h=hd.*windows′；

%上面5個語句可用fir1語句代替：h=fir1(29，0.4，kaiser(30，4.55))；

[h1，w1]=freqz(h，1)；

plot(w1/pi，20*log10(abs(h1)))；

axis([0，1，-80，10])； grid；

xlabel(′歸一化頻率/＼pi′) ylabel(′幅度/dB′)

圖5.2.10給出了同樣階數(shù)采用矩形窗函數(shù)和凱塞窗函數(shù)設(shè)計的濾波器的特性分析結(jié)果。結(jié)果表明，矩形窗設(shè)計的阻帶衰減在28dB左右，不符合指標(biāo)要求；而凱塞窗函數(shù)設(shè)計的阻帶衰減達50dB，符合指標(biāo)要求?？梢?，選擇不同的窗函數(shù)設(shè)計結(jié)果將不一樣。因此，采用窗函數(shù)法設(shè)計時一定要選擇合適的窗函數(shù)。圖5.2.10矩形窗、凱塞窗設(shè)計的濾波器的特性

【例5.2.3】

設(shè)計滿足下面技術(shù)指標(biāo)的數(shù)字帶通濾波器：

低頻端阻帶邊緣：ws1=0.2π，Rr=60dB

低頻端通帶邊緣：wp1=0.35π，Rp=1dB

高頻端通帶邊緣：wp2=0.65π，Rp=1dB

高頻端阻帶邊緣：ws2=0.8π，Rr=60dB

解根據(jù)阻帶衰減60dB選擇blackman窗函數(shù)，程序如下：

wp1=0.35；wp2=0.65；ws1=0.2；ws2=0.8；

deltaw=min((wp1-ws1)，(ws2-wp2))；%求過渡帶的小者

wc1=(ws1+wp1)/2；wc2=(wp2+ws2)/2；

N0=ceil(11/deltaw)%按blackman窗求濾波器長度N0

N=N0+mod(N0+1，2)%使濾波器長度N為奇數(shù)

h=fir1(N，[wc1，wc2]，blackman(N+1))；

freqz(h，1)

程序運行結(jié)果如圖5.2.11所示。圖5.2.11例5.2.3程序運行結(jié)果

5.3頻率采樣法

5.3.1設(shè)計方法

假設(shè)所要設(shè)計的目標(biāo)濾波器為Hd(ejω)，根據(jù)頻率采樣定理可知，如果對Hd(ejω)進行N點等間隔采樣，則可由這些采樣點的值來逼近所要設(shè)計的濾波器。假設(shè)這些采樣點的值為H(k)，經(jīng)過N點IDFT(IFFT)作為所要設(shè)計的濾波器的單位脈沖響應(yīng)h(n)，即(5.3.1)如果要設(shè)計具有線性相位的FIR數(shù)字濾波器，就要考慮線性相位的FIR濾波器必須滿足的條件。根據(jù)5.1節(jié)的分析結(jié)果，下面討論四種類型的線性相位FIR濾波器的約束條件。

類型1：h(n)偶對稱，N為奇數(shù)，其幅度函數(shù)H(ω)應(yīng)具有偶對稱性，即

H(ω)=H(2π－ω)

(5.3.2)

則對H(ejω)=H(ω)e－jω[(N－1)/2]的N點采樣可以表示為H(k)=

Hkejθk，可以導(dǎo)出：(5.3.3)(5.3.4)類型2：h(n)偶對稱，N為偶數(shù)，幅度特性具有奇對稱，即H(ω)=－H(2π－ω)，則可以導(dǎo)出：

Hk=－HN－k

(5.3.5)(5.3.6)類型3：h(n)奇對稱，N為奇數(shù)，幅度特性是奇對稱，則可以導(dǎo)出：(5.3.7)(5.3.8)類型4：h(n)奇對稱，N為偶數(shù)，幅度函數(shù)H(ω)偶對稱，則可以導(dǎo)出：(5.3.9)(5.3.10)5.3.2逼近誤差

下面討論上述設(shè)計過程的逼近誤差。對式(5.3.1)求Z變換，得(5.3.11)令W=e－j2π/N，則有式(5.3.11)是頻率采樣法的理論基礎(chǔ)。下面討論其頻率響應(yīng)，在單位圓上的頻響為(5.3.12)這是一個內(nèi)插公式，式中：(5.3.13)稱為內(nèi)插函數(shù)，當(dāng)ω=(2π/N)i，i=0，1，…，N－1時，有可見，在每個采樣點上，所設(shè)計的濾波器的頻響嚴格地與目標(biāo)濾波器的頻響的采樣值H(k)相等，即H(ejωk)=H(k)，逼近誤差為零。在采樣點之間，頻響由各采樣點的內(nèi)插函數(shù)延伸迭加而形成，因而有一定的逼近誤差。誤差的大小與目標(biāo)濾波器頻率響應(yīng)的曲線形狀有關(guān)，響應(yīng)曲線越平滑，誤差越?。环粗?，誤差越大。在目標(biāo)濾波器頻率響應(yīng)的不連續(xù)點附近，會產(chǎn)生肩峰和波紋。如圖5.3.1所示，梯形采樣值變化比較平緩，產(chǎn)生的肩峰和波紋比矩形采樣小得多。另外，當(dāng)N增大時，采樣點變密，逼近誤差也將減小。圖5.3.1不同頻率采樣的頻率特性

【例5.3.1】

設(shè)計一個FIR數(shù)字低通濾波器，其技術(shù)指標(biāo)為：通帶截止頻率為0.5π，阻帶最小衰減為50dB，N=33，求線性相位。

解以理想的低通濾波器特性為目標(biāo)濾波器模型，則根據(jù)表5.1.1，能設(shè)計低通線性相位數(shù)字濾波器的只有一、二兩種，因N為奇數(shù)，所以只能選擇類型1，即h(n)=h(N－1－n)，幅頻特性偶對稱。根據(jù)指標(biāo)要求，在0～2π內(nèi)有33個取樣點，所以第k點對應(yīng)頻率為(2π/33)k。令(2π/33)kc=ωc=0.5π，則kc=8.25;因為kc必須為整數(shù)，所以取截止頻率0.5π對應(yīng)的采樣點kc=8。當(dāng)幅度采樣點在k=0～8時,采樣值為1；由偶對稱性Hk=HN－k,可得k=25～32時，取樣值也為1，其他取樣值為0，如圖5.3.2所示。

由圖5.3.2可以寫出幅度采樣值的表達式為由式(5.3.4)可得圖5.3.2幅度采樣值所以把H(k)帶入式(5.3.11)得設(shè)計結(jié)果H(z)。下面進行頻譜驗證：(5.3.14)將設(shè)計結(jié)果作頻率特性分析，得如圖5.3.3(b)所示的特性曲線，由圖5.3.3(b)可以看出，其過渡帶寬為一個頻率采樣間隔2π/N=2π/33，而最小阻帶衰減略小于20dB。對大多數(shù)應(yīng)用場合，阻帶衰減如此小的濾波器是不能令人滿意的，因此，必須采取措施增大阻帶衰減。有如下三種方法可以增大阻帶衰減：

(1)增加過渡帶采樣點。

通過增加一個或者多個過渡帶采樣點，使目標(biāo)濾波器的特性變?yōu)楸容^平滑，減少通帶和阻帶的紋波，以增大阻帶衰減。但是這種方法會使過渡帶變寬，假設(shè)增加m個過渡點，則過渡帶Δω將變?yōu)?/p>

(5.3.15)圖5.3.3不同過渡點的設(shè)計結(jié)果比較因此，這種方法是通過加寬過渡帶寬，以犧牲過渡帶換取阻帶衰減的增加。

例如：在本例中，可在k=9和k=24處各增加一個過渡帶采樣點H9=H24=0.5，使過渡帶寬增加到二個頻率采樣間隔4π/33，重新計算的H(ejω)見圖5.3.3(c)，其阻帶衰減增加到約－40dB，如圖5.3.3(d)所示。

(2)過渡點的優(yōu)化。

利用線性最優(yōu)化的方法確定過渡帶采樣點的值，得到要求的濾波器的最佳逼近(而不是盲目地設(shè)定一個過渡帶值)。

例如，本例中可以用簡單的梯度搜索法來選擇H9、H24，使通帶或阻帶內(nèi)的最大絕對誤差最小化。經(jīng)計算得H9=0.3904時，對應(yīng)的H(ejω)的幅頻特性阻帶衰減約為－50dB，比H9=0.5時的阻帶衰減明顯增大。如果還要進一步改善阻帶衰減，可以增加第二個甚至第三個過渡點，當(dāng)然也可用線性最優(yōu)化求取這些取樣值。

(3)增大N。

如果要進一步增加阻帶衰減，但又不增加過渡帶寬，可通過增加采樣點數(shù)N來實現(xiàn)。

例如，對于同樣的截止頻率ωc=0.5π，以N=65采樣，并在k=17和k=48插入經(jīng)阻帶衰減最優(yōu)化計算得到的過渡點采樣值H17=H48=0.5886，在k=18和k=47處插入經(jīng)阻帶衰減最優(yōu)化計算得到的過渡點采樣值H18=H47=0.1065，如圖5.3.3(e)所示。這時得到的H(ejω)，過渡帶為6π/65，沒有增加，而阻帶衰減則達到了-60dB以上，如圖5.3.3(f)所示。當(dāng)然，代價是濾波器階數(shù)增加，運算量增加。5.3.3MATLAB輔助設(shè)計

為了更好地理解和掌握頻率采樣設(shè)計法，根據(jù)5.2.1節(jié)的分析和討論，歸納線性相位FIR數(shù)字濾波器的設(shè)計步驟如下：

(1)確定所要設(shè)計的FIR數(shù)字濾波器類型和技術(shù)指標(biāo)。

(2)根據(jù)濾波器類型和技術(shù)指標(biāo)要求確定目標(biāo)濾波器的模型和參數(shù)。通常采用理想濾波器模型作為目標(biāo)濾波器的模型，參數(shù)包括：插入過渡采樣點的數(shù)量和濾波器階數(shù)N。

(3)求頻譜采樣值H(k)。根據(jù)濾波器類型和線性相位的約束條件，計算幅度采樣值和相位采樣值。

(4)利用IFFT計算得到h(n)——使用MATLAB的IFFT函數(shù)來計算。

(5)性能分析和優(yōu)化(MATLAB編程實現(xiàn))。

【例5.3.2】

設(shè)計一個線性相位FIR數(shù)字低通濾波器，其技術(shù)指標(biāo)為：通帶截止頻率為0.4π，阻帶截止頻率為0.6π，通帶衰減為2dB，阻帶最小衰減為50dB。

解法一：

(1)確定技術(shù)指標(biāo)。

通帶截止頻率為0.4π，阻帶截止頻率為0.6π，通帶衰減為2dB，阻帶最小衰減為50dB。

(2)根據(jù)指標(biāo)要求確定插入過渡點數(shù)量和濾波器階數(shù)N。

選擇插入1個過渡點，根據(jù)式(5.3.15)得過渡帶0.2π=4π/N

，因此可求出N=20。

(3)根據(jù)濾波器類型和約束條件，計算幅度采樣值和相位采樣值。由于要設(shè)計的濾波器類型為低通，選擇情況1設(shè)計，取N=21，幅度函數(shù)采樣Hk具有偶對稱性，即有Hk=HN－k。令ωp=0.4π=(2π/N)kp，則kp=4.2，因此，截止頻率點取kp=4，則幅度函數(shù)采樣為相位采樣為步驟(4)、(5)采用MATLAB編程實現(xiàn)，具體代碼如下：

N=21；wp=0.4*pi；dw=2*pi/N；hg1=0；

wk=round(wp/dw)；%通帶截止頻率的采樣位置點

fori=0：wk

hg1(i+1)=1；

end

hg2=0.3904；

fori=(wk+2)：(N-(wk+2))

hg3(i-wk-1)=0；

end

hg=[hg1，hg2，hg3，hg2，hg1]%頻率采樣的幅度值

i=0：N；

subplot(3，1，1)；

stem(i，hg)；ylabel(′Hk′)；

xlabel(′k′)；

fort=1：N；

ct=-(N-1)/N*pi*(t-1)；%頻率采樣的相位值

ck=exp(j*ct)；

hk(t)=hg(t)*ck；%頻率采樣的頻譜值即H(k)

end；

hn=ifft(hk)%IFFT

[h，w]=freqz(hn，1，512)；

amp=abs(h)；ampdb=20*(log10(amp))；

pha=angle(h)；

%figure；

subplot(3，1，2)；

plot(w/pi，ampdb)；

line([0，1]，[-2，-2])；

ylabel(′幅度(dB)′)；xlabel(′頻率(w/pi)′)；

axis([01-605])；line([0，1]，[-50，-50])；

subplot(3，1，3)；

plot(w/pi，pha)；

ylabel(′相位(ard)′)；xlabel(′頻率(w/pi)′)；

程序運行結(jié)果：

h(n)=｛0.0101，0.0020，-0.0192，-0.0186，0.0259，0.0440，-0.0306，-0.0941，0.0333，0.3143，0.4658，0.3143，0.0333，-0.0941，-0.0306，0.0440，0.0259，-0.0186，-0.0192，0.0020，0.0101｝

特性分析結(jié)果如圖5.3.4所示，由圖可見阻帶衰減不滿足要求，因此需要優(yōu)化。采用插入2個過渡點進行優(yōu)化，則0.2π=6π/N，求得N=30。取N=31，并取兩個過渡點值分別為0.5886、0.1065。優(yōu)化后特性分析結(jié)果如圖5.3.5所示，由圖可見阻帶衰減超出指標(biāo)要求。因此，設(shè)計結(jié)果不是最優(yōu)的。圖5.3.4插入一個過渡點的低通濾波器特性對于高通、帶通、帶阻或者是多帶線性相位FIR數(shù)字濾波器的設(shè)計，方法、步驟完全相同。圖5.3.5插入優(yōu)化二個過渡點的低通濾波器特性解法二：用MATLAB的fir2函數(shù)設(shè)計。

MATLAB提供頻率樣本法的設(shè)計函數(shù)fir2。它的典型調(diào)用方法為

h=fir2(M，f，A)

其中，M為FIR數(shù)字濾波器的階數(shù)(濾波器的長度為N＝M＋1)，數(shù)組f和A給出它的預(yù)期頻率響應(yīng)。算出的濾波器系數(shù)為h，其長度為N。矢量f包含各邊緣頻率，單位為π，即0.0≤f≤1.0。矢量A為各指定頻率上預(yù)期的幅度響應(yīng)。fir2函數(shù)缺省地使用hamming窗函數(shù)。

如果采樣其他窗函數(shù)(包括Boxcar、Hann、Bartlett、Blackman、Kaiser和Chebwin等)，則調(diào)用方法為

h=fir2(M，f，A，window(M+1))本例要求阻帶衰減為50dB，可采用hamming窗設(shè)計，由過渡帶0.2π可計算出N=6.6π/(0.2π)33，程序如下：

f=[0，0.4，0.6，1]；

A=[1，10^(-2/20)，10^(-50/20)，0]；

h=fir2(32，f，A)

freqz(h，1)；

程序運行結(jié)果如圖5.3.6所示。

h=[

0.00000.00030.0000-0.00050.00020.00060.0000

0.00180.0007

-0.00740.00100.02590.0015-0.0569

0.02310.28180.43830.28180.0231-0.05690.0015

0.02590.0010-0.00740.00070.00180.00000.0006

0.0002-0.00050.00000.00030.0000圖5.3.6采樣FIR2設(shè)計結(jié)果

5.4FIR數(shù)字濾波器的最優(yōu)化設(shè)計

5.4.1最優(yōu)設(shè)計準則

最優(yōu)化設(shè)計的前提是最優(yōu)設(shè)計準則的確定。在FIR數(shù)字濾波器的設(shè)計中常用的準則有：均方誤差最小化準則、最大誤差最小化準則。

(1)均方誤差最小化準則。

若以E(ejω)表示逼近誤差，則定義

E(ejω)=Hd(ejω)－H(ejω)

(5.4.1)

式中：Hd(ejω)為所要逼近的目標(biāo)濾波器的頻率響應(yīng)；H(ejω)為所設(shè)計的濾波器的頻率響應(yīng)。那么均方誤差為

(5.4.2)均方誤差最小化準則就是選擇一組采樣值，以使均方誤差最小。這一方法注重的是在整個-π～π頻率區(qū)間內(nèi)總誤差的全局最小，但不能保證局部頻率點的性能，有些頻率點可能會有較大的誤差。例如：矩形窗函數(shù)法FIR濾波器設(shè)計，因采用有限項的h(n)逼近理想的目標(biāo)濾波器hd(n)，所以其逼近誤差為式中因此可得(5.4.3)

由式(5.4.3)可知，誤差只取決于給定的hd(n)，與設(shè)計值h(n)無關(guān)，故是一個常數(shù)。可以證明，這是一個最小均方誤差。所以，矩形窗函數(shù)設(shè)計法是一個最小均方誤差FIR設(shè)計，該方法保證在整個-π～π頻率區(qū)間內(nèi)總誤差的全局最小，但不能保證局部頻率點的性能，如在間斷點處就會出現(xiàn)較大的過沖(Gibbs現(xiàn)象)誤差。在MATLAB中，采樣均方誤差最小化設(shè)計的函數(shù)是firls，使用格式：

h=FIRLS(M，F(xiàn)，A)

其中：M為FIR數(shù)字濾波器的階數(shù)(濾波器的長度為N＝M＋1)，數(shù)組f和A給出它的預(yù)期頻率響應(yīng)。算出的濾波器系數(shù)為h，其長度為N。矢量f包含各邊緣頻率，單位為π，即0.0≤f≤1.0。矢量A為各指定頻率上預(yù)期的幅度響應(yīng)。

【例5.4.1】

采樣均方誤差最小化實現(xiàn)例5.3.2指標(biāo)，即設(shè)計一個線性相位FIR數(shù)字低通濾波器，其技術(shù)指標(biāo)為：通帶截止頻率為0.4π，阻帶截止頻率為0.6π，通帶衰減為

2dB，阻帶最小衰減為50dB。

解本例要求阻帶衰減50dB，可采樣hamming窗計算濾波器階數(shù)N=33，程序如下：

f=[0，0.4，0.6，1]；

A=[1，1，0，0]；

h=firls(32，f，A)；

freqz(h，1)；

程序運行結(jié)果如圖5.4.1所示，滿足設(shè)計要求。圖5.4.1例5.4.1設(shè)計結(jié)果

(2)最大誤差最小化準則(也叫最佳一致逼近準則)。

最大誤差最小化準則的設(shè)計思想是使逼近誤差的最大值最小化，即

max|E(ejω)|=min

ω∈F

(5.4.4)

其中，F(xiàn)是根據(jù)要求預(yù)先給定的一個頻率取值范圍，可以是濾波器的通帶，也可以是阻帶。設(shè)計時可選擇N個頻率采樣值(或時域h(n)值)，在給定頻帶范圍內(nèi)使頻響的最大逼近誤差達到最小，通常達到最小時各個采樣點的逼近誤差也相等，因此也叫等波紋逼近。5.4.2切比雪夫最佳一致逼近

1.切比雪夫最佳一致逼近準則

假設(shè)希望設(shè)計的線性相位FIR數(shù)字濾波器的幅度特性為Hd(ω)，實際設(shè)計的FIR數(shù)字濾波器的幅度特性為H(ω)，其加權(quán)誤差E(ω)定義為

E(ω)=W(ω)[Hd(ω)－H(ω)]

(5.4.5)

式中，W(ω)為誤差加權(quán)函數(shù)，它與通帶或者阻帶逼近精度有關(guān)。當(dāng)要求的逼近精度高時，W(ω)取值大；當(dāng)要求的逼近精度低時，W(ω)取值小。設(shè)計時W(ω)可由設(shè)計的技術(shù)指標(biāo)確定。圖5.4.2低通濾波器的誤差分配例如：假設(shè)要設(shè)計的線性相位FIR數(shù)字濾波器是如圖5.4.2所

示的誤差容限的低通濾波器，其通帶最大波動為δ1，阻帶最大波動為δ2，則定義(5.4.6)式中k滿足：δ1=kδ2。式(5.4.5)中，Hd(ω)為理想濾波器的幅頻特性，對于理想低通有(5.4.7)式(5.4.5)中，H(ω)為所設(shè)計濾波器的幅頻特性，由表5.1.1中的第一種類型濾波器，可得：則(5.4.8)因此，切比雪夫逼近問題變?yōu)椋簩で笠唤M系數(shù)a(n)，n=0，1，…，M，使逼近誤差的最大值達到最小，即

max|E(ω)|=minω∈A

(5.4.9)式中A表示研究的通帶或阻帶?？梢姡?5.4.8)中的是一個M次多項式。因此，最佳一致逼近問題實際上是如何構(gòu)造一個由M次多項式，使其按式(5.4.9)準則逼近一個連續(xù)函數(shù)的問題。

切比雪夫理論解決了這個多項式的存在性、惟一性及如何構(gòu)造等一系列問題。構(gòu)造該多項式的方法是“交錯點組定理”。該定理表明：式(5.4.8)滿足最佳一致逼近的充分必要條件為：E(ω)在A上至少呈現(xiàn)M+2個“交錯”，并使得：

E(ωi)=－E(ωi－1)=max|E(ω)|

式中：ω0≤ω1≤ω2≤…≤ωM+1。

2.等波紋逼近法設(shè)計線性相位FIR濾波器

假設(shè)希望設(shè)計的濾波器為如圖5.4.2所示的線性相位低通濾波器，其逼近誤差函數(shù)為式(5.4.8)?？疾斓屯V波器的閉子集A包括區(qū)間0≤ω≤ωc和ωr≤ω≤π，因為濾波器頻響是逐段恒定的，所以對應(yīng)于誤差函數(shù)各峰值點的頻率也對應(yīng)于恰好滿足誤差容限時的頻率。由于cos(ωn)可以采用其冪級數(shù)展開的形式表示，因此，濾波器的幅頻特性H(ω)可以變換為(5.4.10)式中，ak為與單位脈沖響應(yīng)有關(guān)的常數(shù)?？梢?，H(ω)是M階余弦多項式，因此它在區(qū)間0<ω<π內(nèi)有M－1個最大值和最小值。對式(5.4.10)求導(dǎo)，得(5.4.11)可見，ω=0，ω=π也是它的極值點，因此，在區(qū)間0≤ω≤π內(nèi)H(ω)最多有M+1個極值點。如果考慮邊界頻率點ωc、ωr，使|H(ωc)|=1－δ1，|H(ωr)|=δ2，則該誤差函數(shù)具有M+2個交錯，滿足“交錯點組定理”。假設(shè)A

上的M+2個交錯點頻率ω0，ω1，...，ωM+1已知，則由(5.4.8)式可得到M+2個方程：

(5.4.12)式中：ρ=max|E(ω)|為極值點頻率對應(yīng)的誤差函數(shù)值，{ωi}(i=0，1，…，M+1)為極值點頻率。必須指出：由于ωc和ωr固定為交錯點，因而必為這些極值頻率中的一個，如果假設(shè)ωc=ωl(0<l<M+1)，則ωr=ωl+1。將式(5.4.12)寫成矩陣形式，得(5.4.13)求解上述方程組可得到全部系數(shù)a(0)，a(1)，…，a(M)及誤差ρ，再由a(n)求出濾波器的h(n)。但實際情況下，M+2個極值點頻率是未知的，并且直接求解上述非線性方程組也是比較困難的。利用數(shù)值分析的雷米茲(Remez)迭代算法可以求解該方程組，其算法步驟如下：

(1)在頻率子集A上等間隔地選取M+2個極值點頻率ω0，ω1，…，ωM+1作為初值，并計算ρ：(5.4.14)式中

由于初始值ωi并不剛好就是極值點，因此ρ并不是最佳估計誤差，而是相對于初始值產(chǎn)生的偏差。

(2)由{ωi}(i=0，1，…，M+1)求H(ω)和E(ω)。利用重心形式的拉格朗日插值公式計算H(ω)：(5.4.15)其中：再將H(ω)代入式(5.4.5)，求得E(ω)。

如果在頻帶A上，所有頻率都有|E(ω)|≤ρ，則ρ為所求，ω0，ω1，…，ωM+1為極值點頻率，得到設(shè)計結(jié)果。否則進入第(3)步。

(3)對上次確定的極值點頻率中的每一點，在其附近檢查是否在某一頻率處有|E(ω)|>ρ，如有，則以該頻率點作為新的局部極值點。對M+2個極值點頻率依次進行檢查，得到一組新的極值點頻率，重復(fù)上述步驟求出ρ、H(ω)、E(ω)。直到ρ的值改變很小，迭代結(jié)束，由最后一組極值點頻率求出h(n)。上述迭代算法直接計算相當(dāng)困難，一般利用計算機進行計算。這里，可采用MATLAB工具軟件進行輔助設(shè)計。MATLAB中有兩個函數(shù)：remezord和remez，用于雷米茲交替算法線性相位FIR數(shù)字濾波器設(shè)計。

remezord用于確定濾波器的階數(shù)，其函數(shù)形式為：

[N，fpts，mag，wt]=remezord(fedge，mval，dev，fs)；

其中函數(shù)的形式參數(shù)fedge、mval、dev、fs分別為邊界頻率、幅度值、偏差值(與幅度值的偏差)、采樣頻率；返回參數(shù)中N為濾波器的階數(shù)，fpts為邊界頻率點，mag為期望幅度值，wt為權(quán)向量。

remez用于濾波器的設(shè)計，其函數(shù)形式為：

b=remez(N，fpts，mag，wt)；

【例5.4.2】

利用雷米茲交替算法，設(shè)計一個線性相位低通FIR數(shù)字濾波器，其指標(biāo)為：通帶邊界頻率fc=800Hz，阻帶邊界fr=1000Hz，通帶波動0.5dB，阻帶最小衰減At=

40dB，采樣頻率fs=4000Hz。

解根據(jù)題目所給的技術(shù)指標(biāo)可直接編寫以下代碼進行設(shè)計:

fedge=[8001000]；

mval=[10]；

dev1=1-1/10^(0.5/20)%根據(jù)通帶波動0.5dB計算偏差值

dev2=1/10^(40/20)-0

%根據(jù)阻帶衰減40dB計算偏差值

dev=[dev1，dev2]

%邊界偏差值

%dev=[0.05590.01]；

fs=4000；

[N，fpts，mag，wt]=remezord(fedge，mval，dev，fs)

b=remez(N，fpts，mag，wt)

[h，w]=freqz(b，1，256)；

subplot(2，1，1)；

plot(w*2000/pi，20*log10(abs(h)))；

grid；

axis([0，fs/2，-70，5])；

xlabel(′頻率/Hz′)；

ylabel(′幅度/dB′)；

subplot(2，1，2)；

plot(w*2000/pi，angle(h))；

xlabel(′頻率/Hz′)；

ylabel(′相位/ard′)；設(shè)計結(jié)果得到：N=28，fpts=[00.40000.50001.0000]，mag=[1100]，wt=[1.0000，5.5900]。濾波器的特性如圖5.4.3所示。由圖可見，濾波器的技術(shù)指標(biāo)阻帶衰減不滿足要求，需要進一步的優(yōu)化。辦法是將濾波器的階數(shù)增大，先取N+1，重新設(shè)計看得到的設(shè)計結(jié)果符合不符合要求，如果不符合，則再增大N，直到設(shè)計指標(biāo)滿足要求為止。圖5.4.3雷米茲交替算法設(shè)計本題取濾波器的階數(shù)N+2=30階時，設(shè)計的濾波器的特性符合技術(shù)指標(biāo)的要求，如圖5.4.4所示。所得到的濾波器的單位脈沖響應(yīng)為：

b=[0.0073-0.0041-0.0203-0.01970.00530.0209

-0.0023-0.0320

-0.01130.04160.0359-0.0501

-0.08800.05580.31200.44220.31200.0558

-0.0880-0.05010.03590.0416-0.0113-0.0320

-0.00230.02090.0053-0.0197-0.0203-0.0041

0.0073]

雷米茲交替算法的優(yōu)點是：ωc和ωr可準確確定；H(ω)逼近誤差均勻分布，相同指標(biāo)下，濾波器所需階數(shù)低。圖5.4.4雷米茲交替算法設(shè)計的優(yōu)化結(jié)果

5.5IIR與FIR數(shù)字濾波器的比較

對于同樣技術(shù)指標(biāo)要求的數(shù)字濾波器，我們既可以采用IIR數(shù)字濾波器實現(xiàn)，也可以采用FIR數(shù)字濾波器實現(xiàn)，那么，設(shè)計時應(yīng)選擇哪一種呢？下面，對IIR與FIR數(shù)字濾波器進行比較。首先，從性能上來說，IIR濾波器的系統(tǒng)函數(shù)的極點可位于z平面單位圓內(nèi)的任何地方，因此可用較低的階數(shù)獲得高的選擇性，所用的存儲單元少，所以經(jīng)濟而效率高，但是這個高效率是以相位的非線性為代價的。其次，從結(jié)構(gòu)上看，IIR濾波器必須采用遞歸結(jié)構(gòu)，極點位置必須在單位圓內(nèi)，否則系統(tǒng)將不穩(wěn)定。第三，從設(shè)計方法上看，IIR濾波器可以借助于模擬濾波器的成果，因此一般都有有效的封閉形式的設(shè)計公式可供準確計算，計算工作量比較小，對計算工具的要求不高。可見，IIR與FIR數(shù)字濾波器各有優(yōu)缺點，應(yīng)根據(jù)具體使用情況選擇。為了方便實際應(yīng)用，表5.5.1給出了IIR與FIR數(shù)字濾波器的性能比較。

5.6MATLAB中的濾波器設(shè)計工具FDATool

5.6.1FDATool使用環(huán)境介紹

在MATLAB6.0以上的版本中，為使用者提供了一個圖形化的濾波器設(shè)計與分析工具——FDATool。不同的版本其工作界面略有差別，設(shè)計結(jié)果也不盡相同，下面以MATLAB2006a版本為例進行介紹。

利用FDATool這一工具，可以進行FIR和IIR數(shù)字濾波器的設(shè)計，并且能夠顯示數(shù)字濾波器的幅頻響應(yīng)、相位響應(yīng)以及零極點分布圖等；產(chǎn)生的數(shù)字濾波器系數(shù)在存儲為文件后，可以直接提供給DSP程序代碼調(diào)試工具——CCS或DSP存儲器，以完成實際的數(shù)字濾波器的程序調(diào)試，從而實現(xiàn)實際的濾波器。由于本課程未涉及數(shù)字信號處理芯片部分，因而僅介紹由FDATool生成數(shù)據(jù)文件的方法。在MATLAB命令窗中輸入命令“fdatool”，打開如圖5.6.1所示的FDATool工作界面。

工作界面的第一行為主菜單，分別為File(文件)、Edit(編輯)、Analysis(分析)、Targets(目標(biāo))、View(觀察)、Window(窗口)、Help(幫助)菜單項，每一菜單項又分為多個二級菜單項。工作界面的第二行為圖形按鈕菜單，是由主菜單中常用的菜單構(gòu)成的，方便用戶操作。

下方是工作主界面，它大體分為上、下兩個小界面，為了方便，我們分別稱之為上半界面、下半界面。最下面正中有一個【DesignFilter】按鈕，那是在設(shè)定所設(shè)計濾波器的全部參數(shù)后，指揮計算機進行設(shè)計的確認按鈕。圖5.6.1FDATool工作界面上半界面用來顯示濾波器設(shè)計結(jié)果信息。其中，右半部顯示當(dāng)前濾波器的圖形信息，它的顯示內(nèi)容受主菜單【Analysis】的二級子菜單控制；左半部顯示當(dāng)前濾波器的結(jié)構(gòu)信息，它的顯示內(nèi)容受主菜單【Edit】的二級子菜單【ConvertStructure】的控制，使用戶能得到不同濾波器結(jié)構(gòu)的系數(shù)。這些濾波器系數(shù)可以由主菜單【File】的二級子菜單【Export】導(dǎo)出，輸出變量的名稱可以指定，它們的目標(biāo)位置可以是MATLAB工作空間，或者是其他指定的文件。下半界面用來設(shè)置濾波器設(shè)計指標(biāo)及參數(shù)，顯示的內(nèi)容由最左邊的7個按鈕控制。默認顯示內(nèi)容由第7個【DesignFilter】按鈕控制，如圖5.6.1所示，用來指定濾波器設(shè)計指標(biāo)，可設(shè)置濾波器類型、濾波器階數(shù)、頻率參數(shù)、幅度參數(shù)等。其他6個按鈕從上到下分別是：【CreateaMultirateFilter】，用于建立多采樣率濾波器，如圖5.6.2所示；【Transformfilter】，用于濾波器的頻率轉(zhuǎn)換設(shè)計，如圖5.6.3所示；【SetQuantizationParameters】，用于設(shè)置量化參數(shù)，如圖5.6.4所示；【RealizeModel】，用于實現(xiàn)濾波器模型，如圖5.6.5所示；【Pole/ZeroEditor】，用于增加、刪除和修改濾波器的零點/極點，如圖5.6.6所示；【ImportFilterFromWorkspace】，用于導(dǎo)入濾波器，如圖5.6.7所示。圖5.6.2建立多采樣率濾波器界面圖5.6.3濾波器的頻率轉(zhuǎn)換設(shè)計界面圖5.6.4設(shè)置量化參數(shù)界面圖5.6.5實現(xiàn)濾波器模型界面圖5.6.6濾波器的零點/極點編輯界面圖5.6.7導(dǎo)入外部濾波器界面5.6.2利用FDATool設(shè)計數(shù)字濾波器

【例5.6.1】

利用FDATool工具，設(shè)計一個FIR帶通數(shù)字濾波器。指標(biāo)如下：采樣頻率為20kHz，低端通帶截止頻率為2.5kHz，高端通帶截止頻率為5.5kHz，通帶范圍內(nèi)波動小于1dB；低端阻帶截止頻率為2kHz，高端阻帶截止頻率為6kHz，阻帶衰減大于60dB。

解本題設(shè)計的參數(shù)設(shè)置界面如圖5.6.8所示。

利用FDATool設(shè)計工具進行FIR數(shù)字濾波器設(shè)計,步驟如下：

(1)確定濾波器種類、設(shè)計方法。根據(jù)任務(wù)，“ResponseType”選項框選擇“Bandpass”；“DesignMethod”選項框選擇“FIR”→“Window”，并在“Options”選項框中選擇“Kaiser”窗。圖5.6.8例5.6.1的參數(shù)設(shè)置界面

(2)確定濾波器階數(shù)。如果設(shè)計指標(biāo)中給定了濾波器的階數(shù)，則“FilterOrder”選項框應(yīng)選擇“Specifyorder”，并輸