傅里葉變換的本質(zhì)傅里葉變換的公式為F(o)=jf(t)e-網(wǎng)出-s可以把傅里葉變換也成另外一種形式：F(o)=—Yf(t),ejot>2?？梢钥闯觯道锶~變換的本質(zhì)是內(nèi)積，三角函數(shù)是完備的正交函數(shù)集，不同頻率的三角函數(shù)的之間的內(nèi)積為0，只有頻率相等的三角函數(shù)做內(nèi)積時，才不為0。YejQ/,e02A=jejq-q2)tdt=2兀8(?！?。)下面從公式解釋下傅里葉變換的意義因為傅里葉變換的本質(zhì)是內(nèi)積，所以f(t)和ejo求內(nèi)積的時候，只有f(t)中頻率為O的分量才會有內(nèi)積的結(jié)果，其余分量的內(nèi)積為0?？梢岳斫鉃閒(t)在ejot上的投影，積分值是時間從負無窮到正無窮的積分，就是把信號每個時間在o的分量疊加起來，可以理解為f(t)在ejot上的投影的疊加，疊加的結(jié)果就是頻率為o的分量，也就形成了頻譜。傅里葉逆變換的公式為1+sf(t)= jF(o)ejodo2?！猻下面從公式分析下傅里葉逆變換的意義傅里葉逆變換就是傅里葉變換的逆過程，在F(o)和ej求內(nèi)積的時候，F(xiàn)(o)只有t時刻的分量內(nèi)積才會有結(jié)果，其余時間分量內(nèi)積結(jié)果為0，同樣積分值是頻率從負無窮到正無窮的積分，就是把信號在每個頻率在t時刻上的分量疊加起來，疊加的結(jié)果就是f(t)在t時刻的值，這就回到了我們觀察信號最初的時域。離散付立葉變換的理解FFT是離散傅立葉變換的快速算法，可以將一個信號變換到頻域。有些信號在時域上是很難看出什么特征的，但是如果變換到頻域之后，就很容易看出特征了。這就是很多信號分析采^FT變換的原因。另外，F(xiàn)FT可以將一個信號的頻譜提取出來，這在頻譜分析方面也是經(jīng)常用的。雖然很多人都知道FFT是什么，可以用來做什么，怎么去做，但是卻不知^FT之后的結(jié)果是什意思、如何決定要使用多少點來做FFT?，F(xiàn)在就根據(jù)實際經(jīng)驗來說說FFT結(jié)果的具體物理意義。一個模擬信號，經(jīng)過ADC采樣之后，就變成了數(shù)字信號。采樣定理告訴我們，采樣頻率要大于信號頻率的兩倍，這些我就不在此啰嗦了。采樣得到的數(shù)字信號，就可以做FFT變換了。N個采樣點，經(jīng)過FFT之后，就可以得到N個點的FFT結(jié)果。為了方便進行FFT運算，通常N取2的整數(shù)次方。假設(shè)采樣頻率為Fs，信號頻率F，采樣點數(shù)為N。那么FFT之后結(jié)果就是一個為N點的復(fù)數(shù)。每一個點就對應(yīng)著一個頻率點。這個點的模值，就是該頻率值下的幅度特性。具體跟原始信號的幅度有什么關(guān)系呢？假設(shè)原始信號的峰值為A，那么FFT的結(jié)果的每個點(除了第一個點直流分量之外)的模值就是A的N/2倍。而第一個點就是直流分量，它的模值就是直流分量的N倍。而每個點的相位呢，就是在該頻率下的信號的相位。第一個點表示直流分量(即0日2)，而最后一個點N的再下一個點(實際上這個點是不存在的，這里是假設(shè)的第N+1個點，也可以看做是將第一個點分做兩半分，另一半移到最后)則表示采樣頻率Fs，這中間被N-1個點平均分成N等份，每個點的頻率依次增加。例如某點n所表示的頻率為：Fn=(n-1)*Fs/N。由上面的公式可以看出，F(xiàn)n所能分辨到頻率為為Fs/N，如果采樣頻率Fs為1024Hz，采樣點數(shù)為1024點，則可以分辨到1Hz。1024Hz的采樣率采樣1024點，剛好是1秒，也就是說，采樣1秒時間的信號并做FFT，則結(jié)果可以分析到1Hz，如果采樣2秒時間的信號并做FFT，則結(jié)果可以分析到0.5Hz。如果要提高頻率分辨力，則必須增加采樣點數(shù)，也即采樣時間。頻率分辨率和采樣時間是倒數(shù)關(guān)系。這一部分的描述很不清晰這部分的分析很關(guān)鍵，有利于理解excel假設(shè)FFT之后某點n用復(fù)數(shù)a+bi表示，那么這個復(fù)數(shù)的模就是An=根號a*a+b*b，相位就是Pn=atan2(b,a)。根據(jù)以上的結(jié)果，就可以計算出n點(n手1且n<=N/2)對應(yīng)的信號的表達式為：An/(N/2)*cos(2*pi*Fn*t+Pn)，即2*An/N*cos(2*pi*Fn*t+Pn)。對于n=1點的信號，是直流分量，幅度即為A1/N。由于FFT結(jié)果的對稱性，通常我們只使用前半部分的結(jié)果，即小于采樣頻率一半的結(jié)果。下面以一個實際的信號來做說明。假設(shè)我們有一個信號，它含有2V的直流分量，頻率為50Hz、相位為-30度、幅度為3V的交流信號，以及一個頻率為75Hz、相位為90度、幅度為1.5V的交流信號。用數(shù)學(xué)表達式就是如下：S=2+3*cos(2*pi*50*t-pi*30/180)+1.5*cos(2*pi*75*t+pi*90/180)。式中cos參數(shù)為弧度，所以-30度和90度要分別換算成弧度。我們以256Hz的采樣率對這個信號進行采樣，總共采樣256點。按照我們上面的分析，F(xiàn)n=(n-1)*Fs/N，我們可以知道，每兩個點之間的間距就是1Hz，第n個點的頻率就是n-1。我們的信號有3個頻率：0Hz、50Hz、75Hz，應(yīng)該分別在第1個點、第51個點、第76個點上出現(xiàn)峰值，其它各點應(yīng)該接近0。實際情況如何呢？我們來看看FFT的結(jié)果的模值如圖所示。從圖中我們可以看到，在第1點、第51點、和第76點附近有比較大的值。我們分別將這三個點附近的數(shù)據(jù)拿上來細看：點：512+0i點：-2.6195E-14-1.4162E-13i點：-2.8586E-14-1.1898E-13i點：-6.2076E-13-2.1713E-12i點:332.55-192i點：-1.6707E-12-1.5241E-12i點:-2.2199E-13-1.0076E-12i點：3.4315E-12+192i點:-3.0263E-14+7.5609E-13i很明顯，1點、51點、76點的值都比較大，它附近的點值都很小，可以認為是0,即在那些頻率點上的信號幅度為0。接著，我們來計算各點的幅度值。分別計算這三個點的模值，結(jié)果如下：1點：51251點:38476點：192按照公式，可以計算出直流分量為：512/N=512/256=2；50Hz信號的幅度為：384/(N/2)=384/(256/2)=3；75Hz信號的幅度為192/(N/2)=192/(256/2)=1.5?？梢?，從頻譜分析出來的幅度是正確的。然后再來計算相位信息。直流信號沒有相位可言，不用管它。先計算50Hz信號的相位，atan2(-192,332.55)=-0.5236,結(jié)果是孤度，換算為角度就是180*(-0.5236)/pi=-30.0001。再計算75Hz信號的相位，atan2(192,3.4315E-12)=1.5708孤度，換算成角度就是180*1.5708/pi=90.0002?？梢?，相位也是對的。根據(jù)FFT結(jié)果以及上面的分析計算，我們就可以寫出信號的表達式了，它就是我們開始提供的信號。利用Excel進行FFT和Fourier分析的基本步驟杭州市2000人口分布密度［根據(jù)2000年人口普查的街道數(shù)據(jù)經(jīng)環(huán)帶(rings)平均計算得到的結(jié)果，數(shù)據(jù)由馮健博士處理］。下面的變換實質(zhì)是一種空間自相關(guān)的分析過程。第一步，錄入數(shù)據(jù)在Excel中錄入數(shù)據(jù)不贅述(見表1)。表1原始數(shù)據(jù)序列 表2補充后的數(shù)據(jù)序列

