数学建模优秀论文-空中飞行器无源定位

1、2012年军队院校军事建模竞赛承诺书我们仔细阅读了军队院校军事建模竞赛的竞赛规则。我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网 上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的，如果引用别人的成果或其他公开的 资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在止文引用处和参 考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如冇违反竞赛规 则的行为，我们将受到严肃处理。我们参赛选择的题号是（从a/b/c/d中选择一项填写）：a我们的参赛队号为： 0512078所属学校（

2、请填写完整的全名）：信息工程大学参赛队员（打印并签名）：1. 何杰2. 张洋3. 赵永胜指导教师（打印并签名）:日期：2012年6月25日2012年军队院校军事建模竞赛编号专用页评阅编号:评阅记录:评阅人评分备注空中飞行器无源定位摘要本文针对空中飞行器无源定位问题，采用最小二乘、遗传算法、matlab仿真等方法, 得到了飞行器静止和运动时的定位、可靠性分析、卫星选择策略。针对问题一，根据基于测向夹角的飞行器无源定位方法，得出三颗星可以定位。本 文采用了两个模型，第一个使用遗传算法。对于84组数据的处理，我们采用了先删除误 差点大的数据，再利用遗传算法去逼近实际位置，求出飞行器的位置为 (-51

3、40.2,6512.5,3277.7) (km),距离地球表面25563.6 km；第二个角度误差平方和最 小原理，用最小二乘法，用方程组得出的角度去逼近角度的实际测量值，进而得出比较 接近实际值的数据，求出飞行器的位置为(-5201.2,6603.5,3125.6) (km),距离地球表面 2600.5 km,两种方法确定的飞行器的位置相差不大。针对问题二，我们利用20$时刻的位置和速度矢量表示出其余4个时刻的位置。 利用附表给出的5个时刻的测量数据，建立33个方程，并联立方程。利用最小二乘法， 使得各个时刻的角度与其实际测量值之间的误差平方和最小。从而逼近飞行器的真实位 置，获得飞行器20

4、$吋刻的位置为(-41828.8,4627.0,6190.6) (km),速度矢量 (4.0,2.9,-0.2)伙加/$)。然后利用求得的/ = os吋刻的位置和速度矢量，计算出t = 70s时 飞行器的位置为(-3902.8,4830,6176.6)(血)。对于t = 70s时位置的可靠性的分析，影 响可靠性的因素有两个，一是角度的测量误差，二是最小二乘法解方程组求岀的解与真 实解间存在偏差。分析出这两种谋差后，将其作为噪声加入到附表小的角度测量值，然 后利用这些加入噪声后的角度进行仿真定位实验，统计定位位置的变化。得岀在角度测 量的误差服从n(0,0.1)分布时x坐标的均值ux的置信度为9

5、5%置信区间为： -3871.6,-3943.5； y坐标的均值叭的置信度为95%置信区间为：4802.4,4857.5； z坐 标的均值冬的置信度为95%置信区间为:6109.4,6253.3 o针对问题三，本文以几何精度因子gdop作为卫星分布对定位精度影响的指标。通 过分析卫星数口及一颗卫星对gdop的影响，可知gdop随卫星数口增加单调递减，但 递减幅度变小。综合考虑卫星定位精度和定位效率，并根据不同情况卜对于定位精度和 定位效率的不同要求，给出了偏重精度、精度效率折中和偏重效率三种情况下，选星数 目分别为5颗、6颗和7颗三种数目的优选方案。并着重分析精度和效率折中的6星组 合方案。遍

6、历所有的6星组合方案，找出其gdop最小的卫星组合作为最终的优选方 案。在测量角度存在0.1°误差限的情况下，我们对附表中的数据加入最大为0.1°高斯噪 声，然后进行大量的仿真实验，观察加入噪声后产生的定位误并。统计定位误差，发现 加入噪声后有95.9%的定位结果误差在600如2以内，即定位精度可以认为600如2。在本文的最后，针对每个问题对其结果进行了分析、对每个问题解决方法的优缺点 进行了分析，并提出了相应的改进方案。关键词：基于测向夹角的飞行器无源定位；最小二乘法；遗传算法；gdop,精度分 析1. 问题的重述目标定位技术是导航与制导技术的重要基础。在现有的导航与制导

7、技术中，卫星定 位技术是精度最高的，也是较为理想的导航与制导技术。目前，较为成熟的卫星导航系 统有gps系统、galileo系统等。卫星定位的基木原理是目标接收机通过接收多颗卫星 的信号测量出目标距各卫星的距离（伪距），再通过一定的计算确定出目标的位置。对于空中飞行器，在其飞行过程中很容易接收到太空卫星的信号。现在考虑通过测 量飞行器与地球同步卫星的方向角來实现空中飞行器的自定位。在球心坐标系下，空中 飞行器p的空间坐标记为（齐不妨设它同时能接收到n颗同步卫星的信号，其n 颗同步卫星&的空间坐标分别记为仗略网）（1 = 1為n）。为了方便检测与同步卫星的 方向角，在空中飞行器上固定安装

