1、2014 高教社杯全国大学生数学建模竞赛 承 诺 书 我们仔细阅读了 全国大学生数学建模竞赛章程和 全国大学生数学建模竞赛 参赛规则（ 以下简称为“竞赛章程和参赛规则”，可从全国大学生数学建模竞赛网站下载 ） 。 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。 我们知道，抄袭别人的成果是违反 竞赛 章程和参赛 规则的 ， 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用 处 和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛 章程和参赛 规则，以保证
2、竞赛的公正、公平性。如有违反竞赛 章程和参赛 规则的行为，我们将 受到严肃处理 。 我们授权 全国大学生数学建模竞赛组委会 ，可 将 我们的 论文 以任何形式 进行 公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式 发表 等） 。 我们参赛选择的题号是（从 A/B/C/D 中选择一项填写）： A 我们的报名参赛队号为（ 8 位数字组成的编号） ： 11102007 所属学校（请填写完整的全名）： 杭州电子科技大学 参赛队员 (打印并 签名 ) ： 1. 傅凌燕 2. 王森贤 3. 吴晓雷 指导教师 或 指导教师 组负责人 (打印并 签名 )： 数模组 （ 论文纸质版与电子版中
3、的以上信息必须一致，只是电子版中无需签名。 以上 内容请仔细核对，提交后将不再允许做任何修改 。如填写错误，论文可能被取消评奖资格。 ） 日期： 2014 年 9 月 15 日 赛区评阅编号（由赛区组委会评阅前进行编号）： 2014 高教社杯全国大学生数学建模竞赛 编 号 专 用 页 赛区评阅编号（由赛区组委会评阅前进行编号）： 赛区评阅记录（可供赛区评阅时使用）： 评 阅 人 评 分 备 注 全国统一编号（由赛区组委会送交全国前编号）： 全国评阅编号（由全国组委会评阅前进行编号）： 1 嫦娥 三 号软着陆轨道 设计与控制策略 摘要 嫦娥三号是我国首次对地外天体进行软着陆的探测器。软着陆时对月
4、速度较小，着陆后可继续在月面进行考察，具有较大实际意义。本文对其软着陆轨道设计和控制策略做一个探究。 在软着陆轨道的设计上，首先，我们 以月心为原点 ，建立空间直角坐标系。由开普勒定律和机械能守恒，得到远日点的速度为 1667.19m/s，近日点的速度为 1748.04m/s。又考虑到轨道平面经过月心和着陆点，我们认为在轨道平面与赤道平面的二面角等于着陆点的月心连线与赤道平面的夹角时，轨道平面最优。根据着陆点，月心，近日点三点共面原理得出轨道平面方程为 14.20x-5.04y+15.55z=0；根据速度大小和平面方程，得到近日点速度方向为（ -1192,394,1216），远日点的速度方向为
5、（ 1136,-375,-1160）。再对这些数据结果进行分析，从而确定轨迹图。 在主减速阶段的控制策略上，考虑燃耗最少为目标函数建立多元线性规划模型。然而由 于该模型计算量庞大以及实时数据监测的困难性，很难得到解析解，所以我们参考了一种显式制导方法。首先在一些合理的假设基础上对模型进行简化处理，然后通过查阅文献，确定一些模型的参数。在实现这两者的基础上，进行仿真分析。发现各参数关于时间的曲线都比较平缓，不存在不稳定现象，从而说明采用显式制导律可以得到一个燃耗最优的减速轨迹，且具有一个较好的动态特征。 由于测量误差和系统参数偏差是影响嫦娥三号软着陆制导精度的主要原因，所以我们接着对减速阶段建立
6、的模型进行误差和敏感性分析。选取发动机推力为偏差参数，考察其大小的改变对 制导过程的影响。我们用 MATLAB 进行绘图分析，比较发现：主减速过程中，该轨迹对推力 F 误差敏感性较弱。从而说明显式制导律只和当前的状态信息、终端约束有关，具有较好的鲁棒性。 在避障阶段的控制策略上，我们通过 MATLAB 对已给的高程图进行数据的读取， 得到 2300*2300和 1000*1000两个矩阵。为了减小计算量，我们将两个矩阵分为 10*10个区域，计算每个区域的方差。以区域中心为起点，方差大小来作为每个区域的平坦程度的指标， 当区域方差大于平均方差时，认为该区域是危险区域，需要规避，否则为安全着陆区
7、域。同时 以燃耗最优为目标采用螺旋搜索算法进行搜索寻找。随之用 MATLAB 对该算法进行编程，绘出安全着陆区域的位置。 关键字：轨道设计 燃耗最优 显式制导律 仿真分析 螺旋搜索算法 2 1.问题重述 1.1问题由来 嫦娥三号于 2013年 12月 2日 1时 30分成功发射， 12月 6日 抵达 月球轨道 。嫦娥三号在着陆准备轨道上的运行质量为 2.4t，其安装在下部的主减速发动机能够产生1500N到 7500N的可调节推力，其比冲（即 单位质量的推进剂产生的推力 ）为 2940m/s，可以 满足调整 速度 的 控制 要求。在四周安装有姿态调整发动机，在给定主减速发动机的推力方向后，能够自
8、动通过多 个发动机的脉冲组合实现 各种 姿态 的调整控制。嫦娥三号的 预定着陆点 为 19.51W， 44.12N， 海拔 为 -2641m 1.2问题要求 嫦娥三号在高速飞行的情况下，要保证准确地在月球预定区域内实现软着陆，关键问题是着陆轨道与控制策略的设计。其着陆轨道设计的基本要求：着陆准备轨道为近月点 15km，远月点 100km的椭圆形轨道；着陆轨道为从近月点至着陆点，其软着陆过程共分为 6个阶段要求满足每个阶段在关键点所处的状态；尽量 减少软着陆过程的燃料消耗 。 根据上述的基本要求，请你们建立数学模型解决下面的问题： （ 1） 确定着陆 准备轨道近月点 和远月点的 位置 ，以及嫦娥
9、三号相应速度的大小与方向。 （ 2） 确定嫦娥三号的着陆轨道和在 6个阶段的最 优控制策略。 （ 3）对于你们设计的着陆轨道和控制策略做相应的误差分析和敏感性分析。 2.问题分析 在软着陆轨道的设计上，我们查阅文献，得知飞行器在减速阶段的燃耗会随着其初始速度的减小而减少。为了使近地点的速度最小，着陆轨道平面与月球赤道平面的夹角应最小。又考虑到轨道平面经过月心和着陆点，我们认为在轨道平面与赤道平面的二面角等于着陆点的月心连线与赤道平面的夹角时，轨道平面最优。为了计算方便，我们建立月心坐标系，求出软着陆轨道平面的方程和近日点，远日点的坐标。最后结合开普勒定律，得到嫦娥三号在近日点和远日点的速度。
10、在主减速阶段的控制策略上，我们 以燃耗最优性为出发点建立多元线性目标规划。然而由于计算量庞大以及实时数据监测的困难性，该模型很难得到解析解。所以我们参考了一种显式制导方法。通过仿真计算和敏感性分析，发现这种制导方法可以实现燃耗最优的软着陆，且与推力大小等既定参数关联较小，只与初始状态有关，因此具有较好的鲁棒性。 在避障阶段的控制策略上，我们通过 MATLAB 对已给的高程图的数据进行读取并划分区域。计算每个区域的方差，并与平均方差作比较。当区域方差大于平均方差时，认为该区域是危险区域，需要规避，否则为安全着陆区域 。同时以燃耗最优为目标采用螺旋搜索算法，从区域中心出发进行螺旋式搜索，直到找到安
11、全着陆区域。 3 3.问题假设和符号说明 3.1问题假设 1) 着陆准备轨道的选取不考虑太空漂浮物的影响 ； 2) 软着陆阶段忽略辐射，粒子撞击对嫦娥三号运动速度与方向的影响 ； 3) 忽略地球等其他星体引力对嫦娥三号的影响 ； 4) 假设嫦娥三号在轨道运行中不存在内在故障或在外在事故； 5) 考虑嫦娥三号在姿态调整时，对其主减速过程的轨道不产生影响，将嫦娥三号看做一个质点。 3.2符号说明 G 万有引力常数 嫦娥三号质量 月球质量 d 嫦娥三号在近月点距月心的距离 v 嫦娥三号在近月点的速度 角速度 水平速度 径向速度 径向加速度 水平加速度 推力 F 的径向夹角 比冲 4.模型建立 嫦娥三
12、号是我国首次地外天体软着落任务。所谓软着陆，是指着陆器在制动系统作用下以很小的速度准确降落到月球指定区域，以保证实验设备的安全。 现在国际上，总的来说，月球软着陆有两种形式：一是自地月转移 轨道直接实现软着陆；二是自月球停泊轨道变轨至近月点然后实现软着陆。同第一种直接软着陆相比较，自环月轨道开始的软着陆方案具有很长的软着陆准备时间、可选择更大的着陆区域、可减少着陆舱部分的燃料消耗等优点。嫦娥三号选用的就是这样着陆方式。 4.1软着陆轨道 分析 从月球停泊轨道下降到月面，软着陆过程大致可成为 4 个阶段 1，如图 1-1,1-2 所示。这四个阶段中的每一个阶段都不是独立地，前一个阶段要完成的动作
13、，要考虑到后面几个阶段的技术要求。 4 图 1-1 软着陆阶段轨道示意图 图 1-2 软着陆各阶段分解图 如上述两幅图所示，嫦娥三号的 4 个阶段大致是这样的： 1) 离轨段。嫦娥三号到达月球附近时，通过变轨进入月球停泊轨道，它是一条近圆的卫星轨道。 2) 自由下降段，也称霍曼转移段 2。嫦娥三号脱离原来的停泊轨道，转入过 渡轨道。过渡轨道是一条近月点在所选着陆点附近的椭圆轨道。 3) 动力下降段，也称主动制段。嫦娥三号沿着过渡轨道下降到距离月面一定高度时，制动发动机点火工作，开始动力下降。 4) 避障阶段（包括姿态调整段），也称最终着陆段。嫦娥三号水平方向为 0，为在选定区域内找寻平坦着陆点
14、，进行粗细避障侦查。 5 4.2 霍曼转移轨道探测器最优 轨道平面 模型 嫦娥三号探测器在制动冲量结束后，脱离原来的停泊轨道，转入一条远月点高度约为 100km，近月点高度约为 15km 的椭圆运行轨道。则采用这个无动力下降段一般比从停泊轨道连续采用有动力下降的方式更节省能量。 项仁寿教授 2已经说明探测器沿着轨道相切的椭圆运行时，当所需出发速度越小，则近日点速度越小，才能保证下一阶段 主动制段所需要消耗的燃料越小，从而达到整个软着陆过程消耗的燃料越小。 4.2.1基于万有引力的模型分析 嫦娥三号在轨道上运动，当其运动到近月点时，若它和月心 的距离为 d，则该点处嫦娥三号受到月球的引力和其由于
15、运动产生的向心力分别为： (2.1) (2.2) (G 为万有引力常数； 为嫦娥三号质量； 为月球质量； d 为嫦娥三号在近月点距月心的距离； v 为嫦娥三号在近月点的速度 ) 由于在该点处 = ，则有 (2.3) 根据 (2-3)我们可得出结论：嫦娥三号在近日点的速度 v 随着它与月心距离 d 的增大而减小。我们已知近日点到月表距离为 15km 已恒定，所以要使模型达到最优，需要控制同一水平线上，月心到月表面的距离。 4.2.2月心坐标系的定义 为了方便计算，我们建立月心坐标系 O,x,y,z：以原点 O 为月心， 轴垂直月球赤道面指向月球北极， 轴在月球赤道面内， 轴和与其他两轴垂直并构成
16、右手坐标系。如图 2-1 所示。 图 2-1 月心坐标系 O,x,y,z 6 4.2.3软着陆轨道平面确定 根据图 2-1 我们可以做以下分析：假定月心 0，着陆点 P（ 1175.376， -416.453,1209.24）和近日点 B 三点控制的轨道平面为 ，赤道平面为 ， 与 之间的夹角为 。在满足以上条件下，为了使月心 0 到近日点 B 距离 d 最大，则月表距月心最远点 A 应该与 OB共线，即 应该达到最小值。 综上分析：嫦娥三号着陆点、月心和表面距月心最远点三点可控制一个平面，则近日点只有在该平面上，且在表面距离月心最远点和月心的延展线上时，才能保证嫦娥三号在近日点达到速度最优。
17、 由此，我们建立模型如下： max cos = (2.4) (2.5) 图 2-2 与 二面角示意图 如图 2-2 所示： n 为平面 法向量， z 轴为平面 法向量。根据二面角原理以及外积的几何意义 3,得到算法如下： (2.6) 综合式 (2.5)和式 (2.6),得到轨道平面方程为 14.20x-5.04y+15.55z=0 (2.7) 4.2.3基于开普勒公式确定两点速度大小 嫦娥三号绕着月球做椭圆轨道运动，不妨将它看做是常规一个天体运动。对于行星运动定律，开普勒第三定律证实 7 近 远远 近(2.8) (T 为嫦娥三号绕月球椭圆运转的周期； d 为行星公转轨道半轴长； k 为只用中心
18、天体质量有关的常数 ) 由于嫦娥三号在椭圆轨道运行时没有大气摩擦的干扰并且不考虑粒子、辐射等其他因素对其的影响，所以其在这无动力下降过程中保持了机械能守恒 近 远 (2.9) 结合式 (2.8)和式 (2.9)，得出近日点和远日点速度大小分别为： 近 1748.04m/s 远 1667.19m/s 4.2.4结果分 析 综上分析，在所建立得月心坐标系中，我们得到结果如表 1 表 2-1 近日点与远日点结果分析表 速度大小 速度方向 坐标 近日点 1748.04m/s （ -1192,394,1216） (0,-1666.7,540.2) 远日点 1667.19m/s （ 1136,-375,-
19、1160） (0,1733.4,-561.8) 嫦娥三号到达此位置和速度时，才能保证在下一阶段主动制段消耗的燃料最少。 4.3主减速阶段最优策略模型 嫦娥三号从霍曼转移轨道离轨进 入主减速阶段后，自近月点开始下降，发动进入动力下降段，也就是主动制段。在该阶段，其发动机点火工作，开始动力下降段。这个阶段的主要任务在于消除探测器速度的水平分量。由于探测器在月面上软着陆只能依靠制动发动机的能量来实现，因此，在这一过程中，如何使能量最省成为所研究的关键问题。 4.3.1燃耗最优轨迹模型建立 由于姿态控制系统的对称性，在姿态调整时，只产生力偶，对着探测器的轨道不产生影响 4。因此，我们在建立最优轨迹时候
20、，忽略探测器姿态影响，将其看做一个质点。只以考虑燃耗最少为目标函数建立多元线性规划模型。 (3.1) （ 表示角速度； 表示水平速度； 表示径向速度； 表示径向加速度； 表示水平加速度； 表示推力 F 的径向夹角； 表示比冲） 结合已给出的坐标系，化为空间形式为 8 (3.2) 实现燃耗最优着陆定义如下目标函数 (3.3) 然而由于计算量庞大，一般来说，式（ 3.3）很难得到解析解。所以在求解模型的时候我们不考虑运用极大值原理得到的制导律，而是采用显式制导律。 4.3.2显式制导律 显式制导律，即把控制量由各种可测得的状态以及时间表示出来这种制导方法无需迭代计算，满足实时性要求。 假设制动发动
21、机工作在全推力状态，此时，燃耗最优问题就是时间最优问题。由式(3.1)可以看出，软着陆制导控制系统是一个强非线性系统。为求得显式制导律，需要在一些合理的假设基础上对其进行简化处理。 这种显式制导律虽然需要探测器位置、速度、推力加速度的信息，这些信息量需要雷达高度表、多普勒雷达等多种先进仪器测得，但我们可以查阅文献，确定一些参数，并在模型上进行一定的简化。 表示初始质量，根据文献 5燃耗最优的推力方向角可以表示为时间的一次多项式 (3.4) 其中 是只满足终端速度要求的最优方向 、 、 均为常值。 根据文献 6，在垂直方向的最优着陆轨迹可由一关于时间 t 的 3 次多项式完全表示 (3.5) 假
22、设，在燃耗最优的条件下，制动推力主要用来抵消着探测器的水平初始速度。设 为推力加速度的水平分量 (3.6) 由于 、 为定值， 变化区域在 01 之间，所以 在该阶段变化很小。我们不妨假设 为定值，则嫦娥三号在减速阶段水平方向做匀减速运动。 4.3.3仿真分析 根据文献和题目，仿真参数设置如下： 月球半径： 1737.013km； 引力常数： =4.9 ； 初始质量： ； 发动机比冲： ； 发动机推力： F=7500N 9 表 3-1 着陆边界约束 始端约束终端约束 终端约束 1737013m 174000m 0 自由 1692m/s 0 0 -1.5m/s 仿真结果如下图所示： a）径向速度
23、 时间 t b）径向位移 时间 t C）水平速度 时间 t d）水平位移 时间 t 0 100 200 300 400 500 600- 3 5- 3 0- 2 5- 2 0- 1 5- 1 0-50时间径向速度0 100 200 300 400 500 6001 . 7 41 . 7 4 21 . 7 4 41 . 7 4 61 . 7 4 81 . 7 51 . 7 5 2x 1 06时间径向位移0 100 200 300 400 500 600020040060080010001200140016001800时间水平速度0 100 200 300 400 500 60000 . 511
24、. 522 . 533 . 544 . 55x 1 05时间水平位移10 e）推力方向角 t 图 3-1 显式制导方法软着陆曲线 从上述图得到如下结果：采用显示制导方法律径向速度比较平缓，说明加速度变化小，这样有利于嫦娥三号各部件的正常工作；水平速度随着时间的增大而匀速变化，保证其在规定的时间内平稳地将水平速度降为 0；同时，推力方向角变化的平缓，有利于降低探测器控制系统的负担。 综上，应用显式制导律可以得到一个燃耗最优的减速轨迹，且具有一个较好的动态特征。 4.3.4误差 和敏感性 分析 测量误差和系统参数偏差是影响嫦娥三号软着陆制导精度的主要原因，所以我们要对所建立的模型进行误差和敏感性分
25、析。 模型中很多初始数据，例如：探测器初始速度，探测器初始高度等虽然都是在一定假设条件下得到数据，可能和实际数据存在一定误差，然而，这些数据我们都是通过大量文献查阅得到，具有一定的权威性和准确性，在建立优化模型中，其误差可以忽略不计。 系统参数的敏感性分析 系统参数如制动发动机的推力，比冲以及着陆舱的质量都是不可测得，它们 是 发射以前在地面标定给出的，但在飞行 过程中参数会由于一定的影响而产生一定的偏差。所以考察这些偏差对制导过程的影响就显得十分重要。 下面，我们就推力这一参数进行敏感性分析，得到结果如下： 0 100 200 300 400 500 6000 . 20 . 40 . 60
26、. 811 . 21 . 41 . 6时间推力方向角11 a） 径向速度 时间 t b）径向位移 时间 t 0 100 200 300 400 500 600- 4 0- 3 5- 3 0- 2 5- 2 0- 1 5- 1 0-50时间径向速度7600N7400N0 100 200 300 400 500 6001 . 7 3 81 . 7 41 . 7 4 21 . 7 4 41 . 7 4 61 . 7 4 81 . 7 51 . 7 5 2x 1 06时间径向位移7400N7600N12 C）水平速度 时间 t d）水平位移 时间 t 图 3-2 主减速过程各参数对推力 F 的敏感度分
27、析 如上图反映：显式制导律只和当前的状态信息、终端约束有关。而作为既定参数之一的推力 F 对其着陆效果并没有实质上的影响，即主减速过程中，该轨迹对推力 F 误差敏感性较弱。同时其满足终端约束条件，因此具有较好的鲁棒性。 4.4避障阶段的策略 优化 嫦娥三号经过上述两个阶段后水平方向为 0，在经历短暂的调准姿态阶段后，逐进行垂直下落。 已有的月球图像和高程数据表明 ,月球表面分布着各种高山壑谷 ,即使在相对平坦的月海地区也遍布着大小不一的岩石和陨石坑 .这种地形、地貌以及石块和陨石坑会给着陆器安全软着陆带来较大风险 .只有着陆器具有发现和识别障碍并进行机动避障的能力 ,才0 100 200 30
28、0 400 500 600020040060080010001200140016001800时间水平速度7400N7600N0 100 200 300 400 500 60000 . 511 . 522 . 533 . 544 . 55x 1 05时间水平位移7400N7600N13 能保证软着陆的高安全和高可靠 7. 嫦娥三号的核心任务是实施高可靠高安全的月面软着陆 ，是 世界上首次成功实现了利用机器视觉的地外天体软着陆自主避障 。 下文，我们给出了 针对嫦娥三号软着陆任务提出的粗 -精接力自主避障软着陆控制技术以及在轨飞行结果 。 4.4.1 粗精避障的不同目的 粗避障主要目的是在较大着陆
29、范围内剔除明显危及着陆安全的大尺度障碍 ,为精避障提供较好的安全点选择区域 ,避免出现近距离精避障避无可避的风险 。 此外 ,还要综合 考虑燃料的 消耗来选择最优位置 。 精避障主要目的是在粗避障选取的较安全区域内进行精确的障碍检测 ,务必识别并剔除危及安全的小尺度障碍 ,确保了落点安全 . 4.4.2螺旋搜索算法 用 MATLAB 读取分别在粗避障阶段和精避障阶段的高程图，得到 2300*2300 和1000*1000两个矩阵。为了减小计算量，我们将两个矩阵分为 10*10个区域，以方差 大小来作为 每个区域的平坦 程度的指标 。 探测器 在粗避障阶段垂直下落的过程中，我们采用从探测器正下方
30、即预定着陆区域中心，考虑逆时针螺旋前进搜索的方法。如图 4-1，图 4-2所示 图 4-1 螺旋搜索算法示意图 14 图 4-2 螺旋搜索算法流程图 4.4.3MATLAB编程算法找寻安全区域 我们用 MATLAB 对螺旋搜索算法进行编程，从而来选取着陆安全区域。找寻的结果如图 4-3所示： 15 图 4-3 探测点着陆安全区域图 4.4.4结果分析 如图 4-3所示，我们采用螺旋搜索算法找到的安全着陆区域位于高程图所示区域中心附近，且地势平担，无 高山壑谷 影响探测器着陆安全。同时，和传统的在高空先选定下落点再定位降落的方法相比较，螺旋搜索算法可以有效减少探测器燃料的损耗，从而达到选址安全，
31、燃耗最优的双效果。同时，该算法计算量小，适合月球表面这样的复杂地形的选址。 5.模型评价与分析 本文根据显式制导律， 把控制量由各种可测得的状态以及时间表示出来 ， 无需迭代计算，满足实时性要求 。在粗避障和细避障时，根据燃耗最优目标，选择螺旋搜索算法。在寻找安全着 陆区时，该算法计算量小，计算时间少，有利于对复杂地形的选址。 本文的不足之处在于软着陆时没有考虑月球科氏力和地球引力等因素对探测器着陆轨道的影响，从而产生误差。 6.参考文献 1 郗晓宁 ,曾国强 ,任营 ,赵汉元 .月球探测器轨道设计 .国防工业出版社 ,2001:112122 2 崔枯涛等 .附着小行星的视线制导律 .空间科学
32、学报 ,2002,(3):4749 3 Dayi Wang, Tieshou Li, Xingrui Ma. Optimal guidance for lunar gravity-turn descent. Acta Automatica Sinica,May, 2002, 28(3): 385390 4 任萱 .人造地球卫星轨道力学 M.长沙 :国防科技大学出版社 ,1988. 5 韩祝斋 . 用于大型运载火箭的迭代制导方法 . 宇航学报 , 1983, (1): 920 6 Ueno S, Yamaguchi Y. Near-Minimum Fuel Guidance Law of a L
33、unar Landing Module. Proc. Of 14 th IFAC Symposium on Automatic Control in Aerospace,Seoul, 1998: 337 16 附录 1（近日点远日点速度求值） clear syms x y z; x,y,z=solve(14.2*x-5.04*y+15.55*z=0,x2+y2+z2=1748.0442,-1666.7*y+540.2*z=0,x,y,z) 附录 2（数值仿真） figure b2=-167/676; b3=321/703040; t=0:1:520; vr=b2*t+b3*t.2; plot(
34、t,vr); xlabel(时间 ); ylabel(径向速度 ); figure r=1752000+b2/2*t.2+b3/3*t.3; plot(t,r); xlabel(时间 ); ylabel(径向位移 ); figure k=0.00183; c=k.*t; h=acos(c); plot(t,h); xlabel(时间 ); ylabel(推力方向角 ); figure h=3; vx=1692-h.*t; plot(t,vx); xlabel(时间 ); ylabel(水平速度 ); figure x=1692.*t-h/2.*t.2; plot(t,x); xlabel(时间
35、 ); ylabel(水平位移 ); 17 附录 3（粗避障） clear pic=imread(D:03.tif); pic=double(pic); a=meshgrid(1:1:10); b=meshgrid(1:1:230); for i=1:1:10; for j=1:1:10; for m=1:1:230; for n=1:1:230; b(m,n)=pic(230*i+m-230,230*j+n-230); end end a(i,j)=var(b(:); end end i=6; j=5; var0=mean(a(:); for m=1:1:100 if a(i,j)var0 if i+j+1=10 end if i+j+1=10 end if i+j+1=xR(1) end if i+j+1=10 end if i+j+1=xR(1) z2=pic; z2(z1(:)=1)=6; mesh(x,y,pic,z2);
