Improving BDS broadcast ephemeris accuracy using ground-satellite-link observations

The BeiDou Navigation Satellite System (BDS) employs a hybrid constellation including GEO (Geosynchronous Earth Orbit), IGSO (Inclined Geosynchronous Orbit), and MEO (Medium Earth Orbit) satellites, where the GEO and IGSO satellites are critical to providing continuous and reliable Positioning, Navigation, and Timing (PNT) services in the Asia–Pacific region. To handle the inconsistency between the satellite orbits and clocks in the broadcast ephemeris, which are determined by the Orbit Determination and Time Synchronization (ODTS) and the Two-way Satellite Time Frequency Transfer (TWSTFT) technique, respectively, we present the strategies using ground-satellite-link observations to improve the accuracy of broadcast ephemeris. The clock differences between the ODTS and TWSTFT techniques are used for correcting the radial orbit component to derive the refined orbits, which are used to generate the refined broadcast ephemeris. The test results show the precision of the refined orbits is improved by 50–60% in the 3-h to 12-h predicted arcs for the GEO satellites, and by 40–50% for the IGSO satellites. Moreover, the validation using satellite laser ranging observations shows the mean precision of the refined broadcast ephemeris is improved by 27% compared to the original one. Applying the proposed strategies in the BDS Operational Control Segment (OCS), the time evolution of BDS Single Point Positioning (SPP) in the period from Jan. 2016 to April 2021 is evaluated. The SPP accuracy is improved from 1.94, 2.06 and 3.29 m to 1.39, 1.85, and 2.39 m in the north, east, and up components, respectively. Further update with the inclusion of BDS-3 satellites improve the corresponding SPP precision to 0.68, 0.70 and 1.91 m.


Introduction
The BeiDou Satellite Navigation System (BDS) has provided regional Positioning, Navigation, and Timing (PNT) service in the Asia-Pacific region since Dec. 2012, and the global service since Dec. 2018 (Lu et al., 2020;Yang et al., 2020). A hybrid constellation is adopted in BDS including the Geosynchronous Earth Orbit (GEO), Inclined Geosynchronous Orbit (IGSO), and Medium Earth Orbit (MEO) satellites. The IGSO and GEO satellites play important roles for the continuous and reliable PNT service in the Asia-Pacific region as they are most visible in this region and provide augmentation information (Chen et al., 2022).
As the key mission for the ground control segment of the Global Navigation Satellite System (GNSS), the determination of satellite orbits and clocks provides the fundamental spatiotemporal reference for its PNT services. The BDS satellite orbit determination and time synchronization processing are based on two different techniques, i.e., the Precise Orbit Determination and Time Synchronization (ODTS) process for deriving satellite orbits and the Two-way Satellite Time Frequency Transfer (TWSTFT) technique for obtaining satellite