8、了两个相互垂直的测向阵列，它们的指向分别为 血ss氐切血）和鬲（抵皿吋如）。地球同步卫星堆与空中飞行器p的位置关系示意 图如图所示，,氏分别裘示空中飞行器p的测向阵列方向右，血与地球同步卫星 （1=1x-n）的夹角。现在请你们建立数学模型研究解决下面的问题：（1）通过测量空中飞行器测向阵列方向心和血与多颗地球同步卫星的夹角阿和何， 建立空中飞行器定位的数学模型；对于附表1所给出的9颗同步卫星的数据，试确定空 中飞行器p的位置参数。（2）在某些特殊情况下，空屮飞行器能直接检测到的同步卫星数量较少，可以利 用空屮飞行器在匀速飞行过程中多次检测的结杲来实现定位。针对这种情况，试建立空 中飞行器定位的

9、数学模型；对附表2中给出的3颗同步卫星的检测数据，确定空中飞行 器p在第70秒吋的位置参数，并分析其可靠性。（3）当可用同步卫星数量较多时，为了提高定位精度和定位效率，需要对可用的同 步卫星进行一定的优选。试研究具体的优选策略，并通过仿真，分析在检测方向角误差 限为0.1。时空中飞行器的定位方法和精度。2. 问题的分析木题旨在研究分析基于测向夹角的飞行器无源定位。根据测向阵列方向和地球同步 卫星夹角对空中飞行器无源定位，需耍分析研究飞行器静止和运动吋的定位、卫星优选 等问题。针对题目中的问题，我们进行了如下的分析过程。2.1飞行器静止时的定位分析由于附件1列岀了某时刻地球同步卫星的经度和测向阵

10、列方向和地球同步卫星夹角 等信息，为了更好地描述地球同步卫星和飞行器的位置以及测向阵列方向的关系，所以, 我们建立了以地球球心为中心的大地坐标系。在统一的坐标系下，需要分析测向阵列方向和地球同步卫星夹角与飞行器位置的 关系，此吋有三个矢量：飞行器和卫星连线矢量和两个测向阵列方向矢量，9个未知数: 飞行器的位置、两个测向阵列的方向，需要9个方程。两个测向阵列本身相互垂直，模 为1可以确定3个方程，另外，已知了矢量间的两个夹角，利用矢量点乘的关系构造方 程（每颗卫星可以构造两个方程，3颗卫星确定6个方程），因此3颗卫星总共可以确定 9个独立方程，即3颗卫星可以确定飞行器的一个位置。但是考虑到3颗卫

11、星定位的误 差非常大，精度难以满足要求，因此应该联立其他卫星的测量数据，共同确定飞行器位 置。题口给出了 9颗卫星，对于这9颗卫星测量数据的利用，我们有以下两种方法：在3颗卫星和测向阵列方向的约束关系的基础上建立9个方程，然后利用最小二乘 法，求出飞行器的一个位置，由于卫星共有84种组合，所以共有84个位置。但是 怎样利用84组数据來逼近飞行器的实际位置呢？我们的求解思路如卜stepl:求出位置的平均位置元、y> iostep2:算岀每个点到平均位置的距离,去除距离最大的那一点，求出剩下数据的 ©平方和的平均值。step3:重复第一步和第二步，在距离平方的平均值的一组数据屮找出

12、最小值，然后 把最小值z前数据点全部删除，接着求出剩下数据点的平均值和变化区间。step4:把第三步求出的区间当作飞行器位置的取值区间，用遗传算法寻找此位置区 间中方程误差平方和最小吋的位置，我们认为此位置就是飞行器位置的最优解。（2）利用方程算出的角度，与实际测量角度值的误差平方和作为冃标函数，则当9 颗卫星的所有角度与其实际测量值的谋差平方和最小时（用最小二乘法），可逼近飞行 器的真实位置。2. 2飞行器运动时的定位分析当同步卫星数量较少时，可利用飞行器在匀速飞行过程屮多次检测的结果来提高定 位的精度并确定速度参数。由问题（一）的分析可知，三颗卫星即町确定飞行器的位置，利用五组数据，可以

13、分别确定五个时刻的位置，然后除以时间，即可求出每一段的速度。但是，由于三星定 位误差过大，相比丁 10秒的位移量，三星定位的误差是难以容忍的。所以我们冇以下 思路：假设飞行器是匀速直线运动，设岀飞行器在t = os时位置和速度，t = os后的位置 可以由t=qs时的位置和速度表示。利用5个时刻的数据，建立方程组。利用最小二乘 法，使角度与其测量值z间的误差平方和最小，从而求出最小二乘解，即飞行器的位置 参数和速度参数。然后根据t二os时的位置和速度，即可预测t=70s时的位置。对于预测位置的可靠性的分析，影响可靠性的因素主要有两个:一是角度的测量误 井；二是最小二乘法解方程组求出的解与真实解

