北交大數字信號處理研究性學習報告

1、數字信號處理課程研究性學習報告姓名 學號 同組成員 指導教師 申艷 時間 2015年6月10日星期三 題目一：基本概念和技能學習報告【目的】(1) 掌握離散信號和系統(tǒng)時域、頻域和z域分析中的基本方法和概念；(2) 學會用計算機進行離散信號和系統(tǒng)時域、頻域和z域分析。(3) 培養(yǎng)學生自主學習能力，以及發(fā)現問題、分析問題和解決問題的能力。M1-1已知，以抽樣頻率對上述三個信號進行抽樣，在同一張圖上畫出，和及其抽樣點，對所得結果進行討論?！绢}目目的】1. 掌握抽樣的基本概念；2. 學會MATLAB中對信號抽樣的方法。【仿真結果】【結果分析】通過仿真抽樣前后的信號圖像我們可知，g1=cos(6*pi*

2、t),g2=cos(14*pi*t),g3=cos(26*pi*t)雖然它們信號時域表達式不同，但在以抽樣頻率fm=10Hz抽樣后，各離散信號具有相同的時域表達式，所以當抽樣頻率過小時，所得離散信號可能存在較大失真，無法反映原信號的特征，所以在給定一個信號后，我們應給出合適的抽樣頻率，抽樣頻率過小時失真較多，抽樣頻率過大時又造成了浪費，抽樣信號符合我們給定的要求即可，這樣得到的信號才能在最大程度保留原信號的特征的基礎上節(jié)約資源?！痉抡娉绦颉縦1=0:0.01:2;k2=0:0.1:2;g1=cos(6*pi*k1);g2=cos(6*pi*k2);subplot(3,2,1);plot(k1,

3、g1);axis(0,2,-1.5,1.5);title(g1=cos(6*pi*t) 原信號);xlabel(t);ylabel(g1(t);subplot(3,2,2);stem(k2,g2);axis(0,2,-1.5,1.5);title(g1=cos(6*pi*t) 抽樣頻率fm=10Hz);xlabel(k);ylabel(g1k); k3=0:0.01:2;k4=0:0.1:2;g3=cos(14*pi*k3);g4=cos(14*pi*k4);subplot(3,2,3);plot(k3,g3);axis(0,2,-1.5,1.5);title(g2=cos(14*pi*t)

4、原信號);xlabel(t);ylabel(g2(t);subplot(3,2,4);stem(k4,g4);axis(0,2,-1.5,1.5);title(g2=cos(14*pi*t) 抽樣頻率fm=10Hz);xlabel(k);ylabel(g2k); k5=0:0.001:2;k6=0:0.1:2;g5=cos(26*pi*k5);g6=cos(26*pi*k6);subplot(3,2,5);plot(k5,g5);axis(0,2,-1.5,1.5);title(g3=cos(26*pi*t) 原信號);xlabel(t);ylabel(g3(t);subplot(3,2,6)

5、;stem(k6,g6);axis(0,2,-1.5,1.5);title(g3=cos(26*pi*t) 抽樣頻率fm=10Hz);xlabel(k);ylabel(g3k);M1-2利用MATLAB的filter函數，求出下列系統(tǒng)的單位脈沖響應，并判斷系統(tǒng)是否穩(wěn)定。討論本題所獲得的結果?！绢}目目的】1. 掌握LTI系統(tǒng)單位脈沖響應的基本概念、系統(tǒng)穩(wěn)定性與單位脈沖響應的關系；2. 學會filter函數的使用方法及用filter函數計算系統(tǒng)單位脈沖響應；3. 體驗有限字長對系統(tǒng)特性的影響?！痉抡娼Y果】下面是兩個系統(tǒng)的單位脈沖響應：下面是兩個系統(tǒng)的零極點分布圖：極點1 0.9430 0.9020

6、極點2 1.0000 0.8500【結果分析】根據以上仿真結果可知，兩個系統(tǒng)函數H1，H2雖然只是分母系數差了一點，但最終實現時系統(tǒng)的單位脈沖響應會出現很大的不同，系統(tǒng)的極點也發(fā)生了不小的變化，這就造成了系統(tǒng)的失真。比如，完成數字濾波器設計后，具體實現時，可能會出現較大的誤差，甚至使設計的穩(wěn)定系統(tǒng)變成不穩(wěn)定系統(tǒng)，從而無法達到設計的要求。引起這些誤差的根本原因在于寄存器（存儲單元）的字長有限。誤差的特性與系統(tǒng)的類型、結構形式、數字的表示法、運算方式及字的長短有關。在通用計算機上，字長較長，量化步很小，量化誤差不大。但在專用硬件，如FPGA，實現數字系統(tǒng)時，其字長較短，就必須考慮有限字長效應了?！?/p>

7、問題探究】已知LTI系統(tǒng)的系統(tǒng)函數，有哪些計算系統(tǒng)單位脈沖響應方法，比較這些方法的優(yōu)缺點。 filter函數，可計算出差分方程的零狀態(tài)響應，既可以用來求yk,也可以求出hk，是最基本的方法，也是最復雜的方法;impulse函數，只是用來實現沖擊響應的，是最方便的方法；conv函數，是用來計算卷積的，可以用來求yk，這種方法也較麻煩，間接求出單位脈沖響應，不能直接求出單位脈沖響應?！痉抡娉绦颉縝1=1;a1=1,-1.845,0.850586;b2=1;a2=1,-1.85,0.85;x=1,zeros(1,49);y1=filter(b1,a1,x);subplot(2,1,1);stem(y

8、1);xlabel(x);ylabel(h1k);axis(0 50 0 8);y2=filter(b2,a2,x);subplot(2,1,2);stem(y2);xlabel(x);ylabel(h2k);subplot(2,1,1);r1,p1,m1=residuez(b1,a1); disp(極點1);disp(p1);zplane(b1,a1);title(H1);subplot(2,1,2);r2,p2,m2=residuez(b2,a2); disp(極點2);disp(p2);zplane(b2,a2);title(H2);M1-3(1)利用MATLAB語句x=firls(51

9、1,0 0.4 0.404 1,1 1 0 0）產生一個長度為512的序列xk，并畫出該序列的幅度頻譜。(2) 已知序列，分別畫出時序列yk的幅度頻譜。解釋所得到的結果。【題目目的】1. 學會用MATLAB函數freqz計算序列頻譜；2. 掌握序列頻譜的基本特性及分析方法。【溫磬提示】只需知道MATLAB語句x=firls(511,0 0.4 0.404 1,1 1 0 0）產生一個長度為512的序列xk，該序列滿足不需知道其他細節(jié)。用函數freqz計算該序列的頻譜，在畫幅度頻譜時，建議用歸一化頻率?！痉抡娼Y果】改進后：【結果分析】yk=對信號進行Fourier變換：所以，yk的幅度頻譜是xk

10、幅度頻譜左移與右移之和的一半，在平移的過程中，會出現混疊現象，然而這樣的理論分析與實際的結果不同，比如w0=0.8pi時，本來值較高的地方實際的值卻很低，原因如下：在xk變化較劇烈的點的值處與cos函數值為零的值相乘后消去，所以將高頻分量丟失，所以頻譜會出現理論與實際結果不同的現象。所以，是512的偶數倍的點均會出現這種現象，所以可以選擇將xk變成3000個點的序列，即可解決這種現象?！締栴}探究】有部分的計算結果可能與理論分析的結果不一致，分析出現該現象的原因，給出解決問題方法并進行仿真實驗?！痉抡娉绦颉縳=firls(511,0 0.4 0.404 1,1 1 0 0);b=x;a=1;w=

11、linspace(-pi,pi,512);y=freqz(b,a,w);subplot(3,2,1:2);plot(w/pi,abs(y);axis(-2,2,-1,2);title(生成序列xk的幅度頻譜);k=0:511;y1=x.*cos(0.4*pi*k);g1=freqz(y1,a,w);subplot(3,2,3);plot(w/pi,abs(g1);xlabel(w0=0.4pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y2=x.*cos(0.8*pi*k);g2=freqz(y2,a,w);subplot(3,2,4);plot(w/pi,abs(g2);xla

12、bel(w0=0.8pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y3=x.*cos(0.9*pi*k);g3=freqz(y3,a,w);subplot(3,2,5);plot(w/pi,abs(g3);xlabel(w0=0.9pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y4=x.*cos(pi*k);g4=freqz(y4,a,w);subplot(3,2,6);plot(w/pi,abs(g4);xlabel(w0=pi);ylabel(yk=xkcos(w0k)的幅度頻譜);改進后的程序：x=firls(2999,0 0.4 0.404 1,1 1

13、0 0);b=x;a=1;w=linspace(-pi,pi,512);y=freqz(b,a,w);subplot(3,2,1);plot(w/pi,abs(y);axis(-2,2,-1,2);title(生成序列xk的幅度頻譜);k=0:2999;y1=x.*cos(0.4*pi*k);subplot(3,2,2);plot(k,x);title(xk的時域圖像);g1=freqz(y1,a,w);subplot(3,2,3);plot(w/pi,abs(g1);xlabel(w0=0.4pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y2=x.*cos(0.8*pi*k)

14、;g2=freqz(y2,a,w);subplot(3,2,4);plot(w/pi,abs(g2);xlabel(w0=0.8pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y3=x.*cos(0.9*pi*k);g3=freqz(y3,a,w);subplot(3,2,5);plot(w/pi,abs(g3);xlabel(w0=0.9pi);ylabel(yk=xkcos(w0k)的幅度頻譜);y4=x.*cos(pi*k);g4=freqz(y4,a,w);subplot(3,2,6);plot(w/pi,abs(g4);xlabel(w0=pi);ylabel(yk=x

15、kcos(w0k)的幅度頻譜);M1-4已知，當，比較的幅度響應。從中你能得出什么結論？ 【題目目的】1. 掌握低通濾波器、高通濾波器、FIR濾波器和IIR濾波器的概念；2. 學會使用freqz函數?！痉抡娼Y果】【結果分析】由H1k的幅度響應圖線可知，H1代表的是數字低通FIR濾波器，當=0.8時，H2代表的是數字低通IIR濾波器，當=-0.8時，H2代表的系統(tǒng)是數字IIR濾波器，并且高頻分量幅度響應很大。當=0.8時，H3代表的是數字低通IIR濾波器，是H1和H2級聯的結果，當=-0.8時，H3代表的是數字全通IIR濾波器，幅度響應近似為1，是H1和H2級聯的結果。由此可知，當參數不同時，系

16、統(tǒng)會有較大的不同?！痉抡娉绦颉縝1=1,1;a1=2;w=linspace(-pi,pi,512);y1=freqz(b1,a1,w);subplot(3,1,1);plot(w/pi,abs(y1);axis(-1.5,1.5,-0.5,1.5);title(H1(z)的幅度響應);xlabel(歸一化頻率w/pi);ylabel(|H1|);subplot(3,1,2);b21=0.2;a21=1 -0.8;b22=1.8;a22=1 0.8;y21=freqz(b21,a21,w);y22=freqz(b22,a22,w);plot(w/pi,abs(y21),r,w/pi,abs(y2

17、2);title(H2(z)的幅度響應);xlabel(歸一化頻率w/pi);ylabel(|H2|);axis(-1.5,1.5,-2,12);b31=conv(b1,b21);a31=conv(a1,a21);b32=conv(b1,b22);a32=conv(a1,a22);y31=freqz(b31,a31,w);y32=freqz(b32,a32,w);subplot(3,1,3);plot(w/pi,abs(y31),r,w/pi,abs(y32);title(H3(z)的幅度響應);xlabel(歸一化頻率w/pi);ylabel(|H3|);axis(-1.5,1.5,-0.5

18、,1.5);M1-6已知一因果系統(tǒng)的H(z) 為試求出和H(z)具有相同幅度響應的最小相位系統(tǒng)Hmin(z)和最大相位系統(tǒng)Hmax(z)。畫出并比較H(z)、Hmin(z) 和Hmax(z)的相位響應?！绢}目目的】1. 掌握全通濾波器的基本特征和特性；2. 學會計算具有相同幅度響應的最小相位系統(tǒng)Hmin(z)和最大相位系統(tǒng)Hmax(z)?！緶仨嗵崾尽吭诋嬒到y(tǒng)的相位響應時，要把三個系統(tǒng)的相位響應畫在一張圖上，比較其解卷繞(unwrap)后的相位特點。【仿真結果】手算結果：Hmin(z)=，Hmax(z)=計算機自動求解的結果：【結果分析】根據圖像可以看出，在具有相同幅度響應的因果穩(wěn)定系統(tǒng)中，最小

19、相位系統(tǒng)具有最具有最小的相位滯后，最大相位系統(tǒng)具有最大的相位滯后，原系統(tǒng)的相位滯后介于兩者之間?！締栴}探究】手算時，如何找出最小相位系統(tǒng)Hmin(z)和最大相位系統(tǒng)Hmax(z)？用計算機自動求解時，如何找出最小相位系統(tǒng)Hmin(z)和最大相位系統(tǒng)Hmax(z)？ 你所用的算法是一樣的嗎？ 原理類似?！痉抡娉绦颉渴炙愫螅篵1=1.1,1.99,0.9;a1=1,1.55,0.6;b2=0.9,1.99,1.1;a2=a1;b3=1,2,0.99;a3=a2;w=linspace(0,pi,1000);h1=freqz(b1,a1,w);h2=freqz(b2,a2,w);h3=freqz(b3

20、,a3,w);plot(w,angle(h1),w,unwrap(angle(h2),r,w,angle(h3),g);axis(0,pi,-8,2);xlabel(w);ylabel(相位響應);計算機自動求解：b1=1,2,0.99;a1=1,1.55,0.6;w=linspace(0,pi,1000);h1=freqz(b1,a1,w);j1=roots(b1);j2=j1;j3=j1;if j1(1)1|j1(1)1|j1(2)N時則將原序列截短為N點序列再進行計算，若NM則將序列補零后在做N點的DFT。【仿真結果】L=32L=64L=128L=256L=512L=32L=128【結果

21、分析】對序列補零后再做DFT相當于增加了DFT的點數，即頻域抽樣的點數，而原離散序列沒有改變，其傅里葉變換結果也沒有改變，同時若DFT點數太少則獲得的頻譜信息過少，有可能會使得重要的頻率信息丟失。由結果可知，DFT點數越多，產生的離散譜中含有的信息也就越多，得到的頻譜能更好的反應原連續(xù)譜中的信息。在對離散序列用DFT做譜分析時，應當適當增加DFT的點數，以減小柵欄效應。DFT點數越多，則L越大，即fsam/N越小，那么顯示分辨率就越高。 【自主學習內容】已知幅度譜，用matlab求譜峰所對應的頻率值。利用：a,b=max(X);其中X是一個向量，返回的b值對應的是向量X中的最大值，返回的a值為

22、最大值b對應的索引?！鹃喿x文獻】1陳后金，薛健，胡健. 數字信號處理M. 北京：高等教育出版社，2006 .【發(fā)現問題】 (專題研討或相關知識點學習中發(fā)現的問題)：隨著DFT點數的增加，所獲得的頻譜圖像與實際值誤差減小，而譜峰值更精確，且其對應的譜峰頻率變大。而在用matlab模擬近似頻譜時，由于是離散的點構成的，無法正常表現出真實頻譜的圖像【問題探究】【仿真程序】 k=0:31; x=sin(0.2*pi*k); L=0:31; X=fft(x,32);%對函數進行32點的DFT抽樣，計算頻譜 plot(2*L/32,abs(X)%abs()是求絕對值 title(L=32時序列的頻譜) g

23、rid on %加網格線 Xmax=max(abs(X);%找到譜峰 y1=find(abs(X)=Xmax)-1;%找到譜峰對應的點,因為點從0開始所以要減1 m1=y1*2/32; disp(Xmax);14.5561 disp(m1);0.1875 1.8125k=0:31;x=sin(0.2*pi*k);L=0:63;X=fft(x,64);plot(2*L/64,abs(X)title(L=64時序列的頻譜)grid on Xmax=max(abs(X);%找到譜峰 y1=find(abs(X)=Xmax)-1;%找到譜峰對應的點,因為點從0開始所以要減一 m1=y1*2/64; d

24、isp(Xmax); 14.5561 disp(m1); 0.0938 0.9063 k=0:31;x=sin(0.2*pi*k);L=0:127;X=fft(x,128);plot(2*L/128,abs(X)title(L=128時序列的頻譜)grid on Xmax=max(abs(X);%找到譜峰 y1=find(abs(X)=Xmax)-1;%找到譜峰對應的點,因為點從0開始所以要減一 m1=y1*2/128;disp(Xmax); 14.5561 disp(m1); 0.0469 0.4531 k=0:31;x=sin(0.2*pi*k);L=0:255;X=fft(x,256);

25、plot(2*L/256,abs(X)title(L=256時序列的頻譜)grid on Xmax=max(abs(X);%找到譜峰 y1=find(abs(X)=Xmax)-1;%找到譜峰對應的點,因為點從0開始所以要減一 m1=y1*2/256;disp(Xmax); 14.5561 disp(m1); 0.0234 0.2266 k=0:31;x=sin(0.2*pi*k);L=0:511;X=fft(x,512);plot(2*L/512,abs(X)title(L=512時序列的頻譜)grid on Xmax=max(abs(X);%找到譜峰 y1=find(abs(X)=Xmax)

26、-1;%找到譜峰對應的點,因為點從0開始所以要減一 m1=y1*2/512;disp(Xmax); 14.5561 disp(m1); 0.0117 0.1133M2-3 已知一離散序列為 x k=AcosW0k+Bcos ( (W0+DW)k)。用長度N=64的哈明窗對信號截短后近似計算其頻譜。試用不同的A和B的取值，確定用哈明窗能分辨的最小的譜峰間隔中c的值?！绢}目分析】本題討論用哈明窗計算序列頻譜時的頻率分辨率。Hamming窗函數的幅值有中心向兩端逐漸減弱，因而其高頻分量明顯減小，頻譜中旁瓣的幅度較小，主瓣峰值與第一個旁瓣峰值相對衰減很大，hamming窗以增加主瓣寬度來降低旁瓣能量，

27、用hamming窗極端頻譜時要求能分辨的譜峰的間隔c/Tp=c*fs/N。針對不同的AB值，調整頻率的差值來比較哈明窗所分辨的最小的譜峰。利用控制變量的思想，首先固定AB值，比較不同頻率所造成的影響。然后固定頻率，比較不同的AB值所造成的影響。解題思路：XK加窗截斷，得到XnK，對其進行DFT計算，然后畫出頻譜。【仿真結果】Dw1=pi/32,C=1Dw2=pi/16,C=2Dw3=3*pi/32,C=3 dw4=pi/8,C=4 譜峰間隔相等時AB相差較大與AB近似相等兩種情況：A：B=1:10A:B=1：1【結果分析】將實驗結果與教材中定義的窗函數的有效寬度相比較，發(fā)表你的看法。（1） 當

28、AB近似時，N=3時就能有效的區(qū)別能分辨出譜峰間隔。C=2是教材中定義的哈明窗的有效寬度，在實際中有時候是區(qū)別分開的。從圖像可以看出C=3的時候就能完全分辨出來，C值越大，顯示的額頻譜越接近真實值。C=4是哈明窗的主瓣寬度，就能完全區(qū)別譜峰間隔。（2） 當AB相差較大的時候，采用相同的分析方法。我們可以確定大約C=3能夠有效的分辨出譜峰間隔?！咀灾鲗W習內容】幾種常見的窗函數的特點與應用。有關DFT頻譜分析，對連續(xù)離散信號的分析處理，頻譜分辨率，譜峰間隔的知識。【閱讀文獻】1陳后金，薛健，胡健. 數字信號處理M. 北京：高等教育出版社，2006 .【問題探究】 在離散序列頻譜計算中為何要用窗函數

29、？用不同的窗函數對計算結果有何影響？與矩形窗相比哈明窗有何特點?如何選擇窗函數？【仿真程序】A=B時，比較不同的譜峰間隔 N=64; k=0:63; L=0:511; A=1;B=1; dw1=pi/32; x1=A*cos(pi*k./4)+B*cos(pi*k./4+dw1*k);wh=(hamming(N); x1=x1.*wh; y1=fft(x1,512); figure;plot(L/512,abs(y1);xlabel(歸一化頻率); ylabel(幅值); title(dw1=pi/32,C=1) dw2=pi/16;x2=A*cos(pi*k./4)+B*cos(pi*k./

30、4+dw2*k);wh=(hamming(N);x2=x2.*wh;y2=fft(x2,512);figure;plot(L/512,abs(y2); xlabel(歸一化頻率);ylabel(幅值);title(dw2=pi/16,C=2)dw3=3*pi/32;x3=A*cos(pi*k./4)+B*cos(pi*k./4+dw3*k);wh=(hamming(N);x3=x3.*wh;y3=fft(x3,512);figure;plot(L/512,abs(y3); xlabel(歸一化頻率);ylabel(幅值);title(dw3=3*pi/32,C=3) dw4=pi/8;x4=A

31、*cos(pi*k./4)+B*cos(pi*k./4+dw4*k);wh=(hamming(N);x4=x4.*wh;y4=fft(x4,512);figure;plot(L/512,abs(y4); xlabel(歸一化頻率);ylabel(幅值);title(dw4=pi/8,C=4)A與B相差較大時 N=64;k=0:63;L=0:511;A=1;B=10;dw1=pi/32;x1=A*cos(pi*k./4)+B*cos(pi*k./4+dw1*k);wh=(hamming(N); x1=x1.*wh;y1=fft(x1,512);figure;plot(L/512,abs(y1);

32、xlabel(歸一化頻率); ylabel(幅值);title(dw1=pi/32,C=1);M2-4 已知一離散序列為 xk=cos(W0k)+0.75cos(W1k), 0 k 63 其中W0=0.4p, W1=W0+p/64(1) 對xk做64點FFT, 畫出此時信號的頻譜。(2) 如果(1)中顯示的譜不能分辨兩個譜峰，是否可對(1)中的64點信號補零而分辨出兩個譜峰。通過編程進行證實，并解釋其原因 。(3) 給出一種能分辨出信號中兩個譜峰的計算方案，并進行仿真實驗。【題目分析】分析影響譜峰分辨率的主要因素，進一步認識補零在在頻譜計算中的作用。分析影響譜峰分辨率的主要因數，進一步認識補領

33、在頻譜計算中的作用。影響譜峰分辨率的主要因數是序列的頻率，即連續(xù)時間的抽樣點數?！痉抡娼Y果】（1）（3） L=128L=512(3)64點信號改為128點信號，即能分辨出來。L=128L=256【結果分析】由結果可知，對信號增加零點不能提高頻譜的分辨率；為了提高頻譜的分辨率，可以增加序列的抽樣點數，比如開始的時候，K的范圍是0，64，后來K的范圍是0，128，然后是0，256隨著抽樣點數的增加，頻譜的分辨率逐漸提高。【自主學習內容】如何增加頻譜分辨率【閱讀文獻】1陳后金，薛健，胡健. 數字信號處理M. 北京：高等教育出版社，2006 .【問題探究】離散信號的頻譜具有周期性，會因抽樣頻率的不適而

34、產生混疊，導致失真?！痉抡娉绦颉浚?） 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=fftshift(fft(x); f=-pi:2*pi/L:pi-2*pi/L; plot(f/pi,abs(X) grid（2） k=0:63;L=128;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

35、*pi/L;plot(f/pi,abs(X)grid(3)k=0:127;L=128;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); title(128點DFT);gridk=0:255;L=256;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

36、=-pi:2*pi/L:pi-2*pi/L;plot(f/pi,abs(X);title(256點DFT);gridM2-6 試用DFT近似計算高斯信號的頻譜抽樣值。高斯信號頻譜的理論值為通過與理論值比較，討論信號的時域截取長度和抽樣頻率對計算誤差的影響?！绢}目分析】連續(xù)非周期信號頻譜計算的基本方法。計算中出現誤差的主要原因及減小誤差的方法。利用連續(xù)非周期信號頻譜與離散Fourier變換的關系。對于連續(xù)非周期信號，要對其進行頻譜計算，需要先經過抽樣，變?yōu)殡x散信號，再對離散信號的頻譜在一個周期內抽樣，即做DFT。由于題目所給信號時域無線長，所以為了在抽樣后頻譜不混疊，對信號進行截短即加抗混疊濾波

37、器后在進行頻譜分析，因為加濾波器后照成的誤差遠小于頻譜混疊照成的誤差。在對信號抽樣時還需要考慮抽樣頻率，以使抽樣后所得頻譜不會發(fā)生混疊。本題以d=0.05為例，討論fs、N變化帶來的影響；再討論d不同帶來的影響。【仿真結果】fs=1,N=4,d=0.05fs=1,N=6,d=0.05fs=1,N=8,d=0.05fs變化：fs=0.1;N=128;d=0.05fs=0.2;N=128;d=0.05fs=0.4;N=128;d=0.05d變化時：fs=1;N=128;d=0.05fs=1;N=128;d=0.5fs=1;N=128;d=2【結果分析】由于信號及頻譜都有理論表達式，在進行誤差分析時

38、希望給出一些定量的結果。當N逐漸變大時，譜線逐漸變?yōu)樽V峰，而且譜峰的寬度隨在一定范圍內變化時逐漸變窄；當f變化時，對DFT分析信號頻譜的精度影響較大，因為其直接影響頻譜混疊的精度；當d逐漸變大時，N合適大小時，需要改變抽樣頻率來使的近似值逼近理論曲線?！咀灾鲗W習內容】利用DFT分析連續(xù)非周期信號的頻譜的相關概念及分析過程中出現的現象。【閱讀文獻】1陳后金，薛健，胡健. 數字信號處理M. 北京：高等教育出版社，2006 .【發(fā)現問題】 (專題研討或相關知識點學習中發(fā)現的問題)：信號關于t=0偶對稱，在計算過程中如何處理這個問題。【問題探究】N，fs，d的取值問題【仿真程序】 fs=1; N=4;

39、 T=1/fs;ws=2*pi*fs; d=0.05; L=512; g=exp(-d*(-N+1:N-1)*T.*(-N+1:N-1)*T); X=fftshift(fft(g,L); w=(-ws/2+(0:L-1)*ws/L)/(2*pi); wt=linspace(-2*pi,2*pi,2048); G=1/T*sqrt(pi/d).*exp(-wt.*wt)/(4*d); plot(w,abs(X),wt/(2*pi),abs(G); legend(近似值,理論值) grid擴展題M2-7 本題研究連續(xù)周期信號頻譜的近似計算問題。 周期為T0的連續(xù)時間周期信號x(t)可用Fourie

40、r級數表示為其中X(nw0)稱為連續(xù)時間周期信號x(t)的頻譜函數。稱為信號的基頻（基波），稱為信號的諧波。如果信號x(t)函數表達式已知，則可由積分得出信號的頻譜。如果信號x(t)函數表達式未知，或者x(t)函數表達式非常復雜，則很難由積分得信號的頻譜。本題的目的就是研究如何利用DFT近似計算連續(xù)時間周期信號的頻譜。(1)若在信號x(t)的一個周期T0內抽樣N個點，即， T為抽樣周期(間隔)，可獲得序列xk試分析序列xk的DFT與連續(xù)時間周期信號x(t)的頻譜X(nw0)的關系；(2)由(1)的結論，給出由DFT近似計算周期信號頻譜X(nw0)的方案；(3)周期信號x(t)的周期T0=1，x

41、(t)在區(qū)間0，1的表達式為x(t)=20t2(1-t)4cos(12pt)(a)試畫出信號x(t)在區(qū)間0，1的波形；(b)若要用10次以內的諧波近似表示x(t)，試給出計算方案，并計算出近似表示的誤差。討論出現誤差的原因及減小誤差的方法。 【理論推導】DFT計算所得結果Xm與連續(xù)周期信號頻譜X(nw0)的關系?！居嬎惴桨浮扛鶕碚撏茖ЫY果設計近似計算方案。分析產生誤差的主要原因。Xm是連續(xù)周期信號頻譜X(n0)的周期化，在乘上比例因子N所得，由于對于一般連續(xù)周期信號的頻譜，當n時，X（n0）0,所以可以用Xm來近似較為低次的諧波分量，但仍會有一些高頻諧波分量在周期化是發(fā)生混疊。如果信號是帶

42、限信號，但抽樣時不滿足抽樣定理，也同樣會產生混疊誤差。用矩形窗截取了一個周期長度上的信號抽樣值，因此會產生頻率泄露。如果信號不是帶限信號，在分析其頻譜的時候，要截斷成時限的譜，此時會產生截斷誤差。在對信號做DFT時，相當于在時域做了周期化，則對應的頻譜進行離散化，就會產生柵欄效應。用有限次諧波合成信號，丟失的諧波有一定的誤差在信號的重建是，諧波分量選取的次數多少也對誤差有著較大的影響。對于一般的周期信號來說，通常是不帶限的，而我們只能選擇次數較低的諧波分量（占原有信號的主要能量）進行信號重建，這是不可避免的會發(fā)生泄漏，產生誤差?！緮U展分析】如果周期信號x(t)是帶限信號，即信號的最高頻率分量為

43、Mw0(是正整數)，試確定在一個周期內的最少抽樣點N，使得在頻譜的計算過程當中不存在混疊誤差。與抽樣定理給出的結論比較，發(fā)表你的看法。由上面分析可知，頻譜以N0周期化，原信號x(t)帶限，最高頻率分量為M0，當N2M時，頻譜周期化不會發(fā)生混疊，可以看出這個關系式與抽樣定理十分相像。周期信號同樣遵循抽樣定理，當抽樣頻率fsam2fm時，其中fsam=N/t0,fm=M0/2=M/T0，于是有N2M為不發(fā)生混疊的條件?！痉抡娼Y果】【結果分析】討論DFT點數對近似計算的影響，討論所取諧波項的多少對近似計算的影響。誤差分析要給出定量的結果，如平均誤差，最大誤差等?！痉抡娉绦颉縉=32; T0=1;n=

44、10; T=1/N; ts=(0:N-1)*T; xs=20.*ts.*ts.*(1-ts).4).*cos(12*pi*ts); X=fft(xs,N);C=X/N;t=linspace(0,1,1024);xt=20.*t.*t.*(1-t).4).*cos(12*pi*t);figure(1);plot(t,xt);xlabel(t)ylabel(xt)title(xt的時域波形);grid onfigure(2)k=0:N-1;stem(k,abs(X);xlabel(k)ylabel(X)title(xs的N點DFT);grid onfigure(3);xe=C(1)*ones(1,length(t); for k=1:n-1 xe=xe+C(k+1)*exp(j*k*2*pi/