RTKPPP定位算法流程.pdf

《RTKPPP定位算法流程.pdf》由会员分享，可在线阅读，更多相关《RTKPPP定位算法流程.pdf（38页珍藏版）》请在三一文库上搜索。

1、实用标准文案 文档 1基础知识 1.1GPS 精密单点定位的基本原理 GPS 精密单点定位一般采用单台双频GPS 接收机，利用 IGS提供的精密星历和卫星钟差， 基于载波相位观测值进行的高精度定位。观测值中的电离层延迟误差通过双频信号组合消 除，对流层延迟误差通过引入未知参数进行估计。 1.2时间系统 RTKLIB内部使用 GPST(GPST 时间 ) 用于 GNSS 的数据处理和定位算法。数据在 RTKLIB内部处 理之前， 需要转换成 GPST 时间。 使用 GPST 的原因是避免处理润秒。RTKLIB使用以下结构体表 示时间： typedef struct time_t time; /*

2、 time(s) expressed by standard time_t */ double sec; /* fraction of second under 1 s */ gtime_t; 1.2.1GPST 和 UTC(Universal Time Coordinated) 关系参考【图1】，参考【图 2】： 图1 转换关系公式 实用标准文案 文档 图 2 通过使用 GPS 导航信息中的 UTC 参数， GPST 到UTC 或者 UTC 到GPST 之前的转换可以用更准确 的表达方式，如【图3】。 图 3 这些参数是由 GPS 导航消息提供的。 1.2.2BDT （北斗导航卫星系统时间）

3、 BDT （北斗导航卫星系统时间）是一个连续的时间系统，没有润秒。开始历元的时间是 【UTC 2006年1月1号00:00:00 】。 北斗时间计算公式【图4】： 图 4 UTC 和GPST 时间转换同上面的GPS 一样，只不过 UTC 参数来自与北斗导航信息中。 1.3坐标系统 接收机和卫星的位置在RTKLIB中表示为在 ECEF(地心地固坐标系) 坐标系统中的X, Y, Z 组件。 1.3.1大地坐标到 ECEF 坐标的转换 转换公式如【图5】。 实用标准文案 文档 第三个公式最后一行有错，应该为: (v(1 e2)+h)sin 图 5 参数说明：a :地球参考椭球的长半径 f : 地球参

4、考椭球的扁平率 h: 椭球高度 :纬度 : 经度 当前版本的 RTKLIB使用的值为【图6】： 图 6 图7 参考椭球体 实用标准文案 文档 1.3.2ECEF 坐系到大地坐标的转换 转换公式如【图8】 图 8 1.3.3本地坐标到 ECEF 坐标的转换 在接收机位置的本地坐标，也被称为ENU 坐标，通常使用在GNSS 导航处理。 ECEF 坐标到 本地坐标转换的旋转矩阵表示为【图9】。 图 9 Er旋转矩阵 参数说明：接收机位置的纬度 ：接收机位置的经度 通过使用Er和接收机的坐标rr【ECEF 】，坐标recef【 ECEF 】可以被转换到本地坐标的坐 标rlocal，公式如【图10】。

5、图 10 实用标准文案 文档 2RTKPPP 定位算法 2.1单点定位（ pntpos ） 1:satposs 2: estpos 3: estvel 1.计算计算卫星位置、速度和时钟（satposs ） a)通过广播星历计算卫星钟差（ephclk ） 卫星编号到卫星系统的转换（satsys ） 根据卫星的编号，获取到对应的卫星导航系统。 选择星历（ seleph ） 1、传入信号传输时间，卫星编号，导航数据等参数。 2、遍历导航数据，遍历导航数据里面的星历数据，判断星历数据的卫星 编号是否和传入的卫星编号相等。 3、如果星历数据的卫星编号和传入的卫星编号相等，就计算星历参考时 间(toe)和