14、间的偏差。通过定量分析出这两种误井, 将其作为噪声加到附表中的角度测量值上，然后利用这些加入噪声后的角度进行仿真定 位实验，接着对定位位置进行统计，从而可以给出以置信度为0. 95的置信区间。2. 3卫星选择策略分析当可用同步卫星数量较多时，为了提高定位精度和定位效率，需要对可用的同步卫 星进行一定的优选，因此需要建立具体的优选策略，也就是如何选，选儿颗，选哪儿颗。 而选择几何分布好的卫星组合可以提高定位精度和定位效率。通常我们用几何精度因子 gdop來表征用户和可见卫星在空间几何分布的好坏。从m颗卫星中选择n颗卫星， 共冇c：种选择方案。从可见卫星屮选择参加导航定位计算的卫星数目不同，gdo

15、p的 取值也不相同，gdop与卫星数目之间有一定的变化规律，分析卫星数目和gdop之间 的关系，从而找到一个合适的卫星数目。同时，卫星数目越多，意味着方程数量越多， 计算量越大，定位效率也就越低。因此，为了兼顾定位的效率，应该在保证定位精度的 前提下，尽可能减小定位卫星的数量。根据对定位精度和定位效率的要求，确定一个合 适的卫星数目。确定一卫星数目为n之后，还需对一卫星的几何布局进行优选，遍历从所冇 卫星中选择n颗卫星的所有方案，从中选择gdop最小的卫星组合，即为优选的结杲。确定定位卫星后，给测量角度加上最大为0.1。的高斯噪声，进行多次的仿真实验，统 计定位精度。3. 模型的假设与符号说明

16、3. 1模型的假设（1）假设附表中所有的数据都是真实可靠的；（2）假设飞行器是一个质点，不考虑飞机的飞行姿态；（3）卫星信号强度都足够的强；（4）不考虑地球的实际曲率变化，认为地球是一个均匀球体。3. 2符号说明同步卫星x,的位置e行器的位置飞行器指向同步卫星x,方向的矢量航与测向阵列瓦住“皿切俎）的夹角就与测向阵列屁（抵禺、虫3）的夹角飞行器的速度矢最 每个点到平均位置的距离%a与实际测量值之间的误差 a与实际测量值z间的误丼&dld2%2测向阵列d|,2相互垂直关系确定的方程的解的误差测向阵列模为1确定的方程的解的误差测向阵列 2模为1确定的方程的解的谋差4. 模型的建立（与求解）

17、4. 1问题（一）的模型建立与求解4. 1. 1坐标系的建立建立以地球中心为原点，以赤道平面为基木平面，x轴在基本平面内由地心向外指 向格林威治子午圈与赤道的交点，z轴垂直基本平面，与地球自转轴重合，y轴与x, z轴 组成右手系的三维直角坐标系，如图（1）所示。p(“,y,z)图（1）建立坐标系在此坐标系下，同步卫星的位置是凤仗“幵石）（1 = 1,29）,飞行器的位置是 p（x,yrz）,测向阵列方向矢量为瓦（仏山切cf“）和石s比覘严血j,矢量忒和可、石的 夹角分别为w、仇。则利用矢量点乘关系a 巨=间|司可以得到如卜非线性方程组：（阿石=i醱ii石i1 j 丄斗 q （1 =切4遍=|也

18、|亦隔除由两个测向阵列本身相互垂直、模为1,得国| - 1 =式中,srfltfs >盹i、%2分别为计算误差。4.1.2基于遗传算法的模型的建立由上面的分析可知，矢量点乘关系px 一co<xt = aafi込补石一 i丽石|c“a =瞅务、呗为计算课差，联立3颗卫星和测向阵列关系，建立9个方程组，使z满足如卜函数优化问题：min吝=叱+嘱）+雉倔+ &紅+碓2 t=l気t x2 + ya + s2 > r2式中r为地球半径，表示飞行器到地球质心的距离大于地球半径，利用最小二乘法就可以求解出飞行器位置p(yzx),测向阵列忑(爲屮圧切和忑(怠罔尸爲m)这9个 未知参数

19、。所以存在84个位置 得到一个位置数据比 曲于目标函数是一个 为此需耍把此区域当由于毎3颗卫星就可以确定飞行器的一个位置，飞行器共有c= 84种不同的位置。 由于在实际测量时，测向阵列和地球同步卫星夹角存在测量误差， 是合理的。为了得到精确的位置，我们采用删除误差位置的方法， 较集中的区域，因为遗传算法寻找方程的解时需要一个初始区间， 多元非线性特别复杂的方程，没有全局最优解，只有局部最优解， 作飞行器位置初始区域，再用遗传算法寻找飞行器位置的局部最优解，具体的步骤如下:4.1.2. 1删除部分误差较大的数据stepl： 设定初值k二84.step2：利用平均公式求出位置的平均值壬、予、云k牙

20、二 一kk»=1k(1)式中，k为剩余的位置个数。step3：按式(1)算出每个点到平均位置的距离兮，删除距离最大的那一点，令k=k-l,按式(2)求出距离均方误差酉。xu,-x)2+(；.-刃2 +忆 _汀0=7kk玄二上(2)k式屮，k为剩余的位置个数。step4：重复step2和step3,找出距离平均值最小的那一，点，把均值最小点z前的 位置全部去除掉，算出均值最小点之后的平均位置和位置所在的区间。鮒址口牡ms0.81020304050数据点删除个数60图（2）删除部分误差较大的数据按如上步骤，画出了以位置剩余点数为横轴，距离均方误差为纵轴的趋势图，如图 （2）所示，由图可以

21、看出红点（去除了 25个位置时）所在处距离平均值最小，说明去 除25个位置后，误差较大的点已经删除，数据更加集屮。算出位置所在的区间范围， 如表（1）所示。表（1）删除误差较大点后的位置范围坐标轴下限（单位：km）上限（单位：km）x-5227. 24-5071.48y6534.756630. 44z2717.053336. 544. 1.2. 2利用遗传算法搜寻最优解群体大小：800 参数范围：见表（1）交叉概率:0.8变异概率：0.2终止代数：200目标函数:式屮，勺为计算误差。由于目标函数是一个多元非线性特别复朵的方程，没冇全局 最优解，只冇局部最优解，为了能更好的求出飞行器的最优位置，

22、算法流程如图（3）所示。图(3)遗传算法流程图对应的具体步骤如下：stepl：设定x、y、z的初始空间和迭代步长;step2：保持%、% z的上限不变,下限分别以步长逼近上限;step3：调用mat lab遗传算法工具箱求出此区间内的局部最优解并保存，重复step2直 到跳出循环；step4：输岀所冇局部最优解屮，最小目标函数值对应的兄、¥、z坐标。用matlab的遗传算法工具箱，画出目标函数适应度变化趋势如图(4)所示：best fitness mea n fitness0)n-ea ssasiz406080100120140160180200gen eration图(4)遗传算法

