2.8.傅立葉分析

源自http://civil.nju.edu.tw/weng/excel/lectureNote/2.8.htm

2.8.傅立葉分析

2.8.1.簡介

傅立葉轉換(Fourier Transform)是一種目前十分重要而且廣泛應用於各行業的數位訊號分析技術，當儀器測量所得的數位訊號為時間-振幅的數據時，您可以使用傅立葉轉換將此一訊號轉換為頻率-振幅，從而進行此一訊號的頻率特性的分析，反之若您有頻率-振幅的數據，您可使用逆傅立葉轉換，將此一訊號數據轉換為時間-振幅的數據。

Excel 2000的增益集中有一分析工具箱，其中已附有傅立葉轉換的功能，對於資料點數小於4096點的數據，您可以使用此一增益集中的分析工具箱的傅立葉分析工具，雖然其功能没有專業軟體的功能強大，但對於簡易的數位訊號處理，他提供了一個親和介面的解決方案。

傅立葉積分定義為[1]

　

                                            (2.8.1)   

                                                                           

上式中j代表，在某些數學書中也有以i代表，本文中比照參考文獻1所使用的符號規則，仍以j代表。

當公式(2.8.1)中任一點的f為可積分時則H(f)存在，其中h(t)為時間的函數，而H(f)則為頻率的函數。一般而言H(f)為一複數的量，所以H(f)應可改寫為

                                  (2.8.2)

上式中R(f)為傅立葉轉換的實數部，I(f)為傅立葉轉換的虛數部。而絕對值為其振幅或稱為傅立葉頻譜(Fourier Spectrum)，的定義為:

                                                (2.8.3)

而稱為相位角，可表示為

                                                 (2.8.4)                                                                                   

由於傅立葉轉換的特性，當時間域的函數h(t)為實數且對y軸對稱時，則轉換後的頻率域函數H(f)才仍為實數且對y軸對稱，對於常見的數位訊號常為實數且非對稱，因而轉換後的函數H(f)大常為複數。

傅立葉轉換前後的實數-複數特性可參考表2.8.1，其中偶函數(Even)是指對y軸對稱的函數，即f(-t)=f(t)，奇函數為對y軸反對稱的函數，即f(-t)=-f(t)。

時間域函數h(t)		頻率域函數H(t)
實數		實部-偶函數 虛部-奇函數
虛數		實部-奇函數 虛部-偶函數
實部-偶函數 虛部-奇函數		實數
實部-奇函數 虛部-偶函數		虛數
實數偶函數		實數偶函數
實數奇函數		虛數奇函數
虛數偶函數		虛數偶函數
虛數奇函數		實數奇函數
複數偶函數		複數偶函數
複數奇函數		複數奇函數

表2.8.1 Fourier轉前後關係特性表

傅立葉轉換應用於數位取樣器,例如氣象局使用的A900及A900a型的地震儀,其取得的訊號為每秒200點,連續取樣時間90秒,最大取樣點數為N=18000點,取樣點間隔為1/200秒,即=0.005秒.此時所得的為一連串的數據,並為數學函數,此等數據的傅立葉轉換稱離散傅立葉轉換(Discrete Fourier Transform), 離散傅立葉轉換所使用的計算公式與數學式連續函數的相類似,公式為:

                                                                   _(2.8.5)

其逆向離散傅立葉轉換的公式可表示為

                                                                _(2.8.5)

上式中N為資料點的總點數，為資料點的時間問隔，假如應用於上述的A900型地震儀所記錄,全滿90秒的記錄檔，其、N分別為=0.005,N=18000,而n、k則為記數用的變數 。轉換後頻率域的函數為 ，時間率的函數則為。

離散傅立葉轉換所需演算相當耗時，例如有N=2^γ點，則離散傅立葉轉換需要N²次的複數相乘與N(N-1)次的複數相加。因而另一種稱為快速離散傅立葉轉換(Fast Foruier Transform，簡寫為FFT)的演算法因而發展出來，目前常的離散傅立葉轉換演算法均是採用依據Cooley-Tukey Algorithm演算法所發展的FFT程序。同樣的N=2^γ點的數據，依FFT演算法只需Nγ/2次的複數相乘與Nγ次的複數相加運算，若以相乘運算來計算耗時，則N=1024點的傅立葉轉換運算， 傳統的傅立葉轉換矩陣運算與FFT的運算耗時會相差到200:1，由此可見FFT的優點所在。FFT的演算法原理請參閱參考文獻1。

