CN108169774A - 支持rtppp和rtk的多模gnss单频周跳探测与修复方法 - Google Patents

Abstract

本发明提供了一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法，包括构建历元间差分，消除与时间相关的电离层误差、对流层误差、不同系统间钟差时延误差，经整体最小二乘估计，若位置改正数协方差则循环选择5颗卫星，找出可信度最高的无周跳组合，得到位置和钟差漂移改正数，反算出新的残差向量，该残差向量也即历元间的周跳浮点解，再在一个可接受的范围内对周跳浮点解进行“四舍五入”，进而判定周跳是否可修复或为粗差。本发明根据检测量含义的不同设置了三种不同的经验阈值，降低了周跳修复过程中漏判和误判的概率。

Description

支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法

技术领域

本发明涉及一种实时GNSS数据处理技术，特别是一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法。

背景技术

2012年北斗卫星导航系统正式向亚太地区提供导航、定位、授时服务，继GPS/GLONASS组合使用之后，再添加新的可用导航系统。为了得到高精度、高可靠性的定位结果，因多模GNSS具有更多的观测卫星数，使其组合定位成为新的研究热点；另一方面，差分或者精密单点等高精度定位模式通常依靠于载波相位观测值。因为卫星信号被遮挡、外界干扰等可能原因引起初始整周计数产生偏差，即发生了周跳，相应的会转移成模糊度的不连续性，这是使用载波相位观测值必需要处理的问题。

周跳探测的思路是去除掉卫星定位过程中的其他误差，构建只剩下与初始整周模糊度相关的检测量，当周跳发生时，检测量会产生较大的跳变，若超过了阈值，即认为当前历元发生了周跳，具体可分为两类。一类是每个卫星独立探测，这一类方法直接利用不同频率、不同类型的原始观测值进行组合，如高次差法、多项式拟合法、伪距相位组合、MW组合、多普勒积分法、传感器辅助探测法、电离层残差法等。高次差法和多项式拟合法主要用于静态情况下，需要开始阶段的几个历元是干净的观测值，探测效果较差且不适用于实时动态应用；伪距相位组合和MW组合受限于伪距的测量精度，只适用于大周跳；多普勒观测了载波频率瞬间的漂移，但其精度跟接收机内部晶振有关，不同品牌可能不尽相同，同样只适用于大周跳；为了消除伪距测量精度低的影响，得到卫星至接收机间精准的几何距离，借鉴MW组合探测原理，很多学者提出使用惯导等其它外部传感器提供精度较高的位置初值，但是该方法对于惯性传感器的精度有一定要求，并且不适用于普通的单GNSS定位模式；电离层残差法依赖于双频或三频载波观测值，不适用于单频接收机，并且对特定周跳不敏感。因而这类方法多是受限于观测量自身精度差的影响，而可探测出小周跳的方法通常需要双频以上的数据，并且即便可成功探测出出现周跳的卫星，也难以进行修复，再接下来的定位解算中，只能去除该卫星，重新进行模糊度的固定工作，损失了系统的精度。

第一类方法孤立的看待每颗卫星，而完全忽略了定位时卫星间的整体性，尤其是多系统的优势。因而另一类方法是将所有的卫星看成一个整体，前后历元共视的所有卫星进行历元间差分，在短时间内可基本消除时间相关性误差，如对流层误差、电离层误差，在进行整体最小二乘后，对残差进行假设性检验判断是否含有周跳，典型代表如三差法。单个频率的三差法，可以完成单频的周跳探测工作，但因采用星间差分，导致观测量之间不独立，因而当参考卫星出现周跳或粗差时，会引起其他所有观测值均含误差，此时只能进行循环迭代，重新选用未发生周跳的卫星作为参考卫星。对于成功探测出周跳的卫星，目前通常采用LAMBDA方法对周跳的数值修复，或者利用该当前历元该卫星的残差值和之前无周跳历元对应的残差值进行计算，得到周跳值。前者搜索成功的概率受限于利用伪距观测值得到的初值偏差大小，与定位解算阶段进行模糊度的重新固定原理基本一致；而后者方法中，最小二乘会对所有观测值携带的误差再次平衡后才能得到某个观测量的残差，因而成功修复的概率存在诸多不确定因素。

