第6章 MATLAB在信号与系统中的应用.ppt

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,第,6,章,MATLAB,在信号与系统中的应用,6.1,连续信号及其表示,6.2,线性时不变系统,6.3,连续时间系统的时域分析,6.4,连续时间系统的变换域分析,例,6-1,用,MATLAB,命令绘制连续时间信号,MATLAB,程序,t=0:0.01:10;%,设定时间向量,f=3*exp(-0.1*t).*cos(2*t);%,产生信号,f(t,),plot(t,f,)%,绘制,f(t,),xlabel(t),ylabel(f(t,)%,设定坐标轴名,程序运行结果如图：,6.1,连续信号及其表示,MATLAB,中，用等时间间隔点的样值来近似表示连续信号。,例,6-2,用,MATLAB,命令绘制符号函数和阶跃信号。,MATLAB,程序,t=-200:200;%,设定时间向量,f1=,sign(t,);%,产生符号函数,f2=1/2+f1/2;%,产生单位阶跃信号,subplot(1,2,1),plot(t,f1)%,绘制符号函数,axis(-200 200-1.1 1.1)%,设置坐标范围,xlabel(t),ylabel(sgn(t,)%,设置坐标轴名,title(,符号函数,)%,设置图名,subplot(1,2,2),plot(t,f2)%,绘制单位阶跃信号,axis(-200 200-0.1 1.1)%,设置坐标范围,title(,单位阶跃信号,)%,设置图名,xlabel(t),ylabel(u(t,)%,设置坐标轴名,程序运行结果如图,在,MATLAB,中，还用函数,stepfun,(),实现，其调用格式为,stepfun(T,T0),：,T,为时间向量；,T0,为跃变时刻；,当,TT0,时返回,1,，且返回的向量与,T,具有相同的长度,。,t=linspace(-1,4,1000);,x=t1;,y=stepfun(t,1);,subplot(1,2,1),plot(t,x),xlabel(t),ylabel(x,),axis(-1.1 4.1-0.1 1.1),subplot(1,2,2),plot(t,y),xlabel(t),ylabel(y,),axis(-1.1 4.1-0.1 1.1),表,6-1,常用的信号产生函数,函数名,函数功能,函数名,函数功能,square,产生周期方波信号,rectpuls,产生非周期方波信号,sawtooth,产生周期锯齿波、三角波信号,tripuls,产生非周期三角波信号,sinc,产生,sinc,信号,pulstran,产生冲激串,chirp,产生调频余弦信号,diric,产生,Dirichlet,或周期,sinc,函数,gauspuls,产生高斯正弦脉冲信号,gmonopuls,产生高斯单脉冲信号,例,6-3,产生一个频率为,50Hz,的周期锯齿波。,MTALAB,程序,t=0:0.0001:0.2;%,设定时间向量,f=2*sawtooth(2*pi*50*t);%,产生锯齿波信号,plot(t,f,)%,绘制锯齿波信号,title(,锯齿波信号,)%,加图名,6.2,线性时不变系统,参数不随时间改变，且满足叠加性和均匀性的系统称为,线性时不变（,linear time invariant,LTI,）系统。,分析,LTI,系统可采用时域方法或变换域方法如傅立叶变换、拉普拉斯变换和,Z,变换。,LTI,系统可分为连续时间系统和离散时间系统,描述,连续,LTI,系统,方法：用常系数微分方程、系统的传递函数或状态方程。,(2),连续系统的传递函数,函数,功能,调用格式,tf2ss,传递函数模型转换为状态空间模型,A,B,C,D=tf2ss(b,a),tf2zp,传递函数模型转换为零,-,极点增益模型,z,p,k,=tf2zp(b,a),ss2tf,状态空间模型转换为传递函数模型,b,a,=ss2tf(A,B,C,D,iu),iu,为输入量,ss2zp,状态空间模型转换为零,-,极点增益模型,z,p,k,=ss2zp(A,B,C,D,iu),iu,为输入量,zp2tf,零,-,极点增益模型转换为传递函数模型,b,a,=zp2tf(z,p,k),zp2ss,零,-,极点增益模型转换为状态空间模型,A,B,C,D=zp2ss(z,p,k),表,6-2,系统模型转换函数,b=1 5 6;a=2 6 0 4;,A,B,C,D=tf2ss(b,a),A=,-3 0 -2,1 0 0,0 1 0,B=,1,0,0,C=,0.5000 2.5000 3.0000,D=,0,z,p,k,=tf2zp(b,a),z=,-3.0000,-2.0000,p=,-3.1958,0.0979+0.7850i,0.0979-0.7850i,k=,0.5000,例,6-4,某连续系统的,传递函数,为,求该系统的状态空间模型和零,-,极点增益模型。,例,6-5,某系统的传递函数为,求该系统的系数向量和状态空间模型。,z=-1-2;p=-3-4 5;k=2;,b,a,=zp2tf(z,p,k),b=0 2 6 4,a=1 2 -23 -60,A,B,C,D=zp2ss(z,p,k),A=,5.0000 0 0,1.0000 -7.0000 -3.4641,0 3.4641 0,B=,1,0,0,C=2.0000 -8.0000 -5.7735,D=0,A=-2 1;-3 0;B=2 1;C=1 0;D=1;,b,a,=ss2tf(A,B,C,D),b=,1.0000 4.0000,4.0000,a=,1.0000 2.0000 3.0000,z,p,k,=ss2zp(A,B,C,D),z=,-2.0000+0.0000i,-2.0000-0.0000i,p=,-1.0000+1.4142i,-1.0000-1.4142i,k=1,6.3,连续时间系统的时域分析,LTI,连续时间系统在时域中以常系数微分方程来描述：,可利用,dsolve,函数求解常系数微分方程来得到。,解：先求零输入响应，再求零状态响应，然后两部分加在一起。,MATLAB,程序,eq1=D2y+5*Dy+6*y=0;%,设定零输入条件下的微分方程,ic1=y(0)=1,Dy(0)=1;%,输入初始状态,yzi,=dsolve(eq1,ic1);%,求解微分方程，得到零输入响应,yzi,=,simplify(yzi,)%,化简零输入响应,eq2=D2y+5*Dy+6*y=exp(-t)*,Heaviside(t,);,%,设定给定输入条件下的微分方程,ic2=y(-0.01)=0,Dy(-0.01)=0;%,设定初始状态为,0,yzs,=dsolve(eq2,ic2);%,求解微分方程，得到零状态响应,yzs,=,simplify(yzs,)%,化简零状态响应,y=,yzi+yzs,%,全响应,运行程序，得到：,yzi,=,4*exp(-2*t)-3*exp(-3*t),yzs,=,1/2*,heaviside(t,)*(exp(-t)-2*exp(-2*t)+exp(-3*t),y=,4*exp(-2*t)-3*exp(-3*t)+1/2*,heaviside(t,)*(exp(-t)-2*exp(-2*t)+exp(-3*t),MATLAB,程序,eq1=D2y+5*Dy+6*y=,Dirac(t,);,%,设定输入为冲激信号时的微分方程,eq2=D2y+5*Dy+6*y=,Heaviside(t,);,%,设定输入为阶跃信号时的微分方程,ic,=y(-0.01)=0,Dy(-0.01)=0;%,设定初始状态为,0,h=dsolve(eq1,ic);%,求解微分方程，得到冲激响应,h=,simplify(h,)%,化简冲激响应,g=dsolve(eq2,ic);%,求解微分方程，得到阶跃响应,g=,simplify(g,)%,化简阶跃响应,求解例,6-7,中系统的冲激响应和阶跃响应,程序运行结果为：,h=-,heaviside(t,)*(exp(-3*t)-exp(-2*t),g=-1/6*,heaviside(t,)*(-1+3*exp(-2*t)-2*exp(-3*t),。,MATLAB,程序,delta=0.01;%,采样时间间隔,k1=0:delta:5;%,信号,f1,对应的时间向量,k2=0:delta:5;%,信号,f2,对应的时间向量,f1=sin(k1);%,输入信号,f1,f2=exp(-2*k2);%,输入信号,f2,f=ctconv(f1,k1,f2,k2,delta);,%,调用子函数,ctconv,计算,f1,和,f2,的卷积,%,子函数,ctconv,的定义如下,function f=ctconv(f1,k1,f2,k2,delta),f=delta*conv(f1,f2);%,计算卷积，将结果赋给,f,k0=k1(1)+k2(1);%,计算,f,序列非零样值的起始位置,k3=length(f1)+length(f2)-2;%,计算,f,序列的长度,k=k0:delta:k0+k3*delta;%f,序列对应的时间向量,例,6-8,已知,计算,时的卷积,subplot(2,2,1);plot(k1,f1);%,在子图,1,中绘制,f1(t),xlabel(t);ylabel(f1(t);title(f1(t);,subplot(2,2,2);plot(k2,f2);%,在子图,2,中绘制,f2(t),xlabel(t);ylabel(f2(t);title(f2(t);,subplot(2,2,3);plot(k,f);%,在子图,3,中绘制,f(t,),h=,get(gca,position,);,h(3)=2.3*h(3);set(gca,position,h);%,将子图,3,的横坐标范围扩至原来的,2.3,倍,xlabel(t);ylabel(f(t);title(f(t,)=f1(t)*f2(t);,MATLAB,程序,delta=0.01;%,采样时间间隔,k1=0:delta:5;%,设定,f,对应的时间向量,k2=0:delta:4;%,设定,h,对应的时间向量,f=exp(-k1);%,激励信号,f1,a=1 5 6;%,系数向量,b=1;,r,p,k,=,residue(b,a,);%,极点留数法求冲激响应,h=r(1)*exp(p(1)*k2)+r(2)*exp(p(2)*k2);,y=ctconv(f,k1,h,k2,delta);%,调用子函数,ctconv,计算零状态响应,(1)impulse(),：单位冲激响应，调用格式：,impulse(sys,),：绘制系统的单位冲激响应。,impulse(sys,t,),：绘制在指定时间范围内系统的单位冲激响应。,y=,impulse(sys,t,),：参数同上，不绘图，将结果保留在,y,中。,MATLAB,提供 的求解系统响应的函数,。,(2)step(),函数：求阶跃响应，调用格式：,step(sys),：绘制系统的阶跃响应。,step(,sys,t,),：绘制在指定时间范围内系统的阶跃响应。,y=,step(sys,t,),：参数同上，不绘图，将结果保留在,y,中。,(3),lsim,(),函数：求时域响应，调用格式：,lsim(sys,u,t,),：在,t,指定的时间范围内绘制当输入为,u,时系统的零状态响应波形。,y=,lsim(sys,u,t,),：参数同上，不绘图，将结果保留在,y,中。,lsim(sys,u,t,x0),：在,t,指定的时间范围内绘制当输入为,u,时系统的全响应波形。,y=lsim(sys,u,t,x0),：参数同上，不绘图，将结果保留在,y,中。,(4),initial(),函数：求零输入响应调用格式：,initial(sys,x0),：绘制当初始状态为,x0,时系统的零输入响应波形。,initial(sys,x0,t),：在,t,指定的时间范围内绘制当初始状态为,x0,时系统的零输入响应波形。,y=initial(sys,x0,t),：参数同上，不绘图，将结果保留在,y,中。,注意：调用,impulse(),、,step(),和,lsim,(),时的系统模型可以是传递函数模型，也可以是状态空间模型，但调用,initial(),时的系统模型,sys,必须是状态空间模型。,MATLAB,程序,a=2 1 4;%,设定系数向量,a,b=1 1;%,设定系数向量,b,sys=,tf(b,a,);%,得到系统的传递函数模型,A,B,C,D=tf2ss(b,a);%,由传递函数模型求得状态空间模型,sys1=,ss(A,B,C,D,);%,状态空间模型,y0=1 0;%,输入系统的初始条件,t=0:0.01:5;%,设定时间向量,f=exp(-t);%,系统输入,subplot(2,2,1);impulse(sys,t);,%,impulse(b,a,t,),绘制系统的单位冲激响应,subplot(2,2,2);step(sys,t);%,绘制系统的阶跃响应,subplot(2,2,3);initial(sys1,y0,t);,%initial(A,B,C,D,y0,t),绘制系统的零输入响应,subplot(2,2,4);lsim(sys,f,t);%,绘制系统的零状态响应,%,也可以先求出各种响应，再绘制出来，具体如下,y1=,impulse(sys,t,);subplot(2,2,1);plot(t,y1);,y2=,step(sys,t,);subplot(2,2,2);plot(t,y2);,y3=initial(sys1,y0,t);subplot(2,2,3);plot(t,y3);,y4=,lsim(sys,f,t,);subplot(2,2,4);plot(t,y4);,%,直接绘制各种响应,6.4,连续时间系统的变换域分析,MATLAB,程序,T=4;tao=2;%,设定矩形脉冲信号的周期和脉冲宽度,w=2*pi/T;%,基频,a0=,quadl(abs(t,)=1,-2,2)/T;%,求,a0,N=100;%,设定最大谐波数,an=zeros(1,N);bn=zeros(1,N);%,初始化,an,和,bn,for k=1:N,an(k,)=quadl(rectcos,-2,2,k,w)*2/T;%,计算,an,bn(k,)=quadl(rectsin,-2,2,k,w)*2/T;%,计算,bn,end,例,6-11,利用,MATLAB,将重复周期为,4,、脉冲宽度为,2,、脉冲幅度为,1,的矩形波信号展开成傅立叶级数。,x=pulstran(t,-6:4:6,rectpuls,2);%,生成周期矩形脉冲信号,subplot(3,2,1);plot(t,x);%,绘制周期矩形脉冲信号,xlabel(t);ylabel(Amplitude);axis(-4,4,-0.2,1.2);grid on;,title(Rectangular,wave);,%,有限项级数逼近,k_max,=1 10 30 50 100;%,最大谐波数向量,for m=1:length(k_max),f=a0;,for k=1:k_max(m),f=,f+an(k,)*,cos(k,*w*,t)+bn(k,)*,sin(k,*w*t);,end,subplot(3,2,m+1);plot(t,f);grid on;,xlabel(t);ylabel(Partial,Sum);,axis(-4,4,-0.2,1.2);,title(Max,Harmonics=,num2str(k_max(m);,end,%,绘制周期矩形脉冲信号,t=-4:0.01:4;%,设定时间向量,其中用到两个自定义的子函数：,rectcos,和,rectsin,，它们的定义如下：,function y=,rectcos(t,n,w,),y=(,abs(t,)=1).*,cos(n,*w*t);,function y=,rectsin(t,n,w,),y=(,abs(t,)0,phi(k,)=pi;,end,end,end,%,绘制幅度谱和相位谱,n=1:N;,subplot(1,2,1);stem(n*,w,cn,);,xlabel(nw);ylabel(cn,);,axis(0,21,0,0.6);,title(,幅度谱,);,subplot(1,2,2);stem(n*,w,phi,);,xlabel(nw);ylabel(phi,);,axis(0,21,0,3.5);,title(,相位谱,);,例,6-12,中周期矩形波的脉冲宽度不变，将重复周期由原来的,6,变为,10,，重绘其幅度谱，与上图,6-9,相比较，可以看出，当脉冲宽度不变、周期增大时，谱线的幅度随之减小，谱线之间的间隔变小，谱线更密。,例,6-12,中周期矩形波的周期不变，脉冲宽度由原来的,2,变为,1,，其幅度谱如下，与图,6-9,相比较，可以看出，当周期不变、脉冲宽度减小时，谱线的幅度随之减小，信号频带宽度增加。,6.4.2,傅立叶变换,MATLAB,程序,syms,t;%,定义符号变量,t,f=sym(Heaviside(t+1)-sym(Heaviside(t-1);%,输入符号函数,f,F=,fourier(f,);%,求,f,的傅立叶变换,F,ezplot(abs(F,);%,绘制幅度谱,例,6-14,进行逆傅立叶变换,syms,x w u v;,F1=exp(-w2/4);%,符号函数的变量是,w,f1=ifourier(F1)%,返回函数的参量是默认变量,x,f1=,1/pi(1/2)*exp(-x2),F2=exp(-,abs(x,);%,符号函数的变量是,x,f2=ifourier(F2)%,返回函数的参量是变量,t,f2=,1/(1+t2)/pi,F3=2*exp(-abs(w)-1;%,符号函数的变量是,w,f3=ifourier(F3,u)%,返回函数的参量是指定变量,u,f3=,-dirac(u)+2/pi/(1+u2),T=4;%,时间宽度,N=200;%,时域采样点数,delta=T/N;%,时域采样间隔,t=linspace(-T/2,T/2-T/N,N);%,时间向量,f=(t=-1)-(t=1);%,f(t,),f1=f.*exp(-j*10*t);%f1(t),f2=f.*,exp(j,*10*t);%f2(t),W=20*pi;%,频域范围,K=1000;%,频域采样点数,w=linspace(-W/2,W/2-W/K,K);%,频率序列,F=f*exp(-j*t*w)*delta;%,计算,f(t,),的傅立叶变换,F(jw,),F1=,f1,*exp(-j*t*w)*delta;%,计算,f1(t),的傅立叶变换,F1(jw),F2=,f2,*exp(-j*t*w)*delta;%,计算,f2(t),的傅立叶变换,F2(jw),F=,real(F,);,F1=real(F1);,F2=real(F2);,subplot(3,1,1);plot(w,F);,xlabel(w);ylabel(F(jw,);,title(,频谱,F(jw,),subplot(3,1,2);plot(w,F1);,xlabel(w);ylabel(F1(jw);,title(F(jw,),左移到,w=-10,处的频谱,F1(jw),subplot(3,1,3);plot(w,F2);,title(F(jw,),右移到,w=10,处的频谱,F2(jw),6.4.3,连续时间系统的频域分析,MATLAB,提供了求解系统,频率响应的函数,freqs,(),，,其调用格式如下：,(1)h=,freqs(b,a,w,),：返回在指定频率范围,w,内系统的频率响应，其中,b,和,a,是系统对应的常系数微分方程的系数向量。,(2),h,w,=,freqs(b,a,),：返回默认频率范围内,200,个频率点的系统频率响应的样值和这,200,频率点。,(3),h,w,=,freqs(b,a,f,),：返回默认频率范围内,f,个频率点的系统频率响应的样值和这,f,个频率点。,(4),freqs(b,a,),：以对数坐标绘制系统的幅频特性和相频特性。,MATLAB,程序,b=2;a=1 5 6;%,分子、分母多项式系数向量,w=logspace(-1,1);%,设定频率向量,figure(1);freqs(b,a,w);%,绘制系统的频率响应,%,还可以先求出系统的频率响应，再绘制其幅频特性和相频特性,H=,freqs(b,a,w,);%,求频率响应,mag,=,abs(H,);%,取频率响应的模,phase=,angle(H,)*180/pi;%,取频率响应的相位，并转换成角度,figure(2);,subplot(2,1,1),loglog(w,mag)%,以对数坐标绘制幅频特性曲线,xlabel(Frequency(rad/s),ylabel(Magnitude),grid,subplot(2,1,2),semilogx(w,phase)%,以对数坐标绘制相频特性曲线,xlabel(Frequency(rad/s),ylabel(Phase(degrees),grid,%,下面是以线性坐标绘制频响特性曲线,figure(3);,subplot(2,1,1),plot(w,mag)%,以线性坐标绘制幅频特性曲线,xlabel(Frequency(rad/s),ylabel(Magnitude),grid,subplot(2,1,2),plot(w,phase)%,以线性坐标绘制相频特性曲线,xlabel(Frequency(rad/s),ylabel(Phase(degrees),grid,直接绘制系统的,幅频特性和相频特性,(b),先求系统频响，再以对数,坐标绘制其特性曲线,(c),先求系统频响，再以线性坐标绘制其特性曲线,例,6-17,某一,LTI,连续时间系统由 表征，求该系统的频率响应。,MATLAB,程序,b=1 8-6 4;a=10 5;%,设定分子、分母多项式系数向量,w=-20:0.05:20;%,设定频率向量,H=,freqs(b,a,w,);%,求频率响应,mag,=,abs(H,);%,取频率响应的模,phase=,angle(H,)*180/pi;%,取频率响应的相位，并转换成角度,subplot(2,1,1),plot(w,mag)%,以线性坐标绘制幅频特性曲线,xlabel(Frequency(rad/s),ylabel(Magnitude),grid,subplot(2,1,2),plot(w,phase)%,以线性坐标绘制相频特性曲线,xlabel(Frequency(rad/s),ylabel(Phase(degrees),grid,图,6-15,系统的幅频特性和相频特性,6.4.4,拉普拉斯变换,例,6-18,进行,Laplace,变换,syms,t s a z w,f1=t4;%,符号函数的自由变量为,t,L1=laplace(f1)%,返回函数的参量为默认变量,s,L1=24/s5,f2=1/sqrt(s);%,符号函数的自由变量为,s,L2=laplace(f2)%,返回函数的参量为变量,t,L2=(pi/t)(1/2),f3=exp(-a*t);%,符号函数的自由变量为,t,L3=laplace(f3,z)%,返回函数的参量为指定变量,z,L3=1/(z+a),f4=,abs(t,)*exp(-a*,abs(w,);%,符号函数有两个自由变量：,t,和,w,L4=laplace(f4,w,z)%,对参量为,w,的函数进行拉氏变换，,%,返回函数的参量为指定变量,z,L4=,abs(t)/(z+a,),MATLAB,程序,syms,s;,F=sym(2*s/(s2+1)2);,f=,ilaplace(F,);,subplot(121);ezplot(F);,subplot(122);ezplot(f);,6.4.5,连续时间系统的复频域分析,MATLAB,程序,b=1-4 5;a=3 2-4 6;%,分子、分母多项式系数向量,%,求零极点,zeros=,roots(b,);%,求系统函数的零点,poles=,roots(a,);%,求系统函数的极点,%,利用,plot,绘制系统函数的零极点图,x=,max(max(abs(zeros),max(abs(poles,);,%,确定离坐标原点最远的零极点的距离,x=x+0.1;%,确定横坐标范围,y=x;%,确定纵坐标范围,axis(-x x-y y);%,设置坐标轴的显示范围,figure(1),hold on,plot(-x x,0 0);%,绘制横坐标轴,plot(0 0,-y y);%,绘制纵坐标轴,plot(real(zeros),imag(zeros),o,);%,绘制零点,plot(real(poles),imag(poles),X,);%,绘制极点,title(Pole,-Zero Map);%,加图题,text(0.1,x-0.2,Imaginary Axis);%,标注虚轴,text(y-0.5,0.2,Real Axis);%,标注实轴,%,直接利用,pzmap,函数绘制系统函数的零极点图,sys=,tf(b,a,);%,系统函数,figure(2),pzmap(sys);%,绘制系统函数的零极点图,%,绘制系统的幅频特性和相频特性,figure(3),freqs(b,a);,运行程序，系统函数的零点为：,zeros=2.0000+1.0000i 2.0000-1.0000i,极点为：,poles=-1.9115 0.6224+0.8117i 0.6224-0.8117i,(a),用,plot,绘制的系统函数零极点图,(b),用,pzmap,函数绘制的系统函数零极点图,MATLAB,程序,b=1 4;,a=1 6 11 6;,r,p,k,=,residue(b,a,),程序运行结果为：,r=,0.5000,-2.0000,1.5000,p=,-3.0000,-2.0000,-1.0000,k=,