1、1快速傅里葉變換（FFT）1.1 葉變換簡介快速傅里有限長序列可以通過離散傅里葉變換(DFT)將其頻域也離散化成有限長序列.但其計算量太大,很難實時地處理問題,因此引出了快速傅里葉變換(FFT). 1965 年，Cooley 和Tukey 提出了計算離散傅里葉變換（DFT）的快速算法，將DFT 的運算量減少了幾個數(shù)量級。從此，對快速傅里葉變換（ FFT）算法的研究便不斷深入，數(shù)字信號處理這門新興學(xué)科也隨FFT 的出現(xiàn)和發(fā)展而迅速發(fā)展。根據(jù)對序列分解與選取方法的不同而產(chǎn)生了FFT 的多種算法，基本算法是基-2DIT 和基-2DIF。FFT 在離散傅里葉反變換、線性卷積和線性相關(guān)等方面也有重要應(yīng)用

2、。快速傅氏變換（FFT），是離散傅氏變換的快速算法，它是根據(jù)離散傅氏變換的奇、偶、虛、實等特性，對離散傅立葉變換的算法進(jìn)行改進(jìn)獲得的。它對傅氏變換的理論并沒有新的發(fā)現(xiàn)，但是對于在計算機(jī)系統(tǒng)或者說數(shù)字系統(tǒng)中應(yīng)用離散傅立葉變換，可以說是進(jìn)了一大步。快速傅立葉變換作為一種數(shù)學(xué)方法，已經(jīng)廣泛地應(yīng)用在幾乎所有領(lǐng)域的頻譜分析中，而且經(jīng)久不衰，因為信號處理方法沒有先進(jìn)和落后之分，只有經(jīng)典和現(xiàn)代之別，在實際系統(tǒng)中用得最好的方法就是管用的方法。換句話說，信號處理方法與應(yīng)用背景和目的的貼近程度是衡量信號處理方法優(yōu)劣的唯一標(biāo)準(zhǔn)。FFT 是快速傅利葉變換(Fast FourierTransform 簡稱FFT)的英文

3、縮寫,它在當(dāng)今科技世界中的應(yīng)用相當(dāng)活躍,無論是在時間序列分析領(lǐng)域中,還是在我國剛剛興起的生物頻譜治療的研究與應(yīng)用中,都有著重要的作用。同時,它又是軟件實現(xiàn)數(shù)字濾波器的必備組成部分之一。FFT 算法的基本思想：利用DFT 系數(shù)的特性，合并DFT 運算中的某些項，把長序列的DFT>短序列的DFT，從而減少其運算量。FFT 算法分類: 時間抽選法DIT: Decimation-In-Time ； 頻率抽選法DIF:Decimation-In-Frequency。在此以按時間抽選的基-2FFT 算法為例進(jìn)行說明。1.2 算法原理N 為2 的整數(shù)冪的FFT 算法稱基-2FFT 算法。設(shè)N點序列，=

4、0,1,2,N-1按n的奇偶分成兩個長為N/2的序列 =0,1,2, -1由于則的DFT： 為偶數(shù) 為奇數(shù) = = = =0,1,2,，N-1 （1）其中：= （2）= （3）按式（1）計算得到的只是的前一半項數(shù)的結(jié)果，即（01）。由系數(shù)的周期性可推出的后一半值，即（N1）：. （=1）又因為，是以為周期的，可推出： （4）同理可得，= （5）由此可得N點對應(yīng)的DFT的計算式子： =0,1,2, -1 （6）上一式計算的前一半值，下一式計算的后一半值。因此只要求出01區(qū)間內(nèi)各個整數(shù)值所對應(yīng)的和，便可求出0到N-1點全部的值。明顯節(jié)約了運算量。式中因子在復(fù)數(shù)乘法中起一個旋轉(zhuǎn)的作用，稱為旋轉(zhuǎn)因子。