发明内容

本发明提供了一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法，用以解决现有技术中难以探测小周跳、修复概率低等问题。

实现本发明目的的技术方案为：一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法，包括以下步骤：

步骤1，获取游动站接收机观测到的GPS、BDS、GLONASS原始数据，若是RTPPP定位模式，则仅需该数据；若是RTK定位模式，则需同时获取基准站的观测数据；

步骤2，获得相邻历元t_m+1、t_m的载波相位观测数据，获得GPS、BDS、GLONASS系统历元间差分观测方程的残差V

V＝BΔX-L+ΔN (1)

其中，矩阵由历元t_m+1时刻卫星至接收机的三维视线向量构成，

待估未知数表示接受机在t_m+1时刻的位置改正数和接收机的GPS系统钟差漂移，

观测向量由载波相位观测值和历元t_m+1、t_m时刻的近似伪距构成，λ表示载波的波长；表示载波相位观测值在历元间作差；表示历元t_m+1时刻的近似伪距；表示历元t_m时刻的近似伪距；上角标j、k、l依次代表GPS、BDS、GLONASS卫星；

周跳值ΔN表示卫星在历元t_m+1、t_m时刻发生的周跳值；

步骤3，利用最小二乘计算式(1)构成的方程组，得到相应的残差向量，计算出未知数X的协方差矩阵其中σ₀表示单位权中误差，P表示观测值的权矩阵；选择协方差矩阵的前三维位置改正数协方差的迹作为周跳探测的检测量；

步骤4，比较协方差矩阵的迹和预设的经验阈值threshold1；当迹超过预设的经验阈值threshold1，则有周跳发生ΔN≠0，执行步骤5；若协方差矩阵小于预设的检验阈值threshold1，则表示无周跳发生ΔN＝0，直接退出周跳探测程序；

步骤5，在观测的n颗卫星中，循环选用五颗卫星进行最小二乘，采用步骤4的判断条件，在满足的五颗卫星中依次选择四颗构成方程组，解该方程组得到改正数ΔX，比较改正数ΔX互差diff(ΔX)的模||diff(ΔX)||和预设的经验阈值threshold2，若都满足条件||diff(ΔX)||＜threshold2，则记录该5颗卫星的组合，若不满足则放弃该组合；

步骤6，在记录的卫星组合中找出迹最小的组合，标记为最佳组合，得到接收机位置改正数和钟差偏差ΔX，进入步骤7，若没有满足步骤5中条件的组合，则表示t_m+1时刻有周跳，则剔除掉某一个系统的卫星后，重复步骤2至步骤6；

步骤7，步骤6中得出的接收机位置和钟差改正数，重新代入公式(1)中，反算出周跳ΔN的浮点解N_CS。

本发明与现有技术相比，具有以下优点：(1)本发明经过对大量不同品牌接收机的长时间数据测试后，认为不同系统间ISB时延参数在短时间内可认为不变，通过历元间差分时可直接消除，进而统一了多模和单模系统的未知数参数个数，即把传统的三模系统中需要的7颗必要观测数统一成5颗，降低了多模系统中必要观测数的个数，对于城市峡谷等恶劣观测条件下周跳探测与修复尤为重要；(2)本发明设置中根据检测量含义的不同，设置了三种不同精度的经验阈值，并且其中一个实现自适应，根据具体数据环境而变化，降低了周跳修复过程中漏判和误判的概率，提高了周跳修复的成功率和正确率；(3)本发明同时适用于单频RTPPP定位和单频RTK定位模式.

下面结合说明书附图对本发明做进一步描述。

附图说明

图1为本发明的方法流程示意图。

图2为5个站点不同品牌接受机30s间隔的ISB差值统计示意图。

具体实施方式

结合图1，一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法，,包括以下步骤：