假如原來在時間域的資料點有N點，則經過FFT轉換後，在頻率域其數據仍為N點，但己經是複數了。這N點中第一個點為所有其他點的總合，而且前半數的N/2點與後半數的N/2點是共軛對稱的複數，對稱點在中點，例如有256點，則128點為其共軛對稱點。因此在時間域有N點，則以FFT轉換到頻率頻後只有N/2點可用。至於各點頻率差，若在時間域的數據每點間隔則在頻率域的頻率間隔為，例如前述的A900地震儀，其=0.005，取N=1024點做FFT轉換後可用點只有其半數512點，各點頻率間隔為=0.1953125cps.

2.8.2.設定傅立葉分析工具箱

在Excel 2000中的增益集中已附有傅立葉分析工具箱，但在預設的選單中並没有這一項目，要使用傅立葉分析工具箱您必需做簡單的設定，設定方法應先開啟Excel 2000，選取下拉選單的工具(T)/增益集(I)，增益集的對答框如圖2.8.1，勾選分析工具箱選項，按下確定按鈕即完成設定。

圖2.8.1.設定傅立葉分析工具箱增益集

完成增益集設定之後，當您再選取下拉選單的工具(T)項時，這一下拉選單會多出一個資料分析(D)的選項，選取這一選項，圖2.8.2的對答框會跳出。接著選取傅立葉分析，您就可以開始進行您的數位訊號傅立葉分析工作了。

圖2.8.2.資料分析工具箱

2.8.3.傅立葉分析範例

例題 2.8.1.如範例2.8.1數據檔有1956至1987年間，計32年的太陽黑子活動數據，試分析太陽黑子的活動周期。

解:

1.首先以Excel讀取文字檔數據資料，將文字資料讀入Excel中的方法諘參考第1.4節。

2.在Excel中以繪製X-Y散佈圖的圖表繪出太陽黑子活動的振幅-周期關係圖 (Periodogram) 以供參考，如圖2.8.3所示。

圖 2.8.3.太陽黑子活動的振幅-周期關係圖(Periodogram)

3.由下拉選單選取:工具(T)/資料分析(D)-傅立葉分析，傅立葉分析對答框出現如圖2.8.4。

圖2.8.4.傅立葉分析對答框

4.在傅立葉分析對答框中，將數據中的太陽黑子Size資料範圍填入”輸入範圍”，”輸出範圍”則可在工作表中選取一塊與輸入範圍大小相同的空白區域，接著下確定按鈕，即可完成傅立葉轉換。參考圖2.8.5。

圖2.8.5. 傅立葉轉換工作表

5.由於實數經過傅立葉轉換後的數據為複數，為繪圖表示起見，我們可使用IMABS()函數取其複數的振幅值。如圖2.8.5中的第E欄即是， 舉例而言，第儲存格E3的計算公式為

=IMABS(D3)

6.在太陽黑子的Size完成轉換後，接著我們必需計算頻率域各點的頻率間距，如前述其頻率間距為，因為記錄為年所以應為

參考圖2.8.5，第一列數據因為受到扭曲，其數值可可信，所以E2以手動特別改為零，而其對應頻率，即儲存格F2應為0.03125，F3應為F3+0.03125，其餘依此類推，所以只要複製F3儲存格，貼到F4:F17即可。接著我們繪出頻率與振幅的關係圖即可觀察出太陽黑子的活動頻率，這種頻率與振幅的關係圖我們習慣稱為傅氏頻譜(Fourier Spectrum)

假如以多少年會有一次循環來考慮，則應以t = 1/f 表示才較恰當，所以在第G欄填入第F欄的倒數。

至此已完成所需數據的計算，接著我們繪出週期與振幅的關係圖即可觀察出太陽黑子的活動周期如圖2.8.6，由圖中可知，其活動周期約在8年左右。

圖 2.8.6.太陽黑子活動周期分析圖

由於本例題所使用的數據較少，以致其精確度受影響，讀者可試用習題2.8.1的數據重做一次，您會發現它的活動周期已偏到11年。

習題 2.8.1.在習題2.8.1的數據檔數據中有1932至1987年間，計256年的太陽黑子活動數據，試依範例2.8.1的方式再分析太陽黑子的活動周期。

參考文獻：

1.The Fast Fourier Transform And Its Application，E.Oran Brigham，Avantek，1988.