5、公式（6）的運算可以歸結(jié)為兩個復(fù)數(shù)a、b求得復(fù)數(shù)和。用信號流程圖的方法可以簡單的表示為如圖2所示。這樣的運算稱為蝶形運算，在FFT算法中占有核心的地位。圖2 時間抽選蝶形運算流圖 顯然，每個蝶形運算對應(yīng)于一次復(fù)數(shù)乘法和兩次復(fù)數(shù)加法運算。如果用DFT方法直接計算出和，共需次復(fù)數(shù)乘法運算，再作N/2次蝶形運算，又需N/2次復(fù)數(shù)乘法和N次復(fù)數(shù)加法。這樣算出N點DFT共需要（+N）/2次復(fù)數(shù)乘法和(/2)+N次復(fù)數(shù)加法。當(dāng)N較大時，同直接計算N點DFT所需的次復(fù)數(shù)乘加次數(shù)相比，幾乎減少了一半的計算量。 假設(shè)N/2還是偶數(shù)，則和這兩個N/2點的DTF計算，又分別可以通過計算N/4點DFT和蝶形運算而得到

6、。這時蝶形有兩組，每組N/4個，總數(shù)也是N/2個，所以也需要N/2次復(fù)數(shù)乘法和N次復(fù)數(shù)加法。如果N=，則分解過程一直可以進(jìn)行下去，共分解M次，到2點DFT為止。以上所述就是FFT算法的核心思想。可以看出，這種罪基本形式的快速傅里葉算法要求點數(shù)是2的冪次。當(dāng)序列長度不具有的形式時，可以補(bǔ)上一段0，使總長度為。1.3時間抽取過程的流圖表示現(xiàn)在就N=8的情況對算法結(jié)構(gòu)，即計算流圖作詳細(xì)的說明，由此可以直接得到關(guān)于N=的情況下得一些一般性結(jié)論。圖3是8點DFT歐諾個兩個4點DFT和4個蝶形運算實現(xiàn)的圖示，圖4中進(jìn)一步將4點DFT作類似的分解，最后將圖4中的2點DFT也畫成蝶形運算流圖，如圖5所示。整個

7、運算有3輪蝶形運算構(gòu)成，每一輪有4個蝶形，但它們的旋轉(zhuǎn)因子和循環(huán)方式各輪有所不同。 G(0)-1 4點DFT4點DFTG(1)G(2)G(3)H(0)H(1)H(2)H(3)-1-1-1圖3 8點DFT時間抽選分解2點DFT-1-1-1-1-12點DFT2點DFT2點DFT-1-1-1圖4 8點DFT時間抽選分解（續(xù)） -1-1-1-1-1-1-1-1-1-1-1-1圖5 反序輸入順序輸出時間抽選FFT算法流圖圖5那樣的結(jié)構(gòu)對任何N=的序列都是適用的。整個N點DFT的計算可用M=輪蝶形運算構(gòu)成，每輪有N/2個蝶形，因此總的運算就由這N個蝶形所構(gòu)成。因此總的復(fù)數(shù)乘法次數(shù)是N，復(fù)數(shù)加法次數(shù)是N。1

8、.4 運算量下面討論按時間抽取法的運算量，對于任意一個N=的DFT運算，都可以通過M次分解，最后分解成2點DFT運算的組合。這樣的M次分解，就構(gòu)成由到的M級運算過程，每一級運算都有N/2個蝶式運算完成。而每一蝶式運算有一次復(fù)數(shù)乘法和兩次復(fù)數(shù)加法。因此，每一級運算都需要N/2次復(fù)數(shù)乘法和N次復(fù)數(shù)加法。這樣M級運算共需要復(fù)數(shù)乘法： 復(fù)數(shù)加法： 為了比較FFT算法與直接DFT運算的運算量，現(xiàn)將二者在不同N值時的乘法次數(shù)列于表1中，并將直接計算DFT與FFT算法的計算量之比也列于表中。從表中可以看出，當(dāng)N較大時，時間抽取法要比直接計算塊一、二個數(shù)量級。表1 FFT算法與直接計算DFT算法的比較N 24