步骤1：获取游动站接收机观测到的BDS、GPS、GLONASS原始数据，若是实时RTPP定位模式，则仅需该数据；若是RTK定位模式，则需同时获取基准站的观测数据。

步骤2：得到相邻历元t_m+1、t_m的原始数据后，依次按照GPS、BDS、GLONASS系统历元间差分观测方程获得

V＝BΔX-L+ΔN (17)

其中

设计矩阵由历元t_m+1时刻卫星至接收机的三维视线向量构成；

待估未知数表示接受机在t_m+1时刻的位置改正数和接收机的GPS系统钟差漂移；

观测向量由载波相位观测值和历元t_m+1、t_m时刻的近似位置构成；

周跳值ΔN表示卫星在历元t_m+1、t_m时刻发生的周跳值。

步骤3：利用最小二乘计算式(17)构成的方程组，得到相应的残差向量，计算出未知数X的协方差矩阵其中σ₀表示单位权中误差。为了消除残余误差的影响，选择位置改正数的协方差矩阵作为周跳探测的检测量。

步骤4：比较协方差矩阵和预设的经验阈值，当协方差矩阵超过预设的经验阈值threshold1，则有周跳发生ΔN≠0，执行步骤5；若协方差矩阵小于预设的检验阈值threshold 1，则表示无周跳发生ΔN＝0，可直接退出周跳探测程。

步骤5：在观测的n颗卫星中，循环选用五颗卫星进行最小二乘，按照步骤4的判断条件，对无周跳发生的五颗卫星组合，在以上五颗卫星中依次选择四颗，解方程组得到改正数ΔX，比较改正数ΔX的互差，是否都满足条件

Det_ΔX＜threshold2，若满足条件则记录该5颗卫星的组合，若不满足则放弃该组合。

步骤6：若观测了n颗卫星，共有个组合中，找出满足步骤5的条件中残差最小的组合，标记为最佳组合，得到接收机位置改正数和钟差偏差ΔX，进入步骤6，若没有满足步骤5中条件的组合，则表示t_m+1时刻有周跳，且无法修复，针对多模系统，则可剔除掉某一个星座后，重复步骤2至步骤6。

步骤7：步骤6中得出的接收机位置和钟差改正数，重新代入公式(17)中，即可反算出周跳的浮点解N_CS。首先比较步骤6中最佳5颗星组合对应的周跳N_CS是否为零，若不是，则表示周跳修复失败，针对多模系统，则可剔除掉某一个星座后，重复步骤2至步骤6；若为零，则进一步考虑周跳具有的整数特性，周跳修复值可在一定范围内进行四舍五入，即满足公式(16)

|N_CS-(N_CS)_ROUND|＜threshold3 (16)

若满足则认为周跳修复成功，若不满足，则认为这该卫星的周跳难以修复，标记为粗差。

具体的，步骤2中，采用历元间差分的方式，消除了不同系统间时延偏差，进而把多模和单模的必要观测数统一为5颗。

具体的，步骤4中，经整体最小二乘后，仅采用的是位置改正数的协方差矩阵判定是否有周跳，而不是所有的改正数。进一步对于RTPPP而言，周跳判定采用的阈值threshold1跟历元间隔成正比；对于RTK而言，周跳判定采用的阈值threshold1跟历元间隔、基线距离成正比。

具体的步骤5和步骤6，选择出5颗卫星的最可信度最高的组合，需经过三层筛选：

(1)需满足步骤4的检验条件；

(2)任意4颗卫星构成的四元一次方程组解的互差，应小于阈值threshold2,因而threshold2不应较大；

(3)所有的组合中，以满足条件(1)和(2)中位置改正数协方差矩阵的迹最小的一组，作为可信度最高的无周跳组合。

具体的，步骤2中的具体过程为：

对于GPS、BDS、GLONASS在历元t_m时刻L1/B1/L1频率上的载波相位观测方程如下：

