數(shù)字信號處理-快速傅里葉變換課件

1、第四章 快速傅立葉變換 Fast Fourier Transform 主要內(nèi)容：第一節(jié) 直接計算DFT的問題及改進途徑第二節(jié) 改善DFT運算效率的基本途徑第三節(jié) 按時間抽選的基2-FFT算法第四節(jié) 按頻率抽選的基2-FFT算法第一節(jié) 直接計算DFT的問題及改進途徑1、問題的提出 設(shè)有限長序列x(n)，非零值長度為N，若對x(n)進行一次DFT運算，共需多大的運算工作量？計算成本?計算速度?2. DFT的運算量 回憶DFT和IDFT的變換式： 1）x(n)為復數(shù)， 也為復數(shù)。2）DFT與IDFT的計算量相當。注意： 計算機運算時（編程實現(xiàn)）： N次復乘，N-1次復加 N個點 以DFT為例： 復數(shù)

2、乘法復數(shù)加法一個X(k)NN 1N個X(k)(N點DFT)N 2N (N 1)實數(shù)乘法實數(shù)加法一次復乘42一次復加2一個X (k)4N2N+2 (N 1)=2 (2N 1)N個X (k)(N點DFT)4N 22N (2N 1)運算量(a+jb)(c+jd)=(ac-bd)+j(bc+ad)例：計算一個 N點DFT ，共需N2次復乘。以做一次 復乘1s計，若N =4096，所需時間為例：石油勘探，有24個通道的記錄，每通道波形記 錄長度為5秒，若每秒抽樣500點/秒， 1）每道總抽樣點數(shù)：500*5=2500點 2）24道總抽樣點數(shù)：24*2500=6萬點 3）DFT復乘運算時間：N2=(600

3、00)2=36*108次 由于計算量大，且要求相當大的內(nèi)存，難以實現(xiàn)實時處理，限制了DFT的應(yīng)用。長期以來，人們一直在尋求一種能提高DFT運算速度的方法。 FFT便是 Cooley & Tukey 在1965 年提出的的快速算法，它可以使運算速度提高幾百倍，從而使數(shù)字信號處理學科成為一個新興的應(yīng)用學科。第二節(jié) 改善DFT運算效率的基本途徑 1、利用DFT運算的系數(shù) 的固有對稱性和周期 性，改善DFT的運算效率。 1）對稱性 2）周期性 3）可約性2、將長序列DFT利用對稱性和周期性分解為短 序列DFT的思路 因為DFT的運算量與N2成正比的，如果一個大點數(shù)N的DFT能分解為若干小點數(shù)DFT的組

4、合，則顯然可以達到減少運算工作量的效果。N點DFTN/2點DFTN/2點DFTN/4點DFTN/4點DFTN/4點DFTN/4點DFT.復乘： FFT算法的基本思想： 利用DFT系數(shù)的特性，合并DFT運算中的某些項 把長序列DFT短序列DFT，從而減少運算量。FFT算法分類:時間抽選法 DIT: Decimation-In-Time頻率抽選法 DIF: Decimation-In-Frequency第三節(jié) 按時間抽選的基2-FFT算法1、算法原理 設(shè)輸入序列長度為N=2M(M為正整數(shù)，將該序列按時間順序的奇偶分解為越來越短的子序列，稱為基2按時間抽取的FFT算法。也稱為Coolkey-Tuke

5、y算法。 其中基2表示：N=2M，M為整數(shù).若不滿足這個條件，可以人為地加上若干零值（加零補長）使其達到 N=2M。先將x(n)按n的奇偶分為兩組，作變量置換: 當n=偶數(shù)時，令n=2r; 當n=奇數(shù)時，令n=2r+1; 分組，變量置換2、算法步驟得到： 帶入DFT中所以 由于? X1(k)、X2(k)只有N/2個點，以N/2為周期；而X (k)卻有N個點，以N為周期。要用X1(k)、X2(k)表達全部的X (k) 值，還必須利用WN系數(shù)的周期特性。后半部分前半部分又考慮到 的對稱性：有：后半部分前半部分蝶形運算流圖符號說明： (1) 左邊兩路為輸入 (2) 右邊兩路為輸出 (3) 中間以一個

6、小圓表示加、 減運算（右上路為相加 輸出、右下路為相減輸 出)1個蝶形運算需要1次復乘，2次復加復數(shù)乘法復數(shù)加法一個N 點DFTN 2N (N1)一個N / 2點DFT(N / 2)2N / 2 (N / 2 1)兩個N / 2點DFTN 2 / 2N (N / 2 1)一個蝶形12N / 2個蝶形N / 2N總計N2/2 + N/2 N2/2N(N/2-1) + N N2/2運算量減少了近一半 分解后的運算量：先將N=8點的DFT分解成2個4點DFT：可知：時域上：x(0),x(2),x(4),x(6)為偶子序列 x(1),x(3),x(5),x(7)為奇子序列 頻域上：X(0)X(3)，由

