FIR濾波器的等波紋最佳逼近法優(yōu)化設(shè)計(jì)

蔚瑞(渤海大學(xué) 遼寧 錦州 121000)

FIR濾波器的等波紋最佳逼近法優(yōu)化設(shè)計(jì)

蔚瑞
(渤海大學(xué) 遼寧 錦州 121000)

本文針對(duì)窗函數(shù)設(shè)計(jì)法和頻率采樣法在FIR濾波器設(shè)計(jì)中存在較大的資源浪費(fèi)的問題提出采用等波紋最佳逼近法進(jìn)行優(yōu)化，并且采用MATLAB編程對(duì)窗函數(shù)設(shè)計(jì)法和該方法進(jìn)行比較以求能直觀地驗(yàn)證其優(yōu)越性。

FIR濾波器優(yōu)化設(shè)計(jì)；等波紋最佳逼近法；MATLAB仿真

1.引言

在現(xiàn)代信號(hào)處理中，例如圖像處理、數(shù)據(jù)傳輸、雷達(dá)接收以及一些要求較高的系統(tǒng)對(duì)相位要求較為嚴(yán)格所以通常采用有限沖激響應(yīng)濾波器即FIR濾波器。設(shè)計(jì)該濾波器的方法一般有窗函數(shù)設(shè)計(jì)法和頻率采樣法，雖然它們簡(jiǎn)單方便，易于實(shí)現(xiàn)但都有缺點(diǎn)，總的來說是所設(shè)計(jì)的濾波器性價(jià)比低，所以FIR數(shù)字濾波器的最優(yōu)化設(shè)計(jì)就顯得格外重要。

2.FIR數(shù)字濾波器的優(yōu)化設(shè)計(jì)

2.1 優(yōu)化設(shè)計(jì)準(zhǔn)則

優(yōu)化設(shè)計(jì)有以下三種方法：均方誤差最小化準(zhǔn)則、最大誤差最小化準(zhǔn)則和切比雪夫最佳一致逼近。其中均方誤差最小準(zhǔn)則就是選擇一組時(shí)域采樣值，以使均方誤差最小。這一方法注重的是在整個(gè)-π～π頻率區(qū)間內(nèi)總誤差的全局最小，但不能保證局部頻率點(diǎn)的性能，有些頻點(diǎn)可能會(huì)有較大的誤差。最大誤差最小化準(zhǔn)則（也叫最佳一致逼近準(zhǔn)則）是通過改變N個(gè)頻率采樣值（或時(shí)域h（n）值），使頻響誤差在給定頻帶范圍內(nèi)最大逼近誤差達(dá)到最小。

2.2 等波紋最佳逼近法

2.2.1含義及優(yōu)點(diǎn) 使用的最佳化準(zhǔn)則是“最大誤差最小化準(zhǔn)則”，所謂的等波紋是指用此方法設(shè)計(jì)的FIR數(shù)字濾波器的幅頻響應(yīng)在通帶和阻帶都是等波紋的，而且分別控制通帶和阻帶波紋幅度。最佳逼近是指在濾波器長(zhǎng)度給定的條件下，使加權(quán)誤差波紋幅度最小化。優(yōu)點(diǎn)：階數(shù)相同時(shí)，使濾波器的最大逼近誤差最小，也就是通帶最大衰減最小，阻帶最小衰減最大；指標(biāo)相同時(shí)，可使濾波器階數(shù)最低。

2.2.2基本思想 等波紋最佳逼近基于切比雪夫逼近，在通帶和阻帶以|E(ω)|的最大值最小化為準(zhǔn)則，采用Remez多重交換迭代算法求解濾波器系數(shù)h（n）。定義加權(quán)誤差函數(shù)E（ω）為

其中W(ω)稱為誤差加權(quán)函數(shù)，用來控制不同頻段（一般指通帶和阻帶）的逼近精度。表示希望逼近的幅度特性函數(shù)，設(shè)計(jì)線性相位FIR數(shù)字濾波器時(shí)必須滿足線性相位約束條件，表示實(shí)際設(shè)計(jì)的濾波器幅度特性函數(shù)。

