信号与系统实验6拉普拉斯变换

一、目的

（1）掌握连续系统及信号拉普拉斯变换概念

（2）掌握利用MATLAB绘制系统零极点图的方法 （3）掌握利用MATLAB求解拉普拉斯逆变换的方法

二、拉普拉斯变换曲面图的绘制

连续时间信号f(t)的拉普拉斯变换定义为:

F(s)0f(t)estdt （6-1）

其中sj，若以为横坐标（实轴），j为纵坐标（虚轴），复变量s就构成了一个复平面，称为s平面。

显然，F(s)是复变量s的复函数，为了便于理解和分析F(s)随s的变化规律，可以将F(s)写成：

F(s)F(s)ej(s) （6-2）

其中，F(s)称为复信号F(s)的模，而(s)则为F(s)的幅角。

从三维几何空间的角度来看，F(s)和(s)对应着复平面上的两个平面，如果能绘出它们的三维曲面图，就可以直观地分析连续信号的拉普拉斯变换F(s)随复变量s的变化规律。

上述过程可以利用MATLAB的三维绘图功能实现。现在考虑如何利用MATLAB来绘制s平面的有限区域上连续信号f(t)的拉普拉斯变换F(s)的曲面图，现以简单的阶跃信号u(t)为例说明实现过程。

1我们知道，对于阶跃信号f(t)u(t)，其拉普拉斯变换为F(s)。首先，利用两

s个向量来确定绘制曲面图的s平面的横、纵坐标的范围。例如可定义绘制曲面图的横坐标范围向量x1和纵坐标范围向量y1分别为：

x1=-0.2:0.03:0.2; y1=-0.2:0.03:0.2;

然后再调用meshgrid()函数产生矩阵s，并用该矩阵来表示绘制曲面图的复平面区域，对应的MATLAB命令如下：

[x,y]=meshgrid(x1,y1); s=x+i*y;

上述命令产生的矩阵s包含了复平面0.20.2， 0.2j0.2范围内以时间间隔0.03取样的所有样点。

最后再计算出信号拉普拉斯变换在复平面的这些样点上的值，即可用函数mesh()绘出其曲面图，对应命令为：

fs=abs(1./s); mesh(x,y,fs); surf(x,y,fs);

title('单位阶跃信号拉氏变换曲面图'); colormap(hsv);

axis([-0.2,0.2,-0.2,0.2,0.2,60]); rotate3d;

执行上述命令后，绘制的单位阶跃信号拉普拉斯变换曲面图如图6-1所示。

45

图6-1 阶跃信号拉普拉斯变换曲面图

例6-1：已知连续时间信号f(t)sin(t)u(t)，求出该信号的拉普拉斯变换，并利用MATLAB绘制拉普拉斯变换的曲面图。 解：该信号的拉普拉斯变换为：

1 F(s)2s1利用上面介绍的方法来绘制单边正弦信号拉普拉斯变换的曲面图，实现过程如下： %绘制单边正弦信号拉普拉斯变换曲面图程序 图6-2 单边正弦信号拉氏变换曲面图 clf;

a=-0.5:0.08:0.5; b=-1.99:0.08:1.99;

46

[a,b]=meshgrid(a,b); d=ones(size(a)); c=a+i*b; c=c.*c; c=c+d; c=1./c;

c=abs(c); mesh(a,b,c); surf(a,b,c);

axis([-0.5,0.5,-2,2,0,15]);

title('单边正弦信号拉氏变换曲面图'); colormap(hsv);

上述程序运行结果如图6-2所示。

%确定绘制曲面图的复平面区域

%计算拉普拉斯变换的样值 %绘制曲面图

二、由拉普拉斯曲面图观察频域与复频域的关系

如果信号f(t)的拉普拉斯变换F(s)的极点均位于s平面左半平面，则信号f(t)的傅立叶变换F(j)与F(s)存在如下关系：

