為什么需要FFT？

FFT(快速傅立葉變換)是離散傅立葉變換的快速算法.
傅立葉變換的物理意義在哪里？
傅立葉原理表明：任何連續(xù)測量的時序或信號，都可以表示為不同頻率的正弦波信號的無限疊加。而根據(jù)該原理創(chuàng)立的傅立葉變換算法利用直接測量到的原始信號，以累加方式來計(jì)算該信號中不同正弦波信號的頻率、振幅和相位。當(dāng)然這是從數(shù)學(xué)的角度去看傅立葉變換。
那么從物理的角度去看待傅立葉變換，它其實(shí)是幫助我們改變傳統(tǒng)的時間域分析信號的方法轉(zhuǎn)到從頻率域分析問題的思維，下面的一幅立體圖形可以幫助我們更好得理解這種角度的轉(zhuǎn)換：

所以，最前面的時域信號在經(jīng)過傅立葉變換的分解之后，變?yōu)榱瞬煌也ㄐ盘柕寞B加，我們再去分析這些正弦波的頻率，可以將一個信號變換到頻域。有些信號在時域上是很難看出什么特征的，但是如果變換到頻域之后，就很容易看出特征了。這就是很多信號分析采用FFT變換的原因。另外，F(xiàn)FT可以將一個信號的頻譜提取出來，這在頻譜分析方面也是經(jīng)常用的。
傅立葉變換提供給我們這種換一個角度看問題的工具，看問題的角度不同了，問題也許就迎刃而解！

變換是如何進(jìn)行的？

    首先，按照被變換的輸入信號類型不同，傅立葉變換可以分為4種類型：  
     1、 非周期性連續(xù)信號傅立葉變換（Fourier Transform）  
     2、 周期性連續(xù)信號傅立葉級數(shù)(Fourier Series)  
     3、 非周期性離散信號離散時域傅立葉變換（Discrete Time Fourier Transform）  
     4、 周期性離散信號離散傅立葉變換(Discrete Fourier Transform)  

下面是四種原信號圖例：

這里我們要討論是離散信號，對于連續(xù)信號我們不作討論，因?yàn)橛?jì)算機(jī)只能處理離散的數(shù)值信號，我們的最終目的是運(yùn)用計(jì)算機(jī)來處理信號的。所以對于離散信號的變換只有離散傅立葉變換（DFT）才能被適用，對于計(jì)算機(jī)來說只有離散的和有限長度的數(shù)據(jù)才能被處理，對于其它的變換類型只有在數(shù)學(xué)演算中才能用到，在計(jì)算機(jī)面前我們只能用DFT方法，我們要討論的FFT也只不過是DFT的一種快速的算法。
DFT的運(yùn)算過程是這樣的：

其中，
X(k)—頻域值
X(n)—時域采樣點(diǎn)
n—時域采樣點(diǎn)的序列索引
k—頻域值的索引
N—進(jìn)行轉(zhuǎn)換的采樣點(diǎn)數(shù)量

可見，在計(jì)算機(jī)或者示波器上進(jìn)行的DFT，使用的輸入值是數(shù)字示波器經(jīng)過ADC后采集到的采樣值，也就是時域的信號值，輸入采樣點(diǎn)的數(shù)量決定了轉(zhuǎn)換的計(jì)算規(guī)模。變換后的頻譜輸出包含同樣數(shù)量的采樣點(diǎn)，但是其中有一半的值是冗余的，通常不會顯示在頻譜中，所以真正有用的信息是N/2+1個點(diǎn)。
FFT的過程大大簡化了在計(jì)算機(jī)中進(jìn)行DFT的過程，簡單來說，如果原來計(jì)算DFT的復(fù)雜度是N2次運(yùn)算（N代表輸入采樣點(diǎn)的數(shù)量），進(jìn)行FFT的運(yùn)算復(fù)雜度是Nlg10（N），因此，計(jì)算一個1,000采樣點(diǎn)的DFT，使用FFT算法只需要計(jì)算3,000次，而常規(guī)的DFT算法需要計(jì)算1,000,000次！
我們以一個4個點(diǎn)的DFT變換為例來簡單說明FFT是怎樣實(shí)現(xiàn)快速算法的：

計(jì)算得出：

其中的紅色部分在FFT中是必須計(jì)算的分量，其他藍(lán)色部分不需要直接計(jì)算，可以由紅色的分量直接推導(dǎo)得到，比如：
x(1)e-j0 = -1*x(1)e-jπ
x(2)e-j0 = x(2)e-j2π
… …

這樣，已經(jīng)計(jì)算出的紅色分量只需要計(jì)算機(jī)將結(jié)果保存下來用于之后計(jì)算時調(diào)用即可，因此大大減少了DFT的計(jì)算量。

