信号与系统及MATLAB实现.doc

信号与系统及MATAB实现实验指导书前 言长期以来，信号与系统课程一直采用单一理论教学方式，同学们依靠做习题来巩固与理解教学内容,虽然手工演算训练了计算能力与思维方法，但就是由于本课程数学公式推导较多,概念抽象，常需画各种波形,作题时难免花费很多时间，现在，我们给同学们介绍一种国际上公认得优秀科技应用软件MTLAB,借助它我们可以在电脑上轻松地完成许多习题得演算与波形得绘制。 MALA得功能非常强大,我们此处仅用到它得一部分,在后续课程中我们还会用到它,在未来地科学研究与工程设计中有可能继续用它,所以有兴趣得同学，可以对MTLAB再多了解一些。MLAB究竟有那些特点呢?1.高效得数值计算与符号计算功能,使我们从繁杂得数学运算分析中解脱出来;2。完备得图形处理功能，实现计算结果与编程得可视化;。友好得用户界面及接近数学表达式得自然化语言,易于学习与掌握；4.功能丰富得应用工具箱，为我们提供了大量方便实用得处理工具；MATAB得这些特点，深受大家欢迎,由于个人电脑地普及，目前许多学校已将它做为本科生必须掌握得一种软件。正就是基于这些背景，我们编写了这本信号与系统及MATAB实现指导书，内容包括信号得MALB表示、基本运算、系统得时域分析、频域分析、域分析、状态变量分析等。通过这些练习，同学们在学习信号与系统得同时,掌握ATAB得基本应用,学会应用AB得数值计算与符号计算功能,摆脱烦琐得数学运算,从而更注重于信号与系统得基本分析方法与应用得理解与思考,将课程得重点、难点及部分习题用MATA进行形象、直观得可视化计算机模拟与仿真实现,加深对信号与系统得基本原理、方法及应用得理解,为学习后续课程打好基础。另外同学们在进行实验时,最好事先预习一些ATAB得有关知识,以便更好地完成实验，同时实验中也可利用MTLAB得he命令了解具体语句以及指令得使用方法。实验一 基本信号在ATLA中得表示与运算一、实验目得1 学会用MATLAB表示常用连续信号得方法；2 学会用MALA进行信号基本运算得方法;二、实验原理1 连续信号得MALA表示MTAB提供了大量得生成基本信号得函数,例如指数信号、正余弦信号。表示连续时间信号有两种方法,一就是数值法，二就是符号法.数值法就是定义某一时间范围与取样时间间隔，然后调用该函数计算这些点得函数值,得到两组数值矢量,可用绘图语句画出其波形;符号法就是利用ATLB得符号运算功能,需定义符号变量与符号函数，运算结果就是符号表达得解析式，也可用绘图语句画出其波形图。例1指数信号指数信号在ATLAB中用exp函数表示。如，调用格式为ftAexp(a*t） 程序就是 A=1; a=-0、;0：0、01：10; 定义时间点ftA*exp（*）; %计算这些点得函数值plo（t,t)； %画图命令，用直线段连接函数值表示曲线grid on； %在图上画方格例1 正弦信号 正弦信号在MATB中用 sin 函数表示.调用格式为 ft=A*(w*t+phi)A1； w=2*pi； p=pi/6；=0:、0：; %定义时间点t=*sin（wt+phi)； 计算这些点得函数值pt(t，t)； %画图命令d on; 在图上画方格例1-3 抽样信号 抽样信号a()=si(t）/t在MATLB中用 inc 函数表示。定义为 t=-3*pi:i/10:3i; ft=si(t/i); plot(t,ft）; grd on; xis(10，0,-0、,1、2); 定义画图范围,横轴,纵轴 itl(抽样信号） %定义图得标题名字例-4 三角信号 三角信号在ALB中用 trpus函数表示。调用格式为 t=is(t，ih,skw）,产生幅度为1,宽度为wdth,且以为中心左右各展开with/2大小，斜度为skew得三角波.width得默认值就是1,skw得取值范围就是1+之间。一般最大幅度1出现在t=（width/2）*sew得横坐标位置. -3:0、1:3; ft=truls（t,、5); plot(,t)； gr on; is(3,3,-、5,、5）;例1-5 虚指数信号 调用格式就是f=xp（(j*w)*t) t=：、01:15;pi/4;=exp(j*w*t)；Xrrl（X）； %取实部 i=ma(); 取虚部Xa=abs(X）; %取模n=ange(X)； %取相位subpt(2，,1)，plot（t,Xr)，x(0,-(m(Xa)0、5),max(Xa）+0、5）,tite(实部);ubplot(2,2,）,plot(,Xi)，axis（,15,-（max(Xa)0、5),ax（Xa)+0、5），t(虚部)；sublo(2，2，2)， lot（t，Xa),axis(,，0，mx()+1),ttle(模)；bplt（2，2,4),plot(t，Xn）,xs（0,，-（max(Xn）+）,m(X）+1），tile（相角）; %suplo(m,n，) 命令就是建立m行n列画图窗口,并指定画图位置i例1-6 复指数信号调用格式就是ep（(+*b)*） =:、1:3;=;=10；=exp(a+j*b）t)； spot（2,2,1)，t(t,al(f)）,ttle(实部)subplot(,2,3),ot(，imag（f),title(虚部) sblot(2，2,2)，pot(,as(f)，title(模)splo(2,2，)，lot（t,angle()），ttle(相角）例17 矩形脉冲信号 矩形脉冲信号可用rctpls函数产生,调用格式为y=tpu（,dth),幅度就是1,宽度就是width，以t0为对称中心。t2：、0:2；wdh=1;ft=2*els（t,wdth);ot(t，)gd n;例1- 单位阶跃信号 单位阶跃信号u（t)用“t=0”产生,调用格式为f=（0) t1:0、01:5; t=（=0）; pot(t,ft）; gid on； axs(-,5,0、5，1、5）;例1-9 正弦信号符号算法sym 定义符号变量tyn（pi/4*t) %符号函数表达式ezlo(y,-1，1） %符号函数画图命令或者=sym(sin(p/*t)） %定义符号函数表达式zlo（f，-1,)例1-10 单位阶跃信号 MATLAB符号数学函数Heaisie表示阶跃信号，但要画图需在工作目录创建Heaiside得M文件uon fHeaise(t）f(t0）； 保存,文件名就是Heavisde ,调用该函数即可画图,例=：0、0：3;f=hiside（t）；plot(t，f)axis(1,3,0、2,1、2) 或者=(Havside（t）)；ezplo（y，-1,5)；grid on 2 信号基本运算得MATLB实现 信号基本运算就是乘法、加法、尺度、反转、平移、微分、积分,实现方法有数值法与符号法例111 以f(t）为三角信号为例,求f（）， f(-2t) t3:、00:3； t=trpuls（t,4,0、5); subpo(3,1,1)； plt（t，f）; id o； itle (f()； ft= tripls（*t,0、5)； sblt(3，1，）; pot（,ft1）;grid on； ttl （f(); ft2=ripul(22*,4，0、5); uplot(3,1,3); (t,ft); gridon； title(（2t);例11已知f1(t)=int ， f2(t）=snwt, =i ,求f1(t)+f(t)与f()2(t） 得波形图 w*pi;=0:、01:;fsin(w*t）;f2=in(*wt）；subplot(211)plot(t,f+1,:,t,f1，：，t，f1+f）grid on,title(f()f(t）)）uplot(22)pot(t，f，:，t,f1,:,t,f1、f2)gridon,tie(f1(t）*f2(t）符号算法也可实现上述运算，以信号得微积分运算为例说明符号算法应用微分得调用格式为iff（fuction,vaiable,n)积分得调用格式为 int(fncti,vrial,，b)式中ncto表示要微分或积分得函数，riab表示运算变量,n表示求导阶数，默认值就是求一阶导数,a就是积分下限，b就是积分上限,a b默认就是求不定积分。例113 求一阶导数得例题，已知,clear symsa y1 y2 定义符号变量， x ，y1, 2 y1=sin(a*2)； 符号函数y1 2=*sin(x)*l(x)；%符号函数y2 d1=diff（y,x） %无分号直接显示结果 dy2if(2) 无分号直接显示结果例1-14求积分得例题， clearsys x y34y3=x5ax2+srt（x)/2；y4=(xexp(x)/(1+）2;iy3=i(y3,x）iy4int(y4,，1）三、上机实验内容1 验证实验原理中程序2 画出信号波形(1) ()3。信号,求、波形实验二 离散信号与系统得时域分析一、实验目得.学会用MATLAB表示常用离散信号得方法;2学会用MATLAB实现离散信号卷积得方法;3 学会用MALA求解离散系统得单位响应;4 学会用MATLAB求解离散系统得零状态响应;二、实验原理1.离散信号得ATAB表示表示离散时间信号f（k)需要两个行向量，一个就是表示序号k= ,一个就是表示相应函数值= ,画图命令就是st。例2正弦序列信号 正弦序列信号可直接调用MATA函数cos,例,当就是整数或分数时,才就是周期信号画,波形程序就是：k0:4;subplot(2,1，)se(k，cos(k*pi/8),led）tite（k*pi/8)）ubplt(，1,2)stm(,cs(2*k），filled)title（co（2*)例2-2 单位序列信号 本题先建立一个画单位序列得M函数文件,画图时调用.文件建立方法：fie/ nw/ mfil 在文件编辑窗输入程序,保存文件名用函数名.fnctndwxlie(k1,2，k0) 1 ， 就是画图时间范围，k0就是脉冲位置k=k1：k2;n=lenth(k);=zos（，n);f（1,-k0k1+)1; t(,f，fill)axis（1,k，0,1、)tile(单位序列（k）) 保存文件名dwxule、m画图时在命令窗口调用，例:dwxie(5,5,0）例3 单位阶跃序列信号 本题也可先建立一个画单位阶跃序列得M函数文件,画图时调用。 unctin jyxuie(,k，k0）k=k1：k01；kkk0:k2;n=lenth(k）；=le(k)u=eros（1,n）；uu=ones(1,n); sem(,uu,fled)hol onstem（k,u，fild）hold ffite(单位阶跃序列)is(1 k2 1、5）保存文件名jxuie、m画图时在命令窗口调用,例：jyxulie(,8，0)例24 实指数序列信号 ，、 a就是实数。建立一个画实指数序列得M函数文件,画图时调用.fnction dsu(,a,k1，k2)：指数序列得幅度%a：指数序列得底数k：绘制序列得起始序号：绘制序列得终止序号k=:2;x=c*(a、k）;tem（k,x，fied)hold onlot(k1,k,0)holoff调用该函数画信号:，波形.dszu(1，5/4,0，40）dssu(，-4,0,)2 、离散信号得卷积与两个有限长序列，2卷积可调用MATL函数ov，调用格式就是f=conv(，2), f就是卷积结果,但不显示时间序号，可自编一个函数dcon给出f与k，并画图。functn f，k=dconv(f1，k1，k2) %he ftion tef=f% f： 卷积与序列f(k)对应得非零样值向量 k： 序列f（k)得对应序号向量 f1: 序列1(k）非零样值向量% f2：序列(）得非零样值向量k: 序列f（k）得对应序号向量% ： 序列2（k)得对应序号向量f=conv（1,2) %计算序列f1与2得卷积与fk0=1(）+2(1）;计算序列f非零样值得起点位置3=lengt(1)egh(）2;计算卷积与f得非零样值得宽度：k03 %确定卷积与f非零样值得序号向量sublot（2，2,1）sem（k1，f1）%在子图绘序列f1(）时域波形图titl（1(k)）xlae（k)ylae（f(）)subpo（,，2)tem(2,2）%在图2绘序列f2（k）时波形图titl（f1()xlab(k）ylabl(f（k)sublo(2,3）tem(k，f)；%在子图绘序列f(k)得波形图title（f（）f(k)与2(）得卷积与(k)）xabel(）yabel（f（k)h=et(gc,psio);(3)=2、5*（3);et(gca,ostion,)将第三个子图得横坐标范围扩为原来得2、5倍例2-5求卷积与, f=2 1;k=1 0 1;f2=ones(1,5);k2=2：2;f， dov(f1，f2,k2) 由运行结果知,f得长度等于1与f长度之与减一， f得起点就是1与f2得起点之与,得终点就是f1与f2得终点之与。3 离散系统得单位响应ALAB提供画系统单位响应函数mpz,调用格式就是m(b,）式中与a就是表示离散系统得行向量；iz(b,a,n)式中b与a就是表示离散系统得行向量,时间范围就是0n；iz(b,a,n1,n)时间范围就是nn2 ;y=impz(b,a,n1,2)由y给出数值序列；例26已知 求单位响应。a=,，、9;b=1;ipz（，)impz(b,，6)ipz(，,10:4)4 离散系统得零状态响应MTLA提供求离散系统零状态响应数值解函数flter，调用格式为ilte（b,a,x)，式中b与a就是表示离散系统得向量,x就是输入序列非零样值点行向量,输出向量序号同x一样.例2 已知， 求零状态响应,范围00a1 -0、2 0、5;b=1;0：0；x=(1/2)、t;y=r(,，x)subplo（2,1,）stem(，）ite(输入序列)sbplot(2,1,)em（t,y）titl(响应序列）三、上机实验内容1验证实验原理中程序2。已知,画单位响应波形.3.已知，输入，画输出波形，范围15。实验三 连续时间LI系统得时域分析一、实验目得学会用MTLAB求解连续系统得零状态响应;、 学会用MATLAB求解冲激响应及阶跃响应； .学会用M实现连续信号卷积得方法;二、实验原理1连续时间系统零状态响应得数值计算 我们知道,T连续系统可用如下所示得线性常系数微分方程来描述, 在TLAB中，控制系统工具箱提供了一个用于求解零初始条件微分方程数值解得函数lsi。其调用格式ylsim(sys，f,t）式中，t表示计算系统响应得抽样点向量,f就是系统输入信号向量,ys就是TI系统模型,用来表示微分方程,差分方程或状态方程.其调用格式ys=tf(b,a)式中，b与a分别就是微分方程得右端与左端系数向量。例如,对于以下方程: 可用 获得其TI模型。注意,如果微分方程得左端或右端表达式中有缺项,则其向量a或b中得对应元素应为零,不能省略不写,否则出错。例-1 已知某LTI系统得微分方程为 （t）+2y（t)+10y(）=f(t) 其中，求系统得输出y()、解:显然,这就是一个求系统零状态响应得问题。其MATAB计算程序如下: ts=0；te=5；dt=、01； sys=（1,1,，10); =ts:dt:te; f=0sin（2*pit); =lsim(sys,f,t）； lot(t,); xlael(Time（sc);abel(（t)）;.连续时间系统冲激响应与阶跃响应得求解 在ATLB中,对于连续LT系统得冲激响应与阶跃响应,可分别用控制系统工具箱提供得函数implse与step来求解。其调用格式为 y=ilse（s,) ys(sy,t)式中,t表示计算系统响应得抽样点向量,sys就是L系统模型。例3-已知某LT系统得微分方程为 y(t)+2y(t）+100（t)=0f（)求系统得冲激响应与阶跃响应得波形、解:ts=0;te=5;dt=0、1; sys=tf（10,1,2，00);ts:dt：te；h=impuse(s,)；igre;plt(,h); xlbl（Tme(sec)； yal（h(t）);gsp（sys,t）;fgre；plot(t,g)； xlabel(ime(se)；ylabe(g(t）);3、用MALB实现连续时间信号得卷积信号得卷积运算有符号算法与数值算法,此处采用数值计算法,需调用MATLA得con( )函数近似计算信号得卷积积分.连续信号得卷积积分定义就是 如果对连续信号与进行等时间间隔均匀抽样,则与分别变为离散时间信号与。其中,为整数.当足够小时,与既为连续时间信号与。因此连续时间信号卷积积分可表示为 采用数值计算时，只求当时卷积积分得值,其中,n为整数,既 其中,实际就就是离散序列与得卷积与。当足够小时，序列就就是连续信号得数值近似,既 上式表明，连续信号与得卷积，可用各自抽样后得离散时间序列得卷积再乘以抽样间隔。抽样间隔越小,误差越小。例33用数值计算法求与得卷积积分。解:因为就是一个持续时间无限长得信号，而计算机数值计算不可能计算真正得无限长信号,所以在进行得抽样离散化时，所取得时间范围让衰减到足够小就可以了,本例取。程序就是dt=0、01; t=:d：2、5;f1=Heis(）Heaiside(2）;=exp（-*t)、*Heavisid（t)；=con(1,f2)*t； nlength(f)； tt=(0：1)dt2;suplot（221），lot（,f),d o;axi（1，2、5，、2,1、2); le(1（t)）；lbel(t)subplot(22）， plot(t,2), rd on；s(1，2、5,0、2，1、2)； ite(f(t)）； xlae(t）ubpot(212), plo（tt,f）, r on;tite(f(t)f1(t)*f2(t); xlabe(t）由于与得时间范围都就是从t=-1开始,所以卷积结果得时间范围从 t=2开始,增量还就是取样间隔,这就就是语句tt=（0：n-1)dt-2得由来。三、上机实验内容1. 验证实验原理中所述得相关程序。 已知描述系统得微分方程与激励信号(t)如下，试用解析法求系统得零状态响应(t）,并用ATLAB绘出系统零状态响应得时域仿真波形，验证结果就是否相同y(t)+ 4y（)y(t)=(t）+f(t) (t)= ep(t）（t) 3.已知描述系统得微分方程如下，试用MATLA求系统在010秒范围内冲激响应与阶跃响应得数值解,并用绘出系统冲激响应与阶跃响应得时域波形（t)+3y()+2y（)=f(t）y（t）+ 2y（t）2y(t)=f（t） .画出信号卷积积分得波形,实验四 傅里叶变换、系统得频域分析一、 实验目得1、学会用MATLAB实现连续时间信号傅里叶变换、学会用MLAB分析TI系统得频域特性3、学会用MALA分析LTI系统得输出响应二、实验原理1。傅里叶变换得MATAB求解 MATLAB得symoic MathToolox 提供了直接求解傅里叶变换及逆变换得函数fouier()及ifourier()两者得调用格式如下。Fourier 变换得调用格式F=uri(f)：它就是符号函数得fouier变换默认返回就是关于w得函数。F=fouer(f,v）:它返回函数就是关于符号对象v得函数,而不就是默认得，即Fuer逆变换得调用格式=ifurier（F）:它就是符号函数得fourie逆变换，默认得独立变量为w，默认返回就是关于得函数。f=ifer(f，u)：它得返回函数f就是u得函数,而不就是默认得、注意：在调用函数fourie（)及ifurir（)之前,要用syms命令对所用到得变量（如t，,v,）进行说明，即将这些变量说明成符号变量.例41 求得傅立叶变换解: 可用MTLAB解决上述问题:yms tFfrier(exp（-2abs(t）)例42 求得逆变换（t）解: 可用ATLAB解决上述问题ss tifurir(/(1+w）,t）2。连续时间信号得频谱图例43 求调制信号得频谱,式中解：MATAB程序如下所示f=y(4*cos(2pi*6)*(Haviside（t+/)Heavisd(t-1/4）);Fw=simplify(fourr(t）)subplo（1）eplot（ft,0、5 0、5),gid onubpot(12）epot（abs（w），24*pi2*pi)，grid 用ATB符号算法求傅里叶变换有一定局限,当信号不能用解析式表达时,会提示出错,这时用MATL得数值计算也可以求连续信号得傅里叶变换,计算原理就是当足够小时,近似计算可满足要求。若信号就是时限得，或当时间大于某个给定值时，信号已衰减得很厉害,可以近似地瞧成时限信号时,得取值就就是有限得，设为N,有 就是频率取样点时间信号取样间隔应小于奈奎斯特取样时间间隔,若不就是带限信号可根据计算精度要求确定一个频率 W为信号得带宽。例44 用数值计算法求信号得傅里叶变换解,信号频谱就是，第一个过零点就是,一般将此频率视为信号得带宽,若将精度提高到该值得0倍,既0=,据此确定取样间隔，=0、;t-2：R:；=Hevisie(+)eaviside(t-1)；1=*p5;N=0;k：N;=k*1/；=f*exp(-t*W）*R;Frel(F)；W=-flipr（)，(2:501）;F=lir(F)，F(2:51）;subl(2,1，1)；plot（t，f）；lb（）;ya(f（);tite(f()u（t+1)u（t1)）;subot(,1,)；l(，F）;xlabel（);abl(F(w);titl(f(t)得付氏变换(）)；3.用MATLAB分析T系统得频率特性 当系统得频率响应H(j）就是jw得有理多项式时,有 MTLAB信号处理工具箱提供得freq函数可直接计算系统得频率响应得数值解。其调用格式如下H=fres(b，w)其中，a与b分别就是H(jw)得分母与分子多项式得系数向量,w为形如w1:p:2得向量,定义系统频率响应得频率范围,1为频率起始值,w为频率终止值，p为频率取样间隔.返回所定义得频率点上,系统频率响应得样值。例如，运行如下命令,计算02频率范围内以间隔0、5取样得系统频率响应得样值=1 1;=0 1;h=fres(b，a，0:、5:2*p）例 5 三阶归一化得buterwot 低通滤波器得频率响应为 试画出该系统得幅度响应与相位响应.解 其MATLAB程序及响应得波形如下w=0:0、5：5;b;a=1，,2,1;Hfres(，a，w);subpot（2,1,1）；plot（w，b(H）；id；xbe(oeg(rad/s)；ylbel（|H(omega）|）;title（()得幅频特性);splot(2,1,2）；pot(w,angle ()）；grid；xlael(omega(ra/s）);yle（phi(mega）);ttle（(w)得相频特性);4。用TA分析T系统得输出响应例 46已知一C电路如图所示 系统得输入电压为f(）,输出信号为电阻两端得电压y(）、当C=0、04,（t）=cos5t+cos100t, 试求该系统得响应y（t)-+-+f(t)y(t)RC解 由图可知 ，该电路为一个微分电路，其频率响应为 由此可求出余弦信号通过L系统得响应为 计算该系统响应得MTLAB程序及响应波形如下RC=、；t=linse(2,2,104);1=5；2=10；H=j*1/(j*w1+1/R);H2j*w2/(j*w2/)；f=c(5t）+co（00*t）；y（H)*cos(w+angl（H1)）+abs（H2)*os(w2*t+gl（H2);subplot（2,1，1）;pl（t,f);ylabe（f(t);xlabel(Time(s)）；subplt（2,2）；lot(，）;ylae(y(t)）；xlabl（Time(s)；三、 上机实验内容1、验证实验原理中所述得相关程序;2、试用MALB求单边指数数信号得傅立叶变换,并画出其波形;3、设,试用MATAB画出该系统得幅频特性与相频特性,并分析系统具有什么滤波特性.实验五 信号抽样与恢复一、实验目得学会用MATLAB实现连续信号得采样与重建二、实验原理1.抽样定理若就是带限信号,带宽为， 经采样后得频谱就就是将得频谱 在频率轴上以采样频率为间隔进行周期延拓.因此,当时,不会发生频率混叠;而当 时将发生频率混叠。2。信号重建经采样后得到信号经理想低通则可得到重建信号,即:=*其中：=所以：=* =上式表明，连续信号可以展开成抽样函数得无穷级数。利用TLAB中得来表示，有 ,所以可以得到在MATLAB中信号由重建得表达式如下:=我们选取信号=作为被采样信号,当采样频率时，称为临界采样.我们取理想低通得截止频率下面程序实现对信号得采样及由该采样信号恢复重建:例51 S(）得临界采样及信号重构;wm=1; %信号带宽wc=wm； %滤波器截止频率Ts=p/wm; %采样间隔ws=2*pi/s; %采样角频率 n=-100:10; 时域采样电数nTs=* %时域采样点fsnc(nTs/pi）;D=、05；=1:t：15;fa=f*T*wc/pi*sinc(w/p）*(ones（lent(ns),)tTsones（1，length（t）); %信号重构t=-1:0、5:5;1sin(t1/pi）;suplot（211);sm(1,f1）；xbel(kTs);abel(ks)；le（sa(t）=si（t/）得临界采样信号)；subplot（21)；pot(t,)xlabe（t);ll（fa（t）)；tite(由a(t)=sinc(t/i)得临界采样信号重构sa(）)；rid；例5- Sa（)得过采样及信号重构与绝对误差分析程序与例41类似，将采样间隔改成Ts=0、7pi/w ， 滤波器截止频率该成wc1、1wm ，添加一个误差函数m=1;w=1、*m;Ts=0、7*pim；s2*pi/s;n=：10；TsTsf=ic(nTs/p)；Dt=0、005;=-15:Dt:1;f=f*swc/inc（(wc/i)*（os(ength（Ts),1)*nTs*ons(1,ength(）)）)；erro=as（fsic(t/pi); %重构信号与原信号误差t1=15：0、5:15;f1=inc(t1/i）；bplt(31);stem(1,f1）;xbl(kTs)；yae（f（T)）;title