6、信号传输时间的时间差。如果不相等，继续处理下一条星历数 据。 4、判断计算出来的时间差，如果时间差大于了允许的最大时间差，继续 查找下一个星历数据。否则， 判断时间差最小的星历数据，记录星历数据 的位置。 5、返回之前记录出来的星历数据。 使用广播星历计算卫星时钟偏差（eph2clk ） 1、 传入信号发射时刻的时间和星历数据。 2、 计算信号发射时刻的时间和本时段钟差参数参考时间( 星历参数toc) 的时间差。 3、 通过下式计算钟差，这里还没有处理相对论校正项和tgd ：（ 代码中 有个迭代过程，资料上没看到写） 卫星钟差计算出来之后，信号发射时刻的时间还要加上这个钟差。 b)计算卫星在信

7、号发射时刻的位置、速度和时钟（satpos ） 根据星历表选项来选择不同的处理，如下： 广播星历（ EPHOPT_BRDC）: ephpos （广播星历到卫星位置和钟差） 1、 根据公式计算出 tk ； 实用标准文案 文档 2、 根据使用的卫星系统，选择使用的地球引力常数(mu) 和地球的角速度 (omge) 3、 根据公式计算出平近点角M 。 4、 求解开普勒方程，按照以下公式迭代求解。 5、 根据以下公式计算出u( 改正后的纬度幅角), r(改正后的径向), i(改 正后的轨道倾角) 的值。 6、 根据以下公式计算卫星在轨道平面内的坐标。 7、 根据不同的卫星系统，做不同的计算。 GPS

8、计算方式： 北斗计算方式： 实用标准文案 文档 其中 8、 按照公式计算出时间 tc 。 9、 按照以下公式计算出钟差和钟漂。 精密星历（ EPHOPT_PREC）: peph2pos 广播 +SBAS （EPHOPT_SBAS） ：satpos_sbas 广播 +SSR_APC （ EPHOPT_SSRAPC）: satpos_ssr 广播 +SSR_COM（ EPHOPT_SSRCOM） ：satpos_ssr QZSS LEX星历（ EPHOPT_LEX） ：lexeph2pos 2.使用伪距估算接收机的位置, 返回估算状态结果（estpos ） a)伪距残差（ rescode ） 把

9、ecef 坐标系转换成大地坐标系（ecef2pos ） 1、 按照以下公式做转换，暂时还没看懂。 计算几何距离和接收机到卫星的单位矢量（geodist ） 1、 用卫星的坐标向量做欧几里德范数，返回值和地球长半轴(WGS84) 比 较。小于地球长半轴(WGS84)，返回 -1 ； 2、 计算卫星坐标和接收机坐标的差值向量。 实用标准文案 文档 3、 用差值向量做欧几里德范数，再用上一步计算出来的差值向量和计算 结果做除法，得到视线向量。按照如下公式： 4、 使用以下公式计算几何距离。 计算卫星方位角/ 仰角（ satazel ） 1、 把接收机 ecef 坐标转换到大地坐标； 2、 判断高度是

10、否大于地球半长轴（WGS84 ）的负数值； 3、 如果高度小于等于地球半长轴（WGS84 ）的负数值，方位角为0，仰角 为PI/2; 4、 如果高度大于地球半长轴（WGS84 ）的负数值，把ECEF 向量转换到局 部坐标。然后对转换出来的坐标做内积。 5、 使用如下公式计算卫星方位角和仰角。 伪距使用编码残差改正（prange ） 暂时没找到对应的文档对应。 电离层改正（ ionocorr） 通 过 广 播 电 离 层 模 型 （ klobuchar模 型 ） 计 算 出 电 离 层 延 迟 （ionmodel ） 1、 校验传入的电离层模型参数，校验失败， 使用默认的电离层模型参数； 2、

