1、第五章 快速傅里葉變換,2,本章目錄,直接計算DFT的問題及改進的途徑,按時間抽取的基2-FFT算法,按頻率抽取的基2-FFT算法,快速傅里葉逆變換(IFFT)算法,Matlab實現(xiàn),3,5.1 引言,DFT在實際應(yīng)用中很重要: 可以計算信號的頻譜、功率譜和線性卷積等。 直接按DFT變換進行計算，當序列長度N很大時，計算量非常大，所需時間會很長。 FFT并不是一種與DFT不同的變換，而是DFT的一種快速計算的算法。,4,5.2 直接計算DFT的問題及改進的途徑,DFT的運算量,設(shè)復(fù)序列x(n) 長度為N點，其DFT為,k=0，N-1,（1）計算一個X(k) 值的運算量,復(fù)數(shù)乘法次數(shù)：,N,復(fù)數(shù)

2、加法次數(shù)：,N1,5,5.2.1 DFT的運算量,（2）計算全部N個X(k) 值的運算量,復(fù)數(shù)乘法次數(shù)：,N2,復(fù)數(shù)加法次數(shù)：,N(N1),（3）對應(yīng)的實數(shù)運算量,6,一次復(fù)數(shù)乘法：,4次實數(shù)乘法,2次實數(shù)加法,一個X(k) ：,4N次實數(shù)乘法,2N+2(N-1)= 2(2N-1)次實數(shù)加法,所以,整個N點DFT運算共需要：,N2(2N-1)= 2N(2N-1),實數(shù)乘法次數(shù)：,4 N2,實數(shù)加法次數(shù)：,7,DFT運算量的結(jié)論,N點DFT的復(fù)數(shù)乘法次數(shù)舉例,結(jié)論：當N很大時，其運算量很大，對實時性很強的信號處理來說，要求計算速度快，因此需要改進DFT的計算方法，以大大減少運算次數(shù)。,8,5.2

3、.2 減少運算工作量的途徑,主要原理是利用系數(shù) 的以下特性對DFT進行分解：,（1）對稱性,（2）周期性,（3）可約性,另外，,9,5.3 按時間抽取的基2-FFT算法,算法原理 按時間抽取基-2FFT算法與直接計算DFT運算量的比較 按時間抽取的FFT算法的特點 按時間抽取FFT算法的其它形式流程圖,10,5.3.1 算法原理,設(shè)N2L，將x(n)按 n 的奇偶分為兩組：,r =0，1，,則,11,式中，X1(k)和X2(k)分別是x1(n)和x2(n)的N/2的DFT。,另外，式中k的取值范圍是：0，1， ，N/21 。,12,因此， 只能計算出X(k)的前一半值。,后一半X(k) 值，

4、N/2 ， N/2 1， ，N ？,利用,可得到,同理可得,13,考慮到,因此可得后半部分X(k),及前半部分X(k),k=0，1， ，N/21,k=0，1， ，N/21,14,蝶形運算,蝶形運算式,蝶形運算信號流圖符號,因此，只要求出2個N/2點的DFT，即X1(k)和X2(k)，再經(jīng)過蝶形運算就可求出全部X(k)的值，運算量大大減少。,15,以8點為例第一次按奇偶分解,以N=8為例，分解為2個4點的DFT，然后做8/2=4次蝶形運算即可求出所有8點X(k)的值。,16,蝶形運算量比較,復(fù)數(shù)乘法次數(shù)：,N2,復(fù)數(shù)加法次數(shù)：,N(N1),復(fù)數(shù)乘法次數(shù)：,2*(N/2)2+N/2=N2/2+N/

5、2,復(fù)數(shù)加法次數(shù)：,2*(N/2)(N/21)+2*N/2=N2/2,N點DFT的運算量,分解一次后所需的運算量2個N/2的DFTN/2蝶形：,因此通過一次分解后，運算工作量減少了差不多一半。,17,進一步按奇偶分解,由于N2L，因而N/2仍是偶數(shù) ，可以進一步把每個N/2點 子序列再按其奇偶部分分解為兩個N/4點的子序列。,以N/2點序列x1(r)為例,則有,k=0,1,18,且,k=0,1,由此可見，一個N/2點DFT可分解成兩個N/4點DFT。,同理，也可對x2(n)進行同樣的分解，求出X2(k)。,19,以8點為例第二次按奇偶分解,20,算法原理,對此例N=8，最后剩下的是4個N/4=

6、 2點的DFT，2點DFT也可以由蝶形運算來完成。以X3(k)為例。,k=0, 1,即,這說明，N=2M的DFT可全部由蝶形運算來完成。,21,以8點為例第三次按奇偶分解,N=8按時間抽取法FFT信號流圖,22,5.3.2 按時間抽取基2-FFT算法與直接計算DFT運算量的比較,由按時間抽取法FFT的信號流圖可知，當N=2L時，共有 級 蝶形運算；每級都由 個蝶形運算組成，而每個蝶形有 次復(fù)乘、 次復(fù)加，因此每級運算都需 次復(fù)乘和 次復(fù)加。,L,N/2,N/2,1,2,N,23,這樣 級運算總共需要：,L,復(fù)數(shù)乘法：,復(fù)數(shù)加法：,直接DFT算法運算量,復(fù)數(shù)乘法：,復(fù)數(shù)加法：,N2,N(N1),

