第9次课-信号处理算法Matlab仿真案例.ppt

1、第9讲 MATLAB知识回顾,数字信号处理与Matlb仿真,题型,一、填空：10个题，每空2分，共20分； 二、选择：5个题，每题2分，共10分； 三、分析与计算题（共4题, 每题6分， 共24分） 三、综合编程题（共4题, 共46分）,引言,MATLAB软件介绍： MATLAB是矩阵实验室（Matrix Laboratory）的简称，是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。,4,版本更新,5,Matlab 的特点与功能,Matlab 具有很强的数值计算功能,M

2、atlab 以矩阵作为数据操作的基本单位，但无需预先指定 矩阵维数（动态定维）,按照 IEEE 的数值计算标准进行计算,提供十分丰富的数值计算函数，方便计算，提高效率,Matlab 命令与数学中的符号、公式非常接近，可读性强， 容易掌握,6,Matlab 的特点与功能,7,Matlab 丰富的工具箱（toolbox）,根据专门领域中的特殊需要而设计的各种可选工具箱,Matlab 的特点与功能,8,主要应用,数值分析 数值和符号计算 工程与科学绘图 控制系统的设计与仿真 数字图像处理 数字信号处理 通讯系统设计与仿真 财务与金融工程,9,10,界面介绍,11,11,一、命令窗的常用控制指令,二、

3、MATLAB基本操作,变量与矩阵 矩阵及元素的赋值 矩阵运算,13,1、变量-命名规则,变量名是对字母大小写敏感的。 如:变量myvar和MyVar是不同变量。sin是MATLAB定义的正弦名，但SIN、Sin等都不是。 变量名的第一个字母必须是英文字母，最多可包含31个字符（英文、数字和下划线）。 如：myvar201合法。 变量名中不得含空格、标点，但可以有下连符号 如:变量名my_var_201合法的，且读起来更方便。而my,var201不合法。,14,2. 矩阵及元素的赋值,Matlab中的变量或常量都代表矩阵，标量看作11阶的矩阵 赋值语句：变量=表达式（或数） a = 1 2 3;

4、4 5 6;7 8 9 元素的输入可以用表达式代替 x = -1.3, sqrt(3), (1+2+3)/5*4 矩阵的值放在方括号中，同一行中各元素之间用逗号或空格分开，不同的行则以分号隔开，此时会立即显示运算结果。如果不希望显示结果，以分号结尾，此时结果在workspace可以看到。,15,变量的元素用( )中的数字标明，一维矩阵中的元素用一个下标表示，二维矩阵可有两个下标数，以逗号隔开。三维或更高维的矩阵可有三个或更多下标。 可单独给元素赋值 x(2) = 1.7321 a(2,3) = 6 linspace(a,b,n) %在a和b之间均匀,3. 矩阵运算,（一）预定义变量,（二）运算

5、符,+ - * / sqrt = “点乘”与普通乘法的区别； .* and * 2. 运算优先级：指数运算乘除加减； 3. “=” 和其他运算符两侧，允许有空格，以增加程序可读性。,（三）矩阵（数组）生成函数,19,随机矩阵,为了进行信号分析与处理，经常需要对接收信号进行仿真，而在信号仿真中离不开随机数的产生，MATLAB提供的rand和randn函数可分别产生均匀分布和正态分布的随机数。例如，要产生0，1之间均匀分布的随机向量r1(23)，可输入： r1=rand(2, 3)产生0，1之间均匀分布的随机矩阵 r1 = 0.4565 0.8214 0.6154 0.0185 0.4447 0.

