MatLab_基础篇

1、MatLab编程基础编程基础 MatLab编程基础MatLab编程基础编程基础 MatLab编程基础 Mathematica或Maple的使用MatLab编程基础参考文献参考文献理论力学计算机模拟，清华大学出版社，彭芳麟等MatLab工程数学应用，许波，刘征编著Mathematica 4.0使用教程，刘元高，刘耀儒Maple计算机代数系统应用及程序设计，李世奇，杜慧琴MatLab编程基础第一章第一章 MatLab简介简介第一节 预备知识第二节 矩阵与表达式第三节 编程第四节 作图与动画第五节 常用的计算指令第六节 工具箱MatLab编程基础第一节 预备知识MatLab是矩阵工作室（Matrix

2、 Laboratory）的缩写MatLab的特点n极强的数值计算功能和作图功能n图形窗口式的操作界面n用简单的指令可以完成大量的计算和作图n有大量的不同领域的专用工具箱，用户也可以开发自己的专用工具箱MatLab编程基础第一节 预备知识操作界面 Desktop;指令窗口指令窗口指令记录指令记录当前工作当前工作目录目录当前内存空间浏览器当前内存空间浏览器矩阵矩阵双精度双精度大小写大小写使用方法使用方法MatLab编程基础第一节 预备知识例例 1+2 2*3 4/5ans = 3.0000 6.0000 0.8000 2/3ans = 0.6667 format long 2/3ans = 0.6

3、6666666666667MatLab编程基础第一节 预备知识在线帮助在指令窗口中直接使用在指令窗口中直接使用：help+指令名称指令名称MatLab编程基础第一节 预备知识使用菜单中的使用菜单中的helpFull Product Family HelpMATLAB Help F1MatLab编程基础第二节 矩阵与表达式数据格式算符表达式变量函数矩阵其它数据结构MatLab编程基础第二节 矩阵与表达式数据格式使用常规的十进制表达法，小数位数不限，可以使用常规的十进制表达法，小数位数不限，可以用加号和减号表示正负数用加号和减号表示正负数 1010的幂次表示为的幂次表示为e e或或E E加上正负数

4、表示加上正负数表示默认的虚数单位为默认的虚数单位为i i和和j j，数字后面直接写，数字后面直接写i i或或j j表表示虚数，也可以示虚数，也可以* *i i或或j jMatLab编程基础第二节 矩阵与表达式算符加减乘除幂括号+-*/()算术运算算符算术运算算符关系运算算符关系运算算符逻辑运算算符逻辑运算算符与或非&|小于 小于等于 等于 大于 大于等于 不等于=MatLab编程基础第二节 矩阵与表达式表达式用算符把数、变量和函数组合在一起，就建立了用算符把数、变量和函数组合在一起，就建立了一个表达式。一个表达式。 可以对整个矩阵进行运算。可以对整个矩阵进行运算。 表达式中运算常规的优

5、先级自左至右执行，优先表达式中运算常规的优先级自左至右执行，优先级的顺序为指数运算最优先，其次是乘除，最后级的顺序为指数运算最优先，其次是乘除，最后是加减，可以用括号改变运算顺序。是加减，可以用括号改变运算顺序。 MatLab编程基础第二节 矩阵与表达式变量变量名用字母打头，后面可以跟字母、数字、下变量名用字母打头，后面可以跟字母、数字、下划线，数目不限。但是只有变量名的前面划线，数目不限。但是只有变量名的前面3131个符个符号有效。号有效。 变量名可以用来表示具体的数据（如矩阵），这变量名可以用来表示具体的数据（如矩阵），这称之为对一个变量赋值。默认为双精度的复数。称之为对一个变量赋值。默认

6、为双精度的复数。清除所有变量用指令清除所有变量用指令 clear allclear allMatLab编程基础第二节 矩阵与表达式常用的常量pi 3.14159265358979pi 3.14159265358979i,j i,j 复数单位复数单位Inf Inf 无穷大无穷大epseps 2(-52) 2(-52)2.2204e-0162.2204e-016Nan Nan 非数非数realminrealmin 最小浮点数最小浮点数 2(-1022)2(-1022) 2.22512.2251e-308e-308realmaxrealmax 最大浮点数最大浮点数 2(1024)2(1024) 1.

7、7977 1.7977e+308e+308MatLab编程基础第二节 矩阵与表达式函数基本函数列表基本函数列表 help elfunhelp elfun特殊函数列表特殊函数列表 help specfunhelp specfun使用菜单中的使用菜单中的 helphelpMatLab编程基础第二节 矩阵与表达式矩阵定义矩阵：定义矩阵：A=2 2 3 6;2 12 5 6;23 6 8 11A= 2 2 3 6 2 12 5 6 23 6 8 11x=0:pi/5:pix= 0 0.6283 1.2566 1.8850 2.5133 3.1416MatLab编程基础第二节 矩阵与表达式常用的生成矩阵