式中，角标j、k、l依次表示GPS、BDS、GLONASS系统；是载波相位观测值，若是RTK模式，则为载波相位站间差值，若是RTPPP模式，则为载波相位观测值；λ为卫星对应的载波波长；ρ为站星间距离；t_G、t_B、t_R、t^j、t^k、t^l依次为GPS、BDS、GLONASS系统的接受机钟差和卫星钟差漂移量；I、T为电离层误差、对流层误差；N表示初始整周模糊度；IFB为GLONASS系统不同卫星之间的频率偏差；ε包含多路径和测量噪声偏差。

利用伪距单点定位，计算出t_m历元站星距的近似值ρ₀，在相邻历元t_m+1、t_m之间求差，可得：

式(4)(5)(6)中，Δ表示历元间差分算子，并且(3)式中IFB在较长的一段时间内保持稳定，历元间差分可直接消除。若是RTK模式，利用广播星历计算出卫星轨道、卫星钟差，并且站间单差可消除其影响，若是RTPPP模式，可得到精密轨道和精密钟差；电离层误差、对流层误差等先利用模型改正削弱，又因具有较强的时间相关性，因而在短时间内忽略此项。

现考虑式(5)、(6)中BDS、GLONASS系统的接收机钟差Δt_B、Δt_GL与GPS接收机钟差存在如下关系：

Δt_B＝Δt_G+ΔISB_B-G (7)

Δt_R＝Δt_G+ΔISB_R-G (8)

式中，ISB_B-G、ISB_R-G依次为BDS、GLONASS与GPS的系统间时延偏差。为了验证其稳定性，选择不同品牌的5台接收机，在2017年1月5日观测24h，并且使用Gamit软件同IGS站联网解算，获取站点精确坐标，结合精密星历反算出BDS、GLONASS与GPS的30s系统间时延偏差，如图2所示。

图2中可看出，BDS、GLONASS与GPS的30s间隔接受机钟差偏差的均值在为0.008ns，中误差约为0.04ns，其等效距离分别约为2.4mm、12mm，因而若在短时间内，ISB差值将远小于第一频率的波长，故在历元间差分方程中，可忽略此项，即多系统的接收机钟差可统一为单一系统钟差，因此式(4)～(6)，可写成：

式(9)～(11)中的Δρ^j可表示成：

式中，分别表示历元t_m+1、t_m时刻的近似距离，代入式(9)～(11)中，可得：

当无周跳发生时，ΔN为零，则式(13)～(15)可进一步整理成：

V＝BΔX-L+ΔN (17)

其中：

在对式(17)进行最小二乘时，近似距离引入的误差等其它未被完全消除的误差会主要被接收机钟差偏移Δt_G吸收，因而残差基本属于正态分布。若发生了1周的跳变，则等价于出现与波长长度相等的粗差。至此，得到了一种实时动态多模GNSS单频周跳探测量。

因为未知量X受到历元间电离层残差、对流层残差大小等未被消除的误差的影响，其残余误差跟历元间隔长短有关，因而在RTPPP模式下经验阈值threshold1与历元间隔长短成正比；RTK模式下，残余误差量还跟基线距离相关，故经验阈值threshold1与历元间隔、基线距离成正比。

本发明中通过大量的实验分析，得出多模系统之间的钟差时延ISB经历元间差分后可忽略不计，是得未知数在多模系统和单模系统得到了统一，这对于在城市峡谷复杂环境下可用卫星数较少的情况有显著作用。在观测的n颗卫星中，循环选用五颗卫星进行最小二乘，与预设的经验阈值threshold 1进行比较，若五颗卫星的位置改正数协方差值超过threshold1，则放弃该组合；若小于经验阈值threshold1，对于最优的五颗卫星中，任意的四颗卫星构成的四元一次方程组的解的互差都应在预设阈值threshold2之内，若满足则记录下该五颗卫星的组合，若不满足则同样放弃该组合。

经验阈值threshold2是为了测试所选五颗卫星是否为最优组合，即在无多余观测条件下，任意四颗卫星构成的四元一次方程组的解可认为是无差别的，因而该数值不应太大。