Open Access
Satellite Navigation https://satellite-navigation.springeropen.com/ *Correspondence: junping@shao.ac.cn clocks. The ODTS processing for GEO and IGSO satellites is much more difficult than for MEO satellites due to their regional visibility, nearly-static characters, and poor observation geometry (Steigenberger et al., 2013;Zhao et al., 2013), resulting in a high correlation between orbit radial component and satellite clock parameter, which further deteriorate the orbit accuracy, especially for GEO satellites.
Considering the challenges in the Precise Orbit Determination (POD) processing for BeiDou-2 Navigation Satellite System (BDS-2) satellites, including the limited number of regional tracking stations, the relatively poor orbit modeling due to the insufficient knowledge of the satellite physical property, and the weak observation geometry caused by the stationary GEO satellites, the satellite clocks estimated in the ODTS processing are not used for the production of BDS broadcast ephemeris. Instead, satellite clock parameters in BDS broadcast ephemeris are produced based on the TWSTFT real-time tracking system. TWSTFT utilizes the ground-satellitelink range observations between the ground stations and satellites. By differencing between uplink and downlink observations, the derived TWSTFT satellite clocks are free of the ionosphere and troposphere delay modeling errors, satellite orbit errors, etc., and have an accuracy of better than 0.06 m (Tang et al., 2016). For BDS-3 satellites, the TWSTFT observations are combined with the Ka-band Inter-Satellite-Link (ISL) observations, which expands the system tracking capability from regional to almost global coverage.
As the satellite orbits and clocks in the BDS broadcast ephemeris are generated based on the ODTS and TWSTFT techniques, respectively, the inconsistency between the BDS ODTS and TWSTFT techniques was reported (e.g., He et al., 2014;Zhou et al., 2016;Chen et al., 2017;Chen et al., 2020). As explained in He et al., 2014 andChen et al., 2017, the inconsistency between radial orbits and satellite clocks is because the satellite clock parameters estimated through the ODTS process are highly correlated with the orbit parameters and thus partly absorb the radial orbital parameters. On the other hand, the clocks from TWSTFT are free of any orbital errors. The inconsistency between orbits and clocks affects the overall Signal-in-Space Ranging Error (SISRE) of BDS broadcast ephemeris.
Noticing the between-day similarity of the ODTS-TWSTFT clock differences, He et al. (2014) proposed to fit the differences of the previous 24-h with a cubic spline function, and then shift the fitted function by 24-h as the predicted clock differences, which can be directly used for correcting the radial component of the predicted orbits. However, the clock differences between consecutive 24-h arcs do not always show high degree similarities, and are sometimes even quite different (e.g., 0.2 m for C04 in Fig. 4, He et al., 2014). In other words, this mathematical approximation without any physical background may not suit the real satellite orbits, which are usually described by the dynamic function, including the Kepler orbital elements and other parameters such as the Solar Radiation Pressure (SRP) model. Tang et al. (2016) further proposed to directly integrate the TWSTFT observations in the orbit determination processing, i.e. improving the orbit accuracy by additional constraints to the satellite clocks. Despite of much better SISRE attributed to the improved consistency between orbits and clocks, the orbit quality improvement is not optimal, and the 1-day orbit overlap deteriorates. The major explanations include that the TWSTFT observations are noisier than the carrier-phase observations, and systematic biases are not properly handled.
In this paper, we present a new method to mitigate the orbit and clock inconsistency in the BDS broadcast ephemeris. The clock differences between ODTS and TWSTFT are directly used for correcting the radial orbit components derived from the ODTS, and then the refined orbits are fitted with a dynamic model, i.e. the orbital elements and the SRP parameters, which are further used to predict the satellite orbits. Section 2 briefly describes the BDS broadcast ephemeris processing procedure. Section 3 analyzes the satellite clock differences between ODTS and TWSTFT techniques and presents a new approach for the production of the refined BDS broadcast ephemeris. Section 4 evaluates the performance of the new approach in terms of orbit precision, SISRE, and user positioning. Section 5 gives the main conclusion and discussions.

Procedure of the BDS broadcast ephemeris production
The two main types of parameters in the BDS broadcast ephemeris, i.e., satellite orbits and clocks, provide the spatiotemporal reference for BDS users. Figure 1 illustrates the main procedure of the BDS broadcast ephemeris production.

Precise orbit determination and time synchronization
To obtain the precise satellite orbits, the ODTS approach is implemented in the BDS Operational Control Segment (OCS). The ODTS is performed on an hourly basis, which processes the 3 days' observations of the regional tracking network for each orbit determination arc. Ionosphere-Free (IF) pseudo-range and carrier-phase observations serve as the basic observations, whereas for BDS-3 the ISL range observations between BDS-3 satellites are also utilized. The orbital elements of each satellite together with the additional orbit-dynamic parameters of a reference epoch are then estimated in the ODTS process. For BeiDou-3 Navigation Satellite System (BDS-3) satellites, the difficulties in the orbit determination of GEO and IGSO satellites are partly resolved by combing the observations of the ground tracking network and the Ka-band ISL observations (e.g. Tang et al., 2018;Yang et al., 2020;Ruan et al., 2020). The estimated orbit parameters are further used in the orbit integration process to derive the precise orbits covering the data span and a given period of predicted orbits beyond data coverage. Finally, the integrated satellite orbits covering both the data span and prediction period are used to produce the typical Keplerian orbital model (16 parameters and 18 parameters) via a least-squares estimator through the orbit fitting process, which provides the best trajectory fit in an Earth-Centered, Earth-Fixed (ECEF) coordinate system for each specific fit interval. Note that the satellite clock estimates simultaneously determined with the satellite orbits in the ODTS process are not used to generate the broadcast ephemeris.