8、指令：常用的生成矩阵指令：zeros 零矩阵ones 全部元素为1的矩阵eye 单位矩阵rand 均匀分布的随机数矩阵MatLab编程基础第二节 矩阵与表达式例：例：exp 与与 expm 的区别的区别借助于借助于help exp来查询来查询MatLab编程基础第二节 矩阵与表达式矩阵的识别与修改：矩阵的识别与修改：矩阵操作的指令和算符：矩阵操作的指令和算符：共轭转置加减乘右除 左除幂矩阵运算算符A+-*/数组运算算符A+-.*./.标量作用指令；矢量作用指令；矩阵作用指令行标与列标；合并；删去行与列MatLab编程基础第二节 矩阵与表达式其它数据结构借助于借助于help datatypes来

9、查询来查询列阵：行、列、层、块、堆、列阵：行、列、层、块、堆、三维列阵、四维列阵、五维列阵、三维列阵、四维列阵、五维列阵、cat(dim,A1,A2,A3,) repmat(A,m,n,p,) MatLab编程基础数据网格数据网格第二节 矩阵与表达式X,Y=meshgrid(x,y) X,Y=meshgrid(x) X,Y,Z=meshgrid(x,y,z) 生成二维或三维数据生成二维或三维数据生成生成n维数据维数据X1,X2,X3,=ndgrid(x1,x2,x3,)X1,X2,X3,=ndgrid(x)MatLab编程基础第二节 矩阵与表达式x=0 1 2 3 4;X,Y=meshgrid

10、(x) 二维数据网格示例二维数据网格示例X = 0 1 2 3 4 0 1 2 3 4 0 1 2 3 4 0 1 2 3 4 0 1 2 3 4Y = 0 0 0 0 0 1 1 1 1 1 2 2 2 2 2 3 3 3 3 3 4 4 4 4 4MatLab编程基础第二节 矩阵与表达式x=0 1 2 3 4;X,Y,Z=meshgrid(x) X,Y,Z=ndgrid(x)三维数据网格示例三维数据网格示例MatLab编程基础第二节 矩阵与表达式基元列阵基元列阵可以将不同类型的数据按照与矩阵可以将不同类型的数据按照与矩阵形式相似的结构组织起来加以应用形式相似的结构组织起来加以应用结构数据结

11、构数据比基元列阵更复杂的数据结构，不比基元列阵更复杂的数据结构，不仅可以将不同类型的数据组织在一仅可以将不同类型的数据组织在一起，还可以为它们赋予名称，是一起，还可以为它们赋予名称，是一种有域名的数据种有域名的数据S=struct(field1,values1, field2,values2,)MatLab编程基础第三节 编程程序文件的编辑与调试函数类文件流程控制数据输入与输出MatLab编程基础第三节 编程-程序文件的编辑与调试程序文件的编辑与调试程序编辑器编辑M文件程序运行及调试设置路径MatLab编程基础程序编辑器第三节 编程-程序文件的编辑与调试% 路径及文件名路径及文件名文件名文件名

12、程序注释程序注释MatLab编程基础第三节 编程-程序文件的编辑与调试% clear allA=1 2 12 3 3 5 8 9 11 2 7 9 15 9 11 8;inv(A)det(A)路径及文件名路径及文件名程序编辑器窗口程序编辑器窗口程序注释程序注释编辑M文件文件名文件名MatLab编程基础第三节 编程-程序文件的编辑与调试程序运行及调试% clear allA=1 2 12 3 3 5 8 9 11 2 7 9 15 9 11 8;inv(A)det(A)路径及文件名路径及文件名文件名文件名程序编辑器窗口程序编辑器窗口运行运行文件名文件名 回车回车ansans = = -0.010

13、6 -0.0946 0.0614 0.0412 -0.0106 -0.0946 0.0614 0.0412 -0.0647 0.0942 -0.1716 0.1113 -0.0647 0.0942 -0.1716 0.1113 0.1095 -0.0391 0.0061 -0.0040 0.1095 -0.0391 0.0061 -0.0040 -0.0578 0.1250 0.0694 -0.0720 -0.0578 0.1250 0.0694 -0.0720ansans = = 6039 6039MatLab编程基础第三节 编程-程序文件的编辑与调试设置路径File Set path Ma

