信號(hào)與系統(tǒng)及MATLAB實(shí)現(xiàn)

1、信號(hào)與系統(tǒng)及MATLAB實(shí)現(xiàn)實(shí)驗(yàn)指導(dǎo)書前 言長(zhǎng)期以來(lái)，信號(hào)與系統(tǒng)課程一直采用單一理論教學(xué)方式，同學(xué)們依靠做習(xí)題來(lái)鞏固和理解教學(xué)內(nèi)容，雖然手工演算訓(xùn)練了計(jì)算能力和思維方法，但是由于本課程數(shù)學(xué)公式推導(dǎo)較多，概念抽象，常需畫各種波形，作題時(shí)難免花費(fèi)很多時(shí)間，現(xiàn)在，我們給同學(xué)們介紹一種國(guó)際上公認(rèn)的優(yōu)秀科技應(yīng)用軟件MATLAB，借助它我們可以在電腦上輕松地完成許多習(xí)題的演算和波形的繪制。 MATLAB的功能非常強(qiáng)大，我們此處僅用到它的一部分，在后續(xù)課程中我們還會(huì)用到它，在未來(lái)地科學(xué)研究和工程設(shè)計(jì)中有可能繼續(xù)用它，所以有興趣的同學(xué)，可以對(duì)MATLAB再多了解一些。MATLAB究竟有那些特點(diǎn)呢？1高效的數(shù)值

2、計(jì)算和符號(hào)計(jì)算功能，使我們從繁雜的數(shù)學(xué)運(yùn)算分析中解脫出來(lái)；2完備的圖形處理功能，實(shí)現(xiàn)計(jì)算結(jié)果和編程的可視化；3友好的用戶界面及接近數(shù)學(xué)表達(dá)式的自然化語(yǔ)言，易于學(xué)習(xí)和掌握；4功能豐富的應(yīng)用工具箱，為我們提供了大量方便實(shí)用的處理工具；MATLAB的這些特點(diǎn)，深受大家歡迎，由于個(gè)人電腦地普及，目前許多學(xué)校已將它做為本科生必須掌握的一種軟件。正是基于這些背景，我們編寫了這本信號(hào)與系統(tǒng)及MATLAB實(shí)現(xiàn)指導(dǎo)書，內(nèi)容包括信號(hào)的MATLAB表示、基本運(yùn)算、系統(tǒng)的時(shí)域分析、頻域分析、S域分析、狀態(tài)變量分析等。通過(guò)這些練習(xí)，同學(xué)們?cè)趯W(xué)習(xí)信號(hào)與系統(tǒng)的同時(shí)，掌握MATLAB的基本應(yīng)用，學(xué)會(huì)應(yīng)用MATLAB的數(shù)值計(jì)算

3、和符號(hào)計(jì)算功能，擺脫煩瑣的數(shù)學(xué)運(yùn)算，從而更注重于信號(hào)與系統(tǒng)的基本分析方法和應(yīng)用的理解與思考，將課程的重點(diǎn)、難點(diǎn)及部分習(xí)題用MATLAB進(jìn)行形象、直觀的可視化計(jì)算機(jī)模擬與仿真實(shí)現(xiàn)，加深對(duì)信號(hào)與系統(tǒng)的基本原理、方法及應(yīng)用的理解，為學(xué)習(xí)后續(xù)課程打好基礎(chǔ)。另外同學(xué)們?cè)谶M(jìn)行實(shí)驗(yàn)時(shí)，最好事先預(yù)習(xí)一些MATLAB的有關(guān)知識(shí)，以便更好地完成實(shí)驗(yàn)，同時(shí)實(shí)驗(yàn)中也可利用MATLAB的help命令了解具體語(yǔ)句以及指令的使用方法。實(shí)驗(yàn)一 基本信號(hào)在MATLAB中的表示和運(yùn)算一、實(shí)驗(yàn)?zāi)康? 學(xué)會(huì)用MATLAB表示常用連續(xù)信號(hào)的方法；2 學(xué)會(huì)用MATLAB進(jìn)行信號(hào)基本運(yùn)算的方法；二、實(shí)驗(yàn)原理1 連續(xù)信號(hào)的MATLAB表示M

4、ATLAB提供了大量的生成基本信號(hào)的函數(shù)，例如指數(shù)信號(hào)、正余弦信號(hào)。表示連續(xù)時(shí)間信號(hào)有兩種方法，一是數(shù)值法，二是符號(hào)法。數(shù)值法是定義某一時(shí)間范圍和取樣時(shí)間間隔，然后調(diào)用該函數(shù)計(jì)算這些點(diǎn)的函數(shù)值，得到兩組數(shù)值矢量，可用繪圖語(yǔ)句畫出其波形；符號(hào)法是利用MATLAB的符號(hào)運(yùn)算功能，需定義符號(hào)變量和符號(hào)函數(shù)，運(yùn)算結(jié)果是符號(hào)表達(dá)的解析式，也可用繪圖語(yǔ)句畫出其波形圖。例1-1指數(shù)信號(hào) 指數(shù)信號(hào)在MATLAB中用exp函數(shù)表示。如，調(diào)用格式為 ft=A*exp(a*t) 程序是 A=1; a=-0.4;t=0:0.01:10; %定義時(shí)間點(diǎn)ft=A*exp(a*t); %計(jì)算這些點(diǎn)的函數(shù)值plot(t,ft

5、); %畫圖命令，用直線段連接函數(shù)值表示曲線grid on; %在圖上畫方格例1-2 正弦信號(hào) 正弦信號(hào)在MATLAB中用 sin 函數(shù)表示。調(diào)用格式為 ft=A*sin(w*t+phi)A=1; w=2*pi; phi=pi/6;t=0:0.01:8; %定義時(shí)間點(diǎn)ft=A*sin(w*t+phi); %計(jì)算這些點(diǎn)的函數(shù)值plot(t,ft); %畫圖命令grid on; %在圖上畫方格例1-3 抽樣信號(hào) 抽樣信號(hào)Sa(t)=sin(t)/t在MATLAB中用 sinc 函數(shù)表示。定義為 t=-3*pi:pi/100:3*pi; ft=sinc(t/pi); plot(t,ft); grid

