數(shù)字信號處理實驗報告

1、實 驗 報 告實驗名稱 _課程名稱 院 系 部: 專業(yè)班級： 學(xué)生姓名： 學(xué) 號: 同 組 人: 實驗臺號: 指導(dǎo)教師： 成 績： 實驗日期: 華北電力大學(xué)（北京）一：實驗?zāi)康募耙?深刻理解離散時間系統(tǒng)的系統(tǒng)函數(shù)在分析離散系統(tǒng)的時域特性、頻域特性以及穩(wěn)定性中的重要作用及意義，熟練掌握利用MATLAB分析離散系統(tǒng)的時域響應(yīng)、頻響特性和零極點的方法。掌握利用DTFT和DFT確定系統(tǒng)特性的原理和方法。 二：實驗儀器儀器名稱規(guī)格/型號數(shù)量備注計算機一臺裝有Matlab三：實驗原理 MATLAB提供了許多可用于分析線性時不變連續(xù)系統(tǒng)的函數(shù)，主要包含有系統(tǒng)函數(shù)、系統(tǒng)時域響應(yīng)、系統(tǒng)頻域響應(yīng)等分析函數(shù)。1.

2、 離散系統(tǒng)的時域響應(yīng)離散時間LTI系統(tǒng)可用如下的線性常系數(shù)差分方程來描述： 已知輸入信號xk以及系統(tǒng)初始狀態(tài)y-1，y-2 .，就可以求出系統(tǒng)的響應(yīng)。MATLAB提供了計算差分方程的數(shù)值解的函數(shù)，可以計算上述差分方程描述的離散系統(tǒng)的單位脈沖響應(yīng)、零狀態(tài)響應(yīng)和完全響應(yīng)等。在調(diào)用MATLAB函數(shù)時，需要利用描述該離散系統(tǒng)的系數(shù)函數(shù)。對差分方程進行Z變換即可得系統(tǒng)函數(shù)： 在MATLAB中可使用向量a和向量b分別保存分母多項式和分子多項式的系數(shù)： 這些系數(shù)均從z0按z的降冪排列。離散系統(tǒng)的單位沖激響應(yīng)hk的計算 h,k = impz(b, a)：計算系統(tǒng)的單位脈沖響應(yīng)hk和相應(yīng)的時間向量k；也可簡寫為

3、：h = impz(b, a)。其中： h,k = impz(b, a, n)：計算n點單位脈沖響應(yīng)hk； 也可簡寫為：h = impz(b, a, n)。impz(b, a)：繪制單位脈沖響應(yīng)hk的圖形離散系統(tǒng)響應(yīng)yk的計算 y = filter(b, a, x)： 計算系統(tǒng)在輸入x作用下的零狀態(tài)響應(yīng)yk； y = filter(b, a, x, zi)：計算系統(tǒng)在輸入x和初始狀態(tài)作用下的完全響應(yīng)yk。zi是由系統(tǒng)的初始狀態(tài)經(jīng)過filtic函數(shù)轉(zhuǎn)換而得到的初始條件：zi= filtic (b, a, Y0) ，Y0為系統(tǒng)的初始狀態(tài)，Y0= y-1，y-2，y-3，.。2離散系統(tǒng)的系統(tǒng)函數(shù)零極

4、點分析離散LTI系統(tǒng)的系統(tǒng)函數(shù)H(z)可以表示為零極點形式：使用MATLAB提供的roots函數(shù)計算離散系統(tǒng)的零極點；使用zplane函數(shù)繪制離散系統(tǒng)的零極點分布圖。注意：在利用這些函數(shù)時，要求H(z)的分子多項式和分母多項式的系數(shù)的個數(shù)相等，若不等則需要補零。3離散系統(tǒng)的頻率響應(yīng) 當(dāng)離散因果LTI系統(tǒng)的系統(tǒng)函數(shù)H(z)的極點全部位于z平面單位圓內(nèi)時，系統(tǒng)的頻率響應(yīng)可由H(z)求出，即H, w=freqz(b, a, n)： 計算系統(tǒng)的n點頻率響應(yīng)H，w為頻率點向量。H=freqz(b, a, w) ：計算系統(tǒng)在指定頻率點向量w上的頻響；freqz(b,a)： 繪制頻率響應(yīng)曲線。其中：b和a分

