matlab_拉普拉斯变换.doc

实验六 拉普拉斯变换及其逆变换

一、目的
（1）掌握连续系统及信号拉普拉斯变换概念
（2）掌握利用MATLAB绘制系统零极点图的方法
（3）掌握利用MATLAB求解拉普拉斯逆变换的方法
二、拉普拉斯变换曲面图的绘制
连续时间信号的拉普拉斯变换定义为:
（6-1）
其中，若以为横坐标（实轴），为纵坐标（虚轴），复变量就构成了一个复平面，称为平面。
显然，是复变量的复函数，为了便于理解和分析随的变化规律，可以将写成：
（6-2）
其中，称为复信号的模，而则为的幅角。
从三维几何空间的角度来看，和对应着复平面上的两个平面，如果能绘出它们的三维曲面图，就可以直观地分析连续信号的拉普拉斯变换随复变量的变化规律。
上述过程可以利用MATLAB的三维绘图功能实现。现在考虑如何利用MATLAB来绘制平面的有限区域上连续信号的拉普拉斯变换的曲面图，现以简单的阶跃信号为例说明实现过程。
我们知道，对于阶跃信号，其拉普拉斯变换为。首先，利用两个向量来确定绘制曲面图的平面的横、纵坐标的范围。例如可定义绘制曲面图的横坐标范围向量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;
上述命令产生的矩阵包含了复平面， 范围内以时间间隔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所示。
图6-1 阶跃信号拉普拉斯变换曲面图

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

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