6、7919,20, r2=5-10*rand(2, 3)产生-5，5之间均匀分布的随机矩阵 r2 = 0.8973 4.4211 -3.1317 -3.9365 1.4713 4.9014 r3=randn(2,3)产生均值为0、方差为1的标准正态分布的随机矩阵 r3 = 1.1892 0.3273 -0.1867 -0.0376 0.1746 0.7258,21,通过下标引用矩阵的元素，例如A(3,2)=200 矩阵元素按列存储，例如 A=1,2,3;4,5,6; A(3) ans = 2 显然，序号(Index)与下标(Subscript )是一一对应的，以mn矩阵A为例，矩阵元素A(i,j

7、)的序号为(j-1)*m+i。其相互转换关系也可利用sub2ind和ind2sub函数求得。,矩阵元素,（四）矩阵操作,A (r,c)A第r行第c列元素 A (r,: )A第r行所有列元素 A ( :, c)A所有行第c列元素 A(:,k:m) A的第 k 到第 m 列组成的子矩阵 A ( : )矩阵所有元素排列为一维向量（从左到右按列连接） A (s)按序号寻访元素 A (r,c) = Kr行c列元素赋值为K A( : ) = D( : )A、D两矩阵总元素数相同（行列不一定），将D每个元素赋给A中每个元素,A(:) 与 A(:,:) 的区别 ?,如何获得由 A 的第一、三行和第一、二列组成

