1、第4章 快速傅里葉變換 4.1 引言 4.2 直接計算DFT的問題及改進的途徑 4.3 按時間抽選的基2-FFT算法（Cooley-Turkey）4.4 按頻率抽取的基2-FFT算法（Sande-Turkey）4.5 離散傅立葉反變換的快速計算方法 4.6 線性調頻z變換算法 4.7 FFT實例分析 4.1 引 言 注意：快速傅里葉變換（FFT）不是一種新的變換，是DFT的一種快速算法。 1. DFT在時域和頻域都是的，可以采用計算機進行運算；2. 直接計算DFT的運算量很大，即使采用計算機運算，也不能解決實時性問題，影響其實際應用；3. 1965年首次提出了DFT運算的一種快速算法，并發展和

2、形成了一套高速有效的運算方法， 統稱為快速傅里葉變換（FFT）的算法。4.2 直接計算DFT的問題及改進的途徑 4.2.1 DFT的運算量k=0, 1, , N-1 （4.1） 反變換（IDFT）為 n=0, 1, , N-1 （4.2） 設x(n)為N點有限長序列，其DFT為 二者的差別: 1) WN的指數符號不同， 2) 差一個常數乘因子1/N，IDFT與DFT具有相同的運算量， 可以只討論DFT的運算量。 x (n)和WNnk都是復數，X(k)也是復數，完成整個DFT運算共需要：N2次復數乘法 N(N-1)次復數加法。 X (k)共有N個值 （k=0，1，N-1）每計算一個X(k)值，需

3、要：N次復數乘法 N-1次復數加法。（4.3） 一次復數乘法：4次實數乘法 2次實數加法 一次復數加法：2次實數加法 復數運算實際上是由實數運算來完成的，可將DFT運算式寫成 整個DFT運算共需要：4N 2次實數乘法 2N (2N-1)次實數加法。 一次復乘：4次實數乘法 2次實數加法 一次復加：2次實數加法 一個X(k)值：N次復數乘法 N-1次復數加法每計算一個X(k)需要：4N 次 實數乘法 2N+2(N-1)=2(2N-1)次 實數加法上述統計與實際需要的運算次數稍有出入，某些WNnk可能是1或 j，如W0N=1，WN= -1，WNN/4=-j等 但是為了便于和其他運算方法作比較， 一

4、般都不考慮這些特殊情況，而是把WNnk都看成復數，當N很大時，這種特例的影響很小。因此，直接計算DFT，乘法次數和加法次數都和N2成正比， 當N很大時，運算量是很可觀的。 例 根據式（4.1），對一幅NN點的二維圖像進行DFT變 換，如用每秒可做10萬次復數乘法的計算機，當N=1024時，問需要多少時間（不考慮加法運算時間）？解: 直接計算DFT所需復乘次數為（N 2)21012次， 用每秒可做10萬次復數乘法的計算機，需要近3000小時。 對實時性很強的信號處理，改進方法： 1）提高計算速度，（這樣，對計算速度要求太高了）； 2) 改進DFT的計算方法，以大大減少運算次數。 4.2.2 減少