clf;
a=-0.5:0.08:0.5;
b=-1.99:0.08:1.99;
[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所示。
二、由拉普拉斯曲面图观察频域与复频域的关系
如果信号的拉普拉斯变换的极点均位于平面左半平面，则信号的傅立叶变换与存在如下关系：
（6-3）
即在信号的拉普拉斯变换中令，就可得到信号的傅立叶变换。从三维几何空间角度来看，信号的傅立叶变换就是其拉普拉斯变换曲面图中虚轴所对应的曲线。可以通过将曲面图在虚轴上进行剖面来直观的观察信号拉普拉斯变换与其傅立叶变换的对应关系。
例6-2：试利用MATLAB绘制信号的拉普拉斯变换的曲面图，观察曲面图在虚轴剖面上的曲线，并将其与信号傅立叶变换绘制的幅度频谱相比较。
解：根据拉普拉斯变换和傅立叶变换定义和性质，可求得该信号的拉普拉斯变换和傅立叶变换如下：

利用前面介绍的方法绘制拉普拉斯变换曲面图。为了更好地观察曲面图在虚轴剖面上的曲线，定义绘制曲面图的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所示。从该曲面图可以明显地观察到在虚轴剖面上曲线变化情况。

图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对比可直观地观察到拉普拉斯变换与傅立叶变换的对应关系。
三、拉普拉斯变换零极点分布对曲面图的影响
从单位阶跃信号和单边正弦信号的拉普拉斯变换曲面图可以看出，曲面图中均有突出的尖峰，仔细观察便可得出，这些峰值点在S平面的对应点就是信号拉普拉斯变换的极点位置。
我们再来看拉普拉斯变换零极点对曲面图的影响，考虑如下信号：

该信号的零点为，极点为，。利用如下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可明显看出，曲面在和处有三个峰点，对应着拉普拉斯变换的极点位置，而在处有两个谷点，对应着拉普拉斯变换的零点位置。因此，信号的拉普拉斯变换的零极点位置，决定了其拉氏变换曲面图的峰点和谷点位置。
四、连续系统零极点图的绘制
线性时不变系统可用如下所示的线性常系数微分方程来描述：
（6-4）
其中，为系统输出信号，为输入信号。
将上式两边进行拉普拉斯变换，则该系统的系统函数为：
（6-5）
将式（6-5）因式分解后有：
（6-6）
其中为常数为系统的零点，为系统的极点。
可见，若连续系统函数的零极点已知，系统函数便可确定下来。即系统函数的零极点分布完全决定了系统的特性。
因此，在连续系统的分析中，系统函数的零极点分布具有非常重要的意义。通过对系统函数零极点的分析，我们可以分析连续系统以下几方面的特性：
l 系统冲激响应的时域特性；
l 判断系统的稳定性；
l 分析系统的频率特性（幅频响应和相频响应）。
通过系统函数零极点分布来分析系统特性，首先就要求出系统函数的零极点，然后绘制系统零极点图。下面介绍如何利用MATLAB实现这一过程。
设连续系统的系统函数为：

则系统函数的零极点位置可用MATLAB的多项式求根函数roots()来求得，调用函数roots()的命令格式为：
p=roots(A)
其中A为待求根的关于s的多项式的系数构成的行向量，返回向量p则是包含该多项式所有根位置的列向量。如多项式为：

则求该多项式根的MATLAB命令为：
A=[1 3 4];
p=roots(A)
运行结果为：
p =
-1.5000 + 1.3229i
-1.5000 - 1.3229i
需要注意的是，系数向量A的元素一定要由多项式最高次幂开始直到常熟项，缺项要用0补齐。如多项式为：

则表示该多项式的系数向量为：
A=[1 0 3 0 2 1 -4];
用roots()函数求得系统函数的零极点后，就可以绘制零极点图，下面是求连续系统的系统函数零极点，并绘制其零极点图的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绘出系统零极点图。

解：MATLAB处理命令如下：
a=[1 2 -3 2 1];
b=[1 0 -4];
sjdt(a,b)
运行结果如图6-6所示。
五、拉普拉斯逆变换及MATLAB实现
连续信号的拉普拉斯变换具有如下一般形式：

若，则可以分解为有理多项式与真分式之和，即

其中，是关于的多项式，其逆变换可直接求得（冲激信号及其各阶导数），为关于的有理真分式，即满足。以下进讨论的情况。
设连续信号的拉普拉斯变换为，则

图6-6 系统零极点图
在满足情况下，有以下几种情况：

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

其中，称为有理函数的留数。则的拉普拉斯逆变换为：

（2）有重极点，设为，则部分分式展开为

可用下式求得

则的拉普拉斯逆变换为：

（3）有共轭极点

设有一对共轭极点，则

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

从以上分析可以看出，只要求出部分分式展开的系数（留数），就可直接求出的逆变换。
上述求解过程，可以利用MATLAB的residue()函数来实现。令A和B分别为的分子和分母多项式构成的系数向量，则函数：
[r,p,k]=residue(B,A)
将产生三个向量r、p和k，其中p为包含所有极点的列向量，r为包含部分分式展开系数的列向量，k为包含部分分式展开的多项式的系数行向量，若，则k为空。
例6-4：已知连续信号的拉普拉斯变换为：

试用MATLAB求其拉普拉斯逆变换。
解：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 =
[]
由上述结果可以看出，有三个极点，，为了求得共轭极点对应的信号分量，可用abs()和angle()分别求出部分分式展开系数的模和幅角，命令如下：
abs(r)
ans =
0.7071
0.7071
1.0000
angle(r)/pi
ans =
-0.7500
0.7500
0
由上述结果可得。
例6-5：求下式函数的逆变换

解：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 =
[]
则，对应的逆变换为。
六、实验内容
1、求解下述信号的拉普拉斯变换，并利用MATLAB绘制拉普拉斯变换的曲面图：
（1）
（2）
2、已知连续时间信号，试求出该信号的拉普拉斯变换和傅立叶变换，用MATLAB绘出拉普拉斯变换曲面图及幅频曲线，观察曲面图在虚轴上的剖面图，并将其与幅频曲线相比较，分析频域与复频域的对应关系。
3、已知信号的拉普拉斯变换如下所示，试用MATLAB绘制曲面图，观察拉普拉斯变换零极点分布对曲面图的影响。
（1）
（2）
4、试用MATLAB求下列信号的拉普拉斯逆变换
（1）
（2）
54