1、第十章上機(jī)實(shí)驗(yàn)第十章上機(jī)實(shí)驗(yàn)數(shù)字信號處理是一門理論和實(shí)際密切結(jié)合的課程，為深入掌握課程內(nèi)容，最好在學(xué)習(xí)理論的同時，做習(xí)題和上機(jī)實(shí)驗(yàn)。 上機(jī)實(shí)驗(yàn)不僅可以幫助讀者深入的理解和消化基本理論，而且能鍛煉初學(xué)者的獨(dú)立解決問題的能力。 本章在第二版的基礎(chǔ)上編寫了六個實(shí)驗(yàn)， 前五個實(shí) 驗(yàn)屬基礎(chǔ)理論實(shí)驗(yàn)，第六個屬應(yīng)用綜合實(shí)驗(yàn)。實(shí)驗(yàn)一系統(tǒng)響應(yīng)及系統(tǒng)穩(wěn)定性。實(shí)驗(yàn)二時域采樣與頻域采樣。實(shí)驗(yàn)三用FFT對信號作頻譜分析。實(shí)驗(yàn)四IIR數(shù)字濾波器設(shè)計及軟件實(shí)現(xiàn)。實(shí)驗(yàn)五FIR數(shù)字濾波器設(shè)計與軟件實(shí)現(xiàn)實(shí)驗(yàn)六應(yīng)用實(shí)驗(yàn)一一數(shù)字信號處理在雙音多頻撥號系統(tǒng)中的應(yīng)用任課教師根據(jù)教學(xué)進(jìn)度，安排學(xué)生上機(jī)進(jìn)行實(shí)驗(yàn)。建議自學(xué)的讀者在學(xué)習(xí)完第一章

2、后作實(shí)驗(yàn)一；在學(xué)習(xí)完第三、四章后作實(shí)驗(yàn)二和實(shí)驗(yàn)三；實(shí)驗(yàn)四IIR數(shù)字濾波器設(shè)計及軟件實(shí)現(xiàn) 在。學(xué)習(xí)完第六章進(jìn)行； 實(shí)驗(yàn)五在學(xué)習(xí)完第七章后進(jìn)行。 實(shí)驗(yàn)六綜合實(shí)驗(yàn)在學(xué)習(xí)完第七章或 者再后些進(jìn)行；實(shí)驗(yàn)六為綜合實(shí)驗(yàn)，在學(xué)習(xí)完本課程后再進(jìn)行。10.1 實(shí)驗(yàn)一:系統(tǒng)響應(yīng)及系統(tǒng)穩(wěn)定性1 .實(shí)驗(yàn)?zāi)康?1)掌握求系統(tǒng)響應(yīng)的方法。(2)掌握時域離散系統(tǒng)的時域特性。(3)分析、觀察及檢驗(yàn)系統(tǒng)的穩(wěn)定性。2 .實(shí)驗(yàn)原理與方法在時域中，描寫系統(tǒng)特性的方法是差分方程和單位脈沖響應(yīng)，在頻域可以用系統(tǒng)函數(shù)描述系統(tǒng)特性。已知輸入信號可以由差分方程、單位脈沖響應(yīng)或系統(tǒng)函數(shù)求出系統(tǒng)對于該輸 入信號的響應(yīng)，本實(shí)驗(yàn)僅在時域求解。 在計算機(jī)

3、上適合用遞推法求差分方程的解，最簡單的方法是采用MATLAB語言的工具箱函數(shù) filter函數(shù)。也可以用MATLAB語言的工具箱函數(shù) conv函數(shù)計算輸入信號和系統(tǒng)的單位脈沖響應(yīng)的線性卷積，求出系統(tǒng)的響應(yīng)。系統(tǒng)的時域特性指的是系統(tǒng)的線性時不變性質(zhì)、因果性和穩(wěn)定性。重點(diǎn)分析實(shí)驗(yàn)系統(tǒng) 的穩(wěn)定性，包括觀察系統(tǒng)的暫態(tài)響應(yīng)和穩(wěn)定響應(yīng)。系統(tǒng)的穩(wěn)定性是指對任意有界的輸入信號，系統(tǒng)都能得到有界的系統(tǒng)響應(yīng)。或者系統(tǒng) 的單位脈沖響應(yīng)滿足絕對可和的條件。系統(tǒng)的穩(wěn)定性由其差分方程的系數(shù)決定。實(shí)際中檢查系統(tǒng)是否穩(wěn)定，不可能檢查系統(tǒng)對所有有界的輸入信號，輸出是否都是有界輸出，或者檢查系統(tǒng)的單位脈沖響應(yīng)滿足絕對可和的條件。