在全部观测卫星个组合中，选择出满足threshold2阈值且位置改正数协方差值最小的组合，即认为该五颗卫星是可信度最高的无周跳组合，得到接收机位置改正数和钟差偏差X；若没有满足threshold2阈值的组合，则表示t_m+1时刻有周跳，且无法修复。

经验阈值threshold3是为了约束修复出的周跳符合整数特性，本发明中取为0.2周。

步骤7反算出的残差向量作为周跳的浮点解，若能在阈值threshold3的范围内进行“四舍五入”，则认为该颗卫星的周跳可成功修复；若不能再阈值threshold3的范围内进行“四舍五入”，则认为该颗卫星有粗差。此处的阈值threshold3用于约束周跳的整数性，本发明取为0.2。

Claims (9)

1.一种支持RTPPP和RTK的多模GNSS单频周跳探测与修复方法，其特征在于，包括以下步骤：

步骤1，获取游动站接收机观测到的GPS、BDS、GLONASS原始数据；

步骤2，获得相邻历元t_m+1、t_m的载波相位观测数据，获得GPS、BDS、GLONASS系统历元间差分观测方程的残差V

V＝BΔX-L+ΔN (1)

其中，矩阵由历元t_m+1时刻卫星至接收机的三维视线向量构成，

待估未知数表示接受机在t_m+1时刻的位置改正数和接收机的GPS系统钟差漂移，

观测向量由载波相位观测值和历元t_m+1、t_m时刻的近似伪距构成，λ表示载波的波长；表示载波相位观测值在历元间作差；表示历元t_m+1时刻的近似伪距；表示历元t_m时刻的近似伪距；上角标j、k、l依次代表GPS、BDS、GLONASS卫星；

周跳值ΔN表示卫星在历元t_m+1、t_m时刻发生的周跳值；

步骤3，利用最小二乘计算式(1)构成的方程组，得到相应的残差向量，计算出未知数X的协方差矩阵其中σ₀表示单位权中误差，P表示观测值的权矩阵；选择协方差矩阵的前三维位置改正数协方差的迹作为周跳探测的检测量；

步骤4，比较协方差矩阵的迹和预设的经验阈值threshold1；当迹超过预设的经验阈值threshold1，则有周跳发生ΔN≠0，执行步骤5；若协方差矩阵小于预设的检验阈值threshold1，则表示无周跳发生ΔN＝0，直接退出周跳探测程序；

步骤5，在观测的n颗卫星中，循环选用五颗卫星进行最小二乘，采用步骤4的判断条件，在满足的五颗卫星中依次选择四颗构成方程组，解该方程组得到改正数ΔX，比较改正数ΔX互差diff(ΔX)的模||diff(ΔX)||和预设的经验阈值threshold2，若都满足条件||diff(ΔX)||＜threshold2，则记录该5颗卫星的组合，若不满足则放弃该组合；

步骤6，在记录的卫星组合中找出迹最小的组合，标记为最佳组合，得到接收机位置改正数和钟差偏差ΔX，进入步骤7，若没有满足步骤5中条件的组合，则表示t_m+1时刻有周跳，则剔除掉某一个系统的卫星后，重复步骤2至步骤6；

步骤7，步骤6中得出的接收机位置和钟差改正数，重新代入公式(1)中，反算出周跳ΔN的浮点解N_CS。

2.根据权利要求1所述的方法，其特征在于，步骤2的具体过程为：

步骤2.1，获得历元t_m时刻GPS、BDS、GLONASS分别在L1、B1、L1频率上的载波相位观测数据

式中，角标j、k、l分别表示GPS、BDS、GLONASS系统，

是载波相位观测值，若是RTK模式则为载波相位站间差值；若是RTPPP模式则为载波相位观测值，

λ为卫星对应的载波波长，ρ为站星间距离，

t_G、t_B、t_R、t^j、t^k、t^l依次为GPS、BDS、GLONASS系统的接受机钟差和卫星钟差漂移量，

I、T分别为电离层误差、对流层误差，

N表示初始整周模糊度，

