2.4 连续时间系统的时域分析
MATLAB信号处理超级学习手册
2.4.1 连续时间系统求解
连续时间线性非时变系统（LTI）可以用如下的线性常系数微分方程来描述：
系统的响应一般包括两个部分，即由当前输入所产生的响应（零状态响应）和由历史输入（初始状态）所产生的响应（零输入响应）。
对于低阶系统，一般可以通过解析的方法得到响应。但是，对于高阶系统，手工计算就比较困难，这时MATLAB强大的计算功能就比较容易确定系统的各种响应，如冲激响应、阶跃、零状态响应、全响应等。
涉及到的MATLAB函数有：impulse（冲激响应）、step（阶跃）、roots（零输入响应）、lsim（零状态响应）等。在MATLAB中，要求以系统向量的形式输入系统的微分方程，因此，在使用前必须对系统的微分方程进行变换，得到其传递函数。其分别用向量a和b表示分母多项式和分子多项式的系数（按照s的降幂排列）。
根据系统的单位冲激响应，利用卷积计算的方法，也可以计算任意输入状态下系统的零状态响应。设一个线性零状态系统，已知系统的单位冲激响应为h(t)，当系统的激励信号为f(t)时，系统的零状态响应为：
式中y zs(k)、f(k)、h(k)分别对应以T为时间间隔对连续时间信号y zs(t)、f(t)和h(t)进行采样得到的离散序列。
2.4.2 连续时间系统数值求解
在MATLAB中，控制系统工具箱提供了一个用于求解零初始条件微分方程数值解的函数lsim。其调用格式：
y=lsim(sys,f,t)
式中，t表示计算系统响应的抽样点向量，f是系统输入信号向量，sys是LTI系统模型，用来表示微分方程、差分方程或状态方程。其调用格式：
sys=tf(b,a)
式中，b和a分别是微分方程的右端和左端系数向量。例如，对于以下方程：
可用a = [a_3 ,a_2 ,a_1 ,a_0 ];b = [b_3 ,b_2 ,b_1 ,b_0 ]; sys = tf(b,a) 获得其LTI模型。注意，如果微分方程的左端或右端表达式中有缺项，则其向量a或b中的对应元素应为零，不能省略不写，否则出错。
【例2-23】有一物理学系统，用微分方程描述为 y^{''} (t) + 2y^' (t) + 100y(t) = 10sin 2pi t ，求系统的零状态响应。运行程序如下：
clear
ts=0;te=5;dt=0.01;
sys=tf([1],[1 2 100]);
t=ts:dt:te;
f=10*sin(2*pi*t);
y=lsim(sys,f,t);
plot(t,y);
xlabel('t(s)');ylabel('y(t)');
title('零状态响应')
grid on;
运行结果如图2-23所示。
在MATLAB中，求解系统冲激响应可应用控制系统工具箱个提供的函数impulse，求解阶跃响应可利用函数step，其调用形式为：
y=impulse(sys,t)
y=step(sys,t)
式中，t表示计算系统响应的抽样点向量，sys是LTI系统模型。
【例2-24】计算下述系统在冲激、阶跃、斜坡、正弦激励下的零状态响应：
运行程序如下：
b=[-0.475 -0.248 -0.1189 -0.0564];a=[1 0.6363 0.9396 0.5123 0.0037];
sys=tf(b,a);
T=1000;
t=0:1/T:10;t1=-5:1/T:5;
f1=stepfun(t1,-1/T)-stepfun(t1,1/T);
f2=stepfun(t1,0);
f3=t;
f4=sin(t);
y1=lsim(sys,f1,t);
y2=lsim(sys,f2,t);
y3=lsim(sys,f3,t);
y4=lsim(sys,f4,t);
subplot(221);
plot(t,y1);
xlabel('t');ylabel('y1(t)');
title('冲激激励下的零状态响应');
grid on;axis([0 10 -1.2 1.2]);
subplot(222);
plot(t,y2);
xlabel('t');ylabel('y2(t)');
title('阶跃激励下的零状态响应');
grid on;axis([0 10 -1.2 1.2]);
subplot(223);
plot(t,y3);
xlabel('t');ylabel('y3(t)');
title('斜坡激励下的零状态响应');
grid on;axis([0 10 -5 0.5]);
subplot(224);
plot(t,y4);
xlabel('t');ylabel('y4(t)');
title('正弦激励下的零状态响应');
grid on;axis([0 10 -1.5 1.2]);
运行程序如图2-24所示。
2.4.3 连续时间系统符号求解
连续时间系统可以使用常系数微分方程来描述，其完全响应由零输入响应和零状态响应组成。MATLAB符号工具箱提供了dsolve函数，可以实现对常系数微分方程的符号求解，其调用格式为：
dsolve('eq1,eq2…','cond1,cond2,…','v')
其中，参数eq表示各个微分方程，它与MATLAB符号表达式的输入基本相同，微分和导数的输入是使用Dy，D2y，D3y来表示y的一价导数，二阶导数，三阶导数；参数cond表示初始条件或者起始条件；参数v表示自变量，默认是变量t。
通过使用dsolve函数可以求出系统微分方程的零输入响应和零状态响应，进而求出完全响应。
【例2-25】已知某线性时不变系统的动态方程为：y''(t) + 4y'(t) + 4y(t) = 2f'(t) + 3f(t), t>0系统的初始状态为y(0) = 0,y'(0) = 1, 求系统的零输入响应y x(t)。
运行程序如下：
eq=’D2y+4*Dy+4*y=0’;
cond=’y(0)=0,Dy(0)=1’;
yx=dsolve(eq,cond);
yx=simplify(yx);
ezplot(yx,[0,10]);
xlabel('t');ylabel('yx(t)');
title('系统的零输入响应');
grid on;
运行结果如图2-25所示。
2.4.4 连续时间系统卷积求解
信号的卷积运算有符号算法和数值算法，此处采用数值计算法，需调用MATLAB的conv函数近似计算信号的卷积积分。连续信号的卷积积分定义是：
上式表明，连续信号f 1(t)和f 2(t)的卷积，可用各自抽样后的离散时间序列的卷积再乘以抽样间隔Δ。抽样间隔Delta 越小，误差越小。
【例2-26】用数值计算法求f_1 (t) = u(t) - u(t - 2) 与f_2 (t) = e^{ - 3t} u(t) 的卷积积分。
运行程序如下：
dt=0.01; t=-1:dt:2.5;
f1=heaviside(t)-heaviside(t-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'); ylabel('f1(t)');
subplot(222);
plot(t,f2);
grid on;
axis([-1,2.5,-0.2,1.2]);
title('f2(t)');
xlabel('t'); ylabel('f2(t)');
subplot(212);
plot(tt,f);
grid on;
title('f(t)=f1(t)*f2(t)');
xlabel('t'); ylabel('f3(t)');
运行结果如图2-26所示。