《信号与系统(第5版)MATLAB仿真举例.docx》由会员分享，可在线阅读，更多相关《信号与系统(第5版)MATLAB仿真举例.docx（34页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、最好的沉淀整理MATLAB 仿真举例一、MATLAB 应用概要1. MATLAB 的知识与特点MATLAB 是由美国 Mathworks 公司于 1984 年推出的一种功能强大的科学计算软件，主要适用于矩阵运算、信息处理和控制等领域的分析与设计。 MATLAB 是“矩阵实验室 ”（MATrix LABoratoy）的缩写。它以矩阵运算为基础，以简洁便用的“科学便笺式”编程语言， 可交互的集成环境，强大的图形可视化功能，广泛应用于科学、工程和绘图的各种领域，受到科学界和教育界的高度好评。其特点是：(1) 数值计算功能强，每个数值元素都视为复数，用双精度（64 位）一种格式，编程简化；(2) 具有
2、符号计算功能；(3) 简单而强大的作图功能；(4) 人机界面适合科技人员；(5) 功能多，可扩展性强，有多种工具箱作专门处理。2. MATLAB 的窗口与应用常识(1) 命令窗口以 MATLAB6.x 为例，其集成环境有桌面平台及8 个组成部分：指令窗口、历史指令窗口、工作台及工具箱窗口、工作目录窗口、工作空间窗口、矩阵编辑器、程序编辑器和帮助浏览器。点击桌面上的MATLAB 图标进入 MATLAB后，即可看到命令窗口（Command Window）,它是主窗口，当显示“ ”时，表示系统处于准备接受命令的状态。(2) 图形窗口一般，只要执行任一种绘图命令，就会自动产生图形窗口,以后的绘图都在这
3、一个图形窗口中进行。(3) 文本编辑窗口用 MATLAB 计算有两种方式，一种是简单计算，可在命令窗口一行一行地输入各命令， 另一种是对复杂计算时，把多行命令组成一个M 文件，让 MATLAB 自动执行。(4) M 文件采用 MATLAB 所规定的一套语言及语法编写的源文件记作.m 文件。其文件名不能以数字开头，也不能包含文字，但可以用下划线“_”。(5) Help 命令Help 命令是查询函数相关信息的最基本方式，信息会直接显示在命令窗口中。(6) 应用常识 MATLAB 使用双精度数据，所有系统命令都是小写形式。 矩阵是 MATLAB 进行数据运算的基本元素，矩阵中的下标从1 开始，而不是
4、0，标量是作为 11 的矩阵来处理的。 语句或命令的结尾的分号“；”会屏蔽当前结果的显示。 注释位于%之后，不被执行。 使用上下箭头实现命令的滚动显示，可用于再编辑和再执行。 变量名必须以字母开头，可以由字母、数字和下划线组成。 要画出一条平滑的连续曲线，最少需要200 个数据点。 即使预先对 t 做了定义，也应该使用数组操作，如 f(t)=(2.(-t).*sin(2*t); 若写作2(-t)*sin(2*t)时会提示出错。 要想在同一个图形坐标上绘制多条曲线，应使用plot(t,f)，hold on，plot(t1,f1)，hold off 或 plot(t,f,t1,f1)均可。 为了在
5、图形中标出w ，j ，p 等特殊符号，应采用MATLAB 提供的专门字母表示。如表 1 所示。表 1特殊符号的表示法字符串alphabeta符号ab字符串thetalambda符号ql字符串phiPhi符号jfgammagtaotpipdeltadintinftyepsiloneomegawOmegaW3. MATLAB 的常用函数(1) 基本数学函数常用的基本数学函数如表 2 所示。表 2基本数学函数正弦sin平方根sqrt余弦cos求数组的最大值max正切tan求数组的最小值min指数exp求数组的平均值mean自然对数log求和sum矩阵指数函数expm求矩阵的特征值eig绝对值abs求
6、复数的相角angle求实部real求共轭conj求虚部imag辛格函数sinc(2) 基本作图函数 作图函数如表 3 所示表 3常用的作图函数plot绘制连续波形title为图形加标题stem绘制离散波形grid画网格线loglog绘双对数坐标图xlabel为 X 轴加上轴标subplot分割图形窗口ylabel为 Y 轴加上轴标hold保留当前曲线text在图上加文字说明figure定义图形窗口gtext用鼠标在图上加说明line画直线axis zzplot定义 x,y 坐标轴标度画符号函数的图形4 矩阵的创建在 MATLAB 中，矩阵是运算的基本单元，所以创建矩阵非常重要。创建矩阵有三种方
7、法：直接输入法、内部函数创建法和利用矩阵编辑器和修改矩阵。(1) 直接输入法实数矩阵：对于维数较小的实数矩阵，可以从键盘直接输入。具体方法是：将矩阵的元素用方括号“”括起来，按行输入各元素，各行用分号“；”或 Enter 隔开。 123 例如矩阵 a= 456 789 可写成指令：a=1 2 3; 4 5 6; 7 8 9复数矩阵：矩阵的各元素为复数时，与上述实数矩阵类似。不过写法有规定，如 6+4j不可写成 6+j4。应该写为 6+j*4 或 6+i*4。(2) 内部函数创建法利用 MATLAB 提供的许多函数可以创建一些特殊矩阵。例如： ones(m,n)生成一个元素全部为 1 的 m 行
8、 n 列矩阵eye(m,n)生成一个主对角线元素为 1 的其余为 0 的矩阵pascal(n)生成一个 n 维帕斯卡矩阵(3) 向量的创建创建向量最常用的方法是利用冒号“：”运算符生成。例如：a=m:p:n是生成起始值为 m，终值为 n，步长为 p 的均匀等分行向量 a。例如：a=1:-0.1:0则输出结果为：a=10.90.80.70.60.50.40.30.20.105. 信号波形的产生和运算函数这里列出常用的波形产生和运算的函数如下：sawtooth产生周期锯齿波或三角波cumsum信号累加square产生周期方波sum信号求和pulstran脉冲串diff信号微分sinc辛格函数qua
9、d信号积分fliplr信号翻转conv信号卷积6. 连续系统时域分析的函数在时域分析中，常用的符号形式如表 4 所示表 4时域分析符号名称阶跃函数e (t)冲激函数d (t) 定义符号变量化简微 分 二阶导数定积分求一阶微分方程全解求二阶微分方程全解符号形式Heaviside(t) Dirac(t) symssimple(f) diff(f) diff(f,2)int(f,a,b)dsolve(Dy+a*y=b,y(0)=k) dsolve(D2y+Dy)说明使之包含的字符最少f 被积函数，a,b 为积分的上下限a,b,k 均为数值向量求冲激响应impulse(b,a)a = a0, a ,
10、aL1N, b = b , b12, b LM求阶跃响应step(b,a)同上7. 连续信号频域分析函数进行周期信号的频谱分析，其符号形式为：int(f,t,a,b)其中 f 为符号表达式，t 为积分变量，a 和 b 分别为积分的下限和上限例如，计算傅里叶系数的被积函数 Ae - jnw1t ，则有：f=A*exp(-j*n*2*pi/T*t); Fn=int(f,t,-tao/2,tao/2)/T;求非周期信号的傅里叶变换时，调用格式为：F=fourier(f,t,w)其中 f 为信号表达式，t 为积分变量，w 为角频率，F 表示 f(t)的傅氏变换。傅里叶反变换指令格式为：f=ifouri
11、er(F,w,t) f=ifourier(F)例如，对 f (t) = 3e-2t e (t) ，则有f=3*exp(-2*t)*sym(Heaviside(t); F=fourier(f);分析系统的频率响应，调用格式为：H=fregs(b,a,w)式中，b 为 H(w)中的分子多项式的系数向量，a 为分母多项式系数向量，w 为计算 H(w)的取样点数。8. 系统的复频域分析函数 信号的拉氏变换指令为：Fs=laplace(ft,t ,s)式中，ft 为信号 f(t)的符号表达式，t 为积分变量，s 为复频率，Fs 为拉氏变换 F(s)。 信号的拉氏反变换格式为：ft=ilaplace(fs
12、,s,t)或者ft=ilaplace(fs)若利用部分分式展开求反变换，则展开式的指令格式：r,p,k=residue(b,a)其中，b = b ,b ,L,b 为 F(s)的分子多项式系数向量；a = a , a ,L, a 为分母多项式的01m01n系数向量， r = r , r ,L, r T 为留数列向量，即展开的各分式的系数；p 为各极点列向量，1 2nk 为直接项系数行向量。 由 H(S)确定零、极点，指令格式为：sys=tf(b,a)指定系统函数sys=zpk(z,p,k)由零、极点构成系统函数z=tzero(sys),p=pole(sys)求系统的零、极点pzmap(sys)绘
13、制零、极点图impulse(sys)由 H(s)绘制冲激响应曲线step(sys)由 H(s)绘制阶跃响应曲线9. 离散系统分析的基本函数 求差分方程零状态响应时，其调用格式为：y=filter(b,a,x)其中，x 为输入信号行向量，； a = a0, a , aL1N为差分方程输出变量的系数向量，b = b ,b ,L,b 为输入变量的系数向量。01M求单位响应 h(n)是调用格式为：h=impz(b,a) 求离散卷积和时，调用格式为：y=conv(h,x)其中，h 为单位响应序列值，x 为输入序列值，分别存储 h(n)和 x(n)的行向量。 z 变换，调用格式为：Fz=ztrans(fn
14、,n,z)其中，fn 为 f(n)的符号表达式；n 为序号；z 为复频率。用长除法进行 z 反变换的指令为：impz(b,a,N)， N 为样点的长度也可以用部分分式展开法，即r,p,k=residue(b,a) 计算频率响应，调用格式为：freqz(b,a) freqz(b,a,n)其中，b 和 a 分别为系统函数 H(z)的分子、分母之系数向量，n 为频率的计算点数，常取 z 的整数次幂；绘制频率特性的横坐标W （即w T）的范围为 0 到p 。二、MATLAB 方法实现信号波形例 1 应用 MATLAB 方法实现单位阶跃信号和矩形脉冲。解：对于阶跃函数，MATLAB 中有专门的 stai
15、rs 绘图命令。例如，实现e (t) 和矩形脉冲的程序如下：t=-1:2;x=(t=0);% 定义时间范围向量 tsubplot(1,2,1),stairs(t,x);axis(-1,2,-0.1,1.2);grid on% 绘制单位阶跃信号波形t=-1:0.001:1;% 定义时间范围向量 tg=(t=(-1/2)-(t=(1/2);subplot(1,2,2),stairs(t,g);axis(-1,1,-0.1,1.2);运行结果如图 1 所示。grid on% 绘制矩形脉冲波形图 1例 1 图例 2 应用 MATLAB 方法生成信号 f (t) = sin c(t) 和 f (t) =
16、 S (t) 的波形。a解：为生成函数sin c(t) = sin ptpt可直接调用 MATLAB 中的专门命令，程序如下：t=-5:0.01:5;% 定义时间范围向量 tf=sinc(t);% 计算 Sa(t)函数plot(t,f); grid on% 绘制 Sa(t)的波形运行结果如图 2 所示。图 2例 2 程序运行结果一S (t) 和sin c(t) 的关系如下：asin(p t )p f (t) = S(t) = sin t =p= sin(t ) = sin c(t)atptptp生成信号 f (t) = Sa(t) 波形的 MATLAB 程序如下：t=-3*pi:0.01*pi
17、:3*pi;% 定义时间范围向量 t f=sinc(t/pi);% 计算 Sa(t)函数plot(t,f); grid on% 绘制 Sa(t)的波形运行结果如图 3 所示。图 3例 2 程序运行结果二例 3 应用 MATLAB 方法生成相加信号 f (t) = cos18pt + cos 20pt 和相乘信号f (t) = sin c(t) cos(20p t) 的波形。解：对相加信号 f (t) = cos18pt + cos 20pt ，程序如下：syms t;% 定义符号变量 tf=cos(18*pi*t)+cos(20*pi*t);% 计算符号函数 f(t)=cos(18*pi*t)
18、+cos(20*pi*t) ezplot(f,0 pi); grid on% 绘制 f(t)的波形运行结果如图 4 所示。图 4例 3 程序运行结果一对相乘信号 f (t) = sin c(t) cos(20p t) ，程序如下：t=-5:0.01:5;% 定义时间范围向量f=sinc(t).*cos(20*pi*t);% 计算函数 f(t)=sinc(t)*cos(20*pi*t)plot(t,f);% 绘制 f(t)的波形title(sinc(t)*cos(20*pi*t); grid on% 加注波形标题运行结果如图 5 所示。图 5例 3 程序运行结果二例 4 应用 MATLAB 方法
19、生成调制信号 f (t) = (2 + 2 sin 4p ) cos 50pt 的波形。解：对调制信号 f (t) = (2 + 2 sin 4p ) cos 50pt ，程序如下：syms t;% 定义符号变量 tf=(2+2*sin(4*pi*t)*cos(50*pi*t);% 计算符号函数 f(t)=(2+2*sin(4*pi*t)*cos(50*pi*t) ezplot(f,0 pi); grid on% 绘制 f(t)的波形运行结果如图 6 所示。图 6例 4 图三、MATLAB 方法用于时域分析例 5 设方程 y (t) + 5 y(t) + 6 y(t) = 2e-t e (t)
20、 ，试求零状态响应 y(t) 。解：程序如下：yzs=dsolve(D2y+5*Dy+6*y=2*exp(-t),y(0)=0,Dy(0)=0)% 利用 dslove 命令求解零状态响应ezplot(yzs,0 8); grid on% 绘制零状态响应曲线运行结果：即：图形如图 7 所示。yzs =exp(-t)+exp(-3*t)-2*exp(-2*t)y(t) = (e-t - 2e-2t + e-3t )e (t)图 7例 5 图例 6 已知二阶系统方程u (t) + R u (t) +1 u (t) = 1 d (t)cL cLC cLC对下列情况分别求h(t) ，并画出其波形。a.
21、R = 4W, L = 1H , C = 1/ 3Fb. R = 2W, L = 1H , C = 1Fc. R = 1W, L = 1H , C = 1Fd. R = 0W, L = 1H , C = 1F解：程序如下：R=input(电阻 R=); % 以交互方式输入电阻R 的值L=input(电感 L=); % 以交互方式输入电阻L 的值C=input(电容 C=); % 以交互方式输入电阻C 的值b=1/(L*C);a=1 R/L 1/(L*C);impulse(b,a);% 绘制脉冲响应h(t)的波形运行结果如图 8 至图 11 所示。a. 电阻 R=4 电感L=1 电容 C=1/3
22、图 8例 6 程序运行结果一b. 电阻 R=2 电感 L=1 电容 C=1图 9例 6 程序运行结果二c. 电阻 R=1 电感L=1 电容 C=1图 10例 6 程序运行结果三d. 电阻 R=0 电感 L=1 电容 C=1图 11例 6 程序运行结果四例 7 实现卷积 f (t) * h(t) ，其中： f (t) = 2e (t) - e (t - 1), h(t) = e (t) - e (t - 2) 。解：主程序如下：p=0.01;% 取样时间间隔nf=0:p:1;% f(t)对应的时间向量f=2*(nf=0)-(nf=1);% 序列 f(n)的值nh=0:p:2;% h(t)对应的时
23、间向量h=(nh=0)-(nh=2);% 序列 h(n)的值y,k=sconv(f,h,nf,nh,p);% 计算 y(t)=f(t)*h(t) subplot(3,1,1),stairs(nf,f); grid on% 绘制 f(t)的波形title(f(t);axis(0 3 0 2.1);subplot(3,1,2),stairs(nh,h); grid on% 绘制 h(t)的波形title(h(t);axis(0 3 0 2.1);subplot(3,1,3),plot(k,y); grid on% 绘制 y(t)=f(t)*h(t)的波形title(y(t)=f(t)*h(t);a
24、xis(0 3 0 2.1);子程序 sconv 如下：% 此函数用于计算连续信号的卷积 y(t)=f(t)*h(t) function y,k=sconv(f,h,nf,nh,p)% y:卷积积分 y(t)对应的非零样值向量% k:y(t)对应的时间向量% f:f(t)对应的非零样值向量% nf:f(t)对应的时间向量% h:h(t)对应的非零样值向量% nh:h(t)对应的时间向量% p:取样时间间隔y=conv(f,h);% 计算序列 f(n)与 h(n)的卷积和 y(n)y=y*p;% y(n)变成 y(t)left=nf(1)+nh(1);% 计算序列 y(n)非零样值的起点位置ri
25、ght=length(nf)+length(nh)-2;% 计算序列 y(n)非零样值的终点位置k=p*(left:right);% 确定卷积和 y(n)非零样值的时间向量运行结果如图 12 所示。图 12例 7 图例 8 实现卷积 f (t) * h(t) ，其中： f (t) = 2e (t) - e (t - 2), h(t) = e-t e (t) 。解：主程序如下：p=0.01;% 取样时间间隔nf=0:p:2;% f(t)对应的时间向量f=2*(nf=0)-(nf=2);% 序列 f(n)的值nh=0:p:4;% h(t)对应的时间向量h=exp(-nh);% 序列 h(n)的值y
26、,k=sconv(f,h,nf,nh,p);% 计算 y(t)=f(t)*h(t) subplot(3,1,1),stairs(nf,f); grid on% 绘制 f(t)的波形title(f(t);axis(0 6 0 2.1);subplot(3,1,2),plot(nh,h); grid on% 绘制 h(t)的波形title(h(t);axis(0 6 0 2.1);subplot(3,1,3),plot(k,y); grid on% 绘制 y(t)=f(t)*h(t)的波形title(y(t)=f(t)*h(t);axis(0 6 0 2.1);子程序 sconv 同例 7。运行结
27、果如图 13 所示。图 13例 8 图四、MATLAB 方法用于频域分析例 9 求图 14(a)所示周期矩形脉冲信号的傅里叶级数表示式，并用 MATLAB 方法求出N=7 和 N=21 时的合成图。解：该信号的系数:F = 0.5Sa( np )n2前 N 项的合成表达式为：f (t) = 0.5 + NNn=1Sa( np ) cos nptn 为奇数2利用 MATLAB 工具分析的程序如下：t=-3:0.001:3;% 定义时间范围向量 tN=input(N=);% 以交互方式输入 N 的值F0=0.5;fN=F0*ones(1,length(t); for n=1:2:NfN=fN+co
28、s(pi*n*t)*sinc(n/2);endplot(t,fN);% 绘制 fN 的波形title(N= num2str(N);axis(-3 3 -0.2 1.2); grid on运行结果如图 14(b)和图 14(c)所示。f(t)1-2.5 -2-1.5 -0.5 0 0.5 1.5 2 2.5t(a)(a)(b)(c)图 14周期矩形脉冲的合成例 10 如图 15 所示周期矩形脉冲 f (t) ，试绘出其频谱图。1f(t)t-4.5 -4 -3.5-0.5 0 0.53.5 4 4.5图 15例 10 图解：程序如下：clear allsyms t n T tao A T=4;A=
29、1;tao=1;f=A*exp(-j*n*2*pi/T*t);fn=int(f,t,-tao/2,tao/2)/T;% 计算傅里叶系数fn=simple(fn);% 化简n=-20:-1,eps,1:20;% 给定频谱的整数自变量,eps 代表 0fn=subs(fn,n,n);% 计算傅里叶系数对应各个 n 的值subplot(2,1,1),stem(n,fn,filled);% 绘制频谱line(-20 20,0 0);% 在图形中添加坐标线title( 周 期 矩 形 脉 冲 的 频 谱 ); grid on subplot(2,1,2),stem(n,abs(fn),filled);%
30、 绘制频谱title(周期矩形脉冲的幅度谱); grid on axis(-20 20 0 0.3);运行结果如图 15 所示。图 15例 10 频谱图例 11 如图 16(a)所示三角波信号，即： f (t) = 1 -1f(t),-2 t 2 ，试求其频谱 F (w ) 。t2t-202(a) (b)图 16例 11 图解：程序如下：syms t w f ft;% 定义符号变量f=(1-(abs(t)/2);% 三角波信号ft=f*exp(-j*w*t);% 计算被积函数F=int(ft,t,-2,2);% 计算傅里叶变换 F(w)F=simple(F);F% 化简axis(-3 3 0
31、1.1);title(三角波信号);ezplot(abs(F),-8:0.01:8); grid on% 绘制三角波信号的频谱title(三角波信号的频谱);运行结果：F =-(cos(2*w)-1)/w2即：1 - cos(2w)2 sin 2 (w)F (w) = = = 2S 2 (w)w 2频谱如图 16(b)所示。w 2a例 12 二阶低通滤波器特性为H (w) =1 w 21 w 1 - w + j Q w 00即：11 w Q w1 - w + w 2 1 w 200QwH (w) =和j(w) = - arctan0 w 2 令Q =12和 1 时，分别求幅频特性和相频特性。
32、1 - w0 解：程序如下：Q=input(输入 Q=);% 以交互方式输入 Q normalizedw=linspace(0.1,10,100);H=1./(1-normalizedw.2+j*normalizedw/Q);% 二阶低通滤波器的频率特性表达式subplot(1,2,1),plot(normalizedw,abs(H);% 绘制幅频特性曲线title(幅频特性曲线);gridsubplot(1,2,2),plot(normalizedw,angle(H);% 绘制相频特性曲线title(相频特性曲线);grid运行结果如图 17 和图 18 所示。输入 Q=1/sqrt(2)输
33、入 Q=1图 17例 12 程序运行结果一图 18例 12 程序运行结果二例 13 三阶低通滤波器特性为H (w) =1( jw)3 + 3( jw)2 + 2( jw) + 1试求：a. 该系统的幅频特性H (w ) 和相频特性j(w) ，b. 该系统的冲激响应h(t) 。解：求幅频特性 H (w )和相频特性j(w) 的程序如下：w=0:0.01:5;H=1./(j*w).3+3*(j*w).2+2*j*w+1);% 三阶低通滤波器的频率特性表达式subplot(1,2,1),plot(w,abs(H);% 绘制幅频特性曲线title( 幅 频 特 性 曲 线 );grid;axis ti
34、ght; subplot(1,2,2),plot(w,angle(H);% 绘制相频特性曲线title(相频特性曲线);grid;axis tight;图 19例 13 程序运行结果一求该系统的冲激响应h(t) 的程序如下：b=1;% 分子多项式系数运行结果如图 19 所示。a=1 3 2 1;% 分母多项式系数impulse(b,a);% 冲激响应 h(t)运行结果如图 20 所示。图 20例 13 程序运行结果二例 14 应用 MATLAB 方法生成信号 f (t) = S (t) p(t)，其中 p(t) 的波形如图 21 所示。a1p(t)t-1.4-1 -0.8 -0.4 -0.2
35、0 0.2 0.4 0.8 11.4图 21例 14 图解：程序如下：t=-3*pi:0.01:3*pi;% 定义时间范围向量s=sinc(t/pi);% 计算 Sa(t)函数subplot(3,1,1),plot(t,s); grid on% 绘制 Sa(t)的波形p=zeros(1,length(t);% 预定义 p(t)的初始值为 0 for i=16:-1:-16p=p+rectpuls(t+0.6*i,0.4);% 利用矩形脉冲函数 rectpuls 的平移来产生宽度为 0.4,幅度为 1 的矩形脉冲序列 p(t)endsubplot(3,1,2),stairs(t,p);% 用阶梯
36、图形表示矩形脉冲axis(-10 10 0 1.2); grid on f=s.*p;subplot(3,1,3),plot(t,f); grid on% 绘制 f(t)=Sa(t)*p(t)的波形运行结果如图22所示。图 22例 14 程序运行结果例 15 分析如图 23 所示三角波信号 f (t) 的取样过程，并画出 f (t), y (t) 和 y (t) 的频谱12图。f (t)H (w )w低通p= 4y (t)1cy(t)d (t)T(a)d (t)Tf(t)1-0.500.5tf (t) = 1- 2 | t |-Ts0Tst Ts=0.2s(b) (c)图 23例 15 图解：
37、程序如下：syms t w f;% 定义符号变量f=(1-2*abs(t)*exp(-j*w*t);% 计算被积函数F=int(f,t,-1/2,1/2);% 计算傅里叶系数 F(w)F=simple(F);F% 化简subplot(3,1,1),% 绘制三角波的幅频特性曲线 F(w)low=-26*pi;high=-low;% 设置 w 的上界和下界ezplot(abs(F),low:0.01:high); grid onaxis(low high -0.1 0.5); xlabel(omega); title(三角波的频谱);subplot(3,1,2),% 绘制经过截止频率为 4*pi
38、低通滤波器后的频谱 Y1(w)ezplot(abs(F),-4*pi:0.01:4*pi); grid onaxis(low high -0.1 0.5); xlabel(omega); title(低通滤波后的频谱);% 取样信号的频谱是原信号频谱的周期延拓，延拓周期为(2*pi)/Ts% 利用频移特性 Ff(t)*exp(-j*w0*t)=F(w+w0)来实现subplot(3,1,3);% 绘制取样后的频谱 Y(w)Ts=0.2;% 取样信号的周期w0=(2*pi)/Ts;% 延拓周期 10*pifor k=-2:2ft=f*exp(-j*w0*k*t);FT=int(ft,t,-1/2
39、,1/2);ezplot(1/Ts)*abs(FT),(-4*pi-k*w0):0.01:(4*pi-k*w0); hold onendaxis(low high -0.1 2.5); xlabel(omega); grid on title(取样后的频谱);运行结果如图 24 所示。图 24例 15 各频谱图五、MATLAB 方法用于复频域分析2s + 1例 16部分分式展开： F (s) = s 3 + 2s 2 + 5s解：程序如下：b=2 1;% 分子多项式系数a=1 2 5 0;% 分母多项式系数r p k=residue(b,a)% 部分分式展开运行结果：r =-0.1000 -
40、0.4500i-0.1000 + 0.4500i0.2000p =-1.0000 + 2.0000i-1.0000 - 2.0000i0k = - 0.1 - 0.45 j- 0.1 + 0.45 j0.2故 F (s) = + + s - (-1 + 2 j)s - (-1 - 2 j)s例 17 求拉氏变换：a. f (t) = e-t cosw tb. f (t) = 3e-2t e (t)解：程序如下：syms t w% 指定 t 和 w 为符号变量fat=exp(-t)*cos(w*t);fbt=3*exp(-2*t);fas=laplace(fat)% 拉氏变换fbs=laplac
41、e(fbt)% 拉氏变换运行结果：fas =(s+1)/(s+1)2+w2) fbs =3/(s+2)即：s + 13F (s) = ， F (s) = a(s + 1)2 + w 2bs + 2例 18 求拉氏反变换：2s + 1s 2a. F (s) = b. F (s) = s 2 + 7s + 10s 2 + 3s + 2解：程序如下：syms s% 指定 s 为符号变量fas=(2*s+1)/(s2+7*s+10); fbs=s2/(s2+3*s+2);fat=ilaplace(fas)% 拉氏反变换fbt=ilaplace(fbs)% 拉氏反变换运行结果：fat =3*exp(-5*t)-exp(-2*t)fbt =Dirac(t)-4*exp(-2*t)+exp(-t)即：f (t) = (3e-5t - e-2t )e (t)fab(t) = d (t) + (-4e-2t + e-t )e (t)例 19 零极点分析s + 2a. H (s) = ，求零极点并画出零极点图，并求阶跃响