第一章 MATLAB及其操作简述
§1.1 MATLAB简介
MATLAB是由美国Math Works公司推出的软件产品。MATLAB是“Matric Laboratory”的缩写，意及“矩阵实验室”。MATLAB是一完整的并可扩展的计算机环境，是一种进行科学和工程计算的交互式程序语言。它的基本数据单元是不需要指定维数的矩阵，它可直接用于表达数学的算式和技术概念，而普通的高级语言只能对一个个具体的数据单元进行操作。因此，解决同样的数值计算问题，使用MATLAB要比使用Basic、Fortran和C语言等提高效率许多倍。许多人赞誉它为万能的数学“演算纸”。MATLAB采用开放式的环境，你可以读到它的算法，并能改变当前的函数或增添你自己编写的函数。在欧美的大学和研究机构中，MATLAB是一种非常流行的计算机语言，许多重要的学术刊物上发表的论文均是用MATLAB来分析计算以及绘制出各种图形。它还是一种有利的教学工具，它在大学的线性代数课程以及其它领域的高一级课程的教学中，已成为标准的教学工具。
最初的MATLAB是用FORTRAN编写的，在DOS环境下运行。新版的MATLAB 是C语言编写的高度集成系统。它在几乎所有流行的计算机机种，诸如PC、MACINTOSH、SUN、VAX上都有相应的MATLAB版本。新版的MATLAB增强了图形处理功能，并在WINDOWS环境下运行。现今，MATLAB的发展已大大超出了“矩阵实验室”的范围，在许多国际一流专家学者的支持下，Maths Works公司还为MATLAB配备了涉及到自动控制、信息处理、计算机仿真等种类繁多的工具箱（Tool Box），这些工具箱有数理统计、信号处理、系统辨识、最优化、稳健等等。近年来一些新兴的学科方向，Maths Works公司也很快地开发了相应的工具箱，例如：神经网络、模糊逻辑等。
§1.2 MATLAB操作说明
一．启动MATLAB
启动MATLAB有两种方法：
1．在WINDOWS下，点D:\\matlab\\bin\\matlab.exe进入MATLAB环境
2．在DOS下，键入matlab，这一命令会自动执行WINDOWS，并启动MATLAB 以上两种操作的结果都会出现一个MATLAB的命令窗口。 二．MATLAB的一些基本操作及命令函数 1．矩阵的输入
矩阵可以用几种不同的方法输入到MATLAB语言中： ⑴以直接列出元素的形式输入 ⑵通过语句和函数产生 ⑶建立在M文件中
⑷从外部的数据文件中装入
在MATLAB语言中不必描述矩阵的维数和类型，它们是由输入的格式和内容来确定的。
输入小矩阵最简单的方法是使用直接排列的形式，把矩阵的元素直接排列到方括号中，每行内的元素用空格或逗号分开，行与行的内容用分号格开。例如输入
A=[1 2 3；4 5 6；7 8 9] 或 A=[1，2，3；4，5，6；7，8，9] 都将得到输出结果：
A=
1 2 3 4 5 6 7 8 9
大的矩阵可以分行输入，用回车号代替分号。输入后矩阵A将一直保存在工作空间中，除非被替代和清除，A矩阵可以随时被调出来。若在命令末尾加上“；”号，则表示结果不显示，除非再次调用。
其余输入在实验中再做说明。 2．矩阵的运算
如果一个矩阵A有n行、m列元素，则称A矩阵为n×m矩阵，如果n=m， 则称矩阵A又称为方阵。MATLAB定义了下面各种矩阵的基本运算：
⑴矩阵转置
用符号“’”来表示矩阵的转置。如输入：
x=[-1 0 2]’
则输出为： x = -1 0 2 ⑵矩阵加、减
矩阵的加、减由符号“+”、“-”表示，它有两种格式：
a. 两种矩阵进行加减运算，其对应的元素进行加减，得到一新矩阵。如
输入：
A=[1 2 3；4 5 6；7 8 9]； B=[1 4 7；8 9 10；11 12 13]； C=A+B 则输出为： C= 2 6 10 12 14 16 18 20 22
b. 矩阵与标量进行加减运算，则矩阵中每个元素都与标量进行加减运 算，如输入： x=[-1，0，2]’； y=x-1 则输出为： y= -2 -1 1
⑶矩阵乘法，以符号“*”表示 a. 两矩阵相乘，如输入 x=[2 3 4 5；1 2 2 1];
y=[0 1 1；1 1 0；0 0 1；1 0 0]；
z=xy 则输出为： z=
8 5 6 3 3 3
b. 矩阵与标量相乘，如输入 A=[2，3，4，5；1，2，2，1]； B=A2 则输出为： B=
3 6 8 10 2 4 4 2
⑷矩阵的求逆
以“inv”来表示，如inv(A)表示A的逆。如输入： A=[1 2 3；4 2 6；7 4 9]； inv(A) 则输出为：
ans= -0.5000 -0.5000 0.5000
0.5000 -1.0000 0.5000 0.1667 0.8333 -0.5000
⑸矩阵的乘方，以符号“^”表示
a^p表示a的p次方，即a自乘p次。当p为矩阵时，运算会出错。 3．MATLAB定义的点运算
MATLAB中定义了一种特殊的运算，即所谓的点运算。两个矩阵之间的点运算是该矩阵对应元素的直接运算。注意：点乘积运算要求参与运算的两个矩阵维数相同，其具体操作如下所示： ⑴乘除运算
在MATLAB中，符号“.”表示乘法运算，相乘的数组要有相同的维 数，而符号“./”表示数组除法运算，且为对应元素进行乘除。 如输入：
A=[1，2，3]； B=[4，5，6]；
C=A. B 输出为： C=
4 10 18
输入： D=B./A 输出为： D=
4.00 2.5000 2.0000
⑶乘方运算，以符号“.^”表示
a. 当x、y均为向量时，z=x.^y表示对应元素的乘方。 b. 当x为向量，y为标量时，z=x.^y表示z(i)=x(i)y。 c. 当x为标量，y为向量时，z=x.^y表示z(i)=xy(i)。 该运算也可用于任意维矩阵。
4．在MATLAB中，复数序列用“i”或“j”表示sqrt(-1)，即5．关系运算
MATLAB中共有六种关系，分别为： < 小于 > 大于 == 等于 <= 小于等于 >= 大于等于 ~= 不等于
对两矩阵的对应元素进行比较，若关系成立则为1，否则为0。 6．常用命令函数
！ 执行操作系统命令 abs 绝对值函数 angle 相角函数 axis 坐标轴标度设定
cla 清除当前坐标轴 clc 清除命令窗口显示 clf 清除当前图形窗口 close 关闭图形窗口 delete 删除文件
demo 运行MATLAB演示程序 function grid gtext help hold imag length linspace logspace pi plot quit real size sqrt stem subplot title xlabel ylabel MATLAB函数表达式的引导符给图形加网格线
在鼠标指定的位置加文字说明 启动联机帮助文件显示 当前图形保护模式 求取虚部函数 查询向量的维数 构造线性分布的向量 构造等对数分布的向量 圆周率π
线性坐标图形绘制 退出MATLAB环境 求取实部函数 查询矩阵的维数 平方根函数
函数序列柄状图形绘制 将图形窗口分成若干个区域 给图形加标题 给图形加x坐标说明 给图形加y坐标说明
实验一、熟悉MATLAB环境
[实验目的]
1． 熟悉MATLAB主界面，并学会简单的菜单操作。 2． 学会简单的矩阵输入与信号输入。 3． 掌握部分绘图函数。 [实验原理]
MATLAB 是以复杂矩阵作为基本编程单元的一种程序设计语言。它提供了各种矩阵的运算与操作，并有较强的绘图功能。
用户第一次使用MATLAB时，建议首先在屏幕上键入DEMO命令，它将启动MATLAB的演试程序，用户可在此演示程序中领略MATLAB所提供的强大的运算与绘图功能。也可以键入HELP进行进一步了解。
MATLAB命令窗口的界面下有一个菜单条，如图所示：
利用FILE菜单可方便对文件或窗口进行管理。其中FILE|NEW的各子菜单，M-FILE（M文件）、FIGURE（图形窗口）、或MODEL（SIMULINK编辑界面）分别可创建对应文件或模块。EDIT菜单允许用户和WINDOWS的剪切板交互信息。
MATLAB语言最基本的赋值语句结构为：变量名列表=表达。表达式由操作符或其它字符，函数和变量名组成，表达式的结果为一个矩阵，显示在屏幕上，同时输送到一个变量中并存放于工作空间中以备调用。如果变量名和“=”省略，则ANS变量将自动建立，例如键入：1900/81，得到输出结果：ans =23.4568。
在MATLAB中把数据绘成图形可有多种命令一供选择。下面列出了这些命令：
绘图命令 Plot 线性X-Y坐标图 Loglog 双对数坐标图 Semilogx X轴对数半对数坐标图 Semilogy Y轴对数半对数坐标图 Polar 极坐标图 Mesh 三维消隐图 Contour 等高线图 Bar 条形图 Stairs 阶梯图 除了可以在屏幕上显出图形外，还可以对屏幕上已有的图形加注释、题头或坐标网格。
图形加注 Title 画题头 Xlabel x轴标注 Ylabel y轴标注 Text 任意定位的标注 Gtext 鼠标定位标注 Dgrid 网格 关于坐标轴尺寸的选择和图形处理等控制命令：
图形控制命令 Axis 人工选择坐标轴尺寸 Clg 清除图形窗口 Ginput 利用鼠标的十字准线输入 Hold 保持图形 Shg 显示图形窗口 Subplot 将图形窗口分成N块子窗 [实验内容]
一．熟悉简单的矩阵输入
1． 从屏幕上输入矩阵A=[1 2 3；4 5 6；7 8 9] 或A=[1，2，3；4，5，6；7，8，9] 观察输出结果。
2． 试用回车代替分号，观察输出结果。 3． 输入矩阵B=[9，8，7；6，5，4；3，2，1]
C=[4，5，6；7，8，9；1，2，3], 键入A B C 观察结果 4． 选择File|new菜单中的M-file，输入B=[9 ，8，7；6，5，4；
3，2，1]，保存为B.M文件，退出编辑环境。此时在工作环境中使用B命令就可调出B矩阵。
[注]4.2版的MATLAB不能直接存为.m文件，而存为.txt文件， 需在DOS下改为M文件，即在工作环境下键入 !rename B.txt B.m 既可。
5． 再试着输入一些矩阵，矩阵中的元素可为任意表达式，但注意 矩阵中各行各列的元素个数需分别相等，否则会给出出错信息。 6． 输入who和whos观察结果，了解其作用。 二．基本序列运算
1． 数组的加减乘除和乘方运算
输入A=[1 2 3]，B=[4 5 6]，求C=A+B，D=A-B，E=A.*B，F=A./B，G=A.^B,并用stem画出A，B，C，D，E，F，G。再输入一些数组，进行类似运算。 2．粗略描绘下列各函数的波形（其中对于连续信号可取时间间隔为0.001）
a.
b. c. d.
e. f. g.
实验二 离散时间信号与系统
[实验目的]
1．了解信号处理的基本操作 2．熟悉一些常用的序列及其应用 [实验原理]
我们所接触的信号大多为连续信号，为使之便于处理，往往要对其进行采样，对信号抽样并保证其能完全恢复，对抽样频率有一定的限制。
基本的离散序列的定义如下： 1．单位采样序列
2．单位阶跃序列
3．实指数序列
，
4．复数指数序列
5．正余弦序列
，
6．周期序列
，
[实验内容]
； a为实数
，
1．用MATLAB实现函数impseq(n0,n1,n2)，使函数实现该函数的格式为：
Function [x,n]=impseq(n0,n1,n2)
% Generate x(n)=delta(n-n0);n1<=n<=n2 % [x,n]=impseq(n0,n1,n2)
，。
2．用MATLAB实现函数stepseq(n0,n1,n2)，使函数实现u(n-n0)，该函数的格式为：
function [x,n]=stepseq(n0,n1,n2) %Generate x(n)=u(n-n0);n1<=n<=n2 %[x,n]=stepseq(n0,n1,n2) 3．用MATLAB实现下列序列 a. b. c.
d. 将c中的
扩展为
。
，周期数为4
和
的长度必须相
4．MATLAB中可用算术运算符“+”实现信号相加，但
等。如果序列长度不等，或者长度虽然相等但采样的位置不同，就不能运用“+”了。试用MATLAB写出任意序列相加的函数sigadd,其格式如下： Function % 实现
%
% y=在包括n1和n2的n上求和序列
% x1=长为n1的第一个序列，x2=长为n2的第二个序列（n2 可与n1不等） 5．与sigadd相仿，建立一个信号相乘sigadd函数 6．建立一个函数sigshift，实现fuction % 实现%
，函数格式如下：
7．建立一个函数sigfold，实现。MATLAB中，这一运算由fliplr(x)函数
实现，而对采样位置则由-fliplr(n)得到。格式与上类同。 8．用 MATLAB产生并画出（用stem 函数）下列序列的样本：
a)
b) c) d)
(其中
是一个在[0，1]之间均匀分布的随机序列，用rand(1,N)实现，其中
N表示长度)
e)
，画出五个序列
9．令x(n)=[1,-2,4,6,-5,8,10]，产生并画出下列序列的样本
a)
b)
10．将题9中的序列分解为偶和奇分量。用stem画出这些分量 其中偶部： 奇部：
创建函数evenodd，实现奇偶分量 11．考虑模拟信号
。分别用Ts=0.5,0.25,0.1秒时的采样
间隔对它采样以获得x(n)，对每个Ts，画出x(n)，讨论所得结果。 12．信号的扩展（或抽取，或降低采样频率）定义为：
其中x(n)的采样频率被降低了整数因子M。 a. 开发一个MATLAB函数dnsample，其格式为 Function y=dnsample(x,M)
用以实现上述运算。在应用MATLAB的下标功能时要特别注意时间轴的原点n=0。 b.
。频率降低因子为4，求y(n)。用subplot
函数分别画出x(n)和y(n)，并对结果进行讨论。 c. 用
信号的影响。
重复上题，定性地讨论降低采样频率对
[思考题]
1．任意复值序列x(n)均可分解为： 其中
和
a. 修改evenodd函数，使它能接受任意序列并把它分解为上式表示的分量。 b. 分解下列序列性。
2．对11题中的序列用sinc内插（取并从图中求出
），用
样本集重构模拟信号
,
，画出它的实部和虚部，验证共轭对称
的频率，忽略尾部效应。
。
用三次样条内插，重构
实验三 卷积实验
[实验目的]
1．熟悉并验证卷积的性质
2．利用卷积生成新的波形，建立波形间的联系 3．验证卷积定理
4．分析卷积与圆周卷积的区别与联系 [实验原理]
信号的卷积是针对时域信号处理的一种分析方法。信号的卷积一般用于求取信号通过某系统后的响应。在信号与系统中，我们通常求取某系统的单位冲激响应，所求得的h(k)可作为系统的时域表征。任意系统的系统响应可用卷积的方法求得：
y(k)?x(k)?h(k)
圆周卷积是数字信号处理的重要内容，由于它也是卷积的一种，必须了解其定义，圆周卷积的条件为卷积长度必须大于每个序列长度和。 [实验内容]
1．MATLAB提供了一个内部函数conv来计算两个有限长序列的卷积。conv函数假定两个序列都从n=0开始。给出序列x=[3,11,7,0,-1,4,2]； h=[2,3,0,-5,2,1]；求两者的卷积y 。
将函数conv稍加扩展为函数conv_m，它可以对任意基底的序列求卷积。格式如下：
function [y,ny]=conv_m(x,nx,h,nh) % 信号处理的改进卷积程序 % [y,ny]=conv_m(x,nx,h,nh) % [y,ny]=卷积结果 % [x,nx]=第一个信号 % [h,nh]=第二个信号
2．创建函数circonv，来实现序列的圆周卷积，格式如下： function y=circonv(x1,x2,N)
3．对下面三个序列，用conv_m函数，验证卷积特性（交换律、结合律、分配 律、同一律）
x1(n)?x2(n)?x2(n)?x1(n) 交换律 [x1(n)?x2(n)]?x3(n)?x1(n)?[x2(n)?x3(n)] 结合律 x1(n)?[x2(n)?x3(n)]?x1(n)?x2(n)?x1(n)?x3(n) 分配律 x1(n)??(n?n0)?x1(n?n0) 同一律 其中：x1(n)=n[u(n+10)-u(n-20)] x2(n)=cos(0.1?n)[u(n)-u(n-30)] x3(n)=(1.2)n[u(n+5)-u(n-10)]
4．求出下列序列的自相关序列rxx(l)和互相关函数rxy(l)
nx(n)?0.9 0?n?20 ?n y(n)?0.8 ?20?n?0
你能观察什么结果?
1,2,2?，x2(n)??1,2,3,4?，分别计算8点、5点、6点循环卷积，对 5．设x1(n)?? 所得结果进行讨论。
?6．令x(n)?3cos(0.5?n?60)?2sin(0.3?n) |n| h1(n)?0.9
h2(n)?sa(0.2n)[u(n?20)?u(n?20)]，其中sa(0)?1
（MATLAB中sa函数用sinc表示）
nnh(n)?(0.5?0.4)u(n) 3 对每一种情况求出其输出y(n) [思考题]
1． MATLAB提供了一个称为toeplitz的函数，可根据第一行和第一列生成 toeplitz矩阵。用此函数开发另一个MATLAB函数来执行线性卷积，此函数的规范格式为： Function [y,H]=conv_tp(h,x) %用toeplitz矩阵的线性卷积 % [y,H]=conv_tp(h,x) %y=列向量形式的输出序列
%H=对应于序列h的toeplitz矩阵，因而y=Hx %h=列向量形式的脉冲响应序列
%x=列向量形式的输入序列
2． 设x(n)=n+1，0?n?9，h(n)={1,0,-1}，按N=6用重叠保留法计算 y(n)?x(n)?h(n)
实验四 零极点实验及其频响
[实验目的]
1．掌握系统函数零极点定义 2．零极点与频响的关系 3．零极点与系统稳定性的关系 4．状态方程含义 5．使用zplane函数 [实验原理]
该实验用MATLAB中库函数，如tf2zp(b,a)，ss2zp(A,B,C,D)，zplane(z,p)，freqz(b,a)等。例如：
s2?0.5s?2H(s)?2s?0.4s?1，求其零极点图。 1．传递函数为
程序如下：
num=[1 0.5 2]； 分子系数，按降幂顺序排列 den=[1 0.4 1]； 分母系数，按降幂顺序排列 [z,p]=tf2zp(num,den)； 用tf2zp函数求出其零点z和极点p zplane(z,p) 作出零极点图
2．若给出的是滤波器的输入与输出的状态方程，如：
?10??1??1?x???x?uy??1?x??0?u??0??1?34??????，，求其零极点图。 程序如下： A=[1，0；1，-3]； B=[1；0]；
1 C=[-4，1]；
D=0；
[z,p]=ss2zp(A,B,C,D)；求出其零极点z,p zplane(z,p)
在连续时间系统中，当极点在虚轴的右半平面时，系统不稳定，在虚轴上，
单阶极点系统稳定；若零点均处于左半平面内，则系统为最小相位系统。在离散系统中，极点在单位圆外系统不稳定，在单位圆上，单阶极点系统稳定；零点在单位圆内，系统为最小相位系统。
对一滤波器，我们不仅要知道它的零点和极点，还要了解它的频率特性，本实验可求其频率特性。
对模拟滤波器，可用freqs函数求得其频率特性，对数字滤波器，则用freqz函数求得。
0.2s2?0.3s?1H(s)?2s?0.4s?1，求其频率特性。 例3．已知模拟滤波器的传递函数为
程序如下：
num=[0.2 0.3 1]; den=[1 0.4 1];
w=logspace(-1,1); 频率范围 freqs(num,den,w)
0.2?0.3z?1?z?2H(z)?1?0.4z?1?z?2，取样点数为128点，求其频率特性。 例4．数字滤波器
程序如下：
num=[0.2 0.3 1]; den=[1 0.4 1]; freqz(num,den,128) [实验内容]
1．已知下列传递函数H(s)或H(z)，求其零极点，并画出零极点图。
H(s)?H(s)?a. b.
3(s?1)(s?2)(s?1)(s?2) 1s
s2?1H(s)?2s?2s?5 c.
3z3?2z2?2z?5H(z)?0.6?3z?3z2?2z?1 d.
10??0?0??x??0?ux???001????????6?11?6???1?? y??451?x e.
2．求出下列系统的零极点，分析其稳定性，并判断它们是否为最小相位系统。
?1?5210???2???x???0460?x???u?3?????0?3?6?1???4? y??1252?x a.
?221??3??x??3?ux???131???????122???4?? y?x b.
c.
s4?35s3?291s2?1093s?1700G1(s)?9s?9s8?66s7?294s6?1029s5?2541s4?4684s3?5856s2?4629s?1700
d.
G2(s)?15(s?3)(s?1)(s?5)(s?15)
100s(s?2)2(s2?3s?2)2G3(s)?(s?1)(s?1)(s3?3s2?5s?2)((s2?1)2?3)2 e.
1?1.414z?1?z?2z?1H(z)??1?0.9z?1?0.81z?2z?0.3 f.
3．已知下列H(s)或H(z)，求其频响。 a. b.
H(z)?11?z?1
H(s)?2ss2?2s?1
(1?z?1)2H(z)?1?0.61z?2 c.
d.
H(s)?3(s?1)(s?2)(s?1)(s?2)
?10??1??1?x???x?ey??1?x??0??1?3?4? ???? e. [思考题]
1?2z?1?3z?27?8z?1?9z?2H(z)???1?24?5z?6z10?11z?1?12z?2，求其零极点。 1．已知
（提示：令H(z)?H1(z)?H2(z)，分别求出零极点，再合并） 2．已知数字滤波器的状态方程为：
?10??1??1?x???x?e,y??1?x?0?e?????4??1?3??0? 求其频响特性。 3．设有一模拟滤波器
Ha(s)?1s2?s?1，采样周期为T=2，试用双线性变换法
将它转变为数字滤波器H(z)，并求原模拟滤波器与转换后滤波器的频响特性 （已知双线性变换法的函数为bilinear，其使用方法为：
[numd,dend]=bilinear(num,den,Fs)，其中num,numd,den,dend分别为转变前 后传递函数的分子、分母系数，Fs为采样频率）改变采样周期，观察输出 信号有何不同？
实验九 综合实验
[实验目的]
1．掌握综合实验的能力
2． 复习所学内容，综合所学知识
[实验原理] 1．
H(s)?112??s?1s?3s?4，求其零极点图，并判断系统的稳定性
2． 设计一5阶低通Bessel（贝塞尔）模拟滤波器，截止频率为1000rad/s，求其零
极点及频响特性。
3． 设计一6阶带通Bessel模拟滤波器，通带为2000???10000，求其零极点及频响
特性。
4． 设采样频率为fs?1.2kHz，用脉冲响应不变法设计一三阶巴特沃兹数字底通滤波
器，截止频率fc?400Hz。
5． 设有两序列，试比较X(enj?)和Y(ej?)
6． 用式x(n)?10?0.8，0?n?10，验证DFT的Parseval定理：
1N?1|x(n)|??|X(k)|2?Nk?0 n?0
2N?17． 一因果系统LTI系统
y(n)?0.81y(n?2)?x(n)?x(n?2)
j?求1)H(z) 2)冲激响应h(n) 3)单位阶跃响应u(n) 4)H(e)，并绘出幅频和
相频特性。 8．求解差分方程
1y(n)?[x(n)?x(n?1)?x(n?2)]?0.95y(n?1)?0.9025y(n?2)3 n?0
其中
x(n)?cos(?3n),y(?1)??2,y(?2)??3,x(?1)?1,x(?2)?1