fft演算法的基本思想
Ⅰ 一維實序列的快速傅里葉變換(FFT)
通過前面的分析,我們認識到傅里葉變換本身是復數運算,地球物理獲取的數據大多數是實數,對於實數的變換原則上可直接套用復序列的FFT演算法,但那樣是把實數序列當作虛部為零的復數對待,顯然需要存儲虛部的零並進行無功的運算,既浪費了一倍的計算內存,又降低了約一半的運算速度。
為了不浪費不可不設的虛部內存和必然出現的復數運算,可否將一個實序列分為兩個子實序列,分別作為實部與虛部構成一個復數序列,然後用復序列的FFT演算法求其頻譜,對合成的復序列頻譜進行分離和加工得到原實序列的頻譜呢?答案是肯定的,實現這一過程思路就是實序列FFT演算法的基本思想。
1.實序列的傅里葉變換性質
對於一個N個樣本的實序列x(k),其頻譜為X(j),用Xr(j)和Xi(j)表示X(j)的實部和虛部, 表示X(j)的共軛,則
證明:已知 則
地球物理數據處理基礎
上式兩端取共軛,並注意到x(k)是實序列,則
地球物理數據處理基礎
這就是實序列的傅里葉變換具有復共軛性。
其同樣具有周期性,即
地球物理數據處理基礎
2.一維實序列的FFT演算法
(1)同時計算兩個實序列的FFT演算法
已知兩個實序列h(k),g(k)(k=0,1,…,N-1),例如重磁異常平面數據中的兩條剖面,或地震勘探中的兩道地震記錄,可以人為地構成一個復序列:
地球物理數據處理基礎
設h(k)的頻譜為H(j)=Hr(j)+iHi(j)
g(k)的頻譜為G(j)=Gr(j)+iGi(j)
y(k)的頻譜為Y(j)=Yr(j)+i Yi(j)
利用上節的復序列FFT演算法,求得Y(j),即Yr(j)和Yi(j)已知,來尋找Hr(j),Hi(j),Gr(j),Gi(j)與Yr(j),Yi(j)之間的關系。
對式(8-22)作傅里葉變換:
地球物理數據處理基礎
由於H(j),G(j)本身是復序列,所以不能僅從上式分離出H(j)和G(j)。應用Y(j)的周期性,容易得到
Y(N-j)=H(-j)+iG(-j)
上式取共軛:
地球物理數據處理基礎
由於h(k),g(k)為實序列,對上式右端應用復共軛定理,得
地球物理數據處理基礎
對式(8-23)展開,得
地球物理數據處理基礎
對式(8-24)展開,並應用共軛關系,得
地球物理數據處理基礎
把式(8-25)和式(8-26)與Y(j)=Yr(j)+iYi(j)進行對比,有
地球物理數據處理基礎
整理得
地球物理數據處理基礎
因此,對於兩個實序列,通過構造一個復序列,應用復序列的FFT演算法和式(8-28)的分離加工,即可得到兩個實序列的頻譜。
(2)計算2 N個數據點的實序列FFT演算法
設有2N點的實序列u(k)(k=0,1,…,2N-1),首先按k的偶、奇分成兩個子實序列,並構成復序列,即
地球物理數據處理基礎
通過調用復序列FFT演算法,求得y(k)的頻譜為Y(j)。另記h(k),g(k)的頻譜為H(j)和G(j)。
利用前面式(8-23)和式(8-24),容易求得
地球物理數據處理基礎
下面分析用H(j),G(j)形成u(k)頻譜的問題。記u(k)(k=0,1,…,2 N-1)的頻譜為V(j),分析V(j),H(j),G(j)之間的關系,根據定義
地球物理數據處理基礎
利用式(8-31)和式(8-34)可換算出u(k)的前N個頻譜V(j)(j=0,1,…,N-1),還要設法求u(k)的後N個頻譜V(N+j)(j=0,1,…,N-1)。利用實序列其頻譜的復共軛和周期性:
(1)H(N)=H(0),G(N)=G(0),WN1=-1,得
地球物理數據處理基礎
(2)由於u(k)(k=0,1,…,2N-1)是實序列,同樣利用實序列其頻譜的復共軛和周期性,用已求出的前N個頻譜V(j)表示出後面的N-1個頻譜V(N+j):
地球物理數據處理基礎
由於0<2N-j<N,所以可從V(j)(j=0,1,…,N-1)中選出V(2N-j)(j=N+1,N+2,…,2 N-1),並直接取其共軛 即可得到V(N+1)~V(2 N-1),從而完成整個實序列頻譜的計算。
總結以上敘述,一維實序列u(k)(k=0,1,…,2N-1)的FFT計算編程步驟如下:
(1)按偶、奇拆分實序列u(k),並構造復序列:
地球物理數據處理基礎
(2)調用復序列的FFT計算y(k)的頻譜Y(j)(j=0,1,…,N-1);
(3)用下式計算形成h(k),g(k)的頻譜H(j)和G(j);
地球物理數據處理基礎
(4)用下式換算實序列u(k)的頻譜V(j)(j=0,1,…,2 N-1):
地球物理數據處理基礎
[例3]求實序列樣本u(k)={1,2,1,1,3,2,1,2}(k=0,1,…,7)的頻譜。
解:按偶、奇拆分實序列u(k),按式(8-37)構造復序列c(j)(j=0,1,2,3),即
c(0)=1+2i; c(1)=1+i; c(2)=3+2i; c(3)=1+2i。
(1)調用復序列FFT求c(j)(j=0,1,2,3)的頻譜Z(k)(k=0,1,2,3),得
Z(0)=6+7i; Z(1)=-3; Z(2)=2+i; Z(3)=-1。
地球物理數據處理基礎
(3)運用公式(8-38)計算H(j),G(j):
地球物理數據處理基礎
(4)根據式(8-39)求出u(k)(k=0,1,…,7)的8個頻譜V(j)(j=0,1,…,7):
地球物理數據處理基礎
地球物理數據處理基礎
由上例可見,完成全部2 N個實序列頻譜的計算只需做N次FFT計算,相比直接用復序列的FFT演算法節省了約一半的計算量。
Ⅱ FFT原理的FFT基本原理
FFT是一種DFT的高效演算法,稱為快速傅立葉變換(fast Fourier transform)。FFT演算法可分為按時間抽取演算法和按頻率抽取演算法,先簡要介紹FFT的基本原理。從DFT運算開始,說明FFT的基本原理。
DFT的運算為:
式中
由這種方法計算DFT對於X(K)的每個K值,需要進行4N次實數相乘和(4N-2)次相加,對於N個k值,共需N*N乘和N(4N-2)次實數相加。改進DFT演算法,減小它的運算量,利用DFT中
的周期性和對稱性,使整個DFT的計算變成一系列迭代運算,可大幅度提高運算過程和運算量,這就是FFT的基本思想。
FFT基本上可分為兩類,時間抽取法和頻率抽取法,而一般的時間抽取法和頻率抽取法只能處理長度N=2^M的情況,另外還有組合數基四FFT來處理一般長度的FFT 設N點序列x(n),,將x(n)按奇偶分組,公式如下圖
改寫為:
一個N點DFT分解為兩個 N/2點的DFT,繼續分解,迭代下去,其運算量約為
其演算法有如下規律
兩個4點組成的8點DFT
四個2點組成的8點DFT
按時間抽取的8點DFT
原位計算
當數據輸入到存儲器中以後,每一級運算的結果仍然儲存在同一組存儲器中,直到最後輸出,中間無需其它存儲器
序數重排
對按時間抽取FFT的原位運算結構,當運算完畢時,這種結構存儲單元A(1)、A(2),…,A(8)中正好順序存放著X(0),X(1),X(2),…,X(7),因此可直接按順序輸出,但這種原位運算的輸入x(n)卻不能按這種自然順序存入存儲單元中,而是按X(0),X(4),X(2),X(6),…,X(7)的順序存入存儲單元,這種順序看起來相當雜亂,然而它也是有規律的。當用二進製表示這個順序時,它正好是「碼位倒置」的順序。
蝶形類型隨迭代次數成倍增加
每次迭代的蝶形類型比上一次蝶代增加一倍,數據點間隔也增大一倍 頻率抽取2FFT演算法是按頻率進行抽取的演算法。
設N=2^M,將x(n)按前後兩部分進行分解,
按K的奇偶分為兩組,即
得到兩個N/2 點的DFT運算。如此分解,並迭代,總的計算量和時間抽取(DIT)基2FFT演算法相同。
演算法規律如下:
蝶形結構和時間抽取不一樣但是蝶形個數一樣,同樣具有原位計算規律,其迭代次數成倍減小 時,可採取補零使其成為
,或者先分解為兩個p,q的序列,其中p*q=N,然後進行計算。 前面介紹,採用FFT演算法可以很快算出全部N點DFT值,即z變換X(z)在z平面單位圓上的全部等間隔取樣值。實際中也許①不需要計算整個單位圓上z變換的取樣,如對於窄帶信號,只需要對信號所在的一段頻帶進行分析,這時希望頻譜的采樣集中在這一頻帶內,以獲得較高的解析度,而頻帶以外的部分可不考慮,②或者對其它圍線上的z變換取樣感興趣,例如語音信號處理中,需要知道z變換的極點所在頻率,如極點位置離單位圓較遠,則其單位圓上的頻譜就很平滑,這時很難從中識別出極點所在的頻率,如果采樣不是沿單位圓而是沿一條接近這些極點的弧線進行,則在極點所在頻率上的頻譜將出現明顯的尖峰,由此可較准確地測定極點頻率。③或者要求能有效地計算當N是素數時序列的DFT,因此提高DFT計算的靈活性非常有意義。
螺旋線采樣是一種適合於這種需要的變換,且可以採用FFT來快速計算,這種變換也稱作Chirp-z變換。
Ⅲ FFT的公式是什麼和演算法是怎樣實現
二維FFT相當於對行和列分別進行一維FFT運算。具體的實現辦法如下:
先對各行逐一進行一維FFT,然後再對變換後的新矩陣的各列逐一進行一維FFT。相應的偽代碼如下所示:
for (int i=0; i<M; i++)
FFT_1D(ROW[i],N);
for (int j=0; j<N; j++)
FFT_1D(COL[j],M);
其中,ROW[i]表示矩陣的第i行。注意這只是一個簡單的記法,並不能完全照抄。還需要通過一些語句來生成各行的數據。同理,COL[i]是對矩陣的第i列的一種簡單表示方法。
所以,關鍵是一維FFT演算法的實現。下面討論一維FFT的演算法原理。
【1D-FFT的演算法實現】
設序列h(n)長度為N,將其按下標的奇偶性分成兩組,即he和ho序列,它們的長度都是N/2。這樣,可以將h(n)的FFT計算公式改寫如下 :
(A)
由於
所以,(A)式可以改寫成下面的形式:
按照FFT的定義,上面的式子實際上是:
其中,k的取值范圍是 0~N-1。
我們注意到He(k)和Ho(k)是N/2點的DFT,其周期是N/2。因此,H(k)DFT的前N/2點和後N/2點都可以用He(k)和Ho(k)來表示
Ⅳ FFT演算法分幾種
FFT演算法分析FFT演算法的基本原理是把長序列的DFT逐次分解為較短序列的DFT。按照抽取方式的不同可分為DIT-FFT(按時間抽取)和DIF-FFT(按頻率抽取)演算法。按照蝶形運算的構成不同可分為基2、基4、基8以及任意因子(2n,n為大於1的整數),基2、基4演算法較為常用。 網上有幫助文檔: http://www.5doc.com/doc/123035(右上角有點擊下載)
Ⅳ 什麼是FFT演算法DSP是什麼
FFT是快速傅里葉變換( Fast Fourier Transform )
DSP是數字信號處理 ( Digital Signal Processing )