變換前后信號有何種對應(yīng)關(guān)系？

我們以一個實(shí)際的信號為例來說明：
示波器采樣得到的數(shù)字信號，就可以做FFT變換了。N個采樣點(diǎn)，經(jīng)過FFT之后，就可以得到N個點(diǎn)的FFT結(jié)果。為了方便進(jìn)行FFT運(yùn)算，通常N取2的整數(shù)次方。

假設(shè)采樣頻率為Fs，信號頻率F，采樣點(diǎn)數(shù)為N。那么FFT之后結(jié)果就是一個為N點(diǎn)的復(fù)數(shù)。每一個點(diǎn)就對應(yīng)著一個頻率點(diǎn)。這個點(diǎn)的模值，就是該頻率值下的幅度特性。具體跟原始信號的幅度有什么關(guān)系呢？假設(shè)原始信號的峰值為A，那么FFT的結(jié)果的每個點(diǎn)（除了第一個點(diǎn)直流分量之外）的模值就是A的N/2倍。而第一個點(diǎn)就是直流分量，它的模值就是直流分量的N倍。而每個點(diǎn)的相位呢，就是在該頻率下的信號的相位。第一個點(diǎn)表示直流分量（即0Hz），而最后一個點(diǎn)N的再下一個點(diǎn)（實(shí)際上這個點(diǎn)是不存在的，這里是假設(shè)的第N+1個點(diǎn)，也可以看做是將第一個點(diǎn)分做兩半分，另一半移到最后）則表示采樣頻率Fs，這中間被N-1個點(diǎn)平均分成N等份，每個點(diǎn)的頻率依次增加。例如某點(diǎn)n所表示的頻率為：Fn=(n-1)*Fs/N。由上面的公式可以看出，F(xiàn)n所能分辨到頻率為為Fs/N，如果采樣頻率Fs為1024Hz，采樣點(diǎn)數(shù)為1024點(diǎn)，則可以分辨到1Hz。1024Hz的采樣率采樣1024點(diǎn)，剛好是1秒，也就是說，采樣1秒時間的信號并做FFT，則結(jié)果可以分析精確到1Hz，如果采樣2秒時間的信號并做FFT，則結(jié)果可以分析精確到0.5Hz。如果要提高頻率分辨率，則必須增加采樣點(diǎn)數(shù)，也即采樣時間。頻率分辨率和采樣時間是倒數(shù)關(guān)系。

下面這幅圖更能夠清晰地表示這種對應(yīng)關(guān)系：

變換之后的頻譜的寬度（Frequency Span）與原始信號也存在一定的對應(yīng)關(guān)系。根據(jù)Nyquist采樣定理，F(xiàn)FT之后的頻譜寬度（Frequency Span）最大只能是原始信號采樣率的1/2，如果原始信號采樣率是4GS/s，那么FFT之后的頻寬最多只能是2GHz。時域信號采樣周期（Sample Period）的倒數(shù)，即采樣率（Sample Rate）乘上一個固定的系數(shù)即是變換之后頻譜的寬度，即 Frequency Span = K*（1/ΔT），其中ΔT為采樣周期，K值取決于我們在進(jìn)行FFT之前是否對原始信號進(jìn)行降采樣（抽點(diǎn)），因?yàn)檫@樣可以降低FFT的運(yùn)算量。如下圖所示：

可見，更高的頻譜分辨率要求有更長的采樣時間，更寬的頻譜分布需要提高對于原始信號的采樣率，當(dāng)然我們希望頻譜更寬，分辨率更精確，那么示波器的長存儲就是必要的！它能提供在高采樣率下采集更長時間信號的能力。

頻譜泄露

所謂頻譜泄露，就是信號頻譜中各譜線之間相互干擾，使測量的結(jié)果偏離實(shí)際值，同時在真實(shí)譜線的兩側(cè)的其它頻率點(diǎn)上出現(xiàn)一些幅值較小的假譜。產(chǎn)生頻譜泄露的主要原因是采樣頻率和原始信號頻率不同步，造成周期的采樣信號的相位在始端和終端不連續(xù)。簡單來說就是因?yàn)橛?jì)算機(jī)的FFT運(yùn)算能力有限，只能處理有限點(diǎn)數(shù)的FFT，所以在截取時域的周期信號時，沒有能夠截取整數(shù)倍的周期。信號分析時不可能取無限大的樣本。只要有截?cái)嗖煌骄蜁行孤?。如下圖所示：

圖中被測信號的開始端相位和截止端相位相同，表示在采集時間內(nèi)有整數(shù)倍周期的信號被采集到，所以此時經(jīng)行FFT運(yùn)算后得出的頻譜不會出現(xiàn)泄露。