4、可行的方法是在系統(tǒng)的輸入端加入單位階躍序列，如果系統(tǒng)的輸出趨近一個常數(shù) (包括零)，就可以斷定系統(tǒng)是穩(wěn)定的19 o 系統(tǒng)的穩(wěn)態(tài)輸出是指當(dāng) n 時，系統(tǒng)的輸出。如果系統(tǒng)穩(wěn)定，信號加入系統(tǒng)后，系統(tǒng)輸 出的開始一段稱為暫態(tài)效應(yīng)，隨n的加大，幅度趨于穩(wěn)定，達(dá)到穩(wěn)態(tài)輸出。注意在以下實(shí)驗(yàn)中均假設(shè)系統(tǒng)的初始狀態(tài)為零。3 .實(shí)驗(yàn)內(nèi)容及步驟4 1)編制程序，包括產(chǎn)生輸入信號、 單位脈沖響應(yīng)序列的子程序，用filter函數(shù)或conv函數(shù)求解系統(tǒng)輸出響應(yīng)的主程序。程序中要有繪制信號波形的功能。5 2)給定一個低通濾波器的差分方程為y(n) 0.05x(n) 0.05x(n 1) 0.9y(n 1)輸入信號 xi(n

5、)R8(n)X2(n) u(n)a)分別求出系統(tǒng)對xi(n)Ra(n)和X2(n) u(n)的響應(yīng)序列，并畫出其波形。b)求出系統(tǒng)的單位沖響應(yīng)，畫出其波形。6 3)給定系統(tǒng)的單位脈沖響應(yīng)為hi(n)Ri0(n)h2(n)(n) 2.5 (n 1) 2.5 (n 2) (n 3)用線性卷積法分別求系統(tǒng)h1(n)和h2(n)對x(n) Rg(n)的輸出響應(yīng)，并畫出波形。7 4)給定一諧振器的差分方程為y(n) 1.8237 y(n 1) 0.9801y(n 2) b0x(n) b0x(n 2)令b01/100.49 ,諧振器的諧振頻率為 0.4rad。a)用實(shí)驗(yàn)方法檢查系統(tǒng)是否穩(wěn)定。輸入信號為u(

6、n)時，畫出系統(tǒng)輸出波形。b)給定輸入信號為x(n) sin(0.014n) sin(0.4n)求出系統(tǒng)的輸出響應(yīng)，并畫出其波形。4 .思考題(1)如果輸入信號為無限長序列，系統(tǒng)的單位脈沖響應(yīng)是有限長序列，可否用線性卷積 法求系統(tǒng)的響應(yīng)？如何求？(2)如果信號經(jīng)過低通濾波器，把信號的高頻分量濾掉，時域信號會有何變化，用前面第一個實(shí)驗(yàn)結(jié)果進(jìn)行分析說明。5.實(shí)驗(yàn)報告要求(1)簡述在時域求系統(tǒng)響應(yīng)的方法。(2)簡述通過實(shí)驗(yàn)判斷系統(tǒng)穩(wěn)定性的方法。分析上面第三個實(shí)驗(yàn)的穩(wěn)定輸出的波形。(3)對各實(shí)驗(yàn)所得結(jié)果進(jìn)行簡單分析和解釋。(4)簡要回答思考題。(5)打印程序清單和要求的各信號波形。10.1.2 實(shí)驗(yàn)程

