實(shí)驗(yàn)二 用FFT對(duì)信號(hào)進(jìn)行頻譜分析

1、東莞理工學(xué)院實(shí)驗(yàn)報(bào)告 課程名稱： 數(shù)字信號(hào)處理 實(shí)驗(yàn)室名稱： 實(shí)驗(yàn)名稱：實(shí)驗(yàn)二 用FFT對(duì)信號(hào)進(jìn)行頻譜分析 指導(dǎo)老師： 所在院系： 專業(yè)班級(jí)： 姓名： 學(xué)號(hào)： 日期： 成績(jī)： 1、實(shí)驗(yàn)?zāi)康膶W(xué)習(xí)用FFT對(duì)連續(xù)信號(hào)和時(shí)域離散信號(hào)進(jìn)行譜分析的方法，了解可能出現(xiàn)的分析誤差及其產(chǎn)生原因，以便正確應(yīng)用FFT。 2、實(shí)驗(yàn)原理與方法用FFT對(duì)信號(hào)作頻譜分析是學(xué)習(xí)數(shù)字信號(hào)處理的重要內(nèi)容。經(jīng)常需要進(jìn)行譜分析的信號(hào)是模擬信號(hào)和時(shí)域離散信號(hào)。對(duì)信號(hào)進(jìn)行譜分析的重要問(wèn)題是頻譜分辨率F和分析誤差。頻譜分辨率直接和FFT的變換區(qū)間N有關(guān)，因?yàn)镕FT能夠?qū)崿F(xiàn)的頻率分辨率是2/N，因此要求2/NF?？梢愿鶕?jù)此式選擇FFT的變換

2、區(qū)間N。誤差主要來(lái)自于用FFT作頻譜分析時(shí)得到的是離散譜，而信號(hào)（周期信號(hào)除外）是連續(xù)譜，只有當(dāng)N較大時(shí)離散譜的包絡(luò)才能逼近于連續(xù)譜，因此N要適當(dāng)選擇大一些。對(duì)于長(zhǎng)度為M的有限長(zhǎng)序列x(n)，其N點(diǎn)DFT（NM）X(k)就是x(n)的FT即在0,2內(nèi)的N點(diǎn)等間隔采樣，頻譜分辨率就是采樣間隔2/N。對(duì)于周期為N的周期序列，其頻譜是離散譜：其中，因此周期序列的頻譜結(jié)構(gòu)也可以用離散傅立葉級(jí)數(shù)系數(shù)表示。截取的主值序列，并進(jìn)行N點(diǎn)DFT，得到：因此也可以用表示的頻譜結(jié)構(gòu)。如果截取長(zhǎng)度為的整數(shù)個(gè)周期mN，m為整數(shù)，即，那么于是，也可以表示的頻譜結(jié)構(gòu)。由此可見(jiàn)，對(duì)于周期序列，只有用整數(shù)倍周期的長(zhǎng)度作FFT，

3、得到的離散譜才能代表周期信號(hào)的頻譜。如果不知道信號(hào)周期，可以盡量選擇信號(hào)的觀察時(shí)間長(zhǎng)一些。對(duì)一般的非周期模擬信號(hào)進(jìn)行譜分析時(shí)，首先要按照采樣定理將其變成時(shí)域離散信號(hào)。如果模擬信號(hào)持續(xù)時(shí)間很長(zhǎng)，采樣點(diǎn)數(shù)太多，會(huì)導(dǎo)致無(wú)法存儲(chǔ)和計(jì)算，需截取有限點(diǎn)進(jìn)行DFT。因此對(duì)模擬信號(hào)進(jìn)行頻譜分析必然是近似的，其近似程度與信號(hào)帶寬、采樣頻率和截取長(zhǎng)度有關(guān)。工程上，濾除幅度很小的高頻成分和截去幅度很小的部分信號(hào)是允許的。如果是模擬信號(hào)為周期信號(hào)，也應(yīng)該選取整數(shù)倍周期的長(zhǎng)度，經(jīng)過(guò)采樣后形成周期序列，按照周期序列的譜分析進(jìn)行。3、實(shí)驗(yàn)內(nèi)容及步驟（1）對(duì)以下序列進(jìn)行譜分析： 先繪制各序列的FT頻譜圖，然后選擇FFT的變換

4、區(qū)間N為8和16 兩種情況，通過(guò)序列的DFT分析其頻譜。分別打印其幅頻特性曲線，并進(jìn)行對(duì)比、分析和討論。（2）對(duì)以下周期序列進(jìn)行譜分析： 選擇FFT的變換區(qū)間N為8和16 兩種情況分別對(duì)以上序列進(jìn)行頻譜分析。分別打印其幅頻特性曲線。并進(jìn)行對(duì)比、分析和討論。（3）對(duì)模擬周期信號(hào)進(jìn)行譜分析 選擇采樣頻率，變換區(qū)間N=16，32，64 三種情況進(jìn)行譜分析。分別打印其幅頻特性，并進(jìn)行分析和討論。4、實(shí)驗(yàn)程序清單clear all;close all%*內(nèi)容一：有限長(zhǎng)序列的頻譜分析*%產(chǎn)生三個(gè)序列x1n=ones(1,4);M=8;xa=1:(M/2); xb=(M/2):-1:1;x2n=xa,xb;