上圖的信號頻率為2.1MHz，采集時間內(nèi)沒有截取整數(shù)倍周期的信號，F(xiàn)FT運(yùn)算之后譜線的泄露現(xiàn)象嚴(yán)重，可以看到能量較低的譜線很容易被臨近的能量較高的譜線的泄露給淹沒住。
因此，避免頻譜泄露的方法除了盡量使采集速率與信號頻率同步之外，還可以采用適當(dāng)?shù)拇昂瘮?shù)。
另外一個方法是采集信號時間足夠長，基本上可以覆蓋到整個有效信號的時間跨度。這種方法經(jīng)常在瞬態(tài)捕捉中被使用到，比如說沖擊試驗(yàn)，如果捕捉的時間夠長，捕捉到的信號可以一直包括了振動衰減為零的時刻。在這種情況下，可以不加窗函數(shù)。
窗函數(shù)其實(shí)就是一個加權(quán)函數(shù)，它在截取的信號時間段內(nèi)有值，時間段之外值為0:，記為：
w(t)=g(t) -T/2<t<T/2
w(t)=0 其它
加窗在時域上表現(xiàn)的是點(diǎn)乘，因此在頻域上則表現(xiàn)為卷積。卷積可以被看成是一個平滑的過程。這個平滑過程可以被看出是由一組具有特定函數(shù)形狀的濾波器，因此，原始信號中在某一頻率點(diǎn)上的能量會結(jié)合濾波器的形狀表現(xiàn)出來，從而減小泄漏?；谶@個原理，人們通常在時域上直接加窗。
大多數(shù)的信號分析儀一般使用矩形窗（rectangular），漢寧（hann），flattop和其它的一些窗函數(shù)。

不同的窗函數(shù)對頻譜譜線的影響不同，基本形狀可以參看下圖：

可以看到，不同的窗函數(shù)的主瓣寬度和旁瓣的衰減速度都不一樣，所以對于不同信號的頻譜應(yīng)該使用適當(dāng)?shù)拇昂瘮?shù)進(jìn)行處理。
矩形窗(Rectangular)：加矩形窗等于不加窗，因?yàn)樵诮厝r域信號時本身就是采用矩形截取，所以矩形窗適用于瞬態(tài)變化的信號，只要采集的時間足夠長，信號寬度基本可以覆蓋整個有效的瞬態(tài)部分。
漢寧窗(Von Hann)：如果測試信號有多個頻率分量，頻譜表現(xiàn)的十分復(fù)雜，且測試的目的更多關(guān)注頻率點(diǎn)而非能量的大小。在這種情況下，需要選擇一個主瓣夠窄的窗函數(shù)，漢寧窗是一個很好的選擇。
flattop窗：如果測試的目的更多的關(guān)注某周期信號頻率點(diǎn)的能量值，比如，更關(guān)心其EUpeak,EUpeak-peak,EUrms，那么其幅度的準(zhǔn)確性則更加的重要，可以選擇一個主瓣稍寬的窗，flattop窗在這樣的情況下經(jīng)常被使用。

好了，說了半天，看著公式也暈，下面以一個實(shí)際的信號來做說明。
假設(shè)我們有一個信號，它含有2V的直流分量，頻率為50Hz、相位為-30度、幅度為3V的交流信號，以及一個頻率為75Hz、相位為90度、幅度為1.5V的交流信號。用數(shù)學(xué)表達(dá)式就是如下：

S=2+3cos(2pi50t-pi30/180)+1.5cos(2pi75t+pi90/180)

式中cos參數(shù)為弧度，所以-30度和90度要分別換算成弧度。我們以256Hz的采樣率對這個信號進(jìn)行采樣，總共采樣256點(diǎn)。按照我們上面的分析，F(xiàn)n=(n-1)*Fs/N，我們可以知道，每兩個點(diǎn)之間的間距就是1Hz，第n個點(diǎn)的頻率就是n-1。我們的信號有3個頻率：0Hz、50Hz、75Hz，應(yīng)該分別在第1個點(diǎn)、第51個點(diǎn)、第76個點(diǎn)上出現(xiàn)峰值，其它各點(diǎn)應(yīng)該接近0。實(shí)際情況如何呢？我們來看看FFT的結(jié)果的模值如圖所示：

從圖中我們可以看到，在第1點(diǎn)、第51點(diǎn)、和第76點(diǎn)附近有比較大的值。我們分別將這三個點(diǎn)附近的數(shù)據(jù)拿上來細(xì)看：
1點(diǎn)： 512+0i
2點(diǎn)： -2.6195E-14 - 1.4162E-13i
3點(diǎn)： -2.8586E-14 - 1.1898E-13i