F(j)F(s)sj （6-3） 即在信号的拉普拉斯变换F(s)中令0，就可得到信号的傅立叶变换。从三维几何空间角度来看，信号f(t)的傅立叶变换F(j)就是其拉普拉斯变换曲面图中虚轴所对应的曲线。可以通过将F(s)曲面图在虚轴上进行剖面来直观的观察信号拉普拉斯变换与其傅立叶变换的对应关系。

例6-2：试利用MATLAB绘制信号f(t)etsin(t)u(t)的拉普拉斯变换的曲面图，观察曲面图在虚轴剖面上的曲线，并将其与信号傅立叶变换F(j)绘制的幅度频谱相比较。 解：根据拉普拉斯变换和傅立叶变换定义和性质，可求得该信号的拉普拉斯变换和傅立叶变换如下：

11F(s)F(j) 22(s1)1(j1)1利用前面介绍的方法绘制拉普拉斯变换曲面图。为了更好地观察曲面图在虚轴剖面上的曲线，定义绘制曲面图的S平面实轴范围从0开始，并用view函数来调整观察视角。实现命令如下：

clf;

a=-0:0.1:5;

b=-20:0.1:20; [a,b]=meshgrid(a,b); c=a+i*b; %确定绘图区域 c=1./((c+1).*(c+1)+1); c=abs(c); %计算拉普拉斯变换 mesh(a,b,c); %绘制曲面图 surf(a,b,c);

view(-60,20) %调整观察视角 axis([-0,5,-20,20,0,0.5]);

title('拉普拉斯变换（S域像函数）'); colormap(hsv);

上述程序绘制的拉普拉斯变换的曲面如图6-3所示。从该曲面图可以明显地观察到F(s)在虚轴剖面上曲线变化情况。

47

图6-3 指数衰减正弦信号拉氏变换曲面图

图6-4 指数衰减正弦信号傅氏变换曲幅频图

利用MATLAB绘制该信号的傅立叶变换幅频曲线命令如下： w=-20:0.1:20; %确定频率范围 Fw=1./((i*w+1).*(i*w+1)+1); %计算傅里叶变换 plot(w,abs(Fw)) %绘制信号振幅频谱曲线 title('傅里叶变换（振幅频谱曲线）') xlabel('频率w')

运行结果如图6-4所示。通过图6-3和图6-4对比可直观地观察到拉普拉斯变换与傅立叶变换的对应关系。

48

三、拉普拉斯变换零极点分布对曲面图的影响

从单位阶跃信号和单边正弦信号的拉普拉斯变换曲面图可以看出，曲面图中均有突出的尖峰，仔细观察便可得出，这些峰值点在S平面的对应点就是信号拉普拉斯变换的极点位置。

我们再来看拉普拉斯变换零极点对曲面图的影响，考虑如下信号：

2(s3)(s3)F(s)

(s5)(s210)该信号的零点为z1,23，极点为p1,2j3.1623，p35。利用如下MATLAB命令绘制出的曲面图如图6-5所示。

clf;

a=-6:0.48:6; b=-6:0.48:6;

[a,b]=meshgrid(a,b); c=a+i*b;

d=2*(c-3).*(c+3); e=(c.*c+10).*(c-5); c=d./e; c=abs(c); mesh(a,b,c); surf(a,b,c);

axis([-6,6,-6,6,0,4.5]);

title('拉普拉斯变换曲面图'); colormap(hsv); view(-25,30)

图6-5 拉氏变换零极点分布曲面图

从图6-5可明显看出，曲面在sj3.1623和s5处有三个峰点，对应着拉普拉斯变换的极点位置，而在s3处有两个谷点，对应着拉普拉斯变换的零点位置。因此，信号的拉普拉斯变换的零极点位置，决定了其拉氏变换曲面图的峰点和谷点位置。

四、连续系统零极点图的绘制

49

线性时不变系统可用如下所示的线性常系数微分方程来描述：

ay(t)bf(t) （6-4）

(i)(j)i0ij0jNM其中，y(t)为系统输出信号，f(t)为输入信号。