5、別為系統(tǒng)函數(shù)H(z)的分子分母系數(shù)矩陣；4利用DTFT和DFT確定離散系統(tǒng)的特性 在很多情況下，需要根據(jù)LTI系統(tǒng)的輸入和輸出對系統(tǒng)進行辨識，即通過測量系統(tǒng)在已知輸入xk激勵下的響應(yīng)yk來確定系統(tǒng)的特性。 若系統(tǒng)的脈沖響應(yīng)為hk，由于存在yk=xk* hk，因而可在時域通過解卷積方法求解hk。 在實際應(yīng)用中，進行信號解卷積比較困難。因此，通常從頻域來分析系統(tǒng)，這樣就可以將時域的卷積轉(zhuǎn)變?yōu)轭l域的乘積，從而通過分析系統(tǒng)輸入序列和輸出序列的頻率特性來確定系統(tǒng)的頻率特性 ，再由 得到系統(tǒng)的脈沖響應(yīng)hk。四：實驗步驟第一：實驗內(nèi)容1.已知某LTI系統(tǒng)的差分方程為：（1）初始狀態(tài) ，輸入 計算系統(tǒng)的完全響

6、應(yīng)。（2）當(dāng)以下三個信號分別通過系統(tǒng)時，分別計算離散系統(tǒng)的零狀態(tài)響應(yīng)：（3）該系統(tǒng)具有什么特性？2.已知某因果LTI系統(tǒng)的系統(tǒng)函數(shù)為：（1）計算系統(tǒng)的單位脈沖響應(yīng)。（2）當(dāng)信號 通過系統(tǒng)時，計算系統(tǒng)的零狀態(tài)響應(yīng)。3.已知LTI系統(tǒng)的輸入輸出序列分別為（1）利用系統(tǒng)辨識原理從頻域分別求解系統(tǒng)的單位脈沖響應(yīng) 。（2）利用解析方法確定 ，并求出系統(tǒng)的單位脈沖響應(yīng) 。比較解析方法與系統(tǒng)辨識方法得到的系統(tǒng)單位脈沖響應(yīng)，分析誤差原因。4.已知某離散系統(tǒng)的輸入輸出序列。輸入序列：2，0.8333，0.3611，0.162，0.0748，0.0354，0.017，0.0083，0.0041，0.002，0.0

7、01，0.0005，0.0002，0.0001，0.0001，后面的數(shù)值均趨于0；輸出序列：0.0056，-0.0259，0.073，-0.1593，0.297， -0.4974，0.7711，-1.1267，1.5702，-2.1037，2.724，-3.4207，4.174，-4.9528，5.7117，-6.3889，6.9034，-7.1528，7.012，-6.3322，4.9416，-2.648，-0.7564，5.4872，-11.7557，19.7533，-29.6298，41.4666，-55.2433，70.7979，-87.7810（1）繪出輸入輸出信號的波形。（2）計算

8、該系統(tǒng)的頻率響應(yīng) ，并繪出其幅頻特性。（3）計算該系統(tǒng)的單位脈沖響應(yīng) ，并繪出其波形。第二：實驗結(jié)果，波形：第一題（1）%xk=uk求全響應(yīng)N=100; b=1,-1.143,0.412; a=0.0675,0.1349,0.0675;x=ones(1,N);zi=filtic(b,a,1,2);y=filter(b,a,x,zi)plot(y)xlabel('k');ylabel('幅值');title('N=100 xk=uk時的全響應(yīng)')（2）X1 作用時%N=100時在X1k作用下的響應(yīng)N=100; k=0:N-1;b=1,-1.143,

9、0.412; a=0.0675,0.1349,0.0675;x1=cos(pi*k/10);zi=filtic(b,a,1,2);y=filter(b,a,x1,zi)plot(y)xlabel('k');ylabel('幅值');title('N=100時 在X1k作用下的響應(yīng)')X2 作用時%N=100時在X2k作用下的響應(yīng)N=100; k=0:N-1;b=1,-1.143,0.412; a=0.0675,0.1349,0.0675;x2=cos(pi*k/5);zi=filtic(b,a,1,2);y=filter(b,a,x2,zi)pl