11、使用以下公式计算出地球为中心的角度（半圆）； 3、 使用以下公式计算子的电离层的纬度/ 经度（半圆）； 实用标准文案 文档 4、 使用以下公式计算地磁纬度。 5、 计算本地时间，返回值按这个公式 【tt-=floor(tt/86400.0)*86400.0;】处理，保证 tt 的范围 （0=tt86400 ）。 6、 计算倾斜因子。 7、 电离层延迟计算，公式如下。 对流层改正（ tropcorr） 通过标准大气压和saastamoinen模型计算对流层延迟（tropmodel ） 1、 使用以下公式计算总气压。 2、 使用以下公式计算绝对温度。 3、 使用以下公式计算水蒸汽的分压。代码中没有

12、使用。 4、使用以下公式计算【Saastamoinen 模型】。 实用标准文案 文档 伪距残差 残差值 =伪距 -( 物理距离 +dtr- 光速 * 时钟偏差 +电离层误差 +对流层误差 ) 时钟系统和接收器的偏置补偿 伪距测量误差方差（varerr ） b)方差权重值 (weight by variance) c)最小二乘估计（lsq ） 最小二乘估计通过求解正规方程（X =（ A * A） - 1* A * Y） 1、 计算矩阵 A*Y的结果，保存到矩阵Ay； 2、 计算矩阵 A*A的结果，保存到矩阵Q ； 3、 求的矩阵 Q 的逆矩阵 , 结果保存到矩阵Q； 4、 最后 x=Q*Ay d

13、)欧几里德范数（norm） 公式 ： 1、 通过最小二乘法估算出来的参数dx 【长度为 4】 ，分别加到位置向量上x 【长 度为 4】。 2、 对这个速度向量做欧几里德范数，返回值同1E-4比较，小于这个值时， 就得到估算出的接收机的位置为向量x的值。 e)验证求解 (valsol) 1、 对伪距残差值做内积，然后同卡方分布的自由度分布值(alpha=0.001)做 比较，大于卡方分布的自由度分布值的数据，无效。 2、 计算 dops； 暂时还没没找到资料。 3、 用计算出来的dops和配置的最大dops阈值比较， 大于配置的最大dops阈值 的数据无效。 3.使用多普勒估算接收机速度（est

14、vel ） a)多普勒残差（resdop ） 把 ecef 坐标系转换成大地坐标系（ecef2pos ） 把 xyz 坐标转换成enu 坐标（ xyz2enu ） 计算 ecef 中的瞄准线向量 计算相对于接收机在ECEF中的卫星速度 多普勒残差计算 1、速率公式： 实用标准文案 文档 b)最小二乘估计（lsq ） 最小二乘估计通过求解正规方程（X =（ A * A） - 1* A * Y） 1、 计算矩阵 A*Y 的结果，保存到矩阵Ay； 2、 计算矩阵 A*A的结果，保存到矩阵Q ； 3、 求的矩阵 Q的逆矩阵 , 结果保存到矩阵Q ； 4、 最后 x=Q*Ay c)欧几里德范数（norm

15、） 公式 ： 3、 通过最小二乘法估算出来的参数dx 【长度为 4】 ，分别加到速度向量上x 【长 度为 4】。 4、 对这个速度向量做欧几里德范数，返回值同1E-6比较，小于这个值时， 就得到估算出的接收机速度为向量x的值。 2.2精确定位（pppos） 1:udstate_ppp 2: satposs 3: testeclipse 4: res_ppp 5: res_ppp 6: filter 7: res_ppp 1.暂时更新状态（udstate_ppp ） a)位置更新 b)时钟更新 c)对流层参数更新 d)相位偏差更新 1、 通过 LLI 检测周跳； 2、 通过 Geometry-F

16、ree 相位跳变检测周跳( 如果双频测量值可用) ； 2.计算卫星的位置和速度和时钟（satposs ） c)通过广播星历计算卫星钟差（ephclk ） 卫星编号到卫星系统的转换（satsys ） 根据卫星的编号，获取到对应的卫星导航系统。 选择星历（ seleph ） 1、传入信号传输时间，卫星编号，导航数据等参数。 2、遍历导航数据，遍历导航数据里面的星历数据，判断星历数据的卫星 编号是否和传入的卫星编号相等。 3、如果星历数据的卫星编号和传入的卫星编号相等，就计算星历参考时 实用标准文案 文档 间(toe)和信号传输时间的时间差。如果不相等，继续处理下一条星历数 据。 4、判断计算出来的