9、14.041644.0864125.416256328.03210248012.864409619221.41281638444836.6256262144102464.051210485762304113.81024419430411264204.81.5按時間抽取的FFT算法特點（1）原位運算（同址運算）由所述算法原理及圖5的N=8的信號流圖，F(xiàn)FT的每級（列）計算都是由N個復(fù)數(shù)據(jù)經(jīng)N/2個蝶形運算變成了另外N個復(fù)數(shù)據(jù)。每一個蝶形運算結(jié)構(gòu)完成下述基本迭代運算 （7） （8）式中表示第列迭代，、為數(shù)據(jù)所在的行數(shù)。上式的運算如圖6所示，由一次復(fù)數(shù)乘法和一次復(fù)數(shù)加減法組成。任何兩個節(jié)點和的節(jié)點變

10、量進(jìn)行蝶形運算后，其結(jié)果為下一列的和兩點的節(jié)點變量，而和其他節(jié)點變量無關(guān)。-1圖6 時間抽取蝶形運算結(jié)構(gòu)（2）輸入序列的序號及整序規(guī)律 由圖5可見，當(dāng)按原位進(jìn)行計算時，F(xiàn)FT輸出端的的次序正好是順序排列的。但這時輸入?yún)s不能按自然順序存入存儲單元中，而是按，的順序存入存儲單元，因而是亂序的。這就使得運算時取數(shù)據(jù)的地址編排“混亂無序”。但實際上是由規(guī)律可尋的，我們稱為倒位序。亂序的原因是由按時間抽取進(jìn)行的FFT運算的原理造成的，是由于輸入按標(biāo)號n的奇偶不斷分組造成的。如果你用二進(jìn)制數(shù)表示為（當(dāng)N=8=時，用三位二進(jìn)制符號表示），第一次分組n為偶數(shù)，在上半部分，n為奇數(shù)在下半部分，這可以觀察n的二進(jìn)

11、制的最低位，=0則序列值對應(yīng)于偶數(shù)抽樣，=1則序列值對應(yīng)于奇數(shù)抽樣。下一次則根據(jù)次最低位的0，1來分奇偶（而不管原來的序列的子序列是偶序列還是奇序列）。在實際運算中，直接將輸入數(shù)據(jù)按原位運算要求的“亂序”存放是很不方便的。因此總是先按自然順序?qū)⑤斎胄蛄写嫒氪鎯卧偻ㄟ^變址運算將自然順序變換成按時間抽取的FFT算法要求的順序。變址的過程可以用程序安排加以實現(xiàn)，稱為“整序”或“重排”。 整序規(guī)律：如果輸入序列的序號用二進(jìn)制數(shù)表示，如，若其反序二進(jìn)制用來表示，則原來在自然順序時應(yīng)該放的存儲單元，現(xiàn)在放著的是。表2列出了N=8時的自然順序二進(jìn)制數(shù)以及相應(yīng)的倒位序二進(jìn)制數(shù)。表2 順序和倒位序二進(jìn)制數(shù)

12、自然順序n二進(jìn)制數(shù)倒位序二進(jìn)制數(shù)倒位序順序0000000010011004201001023011110641000011510110156110011371111117將按自然順序存放在存儲單元中的數(shù)據(jù)，調(diào)換成FFT原位運算所要求的倒位序的變址功能，如圖8所示。當(dāng)n=時，數(shù)據(jù)不需要調(diào)換，當(dāng)n時，必須將原來存放數(shù)據(jù)的存單元內(nèi)調(diào)入數(shù)據(jù)，而將存放的單元內(nèi)調(diào)入。為了避免再次考慮前面已調(diào)換過得數(shù)據(jù)，保證調(diào)換只進(jìn)行一次，只需檢查一下是否比n小，假若比n小，則意味著在前面已經(jīng)與互相調(diào)換過。只有當(dāng)比n大時，才將原存放及存放的存儲單元的內(nèi)容互相調(diào)換。經(jīng)過上述變址運算以后，所得到的數(shù)據(jù)順序正是輸入所要求的順序。

