MATLAB在地学中的应用

1、The first part:nMatlab在地在地学中的应用学中的应用一一MATLABMATLAB入门入门nMATLABMATLAB是是MathWorksMathWorks公司的产品，公司的产品， MATLABMATLAB名字由名字由MATrixMATrix和和 LABoratory LABoratory 两词的前三个字母组合而成。两词的前三个字母组合而成。 MATLABMATLAB是一种交互式、面向对象的程序设计语是一种交互式、面向对象的程序设计语言，广泛应用于工业界与学术界，主要言，广泛应用于工业界与学术界，主要用于矩阵运算，同时在数值分析、自动用于矩阵运算，同时在数值分析、自动控制模拟

2、、数字信号处理、动态分析和控制模拟、数字信号处理、动态分析和绘图等方面也具有强大的功能。绘图等方面也具有强大的功能。MATLABMATLAB能做什么能做什么nMATLABMATLAB程序设计语言结构完整，且具程序设计语言结构完整，且具有优良的移植性。它可以高效率地解决有优良的移植性。它可以高效率地解决工业计算问题，特别是关于矩阵和向量工业计算问题，特别是关于矩阵和向量的计算。的计算。MATLABMATLAB与与C C语言和语言和FORTRANFORTRAN语言相比更容易被掌握，寥语言相比更容易被掌握，寥寥几行寥几行MATLABMATLAB程序就可以实现意想不程序就可以实现意想不到的功能。到的功

3、能。MATLABMATLAB能做什么能做什么n另外，另外，MATLABMATLAB提供了一种特殊的工具：提供了一种特殊的工具：工具箱工具箱(Toolbox)(Toolbox)。这些工具箱主要包括：。这些工具箱主要包括：信号处理信号处理(Signal Processing)(Signal Processing)、控制系、控制系统统(Control System)(Control System)、神经网络、神经网络(Neural (Neural Network)Network)、虚拟现实、虚拟现实(Virtual Reality)(Virtual Reality)、小波小波(Wavelets)(W

4、avelets)和和地图制作地图制作(Mapping)(Mapping)等。等。MathWorksMathWorks公司和各个领域的专家还在公司和各个领域的专家还在积极地完善现有工具箱的功能，并不断积极地完善现有工具箱的功能，并不断推出新的工具箱。推出新的工具箱。MATLABMATLAB能做什么能做什么nMATLABMATLAB是一个可视化的计算程序，被是一个可视化的计算程序，被广泛地使用于从个人计算机到超级计算广泛地使用于从个人计算机到超级计算机范围内的各种计算机上。机范围内的各种计算机上。nMATLABMATLAB包括命令控制、可编程，有上包括命令控制、可编程，有上百个预先定义好的命令和函

5、数。百个预先定义好的命令和函数。nMATLABMATLAB有强有力的二维、三维图形工有强有力的二维、三维图形工具。具。MATLABMATLAB能与其他程序一起使用。能与其他程序一起使用。n不同的不同的MATLABMATLAB工具箱工具箱(Toolbox)(Toolbox)可应用可应用于特殊的应用领域。于特殊的应用领域。 初识初识MATLABMATLABnMATLABMATLAB启动后的主界面中包含一个命令窗口启动后的主界面中包含一个命令窗口（Command WindowCommand Window），在其中可以运行各），在其中可以运行各种种MATLABMATLAB命令。命令。n例例1 1：求：