HAHAB□1距離人」密度20.3281S430.92682141.52462152.12317662.71891013.31960183.91694594.51032^105.17282115.T6200126.3S&44136.94297147.53306158.13153168.T2683179.32354189.920281910.51S282011.11S512111.715812212.314902312.914您2413.512782514.110332614.79582715.3882183.91694F94.51082^105.17282115.T6200126.35&44136.94297147.53806158.13153168.72683179.323541S9.920281910.513282011.11&512111.715812212.314902312.914652413.512782514.110332614.795S2715.38822S029030031032033034第二步，補充數(shù)據(jù)由于Fourier變換（FT）一般是借助快速Fourier變換（FastFourierTransformation,FFT）算法，而這種算法的技術(shù)過程涉及到對稱處理，故數(shù)據(jù)序列的長度必須是2n（N=1,2,3,???，）。如果數(shù)據(jù)序列長度不是2n,就必須對數(shù)據(jù)進行補充或者裁減?，F(xiàn)在數(shù)據(jù)長度是26,介于24=16到25=32之間，而26到32更近一些，如果裁減數(shù)據(jù)，就會損失許多信息。因此，采用補充數(shù)據(jù)的方式。補充的方法非常簡單，在數(shù)據(jù)序列后面加0，直到序列長度為32=25為止（表2）。當然，延續(xù)到64=26也可以，總之必須是2的整數(shù)倍。不過，補充的“虛擬數(shù)據(jù)”越多，變換結(jié)果的誤差也就越大。第三步，F(xiàn)ourier變換的選項設(shè)置沿著工具（Tools）一數(shù)據(jù)分析（DataAnalysis）的路徑打開數(shù)據(jù)分析復(fù)選框（圖1）。

