信号与系统(陈后金)_MATLAB

信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLABnMATLAB简介简介n信号的信号的MATLAB表示表示n利用利用MATLAB进行系统的时域分析进行系统的时域分析n利用利用MATLAB进行信号的频域分析进行信号的频域分析n利用利用MATLAB分析系统的频率特性分析系统的频率特性n利用利用MATLAB进行连续系统的进行连续系统的s域分析域分析n利用利用MATLAB进行离散系统的进行离散系统的z域分析域分析n利用利用MATLAB进行系统的状态变量分析进行系统的状态变量分析信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB（1）窗口命令方式（2）运行以 .M 为扩展名磁盘文件信号与系统(陈后金)_MATLAB%用plot函数画一个方波t=-1 0 0 1 1 3;x=0 0 1 1 0 0;plot(t,x);xlabel(t);ylabel(x(t);axis(-1 3 0 2);直接在命令窗口输入以上命令建一个名字为my_file.M的文件，然后在命令窗口输入文件名及回车。信号与系统(陈后金)_MATLAB命令窗口输入： help+函数名例如 help plot信号与系统(陈后金)_MATLAB不需要变量的类型说明变量名的第一个字符必须是字母变量名长度：不超过31个字符大写和小写的字母视为不同的字符 例如：num_students = 25特殊变量： pi 表示圆周率，inf 表示无穷大，NaN(Not a Number)表示不定量，如0/0。信号与系统(陈后金)_MATLAB MATLAB用常规的十进制表示数值 用i或j作为后缀来表示复数的虚部 例 1.235e5表示1.235105，x=2+3jabs(x) 求复数x的模angle(x) 求复数x的相角(弧度)real(x) 求复数x的实部imag(x) 求复数x的虚部conj(x) 求复数x的共轭信号与系统(陈后金)_MATLABu + 加u - 减u * 乘u / 除u 乘方u 矩阵的复共轭转置信号与系统(陈后金)_MATLABu A & B逻辑与(and)u A | B逻辑或(or)u A 逻辑非(not) 值为0时表示逻辑假(F)，其它任何非零值表示逻辑真。信号与系统(陈后金)_MATLABu A B 大于uA = B 大于等于uA = B 等于uA = B不等于信号与系统(陈后金)_MATLABu表达式 1:10 表示产生一个行向量，它的值为 1 2 3 4 5 6 7 8 9 10u表达式 10:-2:1 表示产生一个递减的行向量，它的值为 10 8 6 4 2 信号与系统(陈后金)_MATLAB用冒号:产生数组例 x=2:5 产生一个数组，它的值为x(1)=2, x(2)=3, x(3)=4, x(4)=5x=linspace(0,2,11)将区间0，2均匀抽样11点作为数组x 给2维数组赋值时，用分号表示一行的结束，如：z=1 2; 3 4。用linspace产生数组信号与系统(陈后金)_MATLABMATLAB 提供了一些产生基本矩阵的函数 zeros 产生矩阵元素全为0的矩阵 ones 产生矩阵元素全为1的矩阵 rand 产生(0,1)均匀分布随机数矩阵 randn 产生正态分布随机数矩阵信号与系统(陈后金)_MATLABsize(A) 返回值数组A的行数和列数（二维） length(B) 确定数组B的元素个数（一维）信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB例：编一个绘制图示波形的函数。function y=tri(t)y= abs(t) eps s=s+1/(n*n); n=n+1; end+222312111s例例: 计算 的值，且误差小于10-6fprintf(s=%.5fn,s)信号与系统(陈后金)_MATLABt=linspace(0,4*pi,512);plot(t,sin(t),t,cos(t),-.);title(my figure);xlabel(t);legend(sin(t),cos(t);信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLABk=0:39;stem(k,cos(0.9*pi*k);title(cos(0.9pik);信号与系统(陈后金)_MATLAB0510152025303540-1-0.8-0.6-0.4-0.200.20.40.60.81cos(0.9k)cos(0.9k)波形波形信号与系统(陈后金)_MATLAB指数信号Aeat 、指数序列ak 、抽样函数Sa(t)、 正弦型信号、矩形脉冲信号、三角脉冲信号尺度变换、翻转、时移、尺度变换、翻转、时移、 相加、相乘、相加、相乘、 差分与求和、微分与积分差分与求和、微分与积分信号与系统(陈后金)_MATLAB指数信号Aeat y = A*exp(a*t);指数序列ak幂运算a.k实现正弦型信号内部函数cos( ) 和sin( )抽样函数Sa(t) sinc(t) 矩形脉冲信号 y = rectpuls(t,width)三角波脉冲信号y = tripuls(t, width,skew) 信号与系统(陈后金)_MATLAB%decaying exponentialt=0:001:10;A=1;a=-0.4;ft=A*exp(a*t);plot(t,ft)01234567891000.10.20.30.40.50.60.70.80.91t=0:0.1:10;A=1;a=-0.4;ft=A*exp(a*t);stem(t,ft)01234567891000.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB% rectpulst=0:0.001:4;T=1;ft=rectpuls(t-2*T,T);plot(t,ft)axis(0,4,-0.5,1.5)00.511.522.533.54-0.500.511.5信号与系统(陈后金)_MATLAB-3-2-1012300.10.20.30.40.50.60.70.80.91-3-2-1012300.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB-50-40-30-20-100102030405000.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB-50-40-30-20-100102030405000.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB已知三角波x(t)，用MATLAB画出的x(2t)和x(2-2t) 波形-3-2-1012300.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB相加相加用算术运算符用算术运算符“+”实现实现相乘相乘用数组运算符用数组运算符“.*”实现实现 例:画信号Aeatcos(w0t+f)的波形 t=0:0.001:8; A=1; a=-0.4; w0=2*pi;phi=0; ft1=A*exp(a*t).*sin(w0*t+phi); plot(t,ft1)012345678-0.8-0.6-0.4-0.200.20.40.60.81信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB:已知三角波x(t)，画出其微分与积分的波形%differentiation h=0.001;t= -3:h:3;y1=diff(f2_2(t)*1/h;plot(t(1:length(t)-1),y1)%integrationt= -3:0.1:3;for x=1:length(t) y2(x)=quad(f2_2, -3,t(x);endplot(t,y2)信号与系统(陈后金)_MATLAB三角波x(t)微分与积分的波形-3-2-1012300.10.20.30.40.50.60.70.80.91信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLABt 表示计算系统响应的抽样点向量)()( )( )( )()( )( )( 01230123txbtxbtxbtxbtyatyatyatya+a=a3, a2, a1, a0; b=b3, b2, b1, b0; sys=tf(b,a)y=lsim(sys,x,t)sys=tf(b,a)b和a分别为微分方程右端和左端各项的系数向量x 是系统输入信号向量，sys 是LTI系统模型，借助tf函数获得信号与系统(陈后金)_MATLAB 连续时间系统冲激响应可用impulse函数直接求出，其调用形式为y=impulse(sys, t) 连续时间系统阶跃响应可用step函数直接求出，其调用形式为y=step(sys, t)t 表示计算系统响应的抽样点向量sys 是LTI系统模型信号与系统(陈后金)_MATLABb , a 分别是差分方程左、右端的系数向量00jkxbikyajmjini-b=b0,b1,b2,bM;a=a0,a1,a2, ,aN;可用MATLAB表示为y=filter(b,a,x)x 表示输入序列, y 表示输出序列信号与系统(陈后金)_MATLABb, a 分别是差分方程左、右端的系数向量k 表示输出序列的取值范围h 就是单位脉冲响应h=impz(b,a,k)信号与系统(陈后金)_MATLAB例：(s3+2s+3)(s2+3s+2)可用下面MATLAB语句求出 a =1,0,2,3; b =1,3,2; c=conv(a,b)c=conv(a,b)式中a,b为待卷积两序列的向量表示，c是卷积结果。conv函数也可用于计算两个多项式的积信号与系统(陈后金)_MATLAB%program3_1微分方程求解ts=0;te=5;dt=0.01;sys=tf(10,1 2 100);t=ts:dt:te;x=sin(2*pi*t);y=lsim(sys,x,t);plot(t,y);xlabel(Time(sec)ylabel(y(t)012345-0.2-0.100.10.2Time(sec)y(t)信号与系统(陈后金)_MATLAB%program3_2连续时间系统的冲激响应ts=0;te=5;dt=0.01;sys=tf(10,1 2 100);t=ts:dt:te;y=impulse(sys,t);plot(t,y);xlabel(Time(sec)ylabel(h(t)012345-1-0.500.51Time(sec)h(t)信号与系统(陈后金)_MATLABR =51 ; d = rand(1,R) - 0.5;k=0:R-1;s=2*k.*(0.9.k); x=s+d;figure(1); plot(k,d,r-.,k,s,b-,k,x,g-);M =5; b = ones(M,1)/M; a = 1;y = filter(b,a,x);figure(2); plot(k,s,b-,k,y,r-);110nkxMkyMn-信号与系统(陈后金)_MATLAB例例3 3 分析噪声干扰的信号xk=sk+dk通过M点滑动平均系统的响应，其中sk=(2k)0.9k是原始信号，dk是噪声。110nkxMkyMn- 0 10 20 30 40 50 -2 0 2 4 6 8 Time index k dk sk xk 01020304050-202468Time index kskyk噪声干扰信号xk=sk+dk通过M点滑动平均系统的响应信号与系统(陈后金)_MATLAB% program 3_4 离散系统的单位脉冲响应k=0:10;a=1 3 2;b=10; h=impz(b,a,k);stem(k,h)012345678910-2000-10000100020003000单位脉冲响应信号与系统(陈后金)_MATLAB例例5 计算xk* yk并画出卷积结果，已知xk=1,2,3,4; k=0,1,2,3，yk=1,1,1,1,1; k=0,1,2,3,4 。% program 3_5x=1,2,3,4; y=1,1,1,1,1; z=conv(x,y);N=length(z);stem(0:N-1,z);012345670246810信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB 频谱Cn一般为复数，可分别利用abs和angle函数获得其幅度频谱和相位频谱。其调用格式分别为x=abs(Cn) y=angle(Cn) 周期信号的频谱Cn 为离散信号，可以用stem画出其频谱图。信号与系统(陈后金)_MATLAB例例1 试用MATLAB画出图示周期三角波信号的频谱。解：解：周期信号的频谱为-000 )2sin(j422nnnnCn信号与系统(陈后金)_MATLABN=8;n1= -N:-1; %计算n=-N到-1的Fourier系数c1= -4*j*sin(n1*pi/2)/pi2./n1.2;c0=0; %计算n=0时的Fourier系数n2=1:N; %计算n=1到N的Fourier系数c2= -4*j*sin(n2*pi/2)/pi2./n2.2;cn=c1 c0 c2;n= -N:N;subplot(2,1,1);stem(n,abs(cn);ylabel(Cn的幅度);subplot(2,1,2);stem(n,angle(cn);ylabel(Cn的相位);xlabel(omega/omega0);信号与系统(陈后金)_MATLAB-8-6-4-202468-2-1012Cn的相位w /w 0-8-6-4-20246800.10.20.30.40.5Cn的幅度信号与系统(陈后金)_MATLAB例例2 求周期矩形脉冲的Fourier级数表示式。并用MATLAB求出由前N项Fourier级数系数得出的信号近似波形。)2/(Sa0wnTACntnnTATACtxnntnnT001jcos)2(Sa2e)(0www-+取A=1, T=2, =1, w0=)cos()2(Sa5 . 0)(1tnntxnT+信号与系统(陈后金)_MATLABt=-2:0.001:2;N=input(Number of harmonics= );c0=0.5;xN=c0*ones(1,length(t); %dc componentfor n=0:2:N % even harmonics are zero xN=xN+cos(pi*n*t)*sinc(n/2);endplot(t,xN);信号与系统(陈后金)_MATLAB-2-1.5-1-0.500.511.52-0.200.20.40.60.811.2-2-1.5-1-0.500.511.52-0.200.20.40.60.811.2-2-1.5-1-0.500.511.52-0.200.20.40.60.811.2-2-1.5-1-0.500.511.52-0.200.20.40.60.811.2N=5N=15N=50N=500信号与系统(陈后金)_MATLAB数值函数积分quad8可用来计算非周期信号频谱F 是一个字符串，它表示被积函数的文件名。a,b 分别表示定积分的下限和上限quad8的返回是用自适应Simpson算法得出的积分值y = quad8(F,a,b) 信号与系统(陈后金)_MATLAB例例3 试用数值方法近似计算三角波信号的频谱)()1 ()(21tpttx-解解: 图示三角波可表示为三角波信号频谱的理论值为信号与系统(陈后金)_MATLAB例例3 试用数值方法近似计算三角波信号的频谱function y=sf1(t,w);y=(t=-1 & t=1).*(1-abs(t).*exp(-j*w*t);w=linspace(-6*pi,6*pi,512);N=length(w);X=zeros(1,N);for k=1:N X(k)=quad8(sf1,-1,1,w(k);endfigure(1);plot(w,real(X);title()xlabel(omega);ylabel(X(jomega);figure(2);plot(w,real(X)-sinc(w/2/pi).2);xlabel(omega);title(计算误差);信号与系统(陈后金)_MATLAB例例3 试用数值方法近似计算三角波信号的频谱运行结果运行结果-20-15-10-505101520-20246x 10-9w计算误差信号与系统(陈后金)_MATLAB函数fft可用来计算离散周期信号频谱x 是离散周期信号0N-1 一个周期的序列值Xm = fft(x) Xm 是离散周期信号频谱在0N-1 的值函数fft还可用来计算离散非周期信号频谱、连续周期信号和连续非周期信号的频谱。信号与系统(陈后金)_MATLAB例4：计算下图所示周期矩形序列的频谱xkk0M-M%Program 4_4 计算离散周期矩形序列的频谱N=32; M=4; %定义周期矩形序列的参数x=ones(1,M+1) zeros(1,N-2*M-1) ones(1,M); %产生周期矩形序列X=fft(x); %计算DFS系数m=0:N-1; stem(m,real(X); %画出频谱X的实部title(Xm的实部);xlabel(m)figurestem(m,imag(X); %画出频谱X的虚部title(Xm的虚部);xlabel(m)信号与系统(陈后金)_MATLAB例4：计算下图所示周期矩形序列的频谱信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB) 1()j)(2()j)(1 () 1()j)(2()j)(1 ()()()j (11+-NaaaMbbbABHNNMMwwwwwwwb 分子多项式系数 a 分母多项式系数 w 需计算的H(jw)的抽样点 (数组w中少需包含两个w的抽样点)。计算频响的MATLAB函数H=freqs(b,a,w)信号与系统(陈后金)_MATLAB例例1 三阶归一化的Butterworth低通滤波器的系统函数为w=linspace(0,5,200);b=1;a=1 2 2 1;h=freqs(b,a,w);subplot(2,1,1);plot(w,abs(h);subplot(2,1,2);plot(w,angle(h);1)j (2)j (2)j (1)j (23+wwwwH 试画出|H(jw)| 和(w)。信号与系统(陈后金)_MATLAB三阶Butterworth低通滤波器的幅度响应和相位响应信号与系统(陈后金)_MATLAB)(cosj11)2(Sa21001www+tnRCnnnCRCHwwwj/1j/1)j (+RCwj11+tnnnnHCty0j0)e(j)(ww-2Sa0wnTACn2Sa21n信号与系统(陈后金)_MATLAB%p5_2 Periodic signal pass LTI systemT=4;w0=2*pi/T;RC=0.1;t= -6:0.01:6;N=51;c0=0.5;xN=c0*ones(1,length(t); %dc for n=1:2:N % even harmonics are zero H=abs(1/(1+j*RC*w0*n); phi=angle(1/(1+j*RC*w0*n); xN=xN+H*cos(w0*n*t+phi)*sinc(n*0.5);endplot(t,xN);xlabel(time RC=,num2str(RC);grid;set(gca,xtick,-5 -3 -1 0 1 3 5);信号与系统(陈后金)_MATLAB-5-3-113500.51RC=1t(秒 )y(t)-5-3-113500.51RC=0.1y(t)t(秒 )信号与系统(陈后金)_MATLABNNMMaaabbbABHjj10jj10jjjeeee)e ()e ()e (-+计算频率响应的MATLAB函数b 分子的系数 a 分母系数w 抽样的频率点(至少2点)， w在02之间幅度响应: abs, 相位响应: angleh = freqz(b,a,w)信号与系统(陈后金)_MATLABb=1;a1=1 -0.9; a2=1 0.9;w=linspace(0,2*pi,512);h1=freqz(b,a1,w);h2=freqz(b,a2,w);plot(w/pi,abs(h1),w/pi,abs(h2),:);legend(alpha=0.9,alpha=-0.9);的幅度响应曲线画例：jje11)e (-H信号与系统(陈后金)_MATLAB的幅度响应曲线画例：jje11)e (-H信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLABr,p,k=residue(num,den) num,den分别为X(s)分子多项式和分母多项式的系数向量。 r为部分分式的系数，p为极点，k为多项式的系数。若为真分式，则k为零。信号与系统(陈后金)_MATLAB计算多项式根roots的函数可用于计算H(s)的零极点。r=roots(D) %计算多项式D(s)的根 H(s)零极点分布图可用pzmap函数画出，调用形式为pzmap(sys)表示画出sys所描述系统的零极点图。信号与系统(陈后金)_MATLAB用部分分式展开法求X(s)的反变换。%program6_1format rat %将结果数据以分数的形式输出num=1 2; den=1 4 3 0; r,p=residue(num,den) sssssX342)(23+r = -1/6 ,-1/2 ,2/3 p = -3 ,-1 ,036/115 . 03/2)(+-+-+ssssX)(e61)(e21)(32)()(31tututusXLtxtt-信号与系统(陈后金)_MATLAB用部分分式展开法求X(s)的反变换。% program6_2num=2 3 0 5;den=conv(1 1,1 1 2); %将因子相乘的形式转换成多项式的形式r,p,k=residue(num,den)magr=abs(r) %求r的模angr=angle(r) %求r的相角)2(1532)(223+ssssssX）（信号与系统(陈后金)_MATLAB用部分分式展开法求X(s)的反变换。)2(1532)(223+ssssssX）（lr =-2.0000 + 1.1339i, -2.0000 - 1.1339i, 3.0000lp =-0.5000 + 1.3229i, -0.5000 - 1.3229i, -1.0000lk =2lmagr =2.299, 2.2991, 3.0000 langr =2.6258, -2.6258, 0 3229. 1 j5 . 0e2991. 23229. 1 j5 . 0e2991. 2132)(6258. 2j6258. 2j-+-ssssX)()6258. 23229. 1cos(e1495. 1)(e3)(2)(5 . 0tuttuttxtt+-信号与系统(陈后金)_MATLAB例3 试画出系统 的零极点分布图，求其单位冲激响应h(t)和频率响应H(jw)，并判断系统是否稳定。num=1;den=1 2 2 1;sys=tf(num,den); poles=roots(den)figure(1);pzmap(sys);t=0:0.02:10;h=impulse(num,den,t);figure(2);plot(t,h)title(Impulse Respone)H,w=freqs(num,den);figure(3);plot(w,abs(H)xlabel(omega)title(Magnitude Respone)1221)(23+ssssH信号与系统(陈后金)_MATLAB运行结果运行结果-1.4-1.2-1-0.8-0.6-0.4-0.20-1-0.8-0.6-0.4-0.200.20.40.60.81Real AxisImag AxisPole-zero map012345678910-0.0500.050.10.150.20.250.30.350.40.45Impulse Responetime (s)01234567891000.10.20.30.40.50.60.70.80.91FrequencywMagnitude Respone信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLABr,p,k=residuez(num,den) num,den分别为X(z)分子多项式和分母多项式的系数向量。 r为部分分式的系数，p为极点，k为多项式的系数。若为真分式，则k为零。信号与系统(陈后金)_MATLAB利用tf2zp函数计算H(z)的零极点，调用形式为z,p,k=tf2zp(b,a) H(z)零极点分布图可用zplane函数画出，调用形式为zplane(b,a)b和a分别为H(z)分子多项式和分母多项式的系数向量。返回值z为零点、p为极点、 k为增益常数。信号与系统(陈后金)_MATLAB例例1 1 将X(z)用部分分式展开。%program7_1num = 18; den = 18 3 -4 -1;r,p,k = residuez(num,den)运行结果为r =0.3600 , 0.2400 , 0.4000 p =0.5000 , -0.3333 , -0.3333 k =故X(z)可展开为321431818)(-+zzzzX2111)3333. 01 (4 . 03333. 0124. 05 . 0136. 0)(-+-zzzzX信号与系统(陈后金)_MATLAB例2 试画出系统 的零极点分布图，求其单位冲激响应hk和频率响应H(ej) 。% program 7_2b =1 2 1;a =1 -0.5 -0.005 0.3;figure(1);zplane(b,a);num=0 1 2 1;den=1 -0.5 -0.005 0.3;h=impz(num,den);figure(2);stem(h)xlabel(k)title(Impulse Respone)H,w=freqz(num,den);figure(3);plot(w/pi,abs(H)xlabel(Frequency omega)title(Magnitude Respone)3213213 . 0005. 05 . 012)(-+-+zzzzzzzH信号与系统(陈后金)_MATLAB运行结果运行结果2Real part-1-0.500.51Imaginary part-1-0.500.5105101520253035-1-0.500.511.522.5kImpulse Respone00.10.20.30.40.50.60.70.80.91012345678FrequencyMagnitude Respone信号与系统(陈后金)_MATLAB信号与系统(陈后金)_MATLAB num,den 分别表示系统函数H(s)的分子和分母多项式； A,B,C,D 分别为状态方程的矩阵。 A,B,C,D= tf2ss(num,den)信号与系统(陈后金)_MATLABA, B, C, D 分别表示状态方程的矩阵。K 表示函数ss2tf计算的与第k个输入相关的系统函数，即H(s)的第k列。num 表示H(s)第k列的m个元素的分子多项式den 表示H(s)公共的分母多项式。num,den= ss2tf (A,B,C,D,k)信号与系统(陈后金)_MATLABsys 由函数ss构造的状态方程模型t 需计算的输出样本点, t=0:dt:Tfinalx(:,k) 系统第k个输入在t上的抽样值q0 系统的初始状态(可缺省)y(:,k) 系统的第k个输出 to 实际计算时所用的样本点； q 系统的状态 获得连续系统状态方程的计算机表示模型sys = ss(A,B,C,D)求解状态方程y,to,q=lsim(sys,x,t,q0)信号与系统(陈后金)_MATLABsys 由函数ss构造的状态方程模型x(:,k)系统第k个输入序列q0 系统的初始状态(可缺省)y(:,k)系统的第k个输出n 序列的下标; q 系统的状态获得离散状态方程的计算机表示模型sys = ss(A,B,C,D, )求解状态方程y,n,q=lsim(sys,x, ,q0)或直接利用y, q=dlsim(A,B,C,D,x,q0)信号与系统(陈后金)_MATLAB例1 写出系统 的状态方程。由A,B,C,D=tf2ss(1,1 5 10)可得)()(10)( 5)( txtytyty+-01105A01B所以系统的状态方程为)(01)()(01105)()(2121txtqtqtqtq+-)()(10)(21tqtqty信号与系统(陈后金)_MATLAB例2 已知某连续系统的状态方程和输出方程为+-)()(0110)()(1032)()(212121txtxtqtqtqtq+-)()(0101)()(1011)()(212121txtxtqtqtyty其初始状态和输入分别为-12)0()0(21qq-)(e)()()(321tututxtxt求该系统的系统函数矩阵H(s) 和输出。信号与系统(陈后金)_MATLAB计算系统函数矩阵H(s)A=2 3;0 -1;B=0 1; 1 0;C=1 1; 0 -1;D=1 0; 1 0;num1,den1=ss2tf(A,B,C,D,1)num2,den2=ss2tf(A,B,C,D,2)运行结果num1 =1 0 -1 1 -2 0den1 =1 -1 -2num2 =0 1 1 0 0 0den2 =1 -1 -2所以系统函数矩阵H(s)为-+-021121)(222sssssssH信号与系统(陈后金)_MATLAB%Program 8_2A=2 3;0 -1;B=0 1; 1 0;C=1 1; 0 -1;D=1 0; 1 0;q0=2 -1;dt=0.01;t=0:dt:2; x(:,1)=ones(length(t),1);x(:,2)=exp(-3*t);sys=ss(A,B,C,D);y=lsim(sys,x,t,q0);subplot(2,1,1);plot(t,y(:,1),r);ylabel(y1(t);xlabel(t);subplot(2,1,2);plot(t,y(:,2);ylabel(y2(t);xlabel(t);信号与系统(陈后金)_MATLAB00.20.40.60.811.21.41.61.82020406080100y1(t)t00.20.40.60.811.21.41.61.8200.511.52y2(t)t信号与系统(陈后金)_MATLAB例3 已知某离散系统的状态方程和输出方程为其初始状态和输入分别为求该系统的输出。+-+103210 1 1212121kxkxkqkqkqkq-12112121kqkqkyky-110021qqkukx信号与系统(陈后金)_MATLAB%Program 8_3A=0 1; -2 3;B=0;1;C=1 1; 2 -1;D=zeros(2,1);q0=1; -1;N=10;x=ones(1,N);y=dlsim(A,B,C,D,x,q0);subplot(2,1,1);y1=y(:,1);stem(0:N-1),y1);xlabel(k);ylabel(y1);subplot(2,1,2);y2=y(:,2);stem(0:N-1),y2);xlabel(k);ylabel(y2);信号与系统(陈后金)_MATLAB0123456789-2000-1500-1000-5000ky10123456789-6-4-2024ky2