17、时间差，如果时间差大于了允许的最大时间差，继续 查找下一个星历数据。否则， 判断时间差最小的星历数据，记录星历数据 的位置。 5、返回之前记录出来的星历数据。 使用广播星历计算卫星时钟偏差（eph2clk ） 1、 传入信号传输时间和星历数据。 2、 传入信号发射时刻的时间和星历数据。 3、 计算信号发射时刻的时间和本时段钟差参数参考时间( 星历参数toc) 的时间差。 4、 通过下式计算钟差，这里还没有处理相对论校正项和tgd ：（ 代码中 有个迭代过程，资料上没看到写） d)计算卫星位置、速度和时钟（satpos ） 根据星历表选项来选择不同的处理，如下： 广播星历（ EPHOPT_BRD

18、C）: ephpos （广播星历到卫星位置和钟差） 1、 根据公式计算出 tk ； 2、 根据使用的卫星系统，选择使用的地球引力常数(mu) 和地球的角速度 (omge) 3、 根据公式计算出平近点角M 。 4、 求解开普勒方程，按照以下公式迭代求解。 5、 根据以下公式计算出u( 改正后的纬度幅角), r(改正后的径向), i(改 正后的轨道倾角) 的值。 实用标准文案 文档 6、 根据以下公式计算卫星在轨道平面内的坐标。 7、 根据不同的卫星系统，做不同的计算。 GPS 计算方式： 北斗计算方式： 其中 8、 按照公式计算出时间 tc 。 9、 按照以下公式计算出钟差和钟漂。 实用标准文案

19、 文档 精密星历（ EPHOPT_PREC）: peph2pos 广播 +SBAS （EPHOPT_SBAS） ：satpos_sbas 广播 +SSR_APC （ EPHOPT_SSRAPC）: satpos_ssr 广播 +SSR_COM（ EPHOPT_SSRCOM） ：satpos_ssr QZSS LEX星历（ EPHOPT_LEX） ：lexeph2pos 3.排除对被遮蔽卫星的测量（testeclipse） a)获得太阳和月亮在ecef 中的位置（ sunmoonpos） b)归一化三维矢量（normv3） c)欧几里德范数的向量（norm） d)计算太阳 - 地球 -卫星角（

20、dot ） 4.相位和编码残差（res_ppp ） a)地球潮汐改正（tidedisp） b)计算几何距离（geodist ） 1、 用卫星的坐标向量做欧几里德范数，返回值和地球长半轴(WGS84) 比 较。小于地球长半轴(WGS84)，返回 -1 ； 2、 计算卫星坐标和接收机坐标的差值向量。 3、 用差值向量做欧几里德范数，再用上一步计算出来的差值向量和计算 结果做除法，得到视线向量。 4、 使用以下公式计算几何距离。 备注：这个公式是地球自转改正公式。 地球自转改正产生地球自转改正的机制在于地固系随地球的自转而旋转, 地 固系是非惯性系统, 用地固系中的坐标计算卫星到接收机几何距离时,

21、由于 卫星位置和接收机位置是两个不同时刻的位置矢量, 而且这两个时刻的地固 系相对于惯性系是变化的, 因此要考虑地球自转引起的距离变化。其改正公式 为: 实用标准文案 文档 其中 , w 为地球自转角速度; C 为光速 ; x si和x ri分别代表卫星位置矢量和 测站位置矢量的分量, i= 1, 2 对应 x , y 分量。地球自转对纬度影响甚小, 对经 度影响最大 , 其次是高度。 同时地球的自转影响还跟测站的经纬度以及测站与卫星 的几何关系有关, 对两极的测站 , 影响为零 ; 对赤道上的测站影响最大; 当卫星 在测站子午面内影响为零; 卫星在测站东方时, 影响为负 , 卫星在测站西方时