13、與圖4的輸入順序?qū)φ眨浑y看出上述整序規(guī)律是正確的。存儲單元 A(1) A(3) A(5) 自然順序輸入 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)圖8 倒位序的變址處理（3）各類蝶形運算兩節(jié)點的“距離”及的變化規(guī)律以8點FFT為例。第一列蝶形運算只有一種類型：系數(shù)為=1，參加蝶形運算的兩個數(shù)據(jù)點的間距為1。第二列有兩種類型的蝶形運算：系數(shù)分別為，參加蝶形運算的兩個數(shù)據(jù)節(jié)點間的間距為2。第三列有4種類型的蝶形運算：系數(shù)分別為，參加蝶形運算的兩個數(shù)據(jù)節(jié)點間的間距為4。可見

14、，每列的蝶形運算類型比前1列增加1倍，參加蝶形運算的兩個數(shù)據(jù)節(jié)點間的間距也增大1倍。由此尅類推，對于N=點FFT，當(dāng)輸入為倒位序，輸出為正常順序時，其第m級運算，參與每個蝶形運算的數(shù)據(jù)節(jié)點間距為。針對具體的第m級蝶形運算，一個DIT蝶形運算的來那個節(jié)點“距離”為，因而式7和式8可寫成 （9） （10）中的可用兩種方法確定。 第1種方法：逆推法。由于蝶形運算最后一級適用了N/2中系數(shù)，分寫為，而前一列使用了它后面一列所用系數(shù)對應(yīng)的偶數(shù)序號的一半，依次類推，可以求出所有的需要適用的系數(shù)。 第2種方法：直接求解。首先把式9或式10中兩節(jié)點中的第1個節(jié)點標(biāo)號值，即值，表示成L位（N=）二進(jìn)制數(shù)。然后把

15、此二進(jìn)制數(shù)乘上，即將此L位二進(jìn)制數(shù)左移位（注意m是第m級的運算），把右邊空出的位置補(bǔ)零，此數(shù)就是所求的二進(jìn)制數(shù)。1.6 快速傅立葉變換的C語言實現(xiàn)由DIT基-2FFT算法原理及特點可知，F(xiàn)FT計算程序主要包括變址和M級遞推計算兩部分。 （1）變址（倒序）運算在實際運算中，一般總是按自然順序?qū)⑤斎胄蛄写嫒氪鎯卧虼藶閷崿F(xiàn)DIT基-2FFT首先必須將輸入序列按二進(jìn)制倒序數(shù)重排。由表2知，自然順序每增加1，相應(yīng)于二進(jìn)制順序數(shù)最低位加1，并向左進(jìn)位。而倒序數(shù)則是在M位二進(jìn)制數(shù)最高位加1，并向右進(jìn)位，用這種反向進(jìn)位加法可以從當(dāng)前任意倒序數(shù)值求得下一個倒序數(shù)值，完成倒序重排。反向進(jìn)位加法實現(xiàn)變址運算的

16、計算流程，如圖9所示，稱為雷德（Rader）算法。設(shè)A(I)表示存放原自然順序輸入數(shù)據(jù)的內(nèi)存單元，A(J)表示存放倒序順序數(shù)據(jù)的內(nèi)存單元。I,J都從0開始。若已知某倒序數(shù)J，要求下一個倒序數(shù)，先判斷J的最高位是否為“0”，這可與K=N/2相比較，因為N/2的二進(jìn)制描述中，最高位為1，其余位為0，因為如果KJ，則J的最高位為“0”，只要該位變?yōu)椤?”（令J與K=N/2相加即可），就得到了下一個倒序數(shù)。如果KJ，則J的最高位為“1”，可將其改為“0”（減去N/2即可）。然后再判斷次高位（與K=N/4比較），若其為“0”，則將它改為“1”（與N/4相加）,此外還需判斷再下一位（與K=N/8比較）。如