将上式两边进行拉普拉斯变换，则该系统的系统函数为：

H(s)将式（6-5）因式分解后有：

Y(s)F(s)bsasi0ij0NjMjiB(s) （6-5） A(s)H(s)C(sz)(sp)i0ij0NjM （6-6）

其中C为常数zj为系统的零点，pi为系统的极点。

可见，若连续系统函数的零极点已知，系统函数便可确定下来。即系统函数H(s)的零极点分布完全决定了系统的特性。

因此，在连续系统的分析中，系统函数的零极点分布具有非常重要的意义。通过对系统函数零极点的分析，我们可以分析连续系统以下几方面的特性：

 系统冲激响应h(t)的时域特性；  判断系统的稳定性；

 分析系统的频率特性H(j)（幅频响应和相频响应）。

通过系统函数零极点分布来分析系统特性，首先就要求出系统函数的零极点，然后绘制系统零极点图。下面介绍如何利用MATLAB实现这一过程。

设连续系统的系统函数为：

B(s)H(s)

A(s)则系统函数的零极点位置可用MATLAB的多项式求根函数roots()来求得，调用函数roots()的命令格式为：

p=roots(A)

其中A为待求根的关于s的多项式的系数构成的行向量，返回向量p则是包含该多项式所有根位置的列向量。如多项式为：

A(s)s23s4

则求该多项式根的MATLAB命令为： A=[1 3 4]; p=roots(A) 运行结果为： p =

-1.5000 + 1.3229i -1.5000 - 1.3229i 需要注意的是，系数向量A的元素一定要由多项式最高次幂开始直到常熟项，缺项要用0补齐。如多项式为：

A(s)s63s42s2s4

则表示该多项式的系数向量为：

A=[1 0 3 0 2 1 -4];

50

用roots()函数求得系统函数H(s)的零极点后，就可以绘制零极点图，下面是求连续系统的系统函数零极点，并绘制其零极点图的MATLAB实用函数sjdt()。

function [p,q]=sjdt(A,B)

%绘制连续系统零极点图程序 %A:系统函数分母多项式系数向量 %B:系统函数分子多项式系数向量

%p:函数返回的系统函数极点位置行向量 %q:函数返回的系统函数零点位置行向量

p=roots(A); %求系统极点 q=roots(B);; %求系统零点 p=p'; %将极点列向量转置为行向量 q=q'; %将零点列向量转置为行向量 x=max(abs([p q])); %确定纵坐标范围 x=x+0.1; y=x; %确定横坐标范围 clf hold on

axis([-x x -y y]); %确定坐标轴显示范围 axis('square')

plot([-x x],[0 0]) %画横坐标轴 plot([0 0],[-y y]) %画纵坐标轴 plot(real(p),imag(p),'x') %画极点 plot(real(q),imag(q),'o') %画零点 title('连续系统零极点图') %标注标题 text(0.2,x-0.2,'虚轴') text(y-0.2,0.2,'实轴')

例6-3：已知连续系统的系统函数如下，试用MATLAB绘出系统零极点图。

s24 F(s)432s2s3s2s1解：MATLAB处理命令如下： a=[1 2 -3 2 1];

b=[1 0 -4]; sjdt(a,b)

运行结果如图6-6所示。

五、拉普拉斯逆变换及MATLAB实现

连续信号f(t)的拉普拉斯变换具有如下一般形式：

F(s)C(s)D(s)csdsi1ij1LjKj

i若KL，则F(s)可以分解为有理多项式与真分式之和，即

F(s)P(s)R(s)P(s)B(s)P(s)A(s)bsasi1ij1NjMj

i其中，P(s)是关于s的多项式，其逆变换可直接求得（冲激信号及其各阶导数），R(s)为

51

关于s的有理真分式，即满足MN。以下进讨论MN的情况。

设连续信号f(t)的拉普拉斯变换为F(s)，则

B(s)B(s)F(s)N

A(s)(sp)i1i在满足MN情况下，有以下几种情况：

图6-6 系统零极点图

