北京交通大學(xué)DSP研究性學(xué)習(xí)報(bào)告頻譜計(jì)算

1、數(shù)字信號(hào)處理課程研究性學(xué)習(xí)報(bào)告DFT近似計(jì)算信號(hào)頻譜專題研討姓 名 李 帆 學(xué) 號(hào) 11214008 同組成員 張 靜 11214028 林 恒 11214068 王亞君 11214025 李亞偉 11214009 指導(dǎo)教師 薛 健 時(shí) 間 2013年5月8日 利用DFT近似計(jì)算信號(hào)頻譜專題研討【目的】(1) 掌握利用DFT近似計(jì)算不同類型信號(hào)頻譜的原理和方法；(2) 理解誤差產(chǎn)生的原因及減小誤差的方法；(3)研究用DFT近似計(jì)算連續(xù)周期信號(hào)的方法；(4) 培養(yǎng)學(xué)生自主學(xué)習(xí)能力，以及發(fā)現(xiàn)問(wèn)題、分析問(wèn)題和解決問(wèn)題的能力?！狙杏憙?nèi)容】 基本題基本題是課程的基本要求，所有的人都需完成。問(wèn)題一已知某離

2、散序列為 (1)用L=32點(diǎn)DFT計(jì)算該序列的頻譜，求出頻譜中譜峰的頻率；(2)對(duì)序列進(jìn)行補(bǔ)零，然后分別用L=64、128、256、512點(diǎn)DFT計(jì)算該序列的頻譜，求出頻譜中譜峰的頻率；(3)討論所獲得的結(jié)果，從中你能得到了什么結(jié)論？該結(jié)論對(duì)序列的頻譜計(jì)算有何指導(dǎo)意義？【題目分析】本題討論補(bǔ)零對(duì)離散序列頻譜計(jì)算的影響。 補(bǔ)零可以使DFT計(jì)算得出的頻譜更加細(xì)致，但是不能改變序列的DTFT【溫磬提示】在計(jì)算離散非周期序列頻譜時(shí)常用W/p作為橫坐標(biāo)，稱W/p為歸一化頻率(normalized frequency)。在畫頻譜時(shí)需給出橫坐標(biāo)。每幅圖下都需給出簡(jiǎn)要的文字說(shuō)明。由于離散非周期序列頻譜是周期的

3、，所以在計(jì)算時(shí)不必用fftshift 函數(shù)對(duì)fft計(jì)算的結(jié)果進(jìn)行重新排列?！拘蛄蓄l譜計(jì)算的基本方法】在MATLAB中，用函數(shù)fft（x，N）可以計(jì)算Xk序列的N點(diǎn)DFT【仿真結(jié)果】【結(jié)果分析】 通過(guò)對(duì)序列補(bǔ)零，使DFT在計(jì)算頻譜時(shí)，頻譜更加清晰，容易觀察，隨著點(diǎn)數(shù)的增加，頻譜的很多細(xì)節(jié)之處都顯示出來(lái)，頻譜也越來(lái)越精確。但是當(dāng)點(diǎn)數(shù)增加到一定范圍，頻譜基本不再變化，誤差也沒(méi)有減小，因此，在用DFT計(jì)算離散序列的頻譜時(shí)，點(diǎn)數(shù)合適即可， 不宜過(guò)少，也不宜過(guò)多。同時(shí)，不管取得的點(diǎn)數(shù)是多少，頻譜的譜峰所對(duì)應(yīng)的數(shù)值都是0.2，可見(jiàn)，對(duì)序列補(bǔ)零并不能改變序列的DTFT.【閱讀文獻(xiàn)】1、 數(shù)字信號(hào)處理2、 老師

4、的課件【發(fā)現(xiàn)問(wèn)題】 (專題研討或相關(guān)知識(shí)點(diǎn)學(xué)習(xí)中發(fā)現(xiàn)的問(wèn)題)：對(duì)序列后面補(bǔ)零能不能提高頻譜的分辨率?！締?wèn)題探究】 對(duì)序列后面補(bǔ)零不能提高頻譜的分辨率。DFT是對(duì)信號(hào)fourier變換的離散化處理，對(duì)序列后面補(bǔ)零，只是增加了信號(hào)fourier變換后的離散抽樣點(diǎn)，并不能改變信號(hào)本身的采樣點(diǎn)，故不能提高頻譜的分辨率?！痉抡娉绦颉縦=0:31x=sin(0.2*pi*k)k1=0:31x_32=fft(x,32)subplot(5,1,1)plot(2*k1/32,abs(x_32),'g')title('L=32') k2=0:63x_64=fft(x,64)subp

