數(shù)字信號處理 第4章 快速傅里葉變換(FFT)-

第4章快速傅里葉變換(FFT)4.1引言4.2基2FFT算法問題的提出解決問題的思路與方法基2時間抽取FFT算法基2頻率抽取FFT算法4.3進一步減少運算量的措施

4.1引言

DFT是信號分析與處理中的一種重要變換。因直接計算DFT的計算量與變換區(qū)間長度N的平方成正比，當N較大時，計算量太大，所以在快速傅里葉變換(簡稱FFT)出現(xiàn)以前，直接用DFT算法進行譜分析和信號的實時處理是不切實際的。直到1965年發(fā)現(xiàn)了DFT的一種快速算法以后，情況才發(fā)生了根本的變化。4.2基2FFT算法

問題的提出4點序列{2，3，3，2}DFT的計算復雜度復數(shù)加法N(N-1)復數(shù)乘法N

2如何提高DFT的運算效率?時間抽取FFT解決問題的思路1.將長序列DFT分解為短序列的DFT2.利用旋轉(zhuǎn)因子的周期性、對稱性、可約性。時間抽取FFT旋轉(zhuǎn)因子的性質(zhì)(1)周期性(2)對稱性(3)可約性時間抽取FFT解決問題的方法將時域序列逐次分解為一組子序列，利用旋轉(zhuǎn)因子的特性，由子序列的DFT來實現(xiàn)整個序列的DFT?；?時間抽取(Decimationintime)FFT算法基2頻率抽取(Decimationinfrequency)FFT算法時間抽取FFT基2時間抽取FFT算法基2時間抽取FFT算法推導基2時間抽取FFT算法流圖基2時間抽取FFT算法的計算復雜度基2時間抽取FFT算法流圖規(guī)律時間抽取FFT基2時間抽取FFT算法推導時間抽取FFT基2時間抽取FFT算法推導因此有:時間抽取FFT基2時間抽取FFT算法流圖N=2x[k]={x[0],x[1]}時間抽取FFT4點基2時間抽取FFT算法流圖x[0]x[2]x[1]x[3]X1[0]X1[1]X2[0]X2[1]2點DFT2點DFT-1-1-1-1X

[0]X

[1]X

[2]X

[3]時間抽取FFT4點基2時間抽取FFT算法流圖時間抽取FFT8點基2時間抽取FFT算法流圖4點DFT4點DFTx[0]x[2]x[4]x[6]x[1]x[3]x[5]x[7]X1[0]X1[1]X1[2]X1[3]X2[0]X2[1]X2[2]X2[3]X

[0]X

[1]X

[2]X

[3]X

[4]X

[5]X

[6]X

[7]-1-1-1-1時間抽取FFT4點DFT4點DFTx[0]x[2]x[4]x[6]x[1]x[3]x[5]x[7]X1[0]X1[1]X1[2]X1[3]X2[0]X2[1]X2[2]X2[3]X

[0]X

[1]X

[2]X

[3]X

[4]X

[5]X

[6]X

