實驗二-利用DFT分析離散信號頻譜

實驗二利用DFT分析離散信號頻譜1.利用FFT分析信號的頻譜；(1)確定DFT計算的參數(shù)；(2)進行理論值與計算值比擬，討論信號頻譜分析過程中誤差原因及改善方法。答：〔1〕信號x[k]的周期N=16，角頻率w0=2*pi/N=pi/8;作fft的點數(shù)N=16?！?〕計算程序：N=16;k=0:N-1;x=cos(3*pi/8*k);X=fft(x,N);subplot(2,1,1);stem(k-N/2,abs(fftshift(X)));ylabel('幅度譜')xlabel('m')subplot(2,1,2);stem(k-N/2,angle(fftshift(X)));ylabel('相位譜')xlabel('m')計算結(jié)果：我們還可以利用ifft函數(shù)進行相應(yīng)的傅里葉反變換，以得到時域上的x[k]，分析程序如下：N=16;k=0:N-1;x1=cos(3*pi/8*k);X=fft(x1,N);subplot(2,1,1)stem(k,x1)ylabel('傅里葉變換前的時域波形')xlabel('k')x2=ifft(X,N);subplot(2,1,2)stem(k,x2)ylabel('傅里葉反變換后的時域波形')xlabel('k')結(jié)果如下：可以看出：經(jīng)過傅里葉變換和相應(yīng)反變換后的信號與原信號的時域波形完全相同。2.利用FFT分析信號的頻譜；(1)確定DFT計算的參數(shù)；(2)進行理論值與計算值比擬，討論信號頻譜分析過程中誤差原因及改善方法。答：〔1〕信號無限長，因此需要對其進行截短。該序列單調(diào)衰減，當(dāng)k>=10時，序列已幾乎衰減為0，因此只取序列在區(qū)間[0，10]上的數(shù)值進行分析?！?〕計算程序：k=0:10;x=(1/2).^k;subplot(2,1,1);stem(k,x);subplot(2,1,2);w=k-5;plot(w,abs(fftshift(fft(x))));ylabel('幅度譜')xlabel('w')計算結(jié)果：3.有限長脈沖序列，利用FFT分析其頻譜，并繪出其幅度譜與相位譜。答：計算程序：k=0:5;x=[233105];w=k-3;subplot(2,1,1);stem(w,abs(fftshift(fft(x,6))));ylabel('幅度譜')xlabel('w')subplot(2,1,2);stem(w,angle(fftshift(fft(x,6))));ylabel('相位譜')xlabel('w')計算結(jié)果：4.某周期序列由3個頻率組成：，利用FFT分析其頻譜。如何選取FFT的點數(shù)N？此3個頻率分別對應(yīng)FFT計算結(jié)果X[m]中的哪些點？假設(shè)選取的N不適宜，F(xiàn)FT計算出的頻譜X[m]會出現(xiàn)什么情況？答：序列x[k]的周期為N=32，取FFT的點數(shù)N=32進行頻譜分析，程序如下：N=32;k=0:N-1;x=cos(7*pi/16*k)+cos(9*pi/16*k)+cos(pi/2*k);X=fft(x,N);subplot(2,1,1);stem(k-N/2,abs(fftshift(X)));ylabel('幅度譜')xlabel('m')subplot(2,1,2);stem(k-N/2,angle(fftshift(X)));ylabel('相位譜')xlabel('m’)分析結(jié)果：信號x[k]的三個頻率分別對應(yīng)FFT計算結(jié)果X[m]中的m=7,m=9,m=8對于離散周期序列進行FFT頻譜分析時取計算的點數(shù)為序列的周期N，如果計算的點數(shù)不是此離散周期序列的周期時，X[m]會發(fā)生變化，此時并不表示序列x[k]的頻譜，而相當(dāng)于將序列x[k]的周期擴大或者縮減后得到新信號的頻譜。在本例中，假設(shè)取N=20<32，那么可得到以下圖所示的幅度譜和相位譜：同理，假設(shè)取N=40>32，那么可以得到另外不同的幅度譜和相位譜，這里不再詳細(xì)畫出。5.某離散序列由3個頻率組成:利用FFT分析其頻譜。(1)對x[k]做64點FFT，繪出信號頻譜，能分辨出其中的兩個頻率嗎？(2)對x[k]補零到256點后計算FFT，能分辨出其中的兩個頻率嗎？(3)選用非矩形窗計算FFT，能夠分辨出其中的兩個頻率嗎？(4)假設(shè)不能夠很好地分辨出其中的兩個頻譜，應(yīng)采取哪些措施？答：〔1〕對x[k]做64點FFT,程序如下：k=0:63;x=cos(2*pi/15*k)+0.75*cos(2.3*pi/15*k);X=fft(x,64);stem(k-32,abs(fftshift(X)))ylabel('幅度譜')xlabel('m')可以得到幅度譜和相位譜如下：很明顯，根據(jù)上面的幅度譜不能分辨出信號x[k]中的兩個頻率?！?〕對x[k]補零后做256點FFT：k=0:63;x=cos(2*pi/15*k)+0.75*cos(2.3*pi/15*k);X=fft(x,256);n=-128:127stem(n,abs(fftshift(X)))運行結(jié)果如下：根據(jù)上面的幅度譜不能分辨出信號x[k]的兩個頻率?！?〕選用非矩形窗(主要有漢寧窗、哈明窗、布拉克曼窗和凱塞窗)計算fft，能夠分辨其中的兩個頻率，具體實現(xiàn)過程我將在下學(xué)期和數(shù)字信號處理這門課一起學(xué)習(xí)。6.序列利用FFT分析以下信號的幅頻特性，頻率范圍為，N=500點。(1)(2)(3)假設(shè)將上述x[k]乘以cos(pk/2)，重做(1)和(2)。答：序列x[k]是關(guān)于原點左右對稱的，因此我只考慮原點右側(cè)信號的幅頻特性。程序如下：k=0:50;x=exp(-0.5*(0.1*k).^2);X=fft(x,500);n=-250:249;subplot(3,1,1)stem(n,abs(fftshift(X)))ylabel('x[k]的傅里葉變換')xlabel('n')k=0:25;x=exp(-0.5*(0.2*k).^2);X=fft(x,500);n=-250:249;subplot(3,1,2)stem(n,abs(fftshift(X)))ylabel('x[2k]的傅里葉變換')xlabel('n')k=0:12;x=exp(-0.5*(0.4*k).^2);X=fft(x,500);n=-250:249;subplot(3,1,3)stem(n,abs(fftshift(X)))ylabel('x[4k]的傅里葉變換')xlabel('n')根據(jù)上述程序，我得到的上述三個信號的幅頻特性如以下圖所示：可以看出：信號在時域上進行壓縮，對應(yīng)那么在頻域上進行擴展，時域和頻域上的變換存在著一種相反的關(guān)系?！?〕假設(shè)將上述x[k]乘以cos(pk/2)，取p=3，那么與前問相同進行頻譜分析，程序如下：symspp=3;k=0:50;x=exp(-0.5*(0.1*k).^2).*cos(0.5*p*k);X=fft(x,500);n=-250:249;subplot(3,1,1)stem(n,abs(fftshift(X)))ylabel('x[k]的傅里葉變換')xlabel('n')k=0:25;x=exp(-0.5*(0.2*k).^2).*cos(p*k);X=fft(x,500);n=-250:249;subplot(3,1,2)stem(n,abs(fftshift(X)))ylabel('x[2k]的傅里葉變換')xlabel('n')k=0:12;x=exp(-0.5*(0.4*k).^2).*cos(2*p*k);X=fft(x,500);n=-250:249;subplot(3,1,3)stem(n,abs(fftshift(X)))ylabel('x[4k]的傅里葉變換')xlabel('n')運行結(jié)果如以下圖：四.實驗思考題1.既然可直接由DTFT定義計算序列DTFT，為何利用DFT分析序列的頻譜？答：在時域，兩者都是離散時間序列，在頻域，DTFT是連續(xù)周期曲線X(jw)，而DFT那么是與時間序列同長的離散序列頻譜X(k)。matlab無法計算連續(xù)變量w，只能在范圍內(nèi)把w賦值為很密的、長度很長的向量k*dw來近似連續(xù)變量。也就是說，能夠用計算機以很高效率計算的，只有離散傅里葉變換DFT,所以利用DFT分析序列的頻譜。2.假設(shè)序列持續(xù)時間無限長，且無解析表達式，如何利用DFT分析其頻譜？答：當(dāng)序列為無限長時，需要根據(jù)能量分布，進行截短。確定作FFT的點數(shù)N，然后使用fft函數(shù)作N點FFT計算X[m]。3.在利用DFT分析離散信號頻譜時，會出現(xiàn)哪些誤差？如何克服或改善？答：DFT是DTFT在主值區(qū)間的等間隔采樣，DTFT那么是DFT在密集的頻點上的插值，并作周期延拓。插值肯定存在精度的問題，因此用DFT分析離散信號頻譜時會產(chǎn)生誤差，我認(rèn)為可以通過增加作fft的點數(shù)N來減小插值間隔，從而克服或改善因插值帶來的誤差。4.在利用DFT分析離散信號頻譜時，如何選擇窗函數(shù)？答：(1)假設(shè)序列為離散周期序列，那么先確定離散周期序列x[k]的根本周期N，然后利用fft函數(shù)求其一個周期的DFT，得到X[m]；(2)假設(shè)序列為離散非周期序列，當(dāng)序列為無限長時，需要根據(jù)能量分布，首先進行截短。然后確定fft的點數(shù)N，使用fft函數(shù)作N點FFT計算X[m]；(3)如果序列中有兩個非常接近的頻率，用矩形窗進行fft計算難以分辨時，可以通過選擇適宜的非矩形窗〔漢寧窗