FFT濾波誤差分析

楊 君， 呂鏡清

（①裝備指揮技術學院航天測控工程研究中心，北京 101416；②信息綜合控制國家重點實驗室，四川 成都 610036；③中國西南電子設備研究所，四川 成都 610036）

0 引言

FFT作為DFT的快速算法一直以來都是頻譜分析的重要工具，由于DFT運算與復調制濾波器組運算的相似性[1-2]，FFT的濾波功能也逐漸引起了重視，作為復調制濾波器組的快速算法近年來FFT在圖像傳輸、信道化處理[3]、信號盲識別[4]等領域得到了應用。但是在實際工程設計中FFT的濾波輸出常常與預期存在差異，主要表現在時域波形畸變、起始部分數據丟失等，這些差異可能造成后續數據處理產生嚴重誤差，而對于FFT運算與濾波運算之間的具體差異相關文獻并沒有進一步論述，工程中往往只能通過反復調整其低通原型濾波器來接近預期效果，使得設計效率大大降低且濾波性能仍未得到根本保證，為此本文對FFT的濾波原理進行了推證，找出了存在濾波誤差的原因，并給出了消除誤差的方法。

1 FFT的濾波原理

根據卷積和的定義，一個時域離散信號 x(n)(長度為 L)通過一個FIR線性相位濾波器h(n)(階數為M，長度為M+1)的計算過程如圖1所示。

圖1 濾波過程(y(n)為濾波輸出)

如圖1所示，第n時刻的濾波輸出實際是x(n)及其之前的M個點(將這M+1個點定義為xn)與反向后的h(n)相乘再相加的結果，因此濾波結果完全可以等價表示為如下這種滑動相關的形式：

其中xn(i)=x(n-M+i)。在式(1)基礎上定義一個新的函數：

將求和式展開后比較式(1)和式(2)，有如下結論：

因此 y’(n)是 y(n+1)的近似，而近似程度取決于 h(0)和xn+1(M)的取值。也就是說，當h(0) xn+1(M)相對很小時，信號x(n)通過濾波器h(n)第n+1時刻的濾波輸出近似等于進而有：

以M階梳狀濾波器[5](即系數全為1的濾波器)為低通原型，構造復調制濾波器組如下：

則根據式(4)，信號x(n)通過該濾波器組時，第k個子信道的濾波輸出為：

顯然式(6)等號右邊正是對xn-1做FFT的表達式。而前面定義xn(i)=x(n-M+i)，因此FFT的濾波原理可以表述為：將第n-1時刻的輸入數據及其之前的M個數據作FFT其結果近似等于以M階梳狀濾波器為低通原型的M信道(通道)復調制濾波器組第n時刻的濾波輸出，而濾波誤差為h(0)xn(M)即x(n)。

不難證明用 FFT實現低通原型為其它濾波器(記為h0(n)，且階數為M)的復調制濾波器組的表達式為：

即先加窗再做FFT，所加的窗函數為反向移位后的低通原型h0(n)，同時第n時刻濾波輸出的誤差為h0(0) xn(M)(即h0(0)x(n))。

2 頻率響應誤差

由于第n時刻FFT的濾波輸出中缺少了h0(0) x(n)項(根據式(2))，相當于FFT等效的實現了以(n)為低通原型的復調制濾波器組，其中

顯然 h0(n) 的頻率響應與(n)存在差異，下面進行詳細論述。

2.1 誤差分析

通常采用I類線性相位FIR濾波器（長度為奇數，系數偶對稱）作為低通原型，因此h0(n)的頻率響應為[2]：

在式(10)中令：

相頻響應為：

其中()ωσω滿足：

另一方面考查h0(n)的幅頻、相頻響應，根據式(9)有：

在通帶內 A (ω)遠大于 h0(M ),上式約等于 - 2A( ω )h0(M )·cos(ωM/2)，此時如果的通帶波紋與cos(ωM/2)有著相似的變化規律(比如用等波紋法設計出的h0(n) )，則會造成在取最大值時比大，在取最小值時比小，換句話說(n)增大了h0(n)的通帶波紋，

其增大的幅度由 A (ω)在峰值頻點上的取值以及 h0(M )的大小決定。另一方面在阻帶內 A (ω)與 h0(M )相當，甚至更小，因此式(13)在阻帶內取值會很小，比如用等波紋法設計出128階的h0(n)其阻帶衰減為-140 dB，而此時(n)的阻帶衰減與其最大相差約6 dB，因此(n)相對于(n)在阻帶上的變化是不明顯的。

① 相頻響應非線性；

② 可能造成通帶波紋增大。

這里舉例驗證此結論。取歸一化通帶截止頻率為 1/8，過渡帶寬為1/80，在 Matlab的Filter Design & Analysis Tool中利用Equiripple法設計出不同階數的低通濾波h0(n)，其中最小階數為16最大為1 024，階數增加的步長取15，然后根據式(8)構造出相應的(n)，對(n)通帶內的頻率響應進行誤差分析。各頻點的頻率響應由DTFT(離散時間傅里葉變換)求得，歸一化頻率分辨率取10-5，各頻點導數用斜率代替。仿真結果如圖2所示。

圖2 ()n通帶內的頻率響應誤差

仿真結果表明頻率響應誤差與 h0(M ) 存在密切關系。在階數較小時 h0(M ) 較大，(n)通帶內的頻率響應誤差較大。隨著階數的增加h0(n)本身性能提高導致 h0(M ) 逐漸減小，從而使得(n)通帶內的頻率響應誤差逐漸減小。