Two-way satellite time frequency transfer
To obtain precise satellite clocks, several TWSTFT facilities are deployed in the BDS OCS. Each TWSTFT ground station is equipped with a L-band transmitter generating ranging signals, which are received by the BDS satellites, and the so-called up-link range observations at sampling epochs are generated. The up-link range observations are sent to the OCS via the communication link of the satellite. At the same time, the TWSTFT ground station receives the L-band ranging signals of the same satellite and generates the down-link range observations. The time difference between the TWSTFT ground station and the satellite is derived as the difference between uplink and down-link observations at the same epoch. Considering the clock frequency of the atomic clock of the TWSTFT stations and onboard satellites, the time tag of the up-link and down-link observations differs by a maximum of 1 × 10 −3 s, which may introduce a few ps errors in range observations. Following the same procedure, several time difference pairs between TWSTFT stations and satellites are obtained at the same epoch, and used to generate the clocks of each station and satellite under the same reference time frame, i.e., the BDS Time (BDT). The derived satellite clocks are then used to estimate the clock model parameters in broadcast ephemeris through a clock fitting process, shown in Fig. 1.

Clock differences between ODTS and TWSTFT
In the BDS ODTS processing, the reference clock is defined by the Hydron-maser at the master monitoring station. In the TWSTFT processing, clock reference is also based on the input of a Pulse Per Second (1PPS) of the Hydro-maser, which is connected to the master monitoring station via optical fibers. Therefore, the reference time frame of both the ODTS and TWSTFT techniques is the same, and satellite clocks based on these two techniques should theoretically be the same. Figure 2 shows the satellite clock differences between the clock estimates from ODTS and from TWSTFT for GEO satellites (C01-C05) and IGSO satellites (C06-C10) in three days in 2017. As shown in the figure, apparent bias and periodical signals exist for all the satellites. The constant biases of a few decimeters are due to the inconsistency between the calibrated hardware delays of TWSTFT and ODTS equipment. In addition to the constant biases, periodical variations are also observed, which mostly share a similar period with the satellite orbit cycle, especially for GEO satellites. The satellite clock differences with biases removed between ODTS and TWSTFT techniques are defined as ODTS clock residuals.
In Fig. 2, remarkable differences also exist among all the satellites. In Fig. 3, the three-day time series is divided into three 24-h batches, where each daily time series is shifted with a clock difference at the first epoch re-aligned to 0. For all the satellites, a clear signal related to the satellite revolution period exists in the daily satellite clock difference time series, where the amplitude of GEO satellites is bigger than that of IGSO satellites. The daily clock differences for IGSO show breaks when the As mentioned in previous studies, the TWSTFT clocks are usually considered free of systematic bias, whereas the ODTS clocks are highly correlated with the orbit especially in the radial component. Therefore, the clock differences between ODTS and TWSTFT products presented above are attributed to the ODTS processing, which was also explained by Zhou et al., (2011Zhou et al., ( , 2016, He et al. (2014), and Tang et al. (2016).

Refined orbits using TWSTFT clocks
As the TWSTFT derived satellite clocks are free of satellite orbit errors, the periodical signals in the clock residual time series come from the ODTS satellite clocks. Consequently, the orbit-related periodical errors of the satellite clocks will be compensated by the satellite orbits in the ODTS process, thus maintaining the minimum of the SISRE in the ODTS process. However, the legacy strategy of broadcast ephemeris determination is that satellite clocks are based on the TWSTFT process whereas satellite orbits are from ODTS. Therefore, inconsistency issue in the broadcast ephemeris exists and introduces the errors in the overall SISRE of GEO and IGSO satellites.
As explained in Chen et al. (2017), because the satellite clock parameters estimated from the ODTS process are highly correlated with the orbit parameters, satellite clocks partly absorb radial orbital parameters, resulting in similar periodic signals but in opposite signs for the two parameters. The periodic signals of radial satellite orbits can be evaluated by using Satellite Laser Ranging (SLR) observations, as the SLR residuals mainly reflect orbit errors in the radial component. Using SLR observations, Zhou et al. (2011) showed that the orbit-related periodical errors in the ODTS satellite clocks are compensated by the satellite orbits in the ODTS process.
The orbit-related errors contained in the ODTS clock residuals are used to correct ODTS orbits in this paper. The basic procedure is as follows: the satellite clocks from ODTS and TWSTFT process are differenced to derive residual time series; then the constant bias is removed to obtain the ODTS residuals containing periodical signals only; and finally the ODTS clock residuals are used to correct ODTS satellite orbits to derive refined satellite orbits. As satellite clock is highly correlated with the orbit radial component and they can be compensated by each other, the orbit correction model at each epoch is expressed as: where x y z T M is the modified satellite orbit, x y z T is the orbit from the ODTS process, and M is the conversion matrix from the satellite body frame (R,T,N) to the ECEF frame (X,Y,Z). dR dT dN T are the orbit correction vector of the radial, along, and cross components, which are from the clock residual time series and have the form: where d clk is the original epoch-wise clock difference between ODTS and TWSTFT, and d mean clk is the monitored constant bias term, and c is the light speed.
As ODTS is performed using the three-day observations, the calculation using eqution (1) and (2) is also performed in each three-day batch. After correcting the radial component of the ODTS orbits with the clock differences to generate the refined orbits, the refined orbits are further dynamically fitted to calculate the orbital elements and the SRP parameters, which are used for the orbit integration and prediction.
As shown in Eq.
(1), the basic idea of the orbit refining process is to correct ODTS satellite orbits using the clock difference between ODTS and TWSTFT. Applying the method, the whole procedure of the broadcast  ephemeris production is illustrated in Fig. 4. Compared with the original procedure in Fig. 1, the input satellite orbits for the orbit fitting process are the refined orbits using the proposed correction model. By refining ODTS satellite orbits using the developed algorithm, the consistency between the ODTS orbits and TWSTFT clocks, and consequently the consistency between the satellite orbit and clock parameters in the broadcast ephemeris, is improved. Unlike IGSO and GEO satellites, the non-visible time of a MEO satellite can reach more than 3/4 of its revolution. Therefore, MEO satellites are not considered in the refining process. For IGSO satellites, the orbits during the non-visible periods are not used, which accounts for a small fraction of the whole period, but there is still sufficient data for each processing batch.
For the numerical integration orbit, we adopted the Earth Gravitational Models (EGM09) gravity field with the order of 12 × 12, and applied the planet forces using the NASA (National Aeronautics and Space Administration) JPL (Jet Propulsion Laboratory) Development Ephemerides (DE405) planet ephemeris. The International Earth Rotation and Reference Systems Service (IERS) Conventions 2010 is followed for the effect of solid Earth tides, ocean tides, and pole tide. For solar radiation pressure modeling, the 5-parameter Empirical CODE (Center for Orbit Determination in Europe) Orbit Model (ECOM) is applied. The earth radiation pressure and antenna thrust effects are, however, ignored due to their relatively small effects, especially to the GEO satellites.

Performance evaluation
To evaluate the performance of the refined broadcast ephemeris, ODTS and TWSTFT data of DOY 147-156 in 2017 are used to generate the refined BDS orbits and broadcast ephemeris. In the process, every 3-day arc data form a batch and 8 batches for the 10 days are constructed, where 2-day arc data are overlapped between the consecutive batches.
For BDS users, the satellite clocks are based on the TWSTFT technique. Therefore, the orbits should be more consistent with TWSTFT clocks. Following the orbit refining process, the refined orbits fit better to the TWSTFT satellite clocks than the ODTS orbits, which results in smaller SISRE and better PNT results of BDS satellites. To assess the performance of the refined broadcast ephemeris, the differences between original ODTS orbits and the refined orbits are calculated, and the orbitonly SISRE of the differenced orbits is calculated with the following equation (Montenbruck et al., 2015): where R, T , N are the orbit differences in the R, T, N frame, and α, β are weighting factors with the values of 0.99 and 1/126, respectively, for the BDS GEO and IGSO satellites.

Refined orbits and the predictions
In the upper plot of Fig. 5, the difference between original ODTS orbits and refined orbits for satellite C01 in a 3-day batch is presented. The apparent periodical characters of the radial component are observed for the orbit differences, as drived based on Eq. (2). The same periods and similar amplitudes of the clock residuals (in Fig. 2) and orbit differences (in Fig. 5) are observed, demonstrating that the ODTS clock residuals are mostly absorbed into the refined orbits. The orbit-only SISRE values using the refined and ODTS orbits are presented in the lower plot of Fig. 5, showing the impact of orbits corrections on orbit-only SISRE.
We further investigate the performance of the predicted orbits. Data of a 2-day arc are used for the orbit refining process, and the refined orbits are predicted for another 24-h arc. Both the 2-day refined orbits and predicted orbits are compared with the refined orbits using 3-day arc data as the reference. The results are also plotted in Fig. 5, where we see that the orbit differences between the 2-day arc and the 3-day arc are small with the orbit-only SISRE differences of a few cm. During the prediction periods, the orbit-only SISRE remains small for the first few hours and gradually increases after 12 h.
Performing the same process for the other 7 batches, Fig. 6 shows for each satellite the orbit-only SISRE of the  where POD refers to ODTS orbits and FIT refer to refined orbits and the predictions original ODTS orbits, refined orbits of 2-day arc, 3-h, 6-h, and 12-h predictions, respectively. For all the GEO satellites and most IGSO satellites, the orbit-only SISRE of refined orbits shows significant improvement over the original ODTS orbits. Despite of the overall improvement for most of the satellites in most of the predicted arcs, the IGSO satellite C06 slightly deteriorated the orbit-only SISRE in the 3-h and 12-h prediction arcs by less than 0.05 m. The reason is that the along and cross components from the ODTS solution are treated as the reference "true" orbit, and any deviations are the errors and mapped into the orbit-only SISRE values. As the dynamic orbit fitting inevitably causes the differences in not only the radial but also the cross and along components. Nevertheless, for C06 the insignificant increase in small orbit-only SISRE does not necessarily impair the system performance.

SLR validation
SLR observations are used to further assess the accuracy of the radial components of the refined orbits. Although all BDS-2 satellites are equipped with laser retroreflectors, only three of the GEO and IGSO satellites (C01,   Table 1 summarizes the SLR validation results of these three satellites. During the selected period, the number of Normal Points (NPT) is 182,148 and 140 for C01, C08 and C10, respectively. The offset and STD of the SLR residuals are listed in Table 1. Apparent SLR offsets can be observed between the original and refined orbits of all satellites, which shows that the bias between ODTS and TWSTFT clocks are absorbed by refined satellite orbits, mainly in radial component. A shift in offsets after the refinement processing ensures the minimum SISRE of the broadcast ephemeris. Table 1 also shows that the STDs of refined obits are reduced with maximum rate 37% and average rate 27%, reflecting smaller fluctuation of the refined orbits. The reduction in STD is because the periodical terms between ODTS and TWSTFT clocks are also corrected. Figure 7 shows the SLR validation residuals of these three satellites, the mean offset of each satellite is removed in the figure in order to better illustrate the scatters of the two orbits.

SPP validation
The proposed strategy is implemented in the BDS OCS to generate the improved broadcast ephemeris. We further present the performance improvement of the broadcast ephemeris with Single Point Positioning (SPP) using dual-frequency IF observations. The BDS observations with sampling interval of 30 s on DOY 148 in 2017 at station KARR (Australia) are used. Figure 8 shows the daily time series of each coordinate component of SPP using GEO and IGSO satellites. Using the refined orbits, the mean SPP errors are reduced from 1. 6, 2.22, and 4.24 m to 1.46, 2.15, and 3.74 m in the north, east, and height components, respectively.
Despite of the short period data used in the above demonstration and validation, SPP using BDS observations from Jan. 2016 to April 2021 at station KARR is performed. Figure 9 presents the daily RMS time series for more than 5.4 years.
The time series in Fig. 9 are grouped into 4 periods according to the overall SPP accuracy, and the mean SPP position RMS during the four periods is summarized in Table 2. The first apparent improvement occurs on Jan. 2017, i.e., epoch 1 shown in Fig. 9, when the TWSTFT clocks are used to refine ODTS orbits in BDS OCS. Applying the refined approach, the mean SPP error is reduced from 1.94, 2.06, and 3.29 m to 1.39, 1.85, and 2.39 m in the north (N), east (E), and height (U) components, respectively, which corresponds to an improvement of 28%, 10%, and 27%. It is worth mentioning that the orbit refining method is only applied to the GEO and IGSO satellites, whereas the MEO satellite orbits remain the same. The SPP validation involves all the available BDS-2 satellites during this period, including the three MEO satellites. However, the only difference at the switching epoch "epoch 1" is whether the orbit refining strategy is implemented or not, so it is obvious that the orbit refining strategy is the only contribution to the SPP improvement. The second upgrade occurs in July 2017, the epoch 2 shown in Fig. 9, when the TGD (Timing Group Delay) parameter strategy in BDS OCS was upgraded (Yu et al., 2019). The mean SPP error is further reduced to 1.00 m, 0.72 m, and 2.42 m in the three components.
The third apparent coordinate accuracy improvement occurs since July. 2020, the epoch 3 shown in Fig. 9, when there is a receiver firmware update at the station and BDS-3 satellites can be observed. And the proposed approach in this paper is also applied to BDS-3 satellites. Although BDS satellites with Pseudo-Random Noise (PRN) code greater than 30 cannot be tracked due to the limitation of the receiver channels, another accuracy improvement of 0.32 m (32%), 0.02 m (3%), and 0.51 m  (21%) is observed. The coordinate improvement after the firmware updated is due to more visible satellites and the better orbits quanlity of BDS-3 satellites. The precision of ODTS orbits is improved by 50% for BDS-3 MEO satellites with the contribution of ISL observations (Yang et al., 2020). Furthermore, we also notice the large noise in the SPP time series, which is the result of pseudorange observation noise.

Conclusions
Adopting the hybrid constellation, BDS provides continuous and reliable services in the Asia-Pacific regions thanks to the GEO and IGSO satellites. However, the weak geometry and uneven tracking stations make the POD solution difficult, resulting in a high correlation between the derived satellite orbits and clocks. BDS uses two different techniques for the broadcast ephemeris generation, with the satellite orbits determined by the ODTS techniqueand the satellite clocks determined by the TWSTFT technique. Consequently, the broadcast ephemeris has inconsistency between the orbits and clocks, which degrades the users positioning performance. The refinement of BDS broadcast ephemeris can be achieved by combining the ground-satellite-link and inter-satellite-link observations. The inter-satellite-link observations are available only between BDS-3 satellites, which were discussed in previous studies, e.g., Tang et al. (2018), Yang et al. (2019), and Ruan et al. (2020). In this paper, we focus on the contribution of ground-satellitelink observations and present an algorithm to refine the orbits by correcting the orbit radial components with the clock differences between ODTS and TWSTFT products. The refined orbits are then used for dynamic orbit fitting to generate satellite orbit parameters in the broadcast ephemeris. The refined orbits are therefore consistent with the TWSTFT clocks, and thus the discrepancy between orbits and clocks in the broadcast ephemeris is mitigated.
We demonstrate that the accuracy of satellite broadcast ephemeris is significantly improved, with the predicted SISRE value reduced by 50-60% for the GEO satellites and by 40-50% for the IGSO satellites. The orbit quality improvement is further demonstrated by the satellite laser ranging validation with an average improvement of 27% for the C01, C08, and C09 satellites. Moreover, the time evolution (from Jan. 2016 to April 2021) of BDS single point positioning is presented, which reflects the impact of the implementation of the refined strategies in broadcast ephemeris generation. The overall SPP coordinate error series can be divided into four sections with different upgrade stages of BDS OCS: ODTS orbits refinement using TWSTFT clocks, TGD parameter determination, and BDS-3 satellites included. With the above upgrades, the overall SPP accuracy is improved from 1.94 m, 2.06 m, and 3.29 m to 0.68 m, 0.70 m, and 1.91 m in the north, east, and up components, respectively, corresponding to an improvement of 65%, 66%, and 42%.