RTKPPP定位算法流程

《RTKPPP定位算法流程》由会员分享，可在线阅读，更多相关《RTKPPP定位算法流程（38页珍藏版）》请在装配图网上搜索。

1、1 基础知识1.1 GS精密单点定位的基本原理GPS精密单点定位一般采用单台双频GPS接收机，利用IS提供的精密星历和卫星钟差，基于载波相位观测值进行的高精度定位。观测值中的电离层延迟误差通过双频信号组合消除，对流层延迟误差通过引入未知参数进行估计.1.2 时间系统RTLI内部使用GPT（GPST时间)用于GNSS的数据处理和定位算法。数据在RTLIB内部处理之前，需要转换成GPT时间。使用GPT的原因是避免处理润秒。RTKLI使用以下结构体表示时间：typedeftructtime_t tme； tme（） expresed by stndd t_t doube ec; / facioof

2、scond under s*/ gimet;1.2.1 GST和UTC（Unieral TieCoodina）关系参考【图1】，参考【图2】：图1 转换关系公式图 Error! Bookmark not defined.通过使用GPS导航信息中的T参数，GPST到TC或者UTC到PST之前的转换可以用更准确的表达方式,如【图】。图 Error! Bookmark not defined.这些参数是由GPS导航消息提供的。1.2.2 BDT（北斗导航卫星系统时间）DT(北斗导航卫星系统时间）是一个连续的时间系统，没有润秒。开始历元的时间是【UC 200年1月1号 00：00:00】。北斗时间计算

3、公式【图4】：图 Error! Bookmark not defined.T和GPST时间转换同上面的PS一样，只不过UTC参数来自与北斗导航信息中.1.3 坐标系统接收机和卫星的位置在TK中表示为在ECE（地心地固坐标系）坐标系统中的X， Y,Z组件。1.3.1 大地坐标到ECEF坐标的转换转换公式如【图5】。第三个公式最后一行有错有错,应该为: (1 e2)+h)sin图 1参数说明： a:地球参考椭球的长半径 f： 地球参考椭球的扁平率: 椭球高度 : 纬度 : 经度当前版本的RTKLIB使用的值为【图】:图 Error! Bookmark not defined.图Error! Boo

4、kmark not defined.参考椭球体1.3.2 EF坐系到大地坐标的转换转换公式如【图8】图 Error! Bookmark not defined.1.3.3 本地坐标到EF坐标的转换在接收机位置的本地坐标,也被称为NU坐标,通常使用在SS导航处理.ECEF坐标到本地坐标转换的旋转矩阵表示为【图9】.图 Error! Bookmark not defined. r旋转矩阵参数说明：:接收机位置的纬度：接收机位置的经度通过使用Er和接收机的坐标【】,坐标ref【C】可以被转换到本地坐标的坐标rlocal，公式如【图0】。图 Error! Bookmark not defined.2