6、 on; axis(-10,10,-0.5,1.2); %定義畫圖范圍，橫軸，縱軸 title(抽樣信號(hào)) %定義圖的標(biāo)題名字例1-4 三角信號(hào) 三角信號(hào)在MATLAB中用 tripuls 函數(shù)表示。調(diào)用格式為 ft=tripuls(t,width,skew)，產(chǎn)生幅度為1，寬度為width，且以0為中心左右各展開width/2大小，斜度為skew的三角波。width的默認(rèn)值是1，skew的取值范圍是-1+1之間。一般最大幅度1出現(xiàn)在t=(width/2)*skew的橫坐標(biāo)位置。 t=-3:0.01:3; ft=tripuls(t,4,0.5); plot(t,ft); grid on; ax

7、is(-3,3,-0.5,1.5);例1-5 虛指數(shù)信號(hào) 調(diào)用格式是f=exp(j*w)*t) t=0:0.01:15;w=pi/4;X=exp(j*w*t);Xr=real(X); %取實(shí)部 Xi=imag(X); %取虛部Xa=abs(X); %取模Xn=angle(X); %取相位subplot(2,2,1),plot(t,Xr),axis(0,15,-(max(Xa)+0.5),max(Xa)+0.5),title(實(shí)部);subplot(2,2,3),plot(t,Xi),axis(0,15,-(max(Xa)+0.5),max(Xa)+0.5),title(虛部);subplot(

8、2,2,2), plot(t,Xa),axis(0,15,0,max(Xa)+1),title(模);subplot(2,2,4),plot(t,Xn),axis(0,15,-(max(Xn)+1),max(Xn)+1),title(相角); %subplot(m,n,i) 命令是建立m行n列畫圖窗口，并指定畫圖位置i例1-6 復(fù)指數(shù)信號(hào) 調(diào)用格式是f=exp(a+j*b)*t) t=0:0.01:3;a=-1;b=10;f=exp(a+j*b)*t); subplot(2,2,1),plot(t,real(f),title(實(shí)部)subplot(2,2,3),plot(t,imag(f),t

9、itle(虛部) subplot(2,2,2),plot(t,abs(f),title(模)subplot(2,2,4),plot(t,angle(f),title(相角)例1-7 矩形脈沖信號(hào) 矩形脈沖信號(hào)可用rectpuls函數(shù)產(chǎn)生，調(diào)用格式為y=rectpuls(t,width),幅度是1，寬度是width，以t=0為對(duì)稱中心。t=-2:0.01:2;width=1;ft=2*rectpuls(t,width);plot(t,ft)grid on;例1-8 單位階躍信號(hào) 單位階躍信號(hào)u(t)用“t=0”產(chǎn)生，調(diào)用格式為ft=(t=0) t=-1:0.01:5; ft=(t=0); plo

10、t(t,ft); grid on; axis(-1,5,-0.5,1.5);例1-9 正弦信號(hào)符號(hào)算法syms t %定義符號(hào)變量ty=sin(pi/4*t) %符號(hào)函數(shù)表達(dá)式ezplot(y,-16,16) %符號(hào)函數(shù)畫圖命令或者f=sym(sin(pi/4*t) %定義符號(hào)函數(shù)表達(dá)式ezplot(f,-16,16)例1-10 單位階躍信號(hào) MATTLAB符號(hào)數(shù)學(xué)函數(shù)Heaviside表示階躍信號(hào)，但要畫圖需在工作目錄創(chuàng)建Heaviside的M文件function f=Heaviside(t)f=(t0); 保存，文件名是Heaviside ，調(diào)用該函數(shù)即可畫圖，例t=-1:0.01:3;f

11、=heaviside(t);plot(t,f)axis(-1,3,-0.2,1.2) 或者y=sym(Heaviside(t);ezplot(y,-1,5);grid on 2 信號(hào)基本運(yùn)算的MATLAB實(shí)現(xiàn) 信號(hào)基本運(yùn)算是乘法、加法、尺度、反轉(zhuǎn)、平移、微分、積分，實(shí)現(xiàn)方法有數(shù)值法和符號(hào)法例1-11 以f(t)為三角信號(hào)為例，求f(2t) , f(2-2t) t=-3:0.001:3; ft=tripuls(t,4,0.5); subplot(3,1,1); plot(t,ft); grid on; title (f(t); ft1= tripuls(2*t,4,0.5); subplot(3

12、,1,2); plot(t,ft1); grid on; title (f(2t); ft2= tripuls(2-2*t,4,0.5); subplot(3,1,3); plot(t,ft2); grid on; title (f(2-2t);例1-12 已知f1(t)=sinwt , f2(t)=sin8wt , w=2pi , 求f1(t)+f2(t)和f1(t)f2(t) 的波形圖 w=2*pi;t=0:0.01:3;f1=sin(w*t);f2=sin(8*w*t);subplot(211)plot(t,f1+1,:,t,f1-1,:,t,f1+f2)grid on,title(f1

13、(t)+f2(t)subplot(212)plot(t,f1,:,t,-f1,:,t,f1.*f2)grid on,title(f1(t)*f2(t)符號(hào)算法也可實(shí)現(xiàn)上述運(yùn)算，以信號(hào)的微積分運(yùn)算為例說(shuō)明符號(hào)算法應(yīng)用微分的調(diào)用格式為 diff(function,variable,n)積分的調(diào)用格式為 int(function,variable,a,b)式中function表示要微分或積分的函數(shù)，variable表示運(yùn)算變量，n表示求導(dǎo)階數(shù)，默認(rèn)值是求一階導(dǎo)數(shù)，a是積分下限，b是積分上限，a b默認(rèn)是求不定積分。例1-13 求一階導(dǎo)數(shù)的例題，已知，clear syms a x y1 y2 %定義符

14、號(hào)變量a， x ，y1， y2 y1=sin(a*x2); %符號(hào)函數(shù)y1 y2=x*sin(x)*log(x); %符號(hào)函數(shù)y2 dy1=diff(y1,x) %無(wú)分號(hào)直接顯示結(jié)果 dy2=diff(y2) %無(wú)分號(hào)直接顯示結(jié)果例1-14 求積分的例題， clearsyms a x y3 y4y3=x5-a*x2+sqrt(x)/2;y4=(x*exp(x)/(1+x)2;iy3=int(y3,x)iy4=int(y4,0,1)三、上機(jī)實(shí)驗(yàn)內(nèi)容1 驗(yàn)證實(shí)驗(yàn)原理中程序2 畫出信號(hào)波形（1） （2）3信號(hào)，求、波形實(shí)驗(yàn)二 離散信號(hào)與系統(tǒng)的時(shí)域分析一、實(shí)驗(yàn)?zāi)康?學(xué)會(huì)用MATLAB表示常用離散信號(hào)的

15、方法；2學(xué)會(huì)用MATLAB實(shí)現(xiàn)離散信號(hào)卷積的方法；3 學(xué)會(huì)用MATLAB求解離散系統(tǒng)的單位響應(yīng)；4 學(xué)會(huì)用MATLAB求解離散系統(tǒng)的零狀態(tài)響應(yīng)；二、實(shí)驗(yàn)原理1離散信號(hào)的MATLAB表示表示離散時(shí)間信號(hào)f(k)需要兩個(gè)行向量，一個(gè)是表示序號(hào)k= ，一個(gè)是表示相應(yīng)函數(shù)值f= ，畫圖命令是stem。例2-1正弦序列信號(hào) 正弦序列信號(hào)可直接調(diào)用MATLAB函數(shù)cos，例，當(dāng)是整數(shù)或分?jǐn)?shù)時(shí)，才是周期信號(hào)。畫，波形程序是：k=0:40;subplot(2,1,1)stem(k,cos(k*pi/8),filled)title(cos(k*pi/8)subplot(2,1,2)stem(k,cos(2*k)

16、,filled)title(cos(2*k) 例2-2 單位序列信號(hào) 本題先建立一個(gè)畫單位序列的M函數(shù)文件，畫圖時(shí)調(diào)用。M文件建立方法：file / new / m-file 在文件編輯窗輸入程序，保存文件名用函數(shù)名。function dwxulie(k1,k2,k0) % k1 , k2 是畫圖時(shí)間范圍，k0是脈沖位置k=k1:k2;n=length(k);f=zeros(1,n);f(1,-k0-k1+1)=1; stem(k,f,filled) axis(k1,k2,0,1.5)title(單位序列(k) 保存文件名dwxulie.m畫圖時(shí)在命令窗口調(diào)用，例：dwxulie(-5,5,0

17、)例2-3 單位階躍序列信號(hào) 本題也可先建立一個(gè)畫單位階躍序列的M函數(shù)文件，畫圖時(shí)調(diào)用。 function jyxulie(k1,k2,k0)k=k1:-k0-1;kk=-k0:k2;n=length(k);nn=length(kk)u=zeros(1,n);uu=ones(1,nn); stem(kk,uu,filled)hold onstem(k,u,filled)hold offtitle(單位階躍序列)axis(k1 k2 0 1.5)保存文件名jyxulie.m畫圖時(shí)在命令窗口調(diào)用，例：jyxulie(-3,8,0)例2-4 實(shí)指數(shù)序列信號(hào) ，c、 a是實(shí)數(shù)。建立一個(gè)畫實(shí)指數(shù)序列的M

18、函數(shù)文件，畫圖時(shí)調(diào)用。function dszsu(c,a,k1,k2)%c：指數(shù)序列的幅度%a：指數(shù)序列的底數(shù)%k1：繪制序列的起始序號(hào)%k2：繪制序列的終止序號(hào)k=k1:k2;x=c*(a.k);stem(k,x,filled)hold onplot(k1,k2,0,0)hold off調(diào)用該函數(shù)畫信號(hào)： ，波形。dszsu(1,5/4,0,40)dszsu(1,-3/4,0,40)2 .離散信號(hào)的卷積和兩個(gè)有限長(zhǎng)序列f1，f2卷積可調(diào)用MATLAB函數(shù)conv，調(diào)用格式是f=conv(f1,f2), f是卷積結(jié)果，但不顯示時(shí)間序號(hào)，可自編一個(gè)函數(shù)dconv給出f和k，并畫圖。functi

19、on f,k=dconv(f1,f2,k1,k2) %The function of compute f=f1*f2% f: 卷積和序列f(k)對(duì)應(yīng)的非零樣值向量% k： 序列f(k)的對(duì)應(yīng)序號(hào)向量% f1: 序列f1(k)非零樣值向量% f2: 序列f2(k)的非零樣值向量% k1: 序列f1(k)的對(duì)應(yīng)序號(hào)向量% k2: 序列f2(k)的對(duì)應(yīng)序號(hào)向量f=conv(f1,f2) %計(jì)算序列f1與f2的卷積和fk0=k1(1)+k2(1);%計(jì)算序列f非零樣值的起點(diǎn)位置k3=length(f1)+length(f2)-2;%計(jì)算卷積和f的非零樣值的寬度k=k0:k0+k3 %確定卷積和f非零樣

20、值的序號(hào)向量subplot(2,2,1)stem(k1,f1)%在子圖1繪序列f1(k)時(shí)域波形圖title(f1(k)xlabel(k)ylabel(f1(k)subplot(2,2,2)stem(k2,f2)%在圖2繪序列f2(k)時(shí)波形圖title(f1(k)xlabel(k)ylabel(f2(k)subplot(2,2,3)stem(k,f);%在子圖3繪序列f(k)的波形圖title(f(k)f1(k)與f2(k)的卷積和f(k)xlabel(k)ylabel(f(k)h=get(gca,position);h(3)=2.5*h(3);set(gca,position,h)%將第三

21、個(gè)子圖的橫坐標(biāo)范圍擴(kuò)為原來(lái)的2.5倍例2-5求卷積和， f1=1 2 1;k1=-1 0 1;f2=ones(1,5);k2=-2:2;f, k=dconv(f1,f2,k1,k2) 由運(yùn)行結(jié)果知，f的長(zhǎng)度等于f1和f2長(zhǎng)度之和減一， f的起點(diǎn)是f1和f2的起點(diǎn)之和，f的終點(diǎn)是f1和f2的終點(diǎn)之和。3 離散系統(tǒng)的單位響應(yīng)MATLAB提供畫系統(tǒng)單位響應(yīng)函數(shù)impz，調(diào)用格式是impz(b,a) 式中b和a是表示離散系統(tǒng)的行向量；impz(b,a,n) 式中b和a是表示離散系統(tǒng)的行向量，時(shí)間范圍是0n；impz(b,a,n1,n2) 時(shí)間范圍是n1n2 ；y=impz(b,a,n1,n2) 由y給

22、出數(shù)值序列；例2-6已知 求單位響應(yīng)。a=1,-1,0.9;b=1;impz(b,a)impz(b,a,60)impz(b,a,-10:40)4 離散系統(tǒng)的零狀態(tài)響應(yīng)MATLAB提供求離散系統(tǒng)零狀態(tài)響應(yīng)數(shù)值解函數(shù)filter，調(diào)用格式為filter(b,a,x),式中b和a是表示離散系統(tǒng)的向量，x是輸入序列非零樣值點(diǎn)行向量，輸出向量序號(hào)同x一樣。例2-7 已知 ， 求零狀態(tài)響應(yīng)， 范圍020。a=1 -0.25 0.5;b=1 1;t=0:20;x=(1/2).t;y=filter(b,a,x)subplot(2,1,1)stem(t,x)title(輸入序列)subplot(2,1,2)st

23、em(t,y)title(響應(yīng)序列)三、上機(jī)實(shí)驗(yàn)內(nèi)容1驗(yàn)證實(shí)驗(yàn)原理中程序2已知，畫單位響應(yīng)波形。3已知，輸入，畫輸出波形，范圍015。實(shí)驗(yàn)三 連續(xù)時(shí)間LTI系統(tǒng)的時(shí)域分析一、實(shí)驗(yàn)?zāi)康?學(xué)會(huì)用MATLAB求解連續(xù)系統(tǒng)的零狀態(tài)響應(yīng)；2. 學(xué)會(huì)用MATLAB求解沖激響應(yīng)及階躍響應(yīng)； 3學(xué)會(huì)用MATLAB實(shí)現(xiàn)連續(xù)信號(hào)卷積的方法；二、實(shí)驗(yàn)原理1連續(xù)時(shí)間系統(tǒng)零狀態(tài)響應(yīng)的數(shù)值計(jì)算 我們知道，LTI連續(xù)系統(tǒng)可用如下所示的線性常系數(shù)微分方程來(lái)描述， 在MATLAB中，控制系統(tǒng)工具箱提供了一個(gè)用于求解零初始條件微分方程數(shù)值解的函數(shù)lsim。其調(diào)用格式y(tǒng)=lsim(sys,f,t)式中，t表示計(jì)算系統(tǒng)響應(yīng)的抽樣點(diǎn)向

24、量，f是系統(tǒng)輸入信號(hào)向量，sys是LTI系統(tǒng)模型，用來(lái)表示微分方程，差分方程或狀態(tài)方程。其調(diào)用格式sys=tf(b,a)式中，b和a分別是微分方程的右端和左端系數(shù)向量。例如，對(duì)于以下方程: 可用 獲得其LTI模型。注意，如果微分方程的左端或右端表達(dá)式中有缺項(xiàng)，則其向量a或b中的對(duì)應(yīng)元素應(yīng)為零，不能省略不寫，否則出錯(cuò)。例3-1 已知某LTI系統(tǒng)的微分方程為 y(t)+ 2y(t)+100y(t)=f(t) 其中，求系統(tǒng)的輸出y(t).解：顯然，這是一個(gè)求系統(tǒng)零狀態(tài)響應(yīng)的問(wèn)題。其MATLAB計(jì)算程序如下： ts=0;te=5;dt=0.01; sys=tf(1,1,2,100); t=ts:dt:

25、te; f=10*sin(2*pi*t); y=lsim(sys,f,t); plot(t,y); xlabel(Time(sec); ylabel(y(t);2連續(xù)時(shí)間系統(tǒng)沖激響應(yīng)和階躍響應(yīng)的求解 在MATLAB中，對(duì)于連續(xù)LTI系統(tǒng)的沖激響應(yīng)和階躍響應(yīng)，可分別用控制系統(tǒng)工具箱提供的函數(shù)impluse和step來(lái)求解。其調(diào)用格式為 y=impluse(sys,t) y=step(sys,t)式中，t表示計(jì)算系統(tǒng)響應(yīng)的抽樣點(diǎn)向量，sys是LTI系統(tǒng)模型。例3-2已知某LTI系統(tǒng)的微分方程為 y(t)+ 2y(t)+100y(t)=10f(t)求系統(tǒng)的沖激響應(yīng)和階躍響應(yīng)的波形.解：ts=0;te

26、=5;dt=0.01; sys=tf(10,1,2,100);t=ts:dt:te;h=impulse(sys,t);figure;plot(t,h); xlabel(Time(sec); ylabel(h(t);g=step(sys,t);figure;plot(t,g); xlabel(Time(sec);ylabel(g(t);3. 用MATLAB實(shí)現(xiàn)連續(xù)時(shí)間信號(hào)的卷積信號(hào)的卷積運(yùn)算有符號(hào)算法和數(shù)值算法，此處采用數(shù)值計(jì)算法，需調(diào)用MATLAB的conv( )函數(shù)近似計(jì)算信號(hào)的卷積積分。連續(xù)信號(hào)的卷積積分定義是 如果對(duì)連續(xù)信號(hào)和進(jìn)行等時(shí)間間隔均勻抽樣，則和分別變?yōu)殡x散時(shí)間信號(hào)和。其中，為整

27、數(shù)。當(dāng)足夠小時(shí)，和既為連續(xù)時(shí)間信號(hào)和。因此連續(xù)時(shí)間信號(hào)卷積積分可表示為 采用數(shù)值計(jì)算時(shí)，只求當(dāng)時(shí)卷積積分的值，其中，n為整數(shù)，既 其中，實(shí)際就是離散序列和的卷積和。當(dāng)足夠小時(shí)，序列就是連續(xù)信號(hào)的數(shù)值近似，既 上式表明，連續(xù)信號(hào)和的卷積，可用各自抽樣后的離散時(shí)間序列的卷積再乘以抽樣間隔。抽樣間隔越小，誤差越小。例3-3用數(shù)值計(jì)算法求與的卷積積分。解：因?yàn)槭且粋€(gè)持續(xù)時(shí)間無(wú)限長(zhǎng)的信號(hào)，而計(jì)算機(jī)數(shù)值計(jì)算不可能計(jì)算真正的無(wú)限長(zhǎng)信號(hào)，所以在進(jìn)行的抽樣離散化時(shí)，所取的時(shí)間范圍讓衰減到足夠小就可以了，本例取。程序是dt=0.01; t=-1:dt:2.5;f1=Heaviside(t)-Heaviside(t

28、-2);f2=exp(-3*t).*Heaviside(t);f=conv(f1,f2)*dt; n=length(f); tt=(0:n-1)*dt-2;subplot(221), plot(t,f1), grid on;axis(-1,2.5,-0.2,1.2); title(f1(t); xlabel(t)subplot(222), plot(t,f2), grid on;axis(-1,2.5,-0.2,1.2); title(f2(t); xlabel(t)subplot(212), plot(tt,f), grid on;title(f(t)=f1(t)*f2(t); xlabel

29、(t)由于和的時(shí)間范圍都是從t=-1開始，所以卷積結(jié)果的時(shí)間范圍從 t=-2開始，增量還是取樣間隔，這就是語(yǔ)句tt=(0:n-1)*dt-2的由來(lái)。三、上機(jī)實(shí)驗(yàn)內(nèi)容1 驗(yàn)證實(shí)驗(yàn)原理中所述的相關(guān)程序2 已知描述系統(tǒng)的微分方程和激勵(lì)信號(hào)f(t)如下，試用解析法求系統(tǒng)的零狀態(tài)響應(yīng)y(t)，并用MATLAB繪出系統(tǒng)零狀態(tài)響應(yīng)的時(shí)域仿真波形，驗(yàn)證結(jié)果是否相同y(t)+ 4y(t)+4y(t)=f(t)+3f(t) f(t)= exp(-t)u(t) 3已知描述系統(tǒng)的微分方程如下，試用MATLAB求系統(tǒng)在010秒范圍內(nèi)沖激響應(yīng)和階躍響應(yīng)的數(shù)值解，并用繪出系統(tǒng)沖激響應(yīng)和階躍響應(yīng)的時(shí)域波形y(t)+3y(t)

30、+2y(t)=f(t)y(t)+ 2y(t)+2y(t)=f(t) 4畫出信號(hào)卷積積分的波形，實(shí)驗(yàn)四 傅里葉變換、系統(tǒng)的頻域分析一、 實(shí)驗(yàn)?zāi)康?、學(xué)會(huì)用MATLAB實(shí)現(xiàn)連續(xù)時(shí)間信號(hào)傅里葉變換2、學(xué)會(huì)用MATLAB分析LTI系統(tǒng)的頻域特性3、學(xué)會(huì)用MATLAB分析LTI系統(tǒng)的輸出響應(yīng)二、實(shí)驗(yàn)原理1傅里葉變換的MATLAB求解 MATLAB的symbolic Math Toolbox 提供了直接求解傅里葉變換及逆變換的函數(shù)fourier()及ifourier()兩者的調(diào)用格式如下。Fourier 變換的調(diào)用格式F=fourier(f)：它是符號(hào)函數(shù)f的fourier變換默認(rèn)返回是關(guān)于w的函數(shù)。F=

31、fourier(f，v):它返回函數(shù)F是關(guān)于符號(hào)對(duì)象v的函數(shù)，而不是默認(rèn)的w，即 Fourier逆變換的調(diào)用格式f=ifourier(F):它是符號(hào)函數(shù)F的fourier逆變換，默認(rèn)的獨(dú)立變量為w，默認(rèn)返回是關(guān)于x的函數(shù)。f=ifourier(f,u):它的返回函數(shù)f是u的函數(shù)，而不是默認(rèn)的x.注意：在調(diào)用函數(shù)fourier()及ifourier()之前，要用syms命令對(duì)所用到的變量（如t,u,v,w）進(jìn)行說(shuō)明,即將這些變量說(shuō)明成符號(hào)變量。例4-1 求的傅立葉變換解: 可用MATLAB解決上述問(wèn)題：syms tFw=fourier(exp(-2*abs(t)例4-2 求的逆變換f(t)解：

32、可用MATLAB解決上述問(wèn)題syms t wft=ifourier(1/(1+w2),t)2連續(xù)時(shí)間信號(hào)的頻譜圖例4-3 求調(diào)制信號(hào)的頻譜，式中解：MATLAB程序如下所示ft=sym(4*cos(2*pi*6*t)*(Heaviside(t+1/4)-Heaviside(t-1/4);Fw=simplify(fourier(ft)subplot(121)ezplot(ft,-0.5 0.5),grid onsubplot(122)ezplot(abs(Fw),-24*pi 24*pi),grid 用MATLAB符號(hào)算法求傅里葉變換有一定局限，當(dāng)信號(hào)不能用解析式表達(dá)時(shí)，會(huì)提示出錯(cuò)，這時(shí)用MAT

33、LAB的數(shù)值計(jì)算也可以求連續(xù)信號(hào)的傅里葉變換，計(jì)算原理是當(dāng)足夠小時(shí)，近似計(jì)算可滿足要求。若信號(hào)是時(shí)限的，或當(dāng)時(shí)間大于某個(gè)給定值時(shí)，信號(hào)已衰減的很厲害，可以近似地看成時(shí)限信號(hào)時(shí)，n的取值就是有限的，設(shè)為N，有 是頻率取樣點(diǎn)時(shí)間信號(hào)取樣間隔應(yīng)小于奈奎斯特取樣時(shí)間間隔，若不是帶限信號(hào)可根據(jù)計(jì)算精度要求確定一個(gè)頻率 W0為信號(hào)的帶寬。例4-4 用數(shù)值計(jì)算法求信號(hào)的傅里葉變換解，信號(hào)頻譜是，第一個(gè)過(guò)零點(diǎn)是，一般將此頻率視為信號(hào)的帶寬，若將精度提高到該值的50倍，既W0=50，據(jù)此確定取樣間隔，R=0.02;t=-2:R:2;f=Heaviside(t+1)-Heaviside(t-1);W1=2*pi*

34、5;N=500;k=0:N;W=k*W1/N;F=f*exp(-j*t*W)*R;F=real(F);W=-fliplr(W),W(2:501);F=fliplr(F),F(2:501);subplot(2,1,1);plot(t,f);xlabel(t);ylabel(f(t);title(f(t)=u(t+1)-u(t-1);subplot(2,1,2);plot(W,F);xlabel(w);ylabel(F(w);title(f(t)的付氏變換F(w);3用MATLAB分析LTI系統(tǒng)的頻率特性 當(dāng)系統(tǒng)的頻率響應(yīng)H（jw）是jw的有理多項(xiàng)式時(shí)，有 MATLAB信號(hào)處理工具箱提供的freq

35、s函數(shù)可直接計(jì)算系統(tǒng)的頻率響應(yīng)的數(shù)值解。其調(diào)用格式如下H=freqs(b,a,w)其中，a和b分別是H(jw)的分母和分子多項(xiàng)式的系數(shù)向量，w為形如w1:p:w2的向量，定義系統(tǒng)頻率響應(yīng)的頻率范圍，w1為頻率起始值，w2為頻率終止值,p為頻率取樣間隔。H返回w所定義的頻率點(diǎn)上，系統(tǒng)頻率響應(yīng)的樣值。例如，運(yùn)行如下命令，計(jì)算02pi頻率范圍內(nèi)以間隔0.5取樣的系統(tǒng)頻率響應(yīng)的樣值a=1 2 1;b=0 1;h=freqs(b,a,0:0.5:2*pi)例 4-5 三階歸一化的butterworth 低通濾波器的頻率響應(yīng)為 試畫出該系統(tǒng)的幅度響應(yīng)和相位響應(yīng)。解 其MATLAB程序及響應(yīng)的波形如下w=0

36、:0.025:5;b=1;a=1,2,2,1;H=freqs(b,a,w);subplot(2,1,1);plot(w,abs(H);grid;xlabel(omega(rad/s);ylabel(|H(jomega)|);title(H(jw)的幅頻特性);subplot(2,1,2);plot(w,angle (H);grid;xlabel(omega(rad/s);ylabel(phi(omega);title(H(jw)的相頻特性);4用MATLAB分析LTI系統(tǒng)的輸出響應(yīng)例 4-6已知一RC電路如圖所示 系統(tǒng)的輸入電壓為f(t)，輸出信號(hào)為電阻兩端的電壓y(t).當(dāng)RC=0.04，f

37、(t)=cos5t+cos100t, 試求該系統(tǒng)的響應(yīng)y(t)-+-+f(t)y(t)RC解 由圖可知 ，該電路為一個(gè)微分電路，其頻率響應(yīng)為 由此可求出余弦信號(hào)通過(guò)LTI系統(tǒng)的響應(yīng)為 計(jì)算該系統(tǒng)響應(yīng)的MATLAB程序及響應(yīng)波形如下RC=0.04;t=linspace(-2,2,1024);w1=5;w2=100;H1=j*w1/(j*w1+1/RC);H2=j*w2/(j*w2+1/RC);f=cos(5*t)+cos(100*t);y=abs(H1)*cos(w1*t+angle(H1)+ abs(H2)*cos(w2*t+angle(H2);subplot(2,1,1);plot(t,f)

38、;ylabel(f(t);xlabel(Time(s);subplot(2,1,2);plot(t,y);ylabel(y(t);xlabel(Time(s);三、 上機(jī)實(shí)驗(yàn)內(nèi)容1.驗(yàn)證實(shí)驗(yàn)原理中所述的相關(guān)程序；2.試用MATLAB求單邊指數(shù)數(shù)信號(hào)的傅立葉變換，并畫出其波形；3.設(shè)，試用MATLAB畫出該系統(tǒng)的幅頻特性和相頻特性，并分析系統(tǒng)具有什么濾波特性。實(shí)驗(yàn)五 信號(hào)抽樣與恢復(fù)一、實(shí)驗(yàn)?zāi)康膶W(xué)會(huì)用MATLAB實(shí)現(xiàn)連續(xù)信號(hào)的采樣和重建二、實(shí)驗(yàn)原理 1抽樣定理若是帶限信號(hào)，帶寬為, 經(jīng)采樣后的頻譜就是將的頻譜 在頻率軸上以采樣頻率為間隔進(jìn)行周期延拓。因此，當(dāng)時(shí)，不會(huì)發(fā)生頻率混疊；而當(dāng) 時(shí)將發(fā)生頻率

39、混疊。2信號(hào)重建經(jīng)采樣后得到信號(hào)經(jīng)理想低通則可得到重建信號(hào)，即：=*其中：=所以：=*=* =上式表明，連續(xù)信號(hào)可以展開成抽樣函數(shù)的無(wú)窮級(jí)數(shù)。利用MATLAB中的來(lái)表示，有 ，所以可以得到在MATLAB中信號(hào)由重建的表達(dá)式如下：=我們選取信號(hào)=作為被采樣信號(hào)，當(dāng)采樣頻率=2時(shí)，稱為臨界采樣。我們?nèi)±硐氲屯ǖ慕刂诡l率=。下面程序?qū)崿F(xiàn)對(duì)信號(hào)=的采樣及由該采樣信號(hào)恢復(fù)重建：例5-1 Sa(t)的臨界采樣及信號(hào)重構(gòu)；wm=1; %信號(hào)帶寬wc=wm; %濾波器截止頻率Ts=pi/wm; %采樣間隔ws=2*pi/Ts; %采樣角頻率 n=-100:100; %時(shí)域采樣電數(shù)nTs=n*Ts %時(shí)域采樣點(diǎn)

40、f=sinc(nTs/pi);Dt=0.005;t=-15:Dt:15;fa=f*Ts*wc/pi*sinc(wc/pi)*(ones(length(nTs),1)*t-nTs*ones(1,length(t); %信號(hào)重構(gòu)t1=-15:0.5:15;f1=sinc(t1/pi);subplot(211);stem(t1,f1);xlabel(kTs);ylabel(f(kTs);title(sa(t)=sinc(t/pi)的臨界采樣信號(hào));subplot(212);plot(t,fa)xlabel(t);ylabel(fa(t);title(由sa(t)=sinc(t/pi)的臨界采樣信號(hào)重

41、構(gòu)sa(t);grid;例5-2 Sa(t)的過(guò)采樣及信號(hào)重構(gòu)和絕對(duì)誤差分析程序和例4-1類似，將采樣間隔改成Ts=0.7*pi/wm , 濾波器截止頻率該成wc=1.1*wm ，添加一個(gè)誤差函數(shù)wm=1;wc=1.1*wm; Ts=0.7*pi/wm;ws=2*pi/Ts;n=-100:100;nTs=n*Tsf=sinc(nTs/pi);Dt=0.005;t=-15:Dt:15;fa=f*Ts*wc/pi*sinc(wc/pi)*(ones(length(nTs),1)*t-nTs*ones(1,length(t);error=abs(fa-sinc(t/pi); %重構(gòu)信號(hào)與原信號(hào)誤差t

42、1=-15:0.5:15;f1=sinc(t1/pi);subplot(311);stem(t1,f1);xlabel(kTs);ylabel(f(kTs);title(sa(t)=sinc(t/pi)的采樣信號(hào));subplot(312);plot(t,fa)xlabel(t);ylabel(fa(t);title(由sa(t)=sinc(t/pi)的過(guò)采樣信號(hào)重構(gòu)sa(t);grid;subplot(313);plot(t,error);xlabel(t);ylabel(error(t);title(過(guò)采樣信號(hào)與原信號(hào)的誤差error(t);例5-3 Sa(t)的欠采樣及信號(hào)重構(gòu)和絕對(duì)誤差

43、分析程序和例4-2類似，將采樣間隔改成Ts=1.5*pi/wm , 濾波器截止頻率該成wc=wm=1三、上機(jī)實(shí)驗(yàn)內(nèi)容1驗(yàn)證實(shí)驗(yàn)原理中所述的相關(guān)程序；2設(shè)f(t)=0.5*(1+cost)*(u(t+pi)-u(t-pi) ，由于不是嚴(yán)格的頻帶有限信號(hào)，但其頻譜大部分集中在0，2之間，帶寬wm可根據(jù)一定的精度要求做一些近似。試根據(jù)以下兩種情況用 MATLAB實(shí)現(xiàn)由f(t)的抽樣信號(hào)fs(t)重建f(t) 并求兩者誤差，分析兩種情況下的結(jié)果。(1) wm=2 , wc=1.2wm , Ts=1;(2) wm=2 , wc=2 , Ts=2.5 實(shí)驗(yàn)六 信號(hào)與系統(tǒng)復(fù)頻域分析一、實(shí)驗(yàn)?zāi)康?.學(xué)會(huì)用MA

44、TLAB進(jìn)行部分分式展開；2.學(xué)會(huì)用MATLAB分析LTI系統(tǒng)的特性；3.學(xué)會(huì)用MATLAB進(jìn)行Laplace正、反變換。4.學(xué)會(huì)用MATLAB畫離散系統(tǒng)零極點(diǎn)圖；5.學(xué)會(huì)用MATLAB分析離散系統(tǒng)的頻率特性；二、實(shí)驗(yàn)原理及內(nèi)容1用MATLAB進(jìn)行部分分式展開用MATLAB函數(shù)residue可以得到復(fù)雜有理分式F(s)的部分分式展開式，其調(diào)用格式為 其中，num,den分別為F(s)的分子和分母多項(xiàng)式的系數(shù)向量，r為部分分式的系數(shù)，p為極點(diǎn)，k為F(s)中整式部分的系數(shù)，若F(s)為有理真分式，則k為零。例6-1 用部分分式展開法求F(s)的反變換 解:其MATLAB程序?yàn)閒ormat rat

45、;num=1,2;den=1,4,3,0;r,p=residue(num,den)程序中format rat是將結(jié)果數(shù)據(jù)以分?jǐn)?shù)形式顯示F(s)可展開為 所以，F(xiàn)(s)的反變換為 2用MATLAB分析LTI系統(tǒng)的特性系統(tǒng)函數(shù)H（s）通常是一個(gè)有理分式，其分子和分母均為多項(xiàng)式。計(jì)算H（s）的零極點(diǎn)可以應(yīng)用MATLAB中的roots函數(shù)，求出分子和分母多項(xiàng)式的根，然后用plot命令畫圖。在MATLAB中還有一種更簡(jiǎn)便的方法畫系統(tǒng)函數(shù)H（s）的零極點(diǎn)分布圖，即用pzmap函數(shù)畫圖。其調(diào)用格式為 pzmap(sys)sys表示LTI系統(tǒng)的模型，要借助tf函數(shù)獲得，其調(diào)用格式為sys=tf(b,a)式中，

46、b和a分別為系統(tǒng)函數(shù)H（s）的分子和分母多項(xiàng)式的系數(shù)向量。如果已知系統(tǒng)函數(shù)H（s），求系統(tǒng)的單位沖激響應(yīng)h(t)和頻率響應(yīng)可以用以前介紹過(guò)的impulse和freqs函數(shù)。例6-2 已知系統(tǒng)函數(shù)為 試畫出其零極點(diǎn)分布圖，求系統(tǒng)的單位沖激響應(yīng)h(t)和頻率響應(yīng)，并判斷系統(tǒng)是否穩(wěn)定。解：其MATLAB程序如下：num=1;den=1,2,2,1;sys=tf(num,den);figure(1);pzmap(sys);t=0:0.02:10;h=impulse(num,den,t);figure(2);plot(t,h)title(Impulse Response)H,w=freqs(num,de

47、n);figure(3);plot(w,abs(H)xlabel(omega)title(Magnitude Response) 3用MATLAB進(jìn)行Laplace正、反變換MATLAB的符號(hào)數(shù)學(xué)工具箱提供了計(jì)算Laplace正、反變換的函數(shù)Laplace和ilaplace,其調(diào)用格式為 上述兩式右端的f和F分別為時(shí)域表示式和s域表示式的符號(hào)表示，可以應(yīng)用函數(shù)sym實(shí)現(xiàn)，其調(diào)用格式為 S=sym(A)式中，A為待分析表示式的字符串，S為符號(hào)數(shù)字或變量。例6-3 試分別用Laplace和ilaplace函數(shù)求（1）的Laplace變換；（2）的Laplace反變換。解：（1）其程序?yàn)閒=sym(

48、exp(-t)*sin(a*t);F=laplace(f)或syms a tF=laplace(exp(-t)*sin(a*t)（2）其程序?yàn)镕=sym(s2/(s2+1);ft=ilaplace(F)或syms sft= ilaplace(s2/(s2+1)4離散系統(tǒng)零極點(diǎn)圖離散系統(tǒng)可以用下述差分方程描述：Z變換后可得系統(tǒng)函數(shù)：用MATLAB提供的root函數(shù)可分別求零點(diǎn)和極點(diǎn)，調(diào)用格式是p=a0,a1an,q=b0,b1bm,0,00, 補(bǔ)0使二者維數(shù)一樣。畫零極點(diǎn)圖的方法有多種，可以用MATLAB函數(shù)z,p,k=tf2zp(b,a)和zplane(q,p)，也可用plot命令自編一函數(shù)ljdt.m,畫圖時(shí)調(diào)用。function ljdt(A,B)% The function to draw the pole-zero diagram for discrete systemp=roots(A)； %求系統(tǒng)極點(diǎn)q=roots(B)； %求系統(tǒng)零點(diǎn)p=p； %將極點(diǎn)列向量轉(zhuǎn)置為行向量q=q； %將零點(diǎn)列向量轉(zhuǎn)置為行向量x=max(abs(p q 1)；%確定縱坐標(biāo)范圍x=x+0.1；y=x；%確定橫坐標(biāo)范圍clfhold onaxis(-x x -y y)%確定坐標(biāo)軸顯示范圍w=0:pi/300:2*pi；