在數(shù)據(jù)分析選項框中選擇傅立葉分析(在數(shù)據(jù)分析選項框中選擇傅立葉分析(FourierAnalysis)(圖2)。陞］文件(已編撤日視圖(心插入①格式。工具⑴數(shù)據(jù)E窗口(如幫助皿□芽19藥?B爭龜窿9*拳■排寫⑤,.， F7二自動保存心.，B31v S0ABC共享工作懣(E，,,保護(E) ，聯(lián)機協(xié)作(11) 卜F1距離」人口密度20.32818430.926321規(guī)劃求解毋?,41.52462152.123176加載宏①,,.自定義?心選項(皿,，62.71391073.319601E3.916945數(shù)據(jù)嘗折回,,g4.510829105.17282圖1數(shù)據(jù)分析(DataAnalysis)的路徑圖2數(shù)據(jù)分析(DataAnalysis)在Fourier分析對話框中進行如下設(shè)置：在 輸入?yún)^(qū)域中輸入數(shù)據(jù)序列的單元格范圍“$B$1:$B$33”；選中“標志位于第一行(L)”；將輸出區(qū)域設(shè)為“$C$2”或者“$C$2:$C$33”圖3a)。a

b圖3傅立葉分析(FourierAnalysis)注意：如果“輸入?yún)^(qū)域”設(shè)為“$B$2:$B$33”，則不選“標志位于第一行(L)(圖3b)。表3FFT的結(jié)果

A■_Bc1距離大口密度FFT結(jié)果20.328184218701.857013486 30.926821104459.634221S27-103400.5376153S6141.52462137684.6丑4056542-39512.0012034132152.12317614197.250796317-55505.05387560S8162.71E91011894.2167827967-29641.7572S0&S23173.31960118938.3527863987-20821.94263532651E3.91694520635.720934404-23004.8596457225194.51082914487.9855812052-23932.2738&S&426I105.172829835.411667338-16677.82583265551115.7S20011919.46805S0233-7585.65666S7729S1126.3564416112.0959220229-E044.183634845861136.9429717113.44821S299-7100.569507992161147.5380616598.4374913769-7622.17758&S9195I158.1315315868.1982894284-7046.263462401091IS8.726S313454.062810094-5175.724568&12521179.3235412259.0748348707-1721.407440893981189.9202812261.30462533861910.5182812259.074834870S+1721.4074403940312011.1165113454.062810094+5175.7245686125212111.715S115868.1982894284+7046.26346240109122IS.3149016598.4374913769+7622.1775868919412312.9146517113.448218299+7100.5695079921412413.5127816112.0959220223+5044.1836348458412514.1103311919.4680580232+7585.6566687729912614.79589S35.41166733S02+16677.825832&555i2715.388214487.9855812052+23932.2738636426123020635.7209344041+23004.859645T225129018938.3527863987+20S21.9426353265130011894.216782796S+29641.7572S0SS23i31014197.2507963171+55505.0538T560S91MM*|R原館數(shù)—FF；/業(yè)韭"/第四步，輸出FFT結(jié)果選項設(shè)置完畢以后，確定(OK),立即得到FFT結(jié)果(表3)。顯然，表3給出的都是復(fù)數(shù)(complexnumbers)o假定一個數(shù)據(jù)序列表為f(t)，則理論上Fourier變換的結(jié)果為F(w)=Mfq)e-j^tdt=Ff(t)],(一8<W<8)-3表3中給出的正是相應(yīng)于F(/)的復(fù)數(shù)，這里①為角頻率。第五步，計算功率譜Excel好像不能自動計算功率譜，這需要我們利用有關(guān)函數(shù)進行計算。計算公式為P(w)=T\f(w)|2=*(A2+B2)式中A為復(fù)數(shù)的實部(realnumber),B為虛部(imaginarynumber),T為假設(shè)的周期長度，實則補充后的數(shù)據(jù)序列長度。對于本例，T=32。注意復(fù)數(shù)的平方乃是一個復(fù)數(shù)與其共軛(conjugate)復(fù)數(shù)的乘積，若F(g+bj，則IF(刃)l2=(a+bj)*(a-bj)=a2+b2。這樣，根據(jù)表3中的FFT結(jié)果，我們有(218701.8572+02)/32=1494703196(104459.6342+103400.5382)/32=675108949其余依次類推。顯然，這樣計算非常繁瑣。一個簡單的辦法是調(diào)用Excel的模數(shù)(modulus)計算函數(shù)ImAbs，方法是在函數(shù)類別中找“其他”，在其他類中找“工程”類，在工程類中容易找到ImAbs函數(shù)(圖4)。確定以后，彈出一個選項框，選中第一個FFT結(jié)果，確定，得到218701.857(圖5)。我們知道，復(fù)數(shù)的模數(shù)計算公式為M=(A2+B2)1/2圖4模數(shù)計算函數(shù)對于第一個FFT結(jié)果，由于虛部為0，模數(shù)就是其自身，即(218701.8572+02)1/2=218701.857但對于后面真正的復(fù)數(shù)，就不一樣了。抓住第一個模數(shù)所在的單元格的右下角往下一拉，或者用鼠標雙擊該單元格的右下角，立即得到全部模數(shù)。ABCD距離0.3人口密度2S1S4FFT結(jié)果218701.857013486 :功率譜密度=ImAbs（C2）IniAl'EiRumberKJ="218701.85701348(=218701.B57FieiwnstheEtljEolu+ev：alne(jnodulus.lofaconiplextluhiIier氐sber回 計算結(jié)果=218T01.857 |確定|取稍|圖5計算模數(shù)最后，用模數(shù)的2次方除以數(shù)據(jù)長度32立即得到全部功率譜密度結(jié)果（表4）。表4功率譜密度D E FFFT結(jié)果 模數(shù) 功率請密度阪血M的結(jié)果218701.857013486218701.85714947031961495000000104459.634221S27-103400.53761533611469S1.244967S10894967510000037634.6714056542-39512.0012034132197121.227432銅76665129430000014197.250796317-55505.053S75608Si57291.9971410257415410260000011894.2167827967-29641.7572806823131939.1009231878317.73188000013933.3527363937-20321.9426353265i2S146.3052924756703.22476000020635.720934404-23004.8596457225i30903.9891629S45517.12985000014487.9855812052-23932.2738686426i27975.9800324457983.1244600009835.411667338-16677.8258326555i19361.9523111715162.41172000011919.4680580233-7585.65666877298i14128.549326237997.06623800016112.0959220229-5044.lS3634845S6i16SS3.229068907606.99890800017113.448218299-7100.56950799216118528.0381410727756.2107300001&598.4374913769-7622.17758&8919511S264.S766310425178.71042000015868.1982894284-7046.26346240109117362.302439420298.3942000013454.062310094-5175.72456361252i14415.26729649374T.84649400012259.0748348707-1721.40744089398i12379.344064789004.98479000012261.304625338612261.3046346^3112.224697000下表是利用Mathcad2000計算的功率譜密度（表5）。利用Mathcad進行FFT，過程要簡單得多，只要調(diào)用FFT命令，可以直接給出各種結(jié)果（包括圖表）。但Mathcad的計算不求精度，有一定誤差。將Mathcad的變換結(jié)果copy到Excel中進行比較，可以看到，如果不計誤差，二者是一致的（表4）。表5借助Mathcad2000進行FFT的結(jié)果

0101.495?0916.751?0822.948?0831.026?0843.188?0752.476?0762.985?0772.446?0781.172?0796.238?06108.908?06111.073?07121.042?07139.42?06146.494?06154.79?06164.697?06Power第六步，功率譜分析功率譜分析目前主要用于兩個方面，一是偵測系統(tǒng)變化的某種周期或者節(jié)律，據(jù)此尋找因果關(guān)系（解釋）或者進行某種發(fā)展預(yù)測（應(yīng)用）；二是尋找周期以外的某些規(guī)律，據(jù)此對系統(tǒng)的時空結(jié)構(gòu)特征進行解釋。表6以對稱點（f=0.5）為界，從完整的數(shù)據(jù)序列中截取一半1自然數(shù)列頻率助率譜密度百然數(shù)列頻率助率借輯技2001494703196001494703196310.0312567510894910.03125675108949420.062529476665120.0625294766651530.0937510257415430.09375102574154640.12531878317.740.12531878317.7750.1562S24756703.250.1562524756703.2860.187529845517.160.187529845517.1970.2187524457983.170.21875244579S3.110S0.2511715162.4S0.2511715162.41190.2S1256237997.0690.281256237997.0612100.31258907606.99100.31253907606.9913110.3437510727756.2110.3437510727756.214120.37510425178.7120.37510425178.715130.406259420298.3130.406259420293.316140.43756493747.84140.43756493747.8417150.4&S754789004.98150.468754739004.9318160.54698112.22160.54698112.2219170.531254789004.98上面基于杭州人口密度數(shù)據(jù)的FFT，實際上是一種空間自相關(guān)分析過程，屬于FT的第二類

應(yīng)用。這種過程不以尋找周期為目標，實際上也不存在任何周期。不論目標是什么，都必須借助頻譜圖（頻率一功率譜密度圖）進行分析和解釋。下面第一步就是繪制頻譜圖。首先要計算頻率，線頻或角頻都可以，因為二者相差常數(shù)倍（2n）o一個簡單的辦法是，用0到T=32的自然數(shù)列除以T=32（表6）。如果采用的頻率變化范圍0?1，則繪制的頻譜圖是對稱的（圖6）。實際上，另一半是多余的，Mathcad2000自動生成的頻譜圖就沒有考慮另外一半兒（圖7）。因此，我們可以以對稱點上0.5為界，截取前面一半的數(shù)據(jù)，在Excel上繪制頻譜圖（圖8）。度密譜率功1600000000140000000012000000001000000000800000000600000000