離散時間系統(tǒng)時域分析與仿真

1、鄭州航空工業(yè)管理學院電子信息系統(tǒng)仿真課程設(shè)計 2014 級 電子信息工程 專業(yè) 1413083 班級題 目 離散時間系統(tǒng)時域分析與仿真 姓 名 學號 二一六年十一月月二十五日一、 MATLAB軟件簡介 MATLAB是美國MathWorks公司出品的商業(yè)數(shù)學軟件，用于算法開發(fā)、數(shù)據(jù)可視化、數(shù)據(jù)分析以及數(shù)值計算的高級技術(shù)計算語言和交互式環(huán)境，主要包括MATLAB和Simulink兩大部分。二、 理論分析1、實驗原理離散時間系統(tǒng)其輸入、輸出關(guān)系可用以下差分方程描述： 當輸入信號為沖激信號時，系統(tǒng)的輸出記為系統(tǒng)單位沖激響應(yīng)，則系統(tǒng)響應(yīng)為如下的卷積計算式： 當hn是有限長度的（n：0，M）時，稱系統(tǒng)為F

2、IR系統(tǒng)；反之，稱系統(tǒng)為IIR系統(tǒng)。在MATLAB中，可以用函數(shù)y=Filter(p,d,x) 求解差分方程，也可以用函數(shù) y=Conv(x,h)計算卷積。2、時不變系統(tǒng)就是系統(tǒng)的參數(shù)不隨時間而變化，即不管輸入信號作用的時間先后，輸出信號響應(yīng)的形狀均相同，僅是從出現(xiàn)的時間不同。用數(shù)學表示為Tx(n)=yn則 Tx(n-n0)=yn-n0，這說明序列x(n)先移位后進行變換與它先進行變換后再移位是等效的。3、線性時不變系統(tǒng)既滿足疊加原理又具有時不變特性，它可以用單位脈沖響應(yīng)來表示。單位脈沖響應(yīng)是輸入端為單位脈沖序列時的系統(tǒng)輸出，一般表示為h(n)，即h(n)=T(n)。任一輸入序列x(n)的響應(yīng)

3、y(n)=Tx(n)=T (n-k)；由于系統(tǒng)是線性的，所以上式可以寫成y(n)=T(n-k)；又由于系統(tǒng)是時不變的，即有T(n-k)=h(n-k)；從而得y(n)=h(n-k)=x(n)*h(n)；這個公式稱為離散卷積，用“*”表示。4、線性時不變系統(tǒng)的性質(zhì)（1）、齊次性若激勵f(t)產(chǎn)生的響應(yīng)為y(t)，則激勵A(yù)f(t)產(chǎn)生的響應(yīng)即為Ay(t)，此性質(zhì)即為齊次性。其中A為任意常數(shù)。f(t)系統(tǒng)y(t)，Af(t)系統(tǒng)Ay(t)（2）、疊加性若激勵f1(t)與f2(t)產(chǎn)生的響應(yīng)分別為y1(t)， y2(t)，則激勵f1(t)+f2(t)產(chǎn)生的響應(yīng)即為y1(t)+y2(t)，此性質(zhì)稱為疊加性

4、。（3）、線性若激勵f1(t)與f2(t)產(chǎn)生的響應(yīng)分別為y1(t)， y2(t)，則激勵A(yù)1f1(t)+A2f2(t)產(chǎn)生的響應(yīng)即為A1y1(t)+A2y2(t)，此性質(zhì)稱為線性。（4）、時不變性若激勵f(t)產(chǎn)生的響應(yīng)為y(t)，則激勵f(t-t0)產(chǎn)生的響應(yīng)即為y(t-t0)，此性質(zhì)稱為不變性，也稱定常性或延遲性。它說明，當激勵f(t)延遲時間t0時，其響應(yīng)y(t)也延遲時間t0，且波形不變。（5）、微分性若激勵f(t)產(chǎn)生的響應(yīng)為y(t)，則激勵f'(t)產(chǎn)生的響應(yīng)即y(t),此性質(zhì)即為微分性。（6）、積分性若激勵f(t)產(chǎn)生的響應(yīng)為y(t)，則激勵f(t)的積分產(chǎn)生的響應(yīng)即為