IFB为GLONASS系统不同卫星之间的频率偏差，

ε包含多路径和测量噪声偏差；

步骤2.2，利用伪距单点定位，计算出t_m历元站星距的近似值ρ₀；

步骤2.3，在相邻历元t_m+1、t_m之间的的载波相位观测数据求差，获得

其中，Δ表示历元间差分算子；

步骤2.4，考虑BDS、GLONASS系统的接收机钟差Δt_B、Δt_R与GPS接收机钟差存在如下关系：

Δt_B＝Δt_G+ΔISB_B-G (8)

Δt_R＝Δt_G+ΔISB_R-G (9)

其中，ISB_B-G、ISB_R-G分别为BDS、GLONASS与GPS的系统间时延偏差；

步骤2.5，不考虑卫星钟差漂移量、电离层误差、对流层误差和系统间时延偏差对数据测量的影响，将多系统的接收机钟差统一为单一系统钟差，即

步骤2.6，将式(13)代入式(10)～(12)中

Δρ＝1_xΔx+1_yΔy+1_zΔz+ρ_tm+1,0-ρ_tm,0 (13)

得

其中ρ_tm+1,0、ρ_tm,0分别表示历元t_m+1、t_m时刻的接收机和卫星之间的近似距离，(1_x,1_y,1_z)为历元t_m+1时刻卫星至接收机的三维视线向量；

步骤2.7，将式(14)～(16)进一步整理为：

V＝BΔX-L+ΔN (1)

其中，矩阵由历元t_m+1时刻卫星至接收机的三维视线向量构成，

待估未知数表示接受机在t_m+1时刻的位置改正数和接收机的GPS系统钟差漂移，

观测向量由载波相位观测值和历元t_m+1、t_m时刻的近似伪距构成，λ表示载波的波长；表示载波相位观测值在历元间作差；表示历元t_m+1时刻的近似伪距；表示历元t_m时刻的近似伪距；上角标j、k、l依次代表GPS、BDS、GLONASS卫星，

周跳值ΔN表示卫星在历元t_m+1、t_m时刻发生的周跳值；

当无周跳发生时，ΔN为零。

3.根据权利要求1所述的方法，其特征在于，步骤7的具体过程为：

步骤7.1，比较步骤6中最佳5颗卫星组合对应的周跳N_CS是否为零，若不是，则表示周跳修复失败，针对多模系统，可剔除掉某一个星座后，重复步骤2至步骤6；若为零，转步骤7.2；

步骤7.2，一步考虑周跳具有的整数特性，周跳修复值按照式(16)进行四舍五入

|N_CS-(N_CS)_ROUND|＜threshold3 (16)

其中，(*)_ROUND表示四舍五入取整

若满足则认为周跳修复成功，若不满足，则认为这该卫星的周跳难以修复，标记为粗差。

4.根据权利要求1所述的方法，其特征在于，步骤1中，若是RTK模式，同时获取基准站和游动站数据；若是RTPP模式，仅需游动站数据。

5.根据权利要求1所述的方法，其特征在于，步骤2中，采用历元间差分的方式，消除了不同系统间时延偏差，进而把多模和单模的必要观测数统一为5颗。

6.根据权利要求1所述的方法，其特征在于，步骤4中，经整体最小二乘后，仅采用的是位置改正数的协方差矩阵判定是否有周跳，而不是所有的改正数。

7.根据权利要求6所述的方法，其特征在于，对于RTPPP而言，周跳判定采用的阈值threshold1跟历元间隔成正比；对于RTK而言，周跳判定采用的阈值threshold1跟历元间隔、基线距离成正比。

8.根据权利要求3所述的方法，其特征在于，步骤7反算出的残差向量作为周跳的浮点解，若能在阈值threshold3的范围内进行“四舍五入”，则认为该颗卫星的周跳可成功修复；若不能再阈值threshold3的范围内进行“四舍五入”，则认为该颗卫星有粗差。

9.根据权利要求8所述的方法，其特征在于，阈值threshold3用于约束周跳的整数性，取0.2。