10、ot(y)xlabel('k');ylabel('幅值');title('N=100時 在X2k作用下的響應(yīng)')X3 作用時%N=100時在X3k作用下的響應(yīng)N=100; k=0:N-1;b=1,-1.143,0.412; a=0.0675,0.1349,0.0675;x3=cos(7*pi*k/10);zi=filtic(b,a,1,2);y=filter(b,a,x3,zi)plot(y)xlabel('k');ylabel('幅值');title('N=100時 在X3k作用下的波形)第二題(1)%單

11、位脈沖響應(yīng)N=40;k=0:N-1;a=1,-1.035,0.8246,-0.2605,0.04033;b=0.03571,0.1428,0.2143,0.1428,0.3571; y=impz(b,a,N); stem(y)xlabel('k');ylabel('幅值');title('N=40 單位脈沖響應(yīng)hk')波形圖為：(2)%零狀態(tài)響應(yīng)N=40;k=0:N-1;a=1,-1.035,0.8246,-0.2605,0.04033;b=0.03571,0.1428,0.2143,0.1428,0.3571; x=ones(N); y=fil

12、ter(b,a,x) plot(y);xlabel('k');ylabel('幅值');grid;title('N=40 零狀態(tài)響應(yīng)y0k曲線')波形圖為：第三題%3.1(a)求解單位脈沖響應(yīng)N=50; k=0:N-1;x=(0.5).k;y=0.25*(0.5).k+(0.5).k;X=fft(x,N);Y=fft(y,N);H=Y./X;h=ifft(H,N);plot(h,k);xlabel('k');ylabel('幅值Y');title('單位脈沖響應(yīng)3.1(a)')（2）%3.1(b)求

13、解單位脈沖響應(yīng)N=50; k=0:N-1;x=(0.25).k;y=(0.25).k-(0.25).(k-1);X=fft(x,N);Y=fft(y,N);H=Y./X;h=ifft(H,N);plot(h);xlabel('k');ylabel('幅值Y');title('單位脈沖響應(yīng)3.1(b)')第四題(1)%第四題 繪制輸入信號波形clearN=30;k=0:N;x=2.0000,0.8333,0.3611,0.1620,0.0748,0.0354,0.0170,0.0083,0.0041,0.0020,0.0010,0.0005,0.0

14、002,0.0001,0.0001,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000;plot(x);xlabel('k');ylabel('幅值');title('輸入信號波形');波形圖為：(2) %第四題 系統(tǒng)頻域響應(yīng)的波形%系統(tǒng)頻域響應(yīng)的波形clearN=30;k=0:N;x= 2.0000,0.8333,0.3611,0.1620,0.0748,0.0354,0.0

15、170,0.0083, 0.0041,0.0020,0.0010,0.0005,0.0002,0.0001,0.0001, 0.0000, 0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000;y=0.0056,-0.0259,0.073,-0.1593,0.297,-0.4974,0.7711,-1.1267,1.5702,-2.1037,2.724,-3.4207,4.174,-4.9528,5.7117,-6.3889,6.9034,

16、-7.1528,7.012,-6.3322,4.9416,-2.648,-0.7564,5.4872,-11.7557,19.7533,-29.6298,41.4666,-55.2433,70.7979,-87.7810;X=fft(x,N);Y=fft(y,N);H=Y./X;h=ifft(H,N);subplot(1,2,1);stem(H);xlabel('w');ylabel('幅值');title('系統(tǒng)頻域響應(yīng)的波形H');subplot(1,2,2);plot(abs(H);xlabel('角頻率(rad)');yl

17、abel('幅值');title('幅頻響應(yīng)曲線');%單位脈沖響應(yīng)的波形clearN=30;k=0:N;x= 2.0000,0.8333,0.3611,0.1620,0.0748,0.0354,0.0170,0.0083, 0.0041,0.0020,0.0010,0.0005,0.0002,0.0001,0.0001, 0.0000, 0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000,0.0000;y=0.0056,-0.0259,0.073,-0.1593,0.297,-0.4974,0.7711,-1.1267,1.5702,-2.1037,2.724,-3.4207,4.174,-4.9528,5.7117,-6.3889,6.9034,-7.1528,7.012,-6.3322,4.9416,-2.648,-0.7564,5.4872,-11.7557,19.7533,-29.6298,41.4666,-55.2433,70.7979,-87.7810;X=fft(x,N);Y=fft(y,N);H=Y./X;h