數(shù)字信號處理研討_DFT近似計算信號頻譜

1、beijing jiaotong university數(shù)字信號處理研討dft近似計算信號頻譜學(xué) 院：電子信息工程學(xué)院小組成員：指導(dǎo)教師:2013.06.24dft近似計算信號頻譜專題研討【目的】(1) 掌握利用dft近似計算不同類型信號頻譜的原理和方法。(2) 理解誤差產(chǎn)生的原因及減小誤差的方法。(3) 培養(yǎng)學(xué)生自主學(xué)習(xí)能力，以及發(fā)現(xiàn)問題、分析問題和解決問題的能力?！狙杏戭}目】基本題1. 己知一離散序列為xk = sin(0.2 兀眾)，眾=0，1，31 用l=32點dft計算該序列的頻譜，求出頻譜中譜峰的頻率；(2) 對序列進(jìn)行補零，然后分別用l=64、128、256、512點dft u算該

2、序列的頻譜，求出頻譜中譜峰 的頻率；(3) 討論所獲得的結(jié)果，給出你的結(jié)論。該結(jié)論對序列的頻譜計算有何指導(dǎo)意義？【題目分析】本題討論補零對離散序列頻譜計算的影響?！緶仨嗵崾尽吭谟嬎汶x散非周期序列頻譜吋常用x2/7t作為橫坐標(biāo)，稱/2/k為歸一化頻率(normalized frequency)。 在畫頻譜時需給出橫坐標(biāo)。每幅圖下都需給出簡要的文字說明。由于離散非周期序列頻譜是周期的，所以在計算時不必用ffishift函數(shù)對fft計算的結(jié)果進(jìn)行重 新排列。【序列頻譜計算的基本方法】連續(xù)信號，通過其抽樣的離散信號，和離散信號的dft變換存在如下關(guān)系：n-x(ejq) n =dft氏) = x aot

3、cq=vzwa-=o通過如上關(guān)系，我們就可以通過dft來求信號的頻譜?！痉抡娼Y(jié)果】 dd k 、s os 矚0fm =0.9375fm =0.9219fm =0.9141fm =0.9102fm =0.9082【結(jié)果分析】增加dft的點數(shù)可以使頻譜更容易觀察，即減輕了柵欄效應(yīng)帶來的影響。頻譜的橫坐標(biāo)為歸一 化頻率，所以原信號的峰值第一次應(yīng)該出現(xiàn)在0.2處，隨著dft點數(shù)的增大，頻譜表示也越來越精 確。【自主學(xué)習(xí)內(nèi)容】1. 歸一化頻率相關(guān)知識。2. 通過matlab計算dft和matlab的繪圖操作?！鹃喿x文獻(xiàn)】1. 數(shù)字信號處理2. 補零對有限長序列頻譜及dft的影響。【發(fā)現(xiàn)問題】(專題研討或

4、相關(guān)知識點學(xué)習(xí)屮發(fā)現(xiàn)的問題)：dft點數(shù)的增多是否能提高頻譜分辨率？【問題探究】雖然dft點數(shù)增加使圖像更加細(xì)致，但是因為不論dft點數(shù)是多少，抽樣點數(shù)都是相同的，所 以每個頻譜所包含的信息相同，頻譜分辨率只與抽樣點數(shù)相關(guān)，與點數(shù)無關(guān)，頻譜分辨率相同。 所以不能通過增大dft點數(shù)而減少信息損失。dft點數(shù)的增多不能提過頻率分辨率。【仿真程序】 k = 0:31;n = 16;x = sin (0.2*pi*k);for i = 1:5n = 2*n;x = f ft (x, n); k = 0:n-1; subplot (5, 1, i); plot (k/n,abs(x); xlabel (

5、f/hz; ylabel (丨magnitude *); fm = find( x = max(x)/n; fmend2已知一離散序列為x k=acosfi)k+bcos (z)+a/2)z:)o用長度ams4的哈明窗對信號截短后近似計算其頻譜。試川不同的和z?的值(如和打近似相等，4和z?差距較大)，確定用哈明窗能分辯的最小的譜峰間隔a/2w =十c的值。n【題目分析】本題討論用哈明窗h算序列頻譜時的頻率分辨率?！痉抡娼Y(jié)果】c=4 a/b=16c=4 ab=420100201002010c=4 >vb=120100201002010020100152025c=3 ab=12c=4 ab

