Precision analysis of BDS-3 multi-frequency undifferenced and uncombined precision orbit determination
-
摘要: 针对北斗三号卫星导航系统(BDS-3)五频点观测数据和非差非组合精密定轨理论，介绍了非差非组合观测模型和参数估计方法，提出了利用K均值聚类算法(K-means)进行测站选取的策略，分析并讨论了非差非组合方法的优势. 通过K-means和人工经验选取两种测站选取方案，分别使用BDS-3五频，B1C+B2a、B1I+B3I三种频率选择方式，利用30个观测站，对BDS-3中轨道地球卫星(MEO)和倾斜地球同步轨道卫星(IGSO)进行精密定轨处理. 结果表明：当接收B1C+B2a频点观测数据测站不足时，非差非组合方法可以通过利用五频观测数据增加观测数据数量、优化测站布局，提高定轨精度，与B1C+B2a频率组合相比，五频定轨结果切向(A)、法向(C)、径向(R)和三维(3D)方向均方根(RMS)月均值分别提升0.003 m、0.004 m、0.003 m和0.007 m；K-means算法选取的测站与人工经验选取相比，分布更加合理，定轨精度更高，三种频率选择方案MEO卫星3D RMS月均值精度分别提升0.009 m、0.017 m和0.009 m.
-
关键词:
- 北斗卫星导航系统(BDS) /
- 多频 /
- 非差非组合 /
- K均值聚类算法(K-means) /
- 精密定轨
Abstract: In view of the observation data of BDS-3 at 5 frequencies and the theory of un-differenced and un-combined (UDUC) precise orbit determination, this paper introduces the UDUC model and parameter estimation method, proposes the strategy of station selection using the K-means algorithm, and analyzes the advantages of UDUC method. Through two kinds of station selection schemes, manual experience selection and K-means, 3 frequency selection methods of BDS-3 5-frequency, B1C+B2a, B1I+B3I are used respectively, and 30 IGS observation stations are used to carry out precise orbit determination for BDS-3 MEO and IGSO satellites. The experimental results show that when the stations which can receive B1C+B2a frequency observation data are insufficient, the UDUC method can increase the number of observation data and optimize the station layout by using the 5-frequency observation data so that the orbit determination accuracy can be improved. Compared with B1C+B2a, the monthly average RMS of 5-frequency in A, C and R directions increase by 0.003 m, 0.004 m, 0.003 m respectively, and that of 3D RMS increases by about 0.007 m. The stations selected through the K-means algorithm are more reasonably distributed and have higher orbit determination accuracy than manual experience selection scheme. With the 3 frequency selection methods, the monthly average RMS of MEO satellites in A, C and R directions are improved by 0.009 m, 0.017 m, 0.009 m respectively. -
0. 引 言
用户进行精确导航、定位和授时(PNT)所需的精密轨道产品主要通过国际GNSS服务(IGS)组织等获取. 目前，精密轨道确定主要采用双频无电离层组合模型(IF)[1]. 随着全球卫星导航系统（GNSS）的不断发展，越来越多的导航系统能够提供更多频点的信号[2-3]. 北斗三号全球卫星导航系统(BDS-3)于2020年7月31日开通服务[4]，能够发射五个频点信号，分别是B1I(1561.098 MHz)、B3I(1 268.52 MHz)、B1C(1575.42 MHz)、B2a(1176.45 MHz)、B2b(1207.14 MHz). 多频IF组合模型存在频率组合繁杂、可扩展性差的问题，难以适应日趋丰富的多频点信号资源[1]. 文献[5]证明加入新的频点对定位精度有所改善，而对于定轨精度的影响，研究相对较少.
非差非组合因其能够充分利用多频点观测数据，观测模型简洁、表达形式统一、易扩展且不放大噪声等特点，成为近十年来的研究热点. 文献[3]对多系统非差非组合精密单点定位(PPP)理论方法进行了系统的阐述和研究. 文献[6]通过非差非组合观测模型实现了低轨卫星的精密定轨. 文献[7]使用非差非组合模型对GPS、Galileo、GLONASS三个系统进行了精密定轨处理，对电离层和模糊度处理进行了论述，通过与精密产品对比，证明其轨道产品精度与IGS各分析中心产品的精度相当. 文献[8]通过使用三频信号对GPS、Galileo、北斗二号卫星导航系统（BDS-2）进行了精密定轨，分析了第三频点信号对定轨精度的影响. 目前国内外对BDS-3非差非组合精密定轨研究较少，尤其是BDS-3的多频信号对定轨精度的影响、新体制信号B1C和B2a的定轨精度，以及不同测站分布不同频率选择时BDS-3定轨精度的分析和研究都相对较少.
因此，本文对非差非组合观测模型和参数估计方法进行了简要介绍，通过K均值聚类算法(K-means)优化测站选取策略，对数据处理过程中的误差和模型进行描述，通过两种实验方案，采用五频、B1C+B2a、B1I+B3I三种频率选择方式进行精密定轨，并对BDS-3不同测站分布、不同频率选择的定轨精度进行了对比分析.
1. 非差非组合定轨方法
1.1 观测模型
标准单系统GNSS观测模型的伪距和载波相位观测方程可表示为：
$$ \begin{split} P_{{\rm{r}},{{fre}}}^{\rm{s}} =& \rho _{\rm{r}}^{\rm{s}} + c(\delta {t_{\rm{r}}} - \delta {t^{\rm{s}}}) + T_{\rm{r}}^{\rm{s}} + I_{{\rm{r}},{{fre}}}^{\rm{s}} + \\ & {B_{{\rm{r}},{{fre}}}} + B_{{{fre}}}^{\rm{s}} + e_{{\rm{r}},{{fre}}}^{\rm{s}} , \end{split} $$ (1) $$ \begin{split} L_{{\rm{r}},{{fre}}}^{\rm{s}} =& \rho _{\rm{r}}^{\rm{s}} + c(\delta {t_{\rm{r}}} - \delta {t^{\rm{s}}}) + T_{\rm{r}}^{\rm{s}} - I_{{\rm{r}},{{fre}}}^{\rm{s}} + \\ & {\lambda _{{{fre}}}} \cdot ({b_{{\rm{r}},{{fre}}}} + b_{{{fre}}}^{\rm{s}} + {N_{{{fre}}}}) + e_{{\rm{r}},{{fre}}}^{\rm{s}} \end{split}. $$ (2) 式中：下标
$fre$ 为所选频点；$\;\rho _{\rm{r}}^{\rm{s}}$ 为站星几何距离，指卫星天线相位中心到接收机天线相位中心的几何距离，单位为m；$c$ 为真空中的光速；$\delta {t_{\rm{r}}}$ 和$\delta {t^{\rm{s}}}$ 为接收机和卫星钟差；$T_{\rm{r}}^{\rm{s}}$ 和$I_{{\rm{r}},{{fre}}}^{\rm{s}}$ 为对流层和电离层延迟；${B_{{\rm{r}},{{fre}}}}$ 和$B_{{{fre}}}^{\rm{s}}$ 为接收机和卫星伪距硬件延迟；$e_{{\rm{r}},{{fre}}}^{\rm{s}}$ 为残差；${\lambda _{{{fre}}}}$ 为波长；${b_{{\rm{r}},{{fre}}}}$ 和$b_{{{fre}}}^{\rm{s}}$ 为接收机和卫星载波相位延迟；${N_{{{fre}}}}$ 为载波相位整周模糊度.在精密定轨数据处理中，公式(1)~(2)中的观测模型需要通过泰勒级数展开进行线性化处理，重组后的观测方程为：
$$ \begin{split} P_{{\rm{r}},{{fre}}}^{\rm{s}} =& {\boldsymbol{u}}_{\rm{r}}^{\rm{s}}({\boldsymbol{\varPhi}} ({t_0},t){\boldsymbol{x}}_0^{\rm{s}} - {{\boldsymbol{x}}_{\rm{r}}}) + c(\delta {t_{\rm{r}}} - \delta {t^{\rm{s}}}) + T_{\rm{r}}^{\rm{s}}+ \\ & I_{{\rm{r}},{{fre}}}^{\rm{s}} + {B_{{\rm{r}},{{fre}}}} + B_{{{fre}}}^{\rm{s}} + e_{{\rm{r}},{{fre}}}^{\rm{s}} , \end{split} $$ (3) $$ \begin{split} L_{{\rm{r}},{{fre}}}^{\rm{s}} =&{\boldsymbol{u}}_{\rm{r}}^{\rm{s}}({\boldsymbol{\varPhi}} ({t_0},t){\boldsymbol{x}}_0^{\rm{s}} - {{\boldsymbol{x}}_{\rm{r}}}) + c(\delta {t_{\rm{r}}} - \delta {t^{\rm{s}}}) + T_{\rm{r}}^{\rm{s}}- \\ & I_{{\rm{r}},{{fre}}}^{\rm{s}} + {\lambda _{{{fre}}}} \cdot ({b_{{\rm{r}},{{fre}}}} + b_{{{fre}}}^{\rm{s}} + {N_{{{fre}}}}) + e_{{\rm{r}},{{fre}}}^{\rm{s}} \\ \end{split}. $$ (4) 式中：
${\boldsymbol{u}}_{\rm{r}}^{\rm{s}}$ 为站星几何距离的一阶偏导数；${\boldsymbol{\varPhi}} ({t_0},t)$ 为卫星轨道状态转移矩阵；${\boldsymbol{x}}_0^{\rm{s}}$ 为卫星参考历元轨道参数；${{\boldsymbol{x}}_{\rm{r}}}$ 为接收机质心位置.1.2 参数估计方法
非差非组合方法直接利用观测方程构建方程组进行参数估计、参数估计是通过函数模型利用获取的观测值确定未知参数的过程. 在GNSS数据处理中，函数模型即为观测方程，参数估计方法一般采用最小二乘平差方法，具体过程如下.
首先，将公式(3)~(4)的函数关系简化为
$$ {\boldsymbol{L}} = f({\boldsymbol{X}}) + {\boldsymbol{e}} .$$ (5) 式中：
${\boldsymbol{L}}$ 为伪距和载波相位观测值；$f({\boldsymbol{X}})$ 为观测量和未知参数之间的函数关系，${\boldsymbol{X}}$ 为未知参数；${\boldsymbol{e}}$ 为残差.通过将观测方程进行线性化，将其在初始状态处进行泰勒展开，截断高阶项，可得如下简化的线性方程：
$$ {\boldsymbol{L}} - f({{\boldsymbol{X}}_0}) = {\boldsymbol{H}}({\boldsymbol{X}} - {{\boldsymbol{X}}_0}) + {\boldsymbol{e}}. $$ (6) 式中，
${\boldsymbol{H}} = {\left. {\displaystyle\frac{{\partial F}}{{\partial X}}} \right|_{{X_0}}}$ ；${{\boldsymbol{X}}_0}$ 为初始状态值.令
$\Delta {\boldsymbol{L}} = {\boldsymbol{L}} - f({{\boldsymbol{X}}_0})$ ，$\Delta {\boldsymbol{X}} = {\boldsymbol{X}} - {{\boldsymbol{X}}_0}$ ，则式(6)可化简为$$ \Delta {\boldsymbol{L}} = {\boldsymbol{H}}\Delta {\boldsymbol{X}} + {\boldsymbol{e}}. $$ (7) 观测值的协方差矩阵为
${{\boldsymbol{\varSigma }}}_L$ . 将矩阵求逆可得权矩阵${\boldsymbol{P}} = {\sigma ^2}{{\boldsymbol{\varSigma}}}_L^{ - 1}$ ，${\sigma ^2}$ 为先验方差因子，初始值可设置为1. 则加权的最小化残差平方和为$$ {{\boldsymbol{e}}^{\rm{T}}}{\boldsymbol{P}}{\boldsymbol{e}} = {(\Delta {\boldsymbol{L}} - {\boldsymbol{H}}\Delta {\boldsymbol{X}})^{\rm{T}}}{\boldsymbol{P}}(\Delta {\boldsymbol{L}} - {\boldsymbol{H}}\Delta {\boldsymbol{X}}), $$ (8) 使式(8)最小化可得法方程：
$$ {{\boldsymbol{H}}^{\rm{T}}}{\boldsymbol{PH}}\Delta {\boldsymbol{X}} = {{\boldsymbol{H}}^{\rm{T}}}{\boldsymbol{P}}\Delta {\boldsymbol{L}} ,$$ (9) 求解式（9）可得参数修正的估计值
$\Delta \hat {\boldsymbol{X}} = {({{\boldsymbol{H}}^{\rm{T}}}{\boldsymbol{PH}})^{ - 1}} {{\boldsymbol{H}}^{\rm{T}}}{\boldsymbol{P}}\Delta {\boldsymbol{L}}$ ，则未知参数的估计值为$$ \hat {\boldsymbol{X}} = {{\boldsymbol{X}}_0} + \Delta \hat {\boldsymbol{X}}. $$ (10) 在求解非线性最小二乘平差时，一般通过不断迭代提高参数估计精度，最终获得参数估计解.
2. 处理策略
2.1 K-means算法测站选取策略
将数据样本按照数据之间的内在关系划分为若干类别的算法称为聚类算法，K-means是其中最基础的聚类算法[9]，是一种将向量距离作为相似性评价指标，进行类簇划分的方法. 向量距离即为每一个样本到每一个聚类中心的距离，一般采用欧式距离[10]，公式如下：
$$ D({X_i},{C_j}) = \sqrt {{{({x_{i,{\rm{lon}}}} - {c_{j,{\rm{lon}}}})}^2} + {{({x_{j,{\rm{lat}}}} - {c_{j,{\rm{lat}}}})}^2}}. $$ (11) 式中：
$D({X_i},{C_j})$ 为第$i$ 个样本到第$j$ 个簇类中心的距离；${x_{i,{\rm{lon}}}}$ ，${c_{j,{\rm{lon}}}}$ 和${x_{i,{\rm{lat}}}}$ ，${c_{j,{\rm{lat}}}}$ 分别为第$i$ 个样本和第$j$ 个簇类中心的相关属性.K值代表其需要划分的聚类个数，means代表取聚类中所有样本在各个属性维度的均值作为该簇的中心[10]，公式如下：
$$ {C_j} = \frac{1}{{{n_j}}}\sum\limits_{X_i \in {S_j}} {{X_i}}. $$ (12) 式中：
${C_j}$ 为第$j$ 个聚类的中心；${n_j}$ 为第$j$ 个聚类的样本个数；${S_j}$ 为第$j$ 个聚类；${X_i}$ 为第$j$ 个聚类中第$i$ 个样本.通过计算样本与簇类中心的距离，将样本划分到距离最近的簇类中，再更新簇类中心，不断迭代，最终完成分类.
具体算法流程如下：
1)从样本中随机选取K个对象作为初始簇中心；
2)计算每个样本与K个簇中心的距离；
3)根据距离结果，将所有样本分类到离它最近的簇中心；
4)对重新划分好的聚类重新计算簇中心；
5)重复2~4步骤，直到簇中心不再发生变化.
测站分类计算中，首先获取测站数据，其中包括测站站名和经纬度信息；然后通过K-means算法进行测站分类后，输出结果文件，其中每个簇类中的测站会按照离簇类中心距离由近及远进行排序；最后，按照需要从每个簇类中选取测站，形成测站使用列表.
图1为当K=30，通过K-means算法对全球支持北斗的测站(IGS官网数据)进行簇类分类的结果. 图中圆圈为测站标识，不同颜色表示为不同的簇类，K=30，即分为30个簇类，因此图中共有30种颜色. 从图1可以看出分类效果较好，簇类分布较为均匀.
2.2 误差处理流程
由观测模型可知，GNSS数据处理中涉及卫星端、接收机端、传播路径上等各类误差，需要充分考虑相关误差处理和模型改正. 如信号传播过程中会受到地球引力场引起的时空曲率影响，相对论效应也会对卫星的时钟标称频率产生影响[7]，文献[11]提供了相关改正公式. 同样，天线相位中心偏差(PCO)、相位缠绕和潮汐改正(包括地球固体潮、海潮和极潮)等影响也需要通过相关模型进行改正[1].
对流层延迟改正模型函数为[12]
$$ \begin{gathered} T_{\rm{r}}^{\rm{s}} = {m_{\rm{h}}}(e){D_{{\rm{zh}}}} + {m_{\rm{w}}}(e){D_{{\rm{zw}}}} + {m_{\rm{g}}}(e) ({G_{\rm{N}}}\cos\;a + {G_{\rm{E}}}\sin\;a) . \end{gathered} $$ (13) 式中：
${m_{\rm{h}}}$ 、${m_{\rm{w}}}$ 和${m_{\rm{g}}}$ 分别为静力学、湿延迟和水平延迟梯度的映射函数；${D_{{\rm{zh}}}}$ 为天顶静力学延迟(ZHD)；${D_{{\rm{zw}}}}$ 为天顶湿延迟(ZWD)[13]；$e$ 为俯仰角；$a$ 为方位角；${G_{\rm{N}}}$ 和${G_{\rm{E}}}$ 分别为北、东分量的水平延迟梯度. 在GNSS处理中，ZHD建模较充分，用模型给定初值即可；而ZWD模型改正不足，需要对各个测站的剩余ZWD参数进行估计[14].GNSS信号受电离层影响，会产生测距误差，电离层延迟产生的测距误差是卫星导航中不可忽略的 重要误差[15-16]. 电离层延迟误差一般包括一阶项和其他高阶项误差，文献[17]对电离层改正进行了详细分析，并提出了高阶项误差的改正公式. IF组合模型可以通过不同频点观测数据的线性组合消除一阶电离层误差[1]；非差非组合不对观测值进行线性组合，通过最小二乘平差，对单个历元卫星和接收机所有伪距和相位观测的倾斜总电子含量(STEC)参数进行估计. 而电离层函数的各项误差，均可表示为与STEC相关的函数[7]. 因此，通过第一次参数估计得到的STEC参数值在随后的迭代过程中可以直接参与电离层改正. 同时，估计的电离层参数还可以为差分码偏差（DCB）估计提供观测量[1].
由于相位偏差的存在，模糊度参数不具备整数特性. 本文通过设置相位偏差和整周模糊度参数[7]，估计得到整周模糊度的浮点解，然后通过Lambda方法[18]及改进的Lambda方法[19]进行整周模糊度的固定.
2.3 数据处理流程
通过变分方程和相关模型对轨道进行预处理，地球重力场模型采用GOCO06s模型，其他模型见表1，轨道弧段为24 h，采样率为60 s. 通过估计初始状态和ECOM模型(ECOM)参数，对轨道进行拟合处理.
表 1 轨道积分的力模型项目 模型 地球重力场 GOCO06s N体引力 JPL DE432 行星星历 地球固体潮 IERS Conventions 2010 海潮 FES2014b 全球海潮模型 极潮 IERS Conventions 2010 大气潮汐 AOD1B RL06 大气和海洋质量变化 AOD1B RL06 相对论效应 IERS Conventions 2010 太阳辐射压 Box-wing 模型 地球辐射压 Box-wing 模型 观测数据通过MW组合和GF组合进行周跳探测，使用方差分量估计确定观测值权重.
测站坐标因地球非刚体形变而发生位移，因此，需要对其进行改正. 根据国际地球自转服务(IERS)对地球固体潮、极潮建模，海潮通过全球海潮模型[20]建立. 大气和海洋去混叠1B级(AOD1B)06版模型用于大气潮汐、大气和海洋质量变化建模[21].
最后，使用采样率为30 s的IGS观测站观测数据建立原始观测方程，使用参数估计方法进行迭代计算，最终获得收敛后的定轨结果.
3. 实验与分析
为分析BDS-3多频非差非组合精密定轨精度，验证K-means算法测站选取的优势，实验选取可发射五频信号的BDS-3 MEO和IGSO共27颗卫星，测站选取方案为：
方案1：根据人工经验选取测站（即通过以往进行精密定轨实验的经验，选取所需数量观测数据质量较好的测站）；
方案2：通过K-means算法选取测站.
频率选择方案为：
1) S1：B1I+B3I；
2) S2：B1C+B2a；
3) S3：B1I+B1C+B2a+B3I+B2b.
图2和图3分别为测站选取方案1和方案2使用的测站分布图，图中测站用红色三角表示.
实验定轨时段选取2022年年积日第060—090共一个月的时间，将定轨结果与武汉大学IGS分析中心的精密轨道产品进行对比.
3.1 多频A、C、R三个方向RMS月均值对比分析
通过统计实验数据，首先计算单个年积日中每颗卫星的A、C和R方向RMS值，然后再计算每个年积日A、C和R方向所有卫星的RMS均值，最终形成轨道精度对比图.
图4为方案1中每个年积日A、C、R方向所有卫星的RMS均值对比图. 从图中可以看出，在多数时间内，S2定轨精度要低于其他两个频率选择方案.通过进一步统计计算A、C 、R方向RMS 月平均值，可知， S1定轨结果中A、C、R方向RMS月均值分别为0.052 m、0.044 m、0.044 m；S2定轨结果中A、C、R方向RMS月均值分别为0.055 m、0.048 m、0.047 m；S3定轨结果中A、C、R方向RMS月均值分别为0.052 m、0.044 m、0.044 m. 从统计结果可以看出，S1和S3 A、C、R方向RMS月均值相较于S2，精度均提升约0.003 m、0.004 m和0.003 m.
图5为方案2中每天在A、C、R方向上所有卫星RMS月均值对比图.同样，通过进一步统计计算三个方向RMS 月平均值，可知， 方案2的S1定轨结果中A、C、R方向RMS月均值分别为0.045 m、0.038 m、0.039 m；S2定轨结果中A、C、R方向RMS月均值分别为0.047 m、0.040 m、0.039 m；S3定轨结果中A、C、R方向RMS月均值分别为0.046 m、0.039 m、0.040 m. 可以看出，三种频率选择方案的定轨精度基本一致，而相较于方案1，定轨精度整体都有提升.
3.2 3D RMS月均值对比分析
通过计算得到每颗卫星每天的3D RMS值，然后计算月均值并进行对比.
图6为三种频率选择下两种测站选取方案中BDS-3 MEO和IGSO卫星的3D RMS月均值对比. 可以看出，两种方案的MEO和IGSO定轨精度中，方案2的3D RMS月均值除S2中C38、C39和S3中C38卫星外，其他均低于方案1，很明显，方案2的定轨精度更高.
表2和表3详细统计了两种方案中每颗卫星不同频率选择的3D RMS月均值和所有卫星的总体均值. 从表2可以看出，方案1定轨结果中MEO卫星3D RMS月均值最小值为0.057 m、最大值为0.100 m，IGSO卫星3D RMS月均值最小值为0.127 m、最大值为0.219 m. 从表3可以看出，方案2定轨结果中MEO卫星3D RMS月均值最小值为0.045 m、最大值为0.093 m，IGSO卫星的3D RMS月均值最小值为0.118 m、最大值为0.213 m. 可以看出MEO和IGSO卫星定轨结果中方案2的定轨精度相比于方案1有一定提升.
表 2 方案1 不同频率选择方案的3D RMS月均值m PRN 卫星
类型3D RMS月均值 PRN 卫星
类型3D RMS月均值 S1 S2 S3 S1 S2 S3 C19 MEO 0.061 0.071 0.061 C34 MEO 0.059 0.084 0.063 C20 0.061 0.075 0.061 C36 0.063 0.073 0.067 C21 0.057 0.067 0.057 C37 0.062 0.066 0.061 C22 0.064 0.086 0.064 C41 0.060 0.066 0.057 C23 0.064 0.067 0.066 C42 0.062 0.068 0.061 C24 0.063 0.072 0.067 C43 0.098 0.086 0.100 C25 0.089 0.092 0.088 C44 0.079 0.082 0.082 C26 0.097 0.098 0.095 C45 0.073 0.084 0.071 C27 0.066 0.076 0.068 C46 0.071 0.070 0.074 C28 0.069 0.076 0.070 总体
均值0.070 0.077 0.070 C29 0.070 0.089 0.074 C38 IGSO 0.219 0.188 0.174 C30 0.072 0.079 0.071 C39 0.207 0.194 0.203 C32 0.081 0.076 0.065 C40 0.127 0.129 0.166 C33 0.063 0.076 0.066 总体
均值0.184 0.170 0.181 C35 0.067 0.076 0.065 表 3 方案2 不同频率选择方案的3D RMS月均值m PRN 卫星
类型3D RMS月均值 PRN 卫星
类型3D RMS月均值 S1 S2 S3 S1 S2 S3 C19 MEO 0.055 0.054 0.055 C35 MEO 0.058 0.053 0.054 C20 0.053 0.052 0.053 C36 0.054 0.052 0.054 C21 0.045 0.045 0.046 C37 0.057 0.050 0.052 C22 0.055 0.054 0.056 C41 0.050 0.049 0.050 C23 0.058 0.055 0.057 C42 0.051 0.051 0.052 C24 0.061 0.061 0.063 C43 0.086 0.085 0.065 C25 0.087 0.085 0.086 C44 0.067 0.064 0.080 C26 0.093 0.089 0.091 C45 0.061 0.058 0.059 C27 0.065 0.063 0.065 C46 0.059 0.056 0.058 C28 0.063 0.065 0.062 总体
均值0.061 0.060 0.061 C29 0.067 0.063 0.065 C38 IGSO 0.213 0.204 0.208 C30 0.063 0.061 0.062 C39 0.201 0.196 0.200 C32 0.057 0.055 0.056 C40 0.127 0.118 0.126 C33 0.054 0.056 0.057 总体
均值0.180 0.173 0.178 C34 0.055 0.052 0.058 通过表2和表3对比MEO和IGSO所有卫星的3D RMS月均值可以看出：方案1中S2大部分卫星3D RMS月均值相较于S1和S3精度均偏高，其中MEO卫星总体均值为0.077 m；而其他两种频率选择方案的月均值均为0.070 m，相差0.007 m. 非差非组合方法通过纳入更多观测数据使精度提升约9.1%. 与A、C、R方向RMS值对比结果一致，选用非差非组合方法利用更多频点观测数据进行精密定轨，能够提高定轨精度.
三种频率选择下，MEO 3D RMS方案1均值分别为0.070 m、0.077 m、0.070 m，方案2 均值分别为0.061 m、0.060 m、0.061 m，IGSO 3D RMS方案1均值分别为0.184 m、0.170 m、0.181 m，方案2均值分别为0.180 m、0.173 m、0.178 m. 通过对比可以发现，本文提出的K-means算法测站选取策略，所得定轨精度与人工经验选取方法相比，MEO卫星在S1、S2和S3频率选择下，精度分别提升0.009 m、0.017 m和0.009 m，而IGSO卫星的定轨精度基本一致，各有优劣. 分析认为，由于IGSO卫星轨道的特殊性，覆盖范围为区域覆盖，其他区域的测站不能接收到IGSO的观测信号，虽然方案2的测站选取相较于方案1有明显的优化，但是对于IGSO卫星的测站布局并没有改善，两种方案对比中接收IGSO卫星信号的测站数量和布局没有明显优势，因此其定轨精度相当. 可知，优化全球测站分布的方法对IGSO卫星精度提升不明显. 因其区域性的特点，需要采用中国境内或区域跟踪站优化测站布局改善IGSO卫星定轨精度[22].
从表3中的MEO和IGSO卫星3D RMS月均值可以看出，S2的定轨精度均高于其他两种频率选择方案，3D RMS月均值提升约1 mm. 初步推断，在测站分布较为合理，观测数据充足的情况下，优先选择新体制信号进行精密定轨即能得到较好的定轨结果.
3.3 非差非组合精密定轨优势讨论
表4为不同测站选择方案中，各频率选择方案测站使用数量统计，其中，方案1和方案2可用de 测站数量均为29个，支持B1C+B2a信号的测站分别为20和27个.
表 4 测站使用数量统计测站选取方案 频率选择方案 S1 S2 S3 方案1 29 20 29 方案2 29 27 29 从表4可以看出，方案1中，接收S2频点的测站相比S1和S3使用数量较少，导致测站分布几何构型较差，对比定轨结果，S2定轨精度较低，而S1和S3定轨精度基本一致；方案2中，S2使用测站数量相比S1和S3减少两个，三种方案定轨精度基本一致. 分析认为，方案1中，在不确定所选测站是否支持所选频率的情况下，非差非组合方法纳入了更多频点观测数据，提高了测站使用数量，获得了更优化的测站几何布局，提高了定轨精度. 对比方案1中S1和S3以及方案2的定轨结果，当测站数量相同，测站分布一致、几何构型合理的情况下，纳入更多频点观测数据对定轨精度提升不明显. 通过对比同种频率选择、不同测站选取方案的定轨精度，方案2定轨精度整体优于方案1，因此，测站几何构型相比测站数量和观测数据数量，对定轨精度的提升贡献更大.
从以上对比分析可以看出，影响定轨精度的主要因素为测站数量和测站几何观测信息，而测站几何观测信息对定轨精度影响最大. 非差非组合方法使用更多频点观测数据，通过提高测站使用数量、优化测站几何观测信息，从而提高定轨精度.
非差非组合方法相较于IF组合方法的定轨精度对比中，文献[23]对Galileo双频、三频IF组合和非差非组合的定轨精度进行了评估，两种定轨精度差异在1 mm内；文献[1]示例性地给出了部分GPS和Galileo卫星两种定轨模型的定轨结果，结果显示定轨精度基本相当；文献[8]评定了GPS、Galileo和BDS-2的三频观测量对定轨结果的贡献，结果表明第3频点观测量对精密定轨贡献很小，IF组合和非差非组合的定轨精度结果基本一致.
在两种定轨模型结果基本一致的情况下，IF组合通过不同频点观测值线性组合，消除电离层一阶项误差，会丢失观测值中的一些可用信息，并且放大了观测值的噪声[1]. 与IF组合相比，虽然非差非组合方法增加了待估参数数量，但是在以下四方面具有明显优势：1)模型扩展性强，更加容易纳入多频点观测数据，充分利用BDS-3多频资源；2)对电离层参数进行估计，获得电离层信息，可以进行电离层建模与反演研究；3)根据电离层信息，估计DCB，进一步估计观测量的码偏差（OSB）值，提供丰富的多频相关产品；4)避免IF组合中放大观测噪声的问题.
通过实验分析，引入更多频点观测数据，虽然定轨精度提高存在局限性，但是对于OSB估计、电离层建模等问题具有优势，在计算高精度轨道产品的同时，能够为用户提供更加丰富的多频点OSB等产品. 因此，在多频数据时代，非差非组合模型优势更加明显.
图7显示了通过非差非组合方法利用五频观测数据，估计得到的BDS-3 2022年年积日第060—091伪距OSB序列. 伪距OSB估计与DCB相关，因此其特性与DCB一致. 可以看出，BDS-3伪距OSB序列十分稳定，所有码类型中，C33卫星伪距OSB值较大. 同一频率、不同码类型的伪距OSB数值比较接近，如C1D和C1P，而不同频率观测值的伪距OSB相差较大.
4. 结 论
本文介绍了非差非组合方法的观测模型和参数估计方法，提出了通过K-means算法进行测站选取方法策略. 在此基础上，通过不同测站选取方案和频率选择方案，对BDS-3 MEO和IGSO共27颗卫星进行了精密定轨处理，通过以武汉大学IGS分析中心精密轨道产品为基准，以A、C、R三方向RMS月均值和3D RMS月均值为指标，对定轨精度进行了分析，得出如下结论：
1)选取测站进行精密定轨时，可能存在测站不接收所选频点观测数据的情况，此时，非差非组合方法可以通过纳入更多频点观测数据提高定轨精度.
2) K-means算法测站选取策略，能够提供分布更加均匀的测站，可以得到更好的定轨结果. 当观测数据达到一定数量、测站分布又相同时，加入更多频点的观测数据对定轨精度提升不明显.
3)非差非组合方法灵活性高、可扩展性强，可以充分利用多频观测数据，在电离层特性研究、全频点OSB估计等方面具有明显优势.
致谢：感谢iGMAS分析中心(LSN)提供的硬件支持和技术帮助. 感谢IGS分析中心提供的数据支持. 感谢Graz University of Technology提供的软件支持.
-
表 1 轨道积分的力模型
项目 模型 地球重力场 GOCO06s N体引力 JPL DE432 行星星历 地球固体潮 IERS Conventions 2010 海潮 FES2014b 全球海潮模型 极潮 IERS Conventions 2010 大气潮汐 AOD1B RL06 大气和海洋质量变化 AOD1B RL06 相对论效应 IERS Conventions 2010 太阳辐射压 Box-wing 模型 地球辐射压 Box-wing 模型 
下载: 导出CSV
表 2 方案1 不同频率选择方案的3D RMS月均值
m PRN 卫星
类型3D RMS月均值 PRN 卫星
类型3D RMS月均值 S1 S2 S3 S1 S2 S3 C19 MEO 0.061 0.071 0.061 C34 MEO 0.059 0.084 0.063 C20 0.061 0.075 0.061 C36 0.063 0.073 0.067 C21 0.057 0.067 0.057 C37 0.062 0.066 0.061 C22 0.064 0.086 0.064 C41 0.060 0.066 0.057 C23 0.064 0.067 0.066 C42 0.062 0.068 0.061 C24 0.063 0.072 0.067 C43 0.098 0.086 0.100 C25 0.089 0.092 0.088 C44 0.079 0.082 0.082 C26 0.097 0.098 0.095 C45 0.073 0.084 0.071 C27 0.066 0.076 0.068 C46 0.071 0.070 0.074 C28 0.069 0.076 0.070 总体
均值0.070 0.077 0.070 C29 0.070 0.089 0.074 C38 IGSO 0.219 0.188 0.174 C30 0.072 0.079 0.071 C39 0.207 0.194 0.203 C32 0.081 0.076 0.065 C40 0.127 0.129 0.166 C33 0.063 0.076 0.066 总体
均值0.184 0.170 0.181 C35 0.067 0.076 0.065
下载: 导出CSV
表 3 方案2 不同频率选择方案的3D RMS月均值
m PRN 卫星
类型3D RMS月均值 PRN 卫星
类型3D RMS月均值 S1 S2 S3 S1 S2 S3 C19 MEO 0.055 0.054 0.055 C35 MEO 0.058 0.053 0.054 C20 0.053 0.052 0.053 C36 0.054 0.052 0.054 C21 0.045 0.045 0.046 C37 0.057 0.050 0.052 C22 0.055 0.054 0.056 C41 0.050 0.049 0.050 C23 0.058 0.055 0.057 C42 0.051 0.051 0.052 C24 0.061 0.061 0.063 C43 0.086 0.085 0.065 C25 0.087 0.085 0.086 C44 0.067 0.064 0.080 C26 0.093 0.089 0.091 C45 0.061 0.058 0.059 C27 0.065 0.063 0.065 C46 0.059 0.056 0.058 C28 0.063 0.065 0.062 总体
均值0.061 0.060 0.061 C29 0.067 0.063 0.065 C38 IGSO 0.213 0.204 0.208 C30 0.063 0.061 0.062 C39 0.201 0.196 0.200 C32 0.057 0.055 0.056 C40 0.127 0.118 0.126 C33 0.054 0.056 0.057 总体
均值0.180 0.173 0.178 C34 0.055 0.052 0.058
下载: 导出CSV
-
[1] 李星星, 张伟, 袁勇强, 等. GNSS卫星精密定轨综述: 现状, 挑战与机遇[J]. 测绘学报, 2022, 51(7): 1271-1293. DOI: 10.11947/j.issn.1001-1595.2022.7.chxb202207018 [2] 辜声峰. 多频GNSS非差非组合精密数据处理理论及其应用[D]. 武汉: 武汉大学, 2013. [3] 周锋. 多系统GNSS非差非组合精密单点定位相关理论和方法研究[D]. 上海华东师范大学, 2018. [4] 景鑫, 车通宇, 周舒涵, 等. 北斗三号系统开通前后广播星历精度对比分析[J]. 全球定位系统, 2022, 47(4): 23-30. DOI: 10.12265/j.gnss.2022051 [5] 张小红, 柳根, 郭斐, 等. 北斗三频精密单点定位模型比较及定位性能分析[J]. 武汉大学学报(信息科学版), 2018, 43(12): 2124-2130. [6] ZEHENTNER N, MAYER-GURR T. Precise orbit determination based on raw GPS measurements[J]. Journal of geodesy, 2015, 90(3): 275-286. DOI: 10.1007/s00190-015-0872-7
[7] STRASSER S, MAYER-GURR T, ZEHENTNER N. Processing of GNSS constellations and ground station networks using the raw observation approach[J]. Journal of Geodesy, 2018, 93(8): 1045-1057. DOI: 10.1007/s00190-018-1223-2
[8] 曾添, 隋立芬, 阮仁桂, 等. 三频非组合模型的GPS/BDS/Galileo精密定轨[J]. 测绘学报, 2021, 50(2): 169-180. [9] 彭欢欢, 石教坤, 蔡鑫. 基于单通道的Kmeans算法的图像分割[J]. 信息技术与信息化, 2022(1): 39-42. DOI: 10.3969/j.issn.1672-9528.2022.01.011 [10] 彭海驹, 严科文, 林松, 等. 融合kmeans聚类与Hausdorff距离的点云精简算法改进[J]. 地理空间信息, 2022, 20(8): 59-63. DOI: 10.3969/j.issn.1672-4623.2022.08.014 [11] HOFMANN-WELLENHOF B, LICHTENEGGER H, WASLE E. GNSS-global navigation satellite systems: GPS, GLONASS, Galileo, and more[M]. Wien; New York: Springer, 2007.
[12] PETIT G, LUZUM B. IERS Conventions(2010), IERS Technical Note 36[R].
[13] 谢劭峰, 潘清莹, 黄良珂, 等. 中国区域ZTD, ZWD高程缩放因子的时空特性分析[J]. 大地测量与地球动力学, 2021, 41(12): 1211-1215. [14] 章迪. GNSS对流层天顶延迟模型及映射函数研究[D]. 武汉: 武汉大学, 2017. [15] 李涌涛, 李建文, 庞鹏, 等. 2017年9月磁暴期间电离层TEC变化分析[J]. 全球定位系统, 2018, 43(4): 42-47. [16] 李涌涛, 李建文, 顾晨钟, 等. Linux Shell在电离层TEC格网数据提取和分析中的应用[J]. 全球定位系统, 2019, 44(3): 81-87. [17] HOQUE M M, JAKOWSKI N. Estimate of higher order ionospheric errors in GNSS positioning[J]. Radio science, 2008, 43(5): 1-15. DOI: 10.1029/2007RS003817
[18] TEUNISSEN P J G. The least-squares ambiguity decorrelation adjustment: a method for fast GPS integer ambiguity estimation[J]. Journal of geodesy, 1995, (70): 65-82. DOI: 10.1007/BF00863419
[19] CHANG X-W, YANG X, ZHOU T. MLAMBDA: a modified LAMBDA method for integer least-squares estimation[J]. Journal of geodesy, 2005, 79(9): 552-565. DOI: 10.1007/s00190-005-0004-x
[20] CARRERE L, LYARD F, CANCET M, et al. FES 2014, a new tidal model on the global ocean with enhanced accuracy in shallow seas and in the Arctic region[A]. Egu General Assembly Conference[C]. 2015.
[21] DOBSLAW H, BERGMANN-WOLF I, DILL R, et al. A new high-resolution model of non-tidal atmosphere and ocean mass variability for de-aliasing of satellite gravity observations: AOD1B RL06[J]. Geophysical journal international, 2017, 211(1): 263-269. DOI: 10.1093/gji/ggx302
[22] 张龙平, 党亚民, 成英燕, 等. 北斗GEO/IGSO/MEO卫星定轨地面站构型影响分析及其优化[J]. 测绘学报, 2016, 45(S2): 82-92. [23] 李星星, 黄健德, 袁勇强, 等. Galileo三频非组合精密定轨模型及精度评估[J]. 测绘学报, 2020, 49(9): 1120-1130. DOI: 10.11947/j.AGCS.2020.20200320 -
期刊类型引用(1)
1. 张岩，周康，袁勇强. 基于RINEX4.0广播ERP参数的BDS星历误差实时补偿方法. 全球定位系统. 2024(06): 47-54 . 
本站查看
其他类型引用(0)
计量
- 文章访问数: 290
- HTML全文浏览量: 77
- PDF下载量: 19
- 被引次数: 1