14、tLab编程基础第三节 编程-函数类文件函数类文件函数M文件的结构函数M文件的调用子函数文件全局变量MatLab编程基础第三节 编程-函数类文件函数M文件的结构function 输出变量输出变量Q1,Q2,=函数文件名函数文件名(输入变量输入变量P1,P2,)定义函数行；H1行；函数帮助文本；函数体；注释必须是文件名MatLab编程基础第三节 编程-函数类文件：举例例例1：定义一个二元函数定义一个二元函数% define a function of two variablesfunction yy=fun_1(x,y)yy=2*exp(-(x-2).2-y.2)+exp(-(x+2).2-y.

15、2)-3*exp(-x.2-(y+1).2);例例2：定义定义44的矩阵的矩阵% define a 4-by-4 matrixfunction yy=fun_2yy=1 2 12 3;3 5 8 9;11 2 7 9;15 9 11 8;MatLab编程基础第三节 编程-函数类文件函数M文件的调用输出变量输出变量Q1,Q2,=函数文件名函数文件名(输入变量输入变量P1,P2,)例：例：inv(fun_2)ans = -0.0106 -0.0946 0.0614 0.0412 -0.0647 0.0942 -0.1716 0.1113 0.1095 -0.0391 0.0061 -0.0040

16、-0.0578 0.1250 0.0694 -0.0720fun_1(1,1)ans = 0.2505det(fun_2)ans = 6039MatLab编程基础第三节 编程-函数类文件子函数文件在同一个函数文件中可以建立多个函数，第一个叫主函数，其余的叫子函数。子函数只能被同一个文件中的主函数和其它子函数调用。MatLab编程基础第三节 编程-函数类文件全局变量如果几个函数文件要共用一个变量，那么要在这些函数文件中都定义这个变量是全局变量。如果要在指令窗口中也使用这个变量，还要在指令窗口中定义这个变量为全局变量。必须在使用函数之前定义全局变量。global 变量1 变量2 MatLab编程基