7、直接計算DFT與FFT算法的計算量之比為M,24,FFT算法與直接DFT算法運算量的比較,25,5.3.3 按時間抽取的FFT算法的特點,序列的逆序排列 同址運算（原位運算） 蝶形運算兩節(jié)點間的距離 的確定,26,序列的逆序排列,由于 x(n) 被反復(fù)地按奇、偶分組，所以流圖輸入端的 排列不再是順序的，但仍有規(guī)律可循：,因為 N=2M ，,對于任意 n（0n N-1)，可以用M個二進制碼表示為：,n 反復(fù)按奇、偶分解時，即按二進制碼的“0” “1” 分解。,序列的逆序排列,27,倒位序的樹狀圖（N=8）,28,碼位的倒位序(N=8),29,倒位序的變址處理（N=8）,30,同址運算（原位運算）

8、,某一列任何兩個節(jié)點k 和j 的節(jié)點變量進行蝶形運算后，得到結(jié)果為下一列k、j兩節(jié)點的節(jié)點變量，而和其他節(jié)點變量無關(guān)。這種原位運算結(jié)構(gòu)可以節(jié)省存儲單元，降低設(shè)備成本。,運算前,運算后,例,同址運算（原位運算）,31,觀察原位運算規(guī)律,32,蝶形運算兩節(jié)點間的距離,以N=8為例：,第一級蝶形，距離為：,第二級蝶形，距離為：,第三級蝶形，距離為：,規(guī)律：對于共L級的蝶形而言，其m級蝶形運算的節(jié) 點間的距離為,1,2,4,蝶形運算兩節(jié)點間的距離,33,的確定,以N=8為例：,的確定,34,5.4 按頻率抽取的基2-FFT算法,算法原理,再把輸出X(k)按k的奇偶分組,先把輸入按n的順序分成前后兩半,

9、設(shè)序列長度為N=2L，L為整數(shù),前半子序列x(n),后半子序列,0n,0n,35,5.4.1 算法原理,由DFT定義得,k=0，1， ，N,36,由于,所以,則,k=0，1， ，N,37,然后按k的奇偶可將X(k)分為兩部分,r=0，1， ，,則式,可轉(zhuǎn)化為,38,令,n=0，1， ，,代入,r=0，1， ，,可得,為2個N/2點的DFT，合起來正好是N點X(k)的值。,39,蝶形運算,將,稱為蝶形運算,與時間抽選基2FFT算法中的蝶形運算符號略有不同。,40,例 按頻率抽取(N=8),例 按頻率抽取，將N點DFT分解為兩個N/2點DFT的組合(N=8),41,與時間抽取法的推導(dǎo)過程一樣，由于

10、 N=2L，N/2仍然是 一個偶數(shù)，因而可以將每個N/2點DFT的輸出再分解為偶數(shù)組 與奇數(shù)組，這就將N/2點DFT進一步分解為兩個N/4點DFT。,N=8,42,5.4.2 頻率抽取法與時間抽取法的異同,頻率抽取法輸入是自然順序，輸出是倒位序的；時間抽取法正好相反。 頻率抽取法的基本蝶形與時間抽取法的基本蝶形有所不同。 頻率抽取法運算量與時間抽取法相同。 頻率抽取法與時間抽取法的基本蝶形是互為轉(zhuǎn)置的。,43,5.5 快速傅里葉逆變換(IFFT)算法,IDFT公式,DFT公式,比較可以看出，,IDFT多出,M個1/2可分解到M級蝶形運算中。,44,例 頻率抽取IFFT流圖(N=8),45,快速

11、傅里葉逆變換另一種算法,46,5.8 Matlab實現(xiàn),用FFT進行譜分析的Matlab實現(xiàn) 用CZT進行譜分析的Matlab實現(xiàn),在Matlab中使用的線性調(diào)頻z變換函數(shù)為czt，其調(diào)用格式為 X= czt(x, M, W, A) 其中，x是待變換的時域信號x(n)，其長度為N，M是變換的長度，W確定變換的步長，A確定變換的起點。若M= N，A= 1，則CZT變成DFT。,47,5.8.1 用FFT進行譜分析的Matlab實現(xiàn),例5.1 設(shè)模擬信號 ，以 t= 0.01n (n=0: N-1) 進行取樣，試用fft函數(shù)對其做頻譜分析。N分別為：(1) N=45；(2) N=50；(3) N=

12、55；(2) N=60。,程序清單如下,%計算N=45的FFT并繪出其幅頻曲線 N=45;n=0:N-1;t=0.01*n; q=n*2*pi/N; x=2*sin(4*pi*t)+5*cos(8*pi*t); y=fft(x,N); figure(1) subplot(2,2,1) plot(q,abs(y) title(FFT N=45),48,例5.1程序清單,%計算N=50的FFT并繪出其幅頻曲線 N=50;n=0:N-1;t=0.01*n; q=n*2*pi/N; x=2*sin(4*pi*t)+5*cos(8*pi*t); y=fft(x,N); figure(1) subplot

13、(2,2,2) plot(q,abs(y) title(FFT N=50),49,%計算N=55的FFT并繪出其幅頻曲線 N=55;n=0:N-1;t=0.01*n; q=n*2*pi/N; x=2*sin(4*pi*t)+5*cos(8*pi*t); y=fft(x,N); figure(1) subplot(2,2,3) plot(q,abs(y) title(FFT N=55),50,%計算N=60的FFT并繪出其幅頻曲線 N=60;n=0:N-1;t=0.01*n; q=n*2*pi/N; x=2*sin(4*pi*t)+5*cos(8*pi*t); y=fft(x,N); figure(1) subplot(2,2,4) plot(q,abs(y) title(FFT N=60