-
现今全球导航卫星系统(global navigation satellite system，GNSS)技术得到空前发展，各国都在致力于建设自主的导航系统，主要包括美国的GPS、中国的北斗卫星导航系统(BeiDou navigation satellite system，BDS)、俄罗斯的GLONASS、欧盟的Galileo，以及一些区域的导航系统，如日本的天顶卫星系统和印度的区域导航卫星系统[1]。拥有中国自主产权的BDS已于2020年完成全球组网工作。截至2019年9月，中国共发射了48颗BDS卫星，其中第二代BDS导航系统由5颗地球静止轨道(geostationary earth orbit，GEO)卫星、6颗倾斜地球同步轨道(inclined geosynchronous satellite orbit，IGSO)卫星和3颗中轨道(medium earth orbit，MEO)卫星组成，并免费为亚太地区用户提供服务。
信号在产生到发出或在接收到转变的时间延迟称为仪器硬件延迟，差分码偏差指的是不同频率或者不同观测量间的延迟差[2]，分为频内偏差与频间偏差(differential code bias，DCB)，前者是同一频率不同观测量间的差分码偏差，后者是不同频率相同观测量间的差分码偏差[2-4]。频内偏差的估计方法比较成熟且简单，本文只讨论DCB的估计方法。目前确定DCB的策略分为两大类: (1)在设备出厂前，厂家利用仪器对接收机或者卫星信道进行检测并固定DCB[5]; (2)利用GNSS软件估计的方法对仪器DCB进行事后估计[6]。由于设备DCB易受周边环境以及设备内部配置变化的影响，故第二种策略精度较高，是科研工作中的常用方法。由于DCB与电离层观测耦合，所以利用GNSS软件估计DCB又可分为两种方式: (1)在对电离层建模的过程中，把DCB设为一天时间跨度内的常量参数，随电离层模型参数同步估计[4, 7-9]; (2)先使用经验的电离层产品修正电离层误差，以提取“纯净”的DCB参数[10]。作为信号干扰最大的误差源之一，电离层时延一直是导航定位服务关注的研究热点之一。国际卫星导航服务(International GNSS Service，IGS)先后成立了7个电离层分析中心，旨在提高电离层及相关服务的精度和时效性。其中欧洲定轨中心(Center for Orbit Determination in Europe，CODE)、美国喷气推进实验室(Jet Propulsion Laboratory，JPL)和欧洲空间操作中心(European Space Operations Center of European Space Agency，ESOC)采用了第一种DCB估计方法，加泰罗尼亚理工大学(Universidad Politécnica de Cataluña，UPC)和德国宇航中心(Deutsches Zentrum für Luft- und Raumfahrt，DLR)采用第二种DCB估计方法，且DLR以Bias-SINEX格式提供了多系统的频间偏差[2]。已有研究表明，卫星和接收机DCB稳定性较高，年内变化量级达十分之一纳秒级[11]，因此，国际上在估计DCB时，一般将其在一天内表示为一个常量。为了减少电离层模型对穿刺点分布的苛刻要求，并获得精度更高的DCB产品，李子申等[2, 12]采用广义三角级数函数对电离层逐基站建模，并兼顾卫星实际的稳定性，提出一种“拟稳”基准分离卫星和接收机DCB，具有较高的估计精度。为了评估不同DCB产品的精度，曾添等[13-14]根据不同类型DCB间的代数关系，推导了不同形式的单点定位DCB改正形式，结果表明，不同的表达形式对单点定位的影响在厘米量级，其中对精密单点定位(precise point positioning，PPP)影响较大。
在使用DCB产品时，一般是通过DCBC2I-C7I、DCBC2I-C6I作差来确定DCBC7I-C6I，这将降低DCBC7I-C6I的产品精度。本文首先将3种频间偏差闭合差约束引入DCB观测方程; 然后，在此基础上选取北斗MEO和IGSO卫星作为参考卫星，解决卫星和接收机DCB间的共线性问题; 最后，采用附加限制条件的间接平差方法同步估计北斗二代卫星的3种DCB，提高DCB的估计精度。
-
载波平滑伪距技术是如今最为普适的利用GNSS数据提取电离层原始信息的技术。为不失一般性，本文采用载波平滑伪距技术提取精度可靠的电离层原始观测量[15]。理论上，在BDS二代发射的B1、B2和B3信号间存在3种频间偏差，且其代数和为零。多模GNSS实验(multi-GNSS experiment，MGEX)工作组在2015年之后只发布两种BDS频间偏差产品，第三种使用发布的两种DCB产品线性表示。本文将这3种频间偏差的代数和为零作为约束附加到观测方程，同步估计BDS多种频率组合的频间偏差。
-
BDS观测文件中包含有3种频点的伪距和载波观测值，其观测方程为:
$$ \left\{\begin{array}{l}{O}_{f, r}^{s}={\rho }_{r}^{s}+c{d}_{{t}_{r}}-c{d}_{{t}^{s}}+{I}_{f, r}^{s}+{T}_{r}^{s}+c{d}_{f}^{s}+\\ c{d}_{f, r}+{\epsilon }_{P, f, r}^{s}\\ {L}_{f, r}^{s}={\rho }_{r}^{s}+c{d}_{{t}_{r}}-c{d}_{{t}^{s}}-{I}_{f, r}^{s}+{T}_{r}^{s}-{\lambda }_{f}({N}_{f, r}^{s}+\\ {b}_{f, r}^{s})+{\epsilon }_{L, f, r}^{s}\end{array}\right. $$ (1) 式中，上标$ s $、下标$ r $和$ f $分别表示特定的卫星、接收机和卫星信号频段; $ O $和$ L $分别表示北斗卫星导航系统伪距和载波观测值; $ \rho $表示接收机与北斗卫星间的真实几何距离; $ c $表示真空中的光速; $ {d}_{{t}_{r}} $和$ {d}_{{t}^{s}} $分别表示接收机和卫星钟差; $ {I}_{f, r}^{s} $和$ {T}_{r}^{s} $分别表示电离层和对流层时延等效距离误差，其中前者与信号频率相关; $ {\lambda }_{f} $表示观测信号的波长; $ {N}_{f, r}^{s} $表示载波相位观测的整周模糊度，是制约PPP精度和收敛时间的主要因素之一; $ {b}_{f, r}^{s} $表示接收机端与卫星端的相位硬件延迟之和; $ {d}_{f}^{s} $和$ {d}_{f, r} $分别表示卫星端与接收机端的码偏差; $ {\epsilon }_{L, f, r}^{s} $和$ {\epsilon }_{P, f, r}^{s} $分别为相位和伪距观测误差。
根据北斗二代卫星发射的3种频率的信号，消除与频率无关项可组成两种无几何组合观测量，计算式为:
$$ \left\{\begin{array}{l}{O}_{4, i, j}={O}_{i, r}^{s}-{O}_{j, r}^{s}=({I}_{i, r}^{s}-{I}_{j, r}^{s})+c{D}_{i, j}^{s}+c{D}_{i, j, r}\\ {L}_{4, i, j}={L}_{i, r}^{s}-{L}_{j, r}^{s}=-({I}_{i, r}^{s}-{I}_{j, r}^{s})-{\lambda }_{i, j}({b}_{i, r}^{s}-{b}_{j, r}^{s})-\\ {\lambda }_{i, j}({N}_{i, r}^{s}-{N}_{j, r}^{s})\end{array}\right. $$ (2) 式中，$ {O}_{4, i, j} $和$ {L}_{4, i, j} $分别表示北斗卫星伪距和载波无几何观测量; $ {D}_{i, j}^{s} $和$ {D}_{i, j, r} $分别表示北斗卫星端和接收机端的频间偏差，$ {D}_{i, j}^{s}={d}_{i}^{s}-{d}_{j}^{s} $，$ {D}_{i, j, r}={d}_{i, r}-{d}_{j, r}(i=\mathrm{1, 2};j=\mathrm{2, 3};i\ne j) $。
从式(2)可以看出，伪距无几何观测量可以作为电离层的直接观测量，但噪声较大; 载波无几何观测量虽然精度较高，但很难直接解算组合观测的模糊度[16]。所以，本文假设在连续弧段内，模糊度和DCB的值是固定的，采用载波平滑伪距技术来平滑伪距观测量，计算式为:
$$ \left\{\begin{array}{l}{O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}={\omega }_{t}{O}_{4, i, j}\left(t\right)+(1-{\omega }_{t}){O}_{4, i, {j}_{\mathrm{p}\mathrm{r}\mathrm{d}}}\left(t\right)\\ {O}_{4, i, {j}_{\mathrm{p}\mathrm{r}\mathrm{d}}}\left(t\right)={O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}(t-1)+\left[{L}_{4, i, j}\right(t)-\\ {L}_{4, i, j}(t-1)], t>1\end{array}\right. $$ (3) 式中，$ {O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}} $表示平滑后的伪距无几何观测量; $ {O}_{4, i, {j}_{\mathrm{p}\mathrm{r}\mathrm{d}}} $表示载波和伪距的组合观测量; $ t $表示观测时间; $ {\omega }_{t} $表示与时间相关的权重因子，$ {\omega }_{t}=1/t $。当$ t=1 $时，$ {O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}={O}_{4, i, j} $。
在平滑之前，必须探测并剔除含有周跳的历元，本文采用TurboEdit方法探测周跳，并剔除发生周跳的历元[17]。根据电离层对无线电信号的作用机理，电离层时延等效距离误差的计算式为:
$$ {I}_{f, r}^{s}=\frac{40.3}{{{f}_{i}}^{2}}{S}_{\mathrm{T}\mathrm{E}\mathrm{C}} $$ (4) 式中，$ f $表示信号频率; $ {S}_{\mathrm{T}\mathrm{E}\mathrm{C}} $表示沿传播路径底面积为1的圆柱体中所含电子总量。将式(4)代入式(3)，可得:
$$ {O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}=40.3(\frac{1}{{f}_{i}^{2}}-\frac{1}{{f}_{j}^{2}}){S}_{\mathrm{T}\mathrm{E}\mathrm{C}}+c{D}_{i, j}^{s}+c{D}_{i, j, r} $$ (5) -
根据式(5)$ \mathrm{可}\mathrm{得}{S}_{\mathrm{T}\mathrm{E}\mathrm{C}} $计算式为:
$$ {S}_{\mathrm{T}\mathrm{E}\mathrm{C}}=-\frac{{f}_{i}^{2}{f}_{j}^{2}}{40.3({f}_{i}^{2}-{f}_{j}^{2})}({O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}-c{D}_{i, j}^{s}-c{D}_{i, j, r}) $$ (6) 为了简化电离层建模问题，常常假设电离层电子聚集在地球上空60~1 000 km的范围内的一个无限薄的薄层上，即薄层假设。视线与薄层的交点称为穿刺点(ionospheric pierce point，IPP)。为了统一电离层模型，常设计一个投影函数将视线方向的电离层电子含量投影到穿刺点天顶方向，计算式为:
$$ \left\{\begin{array}{l}{V}_{\mathrm{T}\mathrm{E}\mathrm{C}}=M\left(z\right){S}_{\mathrm{T}\mathrm{E}\mathrm{C}}\\ M\left(z\right)=\mathrm{c}\mathrm{o}\mathrm{s}\left\{\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}\right[\frac{R}{R+H}\mathrm{s}\mathrm{i}\mathrm{n}\left(\alpha z\right)\left]\right\}\end{array}\right. $$ (7) 式中，$ {V}_{\mathrm{T}\mathrm{E}\mathrm{C}} $表示穿刺点天顶方向的电离层电子含量; $ M\left(z\right) $表示电离层穿刺点处的投影函数; z表示卫星相对于接收机的天顶距，为卫星高度角的余角; R为地球半径; H为电离层薄层所在高度; $ \alpha $为比例因子。本文选取H=506.7 km，$ \alpha $=0.978 2。用球谐函数来表示穿刺点天顶方向的电离层电子含量的时空变化，计算式为:
$$ {V}_{\mathrm{T}\mathrm{E}\mathrm{C}}=\sum\limits_{n=0}^{{n}_{\mathrm{m}\mathrm{a}\mathrm{x}}}\sum \limits_{m=0}^{n}{\tilde{P}}_{nm}\left(\mathrm{s}\mathrm{i}\mathrm{n}\beta \right)\left[{a}_{nm}\mathrm{c}\mathrm{o}\mathrm{s}\right(ms)+\\ {b}_{nm}\mathrm{s}\mathrm{i}\mathrm{n}(ms\left)\right] $$ (8) 式中，$ \beta $表示穿刺点的地心纬度; $ s=\lambda -{\lambda }_{0} $表示穿刺点相对经度; $ \lambda $和$ {\lambda }_{0} $分别表示穿刺点的地心经度和当前太阳时经度; $ {a}_{nm} $和$ {b}_{nm} $分别表示待求球谐模型系数; n和m分别表示球谐函数的度数和阶数; nmax表示球谐模型最大度数; $ {\tilde{P}}_{nm} $表示n度m阶的归化勒让德函数，计算式为:
$$ {\tilde{P}}_{nm}\left(\mathrm{s}\mathrm{i}\mathrm{n}\beta \right)={M}_{C}(n, m){P}_{nm}\left(\mathrm{s}\mathrm{i}\mathrm{n}\beta \right) $$ (9) 式中，$ {P}_{nm}\left(\mathrm{s}\mathrm{i}\mathrm{n}\beta \right) $为勒让德函数; $ {M}_{C}(n, m) $为归化函数，$ {M}_{C}(n, m)=\sqrt{2\times \frac{2n+1}{1+{\delta }_{0m}}\frac{(n-m)!}{(n+m)!}} $; $ {\delta }_{0m} $为Kronecker函数，当m=0时，$ {\delta }_{0m} $=1，否则$ {\delta }_{0m} $=0。
将式(7)和式(8)代入式(6)，可得:
$$ \begin{array}{l}M\left(z\right)\frac{{f}_{i}^{2}{f}_{j}^{2}}{40.3({f}_{i}^{2}-{f}_{j}^{2})}{O}_{4, i, {j}_{\mathrm{s}\mathrm{m}}}=-\sum\limits_{n=0}^{{n}_{\mathrm{m}\mathrm{a}\mathrm{x}}}\sum\limits_{m=0}^{n}{\tilde{P}}_{nm}\left(\mathrm{s}\mathrm{i}\mathrm{n}\beta \right)\left[{a}_{nm}\mathrm{c}\mathrm{o}\mathrm{s}\right(ms)+{b}_{nm}\mathrm{s}\mathrm{i}\mathrm{n}(ms\left)\right]+\\ M\left(z\right)\frac{{f}_{i}^{2}{f}_{j}^{2}}{40.3({f}_{i}^{2}-{f}_{j}^{2})}(c{D}_{i, j}^{s}+c{D}_{i, j, r})\end{array} $$ (10) 式(10)的矩阵形式为:
$$ \left\{\begin{array}{l}\boldsymbol{L}=\boldsymbol{B}\widehat{\boldsymbol{X}}\\ \widehat{\boldsymbol{X}}={\left[{\widehat{\boldsymbol{K}}}^{\mathrm{T}}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{1, 2}, r}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{1, 2}, s}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{1, 3}, r}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{1, 3}, s}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{2, 3}, r}, {{\widehat{\boldsymbol{X}}}^{\mathrm{T}}}_{\mathrm{2, 3}, s}\right]}^{\mathrm{T}}\end{array}\right. $$ (11) 式中，$ \boldsymbol{L} $表示平滑后的无几何组合观测量; $ \boldsymbol{B} $表示设计矩阵; $ \widehat{\boldsymbol{X}} $表示包含有球谐函数参数和3种接收机和北斗卫星DCB的待估参数向量。由式(10)可知，方程列秩亏，秩亏数为1。本文考虑到GEO卫星服役时间较长，稳定性较差，以MEO和IGSO卫星作为参考卫星分离接收机和卫星DCB。加上3种DCB的代数和为零的约束，即$ {\widehat{\boldsymbol{X}}}_{{}_{\mathrm{1, 2}, s}}^{\mathrm{T}}-{\widehat{\boldsymbol{X}}}_{{}_{\mathrm{1, 3}, s}}^{\mathrm{T}}+{\widehat{\boldsymbol{X}}}_{{}_{\mathrm{2, 3}, s}}^{\mathrm{T}}=0 $，构成约束方程的计算式为:
$$ \boldsymbol{S}\widehat{\boldsymbol{X}}=\left[\begin{array}{l}{\boldsymbol{S}}_{1}\\ {\boldsymbol{S}}_{2}\end{array}\right]\widehat{\boldsymbol{X}}=0 $$ (12) 其中，
$$ \begin{array}{l}{\boldsymbol{S}}_{1}=\left[\begin{array}{ccccccc}0& 0& \underset{1\times {u}_{1}}{{\boldsymbol{e}}_{\mathrm{1, 2}}}& 0& 0& 0& 0\\ 0& 0& 0& 0& \underset{1\times {u}_{2}}{{\boldsymbol{e}}_{{}_{\mathrm{1, 3}}}}& 0& 0\\ 0& 0& 0& 0& 0& 0& \underset{1\times {u}_{3}}{{\boldsymbol{e}}_{{}_{\mathrm{2, 3}}}}\end{array}\right]\\ {\boldsymbol{S}}_{2}=\left[\begin{array}{ccccccc}0& 0& \underset{u\times u}{{\boldsymbol{I}}_{\mathrm{1, 2}}}& 0& \underset{u\times u}{{\boldsymbol{I}}_{\mathrm{1, 3}}}& 0& \underset{u\times u}{{\boldsymbol{I}}_{\mathrm{1, 2}}}\end{array}\right]\end{array} $$ 式中，e是元素全为1的行向量; I为单位矩阵; $ {u}_{1}\mathrm{、}{u}_{2}\mathrm{和}{u}_{3} $分别为3种观测量的MEO和IGSO卫星总数; $ u $表示共有3种观测量的卫星数。
综合式(11)和式(12)，根据附有限制条件的间接平差，可得:
$$ \widehat{\boldsymbol{X}}=({\boldsymbol{B}}^{\mathrm{T}}\boldsymbol{R}\boldsymbol{B}+{\boldsymbol{S}}^{\mathrm{T}}{\boldsymbol{S})}^{-1}{\boldsymbol{B}}^{\mathrm{T}}\boldsymbol{R}\boldsymbol{L} $$ (13) 其中，
$$ \left\{\begin{array}{l}\boldsymbol{R}=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}\left(\frac{1}{{\sigma }_{1}^{2}}, \frac{1}{{\sigma }_{2}^{2}}, \frac{1}{{\sigma }_{3}^{2}}\cdots \frac{1}{{\sigma }_{n}^{2}}\right)\\ {\sigma }_{i}^{2}=\frac{1}{N}{\sigma }_{P, i}^{2}+\frac{N-1}{N}{\sigma }_{L, i}^{2}\end{array}\right. $$ 式中，$ \boldsymbol{R} $为观测量权阵; $ {\sigma }_{P, i}^{2} $、$ {\sigma }_{L, i}^{2} $分别为伪距和载波无几何观测量的方差，一般取$ {\sigma }_{L, i}={\sigma }_{{}_{P, i}}/100 $; N表示该卫星弧段的平滑历元数。
-
本文使用Matlab程序工具，把MGEX发布的RINEX观测文件和德国地学中心(Helmholtz-Centre Potsdam-German Research Centre for Geosciences，GFZ)发布的北斗卫星精密星历文件作为输入文件，输出接收机与卫星的DCB和电离层球谐模型的参数。具体流程如图 1所示。
图 1 DCB估计程序流程图
Figure 1. Flowchart of DCB Solver
-
为了验证本文方法的有效性，本文选取了19个能同步接收北斗二代卫星3种频率观测数据的MGEX跟踪站，其中大多分布在中低纬度地区，跟踪站的具体信息见表 1。
表 1 MGEX基准站基本信息
Table 1. Basic Information of MGEX Tracking Stations
基准站分布 MGEX站 天线类型 低纬 DARW ASH700936D_M KAT1、SIN1 LEIAR25.R3 LMMF、SPTU TRM57971.00 中纬 AUCK、BRST、CHTI、MQZG、REUN、UCAL、WARK、WGTN TRM57971.00 CEDU、HOB2 AOAD/M_T CUT0、KZN2、TLSG TRM59800.00 高纬 DAV1 LEIAR25.R3 本文选取8阶球谐函数模拟电离层精细变化，并将卫星截止高度角设为10°，假设DCB在一天内为一个常数，与时间分辨率为2 h的球谐函数模型参数同步估计，得到2018年1月1日—30日共30 d的北斗二代卫星DCB，结果如图 2所示。
图 2 DCB估值时间序列
Figure 2. Time Series of Estimated DCBs
将本文方法估计的BDS多频DCB与中国科学院(Chinese Academy of Sciences，CAS)发布的BDS多频DCB(DCBC7I-C6I由发布的两种DCB线性表示)取差值，3种不同轨道卫星的DCB差值对比如图 3所示，30 d的平均偏差与均方根误差(root mean squared error，RMSE)见表 2。
图 3 本文方法估计的DCB与CAS发布的DCB偏差
Figure 3. Difference Between DCB Estimated by This Paper and DCB Provided by CAS
表 2 DCB平均偏差与均方根误差/ns
Table 2. Mean Bias and RMSE of DCB Estimated by This Paper and DCB Provided by CAS /ns
卫星号 DCBC2I-C7I DCBC2I-C6I DCBC7I-C6I 平均偏差 RMSE 平均偏差 RMSE 平均偏差 RMSE C01 -0.006 0.097 0.039 0.142 0.055 0.122 C02 -0.216 0.187 0.035 0.355 0.239 0.252 C03 0.167 0.114 -0.271 0.140 -0.210 0.157 C04 -0.281 0.158 -0.560 0.133 -0.089 0.150 C05 -0.105 0.300 0.361 0.372 0.123 0.312 C06 0.034 0.114 -0.042 0.140 0.029 0.180 C07 0.480 0.195 -0.042 0.235 -0.360 0.203 C08 0.205 0.119 0.234 0.125 -0.066 0.156 C09 -0.126 0.102 0.072 0.142 0.066 0.115 C10 0.222 0.105 -0.003 0.169 0.045 0.182 C11 0.000 0.208 0.223 0.272 0.447 0.697 C12 -0.201 0.151 0.017 0.215 0.080 0.487 C13 -0.125 0.098 -0.293 0.176 -0.045 0.147 C14 -0.072 0.194 0.200 0.200 0.435 0.848 从图 2可以看出，本文方法估计的DCB在30 d的时间跨度内比较稳定，仅有轻微的波动。从图 3可以看出，本文方法估计的北斗二代卫星DCB和CAS发布的DCB基本吻合。除了部分卫星(C11、C12和C14估计过程中数据量不足)偏差较大，其他卫星的DCB偏差大多低于0.5 ns，另外，IGSO与MEO卫星的DCB稳定性要高于GEO卫星。
为了验证本文方法的科学性，使用同一组数据，分别采用不加约束的独立求解方法和本文方法估计3种DCB结果，闭合差误差线对比如图 4所示。基于闭合差约束，本文方法估计的DCB闭合差基本控制在0.5 ns内，小于独立求解的DCB闭合差。由于数据缺失，C11、C12和C14号3颗卫星的闭合差和方差比较大，但均控制在0.4 ns和0.8 ns内。需要特别说明的是，由于所选跟踪站分布不均，制约了球谐模型的精度，导致本文方法估计的DCB闭合差不优于CAS在2015年前发布的3种独立估计的DCB闭合差。但可以预见，随着BDS能同时接收3频数据的跟踪站数量的增加，本文方法具有一定的实用性和科学性。
图 4 独立求解与本文方法估计的3种DCB平均闭合差
Figure 4. Mean Misclosure of Three DCBs Calculated by Independent Solution and the Proposed Method
独立求解与本文方法估计的3种DCB平均闭合差
-
目前IGS电离层分析中心DCB的研究主要集中在GPS和GLONASS系统，只有DLR和CAS发布了多系统的DCB产品。本文以DLR和CAS发布的北斗DCB作为参考，分析本文方法估计的频间偏差的日变化特性。DCB的日变化是DCB日估计相对于平均值的偏差情况，一定程度上能反映DCB的估计精度，计算式为:
$$ {B}^{n}=\sqrt{\frac{\sum\limits_{d=1}^{D}({\widehat{b}}^{d, n}-{\overline{b}}^{n}{)}^{2}}{D-1}} $$ (14) 式中，$ {B}^{n} $表示第n号卫星DCB估值的日变化; D表示计算天数; $ {\widehat{b}}^{d, n} $表示第d天第n号卫星的DCB估值; $ {\overline{b}}^{n} $表示第n号卫星的DCB月平均值。
本文方法估计的DCB与DLR、CAS发布的DCB产品的日变化对比如图 5所示。由于DLR和CAS只发布了两种DCB，所以在此只分析DCBC2I-C7I和DCBC2I-C6I的变化特征。
图 5 本文方法估计的DCB与DLR、CAS发布的DCB日变化对比图
Figure 5. Daily Variation Comparison of DCB Estimated by This Paper and DCBs Provided by DLR and CAS
从图 5中可以看出，本文方法估计的DCB与CAS发布的DCB产品日变化特征吻合，但其日变化值要大于DLR发布的DCB产品。有研究表明，单站估计和多站估计DCB精度相差可达5~10倍，所以卫星的频间偏差的精度随着跟踪站的数量的减少会明显降低。本文实验所用数据来自19个MGEX站，而CAS是逐基站对电离层建模并估计卫星DCB，DLR是采用全球分布的跟踪站数据估计BDS卫星DCB，前两者的跟踪站数量均远远少于DLR所用跟踪站数量。在2018年年积日001-030天时间范围内，DLR、CAS与本文方法估计的卫星DCBC2I-C7I估值日变化范围分别是[0.08，0.21]ns、[0.15，0.23]ns和[0.14，0.32]ns，DCBC2I-C6I估值日变化范围分别是[0.09，0.21]ns、[0.35，0.43]ns和[0.35，0.49]ns。
为了分析典型卫星的频间偏差稳定性，选取C02、C05、C08、C13四颗卫星的频间偏差时间序列进行对比分析，如图 6所示。从图 6中可以看出，3种DCB的变化趋势基本相同。由于3种DCB的估计基准不一样，而本文估计的DCB已经和CAS产品对齐，所以两者与DLR产品之间存在一个固定的偏差。另外，C02和C05卫星较C08和C13卫星的频间偏差波动性大，这与前两颗卫星在轨服务时间较长有关。图 6中出现的毛刺点可能和卫星环境突变或者数据质量相关，具体原因需进一步探究。
图 6 部分典型卫星的频间偏差估值变化时间序列
Figure 6. Time Series of DCBs of Some Typical Satellites
-
为了分析频间偏差改正对单点定位的影响，选取MGEX的CUT0、HKWS和JFNG跟踪站2018年1月11日的数据进行伪距单点定位实验。MGEX跟踪站数据的时间分辨率为30 s，采样截止高度角设为7°，进行BDS单系统单频标准单点定位(standard point positioning，SPP)实验。对BDS的B1和B2频点码观测量分别采用3种不同的处理策略: 不改正DCB、用本文实验产品改正DCB和用CAS产品改正DCB，分析其定位结果。3个测站DCB对B1和B2频点码观测量定位精度的改正情况如图 7所示。为了进一步量化3种策略的定位精度，对3个测站在东(E)、北(N)和天顶(U)及三维(3D)方向的定位结果进行统计，结果见表 3。
图 7 DCB改正前后SPP定位精度对比
Figure 7. Accuracy Comparison in SPP Experiments Before and After DCB Correction
表 3 不同DCB改正策略下单频SPP定位结果
Table 3. Single-Frequency SPP Positioning Results of Different DCB Correction Strategies
观测站 方向 B1 B2 不改正/m 本文方法改正/m 精度提高/% CAS改正/m 精度提高/% 不改正/m 本文方法改正/m 精度提高/% CAS改正/m 精度提高/% CUT0 E 2.32 0.55 76.3 0.52 77.6 0.69 0.53 23.2 0.54 21.7 N 2.13 0.78 63.4 0.80 62.4 2.23 0.76 65.9 0.74 66.8 U 6.71 4.39 34.6 4.29 36.1 8.61 6.72 22.0 6.60 23.3 3D 7.41 4.49 39.4 4.40 40.6 8.92 6.78 24.0 6.66 25.3 HKWS E 1.94 0.31 84.0 0.38 80.4 0.48 0.30 37.5 0.33 31.3 N 1.09 0.77 29.4 0.76 30.3 2.22 0.92 58.6 0.92 58.6 U 5.43 3.73 31.3 3.62 33.3 5.63 5.54 1.6 5.46 3.0 3D 5.87 3.82 34.9 3.72 36.6 6.07 5.63 7.2 5.55 8.6 JFNG E 1.85 0.32 82.7 0.38 79.5 0.48 0.35 27.1 0.38 20.8 N 1.23 0.82 33.3 0.80 35.0 2.14 0.80 62.6 0.98 54.2 U 6.73 4.67 30.6 4.56 32.2 9.09 7.01 22.9 6.87 24.4 3D 7.08 4.75 32.9 5.02 29.1 9.35 7.09 24.2 6.95 25.7 从图 7中可以看出，改正DCB误差的定位精度提高比较明显，可达分米级。除了在JFNG站B1频点观测量上，其他使用本文方法和CAS估计的DCB改正效果基本相当。由于DCBC2I-C6I估计的精度较低，所以两种产品在JFNG站B1频点上改正效果不一致。另外，B2频点上的误差改正值要略优于B1频点，这是因为B2频点的观测数据质量较差，所以改善相对比较明显。
从表 3可以看出，B1和B2频点在水平方向上改正效果比较明显，B1频点在E方向上改正效果表现更强，B2频点在N方向改正效果表现更强。但由于本文实验未对测站和卫星天线相位误差进行改正，并且其他误差在垂直方向改正残留更大，导致U方向上DCB改正效果并不明显。从表 3中可以看出，本文估计的和CAS发布的DCB在各方向上的改正精度提高基本相当，本文获得的DCB改正效果有时还要优于CAS。由于B2频点的观测数据质量较差及DCBC7I-C6I估计精度较低，B1频点的三维改正精度提高更为明显，达到50%。另外，由于JFNG站的数据质量及周边环境的影响，其定位改正精度要低于CUT0和HKWS两个站。
-
由于缺少官方发布的DCBC7I-C6I产品，DCBC7I-C6I通常使用DCBC2I-C7I、DCBC2I-C6I的线性组合来表示由此产生的误差累积会影响DCBC7I-C6I的精度。本文将北斗二代卫星3种频间偏差闭合差为零约束加入到DCB观测方程，以北斗IGSO和MEO卫星为基准分离卫星和接收机DCB，采用附有限制条件的间接平差方法估计了2018年1月1日—30日的北斗卫星3种频间偏差。为了验证本文方法的有效性和科学性，以CAS发布的DCB产品作为参考，分别分析了本文方法估计的DCB与参考DCB间的差值及日稳定性，结果表明本文估计的DCB与CAS发布的北斗卫星DCB吻合较好，较差在0.5 ns以内，但由于测站数量的限制，二者的日稳定性比DLR的差。将独立求解与本文方法估计的DCB闭合差对比可知，由于考虑了闭合差的约束，本文方法估计的DCB比独立求解的DCB的闭合差小。通过分析典型卫星的DCB时间序列，以及DCB改正对SPP的定位精度的影响，证明本文估计的DCB与CAS发布的DCB对SPP精度改善相当，改正DCB的定位精度提高了50%。随着各大导航系统的升级完善，下一步的研究重点是应用本文方法处理多导航系统的融合数据，获得更可靠的多系统DCB产品，为导航定位服务。
-
摘要: 中国的北斗卫星导航系统(BeiDou navigation satellite system，BDS)是全球导航卫星系统(global navigation satellite system，GNSS)中唯一全星座提供3频信号的卫星导航系统，其信号频率间共存在3种频间偏差(differential code bias，DCB)，分别是DCBC2I-C7I、DCBC2I-C6I、DCBC7I-C6I。理论上，这3种DCB之间存在代数和为零的关系。基于3种频间偏差闭合差约束，加入DCB观测方程，以北斗中轨道(medium earth orbit，MEO)和倾斜地球同步轨道(inclined geosynchronous satellite orbit，IGSO)卫星作为参考卫星，采用附加限制条件的间接平差方法同步估计BDS的3种DCB。选取2018年1月1日—30日多模GNSS实验(multi-GNSS experiment，MGEX)基准站的BDS 3频数据，分别采用附加闭合差约束估计和独立求解两种方法计算北斗二代卫星的3种频间偏差。以中国科学院(Chinese Academy of Sciences，CAS)和德国宇航中心(Deutsches Zentrum für Luft- und Raumfahrt，DLR)的DCB产品作为参考，分析了所提方法估计的DCB精度、稳定性及部分典型卫星的DCB时间序列，验证了所提方法对北斗3频DCB估计的适用性和科学性，并通过BDS单频标准单点定位(standard point positioning，SPP)实验验证了DCB对单点定位精度的影响效果。Abstract:
Objectives BeiDou navigation satellite system (BDS) is the only satellite navigation system that provides triple-frequency signals as a whole constellation of global navigation satellite system (GNSS). There are three kinds of differential code biases (DCBs): DCBC2I-C7I, DCBC2I-C6I and DCBC7I-C6I. Theoretically, there are zero algebraic relationship among the three DCBs. Methods Based on the zero algebraic relationship, BDS triple-frequency data of 2018 DOY 001-030 from multi-GNSS experiment (MGEX) stations are selected, and medium earth orbit (MEO) and inclined geosynchronous satellite orbit (IGSO) satellites are used as reference satellites to estimate three DCBs. Results Taking the DCB products provided by Chinese Academy of Sciences (CAS) and Deutsches Zentrum für Luft und Raumfahrt (DLR) as reference values, we analyze the accuracy and stability of the estimated DCB and the DCB time series of some typical satellites. It suggests that results derived from the proposed method are in good agreement with the reference values. And the experiment of standard point positioning (SPP) in BDS single frequency verifies the effect of DCB on accuracy of SPP. Conclusions All the results demonstrate that the proposed method is applicable and scientific to the estimation of BDS DCBs. -
图 3 本文方法估计的DCB与CAS发布的DCB偏差
Figure 3. Difference Between DCB Estimated by This Paper and DCB Provided by CAS
图 4 独立求解与本文方法估计的3种DCB平均闭合差
Figure 4. Mean Misclosure of Three DCBs Calculated by Independent Solution and the Proposed Method
图 5 本文方法估计的DCB与DLR、CAS发布的DCB日变化对比图
Figure 5. Daily Variation Comparison of DCB Estimated by This Paper and DCBs Provided by DLR and CAS
图 7 DCB改正前后SPP定位精度对比
Figure 7. Accuracy Comparison in SPP Experiments Before and After DCB Correction
表 1 MGEX基准站基本信息
Table 1. Basic Information of MGEX Tracking Stations
基准站分布 MGEX站 天线类型 低纬 DARW ASH700936D_M KAT1、SIN1 LEIAR25.R3 LMMF、SPTU TRM57971.00 中纬 AUCK、BRST、CHTI、MQZG、REUN、UCAL、WARK、WGTN TRM57971.00 CEDU、HOB2 AOAD/M_T CUT0、KZN2、TLSG TRM59800.00 高纬 DAV1 LEIAR25.R3 
下载: 导出CSV
表 2 DCB平均偏差与均方根误差/ns
Table 2. Mean Bias and RMSE of DCB Estimated by This Paper and DCB Provided by CAS /ns
卫星号 DCBC2I-C7I DCBC2I-C6I DCBC7I-C6I 平均偏差 RMSE 平均偏差 RMSE 平均偏差 RMSE C01 -0.006 0.097 0.039 0.142 0.055 0.122 C02 -0.216 0.187 0.035 0.355 0.239 0.252 C03 0.167 0.114 -0.271 0.140 -0.210 0.157 C04 -0.281 0.158 -0.560 0.133 -0.089 0.150 C05 -0.105 0.300 0.361 0.372 0.123 0.312 C06 0.034 0.114 -0.042 0.140 0.029 0.180 C07 0.480 0.195 -0.042 0.235 -0.360 0.203 C08 0.205 0.119 0.234 0.125 -0.066 0.156 C09 -0.126 0.102 0.072 0.142 0.066 0.115 C10 0.222 0.105 -0.003 0.169 0.045 0.182 C11 0.000 0.208 0.223 0.272 0.447 0.697 C12 -0.201 0.151 0.017 0.215 0.080 0.487 C13 -0.125 0.098 -0.293 0.176 -0.045 0.147 C14 -0.072 0.194 0.200 0.200 0.435 0.848
下载: 导出CSV
表 3 不同DCB改正策略下单频SPP定位结果
Table 3. Single-Frequency SPP Positioning Results of Different DCB Correction Strategies
观测站 方向 B1 B2 不改正/m 本文方法改正/m 精度提高/% CAS改正/m 精度提高/% 不改正/m 本文方法改正/m 精度提高/% CAS改正/m 精度提高/% CUT0 E 2.32 0.55 76.3 0.52 77.6 0.69 0.53 23.2 0.54 21.7 N 2.13 0.78 63.4 0.80 62.4 2.23 0.76 65.9 0.74 66.8 U 6.71 4.39 34.6 4.29 36.1 8.61 6.72 22.0 6.60 23.3 3D 7.41 4.49 39.4 4.40 40.6 8.92 6.78 24.0 6.66 25.3 HKWS E 1.94 0.31 84.0 0.38 80.4 0.48 0.30 37.5 0.33 31.3 N 1.09 0.77 29.4 0.76 30.3 2.22 0.92 58.6 0.92 58.6 U 5.43 3.73 31.3 3.62 33.3 5.63 5.54 1.6 5.46 3.0 3D 5.87 3.82 34.9 3.72 36.6 6.07 5.63 7.2 5.55 8.6 JFNG E 1.85 0.32 82.7 0.38 79.5 0.48 0.35 27.1 0.38 20.8 N 1.23 0.82 33.3 0.80 35.0 2.14 0.80 62.6 0.98 54.2 U 6.73 4.67 30.6 4.56 32.2 9.09 7.01 22.9 6.87 24.4 3D 7.08 4.75 32.9 5.02 29.1 9.35 7.09 24.2 6.95 25.7
下载: 导出CSV
-
[1] 徐绍铨. GPS测量原理及应用[M]. 武汉: 武汉大学出版社, 2017 Xu Shaoquan. Principle and Application of GPS Survey[M]. Wuhan: Wuhan University Press, 2017 [2] 李子申. GNSS/Compass电离层时延修正及TEC监测理论与方法研究[D]. 北京: 中国科学院大学, 2012 Li Zishen. Study on the Mitigation of Ionospheric Delay and the Monitoring of Global Ionospheric TEC Based on GNSS/Compass[D]. Beijing: Chinese Academy of Sciences, 2012 [3] 章红平, 韩文慧, 黄玲, 等. 地基GNSS全球电离层延迟建模[J]. 武汉大学学报·信息科学版, 2012, 37(10): 1 186-1 189 http://ch.whu.edu.cn/article/id/355 Zhang Hongping, Han Wenhui, Huang Ling, et al. Modeling Global Ionospheric Delay with IGS Ground-Based GNSS Observations[J]. Geomatics and Information Science of Wuhan University, 2012, 37(10): 1 186-1 189 http://ch.whu.edu.cn/article/id/355 [4] 章红平, 施闯, 唐卫明. 地基GPS区域电离层多项式模型与硬件延迟统一解算分析[J]. 武汉大学学报·信息科学版, 2008, 33(8): 805-809 http://ch.whu.edu.cn/article/id/1686 Zhang Hongping, Shi Chuang, Tang Weiming. United Solution to Polynomial VTEC Modeling and DCB Analysis Using Ground-Based GPS Observations[J]. Geomatics and Information Science of Wuhan University, 2008, 33(8): 805-809 http://ch.whu.edu.cn/article/id/1686 [5] 周东旭. GNSS接收机仪器偏差的时变特性研究[D]. 武汉: 中国科学院测量与地球物理研究所, 2011 Zhou Dongxu. Study on Time-Varying Characteristics of Instrument Deviation in GNSS Receiver[D]. Wuhan: Institute of Geodesy and Geophysics, Chinese Academy of Sciences, 2011 [6] Guo F, Zhang X H, Wang J L. Timing Group Delay and Differential Code Bias Corrections for BeiDou Positioning[J]. Journal of Geodesy, 2015, 89(5): 427-445 doi: 10.1007/s00190-015-0788-2 [7] 安家春, 艾松涛, 王泽民. 极区电离层TEC监测和发布系统[J]. 极地研究, 2010, 22(4): 423-430 https://www.cnki.com.cn/Article/CJFDTOTAL-JDYZ201004011.htm An Jiachun, Ai Songtao, Wang Zemin. Monitoring and Distribution System About TEC of Polar Ionosphere[J]. Chinese Journal of Polar Research, 2010, 22(4): 423-430 https://www.cnki.com.cn/Article/CJFDTOTAL-JDYZ201004011.htm [8] 吕志伟, 郝金明, 曾志林, 等. 利用GPS观测资料确定接收机差分码偏差的算法[J]. 全球定位系统, 2010, 35(2): 14-17 doi: 10.3969/j.issn.1008-9268.2010.02.004 Lü Zhiwei, Hao Jinming, Zeng Zhilin, et al. Algorithm of Receiver Differential Code Bias Determination Using GPS Measurements[J]. GNSS World of China, 2010, 35(2): 14-17 doi: 10.3969/j.issn.1008-9268.2010.02.004 [9] Tang L, Zhang X H, Lin X J, et al. Estimation of Compass Satellite Differential Code Biases Using Combined GPS/Compass Observation[M]//Lecture Notes in Electrical Engineering. Berlin, Heidelberg: Springer, 2012: 365-371 [10] Wang N B, Yuan Y B, Li Z S, et al. Determination of Differential Code Biases with Multi-GNSS Observations[J]. Journal of Geodesy, 2016, 90(3): 209-228 doi: 10.1007/s00190-015-0867-4 [11] Feltens J. The Activities of the Ionosphere Working Group of the International GPS Service (IGS)[J]. GPS Solutions, 2003, 7(1): 41-46 doi: 10.1007/s10291-003-0051-9 [12] Li Z S, Yuan Y B, Li H, et al. Two-Step Method for the Determination of the Differential Code Biases of COMPASS Satellites[J]. Journal of Geodesy, 2012, 86(11): 1 059-1 076 doi: 10.1007/s00190-012-0565-4 [13] 曾添, 隋立芬, 肖国锐, 等. BDS卫星端DCB不同表达对单点定位精度的影响[J]. 测绘科学技术学报, 2016, 33(4): 362-367 doi: 10.3969/j.issn.1673-6338.2016.04.007 Zeng Tian, Sui Lifen, Xiao Guorui, et al. The Influence on Single-Receiver Positioning Accuracy of the Different Expressions of BDS Satellite Differential Code Bias[J]. Journal of Geomatics Science and Technology, 2016, 33(4): 362-367 doi: 10.3969/j.issn.1673-6338.2016.04.007 [14] 曾添, 隋立芬, 鲍亚东, 等. BDS卫星端差分码偏差对定位的影响及改正模型研究[J]. 大地测量与地球动力学, 2017, 37(1): 53-57 https://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201701012.htm Zeng Tian, Sui Lifen, Bao Yadong, et al. The Impact of Satellite Differential Code Bias on BDS Positioning and Correction Model Research[J]. Journal of Geodesy and Geodynamics, 2017, 37(1): 53-57 https://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201701012.htm [15] Jin R, Jin S G, Feng G P. M_DCB: Matlab Code for Estimating GNSS Satellite and Receiver Differential Code Biases[J]. GPS Solutions, 2012, 16(4): 541-548 doi: 10.1007/s10291-012-0279-3 [16] 霍星亮. 基于GNSS的电离层形态监测与延迟模型研究[D]. 武汉: 中国科学院测量与地球物理研究所, 2008 Huo Xingliang. Study on Ionospheric Morphology Monitoring and Delay Model Based on GNSS[D]. Wuhan: Institute of Geodesy and Geophysics, Chinese Academy of Sciences, 2008 [17] 吴继忠, 施闯, 方荣新. TurboEdit单站GPS数据周跳探测方法的改进[J]. 武汉大学学报·信息科学版, 2011, 36(1): 29-33 http://ch.whu.edu.cn/article/id/952 Wu Jizhong, Shi Chuang, Fang Rongxin. Improving the Single Station Data Cycle Slip Detection Approach TurboEdit[J]. Geomatics and Information Science of Wuhan University, 2011, 36(1): 29-33 http://ch.whu.edu.cn/article/id/952 -
点击查看大图
图(7) / 表(3)
计量
- 文章访问数: 643
- HTML全文浏览量: 267
- PDF下载量: 69
- 被引次数: 0
