蔚瑞(渤海大學 遼寧 錦州 121000)
FIR濾波器的等波紋最佳逼近法優化設計
蔚瑞
(渤海大學 遼寧 錦州 121000)
本文針對窗函數設計法和頻率采樣法在FIR濾波器設計中存在較大的資源浪費的問題提出采用等波紋最佳逼近法進行優化,并且采用MATLAB編程對窗函數設計法和該方法進行比較以求能直觀地驗證其優越性。
FIR濾波器優化設計;等波紋最佳逼近法;MATLAB仿真
在現代信號處理中,例如圖像處理、數據傳輸、雷達接收以及一些要求較高的系統對相位要求較為嚴格所以通常采用有限沖激響應濾波器即FIR濾波器。設計該濾波器的方法一般有窗函數設計法和頻率采樣法,雖然它們簡單方便,易于實現但都有缺點,總的來說是所設計的濾波器性價比低,所以FIR數字濾波器的最優化設計就顯得格外重要。
2.1 優化設計準則
優化設計有以下三種方法:均方誤差最小化準則、最大誤差最小化準則和切比雪夫最佳一致逼近。其中均方誤差最小準則就是選擇一組時域采樣值,以使均方誤差最小。這一方法注重的是在整個-π~π頻率區間內總誤差的全局最小,但不能保證局部頻率點的性能,有些頻點可能會有較大的誤差。最大誤差最小化準則(也叫最佳一致逼近準則)是通過改變N個頻率采樣值(或時域h(n)值),使頻響誤差在給定頻帶范圍內最大逼近誤差達到最小。
2.2 等波紋最佳逼近法
2.2.1含義及優點 使用的最佳化準則是“最大誤差最小化準則”,所謂的等波紋是指用此方法設計的FIR數字濾波器的幅頻響應在通帶和阻帶都是等波紋的,而且分別控制通帶和阻帶波紋幅度。最佳逼近是指在濾波器長度給定的條件下,使加權誤差波紋幅度最小化。優點:階數相同時,使濾波器的最大逼近誤差最小,也就是通帶最大衰減最小,阻帶最小衰減最大;指標相同時,可使濾波器階數最低。
2.2.2基本思想 等波紋最佳逼近基于切比雪夫逼近,在通帶和阻帶以|E(ω)|的最大值最小化為準則,采用Remez多重交換迭代算法求解濾波器系數h(n)。定義加權誤差函數E(ω)為

其中W(ω)稱為誤差加權函數,用來控制不同頻段(一般指通帶和阻帶)的逼近精度。表示希望逼近的幅度特性函數,設計線性相位FIR數字濾波器時必須滿足線性相位約束條件,表示實際設計的濾波器幅度特性函數。
在此方法設計中,把數字頻段分為“逼近(或研究)區域”和“無關區域”。逼近區域通常指通帶和阻帶,而無關區域一般指過渡帶。設計過程中只考慮對逼近區域的最佳逼近,但是無關區寬度不能為零,即不能是理想濾波特性。用等波紋最佳逼近法設計FIR數字濾波器的步驟為:
根據給定的逼近指標估算濾波器階數N和誤差加權函數W(ω);
采用remez算法得到濾波器單位脈沖響應h(n)。
3.1 已有的編程工具
MATLAB信號處理工具箱函數中為我們提供了采用經典窗函數法設計線性相位FIR數字濾波器的函數即fir1且具有標準低通、帶通、高通,帶阻等類型。其調用格式及功能如下:
? hn=fir1(M,wc),返回6dB截止頻率為wc的M階(單位脈沖響應h(n)長度N=M+1)FIR低通(wc為標量)濾波器系數向量hn,默認選用漢明窗。濾波器單位脈沖響應h(n)和向量hn的關系為:h(n)=hn(n+1) n=0,1,2,…,M。而且滿足線性相位條件:h(n)=hn(N-1-n),其中wc是對π歸一化的數字頻率,0≤wc≤1。
當wc=[wcl,wcu]時,得到的是帶通濾波器,其-6dB通帶為wcl≤ω≤wcu。
? hn=fir1(M,wc,‘ftype’),可設計高通和帶阻FIR濾波器。當ftype=high時,設計的是高通濾波器;當ftype=stop,且wc=[wcl,wcu]時,設計的是帶阻濾波器。在設計高通和帶阻FIR濾波器時,階數M只能取偶數(h(n)長度N=M+1為奇數)。不過,即使用戶將M設置為奇數時,fir1也會自動對M加1。
hn=fir1(M,wc,window),可以指定窗函數向量window。如果缺省window參數則默認為漢明窗。
hn=fir1(M,wc,‘ftype’,window),通過選擇wc,ftype和window參數可以設計各種加窗濾波器。
而remez和remezord是用來實現線性相位FIR數字濾波器的等波紋最佳逼近設計的MATLAB信號處理工具箱函數。下面介紹這兩種函數的調用格式和功能:
remez
remez函數實現線性相位FIR數字濾波器的等波紋最佳逼近設計。其調用格式有以下幾種分別為:
?hn=remez(M,f,m,w)%最常用的格式
調用結果返回單位脈沖響應向量hn。remez函數的調用參數(M,f,m,w)一般通過調用remezord函數來計算。調用參數含義為:M為FIR數字濾波器階數,hn 長度N=M+1。f是邊界頻率向量,0≤f≤1,要求f為單調增向量(即f(k)< f(k+1),k=1,2,…),且從0開始,以1結束,1對應數字頻率ω=π(模擬頻率,表示時域采樣頻率)。m是與f對應的幅度向量,m與f長度相同,m(k)表示頻點f(k)的幅度響應值。f和m給出希望逼近的幅度特性,w是誤差加權向量。
?hn=remez(M,f,m)
設計一個M階FIR數字濾波器,其頻率響應在數組f 和m中給定,含義與上述情形相似,只是數組w中的值均為1。
?hn=remez(M,f,m,w,ftype)
含義與上述情形均類似,只有當ftype是字符串“hilbert”或“differentiator”時,它相應地設計數字希爾伯特變換器或數字微分器。
對于數字希爾伯特變換器,數組f中的最低頻率不能等于0,最高頻率不能等于1;而對于數字微分器來說,矢量m不給出每個帶中預期的斜率,而是給出預期的幅度。
remezord
remezord函數可根據逼近指標估算等波紋最佳逼近FIR數字濾波器的最低階數M、誤差加權向量w和歸一化邊界頻率向量f,使濾波器在滿足指標的前提下造價最低。其返回參數作為remez函數的調用參數。其調用格式為:
[M,fo,mo,w]=remezord(f,m,rip,Fs)
參數含義說明:f與remez中的類似,這里f可以是模擬頻率(單位是Hz)或歸一化數字頻率,但必須從0開始,到(用歸一化頻率對應1)結束,而且其中省略了0 和兩個頻點。為采樣頻率,缺省時默認=2Hz。但是f的長度(包括省略的0和兩個頻點)是m的兩倍,即m中的每個元素表示f給定的一個逼近頻段上希望逼近的幅度值。rip表示f和m描述的各逼近頻段允許的波紋幅度(幅頻響應最大偏差),f的長度是rip的兩倍。
所以,調用這兩個函數設計線性相位FIR數字濾波器的關鍵在于要根據設計指標求出remezord函數的調用參數f、m、rip和。
3.2 編程驗證
以設計FIR數字低通濾波器為例,相應的指標為:

其中含義如下:
首先對所需指標進行估算可以得到:選擇凱塞窗可以降低濾波器階數。
3.2.1使用凱塞窗設計的MATLAB程序如下:
fp=1500;fs=2500;rs=40;Fs=10000;
wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;
Bt=ws-wp; %計算過渡帶寬度
alph=0.5842*(rs-21)^0.4+0.07886*(rs-21); %計算kaiser窗的控制參數a
M=ceil((rs-8)/2.285/Bt); %計算kaiser窗所需階數M
wc=(wp+ws)/2/pi;
hn=fir1(M,wc,kaiser(M+1,alph));%調用kaiser計算低通FIRDF的h(n)
運行結果如下圖:

可以看到M=23即使用kaiser窗設計的濾波器階數為23。
3.2.2使用等波紋最佳逼近法設計的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,Fs);%邊界頻率為模擬頻率(Hz)時必須加入采樣頻率
M=M+1; %估算的M值達不到要求,加1后滿足要求
hn=remez(M,fo,mo,w);
運行結果如下圖:

可以看到M=15即使用等波紋最佳逼近法設計的濾波器階數為15。
綜上,我們可以從MATLAB編程運行的結果直觀地看到:指標相同時,即使用等波紋最佳逼近法設計可使濾波器階數明顯降低。

3.2.3同時用窗函數和等波紋最佳逼近設計相同階數的FIR數字濾波器
以設計一個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);% 默認采用漢明窗設計濾波器
[h,w1]=freqz(hn1);%計算濾波器的頻率響應
f=[0,0.3,0.5,1];
m=[1,1,0,0];
hn2=remez(M,f,m);%采用remez設計濾波器
[hh,w2]=freqz(hn2);%計算濾波器的頻率響應
figure(1)
plot(w1/pi,abs(h),w2/pi,abs(hh),f,m);grid on%繪制濾波器幅頻響應
legend('fir1','remez','理想特性');%給出圖例
xlabel('w/pi');ylabel('幅度')
figure(2)
plot(w1/pi,20*log10(abs(h)),w2/pi,20*log10(abs(hh))); %繪制損耗函數
grid on
legend('fir1','remez');%給出圖例
xlabel('w/pi');ylabel('幅度/dB')
figure(3)
impz(hn1,n);grid on%繪制單位脈沖響應
xlabel('n');
ylabel('h(n)')
figure(4)
impz(hn2,n);grid on%繪制單位脈沖響應
xlabel('n');
ylabel('h(n)')
程序運行后可得:


幅頻響應損耗函數
結論:從幅頻響應可以直觀地看到用等波紋最佳逼近法設計FIR低通濾波器相較窗函數設計法更接近于理想特性。但是由于remez要求等波紋的特點,其在通帶內的振動幅度較大;由損耗函數可以清楚地看出阻帶的波紋變化曲線,remez接近等波紋,而fir1的振動幅度較大。并且比較窗函數設計法和等波紋最佳逼近法的濾波器指標可以看到,前者的過渡帶寬較寬,通帶最大衰減較小。
通過MATLAB仿真來直觀地驗證了使用等波紋最佳逼近法設計FIR濾波器時可以達到在指標相同時階數最低,而階數相同時濾波器的最大逼近誤差最小,所以等波紋最佳逼近法在FIR濾波器的設計中有重要意義。
[1]高西全,丁美玉.數字信號處理(第三版)[M].西安:西安電子科技大學出版社,2008.8
[2]萬永革.數字信號處理的MATLAB實現[M].北京:科學出版社,2007
[3]張德豐.MATLAB數字信號處理與應用[M].北京:清華大學出版社,2010.1
TP3
A
1009-5624(2016)06-0045-04