5、y(t)的積分。此性質(zhì)稱為積分性。三、 MATLAB數(shù)值運算功能MATLAB 的強大數(shù)值計算功能使其在諸多數(shù)學計算軟件中傲視群雄，它是 MATLAB軟件的基礎(chǔ)。下面將介紹運用 MATLAB 計算卷積和沖激響應(yīng)。(1)卷積求解：由于系統(tǒng)的零狀態(tài)響應(yīng)是激勵與系統(tǒng)的單位取樣響應(yīng)的卷積，因此，卷積運算在離散時間信號處理領(lǐng)域被廣泛應(yīng)用。MATLAB求離散時間信號的卷積和命令為conv，其調(diào)用格式如下：y=conv(x,h)其中：x與h表示離散時間信號值的向量；y為卷積結(jié)果。用MATLAB進行卷積運算時，無法實現(xiàn)無限的累加，只能計算時限信號的卷積。一系統(tǒng)的單位取樣響應(yīng)為，用MATLAB求當激勵信號為x(n

6、)=u(n)-u(n-4)時，系統(tǒng)的零狀態(tài)響應(yīng)。 在MATLAB中可通過卷積求解零狀態(tài)響應(yīng)，即 x(n)* h( n)。上述系統(tǒng)描述 h(n)向量的長度至少為8，描述 x(n) 向量的長度至少為4，因此為了圖形完整美觀，將 h(n)向量和 x(n)加上一些附加的零值。其實現(xiàn)的MATLAB程序代碼如下：clear all;nx=-1:5; %x(n)向量顯示范圍（添加了附加的零值）nh=-2:10; %h(n)向量顯示范圍（添加了附加的零值）x=uDT(nx)-uDT(nx-4);h=0.8.nh.*(uDT(nh)-uDT(nh-8);y=conv(x,h);ny1=nx(1)+nh(1);

7、%卷積結(jié)果起始點%卷積結(jié)果長度為兩序列長度之和減1，即0到（length(nx)+length(nh)-2）%因此卷積結(jié)果的時間范圍是將上述長度加上起始點的偏移值ny=ny1+(0:(length(nx)+length(nh)-2);subplot(3,1,1);stem(nx,x,'fill');grid on;xlabel('n');title('x(n)');axis(-4 16 0 3);subplot(3,1,2);stem(nh,h,'fill');grid on;xlabel('n');title(&

8、#39;h(n)');axis(-4 16 0 3);subplot(3,1,3);stem(ny,y,'fill');grid on;xlabel('n');title('y(n)=x(n)*h(n)');axis(-4 16 0 3);其波形如下圖：卷積結(jié)果（2）沖激響應(yīng)求解在連續(xù)線性時域系統(tǒng)中，沖激響應(yīng)和階躍響應(yīng)是系統(tǒng)特性的描述，對它們的分析 是線性系統(tǒng)中極為重要的問題。輸入為單位沖激響應(yīng)函數(shù)(t) 所引起的零狀態(tài)響應(yīng)稱 為單位沖激響應(yīng)，簡稱為沖激響應(yīng)，用h(t) 表示；輸入為單位階躍函數(shù)u(t)所引起的 零狀態(tài)響應(yīng)稱為單位階躍響應(yīng)

9、，簡稱為階躍響應(yīng)，用g(t)表示。在MATLAB中，對于連續(xù)時域系統(tǒng)的沖激響應(yīng)和階躍響應(yīng)的數(shù)值解，可分別利用控制系統(tǒng)工具箱提供的函數(shù)impulse和step來求解。其調(diào)用格式如下：y=impulse(sys,t)y=step(sys,t)其中：t表示計算系統(tǒng)響應(yīng)的時間抽樣點向量，sys表示線性時域系統(tǒng)模型。已知一連續(xù)時域系統(tǒng)的微分方程如下：下面用MATLAB命令繪出0 t 4 范圍內(nèi)系統(tǒng)的沖激響應(yīng) h(t)和階躍響應(yīng) g(t)。其實現(xiàn)的MATLAB 代碼如下：clear all;t=0:0.001:4;sys=tf(1,16,1,2,32);h=impulse(sys,t);g=step(sy

10、s,t);subplot(2,1,1);plot(t,h);grid on;xlabel('時間/s');ylabel('h(t)');title('沖激響應(yīng)');subplot(2,1,2);plot(t,g);grid on;xlabel('時間/s');ylabel('g(t)');title('階躍響應(yīng)');其仿真波形如下圖：沖激響應(yīng)和階躍響應(yīng)波形四、 離散線性系統(tǒng)時域分析與仿真1、線性離散時間信號的基本運算對離散時間序列實行基本運算可得到新的序列，這些基本運算主要包括加、減、乘、除、移位和

11、反折等。兩個序列的加減乘除是對應(yīng)離散樣點值的加減乘除，因此，可通過MATLAB點乘和點除、序列移位和反折來實現(xiàn)，與連續(xù)時間信號處理方法基本一樣。下面列舉例子來說明。用MATLAB命令繪出下列離散時間序列的波形圖。設(shè)a=0.8,N=8。其實現(xiàn)的MATLAB程序代碼如下：clear all;a=0.8;N=8;n=-12:1:12;x=a.n.*(uDT(n)-uDT(n-N);n1=n;n2=n1+3;n3=n1-2;n4=-n1;subplot(4,1,1);stem(n1,x,'fill');grid on;axis(-15 15 0 1);xlabel('(a)x1

12、(n)');subplot(4,1,2);stem(n2,x,'fill');grid on;axis(-15 15 0 1);xlabel('(b)x2(n)');subplot(4,1,3);stem(n3,x,'fill');grid on;axis(-15 15 0 1);xlabel('(c)x3(n)');subplot(4,1,4);stem(n4,x,'fill');grid on;axis(-15 15 0 1);xlabel('(d)x4(n)');程序在運行過程中，調(diào)用

13、到自定義的沖激序列uDT.m文件，其源代碼如下：function y=uDT(n)y=n>=0;其仿真波形如下圖：離散序列波形圖2、線性離散系統(tǒng)分析(1)零狀態(tài)響應(yīng)離散線性時域系統(tǒng)可用下式表示，即MATLAB中的函數(shù)filter可對上式差分方程在指定時間范圍內(nèi)的輸入序列所產(chǎn)生的響應(yīng)進行求解。函數(shù)filter的調(diào)用格式如：y=filter(b,a,x)其中：x為輸入的離散序列；y為輸出的離散序列；y的長度和x的長度一樣；b和a分別為差分方程右端與左端的系數(shù)向量。一離散線性時域系統(tǒng)的差分方程如下： 下面用MATLAB命令繪出當激勵信號為時，該系統(tǒng)的零狀態(tài)響應(yīng)。其實現(xiàn)的MATLAB程序代碼如下

14、：clear all;a=3 -4 2;b=1 2;n=0:30;x=(1/2).n;y=filter(b,a,x);stem(n,y,'fill');grid on;xlabel('n');title('系統(tǒng)零狀態(tài)響應(yīng)y(n)');其仿真波形如下圖：圖3 零狀態(tài)響應(yīng)(2)單位取樣響應(yīng) 系統(tǒng)的單位取樣響應(yīng)定義為系統(tǒng)在（t）激勵下系統(tǒng)的零狀態(tài)響應(yīng)，用h（n）表示。MATLAB求解單位取樣響應(yīng)可利用函數(shù)filter，另一種求單位取樣響應(yīng)的方法是利用控制系統(tǒng)工具箱提供的函數(shù)impz來實現(xiàn)。有一離散線性系統(tǒng)的差分方程如下： 下面利用MATLAB的impz

15、函數(shù)繪出該系統(tǒng)的單位取樣響應(yīng)。其實現(xiàn)的程序代碼如下：clear all;a=3 -4 2;b=1 2;n=0:30;impz(b,a,30);grid on;title('系統(tǒng)單位取樣響應(yīng)h(n)');其仿真波形如下圖：圖4 單位取樣響應(yīng)(3)卷積和運算由于系統(tǒng)的零狀態(tài)響應(yīng)是激勵與系統(tǒng)的單位取樣響應(yīng)的卷積，因此，卷積運算在離散時間信號處理領(lǐng)域被廣泛應(yīng)用。MATLAB求離散時間信號的卷積和命令為conv，其調(diào)用格式如下：y=conv(x,h)其中：x與h表示離散時間信號值的向量；y為卷積結(jié)果。用MATLAB進行卷積運算時，無法實現(xiàn)無限的累加，只能計算時限信號的卷積。一系統(tǒng)的單位取

16、樣響應(yīng)為，用MATLAB求當激勵信號為 x（n）=u（n）-u（n-4）時，系統(tǒng)的零狀態(tài)響應(yīng)。 在MATLAB中可通過卷積求解零狀態(tài)響應(yīng)，即 x(n ) * h（ n）。上述系統(tǒng)描述h ( n) 向量的長度至少為8，描述x ( n) 向量的長度至少為4，因此為了圖形完整美觀，將 h(n ) 向量和 x(n )加上一些附加的零值。其實現(xiàn)的MATLAB程序代碼如下：clear all;nx=-1:5; %x(n)向量顯示范圍（添加了附加的零值）nh=-2:10; %h(n)向量顯示范圍（添加了附加的零值）x=uDT(nx)-uDT(nx-4);h=0.8.nh.*(uDT(nh)-uDT(nh-8

17、);y=conv(x,h);ny1=nx(1)+nh(1); %卷積結(jié)果起始點%卷積結(jié)果長度為兩序列長度之和減1，即0到（length(nx)+length(nh)-2）%因此卷積結(jié)果的時間范圍是將上述長度加上起始點的偏移值ny=ny1+(0:(length(nx)+length(nh)-2);subplot(3,1,1);stem(nx,x,'fill');grid on;xlabel('n');title('x(n)');axis(-4 16 0 3);subplot(3,1,2);stem(nh,h,'fill');grid

18、 on;xlabel('n');title('h(n)');axis(-4 16 0 3);subplot(3,1,3);stem(ny,y,'fill');grid on;xlabel('n');title('y(n)=x(n)*h(n)');axis(-4 16 0 3);其波形如下圖：卷積和求解波形圖 (4) 信號抽樣信號抽樣是連續(xù)時間信號分析向離散時間信號分析、連續(xù)信號處理向數(shù)字信號處理的第一步，廣泛應(yīng)用于實際的各類系統(tǒng)中。信號抽樣也稱為取樣或采樣，就是利用抽樣脈沖序列p ( t)從連續(xù)信號 f(t ) 中

19、抽取一系列的離散樣值，通過抽樣過程得到的 離散樣值信號稱為抽樣信號，用fs (t ) 表示。從數(shù)學上講，抽樣過程就是抽樣脈沖 p(t ) 和原連續(xù)信號f (t ) 相乘的過程，即:fs(t)=f(t)p(t)因此，可以用傅里葉變換的頻域卷積性質(zhì)來求抽樣信號 fs(t) 的頻譜。信號在時域被抽樣后，它的頻譜是原連續(xù)信號的頻譜以抽樣和角頻率為間隔周期的延拓，即信號在時域抽樣或離散化，相當于頻域周期化。升余弦脈沖信號為 下面用MATLAB編程實現(xiàn)該信號經(jīng)沖激脈沖抽樣后得到的抽樣信號 fs(t)及其頻譜。參數(shù) E =1 , = ,則。當采樣間隔Ts=1時，其實現(xiàn)的MATLAB程序代碼如下：clear

20、all;Ts=1;dt=0.1;t1=-4:dt:4;ft=(1+cos(t1)/2).*(uDT(t1+pi)-uDT(t1-pi);subplot(2,2,1);plot(t1,ft);grid on;axis(-4 4 -0.1 1.1);xlabel('Time(sec)');ylabel('f(t)');title('升余弦脈沖信號');N=500;k=-N:1:N;W=pi*k/(N*dt);Fw=dt*ft*exp(-j*t1'*W);subplot(2,2,2);plot(W,abs(Fw);grid on;axis(-1

21、0 10 -0.2 1.1*pi);title('升余弦脈沖信號的頻譜');t2=-4:Ts:4;fst=(1+cos(t2)/2).*(uDT(t2+pi)-uDT(t2-pi);subplot(2,2,3);plot(t1,ft,':');hold on;stem(t2,fst);grid on;axis(-4 4 -0.1 1.1);xlabel('Time(sec)');ylabel('fs(t)');title('抽樣后的信號');hold off;Fsw=Ts*fst*exp(-j*t2'*W)

22、;subplot(2,2,4);plot(W,abs(Fsw);grid on;axis(-10 10 -0.2 1.1*pi);xlabel('omega');ylabel('Fs(w)');title('抽樣信號的頻譜');其波形如下圖：信號抽樣波形圖3、離散系統(tǒng)的仿真已知某濾波器的差分方程：，求沖激響應(yīng)。系統(tǒng)由差分方程描述： 分別用filter函數(shù)和impz函數(shù)繪出其單位沖激響應(yīng)h(n)的線圖。MATLAB程序如下：clear all;close all;clc;b=1 -0.6;a=1 -1 0.9;n1=0:89;x=1zeros(1,

23、89);h1=filter(b,a,x);h,n=impz(b,a,90);subplot(2,1,1),stem(n1,h1,'.');xlabel('n');ylabel('h(n)');grid;title('filter函數(shù)求出的單位沖激響應(yīng)');axis(-5 95 -1.2 1.2);subplot(2,1,2),stem(n,h,'.');xlabel('n');ylabel('h(n)');grid;title('impz函數(shù)求出的單位沖激響應(yīng)');a

24、xis(-5 95 -1.2 1.2);set(gcf,'color','w'); filter函數(shù)和impz函數(shù)單位沖激響應(yīng)圖4、離散LTI的級聯(lián)在實際應(yīng)用中高階因果線性時不變系統(tǒng)可以用低階因果線性時不變系統(tǒng)級聯(lián)得到，這可簡化系統(tǒng)的設(shè)計與實現(xiàn)。例如，對于四階線性時不變系統(tǒng)，可以用二個二階系統(tǒng)級聯(lián)實現(xiàn)。第一級 第二級 用MATLAB語言編程驗證系統(tǒng)的級聯(lián)。matlab程序B1=1,0.9,0.8;A1=0.2,-0.2,0.4;xn=1,zeros(1,30);hn1=filter(B1,A1,xn);B2=1,0.7,0.85;A2=0.2,-0.5,0.3;

25、hn2=filter(B2,A2,hn1); n2=0:length(hn2)-1;subplot(2,1,1)stem(n2,hn2,'-'),title('級聯(lián)后的響應(yīng)')xlabel('n');ylabel('h2(n)')% the original serialsB3=1,1.6,2.28,1.325,0.68;A3=0.06,-0.19,0.27,-0.26,0.12;xn=1,zeros(1,30);hn=filter(B3,A3,xn);n=0:length(hn)-1;subplot(2,1,2)stem(n,hn,'.'), title('原始序列響應(yīng)')xlabel('n');ylabel('h(n)') 級聯(lián)響應(yīng)5、線性時不變系統(tǒng)的穩(wěn)定性分析若一個線性時不變系統(tǒng)的沖激響應(yīng)是絕對可和，則此系統(tǒng)就是BIBO的穩(wěn)定系統(tǒng)。由此，無限沖激響應(yīng)線性時不變系統(tǒng)穩(wěn)定的必要條件是，隨著輸入序列點的增加，沖激響應(yīng)衰減到零。用