22、, 影 响为正。 c)卫星方位角 / 仰角计算（ satazel ） 1、 把接收机 ecef 坐标转换到大地坐标； 2、 判断高度是否大于地球半长轴（WGS84 ）的负数值； 3、 如果高度小于等于地球半长轴（WGS84 ）的负数值，方位角为0，仰角 为PI/2; 4、 如果高度大于地球半长轴（WGS84 ）的负数值，把ECEF 向量转换到局 部坐标。然后对转换出来的坐标做内积。 5、 使用如下公式计算卫星方位角和仰角。 d)对流层延迟修正： 通过标准大气压和saastamoinen 模型计算对流层延迟（tropmodel ） 1、 使用以下公式计算总气压。 2、 使用以下公式计算绝对温度。

23、 3、 使用以下公式计算水蒸汽的分压。代码中没有使用。 4、使用以下公式计算【Saastamoinen 模型】。 实用标准文案 文档 e)卫星天线模型（satantpcv ） f)接收机天线模型，通过天线相位中心参数计算天线偏移（antmodel ） g)相位缠绕校正（windupcorr ） h)电离层和天线相位校正测量（corrmeas ） i)计算卫星时钟和对流层延迟 5.卡尔曼滤波（ filter） 卡尔曼滤波按照以下的公式更新状态： K = P * H * (H * P * H + R)-1 ; xp = x + K * v ; Pp = (I K * H) * P ; 参数解释：

24、x: 状态向量 (n x 1) p: 状态的协方差矩阵 (n x n) H: 设计矩阵的转置矩阵 (n x m) v: 创新 (测量模型 ) (m x 1) R: 测量误差的协方差矩阵 (m x m) n,m: 状态和测量值的数量 xp：状态更新后的向量 (n x 1) Pp: 状态更新后的协方差矩阵 (n x n) 6.组合残差（res_ppp ） 7.解析 PPP整周模糊度（ pppamb ） a)average LC （average_LC） b)固定宽巷模糊（fix_amb_WL） 计算 LC波长度（ m ） （lam_LC） 宽巷模糊 宽巷模糊度的方差 c)fix narrow-la

25、ne ambiguity 【AR mode: PPP-AR】 d)fix narrow-lane ambiguity by ILS 【AR mode: PPP-AR ILS 】 实用标准文案 文档 3以下是上面一些算法的参考资料，信息来自于Rtklib的 手册 3.1接收机和卫星天线之间的几何距离 实用标准文案 文档 3.2卫星方向的方位角和仰角 实用标准文案 文档 实用标准文案 文档 3.3GPS 、Galileo和 QNSS 广播星历和时钟 实用标准文案 文档 3.4北斗广播星历和时钟 实用标准文案 文档 3.5对流层和电流层模型 3.5.1对流层模型 实用标准文案 文档 3.5.2电离层

26、模型 实用标准文案 文档 实用标准文案 文档 3.6单点定位 3.6.1线性最小二乘估计 实用标准文案 文档 实用标准文案 文档 实用标准文案 文档 3.6.2估算接收机的位置和钟差 实用标准文案 文档 实用标准文案 文档 3.6.3估算接收机的速度和钟漂 实用标准文案 文档 实用标准文案 文档 3.6.4求解验证和 RAIM FDE 实用标准文案 文档 3.7精密定位 3.7.1用于 PPP的 ZD测量模式 3.7.2接收机天线相位中心模型 实用标准文案 文档 3.7.3卫星天线相位中心模型 实用标准文案 文档 实用标准文案 文档 3.7.4由地球潮汐引起的位移 3.7.5相位缠绕校正 实用标准文案 文档 3.7.6通过 PPP估算接收机的位置 实用标准文案 文档