在此方法設(shè)計(jì)中，把數(shù)字頻段分為“逼近（或研究）區(qū)域”和“無(wú)關(guān)區(qū)域”。逼近區(qū)域通常指通帶和阻帶，而無(wú)關(guān)區(qū)域一般指過渡帶。設(shè)計(jì)過程中只考慮對(duì)逼近區(qū)域的最佳逼近，但是無(wú)關(guān)區(qū)寬度不能為零，即不能是理想濾波特性。用等波紋最佳逼近法設(shè)計(jì)FIR數(shù)字濾波器的步驟為：

根據(jù)給定的逼近指標(biāo)估算濾波器階數(shù)N和誤差加權(quán)函數(shù)W(ω)；

采用remez算法得到濾波器單位脈沖響應(yīng)h(n)。

3.計(jì)算機(jī)模擬

3.1 已有的編程工具

MATLAB信號(hào)處理工具箱函數(shù)中為我們提供了采用經(jīng)典窗函數(shù)法設(shè)計(jì)線性相位FIR數(shù)字濾波器的函數(shù)即fir1且具有標(biāo)準(zhǔn)低通、帶通、高通，帶阻等類型。其調(diào)用格式及功能如下：

? hn=fir1(M，wc)，返回6dB截止頻率為wc的M階（單位脈沖響應(yīng)h（n）長(zhǎng)度N=M+1）FIR低通（wc為標(biāo)量）濾波器系數(shù)向量hn，默認(rèn)選用漢明窗。濾波器單位脈沖響應(yīng)h(n)和向量hn的關(guān)系為：h（n）=hn（n+1） n=0，1，2，…，M。而且滿足線性相位條件：h（n）=hn（N-1-n），其中wc是對(duì)π歸一化的數(shù)字頻率，0≤wc≤1。

當(dāng)wc=[wcl，wcu]時(shí)，得到的是帶通濾波器，其-6dB通帶為wcl≤ω≤wcu。

? hn=fir1(M，wc，‘ftype’)，可設(shè)計(jì)高通和帶阻FIR濾波器。當(dāng)ftype=high時(shí)，設(shè)計(jì)的是高通濾波器；當(dāng)ftype=stop，且wc=[wcl，wcu]時(shí)，設(shè)計(jì)的是帶阻濾波器。在設(shè)計(jì)高通和帶阻FIR濾波器時(shí)，階數(shù)M只能取偶數(shù)（h（n）長(zhǎng)度N=M+1為奇數(shù)）。不過，即使用戶將M設(shè)置為奇數(shù)時(shí)，fir1也會(huì)自動(dòng)對(duì)M加1。

hn=fir1(M，wc，window)，可以指定窗函數(shù)向量window。如果缺省window參數(shù)則默認(rèn)為漢明窗。

hn=fir1(M，wc，‘ftype’，window)，通過選擇wc，ftype和window參數(shù)可以設(shè)計(jì)各種加窗濾波器。

而remez和remezord是用來實(shí)現(xiàn)線性相位FIR數(shù)字濾波器的等波紋最佳逼近設(shè)計(jì)的MATLAB信號(hào)處理工具箱函數(shù)。下面介紹這兩種函數(shù)的調(diào)用格式和功能：

remez

remez函數(shù)實(shí)現(xiàn)線性相位FIR數(shù)字濾波器的等波紋最佳逼近設(shè)計(jì)。其調(diào)用格式有以下幾種分別為：

?hn=remez（M，f，m，w）%最常用的格式

調(diào)用結(jié)果返回單位脈沖響應(yīng)向量hn。remez函數(shù)的調(diào)用參數(shù)（M，f，m，w）一般通過調(diào)用remezord函數(shù)來計(jì)算。調(diào)用參數(shù)含義為：M為FIR數(shù)字濾波器階數(shù)，hn 長(zhǎng)度N=M+1。f是邊界頻率向量，0≤f≤1，要求f為單調(diào)增向量（即f（k）< f（k+1），k=1，2，…），且從0開始，以1結(jié)束，1對(duì)應(yīng)數(shù)字頻率ω=π（模擬頻率，表示時(shí)域采樣頻率）。m是與f對(duì)應(yīng)的幅度向量，m與f長(zhǎng)度相同，m（k）表示頻點(diǎn)f（k）的幅度響應(yīng)值。f和m給出希望逼近的幅度特性，w是誤差加權(quán)向量。

