FIR數字濾波器設計（下）
今天給大俠帶來FIR數字濾波器設計，由於篇幅較長，分三篇。今天帶來第三篇，FIR數字濾波器設計，包括窗函式法設計FIR濾波器、頻率取樣法設計FIR濾波器以及基於firls函式和remez函式的最優化方法設計FIR濾波器。話不多說，上貨。
這裡也超連結了上、中兩篇，方便參考學習。
數字濾波器的輸入輸出均為數字訊號，訊號通過數字濾波器後，可以改變頻率成分的相對比例或濾除某些頻率成分。數字濾波器可以分為IIR數字濾波器和FIR數字濾波器。
本篇只介紹FIR數字濾波器的設計，可以根據所給定的頻率特性直接設計FIR數字濾波器。FIR數字濾波器在保證幅度特性滿足要求的同時，能夠做到嚴格的線性特性。
本篇採用了窗函式法、頻率取樣法以及基於firls函式和remez函式的最優化方法設計FIR濾波器。對FIR濾波器進行了詳細的理論分析，並且對應於每種方法都給出了設計例項。通過編寫MATLAB語言程式，執行程式，得到幅頻和相頻特性圖。
對於窗函式和firls函式設計的濾波器，還通過建立Simulink系統模組進行模擬，觀察濾波器濾波情況。
FIR數字濾波器設計
FIR數字濾波器的設計方法有窗函式法、頻率取樣法和基於firls函式和remez函式的最優化方法。MATLAB語言中的數字訊號處理工具箱，提供了一些濾波器的函式，使FIR濾波器的運算更加方便和快捷。
在MATLAB中提供的濾波函式有fir1()，此函式以經典的方法實現加窗線性相位FIR數字濾波器設計，可以設計出低通、高通、帶通和帶阻濾波器；fir2函式設計的FIR濾波器，其濾波的頻率特性由向量f和m決定，f和m分別為濾波器的期望幅頻響應的頻率相量和幅值相量。
Firls()和remez()的基本格式用於設計I型和II型線性相位FIR濾波器，I型和II型的區別是偶函式還是奇函式。freqz()用於求數字濾波器的頻率響應。並且提供了各種窗函式的函式，比如，hamming()是海明窗函式，hanning()是漢寧窗函式，kaiser()是凱澤窗函式，使在設計的過程中，不用自己重新設計窗函式。
Simulink是MATLAB眾多工具包中的一員，對於建模，Simulink提供了一個圖形化的使用者介面（GUI）。Simulink包括一個複雜的由接收器、訊號源、線性和非線性元件以及連線件組成的模組庫。定義完一個模型後，就可以通過Simulink的選單或者在MATLAB的命令視窗輸入命令對它進行模擬。使用Scopes或者其他的顯示模組，可以在執行模擬時觀察到模擬的結果。另外，還可以在模擬時改變引數，並且立即就可看到變化。
一、窗函式法設計FIR濾波器
設我們所要設計的FIR濾波器的傳輸函式是，是與其對應的單位脈衝響應，因此
（3-1）
（3-2)
如果我們能夠在已知的情況下，求出，經過Z變換可得到濾波器的系統函式。通常情況下理想數字濾波器的單位脈衝相應是無限長的，且是非因果序列。獲得有限脈衝響應濾波器的一種可能方法是對擷取一段來近似代替，可是這樣會改變原來的濾波器指標，出現吉布斯效應誤差。
窗函式法就是用被稱為窗函式的有限加權序列w(n)來修正式(1)的傅立葉基數以求得要求的有限脈衝響應序列，即：
（3-3）
w(n)是有限長序列，當n<0或n>N-1時，w(n)=0。
這種方法的重點在於選擇某種合適的窗函式。要求窗函式主瓣寬度儘可能窄，以獲得最小的過渡帶；旁瓣相對值儘可能小，以使得通帶波紋小，並且阻帶衰減大。
下面介紹幾種常用的窗函式：
1、矩形窗(Rectangle Window)
（3-4）
其頻率函式為：
（3-5）
2、三角形窗(Bartlett Window)
（3-6）
其頻率函式為：
（3-7)
3、漢寧(Hanning)窗，又稱升餘弦窗
（3-8）
利用傅立葉變換得到頻率函式為：
（3-9）
當時，，所以窗函式的幅度函式為：
(3-10)
4、漢明(Hamming)窗，又稱改進的升餘弦窗
（3-11）
其幅度函式為：
（3-12）
5、布萊克曼(Blankman)窗，又稱二階升餘弦窗
（3-13)
其幅度函式為：
(3-14）
6、凱澤(Kaiser)窗
（3-15）
其中：β是一個可自由選擇的引數，是第一類修正零階貝塞爾函式。
上述窗函式的基本引數如下表：
窗函式法設計濾波器的步驟：
1）根據技術要求確定等待求值濾波器的單位取樣響應。
2）根據對過渡帶和阻帶衰減的要求，選擇窗函式的形式，並估計視窗長度N。
3）計算濾波器的單位取樣響應：
（3-16)
上式中，是前面所選擇好的窗函式。
4）檢驗技術指標是否滿足要求。根據下式計算：
(3-17)
如果不滿足要求，根據具體情況重複步驟(2)(3)(4)步，直到滿足要求為止。
本篇以一個FIR濾波器的設計為例說明如何使用MATLAB設計數字濾波器。
設計例項：用窗函式法設計線性相位FIR低通數字濾波器，要求通帶截止頻率Wp=0.4*，阻帶截止頻率Ws=0.5*, 通帶衰減不大於3db，阻帶衰減不小於40db。
1、程式如下：
Wp=0.4*pi;Ws=0.5*pi;Wdel=Ws-Wp;N=ceil(8*pi/Wdel);Wn=(0.4+0.5)*pi/2;window=hanning(N+1);b=fir1(N,Wn/pi,window);freqz(b,1,512)
程式執行後所得幅頻和相頻如下圖所示：
2、用Simulink模擬對窗函式設計的低通濾波器進行濾波：
FIR低通濾波器的Simulink模擬模組圖如下所示：
輸入訊號：s(t)=sin(30πt)+sin(56πt)+sin(70πt)
數字濾波器引數設定：已知濾波器的階數n=80。
首先在Filter Type中選擇lowpass；
在Design Method選項中選擇FIR Window，接著在Window選項中選取hann；指定Filter Order項中的Specify order為80；
取樣頻率Fs=100Hz,截止頻率Fc=22.5HZ。
設定完以後點選視窗下方的Design Filter，在視窗上方就會看到所設計濾波器的幅頻響應和相頻響應曲線。
設計介面如下圖所示：
訊號源引數設定：三個訊號源的頻率分別是15Hz、28Hz、35Hz ：
執行模擬，得波形如下圖所示：
15HZ訊號波形如下：
28HZ訊號波形如下：
35HZ訊號波形如下：
濾波前的疊加波形如下：
濾波後波形如下：
由上圖可知，15Hz、28Hz、35Hz疊加後通過截止頻率為22.5HZ的FIR低通濾波器，輸出與15Hz訊號波形相似，其它的高頻被阻隔。
二、頻率取樣法設計FIR濾波器
1、對理想濾波器的系統函式進行頻率取樣以得到系統的理想頻響的等間隔取樣值H(k)。
H(k)實際上是所要求的濾波器的單位取樣響應(h(n))的離散傅立葉變換（DFT），如下關係式：
(3-18)
(3-19)
為了減小H（k）的通帶邊緣由於抽樣點的變化而引起的起伏振盪，可以增加過渡點，加寬過渡帶以減小通帶的起伏。
每一個抽樣值產生一個與成正比，並位移的頻率響應，而H（k）與內插函式的線性組合就是FIR濾波器的頻率響應，增加一點過渡可以使阻帶衰減提高到-44~54dB，二點過渡衰減-65~75，三點過渡衰減-85~95dB。
如果不能使過渡帶太寬，同時要求增大阻帶衰減，可以增加取樣點數N，但這樣會增加計算量、延時和誤差。
頻率取樣型FIR濾波器設計步驟：
(1)給定理想濾波器頻率響應。
(2)根據過渡頻寬和阻帶衰減確定過渡點數和h(n)的長度N。
(3-20)
(3)由IFFT計算IDFT得到：
(3-21)
設計例項：
頻率取樣法設計一個帶通濾波器，
滿足：低阻帶邊緣：w1s=0.2*；低通帶邊緣：w1p=0.35*；高通帶邊緣：w2p=0.65*；高阻帶邊緣：w2s=0.8*。設計過渡帶中的頻率樣本值為t1和t2，取t1=0.109021，t2=0.59417456。
設計程式如下：
M=40;al=(M-1)/2;l=0:M-1;t1=0.109021;t2=0.59417456;Hrs=[zeros(1,5),t1,t2,ones(1,7),t2,t1,zeros(1,9),t1,t2,ones(1,7),t2,t1,zeros(1,4)];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));freqz(h,1,512,1000)
實驗得幅頻相頻特性如下圖所示：
三、最優化方法設計FIR數字濾波器
MATLAB訊號處理工具箱提供了通用的工具函式firls和remez，採用不同的優化方法可以設計最優濾波器和多頻帶濾波器。
firls是fir1和fir2函式的擴充套件，其基本設計準則是利用最小二乘法使期望的頻率響應和實際的頻率響應之間的整體誤差最小。
Remez函式則利用Parks-McClellan演算法設計線性相位FIR濾波器。
該演算法能使期望頻率響應和實際頻率響應之間的最大誤差最小。Firls函式設計的濾波器在頻率響應上出現等波紋，因此也叫等波紋濾波器。
函式firls、remez呼叫的語法規則相同。
1） firls。
其呼叫格式為：
b=firls(n,f,a)
b=firls(n,f,a,w)
b=firls(n,f,a,’ftype’)
b=firls(n,f,a,w,’ftype’)
該函式用於設計n階FIR濾波器，其幅頻特性由f和a向量確定，f是頻率點向量，其範圍為0~1，是遞增向量；a是指定頻率點的幅度響應，與f必須同長度；b為函式返回的濾波器係數，長度為n+1,且具有偶對稱關係：b(k)=-b(n+2-k), k=1,2,…,n+1；
b=firls(n,f,a,w)則使用權係數w給誤差加權。
b=firls(n,f,a,w,’ftype’)呼叫形式中，引數’ftype’用於指定所設計的濾波器型別，ftype=Hilbert，為奇對稱的線性相位濾波器。
返回的濾波器係數滿足b(k)=-b(n+2-k), k=1,2,…,n+1；ftype=differentiator，則採用特殊加權技術，生成奇對稱的線性相位濾波器，使低頻段誤差大大小於高頻段誤差。
2） 函式remez的基本呼叫格式為：b=remez(n,f,a)。
其中各項意義同函式firls所述。
設計例項：
分別使用函式firls和remez函式設計20階帶通濾波器，阻帶頻率為ws1=0.3*，ws2=0.7*，通帶頻率為wp1=0.4*,wp2=0.6*，並畫出幅頻特性曲線。
1、程式如下所示：
f=[0 0.3 0.4 0.6 0.7 1];m=[0 0 1 1 0 0];n=20;b=firls(n,f,m); %firls函式設計FIR數字濾波器[h,w]=freqz(b);bb=remez(n,f,m); %remez函式設計FIR數字濾波器[hh,w]=freqz(bb);axes('position',[0.2 0.2 0.5 0.5]);plot(w/pi,abs(h),'b:',w/pi,abs(hh),'r-');xlabel('frequency');ylabel('magnitude');grid on;
實驗得幅頻相頻特性如下圖所示：
2、用Simulink模擬對用firls函式設計的帶通濾波器進行濾波：
模擬模組如下：
輸入訊號：s(t)=sin(20πt)+sin(50πt)+sin(90πt)
數字濾波器引數設定：已知濾波器的階數n=20。
首先在Filter Type中選擇bandpass；
在Design Method選項中選擇FIR least-squares；
指定Filter Order項中的Specify order為20；
取樣頻率Fs=100Hz,Fstop1=15HZ，Fpass1=20HZ，Fpass2=30HZ，Fstop2=35HZ。
設定完以後點選視窗下方的Design Filter，在視窗上方就會看到所設計濾波器的幅頻響應曲線。
設計介面如下圖所示：
訊號源引數設定：三個訊號源的頻率分別是10Hz、25Hz、45Hz ：
執行模擬，得波形如下圖所示：
（1）10HZ訊號波形如下所示：
（2） 25HZ訊號波形如下所示：
（3）45HZ訊號波形如下所示：
（4）濾波前波形如下所示：
（5）濾波後波形如下所示：
由以上波形可知，10Hz、25Hz、45Hz疊加後通過通頻帶為20Hz~30Hz的FIR帶通濾波器，輸出與25Hz訊號相似，低頻和高頻部分被阻隔。
第三篇就到這裡，本篇在數字濾波器的設計過程中，採用的設計方法是基於MATLAB的數字濾波器的設計。設計過程中充分發揮了MATLAB語言程式設計效率高的特點以及Simulink模擬能力強的特點，完成各種計算與圖形繪製工作，使得數字濾波器的設計變得簡便快捷。通過FIR數字濾波器的設計例項，說明如何利用MATLAB來完成數字濾波器的設計，設計結果的各項效能指標均達到指定要求，Simulink模擬驗證了設計方法的正確性。本篇所述設計過程簡便易行，具有很好的實際運用意義。
END
後續會持續更新，帶來Vivado、 ISE、Quartus II 、candence等安裝相關設計教程，學習資源、專案資源、好文推薦等，希望大俠持續關注。大俠們，江湖偌大，繼續闖蕩，願一切安好，有緣再見！
往期推薦