17、此依次進(jìn)行，直到遇到某位為“0”,這時把這個“0”改為“1”，就得到了下一個倒序數(shù)。求得新的倒序數(shù)J以后，當(dāng)IJ時，進(jìn)行 （2）M級遞推計算由DIT基-2FFT算法的特點，可得出M級遞推計算的流程圖，如圖10所示。整個M輪遞推過程由三個嵌套循環(huán)構(gòu)成。外層循環(huán)控制M輪的順序運算。內(nèi)層的兩個循環(huán)控制同一輪的各個蝶形運算，其中最內(nèi)一層循環(huán)控制同類（指有相同的旋轉(zhuǎn)因子）蝶形的運算，而中間一層循環(huán)則針對不同類的蝶形運算。I和IP代表的是一個蝶形運算的兩個節(jié)點。設(shè)N=，M表示第M輪，M=1,2.,L。在第M輪，有個不同的旋轉(zhuǎn)因子，它們?yōu)椋↘=0,1,., -1）。每輪中，同類旋轉(zhuǎn)因子對應(yīng)的蝶形有個。各蝶形

18、依次相距D=點。最內(nèi)層循環(huán)的循環(huán)變量為I，它用來控制同類蝶形運算。其步進(jìn)值為蝶形間距D=，同類蝶形中參加運算的兩個節(jié)點的間距為D1= 點。第二層循環(huán)的循環(huán)變量J用來控制計算不同類型的蝶形，J的步進(jìn)為1。最內(nèi)層循環(huán)完成每輪每輪的蝶形運算，而中間層（第二層）循環(huán)完成系數(shù)的運算。最外層循環(huán)的循環(huán)變量為M，它用來控制運算的輪數(shù)，M由1變到L，步進(jìn)值為1。當(dāng)L改變時，則D，D1,U都會改變。U，W，J為存放復(fù)數(shù)單元，其乘法為復(fù)數(shù)乘法，旋轉(zhuǎn)因子根據(jù)下面遞推公式計算： = (11) 輸入NNV2N/2NM1N-1JOIOIJ?TA(J)A(J)A(I)A(I)TKNV2K>J?JJ+K YII+1IN

19、M1YNJJ-KKK/2N N 圖9 雷德算法流程圖數(shù)據(jù)倒序Y 1MDD/2D11UExp(-j/D1)WOJJII+d1IPA(IP)*UTA(I)-TA(IP)A(I)+TA(I)I+DII>N-1?U*WUJ+IJJ>D-1?IML>M?N YNY N圖10 基-2DIT-FFT算法流程圖2.用FFT計算相關(guān)函數(shù) 相關(guān)概念很重要，互相關(guān)運算廣泛應(yīng)用于信號分析與統(tǒng)計分析，如通過相關(guān)函數(shù)峰值的檢測測量兩個信號的時延差等。 兩個長為N的實離散時間序列x（n）與y（n）的互相關(guān)函數(shù)定義為 則可以證明，rxy()的離散付里葉變換為 Rxy（k）=X*(k)Y(k) 0kN-1 其

20、中X(k)=DFTx(n),Y(k)=DFTy(n),Rxy(k)=DFTrxy() 證：將x(n)、y(n)的逆離散付里葉變換代入互相關(guān)函數(shù)定義式 因x（n）是實序列，所以x（n）=x*（n），得 因 得 證畢。 當(dāng)x（n）=y（n）時，得到x（n）的自相關(guān)函數(shù)為： 上面的推導(dǎo)表明，互相關(guān)和自相關(guān)函數(shù)的計算可利用FFT實現(xiàn)。由于離散付里葉變換隱含著周期性，所以用FFT計算離散相關(guān)函數(shù)也是對周期序列而言的。直接做N點FFT相當(dāng)于對兩個N點序列x(n)、y(n)作周期延拓，作相關(guān)后再取主值（類似圓周卷積）。而實際一般要求的是兩個有限長序列的線性相關(guān)，為避免混淆，需采用與圓周卷積求線性卷積相類似的