?hn=remez(M，f，m)

設(shè)計(jì)一個(gè)M階FIR數(shù)字濾波器，其頻率響應(yīng)在數(shù)組f 和m中給定，含義與上述情形相似，只是數(shù)組w中的值均為1。

?hn=remez（M，f，m，w，ftype）

含義與上述情形均類似，只有當(dāng)ftype是字符串“hilbert”或“differentiator”時(shí)，它相應(yīng)地設(shè)計(jì)數(shù)字希爾伯特變換器或數(shù)字微分器。

對(duì)于數(shù)字希爾伯特變換器，數(shù)組f中的最低頻率不能等于0，最高頻率不能等于1；而對(duì)于數(shù)字微分器來說，矢量m不給出每個(gè)帶中預(yù)期的斜率，而是給出預(yù)期的幅度。

remezord

remezord函數(shù)可根據(jù)逼近指標(biāo)估算等波紋最佳逼近FIR數(shù)字濾波器的最低階數(shù)M、誤差加權(quán)向量w和歸一化邊界頻率向量f，使濾波器在滿足指標(biāo)的前提下造價(jià)最低。其返回參數(shù)作為remez函數(shù)的調(diào)用參數(shù)。其調(diào)用格式為：

[M，fo，mo，w]=remezord(f，m，rip，F(xiàn)s)

參數(shù)含義說明:f與remez中的類似，這里f可以是模擬頻率（單位是Hz）或歸一化數(shù)字頻率，但必須從0開始，到（用歸一化頻率對(duì)應(yīng)1）結(jié)束，而且其中省略了0 和兩個(gè)頻點(diǎn)。為采樣頻率，缺省時(shí)默認(rèn)=2Hz。但是f的長(zhǎng)度（包括省略的0和兩個(gè)頻點(diǎn)）是m的兩倍，即m中的每個(gè)元素表示f給定的一個(gè)逼近頻段上希望逼近的幅度值。rip表示f和m描述的各逼近頻段允許的波紋幅度（幅頻響應(yīng)最大偏差），f的長(zhǎng)度是rip的兩倍。

所以，調(diào)用這兩個(gè)函數(shù)設(shè)計(jì)線性相位FIR數(shù)字濾波器的關(guān)鍵在于要根據(jù)設(shè)計(jì)指標(biāo)求出remezord函數(shù)的調(diào)用參數(shù)f、m、rip和。

3.2 編程驗(yàn)證

以設(shè)計(jì)FIR數(shù)字低通濾波器為例，相應(yīng)的指標(biāo)為：

其中含義如下：

首先對(duì)所需指標(biāo)進(jìn)行估算可以得到：選擇凱塞窗可以降低濾波器階數(shù)。

3.2.1使用凱塞窗設(shè)計(jì)的MATLAB程序如下：

fp=1500；fs=2500；rs=40；Fs=10000；

wp=2*pi*fp/Fs；ws=2*pi*fs/Fs；

Bt=ws-wp； %計(jì)算過渡帶寬度

alph=0.5842*(rs-21)^0.4+0.07886*(rs-21)； %計(jì)算kaiser窗的控制參數(shù)a

M=ceil((rs-8)/2.285/Bt)； %計(jì)算kaiser窗所需階數(shù)M

wc=(wp+ws)/2/pi；

hn=fir1(M，wc，kaiser(M+1，alph))；%調(diào)用kaiser計(jì)算低通FIRDF的h（n）

運(yùn)行結(jié)果如下圖：

可以看到M=23即使用kaiser窗設(shè)計(jì)的濾波器階數(shù)為23。

3.2.2使用等波紋最佳逼近法設(shè)計(jì)的MATLAB程序如下：

Fs=10000；f=[1500，2500]；m=[1，0]；

rp=1；rs=40；

dat1=(10^(rp/20)-1)/(10^(rp/20)+1)；

dat2=10^(-rs/20)；

rip=[dat1，dat2]；