6、=1220100152025c=22 a/b=1220100152025201002010020100201002010v152025c=3 a/b=16152520101520c=22 ab=16252010201015252025152025152025c=2 a«=12c=2 a/b=1gc=2 ab=41520c=19 a/b=122010020100251520c=19 ab=16152025c:18 ./vb=122010025201002010020100152025c:1861520c:18 abm2010251520c=15 ab=161520c=15，a/b=12

7、2010020100«:20252520100z/15202s2010252015file editinsert took desktopwindow help、【結(jié)果分析】第一排是c=4,第二排c=3,遠(yuǎn)大于教材定義的標(biāo)準(zhǔn)c=2。無論a/b的比值是多少，區(qū)分譜峰毫 無壓力。第3排巾，c=2是教材定義的標(biāo)準(zhǔn)，a與b接近的時候，區(qū)分較為容易，當(dāng)a/b為1.6的時候, 幾乎是無法區(qū)分，而a/b=4的時候，完企無法區(qū)分。所以說在沒計使用哈明窗的時候，如果不同頻 譜的信號所占比例相差較大，哈明窗的c要比2大。從第一列可以看岀，如果a/b象近的時候，即 使c稍小于2,也是可以區(qū)分的。每一列中，

8、a/b的比值是相同的，隨著c的減小，譜峰的區(qū)別能力很容易看出是降低的。【自主學(xué)習(xí)內(nèi)容】不同窗函數(shù)的屏幕分辨率。【閱讀文獻(xiàn)】數(shù)字信號處理【發(fā)現(xiàn)問題】（專題研討或相關(guān)知識點學(xué)習(xí)中發(fā)現(xiàn)的問題）：哈明窗所起的作川？哈明窗與矩形窗誰的分辨率高？在仿真的時候c該選用怎樣的值？【問題探宄】哈明窗主瓣較寬，分辨率比矩形窗差，但是哈明窗可以使能呈集中在主瓣處。所以哈明窗可以 起到減小泄漏的作用。在仿真的時候c該選用怎樣的值？在開始的時候選用的c為1248 16,對比效果不明顯，發(fā)現(xiàn)c的取值過于大，反s嘗試確定在2附近多取值，其他值也不可取得過于大，因為c稍大實驗現(xiàn)象基本相同，所以沒有意義?！痉抡娉绦颉縩 = 6

9、4;k = 0:n-1;c_set = 4 3 2.2 2 1.9 1.8 1.5;a_set = 1 1.2 1.6 4;b = 1;wo = 0.5*pi;for c_index = 1:length(c_set) c = c_set(c_index); dw = c*2*pi/n;for a_index = 1:length(a_set)a = a_set(a_index); x = a*cos(wo*k)+b*cos(wo+dw)*k);xh = fft (x.*hamming(n) 1);subplot(length(c_set),length(a_set),length(a_set

10、)*(c_index-l)+a_index); plot(k,abs(xh); hold on;title ( 1c=1 num2str(c);a/b=1 num2str(a);axis ( 12 25 0 25);endend3 己知一離散序列為 x*=cos(/)+0.75cos(/2|z:)，0<*<63 其屮從二0.471，,=+71/64(1) 對x幻做64點fft,畫山此時信號的頻譜。(2) 如果(1)中顯示的譜不能分辨兩個譜峰，是否可對(1)中的64點信號補零而分辨出兩個譜 峰。通過編程進(jìn)行證實，并解釋其原因。(3) 給出一種能分辨出信號中兩個譜峰的計算方案，并進(jìn)行仿

11、真實驗?！绢}目分析】分析影響譜峰分辨率的主要因數(shù)，進(jìn)一步認(rèn)識補苓在在頻譜計算屮的作用?！痉抡娼Y(jié)果】 figure 1rue toix v»ewioo«5 uewcxopwindow ne<p> jj ji k c © 3 x -q on=4 l=45.380.390.40410 420 430 44l=256401 0 380 39040 410 420 430 44n=128 l=266恥 256 l=256n=256 l=32768【結(jié)果分析】頻譜的橫坐標(biāo)是歸一化頻率，圖示只觀察要區(qū)別的部分。每一個列的抽樣點相同，每一行dft 點相同。但看第一列，

12、雖然增加dft點數(shù)，甚至是增加到32768,也不能夠區(qū)分兩個譜峰，所以要想增 加頻率的分辨率，可以增加抽樣點數(shù)，第二列和第三列抽樣點數(shù)分別是128和256,都可以區(qū)分兩 個譜峰。【自主學(xué)習(xí)內(nèi)容】如何增加頻譜分辨率。【閱讀文獻(xiàn)】數(shù)字信號處理?！景l(fā)現(xiàn)問題】(專題研討或相關(guān)知識點學(xué)習(xí)中發(fā)現(xiàn)的問題)：如何計算所需抽樣點數(shù)？可以通過信號頻率的差值計算?！締栴}探究】1、2、3題討論的是離散信號頻譜的計算問題。與連續(xù)信號頻譜計算問題相比較，其計算誤差 有何不同？連續(xù)信號的頻率是非周期的，離散信號的頻譜是連續(xù)信號頻譜的周期話，可能會有混疊誤差?！痉抡娉绦颉縲0 = 0.4*pi;dw = pi/64;wl =