21、方法，先將序列延長補(bǔ)0后再用上述方法。 利用FFT求兩個有限長序列線性相關(guān)的步驟： 設(shè)x（n）長為N1，y（n）長為N2，求線性相關(guān)。 （1）為了使兩個有限長序列的線性相關(guān)可用其圓周相關(guān)代替而不產(chǎn)生混淆，選擇周期NN1+N2-1，且N=2m，以便使用FFT，將x（n），y（n）補(bǔ)零至長為N。即： （2）用FFT計算X(k)，Y(k)（k=0，1，N-1） （3）R（k）=X*（k）Y（k） （4）對R（k）作IFFT，得到r（n）（n=0，1，N-1）/*fft*/先創(chuàng)建兩個txt格式的文件，分別為x.txt和y.txt，用來存放互相關(guān)變換的兩組數(shù)據(jù)；運行完產(chǎn)生result.txt 和resu

22、lt2.txt和result.txt，用來保存結(jié)果。#include <stdio.h> #include <math.h> #include <stdlib.h>#include <iostream>#include <fstream>using namespace std;#define PI 3.1415926535897932384626433832795028841971 #define N 128 typedef struct double real;/*實部*/ double img;/*虛部*/ complex;com

23、plex xN, yN,RN, *W;/*輸出序列的值*/int size_x=0,size_y=0,length=0;/序列長度void add(complex a,complex b,complex *c) c->real=a.real+b.real; c->img=a.img+b.img; void sub(complex a,complex b,complex *c) c->real=a.real-b.real; c->img=a.img-b.img; void mul(complex a,complex b,complex *c) c->real=a.r