（1）极点均为单重情况下，可对其直接进行部分分式展开得：

rrrF(s)12N

sp1sp2spN其中，ri(spi)F(s)spi(i1,2,,N)称为有理函数F(s)的留数。则F(s)的拉普拉斯逆变换为：

f(t)riepitu(t)

i1N（2）有k重极点，设为p1，则部分分式展开为 K11K12K1kE(s)F(s)

(sp1)k(sp1)k1(sp1)D(s)K1i可用下式求得

1di1(sp1)kF(s) K1ii1(i1)!dssp1则F(s)的拉普拉斯逆变换为：

NK1jkjpitf(t)teu(t)riepitu(t)

j1(kj)!i2k（3）有共轭极点

52

F(s)r1rrr22N sp1sp2sp3spNf2(t)设F(s)有一对共轭极点p1,2j，则

r1(sp1)F(s)sp1r1ej

r2r1

由共轭极点所决定的两项复指数信号可以合并成一项，故有

f2(t)2r1etcos(t)u(t)

*从以上分析可以看出，只要求出F(s)部分分式展开的系数（留数）ri，就可直接求出F(s)的逆变换f(t)。

上述求解过程，可以利用MATLAB的residue()函数来实现。令A和B分别为F(s)的分子和分母多项式构成的系数向量，则函数：

[r,p,k]=residue(B,A)

将产生三个向量r、p和k，其中p为包含F(s)所有极点的列向量，r为包含F(s)部分分式展开系数ri的列向量，k为包含F(s)部分分式展开的多项式的系数行向量，若

MN，则k为空。

例6-4：已知连续信号的拉普拉斯变换为：

F(s)2s4 3s4s试用MATLAB求其拉普拉斯逆变换f(t)。 解：MATLAB命令如下：

a=[1 0 4 0]; b=[2 4];

[r,p,k]=residue(b,a) 运行结果： r =

-0.5000 - 0.5000i -0.5000 + 0.5000i 1.0000 p =

0 + 2.0000i 0 - 2.0000i 0 k = []

由上述结果可以看出，F(s)有三个极点p1,2j2，p30，为了求得共轭极点对应的信号分量，可用abs()和angle()分别求出部分分式展开系数的模和幅角，命令如下：

abs(r) ans =

0.7071 0.7071 1.0000 angle(r)/pi ans =

-0.7500

53

0.7500 0

3由上述结果可得f(t)[12cos(2t)]u(t)。

4例6-5：求下式函数的逆变换

s2F(s) 3s(s1)解：MATLAB程序如下：

a=[1 3 3 1 0]; b=[1 -2];

[r,p,k]=residue(b,a) 运行结果： r =

2.0000 2.0000 3.0000 -2.0000 p =

-1.0000 -1.0000 -1.0000 0 k = []

223232t则F(s)，对应的逆变换为 f(t)[(t2t2)e2]u(t)。23(s1)(s1)(s1)s2六、实验内容

1、求解下述信号的拉普拉斯变换，并利用MATLAB绘制拉普拉斯变换的曲面图：

（1）f(t)cos(2t)u(t)

求解：

a=-0.5:0.08:0.5; b=-3.99:0.08:3.99; [a,b]=meshgrid(a,b); c=a+i*b; syms t;

f = cos(2*t)*heaviside(t); l = laplace(f); x=subs(l,'s',c); y=abs(x);

mesh(a,b,y); surf(a,b,y);

axis([-0.5,0.5,-4,4,0,15]); title('拉氏变换曲面图'); colormap(hsv);

54

（2）f(t)e2tsin(t)u(t)

a=-4:0.08:4; b=-2:0.08:0;

[a,b]=meshgrid(a,b); c=a+i*b; syms t;

f = exp(-2*t)*sin(t)*heaviside(t); l = laplace(f); x=subs(l,'s',c); y=abs(x);

mesh(a,b,y); surf(a,b,y);

axis([-4,4,-2,0,0,15]); title('拉氏变换曲面图'); colormap(hsv);

