信号与系统实验教程(MATLAB)

信 号 与 系 统

实 验 教 程

目 录

实验一：连续时间信号与系统的时域分析—————-----—--—--—---——--—--————---—-——————————-——6

一、实验目的及要求—-——-————-————————-—--——-----——------

—--—-—-—-—--——-—--—-——--—---—---------6

二、实验原理—-—-———---——————-——-—-—-————-———--—-——

-—---——--—-—-—-——-—--———————-———-———-—---—-——6

1、信号的时域表示方法-—-————-————--—--—————-—--———-

——-———-——-——---—--—---—-——-—--—--—--6

2、用MATLAB仿真连续时间信号和离散时间信号-——-—--—--—----——-—-

——--—-——-—-———7

3、LTI系统的时域描述———-—----—--——----——--—-—----———----—--

—-——-—--——--—-———--——--—--11

三、实验步骤及内容——-——-——————-———--—-----—-——---—-———

—-———-——————-——---———-——-—--——--——-—-—15

四、实验报告要求—-——-—-————--——-———-———-—-----—-—-——-—

——-—--——--————-—————---—-——-———-——----—26

实验二:连续时间信号的频域分析————---——---—--—--——-———-————-----———————-—————-——--—————27

一、实验目的及要求--—--—-----———--—--—-—-——-—-—-————-—--—--

———----——-—--———-——-——-—-----———-27

二、实验原理—-—---——--—--—--——-—--———-——————----—-—-—---

————-——————-————--——-—--———--—--—-—--—27

1、连续时间周期信号的傅里叶级数CTFS—--——--—--—--—-——-——-———

—-—--—-——--—————---——27

2、连续时间信号的傅里叶变换CTFT————--—--———-—--—-—-—--—--—

—--———————-—-——————--——28

3、离散时间信号的傅里叶变换DTFT -—----————-—------—-—-------———

—------——-—-—-—---28

4、连续时间周期信号的傅里叶级数CTFS的MATLAB实现-—--——-————--

—-—-—-----—29

5、用MATLAB实现CTFT及其逆变换的计算——--——-——-—-————-—--

———---—-——---———--—33

三、实验步骤及内容--—--——-—----——----—---———--———--—-—-—--——-—--—-———-———-——--——-———-—-—34

四、实验报告要求--—--————---—-————-——————-—————---—--—---——-———--—--——-—----———-—-——-—-——48

实验三：连续时间LTI系统的频域分析--———-———-——-—-—---——--———————---——---—----—-—-—---49

一、实验目的及要求—----——-——--——————--—————--—-——-——

-----————-—-------——-——-—-———-—————----—49

二、实验原理—-—--—-—-—-—-—-—-——-———---————————----—-—

——-—---—---——---——-——-—--————---——--—--——49

1、连续时间LTI系统的频率响应————————--——-—-—---—-—-—--—

—-———-————-———---—-—---——-—-49

2、LTI系统的群延时—--—--—--——-—--——-—-—-——-—------——-—-—---

——-———------———--—-——-----—-50

3、用MATLAB计算系统的频率响应——-—————-—-————---—-—-----

———-—---—--—--———--——-—-50

三、实验步骤及内容-————-—---——-—-—-——-—----——-----—-—-—-----——-—-—-————-—---———————--——-51

四、实验报告要求——————--—-----—--——-——-——--———-——-—-

2

—-——--——-—-—------—--—————--—----—-——58

实验四：调制与解调以及抽样与重建---------—————————-——-————-——--—--—--—--———-—————-——--59

一、实验目的及要求—-—--——-—-——-—----—--—-———--—--——--—-—--—

———-——--———-—--—-————----—---—--—59

二、实验原理—-—--————-—-———-—-—--——-—-——-—-—-————-—

——----—-----—--—-———————-----—----——--—--——59

1、信号的抽样及抽样定理--———---—-—--——-—————-———-—---—---

—-—-———-———--—--——-----—--——-59

2、信号抽样过程中的频谱混叠———--————--——-————-——-—--——---

—-—-—---————--—-———---—-———-62

3、信号重建-————-—-——-—-—-————-- —-——--————--———

----------——---—---——------——--—-—-———-——--62

4、调制与解调--——---—-—----—----—-———---—----————-——---———-——--————----—————---———-----—----—-—

5、通信系统中的调制与解调仿真—-—---—-———-——-———---—---—-———--—

———-—--——--—--—--—-————-66

三、实验步骤及内容————-—————-———--—-——-—---—------——-------

——-———--—--——-—-———-—--—---—-—-66

四、实验报告要求--—---—--——-—---——--—--—-—-—--—-——--———-—-—

——--———--——————---—-——————-—----75

实验五：连续时间LTI系统的复频域分析----—-—--—-——----—-—-——-——————————-—--——----——76

一、实验目的及要求-———-——---—--—-————----—-——————---—--—-

———---———————-—-——-—-——-—————-——-76

二、实验原理--—-—-------—-—-—-—-—-------——----——-----—-—---—-—-—-

—--—--—--———-—---—-----—-——76

1、连续时间LTI系统的复频域描述—-—---——-——-—--—————-—----——

—----—----—-—-—-—---——76

2、系统函数的零极点分布图—--—-———-——-————--—-——-—---——-—---—--———---—-———-—-—--—-----——-——77

3、拉普拉斯变换与傅里叶变换之间的关系-—-—--—-—-——-—-—-—--——---—-—---

3

——--———----——-—-78

4、系统函数的零极点分布与系统稳定性和因果性之间的关系--—-——-——-—---——

—-—-——-—79

5、系统函数的零极点分布与系统的滤波特性—--———--——--—---—---—————

——-—--—--—---—————80

6、拉普拉斯逆变换的计算———-——-———-—-----———----——-———---

—————---——-——-—--———-———-—-—81

三、实验步骤及内容-—-------—-—--—---—--———--————--—-————-——

—---—-————--———--—-—---—-—---——82

四、实验报告要求——--—---———--—-—-----——--——-—--——------—-—--—

-------——--——--—-—--—-—-—-----87

附录：授课方式和考核办法—----—------——-—-—-----—-—--——-—-—-—--—-——----———-—-——-------—-——88

实验一 信号与系统的时域分析

一、实验目的

1、熟悉和掌握常用的用于信号与系统时域仿真分析的MATLAB函数；

2、掌握连续时间和离散时间信号的MATLAB产生，掌握用周期延拓的方法将一个非周期信号进行周期信号延拓形成一个周期信号的MATLAB编程；

3、牢固掌握系统的单位冲激响应的概念，掌握LTI系统的卷积表达式及其物理意义，掌握卷积的计算方法、卷积的基本性质；

4、掌握利用MATLAB计算卷积的编程方法，并利用所编写的MATLAB程序验证卷积的常用基本性质;

掌握MATLAB描述LTI系统的常用方法及有关函数，并学会利用MATLAB求解LTI系统响应，绘制相应曲线。

基本要求：掌握用MATLAB描述连续时间信号和离散时间信号的方法,能够编写MATLAB程序，实现各种信号的时域变换和运算，并且以图形的方式再现各种信号的波形.掌握线性时不变连续系统的时域数学模型用MATLAB描述的方法，掌握卷积运算、线性常系数微分方程的求解编程.

二、实验原理

信号（Signal）一般都是随某一个或某几个变量的变化而变化的，例如,温度、压力、

4

声音,还有股票市场的日收盘指数等，这些信号都是随时间的变化而变化的，还有一些信号,例如在研究地球结构时，地下某处的密度就是随着海拔高度的变化而变化的。一幅图片中的每一个象素点的位置取决于两个坐标轴，即横轴和纵轴，因此,图像信号具有两个或两个以上的变量。

在《信号与系统》课程中，我们只关注这种只有一个变量（Independent variable)的信号，并且把这个变量统称为时间变量(Time variable），不管这个变量是否是时间变量. 在自然界中,大多数信号的时间变量都是连续变化的,因此这种信号被称为连续时间信号（Continuous—Time Signals）或模拟信号（Analog Signals）,例如前面提到的温度、压力和声音信号就是连续时间信号的例子。但是，还有一些信号的时间变量是离散变化的，这种信号称为离散时间信号。前面提到的股票市场的日收盘指数，由于相邻两个交易日的日收盘指数相隔24小时，这意味着日收盘指数的时间变量是不连续的,因此日收盘指数是离散时间信号. 而系统则用于对信号进行运算或处理,或者从信号中提取有用的信息，或者滤出信号中某些无用的成分，如滤波,从而产生人们所希望的新的信号.系统通常是由若干部件或单元组成的一个整体（Entity).系统可分为很多不同的类型，例如，根据系统所处理的信号的不同，系统可分为连续时间系统(Continuous-time system)和离散时间系统（Discrete—time system）,根据系统所具有的不同性质,系统又可分为因果系统（Causal system）和非因果系统（Noncausal system）、稳定系统(Stable system）和不稳定系统(Unstable system）、线性系统(Linear system)和非线性系统（Nonlinear system）、时变系统（Time-variant system）和时不变系统（Time-invariant system）等等。

然而，在信号与系统和数字信号处理中，我们所分析的系统只是所谓的线性时不变系统，这种系统同时满足两个重要的基本性质,那就是线性性和时不变性，通常称为线性时不变（LTI）系统。

1. 信号的时域表示方法

1.1将信号表示成时间变量的函数

例如 x(t)=sin(ωt) 和 x[n］=n（0.5)nu［n]

分别表示一个连续时间信号和一个离散时间信号.在MATLAB中有许多内部函数，可以直接完成信号的这种表达，例如：

sin（）：正弦信号 cos（）：余弦信号 exp（):指数信号

5

1。2用信号的波形图来描述信号

用函数曲线表示一个信号，图1。1就是一个连续时间信号和一个离散时间信号的波形图。

图1.1 连续时间信号与离散时间信号的波形图

1。3将信号用一个数据序列来表示

对于离散时间信号,还可以表示成一个数的序列，例如: x［n］=｛。。。., 0.1, 1.1, -1。2， 0, 1.3， …。｝ ↑n=0

在《信号与系统》和《数字信号处理》课程中，上述三种信号的描述方法是经常要使用的。

2 用MATLAB仿真连续时间信号和离散时间信号

在MATLAB中，无论是连续时间信号还是离散时间信号，MATLAB都是用一个数字序列

来表示信号,这个数字序列在MATLAB中叫做向量（vector）。通常的情况下，需要与时间变量相对应。

如前所述,MATLAB有很多内部数学函数可以用来产生这样的数字序列，例如sin（)、cos(）、exp（)等函数可以直接产生一个按照正弦、余弦或指数规律变化的数字序列。

2。1连续时间信号的仿真

程序Program1_1是用MATLAB对一个正弦信号进行仿真的程序，请仔细阅读该程序，

并在计算机上运行，观察所得图形.

% Program1_1

％ This program is used to generate a sinusoidal signal and draw its plot clear， % Clear all variables

close all, ％ Close all figure windows

dt = 0.01; ％ Specify the step of time variable t = —2:dt:0.2； ％ Specify the interval of time x = sin（2＊pi*t)； ％ Generate the signal

plot（t,x) % Open a figure window and draw the plot of x（t) title(’Sinusoidal signal x(t)’） xlabel('Time t （sec)')

常用的图形控制函数

6

axis(［xmin，xmax，ymin，ymax])：图型显示区域控制函数，其中xmin为横轴的显示起点，xmax为横轴的显示终点，ymin为纵轴的显示起点，ymax为纵轴的显示终点。

有时，为了使图形具有可读性，需要在所绘制的图形中,加上一些网格线来反映信号的幅度大小。MATLAB中的grid on/grid off可以实现在你的图形中加网格线.

grid on：在图形中加网格线。 grid off:取消图形中的网格线.

x = input（‘Type in signal x(t） in closed form：’)

在《信号与系统》课程中，单位阶跃信号u（t） 和单位冲激信号δ(t） 是二个非常有用的信号。它们的定义如下

t(t)dt1t0 1。1（a)

(t)0,1,u(t)0,t0 1.1（b) t0这里分别给出相应的简单的产生单位冲激信号和单位阶跃信号的扩展函数.产生单位冲激信号的扩展函数为：

function y = delta(t) dt = 0.01；

y = (u(t）-u（t—dt)）/dt;

产生单位阶跃信号的扩展函数为： % Unit step function function y = u（t)

y = （t〉=0）; % y = 1 for t > 0， else y = 0

请将这二个MATLAB函数分别以delta 和u为文件名保存在work文件夹中，以后,就可以像教材中的方法使用单位冲激信号δ(t） 和单位阶跃信号u(t）。

2.2离散时间信号的仿真

程序Program1_2用来产生离散时间信号x［n］=sin（0。2πn）。

％ Program1_2

% This program is used to generate a discrete—time sinusoidal signal and draw its plot clear, ％ Clear all variables

close all, % Close all figure windows n = -10:10; ％ Specify the interval of time x = sin（0.2＊pi＊n）； ％ Generate the signal

stem (n，x） ％ Open a figure window and draw the plot of x［n］ title ('Sinusoidal signal x[n]') xlabel （'Time index n')

请仔细阅读该程序,比较程序Program1_1和Program1_2中的不同之处,以便自己编程时能够正确使用这种方法方针连续时间信号和离散时间信号. 程序Program1_3用来仿真下面形式的离散时间信号:

7

x[n］={。.。。， 0.1, 1。1, —1.2， 0, 1.3, …。} ↑n=0

% Program1_3

% This program is used to generate a discrete-time sequence % and draw its plot

clear， % Clear all variables

close all, ％ Close all figure windows n = —5:5； ％ Specify the interval of time， the number of points of n is

11.

x = ［0， 0， 0, 0， 0.1, 1。1， -1。2, 0， 1.3, 0， 0]； % Generate the signal stem（n,x,'.') ％ Open a figure window and draw the plot of x［n] grid on,

title (’A discrete—time sequence x[n]') xlabel （'Time index n'）

由于在程序的stem(n，x,’。’) 语句中加有’.'选项，因此绘制的图形中每根棒条线的顶端是一个实心点.

如果需要在序列的前后补较多的零的话，可以利用函数zeros(）,其语法为：

zeros(1， N)：圆括号中的1和N表示该函数将产生一个一行N列的矩阵，矩阵中的所有元素均为零.利用这个矩阵与序列x[n]进行组合，从而得到一个长度与n相等的向量。

例如,当 x[n]=｛ 0。1, 1.1, -1.2， 0， 1。3} 时,为了得到程序Program1_3中的序列， ↑n=0

可以用这个MATLAB语句x = [zeros（1,4) x zeros(1， 2）] 来实现。用这种方法编写的程序如下：

％ Program1_4

% This program is used to generate a discrete-time sinusoidal signal ％ and draw its plot

clear， ％ Clear all variables close all, % Close all figure windows n = -5:5； ％ Specify the interval of time x = ［zeros（1，4）， 0.1, 1。1， -1.2， 0， 1。3, zeros（1,2）]; ％ Generate the

