用MATLAB计算GPS卫星位置.docx

用MATLAB计算GPS卫星位置.docx

《用MATLAB计算GPS卫星位置.docx》由会员分享，可在线阅读，更多相关《用MATLAB计算GPS卫星位置.docx（11页珍藏版）》请在冰豆网上搜索。

用MATLAB计算GPS卫星位置

用MATLAB计算GPS卫星位置

作者：

罗利娟杨乐

来源：

《科技资讯》2013年第21期

罗利娟杨乐

（西安翻译学院陕西西安710061）

摘要:

本文主要介绍了GPS测量数据的常用格式RINEX标准文件格式,并利用MATLAB工具计算出所观测卫星里的五颗卫星（14、20、29、31和32五颗）在283个历元的瞬时位置,即所观测时间段里五颗卫星在WGS-84坐标下的空间运行轨迹。

关键词:

RINEX标准文件WGS-84下卫星位置MATLAB工具

中图分类号:

P228.4文献标识码:

A文章编号:

1672-3791（2013）07（c）-0005-05

GPS定位的基本原理简单来说就是在WGS-84空间直角坐标系中,确定未知点与GPS卫星的空间几何关系。

因此利用GPS进行导航和测量时,卫星是作为位置已知的高空观测目标。

那么如何精确快速的解算出卫星在空间运行的轨迹即其轨道是实现未知点快速定位的关键。

1标准格式RINEX格式简述

在进行GPS数据处理时,由于接收机出自于不同厂家,所以厂家设计的数据格式也是五花八门的,但是在实际中,很多时候需要把来自不同型号的接收机的数据放在一块进行处理,这就需要数据格式的统一,为了解决这种矛盾,RINEX（英文全称为:

TheReceiverIndependentExchangeFormat）格式则应运而生,该格式存储数据的类型是文本文件,数据记录格式是独立于接收机的出自厂家和具体型号的。

由此可见,其特点是:

由于是通用格式,所以可将不同型号接收机收集的数据进行统一处理,并且大多数大型数据处理软件都能够识别处理,此外也适用于多种型号的接收机联合作业,通用性很强。

RINEX标准文件里不是单一的一个文件,而是包括如下几种类型的文件[1]。

（1）观测数据文件（ssssdddf.yyo）,记录的是GPS观测值信息,（OBServationdata,简写OBS,为接收机记录的伪距、相位观测值;O文件,如XG012191.10O）。

（2）导航电文文件（ssssdddf.yyn）,记录的是GPS卫星星历信息（NAVavigationdata,简写NAV,记录实时发布的广播星历;N文件,如XG012191.10N）。

（3）气象数据文件（ssssdddf.yym）,主要是在测站处所测定的气象数据（METerologicaldata,简写MET,记录气象仪器观测的温、压、湿度状况;M文件,如XG012191.10M）。

（4）GLONASS导航电文文件（ssssdddf.yyg）,记录的是地球同步卫星的导航电文。

由上述可见,RINEX文件的命名规则为ssssdddf.yyt（t指的是数据类型,不同的文件,t所代表的字母不同）,其中文件名前四个字母（ssss）指的是测站名,一般是用字母和数字的组合来定义,方便识别,用户自己定义。

紧跟着的三个字母（ddd）指的是第一组数据的年积日（年积日是仅在一年中使用的连续计算日期的方法,是从当年1月1日起开始计算的天数。

例如:

每年的1月1日为第1日,2月1日为第32日,以此类推）如219表示8月7号,年积日的计算可通过在网上下载软件进行快速计算,也可自己编一个小程序计算。

f指的是观测当天文件的观测序列号,其可以在0～9或A～Z中取值,如果文件序号取值为0,则其意味着今天一天观测的所有数据都放在该文件中,用户则不用在所有文件中一一去找当天的某个文件,只需在此文件中进行查询即可。

对于f的理解,还有一点需要注意,通过下面的例子来进行说明:

在某一天,用GPS测量进行某一项目的时候,共使用了3台GPS接收机,由于各种原因,用户分为三个时间段来进行:

第一时间段,三台接收机全部启用；第二时间段,启用两台接收机；第三时间段,三台接收机又全部工作,那么在第一时间段,三台接收机所接收到的数据文件编号则都为1,在第二时间段,参与工作的两台接收机接收到的观测数据文件编号为2,在第三时间段,又是所有三台接收机参与观测,那么在该时间段,这三台接收机所对应的文件编号就为3。

由此可见,文件序号的编排不是以某台接收机在该天的观测时段为基础定义,而是用用户所进行的整体项目在该天的同步观测时段为基础的。

RINEX文件命名里的yy指的是年的后两位数字,如对于2013年,则yy为13。