50點(diǎn)：-6.2076E-13 - 2.1713E-12i
51點(diǎn)：332.55 - 192i
52點(diǎn)：-1.6707E-12 - 1.5241E-12i

75點(diǎn)：-2.2199E-13 -1.0076E-12i
76點(diǎn)：3.4315E-12 + 192i
77點(diǎn)：-3.0263E-14 +7.5609E-13i

很明顯，1點(diǎn)、51點(diǎn)、76點(diǎn)的值都比較大，它附近的點(diǎn)值都很小，可以認(rèn)為是0，即在那些頻率點(diǎn)上的信號幅度為0。接著，我們來計(jì)算各點(diǎn)的幅度值。分別計(jì)算這三個點(diǎn)的模值，結(jié)果如下：
1點(diǎn)： 512
51點(diǎn)：384
76點(diǎn)：192
按照公式，可以計(jì)算出直流分量為：512/N=512/256=2；50Hz信號的幅度為：384/(N/2)=384/(256/2)=3；75Hz信號的幅度為192/(N/2)=192/(256/2)=1.5?？梢姡瑥念l譜分析出來的幅度是正確的。
然后再來計(jì)算相位信息。直流信號沒有相位可言，不用管它。先計(jì)算50Hz信號的相位，atan2(-192, 332.55)=-0.5236,結(jié)果是弧度，換算為角度就是180(-0.5236)/pi=-30.0001。再計(jì)算75Hz信號的相位，atan2(192, 3.4315E-12)=1.5708弧度，換算成角度就是1801.5708/pi=90.0002?？梢?，相位也是對的。
根據(jù)FFT結(jié)果以及上面的分析計(jì)算，我們就可以寫出信號的表達(dá)式了，它就是我們開始提供的信號。

總結(jié)

假設(shè)采樣頻率為Fs，采樣點(diǎn)數(shù)為N，做FFT之后，某一點(diǎn)n（n從1開始）表示的頻率為：Fn=(n-1)*Fs/N；該點(diǎn)的模值除以N/2就是對應(yīng)該頻率下的信號的幅度（對于直流信號是除以N）；該點(diǎn)的相位即是對應(yīng)該頻率下的信號的相位。相位的計(jì)算可用函數(shù)atan2(b,a)計(jì)算。atan2(b,a)是求坐標(biāo)為(a,b)點(diǎn)的角度值，范圍從-pi到pi。要精確到xHz，則需要采樣長度為1/x秒的信號，并做FFT。要提高頻率分辨率，就需要增加采樣點(diǎn)數(shù)，這在一些實(shí)際的應(yīng)用中是不現(xiàn)實(shí)的，需要在較短的時間內(nèi)完成分析。解決這個問題的方法有頻率細(xì)分法，比較簡單的方法是采樣比較短時間的信號，然后在后面補(bǔ)充一定數(shù)量的0，使其長度達(dá)到需要的點(diǎn)數(shù)，再做FFT，這在一定程度上能夠提高頻率分辨力。具體的頻率細(xì)分法可參考相關(guān)文獻(xiàn)。

% [附錄：本測試數(shù)據(jù)使用的matlab程序]
close all; % 先關(guān)閉所有圖片
Adc=2; % 直流分量幅度
A1=3;   % 頻率F1信號的幅度
A2=1.5; % 頻率F2信號的幅度
F1=50; % 信號1頻率(Hz)
F2=75; % 信號2頻率(Hz)
Fs=256; % 采樣頻率(Hz)
P1=-30; % 信號1相位(度)
P2=90; % 信號相位(度)
N=256; % 采樣點(diǎn)數(shù)
t=[0:1/Fs:N/Fs]; % 采樣時刻

% 信號
S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);
% 顯示原始信號
plot(S);
title('原始信號');

figure;
Y = fft(S,N); % 做FFT變換
Ayy = (abs(Y)); % 取模
plot(Ayy(1:N)); % 顯示原始的FFT模值結(jié)果
title('FFT 模值');

figure;
Ayy=Ayy/(N/2);   % 換算成實(shí)際的幅度
Ayy(1)=Ayy(1)/2;
F=([1:N]-1)*Fs/N; % 換算成實(shí)際的頻率值
plot(F(1:N/2),Ayy(1:N/2));   % 顯示換算后的FFT模值結(jié)果
title('幅度-頻率曲線圖');

figure;
Pyy=[1:N/2];
for i=1:N/2
Pyy(i)=phase(Y(i)); % 計(jì)算相位
Pyy(i)=Pyy(i)*180/pi; % 換算為角度
end;
plot(F(1:N/2),Pyy(1:N/2));   % 顯示相位圖
title('相位-頻率曲線圖');