5、lot(5,1,2)plot(2*k2/64,abs(x_64),'r')title('L=64') k3=0:127x_128=fft(x,128)subplot(5,1,3)plot(2*k3/128,abs(x_128),'m')title('L=128') k4=0:255x_256=fft(x,256)subplot(5,1,4)plot(2*k4/256,abs(x_256),'y')title('L=256') k5=0:511x_512=fft(x,512)subplot(5,1,5

6、)plot(2*k5/512,abs(x_512),'k')title('L=512') figurek=0:31;x=sin(0.2*pi*k)L=0:1023X=fft(x,1024)subplot(2,1,1)plot(2*L/1024,abs(X),'r')hold on;k1=0:31x_32=fft(x,32)stem(2*k1/32,abs(x_32),'x','b')title('L=32') k=0:31x=sin(0.2*pi*k)L=0:1023X=fft(x,1024)sub

7、plot(2,1,2)plot(2*L/1024,abs(X),'r')hold onk3=0:127x_128=fft(x,128)stem(2*k3/128,abs(x_128),'*','b')title('L=1024')問(wèn)題二 某離散序列為 xk=AcosW0k+Bcos ( (W0+DW)k)。用長(zhǎng)度N=64的哈明窗對(duì)信號(hào)截短后近似計(jì)算其頻譜。試用不同的A和B的值(如 A和B近似相等，A和B近差距較大)，確定用哈明窗能分辯的最小的譜峰間隔中c的值?！绢}目分析】本題討論用哈明窗計(jì)算序列頻譜時(shí)的頻率分辨率問(wèn)題。 用哈明窗計(jì)

8、算序列頻譜時(shí)，可以減小序列的泄漏現(xiàn)象，即減少頻譜中出現(xiàn)的多余高頻分量，但是其是在降低頻譜的分辨率的基礎(chǔ)上實(shí)現(xiàn)的，故使用哈明窗會(huì)降低頻譜的分辨率?！痉抡娼Y(jié)果】1A=B=2時(shí)C=2C=1C=0.5C=0.7C=0.8C=0.9C=1.3C=1.5C=4C=32A=4，B=1C=4C=3C=2C=13A=1，B=4C=4C=3C=2C=1【結(jié)果分析】將實(shí)驗(yàn)結(jié)果與教材中定義的哈明窗有效寬度相比較，發(fā)表你的看法。 教材定義的標(biāo)準(zhǔn)是C=2，根據(jù)結(jié)果可以看出，當(dāng)C2時(shí)，隨著C的減少，頻譜譜峰的分辨率隨之降低，但是當(dāng)C由2變?yōu)?時(shí)，譜峰的分辨率逐漸提高，當(dāng)C由1變小時(shí)，分辨率也變低。這個(gè)規(guī)律與AB的取值無(wú)關(guān)。

9、但是AB的取值影響著兩個(gè)譜峰的幅度，當(dāng)A與B 相近時(shí)，兩個(gè)譜峰相當(dāng)，故當(dāng)C=2,譜峰還是能夠區(qū)分的，但是當(dāng)A、B相差很大時(shí)，兩個(gè)譜峰幅度不相當(dāng)，當(dāng)C取2時(shí)，譜峰基本不能分辨?！咀灾鲗W(xué)習(xí)內(nèi)容】不同窗函數(shù)對(duì)頻譜分辨率的影響【閱讀文獻(xiàn)】1、 數(shù)字信號(hào)處理2、 老師的課件【發(fā)現(xiàn)問(wèn)題】 (專題研討或相關(guān)知識(shí)點(diǎn)學(xué)習(xí)中發(fā)現(xiàn)的問(wèn)題)：理論上，隨著C數(shù)值的減小，頻譜分辨率會(huì)逐漸降低，但是C取1時(shí)出現(xiàn)了特殊，頻譜分辨率相對(duì)1周圍的數(shù)值較大?！締?wèn)題探究】在離散序列頻譜計(jì)算中為何要用窗函數(shù)？用不同的窗函數(shù)對(duì)計(jì)算結(jié)果有何影響？與矩形窗相比哈明窗有何特點(diǎn)?如何選擇窗函數(shù)？ 在實(shí)際中的信號(hào)一般都不是帶限信號(hào)，離散化后的信號(hào)