13、 wo+dw;n_set = 64 128 256;l_set = 64 128 256 512 32768;for n_index = 1:length(n_set);n = n_set(n_index); k = 0:n-l;x = cos (w0*k)+0.75*cos(wl*k); for l_index =1:length(l_set)l = l_set(l_index); if l<ncontinueendx = f ft (x, l); m = (o:l-1)*2/l;subplot(length(l一set),length(n_set),(lindex1*length(n

14、set)+nindex); plot(m,abs(x); axis (0.38 0.44 0 40);title(1n=1 num2str(n) l= num2str(l);endend4試用dft近似計算島斯信號以z) = cxp(-6/z2)的頻譜抽樣值。島斯信號頻譜的理論值為 g( j co)=通過與理論伉比較，討論信號的時域截取長度和抽樣頻率對計算誤差的影響?！绢}目分析】連續(xù)非周期信號頻譜計算的基本方法。計算中出現(xiàn)誤差的主要原因及減小誤差的方法。【仿真結(jié)果】filt ea vitwtook desktop window h«4p d k %b© a os cutli

15、m»=0 125cutlim«so 0625cimtime=0 25 fsam=16【結(jié)果分析】如圖所示，橫軸坐標(biāo)為角頻率，左斜線上的閣截取時間相同，垂線上的抽樣頻率相同，很容易能否看出，增加截取時間和抽樣頻率都可以減小誤差【自主學(xué)習(xí)內(nèi)容】連續(xù)非周期信號的頻譜計算方法?！痉抡娉绦颉縟 = pi;fsam_set =1248 16;n_set = 1248;l = 64 ;for fsam_index = 1:length(fsam_set); fsam = fsam_set(fsam_index);tsam = 1/fsam;for nindex = 1:length(ns

16、et);n = n_set(n_index); cuttime = n*tsam; k = -n/2:n/2-l; t = k*tsam; x = exp(-d*t.a2);x = tsam*fftshift(fft(xz l); w = 1inspace(-pi*fsam,pi*fsam,l);subplot(length(n_set),length(fsam_set),fsam_index+length(fsam_set)*(n_index -1);plot(w,abs(x); hold on;g = exp(-w.*w/(4*d); plot(w,g,'r)；title(cutt

17、ime=1 num2str(cuttime) fsam=1 num2str(fsam);endend擴展題5本題研究連續(xù)周期信號頻譜的近似計算w題。周期為7；的連續(xù)時間周期信號;v（0可用fourier級數(shù)表示為zi=oo其屮稱力連續(xù)時間周期信號a;的頻譜函數(shù)。叫=2兀/7； = 2兀乂稱力信號的基頻（基波），n叫稱力 信號的諧波。如果信號x(z)函數(shù)表達(dá)式已知，則可由積分得出信號的頻譜。如果信號x(0函數(shù)表達(dá)式未知，或 者40闌數(shù)表達(dá)式非常復(fù)雜，則很難巾積分得信號的頻譜。本題的fi的就是研宂如何利用dft近似 計算連續(xù)時間周期信號的頻譜。(1)若在信號的一個周期r0閃抽樣tv個點，即tq=

18、nt, r為抽樣周期(間隔)，可獲得序列對%/: = x(f)| z=at; z:二 0,1，tv -1試分析序列對幻的dft與連續(xù)時間周期信號x(0的頻譜汾似的關(guān)系；(2) 由(1)的結(jié)論，給出由dft近似計算周期信號頻譜的方案；(3) 周期信號x的周期7>1，x在區(qū)間0, 1的表達(dá)式為x(0=20r( 1 -z)4cos( 1 2tu)(a) 試畫出信號x(/)在區(qū)間0, 1的波形；(b) 若要用6次以a的諧波近似表示xw，試給出計算方案，并計算出近似表示的誤差。討論出現(xiàn)誤 差的原因及減小誤差的方法。【理論推導(dǎo)】dft計算所得結(jié)果xm與連續(xù)周期信號頻譜x(rtob)的關(guān)系x0= ix

19、cziq,)/=-00l=kt"=-oo2/r由于 = nt為t = -，于是得co.2/rk=nnk今n=m+rn;m=0，l,2,n 1，r為整數(shù)，上式可化為(m+rn)k'v-1<»y4眾=x+ /w)690)-c a，'nt=q /i=oo化簡可得n-«>刺zx(m + raqty0).e vzw=o /!=«>對比idft 2/r1 amaw=exm，e，amkm=0xmoo可得x (m + zw)6y0); m = 0,1，2,/v -1"=-oo【計算方案】根據(jù)理論推導(dǎo)結(jié)果設(shè)計近似計算方案。分析產(chǎn)

20、生誤差的主要原因。dft是連續(xù)信號頻譜的周期化，如果頻譜信號高頻信號較多，就會產(chǎn)生混疊，所以有誤差。如 果主要是低頻信號，周期話后互補產(chǎn)生影響，就沒有誤差?！緮U展分析】如果周期信號x（0是帶限信號，即信號的最高頻率分量為m吻（是正整數(shù)），試確定在一個周期閃 的最少抽樣點m使得在頻譜的計算過程當(dāng)屮不存在混疊誤差。與抽樣定理給出的結(jié)論比較，發(fā)表 你的看法。很容易看出，要想不產(chǎn)生混疊，必須7v22m+1。該結(jié)論與抽樣定理相似。因為周期信號是一種特殊的連續(xù)信號，所以抽樣定理對周期信號也是適用的?！痉抡娼Y(jié)果】n=16k=6diviation:mean:0.07187 max:0.1339diviatio

21、n:n=12k=6diviation:mean:0.1564 max:0.3297diviation:【結(jié)果分析】在n=16，k=6時，誤差最小，最大誤差力0.1229，平均誤差為0.07187。無論是減小n和減小k都很大的增大了誤差。為了分別測試n和k對誤差的影響，分別設(shè)k=6,求不同n吋候的誤差，n=16,不同k對誤差 的影響，分別得到了下圖（在原程序上稍加改動，所以在仿真程序中不再另寫了）在k=6的時候，誤差有明顯的減小?？梢娫瘮?shù)最大的高頻分量為6*0,根據(jù)前面推出的理論，n取6*2+1=13點的吋侯誤差減小。 實驗結(jié)果符合結(jié)論 【自主學(xué)習(xí)內(nèi)容】周期信號的頻譜與dft的關(guān)系?！景l(fā)現(xiàn)問題】原函數(shù)最高頻分量大約是多少？【問題探究】根據(jù)前面推理，最高頻分量為6*w0?！痉抡娉绦颉縩 = 1024;t = 1inspace (0,1,n);n_set = 12 16;h_set = 4 6;to = 1;count = 0;for n_index = 1:length(n_set)n = n_set(n_index); dt = to