5、RPPP定位算法2.1 单点定位（pntpos):saoss 2:epos ： etvl1. 计算计算卫星位置、速度和时钟(satss）a) 通过广播星历计算卫星钟差（ephclk） 卫星编号到卫星系统的转换(satsys）根据卫星的编号，获取到对应的卫星导航系统。 选择星历(seleph）、传入信号传输时间，卫星编号，导航数据等参数。2、遍历导航数据,遍历导航数据里面的星历数据，判断星历数据的卫星编号是否和传入的卫星编号相等.、如果星历数据的卫星编号和传入的卫星编号相等，就计算星历参考时间（e）和信号传输时间的时间差.如果不相等,继续处理下一条星历数据。、判断计算出来的时间差,如果时间差大于

6、了允许的最大时间差,继续查找下一个星历数据。否则，判断时间差最小的星历数据，记录星历数据的位置。5、返回之前记录出来的星历数据。 使用广播星历计算卫星时钟偏差(phclk）1、 传入信号发射时刻的时间和星历数据.2、 计算信号发射时刻的时间和本时段钟差参数参考时间(星历参数to）的时间差。3、 通过下式计算钟差，这里还没有处理相对论校正项和td:（代码中有个迭代过程,资料上没看到写）卫星钟差计算出来之后,信号发射时刻的时间还要加上这个钟差。b) 计算卫星在信号发射时刻的位置、速度和时钟（sapos）根据星历表选项来选择不同的处理，如下: 广播星历（EHOTBRDC）:phs （广播星历到卫星位

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

8、（PHPT_SSRAPC): satos_ssr 广播+SR_CM（EPHOTSSROM)：sap_ssr SS LEX星历（EPHP_LEX）：lexepo2. 使用伪距估算接收机的位置，返回估算状态结果（ets）a) 伪距残差（resco) 把ece坐标系转换成大地坐标系（ece2os）1、 按照以下公式做转换,暂时还没看懂. 计算几何距离和接收机到卫星的单位矢量(odis)1、 用卫星的坐标向量做欧几里德范数,返回值和地球长半轴(WGS4)比较。小于地球长半轴(WS84）,返回;2、 计算卫星坐标和接收机坐标的差值向量。3、 用差值向量做欧几里德范数,再用上一步计算出来的差值向量和计算结

9、果做除法，得到视线向量。按照如下公式:4、 使用以下公式计算几何距离. 计算卫星方位角/仰角（satazel）1、 把接收机cef坐标转换到大地坐标;2、 判断高度是否大于地球半长轴(GS84)的负数值；3、 如果高度小于等于地球半长轴(GS84)的负数值，方位角为，仰角为I/；4、 如果高度大于地球半长轴（S8）的负数值,把ECF向量转换到局部坐标。然后对转换出来的坐标做内积.5、 使用如下公式计算卫星方位角和仰角。 伪距使用编码残差改正（pange)暂时没找到对应的文档对应。 电离层改正（ionoorr）n 通过广播电离层模型（kobuchar模型)计算出电离层延迟（omodel)1、 校

10、验传入的电离层模型参数,校验失败，使用默认的电离层模型参数;2、 使用以下公式计算出地球为中心的角度（半圆）；3、 使用以下公式计算子的电离层的纬度经度（半圆)；4、 使用以下公式计算地磁纬度。5、 计算本地时间,返回值按这个公式【ttflo（t/8640)640.;】处理，保证tt的范围（0=tt840）。6、 计算倾斜因子。7、 电离层延迟计算，公式如下. 对流层改正(tropcorr)n 通过标准大气压和stone模型计算对流层延迟（tropmodel)1、 使用以下公式计算总气压.2、 使用以下公式计算绝对温度。3、 使用以下公式计算水蒸汽的分压。代码中没有使用.4、 使用以下公式计算

11、【Saastoin模型】。 伪距残差残差值=伪距-（物理距离+dtr光速时钟偏差电离层误差对流层误差） 时钟系统和接收器的偏置补偿 伪距测量误差方差（vaerr）b) 方差权重值 （weigh bvariace）c) 最小二乘估计 (lsq) 最小二乘估计通过求解正规方程（X（A ) - 1 A * Y）1、 计算矩阵*Y的结果，保存到矩阵Ay;2、 计算矩阵A的结果，保存到矩阵；3、 求的矩阵的逆矩阵,结果保存到矩阵Q；4、 最后xQAyd) 欧几里德范数(rm)公式:1、 通过最小二乘法估算出来的参数dx【长度为】，分别加到位置向量上x【长度为】.2、 对这个速度向量做欧几里德范数，返回值

12、同1E4比较,小于这个值时,就得到估算出的接收机的位置为向量x的值。e) 验证求解(valso）1、 对伪距残差值做内积,然后同卡方分布的自由度分布值(lpa=0.01）做比较，大于卡方分布的自由度分布值的数据，无效。2、 计算ps；暂时还没没找到资料。3、 用计算出来的dops和配置的最大dos阈值比较，大于配置的最大dps阈值的数据无效。3. 使用多普勒估算接收机速度(estvel）a) 多普勒残差（sp） 把ef坐标系转换成大地坐标系(ef2pos) 把yz坐标转换成eu坐标（xyzenu) 计算ecef中的瞄准线向量 计算相对于接收机在CE中的卫星速度 多普勒残差计算1、速率公式：b)

