本发明涉及卫星定位导航GNSS数据处理领域，尤其涉及一种用于卫星定位导航等应用中单频CDMA卫星系统数据预处理的周跳探测与修复方法。
背景技术：
全球卫星导航系统(GNSS)发展至今，精密定位技术已较为成熟。高精度RTK定位技术可以提供高质量的实时定位服务，其应用范围已经扩展到导航定位、水利、国土、城市规划、国家重大工程建设等各个领域。载波相位观测值的应用是高精度GNSS测量的关键。在载波相位测量过程中，由于卫星信号被障碍物遮挡，大气条件恶劣或无线电信号干扰等原因，会导致计数器无法计数而使观测值发生整周跳变。即使只有一周的周跳也会对定位结果产生20cm左右的影响，因此在用载波相位观测值进行解算之前，必须要对周跳进行处理。
由于低成本终端与精密的测地型仪器相比性能较差，其周跳探测技术更具有挑战性。低成本导航终端一般仅接收单频信号，因而无法像双频、三频接收机，可以通过两个或三个频率间的组合来探测周跳。此外，由于硬件性能较差，低成本终端接收的数据往往具有信噪比低，易频繁中断，卫星失锁多，受多路径效应影响大等问题，这些特点均为单频周跳探测技术带来了困难和挑战。
目前对于单频观测值，一般可采用多项式拟合法，高次差分法，三差法及伪距和相位组合法等周跳探测方法。多项式拟合法这种基于观测值域逐个卫星的拟合方法模型较弱，对于频繁发生的周跳无法有效探测。高次差分法则是利用相位观测值高次差的规律性进行周跳探测。上述两种方法均适用于静态或低动态的GNSS解算。三差法的原理是利用三差观测值(星间、站间和历元间差分)对周跳的敏感性进行周跳探测。经历元间差分后模糊度参数被消去，可利用三差解基线的残差值判断有无周跳发生，或将周跳视为未知数与位置参数一同求解。但是三差法对多余观测量要求较高，当三差观测值个数小于4时该方法就无法使用，且实际应用中一般需要5-6个或更多观测值才能得到较好的探测效果。伪距和相位组合的方法也适用于单频观测值，但该方法的精度受伪距观测值精度的影响较大，难以实现较高精度的周跳探测。
综上所述，本发明提出一种基于位置域曲线约束的单频动态周跳探测与修复方法，克服了现有方法的局限性，具有更高的准确性、稳定性及探测精度。
技术实现要素：
本发明的目的在于提供一种单频动态周跳探测与修复方法，能够处理卫星定位导航等应用中单频CDMA卫星系统数据预处理的周跳问题。
为了达到上述目的，本发明提供了一种基于位置域曲线约束的单频动态周跳探测与修复方法，包括：
步骤L1：获取观测对象的单频实时GNSS观测数据，根据所述实时GNSS观测数据构建三差观测方程；(本领域常规技术)
步骤L2：选取长度为n的时间窗口，设当前历元为i+n，判断所述观测对象的当前累积观测历元数是否小于n，其中，4≤n≤14，i≥1，n及i均为正整数；
若是，进入步骤L3；
若否，进入步骤L4；
步骤L3：根据所述三差观测方程，获取第一周跳参数浮点解；
步骤L4：根据所述时间窗口中各历元的位置参数获取动力学模型，获取所述观测对象在该时间窗口中的运动状态；
步骤L5：根据所述动力学模型及窗口内各个历元的观测信息，获包括位置约束条件的观测方程组；
步骤L6：根据所述包括位置约束条件的观测方程，获取第二周跳参数浮点解。
可选的，在上述单频动态周跳探测与修复方法中，所述步骤L1包括：
获取所述观测对象的流动站与基准站的单频实时GNSS观测数据；
对所述实时GNSS观测数据进行数据预处理；
构建相位与伪距双差观测方程；
获取三差观测方程。
可选的，在上述单频动态周跳探测与修复方法中，对所述实时GNSS观测数据进行数据预处理的步骤包括：
流动站与基准站的单点定位、卫星截止高度角设置、时标校正、大气延迟改正、粗差探测与处理以及卫星和接收机的天线相位中心修正。
可选的，在上述单频动态周跳探测与修复方法中，所述步骤L3包括：
将所述三差观测方程线性化；
联立所述相位与伪距三差观测方程，获取所述第一周跳参数浮点解；
根据所述第一周跳参数浮点解进行单频动态周跳探测与修复。
可选的，在上述单频动态周跳探测与修复方法中，通过最小二乘法获取所述第一周跳参数浮点解。
可选的，在上述单频动态周跳探测与修复方法中，所述步骤L4包括：
根据观测对象的第(i)、(i+1)、...、(i+n-2)以及第(i+n-1)历元的位置参数获取动力学模型，获取所述观测对象在n个历元中的任一历元的运动状态。
可选的，在上述单频动态周跳探测与修复方法中，所述步骤L5还包括：
根据所述动力学模型，获取当前历元位置的预报值，并将其作为约束条件，与所述三差观测方程联立，获取所述包括位置约束条件的观测方程。
可选的，在上述单频动态周跳探测与修复方法中，所述观测对象在n个历元中的动力学模型为一元二次多项式。
可选的，在上述单频动态周跳探测与修复方法中，还包括：
步骤L7：根据所述第二周跳参数浮点解，进行整数周跳固定；
若固定成功，则输出所述整数周跳；
若固定失败，则周跳探测失败。
可选的，在上述单频动态周跳探测与修复方法中，所述步骤L3还包括：根据所述第一周跳参数浮点解，进行整数周跳固定；
若固定成功，则输出所述整数周跳；
若固定失败，则周跳探测失败。
综上所述，本发明采用位置域曲线拟合，通过观测对象的各个历元的位置参数的相关性，将各个历元的观测信息有效利用，增强了模型强度以及方法的可用性。
具体的，与现有技术相比，本发明具有以下优点：
现有的多项式拟合周跳探测方法普遍基于各个卫星观测值域的拟合，在逐个卫星拟合时，拟合结果容易受到各类观测误差的影响，模型强度较弱。且容易受到卫星频繁升降的影响，当周跳发生频率较高时无法有效探测。而基于位置域的曲线拟合方法用全体观测值联合解算，模型强度更强。当相邻历元的可用卫星发生变化或一颗卫星频繁周跳时均都能有效探测。
现有的三差周跳探测方法普遍利用当前历元的三差观测值进行解算，由于需要同时解算当前历元三差基线以及周跳参数，该方法对多余观测量要求较高。一旦当前历元的三差观测值少于4，该方法不可用。而结合了空间域曲线拟合的周跳探测方法，有效的引入了其他历元的观测信息，增加了多余观测值数量，对当前历元的多余观测没有要求。
现有单频周跳探测与修复方法均没有很好的顾及单频低成本导航型接收机的特点，当数据质量较差，信噪比低，卫星失锁频繁，可用卫星数较少时，现有方法均具有一定的局限性。本发明考虑了低成本终端GNSS数据的特点，当观测条件恶劣，遮挡较多，可用卫星较少时具有更高的稳定性及探测精度。
附图说明
图1为本发明一优选实施例中的单频动态周跳探测与修复方法流程示意图；
图2为图1中步骤S1的具体流程示意图；
图3为图1中步骤S3的具体流程示意图；
图4为图1中步骤S5的具体流程示意图.
具体实施方式
下面将结合示意图对本发明的具体实施方式进行更详细的描述。根据下列描述和权利要求书，本发明的优点和特征将更清楚。需说明的是，附图均采用非常简化的形式且均使用非精准的比例，仅用以方便、明晰地辅助说明本发明实施例的目的。
参考图1，本发明一优选实施例中，一种单频动态周跳探测与修复方法包括：
步骤S1：获取观测对象的单频实时GNSS观测数据，根据所述实时GNSS观测数据构建三差观测方程；
具体的，参考图2，步骤S1包括(步骤S1为本领域常规技术)：
步骤S11：获取所述观测对象的流动站与基准站的单频实时GNSS观测数据。
步骤S12：对所述实时GNSS观测数据进行数据预处理。
数据预处理包括但不限于为单点定位、卫星截止高度角设置、时标校正、大气延迟改正、粗差探测与处理以及卫星和接收机的天线相位中心修正，本发明对此不做任何限制。
优选的，采用的高度角定权公式为：
其中，下标i表示第i颗卫星，θ表示其高度角，σ0是验前单位权中误差，σ是当前观测值的验前标准差。
步骤S13：构建相位与伪距双差观测方程。
所述双差观测方程如下：
其中，表示差分算子。p，φ分别为单频伪距和相位观测值，ρ为卫地距，c为光速，λ为载波波长，N为整周模糊度，I为双差电离层延迟，T为对流层延迟残差，ε为随机噪声。经双差后，卫星和接收机的钟差以及初始相位偏差都被消除，电离层和对流层延迟也大大削弱，因此，在短基线情况下可以忽略不计。
步骤S14：获取三差观测方程。
GNSS周跳具有整数性和连续性两个特性。与粗差不同，周跳具有整数特性，并且将自周跳发生历元向后延续。因此，周跳的探测必须建立在观测数据前后历元间差分的基础上。所以在进行单频周跳探测时，需将双差观测值作历元间差分，形成三差观测方程：
其中，t表示第t个历元。此处忽略了对流层和电离层延迟。由于一颗卫星整周模糊度不随时间改变，因此未发生周跳时三差模糊度应为一个极小量，此外三差后的卫地距也为极小量；而一旦发生周跳，相应的三差观测值将发生跳变。因此三差观测值对周跳较为敏感，可利用这一特性来探测周跳。
步骤S2：选取长度为n的时间窗口，设当前历元为i+n，判断所述观测对象的当前累积观测历元数是否小于n，其中，4≤n≤14，n及i均为正整数；
若是，进入步骤S3；
若否，进入步骤S4。
步骤S3：根据所述式(3)及(4)，获取第一周跳参数浮点解。
具体的，参考图3，步骤S3包括：
步骤S31：首先将式(3)及(4)线性化，即将所述三差观测方程线性化，得到：
其中，为地固系中测站到卫星三个双差方向余弦；δX、δY以及δZ为基准站固定后流动站的坐标改正数，常数项为残差值。
步骤S32：联立相位与伪距三差观测方程，并将线性化后的三差观测方程写成矩阵形式，将浮点周跳看作未知参数与三差基线一同求解：
其中，和分别表示位置参数和周跳参数，A和B为系数矩阵。
优选的，本发明使用最小二乘法求解，本发明对此不作任何限制，其法方程为：
进一步可求得第一周跳参数浮点解：
可选的，步骤S3还包括：S33：根据所述第一周跳参数浮点解，进行整数周跳固定；
若固定成功，则输出所述整数周跳；
若固定失败，则周跳探测失败。
可选的，通过LAMBDA方法进行整数周跳固定。
步骤S4：根据所述三差观测方程，分别获取所述观测对象在第i、i+1、...、i+n-2以及第i+n-1历元的位置参数，求得该时间段内描述物体运动状态的动力学模型。
由于动态测量应用中，各历元间的位置参数具有高度的相关性。尤其在高采样率的情况下，可认为物体位置服从低阶多项式模型。因此可用曲线拟合的原理，根据时间窗口内各个历元的位置参数，求得该时间段内可以描述物体运动状态的动力学模型。此处，以一元二次多项式为例来描述物体的运动状态。动力学模型的具体求解方法见步骤S5.
步骤S5：根据所述时间窗口内n个历元的位置参数联合求得的动力学模型，并与当前历元的三差观测方程联立，获得带有动力学约束条件的观测方程组。
具体的，参考图4，本发明一优选实施例中的步骤S5包括：
步骤S51：对所述观测对象的第i、i+1、...、i+n-2以及第i+n-1历元的位置参数进行曲线拟合，获取所述观测对象在该时间窗口内的动力学模型。可选的，所述动力学模型可为低阶多项式，本发明一优选实施例中，所述动力学模型为一元二次函数，表示为：
xi＝a0+a1i+a2i2 (11)，
其中，下标i＝1,2,…，m，m表示历元数，为简化计算，用历元的顺序数代替时间t作为自变量，a0，a1，a2为多项式参数。
步骤S52：则所述n个历元中的任一历元的三差基线可被表示为：
将所述时间窗口内前n-1个历史历元的位置信息代入所确立的动力学模型中。
步骤S53：联立所述时间窗口内n个历元(包括当前历元)的三差观测方程以及所采用的动力学模型，得到所述带有动力学约束的三差观测方程组：
其中，下标n表示第n个历元，yi＝[PiT ΦiT]T，δx为各个历元三差基线参数的集合，z为待估周跳参数。设计矩阵分别为：
A＝blkdiag(A1,…,An)；
B＝[0 Im]T；
以及
a为拟合的多项式系数，
本实施例中，设当前历元三差相位观测值的个数m，则式(13)的第一个方程将窗口内各个历元的观测方程同时联立，此时的未知数个数为(3(n-1)+m)个，引入条件(13)的第二个式子后，可将第二个式子带入第一个式子，未知数个数变为(3×2+m)个，当所选窗口长度大于3时，利用历元间位置的相关性可以减少观测方程中未知数的个数，因此，本发明结合了空间域曲线拟合的周跳探测方法，有效的引入了其他历元的观测信息，增加了多余观测值数量。
可选的，获取所述包括位置约束条件的观测方程的步骤包括：
利用所述时间窗口内各历元的观测信息，即联立时间窗口内各个历元的三差观测方程，并将三差基线满足的曲线关系作为约束条件，与当前历元观测方程联立求解；或
利用时间窗口内各历元的位置信息，即通过时间窗口内各历元的位置信息及动力学模型，获取当前历元的预报位置，将该预报位置作为虚拟观测值，与当前历元的观测方程联立求解，本发明对此不做任何限制。
步骤S6：根据所述三差观测方程组，获取第二周跳参数浮点解。
联立式(13)，得法方程为：
优选的，本发明一优选实施例中，采用最小二乘原理，得第二周跳参数浮点解为：
可选的，还包括步骤S7：根据所述第二周跳参数浮点解，进行整数周跳固定；
若固定成功，则输出所述整数周跳；
若固定失败，则周跳检测失败。
可选的，通过LAMBDA方法进行整数周跳固定。
综上所述，本发明采用位置域曲线拟合，通过观测对象的各个历元的位置参数的相关性，将各个历元的观测信息有效利用，增强了模型强度以及方法的可用性。
具体的，与现有技术相比，本发明具有以下优点：
现有的多项式拟合周跳探测方法普遍基于各个卫星观测值域的拟合，在逐个卫星拟合时，拟合结果容易受到各类观测误差的影响，模型强度较弱。且容易受到卫星频繁升降的影响，当周跳发生频率较高时无法有效探测。而基于位置域的曲线拟合方法用全体观测值联合解算，模型强度更强。当相邻历元的可用卫星发生变化或一颗卫星频繁周跳时均都能有效探测。
现有的三差周跳探测方法普遍利用当前历元的三差观测值进行解算，由于需要同时解算当前历元三差基线以及周跳参数，该方法对多余观测量要求较高。一旦当前历元的三差观测值少于4，该方法不可用。而结合了空间域曲线拟合的周跳探测方法，有效的引入了其他历元的观测信息，增加了多余观测值数量，对当前历元的多余观测没有要求。
现有单频周跳探测与修复方法均没有很好的顾及单频低成
本导航型接收机的特点，当数据质量较差，信噪比低，卫星失锁频繁，可用卫星数较少时，现有方法均具有一定的局限性。本发明考虑了低成本终端GNSS数据的特点，当观测条件恶劣，遮挡较多，可用卫星较少时具有更高的稳定性及探测精度。
上述仅为本发明的优选实施例而已，并不对本发明起到任何限制作用。任何所属技术领域的技术人员，在不脱离本发明的技术方案的范围内，对本发明揭露的技术方案和技术内容做任何形式的等同替换或修改等变动，均属未脱离本发明的技术方案的内容，仍属于本发明的保护范围之内。