5、運算量的途徑（2）WNnk的周期性 DFT運算，利用系數Wnk的固有特性，可減少運算量： 能否減少運算量，從而縮短計算時間呢？（1） WNnk的對稱性 （3） WNnk的可約性 則有 另外 FFT算法基本上可以分成兩大類，按時間抽選法 (DIT ，Decimation-In-Time）按頻率抽選法 (DIF， Decimation-In-Frequency） 利用這些特性，使DFT運算中可以合并有些項；利用WNnk的對稱性、周期性、可約性，使DFT分解為更少點數的DFT運算。前面提到，DFT的運算量與N2 成正比，N 越小DFT的運算量越小。快速傅里葉變換算法的基本思想4.3 按時間抽取（DI

6、T）的基 -2 FFT算法 4.3.1 算法原理設序列 x (n) 長度為N，且滿足N=2L，L為正整數步驟：將 x(n)按 n 先奇后偶分解為兩個N/2點的子序列 基-2 FFT算法則可將DFT化為 （4.5） （4.4） 利用WNnk的可約性 X1(k) 與 X2(k) 分別是 x1(r) 及 x2(r) 的N/2點DFT： （4.6） （4.7） （4.5） 只是X(k)的前一半的結果第4章 快速傅里葉變換X1(k)，X2(k)只有N/2個點，即k=0, 1, , N/2-1。X(k)有N個點，即k=0, 1, , N-1， 要用X1(k)，X2(k) ，Wk來表達全部的X(k)值這樣可

7、得到 （4.8） 同理可得 （4.9） 式（3-8）、式（3-9）說明了后半部分k值（N/2kN-1）所對應的X1(k)，2(k)分別等于前半部分k值（0kN/2-1）所對應的X1(k)，X2(k)。 再考慮到WkN 的以下性質: 這樣，把式（3-8）、式（3-9）、式（3-10）代入式（3-5），就可將X(k)表達為前后兩部分： （4.10） （4.11） （4.12） 圖 4.1 時間抽選法蝶形運算流圖符號 圖 4.2 按時間抽選，將一個N點DFT分解為兩個N/2點DFT（N=8） N點DFT，一次分解后次復數乘法次復數加法將N點DFT，進行一次分解后，運算工作量節省了近一半N點DFT，不

8、進行分解次復數加法次復數乘法由于N=2L，N/2仍是偶數可以進一步把每個N/2點子序列分解為兩個N/4點的子序列（4.13） ， 且 式中 （4.14） （4.15） 圖4.3給出N=8時，將一個N/2點DFT分解成兩個N/4點DFT， 由這兩個N/4點DFT組合成一個N/2點DFT的流圖。 圖4.3 由兩個N/4點DFT組合成一個N/2點DFT X2(k)也可進行同樣的分解： 式中 （4.16） （4.17） 圖4.4 按時間抽選，將一個N點DFT分解為四個N/4點DFT （N=8） 根據上面同樣的分析知道，利用四個N/4點的DFT及兩級蝶形組合運算來計算N點DFT，比只用一次分解蝶形組合方

9、式的計算量又減少了大約一半。 式中， , 故上式不需要乘法。 類似地, 可求出X4(k)，X5(k)，X6(k)。 這種方法的每一步分解，都是按輸入序列在時間上的次序是屬于偶數還是屬于奇數來分解為兩個更短的子序列，所以稱為“按時間抽取法”。 圖 3-5 N=8 按時間抽取的FFT運算流圖4.3.2 按時間抽取的FFT算法與直接計算DFT運算量的比較 由按時間抽取法FFT的流圖可見，當N=2M時，共有M級蝶形， 每級都由N/2個蝶形運算組成，每個蝶形需要一次復乘、 二次復加，因而每級運算都需N/2次復乘和N次復加，這樣M級運算總共需要 復乘數 復加數 式中，數學符號lb=log2。 實際計算量與

10、這個數字稍有不同，因為W0N=1，NN/2=-1，NN/2=-j，這幾個系數都不用乘法運算，但是這些情況在直接計算DFT中也是存在的。此外，當N較大時，這些特例相對而言就很少。所以，為了統一作比較起見，下面都不考慮這些特例。 由于計算機上乘法運算所需的時間比加法運算所需的時間多得多，故以乘法為例，直接DFT復數乘法次數是N2，FFT復數乘法次數是(N/2) lbN。 直接計算DFT與FFT算法的計算量之比為 當N=2048時，這一比值為372.4，即直接計算DFT的運算量是FFT運算量的372.4倍。當點數N越大時，FFT的優點更為明顯。 （4.20） 例4-2 用FFT算法處理一幅NN點的二

11、維圖像，如用每秒可做10萬次復數乘法的計算機，當N=1024時，問需要多少時間（不考慮加法運算時間）？ 解 當N=1024點時，FFT算法處理一幅二維圖像所需復數乘法約為次，僅為直接計算DFT所需時間的10萬分之一。 即原需要3000小時，現在只需要2 分鐘。 4.3.3 按時間抽取的FFT算法的特點及DIT-FFT程序框圖 1. 原位運算（同址運算）這種運算是很有規律，其每級（每列）計算都是由N/2 個蝶形運算構成的，每一個蝶形結構完成下述基本迭代運算： （4.21a） （4.21b） 式中，m表示第m列迭代，k, j為數據所在行數。式（4.21）的蝶形運算如圖4-6所示，由一次復乘和兩次復

12、加（減）組成。 圖 4-6 蝶形運算單元 當數據輸入到存儲器后，每一級運算的結果仍然存儲在這同一組存儲器中，直到最后輸出，中間無需其他的存儲器。每一級運算均可在原位進行，這種原位運算的結構可以節省存儲單元，降低設備成本。2. 倒位序規律按時間抽選FFT 輸出X(k) 按自然順序排列在存儲單元 輸入 x(n) 不是按自然順序排列是由于x(n)按標號n先偶后奇不斷分組造成的倒位序的規律：例如圖4-7 倒位序的形成 表4-1 N=8時的自然順序二進制數和相應的倒位序二進制數 自然順序（I） 二進制數 倒位序二進制數 倒位序（J） 01234567000001010011100101110111000

13、10001011000110101111104261537圖 4-8 N=8 倒位序的變址處理 3. 蝶形運算兩節點的“距離” 以8點FFT為例，其輸入是倒位序的，輸出是自然順序的。 其第一級（第一列）每個蝶形的兩節點間“距離”為1， 第二級每個蝶形的兩節點“距離”為2， 第三級每個蝶形的兩節點“距離”為4。 由此類推得，對N=2M點FFT，當輸入為倒位序， 輸出為正常順序時，其第m級運算，每個蝶形的兩節點“距離”為2m-1。 4. WNr的確定 由于對第m級運算，一個DFT蝶形運算的兩節點“距離”為2m-1， 因而（4.22a） （4.22b） 為了完成上兩式運算，還必須知道系數Nr的變換規

14、律。仔細觀察圖4-5的流圖可以發現r的變換規律是： 把式（4.22）中蝶形運算兩節點中的第一個節點標號值， 即k值，表示成M位（注意N=2M）二進制數； 把此二進制數乘上2M-m，即將此M位二進制數左移M-m位（注意m是第m級運算），把右邊空出的位置補零， 此數即為所求r的二進制數。 圖中，WNr因子最后一列有N/2種，順序為WN0， WN1， ， 其余可類推。 5. DIT-FFT程序框圖 圖 4-9 DIT-FFT運算程序框圖 圖 4-10 倒序程序框圖 4.4 按頻率抽?。―IF）的基 -2 FFT算法 4.4.1 算法原理 仍設序列點數為N=2M，M為正整數。在把輸出X(k)按k的奇偶

15、分組之前，先把輸入序列按前一半、后一半分開（不是按偶數、 奇數分開）， 把N點DFT寫成兩部分。 k=0, 1, , N-1 式中，用的是 ，而不是 ，因而這并不是N/2點DFT。 由于 ，故，可得 k=0, 1, , N-1 （4.23） 當k為偶數時，(-1)k=1；k為奇數時，(-1)k=-1。因此，按k的奇偶可將X(k)分為兩部分: （4.24） （4.25） 式（4.24）為前一半輸入與后一半輸入之和的N/2點DFT， 式（4.25）為前一半輸入與后一半輸入之差再與WNn之積的N/2點DFT。 令:（4.27） 圖 4-12 頻率抽取法蝶形運算單元 這樣，我們就把一個N點DFT按k的

16、奇偶分解為兩個N/2點的DFT了。 N=8時，上述分解過程示于圖4-13。 與時間抽取法的推導過程一樣，由于N=2M，N/2仍是一個偶數，因而可以將每個N/2點DFT的輸出再分解為偶數組與奇數組，這就將N/2點DFT進一步分解為兩個N/4 點DFT。 這兩個N/4點DFT的輸入也是先將N/2點DFT的輸入上下對半分開后通過蝶形運算而形成的，圖4-14示出了這一步分解的過程。 圖4-13 按頻率抽取的第一次分解 圖 4-14 按頻率抽取的第二次分解 這樣的分解可以一直進行到第M次（N=2M），第M次實際上是做兩點DFT，它只有加減運算。 然而，為了有統一的運算結構，仍然用一個系數為WN0的蝶形運

17、算來表示， 這N/2個兩點DFT的N個輸出就是x(n)的N點DFT的結果X(k)。 圖4-15表示一個N=8 完整的按頻率抽取的基-2 FFT運算結構。 圖 4-15 按頻率抽取的FFT（N=8）信號流圖 4.4.2 按頻率抽取法的運算特點 按頻率抽取法的運算特點與時間抽取法基本相同。從圖4-15可以看出，它也是通過(N/2)M個蝶形運算完成的。每一個蝶形結構完成下述基本迭代運算： 式中，m表示m列迭代，k,j為整數所在行數，此式的蝶形運算如圖4-16所示，也需要一次復乘和兩次復加。 圖4-16 頻率抽取法蝶形運算單元 1按時間抽選：輸入是倒位序， 輸出是自然順序； 按頻率抽選：輸入是自然順序

18、， 輸出是倒位序。DIT法與DIF法比較：二、相同點 運算量相同一、不同點2按時間抽選：復數乘法出現在加、減運算之前； 按頻率抽選：復數乘法出現在減法運算之后。4.5 離散傅里葉反變換的快速計算方法 1）DFT中的 IDFT中的 2）IDFT 中有常數區別：按時間抽選的 FFT按頻率抽選的 IFFT按頻率抽選的 FFT按時間抽選的 IFFT稍作改動FFT程序，計算IFFT的方法按時間抽選的 IFFT 流圖（N =8） 完全不改變FFT程序，計算IFFT的方法 4.6* 線性調頻z變換算法 4.6.1 CZT變換算法原理（4.32） 為適應z可以沿z平面更一般的路徑取值，故沿z平面上的一段螺線作

19、等分角的采樣，z的這些采樣點zk為 zk=AW-k k=0, 1, , M-1 （4.33） M為所要分析的復頻率的點數，即采樣點的總數，不一定等于N； A和W都是任意復數，可表示為: 已知x(n)（0nN-1）是有限長序列，其z變換為 將式（ 4.34 ）與式（ 4.35 ）代入式（ 4.33 ）， 可得 （4.34）（4.35）（4.36）因此有 圖4.15 (a) CZT在平面的螺旋線抽樣 采樣點在Z平面上所沿的周線如圖4.15(a)所示。由以上討論和圖4.15(a)可以看出： （1）A0表示起始采樣點z0的矢量半徑長度，通常A01; 否則z0將處于單位圓|z|=1 的外部。 （2）0表

20、示起始采樣點z0的相角，它可以是正值或負值。 （3）0表示兩相鄰采樣點之間的角度差。0為正時，表示zk的路徑是逆時針旋轉的；0為負時，表示zk的路徑是順時針旋轉的。 （4）W0的大小表示螺線的伸展率。W01 時，隨著k的增加螺線內縮; W0M，則只需要求得0nM-1一段N點序列值即可; h(n)也可事先準備好，不必實時分析時計算，因此，可不用考慮其計算量。 同時，h(n)的L點FFT即H(r)也可預先計算好。 （3）計算G(r)，H(r), q(n), 需要二次L點FFT（或IFFT），共需要 次復乘。 （4） 計算Q(r)=G(r)H(r), 需要L次復乘。 （5）計算X(zk)= （0kM

21、-1），需要M次復乘。 綜上所述，CZT總的復數乘法次數為前面說過，直接計算式（4.37）的X(zk)需要NM次復數乘法?？梢钥闯觯擭, M都較大時（例如N, M都大于50時），CZT的FFT算法比直接算法的運算量要小得多。 4.7 FFT實例分析4.7.1 利用FFT分析時域連續信號頻譜 圖 4.17 時域連續信號離散傅里葉分析的處理步驟 一、 基本步驟 前置低通濾波器LPF（預濾波器）的引入，是為了消除或減少時域連續信號轉換成序列時可能出現的頻譜混疊的影響。 在實際工作中，時域離散信號x(n)的時寬是很長的， 甚至是無限長的（例如語音或音樂信號）。由于DFT的需要(實際應用FFT計算)，

22、必須把x(n)限制在一定的時間區間之內，即進行數據截斷。數據的截斷相當于加窗處理（窗函數見6.2節）。 因此， 在計算FFT之前，用一個時域有限的窗函數w(n)加到x(n)上是非常必要的。xc(t)通過A/D變換器轉換(忽略其幅度量化誤差)成采樣序列x(n)，其頻譜用X(ej)表示，它是頻率的周期函數，即 (4.54) 式中，Xc(j)或 為xc(t)的頻譜。 在實際應用中，前置低通濾波器的阻帶不可能是無限衰減的， 故由Xc(j)周期延拓得到的X(ej)有非零重疊，即出現頻譜混疊現象。 由于進行FFT的需要，必須對序列x(n)進行加窗處理，即v(n)=x(n)w(n)，加窗對頻域的影響，用卷積

23、表示如下: 最后是進行FFT運算。 加窗后的DFT是 0kN-1 (4.55) 式中，假設窗函數長度L小于或等于DFT長度N，為進行FFT運算，這里選擇N為2的整數冪次即N=2m。 有限長序列v(n)=x(n)w(n)的DFT相當于v(n)傅里葉變換的等間隔采樣。 V(k)便是sc(t)的離散頻率函數。因為DFT對應的數字域頻率間隔為=2/N，且模擬頻率和數字頻率間的關系為=， 其中=2f。所以，離散的頻率函數第k點對應的模擬頻率為: (4.57) (4.58) (4.56) 由式（4.58）很明顯可看出，數字域頻率間隔=2/N對應的模擬域譜線間距為 (4.59) 譜線間距，又稱頻譜分辨率（單

24、位：Hz）。所謂頻譜分辨率是指可分辨兩頻率的最小間距。它的意思是，如設某頻譜分析的F=5Hz，那么信號中頻率相差小于5Hz的兩個頻率分量在此頻譜圖中就分辨不出來。 長度N=16 的時間信號v(n)=(1.1)nR16(n)的圖形如圖 4.18(a)所示, 其16點的DFT V(k)的示例如圖4.18 (b)所示。 其中T為采樣時間間隔（單位：s）；fs為采樣頻率（單位：Hz）；tp為截取連續時間信號的樣本長度（又稱記錄長度，單位：s）；F為譜線間距，又稱頻譜分辨率（單位：Hz）。 注意：V(k)示例圖給出的頻率間距F及N個頻率點之間的頻率fs為對應的模擬域頻率(單位： Hz)。 圖 4.18

25、時間信號v(n)和V(k)的示意圖 由圖可知： (4.60) (4.61) 在實際應用中, 要根據信號最高頻率fh和頻譜分辨率F的要求， 來確定T、tp和N的大小。 (4.62) （1）首先，由采樣定理，為保證采樣信號不失真，fs2fh(fh為信號頻率的最高頻率分量，也就是前置低通濾波器阻帶的截止頻率)， 即應使采樣周期T滿足 （2） 由頻譜分辨率F和T確定N, (4.63) 為了使用FFT運算，這里選擇N為2的冪次即N=2m，由式(4.61)可知，N大，分辨率好，但會增加樣本記錄時間tp。 （3） 最后由N, T確定最小記錄長度，tp=NT。 例 4-1 有一頻譜分析用的FFT處理器，其采樣

26、點數必須是2的整數冪， 假定沒有采用任何特殊的數據處理措施，已給條件為: 頻率分辨率10 Hz;信號最高頻率4kHz。解 (1) 由分辨率的要求確定最小長度tp 所以記錄長度為 試確定以下參量： (1) 最小記錄長度tp； (2) 最大采樣間隔T（即最小采樣頻率）； (3) 在一個記錄中的最少點數N。 (2) 從信號的最高頻率確定最大可能的采樣間隔T（即最小采樣頻率fs=1/T）。 按采樣定理 即 (3) 最小記錄點數N應滿足 取 如果我們事先不知道信號的最高頻率，可以根據信號的時域波形圖來估計它。例如， 某信號的波形如圖 4.19 所示。 先找出相鄰的波峰與波谷之間的距離，如圖中t1，t2，

27、t3，t4。 然后，選出其中最小的一個如t4。這里, t4可能就是由信號的最高頻率分量形成的。 峰與谷之間的距離就是周期的一半。 因此，最高頻率為 知道 fh 后就能確定采樣頻率 圖 4.19 估算信號最高頻率fh 利用FFT對連續信號進行傅里葉分析時可能造成的誤差如下也就是采樣間隔T滿足 二、 可能出現的誤差 1. 頻譜混疊失真 在圖4.17 的基本步驟中，A/D變換前利用前置低通濾波器進行預濾波，使xc(t)頻譜中最高頻率分量不超過fh。假設A/D變換器的采樣頻率為fs，按照奈奎斯特采樣定理，為了不產生混疊， 必須滿足 一般應取 fs=(2.53.0)fh 如果不滿足fs2fh，就會產生頻

28、譜混疊失真。 對于FFT來說，頻率函數也要采樣，變成離散的序列，其采樣間隔為F(即頻率分辨率)。 由式（4.61）可得 （4.64） （4.65） 從以上和tp兩個公式來看，信號的最高頻率分量fh與頻率分辨率F存在矛盾關系，要想fh增加，則時域采樣間隔T就一定減小，而fs就增加，由式（4.63）可知，此時若是固定N，必然要增加F， 即分辨率下降。 反之，要提高分辨率（減小F），就要增加tp，當N給定時， 必然導致T的增加（fs減?。?。 要不產生混疊失真，則必然會減小高頻容量（信號的最高頻率分量）fh。 （4.66） 這個公式是未采用任何特殊數據處理（例如加窗處理）的情況下，為實現基本FFT算法

29、所必須滿足的最低條件。如果加窗處理，相當于時域相乘，則頻域周期卷積，必然加寬頻譜分量， 頻率分辨率就可能變差，為了保證頻率分辨率不變，則須增加數據長度tp。 要想兼顧高頻容量fh與頻率分辨率F，即一個性能提高而另一個性能不變（或也得以提高）的惟一辦法就是增加記錄長度的點數N， 即要滿足 2. 柵欄效應 利用FFT計算頻譜，只給出離散點k=2k/N或k=2k/(NT)上的頻譜采樣值，而不可能得到連續頻譜函數， 這就像通過一個“柵欄”觀看信號頻譜， 只能在離散點上看到信號頻譜， 稱之為“柵欄效應” 。這時，如果在兩個離散的譜線之間有一個特別大的頻譜分量，就無法檢測出來了。 減小柵欄效應的一個方法就

30、是要使頻域采樣更密，即增加頻域采樣點數N， 在不改變時域數據的情況下，必然是在數據末端添加一些零值點，使一個周期內的點數增加，但并不改變原有的記錄數據。 頻譜采樣為2k/N，N增加，必然使樣點間距更近（單位圓上樣點更多），譜線更密，譜線變密后原來看不到的譜分量就有可能看到了。 必須指出，補零以改變計算FFT的周期時，所用窗函數的寬度不能改變。 換句話說，必須按照數據記錄的原來的實際長度選擇窗函數， 而不能按照補了零值點后的長度來選擇窗函數。 補零不能提高頻率分辨率，這是因為數據的實際長度仍為補零前的數據長度。 3. 頻譜泄漏與譜間干擾 對信號進行FFT計算，首先必須使其變成有限時寬的信號， 這

31、就相當于信號在時域乘一個窗函數如矩形窗，窗內數據并不改變。 時域相乘即v(n)=x(n)w(n)， 加窗對頻域的影響，可用卷積公式卷積的結果，得到的頻譜V(ej)與原來的頻譜X(ej)不同，有失真。這種失真主要會造成頻譜“擴散”（拖尾、 變寬），即“頻譜泄漏”。 頻譜泄漏是由于截取有限長信號所造成的。 對具有單一譜線的正弦波來說，它必須是無限長的。即，如果輸入信號是無限長的，那么FFT就能計算出完全正確的單一線頻譜。 可是實際上，只能取有限長記錄樣本。 如果在該有限長記錄樣本中，正弦信號又不是整數個周期時，就會產生泄漏。 例 周期為N=16的余弦信號 x(n)=cos(6n）（2）截取的長度為

32、13，則其16點FFT的譜圖見 4.20 下圖。 由圖可見，頻譜不再是單一的譜線，能量散布到整個頻譜的各處。 這種能量散布到其他譜線位置的現象即為“頻譜泄漏”。 （1）截取一個周期長度 x1(n)=cos(6n/16)R16(n), 其16點FFT的譜圖見 4.20上圖，圖 4.20 余弦信號頻譜泄漏示例圖 泄漏將會導致頻譜的擴展， 從而使最高頻率有可能超過折疊頻率（fs/2），因此，頻譜泄漏也會造成混疊失真。 泄漏將會降低頻譜的分辨率。由于在主譜線兩邊形成很多旁瓣，引起不同頻率分量間的干擾(簡稱譜間干擾)， 特別是強信號譜的旁瓣可能湮沒弱信號的主譜線，或者把強信號譜的旁瓣誤認為是另一信號的譜

33、線，從而造成假信號，這樣就會使譜分析產生較大偏差。 在進行FFT運算時，要注意的兩點問題：第一，時域截斷是必然的，從而頻譜泄漏和譜間干擾也是不可避免的。為盡量減小泄漏和譜間干擾的影響，需增加窗的時域寬度(頻域主瓣變窄)，但這又導致運算量及存儲量的增加；第二，數據不要突然截斷，也就是不要加矩形窗，而是加各種緩變的窗（例如：三角形窗、升余弦窗、改進的升余弦窗等），使得窗譜的旁瓣能量更小，卷積后造成的泄漏減小，這個問題在FIR濾波器設計一章中會討論。 4.7.2 線性卷積和線性相關的FFT算法 一、線性卷積的FFT算法 以FIR濾波器為例，因為它的輸出等于有限長單位脈沖響應h(n)與有限長輸入信號x

34、(n)的離散線性卷積。 設x(n)為L點，h(n)為M點， 輸出y(n)為 y(n)也是有限長序列，其點數為L+M-1 點。 下面首先討論直接計算線性卷積的運算量。 由于每一個x(n)的輸入值都必須和全部的h(n)值相乘一次，因而總共需要LM次乘法，這就是直接計算的乘法次數，以md表示為 md=LM （4.67） 對于線性相位FIR濾波器，滿足 （4-68）（4-69）其加權系數約減少了一半，因而相乘次數大約可以減少一半，即 用FFT算法也就是用圓周卷積來代替這一線性卷積時，為了不產生混疊，其必要條件是使x(n)，h(n)都補零值點，補到至少N=M+L-1， 即: 0nL-1 LnN-1 0n

35、M-1 MnN-1 然后計算圓周卷積 N用FFT計算y(n)的步驟如下： 求H(k)=DFTh(n)，N點； 求X(k)=DFTx(n), N點； 求y(n)=IDFTY(k)，N點。 計算Y(k)=X(k)H(k); 步驟、都可用FFT來完成。 此時的工作量如下： 三次FFT運算共需要 次相乘，還有步驟的N次相乘，因此共需要相乘次數為 （4.70） 比較直接計算線性卷積(簡稱直接法)和FFT法計算線性卷積(簡稱FFT法)這兩種方法的乘法次數。當與點數差不多時，時，則 例如，設式（4.69）與式（4.70）的比值為Km，這樣可得下表： M=L81632641282565121024204840

36、96Km0.5720.9411.62.785.928.821629.2453.999.9當x(n)的點數很多時，即當LM，通常不允許等x(n)全部采集齊后再進行卷積; 否則，使輸出相對于輸入有較長的延時。此外，若N=L+M-1 太大，h(n)必須補很多個零值點，很不經濟，且FFT的計算時間也要很長。這時FFT法的優點就表現不出來了。因此需要采用分段卷積或稱分段過濾的辦法。即將x(n)分成點數和h(n)相仿的段，分別求出每段的卷積結果，然后用一定方式把它們合在一起，便得到總的輸出，其中每一段的卷積均采用FFT方法處理。有兩種分段卷積的辦法: 重疊相加法和重疊保留法。 2重疊相加法 設h(n)的點

37、數為M，信號x(n)為很長的序列。我們將x(n)分解為很多段，每段為L點，L選擇成和M的數量級相同，用xi(n)表示x(n)的第i段： iLn(i+1)L-1 其他n i=0, 1, （4.72） 則輸入序列可表示成 （4.73）這樣，x(n)和h(n)的線性卷積等于各xi(n)與h(n)的線性卷積之和，即 （4.74） 每一個xi(n)*h(n)都可用上面討論的快速卷積辦法來運算。 由于xi(n)*h(n)為L+M-1 點，故先對xi(n)及h(n)補零值點，補到N點。 為便于利用基-2 FFT算法，一般取N=2mL+M-1，然后作N點的圓周卷積： N 由于xi(n)為L點，而yi(n)為（L+M-1）點（設N=L+M-1）， 故相鄰兩段輸出序列必然有（M-1）個點發生重疊，即前一段的后（M-1）個點和后一段的前（M-1）個點相重疊，如圖4.21 所示。按照式（4.74），應該將重疊部分相加再和不重疊的部分共同組成輸出y(n)。 圖 4.21 重疊相加法圖形 圖 4.21重疊相加法圖形