快速付里葉變換FFT2

1、2.3 快速傅立葉變換FFT產(chǎn)生1965年，庫(kù)利-圖基在計(jì)算數(shù)學(xué)（Mathematic of Computation）雜志上發(fā)表了著名的“機(jī)器計(jì)算傅里級(jí)數(shù)的一種算法”文章，提出一種快速計(jì)算DFT的方法和計(jì)算機(jī)程序-揭開了FFT發(fā)展史上的第一頁(yè)。本節(jié)主要內(nèi)容 直接計(jì)算DFT算法存在的問(wèn)題及改進(jìn)途徑。 FFT算法(基2時(shí)間抽取算法） FFT的應(yīng)用直接計(jì)算DFT計(jì)算量問(wèn)題提出：設(shè)有限長(zhǎng)序列x(n)，非零值長(zhǎng)度為N，計(jì)算對(duì)x(n)進(jìn)行一次DFT運(yùn)算，共需多大的運(yùn)算工作量?一、直接計(jì)算DFT算法存在的問(wèn)題及改進(jìn)途徑(1) 以DFT為例，計(jì)算DFT復(fù)數(shù)運(yùn)算量計(jì)算一個(gè)X(k)(一個(gè)頻率成分)值，運(yùn)算量為 例

2、（k=1）則 10)() 1 (NnnNWnxX要進(jìn)行要進(jìn)行N次復(fù)數(shù)乘法次復(fù)數(shù)乘法+(N-1)次復(fù)數(shù)加法次復(fù)數(shù)加法 所以，要完成整個(gè)所以，要完成整個(gè)DFT運(yùn)算，其計(jì)算量為：運(yùn)算，其計(jì)算量為： N N* *N N次復(fù)數(shù)相乘次復(fù)數(shù)相乘+ +N N* *(N-1)(N-1)次復(fù)數(shù)加法次復(fù)數(shù)加法(2) 一次復(fù)數(shù)乘法換算成實(shí)數(shù)運(yùn)算量復(fù)數(shù)運(yùn)算要比加法運(yùn)算復(fù)雜，需要的運(yùn)算時(shí)間長(zhǎng)。一個(gè)復(fù)數(shù)乘法包括4個(gè)實(shí)數(shù)乘法和個(gè)實(shí)數(shù)乘法和2個(gè)實(shí)數(shù)相個(gè)實(shí)數(shù)相法法。(a+jb)(c+jd)=(ac-bd)+j(bc+ad) 4次實(shí)數(shù)乘法2次實(shí)數(shù)加法(3) 計(jì)算DFT需要的實(shí)數(shù)運(yùn)算量 每運(yùn)算一個(gè)X(k)的值，需要進(jìn)行 4N次實(shí)數(shù)相

3、乘和 2N+2(N-1)=2(2N-1)次實(shí)數(shù)相加. 整個(gè)DFT運(yùn)算量為：4N2次實(shí)數(shù)相乘和次實(shí)數(shù)相乘和2N(2N-1)次實(shí)數(shù)相次實(shí)數(shù)相加加. 由此看出：直接計(jì)算DFT時(shí)，乘法次數(shù)與加法次數(shù)都是和N2成比例的。當(dāng)N很大時(shí)，所需工作量非?？捎^。)Re)(ImIm)(Re) Im)(ImRe)(Re)(10knNknNNnknNknNWnxWnxjWnxWnxkX當(dāng)N=1024點(diǎn)時(shí)，直接計(jì)算DFT需要： N2=1048576次，即一百多萬(wàn)次的復(fù)乘運(yùn)算。 這對(duì)實(shí)時(shí)性很強(qiáng)的信號(hào)處理(如雷達(dá)信號(hào)處理)來(lái)講，它對(duì)計(jì)算速度有十分苛刻的要求 迫切需要改進(jìn)DFT的計(jì)算方法，以減少總的運(yùn)算次數(shù)。(4) 比較DFT

4、與IDFT之間的運(yùn)算量10)()()(NnknNDFTWnxkXnx1, 1 , 0Nk101( )( )( )NIDFTknNkX kx nX k WN 1, 1 , 0Nn其中x(n)為復(fù)數(shù)， 也為復(fù)數(shù)。所以DFT與IDFT二者計(jì)算量相同。knNjknNeW2 2. 改善DFT運(yùn)算效率的基本途徑 利用DFT運(yùn)算的系數(shù) 的固有對(duì)稱性和周期性，改善DFT的運(yùn)算效率。 合并法：合并DFT運(yùn)算中的某些項(xiàng)。 分解法：將長(zhǎng)序列DFT利用對(duì)稱性和周 期性，分解為短序列DFT。knNW1)()(1)()(22222jNNjNNjNNjNNeeWeeWrNrNNkNkNjkNNNkNNWWWWeWWW)()