[7]-1-1-1-18點基2時間抽取FFT算法流圖時間抽取FFT第一級第二級第三級8點基2時間抽取FFT算法流圖時間抽取FFT算法的計算復雜度復乘次數(shù)復乘次數(shù)NN2時間抽取FFT基2時間抽取FFT算法流圖第一級第二級第三級時間抽取FFTFFT算法流圖旋轉(zhuǎn)因子規(guī)律第二級的蝶形系數(shù)為，蝶形節(jié)點的距離為2。第一級的蝶形系數(shù)均為，蝶形節(jié)點的距離為1。第三級的蝶形系數(shù)為，蝶形節(jié)點的距離為4。第M級的蝶形系數(shù)為，蝶形節(jié)點的距離為N/2。時間抽取FFT倒序k0k1k2x[k2k1k0]x[000]x[100]x[010]01011]12x[kk0]x[k2k101x[110]x[001]x[101]x[011]x[111]01010101時間抽取FFT例：試利用N=4基2時間抽取的FFT流圖計算8點序列x[k]={1,-1,1,-1,2,1,1,2}的DFT。解：根據(jù)基2時間抽取FFT算法原理，8點序列的DFTX[m]可由兩個4點序列的DFTX1[m]和X2[m]表達。如果按照序列x[k]序號的奇偶分解為x1[k]和x2[k]，則存在其中x1[k]={1,1,2,1}，x2[k]={-1,-1,1,2}，X1[m]和X2[m]可通過4點的FFT來計算。時間抽取FFT例：試利用N=4基2時間抽取的FFT流圖計算8點序列x[k]={1,-1,1,-1,2,1,1,2}的DFT。解：

X1[m]={5,-1,1,-1},X2[m]={1,-2+3j,1,-2-3j}利用上述公式，可得序列x[k]的DFTX[m]為X[m]={6,-0.293+3.535j,1+j,-1.707+3.535j,4,-1.707-3.535j,1-j,-0.293-3.535j}時間抽取FFT則x(n)的DFT為由于所以序列x(n)的長度為N,M為自然數(shù)按n的奇偶把x(n)分解為兩個N/2點的子序列N/2點DFTX(K)前半部分X(K)后半部分基2時間抽取FFT算法流圖N=2x[k]={x[0],x[1]}圖4.2.2N點DFT的一次時域抽取分解圖(N=8)與第一次分解相同，將x1(r)按奇偶分解成兩個N/4長的子序列x3(l)和x4(l)，即那么，X1(k)又可表示為(4.2.9)式中同理，由X3(k)和X4(k)的周期性和WmN/2的對稱性Wk+N/4

N/2=-WkN/2最后得到：(4.2.10)用同樣的方法可計算出(4.2.11)其中圖4.2.3N點DFT的第二次時域抽取分解圖(N=8)圖4.2.4N點DIT―FFT運算流圖(N=8)4.2.3DIT―FFT算法與直接計算DFT運算量的比較每一級運算都需要N/2次復數(shù)乘和N次復數(shù)加(每個蝶形需要兩次復數(shù)加法)。所以，M級運算總共需要的復數(shù)乘次數(shù)為復數(shù)加次數(shù)為例如，N=210=1024時圖4.2.5FFT算法與直接計算DFT所需乘法次數(shù)的比較曲線4.2.4DIT―FFT的運算規(guī)律及編程思想FFT的每級（列）計算都是由N個復數(shù)數(shù)據(jù)（輸入）兩兩構(gòu)成一個蝶型（共N/2個蝶形）運算而得到另外N個復數(shù)數(shù)據(jù)（輸出）。當數(shù)據(jù)輸入到存儲器以后，每一組運算的結(jié)果，仍然存放在這同一組存儲器中直到最后輸出。例：將x(0)放在單元A(0)中，將x(4)放在單元A(1)中，W80

放在一個暫存器中。將x(0)+W80x(4)→送回A(0)單元將x(0)-W80x(4)→送回A(1)單元X3(0)X3(1)x(0)x(4)1）

原位運算(亦稱同址計算)x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7)X(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)回顧：N點DIT―FFT運算流圖(N=8)輸入數(shù)據(jù)、中間運算結(jié)果和最后輸出均用同一存儲器。如上所述，N點DIT―FFT運算流圖中，每級都有N/2個蝶形。每個蝶形都要乘以因子WNP，稱其為旋轉(zhuǎn)因子，p稱為旋轉(zhuǎn)因子的指數(shù)。2)旋轉(zhuǎn)因子的變化規(guī)律觀察FFT運算流圖發(fā)現(xiàn)，第L級共有2L-1個不同的旋轉(zhuǎn)因子。N=23=8時的各級旋轉(zhuǎn)因子表示如下：L=1時，WNp=WN/4J，N/4=21=2L，J=0L=2時，WNp=WN/2J，N/2=22=2L，J=0，1L=3時，WNp=WNJ，N=23=2L，J=0，1，2，3對N=2M的一般情況，第L級的旋轉(zhuǎn)因子為：設(shè)序列x(n)經(jīng)時域抽選(倒序)后，存入數(shù)組X中。如果蝶形運算的兩個輸入數(shù)據(jù)相距B個點，應(yīng)用原位計算，則蝶形運算可表示成如下形式：下標L表示第L級運算，XL(J)則表示第L級運算后數(shù)組元素X(J)的值。蝶形運算兩節(jié)點的距離B=2L-1,L=1…,M例如N=8=23,第一級(列)距離為21-1=1,第二級(列)距離為22-1=2，