10、xk也不是帶限信號(hào)，無(wú)法使用DFT進(jìn)行頻譜分析，故而需要對(duì)其進(jìn)行加窗截短使之成為有限長(zhǎng)序列。頻譜中出現(xiàn)的高頻分量不同，即頻率泄露現(xiàn)象不同；頻譜的分辨率不同。與矩形窗函數(shù)相比，哈明窗減小了頻率泄露現(xiàn)象，因?yàn)楣鞔昂瘮?shù)突然截?cái)嗟默F(xiàn)象不是很明顯，故頻譜中出現(xiàn)的高頻分量比較少。當(dāng)要求頻譜的分辨率較高時(shí)，應(yīng)該選擇矩形窗，當(dāng)要求頻率泄漏現(xiàn)象不能太明顯時(shí)，應(yīng)該選擇哈明窗?！痉抡娉绦颉縦=0:63;A=input('A=');B=input('B=');c=input('c=');N=64;L=512;x=A*cos(0.1*pi*k)+B*cos(0.1*pi

11、+c*2*pi/N)*k);wh=(hamming(N)'x=x.*wh;X=fftshift(fft(x,L);fs=30;ws=2*pi*fsm=(-ws/2+(0:L-1)*ws/L)/(2*pi);plot(m,abs(X);grid on; 問(wèn)題三已知一離散序列為 xk=cos(W0k)+0.75cos(W1k), 0£ k £ 63 其中W0=0.4p, W1=W0+p/64(1) 對(duì)xk做64點(diǎn)FFT, 畫出此時(shí)信號(hào)的頻譜。(2) 如果(1)中顯示的譜不能分辨兩個(gè)譜峰，是否可對(duì)(1)中的64點(diǎn)信號(hào)補(bǔ)零而分辨出兩個(gè)譜峰。通過(guò)編程進(jìn)行證實(shí)，并解釋其原因 。

12、(3) 給出一種能分辨出信號(hào)中兩個(gè)譜峰的計(jì)算方案，并進(jìn)行仿真實(shí)驗(yàn)?！绢}目分析】分析影響譜峰分辨率的主要因數(shù)，進(jìn)一步認(rèn)識(shí)補(bǔ)零在在頻譜計(jì)算中的作用。 影響譜峰分辨率的主要因數(shù)是序列的頻率，即連續(xù)時(shí)間的抽樣點(diǎn)數(shù)?！痉抡娼Y(jié)果】（1）（2）L=128L=256L=512L=1024（3）0<=k<=127;L=128L=2560<=k<=256,L=256L=512L=1024【結(jié)果分析】 由結(jié)果可知，對(duì)信號(hào)增加零點(diǎn)不能提高頻譜的分辨率；為了提高頻譜的分辨率，可以增加序列的抽樣點(diǎn)數(shù)，比如開(kāi)始的時(shí)候，k的范圍是0 64，后來(lái)k的范圍是0 128,然后是0 256隨著抽樣點(diǎn)數(shù)的增加，