13、 最小二乘估计 （） 最小二乘估计通过求解正规方程（X =(A A） A* Y）1、 计算矩阵A*Y的结果，保存到矩阵Ay;2、 计算矩阵A的结果，保存到矩阵Q;3、 求的矩阵Q的逆矩阵，结果保存到矩阵;4、 最后x=QAc) 欧几里德范数（nor）公式：3、 通过最小二乘法估算出来的参数dx【长度为4】，分别加到速度向量上x【长度为4】。4、 对这个速度向量做欧几里德范数，返回值同1E6比较，小于这个值时，就得到估算出的接收机速度为向量x的值。2.2 精确定位 （pppos)1:dstate_ppp 2: atos 3: eecipse ： rpp 5： r_pp : flter 7： es

14、_ppp1. 暂时更新状态（usat_pp）a) 位置更新b) 时钟更新c) 对流层参数更新d) 相位偏差更新1、 通过LL检测周跳；2、 通过Geomet-ee相位跳变检测周跳（如果双频测量值可用）；2. 计算卫星的位置和速度和时钟（satposs)c) 通过广播星历计算卫星钟差（ehclk） 卫星编号到卫星系统的转换(satsy)根据卫星的编号，获取到对应的卫星导航系统。 选择星历（elph）1、传入信号传输时间,卫星编号，导航数据等参数。2、遍历导航数据,遍历导航数据里面的星历数据，判断星历数据的卫星编号是否和传入的卫星编号相等.3、如果星历数据的卫星编号和传入的卫星编号相等,就计算星历

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

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

17、pos 广播SBAS(EHOPTBA)：sos_sbas 广播+SR_PC（EPHOP_SSRA)： satpo_sr 广播SRCOM（PHPT_SROM）：saosssr QZSS EX星历（EHOPT_LEX）：lexeph2po3. 排除对被遮蔽卫星的测量（tetecipse)a) 获得太阳和月亮在cef中的位置（snons）b) 归一化三维矢量(norv3）c) 欧几里德范数的向量（nm）d) 计算太阳-地球-卫星角（dot)4. 相位和编码残差(r_p)a) 地球潮汐改正(tidedisp）b) 计算几何距离(eodist)1、 用卫星的坐标向量做欧几里德范数，返回值和地球长半轴（G

18、8）比较。小于地球长半轴(WGS8),返回-1；2、 计算卫星坐标和接收机坐标的差值向量。3、 用差值向量做欧几里德范数，再用上一步计算出来的差值向量和计算结果做除法,得到视线向量。4、 使用以下公式计算几何距离。备注:这个公式是地球自转改正公式.地球自转改正产生地球自转改正的机制在于地固系随地球的自转而旋转， 地固系是非惯性系统， 用地固系中的坐标计算卫星到接收机几何距离时, 由于卫星位置和接收机位置是两个不同时刻的位置矢量, 而且这两个时刻的地固系相对于惯性系是变化的, 因此要考虑地球自转引起的距离变化。其改正公式为：其中， w 为地球自转角速度； C 为光速; xsi和x ri分别代表卫

19、星位置矢量和测站位置矢量的分量， i= ， 2 对应x， y 分量。地球自转对纬度影响甚小,对经度影响最大， 其次是高度。同时地球的自转影响还跟测站的经纬度以及测站与卫星的几何关系有关， 对两极的测站， 影响为零; 对赤道上的测站影响最大； 当卫星在测站子午面内影响为零；卫星在测站东方时,影响为负, 卫星在测站西方时，影响为正。c) 卫星方位角仰角计算（aazel)1、 把接收机ecef坐标转换到大地坐标；2、 判断高度是否大于地球半长轴（S84)的负数值；3、 如果高度小于等于地球半长轴(WGS4)的负数值，方位角为0,仰角为I/；4、 如果高度大于地球半长轴（WGS84）的负数值，把CF向