2.2 誤差影響

頻率響應誤差導致'

0()hn的濾波性能較 h0(n)有所下降，主要表現在：

①導致信號時域波形畸變；

②影響和信道濾波器的通帶平坦度，降低和信道濾波器的濾波性能。

下面著重說明第二點。在圖像傳輸以及信道化處理中經常需要將相鄰的子信道濾波器相加來增加濾波帶寬[6-7](將相加后得到濾波器稱之為和信道濾波器)，為了保證濾波質量，要求和信道濾波器通帶平坦。假設以h0(n)為低通原型的復調制濾波器組滿足和信道通帶平坦的要求，當用FFT實現此濾波器組時其低通原型變為(n)，其和信道濾波器的平坦度必然受到影響。下面舉例說明。

取歸一化通帶截止頻率為0.10375，過渡帶寬為0.0425，在 Matlab的 Filter Design & Analysis Tool中利用Equiripple(等波紋法)設計出128階低通濾波h0(n)，根據式(8)得到(n)，分別以h0(n)和(n)為低通原型構造8信道濾波器組，將第 2,3,4子信道濾波器合并，得到和信道濾波器通帶平坦度對比如下。

圖3 和信道通帶對比

2.3 誤差消除

根據式(11c)和式(12)可知，若 h0(M )為零則可以完全消除(n)的相頻特性非線性以及幅頻特性誤差，因此在設計低通原型h0(n)時應注意使首末系數為零或盡量接近于零。利用窗函數法就可以很容易的設計出首末系數為零的低通濾波器。窗函數法設計出的低通濾波器系數具有如下特點：

其中cω為通帶截至頻率，M為濾波器階數，W(n)為選擇的窗函數(如 Gaussion窗、Kaiser窗、Chebyshev窗等)。首末系數為：

顯然只要ωcM/2為π的整數倍就可實現首末系數為零。當需要實現2D個信道的濾波器組時，低通原型的通帶截至頻率應設置為π/2D，此時ωcM/2=πM/4D ，因此只要階數M為4D的整數倍就可實現首末系數為零。在設計濾波器時這一點是非常容易做到的。

事實上更為簡單的方法是在已有濾波器的首末各補一個零,這樣既不影響濾波器性能又消除了FFT時的頻率響應誤差。

3 相位超前現象

當h0(M)=0時可以避免頻率響應誤差，但是仍會產生相位超前現象。

考查式(10)，當h0(M)=0時'()hn的幅頻響應為：

如果在用 FFT進行濾波運算時仍然認為群延遲為 M/2個采樣點就會產生一個采樣周期的誤差。這個誤差隨采樣周期的增大而增大，在一些場合中會影響數據處理精度，例如統一測控系統中的側音測距，其利用收側音相對于發側音的傳輸延遲來推算出目標距離，公式可以簡單的表示為R=Δt×C/2，其中R為目標距離，Δt為傳輸延時，C為光速3×108m/s，設中頻欠采樣頻率為6.5 MHz，則算出的距離就會增加約23 m的誤差。

另一種超前現象主要在事后數據處理時發生。在事后處理中由于數據都已獲得因此N點FFT直接從前N個數據開始，等效于圖1中從n=M,(N=M+1)開始向后滑動，濾波輸出產生了M點的超前，這樣起始M個點的濾波結果就會丟失，減去(M-2)/2個點的群延時，則最終造成起始的(M+2)/2個點數據丟失。而在實時處理時采樣數據逐個進入數據緩沖區，而數據未到達的緩沖區數據位為零，等效于數據處理從圖 1中n=0開始，因此實時處理不會出現此種超前現象。

兩種超前現象都會造成濾波誤差。對于第一種超前現象只需要在扣除濾波延時時用(M-2)/2代替M/2即可消除誤差；對于第二種超前現象，在數據前補M個零即可等效于圖 1中從n=0開始濾波，從而避免數據丟失。

4 結語

本文詳細分析了 FFT運算與復調制濾波器組運算之間的差異，推導出了頻率響應的誤差項，并對該誤差造成的濾波性能下降進行了仿真分析，針對誤差產生原因提出了首末置零的簡單方法以完全消除誤差；另一方面文章分析指出了在實時及事后數據處理中 FFT相對于濾波器組存在的固有的相位超前現象，給出具體超前量和有效解決方法。文章為精確控制FFT的濾波性能及進一步提高其濾波質量提供了有效參考。

[1] Harris F J.Time Domain Signal Processing with the DFT，Ch.8 of Elliot,D.F.，Editor，Handbook of Digital Signal Engineering Applications[M]．San Diego．CA:Academic Press,Inc,1987.

[2] 李素芝,萬建偉.時域離散信號處理[M]．長沙：國防科技大學出版社，2000:83-133．

[3] 王甲峰,葛曉碕.無盲區數字信道化實現方法[J].通信技術,2009,42(03):8-10.

[4] 楊偉超,張忠,丁群.低信噪比數字通信信號識別算法研究[J].通信技術,2009,42(01):68-71.

[5] Proakis J G, Manolakis D G.數字信號處理[M].方艷梅，劉永清譯．北京:電子工業出版社,2007:367-372.

[6] Tsui J．寬帶數字接收機[M]．北京：電子工業出版社,2002：242-260．

[7] 陶然,張惠云,王越.多抽樣率數字信號處理理論及其應用[M].北京：清華大學出版社,2008:122-130.