23、结果求出飞行器的位置为（一 514山2£5口5吕（km）。4.1.3基于角度误差平方和最小模型由丁误差主要是由夹角测量不准和求解是带来的，上述模型不能直接分析误差与测 量夹角偏差的关系，因此，基于角度误差平方和最小模型，从上面矢量点乘的分析可知,at 一 arc cos兀-x)%+(x-ym“+(&-z)d)"+ (x - y) + (召 - z) x jdj+di'+d 兀一 x)+ b 一 y)d2y + (& - z)d 土2)，pi - arc cos j(兀一汀 + ( x - y ) + ( « - x (2+2)，2 + d

24、2?久、切为计算误差，同理联立9颗卫星和测向阵列关系，建立21个方程，使之满 足如下函数优化问题：9+嗽）+刃紅+疡+血他t=is.t. x2 + y2 + z2 > r2式中r为地球半径，表示飞行器到地球质心的距离人于地球半径，利用最小二乘法就可以求解出飞行器位置巩爲刃耳，测向阵列忑（瓦疋炉和忑（転心尸爲m）这9个 未知参数。求岀飞行器的位置为（一 5zm2q6w3久25q （km）。9o1best: 367578.8159 mean: 381861.84854. 2问题(二)的模型建立与求解4.2.1飞行器运动时的定位模型图(5)匀速运动的飞行器各时刻位置如图(5)所示，假设飞行器在

25、(=0时刻的位置为(x,y,z),飞行器作匀速直线运动, 设其速度矢量为v = (vx,vv,v.),则飞行器在各个时间的坐标可分别表示为：(x,y,z), (x+vxt,y+vyt,z + v:t) ,(x+2v/, y+2vyt, z + 2v j),(x+3 匕/,y+3vj,z + 3vj),(兀+4叮,)，+4叩採+ 4叩)。利用这5个时刻的数据，根据飞行器指向卫星的矢量与测向阵 列方向的点乘关系，建立下列方程：arccosarccos（21,2,3,4,5）（21,2,3,4,5）（兀一兀一咕，x 一 v 一 m 心一 z - f） （d“，卜,dj、j（兀 _ 兀 _）2 +（x

26、 _ y _ v/ ）2 + （z, _ z _ 匕厶）2 x jd j + d j）丿（兀 一 - 匕人，x y 忖，& 一 z 叩j （2* 性.,2 jj（xj_x_匕右）2 + （开 _ y _片石）2 + （& _ z_匕）2 x風：+ x； + df）联立这10个方程，利用最小二乘法，使各个时刻的角度与其实际测量值z间的谋差平 方和最小，建立目标函数：5min 心工(廟 + 磅)+昉 + /2 + 4.2(3)/=!利用最小二乘法，获得飞行器的t二os时刻的位置为(-41828.8,4627.0,6190.6) (km), 速度矢量(4.0,2.9,-0.2) (k

27、m/s)o因此我们可以计算在t二70s时，飞行器的位置为(-3902.&4830,6176.6) (km)。4. 2. 2飞行器运动时的定位模型的可靠性分析对于预测的"70s吋位置参数，影响其可靠性的因素有：(1) 最小二乘法解方程组给出的解与真实解间存在偏差：角度的测量值与真实值z间存在误差。用最小二乘法求解式(3)时，最后的角度 误差平方和为4.6459x10°,则每个角度的误并平方最大为4.6459x10"°,误并为 ±2.1554x10_5o(2) 角度的测量误差随机产生，我们可假设角度的测量误差符合高斯分布 (0,0.1) 0