17、础第三节 编程-函数类文件例：例： ex_5.mfunction yy=fun_3yy=2 3;3 4;function yy=fun_4global lambdayy=lambda*eye(size(fun_3)-fun_3;clear allglobal lambdalambda=3;D=det(fun_4)MatLab编程基础第三节 编程-流程控制流程控制循环结构分支结构其它指令MatLab编程基础第三节 编程-流程控制循环结构for循环结构while循环结构for 循环变量=起始值:步长:终止值 循环体endwhile 表达式 循环体end用于不能确定循环次数的循环。当表达式为真时，执

18、行循环体；执行后在判断表达式是否为真，是真则执行循环体，否则推出。用于能确定循环次数的循环。步长的默认值为1，可以是任意的正实数或负实数。MatLab编程基础第三节 编程-流程控制：举例例1 ex_1.m：对角矩阵的赋值clear allm=5;for i=1:m a(i,i)=i;endaa = 1 0 0 0 0 0 2 0 0 0 0 0 3 0 0 0 0 0 4 0 0 0 0 0 5MatLab编程基础第三节 编程-流程控制：举例例2 ex_2.m：计算m的阶乘clear allm=10;a=1;for i=1:m a=i*a;endaa = 3628800MatLab编程基础第三

19、节 编程-流程控制：举例例3 ex_3.m：矩阵的赋值a = 1.0000 0.5000 0.3333 0.5000 0.3333 0.2500 0.3333 0.2500 0.2000clear allm=3;n=3;for i=1:m for j=1:n a(i,j)=1/(i+j-1); endendaMatLab编程基础第三节 编程-流程控制：举例例4 ex_4.m：矩阵指数的计算e = 162.7871 224.6754 224.6754 312.5707clear alle=zeros(size(fun_3);f=eye(size(fun_3);k=1;while norm(f,1

20、)0 e=e+f; f=fun_3*f/k; k=k+1;endefunction yy=fun_3yy=2 3;3 4;注：对比expm指令MatLab编程基础第三节 编程-流程控制分支结构if分支结构switch分支结构if 逻辑表达式1 语句体1elseif 逻辑表达式2 语句体2else 语句体3endswitch 表达式(标量或字符串) case 值1 语句1 case 值2 语句2 otherwise 语句3end在语句执行时，若逻辑表达式1为真，则执行语句体1；否则，判断逻辑表达式2是否为真，若真执行语句体2，否则执行语句体3。若表达式与值1相符，则执行语句1；若与值2相符，则执

21、行语句2；依此类推。如果都不符合则执行语句3。MatLab编程基础第三节 编程-流程控制：举例例1 pdbsline.m：B样条函数pdbsline(-1)ans = 0pdbsline(0.6)ans = 0.6000pdbsline(1.36)ans = 0.6400pdbsline(2.6)ans = 0function f=pdbsline(x)if x0 f=0;elseif x1 f=x;elseif x2 f=2-x;else f=0;endMatLab编程基础第三节 编程-流程控制：举例例2 ex_6.m：单位换算clear allx=2.7;units=in;switch u

22、nits case inch,in y=x*2.54; case feet,ft y=x*2.54/12; case meter,m y=x/100; y = 6.8580 case millimeter,mm y=x*10; case centimeter,cm y=x; otherwise disp(unknown units: units) y=nan;endyMatLab编程基础第三节 编程-流程控制其它指令break用于for,while,if语句的终止。如果循环是嵌套的，则只能从最内的一层退出。pause使程序运行暂停，按任意键恢复运行，而pause(n)则是暂停n秒。input将

23、用户从键盘输入的数字、字符串或表达式赋予指定的变量。disp()在屏幕上显示引号中的内容。keyboard暂时终止程序的运行，等候键盘输入的指令。在执行完输入的指令后，键入return，就可恢复程序的运行。MatLab编程基础第三节 编程-数据输入与输出数据输入与输出导出save 路径文件名.dat 内存变量 -ascii 导入load 路径文件名.datMatLab编程基础第四节 作图与动画二维图形三维图形四维数据的表现动画MatLab编程基础第四节 作图与动画-二维图形二维图形plot(x,y)基本作图命令格式基本作图命令格式指定颜色、线型和数据点标志的命令指定颜色、线型和数据点标志的命令

24、subplot(m,n,p)分区作图分区作图plot(x,y,color style marker)MatLab编程基础第四节 作图与动画-二维图形数据点标志数据点标志:颜色颜色:青青c；洋红；洋红m ；黄；黄y ；红；红r ；绿；绿g ；蓝；蓝b ；白；白w ；黑；黑k线型线型:实线实线-；短划线；短划线- ；点线；点线: ；点划；点划线线.- 点点 .；加号；加号+ ；圆圈；圆圈o ；星号；星号*；叉号；叉号 ；四方形；四方形s ；钻石型；钻石型d ；五角形；五角形p；六角形六角形 h；三角形（向上；三角形（向上 ；向下；向下 ；向左向左）MatLab编程基础第四节 作图与动画-二维图形:

25、 举例例例1t=0:pi/100:2*pi;y=sin(t);plot(t,y)t=0:0.1:2*pi;y1=sin(t-0.25);y2=sin(t-1);plot(t,y1,b+,t,y2,ro)0246-1-0.500.510246-1-0.500.51例例2MatLab编程基础t=0:0.1:2*pi;y1=sin(t);y2=cos(t);plot(t,y1,b+,t,y2,ro)legend(sine,cosine)例例3第四节 作图与动画-二维图形: 举例0246-1-0.500.51sinecosineMatLab编程基础第四节 作图与动画-二维图形: 举例t=0:0.1:2

26、*pi;y1=sin(t);y2=cos(t);subplot(2,1,1)plot(t,y1,ro)subplot(2,1,2)plot(t,y2,r*)例例4MatLab编程基础例例5 复数作图第四节 作图与动画-二维图形: 举例t=0:pi/10:2*pi;plot(exp(i*t),-o)axis equal-1-0.500.51-1-0.500.51MatLab编程基础第四节 作图与动画-二维图形: 举例例例6 雪花分形图雪花分形图的生成：将等边三角形的每一条边三等分，将中间的一段去掉，代之以一个更小的等边三角形的两条边，依次做下去，即可得到雪花分形图MatLab编程基础第四节 作图

27、与动画-二维图形: 举例snow.mclear allnew=0.5+sqrt(3)/2*i,-0.5+sqrt(3)/2*i,0,0.5+sqrt(3)/2*i;subplot(2,3,1); plot(new); axis equalfor k=1:5 old=new; n=length(old)-1; for jj=0:n-1 diff=(old(jj+2)-old(jj+1)/3; new(4*jj+1)=old(jj+1); new(4*jj+2)=old(jj+1)+diff; new(4*jj+3)=new(4*jj+2)+diff*(1-sqrt(3)*i)/2); new(4

28、*jj+4)=old(jj+1)+2*diff; end new(4*n+1)=old(n+1); subplot(2,3,k+1); plot(new); axis equalendMatLab编程基础-0.500.500.51-0.500.500.51-0.500.500.51-0.500.500.51-0.500.500.51-0.500.500.51第四节 作图与动画-二维图形: 举例MatLab编程基础第四节 作图与动画-二维图形: 举例例例7 树形分形图树形分形图的生成：将一条线段三等分，在等分点上各画出一条长度为原线段长度三分之一的线段，该线段与原线段成固定的夹角，依次做下去，即

29、可得到树形分形图MatLab编程基础第四节 作图与动画-二维图形: 举例clear alltheta=pi/6; u=0,0;0,1; subplot(2,3,1); plot(u(:,1),u(:,2); axis(-0.5,0.5,0,1)for n=1:5 uuu=; for i=0:length(u)/2-1 p1=(u(2*i+1,:)*2+u(2*i+2,:)/3; p2=(u(2*i+1,:)+u(2*i+2,:)*2)/3; lp=cos(theta),-sin(theta);sin(theta),cos(theta)*. (u(2*i+2,1)-u(2*i+1,1);u(2*

30、i+2,2)-u(2*i+1,2)/3; lp=p1+lp; rp=cos(theta),sin(theta);-sin(theta),cos(theta)*. (u(2*i+2,1)-u(2*i+1,1);u(2*i+2,2)-u(2*i+1,2)/3; rp=p2+rp; uu=u(2*i+1,:);p1;p1;lp;p1;p2;p2;rp;p2;u(2*i+2,:); uuu=uuu;uu; end u=uuu; subplot(2,3,n+1); plot(u(:,1),u(:,2); axis(-0.5,0.5,0,1)endtree.mMatLab编程基础-0.500.500.51

31、-0.500.500.51-0.500.500.51-0.500.500.51-0.500.500.51-0.500.500.51第四节 作图与动画-二维图形: 举例MatLab编程基础第四节 作图与动画-二维图形: 举例例例8 虫口模型-Logistic模型)(21nnnxxx其中是与虫口增长有关的参数，一般称为控制参数。该模型是生物学家May于1976年提出的。设某昆虫第n代的虫口数为xn ，第n+1代的虫口数为xn+1，则这种昆虫的繁殖规律可表为R.May, Nature 261, 459 (1976)MatLab编程基础clear alla=4;m=400;h=a/m;nu=0:h:a

32、;for k=1:m+1; n=1000;x=0.1; for i=1:n+1 x=nu(k)*x*(1-x); end y(1)=x; for j=1:100 y(j+1)=nu(k)*y(j)*(1-y(j); end hold on plot(nu(k),y,.,MarkerSize,5)end第四节 作图与动画-二维图形: 举例bug.mMatLab编程基础三维图形第四节 作图与动画-三维图形三维空间曲线三维空间曲线plot3(x,y,z,s)三维空间曲面三维空间曲面mesh(x,y,u)； surf(x,y,u)空间曲面的等高线空间曲面的等高线contour(x,y,u)MatLab

33、编程基础第四节 作图与动画-三维图形: 举例例例1t=0:pi/100:10*pi;x=t.*sin(t);y=cos(t);z=t;plot3(x,y,z)grid on-20-1001020-10105101520MatLab编程基础第四节 作图与动画-三维图形: 举例例例2x=-5:0.05:5; y=-5:0.05:5;X,Y=meshgrid(x,y);Z=2*exp(-(X-2).2-Y.2)+exp(-(X+2).2-Y.2)-3*exp(-X.2-(Y+1).2);mesh(X,Y,Z)surf(X,Y,Z)MatLab编程基础第四节 作图与动画-三维图形: 举例例例3x=-5

34、:0.05:5;y=-5:0.05:5;X,Y=meshgrid(x,y);Z=2*exp(-(X-2).2-Y.2) +exp(-(X+2).2-Y.2) -3*exp(-X.2-(Y+1).2);contour(X,Y,Z,50)MatLab编程基础四维数据的表现第四节 作图与动画-四维数据的表现四维数据的表现四维数据的表现slice(x,y,z,u,m,n,p)四维等面图四维等面图isosurface(x,y,z,u,isovalue)片组合图片组合图patch(FV,PropertyName,PropertyValue,)FV:isocaps(x,y,z,u,isovalue)isos

35、urface(x,y,z,u,isovalue)MatLab编程基础第四节 作图与动画-四维数据的表现例例x,y,z=meshgrid(-3:0.1:3,-3:0.1:3,-3:0.1:3);phi=sqrt(x.2+y.2).*exp(-sqrt(x.2+y.2+z.2);slice(x,y,z,phi,0,0,0)isosurface(x,y,z,phi,0.25)ex_10.mMatLab编程基础第四节 作图与动画-四维数据的表现x,y,z=meshgrid(-3:0.1:3,-3:0.1:3,-3:0.1:3);phi=sqrt(x.2+y.2).*exp(-sqrt(x.2+y.2+

36、z.2);x(:,:,31:61) = ;y(:,:,31:61) = ;z(:,:,31:61) = ;phi(:,:,31:61)=;patch(isosurface(x,y,z,phi,0.25),FaceColor,blue,EdgeColor,none);patch(isocaps(x,y,z,phi,0.25),FaceColor,interp,EdgeColor,none);view(3); axis tight; daspect(1,1,.4)lighting phonglight(Position,0 0 1,Style,infinite);ex_11.mMatLab编程基础

37、第四节 作图与动画-四维数据的表现MatLab编程基础第四节 作图与动画-四维数据的表现ex_12.mMatLab编程基础动画第四节 作图与动画-动画影片动画影片动画M=moviein(n)定义一个n列的矩阵M,以存储n幅画面，每列对应一幅M(:,j)=getframe 依次将n幅画面存入矩阵M，通常用for语句来完成Movie(M,k,FPS) 将n幅画面演示k次，每秒FPS幅，k默认1次，FPS默认12次MatLab编程基础第四节 作图与动画-动画例例1 ex_13.mclear allt=0:pi/20:4*pi;x=0:0.1*pi:4*pi;m=moviein(40)for i=1:

38、40 y=sin(x-t(i); plot(x,y); axis(0,12,-1,1); m(:,i)=getframe;endmovie(m,3,30)MatLab编程基础实时动画实时动画第四节 作图与动画-动画通过设置图形的属性来实现在保持背景图案的条件下擦除旧对象，显示新图像。当新对象的属性修改后，应该用指令drawnow刷新屏幕，使新对象显示出来set(句柄句柄,XData,x,YData,y)set(句柄句柄,EraseMode,属性属性)xor:擦除旧对象的点并绘制新对象的点none:不做任何擦除，直接在原图上绘制 MatLab编程基础第四节 作图与动画-动画例例1 ex_14.m

39、clear allt=0:pi/20:4*pi;x=0:0.1*pi:4*pi;y=sin(x);h=plot(x,y);axis(0,12,-1,1);set(h,EraseMode,xor)for i=2:60 y=sin(x-t(i); set(h,XData,x,YData,y); drawnow; pause(1)endMatLab编程基础第五节 常用的计算指令插值与曲线拟合微分积分求极值和零点快速Fourier变换解代数方程解常微分方程MatLab编程基础第五节 常用指令-插值与曲线拟合插值与曲线拟合插值当从一组不完全的或离散的数据中获取某些局部的信息时，通常使用函数的插值曲线拟合

40、如果要获取一个给定的离散集合的整体行为，需要采用曲线拟合，通常使用的是最小二乘法的多项式拟合MatLab编程基础第五节 常用指令-插值与曲线拟合插值一维插值yi=interp1(x,y,xi)yi=interp1(y,xi)yi=interp1(x,y,xi,method)nearest 线性最近项插值linear 线性插值spline 三次样条插值cubic 立方插值二维插值zi=interp2(x,y,z,xi,yi)zi=interp2(y,z,xi,yi)zi=interp2(z,ntimes)zi=interp2(x,y,z,xi,yi,method)MatLab编程基础第五节 常用

41、指令-插值与曲线拟合例例1 chazhi.mclear allx=0:10;y=sin(x);xi=0:0.25:10;y1=interp1(x,y,xi,nearest);y2=interp1(x,y,xi,linear);y3=interp1(x,y,xi,spline);y4=interp1(x,y,xi,cubic);subplot(2,2,1)plot(xi,y1,x,y,o)subplot(2,2,2)plot(xi,y2,x,y,o)subplot(2,2,3)plot(xi,y3,x,y,o)subplot(2,2,4)plot(xi,y4,x,y,o/p>

42、-1010510-1010510-101MatLab编程基础第五节 常用指令-插值与曲线拟合例例2 chazhi_1.mclear allt=0:0.2:4;x,y=meshgrid(t,t);z=sin(x.2+y.2);ti=0:0.02:4;xi,yi=meshgrid(ti,ti);zi=interp2(x,y,z,xi,yi,nearest);mesh(xi,yi,zi)figureze=sin(xi.2+yi.2);mesh(xi,yi,ze)nearest linearspline cubicExact resultMatLab编程基础曲线拟合第五节 常用指令-插值与曲线拟合P=

43、polyfit(x,y,N)其中x,y是样本点，N是拟合多项式的的最高幂次，返回的是拟合多项式的系数，按降幂排列。通常结合指令polyval来给出拟合多项式polyval(P,x)MatLab编程基础第五节 常用指令-插值与曲线拟合例例1 nh.mclear allx=0:0.5:5;yy=2*x+0.4*x.2;dy=0.5 -0.5 0.1 -0.4 0.8 . 0.8 -0.2 -0.6 0.4 -0.7 0.6;y=yy+dy;P=polyfit(x,y,3);xx=0:0.1:5;yp=polyval(P,xx);plot(x,y,o,xx,yp)0123450510152025Ma

44、tLab编程基础第五节 常用指令-微分微分差分梯度MatLab编程基础第五节 常用指令-微分例例 wf_1.m2sin()yxclear allh=0.1;x=0:h:pi;y1=diff(sin(x.2)/h;y2=2*cos(x.2).*x;plot(x,y2,x(1:end-1),y1,-.)legend(exact,numerical)差分1nnxxx diff(x)MatLab编程基础第五节 常用指令-微分梯度,dfdfdx dyfx,fy=gradient(f,hx,hy)fx,fy,fz=gradient(f,hx,hy,hz),dfdfdfdx dy dz通常伴随着quiver

45、指令：quiver(x,y,fx,fy)MatLab编程基础第五节 常用指令-微分例例 wf_2.m22xyzxeclear allh=0.2;x=-2:h:2;X,Y=meshgrid(x);Z=X.*exp(-X.2-Y.2);Px,Py=gradient(Z,h);contour(X,Y,Z,. -.4,-.3,-.2,-.1,0,.1,.2,.3,.4,b);hold onquiver(X,Y,Px,Py,k)MatLab编程基础第五节 常用指令-积分积分梯形积分累计梯形积分函数积分二重积分MatLab编程基础第五节 常用指令-积分梯形积分trapz(x,y)累计梯形积分cumtrap

46、z(x,y)注：前者是定积分； 后者是原函数另：sumMatLab编程基础clear allx=0:0.01:2;y=exp(x);z=cumtrapz(x,y);% exact expressionze=exp(x)-1;plot(x,z,x,ze,.)legend(numerical,exact)第五节 常用指令-积分例例 jf_1.m和和jf_2.m2201xe dxe01xtxe dteclear allx=0:0.01:2;y=exp(x);z=trapz(x,y)z = 6.3891MatLab编程基础函数积分第五节 常用指令-积分quad(fun,a,b)quad(fun,a,b

47、,tol,trace,p1,p2,)例例 jf_3.mclear allqx=inline(x.3-a*x.2+x+10,x,a);A=quad(qx,-2,7,8)A = -227.2500与函数与函数M文件的区别，见文件的区别，见jf_31.mMatLab编程基础第五节 常用指令-积分二重积分dblquad(fun,inmin,inmax,outmin,outmax)例例 jf_4.mclear alltol=10(-5);f2=inline(y*sin(x)+x*cos(y),x,y);A=dblquad(f2,pi,2*pi,0,pi,tol)A = -9.8696dblquad(fu

48、n,inmin,inmax,outmin,outmax,tol)MatLab编程基础第五节 常用指令-求极值和零点求极值单变量函数的极值多变量函数的极值求零点单变量函数的零点对多项式求零点MatLab编程基础第五节 常用指令-求极值和零点求极值多变量函数的极值fminsearch(fun,x0)fminsearch(fun,x0,options)x,fval=fminsearch()单变量函数的极值fminbnd(fun,x1,x2)fminbnd(fun,x1,x2,options)x,fval=fminbnd()min, maxMatLab编程基础第五节 常用指令-求极值和零点例例 jz_

49、1.m32810 2 8xxx求在，内的极小值clear allqx=inline(x.3-8*x.2+x+10);fplot(qx,-2,8)x,f=fminbnd(qx,-2,8)x = 5.2701f = -60.5500MatLab编程基础第五节 常用指令-求极值和零点例例 jz_2.m22(21)xyxe求在，邻近的极小值clear allqx=inline(x(1).*exp(-x(1).2-x(2).2);ezmesh(x.*exp(-x.2-y.2),-2,2,-2,2)x,f=fminsearch(qx,-2,-1)x = -0.7071 -0.0000f = -0.4289

50、MatLab编程基础第五节 常用指令-求极值和零点求零点多项式的零点roots(c)单变量函数的零点fzero(fun,x0)fzero(fun,x0,options)x,fval=fzero()MatLab编程基础第五节 常用指令-求极值和零点clear allqx=inline(x.3-8*x.2+x+10);fplot(qx,-2,8)x,fval=fzero(qx,2)%-c=1 -8 1 10;roots(c)例例 ld_1.m32810 xxx求的零点x =1.2984 fval = -1.7764e-015ans =7.7016 1.2984 -1.0000MatLab编程基础第

51、五节 常用指令-快速Fourier变换快速Fourier变换离散Fourier变换10 ,110NnegNfNkitknkn10 ,10NkefgNntinknk设f(t)是定义在时域0,T)上的复值函数TkkhNnTnhtktn2 ,ThNTht2tNhh2频域频域2N/T-2N/T0N/T-N/TTT/20-T/2-T时域MatLab编程基础第五节 常用指令-快速Fourier变换MatLab中的形式12(1)11, 1kNinNnkkxX enNN12(1)1, 1kNinNknnXx ekNN-110-(N-1)2-1(N-1)/N-(N-1)/N0 1/N 2/N-1/N12211,

52、NNNNNN 功率谱2kXN12210,NNNNNN1频域频域负频域负频域用于不同频率处功率的测量MatLab编程基础第五节 常用指令-快速Fourier变换clear alltict=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t);xr=x+2*randn(size(t);plot(1000*t(1:100),xr(1:100)N=512;for k=1:N X(k)=0; for n=1:N X(k)=X(k)+xr(n)*exp(-i*2*pi*(k-1)*(n-1)/N); endendfreq=1000*(0:N/2)/N;PX=X.*con

53、j(X)/N;figureplot(freq,PX(1:N/2+1)toc例例 fly_1.m 直接计算方法直接计算方法计算频率计算频率和功率谱和功率谱MatLab编程基础第五节 常用指令-快速Fourier变换elapsed_time=4.2820MatLab编程基础第五节 常用指令-快速Fourier变换MatLabMatLab的的FFTFFT指令指令一维FFT及其逆变换fft(x), fft(x,n)ifft(x), ifft(x,n)二维FFT及其逆变换fft2(x), fft2(x,m,n)ifft2(x), ifft2(x,m,n)通常伴随着fftshift指令MatLab编程基础

54、第五节 常用指令-快速Fourier变换clear alltict=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t);xr=x+2*randn(size(t);plot(1000*t(1:100),xr(1:100) )N=512; X=fft(xr,N);freq=1000*(0:N/2)/N; PX=X.*conj(X)/N;figureplot(freq,PX(1:N/2+1)toc计算频率计算频率计算功率谱计算功率谱例例 fly_2.m 利用利用FFTelapsed_time=0.2190MatLab编程基础第五节 常用指令-快速Fourier变

55、换clear allf=zeros(30,30);f(5:24,13:17)=1;figure(1)imshow(f,notruesize)colormap(jet)F=fft2(f,256,256);FF=fftshift(F);figure(2)imshow(abs(FF),0,5,notruesize)colormap(jet)colorbarF2=log(abs(FF);figure(3)imshow(F2,-1,5,notruesize)colormap(jet)colorbar例例 fly_3.m，二维二维FFTMatLab编程基础第五节 常用指令-快速Fourier变换Image

56、 of log(abs(FF)Image of abs(FF)MatLab编程基础第五节 常用指令-快速Fourier变换例例 fly_4.m和和fly_5.m，与连续与连续FourierFourier变换的比较变换的比较连续连续FourierFourier变换变换( )( )i xFf x edx1( )( )2i xf xFed32,0( )0, 0 xexf xx2( )3Fi1, 1,1( )0, otherxf x 2sin( )FMatLab编程基础第五节 常用指令-快速Fourier变换例例 fly_4.mclear allN=2048;hx=0.04;x=(-N/2:N/2-1

57、)*hx;hw=2*pi/(N*hx);w=fftshift(-N/2:N/2-1)*hw);for j=1:N if x(j)=0 f(j)=2*exp(-3*x(j); else f(j)=0; endendF=fft(f);Fc=fft(f)*hx;ww=fftshift(w);Fa=2./(3+i*ww);plot(ww,abs(Fa),w,abs(Fc),r.)legend(exact F(omega),DFT FFT)ff=ifft(F);figureplot(x,f,.,x,real(ff),o)legend(exact f(x),DFT FFT and IFFT)MatLab编

58、程基础第五节 常用指令-快速Fourier变换频域比较频域比较时域比较时域比较高频部分高频部分MatLab编程基础第五节 常用指令-快速Fourier变换例例 fly_5.mclear allN=1024;hx=0.04;x=(-N/2:N/2-1)*hx;hw=2*pi/(N*hx);w=fftshift(-N/2:N/2-1)*hw);for j=1:N if x(j)=-1 & x(j)=1 f(j)=1; else f(j)=0; endendF=fft(f);Fc=fft(f)*hx;ww=fftshift(w);Fa=2*sin(ww+eps)./(ww+eps);plot(ww,abs(Fa),w,abs(Fc),r.)legend(exact F(omega),DFT FFT)ff=ifft(F);figureplot(x,f,.,x,real(ff),o)legend(exact f(x),DFT FFT and IFFT)MatLab编程基础第五节 常用指令-快速Fourier变换频域比较频域比较时域比较时域比较高频部分高频部分MatLab编程基础第五节 常用指令-解代数方程