第三級(列)距離為23-1=4。3)編程思想及流程圖開始送入x(n)和N=2M調(diào)整輸入x(n)的順序for(L=1;L<=M;L++)B=2L-1for(J=0;J<=B-1;J++)p=J·2M-Lfor(k=J;k<=N-1;k=k+2L)輸出結(jié)果結(jié)束4）碼位倒序由N=8蝶形圖看出：原位計算時，F(xiàn)FT輸出的X(k)的次序正好是順序排列的，即X(0)…X(7),但輸入x(n)都不能按自然順序存入到存儲單元中，而是按x(0),x(4),x(2),x(6),x(1),x(5),x(3),x(7)的順序存入存儲單元，即為亂序輸入，順序輸出。這種順序看起來相當雜亂，然而它是有規(guī)律的。即碼位倒讀規(guī)則。n=00n=10n=01n=11n=01n=1101010101

(n2)(000)0(100)4(010)2(110)6(001)1(101)5(011)3(111)7（偶）（奇）碼位倒讀規(guī)則由奇偶分組造成的,以N=8為例說明如下:輸入序列的序號n二進制為(n2n1n0)碼位倒序二進制為(n0n1n2)

00

0

00

0

00

10

0

11

0

04

20

1

00

1

02

30

1

11

1

06

41

0

00

0

11

51

0

11

0

15

61

1

00

1

13

71

1

11

1

17自然順序n二進制nnn

倒位序二進制nnn倒位順序n^210012自然順序n二進制碼表示碼位倒讀碼位倒置順序n’以N=8為例：0123456700000101001110010111011100010001011000110101111104261537看出：碼位倒讀后的順序剛好是數(shù)據(jù)送入計算機內(nèi)的順序。A(1)A(2)A(3)A(4)A(5)A(6)A(7)A(8)x(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7)x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7)變址處理方法存儲單元自然順序變址倒位序?qū)τ跀?shù)N，在其二進制最高位加1，等于加N/2。若已知某個反序號為J，為求下一個反序號，可先判J的最高位：

1)若為0，則把該位變成1（即加N/2）就得到下一個反序號，2)若為1，則需判斷次高位：

①若次高位為0，則把最高位變0（相當減去N/2）后，再把次高位變1（即加N/4）。

②若次高位為1，則需判斷次次高位……分析：倒序排列算法的流程圖正序序列已在數(shù)組A[]中，輸入NLH=N/2,J=LH,N1=N-2J=J-kk=k/2k=LHJ<kJ=J+kT=A(I)A(I)=A(J)A(J)=Tfor(i=1;i<=N1;i++)i≥JNYYN第四節(jié)按頻率抽選的基2-FFT算法在基2快速算法中，頻域抽取法FFT也是一種常用的快速算法，簡稱DIF―FFT。設(shè)序列x(n)長度為N=2M，首先將x(n)前后對半分開，得到兩個子序列，其DFT可表示為如下形式DIF―FFT一次分解運算流圖(N=8)4點DFT4點DFTx(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7)X(0)X(2)X(4)X(6)X(1)X(3)X(5)X(7)x1(0)x1(1)x1(2)x1(3)x2(0)x2(1)x2(2)x2(3)DIF―FFT二次分解運算流圖(N=8)DIF―FFT運算流圖(N=8)時間抽取算法與頻率抽取算法的比較1)頻率抽選法和時間抽選法總的計算量是相同的復乘：復加：2)頻率抽取法和時間抽取法一樣，都適用于原位運算，即蝶形的輸入和輸出占用同一個存儲單元。3)均存在碼位倒序問題。4)頻率抽選法和時間抽選法一樣，基本運算也是蝶形運算。但兩者的蝶形形式略有不同。第五節(jié)IDFT的快速算法－IFFT上述FFT算法流圖也可以用于離散傅里葉逆變換(InverseDiscreteFourierTransform，簡稱IDFT)。比較DFT和IDFT的運算公式：1)旋轉(zhuǎn)因子：2)系數(shù)：DIT―IFFT運算流圖DIT―IFFT運算流圖(防止溢出)如果希望直接調(diào)用FFT子程序計算IFFT，則可用下面的方法：對上式兩邊同時取共軛，得：例1、如果通用計算機的速度為平均每次復乘需要5s，每次復加需要0.5s，用它來計算512點的DFT[x(n)]，問：1）直接計算需要多少時間？2）用FFT需要多少時間？解：1）用DFT進行運算：復乘：T1=N2×5×10-6=1.31072秒復加：T2=N(N-1)×0.5×10-6=0.130816秒總共：T=T1+T2=1.441536秒2）用FFT進行運算：復乘：T1’=(N/2)log2N×5×