實驗二參考 快速傅立葉變換(FFT)及其應用

1、實驗二快速傅立葉變換（FFT）及其應用 一、實驗目的 1.在理論學習的基礎上，通過本實驗，加深對FFT的理解，熟悉FFT子程序。 2.熟悉應用FFT對典型信號進行頻譜分析的方法 3.了解應用FFT進行信號頻譜分析過程中可能出現(xiàn)的問題以便在實際中正確應用FFT。二、實驗原理 在各種信號序列中，有限長序列信號處理占有很重要地位，對有限長序列，我們可以使用離散Fouier變換(DFT)。這一變換不但可以很好的反映序列的頻譜特性，而且易于用快速算法在計算機上實現(xiàn)，當序列x(n)的長度為N時，它的DFT定義為：，反換為：有限長序列的DFT是其Z變換在單位圓上的等距采樣，或者是序列Fourier變換的等距

2、采樣，因此可以用于序列的譜分析。 FFT并不是與DFT不同的另一種變換，而是為了減少DFT運算次數(shù)的一種快速算法。它是對變換式進行一次次分解，使其成為若干小點數(shù)的組合，從而減少運算量。常用的FFT是以2為基數(shù)的，其長度 N=2L，它的效率高，程序簡單使用非常方便，當要變換的序列長度不等于2的整數(shù)次方時，為了使用以2為基數(shù)的FFT，可以用末位補零的方法，使其長度延長至2的整數(shù)次方。在運用DFT進行頻譜分析的過程中可能產生幾種問題：(1)混疊 序列的頻譜時被采樣信號的周期延拓，當采樣速率不滿足Nyquist定理時，就會發(fā)生頻譜混疊，使得采樣后的信號序列頻譜不能真實的反映原信號的頻譜。 避免混疊現(xiàn)象

3、的唯一方法是保證采樣速率足夠高，使頻譜混疊現(xiàn)象不致出現(xiàn)，即在確定采樣頻率之前，必須對頻譜的性質有所了解，在一般情況下，為了保證高于折疊頻率的分量不會出現(xiàn)，在采樣前，先用低通模擬濾波器對信號進行濾波。(2)泄漏 實際中我們往往用截短的序列來近似很長的甚至是無限長的序列，這樣可以使用較短的DFT來對信號進行頻譜分析，這種截短等價于給原信號序列乘以一個矩形窗函數(shù)，也相當于在頻域將信號的頻譜和矩形窗函數(shù)的頻譜卷積，所得的頻譜是原序列頻譜的擴展。 泄漏不能與混疊完全分開，因為泄漏導致頻譜的擴展，從而造成混疊。為了減少泄漏的影響，可以選擇適當?shù)拇昂瘮?shù)使頻譜的擴散減至最小。DFT是對單位圓上Z變換的均勻采樣

4、，所以它不可能將頻譜視為一個連續(xù)函數(shù)，就一定意義上看，用DFT來觀察頻譜就好像通過一個柵欄來觀看一個圖景一樣，只能在離散點上看到真實的頻譜，這樣就有可能發(fā)生一些頻譜的峰點或谷點被“尖樁的柵欄”所攔住，不能別我們觀察到。 減小柵欄效應的一個方法就是借助于在原序列的末端填補一些零值，從而變動DFT的點數(shù)，這一方法實際上是人為地改變了對真實頻譜采樣的點數(shù)和位置，相當于搬動了每一根“尖樁柵欄”的位置，從而使得頻譜的峰點或谷點暴露出來。 用FFT可以實現(xiàn)兩個序列的圓周卷積。在一定的條件下，可以使圓周卷積等于線性卷積。一般情況，設兩個序列的長度分別為N1和N2，要使圓周卷積等于線性卷積的充要條件是FFT的

5、長度 NN1N2 對于長度不足N的兩個序列，分別將他們補零延長到N。 當兩個序列中有一個序列比較長的時候，我們可以采用分段卷積的方法。有兩種方法： 重疊相加法。將長序列分成與短序列相仿的片段，分別用FFT對它們作線性卷積，再將分段卷積各段重疊的部分相加構成總的卷積輸出。 重疊保留法。這種方法在長序列分段時，段與段之間保留有互相重疊的部分，在構成總的卷積輸出時只需將各段線性卷積部分直接連接起來，省掉了輸出段的直接相加。 (3) 柵欄效應 DFT是對單位圓上z變換的均勻采樣，所以它不可能將頻譜視為一個連續(xù) 函數(shù)，從某種意義上講，用DFT來觀察頻譜就如同通過一個柵欄來觀看景象一 樣，只能在離散點上看

6、到真實的頻譜，這樣一些頻譜的峰點或谷點就可能被尖 樁的柵欄擋住，也就是正好落在兩個離散采樣點之間，不能被觀察到。減小柵欄效應的一個方法是在原序列的末端填補一些零值，從而變動DFT 的點數(shù)，這一方法實際上是人為地改變了對真實頻譜采樣的點數(shù)和位置，相當于 搬動了尖樁柵欄的位置，從而使得頻譜的峰點或谷點暴露出來。(4) DFT的分辨率填補零值可以改變對DTFT的采樣密度，人們常常有一種誤解，認為補零可以提高DFT的頻率分辨率，事實上，DFT的頻率分辨率通常規(guī)定為，這里的N是指信號的有效長度，而不是補零的長度。不同長度的，其DTFT的結果是不同的;而相同長度的盡管補零的長度不同其DTFT的結果應是相同