8、的子矩阵？,23,23,24,查看矩阵的大小：size与length,size(A) 列出矩阵 A 的行数和列数,size(A,1) 返回矩阵 A 的行数,size(A,2) 返回矩阵 A 的列数,length(x) 返回向量 X 的长度,length(A) 等价于 max(size(A),25,改变矩阵的形状：reshape,reshape(A,m,n): 将矩阵元素按 列方向 进行重组,重组后得到的新矩阵的元素个数 必须与原矩阵元素个数相等！,26,矩阵基本运算,矩阵的加减：对应分量进行运算,要求参与加减运算的矩阵具有 相同的维数,矩阵的普通乘法,要求参与运算的矩阵满足线性代数中矩阵相乘的

9、原则,27,矩阵的数组运算,数组运算：对应元素进行运算,点与算术运算符之间不能有空格！,数组运算包括：点乘、点除、点幂,相应的数组运算符为： “.* ” ， “./ ” ， “. ” 和“ . ”,参与运算的对象必须具有相同的形状！,28,线性代数运算的MATLAB命令 MATLAB是矩阵化程序设计语言，所以处理矩阵和向量运算特别方便。关于矩阵和向量的一些基本运算命令已在前面有所介绍，常用的命令和函数还有,29,Matlab中常见数学函数,log 是自然对数，即以 e 为底数 mod(x,y) 结果与 y 同号，rem(x,y) 则与 x 同号 max 等函数的参数是矩阵时，是作用在矩阵各列上

10、,二、Matlab绘图,作为一个功能强大的工具软件，Matlab具有很强的图形处理功能，提供了大量的二维、三维图形函数。 由于系统采用面向对象的技术和丰富的矩阵运算，所以在图形处理方面既方便又高效。,31,1 概述及图形窗口的控制,【1】figure命令 格式：figure 或 figure(N) 功能：创建一个图形窗口。使编号为N的图形窗口成为当前图形窗口，即图形窗口处于可视状态。如果窗口N不存在，则将创建一个句柄为N的图形窗口。,31,32,1 概述及图形窗口的控制,【2】clf命令 格式：clf 功能：清除当前图形窗口中的所有的内容 例 x=-pi:pi/10:pi; y=sin(x);

11、 plot(x,y) clf,32,33,1 概述及图形窗口的控制,【3】close命令 格式一：close 功能： 关闭当前图形窗口 格式二： close (N) 功能： 关闭指定编号N的图形窗口 格式三：close all 功能： 关闭所有图形窗口,34,二维图形,plot函数 函数格式：plot(x,y) 其中x和y为坐标向量 函数功能：以向量x、y为轴，绘制曲线。 例： 在区间0 x2内，绘制正弦曲线y=sin（x），其程序为： x=0:pi/100:2*pi; y=sin(x); plot(x,y),35,线型与颜色 格式：plot(x,y1,cs,.) 其中c表示颜色， s表示线型

12、。,例：x=0:pi/100:2*pi; y1=sin(x); y2=cos(x); plot(x,y1,g-,x,y2,b-.) 其中参数g-和b-.表示图形的颜色和线型。g表示绿色，-表示图形线型为实线；b表示蓝色，-.表示图形线型为点划线。,36,设置曲线的样式属性,37,figure, t=0:pi/20:2*pi; y=sin(t); y1=sin(t-pi/2); y2=sin(t-pi); plot(t,y, c-v ,t,y1,k-s,t,y2,m:o),38,若使用plot函数时不指定线型，则绘制的曲线将不被连接起来 figure, t=0:pi/20:2*pi; y=sin

13、(t); y1=sin(t-pi/2); y2=sin(t-pi); plot(t,y,cv ,t,y1,ks,t,y2,mo),38,39,使用命令grid，可以将图形窗体的轴在有网格线和无网格线之间切换 grid on图形窗体的轴显示坐标网格线 grid off图形窗体的轴不显示坐标网格线,39,40, t=0:pi/20:2*pi; y=sin(t); y1=sin(t-pi/2); y2=sin(t-pi); plot(t,y,cv ,t,y1,ks,t,y2,mo) grid on grid off,41,图形标记 在绘制图形的同时，可以对图形加上一些说明，如图形名称、图形某一部分的

14、含义、坐标说明等，将这些操作称为添加图形标记。 title(加图形标题); xlabel(加X轴标记); ylabel(加Y轴标记);,42,设定坐标轴 用户若对坐标系统不满意，可利用axis命令对其重新设定。 axis(xmin xmax ymin ymax) 设定最大和最小值 axis (auto) 将坐标系统返回到自动缺省状态 axis (square) 将当前图形设置为方形 axis (equal) 两个坐标因子设成相等 axis (off) 关闭坐标系统 axis (on) 显示坐标系统,43,例： 在坐标范围0X2,-2Y2内重新绘制正弦曲线，其程序为： x=linspace(0,

15、2*pi,60); %生成含有60个数据元素的向量X y=sin(x); plot(x,y); axis (0 2*pi -2 2);%设定坐标轴范围,44,加图例 给图形加图例命令为legend。该命令把图例放置在图形空白处，用户还可以通过鼠标移动图例，将其放到希望的位置。 格式:legend(图例说明,图例说明);,例： 为正弦、余弦曲线增加图例，其程序为： x=0:pi/100:2*pi; y1=sin(x); y2=cos(x); plot(x,y1,x,y2, -); legend(sin(x),cos(x);,45,subplot（m,n,p） 该命令将当前图形窗口分成mn个绘图区

16、，即每行n个，共m行，区号按行优先编号，且选定第p个区为当前活动区。,subplot函数,46,x=linspace(0,2*pi,60); y=sin(x); z=cos(x); t=sin(x)./(cos(x)+eps); %eps为系统内部常数 ct=cos(x)./(sin(x)+eps); subplot(2,2,1); plot(x,y); title(sin(x); axis (0 2*pi -1 1); subplot(2,2,2); plot(x,z); title(cos(x); axis (0 2*pi -1 1); subplot(2,2,3); plot(x,t);

17、 title(tangent(x); axis (0 2*pi -40 40); subplot(2,2,4); plot(x,ct); title(cotangent(x); axis (0 2*pi -40 40);,例： 在一个图形窗口中同时绘制正弦、余弦、正切、余切曲线,47,多图形窗口 需要建立多个图形窗口，绘制并保持每一个窗 口的图形，可以使用figure命令。 每执行一次figure命令，就创建一个新的图形窗口，该窗口自动为活动窗口，若需要还可以返回该窗口的识别号码，称该号码为句柄。句柄显示在图形窗口的标题栏中，即图形窗口标题。用户可通过句柄激活或关闭某图形窗口，而axis、xl

18、abel、title等许多命令也只对活动窗口有效。,48,重新绘制上例4个图形，程序变动后如下： x=linspace(0,2*pi,60); y=sin(x); z=cos(x); t=sin(x)./(cos(x)+eps); ct=cos(x)./(sin(x)+eps); figure; plot(x,y); title(sin(x); axis (0 2*pi -1 1); figure; plot(x,z); title(cos(x);,axis (0 2*pi -1 1); figure; plot(x,t); title(tangent(x); axis (0 2*pi -40

19、 40); figure; plot(x,ct); title(cotangent(x); axis (0 2*pi -40 40);,49,hold命令 若在已存在图形窗口中用plot命令继续添加新的图形内容，可使用图形保持命令hold。发出命令hold on后，再执行plot命令，在保持原有图形或曲线的基础上，添加新绘制的图形。,50,阅读如下程序： x=linspace(0,2*pi,60); y=sin(x); z=cos(x); plot(x,y,b); %绘制正弦曲线 hold on; %设置图形保持状态 plot(x,z,r); %保持正弦曲线同时绘制余弦曲线 axis (0 2

20、*pi -1 1); legend(sin,cos); hold off %关闭图形保持,51,52,二维绘图函数小结,plot 二维图形基本函数 fplot f(x)函数曲线绘制 fill 填充二维多边图形 polar 极坐标图 bar 条形图 loglog 双对数坐标图 semilogx X轴为对数的坐标图 semilogy Y轴为对数的坐标图 stairs 阶梯形图 axis 设置坐标轴 clf 清除图形窗口内容 close 关闭图形窗口,figure 创建图形窗口 grid 放置坐标网格线 gtext 用鼠标放置文本 hold 保持当前图形窗口内容 subplot 创建子图 text

21、放置文本 title 放置图形标题 xlabel 放置X轴坐标标记 ylabel 放置Y轴坐标标记,53,三维图形,plot3函数 最基本的三维图形函数为plot3，它是将二维函数plot的有关功能扩展到三维空间，用来绘制三维图形。 函数格式：plot3(x1,y1,z1,c1,x2,y2,z2,c2,) 其中x1,y1,z1表示三维坐标向量，c1，c2表示线形或颜色。 函数功能：以向量x，y，z为坐标，绘制三维曲线。,54,绘制三维网线图,mesh(z)： z为mn的矩阵 Z坐标矩阵中元素的值 X与Y坐标矩阵中元素的下标 X1：n，对应矩阵的列 Y1：m，对应矩阵的行,55,图像的读/写,用

22、imread函数读取图像，调用格式为： x=imread(filename.fmt) filename是图像的文件名，fmt指定图像文件的格式。,图像读、写与显示,图像的读/写,图像读、写与显示,图像的读/写,56,a=imread(jlu.jpg); size(a),运行结果,ans = 208 208 3,如：,57,在MATLAB中，imshow函数用于图像显示 调用格式为： imshow(x),例子： x=imread(lena.jpg); figure, imshow(x);,图像读、写与显示,图像的显示,三、 Matlab程序设计,三、 Matlab程序设计,常用指令,60,流程控

23、制语句,四种流程语句： if语句 switch语句 while语句 for语句,选择结构,循环结构,61,（一）选择结构,选择结构的语句有if语句和switch语句。 1 if语句 格式一： if 条件 语句组 end %必不可少，若没有则表达式为0时找不到 %继续执行的程序入口 格式二： if 条件 语句组1 else 语句组2 end,62,格式三： if 条件1 语句组1 elseif 条件2 语句组2 elseif 条件m 语句组m else 语句组m+1 end,63,例： 输入三角形的三条边，求面积。,A=input(请输入三角形的三条边：); if A(1)+A(2)A(3) e

24、lse disp(不能构成一个三角形。) end 运行： 请输入三角形的三条边：4 5 6 9.9216,64,（二）循环结构,实现循环结构的语句：for语句和while语句。 1. for语句： 格式： for 循环变量=表达式1:表达式2:表达式3 循环体语句 end 表达式1的值为循环变量的初值; 表达式 2的值为步长; 表达式3的值为循环变量的终值。 步长为1时，表达式2可以省略。,65,例：计算1+3+5+99=?,sum = 0; for n = 1:2:99 sum = sum +n; end sum 运行结果：sum = 2500,66,例： 当n=100时，求 的值。,程序如

25、下: tic y=0; n=100; for i=1:n y=y+1/(2*i-1); end toc y,y = 3.2843,67,在实际MATLAB编程中，采用循环语句会 降低其执行速度，所以前面的程序通常由下面的程序来代替：,tic n=100000; k=1:2:(2*n-1); y=sum(1./k); toc y,y = 3.2843,67,68,2. while语句 格式为： while (条件) 循环体语句 end,例：计算110的和 i=1; sum=0; while(i=10) sum = sum+i; i = i+1; end sum,运行结果：55,69,M文件及程序

26、调试,在入门阶段，通常在行命令模式下工作，键入一条命令后，让系统立即执行该命令。该方法程序可读性差且难以存储。对于复杂的问题应编成可存储的程序文件，让Matlab执行该程序文件，这种工作模式成为程序文件模式。 由Matlab语句构成的程序文件称为m文件，用m作为文件的扩展名，可直接阅读并用任何编辑器来建立。 m文件分主程序和子程序 主程序：由用户为解决问题而编写的 子程序：函数文件，必须由其他m文件来调用,70,主程序：用clear、close all等语句开始，清除工作空间中原有的变量和图形，以避免其他已执行的程序残留数据对本程序的影响。 主程序中如有全局变量，即在子程序中与主程序共用的变量

27、，在程序起始部分注明 global 变量名1 变量名2 整个程序按Matlab标识符的要求起名，加上后缀.m,71,函数文件,用来定义子程序 由function起头，后跟的函数名必须与文件名相同 有输入输出变量，可进行变量传递 除非用global声明，子程序中的变量均为局部变量，不保存在工作空间中,72,function y = mean_my(x) %定义函数名，输入和输出变量 m,n=size(x); if m=1 m=n; end y=sum(x)/m m，n，y为局部变量 在主程序中调用：x = a = mean_my(x) ,73,变量的存储和下载 存储：save：把工作空间的变量存

28、入磁盘，.mat save aa %把工作空间的全部变量存入保存在名为aa.mat文件中; save aa a b c %把a，b，c三个变量保存在名为aa.mat文件中; 下载：load，把mat数据文件传回到工作空间中，但不能选择变量 load aa %把aa.mat文件中的全部数据和变量名一起下载到工作空间中,74,四、离散傅立叶变换(DFT),75,定义：有限长序列x(n)(0nN-1)有N个样本值。它的傅立叶变换X(j)在频率区间(02)的N个等间隔分布的点k=2k/N(0kN-1)上也有N个取样值。 这两组同样长的序列之间可以用简单的关系联系起来： 其中：,1、DFT定义,76,它

29、的反变换称为离散傅立叶反变换(IDFT)，可以表为：,1、DFT定义,77,若用矩阵式 X=W*x 及 x= W-1*X求 变换矩阵：W=dftmtx(N)， 反变换矩阵：W-1=conj(dftmtx(N)/N， 但通常可直接调用FFT程序： 正变换：X=fft(x,N) X= fft(x) N默认为x的长度 反变换：X=ifft(x)或X= ifft(x,N),4、用工具箱函数计算DFT,78,FFT程序求频谱的范围规范化为0,2)，对应的频谱下标为k=0,1,N-1。设频率分辨率dw=2/N，则频率向量为 w=k*dw 若要画出乃奎斯特频率范围-,)内的频谱， 需要把,2)范围的频谱平移

30、到-,0)内，对于这一运算，有专门的工具箱函数fftshift来实现。 X1=fftshift(X) X1就是在-,)内的频谱向量。,4、用工具箱函数计算DFT,79,傅里叶变换函数,常用的FFT及反变换函数如下表所示。,80,80,FFT计算频谱举例,例：选择合适的变换区间长度N，用DFT对下列信号进行谱分析，画出幅频、相频特性。,解：(1)为周期序列，取一个周期用DFT谱分析,% (1) 周期序列,截取一个周期作谱分析 N1 = 100; n1 = 0: N1-1; x1n = 2 * cos(0.2 * pi * n1); % 生成序列x1n X1k = fft(x1n, N1); %

31、x1n的N1点DFT subplot(1,2,1);stem(n1/N1 *pi, abs(X1k),.); title(X1k - Amplitude); subplot(1,2,2);stem(n1/N1 *pi, angle(X1k),.); title(X1k - Phrase);,81,例4.1 求x=2,-1,1,1的离散傅立叶变换。并画到图上。 解： hold on Xd=fft(2,-1,1,1) %其DFT只有四点 plot(0:3*2*pi/4,abs(Xd),o) Xd1=fftshift(Xd); % 画出移位后的FFT plot(-2:1*2*pi/4,abs(Xd1),x),4、用工具箱函数计算DFT,82,上述用频域DFT相乘再求反变换方法3可以用下列方框图来表示： 这样求出的是循环卷积。可以证明，如果把N选得足够大，NN1N2-1，进行N点循环卷积，就有可能避免头尾混叠，结果就等于线性卷积。,3.7 用DFT计算线性卷积,83,例1：已知一个8点的时域非周期离散阶跃信号，