7、序清單%實(shí)驗(yàn)1 ：系統(tǒng)響應(yīng)及系統(tǒng)穩(wěn)定性close all;clear all%=內(nèi)容1:調(diào)用filter解差分方程，由系統(tǒng)對 u(n)的響應(yīng)判斷穩(wěn)定性=A=1,-0.9;B=0.05,0.05;%系統(tǒng)差分方程系數(shù)向量B和Ax1n=1 1 1 1 1 1 1 1 zeros(1,50);% 產(chǎn)生信號 x1(n)=R8(n)x2n=ones(1,128);%產(chǎn)生信號 x2(n)=u(n)hn=impz(B,A,58); %求系統(tǒng)單位脈沖響應(yīng)h(n)subplot(2,2,1);y='h(n)'tstem(hn,y); %調(diào)用函數(shù) tstem 繪圖 title('系統(tǒng)單位脈沖

8、響應(yīng)h(n)');box ony1n=filter(B,A,x1n); % 求系統(tǒng)對 x1(n)的響應(yīng) y1(n) subplot(2,2,2);y='y1(n)'tstem(y1n,y);title('(b) 系統(tǒng)對 R8(n)的響應(yīng) y1(n)');box ony2n=filter(B,A,x2n); % 求系統(tǒng)對 x2(n)的響應(yīng) y2(n)subplot(2,2,4);y='y2(n)'tstem(y2n,y);title('系統(tǒng)對 u(n)的響應(yīng) y2(n)');box on%=內(nèi)容 2：調(diào)用 conv 函數(shù)計算

9、卷積=x1n=1 1 1 1 1 1 1 1 ;%產(chǎn)生信號 x1(n)=R8(n)h1n=ones(1,10) zeros(1,10);h2n=1 2.5 2.5 1 zeros(1,10);y21n=conv(h1n,x1n);y22n=conv(h2n,x1n);figure(2)subplot(2,2,1);y='h1(n)'tstem(h1n,y);%調(diào)用函數(shù) tstem 繪圖title('(d) 系統(tǒng)單位脈沖響應(yīng)h1(n)');box onsubplot(2,2,2);y='y21(n)'tstem(y21n,y);title('

10、;(e) h1(n)與 R8(n)的卷積 y21(n)');box onsubplot(2,2,3);y='h2(n)'tstem(h2n,y);%調(diào)用函數(shù) tstem 繪圖title('系統(tǒng)單位脈沖響應(yīng)h2(n)');box onsubplot(2,2,4);y='y22(n)'tstem(y22n,y);title('(g) h2(n)與 R8(n)的卷積 y22(n)');box on%=內(nèi)容 3：諧振器分析=un=ones(1,256);% 產(chǎn)生彳言號 u(n)n=0:255;xsin=sin(0.014*n)+s

11、in(0.4*n); % 產(chǎn)生正弦信號A=1,-1.8237,0.9801;B=1/100.49,0,-1/100.49;%系統(tǒng)差分方程系數(shù)向量B 和 Ay31n=filter(B,A,un);y32n=filter(B,A,xsin);%諧振器對u(n)的響應(yīng)y31(n)%諧振器對u(n)的響應(yīng)y31(n)figure(3)subplot(2,1,1);y='y31(n)'tstem(y31n,y);title('(h)諧振器對 u(n)的響應(yīng) y31(n)');box onsubplot(2,1,2);y='y32(n)'tstem(y32n

12、,y);title('(i) 諧振器對正弦彳t號的響應(yīng)y32(n)');box on10.1.3 實(shí)驗(yàn)程序運(yùn)行結(jié)果及分析討論程序運(yùn)行結(jié)果如圖10.1.1所示。實(shí)驗(yàn)內(nèi)容(2)系統(tǒng)的單位沖響應(yīng)、系統(tǒng)對x1(n)&(2和乂25)u(n)的響應(yīng)序列分別如圖(a)、(b)和(c)所示;實(shí)驗(yàn)內(nèi)容(3)系統(tǒng)h1(n)和h2(n)對x1(n)R8(n)的輸出響應(yīng)分別如圖(e)和(g)所示；實(shí)驗(yàn)內(nèi)容(4)系統(tǒng)對u(n)和x(n) sin(0.014n) sin(0.4n)的響應(yīng)序列分別如圖(h)和(i)所示。由圖(h)可見，系統(tǒng)對u(n)的響應(yīng)逐漸衰減到零，所以系統(tǒng)穩(wěn)定。由圖 (i)可見

13、,系統(tǒng)對x(n) sin(0.014n) sin(0.4n)的穩(wěn)態(tài)響應(yīng)近似為正弦序列sin(0.4n),這一結(jié)論驗(yàn)證了該系統(tǒng)的諧振頻率是0.4 rad。(a)系統(tǒng)單位脈沖響應(yīng)h(n)(b)系統(tǒng)對R8(n)的 響應(yīng)y1(n)(e) h1(n)與 R8(n)的卷積 y21(n)系統(tǒng)單位脈沖響應(yīng)h2(n)3 10(g) h2(n)與 R8(n)的卷積 y22(n)y(i)諧振器對正弦信號的響應(yīng)y32(n)10.5I 1-0.55010015020025030圖 10.1.110.1.4 簡答思考題(1)如果輸入信號為無限長序列，系統(tǒng)的單位脈沖響應(yīng)是有限長序列，可否用線性卷積 法求系統(tǒng)的響應(yīng)。對輸入信

14、號序列分段； 求單位脈沖響應(yīng) h(n)與各段的卷積；將各段 卷積結(jié)果相加。具體實(shí)現(xiàn)方法有第三章介紹的重疊相加法和重疊保留法。(2)如果信號經(jīng)過低通濾波器，把信號的高頻分量濾掉，時域信號的劇烈變化將被平滑，由實(shí)驗(yàn)內(nèi)容(1)結(jié)果圖10.1.1(a)、(b)和(c)可見，經(jīng)過系統(tǒng)低通濾波使輸入信號(n)、xi(n)®(n)和X2(n) u(n)的階躍變化變得緩慢上升與下降。10.2 實(shí)驗(yàn)二時域采樣與頻域采樣10.2.1 實(shí)驗(yàn)指導(dǎo)1 .實(shí)驗(yàn)?zāi)康臅r域采樣理論與頻域采樣理論是數(shù)字信號處理中的重要理論。要求掌握模擬信號采樣 前后頻譜的變化，以及如何選擇采樣頻率才能使采樣后的信號不丟失信息；要求掌握

15、頻率域采樣會引起時域周期化的概念，以及頻率域采樣定理及其對頻域采樣點(diǎn)數(shù)選擇的指導(dǎo)作用。2 .實(shí)驗(yàn)原理與方法時域采樣定理的要點(diǎn)是：a)對模擬信號Xa(t)以間隔T進(jìn)行時域等間隔理想采樣，形成的采樣信號的頻譜火(j )是原模擬信號頻譜 Xa(j )以采樣角頻率 s( s 2 /T)為周 期進(jìn)行周期延拓。公式為：丸(j ) FTXa(t) 1 Xa(j jn s) I nb)采樣頻率s必須大于等于模擬信號最高頻率的兩倍以上，才能使采樣信號的頻譜不產(chǎn)生頻譜混疊。利用計算機(jī)計算上式并不方便，下面我們導(dǎo)出另外一個公式，以便用計算機(jī)上進(jìn)行實(shí)驗(yàn)。理想采樣信號？a(t)和模擬信號xa(t)之間的關(guān)系為：?a(t

16、) Xa(t)(t nT)n對上式進(jìn)行傅立葉變換，得到：E(j ) Xa(t) (t nT)ejtdtn=xa(t) (t nT)e j tdtn在上式的積分號內(nèi)只有當(dāng)t nT時，才有非零值，因此：雙(j )Xa(nT)e j nTn上式中，在數(shù)值上 xa(nT) = x(n),再將 T代入，得至U：雙(j ) x(n)e j n n上式的右邊就是序列的傅立葉變換X(ej ),即只(j ) X(ej ) t上式說明理想采樣信號的傅立葉變換可用相應(yīng)的 采樣序列的傅立葉變換得到，只要將自變量3用 T代替即可。頻域采樣定理的要點(diǎn)是：a)對信號x(n)的頻譜函數(shù) X(ejC0)在0, 2兀上等間隔采樣

17、N點(diǎn)，得到XN(k) X(ej )2 kN , k 0,1,2,L ,N則N點(diǎn)IDFT Xn (k)得到的序列就是原序列 x(n)以N為周期進(jìn)行周期延拓后的主值區(qū) 序列，公式為：xN(n) IDFTXn(Qn x(n iN)R(n)b)由上式可知，頻域采樣點(diǎn)數(shù)N必須大于等于時域離散信號的長度M(即NW),才能使時域不產(chǎn)生混疊，則N點(diǎn)IDFT XN(k)得到的序列xN(n)就是原序列x(n),即xn (n) =x(n)。如果N>M , xN(n)比原序列尾部多 N-M個零點(diǎn)；如果 N<M , z則 xN(n)=IDFT Xn*)發(fā)生了時域混疊失真， 而且xN(n)的長度N也比x(n)

