单频GNSS相位稳定性监测的高频逐历元相位差方法

摘要

本发明公开了一种单频GNSS相位稳定性监测的高频逐历元相位差方法，在测站点上利用单频GNSS接收机进行高频观测，计算相邻历元的相位差，利用最小二乘方法，计算连续历元下相邻历元间隔内的测站三维(NEU)坐标改正数时间序列，通过时间序列统计分析，分别取得NEU坐标的平均值和标准偏差，计算异常比率。通过提出的方法解算连续观测异常比率，以此来评估单频GNSS接收机相位观测的稳定性。本发明利用单频GNSS数据，进行单点实时高频GNSS相位数据质量监测；无需高精度钟差、各种模型改正即能进行GNSS相位稳定性的评估。

主权项

一种单频GNSS相位稳定性监测的高频逐历元相位差方法，其特征在于，包括以下步骤：第一步，用单频GNSS接受机进行高频观测，计算相邻两个历元之间的相位差；利用单频GNSS接收机在某一测站P上进行高频观测，在观测卫星多于3颗的条件下，高频采集相位数据；对某一GNSS卫星j在历元t₁上对应的相位为在GNSS信号不失锁的情况下，在下一个历元t₂GNSS卫星j对应的测站P相位为式(1)和(2)中，ρ_j(t₁)为历元t₁的测站至卫星j的几何距离，即，${\rho }_j\left(t_1\right)=\sqrt{{[x_j\left(t_1\right)-x_P]}^2+{[y_j\left(t_1\right)-y_P]}^2+{[z_j\left(t_1\right)-z_P]}^2}$ρ_j(t₂)为历元t₂时测站P至卫星j的几何距离，即${\rho }_j\left(t_2\right)=\sqrt{{[x_j\left(t_2\right)-x_P]}^2+{[y_j\left(t_2\right)-y_P]}^2+{[z_j\left(t_2\right)-z_P]}^2}$[x_j(t₁),y_j(t₁),z_j(t₁)]为历元t₁的卫星j的三维空间坐标；[x_j(t₂),y_j(t₂),z_j(t₂)]为历元t₂时卫星j的三维坐标；(x_P,y_P,z_P)为测站P在WGS84下的坐标；S_clock(j)和R_clock分别为GNSS卫星j和接收机的钟差；S_hardware(j)和R_hardware分别为GNSS卫星j和接收机硬件延迟误差；D_troposphere(j)和D_ionosphere(j)分别为对流层和电离层延迟误差；N(j)为整周模糊度；E_multipath(j)为多路径效应；E_phasecenter(j)为天线相位中心偏差及其变化影响；ε₁是残余误差；ε₂是残余误差；对高频GNSS信号来说，在相邻的两个历元上，主要误差诸如钟差、大气延迟、硬件延迟、多路径效应、整周模糊度和天线相位中心偏差及其变化影响基本上保持不变，式(1)减去式(2)，可得式中，v(j)＝ε₁‑ε₂；通过式(3)得出，对高频GNSS信号来说，上述主要误差可以通过逐历元相位差移去或有效降低；相邻历元相位差的精度为其中，σ是相位精度，测站P的最初坐标是[x₀,y₀,z₀]，坐标改正是[Δx,Δy,Δz]；将式(3)在[x₀,y₀,z₀]处进行泰勒展开，忽略高阶项，获得线性形式为其中和是测站P到GNSS卫星j在历元t₁和t₂时的初始距离，Δx,Δy,Δz为相邻测站坐标改正数的三个时间序列；第二步，解算Δx,Δy,Δz设在历元t₁和t₂，测站P同时都观测到的GNSS卫星数为m，那么式(5)用矩阵表达为V＝AX+L   (6)式中A为系数矩阵，$A=\left[\begin{array}{ccc}\frac{\left(x^2\left(1\right)-x_0\right)}{{\rho }_0^2\left(1\right)}-\frac{\left(x^1\left(1\right)-x_0\right)}{{\rho }_0^1\left(1\right)}& \frac{\left(y^2\left(1\right)-y_0\right)}{{\rho }_0^2\left(1\right)}-\frac{\left(y^1\left(1\right)-y_0\right)}{{\rho }_0^1\left(1\right)}& \frac{\left(z^2\left(1\right)-z_0\right)}{{\rho }_0^2\left(1\right)}-\frac{\left(z^1\left(1\right)-z_0\right)}{{\rho }_0^1\left(1\right)}\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ \frac{\left(x^2\left(j\right)-x_0\right)}{{\rho }_0^2\left(j\right)}-\frac{\left(x^1\left(j\right)-x_0\right)}{{\rho }_0^1\left(j\right)}& \frac{\left(y^2\left(j\right)-y_0\right)}{{\rho }_0^2\left(j\right)}-\frac{\left(y^1\left(j\right)-y_0\right)}{{\rho }_0^1\left(j\right)}& \frac{\left(z^2\left(j\right)-z_0\right)}{{\rho }_0^2\left(j\right)}-\frac{\left(z^1\left(j\right)-z_0\right)}{{\rho }_0^1\left(j\right)}\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ \frac{\left(x^2\left(m\right)-x_0\right)}{{\rho }_0^2\left(m\right)}-\frac{\left(x^1\left(m\right)-x_0\right)}{{\rho }_0^1\left(m\right)}& \frac{\left(y^2\left(m\right)-y_0\right)}{{\rho }_0^2\left(m\right)}-\frac{\left(y^1\left(m\right)-y_0\right)}{{\rho }_0^1\left(m\right)}& \frac{\left(z^2\left(m\right)-z_0\right)}{{\rho }_0^2\left(m\right)}-\frac{\left(z^1\left(m\right)-z_0\right)}{{\rho }_0^1\left(m\right)}\end{array}\right],$一般情况下，m＞3，因此利用最小二乘法就可以对式(5)进行解算，获得X的最优估值；式(6)的最小二乘解为X＝(A^TPA)^‑1A^TPL   (7)方差‐协方差矩阵为$D_X={\sigma }_0^2{\left(A^{T}PA\right)}^{-1}---\left(8\right)$其中，单位权方差是${\sigma }_0^2=\frac{V^{T}PV}{m-3}---\left(9\right)$式中，P是权矩阵，初始权是P₀＝I，I为单位矩阵，通过式(9)可以得到新的单位权σ₀，$p\left(j\right)=\frac{{\sigma }_0}{\left|v\left(j\right)\right|+c}---\left(10\right)$其中，c是一个大于零的常数，通过公式(7)即可计算X；第三步，利用Δx,Δy,Δz序列评估相位稳定性相邻高频GNSS测站坐标改正数的三个时间序列Δx,Δy和Δz可以通过式(7)构建，从这三个时间序列中，可以得到平均值和对应的标准偏差，即MEAN_Δx,MEAN_Δy,MEAN_Δz,STD_Δx,STD_Δy和STD_Δy，并定义上界值和下界值分别为UB＝MEAN+kSTD   (11)LB＝MEAN‑kSTD   (12)其中，k是一个常数，如果Δx,Δy和Δz的一个解大于对应的上界值或小于对应的下界值，则被认为是非正常值；在一个时段内，通过连续高频观测获得n个异常解；连续观测异常比率AP被定义为$AP=\frac{n}{q-1}\times 100\%---\left(13\right)$其中，q是连续观测中的观测历元数；当这个连续时间段非常长，例如10天或超过10天时，计算连续观测异常比率AP，其结果作为连续观测异常比率的临界值，把临界值作为判断相位稳定性的一个临界点，用来判断平时GNSS测量时某一段时间内的相位稳定性情况；即在一个时间段内，如果连续观测异常比率AP大于这个临界值，则表明相位是不稳定的，反之则相位是稳定的；在一个GNSS测站上，某一段时间内，根据计算出的连续观测异常比率，结合连续观测异常比率的临界值就可以评估单频相位数据在某个时间段内的稳定性情况。