7、X(k)給出 X(4)X(7)，由X(k+N/2)給出例子：求 N=23=8點FFT變換 按N=8N/2=4，做4點的DFT： N=8點的直接DFT的計算量為： 復乘：N2次 = 64次 復加：N(N-1)次 = 87=56次 此外，還有4個蝶形結(jié)，每個蝶形結(jié)需要1次復乘，2次復加。一共是：復乘4次，復加8次。得到X1(k)和X2(k)需要： 復乘：(N/2)2+ (N/2)2次 = 32次 復加：N/2(N/2-1)+N/2(N/2-1) =12+12 =24次用分解的方法得到X (k)需要： 復乘：32+4 = 36次 復加：24+8 = 32次N點DFT的一次時域抽取分解圖(N=8) 4

8、點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)因為4點DFT還是比較麻煩，所以再繼續(xù)分解。 若將N/2(4點)子序列按奇/偶分解成兩個N/4點(2點)子序列。即對將x1(r)和x2(r)分解成奇、偶兩個N/4點(2點)點的子序列。那么，X1(k)又可表示為 X2(k)也可以進行相同的分解： 注意：通常我們會把 寫成 。N點DFT的第二次時域抽取分解圖(N=8) 2點DFT2點DFT2點DFT2點DFTx(0)x(4)

9、x(2)x(6)x(1)x(5)x(3)x(7)X3(0)X3(1)X4(0)X4(1)X5(0)X5(1)X6(0)X6(1)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)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)88X3(0)X3(1)x(0)=x3(0)x(4)=x3(1)N點DITFFT運算流圖(N=8)

10、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) 3、DITFFT算法與直接計算DFT運算量的比較1)、N=2M的DFT運算可分成M級，每一級有N/2個蝶形 ，每個蝶形有一次復乘兩次復加。2)、所以M級共有 次復乘和 次復加。3)、若直接計算DFT，需N2次復乘和N(N-1)次復加。顯然，當N較大時，有：例如，N=210=1024時FFT算法與直接計算DFT所需乘法次數(shù)的比較曲線4*、DITFFT的運算規(guī)律及編程思想 FFT的每級（列）計算都是由N個復數(shù)數(shù)據(jù)（輸入）兩兩構(gòu)成一個蝶型（共

11、N/2個蝶形）運算而得到另外N個復數(shù)數(shù)據(jù)（輸出）。 當數(shù)據(jù)輸入到存儲器以后，每一組運算的結(jié)果，仍然存放在這同一組存儲器中直到最后輸出。例：將x(0)放在單元A(0)中，將x(4)放在單元A(1)中，W80 放在一個暫存器中。將x(0) + W80 x(4) 送回A(0)單元將x(0) - W80 x(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點DITFFT運算流圖(N=8) 如上所

12、述，N點DITFFT運算流圖中，每級都有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)用原位計算，則蝶

13、形運算可表示成如下形式： 下標L表示第L級運算，XL(J)則表示第L級運算后數(shù)組元素X(J)的值。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 = J2M-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)的順序存入存儲單

14、元，即為亂序輸入，順序輸出。 這種順序看起來相當雜亂，然而它是有規(guī)律的。即碼位倒讀規(guī)則。自然順序n二進制碼表示碼位倒讀碼位倒置順序n以N=8為例：0123456700000101001110010111011100010001011000110101111104261537看出：碼位倒讀后的順序剛好是數(shù)據(jù)送入計算機內(nèi)的順序。倒序規(guī)律 對于數(shù)N，在其二進制最高位加1，等于加N/2。 若已知某個反序號為J，為求下一個反序號，可先判J的最高位： 1) 若為0，則把該位變成1（即加N/2）就得到下 一個反序號， 2) 若為1，則需判斷次高位： 若次高位為0，則把最高位變0（相當減去 N/2）后，再把次

15、高位變1（即加N/4）。 若次高位為1，則需判斷次次高位分析：倒序排列算法的流程圖正序序列已在數(shù)組A 中，輸 入 NLH= N/2 , J=LH , N1 = N-2J=J-kk=k/2 k=LHJkJ=J+k T = A(I) A(I) = A(J) A(J) = Tfor(i=1; i=N1; i+) i JNYYN第四節(jié) 按頻率抽選的基2-FFT算法 在基2快速算法中，頻域抽取法FFT也是一種常用的快速算法，簡稱DIFFFT。 設(shè)序列x(n)長度為N=2M，首先將x(n)前后對半分開，得到兩個子序列，其DFT可表示為如下形式DIFFFT一次分解運算流圖(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)