28、这两个因素影响了求出的飞行器位置和速度的可靠性，从而影响飞行器t=70s时位 置参数。给附表2中的测量角上加上上述两种随机误差，然后进行100次的仿真，得到位置 参数和速度参数大量计算值。然后得到100个t二70s时位置参数。如图(6)所示：(0,0.02)图(6)不同测量误差对位置的影响图(6)屮的四幅图分别是方差为m2、005、010和015时的预测位置的空间分布 图，可以看出噪声方差越小时，预测位置的空间分布越集中，即说明了误差比较小时， 对位置的空间分布影响较小；而方差越大时，预测位置的空间分布越分散，即说明了误 茅比较大时，对位置的空间分布影响较大。统计不同噪声下的位置分布，得到叭匕

29、和 冬的置信度为95%置信区间，如表(2)所示表(2)不同噪声下预测位置的置信区间/ = 0.02ff 2 = 0.05c- = 0.1er2 = 0.15-3905.7,-3895.13943.5, -3871.6-3916.1,-3883.4-3948.8,-3855.7wv4815.1, 4826.94802.4, 4857.54795.9, 4825.14664.1,4766.1l6110.1, 6154.3j6169.4, 6293.3j6055.1,6168.56055.1,6053.4表(2)给出了不同噪声卜预测位置的置信区间，从中可以看出，噪声越大，即测量误 差越大，相同置信度

30、下，x, y, z坐标的置信区间的长度也就越大。这意味着位置确定 的可靠性降低。因此，对于模型确定的飞行器t二70s时的位置，随着角度测量课差的不 同，其可靠性也就不同。角度测量误差越小，可靠性越高；反之，可靠性越低。4.3问题(三)的模型建立与求解4. 3.1基于几何精度(gdop )因子的优选模型定义几何精度因子如下:gdop = + £)22 + d33 + »441111式中，2是g = (l h)7的对角线元素。h是观测矩阵，q(心1,2,3n)是沿飞行器接收机指向卫星的直线方向的单位矢量。从上式可以看岀，gdop仅仅是平星和用户集合布局的函数。为了提高定位精度,

31、 应选择可见卫星和用户有较好的空间儿何分布，即gdop值最小。4. 3. 2卫星数目及一颗卫星对gdop的影响:从可见卫星中选择参加导航定位计算的卫星数目不同，gdop的取值也不相同， gdop与卫星数目z间有一定的变化规律。设0为选择加颗卫星定位时的观测矩阵，从加颗卫星中去掉第i颗星 (21,2,3肋，等到加颗卫星的观测矩阵,两者有如下关系：h：h肿(hqth丄十h轨其中，也二冋,勺2，%，1为第，颗星的观测矢量。由 she rm an -mo rrison 公式可得g = (hj= (h：叽-町订' =gm+gmhf(-higmhfyihigm其中，sh：是一标量，记为因此有gdo

32、p；j = traced) = traced) + teace(gmh；/ 5,)trace(gn) + teace(gmh/ s”)= gd0p： +trace(gmhlhigm / su)去掉一颗卫星后的gdopj比gdo畫多出一项，trace(gmhh,gm / su)项，确定多 出这项的正负，即可确定减少一颗卫星对gdop的影响。为确定tracedhtgtn is.)的正负，将h鳥进行奇异值分解。= usvt其中,u和v分别为(m-l)x(/h-l)和4x4的止交矩阵，是(/h-l)x4对角线矩阵。将式带入，并将等式两边进行相同的止交变换，由sherman-morrison公式得出vt

33、(hlhm)v' = vt(usvt)t usvt)v +vth hy = (z + vvr)_,其中 z = sts =/dg(zii，z22，z33，z44)， v = vthf =vj,v2,v3,v4o 由于正交变化不改变矩阵的迹，因此冇tracehj) = traced) = trace(vt + tracekvytz)= trace(h；hj') + trace(kvytz)其中，v. =z_,v ,n vtrace(kv：v )匹丄因此，4. 3. 3优选方案gdopj = gdop； +q> gdop我们对附表1中给出的9颗卫星按照下面的步骤进行优选，流程

34、图如图(7)所示。stepl:从n颗星中任选3颗星，共有c市种卫星组合，计算每种组合的gdop ,及 其平均值。step2:颗数加1后并判断是否大于n,如果是则执行第三步，反z,重复第一步。step3:分析获得的gdop,兼顾精度和效率原则，确定卫星的数目step4:从n颗星中任选n颗星的用种卫星组合中，找出其中对应gdop值最小的卫 星组合，作为最终的优选方案。设定初始选星个数n=3遍丿丿j种组合求gdop n=n+l根据gdop数据确圮选星个数n遍历种组合取最小gdop対应的卩星分布图（7）选星的流程图选择数目不同卫星对应的gdop的如表3所示,表（3）卫星数目和对应的gdop卫星数量平均