18、的長度 M短，因此。xn (n)與x(n)不相同。在數(shù)字信號處理的應(yīng)用中，只要涉及時域或者頻域采樣，都必須服從這兩個采樣理論的 要點(diǎn)。對比上面敘述的時域采樣原理和頻域采樣原理，得到一個有用的結(jié)論，這兩個采樣理論 具有對偶性：“時域采樣頻譜周期延拓，頻域采樣時域信號周期延拓”。因此放在一起進(jìn)行實(shí) 驗(yàn)。3 .實(shí)驗(yàn)內(nèi)容及步驟(1)時域采樣理論的驗(yàn)證。給定模擬信號，xa(t) Ae t sin( 0t)u(t)式中A=444.128,=50，2兀，o=50< u rad/s，它的幅頻特性曲線如圖10.2.1/Hz圖10.2.1xa(t)的幅頻特性曲線現(xiàn)用DFT(FFT)求該模擬信號的幅頻特性，以

19、驗(yàn)證時域采樣理論。安照xa(t)的幅頻特性曲線，選取三種采樣頻率，即Fs=1kHz, 300Hz, 200H乙觀測時間選Tp 50 ms。p為使用DFT ,首先用下面公式產(chǎn)生時域離散信號，對三種采樣頻率，采樣序列按順序用 x1 (n), x2(n) , x3(n)表示。x(n) xa(nT) Ae nT sin( °nT)u(nT)因?yàn)椴蓸宇l率不同，得到的Xi (n), x2(n) , x3(n)的長度不同，長度(點(diǎn)數(shù))用公式N Tp Fs計算。選FFT的變換點(diǎn)數(shù)為 M=64 ,序列長度不夠64的尾部加零。X(k)=FFT x(n) ,k=0,1,2,3,-, M-1 2 式中k代表

20、的頻率為k k oM要求：編寫實(shí)驗(yàn)程序，計算x1(n)、x2(n)和x3(n)的幅度特性，并繪圖顯示。觀察分析頻譜混疊失真。(2)頻域采樣理論的驗(yàn)證。給定信號如下：n 10n 13x(n) 27 n 14 n 260 其它編寫程序分別對頻譜函數(shù)X(ej ) FTx(n)在區(qū)間0,2 上等間隔采樣32和 16 點(diǎn)，得到 X32(k)和X16(k):X32(k) X(ej )2，k 0，1,2,L 31k32X16(k) X(ej )2, k 0,1,2,L 15、-k16再分別對X32(k)和X16(k)進(jìn)彳:亍32點(diǎn)和16點(diǎn)IFFT,得到x32(n)和x16(n)：X32(n)IFFTX32(

21、k)32 , n 0,1,2,L ,31x16(n) IFFTX16(k)，6n 0,1,2,L ,15分別畫出X(ej )、X32%)和X16(k)的幅度譜，并繪圖顯示 x(n)、x32(n)和x16(n)的波形，進(jìn)行對比和分析，驗(yàn)證總結(jié)頻域采樣理論。提示：頻域采樣用以下方法容易變程序?qū)崿F(xiàn)。 直接調(diào)用MATLABi數(shù)ft計算X32(k) FFTx(n)32就得到X(ej )在0,2 的32點(diǎn)頻率域采樣 抽取X32(k)的偶數(shù)點(diǎn)即可得到X(ej )在0,2 的16點(diǎn)頻率域采樣 X16(k),即X16(k) X32(2k) , k 0,1,2,L ,15。3；當(dāng)然也可以按照頻域采樣理論，先將信號

22、x(n)以16為周期進(jìn)行周期延拓， 取其主值 區(qū)(16點(diǎn))，再對其進(jìn)行16點(diǎn)DFT(FFT),得到的就是X(ej )在0,2 的16點(diǎn)頻率域采樣X16(k)o4 .思考題:如果序列x(n)的長度為M,希望得到其頻譜 X(ej )在0,2 上的N點(diǎn)等間隔采樣,當(dāng)N<M時，如何用一次最少點(diǎn)數(shù)的 DFT得到該頻譜采樣？5 .實(shí)驗(yàn)報告及要求a)運(yùn)行程序打印要求顯示的圖形，。b)分析比較實(shí)驗(yàn)結(jié)果，簡述由實(shí)驗(yàn)得到的主要結(jié)論c)簡要回答思考題d)附上程序清單和有關(guān)曲線。10.2.2 實(shí)驗(yàn)程序清單1時域采樣理論的驗(yàn)證程序清單%時域采樣理論驗(yàn)證程序exp2a.mTp=0.06;%觀察時間Tp=64毫秒%產(chǎn)

23、生M長采樣序列x(n)% Fs=1000;T=1/Fs;Fs=1000;T=1/Fs;M=Tp*Fs;n=0:M-1;A=444.128;alph=pi*50*2A0.5;omega=pi*50*2A0.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T);Xk=T*fft(xnt,M); %M 點(diǎn) FFTxnt) yn='xa(nT)'subplot(3,2,1);tstem(xnt,yn); %調(diào)用自編繪圖函數(shù) tstem繪制序列圖box on;title('(a) Fs=1000Hz');k=0:M-1;fk=k/Tp;subplo