5、(22利用DFT運(yùn)算的系數(shù) 的固有對(duì)稱性和周期性，改善DFT的運(yùn)算效率knNW利用DFT運(yùn)算的系數(shù) 的固有對(duì)稱性和周期性，改善DFT的運(yùn)算效率knNWknNW的對(duì)稱性：*)()(knNnkNnNkNWWWknNW的周期性：)()(kNnNkNnNknNWWW1)()(22kjkNNjkNNeeW因?yàn)椋?)()(22njNnNjNnNeeW由此可得出：knNknNkNNnNkNnkNnkNkNNnNkNWWWWWWWW)()(kNNkNWW)(2同理可得出例:1454)54(494WWWW1898178258WWWW利用以上特性,得到改進(jìn)DFT直接算法的方法.將長(zhǎng)序的DFT利用對(duì)稱性和周期性分解

6、為短序列DFTDFT的運(yùn)算量與N2成正比如果一個(gè)大點(diǎn)數(shù)N的DFT能分解為若干小點(diǎn)數(shù)DFT的組合，則顯然可以達(dá)到減少運(yùn)算工作量的效果。方 法22)()()(NNNkXkXkX把N點(diǎn)數(shù)據(jù)分成二半：其運(yùn)算量為：2)2(N2)2(N2N再分二半：44)()(NNkXkX44)()(NNkXkX+=22N2)4(N2)4(N2)4(N2)4(N+=24N這樣一直分下去，剩下兩點(diǎn)的變換。結(jié) 論快速付里時(shí)變換(FFT)就是在此特性基礎(chǔ)上發(fā)展起來(lái)的，并產(chǎn)生了多種FFT算法，其基本上可分成兩大類：按抽取方法分： 時(shí)間抽取法(DIT)；頻率抽取法(DIF)按“基數(shù)”分：基-2FFT算法；基-4FFT算法；混合基F

7、FT算法；分裂基FFT算法其它方法：線性調(diào)頻Z變換(CZT法)二、按時(shí)間抽取的FFT算法Decimation-in-Time(DIT)1 算法原理設(shè)輸入序列長(zhǎng)度為N=2M(M為正整數(shù)，將該序列按時(shí)間順序的奇偶分解為越來(lái)越短的子序列，稱為基2按時(shí)間抽取的FFT算法。也稱為Coolkey-Tukey算法。其中基數(shù)2-N=2M，M為整數(shù)。若不滿足這個(gè)條件，可以人為地加上若干零值（加零補(bǔ)長(zhǎng)）使其達(dá)到 N=2M。例子設(shè)一序列設(shè)一序列x(n)的長(zhǎng)度為的長(zhǎng)度為L(zhǎng)=9，應(yīng)加零補(bǔ)長(zhǎng)，應(yīng)加零補(bǔ)長(zhǎng)為為N=24=16 應(yīng)補(bǔ)應(yīng)補(bǔ)7個(gè)零值。個(gè)零值。 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 1

8、5 nx(n)2 算法推導(dǎo)（1）分組分組DFT變換：10)()(NnknNWnxkX1,0Nk先將x(n)按n的奇偶分為兩 組，作變量置換:當(dāng)n=偶數(shù)時(shí)，令n=2r；當(dāng)n=奇數(shù)時(shí)，令n=2r+1；得到：x(2r)=x1(r); x(2r+1)=x2(r)；r=0,N/2-1；(2) DFT)( 2)( 1) 12 ()2 ()()rxDFTrxDFTrxDFTrxDFTnxDFTkX（12/, 0Nr生成兩個(gè)子序列x(0),x(2)x(2r)偶數(shù)點(diǎn)x(1),x(3)x(2r+1)奇數(shù)點(diǎn)代入DFT變換式：kNrkNNrNrrkNkrNNrNrrkNWWrxWrxWrxWrxkX)( )(2)(

9、1) 12()2()212/012/02)12(12/012/02（(3) 求出子序列的DFTkNkNjkNjkNWeeW2/2/222212/1 , 0)()()()()212/12/0212/02/1NkkXWkXWWrxWrxkXkNkNrkNNrNrrkN（12/, 0Nk12/012/02/2/2212/012/02/2/11) 12()()()2()()(NrNrrkNrkNNrNrrkNrkNWrxWrxkXWrxWrxkX其中：(4) 結(jié)論1利用周期性質(zhì)，一個(gè)N點(diǎn)的DFT被分解為兩個(gè)N/2點(diǎn)DFT。X1(k),X2(k)這兩個(gè)N/2點(diǎn)的DFT按照：12/1 , 0)()()21