20、量转换到局部坐标。然后对转换出来的坐标做内积.5、 使用如下公式计算卫星方位角和仰角。d) 对流层延迟修正: 通过标准大气压和astminn模型计算对流层延迟(romod）1、 使用以下公式计算总气压.2、 使用以下公式计算绝对温度。3、 使用以下公式计算水蒸汽的分压。代码中没有使用.4、 使用以下公式计算【asanen模型】。e) 卫星天线模型(sttv）f) 接收机天线模型,通过天线相位中心参数计算天线偏移 （tmdel）g) 相位缠绕校正（widupcorr）h) 电离层和天线相位校正测量（corrmeas）i) 计算卫星时钟和对流层延迟5. 卡尔曼滤波（er）卡尔曼滤波按照以下的公式更

21、新状态: P H （HP * H+ R）-1 ;xp = x + K v ;p = （I H) ;参数解释:x： 状态向量 （n x 1）p：状态的协方差矩阵 (n x n）H： 设计矩阵的转置矩阵 （n x ）v: 创新(测量模型） （m x ）R： 测量误差的协方差矩阵 (mx m)n，m: 状态和测量值的数量xp：状态更新后的向量(nx 1）P： 状态更新后的协方差矩阵 (n x）6. 组合残差 (res_ppp)7. 解析PPP整周模糊度(pppab)a) arageLC（aveaeC）b) 固定宽巷模糊(fixambL） 计算LC波长度(m）（lam_LC） 宽巷模糊 宽巷模糊度的方

22、差c) fi narrwlane amiguit 【AR mode：PPPAR】d) ix aro-lane ambiguity by ILS 【A o： PP-AR ILS】3 以下是上面一些算法的参考资料,信息来自于Rtlb的手册 3.1 接收机和卫星天线之间的几何距离对应rtklib中的函数extern double geodist(const double *rs, const double *rr, double *e)3.2 卫星方向的方位角和仰角对应rtklib的函数extern double satazel(const double *pos, const double *e,

23、 double *azel)3.3 GPS、il和NSS广播星历和时钟对应rtklib的函数extern void eph2pos(gtime_t time, const eph_t *eph, double *rs, double *dts, double *var)3.4 北斗广播星历和时钟对应rtklib的函数extern void eph2pos(gtime_t time, const eph_t *eph, double *rs, double *dts, double *var) 3.5 对流层和电流层模型3.5.1 对流层模型对应rtklib的函数extern double tr

24、opmodel(gtime_t time, const double *pos, const double *azel,double humi)3.5.2 电离层模型对应rtklib的函数extern double ionmodel(gtime_t t, const double *ion, const double *pos,const double *azel)3.6 单点定位3.6.1 线性最小二乘估计对应rtklib的函数extern int lsq(const double *A, const double *y, int n, int m, double *x,double *Q)

25、3.6.2 估算接收机的位置和钟差static int estpos(const obsd_t *obs, int n, const double *rs, const double *dts,const double *vare, const int *svh, const nav_t *nav, const prcopt_t *opt, sol_t *sol, double *azel, int *vsat, double *resp, char *msg)3.6.3 估算接收机的速度和钟漂static void estvel(const obsd_t *obs, int n, const

26、 double *rs, const double *dts, const nav_t *nav, const prcopt_t *opt, sol_t *sol, const double *azel, const int *vsat)3.6.4 求解验证和RAIM Estatic int valsol(const double *azel, const int *vsat, int n, const prcopt_t *opt, const double *v, int nv, int nx, char *msg)static int raim_fde(const obsd_t *obs, int n, const double *rs, const double *dts, const double *vare, const int *svh, const nav_t *nav, const prcopt_t *opt, sol_t *sol, double *azel, int *vsat, double *resp, char *msg)3.7 精密定位3.7.1 用于P的ZD测量模式3.7.2 接收机天线相位中心模型3.7.3 卫星天线相位中心模型3.7.4 由地球潮汐引起的位移3.7.5 相位缠绕校正3.7.6 通过PPP估算接收机的位置38 / 38