6、求n（1 1）用键盘在）用键盘在MATLABMATLAB指令窗中输入以下内指令窗中输入以下内容容n (12+2 (12+2* *(7-4)/32 (7-4)/32 n（2 2）在上述表达式输入完成后，按【）在上述表达式输入完成后，按【EnterEnter】键，该指令就被执行。键，该指令就被执行。23)47(212初识初识MATLABMATLABn（3 3）在指令执行后，）在指令执行后，MATLABMATLAB指令窗指令窗中将显示以下结果。中将显示以下结果。nans =ans =n 2 2 n（注意：计算结果中（注意：计算结果中“ans”ans”是英文是英文“answer”answer”的缩的缩

7、写，其含义为写，其含义为“运算答案运算答案”，是，是MATLABMATLAB的一个默认的一个默认变量。）变量。） 初识初识MATLABMATLABn例例2 2：简单矩阵的输入步骤。：简单矩阵的输入步骤。n（1 1）在键盘上输入下列内容）在键盘上输入下列内容nA = 1,2,3; 4,5,6; 7,8,9A = 1,2,3; 4,5,6; 7,8,9 n（2 2）按【）按【EnterEnter】键，指令被执行。】键，指令被执行。n（3 3）在指令执行后，）在指令执行后，MATLABMATLAB指令窗中将显示以下指令窗中将显示以下结果：结果：nA =A =n 1 2 3 1 2 3n 4 5 6

8、4 5 6n 7 8 9 7 8 9 初识初识MATLABMATLABn 矩阵的分行输入矩阵的分行输入: :nA=1 2 3A=1 2 3n 4 5 6 4 5 6n 7 8 9 7 8 9 nA =A =n 1 2 3 1 2 3n 4 5 6 4 5 6n 7 8 9 7 8 9 初识初识MATLABMATLABn例例3 3 二维函数（二维函数（MATLABMATLAB能用于计算、并以二维和三能用于计算、并以二维和三维图形显示各种函数。在维图形显示各种函数。在MATLABMATLAB函数中包括了所有函数中包括了所有主要的数学函数和大量的高级函数。）主要的数学函数和大量的高级函数。） n用简

9、短的用简短的MATLABMATLAB命令计算并绘制在命令计算并绘制在0X60X6范围内范围内的的sin(2x)sin(2x)、sinxsinx2 2和和sin(x)sin(x)2 2。nx = linspace( 0, 6 ) ;x = linspace( 0, 6 ) ; % % 创建一个向量创建一个向量x xny1 = sin (2y1 = sin (2* *x);x);ny2 = sin (x.2) ;y2 = sin (x.2) ; ny3 = sin (x).2 ; y3 = sin (x).2 ; 初识初识MATLABMATLABnplot (x,y1); plot (x,y1);

10、 n%绘制向量绘制向量y1y1，y1y1作为向量作为向量x x的一个函数的一个函数nhold on; hold on; n%保存当前图形保存当前图形nplot (x,y2,r);plot (x,y2,r);nplot (x,y3,+);plot (x,y3,+);nhold offhold off进一步了解进一步了解MATLABMATLAB n一些常用命令和功能键：一些常用命令和功能键：nexitexit，quitquit， 结束结束MATLABMATLAB会话会话（如果没有保存，则变量中的数据丢失）（如果没有保存，则变量中的数据丢失）n 中断一个中断一个MATLABMATLAB任务（但会话任

11、务（但会话并没有结束）并没有结束）ncla cla 清除坐标清除坐标nclf clf 清除图形清除图形nclc clc 清除清除Command WindowCommand Window的所有显示的所有显示内容内容nclear clear 清除所有变量清除所有变量进一步了解进一步了解MATLABMATLABndir dir 列出一个目录中的所有文件列出一个目录中的所有文件nls ls 以不同的输出格式列出文件以不同的输出格式列出文件npwd pwd 列出当前的工作目录列出当前的工作目录ndelete filename delete filename 删除文件删除文件ncd cd 改变当前目录改变

12、当前目录ntype filename type filename 显示文件的内容显示文件的内容nedit filename edit filename 编辑文件编辑文件nhelp help 获取帮助信息获取帮助信息课堂练习课堂练习n在在Command WindowCommand Window中练习使用中练习使用plotplot函数函数: :首先输入一个向量首先输入一个向量, ,然后在然后在同一个坐标系内分别绘制此向量的同一个坐标系内分别绘制此向量的正切与余切函数曲线正切与余切函数曲线, ,并分别表示成并分别表示成不同的颜色或线型。不同的颜色或线型。二二MATLABMATLAB程序设计程序设计n

13、通常，通常，MATLABMATLAB运行在单命令驱动模式下，即运行在单命令驱动模式下，即用户可以通过键盘键入一条一条的用户可以通过键盘键入一条一条的MATLABMATLAB命命令，每键入一条命令，令，每键入一条命令，MATLABMATLAB就立即执行该就立即执行该命令并显示出执行结果。但有时单靠一条一条命令并显示出执行结果。但有时单靠一条一条地输入命令难以实现复杂功能，而地输入命令难以实现复杂功能，而MATLABMATLAB还还可以在文件驱动模式下进行工作，驱动文件称可以在文件驱动模式下进行工作，驱动文件称为为MM文件文件，它由一系列，它由一系列MATLABMATLAB语句构成。也语句构成。也

14、就是说，就是说，MATLABMATLAB可以执行包含一系列命令的可以执行包含一系列命令的文件。这两种模式共同构成了文件。这两种模式共同构成了MATLABMATLAB的交互的交互式工作环境。式工作环境。 MATLABMATLAB程序设计程序设计n由由MATLABMATLAB语言编写的磁盘文件称为语言编写的磁盘文件称为MM文件，文件，其扩展名为其扩展名为.m.m。例如，。例如，belsel.mbelsel.m就是一个由一就是一个由一系列系列MATLABMATLAB语句编写的语句编写的MM文件，它可用来产文件，它可用来产生贝塞尔曲线。生贝塞尔曲线。MM文件中还可以调用其它文件中还可以调用其它MM文文

15、件，也就是说件，也就是说MM文件可以互相调用，甚至可以文件可以互相调用，甚至可以调用自己。调用自己。n（MM文件存盘时不能与文件存盘时不能与MATLABMATLAB内部定义的函内部定义的函数重名，也不能保存成纯数字文件名！）数重名，也不能保存成纯数字文件名！） M M文件文件n为建立新的为建立新的MM文件，启动文件，启动MATLABMATLAB文本编辑器文本编辑器有有3 3种方法种方法: :n（1 1） 菜单操作：从菜单操作：从MATLABMATLAB操作桌面的操作桌面的“File”File”菜单中选择菜单中选择“New”New”菜单项，再选择菜单项，再选择“M-File”M-File”命令，

16、屏幕将出现命令，屏幕将出现MATLABMATLAB文本编文本编辑的窗口辑的窗口; ;n（2 2）命令按钮操作：单击）命令按钮操作：单击MATLABMATLAB命令窗口命令窗口工具栏上的新建命令按钮工具栏上的新建命令按钮; ;n（3 3）命令操作：在）命令操作：在MATLABMATLAB命令窗口输入命命令窗口输入命令令“edit”,edit”,按按键后，即可启动键后，即可启动MATLABMATLAB文本编辑器。文本编辑器。M M文件文件nMATLABMATLAB文本编辑器是一个集编辑与调文本编辑器是一个集编辑与调试两种功能于一体的工具环境。利用它试两种功能于一体的工具环境。利用它不仅可以完成基本

17、的文本编辑操作，还不仅可以完成基本的文本编辑操作，还可以对可以对MM文件进行调试。文件进行调试。n打开已有的打开已有的MM文件（方法类似）。文件（方法类似）。M M文件分类文件分类nMM文件从功能上可分为两类：文件从功能上可分为两类：脚本文件和函数脚本文件和函数文件文件。调用脚本文件会自动执行一系列命令直。调用脚本文件会自动执行一系列命令直至输出结果；函数文件则为至输出结果；函数文件则为MATLABMATLAB提供了扩提供了扩充性，通过函数文件，可以产生完成某一特定充性，通过函数文件，可以产生完成某一特定功能的函数。事实上，功能的函数。事实上，MATLABMATLAB提供的标准函提供的标准函数

18、大部分都是由函数文件定义的，这足以说明数大部分都是由函数文件定义的，这足以说明函数文件的重要性。函数文件的目的是扩展函数文件的重要性。函数文件的目的是扩展MATLABMATLAB功能，也就是说，函数文件利用功能，也就是说，函数文件利用MATLABMATLAB语言构造了一个新的语言构造了一个新的MATLABMATLAB函数，函数，而且这个函数的使用方法与而且这个函数的使用方法与MATLABMATLAB本身提供本身提供的库函数一样。的库函数一样。 函数文件的基本结构函数文件的基本结构n函数文件的第一行必须 包 含 关 键 字函数文件的第一行必须 包 含 关 键 字“function”functio

19、n”，其基本结构为：，其基本结构为：n function function 输出形参表输出形参表=函数名（输入函数名（输入形参表）形参表）n %注释语句注释语句n 函数体语句函数体语句脚本文件举例脚本文件举例n例：使用脚本式例：使用脚本式MM文件绘制一个表面比文件绘制一个表面比较平滑的单位球面。较平滑的单位球面。n n 新建一个新建一个MM文件，在文本编辑器中输文件，在文本编辑器中输入以下代码，然后保存成入以下代码，然后保存成spher.mspher.m文件，文件，在命令窗口中输入在命令窗口中输入spherspher即可生成单位球即可生成单位球面。面。脚本文件举例脚本文件举例nclear al

20、lclear alln%清除所有变量及函数等内容清除所有变量及函数等内容nspherespheren%绘制一个单位球面绘制一个单位球面nshading interpshading interpn%对球面进行插值着色对球面进行插值着色naxis equalaxis equaln%使各坐标系方向上的度量单位相同，如果不使各坐标系方向上的度量单位相同，如果不同则象椭球面同则象椭球面函数文件举例函数文件举例n例：使用函数式例：使用函数式M文件编制一个绘图函数，调文件编制一个绘图函数，调用此函数时可以根据输入参数的不同生成单位用此函数时可以根据输入参数的不同生成单位球面或单位柱面。球面或单位柱面。n 在

21、在M文件编辑器中输入以下代码，然后保文件编辑器中输入以下代码，然后保存到存到MATLAB安装目录下的安装目录下的work目录，文件目录，文件名为名为drawsur.m。在命令窗口中输入。在命令窗口中输入drawsur(sphere)生成单位球面；输入生成单位球面；输入drawsur(cylinder)生成单位柱面。生成单位柱面。函数文件举例函数文件举例nfunction drawsur(surface)n switch surfacen case spheren spheren case cylindern cylindern endn shading interpn axis equal编程

22、语句结构编程语句结构n 从理论上讲，只要有顺序、循环和分支从理论上讲，只要有顺序、循环和分支3 3种基本程序结构，种基本程序结构，MATLABMATLAB就可以构就可以构成任何一种程序并完成相应的工作。与成任何一种程序并完成相应的工作。与大多数计算机语言一样，大多数计算机语言一样，MATLABMATLAB有设有设计程序所必须的程序结构，即顺序结构、计程序所必须的程序结构，即顺序结构、循环结构与分支结构。在循环结构与分支结构。在MATLABMATLAB语言语言中，循环由中，循环由forfor和和whilewhile语句实现，分支语句实现，分支主要由主要由if if语句实现。语句实现。forfor

23、循环语句循环语句nMATLABMATLAB提供的提供的forfor循环语句，其功能是重复循环语句，其功能是重复执行一条或一组语句的固定次数执行一条或一组语句的固定次数( (也可以是零也可以是零次次) )，例如：，例如：nclear allclear allnn=5n=5nfor ifor i1 1：n n n x(i) x(i)0 0nendendn语句的执行结果是把语句的执行结果是把0 0分配给分配给x x的前的前n n个元素。个元素。forfor循环语句循环语句n其中值得注意的是：如果其中值得注意的是：如果n n小于小于1 1，语句，语句的描述仍是合法的，只是的描述仍是合法的，只是MATL

24、ABMATLAB不执不执行其中的语句行其中的语句x(i)x(i)0 0，或者说执行零次。，或者说执行零次。如果变量如果变量x x不存在，不存在，MATLABMATLAB将自动生成将自动生成该变量并为其分配存储空间；如果原先该变量并为其分配存储空间；如果原先x x变量的元素个数少于变量的元素个数少于n n个，则个，则MATLABMATLAB也会自动为其分配所需增加的空间。也会自动为其分配所需增加的空间。forfor循环语句循环语句nforfor循环的循环的般格式为：般格式为：nfor v=expressionfor v=expressionn statement statementnendend

25、n其中其中v v为变量名；为变量名；expression(expression(表达式表达式) )实际上实际上是一个矩阵，因为在是一个矩阵，因为在MATLABMATLAB中基本元素就是中基本元素就是矩阵。当矩阵的矩阵。当矩阵的每一列元素每一列元素依次赋给变量依次赋给变量v v时，时，statement(statement(语句语句) )就执行一次。就执行一次。n注意：注意：forfor语句必须有语句必须有endend与之匹配。与之匹配。 whilewhile循环语句循环语句nMATLABMATLAB还提供了还提供了whilewhile循环语句，其循环语句，其功能是重复执行一条或者一组语句，重功

26、能是重复执行一条或者一组语句，重复次数由某个逻辑条件控制，复次数由某个逻辑条件控制，般来说般来说其执行次数是不能预先确定的。下面通其执行次数是不能预先确定的。下面通过一个简单例子来说明过一个简单例子来说明whilewhile语句的功能。语句的功能。n例：求例：求1 1到到100100的所有偶数项之和。的所有偶数项之和。whilewhile循环语句循环语句nclear allclear allnx=0;x=0;nsum=0;sum=0;nwhile x101while x101n sum=sum+x; sum=sum+x;n x=x+2; x=x+2;nendendnsumsumwhilewhi

27、le循环语句循环语句nwhilewhile循环语句的一般格式是：循环语句的一般格式是：nwhile expressionwhile expressionn statement statementnendendn与与forfor语句语句样，样，whilewhile语句也必须有语句也必须有endend与之相匹配。与之相匹配。 ifif分支结构分支结构nMATLABMATLAB中中if if条件语句的功能是：如果满足某条件语句的功能是：如果满足某个条件则执行一条或一组语句，否则就不执行个条件则执行一条或一组语句，否则就不执行其内部语句，或者执行其内部语句，或者执行elseelse后面的语句。后面的语

28、句。nif if条件语句的一般格式有两种条件语句的一般格式有两种: :n格式一：格式一：nif conditionif conditionn statement statementnendendifif分支结构分支结构n格式二：格式二：nif conditionif conditionn statement1 statement1nelseelsen statement2 statement2nendendn这两种格式类似于这两种格式类似于C C语言中的语言中的if if语句，而且语句，而且if if条条件语句可以嵌套使用。件语句可以嵌套使用。但仍应注意，但仍应注意，if if语句也语句也必须

29、和必须和endend语句成对出现。语句成对出现。ifif分支结构分支结构n例：输入一个数，如果小于例：输入一个数，如果小于1010就输出这个数加就输出这个数加1010；否则输出这个数减；否则输出这个数减1010。nclear allclear allnn=input(n=input(Please enter a number,n=Please enter a number,n=); );nif n10if n load coast load coastn whos whosn Name Size Bytes ClassName Size Bytes Classn lat 9589x1 7671

30、2 double array lat 9589x1 76712 double arrayn long 9589x1 76712 double array long 9589x1 76712 double array1 1 矢量数据矢量数据n变量变量latlat和和longlong都是都是coastcoast文件中的矢量，文件中的矢量，它们组成一幅世界海岸线矢量地图。它们组成一幅世界海岸线矢量地图。n（2 2）要查看地图，输入：）要查看地图，输入：n axesm mercator % mercator axesm mercator % mercator投影投影n framem framemn p

31、lotm(lat,long) plotm(lat,long)1 1 矢量数据矢量数据n（3 3）查看矢量数据的前）查看矢量数据的前2020个坐标：个坐标：nlat(1:20) long(1:20)lat(1:20) long(1:20)n（4 4）为了查看这些矢量点表示的海岸线，）为了查看这些矢量点表示的海岸线，输入：输入：nplotm(lat(1:20),long(1:20),r)plotm(lat(1:20),long(1:20),r)2 2 栅格数据栅格数据n还可以将矩阵数据用地图表示，矩阵中还可以将矩阵数据用地图表示，矩阵中每个行、列交叉点的元素对应特定地理每个行、列交叉点的元素对应特

32、定地理区域的矩形面片，另有数据表示它与相区域的矩形面片，另有数据表示它与相邻面片的连通性，这就是栅格数据，栅邻面片的连通性，这就是栅格数据，栅格的值可以用单元内的平均值或中心值格的值可以用单元内的平均值或中心值表示。表示。2 2 栅格数据栅格数据n例：输入如下命令，载入例：输入如下命令，载入topotopo数据集并数据集并显示它的结构。显示它的结构。n load topo load topon whos whosn Name Size Bytes ClassName Size Bytes Classn topo 180 x360 518400 double array topo 180 x36

33、0 518400 double arrayn topolegend 1x3 24 double array topolegend 1x3 24 double arrayn topomap1 64x3 1536 double array topomap1 64x3 1536 double arrayn topomap2 128x3 3072 double array topomap2 128x3 3072 double array3 3 操作矢量数据操作矢量数据n多边形布尔操作可以回答一系列与矢量多边形布尔操作可以回答一系列与矢量数据多边形对象逻辑关系有关的问题。数据多边形对象逻辑关系有关的问题

34、。标准布尔操作包括交、并、差和异或等，标准布尔操作包括交、并、差和异或等，polyboolpolybool函数可以对两个矢量集合进行函数可以对两个矢量集合进行这些操作。这些操作。n例如：例如：3 3 操作矢量数据操作矢量数据n（1 1）创建一个十二边形：）创建一个十二边形：n theta=(0:pi/6:2 theta=(0:pi/6:2* *pi);pi);n lat1=sin(theta); lat1=sin(theta);n lon1=cos(theta); lon1=cos(theta);n（2 2）创建一个三角形叠加到十二边形上：）创建一个三角形叠加到十二边形上：n lat2=0 1

35、 -1 0; lat2=0 1 -1 0;n lon2=0 2 2 0; lon2=0 2 2 0;3 3 操作矢量数据操作矢量数据n（3 3）蓝线和红线将两个形状显示在一起：）蓝线和红线将两个形状显示在一起：n axesm miller %miller axesm miller %miller投影投影n plotm(lat1,lon1,b) plotm(lat1,lon1,b)n plotm(lat2,lon2,r) plotm(lat2,lon2,r)3 3 操作矢量数据操作矢量数据n（4 4）计算多边形的交，绘制成绿色面片：）计算多边形的交，绘制成绿色面片：n lati,loni=pol

36、ybool(intersection,lat1,lati,loni=polybool(intersection,lat1,lon1,lat2,lon2)lon1,lat2,lon2)n patchm(lati,loni,g) patchm(lati,loni,g)3 3 操作矢量数据操作矢量数据n（5 5）计算多边形的并，绘成洋红色面片：）计算多边形的并，绘成洋红色面片：n latu,lonu=polybool(union,lat1,lonlatu,lonu=polybool(union,lat1,lon1,lat2,lon2)1,lat2,lon2)n patchm(latu,lonu,m)

37、 patchm(latu,lonu,m)3 3 操作矢量数据操作矢量数据n（6 6）对多边形进行异或计算，显示成黄）对多边形进行异或计算，显示成黄色面片：色面片：n latx,lonx=polybool(xor,lat1,lon1,lalatx,lonx=polybool(xor,lat1,lon1,lat2,lon2)t2,lon2)n patchm(latx,lonx,y) patchm(latx,lonx,y)3 3 操作矢量数据操作矢量数据n（7 7）求十二边形与三角形的差，绘制成）求十二边形与三角形的差，绘制成白色面片：白色面片：n latm,lonm=polybool(minus,

38、lat1,lolatm,lonm=polybool(minus,lat1,lon1,lat2,lon2)n1,lat2,lon2)n patchm(latm,lonm,w) patchm(latm,lonm,w)3 3 操作矢量数据操作矢量数据n使用使用buffermbufferm函数计算并返回缓冲区点函数计算并返回缓冲区点集矢量：集矢量：n（1 1）创建一个包围意大利的区域基准图，隐）创建一个包围意大利的区域基准图，隐藏其边界藏其边界nclose allclose allnclear allclear allnworldmap(lo,35,50,3,23,lineonly)worldmap(

39、lo,35,50,3,23,lineonly)nhidem(gca)hidem(gca)3 3 操作矢量数据操作矢量数据n（2 2）读入）读入worldloworldlo数据并提取出表示意大利的数据并提取出表示意大利的多边形多边形nload worldloload worldlonilat,ilon=extractm(POpatch,italy);ilat,ilon=extractm(POpatch,italy);n（3 3）用）用buffermbufferm函数处理意大利的多边形，函数处理意大利的多边形，输出一个距离国界输出一个距离国界1.51.5度的缓冲区度的缓冲区nlatb,lonb=b

40、ufferm(ilat,ilon,1.5,out);latb,lonb=bufferm(ilat,ilon,1.5,out);n（4 4）缓冲区绘成黄色，意大利界内绘成绿色）缓冲区绘成黄色，意大利界内绘成绿色npatchesm(latb,lonb,y)patchesm(latb,lonb,y)npatchesm(ilat,ilon,g)patchesm(ilat,ilon,g)4 4 操作栅格数据操作栅格数据n可以将纬度可以将纬度经度矢量数据转换为任意经度矢量数据转换为任意精度的网格数据，然后可以使用该网格精度的网格数据，然后可以使用该网格数据制作栅格图。地图制作工具箱提供数据制作栅格图。地图

41、制作工具箱提供了此类转换的图形用户界面，也可以从了此类转换的图形用户界面，也可以从命令行中进行转换。矢量数据网格化的命令行中进行转换。矢量数据网格化的主要函数是主要函数是vec2mtxvec2mtx。n例如：例如：4 4 操作栅格数据操作栅格数据n（1 1）用）用extractmextractm函数获取瑞士国界的函数获取瑞士国界的面片数据：面片数据：nswlat,swlon=extractm(worldlo(POpswlat,swlon=extractm(worldlo(POpatch),Switzerland);atch),Switzerland);n（2 2）将网格密度设置为）将网格密度设

42、置为4040单元单元/ /度，用度，用vec2mtxvec2mtx函数光栅化国界，并生成它的函数光栅化国界，并生成它的引用向量：引用向量：nden=40;den=40;nswgrid,swrv=vec2mtx(swlat,swlon,swgrid,swrv=vec2mtx(swlat,swlon,den);den);nwhoswhos4 4 操作栅格数据操作栅格数据n（3 3）用对比色绘制数据网格的地图：）用对比色绘制数据网格的地图：naxesm eqdcylin %axesm eqdcylin %等距离柱面投影等距离柱面投影nmeshm(swgrid,swrv)meshm(swgrid,sw

43、rv)ncolormap jet(4)colormap jet(4)n（4 4）设置地图的范围限制：）设置地图的范围限制：nlatlim lonlim=limitm(swgrid,swrv);latlim lonlim=limitm(swgrid,swrv);nsetm(gca,FlatLimit,latlim,FlonLisetm(gca,FlatLimit,latlim,FlonLimit,lonlim)mit,lonlim)ntightmaptightmap4 4 操作栅格数据操作栅格数据n（5 5）为了填充瑞士的内部区域，选择网）为了填充瑞士的内部区域，选择网格中间的行和列，并在调用格

44、中间的行和列，并在调用encodemencodem函函数生成新网格时选择索引值数生成新网格时选择索引值3 3来识别瑞士来识别瑞士的版图：的版图：nswpt=size(swgrid)/2,3swpt=size(swgrid)/2,3nswgrid3=encodem(swgrid,swpt,1)swgrid3=encodem(swgrid,swpt,1)n%最后一个变量识别边界单元的编码。最后一个变量识别边界单元的编码。4 4 操作栅格数据操作栅格数据n（6 6）用填充后的网格清除和重绘地图：）用填充后的网格清除和重绘地图：nmeshm(swgrid3,swrv)meshm(swgrid3,swr

45、v)n（7 7）在网格上绘原向量，比较栅格化的）在网格上绘原向量，比较栅格化的效果：效果：nplotm(swlat,swlon,k)plotm(swlat,swlon,k)显示一个旋转的地球显示一个旋转的地球n因为从各个方向观察球体时不需要重新因为从各个方向观察球体时不需要重新计算投影，所以比较容易地用动画来表计算投影，所以比较容易地用动画来表现一个旋转的地球。如果显示的数据足现一个旋转的地球。如果显示的数据足够简单，这个动画就可以相对较快的速够简单，这个动画就可以相对较快的速度进行绘制。在下面的例子中，首先用度进行绘制。在下面的例子中，首先用一个一个MM文件设置地球从西向东按照一定文件设置地

46、球从西向东按照一定的增量旋转如何重新设置视图，然后用的增量旋转如何重新设置视图，然后用该该MM文件控制地球旋转时的显示效果。文件控制地球旋转时的显示效果。显示一个旋转的地球显示一个旋转的地球n（1 1）创建一个包含下面代码的）创建一个包含下面代码的MM文件：文件：nfor i=360:-5:0for i=360:-5:0n view(i,0); view(i,0);n drawnow drawnownendendn（2 2）保存名称为）保存名称为spin.mspin.m。显示一个旋转的地球显示一个旋转的地球n（3 3）用格线创建球体显示效果：）用格线创建球体显示效果：naxesm(globe,

47、Grid,on,Gcolor,.7,.axesm(globe,Grid,on,Gcolor,.7,.8,.9,GlineStyle,-)8,.9,GlineStyle,-)n（4 4）隐藏图形包围盒的边线，进行透视）隐藏图形包围盒的边线，进行透视投影：投影：nset(gca,Box,off,Projection,perspeset(gca,Box,off,Projection,perspective)ctive)显示一个旋转的地球显示一个旋转的地球n（5 5）用）用MM文件旋转球体：文件旋转球体：n spin spinn（6 6）创建海平面数据网格，使球体不透）创建海平面数据网格，使球体不透明

48、：明：n base=zeros(180,360); base=zeros(180,360);n baseref=1 90 0; baseref=1 90 0;n hs=meshm(base,baseref,size(base);hs=meshm(base,baseref,size(base);n colormap copper colormap copper显示一个旋转的地球显示一个旋转的地球n（7 7）提升网格，使它浮在球面上方）提升网格，使它浮在球面上方2.5%2.5%半径距离的地方：半径距离的地方：n setm(gca,Galtitude,0.025) setm(gca,Galtitud

49、e,0.025)n（8 8）再次旋转地球：）再次旋转地球：n spin spinn现在旋转速度变慢，因为每次旋转都要现在旋转速度变慢，因为每次旋转都要对对180180* *360360的网格着色。的网格着色。显示一个旋转的地球显示一个旋转的地球n（9 9）准备用地形地貌数据替换目前的球）准备用地形地貌数据替换目前的球面数据：面数据：n clmo(hs) clmo(hs)n load topo load topon（1010）将高程放大显示，绘制表面：）将高程放大显示，绘制表面：显示一个旋转的地球显示一个旋转的地球n topo=topo/(almanac(earth,radius)topo=to

50、po/(almanac(earth,radius)* *20);20);n hs=meshm(topo,topolegend,size(tophs=meshm(topo,topolegend,size(topo),topo);o),topo);n demcmap(topo) demcmap(topo)显示一个旋转的地球显示一个旋转的地球n（1111）再次旋转：）再次旋转：n spin spinn（1212）也可以应用光照：）也可以应用光照：n camlight right camlight rightn lighting phong; lighting phong;n material(.7,

51、.9,.8) material(.7,.9,.8)六数值计算六数值计算n解线性方程组解线性方程组nMATLABMATLAB用用 和和/ /运算符求解线性方程组：运算符求解线性方程组：nX=AB X=AB 求矩阵方程求矩阵方程AX=BAX=B的解的解 nX=B/A X=B/A 求矩阵方程求矩阵方程XA=BXA=B的解的解 n若若A A为为mm* *n n矩阵，则有如下三种情况：矩阵，则有如下三种情况：nm=n m=n 平方系统，可以求得唯一解平方系统，可以求得唯一解nmn mn 超定系统，至少能找到一组解超定系统，至少能找到一组解nmn mn 不定系统，解中至多有不定系统，解中至多有mm个非个非

52、 零元素零元素例：求下列线性方程组的解例：求下列线性方程组的解5426255452321321321xxxxxxxxxMATLABMATLAB程序为：程序为：A=2 -5 4;1 5 -2;-1 2 4A=2 -5 4;1 5 -2;-1 2 4；B=5;6;5B=5;6;5；x=ABx=AB多项式多项式n1 1 根根n 找出多项式的根，即多项式为零时的值，找出多项式的根，即多项式为零时的值，是许多学科共同的问题。是许多学科共同的问题。MATLABMATLAB能求解这个能求解这个问题，并提供其它的多项式操作工具。在问题，并提供其它的多项式操作工具。在MATLABMATLAB里，多项式由一个行向

53、量表示，它的里，多项式由一个行向量表示，它的系数是按降序排列。系数是按降序排列。n例如，输入多项式例如，输入多项式x x4 412x12x3 30 x0 x2 225x25x116116n p=1 -12 0 25 116 p=1 -12 0 25 116np =p =n 1 -12 0 25 116 1 -12 0 25 116多项式多项式n 注意，必须包括具有零系数的项。除非注意，必须包括具有零系数的项。除非特别地辨认，否则特别地辨认，否则MATLABMATLAB无法知道哪一项为无法知道哪一项为零。给出这种形式之后，可以用函数零。给出这种形式之后，可以用函数rootsroots求求多项式的

54、根。多项式的根。n r=roots(p) r=roots(p)nr =r =n 11.7473 11.7473 n 2.7028 2.7028 n -1.2251 + 1.4672i -1.2251 + 1.4672in -1.2251 - 1.4672i -1.2251 - 1.4672i多项式多项式n 在在MATLABMATLAB中，无论一个多项式，中，无论一个多项式，还是它的根，都是向量。还是它的根，都是向量。MATLABMATLAB按惯按惯例规定，例规定，多项式是行向量，根是列向量多项式是行向量，根是列向量。给出一个多项式的根，也可以构造相应给出一个多项式的根，也可以构造相应的多项式。

55、在的多项式。在MATLABMATLAB中，命令中，命令polypoly执执行这个任务。行这个任务。n pp=poly(r) pp=poly(r)；多项式多项式npp=real(pp) pp=real(pp) n%throw away spurious imaginary part%throw away spurious imaginary partnpp =pp =n 1.0000 -12.0000 0.0000 1.0000 -12.0000 0.0000 25.0000 116.000025.0000 116.0000n 因为因为MATLABMATLAB无限地处理复数，所以当无限地处理复数

56、，所以当用根重组多项式时，如果一些根有虚部，由于用根重组多项式时，如果一些根有虚部，由于截断误差，则截断误差，则polypoly的结果有一些小的虚部，这的结果有一些小的虚部，这是很正常的。消除虚假的虚部，如上所示，只是很正常的。消除虚假的虚部，如上所示，只要使用函数要使用函数realreal抽取实部即可。抽取实部即可。 多项式多项式n2 2 加、减法加、减法n 对于次数相同的若干个多项式，可对于次数相同的若干个多项式，可直接对多项式系数向量进行加、减运算。直接对多项式系数向量进行加、减运算。如果多项式的次数不同，则应该把低次如果多项式的次数不同，则应该把低次的多项式系数不足的高次项用零补足，的

57、多项式系数不足的高次项用零补足，使相同式中的各多项式具有相同的次数。使相同式中的各多项式具有相同的次数。多项式多项式n3 3 乘法乘法n 函数函数convconv支持多项式乘法支持多项式乘法( (执行两个数组的执行两个数组的卷积卷积) )。考虑两个多项式。考虑两个多项式a(x)=xa(x)=x3 32x2x2 23x3x4 4和和b(x)= xb(x)= x3 34x4x2 29x9x1616的乘积：的乘积：n a=1 2 3 4 ; b=1 4 9 16; a=1 2 3 4 ; b=1 4 9 16;n c=conv(a , b) c=conv(a , b)nc =c =n 1 6 20

58、50 75 84 64 1 6 20 50 75 84 64n 结果是结果是c(x)=xc(x)=x6 66x6x5 520 x20 x4 450 x50 x3 375x75x2 284x84x6464。两个以上的多项式的乘法需要重。两个以上的多项式的乘法需要重复使用复使用convconv。多项式多项式n4 4 除法除法n 在一些特殊情况，一个多项式需要除以另在一些特殊情况，一个多项式需要除以另一个多项式。在一个多项式。在MATLABMATLAB中，这由函数中，这由函数deconvdeconv完成。用上面的多项式完成。用上面的多项式b b和和c c n q , r=deconv(c , b)

59、q , r=deconv(c , b)nq =q =n 1 2 3 4 1 2 3 4nr =r =n 0 0 0 0 0 0 0 0 0 0 0 0 0 0此结果是此结果是c c被被b b除，给出商除，给出商q q和余数和余数r r。现在情况下。现在情况下r r是零，是因为是零，是因为b b和和q q的乘积恰好是的乘积恰好是c c。多项式多项式n5 5 导数导数n 由 于 一 个 多 项 式 的 导 数 表 示 简 单 ，由 于 一 个 多 项 式 的 导 数 表 示 简 单 ，MATLABMATLAB为多项式求导提供了函数为多项式求导提供了函数polyderpolyder。ng=1 6 2

60、0 48 69 72 44g=1 6 20 48 69 72 44ng =g =n 1 6 20 48 69 72 44 1 6 20 48 69 72 44n h=polyder(g) h=polyder(g)nh =h =n 6 30 80 144 138 72 6 30 80 144 138 72多项式多项式n6 6 估值估值n 根据多项式系数的行向量，可对多项式进行加、根据多项式系数的行向量，可对多项式进行加、减、乘、除和求导，也应该能对它们进行估值。在减、乘、除和求导，也应该能对它们进行估值。在MATLABMATLAB中，这由函数中，这由函数polyvalpolyval来完成。来完成