35、gdop卫星数量平均gdop34.5213x10970.1743x10647.1169x10680.1346x10650.7139x10690.1106x10660.2632x106从表屮可以看出随着卫星数量的增加，gdop的值变小，即精度越高，同时相邻卫 星数量的gdop的差值变化越来越小用图展示如图（8）所示：用图表展示如下：6o1x i84567定位卫星数量图(8)定位卫星数量对精度的影响从上图屮可以看岀，当定位t星的数量达到6颗以后，再增加t星数量，对gdop 的改善已经不明显。根据不同应用场合对定位精度和定位效率的不同要求，我们结合仿 真结果，提出下列的三种卫星数目选择方案。表(4)

36、不同需求下的选星方案注重精度兼顾精度与效率注重效率卫星数量(颗)765下面分析其中兼顾精度与效率的选星方案，即选星数目确定为6颗。从9颗卫星选择6颗卫星，共有c：种选择方案，每一种方案对应的gdop如图(9)所 示：5ox i4doooq |0102030405060708090选星方案图(9)定位卫星几何分布对精度的影响我们选择其中对应gdop最小的一种选择方案，也就是选择图屮第12种方案： 选择第 1,2, 3, 5, 6, 8 颗卫星，经度分别为 e76°, £89°, £110°, £130°, £136&

37、#176;, e163。, 该组合对应的gdop为1.1475x10。4.3. 4角度误差影响的仿真对角度的测量值随机加入高斯噪声，误差限为0. 1。用加入噪声后的角度值进行1000 次仿真定位实验，来观察角度误差对定位的影响。仿真实验的结果如下图所示：o70004000-3000、一z 2000、一1000 八6660064006200-5600yx6800-5200-5400-50004800图(10)测量数据加入噪声后定位位置如图(10)所示，红色的位置表示未加入噪声时的定位位置，而蓝色点为加入噪声 后的定位位置，由于角度测量值声加入了噪声，定位的位置也围绕在红色点附近波动。 计算100

38、0次仿真实验的定位位置与未加噪声时的定位位置z间的距离谋差，对统计这 些距离误差，得到下图：25020015010050-5000500100015002000250030003500定位距离误差(km)图(11)距离误差统计上图为经过统计后的定位距离的谋差的分布图，统计结杲显示，有95.9%的数据谋 差落在了 600km以内。因此，认为定位精度为600kmo5. 模型结果的分析与检验5.1问题一的结果分析(1) 由于在实际测量时，测向阵列和地球同步卫星夹角存在测量误差，产生误差 的原因可分为两部分，英一是由于系统设备的精度和测角方法而产生，称为系统误差叭 其二是由于测向阵列与地球同步卫星之间

39、的环境干扰而产生，称为环境误差。因此模型 一屮对84个位置的进行删除误差位置、确定位置初始区域，再用遗传算法在初始区域 寻找位置的最优解，这样处理是合理的。(2) 模型一求岀的飞行器距离地球表而2553.6千米，模型二求岀的飞行器距离地 球表面2600.5千米。5. 2问题二的结果分析问题二，利用5个时刻的测量数据，联立方程组，用最小二乘法使得角度与测量值 之间的误差最小，从而逼近飞行器的真实位置。得到的位置矢量为 (-41828.8,4627.0,6190.6),速度矢量(4.0,2.9,0.2),计算结果显示e行器大约距离地面 2421km,应该是一颗空间飞行器，其速度为4. 9447km

40、/s,与空间飞行器的飞行速度也较 吻合。在进行可靠性分析时，z轴坐标的置信度为95%的置信区间的长度较x轴，y轴坐 标较宽，说明z轴坐标的可靠性较低。其原因可能是因为定位用的同步卫星在坐标系小 的z坐标均为0,影响了对飞行器的z轴坐标的测算。5. 3问题三的结果分析问题三，利用几何精度因子gdop来衡量卫星几何分布的好坏，遍历所有方案后， 发现6星组合定位的平均精度和效率能够取得一个较好的折中。其定位结果与9星定位 的结果很接近，说明6星定位的组合是较为合理的。在给测量值加入最人为0.1。的噪 声后，大量的仿真实验显示，加入的噪声给定位带600km的误差。说明模型建立的 定位方法对于角度测量冇

41、着很高的要求，否则就冇可能造成较大的误差。即说明本模型 的对测量谋差很灵敏，仍需改进。6. 模型的优缺点分析与改进方向61优点：(1) 在解决问题一时，模型一：通过建立非线性方程组，三颗卫星可以求出一个飞行 器位置参数，而9颗卫星存在84种组合，对于84组数据，首先进行误差分析， 去除一些误差较大的位置，再对剩下的数据使用遗传算法求出飞行器最终的位 置。正由于去除了一些误差较大的数据，在一定程度上提高了精确度；模型二： 充分考虑9颗卫星的数据，首先假设飞行器的位置参数，建立21个方程，然后 利用飞行器位置反求与9颗卫星的18个夹角，再和给岀的测量值相比较，通过 最小二乘法，把误井平方和最小时的

