【正文】 出，一個N點DFT已分解為2個N/2點的DFT，它們按式(52)又組合成一新的DFT。但是x1(r)、x2(r)以及X1(k)、X2(k)都是N/2點的序列，即滿足r,k=0,1,…, N/21。而X(k)卻有N點，用式(52)計算可以得到的只是X(k)的前一半項數(shù)的結(jié)果，要用X1(k)、X2(k)來表達(dá)全部X(k)值，還必須應(yīng)用系數(shù)的周期性，即這樣可以得到 (55)同理可得 (56)式(55) 和(56)說明了后半部分k值()所對應(yīng)的X1(k)、X2(k)分別等于前半部分k值()所對應(yīng)的X1(k)、X2(k)。再考慮到的對稱性 (57)把式(55)~ (57)代入到式(52)中，就可將X(k)表達(dá)為下面所列的前后兩個部分： 前半部分X(k) (58) 后半部分X(k) (59)這樣，只要求出0~(N/21)區(qū)間的所有X1(k)和X2(k)的值，即可求出0~N1區(qū)間內(nèi)的所有X(k)值，這就大大節(jié)省了運算量。式(58)和(59)的運算可以用圖51的蝶形信號流圖符號表示。當(dāng)支路上沒有標(biāo)出系數(shù)時，則該支路的傳輸系數(shù)為1。圖51 時間抽取法蝶形信號流圖符號采用這種方法，可將上面討論的分解過程示意于圖52中，此圖表示N=23=8的情況，其中輸出值X(0)~ X(3)是由公式(58)給出的，而X(4)~ X(7)是由公式(59)給出的給出的。圖52按時間抽取，將一個N點DFT分解為兩個N/2點DFT一個N點DFT分解為兩個N/2點DFT，每一個N/2點DFT只需(N/2) 2=N2/4次復(fù)數(shù)乘法，N/2(N/21)次復(fù)數(shù)加法。兩個N/2點DFT共需2(N/2)2=N2/2次復(fù)數(shù)乘法和N(N/21)次復(fù)數(shù)加法。此外，把兩個N/2點DFT合成為N點DFT時，有N/2個蝶形運算，還需要N/2次復(fù)數(shù)乘法及2*N/2=N次復(fù)數(shù)加法。因而通過第一步分解后，總共需要(N2/2)+(N/2)=N(N+1)/2≈N2/2次復(fù)數(shù)乘法和N(N/21)+N=N2/2次復(fù)數(shù)加法。由此可見，通過這樣分解后運算工作量差不多節(jié)省了一半。既然如此，由于N=2L，因而N/2仍為偶數(shù)，可以進一步把每個N/2點子序列再按其奇偶部分分解為2個N/4點的子序列。 (510) ，并且，式中 (511) (512)圖53示意出了N=8時，將一個N/2點DFT分解成兩個N/4點DFT，由此兩個N/4點DFT組合成了一個N/2點DFT的信號流圖。圖53 由2個N/4點DFT組合成1個N/2點DFTX2(k)也可以進行同樣的分解式中 (513) (514)將系數(shù)統(tǒng)一為，則一個N=8點DFT就可以分解為4個N/4=2點DFT，這樣就可以得到如圖54的信號流圖。圖54 按時間抽取，將一個N點DFT分解為4個N/4點DFT(N=8)根據(jù)上面同樣的分析可知，利用4個N/4點的DFT及2個蝶形組合運算來計算N點DFT，比只用一次分解蝶形組合方式的計算量又減少了大約一半。最后，剩下的是2點DFT，對于此例N=8，就是四個N/4=2點DFT，其輸出為X3(k)、X4(k)、X5(k)、X6(k)，k=0,1，由式(511)~ (514)可以將四個式子都計算出來，例如由式(512)計算得：上式中，故上式無需乘法。類似的可以求出X3(k)、X5(k)、X6(k)，這些2點DFT都可以用一個蝶形表示。由此可得出一個按時間抽取運算的完整的8點DFT信號流圖，如圖55所示。圖55 N=8按時間抽取法FFT信號流圖 按時間抽取的FFT算法與直接計算DFT運算量的比較由按時間抽取法FFT的流圖可見，當(dāng)N=2L時，共有L級蝶形，每級都由N/2個蝶形運算組成，每個蝶形需要一次復(fù)乘、二次復(fù)加，因而每級運算都需N/2次復(fù)乘和N次復(fù)加，這樣L級運算總共需要 復(fù)數(shù)乘法：復(fù)數(shù)加法：實際計算量與這個數(shù)字稍有不同，因為，所以這幾個系數(shù)是不需要乘法運算的，但是這些情況是在直接計算DFT中存在。此外，當(dāng)N較大時，這些特例相對而言就很少。所以，為了統(tǒng)一作比較起見，下面都不考慮這些特例。由于計算機上乘法運算所需的時間比加法運算所需的時間多得多，故以乘法為例，直接DFT復(fù)數(shù)乘法次數(shù)是N2，F(xiàn)FT復(fù)數(shù)乘法次數(shù)是(N/2)log2N。 直接計算DFT與FFT算法的計算量之比為：我們可以將此函數(shù)做出曲線圖，如圖56所示，這樣我們就可以很直觀的看出按時間抽取的FFT算法與直接計算DFT運算量的比較。圖56 FFT算法與直接計算DFT所需乘法次數(shù)的比較曲線例如，當(dāng)N=2048時。當(dāng)點數(shù)N越大時，F(xiàn)FT的優(yōu)勢更為明顯。 算法的C++實現(xiàn)我們從上面的分析可以看出，這種運算是很有規(guī)律的，其每級(每列)計算都是由N/2個蝶形運算構(gòu)成的，每一個蝶形結(jié)構(gòu)完成下述基本迭代運算： (515)式中，m表示第m列迭代，k, j為數(shù)據(jù)所在行數(shù)。式(515)的蝶形運算由一次復(fù)乘和兩次復(fù)加(減)組成。某一列的任何兩個節(jié)點k和j的節(jié)點變量進行蝶形運算后，得到結(jié)果為下一列k, j兩節(jié)點的節(jié)點變量，而和其他節(jié)點變量無關(guān)，因而可以采用原位運算[9]。即某一列的N個數(shù)據(jù)送到存儲器后，經(jīng)蝶形運算，其結(jié)果為下一列數(shù)據(jù)，它們以蝶形為單位仍存儲在這同一組存儲器中，直到最后輸出，中間無需其他存儲器。也就是蝶形的兩個輸出值仍放回蝶形的兩個輸入所在的存儲器中。每列的N/2個蝶形運算全部完成后，再開始下一列的蝶形運算。這樣存儲器數(shù)據(jù)只需N個存儲單元。下一級的運算仍采用這種原位方式，只不過進入蝶形結(jié)的組合關(guān)系有所不同。這種原位運算結(jié)構(gòu)可以節(jié)省存儲單元，降低設(shè)備成本。但是，我們從圖55中很容易看出，由于奇偶的不斷分組，最后輸入的序列順序不是順位輸入的，而是倒位序[9]的。下面首先用遞歸的形式進行算法的描述，由于遞歸方法沒有充分利用DIT基2算法的優(yōu)點——原位計算，因此遞歸形式只是為了描述的清晰，并不是使用的程序。void dit2rec(plexdouble*InData,plexdouble*OutData,int length,int flag){plexdouble*EvenData=new plexdouble(length/2)。//偶序列組成一個子序列plexdouble*OddData=new plexdouble(length/2)。//奇序列組成一個子序列plexdouble*EvenResult=new plexdouble(length/2)。//偶序列的變換結(jié)果plexdouble*OddResult=new plexdouble(length/2)。//奇序列的變換結(jié)果int i,j。if(length==1)//單個數(shù)的變換，不變，也是遞歸結(jié)束的條件{ if (flag==1) OutData[0]=InData[0]/length。 return。 }for(i=0。ilength/2。i++)//形成奇序列和偶序列{ EvenData[i]=InData[2*i]。 OddData[i]=InData[2*i+1]。}dit2rec(EvenData,EvenResult,length/2,sign)。//對奇偶序列分別進行變換dit2rec(OddData,EvenResult,length/2,sign)。for(i=0。ilength/2。i++)//由奇偶序列的變換進行組合{OutData[i]=EvenData+OddData*plexdouble(cos(2*pi*i/length),sin(flag*2*pi*i/length))。OutData[i+length/2]=EvenDataOddData*plexdouble(cos(2*pi*i/length),sin(flag*2*pi*i/length))。}delete[] EvenData,OddData,EvenResult,OddResult。 return。}非遞歸實現(xiàn)如下(現(xiàn)不考慮輸入的倒序數(shù)問題)：//data為數(shù)據(jù)指針，返回值相同，Log2N=log2(length）,flag=1時為正變換，flag=1時為逆變換。 plexdouble* dit2(plexdouble*Data,int Log2N,int flag){ int i,j,k。 int length=1Log2N。 plexdouble wn,temp。 int step。 int index0,index1。 for(i=1。i=Log2N。i++) { step=1i。 for(j=0。jstep/2。j++) { wn=plexdouble(cos(2*pi/step*j),sin(flag*2*pi/step*j))。 for(k=0。klength/step。k++) { index0=k*step+j。 index1=k*step+step/2+j。 temp=Data[index0]。 Data[index0]=temp+wn*Data[index1]。 Data[index1]=tempwn*Data[index1]。 } } } if(flag==1)//如果為反變換，要除以序列的長度 for(int i=0。ilength。i++) Data[i]/=length。 return Data。}當(dāng)i=1時，也就是第一次循環(huán)并沒有必要進行復(fù)數(shù)運算，因為j只能取1，wn為實數(shù)，這個時間可以節(jié)省。因此可以改進為：void dit2(plexdouble*Data,int Log2N,int flag){ int i,j,k,step,length。 plexdouble wn,temp,deltawn。 length=1Log2N。 for(i=0。ilength。i+=2) { temp=Data[i]。 Data[i]=Data[i]+Data[i+1]。 Data[i+1]=tempData[i+1]。 } for(i=2。i=Log2N。i++) {wn=1。step=1i。deltawn=plexdouble(cos(2*pi/step),sin(flag*2*pi/step))。； for(j=0。jstep/2。j++) { for(i=0。ilength/step。i++) { temp=Data[i*step+step/2+j]*wn。 Data[i*step+step/2+j]=Data[i*step+j]temp。 Data[i*step+j]=Data[i*step+j]+temp。 } wn=wn*deltawn。 } } if(flag==1) for(i=0。ilength。i++) Data[i]/=length。}6總結(jié)和展望本文從Fourier變換的基本理論出發(fā)，對其定義和基本性質(zhì)做了相關(guān)的總結(jié)和介紹。接著重點對Fourier的幾種變形形式做了深入的探討和研究，從時域和頻域的信號處理角度來分析了各種不同F(xiàn)ourier分析形式的異同，從定義和性質(zhì)兩個方面進行討論。在此理論基礎(chǔ)之上，對現(xiàn)實的Fourier變換實例展開研究。主要對離散Fourier變換(DFT)和Fourier變換的變形——分?jǐn)?shù)階Fourier變換(FRFT)兩個應(yīng)用最廣泛的案例做了分析，應(yīng)用的領(lǐng)域則采取了以信號處理為主、物理光學(xué)為輔的方式。結(jié)合Fourier變換的各種性質(zhì)，通過對其應(yīng)用的方法研究，對其各種性質(zhì)的應(yīng)用方法有了更深入和更直觀的理解。最后，在這些理解的基礎(chǔ)之上，本文對離散Fourier變換的快速計算方法——快速Fourier變換(FFT)進行了算法研究和C++語言的具體實現(xiàn)，從而在計算機仿真模擬這方面有了初步的了解和體驗。本文的研究也只是粗淺和局限性的。并且在仿真模擬這方面，還有著有很大的研究空間，例如還可以用MATLAB對各種應(yīng)用案例進行實現(xiàn)，這比C++的實現(xiàn)更具直觀性。但是通過此次研究，對各種基本方法都有了一定的了解和體會。帶著這樣的基礎(chǔ)，我們還可以就其他的Fourier變換形式以及更多的領(lǐng)域做更深入的探討和研究，從而得到Fourier變換更多更有用的應(yīng)用方法。參考文獻(xiàn)[1] Lin C. C. amp。 Segle L. A. Mathematics Applied to Deterministic Problems in the Natural Sciences [J]. Macmillan Inc., New York, 1974, 28(7): 1560.[2] 劉粵鉗, 姚紅紅. Fourier變換的改進形式——分?jǐn)?shù)傅立葉變換[J]. 安慶師范學(xué)院學(xué)報(自然科學(xué)版), 2005, 11(3): 59.[3] . Cultural malpractice