13、頻譜的分辨率逐漸提高?！咀灾鲗W(xué)習(xí)內(nèi)容】如何增加頻譜分辨率【閱讀文獻(xiàn)】1、 數(shù)字信號(hào)處理；2、 老師的課件?！景l(fā)現(xiàn)問(wèn)題】 (專題研討或相關(guān)知識(shí)點(diǎn)學(xué)習(xí)中發(fā)現(xiàn)的問(wèn)題)：頻譜的分辨率與k范圍的具體關(guān)系【問(wèn)題探究】問(wèn)題一、二、三討論的是離散信號(hào)頻譜的計(jì)算問(wèn)題。與連續(xù)信號(hào)頻譜計(jì)算問(wèn)題相比較，其計(jì)算誤差有何不同？ 離散信號(hào)的頻譜具有周期性，因而在計(jì)算中，可能會(huì)因抽樣頻率的不適而產(chǎn)生混疊，導(dǎo)致失真，誤差加大?！痉抡娉绦颉?1)k=0:63;L=64;f1=0.4*pi;f2=f1+pi/64;x=cos(f1*k)+0.75*cos(f2*k);x=x zeros(1,L-length(x);X=fftshi

14、ft(fft(x);f=-pi:2*pi/L:pi-2*pi/L;plot(f/pi,abs(X)title('64點(diǎn)FFT')(2)L=input('L=');k=0:63;f1=0.4*pi;f2=f1+pi/64;x=cos(f1*k)+0.75*cos(f2*k);x=x zeros(1,L-length(x);X=fftshift(fft(x);f=-pi:2*pi/L:pi-2*pi/L;plot(f/pi,abs(X)(3)k=0:127;256L=128;/256/512f1=0.4*pi;f2=f1+pi/64;x=cos(f1*k)+0.75

15、*cos(f2*k);x=x zeros(1,L-length(x);X=fft(x);f=-pi:2*pi/L:pi-2*pi/L;plot(f/pi,abs(X);title('128點(diǎn)DFT');問(wèn)題四 試用DFT近似計(jì)算高斯信號(hào)的頻譜抽樣值。高斯信號(hào)頻譜的理論值為通過(guò)與理論值比較，討論信號(hào)的時(shí)域截取長(zhǎng)度和抽樣頻率對(duì)計(jì)算誤差的影響。【題目分析】學(xué)會(huì)連續(xù)非周期信號(hào)頻譜計(jì)算的基本方法。分析計(jì)算中出現(xiàn)誤差的主要原因及減小誤差的方法。 計(jì)算連續(xù)非周期信號(hào)的頻譜時(shí)，先對(duì)其進(jìn)行抽樣，使之變?yōu)殡x散信號(hào)，然后對(duì)其進(jìn)行DFT計(jì)算其頻譜。為了減小誤差，有兩個(gè)方法，一是增加信號(hào)的抽樣頻譜；二是

16、對(duì)信號(hào)加窗進(jìn)行截短，使之變?yōu)橛邢揲L(zhǎng)序列?！痉抡娼Y(jié)果】fs=1; N=0.05N=2N=4N=8N=16當(dāng)N=8確定時(shí)，fs對(duì)計(jì)算誤差的影響Fs=2Fs=4Fs=8Fs=16Fs=32【結(jié)果分析】由于信號(hào)在時(shí)域和頻譜都有理論表達(dá)式，在進(jìn)行誤差分析時(shí)希望給出一些定量的結(jié)果。例如，可參考概率論中的對(duì)高斯概率密度函數(shù)的討論，定義信號(hào)時(shí)域的有效寬度，頻域的有效寬度，分析時(shí)域、頻域有效寬度對(duì)計(jì)算誤差的影響。 由結(jié)果可知，時(shí)域截取長(zhǎng)度和抽樣頻率對(duì)計(jì)算誤差都有影響。當(dāng)抽樣頻率一定時(shí)，時(shí)域截取長(zhǎng)度越大，計(jì)算誤差越??；當(dāng)時(shí)域截取長(zhǎng)度一定時(shí)，增大信號(hào)的抽樣頻率，計(jì)算誤差先減小后增大，在抽樣頻率在2左右時(shí)，誤差相對(duì)較

17、小?！咀灾鲗W(xué)習(xí)內(nèi)容】研究時(shí)域截取長(zhǎng)度和抽樣頻率對(duì)計(jì)算頻譜誤差的影響?！鹃喿x文獻(xiàn)】1、 數(shù)字信號(hào)處理；2、 老師的課件；3、 概率論?！景l(fā)現(xiàn)問(wèn)題】 (專題研討或相關(guān)知識(shí)點(diǎn)學(xué)習(xí)中發(fā)現(xiàn)的問(wèn)題)：如何判斷d，N, fs的取值。【仿真程序】d=input('d=');N=input('N=');fs=input('fs=');w=linspace(-2*pi,2*pi,2048);G=sqrt(pi/d).*exp(-w.*w)/(4*d);plot(w/(2*pi),abs(G);axis(-1.5,1.5,0,12)L=512;T=1/fs;t=(0:

18、N-1)*T;ws=2*pi*fs;g=exp(-d*t.*t);X=fftshift(fft(g,L);w1=(0:L-1)*ws/L-ws/2)/(2*pi);hold on;plot(w1,abs(X),'r');axis(-1.5,1.5,0,12);title('紅線表示非理論值，藍(lán)線表示理論值' );擴(kuò)展題擴(kuò)展題有一定難度，不要求每個(gè)同學(xué)都做。鼓勵(lì)大家進(jìn)行一下探索。問(wèn)題五本題研究連續(xù)周期信號(hào)頻譜的近似計(jì)算問(wèn)題。 周期為T0的連續(xù)時(shí)間周期信號(hào)x(t)可用Fourier級(jí)數(shù)表示為其中X(nw0)稱為連續(xù)時(shí)間周期信號(hào)x(t)的頻譜函數(shù)。稱為信號(hào)的基頻。如果

19、信號(hào)x(t)函數(shù)表達(dá)式已知，則可由積分得出信號(hào)的頻譜。如果信號(hào)x(t)函數(shù)表達(dá)式未知，或者x(t)函數(shù)表達(dá)式非常復(fù)雜，則很難由積分得信號(hào)的頻譜。本題的目的就是研究如何利用DFT近似計(jì)算連續(xù)時(shí)間周期信號(hào)的頻譜。(1)若在信號(hào)x(t)的一個(gè)周期T0內(nèi)抽樣N個(gè)點(diǎn)，即， T為抽樣周期(間隔)，可獲得序列xk試推導(dǎo)序列xk的DFT與連續(xù)時(shí)間周期信號(hào)x(t)的頻譜X(nw0)的關(guān)系；(2)由(1)的結(jié)論，給出由DFT近似計(jì)算周期信號(hào)頻譜X(nw0)的方案；(3)周期信號(hào)x(t)的周期T0=1，x(t)在區(qū)間0，1的表達(dá)式為x(t)=20t2(1-t)4cos(12pt)(a)試畫出信號(hào)x(t)在區(qū)間0，1

20、的波形；(b)若要用6次以內(nèi)的諧波近似表示x(t)，試給出計(jì)算方案，并計(jì)算出近似表示的誤差。討論出現(xiàn)誤差的原因及減小誤差的方法。 【題目分析】 為了利用DFT近似計(jì)算連續(xù)周期信號(hào)的頻譜，可先計(jì)算出序列DFT與連續(xù)周期信號(hào)的頻譜X（nw0）之間的關(guān)系。然后進(jìn)一步利用DFT近似計(jì)算連續(xù)周期信號(hào)的頻譜。【理論推導(dǎo)】DFT計(jì)算所得結(jié)果Xm與連續(xù)周期信號(hào)頻譜X(nw0)的關(guān)系。 由于,于是得： 令n=m+rN;m=0,1,2,，為整數(shù)，上式可化為化簡(jiǎn)可得對(duì)比IDFT 可得 【計(jì)算方案】根據(jù)理論推導(dǎo)結(jié)果設(shè)計(jì)近似計(jì)算方案。分析產(chǎn)生誤差的主要原因。 根據(jù)上面推出的關(guān)系式可以得出，Xm可以近似為較低的諧波分量。

21、這一過(guò)程的誤差主要在于低次諧波分量與高次諧波分量的混疊，如果抽樣點(diǎn)數(shù)N較大，則高頻分量的頻率值可以提高，當(dāng)提高到一定值時(shí)，其幅度值趨于0，這時(shí)就可以近似計(jì)算了。同時(shí)，在信號(hào)重建時(shí)，諧波分量的次數(shù)對(duì)誤差也有影響，在選取較低的諧波分量進(jìn)行信號(hào)重建時(shí)，會(huì)產(chǎn)生頻率泄露現(xiàn)象，增加誤差?！緮U(kuò)展分析】如果周期信號(hào)x(t)是帶限信號(hào)，即信號(hào)的最高頻率分量為Mw0(是正整數(shù))，試確定在一個(gè)周期內(nèi)的最少抽樣點(diǎn)N，使得在頻譜的計(jì)算過(guò)程當(dāng)中不存在混疊誤差。與抽樣定理給出的結(jié)論比較，發(fā)表你的看法。 若不發(fā)生混疊，則抽樣頻率fs>=Mw0，又T=T0/N,則N>=2T0Mw0，與抽樣定理基本相同，當(dāng)取到等號(hào)時(shí)，若邊界值不為0，則會(huì)產(chǎn)生混疊，故不能取到等號(hào)，若邊界值為0，則可取到等號(hào)?！痉抡娼Y(jié)果】N=32N=16N=32，n=8N=32，n=6【結(jié)果分析】討論DFT點(diǎn)數(shù)對(duì)近似計(jì)算的影響，討論所取諧波項(xiàng)的多少對(duì)近似計(jì)算的影響。誤差分析要給出定量的結(jié)果，如平均誤差，最大誤差等。與連續(xù)非周期信號(hào)頻譜計(jì)算過(guò)程中存在的誤