10、NkDFTNkXWkXkXkN中的前半部分點(diǎn)合成（再應(yīng)用再應(yīng)用W系數(shù)的周期性，求出用系數(shù)的周期性，求出用X1(k),X2(k)表達(dá)的后半部的表達(dá)的后半部的X(k+N/2)的值。的值。（5）后半部的表示式12/012/022/2)2/(2/2212/012/012/1)2/(2/11)2/(2/2/)()()()2()()()()2(NrNrrkNkNrNNrNrrkNkNrNNkrNrkNkXWrxWrxkNXkXWrxWrxkNXWW后半部的后半部的k值所對(duì)應(yīng)的值所對(duì)應(yīng)的X1(k)，X2(k)則完全重復(fù)了前半部則完全重復(fù)了前半部分的分的k值所對(duì)應(yīng)的值所對(duì)應(yīng)的X1(k), X2(k)的值。的值

11、。)()()(212/)2/kXWkXkXWWWWkNkNkNNNkNN后半部分：又（(6) 結(jié)論2頻域中的N個(gè)點(diǎn)頻率成分為：)()()2/()()()(2121kXWkXNkXkXWkXkXkNkN后半部分：前半部分：結(jié)論：只要求出（結(jié)論：只要求出（0N/2-1）區(qū)間內(nèi)的各個(gè)整區(qū)間內(nèi)的各個(gè)整數(shù)數(shù)k值所對(duì)應(yīng)的值所對(duì)應(yīng)的X1(k)，X2(k)值，即可以求出值，即可以求出(0N-1)整個(gè)區(qū)間內(nèi)全部整個(gè)區(qū)間內(nèi)全部X(k)值，這就是值，這就是FFT能能節(jié)省計(jì)算的關(guān)鍵。節(jié)省計(jì)算的關(guān)鍵。(7）結(jié)論3由于由于N=2M，因此因此N/2仍為偶數(shù)，可以依仍為偶數(shù)，可以依照上面方法進(jìn)一步把每個(gè)照上面方法進(jìn)一步把每個(gè)

12、N/2點(diǎn)子序列，點(diǎn)子序列，再按輸入再按輸入n的奇偶分解為兩個(gè)的奇偶分解為兩個(gè)N/4點(diǎn)的子點(diǎn)的子序列，按這種方法不斷劃分下去，直到序列，按這種方法不斷劃分下去，直到最后剩下的是最后剩下的是2點(diǎn)點(diǎn)DFT，兩點(diǎn)，兩點(diǎn)DFT實(shí)際上實(shí)際上只是只是加減運(yùn)算加減運(yùn)算。3 蝶形結(jié)（基本單元）即蝶式計(jì)算結(jié)構(gòu)也即為蝶式信號(hào)流圖上面頻域頻域中前/后半部分表示式可以用蝶形信號(hào)流圖表示。X1(k)X2(k)kNW)()(21kXWkXkN)()(21kXWkXkN將N=8點(diǎn)分解成2個(gè)4點(diǎn)的DFT的信號(hào)流圖3821282118210821) 3 () 3 (3) 2() 2(2) 1 () 1 (1) 0() 0(0WX

13、XXWXXXWXXXWXXX）（）（）（）（如：4點(diǎn)DFTx(0)x(2)x(4)x(6)4點(diǎn)DFTx(1)x(3)x(5)x(7)08W18W28W38WX(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)X1(0)X1(1)X1(2)X1(3)X2(0)X2(1)X2(2)X2(3)偶數(shù)序列奇數(shù)序列x1(r)x2(r)一個(gè)2點(diǎn)的DFT蝶形流圖2點(diǎn)DFT2點(diǎn)DFTx(0)x(4)x(2)x(6)X3(0)X3(1)X4(0)X4(1)X1(0)X1(1)X1(2)X1(3)04W14W) 1 () 1 () 3 () 0() 0() 2() 1 () 1 () 1 () 0()

14、0() 0(41431404314143140431XWXXXWXXXWXXXWXX其中 另一個(gè)2點(diǎn)的DFT蝶形流圖) 1 () 1 () 3 () 0() 0() 2() 1 () 1 () 1 () 0() 0() 0(61452604526145260452XWXXXWXXXWXXXWXX其中2點(diǎn)DFT2點(diǎn)DFTx(1)x(5)x(3)x(7)X5(0)X5(1)X6(0)X6(1)X2(0)X2(1)X2(2)X2(3)04W14W同理：將N/4（2點(diǎn)）DFT再分解成2個(gè)1點(diǎn)的DFT0210221202123020231021 , 0; 1 , 0)4()0()4()0() 1 ()4

15、()0()4()0()0()()(2WWknWWWWWxxWxxXWxxWxxXWnxkXnknknkNnkNnnkN，其中，則這里用到對(duì)稱性這是一蝶形結(jié)代入上面流圖可知： 最后剩下兩點(diǎn)DFT,它可分解成兩個(gè)一點(diǎn)DFT，但一點(diǎn)但一點(diǎn)DFT就等于輸入信號(hào)本身就等于輸入信號(hào)本身，所以兩點(diǎn)DFT可用一個(gè)蝶形結(jié)表示。取x(0)、x(4)為例。2個(gè)1點(diǎn)的DFT蝶形流圖 1點(diǎn)DFTx(0)1點(diǎn)DFTx(4)X3(0)X3(1)02W進(jìn)一步簡(jiǎn)化為蝶形流圖：02WX3(0)X3(1)x(0)x(4)4()0()4()0() 1 ()4()0()4()0()0(023023xxxWxXxxxWxX其中：一個(gè)完整

16、N=8的按時(shí)間抽取FFT的運(yùn)算流圖 12/0NNNWW其中蝶形因子，共有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)38W28W18W08W08W08W08W08W08W08W28W28Wm=1m=2m=34 FFT算法中一些概念將N 點(diǎn)DFT先分成兩個(gè)N/2點(diǎn)DFT，再是四個(gè)N/4點(diǎn)DFT直至N/2個(gè)兩點(diǎn)DFT.每分一次稱為“一”級(jí)運(yùn)算。因?yàn)镹=2M所以N點(diǎn)DFT可分成M級(jí)如上圖所示依次m=1,m=2.M共M級(jí)（1 1）級(jí)）級(jí)（2 2）組）組 每一級(jí)都有N/2個(gè)蝶形單元，例如：N=8,則每級(jí)都有4個(gè)蝶形單元。

17、每一級(jí)的N/2個(gè)蝶形單元可以分成若干組，每一組具有相同的結(jié)構(gòu)，相同的 因子分布，第m級(jí)的組數(shù)為：rNWmN2例：N=8=23，分3級(jí)。m=1級(jí),分成四組,每組系數(shù)為m=2級(jí),分成二組,每組系數(shù)為m=3級(jí),分成一組,每組系數(shù)為080204/WWWN28082012/02/,WWWWWWNNNN382818083210,WWWWWWWWNNNN(3) 因子的分布rNW121 , 0,3 , 2 , 10310201238281808828081404408022mkkkkkWmWWWWkWmWWWWkWmWWkWmm級(jí)的系數(shù)為看出：第，級(jí)，級(jí)，級(jí)，由上可知：結(jié)論：每由后向前（結(jié)論：每由后向前（m由

18、由M-1級(jí)）推進(jìn)一級(jí)，級(jí)）推進(jìn)一級(jí)，則此系數(shù)為后級(jí)系數(shù)中偶數(shù)序號(hào)的那一半。則此系數(shù)為后級(jí)系數(shù)中偶數(shù)序號(hào)的那一半。N=8點(diǎn)和N=1024點(diǎn)時(shí)直接計(jì)算DFT與用基2-按時(shí)間抽取法FFT的運(yùn)算量。4.102102227.224648loglog1020102222NNNNNNN看出：當(dāng)N較大時(shí)，按時(shí)間抽取法將比直接法快一、二個(gè)數(shù)量級(jí)之多。5 按時(shí)間抽取FFT算法的特點(diǎn)根據(jù)DIT基2-FFT算法原理，能得出任何N=2m點(diǎn)的FFT信號(hào)流圖，并進(jìn)而得出FFT計(jì)算程序流程圖。最后總結(jié)出按時(shí)間抽取法解過(guò)程的規(guī)律。原位運(yùn)算（in-place)碼位倒序規(guī)則（1） 原位運(yùn)算（in-place)原位運(yùn)算的結(jié)構(gòu)，可以節(jié)