42、位置作为最优解。(2) 在解决问题二时，假设飞行器匀速直线运动，本模型不是单独考虑每一时刻的测 量数据，而是通过速度矢量将5个吋刻的测量数据联系起来，利用最小二乘法， 使得与5个测量值的误差达到最小，从而提高了定位精度；为了验证其数据口j靠 性，我们采用定量的分析方法，把角度进行波动，综合考虑角度测量误差和方程 求解时的偏差带來的影响，经过100次仿真实验，对仿真结果进行统计，给出了 置信水平为0.95的置信区间，定量的描述了预测值的可靠性。(3) 在解决问题三时，考虑卫星优选策略时使用了精度几何因子gdop的方法，兼顾 定位的精度和效率，选择了合适的卫星个数和卫星分布，这种方法精度比较高。

43、当测量角度发生偏差时，假设偏差符合高斯分布，通过1000次的仿真，大量的 数据统计使结果更加符合实际。62缺点：(1) 在问题一的方法一中，三颗定出卫星位置，谋差比较大，即使对84组数据进行 优化处理，仍然是建立在谋差比较大的基础z上的。(2) 在问题二中，最小二乘拟合出的直线给出70s的预测值是理论值，并没有考虑环 境噪声和测量误茅带來的影响。(3) 在问题三中，我们仅仅只考虑了卫星个数带来的运算量对定位效率的影响，没 有没有考虑其他因素的影响，这是不符合实际情况的。6.3改进方向:（1）对于问题一，对丁某吋刻计算出的目标位置，受测量条件的影响，按照本文提出 的定位方法解算出的口标位置与实际

44、位置点总是会存在一定的偏差，因此我们可 以把定位结果带入卡尔曼滤波方程，经过滤波处理，提高其定位精度。（2）对于问题二，假设飞行器匀速圆周运动，利用已知数据建立非线性方程，另外还 需要满足方程2 + /4-r=（/?+/7）2,使用多项式最小二乘曲线拟合，确定70s时 刻的飞行器位置参数。（3）对于问题三，当飞行器运动时，对于实时定位要求就比较高，就要求定位时间尽 量少，而选星时间与选星算法有关，因此我们在选星颗数不变的前提下，可以优 化算法，采用快速计算儿何精度因子的改进算法，该算法基于矩阵的qr分解, 避免了传统算法在计算量和存储量较大导致定位解算实时性降低的问题，在计算 gdop值吋，数

45、值稳定性更好，让算复杂度更低，计算效率更高。当分析定位效 率时，定位效率需要考虑选星时间和方程求解时间两方面因素，可以采用定量的 分析这两方面因素对效率的影响，为提高定位时间提供参考价值。参考文献1 丛丽，谈展屮.提高卫星导航定位精度和实时性的选星算法j,系统工程与电子技 术，2008, (6)2 程铭东，刘利姣，黄光明.基于遗传算法的多传感器网络中目标定位算法j中国 科技论文在线3 沈文亮，李艳斌，陈卫东等基于无源测距的快速定位方法研究j电子学报，2009, (6).4 陈小平，腾云龙，康荣雷等儿何精度因子改进算法研究j.电子科技大学学报， 2008, (6)附录：附录1.第一问遗传算法的求