24、t(3,2,2);plot(fk,abs(Xk);title('(a) T*FTxa(nT),Fs=1000Hz');xlabel('f(Hz)');ylabel('幅度');axis(0,Fs,0,1.2*max(abs(Xk) %=% Fs=300Hz和Fs=200Hz的程序與上面 Fs=1000Hz完全相同。2頻域采樣理論的驗(yàn)證程序清單%實(shí)驗(yàn)二(2頻域采樣理論的驗(yàn)證%產(chǎn)生M長三角波序列x(n),畫出它的波形M=27;n=0:M;xa=1:ceil(M/2); xb= floor(M/2):-1:0; xn=xa,xb;figure;stem

25、(n,xn,'o');axis(0,32,0,16);title('x(n)的波形'); xlabel('n'); ylabel('x(n)');%畫出x(n)的傅里葉變換結(jié)果Xk=fft(xn,1024);%1024點(diǎn)FFTx(n),它與序列x(n)的傅里葉變換近似k=0:1023;wk=2*k/1024;figure; plot(wk,abs(Xk);axis(0,1,0,200);title('FTx(n)'); xlabel('omega/pi'); ylabel('|X(eAjAo

26、mega)|')；X32k=fft(xn,32)%32 點(diǎn) FFTx(n)x32n=ifft(X32k);X16k=X32k(1:2:32);x16n=ifft(X16k,16);%對X32(k)進(jìn)彳:亍32點(diǎn)IFFT,得到x32(n)%抽取X32k的偶數(shù)點(diǎn)得到X16(K)%對X16(k)進(jìn)彳:亍16點(diǎn)IFFT,得到x16(n)%畫出16點(diǎn)頻率采樣結(jié)果k=0:15;figure; stem(k,abs(X16k),'o'); axis(0,16,0,200);title('16 點(diǎn)頻域采樣')； xlabel('k');%畫出x16(n)

27、的波形n=0:15; figure; stem(k,x16n,'o');title('X_1_6(n)的波形)；xlabel('n');ylabel('|X_1_6(k)|');axis(0,32,0,16);ylabel('x_1_6(n)');%畫出32點(diǎn)頻率采樣結(jié)果k=0:31; figure; stem(k,abs(X32k),'o'); axis(0,32,0,200);title('32 點(diǎn)頻域采樣'); xlabel('k');ylabel('|X_3_

28、2(k)|');%畫出x32(n)的波形n=0:31; figure; stem(n1,x32n,'o');axis(0,32,0,16);title(' X_3_2(n)的波形');xlabel('n');ylabel('x_3_2(n)');10.2.3實(shí)驗(yàn)程序運(yùn)行結(jié)果1時域采樣理論的驗(yàn)證程序運(yùn)行結(jié)果exp2a.m如圖10.3.2所示。由圖可見，采樣序列的頻譜的確是以采樣頻率為周期對模擬信號頻譜的周期延拓。當(dāng)采樣頻率為1000Hz時頻譜混疊很小；當(dāng)采樣頻率為 300Hz時，在折疊頻率150Hz附近頻譜混疊很嚴(yán)重；當(dāng)采樣

29、頻率 為200Hz時，在折疊頻率110Hz附近頻譜混疊更很嚴(yán)重。 Fs=1000Hz(a) TTTxa(nl).FS=1000Hz匚mm*50ao2060n(b) FS-300H2f(Hz)(切 TFT快版nD.FgpSOOH?luD100 50L o1015n圖 10.2.22時域采樣理論的驗(yàn)證程序exp2b.m運(yùn)行結(jié)果如圖10.3.3所示。44Xa0的波形該圖驗(yàn)證了頻域采樣理論和頻域采樣定理。對信號x(n)的頻譜函數(shù) X(ej3)在0,2兀上等間隔采樣N=16時，N點(diǎn)IDFT XN(k)得到的序列正是原序列x(n)以16為周期進(jìn)行周期延拓后的主值區(qū)序列：XN(n)IDFTXN(k)N x(

30、n iN)R(n)由于N<M ,所以發(fā)生了時域混疊失真，因此。XN(n)與x(n)不相同，如圖圖10.3.3(c)和(d)所示。當(dāng)N=32時，如圖圖10.3.3(c)和(d)所示，由于N>M ,頻域采樣定理，所以不存在時 域混疊失真，因此。xN(n)與x(n)相同。10.2.4簡答思考題先對原序列x(n)以N為周期進(jìn)行周期延拓后取主值區(qū)序列，xn (n) x(n iN)RN(n)再計算N點(diǎn)DFT則得到N點(diǎn)頻域采樣：XN(k) DFTxn (n)N =X(ej )2 k , k 0,1,2,L ,N 1k N10.3實(shí)驗(yàn)三：用FFT對信號作頻譜分析10.3.1 實(shí)驗(yàn)指導(dǎo)1 .實(shí)驗(yàn)?zāi)康?/p>

31、了解可能出現(xiàn)的分析學(xué)習(xí)用FFT對連續(xù)信號和時域離散信號進(jìn)行譜分析的方法,誤差及其原因，以便正確應(yīng)用FFT。2 .實(shí)驗(yàn)原理用FFT對信號作頻譜分析是學(xué)習(xí)數(shù)字信號處理的重要內(nèi)容。經(jīng)常需要進(jìn)行譜分析的信號是模擬信號和時域離散信號。對信號進(jìn)行譜分析的重要問題是頻譜分辨率D和分析誤差。頻譜分辨率直接和 FFT的變換區(qū)間N有關(guān)，因?yàn)镕FT能夠?qū)崿F(xiàn)的頻率分辨率是 2 /N,因 此要求2 /N D。可以根據(jù)此式選擇 FFT的變換區(qū)間No誤差主要來自于用 FFT作頻譜 分析時，得到的是離散譜，而信號(周期信號除外)是連續(xù)譜，只有當(dāng)N較大時離散譜的包絡(luò)才能逼近于連續(xù)譜，因此 N要適當(dāng)選擇大一些。周期信號的頻譜是離

32、散譜，只有用整數(shù)倍周期的長度作FFT,得到的離散譜才能代表周期信號的頻譜。如果不知道信號周期，可以盡量選擇信號的觀察時間長一些。對模擬信號進(jìn)行譜分析時，首先要按照采樣定理將其變成時域離散信號。如果是模擬 周期信號，也應(yīng)該選取整數(shù)倍周期的長度，經(jīng)過采樣后形成周期序列，按照周期序列的譜分析進(jìn)行。3 .實(shí)驗(yàn)步驟及內(nèi)容(1)對以下序列進(jìn)行譜分析。X(n) R4(n)n1,0n3x2(n)8n,4n70 ,其它n4n,0n3X3(n)n3,4n70,其它n選才i FFT的變換區(qū)間N為8和16兩種情況進(jìn)行頻譜分析。分別打印其幅頻特性曲線。 并進(jìn)行對比、分析和討論。(2)對以下周期序列進(jìn)行譜分析。x4(n)

33、 cos n44x5(n) cos( n/4) cos( n/8)選才i FFT的變換區(qū)間N為8和16兩種情況分別對以上序列進(jìn)行頻譜分析。分別打印 其幅頻特性曲線。并進(jìn)行對比、分析和討論。(3)對模擬周期信號進(jìn)行譜分析x6(t) cos8 t cos16 t cos20 t選擇采樣頻率Fs 64Hz ,變換區(qū)間 N=16,32,64三種情況進(jìn)行譜分析。分別打印其幅頻 特性，并進(jìn)行分析和討論。4 .思考題(1)對于周期序列，如果周期不知道，如何用FFT進(jìn)行譜分析? 如何選擇FFT的變換區(qū)間？（包括非周期信號和周期信號）（3）當(dāng)N=8時，X2（n）和X3（n）的幅頻特性會相同嗎？為什么？N=16呢

34、?5 .實(shí)驗(yàn)報告要求（1）完成各個實(shí)驗(yàn)任務(wù)和要求。附上程序清單和有關(guān)曲線。（2）簡要回答思考題。10.3.2實(shí)驗(yàn)程序清單%第10章實(shí)驗(yàn)3程序exp3.m%用FFT對信號作頻譜分析clear all;close all%實(shí)驗(yàn)內(nèi)容（1）=x1n=ones(1,4);%產(chǎn)生序列向量 x1(n)=R4(n)M=8;xa=1:(M/2);xb=(M/2):-1:1; x2n=xa,xb; %產(chǎn)生長度為 8 的三角波序列 x2(n)x3n=xb,xa;X1k8=fft(x1n,8);X1k16=fft(x1n,16);X2k8=fft(x2n,8);X2k16=fft(x2n,16);X3k8=fft(x

35、3n,8);X3k16=fft(x3n,16);%以下繪制幅頻特性曲線%計算x1n的8點(diǎn)DFT%計算x1n的16點(diǎn)DFT%計算x1n的8點(diǎn)DFT%計算x1n的16點(diǎn)DFT%計算x1n的8點(diǎn)DFT%計算x1n的16點(diǎn)DFTsubplot(2,2,1);mstem(X1k8); % 繪制 8 點(diǎn) DFT 的幅頻特性圖 title('(1a) 8 點(diǎn) DFTx_1(n)');xlabel(' 3 /無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X1k8)subplot(2,2,3);mstem(X1k16); % 繪制

36、 16 點(diǎn) DFT 的幅頻特性圖 title('(1b)16 點(diǎn) DFTx_1(n)');xlabel(' 3 / 無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X1k16)figure(2)subplot(2,2,1);mstem(X2k8); % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(2a) 8 點(diǎn) DFTx_2(n)');xlabel(' 3 /無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X2k8)subplot

37、(2,2,2);mstem(X2k16); % 繪制 16 點(diǎn) DFT 的幅頻特性圖 title('(2b)16 點(diǎn) DFTx_2(n)');xlabel(' 3 / 無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X2k16)subplot(2,2,3);mstem(X3k8); % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(3a) 8 點(diǎn) DFTx_3(n)');xlabel(' 3 /無');ylabel('幅度');axis(0,2,0,1.2*max

38、(abs(X3k8)subplot(2,2,4);mstem(X3k16); % 繪制 16 點(diǎn) DFT 的幅頻特性圖title('(3b)16 點(diǎn) DFTx_3(n)');xlabel(' 3 / 無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X3k16)%實(shí)驗(yàn)內(nèi)容周期序列譜分析=N=8;n=0:N-1;%FFT 的變換區(qū)間 N=8x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k8=fft(x4n);%計算 x4n 的 8 點(diǎn) DFTX5k8=fft(x5n);%計

39、算 x5n 的 8 點(diǎn) DFTN=16;n=0:N-1;%FFT 的變換區(qū)間 N=16x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n);%計算 x4n 的 16 點(diǎn) DFTX5k16=fft(x5n);%計算 x5n 的 16 點(diǎn) DFTfigure(3)subplot(2,2,1);mstem(X4k8); % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(4a) 8 點(diǎn) DFTx_4(n)');xlabel(' 3 /無');ylabel('幅度');axis(0,2,0

40、,1.2*max(abs(X4k8)subplot(2,2,3);mstem(X4k16); % 繪制 16 點(diǎn) DFT 的幅頻特性圖title('(4b)16 點(diǎn) DFTx_4(n)');xlabel(' 3 / 無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X4k16)subplot(2,2,2);mstem(X5k8); % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(5a) 8 點(diǎn) DFTx_5(n)');xlabel(' 3 /無');ylabel('幅度&

41、#39;);axis(0,2,0,1.2*max(abs(X5k8)subplot(2,2,4);mstem(X5k16); % 繪制 16 點(diǎn) DFT 的幅頻特性圖title('(5b)16 點(diǎn) DFTx_5(n)');xlabel(' 3 / 無');ylabel('幅度');axis(0,2,0,1.2*max(abs(X5k16)%實(shí)驗(yàn)內(nèi)容模擬周期信號譜分析=figure(4)Fs=64;T=1/Fs;N=16;n=0:N-1;%FFT 的變換區(qū)間 N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*

42、pi*n*T);%對 x6(t)16 點(diǎn)采樣X6k16=fft(x6nT);%計算 x6nT 的 16 點(diǎn) DFTX6k16=fftshift(X6k16);%將零頻率移到頻譜中心Tp=N*T;F=1/Tp;%頻率分辨率 Fk=-N/2:N/2-1;fk=k*F;%產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率(以零頻率為中心)subplot(3,1,1);stem(fk,abs(X6k16)，');box on % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(6a) 16 點(diǎn) |DFTx_6(nT)|');xlabel('f(Hz)');ylabel('

43、; 幅度');axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16)N=32;n=0:N-1;%FFT 的變換區(qū)間 N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%對 x6(t)32 點(diǎn)采樣X6k32=fft(x6nT);%計算 x6nT 的 32 點(diǎn) DFTX6k32=fftshift(X6k32);%將零頻率移到頻譜中心Tp=N*T;F=1/Tp;%頻率分辨率 Fk=-N/2:N/2-1;fk=k*F;%產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率（以零頻率為中心） subplot(3,1,2);stem

44、(fk,abs(X6k32)，');box on % 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(6b) 32 點(diǎn) |DFTx_6(nT)|');xlabel('f(Hz)');ylabel(' 幅度');axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32)N=64;n=0:N-1;%FFT 的變換區(qū)間 N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%對 x6(t)64 點(diǎn)采樣X6k64=fft(x6nT);%計算 x6nT 的 64 點(diǎn) D

45、FTX6k64=fftshift(X6k64);%將零頻率移到頻譜中心Tp=N*T;F=1/Tp;%頻率分辨率 Fk=-N/2:N/2-1;fk=k*F;%產(chǎn)生16點(diǎn)DFT對應(yīng)的采樣點(diǎn)頻率(以零頻率為中心)subplot(3,1,3);stem(fk,abs(X6k64)，'); box on% 繪制 8 點(diǎn) DFT 的幅頻特性圖title('(6a) 64 點(diǎn) |DFTx_6(nT)|');xlabel('f(Hz)');ylabel(' 幅度');axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64)10

46、.3.3 實(shí)驗(yàn)程序運(yùn)行結(jié)果實(shí)驗(yàn)3程序exp3.m運(yùn)行結(jié)果如圖10.3.1所示楮00.511.62ki/iii口制口點(diǎn)口辦(2b)1E DFT|x3 Iihr 6 a>rpbyi6的KF叼時JuO/TT口百口點(diǎn)口FTk3gJGd/T-I1051 10-30-101C2030Lf(H2)miT巧的柘旦13 rH&tm：5 iUDIT網(wǎng)16點(diǎn)DFTm利川1O(64點(diǎn)。曰>(前3Inffl20口1-J0，* * *2030圖 10.3.1程序運(yùn)行結(jié)果分析討論:請讀者注意，用 DFT (或FFT)分析頻譜，繪制頻譜圖時，最好將X(k)的自變量k換算成對應(yīng)的頻率，作為橫坐標(biāo)便于觀察頻譜

47、。k k, k 0,1,2,L ,N 1 N為了便于讀取頻率值，最好關(guān)于兀歸一化，即以 / 作為橫坐標(biāo)。1、實(shí)驗(yàn)內(nèi)容(1)圖(1a)和(1b)說明x1(n)R4 (n)的8點(diǎn)dft和16點(diǎn)dft分別是x1( n)的頻譜函數(shù)的8點(diǎn)和16點(diǎn)采樣；因?yàn)?x3(n) x2(n 3)8 R8(n),所以，x3(n" x2(n)的 8點(diǎn) dft 的模相等,如圖(2a)和(3a)。但是，當(dāng)N=16時，乂3(門)與乂2(門)不滿足循環(huán)移位關(guān)系，所以圖（2b）和（3b）的模不同。2、實(shí)驗(yàn)內(nèi)容（2）,對周期序列譜分析x4（n）cosn的周期為8,所以n=8和N=i6均是其周期的整數(shù)倍，得到正確4的單一頻

48、率正弦波的頻譜，僅在 0.25兀處有1根單一譜線。如圖（4b）和（4b）所示。x5（n）cos（ n/4） cos（ n/8）的周期為16,所以n=8不是其周期的整數(shù)倍，得到的頻譜不正確，如圖（5a）所示。N=16是其一個周期，得到正確的頻譜，僅在0.25兀和0.125兀處有2根單一譜線，如圖（5b）所示。3、實(shí)驗(yàn)內(nèi)容（3）,對模擬周期信號譜分析x6（t） cos8 t cos16 t cos20 tx6（t）有 3個頻率成分，f1 4Hz, f28Hz, f310Hz。所以 x6（t）的周期為0.5s。采樣頻率F 64Hz16 fl8 f26.4 f3。變換區(qū)間N=16 時，觀察S123時間

49、Tp=16T=0.25s,不是x6（t）的整數(shù)倍周期，所以所得頻譜不正確，如圖（ 6a）所示。變換區(qū)間N=32,64時，觀察時間Tp=0.5s, 1s,是x6（t）的整數(shù)周期，所以所得頻譜正確，如圖（6b）和（6c）所示。圖中3根譜線正好位于 4Hz,8Hz,10Hz處。變換區(qū)間N=64時頻譜幅度是變換區(qū)間 N=32時2倍，這種結(jié)果正好驗(yàn)證了用DFT對中期序列譜分析的理論。注意：（1）用DFT （或FFT）對模擬信號分析頻譜時，最好將 X（k）的自變量k換算成對應(yīng)的 模擬頻率fk,作為橫坐標(biāo)繪圖，便于觀察頻譜。這樣，不管變換區(qū)間N取信號周期的幾倍，畫出的頻譜圖中有效離散諧波譜線所在的頻率值不變

50、，如圖（6b）和（6c）所示。fk Mk Nfk T1k, pk 0,1,2,L ,N 1（2）本程序直接畫出采樣序列N點(diǎn)DFT的模值，實(shí)際上分析頻譜時最好畫出歸一化幅度譜，這樣就避免了幅度值隨變換區(qū)間N變化的缺點(diǎn)。本實(shí)驗(yàn)程序這樣繪圖只要是為了驗(yàn)證了用DFT對中期序列譜分析的理論。10.3.4 簡答思考題思考題（1）和（2）的答案請讀者在教材 3.?節(jié)找，思考題（3）的答案在程序運(yùn)行結(jié) 果分析討論已經(jīng)詳細(xì)回答。10.4 實(shí)驗(yàn)四IIR數(shù)字濾波器設(shè)計及軟件實(shí)現(xiàn)10.4.1 實(shí)驗(yàn)指導(dǎo)1 .實(shí)驗(yàn)?zāi)康?1)熟悉用雙線性變換法設(shè)計IIR數(shù)字濾波器的原理與方法；(2)學(xué)會調(diào)用 MATLAB信號處理工具箱中濾

51、波器設(shè)計函數(shù)(或?yàn)V波器設(shè)計分析工具 fdatool )設(shè)計各種IIR數(shù)字濾波器，學(xué)會根據(jù)濾波需求確定濾波器指標(biāo)參數(shù)。(3)掌握IIR數(shù)字濾波器的 MATLA改現(xiàn)方法。(3)通過觀察濾波器輸入輸出信號的時域波形及其頻譜，建立數(shù)字濾波的概念。2 .實(shí)驗(yàn)原理設(shè)計IIR數(shù)字濾波器一般采用間接法(脈沖響應(yīng)不變法和雙線性變換法)，應(yīng)用最廣泛的是雙線性變換法?；驹O(shè)計過程是：先將給定的數(shù)字濾波器的指標(biāo)轉(zhuǎn)換成過渡模擬濾波 器的指標(biāo)；設(shè)計過渡模擬濾波器；將過渡模擬濾波器系統(tǒng)函數(shù)轉(zhuǎn)換成數(shù)字濾波器的系 統(tǒng)函數(shù)。MATLAB言號處理工具有f中的各種IIR數(shù)字濾波器設(shè)計函數(shù)都是采用雙線性變換法。第六章介紹的濾波器設(shè)計函

52、數(shù)butter、cheby1、cheby2和ellip 可以分別被調(diào)用來直接設(shè)計巴特沃斯、切比雪夫 1、切比雪夫2和橢圓模擬和數(shù)字濾波器。本實(shí)驗(yàn)要求讀者調(diào)用如 上函數(shù)直接設(shè)計IIR數(shù)字濾波器。本實(shí)驗(yàn)的數(shù)字濾波器的 MATLA取現(xiàn)是指調(diào)用MATLAB言號處理工具箱函數(shù) filter 對給 定的輸入信號x(n)進(jìn)行濾波，得到濾波后的輸出信號y(n)。3 .實(shí)驗(yàn)內(nèi)容及步驟(1)調(diào)用信號產(chǎn)生函數(shù)mstg產(chǎn)生由三路抑制載波調(diào)幅信號相加構(gòu)成的復(fù)合信號st ,該函數(shù)還會自動繪圖顯示 st的時域波形和幅頻特性曲線，如圖10.4.1所示。由圖可見，三路信號時域混疊無法在時域分離。但頻域是分離的，所以可以通過濾波

53、的方法在頻域分離， 這就是本實(shí)驗(yàn)的目的。即)的波形w 0.50200 叩。600 aoa ioao 120口 uoo ieoo 1300 2000 ffH20圖10.4.1 三路調(diào)幅信號st的時域波形和幅頻特性曲線(2)要求將st中三路調(diào)幅信號分離，通過觀察 st的幅頻特性曲線，分別確定可以分離st中三路抑制載波單頻調(diào)幅信號的三個濾波器(低通濾波器、帶通濾波器、高通濾波器)的通帶截止頻率和阻帶截止頻率。要求濾波器的通帶最大衰減為0.1dB,阻帶最小衰減為60dB。提示：抑制載波單頻調(diào)幅信號的數(shù)學(xué)表示式為1s(t) cos(2 fot)cos(2 fct) cos(2 (% f0)t) cos(

54、2 (i f°)t) 2其中，cos(2 fct)稱為載波，fc為載波頻率，cos(2 fot)稱為單頻調(diào)制信號，f0為調(diào)制正弦波信號頻率，且滿足fc foo由上式可見，所謂抑制載波單頻調(diào)幅信號，就是 2個正弦信號相乘，它有2個頻率成分：和頻fc fo和差頻fc fo ,這2個頻率成分關(guān)于載波頻率fc對稱。所以，1路抑制載波單頻調(diào)幅信號的頻譜圖是關(guān)于載波頻率 fc對稱的2根譜線，其 中沒有載頻成分，故取名為抑制載波單頻調(diào)幅信號。容易看出，圖 10.4.1中三路調(diào)幅信號 的載波頻率分別為 250Hz、500Hz、1000Hz。如果調(diào)制信號 m(t)具有帶限連續(xù)頻譜，無直流 成分，則s(

55、t) m(t)cos(2 fct)就是一般的抑制載波調(diào)幅信號。其頻譜圖是關(guān)于載波頻率fc對稱的2個邊帶(上下邊帶)，在專業(yè)課通信原理中稱為雙邊帶抑制載波(DSB-SC)調(diào)幅信號，簡稱雙邊帶(DSB)信號。如果調(diào)制信號 m(t)有直流成分，則s(t) m(t)cos(2 fct)就 是一般的雙邊帶調(diào)幅信號。其頻譜圖是關(guān)于載波頻率fc對稱的2個邊帶(上下邊帶)，并包含載頻成分。(3)編程序調(diào)用MATLA雎波器設(shè)計函數(shù)ellipord 和ellip 分別設(shè)計這三個橢圓濾波 器，并繪圖顯示其幅頻響應(yīng)特性曲線。(4)調(diào)用濾波器實(shí)現(xiàn)函數(shù) filter ,用三個濾波器分別對信號產(chǎn)生函數(shù)mstg產(chǎn)生的信號st進(jìn)行濾波，分離出st中的三路不同載波頻率的調(diào)幅信號 y1(n)、y2(n)和y3(n),并繪 圖顯示y1(n)、y2(n)和y3(n)的時域波形，觀察分離效果。4 .信號產(chǎn)生函數(shù)mstg