sequence

stem (n，x,’。’） % Open a figure window and draw the plot of x[n］ grid on,

title ('A discrete—time sequence x［n]'） xlabel （'Time index n'）

离散时间单位阶跃信号u[n］定义为 u[n]1,0,n0 1。2 n0离散时间单位阶跃信号u[n］除了也可以直接用前面给出的扩展函数来产生,还可以利用MATLAB内部函数ones（1,N） 来实现。这个函数类似于zeros(1,N)，所不同的是它产生的矩

8

阵的所有元素都为1。

值得注意的是，利用ones（1,N) 来实现的单位阶跃序列并不是真正的单位阶跃序列，而是一个长度为N单位门(Gate)序列，也就是u［n］-u[n—N］。但是在一个有限的图形窗口中，我们看到的还是一个单位阶跃序列。 在绘制信号的波形图时，有时我们需要将若干个图形绘制在图一个图形窗口中,这就需要使用MATLAB的图形分割函数subplot(），其用法是在绘图函数stem或plot之前，使用图形分割函数subplot（n1,n2，n3）,其中的参数n1，n2和n3的含义是，该函数将把一个图形窗口分割成n1xn2个子图，即将绘制的图形将绘制在第n3个子图中。

2.3信号的时域变换 2.3.1 信号的时移

信号的时移可用下面的数学表达式来描述： 设一个连续时间信号为x(t），它的时移y(t） 表示为：

y(t) = x(t - t0） 1。3

其中，t0为位移量。若t0为正数，则y(t）等于将x（t)右移t0秒之后的结果.反之，若t0为负数,则y(t）等于将x(t)左移t0秒之后的结果。

在MATLAB中，时移运算与数学上习惯表达方法完全相同。

—

程序Program1_5对给定一个连续时间信号x（t) = e0.5tu（t)，对它分别左移2秒钟和右移

—（））

2秒钟得到信号x1（t) = e0.5t+2u(t+2）和x2（t） = e-0.5(t-2u（t-2)。

％ Program1_5

% This program is used to implement the time—shift operation

% on a continuous-time signal and to obtain its time-shifted versions ％ and to draw their plots. clear,close all， t = -5：0.01：5；

x = exp（-0.5*t).＊u(t）； % Generate the original signal x(t) x1 = exp(—0。5*（t+2））。＊u（t+2); ％ Shift x(t） to the left by 2 second to get x1(t） x2 = exp（-0。5*(t-2））。＊u（t-2）； % Shift x(t) to the right by 2 second to get x2（t） subplot（311)

plot(t,x） % Plot x（t） grid on,

title （'Original signal x(t）') subplot （312)

plot （t，x1) ％ Plot x1(t) grid on,

title ('Left shifted version of x(t)’) subplot （313)

plot (t，x2) ％ Plot x2（t) grid on，

title (’Right shifted version of x(t）’) xlabel （'Time t (sec）')

9

2.3.2 信号的时域反褶

对一个信号x[n]的反褶运算在数学上表示为

y［n] = x［-n］ 1.4

这种反褶运算，用MATLAB实现起来也是非常简单的。有多种方法可以实现信号的反褶运算.

方法一，修改绘图函数plot(t,x)和stem（n，x)中的时间变量t和n，即用-t和-n替代原来的t和n,这样绘制出来的图形，看起来就是原信号经时域反褶后的版本。

方法二，直接利用原信号与其反褶信号的数学关系式来实现。这种方法最符合信号反褶运算的实际意义。

方法三，使用MATLAB内部函数fliplr(）来实现信号的反褶运算。其用法如下： y = fliplr(x)：其中x为原信号x(t)或x［n],而y则为x的时域反褶。需要说明的是，函数fliplr（)对信号作时域反褶，仅仅将信号中各个元素的次序作了一个反转，这种反转处理是于时间变量t和n的。因此，如果信号与其时间变量能够用一个数学函数来表达的话，那么建议将时间变量t和n的范围指定在一个正负对称的时间区间即可。

2。3.3 信号的时域尺度变换

信号x（t)的时域尺度变换在数学描述为

y(t) = x(at), 1。5

其中a为任意常数.根据a的不同取值，这种时域尺度变换对信号x(t)具有非常不同的影响。 当a = 1时,y（t） = x(t）；

当a = —1时，y(t) = x(—t）,即y（t)可以通过将x(t)反褶运算而得到； 当a 〉 1时，y（t) = x（at)，y（t）是将x(t)在时间轴上的压缩而得到；

当0 〈 a 〈 1时，y(t） = x（at),y（t）是将x（t)在时间轴上的扩展而得到； 当 -1 〈 a 〈 0时，y（t) = x（at），y(t)是将x（t）在时间轴上的扩展同时翻转而得到； 当 a 〈 -1时，y（t) = x(at)，y（t)是将x(t）在时间轴上的压缩同时翻转而得到;

由此可见，信号的时域尺度变换，除了对信号进行时域压缩或扩展外，还可能包括对信号的时域反褶运算。实际上,MATLAB完成式1。5的运算，并不需要特殊的处理,按照数学上的常规方法即能完成。

2。3.4周期信号

在《信号与系统》课程中，周期信号是一类非常重要的信号.给定一个信号x(t)或x［n］，如果满足

x（t） = x(t+kT） 1。6 x[n］ = x[n+kN] 1。7

则该信号叫做周期信号.其中，k为任意整数,T和N为常数，通常称为信号的基本周期或最小周期.

周期信号可以看作是一个时限的非周期信号经过周期延拓之后形成的。在数字信号处理中，

10

周期延拓这一信号处理方法非常重要。

下面的程序段，就是将一个非周期信号x1（t） = e-2t[u（t）-u（t—2）]经过周期延拓之后而得到一个周期信号:

clear， close all; t = —4：0。001:4; T = 2; x = 0； for k = -2：2；

x = x+exp（-2＊(t-k*T）)。＊(u(t-k*T）—u(t-(k+1)＊T))； end

仔细阅读该程序，可以发现其算法就是：

x(t)kx(tkT) 1。8

1由于k无法计算到无穷，而是以有限值加以替代,反映到有限宽度图形窗口中得到的效果完全符合要求。

3 LTI系统的时域描述 3。1线性时不变系统

在分析LTI系统时，有关LTI系统的两个重要的性质是必须首先掌握和理解的。这就是线

性性（Linearity）和时不变性（Time—invariance）.所谓线性性就是指系统同时满足齐次性和叠加性。这可以用下面的方法来描述。

假设系统在输入信号x1(t)作用时的响应信号为y1(t），在输入信号x2(t）作用时的响应信号为y2(t），给定两个常数a和b，如果当输入信号为x(t）时系统的响应信号为y（t），且满足 x（t) = x1（t） + x2（t) 1。9(a） y（t） = y1(t） + y2（t) 1。9(b） 则该系统具有叠加性(Additivity）.如果满足

x(t） = ax1(t） 1。10（a) y（t） = ay1(t) 1。10（b）

则该系统具有齐次性(Homogeneity）。一个系统如果是线性系统的话，那么这个系统必须同时具有叠加性和齐次性.

又假设系统在输入信号x（t)作用时的响应信号为y（t），对一个给定时间常数t0，如果当输入信号为x(t-t0）时,系统的响应信号为y(t—t0）的话，则该系统具有时不变性。

同时具有线性性和时不变性的系统，叫做线性时不变系统，简称LTI系统。LTI系统有连续时间LTI系统和离散时间LTI系统之分。连续时间系统的输入和输出信号都必须是连续时间

11

信号，而离散时间系统的输入和输出信号都必须是离散时间信号。

3.2 LTI系统的单位冲激响应和卷积模型

给定一个连续时间LTI系统，在系统的初始条件为零时,用单位冲激信号δ(t)作用于系统,

此时系统的响应信号称为系统的单位冲激响应（Unit impulse response)，一般用h(t)来表示。需要强调的是，系统的单位冲激响应是在激励信号为δ（t)时的零状态响应(Zero—state response）。 离散时间LTI系统的单位冲激响应的定义与连续时间LTI系统的单位冲激响应相同，只是离散时间单位冲激函数δ[n]的定义有所不同。 系统的单位冲激响应是一个非常重要的概念，对于一个系统,如果我们知道了该系统的单位冲激响应，那么,该系统对任意输入信号的响应信号都可以求得。也就是说，系统的输入信号x（t)、x［n]和输出信号y（t）、y[n］之间的关系可以用一个数学表达式来描述,这个数学表达式为

y(t) y[n]x()h(t)d 1.11(a)

kx[k]h[nk] 1.11（b）

这个表达式就是LTI系统的卷积模型，它是根据系统的线性性和时不变性以及信号可以分解成单位冲激函数经推理得到的。这个表达式实际上告诉了我们一个重要的结论，那就是，任意LTI系统可以完全由它的单位冲激响应h(t)/h［n］来确定.由于系统的单位冲激响应是零状态响应，故按照式1。11求得的系统响应也是零状态响应.式1.11中的积分运算叫做卷积积分,求和运算叫做卷积和，是描述连续时间系统输入输出关系的一个重要表达式。

3.3卷积的计算

卷积的计算通常可按下面的五个步骤进行（以卷积积分为例）：

1. 该换两个信号波形图中的横坐标，由t改为τ，τ变成函数的自变量； 2. 把其中一个信号反褶，如把h（τ)变成h(-τ);

3. 把反褶后的信号做移位，移位量是t，这样t是一个参变量。在τ坐标系中，t > 0时图形右移， t 〈 0时图形左移。

4. 计算两个信号重叠部分的乘积x(τ)h(t-τ)； 5. 完成相乘后图形的积分。

对于两个时限信号（Time—limited signal)，按照上述的五个步骤，作卷积积分运算时，关键是正确确定不同情况下的积分限。只要正确地确定了积分限都能得到正确定积分结果.尽管如此，在时域中计算卷积积分，总体上来说是一项比较困难的工作。

程序convlution_demo用来演示上述作卷积积分运算的五个步骤。本程序较为复杂，不建议读者读懂该程序，只需执行这个程序，观看程序执行过程中有关卷积积分的运算过程，以便于理解这五个步骤。

12

借助MATLAB的内部函数conv(）可以很容易地完成两个信号的卷积积分运算。其语法为：y = conv(x,h)。其中x和h分别是两个作卷积运算的信号，y为卷积结果。

为了正确地运用这个函数计算卷积，这里有必要对conv(x，h)做一个详细说明。conv（x,h)函数实际上是完成两个多项式的乘法运算.例如，两个多项式p1和p2分别为:

p1s32s23s4 和 p24s33s22s1

这两个多项式在MATLAB中是用它们的系数构成一个行向量来表示的，如果用x来表示多项式p1,h表示多项式p2，则x和h分别为 x = ［1 2 3 4] h = [4 3 2 1］ 在MATLAB命令窗口依次键入

〉〉 x = ［1 2 3 4］； 〉> h = [4 3 2 1]； >〉 y=conv（x,h） 在屏幕上得到显示结果:

y = 4 11 20 30 20 11 4 这表明,多项式p1和p2的乘积为：

p34s611s520s430s320s211s4

正如前所述,用MATLAB处理连续时间信号时，时间变量t的变化步长应该是很小的,假定用符号dt表示时间变化步长，那么，用函数conv(）作两个信号的卷积积分时，应该在这个函数之前乘以时间步长方能得到正确的结果。也就是说,正确的语句形式应为:y = dt*conv(x,h）。

对于定义在不同时间段的两个时限信号x(t)，t0 ≤ t ≤ t1，和h（t），t2 ≤ t ≤ t3. 如果用y(t)来表示它们的卷积结果,则y(t）的持续时间范围要比x（t)或h(t）要长,其时间范围为t0+t2 ≤ t ≤ t1+t3.这个特点很重要，利用这个特点,在处理信号在时间上的位置时，可以很容易地将信号的函数值与时间轴的位置和长度关系保持一致性。

根据给定的两个连续时间信号x（t) = t[u(t)-u（t-1)］和h（t） = u(t）-u（t-1)，编写程序，完成这两个信号的卷积运算，并绘制它们的波形图。范例程序如下：

％ Program1_6

％ This program computes the convolution of two continuou—time signals clear;close all;

t0 = —2； t1 = 4; dt = 0.01； t = t0:dt：t1; x = u(t）-u(t-1）;

h = t。＊(u（t)-u（t—1)）;

y = dt＊conv(x,h)； ％ Compute the convolution of x(t） and h(t) subplot（221)

plot(t,x)， grid on， title（'Signal x（t）'）, axis（[t0,t1,-0。2,1。2］) subplot（222）

plot（t，h)， grid on， title（'Signal h（t）’)， axis（[t0,t1,—0。2，1。2])

13

subplot（212）

t = 2*t0：dt：2*t1； % Again specify the time range to be suitable to the ％ convolution of x and h.

plot(t，y)， grid on， title('The convolution of x(t) and h（t）’)， axis([2＊t0，2＊t1，

—0.1,0.6］）,

xlabel(’Time t sec')

在有些时候，做卷积和运算的两个序列中，可能有一个序列或者两个序列都非常长,甚至是无限长，MATLAB处理这样的序列时，总是把它看作是一个有限长序列，具体长度由编程者确定。实际上,在信号与系统分析中所遇到的无限长序列，通常都是满足绝对可和或绝对可积条件的信号。因此,对信号采取这种截短处理尽管存在误差，但是通过选择合理的信号长度,这种误差是能够减小到可以接受的程度的.若这样的一个无限长序列可以用一个数学表达式表示的话，那么，它的长度可以由编程者通过指定时间变量n的范围来确定。

例如，对于一个单边实指数序列x［n］ = 0.5nu［n]，通过指定n的范围为0 ≤n ≤ 100，则对应的x［n]的长度为101点,虽然指定更宽的n的范围,x［n］将与实际情况更相符合，但是，注意到，当n大于某一数时，x[n］之值已经非常接近于0了。对于序列x[n] = 0.5nu［n]，当n = 7时,x[7] = 0。0078，这已经是非常小了。所以，对于这个单边实指数序列，指定更长的n的范围是没有必要的。当然，不同的无限长序列具有不同的特殊性，在指定n的范围时，只要能够反映序列的主要特征就可以了。

3。4 用线性常系数微分方程描述LTI系统

线性常系数微分方程或差分方程是描述LTI系统的另一个时域模型.一个连续时间LTI系统，它的输入信号x(t)输出信号y(t）关系可以用下面的微分方程来表达

dky(t)Mdkx(t) ak 1。12 bkkkdtdtk0k0N式1.12中，max （N, M）定义为系统的阶。式1。12描述了LTI系统输入信号和输出信号的一种隐性关系（Implicit relationship）。为了求得系统响应信号的显式表达式（Explicit expression）,必须对微分方程和差分方程求解。

在MATLAB中，一个LTI系统也可以用系统微分方程的系数来描述，例如,一个LTI连续时间系统的微分方程为

d2y(t)dy(t)32y(t)x(t)

dtdt2MATLAB则用两个系数向量num = ［1]和den = ［1 3 2]来描述该系统，其中num和den

分别表示系统微分方程右边和左边的系数，按照微分运算的降阶排列.

MATLAB的内部函数impulse(）,step（），initial(），lsim（） 可以用来计算并绘制连续时间LTI系统的单位冲激响应，单位阶跃响应，零输入响应和任意信号作用于系统的零状态响应.这些函数的用法描述如下:

14

h= impulse（num， den， T） 和 impulse(num， den， T） s = step(num， den， T） 和 step（num， den， T）

y = lmis(num, den, x, t) 和 lmis(num， den, x， t)

函数impulse（），step()用来计算由num和den表示的LTI系统的单位冲激响应和单位阶跃响应，响应的时间范围为0~T,其中den和num分别为系统微分方程左右两边的系数向量，T为指定的响应的终点时间.h和s的点数默认值为101点.由此可以计算时间步长为dt = T/（101—1).不带返回值的函数如impulse（num， den, T)和step（num, den, T)将直接在屏幕上绘制系统的单位冲激响和单位阶跃响应曲线。带返回值的函数如y = lmis(num, den， x， t)和y = lmis(num， den, x, t）,用来计算由num和den表示的LTI系统在输入信号x作用下的零状态响应。其中t为指定的时间变化范围，x为输入信号，它们的长度应该是相同的。如带返回参数y，则将计算的响应信号保存在y中，若不带返回参数y,则直接在屏幕上绘制输入信号x和响应信号y的波形图.

例如，编写程序，计算并绘制由下面的微分方程表示的系统的单位冲激响应h(t），单位阶跃响应s(t)。

d2y(t)dy(t)32y(t)8x(t)

dtdt2MATLAB范例程序如下:

％ Program1_7

% This program is used to compute the impulse response h(t） and the step response s(t) of a

％ continuous—time LTI system clear， close all;

num = input（’Type in the right coefficient vector of differential equation：'); den = input('Type in the left coefficient vector of differential equation：')； t = 0：0。01:8;

x = input（'Type in the expression of the input signal x（t）：'); subplot(221)， impulse(num,den,8)； subplot(222), step（num,den，8）

三、实验内容及步骤

实验前，必须首先阅读本实验原理，读懂所给出的全部范例程序。实验开始时,先在计算机上运行这些范例程序，观察所得到的信号的波形图.并结合范例程序应该完成的工作，进一步分析程序中各个语句的作用，从而真正理解这些程序.

实验前，一定要针对下面的实验项目做好相应的实验准备工作，包括事先编写好相应的实验程序等事项。

15

Q1-1：修改程序Program1_1，将dt改为0.2,再执行该程序，保存图形，看看所得图形的效果

如何?

dt = 0.01时的信号波形 dt = 0.2时的信号波形

这两幅图形有什么区别,哪一幅图形看起来与实际信号波形更像？ 答:

Q1-2：修改程序Program1_1,并以Q1_2为文件名存盘，产生实指数信号x(t）=e-2t。 要求在

图形中加上网格线，并使用函数axis()控制图形的时间范围在0~2秒之间。然后执行该程序，保存所的图形。

修改Program1_1后得到的程序Q1_2如下：

% Program1_2

％ This program is used to generate a sinusoidal signal and draw its plot clear， % Clear all variables

close all， % Close all figure windows

dt = 0。2； % Specify the step of time variable t = -2：dt:2; ％ Specify the interval of time x = exp（—2*t）； % Generate the signal

plot(t,x) ％ Open a figure window and draw the plot of x（t） axis（[0 2 0 1]) grid on；

title('Sinusoidal signal x（t)'）

xlabel(’Time t (sec）')

—

信号x(t)=e0.5t的波形图

16

Q1-3：修改程序Program1_1，并以Q1_3为文件名存盘,使之能够仿真从键盘上任意输入的一

个连续时间信号,并利用该程序仿真信号x(t)=e-2t。

修改Program1_1后得到的程序Q1_3如下： 信号x（t）=e

—2t

的波形图

% Program1_3 clc

clear, % Clear all variables

close all, ％ Close all figure windows

dt = 0。01； % Specify the step of time variable t = —2：dt：2； ％ Specify the interval of time

x = input('please input a signal,I will draw its plot for you Singnal x=’）； %Input the signal plot(t,x） % Open a figure window and draw the plot of x（t） xlabel（’Time t （sec）'） grid on

%x=exp(-2*t

）

17

此处粘贴图形

Q1—4：将实验原理中所给的单位冲激信号和单位阶跃信号的函数文件在

辑器中编写好,并分别以以文件名delta和u存入work文件夹中以便于使用.

MATLAB文件编

抄写函数文件delta如下： 抄写函数文件u如下：

t(t)dt1t0 u[n](t)0,

1,0,n0 n0Q1—5：修改程序Program1_4，并以Q1_5为文件名存盘，利用axis()函数，将图形窗口的横

坐标范围改为-2≤n≤5，纵坐标范围改为—1.5≤ x ≤1.5。

修改Program1_4后得到的程序Q1_5如下： 信号的波形图

此处粘贴图形 >〉 % Program1_5

％ This program is used to generate a discrete—time sinusoidal signal % and draw its plot

clear， % Clear all variables

close all, ％ Close all figure windows n = -5：5； % Specify the interval of time

x = ［zeros（1,4), 0.1, 1。1， —1。2， 0, 1.3， zeros（1,2)]； % Generate the sequence stem （n，x,’filled’,’r’) ％ Open a figure window and draw the plot of x[n］

18

axis（［—2 5 —1。5 1。5]) grid on,

title (’A discrete-time sequence x［n]’) xlabel (’Time index n’)

Q1-6:仿照前面的示例程序的编写方法，编写一个MATLAB程序，以Q6为文件名存盘，使

之能够在同一个图形窗口中的两个子图中分别绘制信号x［n］=0。5n 和x（t）=cos(2πt）[u（t)—u（t-3）]。要求选择的时间窗能够表现出信号的主要部分（或特征)。

｜｜

编写的程序Q1_6如下：

信号x[n］=0.5|n| 的波形图和信号x(t）=cos(2πt)［u(t)-u（t—3）］的波形图

19

Q1-7：根据示例程序的编程方法，编写一个MATLAB程序,以Q1_7为文件名存盘，由给定

信号

—。

x（t） = e05tu(t) 求信号y(t） = x(1。5t+3)，并绘制出x(t） 和y(t)的图形.

编写的程序Q1_7如下：

信号x(t)的波形图 和 信号y（t） = x(1.5t+3) 的波形图

此处粘贴图形 此处粘贴图形

Q1-8：给定一个离散时间信号x[n］ = u［n] – u［n—8］，仿照示例程序Program1_5，编写

程序Q1_8,产生x[n]的左移序列x1［n] = x[n+6］和右移序列x2[n］ = x[n-6］，并在同一个图形窗口的三个子图中分别绘制这三个序列的图形。

编写的程序Q1_8如下：

20

信号波形图

此处粘贴图形

Q1-9：编写程序

Q1_9，使之能够接受以键盘方式输入的定义在不同时间段的两个不同连续

时间信号并完成卷积运算，分别绘制这两个信号及其卷积的结果的图形,图形按照22分割成四个子图.

编写的程序Q1_9如下：

信号x (t）、h(t)和x （t)＊h(t)的波形图

此处粘贴图形

21

Q1-10：给定两个离散时间序列

x［n] = 0.5n｛u［n]-u[n-8］｝ h[n］ = u[n]-u[n-8］

编写程序Q1_10，计算它们的卷积，并分别绘制x［n]、h[n］和它们的卷积y［n]的图形.

编写的程序Q1_10如下:

信号x[n］、h［n］和y［n]的波形图

此处粘贴图形

22

Q1—11已知一个序列为

x[n]n,0,0n4

otherwise编写MATLAB程序Q1_11，能够将x[n]以N = 8为周期进行周期延拓得到一个周期为N =8的周期序列y［n］，并分别绘制x[n]和y［n］图形。

编写的程序Q1_11如下：

信号x［n]的波形图 信号y[n]的波形图 此处粘贴图形 此处粘贴图形

23

Q1—12 仿照范例程序Program1_7，编写程序Q1_12,计算并绘制由如下微分方程表示的系统

在输入信号为x（t) = (e-2t — e-3t)u（t）时的零状态响应和你手工计算得到的系统零状态响应曲线。

d2y(t)dy(t)32y(t)8x(t)

dtdt2手工计算得到的系统零状态响应的数学表达式是：

编写的程序Q1_12如下： 用MATLAB绘制的手工计算的系统响应

粘帖用MATLAB绘制的手工计算的系统响应

执行程序Q1_12得到的系统响应

此处粘帖执行程序Q1_12得到的系统响应

Q1-13:Q1—13:利用程序Q1_9，验证卷积的相关性质。

(a） 验证性质：x(t)*(t)x(t)

选择信号x（t)的数学表达式为：

24

x（t）、δ（t）和x(t）*δ（t）的波形

验证所得结论是:

(b) 验证性质：x(t)*(tt0)x(tt0)

选择信号x(t）的数学表达式为：

x(t)、δ（t—t0） 和x(t)*(tt0)的波形

验证所得结论是：

（c) 验证性质:x(tt1)*(tt2)x(tt2)*(tt1)x(tt1t2) 25

选择信号x（t）的数学表达式为： 选择的t1 = 秒,t2 = 秒。 执行程序Q1_9，输入信号x（t-t1) 和δ(t—t2) 的数学表达式，得到的信号及其卷积的波形图如

下：

执行程序Q1_9，输入信号x(t-t2) 和δ（t—t1）如下:

验证所得结论是：

26

的数学表达式，得到的信号及其卷积的波形图

(d) 验证性质：x(t)*u(t)tx(t)dt

选择信号x(t）（建议选择一个时限信号)的数学表达式为：



tx(t)dt的数学表达式为：

t手工绘制的x(t)dt波形如下:

执行程序Q1_9，输入信号x(t) 和u（t) 的数学表达式,得到的信号及其卷积的波形图如下：

验证所得结论是:

(e） 验证性质:x(t)*h(tt0)x(tt0)*h(t)

选择信号x（t)的数学表达式为： 选择信号h（t)的数学表达式为： 选择的t0=：

执行程序Q1_9，输入信号x(t) 和h(t-t0） 的数学表达式，得到的信号及其卷积的波形图如下：

27

执行程序Q1_9,输入信号x（t-t0) 和h(t) 的数学表达式，得到的信号及其卷积的波形图如下:

验证所得结论是:

Q1-14：做如下总结：

1、信号与系统分析，就是基于信号的分解，在时域中，信号主要分解成:

2、写出卷积的运算步骤，并谈谈你对卷积的一些基本性质的理解。利用MATLAB计算卷积的函数是什么？如何使用？

28

3、在时域中，描述一个连续时间LTI系统的数学模型有：

4、MATLAB是如何表示一个由微分方程描述的连续时间LTI系统的？求解连续时间LTI系统的单位冲激响应、单位阶跃响应以及系统在某一个输入信号作用下的零状态响应的MATLAB函数有哪些?

四、实验报告要求

1、按要求完整书写你所编写的全部MATLAB程序

2、详细记录实验过程中的有关信号波形图（存于自带的U盘中),图形要有明确的标题.全部的MATLAB图形应该用打印机打印，然后贴在本实验报告中的相应位置,禁止复印件。

3、实事求是地回答相关问题，严禁抄袭。

本实验完成时间： 年 月 日

实验二 连续时间信号的频域分析

一、实验目的

29

1、掌握连续时间周期信号的傅里叶级数的物理意义和分析方法； 2、观察截短傅里叶级数而产生的“Gibbs现象”，了解其特点以及产生的原因； 3、掌握连续时间傅里叶变换的分析方法及其物理意义；

4、掌握各种典型的连续时间非周期信号的频谱特征以及傅里叶变换的主要性质；

5、学习掌握利用MATLAB语言编写计算CTFS、CTFT和DTFT的仿真程序，并能利用这些程序对一些典型信号进行频谱分析，验证CTFT、DTFT的若干重要性质。

基本要求:掌握并深刻理傅里叶变换的物理意义,掌握信号的傅里叶变换的计算方法，掌握利用MATLAB编程完成相关的傅里叶变换的计算。

二、实验原理及方法

1、连续时间周期信号的傅里叶级数CTFS分析

任何一个周期为T1的正弦周期信号，只要满足狄利克利条件，就可以展开成傅里叶级数。 其中三角傅里叶级数为：

x(t)a0[akcos(k0t)bksin(k0t)] 2。1

k1或： x(t)a0ck1kcos(k0tk) 2。2

其中02，称为信号的基本频率(Fundamental frequency），a0，ak，和bk分别是信号T1x(t)的直流分量、余弦分量幅度和正弦分量幅度,ck、k为合并同频率项之后各正弦谐波分量

的幅度和初相位,它们都是频率k0的函数，绘制出它们与k0之间的图像，称为信号的频谱图(简称“频谱”），ck－k0图像为幅度谱，k－k0图像为相位谱.

三角形式傅里叶级数表明，如果一个周期信号x(t），满足狄里克利条件,那么，它就可以被看作是由很多不同频率的互为谐波关系（harmonically related）的正弦信号所组成,其中每一个不同频率的正弦信号称为正弦谐波分量 （Sinusoid component），其幅度(amplitude)为ck。也可以反过来理解三角傅里叶级数：用无限多个正弦谐波分量可以合成一个任意的非正弦周期信号。

指数形式的傅里叶级数为：

x(t)

kaekjk0t 2。3

30

其中,ak为指数形式的傅里叶级数的系数,按如下公式计算:

1 akT1T1/2T1/2jk0tx(t)edt 2.4 指数形式的傅里叶级数告诉我们，如果一个周期信号x（t)，满足狄里克利条件,那么，它

就可以被看作是由很多不同频率的互为谐波关系（harmonically related）的周期复指数信号所组成，其中每一个不同频率的周期复指数信号称为基本频率分量,其复幅度（complex amplitude）为ak。这里“复幅度（complex amplitude）\"指的是ak通常是复数.

上面的傅里叶级数的合成式说明，我们可以用无穷多个不同频率的周期复指数信号来合成任意一个周期信号。然而，用计算机（或任何其它设备)合成一个周期信号,显然不可能做到用无限多个谐波来合成，只能取这些有限个谐波分量来近似合成。

假设谐波项数为N，则上面的和成式为:

x(t)kNaekNjk0t 2。5

显然，N越大，所选项数越多,有限项级数合成的结果越逼近原信号x(t)。本实验可以比较直观地了解傅里叶级数的物理意义，并观察到级数中各频率分量对波形的影响包括“Gibbs”现象：即信号在不连续点附近存在一个幅度大约为9%的过冲,且所选谐波次数越多，过冲点越向不连续点靠近.这一现象在观察周期矩形波信号和周期锯齿波信号时可以看得很清楚。

2、连续时间信号傅里叶变换———-CTFT

傅里叶变换在信号分析中具有非常重要的意义，它主要是用来进行信号的频谱分析的。傅

里叶变换和其逆变换定义如下:

 X(j)jtx(t)edt 2.6 1x(t)2X(j)ejtd 2.7

连续时间傅里叶变换主要用来描述连续时间非周期信号的频谱。按照教材中的说法，任意非周期信号，如果满足狄里克利条件，那么,它可以被看作是由无穷多个不同频率（这些频率都是非常的接近）的周期复指数信号ejt的线性组合构成的，每个频率所对应的周期复指数信号ejt称为频率分量（frequency component），其相对幅度为对应频率的|X(j）｜之值，其相位为对应频率的X（j）的相位.

X(j）通常为关于的复函数，可以按照复数的极坐标表示方法表示为：

（）

X（j）=｜ X(j）|ej Xj

其中,｜ X（j)|称为x(t）的幅度谱，而X（j）则称为x(t）的相位谱。

给定一个连续时间非周期信号x（t)，它的频谱也是连续且非周期的。对于连续时间周期信号，也可以用傅里变换来表示其频谱，其特点是，连续时间周期信号的傅里叶变换时有冲激序列构成的，是离散的——这是连续时间周期信号的傅里叶变换的基本特征。

31

3、离散时间序列的傅里叶变换---DTFT

给定一个非周期离散时间序列x［n］，它的傅里叶变换定义为 X(ej)nx[n]ejn 2。8

1其反变换定义为 x[n]2jjnX(e)ed 2。9  式2。8称为离散时间傅里叶正变换，式2.9称为离散时间傅里叶反变换.由式2。9可以看

出，一个非周期离散时间序列，总是可以被看作是由无穷多个不同频率的加权的基本频率分量ejωn组合而成的。对序列中的频率为ω的频率分量来说,其权为X(ejω）,通常是复数，也可以将它表示为

（）

X（ejω) = ｜ X（ejω)｜ejXejω

｜X(ejω）｜称为序列的幅度谱，而X（ejω)称为序列的相位谱，它们都是频率ω的周期函数。

4、连续周期信号的傅里叶级数CTFS的MATLAB实现 4。1 傅里叶级数的MATLAB计算

设周期信号x(t)的基本周期为T1，且满足狄里克利条件,则其傅里叶级数的系数可由式2.4计算得到。式2.4重写如下:

1akT1基本频率为: 0T1/2T1/2jk0tx(t)edt 2 T1对周期信号进行分析时，我们往往只需对其在一个周期内进行分析即可，通常选择主周期(Principle period）。假定x1（t)是x（t)中的主周期，则

1akT1T1/2T1/2jk0tx(t)edt 1计算机不能计算无穷多个系数，所以我们假设需要计算的谐波次数为N，则总的系数个数

为2N+1个。在确定了时间范围和时间变化的步长即T1和dt之后，对某一个系数，上述系数的积分公式可以近似为：

11jk0tjk0takx(t)edtx(t)edt/T1 1nT1T1/2n [x(t1),x(t2),x(tM)][e

jk0t1T/2,ejk0t2,ejk0tM]dt/T1

32

对于全部需要的2N+1个系数，上面的计算可以按照矩阵运算实现.MATLAB实现系数计算

的程序如下：

dt = 0.01;

T = 2; t = -T/2：dt:T/2； w0 = 2*pi/T；

x1 = input(‘Type in the periodic signal x(t） over one period x1(t）=’）; N = input（‘Type in the number N=’); k = -N:N； L = 2＊N+1； ak = x1＊exp(-j＊k*w0＊t’)*dt/T;

需要强调的是,时间变量的变化步长

dt的大小对傅里叶级数系数的计算精度的影响非常

大，dt越小,精度越高，但是，计算机计算所花的时间越长。

例题2-1：给定一个周期为T1 = 2s的连续时间周期方波信号，如图所示，其一个周期内的数

学表达式为：

1, x1(t)0,0t1

1t2-2 -1 1 x(t) 0 1 2 t 解：首先，我们根据前面所给出的公式,计算该信号的傅里叶级数的系数。

图2.1 周期方波信号

101akT1jk0ejkt11jk0t1jkt0jkted(jkt)x(t)edtedt01j2k00j2k02T/20T1/21001e1ej2k0kj02ekj02ej2k0kj02ksin(0)jk02e2 k0sin(因为:0 = 2π/T1 = π，代入上式得到:ak(j)kkk)2

在MATLAB命令窗口，依次键入：

〉〉 k = —10：10; >> ak = （（-j）。^k).* (sin（（k+eps)＊pi/2）./（（k+eps)＊pi）） ％ The expression of ak ak =

Columns 1 through 4

-0.0000 0 + 0.0354i —0.0000 0 +

0.0455i

Columns 5 through 8

—0.0000 0 + 0.0637i —0.0000 0 + 0。

1061i

Columns 9 through 12

-0。0000 0 + 0。3183i 0。5000 0 -

0.3183i

33

Columns 13 through 16

-0。0000 0 - 0。1061i —0。0000 0 -

0.0637i

Columns 17 through 20

—0.0000 0 - 0。0455i —0。0000 0 —

0.0354i

Column 21

—0。0000

从MATLAB命令窗口，我们得到了该周期信号从a10到a10共21个系数。 紧接着再键入以下命令：

>> subplot（221）

>〉 stem（k，abs(ak),'k。') //abs是求绝对值的

函数

>〉 title(’The Fourier series coefficients’) >〉 xlabel(’Frequency index k'）

就得到一幅如右图所示的描述ak与k之间的关系的

图形。

以上是我们通过手工计算得到的这个周期信号的傅里叶级数表达式及其频谱图，下面给出完成傅里叶级数系数计算的相应MATLAB范例程序。

% Program2_1

％ This program is used to evaluate the Fourier series coefficients ak of a periodic square

wave

clear, close all

T = 2； dt = 0。00001; t = -2：dt:2; x1 = u(t） - u(t—1-dt); x = 0;

for m = —1：1 ％ Periodically extend x1（t) to form a periodic

signal

x = x + u（t—m＊T) — u(t-1—m*T—dt）； end

w0 = 2*pi/T；

N = 10; ％ The number of the harmonic components L = 2＊N+1;

for k = -N： N； ％ Evaluate the Fourier series coefficients ak ak（N+1+k) = (1/T)*x1*exp（—j＊k*w0＊t’）＊dt； end

phi = anglel(ak)； % Evaluate the phase of ak

执行程序Program2_1后，就完成了信号的傅里叶级数的系数的计算,在命令窗口键入

〉> ak

命令窗口就可以显示傅里叶级数的21个系数：

34

ak =

Columns 1 through 4

0.0000 + 0.0000i 0.0000 + 0。0354i 0。0000 - 0.0000i 0。0000 + 0.0455i Columns 5 through 8

0.0000 — 0。0000i 0.0000 + 0.0637i 0.0000 - 0.0000i 0。0000 + 0。1061i Columns 9 through 12

0.0000 - 0。0000i 0。0000 + 0。3183i 0.5000 0.0000 - 0.3183i Columns 13 through 16

0.0000 + 0.0000i 0。0000 - 0。1061i 0.0000 + 0。0000i 0.0000 — 0.0637i Columns 17 through 20

0。0000 + 0。0000i 0。0000 - 0.0455i 0。0000 + 0。0000i 0.0000 — 0.0354i Column 21

0。0000 - 0.0000i

将这里的ak之值同前面手工计算得到的ak比较，可见两者是完全相同的。

再次特别提示:程序中，时间变量的变化步长dt的大小对傅里叶级数系数的计算精度的影

响非常大,dt越小，精度越高,本程序中的dt之所以选择0.00001就是为了提高计算精度.但是，计算机所花的计算时间越长。

在程序Program2_1中添加相应的计算| ak ｜和绘图语句，就可以绘制出信号的幅度谱和相位谱的谱线图。

4.2 周期信号的合成以及Gibbs现象

从傅里叶级数的合成式（Synthesis equation）

x(t)kaekjk0t

可以看出，用无穷多个不同频率和不同振幅的周期复指数信号可以合成一个周期信号.然而，我们无法用计算机实现对无穷多个周期复指数信号的合成。但是，用有限项来合成却是可行的,在实际应用中，多半也就是这么做的。然而,这样做的一个必然结果，就是引入了误差。

如果一个周期信号在一个周期有内断点存在，那么,引入的误差将除了产生纹波之外，还将在断点处产生幅度大约为9%的过冲（Overshot）,这种现象被称为吉伯斯现象（Gibbs phenomenon）.

为了能够观察到合成信号与原信号的不同以及Gibbs现象，我们可以利用前面已经计算出的傅里叶级数的系数，计算出截短的傅里叶级数：

x(t)kNaekNjk0t

这个计算可用L = 2N+1次循环来完成：

x2x2ak(r)ej(r1N)0t

35

其中r作为循环次数，x2在循环之前应先清零。完成这一计算的MATLAB程序为：

x2 = 0; L = 2＊N+1; for r = 1:L；

x2 = x2+ak（r)＊exp（j＊（r—1-N）*w0＊t); end;

完成了所有的计算之后，就可以用绘图函数：plot()和stem()将计算结果包括x1, x2, abs(ak）和angle（ak)以图形的形式给出,便于我们观察。

观察吉伯斯现象的最好的周期信号就是图2-1所示的周期方波信号,这种信号在一个周期内有两个断点，用有限项级数合成这个信号时,吉伯斯现象的特征非常明显，便于观察。

例题2-2：修改程序Program2_1，使之能够用有限项级数合成例题2-1所给的周期方波信号，

并绘制出原始周期信号、合成的周期信号、信号的幅度谱和相位谱。

为此，只要将前述的for循环程序段和绘图程序段添加到程序Program2_1中即可，范例程序如下：

% Program2_2

% This program is used to compute the Fourier series coefficients ak of a periodic square wave clear，close all

T = 2； dt = 0.00001; t = -2:dt：2; x1 = u(t)—u（t-1-dt）； x = 0； for m = —1：1

x = x + u(t—m＊T） - u（t—1-m*T—dt)； ％ Periodically extend x1(t） to form a periodic signal

end

w0 = 2*pi/T;

N = input（'Type in the number of the harmonic components N = ：'）; L = 2＊N+1；

for k = —N：1:N；

ak(N+1+k） = (1/T)*x1＊exp（—j*k*w0*t’）*dt; end

phi = angle(ak）； y=0；

for q = 1：L; % Synthesiz the periodic signal y(t) from the finite Fourier series y = y+ak(q）＊exp(j＊(—（L—1）/2+q—1）*2*pi＊t/T）； end;

subplot（221),

plot（t，x)， title('The original signal x（t)’), axis([-2,2,-0.2,1。2]), subplot（223), plot（t,y), title（'The synthesis signal y(t)'）， axis([—2,2，-0.2，1.2］）， xlabel（'Time t’),

subplot（222）

k=—N:N; stem(k,abs(ak）,’k.'）， title（'The amplitude ｜ak| of x(t）’)， axis([—N,N，-0.1，0。6］）

36

subplot(224） stem（k,phi，'r。’), title(’The phase phi（k) of x(t）'）, axis([-N,N，-2，2］), xlabel(’Index k'）

在用这个程序观察吉伯斯现象时，可以反复执行该程序,每次执行时，输入不同之N值,比较所的图形的区别，由此可以观察到吉伯斯现象的特征。

5 用MATLAB实现CTFT及其逆变换的计算 5。1 用MATLAB实现CTFT的计算

MATLAB进行傅里叶变换有两种方法，一种利用符号运算的方法计算，另一种是数值计算，本实验要求采用数值计算的方法来进行傅里叶变换的计算。严格来说，用数值计算的方法计算连续时间信号的傅里叶变换需要有个限定条件，即信号是时限信号（Time limited signal），也就是当时间|t｜大于某个给定时间时其值衰减为零或接近于零，这个条件与前面提到的为什么不能用无限多个谐波分量来合成周期信号的道理是一样的。计算机只能处理有限大小和有限数量的数。

采用数值计算算法的理论依据是：

 X(j)x(t)eNjtdtlimx(kT)ejkTT

T0k若信号为时限信号，当时间间隔T取得足够小时,上式可演变为：

X(j)TkNx(kT)ejkT

[x(t1),x(t2),,x(t2N1)][ejt1,ejt2,,ejt2N1]T

上式用MATLAB表示为:

X=x*exp（j＊t’*w)*T

其中X为信号x(t）的傅里叶变换,w为频率Ω，T为时间步长。 相应的MATLAB程序：

T = 0.01； dw = 0。1; ％时间和频率变化的步长 t = -10：T：10;

w = —4＊pi:dw:4*pi;

X(j）可以按照下面的矩阵运算来进行：

X=x*exp(-j＊t’*）*T; %傅里叶变换 X1=abs(X）； %计算幅度谱 phai=angle（X）; %计算相位谱

为了使计算结果能够直观地表现出来,还需要用绘图函数将时间信号x(t)，信号的幅度谱｜X(j）|和相位谱 X(j)分别以图形的方式表现出来，并对图形加以适当的标注。

5.2 用MATLAB实现傅里叶逆变换

连续时间傅里叶逆变换可用式2.7进行计算。式2.7重写如下：

37

1x(t)2jtX(j)ed 从定义式可看出,其计算方法与傅里叶变换是一样的，因此可以采用同样的矩阵运算的方法

来计算，即

x(t)=X(j)*exp(j’＊t）*d 具体的MATLAB函数如下：

t = —5：0.01；5; ％ 指定信号的时间范围，此范围应根据信号的持续时间确定。 dw = 0。1; w = —4＊pi:d：4*pi；

X = input(‘Type in the expression of X(jw）’）； x = X* exp（jw’＊t)＊dw;

然后用绘图函数就可以绘制出逆变换得到的时域信号波形图。

三、实验内容和要求

实验前，必须首先阅读本实验原理，读懂所给出的全部范例程序。实验开始时，先在计算机上运行这些范例程序,观察所得到的信号的波形图.并结合范例程序应该完成的工作，进一步分析程序中各个语句的作用,从而真正理解这些程序。 实验前，一定要针对下面的实验项目做好相应的实验准备工作,包括事先编写好相应的实验程序等事项。

Q2—1 编写程序Q2_1，绘制下面的信号的波形图：

111n x(t)cos(0t)cos(30t)cos(50t)sin()cos(n0t)

35n2n1其中，0 = 0.5π，要求将一个图形窗口分割成四个子图,分别绘制cos(0t)、cos(30t）、cos（50t）

和x(t) 的波形图，给图形加title，网格线和x坐标标签，并且程序能够接受从键盘输入的和式中的项数。

抄写程序Q2_1如下：

执行程序Q2_1所得到的图形如下：

38

Q2-2 给程序Program2_1增加适当的语句，并以Q2_2存盘，使之能够计算例题2-1中的周期

方波信号的傅里叶级数的系数，并绘制出信号的幅度谱和相位谱的谱线图。

通过增加适当的语句修改Program2_1而成的程序Q2_2抄写如下:

执行程序Q2_2得到的图形

此处粘帖执行程序Q2_2所得到的图形

39

Q2-3 反复执行程序Program2_2，每次执行该程序时，输入不同的N值，并观察所合成的周

期方波信号。通过观察，你了解的吉伯斯现象的特点是:

1、周期信号的傅里叶级数与GIBBS现象

给定如下两个周期信号:

x1(t)1x2(t)1t2112t

20.20.22

Q2—4 分别手工计算x1（t) 和x2(t) 的傅里叶级数的系数。 信号x1(t) 在其主周期内的数学表达式为： 计算x1(t） 的傅里叶级数的系数的计算过程如下:

通过计算得到的x1（t)的傅里叶级数的系数的数学表达式是: 信号x2（t) 在其主周期内的数学表达式为： 计算x2(t） 的傅里叶级数的系数的计算过程如下:

40

通过计算得到的x1（t）的傅里叶级数的系数的数学表达式是：

用MATLAB帮助你计算出你手工计算的傅里叶级数的系数ak从—10到10共21个系数。从命令窗口上抄写x1（t)的21个系数如下：

从命令窗口上抄写x2（t)的21个系数如下:

Q2-5 仿照程序Program2_1,编写程序Q2_5，以计算x1(t）的傅里叶级数的系数. 程序Q2_5如下：

41

执行程序Q2_5所得到的x1(t)的傅里叶级数的ak从-10到10共21个系数如下:

与你手工计算的ak相比较，是否相同，如有不同，是何原因造成的？ 答：

Q2—6 仿照程序Program2_1，编写程序Q2_6，以计算x2(t) 的傅里叶级数的系数（不绘图）。 程序Q2_6如下：

执行程序Q2_6所得到的x2（t）的傅里叶级数的ak从-10到10共21个系数如下:

与你手工计算的ak相比较,是否相同，如有不同，是何原因造成的？ 答:

Q2—7 仿照程序Program2_2,编写程序Q2_7，计算并绘制出原始信号x1（t） 的波形图，用

有限项级数合成的y1(t) 的波形图，以及x1（t） 的幅度频谱和相位频谱的谱线图。

42

编写程序Q2_7如下：

执行程序Q2_7，输入N = 10所得到的图形如下：

反复执行程序Q2_7,输入不同的N值，观察合成的信号波形中，是否会产生Gibbs现象？为什

么？；

答：

43

2. 连续时间非周期信号的傅里叶变换

给定两个时限信号:

2t1t2,x1(t)1,1t1 x2(t)cos(t)[u(t1)u(t1)]

2t2,1t2Q2—8 利用单位阶跃信号u(t)，将x1(t） 表示成一个数学闭式表达式，并手工绘制x1（t) 和

x2（t) 的时域波形图.

信号x1(t） 的闭式数学表达式为：x1(t) = ：

手工绘制的x1（t）的时域波形图 手工绘制的x2(t)的时域波形图

Q2-9手工计算x1（t） 和x2（t) 的傅里叶变换（如能够用傅里叶变换的性质计算最好）,并手

工绘制出它们的幅度谱和相位谱;

计算x1(t） 的傅里叶变换的过程:

计算得到的x1(t) 的傅里叶变换的数学表达式为：

计算x2（t) 的傅里叶变换的过程：

44

计算得到的x2(t） 的傅里叶变换的数学表达式为:

手工绘制的x1（t）的幅度频谱图 手工绘制的x2（t）的幅度频谱图

Q2-10 编写MATLAB程序Q2_10，能够接受从键盘输入的时域信号表达式，计算并绘制出

信号的时域波形、幅度谱。

程序Q2_10抄写如下

执行程序Q2_10，输入信号x1（t)的数学表达式，得到的信号时域波形、幅度谱和相位谱如下：

45

执行程序Q2_10，输入信号x2（t)的数学表达式,得到的信号时域波形、幅度谱和相位谱如下：

Q2—11 修改程序Q2_10，并以程序Q2_11为文件名存盘，要求能够接受从键盘输入的时域

信号表达式，计算其傅里叶变换，并分别绘制其傅里叶变换的实部、虚部、幅度频谱和相位频谱的图形。

编写的程序Q2_11如下：

选定适当的信号,该信号的时域表达式为:

执行你编写好的MATLAB程序Q2_11，输入你选定的信号的数学表达式，绘制出的该信号的

傅里叶变换的图形如下：

46

Q2-12 修改程序Q2_11,并以Q2_12存盘，要求程序能接受从键盘输入信号的时域表达式，计

算并绘制信号的时域波形、信号的幅度频谱和相位频谱图。

编写的程序Q2_12如下：

Q2—13选择一个时限信号,执行程序Q2_12以验证性质Duality。 选定适当的信号x1(t),该信号的时域表达式为:x1（t） = 执行程序Q2_12,绘制出的该信号的傅里叶变换的图形如下:

47

选定适当的信号x2（t)，要求其数学表达式与x1(t) 的傅里叶变换的数学表达式相同,该信号的

时域表达式为:x2(t) =

再一次执行程序Q2_12，绘制出的信号x2(t）的傅里叶变换的图形如下：

比较这两次执行的结果，简述Duality性质，并说明验证结论。 答:

Q2-14选择一个时限信号,执行程序Q2_12以验证性质Time scaling。

选定适当的信号x3（t),该信号的时域表达式为：x3（t) =: 执行程序Q2_12,绘制出的该信号的傅里叶变换的图形如下：

48

选定信号

x4（t) = x3(2t),该信号的时域表达式为：x4（t） =： ，再一次执行程序Q2_12,绘制出的信号x4（t）的傅里叶变换的图形如下：

比较这两次执行的结果，简述Time scaling性质,并说明验证结论。 答：

Q2—15选择一个时限信号，执行程序Q2_12以验证性质Time shifting。 选定适当的信号

x1（t)，该信号的时域表达式为：x1(t) = ： ，执行程序Q2_12，绘制出的该信号的傅里叶变换的图形如下：

49

选定信号x（ = x（x2(t） = ： ，2t）1t—1),该信号的时域表达式为：

再一次执行程序Q2_12，绘制出的信号x2（t） 的傅里叶变换的图形如下：

比较这两次执行的结果,简述Time shifting性质，并说明验证结论。 答：

Q2-16选择一个合适的信号，执行程序Q2_12以验证性质Frequency shifting（Multiplication by

ej0t）.

选定适当的信号x1(t），该信号的时域表达式为：x1(t) = ： ,执

行程序Q2_12,绘制出的该信号的傅里叶变换的图形如下：

50

选定信号x2（t） = ej0tx1(t）,其中0 = ，再一次执行程序Q2_12,绘制出的信号x2（t)

的傅里叶变换的图形如下：

比较这两次执行的结果，简述Frequency shifting (Multiplication by ej0t)性质，并说明验证结论。 答:

Q2—17:回答如下问题：

1、从信号分解的角度，谈谈你对周期信号的傅里叶级数的理解。

答：

51

2、从信号分解的角度，谈谈你对傅里叶变换及其物理意义的理解，谈谈你对信号频谱概念的理解.

四、实验报告要求

1、按要求完整书写你所编写的全部MATLAB程序

2、详细记录实验过程中的有关信号波形图（存于自带的U盘中），图形要有明确的标题。全部的MATLAB图形应该用打印机打印，然后贴在本实验报告中的相应位置，禁止复印件.

3、实事求是地回答相关问题,严禁抄袭。

本实验完成时间： 年 月 日

52

实验三 连续时间LTI系统的频域分析

一、实验目的

1、掌握系统频率响应特性的概念及其物理意义; 2、掌握系统频率响应特性的计算方法和特性曲线的绘制方法，理解具有不同频率响应特性的滤波器对信号的滤波作用；

3、学习和掌握幅度特性、相位特性以及群延时的物理意义； 4、掌握用MATLAB语言进行系统频响特性分析的方法。

基本要求:掌握LTI连续和离散时间系统的频域数学模型和频域数学模型的MATLAB描述方法,深刻理LTI系统的频率响应特性的物理意义，理解滤波和滤波器的概念，掌握利用MATLAB计算和绘制LTI系统频率响应特性曲线中的编程。

二、实验原理及方法

1 连续时间LTI系统的频率响应

所谓频率特性，也称为频率响应特性，简称频率响应（Frequency response），是指系统在正弦信号激励下的稳态响应随频率变化的情况，包括响应的幅度随频率的变化情况和响应的相位随频率的变化情况两个方面.

x(t)X(j)LTI系统h(t)H(j)y(t)Y(j)连续时间LTI系统的时域及频域分析图 上图中x（t）、y(t）分别为系统的时域激励信号和响应信号，h（t)是系统的单位冲激响应,它们三者之间的关系为：y(t)x(t)*h(t)，由傅里叶变换的时域卷积定理可得到:

53

Y(j)X(j)H(j)

或者: H(j)3。1

Y(j) 3.2

X(j)H(j)为系统的频域数学模型，它实际上就是系统的单位冲激响应h(t）的傅里叶变换。即

 H(j)jth(t)edt 3.3 由于H（j)实际上是系统单位冲激响应h（t）的傅里叶变换，如果h(t）是收敛的，或者说是绝对可积（Absolutly integrabel）的话，那么H(j）一定存在，而且H(j）通常是复数，因此，也可以表示成复数的不同表达形式。在研究系统的频率响应时,更多的是把它表示成极坐标形式:

H(j)H(j)ej() 3。4

上式中,H(j)称为幅度频率相应（Magnitude response),反映信号经过系统之后，信号各频率分量的幅度发生变化的情况，()称为相位特性（Phase response），反映信号经过系统后,信号各频率分量在相位上发生变换的情况。H(j)和()都是频率的函数。

对于一个系统，其频率响应为H（j）,其幅度响应和相位响应分别为H(j)和()，如果作用于系统的信号为x(t)ej0t,则其响应信号为

y(t)H(j0)ej0t3.5

H(j0)ej(0)ej0tH(j0)ej(0t(0))

若输入信号为正弦信号，即x（t) = sin（0t），则系统响应为

y(t)H(j0)sin(0t)|H(j0)|sin(0t(0)) 3。6

可见，系统对某一频率分量的影响表现为两个方面，一是信号的幅度要被H(j)加权，二是信号的相位要被()移相。

由于H(j)和()都是频率的函数，所以,系统对不同频率的频率分量造成的幅度和相位上的影响是不同的。

54

2 LTI系统的群延时

从信号频谱的观点看，信号是由无穷多个不同频率的正弦信号的加权和(Weighted sum）所组成.正如刚才所述,信号经过LTI系统传输与处理时,系统将会对信号中的所有频率分量造成幅度和相位上的不同影响。从相位上来看,系统对各个频率分量造成一定的相位移（Phase shifting),相位移实际上就是延时（Time delay）.群延时（Group delay）的概念能够较好地反映系统对不同频率分量造成的延时。

LTI系统的群延时定义为：

()7

d() 3。d群延时的物理意义：群延时描述的是信号中某一频率分量经过线性时不变系统传输处理后产生的响应信号在时间上造成的延时的时间.

如果系统的相位频率响应特性是线性的，则群延时为常数,也就是说，该系统对于所有的频率分量造成的延时时间都是一样的，因而,系统不会对信号产生相位失真(Phase distortion）.反之，若系统的相位频率响应特性不是线性的,则该系统对于不同频率的频率分量造成的延时时间是不同的，因此，当信号经过系统后，必将产生相位失真。

3 用MATLAB计算系统频率响应

在本实验中，表示系统的方法仍然是用系统函数分子和分母多项式系数行向量来表示。实验中用到的MATLAB函数如下：

[H,w] = freqs（b,a）:b,a分别为连续时间LTI系统的微分方程右边的和左边的系数向量

（Coefficients vector），返回的频率响应在各频率点的样点值(复数）存放在H中,系统默认的样点数目为200点；

Hm = abs(H）：求模数，即进行HmH运算，求得系统的幅度频率响应，返回值存

于Hm之中。

real（H):求H的实部; imag(H):求H的虚部；

phi = atan(-imag(H)。/（real（H）+eps）)：求相位频率相应特性，atan()用来计算反正

切值;或者

phi = angle(H）：求相位频率相应特性；

tao = grpdelay（num,den,w)：计算系统的相位频率响应所对应的群延时. 计算频率响应的函数freqs()的另一种形式是：

H = freqs(b,a，w）:在指定的频率范围内计算系统的频率响应特性。在使用这种形式的

freqs/freqz函数时，要在前面先指定频率变量w的范围。

55

例如在语句H = freqs（b,a,w）之前加上语句:w = 0：2*pi/256:2＊pi。

下面举例说明如何利用上述函数计算并绘制系统频率响应特性曲线的编程方法。 假设给定一个连续时间LTI系统，下面的微分方程描述其输入输出之间的关系

d2y(t)dy(t)32x(t)

dt2dt编写的MATLAB范例程序，绘制系统的幅度响应特性、相位响应特性、频率响应的实部

和频率响应的虚部。程序如下：

% Program3_1

％ This Program is used to compute and draw the plots of the frequency response ％ of a continuous—time system

b = [1］; ％ The coefficient vector of the right side of the differential equation a = ［1 3 2]； ％ The coefficient vector of the left side of the differential equation ［H,w] = freqs（b，a); % Compute the frequency response H Hm = abs(H); ％ Compute the magnitude response Hm phai = angle（H）; % Compute the phase response phai

Hr = real（H）; % Compute the real part of the frequency response Hi = imag（H）； % Compute the imaginary part of the frequency response subplot(221)

plot（w，Hm)， grid on, title(’Magnitude response’)， xlabel（’Frequency in rad/sec') subplot(223） plot(w,phai）， grid on, title（’Phase response’）, xlabel('Frequency in rad/sec’) subplot（222） plot（w，Hr）， grid on， title（’Real part of frequency response')， xlabel('Frequency in rad/sec’） subplot(224）

plot（w,Hi）, grid on, title(’Imaginary part of frequency response'）, xlabel(’Frequency in rad/sec'）

三、实验内容及步骤

实验前，必须首先阅读本实验原理，了解所给的MATLAB相关函数，读懂所给出的全部范例程序。实验开始时，先在计算机上运行这些范例程序，观察所得到的信号的波形图。并结合范例程序所完成的工作,进一步分析程序中各个语句的作用,从而真正理解这些程序。

实验前，一定要针对下面的实验项目做好相应的实验准备工作，包括事先编写好相应的实验程序等事项。

给定三个连续时间LTI系统，它们的微分方程分别为

d2y(t)dy(t)dx(t)125y(t)系统1： Eq。3。1

dtdtdt2系统2： 系统3:

56

dy(t)dx(t)y(t)x(t) Eq。3。2 dtdt

d6y(t)d5y(t)d4y(t)d3y(t)d2y(t)dy(t)1048148306401262y(t)262x(t)65432dtdtdtdtdtdt Eq.3。3

Q3-1 修改程序Program3_1，并以Q3_1存盘，使之能够能够接受键盘方式输入的微分方程系

数向量.并利用该程序计算并绘制由微分方程Eq。3。1、Eq.3。2和Eq。3.3描述的系统的幅度响应特性、相位响应特性、频率响应的实部和频率响应的虚部曲线图。

抄写程序Q3_1如下：

执行程序Q3_1，绘制的系统1的频率响应特性曲线如下：

从系统1的幅度频率响应曲线看，系统1是低通、高通、全通、带通还是带阻滤波器？

57

答：

执行程序Q3_1,绘制的系统2的频率响应特性曲线如下：

从系统2的幅度频率响应曲线看，系统2低通、高通、全通、带通还是带阻滤波器? 答：

执行程序Q3_1，绘制的系统3的频率响应特性曲线如下：

58

从系统3的幅度频率响应曲线看，系统3是低通、高通、全通、带通还是带阻滤波器? 答：

这三个系统的幅度频率响应、相位频率相应、频率响应的实部以及频率响应的虚部分别具有何

种对称关系？请根据傅里叶变换的性质说明为什么会具有这些对称关系？

答：

Q3-2 编写程序Q3_2，使之能够能够接受键盘方式输入的输入信号x（t）的数学表达式，系

统微分方程的系数向量，计算输入信号的幅度频谱，系统的幅度频率响应，系统输出信号y（t）的幅度频谱，系统的单位冲激响应h（t），并按照下面的图Q3-2的布局，绘制出各个信号的时域和频域图形.

59

图Q3-2

你编写的程序Q3_2抄写如下：

执行程序Q3_2，输入信号x(t) = sin（t) + sin（8t），输入由Eq.3.3描述的系统。得到的图形如

下：

此处粘帖执行程序Q3_2所得到的图形

60

请手工绘制出信号x(t） = sin(t） + sin(8t) 的幅度频谱图如下：

你手工绘制的信号x（t) = sin(t） + sin(8t) 的幅度频谱图与执行程序Q3_2得到的x（t) = sin

（t） + sin（8t） 的幅度频谱图是否相同？如不同，是何原因造成的？

答：

执行程序Q3_2得到的x(t) = sin(t） + sin（8t） 的幅度频谱图实际上是另外一个信号x1(t）的

幅度频谱，这个信号的时域数学表达式为 x1（t） = 请利用傅里叶变换的相关性质计算并绘制信号x1（t)的幅度频谱图.

计算过程:

手工绘制的x1(t） 的幅度频谱图如下：

结合所学的有关滤波的知识，根据上面所得到的信号的时域和频域图形，请从时域和频域两个

61

方面解释滤波的概念。

答：

Q3-3 编写程序Q3_3，能够接受从键盘输入的系统微分方程系数向量,并分别绘制所给三个系

统的群延时曲线图。

抄写程序Q3_3如下：

系统Eq。3。1的群延时曲线图 系统Eq.3。3的群延时曲线图

根据上面的群延时曲线图可以看出，对系统Eq。3。1，当频率为5弧度/秒时，群延时为 秒，

当频率为10弧度/秒时，群延时为 秒，如何解释这两个群延时时间？

62

根据上面的群延时曲线图,说明这两个系统是否会造成对信号的相位失真？为什么？

从系统Eq。3。3的群延时曲线图中可以看出,当信号的频率为1弧度/秒时，系统Eq.3。3对这

一频率的信号的延时是 秒。所以，执行程序Q3_2时，当作用于系统Eq.3.3的输入信号为x(t) = sin（t） + sin（8t)时,其输出信号y（t)的数学表达式为：

四、实验报告要求

1、按要求完整书写你所编写的全部MATLAB程序

2、详细记录实验过程中的有关信号波形图（存于自带的U盘中），图形要有明确的标题.全部的MATLAB图形应该用打印机打印，然后贴在本实验报告中的相应位置，禁止复印件.

3、实事求是地回答相关问题，严禁抄袭.

本实验完成时间： 年 月 日

实验四 信号抽样与调制解调

一、实验目的

1、进一步理解信号的抽样及抽样定理； 2、进一步掌握抽样信号的频谱分析；

3、掌握和理解信号抽样以及信号重建的原理； 4、掌握傅里叶变换在信号调制与解调中的应用.

基本要求：掌握并理解“抽样”的概念,理解抽样信号的频谱特征。深刻理解抽样定理及其重要意义。一般理解信号重建的物理过程以及内插公式所描述的信号重建原理.理解频率混叠的概念。理解调制与解调的基本概念，理解信号调制过程中的频谱搬移。掌握利用MATLAB仿真正弦幅度调制与解调的方法.

二、实验原理及方法 1、信号的抽样及抽样定理

63

抽样（Sampling）,就是从连续时间信号中抽取一系列的信号样本，从而，得到一个离散时间序列（Discrete—time sequence），这个离散序列经量化(Quantize)后，就成为所谓的数字信号（Digital Signal）。今天，很多信号在传输与处理时，都是采用数字系统（Digital system）进行的，但是，数字系统只能处理数字信号，不能直接处理连续时间信号或模拟信号（Analog signal）。为了能够处理模拟信号,必须先将模拟信号进行抽样，使之成为数字信号,然后才能使用数字系统进行传输与处理。所以，抽样是将连续时间信号转换成离散时间信号必要过程。模拟信号经抽样、量化、传输和处理之后，其结果仍然是一个数字信号，为了恢复原始连续时间信号,还需要将数字信号经过所谓的重建（Reconstruction）和平滑滤波（Smoothing）。图4。1展示了信号抽样与信号重建的整个过程。 Antialiasing filter xa(t)Sampler/ Holder A/D convertor Digital Processor D/A convertor Antialiasing filter y(t)p(t)图4.1 模拟信号的数字处理过程

图4.2给出了信号理想抽样的原理图：

x(t)p(t)xs(t)X(j)mm(a)

(b)

图4.2 (a) 抽样原理图，(b) 带限信号的频谱

上图中,假设连续时间信号是一个带限信号（Bandlimited Signal），其频率范围为m~m,抽样脉冲为理想单位冲激串(Unit Impulse Train），其数学表达式为：

p(t)(tnTs) 4。1

由图可见，模拟信号x（t)经抽样后，得到已抽样信号(Sampled Signal）xs(t）,且：

xs(t)x(t)p(t) 4.2

将p（t）的数学表达式代入上式得到：

xs(t)x(nTs)(tnTs) 4.3

显然，已抽样信号xs（t） 也是一个冲激串，只是这个冲激串的冲激强度被x(nTs） 加权了。

从频域上来看,p（t） 的频谱也是冲激序列，且为：

F{p(t)}s(ns) 4。4

根据傅里叶变换的频域卷积定理，时域两个信号相乘，对应的积的傅里叶变换等于这两个信号的傅里叶变换之间的卷积。所以，已抽样信号xs（t)的傅里叶变换为:

1Xs(j)TsnX(j(n)) 4。5

s表达式4.5告诉我们，如果信号x(t)的傅里叶变换为X(j），则已抽样信号xs(t) 的傅里叶变换Xs（j)等于无穷多个加权的移位的X(j）之和，或者说，已抽样信号的频谱等于原连续时间信号的频谱以抽样频率s为周期进行周期复制的结果。如图4。3所示:

x(t)X(j)tp(t)MMst1/TssP(j)Tsxs(t)sXs(j)t图4.3 信号抽样及其频谱图



由图可见，如果抽样频率不小于信号带宽的2倍时，xs（t) 的频谱中,X（j)的各个复制品之间没有混叠（Aliasing)，因此，可以用一个理想低通滤波器来恢复原始信号。由抽样信号恢复原来的原始信号的过程称为信号的重建( Reconstruction ）。反之，如果抽样频率小于信号带宽的2倍时，xs(t) 的频谱中，X(j）的各个复制品之间的距离（也就是s）太近，所以必将造成频谱之间的混叠，在这种情况下,是无论如何也无法恢复出原来的连续时间信号的.

由此，我们得出下面的结论:当抽样频率 s 〉 2M 时，将原连续时间信号x(t)抽样而得到的离散时间序列x[n］可以唯一地代表原连续时间信号，或者说，原连续时间信号x(t)可以完

65

全由x［n]唯一地恢复.

以上讨论的是理想抽样的情形，由于理想冲激串是无法实现的，因此，这种理想抽样过程，只能用来在理论上进行抽样过程的分析。在实际抽样中，抽样往往是用一个A/D转换器实现的。一片A/D转换芯片包含有抽样保持电路和量化器。模拟信号经过A/D转换器后，A/D转换器的输出信号就是一个真正意义上的离散时间信号，而不再是冲激串了.

A/D转换器的示意图如图4.4所示。

Sampler x(t)p(t)TsHolder Quantizer x[n]图4.4 A/D转换器示意图

上述的实际抽样过程，很容易用简单的数学公式来描述。设连续时间信号用x(t）表示,抽样周期（Sampling Period）为Ts，抽样频率（Sampling Frequency）为s，则已抽样信号的数学表达式为

x[n]x(t)tnTsx(nTs) 4。6

在MATLAB中,对信号抽样的仿真，实际上就是完成式4。6的计算。下面给出一个例题和相应的范例程序，来实现信号抽样的仿真运算.

例题4.1 设连续时间信号为一个正弦信号 x(t） = cos(0.5πt)，抽样周期为Ts = 1/4秒，编程

序绘制信号x（t)和已抽样信号x[n］的波形图。

范例程序Sampling如下：

％ Sampling clear, close all， t = 0:0。01:10;

Ts = 1/4； % Sampling period

n = 0:Ts：10; % Make the time variable to be discrete x = cos（0.5*pi＊t）;

xn = cos（0.5＊pi*n); ％ Sampling subplot（221)

plot（t，x）, title('A continuous-time signal x(t）'）, xlabel（'Time t'） subplot（222）

stem(n，xn，’.')， title（'The sampled version x[n] of x（t）’)， xlabel（'Time index n’) 执行该程序后,得到的波形图如图4.5所示。

66

图4.5 连续时间信号及其抽样后的离散时间序列

在这个范例程序中，先将连续时间t进行离散化，使之成为以Ts = 1/4秒的离散时间n，然后，将n代入到信号x（t） 的数学表达式中计算，就完成了抽样过程，且得到了抽样后的离散时间序列x[n].

2、 信号抽样过程中的频谱混叠

为了能够观察到已抽样信号的频谱是否会存在混叠现象,或者混叠程度有多么严重，有必要计算并绘制出已抽样信号的傅里叶变换.

根据式4.5可计算出已抽样信号的频谱。下面给出的范例程序Program4_1就是按照式4.5进行计算的。其中，主要利用了一个for循环程序完成周期延拓运算。

％ Program4_1 clear, close all，

tmax = 4; dt = 0。01; t = 0：dt:tmax； Ts = 1/10; ws = 2*pi/Ts;

w0 = 20*pi； dw = 0.1; w = -w0:dw：w0; n = 0：1：tmax/Ts; x = exp（-4*t）。＊u(t)； xn = exp（—4＊n＊Ts）； subplot（221）

plot（t,x), title(’A continuous—time signal x(t)’)， xlabel(’Time t’), axis(［0，tmax，0，1］)， grid on subplot(223） stem（n,xn,'。’）， title（’The sampled version x[n］ of x(t）'), xlabel（'Time index n'), axis（［0，tmax/Ts，0,1］）, grid on Xa = x*exp(—j＊t'*w)＊dt; X = 0;

for k = —8:8;

X = X + x＊exp（—j*t’＊（w—k＊ws)）*dt；

67

end

subplot(222)

plot（w,abs（Xa）)

title（'Magnitude spectrum of x(t)’）, grid on axis（［-60,60,0,1。8＊max(abs(Xa）)]） subplot（224) plot(w，abs(X））

title（'Magnitude spectrum of x［n］')， xlabel('Frequency in radians/s'),grid on axis（[-60,60,0,1.8＊max(abs(Xa））])

本程序可以用来观察在不同的抽样频率条件下,已抽样信号的频谱的混叠程度，从而更加牢固地理解抽样定理。但是,提请注意的是,在for循环程序段中，计算已抽样信号的频谱X 时,没有乘以系数1/Ts，是为了便于比较X与Xa之间的区别，从而方便观察频谱的混叠程度。另外,程序中的时间步长dt的选择应该与抽样周期Ts保持一定的比例关系，建议Ts不应小于10dt，否则，计算得到的已抽样信号的频谱将出现错误。

3、 信号重建

如果满足抽样定理，那么,我们就可以唯一地由已抽样信号x[n］ 恢复出原连续时间信号x(t）。在理想情况下，可以将离散时间序列通过一个理想低通滤波器，图4。6给出了理想情况下信号重建的原理示意图。

x(t)xp(t)Ideal Lowpass Filter xr(t)p(t)图4.6 信号重建原理图

理想低通滤波器也称重建滤波器,它的单位冲激响应

h(t)cTsin(ct) 4。7

ct已抽样信号xp（t)的数学表达式为:xp(t)表达式，我们有

x(nT)(tnT)，根据系统输入输出的卷积

xr(t)xp(t)h(t) 4。8

将xp(t）代入式4.8，有

xr(t)nx(nT)cTsin(c(tnT)) 4。9

c(tnT)这个公式称为内插公式（Interpolation Formula)，该公式的推导详见教材，请注意复习有

68

关内容。须提请注意的是，这里的内插公式是基于重建滤波器为理想低通滤波器的。若重建滤波器不是理想低通滤波器,则不能用这个内插公式。理想低通滤波器的频率响应特性曲线和其单位冲激响应曲线如图4.7所示。

H(j)Th(t)Tccct图4.7 理想低通滤波器的幅度频率响应和单位冲激响应

范例程序程序Program4_2就是根据这个内插公式来重构原始信号。本程序已经做了较为详细的注释,请结合教材中的内插公式仔细阅读本程序，然后执行，以掌握和理解信号重建的基本原理。范例程序Program4_2如下。

％ Program4_2

％ Signal sampling and reconstruction

％ The original signal is the raised cosin signal： x（t） = ［1+cos（pi＊t)］.*［u

（t+1)—u(t-1）］.

clear; close all,

wm = 2＊pi； ％ The highest frequency of x(t) a = input（’Type in the frequency rate ws/wm=:’）； % ws is the sampling frequency wc = wm; ％ The cutoff frequency of the ideal lowpass filter t0 = 2; t = -t0：0.01:t0;

x = （1+cos（pi＊t）)。＊（u（t+1）—u(t-1))；

subplot（221）; % Plot the original signal x(t) plot（t，x）; grid on， axis（［—2，2,—0。5，2。5］）； title('Original signal x(t）’）;xlabel(’Time t’）;

ws = a＊wm; ％ Sampling frequency Ts = 2＊pi/ws； ％ Sampling period

N = fix（t0/Ts）; ％ Determine the number of samplers n = -N:N；

nTs = n*Ts； ％ The discrete time variable

xs = (1+cos（pi＊nTs))。*（u（nTs+1)—u（nTs—1）)； % The sampled version of

x（t)

subplot(2,2，2) % Plot xs stem（n,xs,'.'）； xlabel(’Time index n’); grid on, title（'Sampled version x［n]’）; xr = zeros(1,length(t））; ％ Specify a memory to save the reconstructed signal L = length（—N：N)； xa = xr;

figure(2）; % Open a new figure window to see the demo of signal reconstruction stem（nTs，xs,'.’); xlabel（'Time index n’); grid on；hold on

69

for i = 1:L m = (L—1)/2+1-i；

xa = Ts*(wc）＊xs（i)*sinc((wc)*(t+m＊Ts）/pi）/pi; plot（t，xa，'b：'）；axis（[—2，2,-0。5,2。5］); hold on pause

xr = xr+xa； % Interpolation end

plot(t,xr,'r’）; axis(［-2,2,—0.5，2.5]）； hold on figure(1）; subplot（223)

plot（t，xr,’r')；axis（[-2,2,—0.5,2。5］); xlabel（'Time t')；grid on

title('Reconstructed signal xr（t)');

% Compute the error between the reconstructed signal and the original signal error = abs（xr—x)； subplot(2，2,4）

plot（t，error）;grid on title（’Error’）；xlabel(’Time t’)

4、 调制与解调

在通信系统(Communication system）中，信号在传输之前，往往需要先对它进行调制（Modulation），然后才能发射出去.在接收端，还要进行解调（Demodulation）,才能恢复原信号。在实际应用中，有多种调制方法，在信号与系统中,仅介绍了模拟调制中的正弦幅度调制（Sinusoidal amplitude modulation）。正弦幅度调制就是利用高频正弦信号的幅度携带调制信号x（t），也就是使高频正弦信号的幅度随调制信号的变化而变化。正弦调制的解调分为同步解调（Synchronous demodulation）和异步解调(Asychronous demodulation)，调制与解调的原理框图如图4.8所示。

图中，需要传输的信号称为调制信号（Modulating signal），频率为c的正弦信号称为载波(Carrier），c称为载频,调制器的输出信号称为已调信号(Modulated signal）。

正弦幅度调制的基本原理，就是将调制信号与一个高频正弦载波相乘,从而将调制信号的频

谱搬移到较高的频段上，以利于发射传输。

70

x(t)(a)

y(t)x(t)cos(ct)Lowpass filter x(t)c(t)cos(ct)c(t)cos(ct)

下面,我们回顾一下调制与解调过程中的时域和频域的有关情况。从时域上看，已调信号的数学表达式为

(b)

图4.8 正弦幅度调制与解调 (a) 调制 (b) 同步解调

y(t)x(t)cos(ct) 4.10

调制信号x（t）、载波c(t)和已调信号y(t）的波形如图4。9所示

图4.9 正弦幅度调制中信号的波形

从频域上看，假设调制信号是一个带限信号，其频谱用X（j) 表示,而正弦载波cos（ct) 的频谱C（j) 由两个冲激构成，即

C(j)[(c)(c)] 4.11

根据傅里叶变换的频域卷积定理可知,已调信号的频谱为

71

Y(j)即

1X(j)C(j)] 4。12 21Y(j)[X(j(c))X(j(c)] 4.13

2式4.13说明，已调信号的频谱由两个移位的X（j)构成，位移量为 ±c。图4。10示出了调制过程中各信号的频谱图。

从已调信号的频谱上看，我们发现，调制信号为低通信号（Lowpass signal),其带宽(Bandwidth)为M，而已调信号则变成了一个带宽为2M的带通信号（Bandpass signal）.这表明，通过调制,信号在传输过程中，与不调制而直接传输相比，需要占居更宽的信道（Channel）带宽.

从图4。8可以看出,同步解调的原理类似于调制原理，只是在乘法器(Multiplier）后面接了一个低通滤波器,请同学们参看教材中的相关内容，自行分析同步解调的原理,并绘制出同步解调过程中各信号的频谱图.

X(j)MM5、 通信系统中的调制与解调仿真

C(j)本实验室利用MATALAB对通信系统中的调制与解

cc调、滤波等进行仿真。具体方法简述如下：

1、在命令窗口键入simulink然后按回车键，这时屏幕上将

Y(j)出现一个新界面： 1/2Simulink Library Browser 界面.

2、新建仿真系统图：

第一步:单击左上角的新建按钮，将在屏幕右部出现建

图4.10 调制过程中各信号的频谱图

模窗口;

第二步：建立仿真系统的系统函数。单击continuous

模块，选择系统函数功能框，将它拖放到空白图面上，然后双击该功能框，又出现参数选择对话框,在该框中设定好仿真系统(滤波器）的参数。如果仿真模型中需要多个系统函数功能框图,可重复第二步;

第三步：选择信号源。单击Sources模块，选择需要的信号源拖放到模型图中，然后双击已设定适当的参数;

第四步：选择信号之间的运算单元,如加法器,乘法器（调制器）等。单击Math Operations 模块，选择所需的运算单元，拖放到模型图中并双击加以设定；

第五步：选择显示器，通常选择示波器。单击Sinks模块，并拖放到模型图中;双击加以设定；

第六步：将模型图中的所有元件调整好位置,然后进行连接。

3、选择仿真时间，单击模型图上部的Simulation菜单，选中Simulation Parameters子菜单设定方针的起止时间，如不设定,则系统默认的起止时间为0～10s。

72

4、单击运行按钮开始仿真，双击示波器即可看到仿真结果。

三、实验内容及步骤

实验前，必须首先阅读本实验原理，了解所给的MATLAB相关函数,读懂所给出的全部范例程序。实验开始时，先在计算机上运行这些范例程序，观察所得到的信号的波形图.并结合范例程序所完成的工作,进一步分析程序中各个语句的作用，从而真正理解这些程序的编程算法。

实验前，一定要针对下面的实验项目做好相应的实验准备工作，包括事先编写好相应的实验程序等事项.

Q4-1 给范例程序Program4_1加注释。

Q4—2 范例程序Program4_1中的连续时间信号x(t） 是什么信号？它的数学表达式为：

Q4-3 在1/2—1/10之间选择若干个不同Ts值,反复执行执行范例程序Program4_1，保存执

行程序所得到的图形。

Ts = 1/2时的信号时域波形和频谱图

Ts = 1/4时的信号时域波形和频谱图

73

Ts = 1/8时的信号时域波形和频谱图

根据上面的三幅图形，作一个关于抽样频率是怎样影响已抽样信号频谱的小结。答:

74

程序Program4_1中的连续信号是否是带限信号？如果不是带限信号，是否可以选择一个抽样

频率能够完全消除已抽样信号中的频谱的混叠？如果不是带限信号，那么，这个连续时间信号在抽样前必须滤波,请你选择一个较为合适的频率作为防混叠滤波器的截止频率，你选择的这个截止频率是多少？说明你的理由.

答：

Q4—4 请手工计算升余弦信号x（t） = [1+cos（pi＊t)].＊[u（t+1）—u(t—1）] 的傅里叶变

换的数学表达式，手工绘制其幅度频谱图.

计算过程:

手工绘制的升余弦信号x(t) = ［1+cos(pi*t）］.＊［u（t+1)-u(t-1)］ 的幅度频谱图

75

从上图的幅度频谱上看，升余弦信号是否是带限信号？能否近似将它看作是一个带限信号？如

果可以，那么,估计信号的最高频率大约是多少？

答：

Q4—5 阅读范例程序答：

Program4_2,在这个程序中，选择的信号的最高频率是多少?这个频率

选择得是否恰当？为什么？

Q4-6 在1—8之间选择抽样频率与信号最高频率之比，即程序Program4_2中的a值，反复

执行范例程序Program4_2，观察重建信号与原信号之间的误差，通过对误差的分析，说明对于带限信号而言,抽样频率越高，则频谱混叠是否越小？

答:

Q4—7 图Q4-7是一个RLC串联电路，在有些场合，可以把它用

来作为一个滤波器使用,如果选择不同的位置的信号作为输出信号，该电路具有不同的滤波特性。该电路的输出信号分别为y1（t)、y2（t)和y3（t）时，电路的输入输出微分方程和频率响应的数学表达式分别形如：

•R ••y3(t)••选择y1（t)为输出信号时（可将电路看成系统1）

x(t)L y2(t)••76 C

y1(t)•••图Q4-7 RLC串联电路

d2y1(t)dy1(t)222y(t)x(t) 微分方程为 nn1n2dtdt2n频率响应为 H(j) 22(j)2n(j)n选择y2（t）为输出信号时（可将电路看成系统2）

d2y2(t)dy2(t)d2x(t)22nny2(t)微分方程为 22dtdtdt(j)2频率响应为 H(j) 22(j)2n(j)n选择y3(t）为输出信号时（可将电路看成系统3）

d2y3(t)dy3(t)dx(t)22y(t)2微分方程为 nn3ndtdtdt2频率响应为 H(j)2n(j) 2(j)22n(j)n请写出上面的微分方程和频率响应表达式中的参数、n与R、L、C之间的数学关系.

Q4—8 编写程序Q4_8,能够接受从键盘输入的、n之值，计算并在同一个图形窗口的三个

子图中绘制出这三个频率响应特性曲线，要求每个子图有标题，绘制的频率范围为0—40弧度/秒.图形布置如图Q4-8所示。

77

图Q4-8 图形布置(zeta = ， wn = n）

抄写程序Q4_8如下：

执行程序Q4_8，输入zeta = 0.7，wn = 15，在图形中的空白处,标上zeta 和wn之值，如图

Q4-8所示.保存所得到的图形如下。 zeta = 0。7,wn = 15时的频率响应曲线图

78

根据上面的图形，请说明系统1、系统2和系统3分别具有何种滤波特性？ 答:

固定zeta = 0。7，在2-30之间选择不同的wn值,反复执行程序Q4_8，保存zeta = 0。7，

wn = 5和zeta = 0。7，wn = 20所得到的两幅图形.根据执行程序所得到的系统频率响应的形状，说明wn的不同取值分别对系统1、系统2和系统3的滤波特性（从通频带的带宽、过渡带宽和截止频率等方面作说明）的影响.

zeta = 0。7,wn = 5时的频率响应曲线图

79

zeta = 0。7，wn = 20时的频率响应曲线图

答：

固定wn = 15，在0.2—1之间选择不同的zeta值,反复执行程序Q4_8，保存zeta = 0.4,wn =

15和zeta = 0.8，wn = 15所得到的两幅图形。根据执行程序所得到的系统频率响应的形状,说明zeta的不同取值分别对系统1、系统2和系统3的滤波特性的影响。

zeta = 0.4，wn = 15时的频率响应曲线图

80

zeta = 0。8，wn = 15时的频率响应曲线图

答:

Q4—9 调制与解调仿真实验.设调制信号为单频正弦信号x(t） = sin(t)，其角频率为1 rad/s，

载波为c（t) = cos（10t）,载频为10rad/s。

请按下面的图Q4-9建立仿真模型图：

图有三个信号源，其中：

Sin Wave为调制信号源即调制信号，可设定其频率为1 rad/s ；

Sin Wave1为载波信号，可设定其频率为30rad/s，Band—Limited White Noise为带限白噪声干扰信号，其频率可认为远大于1 rad/s;

Product和Product1分别为调制器和解调器,完成信号的乘法运算。

81

图Q4-9 信号的调制与解调仿真模型图

第一个乘法器之后的Transfer Fun是一个带通滤波器，数学模型可给定为:

H(j)2n(j) 2(j)22n(j)n可以用Q4-7的RLC串联电路构成的系统3实现。

根据调制信号和载波的频率，以及实验结果，你认为图Q4—9中的带通滤波器的参数  和 

n

应该选择为：

wn =  =

第二个乘法器之后的Transfer Fun是一个低通滤波器，设定其系统函数为：

2n H(j)22(j)2n(j)n可以用Q4—7的RLC串联电路构成的系统1实现。

根据调制信号和载波的频率，以及实验结果，你认为图Q4-9中的低通滤波器的参数  和 

n

应该选择为：

wn =  =

四、实验报告要求

1、按要求完整书写你所编写的全部MATLAB程序

2、详细记录实验过程中的有关信号波形图(存于自带的U盘中)，图形要有明确的标题。全部的MATLAB图形应该用打印机打印,然后贴在本实验报告中的相应位置，禁止复印件。

3、实事求是地回答相关问题，严禁抄袭。

82

本实验完成时间： 年 月 日

实验五 连续时间LTI系统的复频域分析

一、实验目的

1、掌握拉普拉斯变换的物理意义、基本性质及应用；

2、掌握用拉普拉斯变换求解连续时间LTI系统的时域响应；

3、掌握系统函数的概念，掌握系统函数的零、极点分布（零、极点图）与系统的稳定性、时域特性等之间的相互关系；

4、掌握用MATLAB对系统进行变换域分析的常用函数及编程方法。

基本要求：掌握拉普拉斯变换及其基本性质，掌握应用拉普拉斯变换求解系统的微分方程,能够自己编写程序完成对系统时域响应的求解。掌握并理解系统函数的概念，掌握系统函数零极点与系统时域和频域特性之间的关系，能够编写程序完成对系统的一些主要特性如稳定性、因果性等的分析.

二、实验原理及方法

83

1、连续时间LTI系统的复频域描述

拉普拉斯变换（The Laplace transform）主要用于系统分析。描述系统的另一种数学模型就是建立在拉普拉斯变换基础上的“系统函数（System Function）\"——H(s）：

H(s)Y(s)系统冲击响应的拉氏变换Ly(t) 5.1

X(s)系统激励信号的拉氏变换Lx(t)系统函数H(s)的实质就是系统单位冲激响应（Impulse Response）h(t)的拉普拉斯变换。因此,系统函数也可以定义为:

 H(s)sth(t)edt 5。2 所以，系统函数H(s)的一些特点是和系统的时域响应h(t)的特点相对应的。在教材中,我们求系统函数的方法,除了按照拉氏变换的定义式的方法之外，更常用的是根据描述系统的线性常系数微分方程(Linear Constant-Coefficient Defrential Equation)，经过拉氏变换之后得到系统函数H(s)。

假设描述一个连续时间LTI系统的线性常系数微分方程为：

dky(t)Mdkx(t) 5。3 akbkkkdtdtk0k0N对式4.3两边做拉普拉斯变换，则有

kasY(s)bskkX(s)

kk0k0NM即 H(s)Y(s)X(s)bskMkak0k0N 5.4

ksk式5。4告诉我们，对于一个能够用线性常系数微分方程描述的连续时间LTI系统，它的系统函数是一个关于复变量s的有理多项式的分式，其分子和分母的多项式系数与系统微分方程左右两端的系数是对应的。根据这一特点，可以很容易的根据微分方程写出系统函数表达式，或者根据系统函数表达式写出系统的微分方程。

系统函数H(s)大多数情况下是复变函数，因此，H(s)可以有多种表示形式:

84

1、直角坐标形式： H(s)Re(s)jIm(s)

k(szj)2、零极点形式： H(s)j1M(sp)ii1N

3、部分分式和形式: H(s)Ak（假设系统的N>M，且无重极点） k0sskN根据我们所要分析的问题的不同,可以采用不同形式的系统函数H(s)表达式.

在MATLAB中，表达系统函数H(s)的方法是给出系统函数的分子多项式和分母多项式的系数向量.由于系统函数的分子和分母的多项式系数与系统微分方程左右两端的系数是对应的，

因此，用MATLAB表示系统函数，就是用系统函数的两个系数向量来表示.

应用拉普拉斯变换分析系统的主要内容有：

1、分析系统的稳定性; 2、分析系统的频率响应.

分析方法主要是通过绘制出系统函数的零极点分布图，根据零极点分布情况，判断系统的稳定性。

MATLAB中有相应的复频域分析函数，下面简要介绍如下：

［z，p,k] = tf2zp(num，den）：求系统函数的零极点，返回值z为零点行向量，p为极

点行向量,k为系统传递函数的零极点形式的增益。num为系统函数分子多项式的系数向量，den为系统函数分母多项式系数向量。

H = freqs（num,den,w）:计算由num，den描述的系统的频率响应特性曲线.返回值H

为频率向量规定的范围内的频率响应向量值。如果不带返回值H,则执行此函数后，将直接在屏幕上给出系统的对数频率响应曲线(包括幅频特性取向和相频特性曲线)。

［x，y] = meshgrid(x1，y1）：用来产生绘制平面图的区域，由x1，y1来确定具体的区

域范围，由此产生s平面区域。

meshgrid(x,y,fs)：绘制系统函数的零极点曲面图. H = impulse(num，den）：求系统的单位冲激响应，不带返回值，则直接绘制响应曲线，

带返回值则将冲激响应值存于向量h之中。

2、系统函数的零极点分布图

系统函数的零极点图（Zero—pole diagram）能够直观地表示系统的零点和极点在s平面上的位置，从而比较容易分析系统函数的收敛域（Regin of convergence）和稳定性(stablity）。

下面给出一个用于绘制连续时间LTI系统的零极点图的扩展函数splane(num,den）:

% splane

％ This function is used to draw the zero—pole plot in the s—plane

85

function splane（num,den） p = roots（den）； ％ Determine the poles q = roots（num）; % Determine the zeros p = p’; q = q’; x = max(abs（[p q]））； ％ Determine the range of real—axis x = x+1;

y = x; ％ Determine the range of imaginary-axis plot（［—x x］,[0 0]，'：'）；hold on; ％ Draw the real-axis plot(［0 0］，［—y y］，’：')；hold on; % Draw the imaginary-axis plot(real（p)，imag（p)，'x');hold on； % Draw the poles plot(real(q)，imag（q），’o'）；hold on； % Draw the zeros title(’zero-pole plot');

xlabel('Real Part’)；ylabel(’Imaginal Part’） axis([—x x -y y]）; % Determine the display-range 对于一个连续时间LTI系统,它的全部特性包括稳定性、因果性（Causality）和它具有何种滤波特性（Frequency-domain aspect）等完全由它的零极点在s平面上的位置所决定。

3、拉普拉斯变换与傅里叶变换之间的关系

根据课堂上所学的知识可知,拉普拉斯变换与傅里叶变换之间的关系可表述为：傅里叶变换

是信号在虚轴上的拉普拉斯变换,也可用下面的数学表达式表示 H(j)H(s)sj 5。5

上式表明，给定一个信号h(t）,如果它的拉普拉斯变换存在的话，它的傅里叶变换不一定存在，只有当它的拉普拉斯变换的收敛域包括了整个虚轴,则表明其傅里叶变换是存在的。下面的程序可以以图形的方式，表现拉普拉斯变换与傅里叶变换的这种关系。

% Relation_ft_lt

％ This program is used to observe the relationship between the Fourier transform ％ and the Laplace transform of a rectangular pulse. clear, close all， a = -0：0。1:5； b = —20:0.1：20；

［a， b］ = meshgrid （a， b）； c = a+i＊b; %确定绘图区域 c = (1—exp (-2＊ （c+eps)）)。/ (c+eps); c = abs （c); ％计算拉普拉斯变换 subplot (211）

mesh (a，b，c); ％绘制曲面图 surf (a，b,c）； view （-60,20) %调整观察视角

86

axis ([-0,5，—20，20，0,2]）;

title （’The Laplace transform of the rectangular pulse’); w = —20：0。1：20；

Fw = (2＊sin（w+eps)。*exp（i＊(w+eps）)）。/(w+eps）； subplot (212)； plot(w，abs（Fw))

title （'The Fourier transform of the rectangular pulse'） xlabel ('frequence w'）

上面的程序不要求完全读懂,重点是能够从所得到的图形中，观察拉和理解普拉斯变换与傅里叶变换之间的相互关系就行。

4、系统函数的极点分布与系统的稳定性和因果性之间的关系

一个稳定的LTI系统,它的单位冲激响应h（t)满足绝对可积条件，即



h(t)dt 5.6

同时,我们还应该记得，一个信号的傅里叶变换的存在条件就是这个信号满足绝对可积条件,所以，如果系统是稳定的话，那么，该系统的频率响应也必然是存在的。又根据傅里叶变换与拉普拉斯变换之间的关系，可进一步推理出，稳定的系统，其系统函数的收敛域必然包括虚轴。稳定的因果系统,其系统函数的全部极点一定位于s平面的左半平面.

所以，对于一个给定的LTI系统，它的稳定性、因果性完全能够从它的零极点分布图上直观地看出。

例题5-1：已知一个因果的LTI系统的微分方程为

d6y(t)d5y(t)d4y(t)d3y(t)d2y(t)dy(t)1048148306401262y(t)262x(t)65432dtdtdtdtdtdt

编写程序，绘制出系统的零极点分布图，并说明它的稳定性如何。

解：这是一个高阶系统，显然手工计算它的极点是很困难的。可以利用前面给出的扩展函数

splane(),来绘制系统的零极点分布图.范例程序如下：

％ Program5_1

％ This program plots the zero—pole diagram of an LTI system described ％ by the linear constant-coefficient differential

equation

clear, close all， b = 262;

a = ［1 10 48 148 306 401 262］; subplot （221） splane (b，a)

87

title （'The zero—pole diagram’)

执行该程序后，得到系统的零极点分布图如图5。1所示.由于已知该系统是因果系统,从零极点分布图上看，它的全部极点都位于s平面的左半平面上，所以系统是稳定的.

然后，直接在命令窗口键入

>> roots（a)

回车后,就得到系统的极点为： 图5.1

ans =

—0.5707 + 2.4716i -0。5707 - 2。4716i —2。7378 + 0。0956i —2.7378 - 0。0956i —1.6915 + 1。6014i —1.6915 - 1.6014i

若题目中没有说明该系统是否是因果的，则需要做详细的分析。从零极点分布图上可以看出，该系统可能的收敛域共有四种可能，另外三种可能如下：

（a) 收敛域为Re｛s}〈 -2.7378，此种情况说明，该系统是一个反因果系统(Anticausal system），由于收敛域不包含虚轴，故此系统是不稳定的。

（b）、（c)收敛域为 —2.7378 < Re{s} 〈 -1。6915和—1。6915 < Re｛s｝ < -0。5707，此两种情况说明该系统是一个单位冲激响应为双边信号的非因果系统，收敛域仍不包含虚轴，所以,系统是不稳定的.

总之，系统的稳定性主要取决于系统函数的收敛域是否包含整个虚轴，而系统的因果性则取决于系统极点位置的分布.

需要特别强调的是，MATLAB

总是把由分子和分母多项式表示任何系统都当作是因果系

统。所以，利用impulse ()函数求得的单位冲激响应总是因果信号。

5、系统函数的零极点分布与系统的滤波特性

系统具有何种滤波特性,主要取决于系统的零极点所处的位置。没有零点的系统,通常是一个低通滤波器.

例题5-2 已知一个系统的系统函数为

H(s)

1 s1图5.2

88

显然，这是一个一阶系统，无零点。为了确定该系统具有何种滤波特性,需要把系统的频率响应

特性曲线绘制出来加以判断。借助实验三中的范例程序Program3_1，可以绘制系统的频率特性曲线如图5。2所示。

通过编程，可以将系统的零极点分布图和系统的频率响应特性以及系统的单位冲激响应特性绘制在一个图形窗口的各个子图中，这样便于观察系统的零极点分布情况与系统的时域和频域之间的关系。如图5。3所示

图5.3 系统的零极点分布与系统的单位冲激响应、频率相应

6、拉普拉斯逆变换的计算

我们已经知道，直接用拉普拉斯逆变换（Inverse transform）的定义公式计算逆变换是很困难的，通常的计算拉普拉斯逆变换的方法是长除法(Long division）和部分分式分解法(Partial fraction expension）。MATLAB的内部函数residue(）可以帮助我们完成拉普拉斯逆变换的计算.

例题5—3 已知某信号的拉普拉斯变换表达式为

X(s)求该信号的时域表达式。

1

s23s2解:由于题目没有指定收敛域，所以必须考虑所有可能的情况。为此,可以先计算出该信号的拉

普拉斯变换表达式的极点.很显然，X(s)有两个极点，分别为 s = —1,s = -2。零极点分布图如

例题图5-3所示。

在MATLAB命令窗口键入：

〉〉 b = 1；

〉> a = ［1 3 2］；

>> ［r， p, k] = residue (b, a)

命令窗口立即给出计算结果为:

r = -1 1 p =

—2 -1 k =

［ ］

根据r、p、k之值，可以写出X（s)的部分分式和的表达式为： X(s)例题图5-3

11 s2s1然后根据不同的ROC,可写出X(s）的时域表达式x（t）。

第一种情况,ROC为 Re｛s} 〈 -2，则x(t)为反因果信号,其数学表达式为 x(t)e2tu(t)etu(t)

第二种情况，ROC为 —2 < Re{s} 〈 -1，则x（t）为双边非因果信号，其数学表达式为 x(t)e2tu(t)etu(t)

第三种情况，ROC为 Re｛s} > -1，则x（t)为因果信号，其数学表达式为 x(t)e2tu(t)etu(t)

在这个例题中，函数residue()仅仅完成了部分分式分解的任务，至于逆变换的数学表达式的结果是什么，还得结合收敛域的不同才能写出。

如果X(s）的分子的阶不小于分母的阶，则k将不等于一个空矩阵，例如，当

s3X(s)2时，我们在命令窗口中键入：

s3s2>> b = ［1 0 0 0］； 〉〉 a = [1 3 2]； 〉〉 [r,p，k］=residue（b,a)

则：

r = 8 -1 p = -2

90

—1 k =

1 -3

这里的k = [1 3]，实际上是将X（s)做了一个长除法后，得到的商的多项式.所以，根据上面的r、p、k之值，可写出X(s)的部分分式和的表达式为：

X(s)s381 s2s1有关函数residue（）的详细用法，可通过在线帮助加以了解。

三、实验内容及步骤

实验前，必须首先阅读本实验原理,了解所给的MATLAB相关函数，读懂所给出的全部范例程序。实验开始时,先在计算机上运行这些范例程序，观察所得到的信号的波形图.并结合范例程序所完成的工作,进一步分析程序中各个语句的作用，从而真正理解这些程序.

实验前,一定要针对下面的实验项目做好相应的实验准备工作，包括事先编写好相应的实验程序等事项。

Q5-1 将绘制零极点图的扩展函数文件splane以splane为文件名存盘. Q5—2 运行程序

Relation_ft_lt,观察拉普拉斯变换与傅里叶变换之间的关系.在点击工具条上

的旋转按钮,再将鼠标放在曲面图上拖动图形旋转，从各个角度观察拉普拉斯曲面图形，并同傅立叶变换的曲线图比较，加深对拉普拉斯变换与傅里叶变换之间关系的理解与记忆.

Q5—3 编写程序

Q5_3，能够接受从键盘输入的系统函数的分子分母多项式系数向量，并绘

制出系统的零极点图、系统的单位冲激响应、系统的幅度频率响应和相位频率相应的图形。各个子图要求按照图5.3布置。

程序Q5_3抄写如下：

91

Q5-4 执行程序编写Q5_3，输入因果的系统函数H(s)s的分子分母系数向量，

(s1)(s2)绘制所得到的图形如下：

执行Q5_3所得到的图形

从上面的图形中可以看出，该系统的零点和极点分别位于：

从时域和零极点分布特征两个方面说明该系统是否是稳定的系统？ 答:

从频率响应特性上看，该系统具有何种滤波特性？ 答：

92

Q5—5 执行程序编写Q5_3，输入因果的系统函数

1s21H(s)as32s22s1 此处a取1，执行程序Q5_3，输入该系统的分子分母系数向量，得到的图形如下：

从上面的图形中可以看出，该系统的零点和极点分别位于：

从时域和零极点分布特征两个方面说明该系统是否是稳定的系统？ 答:

从频率响应特性上看,该系统具有何种滤波特性？ 答：

93

改变系统函数中的a值，分别取0。6、0。8、4、16等不同的值，反复执行程序Q5-3，观察

系统的幅度频率响应特性曲线（带宽、过渡带宽和阻带衰减等)，贴一张a = 4时的图形如下：

观察a取不同的值时系统的幅度频率响应特性曲线的变化（带宽、过渡带宽和阻带衰减等），

请用一段文字说明零点位置对系统滤波特性的这些影响。

答：

12s1Q5-6 对于因果系统H(s)3a2，已知输入信号为x(t)sin(t)sin(8t)，要

s2s2s1求输出信号y(t)Ksin(t)，K为一个不为零的系数，根据Q5—5所得到的不同a值时的幅度频率响应图形,选择一个合适的a值从而使本系统能够实现本题的滤波要求。

你选择的a值为：

94

选择a值的根据是：

试编写一个MATLAB程序Q5_6，仿真这个滤波过程，要求绘制出系统输入信号、系统的单位

冲激响应和系统的输出信号波形。

抄写程序Q5_6如下：

执行程序Q5_6得到的输入输出信号波形图如下：

95

Q5-7 已知一个因果系统的系统函数为H(s)s5，作用于系统的输入信号为

s36s211s6x(t)e4tu(t)，试用MATLAB帮助你求系统的响应信号y(t）的数学表达式。

请在这里抄写你用MATLAB求解的命令（结合必要的文字说明)：

四、实验报告要求

1、按要求完整书写你所编写的全部MATLAB程序

2、详细记录实验过程中的有关信号波形图(存于自带的U盘中）,图形要有明确的标题.全部的MATLAB图形应该用打印机打印，然后贴在本实验报告中的相应位置，禁止复印件.

3、实事求是地回答相关问题，严禁抄袭。

本实验完成时间： 年 月 日

96