24、eal*b.real - a.img*b.img; c->img=a.real*b.img + a.img*b.real; void divi(complex a,complex b,complex *c) c->real=( a.real*b.real+a.img*b.img )/( b.real*b.real+b.img*b.img); c->img=( a.img*b.real-a.real*b.img)/(b.real*b.real+b.img*b.img); void initW(int size_x)/核運算 int i; W=(complex *)malloc(

25、sizeof(complex) * size_x); for(i=0;i<size_x;i+) Wi.real=cos(2*PI/size_x*i); Wi.img=-1*sin(2*PI/size_x*i); void change()/變址運算 complex temp; unsigned short i=0,j=0,k=0; double t; for(i=0;i<size_x;i+) k=i;j=0; t=(log(size_x)/log(2); while( (t-)>0 ) j=j<<1; j|=(k & 1); k=k>>1; i

26、f(j>i) temp=xi; xi=xj; xj=temp; void fftx()/快速傅里葉函數(shù)int i=0,j=0,k=0,l=0; complex up,down,product; change(); for(i=0;i<log(size_x)/log(2);i+) /*一級蝶形運算*/ l=1<<i; for(j=0;j<size_x;j+=2*l) /*一組蝶形運算*/ for(k=0;k<l;k+) /*一個蝶形運算*/ mul(xj+k+l,Wsize_x*k/2/l,&product); add(xj+k,product,&am

27、p;up); sub(xj+k,product,&down); xj+k=up; xj+k+l=down; void ffty()/快速傅里葉函數(shù)int i=0,j=0,k=0,l=0; complex up,down,product; change(); for(i=0;i<log(size_y)/log(2);i+) /*一級蝶形運算*/ l=1<<i; for(j=0;j<size_y;j+=2*l) /*一組蝶形運算*/ for(k=0;k<l;k+) /*一個蝶形運算*/ mul(yj+k+l,Wsize_x*k/2/l,&product

28、); add(yj+k,product,&up); sub(yj+k,product,&down); yj+k=up; yj+k+l=down; void ifft() /傅里葉逆變換 length=size_x+size_y-1; int i=0,j=0,k=0,l=length; complex up,down; for(i=0;i<(int)(log(length)/log(2);i+) /*一級蝶形運算*/ l/=2; for(j=0;j<length;j+=2*l) /*一組蝶形運算*/ for(k=0;k<l;k+) /*一個蝶形運算*/ add(

29、Rj+k,Rj+k+l,&up); up.real/=2;up.img/=2; sub(Rj+k,Rj+k+l,&down); down.real/=2;down.img/=2; divi(down,Wlength*k/2/l,&down); Rj+k=up; Rj+k+l=down; change(); void outputone() /*輸出x結(jié)果*/ int i; printf("nx傅里葉變換結(jié)果n"); for(i=0;i<size_x;i+) if(i%4=0&&i!=0) printf("n"

30、); printf("%.4f",xi.real); if(xi.img>=0.0001) printf("+%.4fj ",xi.img); else if(fabs(xi.img)<0.0001) printf("+0.0000j "); else printf("%.4fj ",xi.img); printf("n"); void outputtwo() /*輸出x結(jié)果*/ int i; printf("ny傅里葉變換結(jié)果n"); for(i=0;i<

31、size_y;i+) if(i%4=0&&i!=0) printf("n"); printf("%.4f",yi.real); if(yi.img>=0.0001) printf("+%.4fj ",yi.img); else if(fabs(yi.img)<0.0001) printf("+0.0000j "); else printf("%.4fj ",yi.img); printf("n"); void outputthree() /*輸出x

32、y互相關(guān)結(jié)果*/ int i; printf("nxy互相關(guān)變換結(jié)果n"); for(i=0;i<length;i+) if(i%4=0&&i!=0) printf("n"); printf("%.4f",Ri.real); if(Ri.img>=0.0001) printf("+%.4fj ",Ri.img); else if(fabs(Ri.img)<0.0001) printf("+0.0000j "); else printf("%.4fj &

33、quot;,Ri.img); printf("n"); void saveone()/保存x傅里葉變換結(jié)果int i;ofstream outfile("result1.txt",ios:out);if(!outfile)cerr<<"open result1.txt erro!"<<endl;exit(1);outfile<<"x傅里葉變換結(jié)果:"<<endl;for(i=0;i<size_x;i+)outfile<<xi.real;if(xi.i

34、mg>=0.0001) outfile<<"+"<<xi.img<<"j"<<" " else if(fabs(xi.img)<0.0001) outfile<<"+"<<"0.0000"<<"j"<<" " else outfile<<xi.img<<"j"<<" " pr

35、intf("n");outfile.close();void savetwo()/保存y傅里葉變換結(jié)果int i;ofstream outfile("result2.txt",ios:out);if(!outfile)cerr<<"open result2.txt erro!"<<endl;exit(1);outfile<<"y傅里葉變換結(jié)果:"<<endl;for(i=0;i<size_y;i+)outfile<<yi.real;if(yi.img

36、>=0.0001) outfile<<"+"<<yi.img<<"j"<<" " else if(fabs(yi.img)<0.0001) outfile<<"+"<<"0.0000"<<"j"<<" " else outfile<<yi.img<<"j"<<" " prin

37、tf("n");outfile.close();void savethree()/保存xy互相關(guān)變換結(jié)果int i;ofstream outfile("result.txt",ios:out);if(!outfile)cerr<<"open result.txt erro!"<<endl;exit(1);outfile<<"xy互相關(guān)變換結(jié)果:"<<endl;for(i=0;i<length;i+)outfile<<Ri.real;if(Ri.img>=0.0001) outfile<&