7、的，它們的DFT只是反映了對相同的DTFT采用了不同的采樣密度?？偨Y一下: 要提高DFT分辨率只有增加信號的截取長度N。三、實驗用到序列 a) 高斯序列 b) 衰減正弦序列 c) 三角波序列d) 反三角波序列四、實驗內容Matlab編程實現(xiàn)FFT實踐及頻譜分析 1.用Matlab產生正弦波,矩形波,以及白噪聲信號，并顯示各自時域波形圖2進行FFT變換，顯示各自頻譜圖，其中采樣率，頻率、數(shù)據(jù)長度自選3做出上述三種信號的均方根圖譜,功率圖譜,以及對數(shù)均方根圖譜4用IFFT傅立葉反變換恢復信號，并顯示恢復的正弦信號時域波形圖源程序%*1.正弦波*%fs=100;%設定采樣頻率N=128;n=0:N-

8、1;t=n/fs;f0=10;%設定正弦信號頻率%生成正弦信號x=sin(2*pi*f0*t);figure(1);subplot(231);plot(t,x);%作正弦信號的時域波形xlabel(t);ylabel(y);title(正弦信號y=2*pi*10t時域波形);grid;%進行FFT變換并做頻譜圖y=fft(x,N);%進行fft變換mag=abs(y);%求幅值f=(0:length(y)-1)*fs/length(y);%進行對應的頻率轉換figure(1);subplot(232);plot(f,mag);%做頻譜圖axis(0,100,0,80);xlabel(頻率(Hz

9、);ylabel(幅值);title(正弦信號y=2*pi*10t幅頻譜圖N=128);grid;%求均方根譜sq=abs(y);figure(1);subplot(233);plot(f,sq);xlabel(頻率(Hz);ylabel(均方根譜);title(正弦信號y=2*pi*10t均方根譜);grid;%求功率譜power=sq.2;figure(1);subplot(234);plot(f,power);xlabel(頻率(Hz);ylabel(功率譜);title(正弦信號y=2*pi*10t功率譜);grid;%求對數(shù)譜ln=log(sq);figure(1);subplot(

10、235);plot(f,ln);xlabel(頻率(Hz);ylabel(對數(shù)譜);title(正弦信號y=2*pi*10t對數(shù)譜);grid;%用IFFT恢復原始信號xifft=ifft(y);magx=real(xifft);ti=0:length(xifft)-1/fs;figure(1);subplot(236);plot(ti,magx);xlabel(t);ylabel(y);title(通過IFFT轉換的正弦信號波形);grid;%*2.矩形波*%fs=10;%設定采樣頻率t=-5:0.1:5;x=rectpuls(t,2);x=x(1:99);figure(2);subplot

11、(231);plot(t(1:99),x);%作矩形波的時域波形xlabel(t);ylabel(y);title(矩形波時域波形);grid;%進行FFT變換并做頻譜圖y=fft(x);%進行fft變換mag=abs(y);%求幅值f=(0:length(y)-1)*fs/length(y);%進行對應的頻率轉換figure(2);subplot(232);plot(f,mag);%做頻譜圖xlabel(頻率(Hz);ylabel(幅值);title(矩形波幅頻譜圖);grid;%求均方根譜sq=abs(y);figure(2);subplot(233);plot(f,sq);xlabel(

12、頻率(Hz);ylabel(均方根譜);title(矩形波均方根譜);grid;%求功率譜power=sq.2;figure(2);subplot(234);plot(f,power);xlabel(頻率(Hz);ylabel(功率譜);title(矩形波功率譜);grid;%求對數(shù)譜ln=log(sq);figure(2);subplot(235);plot(f,ln);xlabel(頻率(Hz);ylabel(對數(shù)譜);title(矩形波對數(shù)譜);grid;%用IFFT恢復原始信號xifft=ifft(y);magx=real(xifft);ti=0:length(xifft)-1/fs;

13、figure(2);subplot(236);plot(ti,magx);xlabel(t);ylabel(y);title(通過IFFT轉換的矩形波波形);grid;%*3.白噪聲*%fs=10;%設定采樣頻率t=-5:0.1:5;x=zeros(1,100);x(50)=100000;figure(3);subplot(231);plot(t(1:100),x);%作白噪聲的時域波形xlabel(t);ylabel(y);title(白噪聲時域波形);grid;%進行FFT變換并做頻譜圖y=fft(x);%進行fft變換mag=abs(y);%求幅值f=(0:length(y)-1)*fs

14、/length(y);%進行對應的頻率轉換figure(3);subplot(232);plot(f,mag);%做頻譜圖xlabel(頻率(Hz);ylabel(幅值);title(白噪聲幅頻譜圖);grid;%求均方根譜sq=abs(y);figure(3);subplot(233);plot(f,sq);xlabel(頻率(Hz);ylabel(均方根譜);title(白噪聲均方根譜);grid;%求功率譜power=sq.2;figure(3);subplot(234);plot(f,power);xlabel(頻率(Hz);ylabel(功率譜);title(白噪聲功率譜);grid;%求對數(shù)譜ln=log(sq);figure(3