46、解%遗传算法erro=inf；chu_f =；chu_x=；chu_y=;chu_z=；for i=l:2for j =1：2for k=l：10options=gaoptimset('generations 1,800, 'stallgenlimit',600,gaplotbestf);x,f=ga(thefirstyouhua,9, , , , , -5400+(i-l)*200；6300+(j-l)*350； (k-1)*400;-l;-l;-l;-l;-l;-l , -5000;7000；4000;1;1;1;1;1;1, zoptions);chu_x=chu

47、_x；x(1)；chu_y=chu_y；x(2)；chu_z=chu_z；x(3)；chu_f=chu_f；f； if (erro>f)x0=x；erro=f；endendendend附录2第一问遗传算法的求解的主函数%遗传算法函数function f=thefirstyouhua (x) %viuo»ie %acaiadqe«oaviaer=6367;h=35800;thetal=142/360*2*pi;theta2=1630/360*2*pi;theta3=172.0/360*2*pi；theta4=136.0/360*2*pi；theta5=1300/360*

48、2*pi；theta6=125.0/360*2*pi；theta7=110.0/360*2*pi;theta8=89-0/360*2*pi；theta9=76.0/360*2*pi；al=100.32/360*2*pi;bl=ll.58/360*2*pi; a2=125.58/360*2*pi;b2=36.04/360*2*pi; a3=135.93/360*2*pi;b3=4634/360*2*pi; a4=92.91/360*2*pi;b4=5.99/360*2*pi; a5=85.46/360*2*pi;b5=6.95/360*2*pi; a6=79.24/360*2*pi;b6=ll.

49、99/360*2*pi; a7=60.73/360*2*pi;b7=29.81/360*2*pi; a8=3582/360*2*pi;b8=5463/360*2*pi; a9=21.33/360*2*pi；b9=69.25/360*2*pi； xm8=(r+h)*cos(theta8)；ym8=(r+h)*sin(theta8)；zm8=0;xml=(r+h)*cos(thetal)；yml xm2=(r+h)*cos(theta2)；ym2 xm3=(r+h)*cos(theta3)；ym3 xm4=(r+h)*cos(theta4)；ym4 xm5=(r+h)*cos(theta5)；ym

50、5 xm6=(r+h)*cos(theta6)；ym6 xm7=(r+h)*cos(theta7)；ym7(r+h)*sin(thetal)；zml=0 (r+h)*sin(theta2)；zm2=0 (r+h)*sin(theta3)；zm3=0 (r+h)*sin(theta4)；zm4=0 (r+h)*sin(theta5);zm5=0 (r+h)*sin(theta6)；zm6=0 (r+h)*sin(theta7)；zm7=0xm9=(r+h)*cos(theta9)；ym9=(r+h)*sin(theta9)；zm9=0；if(x(l)a2+x(2)a2+x(3)a2<636

51、7a2)f=inf;elsef=(xml-x(1)*x(4)+(yml-x(2)*x(5)+(-x(3)*x(6)-(xml-x(1)a2+(yml-x(2)a2+ (zml-x(3)x2)x0.5*cos(al) ) x2+.(xml-x(1)*x(7)+(yml-x(2)*x(8)+(-x(3)*x(9)-(xml-x(1)x2+(yml-x(2)x2+(z ml-x(3)a2)a0.5*cos(bl)x2+.(xm2-x(1)*x(4)+(ym2-x(2)*x(5)+(-x(3)*x(6)-(xm2-x(1)a2+(ym2-x(2)a2+(z m2-x(3)x2)x0.5*cos(a2)

52、x2+.(xm2-x(1)*x(7)+(ym2-x(2)*x(8)+(-x(3)*x(9)-(xm2-x(1)x2+(ym2-x(2)a2+(z m2-x(3)a2)a0.5*cos(b2)a2+.(xm3-x(1)*x(4)+(ym3-x(2)*x(5)+(-x(3)*x(6)-(xm3-x(1)a2+(ym3-x(2)a2+(z m3-x(3)a2)a0.5*cos(a3)a2+.(xm3-x(1)*x(7)+(ym3-x(2)*x(8)+(-x(3)*x(9)-(xm3-x(1)x2+(ym3-x(2)x2+(z m3-x(3)x2)x0.5*cos(b3)x2+.(xm5-x(1)*x

53、(4)+(ym5-x(2)*x(5)+(-x(3)*x(6)-(xm5-x(1)x2+(ym5-x(2)x2+(z m5-x(3)a2)a0.5*cos(a5)a2+.(xm5-x(1)*x(7)+(ym5-x(2)*x(8)+(-x(3)*x(9)-(xm5-x(1)x2+(ym5-x(2)x2+(z m5-x(3)a2)x0.5*cos(b5)x2+.(xm7-x(1)*x(4)+(ym7-x(2)*x(5)+(-x(3)*x(6)-(xm7-x(1)x2+(ym7-x(2)x2+(z m7-x(3)a2)a0.5*cos(a7)a2+.(xm7-x(l)*x(7)+(ym7-x(2)*x

54、(8)+(-x(3)*x(9)-(xm7-x(l)x2+(ym7-x(2)x2+(z m7-x(3)a2)x0.5*cos(b7)x2+.(xm8-x(1)*x(4)+(ym8-x(2)*x(5)+(-x(3)*x(6)-(xm8-x(1)人2+(ym8-x(2)2+(z m8-x(3)a2)a0.5*cos(a8)a2+.(xm8-x(1)*x(7)+(ym8-x(2)*x(8)+(-x(3)*x(9)-(xm8-x(1)人2+(ym8-x(2)2+(z m8-x(3)a2)a0.5*cos(b8) ) a2+.(x(4)*x(7)+x(5)*x(8)+x(6)*x(9)x2+.(x(4)x

55、2+x(5)x2+x(6)x2-l)x2+.(x(7)a2+x(8)a2+x(9)a2-1)a2;endend附录3最小二乘法clc；clear；%令6xx=-6000;-3000;-1000;yy=2000;4000; 7000;zz= 1000; 3000; 5000;vx=-30;0;30;vy= -30;0;30;vz= -30;0;30;xxx=;x0= -6000;2000;1000;l;0;l;2;l;6;30;30;30;mm=10；nn=x0；for i=l:3for j =1：3for k=l:3for ii=l:3for j j =1：3for kk=l:3aal=xx(i) ；bbl=yy(j) ；ccl = zz(k)；vxl=vx(ii) ；vyl=vy(j j) ；vzl=vz(kk)；x0= aal ；bbl ； ccl ；l；0；l；2；l；6 ； vxl ； vyl ； vzl ； %onaoxid options = optimset(optimset,'maxfunevals1000