~~ 如果有什么問題可以在我的個人博客留言 ,我會及時回復。歡迎來訪交流 ~~
FFT_頻譜分析(數字信號處理)
(一)實驗原理
用FFT對信號作頻譜分析是學習數字信號處理的重要內容。經常需要進行譜分析的信號是模擬信號和時域離散信號。對信號進行譜分析的重點在於頻譜分辨率及分析誤差。頻譜分辨率D和頻譜分析的點數N直接相關,其分辨率為2π/N 。因此2π/N≤D,可以據這個公式確定頻率的分辨率。
FFT分析頻譜的誤差在於得到的是離散譜,而信號(非周期信號)是連續譜,只有當N較大時,離散譜的包絡才能逼近於連續譜。因此N要適當選擇大一些。
周期信號的頻譜是離散譜,只有用整數倍周期的長度作FFT,得到的離散譜才能代表周期信號的頻譜。如果不知道信號周期,可以盡量選擇信號的觀察時間長一些。
對模擬信號進行譜分析時,首先要按照采樣定理將其變成時域離散信號。如果是模擬周期信號,也應該選取整數倍周期的長度,經過采樣后形成周期序列,按照周期序列的譜分析進行。
(二)實驗步驟
1.對以下序列進行頻譜分析
為上三角有限長序列,序列長度為8; 為下三角有限長序列,序列長度為8。選擇FFT變換點數N分別為8和16兩種情況進行頻譜分析,打印出頻譜特性曲線,觀察不同N值, 和 的頻譜特性曲線是否相同,進行討論分析並得出結論。整個頻譜分析過程通過Matlab軟件進行程序設計實現。以下為Matlab的程序實現:
function e_6_1
M=8;xa=1:(M/2);xb=(M/2):-1:1;x1n=[xa,xb]; x2n=[xb,xa];%產生長度為M的上下三角序列
subplot(321);stem([0:M-1],x1n),title('X_1(n)');xlim([-1,M]);%打印序列X_1(n)
subplot(322);stem([0:M-1],x2n),title('X_2(n)');xlim([-1,M]);%打印序列X_2(n)
subplot(323);fft_stem(x1n,8,1);%打印N=8點X_1(n)DFT幅頻特性
subplot(324);fft_stem(x2n,8,2);%打印N=8點X_2(n)DFT幅頻特性
subplot(325);fft_stem(x1n,16,1);%打印N=16點X_1(n)DFT幅頻特性
subplot(326);fft_stem(x2n,16,2);%打印N=16點X_2(n)DFT幅頻特性
end
function fft_stem(A,N,name) % 計算FFT並打印子函數
stem([2/N*(0:(N-1))],abs(fft(A,N)));ylabel('幅度');xlabel('w/π');
title([num2str(N),'點DFT[X_',num2str(name),'(n)]']);xlim([-(2/N),2]);%橫坐標范圍
end
2.對下列模擬信號進行頻譜分析
這是一個含有三個頻率成分的模擬信號,頻率分別為選擇采樣頻率Fs=64Hz,對DFT變換點數N分別為:16、32、64這三種情況進行頻譜分析,分別打印頻譜特性曲線。對三種點數的頻譜分析結果進行討論分析。整個頻譜分析過程通過Matlab軟件進行程序設計實現。以下為Matlab的程序實現:
function e_6_2
global T;%全局變量 采樣周期
Fs=64;T=1/Fs;%Fs采樣頻率
n=0:64-1;nT=n*T;%采樣點數
X3n=cos(8*pi*nT)+cos(16*pi*nT)+cos(20*pi*nT);%采樣
subplot(311);fft_stem(X3n,16);%對N=16的頻譜分析
subplot(312);fft_stem(X3n,32);%對N=32的頻譜分析
subplot(313);fft_stem(X3n,64);%對N=64的頻譜分析
end
function fft_stem(A,N) % 計算FFT並打印子函數
global T;Tp=N*T;F=1/Tp;%Tp 觀察時間 F 頻率分辨率
stem([-N/2:N/2-1]*F,abs(fftshift(fft(A,N))));xlabel('f(Hz)');
%[-N/2:N/2-1]*F產生DFT對應采樣點頻率(fftshift()零點居中,關於中心反轉)
ylabel('幅度');title([num2str(N),'點DFT[X_3(n)]']);xlim([-N*F/2,N*F/2-1]);%x范圍
end
三、實驗結果分析
(一)序列結果分析
通過對 和 頻譜特性曲線結果的比較與分析(為了便於觀察頻譜和讀取頻率值對數據進行歸一化處理,即分析以為橫坐標),可以得出當FFT變換點數N為8 時, 和 頻譜特性相同(見圖3中的b1與b2),而N為16時, 和 頻譜特性曲線不相同(見圖3中的c1與c2)。出現這兩種不同情況的原因為 和 為序列為8的有限長序列(如圖3的a1和a2)當取N為8時,將兩個序列做周期延拓后發現 兩個周期序列的波形僅存在相位上的差別,將 向右移動4個單位長度即可得到 ,即,因而兩序列的頻譜特性相等,而當取N為16時,原 和 序列通過補零的方式變為序列長度為16的序列。此時序列波形完全不同,所以頻譜特性必然不同。
(二)模擬信號結果分析
x_3(n) 共有3個頻率成分,f1=4Hz, f2 =8Hz, f3=10Hz。所以x_3(n)的周期為0.5s。采樣的頻率為 Fs=64Hz=16f1=8f2=6.4f3。N為16時,觀察時間Tp=16T=0.25s,T為Fs的倒數。不是x_3(n)周期的整數倍,出現混疊現象,嚴重失真,所以得到的頻譜不正確(如圖4中a3)。當N=32或64時,觀察時間Tp=0.5s或1s, 是周期的整數倍,所以會得到正確的頻譜(如圖4中b3和c3),頻譜中有3根譜線正好對應位於4Hz,8Hz和10Hz。而N為64時的頻譜幅度為40正好為N為32時頻譜幅度為20的兩倍,即驗證DFT對中期序列譜分析的理論。