19、省存儲(chǔ)單元，降低設(shè)備成本。定義：當(dāng)數(shù)據(jù)輸入到存儲(chǔ)器以后，每一組運(yùn)算的結(jié)果，仍然存放在這同一組存儲(chǔ)器中直到最后輸出。02Wx(0)x(4)X3(0)X3(1)中放在一個(gè)暫存器中放在單元中，將放在單元例：將RWAxAx08) 1 ()4()0()0(單元送回單元送回將) 1 ()4()0()0()4()0(0808AxWxAxWx(2) 碼位倒序規(guī)則我們從輸入序列的序號(hào)及整序規(guī)律得到碼位倒讀規(guī)則。由N=8蝶形圖看出：原位計(jì)算時(shí)，F(xiàn)FT輸出的X(k)的次序正好是順序排列的，即X(0)X(7),但輸入x(n)都不能按自然順序存入到存儲(chǔ)單元中，而是按x(0),x(4),x(2), x(6).的順序存入存

20、儲(chǔ)單元即為亂序輸入亂序輸入，順序輸出順序輸出。這種順序看起來(lái)相當(dāng)雜亂，然而它是有規(guī)律的。即碼位倒序規(guī)則。以N=8為例：01234567000001010011100101110111自然順序二進(jìn)制碼表示碼位倒讀碼位倒置順序00010001011000110101111104261537看出：碼位倒讀后的順序剛好是數(shù)據(jù)送入計(jì)算機(jī)內(nèi)的順序。排序子程序具體執(zhí)行時(shí)，只須將1與4對(duì)調(diào)送入，3與6對(duì)調(diào)送入，而0，2，5，7不變，僅需要一個(gè)中間存儲(chǔ)單元。n01234567n04261537 在實(shí)際運(yùn)算時(shí)，先按自然順序?qū)⑤斎胄蛄写嫒氪鎯?chǔ)單元，再通過(guò)變址運(yùn)算將自然順序變換成按時(shí)間抽取的FFT算法要求的順序。變址的過(guò)程可以用程序安排加以實(shí)現(xiàn)-稱為“整序”或“重排”（采用碼位倒讀）且注意：(1)當(dāng)n=n時(shí)，數(shù)據(jù)不必調(diào)換；(2)當(dāng)nn時(shí)，必須將原來(lái)存放數(shù)據(jù)x(n)送入暫存器R，再將x(n)送入x(n),R中內(nèi)容送x(n).進(jìn)行數(shù)據(jù)對(duì)調(diào)。(3)為了避免再次考慮前面已調(diào)換過(guò)的數(shù)據(jù)，保證調(diào)換只進(jìn)行一次，否則又變回原狀。nn時(shí)，調(diào)換。三、按頻率抽取的FFT算法可以將x(k)細(xì)分為越來(lái)越短的子序列，稱為按頻率抽取的FFT算法。1, 2 , 1)()(0NkWnxkXNnnkN已知分為偶數(shù)和