[M，fo，mo，w]=remezord(f，m，rip，F(xiàn)s)；%邊界頻率為模擬頻率（Hz）時(shí)必須加入采樣頻率

M=M+1； %估算的M值達(dá)不到要求，加1后滿足要求

hn=remez(M，fo，mo，w)；

運(yùn)行結(jié)果如下圖：

可以看到M=15即使用等波紋最佳逼近法設(shè)計(jì)的濾波器階數(shù)為15。

綜上，我們可以從MATLAB編程運(yùn)行的結(jié)果直觀地看到：指標(biāo)相同時(shí)，即使用等波紋最佳逼近法設(shè)計(jì)可使濾波器階數(shù)明顯降低。

3.2.3同時(shí)用窗函數(shù)和等波紋最佳逼近設(shè)計(jì)相同階數(shù)的FIR數(shù)字濾波器

以設(shè)計(jì)一個(gè)20階FIR低通濾波器為例，要求如下：

MATLAB程序如下：

fp=1500；fs=2500；Fs=10000；

wp=2*pi*fp/Fs；ws=2*pi*fs/Fs；

M=20；n=M+1；

wc=(wp+ws)/2/pi；

hn1=fir1(M，wc)；% 默認(rèn)采用漢明窗設(shè)計(jì)濾波器

[h，w1]=freqz(hn1)；%計(jì)算濾波器的頻率響應(yīng)

f=[0，0.3，0.5，1]；

m=[1，1，0，0]；

hn2=remez(M，f，m)；%采用remez設(shè)計(jì)濾波器

[hh，w2]=freqz(hn2)；%計(jì)算濾波器的頻率響應(yīng)

figure(1)

plot(w1/pi，abs(h)，w2/pi，abs(hh)，f，m)；grid on%繪制濾波器幅頻響應(yīng)

legend('fir1'，'remez'，'理想特性')；%給出圖例

xlabel('w/pi')；ylabel('幅度')

figure(2)

plot(w1/pi，20*log10(abs(h))，w2/pi，20*log10(abs(hh)))； %繪制損耗函數(shù)

grid on

legend('fir1'，'remez')；%給出圖例

xlabel('w/pi')；ylabel('幅度/dB')

figure(3)

impz(hn1，n)；grid on%繪制單位脈沖響應(yīng)

xlabel('n')；

ylabel('h(n)')

figure(4)

impz(hn2，n)；grid on%繪制單位脈沖響應(yīng)

xlabel('n')；

ylabel('h(n)')

程序運(yùn)行后可得：

幅頻響應(yīng)損耗函數(shù)

結(jié)論：從幅頻響應(yīng)可以直觀地看到用等波紋最佳逼近法設(shè)計(jì)FIR低通濾波器相較窗函數(shù)設(shè)計(jì)法更接近于理想特性。但是由于remez要求等波紋的特點(diǎn)，其在通帶內(nèi)的振動(dòng)幅度較大；由損耗函數(shù)可以清楚地看出阻帶的波紋變化曲線，remez接近等波紋，而fir1的振動(dòng)幅度較大。并且比較窗函數(shù)設(shè)計(jì)法和等波紋最佳逼近法的濾波器指標(biāo)可以看到，前者的過渡帶寬較寬，通帶最大衰減較小。

4.結(jié)束語(yǔ)

通過MATLAB仿真來直觀地驗(yàn)證了使用等波紋最佳逼近法設(shè)計(jì)FIR濾波器時(shí)可以達(dá)到在指標(biāo)相同時(shí)階數(shù)最低，而階數(shù)相同時(shí)濾波器的最大逼近誤差最小，所以等波紋最佳逼近法在FIR濾波器的設(shè)計(jì)中有重要意義。

[1]高西全，丁美玉．?dāng)?shù)字信號(hào)處理(第三版)[M]．西安:西安電子科技大學(xué)出版社，2008.8

[2]萬(wàn)永革.數(shù)字信號(hào)處理的MATLAB實(shí)現(xiàn)[M].北京：科學(xué)出版社，2007

[3]張德豐.MATLAB數(shù)字信號(hào)處理與應(yīng)用[M].北京：清華大學(xué)出版社，2010.1

TP3

A

1009-5624（2016）06-0045-04