RINEX文件命名里后一位o,n,m,g则标识的分别是观测数据文件、导航电文文件、气象数据文件以及GLONASS导航电文文件。

此标准文件是纯ASCll码文本文件,所有类型的文件都由两部分组成,文件头和数据记录。

文件头和数据记录的区分界是ENDOFHEADER,文件头里每条记录占一行,列宽不超过80,并且第61～80列是对前面列数据记录内容的说明,称为标签。

数据记录里如果一行写不下,则在第二行继续数据的记录。

1.1GPS观测数据文件

GPS观测数据文件中储存的数据是与确定GPS整周未知数息息相关的信息,主要包括星历、卫星数、测站概略坐标、伪距信息、载波相位观测量信息等。

下面分别通过表1和表2对GPS观测数据文件的文件头和数据记录里每个记录的含义进行详细说明,以便更清楚的看懂观测数据文件。

GPS观测数据文件的文件头里每个标签的含义如表1所示。

GPS观测数据文件里数据记录里每一行的含义如表2所示。

本文中所使用的算例里观测数据文件XG032191.10O的头文件和部分数据记录如图1与图2所示。

1.2导航数据文件

GPS导航文件也是计算卫星瞬时位置必不可少的文件,因为卫星星历、时钟改正、偏心率、轨道摄动改正、大气折射改正等导航信息是用户用来实时定位和精确导航的必备数据。

导航文件仍然包括头文件和数据记录两部分,数据文件中数据之间是用空格或回车符隔开的,因为有这样的规律,所以用MATLAB编程读取数据的时候就可利用这点进行数据的获取的控制。

GPS导航数据文件的头文件不像观测数据文件那样信息丰富,其仅仅只是对RINEX版本号、观测类型、文件纲要名称、文件机构名称、文件建立日期相关信息进行说明,同观测文件一样,前60列是数据,61～80是标签说明,并以ENDOFHEADER作为文件头的结束标志。

GPS导航数据文件里的数据记录部分相对较复杂,牵扯到的内容较多,其每一行每一项的含义如表3所示。

（如图3,4）

2卫星坐标的计算步骤

由于在GPS定位和导航的时候,用户都是把GPS卫星的位置作为已知量来对待,并且GPS定位所用的坐标系是世界大地坐标系WGS-84。

所以就先必须根据GPS接收机观测的相应星历数据,解算出GPS卫星在WGS-84坐标系中的瞬时位置。

为了后面计算方便,先对广播星历中涉及到的计算卫星坐标的一些轨道参数进行说明,如表4所示。

由于每隔两个小时,GPS接收机收到的广播星历才更新一次,所以用户在根据接收机收到的卫星导航电文汇总的广播星历参数推算GPS的瞬时坐标的时候,一定要选取与GPS卫星的瞬时坐标时刻最相近的那组广播星历数据[2],否则误差将会很大。

首先由已知的GPS接收机接受时刻的钟面时,根据公式换算出GPS卫星发射时刻的钟面时[3]:

其中：

可由观测文件直接读取,为伪距,是观测GPS卫星与观测GPS接收机之间的距离（由于这个距离里含有电离层误差、对流层误差等各种误差,不是GPS卫星与GPS接收机之间的真实距离,所以称为伪距）,本次算例读取的是P2码伪距观测量,c为光速,值为3×108m/s。

下面就是按照公式,读取相应的广播星历参数,计算观测时刻的GPS卫星瞬时坐标:

（1）GPS卫星在空中运行的平均角速度n:

式中:

GM为地球引力常数,;为卫星轨道长半轴的平方根,为摄动改正数,这两个数值都可在导航电文中直接读取。

（2）归化时刻。

对观测时刻进行卫星的钟差改正:

为求得观测时刻所对应的轨道参数,则需从导航电文中给出的对应于参考时刻的GPS卫星轨道参数推算而得。

两个历元的时间差（归化时间）计算如下：

其中：

、a0、a1、a2可直接从导航文件中读取,为已知量。

GPS卫星发射时刻的GPS标准时计算方法如下:

由具体事例可知,RINEX标准文件里记录的观测时刻t是用年、月、日、时、分、秒表示的,然而GPS广播星历中的参考时刻t是用GPS秒表示的,因此我们需要一个转化过程,把时间单位进行统一,于是把民用日中的年月日时分秒根据相应的公式换算为GPS秒,首先将观测时刻里的时分秒化为实数时,即：

然后将民用日的Y、M、D、UT化为儒略日（不记年月,只记日的历法）,即：

式中：

INT表示的是对括号里的值进行取整运算。

对于y、m的取值遵循下面的规则:

如果月的值不大于2,则y=Y-l,m=M+12;否则y=Y,m=M。

最后把JD的值代入下列公式计算GPS秒:

GPS周=INT（（JD-2444244.5）/7）（7）

GPS秒=（JD-2444244.5-GPS周×7）×24×3600

（8）

此外对于计算出来的还应做如下处理:

s时,=-604800s;当s时,=+604800s。

（3）计算观测时刻卫星平近点角。

式中：

为参考时刻的平近点角,在导航电文中直接读取;n为GPS卫星运行的平均角速度,在第一步已算出,在第二步已算出。

（4）计算偏近点角:

式中各个字母的含义如下:

为卫星轨道偏心率,在导航电文中直接读取;在第三步中已算出,和单位均为弧度,从该式中中可以看出,计算偏近点角必须采用迭代法。

偏近点角的初始值可取,代入上式,计算出偏近点角,再代入上式,计算出偏近点角,依次一直迭代,直至两次迭代结果之差时停止迭代,取该值为最终的偏近点角。

（5）相对论效应:

当s时,=-604800s;当s时,=+604800s。

由于相对论效应的影响,规划时间发生变化,于是需要把代到式（9）和式（10）中重新计算后续过程需要的和的值。

（6）真近点角:

在用MATLAB计算真近点角时应注意:

式中的arctan函数应采用atan2函数。

（7）升交距角:

其中:

已在上一步计算出来了,为近地点角距,可在导航电文中直接读取。

（8）计算升交距角的摄动量、卫星失径r的摄动量、轨道倾角i的摄动量:

其中:

均可从导航电文直接读取,为已知量。

（9）计算经过摄动改正的升交距角、卫星的地心距离及轨道倾角:

式中:

为参考时刻的轨道倾角,为轨道倾角变化率,均可从导航电文直接读取,为已知量。

（10）计算卫星在轨道坐标系中的位置:

（11）计算卫星在世界大地坐标系中的坐标修正后的升交点经度:

其中:

为升交点赤经变化率;为参考时刻的升交点赤经;、均可从导航电文直接读取;为已知量,为地球自转速率,值是一个定数,。

（12）最后将卫星在轨道坐标系的坐标经坐标转换,换算出卫星在WGS-84坐标系的瞬时位置:

3卫星坐标的计算实例

为了验证本文的计算卫星坐标的理论,本人经过实地设置观测点XG01和XG03,从8月7日1:

37:

00开始观测,8月7日2:

47:

30结束观测过程。

数据格式为RINEX标准格式,无需转换。

观测数据文件XG012191.10O里观测数据的读取、观测数据按照定义的卫星名存放,XG012191.10N里定义卫星名的相关的导航数据的读取和存放,GPS卫星瞬时位置的计算,以上的这些过程都在MATLAB环境中实现,图5就是在MATLAB环境里计算、绘制出的GPS卫星名为14、20、29、31和32在观测时间段的空间运行轨迹图。

4结论

快速准确计算出GPS卫星在WGS-84坐标系下的瞬时位置是GPS定位里很基础、很重要的问题,由文中描述可知,计算步骤比较繁琐,需要注意的细节也很多,牵扯到的符号也很多。

从自己对实例的计算经历中来看,我觉得在使用卫星计算公式计算卫星瞬时坐标的时候特别要注意以下三个问题。

（1）由于每隔两个小时,GPS接收机收到的广播星历才更新一次,所以用户在根据接收机收到的卫星导航电文汇总的广播星历参数推算GPS的瞬时坐标的时候,一定要选取与GPS卫星的瞬时坐标时刻最相近的那组广播星历数据,否则误差将会很大。

（2）由于GPS定位系统的高空动态性,又因为卫星的发射时刻可唯一表征该时刻GPS卫星的坐标,因此应先由接收机钟面时计算出公式中所需的GPS卫星发射时刻的钟面时。

（3）必须把GPS接收机的观测历元换算成GPS秒,方可进行下一步。

参考文献

[1]张妮,王标标.基于Matlab读取标准RINEX格式的GPS星历数据[J].电子设计工程,2010,18（8）:

23-25.

[2]王猛,张志伟.利用广播星历计算卫星的瞬时坐标[J].城市勘测,2010,2:

89-90.

[3]刘勤志,尹长林,易重海.计算GPS卫星发射时刻的两种方法[J].长沙电力学院学报,2004,19

（1）:

65.

[4]RobertS.Radovanovic.AdjustmentofSatellite-BasedRangingObservati-onsforprecisePositionandDeformationMonitoring,Ph.D.Dissertation,Depart-mentofGeomaticsEnglneering[J].UniversityofCa1Gary,2002.