5、 %產(chǎn)生長(zhǎng)度為8的三角波序列 x3n=xb,xa;%獲得序列x1n的傅立葉變換X1jw=fft(x1n,1024); %用1024個(gè)點(diǎn)的fft獲得傅立葉變換的密集采樣w=2*pi/1024*(0:1023);%采樣頻率位置subplot(3,1,1);plot(w/pi,abs(X1jw); %繪制FT的幅頻特性title('(1a) X_1(ejomega)'); xlabel('omega/pi');ylabel('幅度');X1k8=fft(x1n,8); %8個(gè)點(diǎn)DFTX1k16=fft(x1n,16); %16個(gè)點(diǎn)DFT%補(bǔ)充繪圖8個(gè)點(diǎn)

6、及16個(gè)點(diǎn)的DFT幅頻特性，提示用stem指令繪圖w=2*pi/8*(0:7);subplot(3,1,2);stem(w/pi,X1k8,'.'); %繪制FT的幅頻特性title('(1b) 8點(diǎn)DFTx1(n)'); xlabel('omega/pi');ylabel('幅度');w=2*pi/16*(0:15);subplot(3,1,3);stem(w/pi,X1k16,'.'); %繪制FT的幅頻特性title('(1c) 16點(diǎn)DFTx1(n)'); xlabel('omega

7、/pi');ylabel('幅度');%補(bǔ)充程序?qū)2n及x3n的頻譜進(jìn)行分析clear all;close allM=8;xa=1:(M/2); xb=(M/2):-1:1;x2n=xa,xb; %產(chǎn)生長(zhǎng)度為8的三角波序列 X2jw=fft(x2n,1024); %用1024個(gè)點(diǎn)的fft獲得傅立葉變換的密集采樣w=2*pi/1024*(0:1023);%采樣頻率位置subplot(3,1,1);plot(w/pi,abs(X2jw); %繪制FT的幅頻特性title('(2a) X_2(ejomega)'); xlabel('omega/pi&

8、#39;);ylabel('幅度');X2k8=fft(x2n,8); %8個(gè)點(diǎn)DFTX2k16=fft(x2n,16); %16個(gè)點(diǎn)DFTw=2*pi/8*(0:7);subplot(3,1,2);stem(w/pi,abs(X2k8),'.'); %繪制FT的幅頻特性title('(2b) 8點(diǎn)DFTx2(n)'); xlabel('omega/pi');ylabel('幅度');w=2*pi/16*(0:15);subplot(3,1,3);stem(w/pi,abs(X2k16),'.');

9、 %繪制FT的幅頻特性title('(2c) 16點(diǎn)DFTx2(n)'); xlabel('omega/pi');ylabel('幅度');clear all;close allM=8;xa=1:(M/2); xb=(M/2):-1:1;x3n=xb,xa;X3jw=fft(x3n,1024); %用1024個(gè)點(diǎn)的fft獲得傅立葉變換的密集采樣w=2*pi/1024*(0:1023);%采樣頻率位置subplot(3,1,1);plot(w/pi,abs(X3jw); %繪制FT的幅頻特性title('(3a) X_2(ejomega)&

10、#39;); xlabel('omega/pi');ylabel('幅度');X3k8=fft(x3n,8); %8個(gè)點(diǎn)DFTX3k16=fft(x3n,16); %16個(gè)點(diǎn)DFTw=2*pi/8*(0:7);subplot(3,1,2);stem(w/pi,abs(X3k8),'.'); %繪制FT的幅頻特性title('(3b) 8點(diǎn)DFTx3(n)'); xlabel('omega/pi');ylabel('幅度');w=2*pi/16*(0:15);subplot(3,1,3);stem(w

11、/pi,abs(X3k16),'.'); %繪制FT的幅頻特性title('(3c) 16點(diǎn)DFTx3(n)'); xlabel('omega/pi');ylabel('幅度');%*內(nèi)容二：周期序列的頻譜分析*N=8; %設(shè)置FFT的變換區(qū)間長(zhǎng)度為8 %根據(jù)區(qū)間長(zhǎng)度8，截取序列n=0:N-1;x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4jw=fft(x4n,1024); X5jw=fft(x5n,1024); w=2*pi/8*(0:8);%補(bǔ)充程序，對(duì)兩序列進(jìn)行8點(diǎn)的DFT%補(bǔ)充

12、程序，根據(jù)區(qū)間長(zhǎng)度16，截取序列，并進(jìn)行16個(gè)點(diǎn)的DFT%補(bǔ)充程序，顯示8個(gè)點(diǎn)及16個(gè)點(diǎn)DFT的幅度頻譜X4k8=fft(x4n,8); %8個(gè)點(diǎn)DFTX5k8=fft(x5n,8); %8個(gè)點(diǎn)DFTN=16; %設(shè)置FFT的變換區(qū)間長(zhǎng)度為16 %根據(jù)區(qū)間長(zhǎng)度16，截取序列n=0:N-1;x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n,16); %16個(gè)點(diǎn)DFTX5k16=fft(x5n,16); %16個(gè)點(diǎn)DFTw=2*pi/8*(0:7);subplot(2,2,1);stem(w/pi,abs(X4k8),'.&

13、#39;); %繪制FT的幅頻特性title('(4b) 8點(diǎn)DFTx4(n)'); xlabel('omega/pi');ylabel('幅度');subplot(2,2,3);stem(w/pi,abs(X5k8),'.'); %繪制FT的幅頻特性title('(5b) 8點(diǎn)DFTx5(n)'); xlabel('omega/pi');ylabel('幅度');w=2*pi/16*(0:15);subplot(2,2,2);stem(w/pi,abs(X4k16),'.&

14、#39;); %繪制FT的幅頻特性title('(4c) 16點(diǎn)DFTx4(n)'); xlabel('omega/pi');ylabel('幅度');subplot(2,2,4);stem(w/pi,abs(X5k16),'.'); %繪制FT的幅頻特性title('(5c) 16點(diǎn)DFTx5(n)'); xlabel('omega/pi');ylabel('幅度');%*內(nèi)容三：模擬周期信號(hào)的頻譜分析*Fs=64;T=1/Fs;%設(shè)置采樣頻率和采樣間隔N=16;n=0:N-1;

15、%設(shè)置采樣點(diǎn)數(shù)為16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對(duì)x6(n)進(jìn)行16點(diǎn)的采樣 X6k16=fft(x6nT); %16點(diǎn)DFTX6k16=fftshift(X6k16); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %計(jì)算出頻譜采樣間隔k=-N/2:N/2-1;fk=k*F; %找出采樣頻譜對(duì)應(yīng)的頻率位置%補(bǔ)充程序，繪制16點(diǎn)DFT頻譜圖，通過(guò)采樣序列16點(diǎn)DFT反映模擬信號(hào)頻譜axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16) subplot(3,1,1);stem(fk,a

16、bs(X6k16),'.'); %繪制FT的幅頻特性title('(6a) 16點(diǎn)DFTx6(nT)'); xlabel('omega/pi');ylabel('幅度');%補(bǔ)充程序，設(shè)置采樣點(diǎn)數(shù)為32，通過(guò)采樣序列的32點(diǎn)DFT反映模擬信號(hào)頻譜%補(bǔ)充程序，設(shè)置采樣點(diǎn)數(shù)為64，通過(guò)采樣序列的64點(diǎn)DFT反映模擬信號(hào)頻譜N=32;n=0:N-1; %設(shè)置采樣點(diǎn)數(shù)為32x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對(duì)x6(n)進(jìn)行32點(diǎn)的采樣 X6k32=fft(x6nT); %

17、32點(diǎn)DFTX6k32=fftshift(X6k32); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %計(jì)算出頻譜采樣間隔k=-N/2:N/2-1;fk=k*F; %找出采樣頻譜對(duì)應(yīng)的頻率位置axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32) subplot(3,1,2);stem(fk,abs(X6k32),'.'); %繪制FT的幅頻特性title('(6b) 32點(diǎn)DFTx6(nT)'); xlabel('omega/pi');ylabel('幅度');N=64;n=0:N-1;

18、%設(shè)置采樣點(diǎn)數(shù)為64x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %對(duì)x6(n)進(jìn)行64點(diǎn)的采樣 X6k64=fft(x6nT); %64點(diǎn)DFTX6k64=fftshift(X6k64); %將零頻率移到頻譜中心 Tp=N*T;F=1/Tp; %計(jì)算出頻譜采樣間隔k=-N/2:N/2-1;fk=k*F; %找出采樣頻譜對(duì)應(yīng)的頻率位置axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64) subplot(3,1,3);stem(fk,abs(X6k64),'.'); %繪制FT的幅頻特性titl

19、e('(6c) 64點(diǎn)DFTx6(nT)'); xlabel('omega/pi');ylabel('幅度');5、實(shí)驗(yàn)結(jié)果及分析（1）程序運(yùn)行結(jié)果及討論1）對(duì)以下序列進(jìn)行譜分析： 分析：圖（1b）和（1c）說(shuō)明x1(n)=R4(n)的8點(diǎn)DFT和16點(diǎn)DFT分別是x1(n)的頻譜函數(shù)的8點(diǎn)和16點(diǎn)采樣。因?yàn)閤3(n）=x2(n+3）8R8（n），所以，x3(n)與x2(n)的8點(diǎn)DFT的模相等，如圖（2b）和（2c），但是當(dāng)N=16時(shí)，x2(n）與x3(n）不滿足循環(huán)移位關(guān)系，模值不相等。2）對(duì)以下周期序列進(jìn)行譜分析：分析：對(duì)于x4(n）= cos(p