2、已知连续时间信号f(t)e3tu(t)，试求出该信号的拉普拉斯变换F(s)和傅立叶变换

F(j)，用MATLAB绘出拉普拉斯变换曲面图F(s)及幅频曲线F(j)，观察曲面图

在虚轴上的剖面图，并将其与幅频曲线相比较，分析频域与复频域的对应关系。

求傅里叶变换和拉普拉斯变换：

55

syms t;

f = exp(-3*t)*heaviside(t); F=fourier(f) L=laplace(f) F =

1/(3+i*w) L =

1/(s+3)

绘制拉普拉斯变换曲面图： a=-3:0.04:-1; b=-2:0.04:2;

[a,b]=meshgrid(a,b); c=a+i*b; c=1./(c+3); c=abs(c); mesh(a,b,c); surf(a,b,c); view(-60,20)

axis([-3,-1,-2,2,0,20]);

title('拉普拉斯变换（S域像函数）'); colormap(hsv);

绘制傅里叶变换振幅频谱曲线: w=-20:0.01:20; Fw=1./(i*w+3); plot(w,abs(Fw))

title('傅里叶变换（振幅频谱曲线）') xlabel('频率w')

56

3、已知信号的拉普拉斯变换如下所示，试用MATLAB绘制曲面图，观察拉普拉斯变换零极点分布对曲面图的影响。

function [p,q]=sjdt(A,B)

%绘制连续系统零极点图程序 %A:系统函数分母多项式系数向量 %B:系统函数分子多项式系数向量

%p:函数返回的系统函数极点位置行向量 %q:函数返回的系统函数零点位置行向量

p=roots(A); %求系统极点 q=roots(B); %求系统零点 p=p'; %将极点列向量转置为行向量 q=q'; %将零点列向量转置为行向量 x=max(abs([p q])); %确定纵坐标范围 x=x+0.1; y=x; %确定横坐标范围 clf hold on

axis([-x x -y y]); %确定坐标轴显示范围 axis('square')

plot([-x x],[0 0]) %画横坐标轴 plot([0 0],[-y y]) %画纵坐标轴 plot(real(p),imag(p),'x') %画极点 plot(real(q),imag(q),'o') %画零点 title('连续系统零极点图') %标注标题 text(0.2,x-0.2,'虚轴') text(y-0.2,0.2,'实轴')

（1）F(s)求零极点： A = [1 5 6 0]; B = [1 5 4]; sjdt(A,B);

(s1)(s4)

s(s2)(s3) 57

绘制曲面图： a=-6:0.04:2; b=-2:0.04:2;

[a,b]=meshgrid(a,b); c=a+i*b;

x = (c+1).*(c+4)./(c.*(c+2).*(c+3)); y=abs(x); mesh(a,b,y);

axis([-6,2,-2,2,0,20]); title('拉氏变换曲面图'); colormap(hsv);

（2）F(s)求零极点： A=[1 0 4]; B=[1 0 -4]; sjdt(A,B);

s4 2s42

58

绘制曲面图： a=-4:0.04:4; b=-4:0.04:4;

[a,b]=meshgrid(a,b); c=a+i*b;

x = (c.^2-4)./(c.^2+4); y=abs(x); mesh(a,b,y);

axis([-4,4,-4,4,0,20]); title('拉氏变换曲面图'); colormap(hsv);

4、试用MATLAB求下列信号的拉普拉斯逆变换

s25s4 （1）F(s)3 2s5s6ssyms s;

F = (s^2+5*s+4)/(s^3+5*s^2+6*s); f=ilaplace(F)

59

f =

2/3+exp(-2*t)-2/3*exp(-3*t)

（2）F(s)1

s32s22s1syms s;

F = 1/(s^3+2*s^2+2*s+1); f=ilaplace(F) f =

exp(-t)+1/3*exp(-1/2*t)*(-3*cos(1/2*3^(1/2)*t)+3^(1/2)*sin(1/2*3^(1/2)*t))

60