RTKPPP定位算法流程.docx
《RTKPPP定位算法流程.docx》由会员分享，可在线阅读，更多相关《RTKPPP定位算法流程.docx（40页珍藏版）》请在冰豆网上搜索。
RTKPPP定位算法流程
1基础知识
2
2.1GPS精密单点定位的基本原理
2.2
GPS精密单点定位一般采用单台双频GPS接收机，利用IGS提供的精密星历和卫星钟差，基于载波相位观测值进行的高精度定位。
观测值中的电离层延迟误差通过双频信号组合消除，对流层延迟误差通过引入未知参数进行估计。
2.3时间系统
2.4
RTKLIB内部使用GPST（GPST时间）用于GNSS的数据处理和定位算法。
数据在RTKLIB内部处理之前，需要转换成GPST时间。
使用GPST的原因是避免处理润秒。
RTKLIB使用以下结构体表示时间：
typedefstruct
time_ttime;/*time（s）expressedbystandardtime_t*/
doublesec;/*fractionofsecondunder1s*/
}gtime_t;
2.4.1GPST和UTC（UniversalTimeCoordinated）
2.4.2
关系参考【图1】，参考【图2】：
图1转换关系公式
图2
通过使用GPS导航信息中的UTC参数，GPST到UTC或者UTC到GPST之前的转换可以用更准确的表达方式，如【图3】。
图3
这些参数是由GPS导航消息提供的。
2.4.3BDT（北斗导航卫星系统时间）
2.4.4
BDT（北斗导航卫星系统时间）是一个连续的时间系统，没有润秒。
开始历元的时间是【UTC2006年1月1号00:
00:
00】。
北斗时间计算公式【图4】：
图4
UTC和GPST时间转换同上面的GPS一样，只不过UTC参数来自与北斗导航信息中。
2.5坐标系统
2.6
接收机和卫星的位置在RTKLIB中表示为在ECEF（地心地固坐标系）坐标系统中的X,Y,Z组件。
2.6.1大地坐标到ECEF坐标的转换
2.6.2
转换公式如【图5】。
第三个公式最后一行有错，应该为:
（v（1–e2）+h）sin
图5
参数说明：
a:
地球参考椭球的长半径
f:
地球参考椭球的扁平率
h:
椭球高度
:
纬度
:
经度
当前版本的RTKLIB使用的值为【图6】：
图6
图7参考椭球体
2.6.3ECEF坐系到大地坐标的转换
2.6.4
转换公式如【图8】
图8
2.6.5本地坐标到ECEF坐标的转换
2.6.6
在接收机位置的本地坐标，也被称为ENU坐标，通常使用在GNSS导航处理。
ECEF坐标到本地坐标转换的旋转矩阵表示为【图9】。
图9Er旋转矩阵
参数说明：
：
接收机位置的纬度
：
接收机位置的经度
通过使用Er和接收机的坐标rr【ECEF】，坐标recef【ECEF】可以被转换到本地坐标的坐标rlocal，公式如【图10】。
图10
3RTKPPP定位算法
4
4.1单点定位（pntpos）
4.2
1:
satposs2:
estpos3:
estvel
1.计算计算卫星位置、速度和时钟（satposs）
2.
a）通过广播星历计算卫星钟差（ephclk）
b）
卫星编号到卫星系统的转换（satsys）
根据卫星的编号，获取到对应的卫星导航系统。
选择星历（seleph）
1、传入信号传输时间，卫星编号，导航数据等参数。
2、遍历导航数据，遍历导航数据里面的星历数据，判断星历数据的卫星编号是否和传入的卫星编号相等。
3、如果星历数据的卫星编号和传入的卫星编号相等，就计算星历参考时间（toe）和信号传输时间的时间差。
如果不相等，继续处理下一条星历数据。
4、判断计算出来的时间差，如果时间差大于了允许的最大时间差，继续查找下一个星历数据。
否则，判断时间差最小的星历数据，记录星历数据的位置。
5、返回之前记录出来的星历数据。
使用广播星历计算卫星时钟偏差（eph2clk）
1、传入信号发射时刻的时间和星历数据。
2、
3、计算信号发射时刻的时间和本时段钟差参数参考时间（星历参数toc）的时间差。
4、
5、通过下式计算钟差，这里还没有处理相对论校正项和tgd：
（代码中有个迭代过程，资料上没看到写）
6、
卫星钟差计算出来之后，信号发射时刻的时间还要加上这个钟差。
c）计算卫星在信号发射时刻的位置、速度和时钟（satpos）
d）
根据星历表选项来选择不同的处理，如下：
广播星历（EPHOPT_BRDC）:
ephpos（广播星历到卫星位置和钟差）
1、根据公式计算出tk；
2、
3、根据使用的卫星系统，选择使用的地球引力常数（mu）和地球的角速度（omge）
4、
5、根据公式计算出平近点角M。
6、
7、求解开普勒方程，按照以下公式迭代求解。
8、
9、根据以下公式计算出u（改正后的纬度幅角）,r（改正后的径向）,i（改正后的轨道倾角）的值。
10、
11、根据以下公式计算卫星在轨道平面内的坐标。
12、
13、根据不同的卫星系统，做不同的计算。
14、
GPS计算方式：
北斗计算方式：
其中
15、按照公式计算出时间tc。
16、
17、按照以下公式计算出钟差和钟漂。
18、
精密星历（EPHOPT_PREC）:
peph2pos
广播+SBAS（EPHOPT_SBAS）：
satpos_sbas
广播+SSR_APC（EPHOPT_SSRAPC）:
satpos_ssr
广播+SSR_COM（EPHOPT_SSRCOM）：
satpos_ssr
QZSSLEX星历（EPHOPT_LEX）：
lexeph2pos
3.使用伪距估算接收机的位置,返回估算状态结果（estpos）
4.
a）伪距残差（rescode）
b）
把ecef坐标系转换成大地坐标系（ecef2pos）
1、按照以下公式做转换，暂时还没看懂。
2、
计算几何距离和接收机到卫星的单位矢量（geodist）
1、用卫星的坐标向量做欧几里德范数，返回值和地球长半轴（WGS84）比较。
小于地球长半轴（WGS84），返回-1；
2、
3、计算卫星坐标和接收机坐标的差值向量。
4、
5、用差值向量做欧几里德范数，再用上一步计算出来的差值向量和计算结果做除法，得到视线向量。
按照如下公式：
6、
7、使用以下公式计算几何距离。
8、
计算卫星方位角/仰角（satazel）
1、把接收机ecef坐标转换到大地坐标；
2、
3、判断高度是否大于地球半长轴（WGS84）的负数值；
4、
5、如果高度小于等于地球半长轴（WGS84）的负数值，方位角为0，仰角为PI/2;
6、
7、如果高度大于地球半长轴（WGS84）的负数值，把ECEF向量转换到局部坐标。
然后对转换出来的坐标做内积。
8、
9、使用如下公式计算卫星方位角和仰角。
10、
伪距使用编码残差改正（prange）
暂时没找到对应的文档对应。
电离层改正（ionocorr）
通过广播电离层模型（klobuchar模型）计算出电离层延迟（ionmodel）
1、校验传入的电离层模型参数，校验失败，使用默认的电离层模型参数；
2、
3、使用以下公式计算出地球为中心的角度（半圆）；
4、
5、使用以下公式计算子的电离层的纬度/经度（半圆）；
6、
7、使用以下公式计算地磁纬度。
8、
9、计算本地时间，返回值按这个公式【tt-=floor（tt/86400.0）*86400.0;】处理，保证tt的范围（0<=tt<86400）。
10、
11、计算倾斜因子。
12、
13、电离层延迟计算，公式如下。
14、
对流层改正（tropcorr）
通过标准大气压和saastamoinen模型计算对流层延迟（tropmodel）
1、使用以下公式计算总气压。
2、
3、使用以下公式计算绝对温度。
4、
5、使用以下公式计算水蒸汽的分压。
代码中没有使用。
6、
7、使用以下公式计算【Saastamoinen模型】。
8、
伪距残差
残差值=伪距-（物理距离+dtr-光速*时钟偏差+电离层误差+对流层误差）
时钟系统和接收器的偏置补偿
伪距测量误差方差（varerr）
c）方差权重值（weightbyvariance）
d）
e）最小二乘估计（lsq）
f）
最小二乘估计通过求解正规方程（X=（A*A'）^-1*A*Y）
1、计算矩阵A*Y的结果，保存到矩阵Ay；
2、
3、计算矩阵A*A’的结果，保存到矩阵Q；
4、
5、求的矩阵Q的逆矩阵,结果保存到矩阵Q；
6、
7、最后x=Q*Ay
8、
g）欧几里德范数（norm）
h）
公式：
1、通过最小二乘法估算出来的参数dx【长度为4】，分别加到位置向量上x【长度为4】。
2、
3、对这个速度向量做欧几里德范数，返回值同1E-4比较，小于这个值时，就得到估算出的接收机的位置为向量x的值。
4、
i）验证求解（valsol）
j）
1、对伪距残差值做内积，然后同卡方分布的自由度分布值（alpha=0.001）做比较，大于卡方分布的自由度分布值的数据，无效。
2、
3、计算dops；暂时还没没找到资料。
4、
5、用计算出来的dops和配置的最大dops阈值比较，大于配置的最大dops阈值的数据无效。
6、
5.使用多普勒估算接收机速度（estvel）
6.
a）多普勒残差（resdop）
b）
把ecef坐标系转换成大地坐标系（ecef2pos）
把xyz坐标转换成enu坐标（xyz2enu）
计算ecef中的瞄准线向量
计算相对于接收机在ECEF中的卫星速度
多普勒残差计算
1、速率公式：
c）最小二乘估计（lsq）
d）
最小二乘估计通过求解正规方程（X=（A*A'）^-1*A*Y）
1、计算矩阵A*Y的结果，保存到矩阵Ay；
2、
3、计算矩阵A*A’的结果，保存到矩阵Q；
4、
5、求的矩阵Q的逆矩阵,结果保存到矩阵Q；
6、
7、最后x=Q*Ay
8、
e）欧几里德范数（norm）
f）
公式：
5、通过最小二乘法估算出来的参数dx【长度为4】，分别加到速度向量上x【长度为4】。
6、
7、对这个速度向量做欧几里德范数，返回值同1E-6比较，小于这个值时，就得到估算出的接收机速度为向量x的值。
8、
4.3精确定位（pppos）
4.4
1:
udstate_ppp2:
satposs3:
testeclipse4:
res_ppp5:
res_ppp6:
filter7:
res_ppp
1.暂时更新状态（udstate_ppp）
2.
a）位置更新
b）
c）时钟更新
d）
e）对流层参数更新
f）
g）相位偏差更新
h）
1、通过LLI检测周跳；
2、
3、通过Geometry-Free相位跳变检测周跳（如果双频测量值可用）；
4、
3.计算卫星的位置和速度和时钟（satposs）
4.
e）通过广播星历计算卫星钟差（ephclk）
f）
卫星编号到卫星系统的转换（satsys）
根据卫星的编号，获取到对应的卫星导航系统。
选择星历（seleph）
升级会员 