Estimation of BDS pseudorange biases with high temporal resolution: feasibility, affecting factors, and necessity

A common practice adopted for the pseudorange bias estimation and calibration assumes that Global Navigation Satellite System satellite-dependent pseudorange biases vary gently over time. Whereupon satellite pseudorange biases are routinely estimated and provided as the products with low temporal resolution, e.g., hourly or daily, by the agencies. The story sounds unquestionably perfect under the acquainted assumption. To validate the inadequacy of the above hypothesis we herein present an approach to the estimate the BeiDou Navigation Satellite System (BDS) pseudorange biases with high temporal resolution. Its feasibility, affecting factors, and necessity are discussed. Concretely, the Geometry-Free function models are first constructed to retrieve the linear combination of the pseudorange biases; then the pseudorange Observable-specific Signal Bias (OSB) values with respect to baseline frequencies (e.g., BDS C2I/C6I) are estimated along with the ionosphere modeling; subsequently, all multi-frequency pseudorange OSBs are determined by using the ionospheric information with constraint conditions; finally, the possible Differential Code Bias sets are attainable with the estimated pseudorange OSBs. Using the observation data of four months when the estimated BDS pseudorange biases are stable, their reliability is demonstrated with the stability at the level of sub-nanosecond and the BeiDou-3 Navigation Satellite System (BDS-3) values more stable than that of BeiDou-2 Navigation Satellite System (BDS-2). The comparison between the estimated pseudorange biases and the Chinese Academy of Sciences products reveals that the accuracy of the estimated pseudorange biases is 0.2–0.4 ns. Moreover, the large magnitude of the short-term pseudorange bias variation in the tens of nanoseconds for the BDS-2 and BDS-3 are found in years 2021 and 2022, which are affected by two types of the satellite flex power for the BDS-2 and BDS-3, respectively. We stress that it’s necessary to estimate the BDS pseudorange biases with high temporal resolution in the case of the satellite flex power and the products currently provided by the agencies cannot reflect the true quantity under the circumstance.


Introduction
Global Navigation Satellite System (GNSS) pseudorange observations are affected by the hardware delays at both satellite and receiver ends arising from the transmission and reception chains, and signal generation by the digital and analog equipment (Montenbruck et al., 2014). Independent on the signal modulation type and frequency, pseudorange delay exhibits the distinctive characteristic (Sanz et al., 2017). Proper handling of pseudorange bias is of a critical concern for pseudorange based positioning and timing applications (Guo et al., 2015;Qin et al., 2021). Besides, to retrieve the unbiased ionospheric Total Electron Content (TEC) for ionospheric sensing, reliable pseudorange biases in the receiver and satellite ends need to be carefully considered Zhang, 2016;Zheng et al., 2022). Hence, pseudorange bias handling is a critical issue in precise GNSS data processing.
With respect to the receiver pseudorange biases, the earliest research assumed that the pseudorange biases did not change apparently over time (Banville & Langley, 2011;Håkansson et al., 2017). Later, many works validated that the receiver pseudorange biases varied significantly on a timescale of hours or less and were closely related to the receiver ambient temperature (Coster et al., 2013;Wanninger et al., 2017;Zhang et al., 2018). The receiver pseudorange bias variation, if not properly handled, will affect the performance of the precise GNSS data processing such as the ionospheric modeling. Some efforts have been made to reduce its effect such as the Modified Carrier-to-Code Leveling (MCCL) and Modified Precise Point Positioning (MPPP) methods (Zhang et al., 2018. Related to the satellite pseudorange biases, the terms designated for Broadcast Group Delay (BGD), Time Group Delay (TGD), Inter Signal Correction (ISC), Differential Code Bias (DCB), and pseudorange Observable-specific Signal Bias (OSB) are widely used. Considering the characteristics of respective GNSS, the BRoaDcast Ephemeris (BRD) provides the BGD, TGD, or ISC parameters using the regional or global networks with the specific receivers for the Global Positioning System (GPS), BeiDou Navigation Satellite System (BDS), Galileo navigation satellite system (Galileo) Quasi-Zenith Satellite System (QZSS), and NAVigation with Indian Constellation (NAVIC) systems in real-time services (CSNO, 2021;EU, 2016;ISRO, 2017;JAXA, 2018;Steigenberger et al., 2015). Exhibiting elevation-and frequency-dependent group delay variation and code phase variation for the pseudorange observables, the DCB in the case of applying the GNSS observables for different signals on different or even the same frequencies denotes the time differences of two individual signals on the equipment signal generation and reception chains (Kersten & Schön, 2017;Wanninger et al., 2017). The Center for Orbit Determination in Europe (CODE) has transformed the Differential Signal Bias (DSB) parameterization to OSB parameterization since 2016, which gradually affects the application of pseudorange biases in the International GNSS Service (IGS) community (Dach & Walser, 2015). With the Radio Technical Commission for Maritime Service (RTCM-SC) format and bias Solution INdependent EXchange (SINEX) format, the OSB application is popular (RTCM, 2016;Schaer, 2016;Villiger et al., 2019). Compared with DSB parameterization such as DCB, OSB concept is more flexible and directly corresponds to raw observations. For instance, pseudorange OSB needs a specific datum to transmit the differential bias to absolute bias and denotes the individual pseudorange bias in the undifferenced format, which provides the sufficient flexibility for the pseudorange bias handling. In summary, although both DCB and pseudorange OSB products are popular in GNSS community, the OSB calibration is more convenient for users.
With the observation data from a global network, a set of the GNSS pseudorange biases are estimated and provided by Different Analysis Centers (ACs) e.g., the product of Chinese Academy of Sciences (CAS) using the observations from a global network by applying the ionospheric modeling approach (Wang et al., 2016(Wang et al., , 2020Zhang et al., 2021), the product of German Aerospace Center (DLR) by applying a priori ionospheric information (Montenbruck et al., 2014), and the product of the CODE by combining the ionosphere and clock analysis (Dach et al., 2009), which exhibits better performance than the broadcast one and is commonly applied to evaluate the GNSS BGD quality (Montenbruck et al., 2014;Wang et al., 2016Wang et al., , 2019. The satellite pseudorange biases associated with the instrumental delays are usually regarded as stable on a continuous arc (e.g. one day) based on the common sense that the GNSS satellites transmit the signals with a constant total power and power ratio on each signal component (Zhong et al., 2016). Hence, almost all ACs or literature provide the pseudorange bias products in a daily scale (Deng et al., 2021;Wang et al., 2020). Nevertheless, the so-called 'flex power' ability trips the balance to arise the GNSS pseudorange bias variations by enabling the redistribution of the power between the signal components (Jimenez-Banos et al., 2010;Thoelert et al., 2019). As a part of the GPS modernization plan, the action of flex power on GPS Block IIF satellites for the remedy against jamming has been active since January 2017 (Esenbuğa & Hauschild, 2020). In 2020, two types of regional GPS satellite flex power were also observed (Yang et al., 2022a). The satellite flex power operation affects the GPS signal to noise density ratio (C/N 0 ) measurements as well as the GPS satellite DCB short-term variation .
Aside from GPS, BDS is also a focus issue in the GNSS community (Yang et al., 2020(Yang et al., , 2022b. However, few works present the BDS satellite pseudorange bias variation characteristics. Cui (2022) found that the BeiDou-2 Navigation Satellite System (BDS-2) Inclined Geo-Synchronous Orbit (IGSO) satellites exhibited the flex power operation and a great variation of the corresponding C2I-C6I DCB values was demonstrated in January of 2021. However, the specific time series of satellite DCB were not provided in the study.
Moreover, to further improve the BeiDou-3 Navigation Satellite System (BDS-3) service performance by providing the safe, stable, and reliable service, the in-orbit software of all BDS-3 satellites was upgraded from June to October in 2022 (http:// www. beidou. gov. cn/). During the upgrading progress, the BDS ran continuously and stably, and the service performance was improved steadily. The BDS pseudorange bias characteristics, discontinuities, and variations at this stage need to be fully investigated and clearly explained. Interpreting and understanding the BDS pseudorange bias stability is critical for the ionospheric study and other applications.
With the aforementioned understanding, we herein extend the Geometry-Free (GF) approach to estimate the BDS pseudorange biases with high temporal resolution for the BDS C2I, C6I, C7I, C1P, and C5P signals. The organization of this study proceeds as follow. The general GNSS equation, inter-and intra-frequency GF function models, and the pseudorange bias estimation are introduced in methodology part. The perspectives of BDS pseudorange bias analysis in terms of the feasibility, affecting factors, and necessity are discussed. Finally, the discussion and conclusion are given.

Pseudorange bias estimation
Studying the pseudorange bias characteristics and variations usually involves two sequential tasks, namely the linear combination of the pseudorange bias retrieval and the pseudorange bias estimation. The section presents the intra-and inter-frequency GF function models for the linear combination of the pseudorange bias retrieval and pseudorange bias estimation in detail after introducing the general GNSS observation equation.

GNSS observation equation
Considering the scenario of tracking m satellites with transmitting signals on jth frequency and on yth channel with one receiver at epoch t, the observation equation reads (Leick et al., 2015): where r, s, j, y, and t denote the receiver, satellite, frequency, channel, and epoch indices;p s r,j,y (t) and φ s r,j,y (t) denote the pseudorange and carrier phase observables, respectively; ρ s r (t) denotes the receiver and satellite (1) geometrical range; dt r (t) and dt s (t) denote the receiver and satellite GNSS clock offsets, respectively; T s r (t) denotes the slant tropospheric delay; I s r,1 (t) denotes the GNSS L1-based slant ionospheric delay; µ s j denotes the frequency-dependent multiplier factor; d r,j,y (t) and d s j,y (t) denote the receiver and satellite Uncalibrated Code Delays (UCDs), respectively; d s r,j,y denotes the Signal Distortion Bias (SDB) related to the receiver-satellite type owing to the response of receivers' correlator and front-end designed bandwidth (Hauschild & Montenbruck, 2016); b r,j,y and b s j,y denote the receiver and satellite Uncalibrated Phase Delays (UPDs), respectively; N s r,j,y denotes the integer ambiguity; ε s p,j,y (t) and ε s φ,j,y (t) denote the pseudorange and carrier phase measurement noises, respectively.
Specially, all quantities are labeled in unit of meters. The time-variant parameters are implemented with the epoch index t and the time-invariant parameters don't have the epoch index. The receiver UCDs are assumed to vary freely over time (Zhang et al., 2018). Although the satellite UCDs are generally regarded as a constant over a continuous period such as few hours or even one day (Xue et al., 2016), we herein consider both receiver and satellite UCDs as the time-varying parameters that possibly have the short-term variations.

Intra-frequency GF function model
By differencing the pseudorange observations on same jth frequency with respect to the xth and yth channels, we can obtain the intra-frequency GF function model as follows (Wang et al., 2020): where (·) GF,ij,xy = (·) i,x − (·) j,y denotes the GF combined operation for each variable. The Satellite-Plus-Receiver (SPR) and receiver-satellite related GF pseudorange biases remain because the slant ionospheric observables are eliminated by using the signals modulated on the same frequency.

Inter-frequency GF function model
Constructing the GF measurements with the GNSS observations on ith and jth frequency with respect to the xth and yth channels, the equation reads: (2) p s r,GF,j,xy (t) = d r,GF,j,xy (t) − d s GF,j,xy (t) + d s r,GF,j,xy + ε s p,GF,j,xy (t) (3) p s r,GF,ij,xy (t) = µ s GF,ij · I s r,1 (t) + d r,GF,ij,xy (t) − d s GF,ij,xy (t) + d s r,GF,ij,xy + ε s p,GF,ij,xy (t) φ s r,GF,ij,xy (t) = −µ s GF,ij · I s r,1 (t) + b r,GF,ij,xy − b s GF,ij,xy + N s r,GF,ij,xy + ε s φ,GF,ij,xy (t) The above equation represents a rank-deficient system, implying that the parameters are not estimable. Using the reparameterization operation, the MCCL approach is built by choosing the minimum set of parameters as datum (Teunissen, 1985;Zhang et al., 2018). Thereby, the full-rank linear equation of the MCCL model for the multi-frequency pseudorange and carrier phase observations by selecting the inter-frequency measurements can be written as : with where the top identifier ' ¯ ' denotes the re-parametrized estimable operation. Then, the full-rank equation of the multi-frequency MCCL approach is established. The estimated parameters include the biased slant ionospheric delay, receiver pseudorange bias variation with regard to each GF combination, multi-frequency inter-frequency bias (IFB), and ambiguity parameters.
Then, the k frequency full-rank MCCL function model with m observed satellites can be described as: which o m denotes the m-dimension row vector; z 2 = [0, 1] T ; ξ MCCL = τ , d r,GF,1,2 (t), · · · , d r,GF,1,k (t), b if , a k−1 T denotes the estimated parameters of the MCCL model, in and ε Φ MCCL denote the pseudorange and carrier phase GF measurement variance vectors; n k = [e k−1 , −1 · I k−1 ] T ; q k = diag q 2 1 , q 2 2 , . . . , q 2 k , in which q i denotes the measurement noise ratio on ith frequency; Q r = diag δ 2 p , δ 2 φ denotes the observed pseudorange and carrier phase precision matrix in zenith direction; Q m denotes the cofactor matrix with respect to the elevation diversity of m satellites; ⊗ denotes the Kronecker product operation.

Pseudorange bias estimation
Using the satellite intra-frequency pseudorange biases from the intra-frequency GF function model, biased slant ionospheric delay, and multi-frequency IFB from the inter-frequency GF function model, that are: The GNSS pseudorange biases with respect to each band are solvable. The pure slant ionospheric delay in the biased slant ionospheric delay needs to be removed using the prior ionospheric information or joint estimation. Consequently, the ionospheric delay in (7) is modeled by a spherical harmonic function with the thin-layer ionospheric model namely the modified single-layer mapping function (MSLM), which reads (Dach & Walser, 2015;Liu et al., 2020) where M s r denotes the ionospheric mapping function, R E = 6371 km denotes the mean Earth radius, H ion = 450 km denotes the assumed single layer height, γ = 0.9782 denotes the MSLM model coefficient, ϕ and s are the coordinates of the Ionosphere Pierce Point (IPP) in the sun-fixed geomagnetic frame, n and m denote the order and degree of the spherical harmonic function, P nm (sin ϕ) denotes the normalized associated Legendre function; Ã nm and B nm are the normalized coefficients.
After modeling the slant ionospheric delay, the equation for solving the satellite pseudorange biases is still rank deficient. The specifc constraints or datums are required to solve the rank deficient problem. Hence, three kinds of the constraints including the zero-mean, zero-ionosphere-free, and geometry-free conditions are applied to elminate the rank deficiency and estimate the satellite pseudorange OSBs, which read (Deng et al., 2021): where m ς denotes the number of satellites broadcasting the signal on the ξ th frequency. The zero-mean condition is applied to separate the satellite pseudorange biases from the receiver part. The zero-ionosphere-free and geometry-free conditions are used for the estimation of various types of satellite OSBs with a unified datum.
The GF combination with respect to the baseline frequency is used for the estimation of spherical harmonic parameters and pseudorange biases. With the solvable ionospheric information, the counterparts in other GF combination can be eliminated. After solving for the (8) satellite pseudorange OSB, the corresponding DCB is available consequently. For the better understanding of the entire processing procedure, Fig. 1 depicts the flowchart of the GNSS pseudorange bias estimation. The flowchart of the pseudorange bias estimation is divided into two steps: the linear combination of the pseudorange bias retrieval and the pseudorange bias estimation.

Pseudorange bias analysis
This section depicts the receiver network for the pseudorange bias estimation, BDS pseudorange bias analysis, and the effect of the satellite flex power on BDS pseudorange biases. The estimated pseudorange biases include both pseudorange OSB and DCB.

Receiver network for pseudorange bias estimation
The possible BDS B1I, B3I, B2I, B1C, and B2a pseudorange OSBs are estimated using the global BDS observation data from the International GNSS Service (IGS) Multi-GNSS Experimental (MGEX) network. Then, the BDS C2I-C6I, C2I-C7I, C1P-C6I, and C1P-C5P DCB sets are constructed consequently. In total, five types of the BDS pseudorange OSB and four types of BDS DCB sets are estimated. Considering that no agency provides the pseudorange bias product of high temporal resolution, the estimated pseudorange biases are compared with BRD as well as CAS products to validate the performance. Although the DCB can be regarded as the linear combination of two OSBs, there are two reasons why the DCB sets are also given here. First, the estimated pseudorange biases can be compared with the reference values, especially for the BRD values. Second, the pseudorange bias variation possibly exists only in the pseudorange OSB not in the DCB. To better explain the mechanism, various types of the pseudorange biases are depicted. Figure 2 shows the network of 125 Septentrio PolaRx5 (TR) receivers and the ground tracks of BDS satellites on January 1, 2022. The different colors denote the different tracked BDS satellites. Considering the receiver front-end characteristic and the correlator design, the SDB may arise and act as the systematic offsets for the GNSS network when using mixed types of receivers (Hauschild & Montenbruck, 2016). Hence, the receivers of the same types are applied to reduce the sensitivity of pseudorange biases in different receiver types and the SDB effects can be eliminated. Table 1 summarizes the information of the operational receivers supporting the BDS multi-frequency signals, including the receiver manufacturers, stations, and estimated pseudorange bias types. Table 2 introduces the processing strategy of the pseudorange bias estimation in two steps. Moreover, for the convenience of identification, the estimated BDS pseudorange bias using the Multi-GNSS Positioning and Analysis System (MGPAS) software is generalized the 'MGP' in this study Su et al., 2022).

BDS pseudorange bias analysis
For the validation of the estimated BDS pseudorange biases, the observation data at 125 IGS MGEX stations during the Day of Year (DOY) 1-120 in 2022 are used to estimate the BDS pseudorange biases with high temporal resolution and analyze their performance. To evaluate the accuracy of the estimated pseudorange biases, the BRD and CAS pseudorange bias products are used as the reference values, which are available at ftp:// gdc. cddis. eosdis. nasa. gov/ pub/ gnss/ data/ daily/ and ftp:// ftp. gipp. org. cn/ produ ct/ dcb/ mgex/, respectively. Thereafter, the estimated pseudorange biases are analyzed in terms of consistency and stability.

GNSS pseudorange and carrier phase observations GF function models
Global ionosphere modeling
Moreover, the multi-GNSS daily pseudorange bias products (i.e. pseudorange OSB and DCB) for BDS, GPS, GLONASS, Galileo are routinely provided by the CAS agency. To validate the reliability of the estimated pseudorange biases, Fig. 4 depicts the time series of the CAS and MGP pseudorange biases from January to April 2022, in which the BDS pseudorange OSB and DCB sets of all possible types are shown. Similar to those of the BDS TGD, the time series of the BDS pseudorange biases in CAS are continuous and quite stable. The mean Standard Deviation (STD) range values of the BDS pseudorange biases are 0.28, 0.43, 0.50, 0.28, 0.36, 0.15, 0.16, 0.12 and 0.13 ns for the C2I, C6I, C7I, C1P, and C5P pseudorange OSB and C2I-C6I, C2I-C7I, C1P-C6I, and C1P-C5P DCB values, respectively, indicating that all the BDS pseudorange biases in CAS products are stable within four months as well. Thereby, we assume that the magnitude of the BDS pseudorange bias variation is small during the period.
Besides, the mean values and STDs of the estimated BDS pseudorange biases for all types are illustrated in Figs. 5 and 6 to depict the performance in terms of the magnitude and stability. The statistic of the mean STD values of all the pseudorange biases is depicted in Fig. 6 as well. It can be seen that the BDS pseudorange biases are extensive and range from − 200 to 100 ns. The BDS-3 pseudorange biases are generally larger than those of the BDS-2 values. The magnitude of the BDS pseudorange biases depends on the BDS satellite types and the differences are significant. The estimated C2I and C6I pseudorange OSBs for a satellite always have the same sign, and the absolute magnitude of the C2I pseudorange OSB is smaller than that of the C6I value. The phenomenon is normal owing to the zero-ionosphere-free and geometryfree conditions applied. Regarding the BDS pseudorange bias stability, the BDS DCB exhibits the better stability than the BDS pseudorange OSB. The BDS C2I and C1P pseudorange OSBs exhibit the best stability among all BDS pseudorange OSB types.
To better understand the accuracy of the pseudorange biases, the pseudorange biases provided by different agencies are generally compared with each other . The pseudorange biases possibly exhibit the systematic errors due to different processing strategies, which need to be removed by realigning the bias series with the common fixed satellite set and datum (Schaer and Dach, 2010). With the common satellite sets applied, the Root Mean Square (RMS) values of all differences are calculated. Herein, we compare the BDS pseudorange biases of MGP, CAS, and BRD products to analyze the bias consistency. Figure 7 shows the RMS of the differences of the pseudorange biases among MGP, BRD, and CAS. Specially, the TGDs are transformed to the corresponding DCB for the comparison. The pseudorange biases for the C2I and C6I signals are compared for all BDS satellites. Also, the pseudorange biases related to the C7I signal are compared for the BDS-2 satellites and the remaining pseudorange biases are compared for the BDS-3 satellites. In the case that the realigned TGD 1 and CAS or MGP C2I-C6I DCB values are directly compared for the selected four months, the RMS values reach approximately 1 ns. As shown in the bottom panel of Fig. 7, the estimated BDS pseudorange biases from MGP and CAS show a great agreement. The bias RMS differences between MGP and CAS DCB are in the order of 0.2-0.4 ns. Particularly, we analyzed the consistency of the pseudorange biases for different constellations in BDS-2, i.e., GEO, IGSO, and MEO, and the results are depicted in each panel. The results show that the consistency of the MEO satellite pseudorange biases is slightly better than that of the IGSO and the GEO satellite exhibits the worse performance.

Effect of the satellite flex power on BDS pseudorange biases
The GNSS C/N 0 observations in geodetic receivers can detect the operation of the satellite flex power. Owing to the antenna gains, low noise amplifiers, and estimated approach in different receivers, the corresponding observed C/N 0 values are different (Falletti et al., 2011). Although the C/N 0 measurements of different receivers are not directly comparable, the variation at a specific station can be used as an indicator to sense the drastic variations of the transmit power.
We firstly focus on an exemplary case for the BDS-2. During the July 1-10, 2022 (DOY 182-191, 2022), the BDS-2 C/N 0 variations were observed at the GNSS stations. Figure 8 illustrates the time series of the C/N 0 measurements of the S2I, S6I, and S7I at the IGS stations COCO, JFNG, WUH2, and URUM, where the power variations of the B1I, B3I, and B2I observations were observed. All receivers consistently detected the power adjustment with the C/N 0 variations for the S6I observations. The power of the B1I and B2I observations exhibits the normal variation characteristics. The S6I time series increase and decrease by 2-10 dB on DOY 183 and DOY187, respectively. The increased C/N 0 values last for 4 days and then return to the previous status. The C/N 0 variation characteristics of the IGSO and MEO satellites is different from that of GEO satellite in that the C/N 0 measurements vary with the elevation. Hence, we can assume that the variation of the BDS-2 C/ The satellite flex power means the change of the modulation mode of a GNSS signal, possibly manifested by pseudorange bias variation. To validate the assumption and analyze the effect of the satellite flex power on BDS pseudorange biases, the pseudorange bias variation is analyzed. Consistent with the above period, Fig. 9 depicts the time series of the corresponding BDS-2 pseudorange biases from BRD, CAS, and MGP. We can see that the three types of the pseudorange biases from different agencies reflect the variation with large magnitude. The large bias changes related to the C6I signal are observed in the case of the satellite flex power, including the C2I, C6I, and C7I pseudorange OSB and C2I-C6I DCB. All types of pseudorange OSBs exhibit the large variation though only S6I observations exhibit the satellite flex power. It's not surprising because the OSB variations are determined by the chosen datum and their relationships cause the variation tendency to some extent. The maximum variation magnitudes of the pseudorange biases from DOY 183 to 192 are summarized in Table 3. The results show that the DCB variation ranges from 4 to 13 ns and the pseudorange OSBs range within the 40 ns. The pseudorange biases with the sampling rate of 30 s can rapidly reflect the practical situation whereas the low-resolution BRD and CAS products have the time delay with varying degrees.
To better analyze the pseudorange bias variation, Fig. 10 depicts the all observed cases of BDS-2 C2I-C6I DCB variations in years 2021 and 2022. Considering the space limitation of the article, all pseudorange OSBs are not shown here for they can be regarded as a linear function of the C2I-C6I DCB. The C/N 0 values observed at the randomly selected station JFNG or COCO are also shown in the figure. Eight periods of the pseudorange bias variation are detected and shown here. Based on the actual condition, the DCB variation time scale in 2022 is larger than that in 2021. A significant jump of the pseudorange biases for nearly all BDS-2 satellites can be observed when the S6I in the signal-to-noise power density ratio in the case that satellite flex power is active. The pseudorange biases in BRD and CAS products with low temporal resolution cannot reflect the actual situation of the pseudorange bias variation. Some sudden jumps can also be observed in BRD TGD, which may be attributed to the gross errors in BRD product. The duration of the pseudorange bias variation in MGP is longer than those by other agencies owing to the processing strategy in terms of the filter noise set and the utilization of the observation duration. Hence, the pseudorange biases with high temporal resolution are necessary in the case of the satellite flex power, otherwise the performance of satellite clock offset, positioning, timing, and ionosphere sensing will be affected.
After discussing the BDS-2 pseudorange bias variations, their impact on BDS-3 pseudorange biases is analyzed subsequently. Herein, the time series of the BDS-2 and BDS-3 C2I-C6I DCB from May to October in 2022 from BRD, CAS, and MGP products are depicted in Fig. 11. We can see that the larger variation of the BDS-2   C01  C02  C03  C04  C05  C06  C07  C08  C09  C10  C11  C12  C13  C14  C16  C19  C20  C21  C22  C23  C24  C25  C26  C27  C28  C29  C30  C32  C33  C34  C35  C36  C37  C38  C39  C40  C41  C42  C43  C44  C45  The mean STDs of each type for all BDS satellites are also shown pseudorange biases result in the BDS-3 pseudorange bias variation in CAS product. The phenomenon is explainable because a common zero-mean condition is applied for BDS-2 and BDS-3 in CAS product. The effect will be eliminated if two independent datums are applied for BDS-2 and BDS-3, and the assumption is also validated in BRD product. In MGP product, the BDS-2 and BDS-3 C2I-C6I DCB sets are estimated with one joint and two individual datums, in which one and two zero-mean conditions are applied for the pseudorange bias estimation in two schemes, respectively. Consistent with the CAS product, the BDS-3 pseudorange biases are seriously affected by the BDS-2 bias variation when a common datum is applied. Based on the same principle, the case of the BDS-3 pseudorange bias variation will also affect the BDS-2 pseudorange biases to a certain extent such as the DOY 258 and 263 in 2022. When using the two individual datums, the pseudorange bias effect is only within each constellation. Hence, it would be much safer to independently estimate BDS-2 and BDS-3 pseudorange biases for the GNSS users. Next, we focus only on the BDS-3 pseudorange biases. Figure 12 depicts the time series of the BDS-3 pseudorange biases from May to October 2022 in BRD, CAS, and MGP products. The time series with one joint and two separate datums are also given. The possible constructed pseudorange OSB and DCB sets are all given, including the C2I, C6I, C1P, and C5P pseudorange OSB and C2I-C6I, C1P-C6I, C1P-C5P, C2I-C1P, and C5P-C6I DCB. The BDS-3 pseudorange biases vary significantly in a disorderly fashion due to the zero-mean condition and the variation of some specific satellites. Compared with the BDS-2, the BDS-3 pseudorange bias variation is relatively smaller except for the two cases of satellites C45 and C46 on DOY 258 and 263 in 2022. When estimating the BDS-3 pseudorange biases with the two individual datums, some BDS-3 pseudorange biases do not vary. Nearly all the BDS-3 pseudorange biases obviously have long-term trends except for the C2I-C1P and C5P-C6I DCB. The reason is that the characteristic of the C2I and C1P or C5P and C6I is similar, which is attributed to Dual-frequency Constant Envelope Multiplexing (DCEM) technique and similar frequency values.
Aside from the BDS-3 pseduorange bias variation caused by the BDS-2, BDS-3 pseudorange biases also have the inherent variation. Taking the C45 and C46 satellites as an example, Fig. 13 depicts the time series of the BDS-3 C45   (Esenbuğa & Hauschild, 2020). We can treat it as another type of the satellite flex power.

Discussion and summary
In summary, this section answers the following questions for a better understanding of the article.
(1) Q: What are the affecting factors of the BDS pseudorange bias variation?
A: The redistribution of the transmit power among the GNSS signal components is known as the satellite flex power. The GPS satellite flex power modes change frequency after 2020 (Esenbuğa et al., 2023). The GPS satellite flex power affects the magnitude of the satellite pseudorange biases for different satellite and type modes. The same also applies to BDS. Based on the analysis, the BDS-2 and BDS-3 are affected by two types of satellite flex power. BDS-2 redistributes the power of S6I signal.
For BDS-3, it operates to simultaneously increase the transmit power on individual signal. The frequency of the BDS-2 satellite flex power is significantly higher than that of BDS-3.

(B)
Q: At which stage we need to estimate the pseudorange biases with high temporal resolution?
A: The large variation of the BDS pseudorange biases has been observed in the case of the satellite flex power. At this stage, BDS pseudorange biases with high temporal resolution are necessary for more reliable services. Therefore, we can provide the "normal" or "flex" pseudorange biases for GNSS users with different sampling rate. This is particularly worthy of attention for some ACs pseudorange bias products such as CAS and DLR.

(C)
Q: Which constraint condition is more reasonable for the BDS pseudorange bias estimation of high temporal resolution in the case of satellite flex power?
A: Using the zero-mean condition, the pseudorange bias stability of some satellites will be affected in the case of the satellite replacement or pseudorange bias variation no matter one or two zero-mean conditions are applied. However, the accuracy of the estimated pseudorange biases will not be deteriorated. Using two zero-mean conditions can reduce pseudorange bias variation, but doesn't affect the correction of the pseudorange observation with the reliable and high-accuracy   A: The ACs release the pseudorange bias products including daily real-time product as well as the postprocessing product. For real-time pseudorange bias product, the ACs can process the GNSS observations consecutively for routine product with the unidirectional filter. For postprocessing pseudorange bias product, the bidirectional filter can be used for the stability of the product.
(E) Q: Which aspects does the pseudorange bias variation affect?
A: Precise estimation of pseudorange biases plays a crucial role in precise ionospheric modeling, positioning, navigation, and timing services. The performance of the satellite clock offset is also affected. Considering the space limitation of the article, the paper only discusses the feasibility, affecting factors, and necessity of the BDS pseudorange bias estimation with high temporal resolution. Future work related to the effect of the pseudorange bias variation on the satellite clock estimation, positioning, timing, and ionosphere sensing will be investigated.

Conclusion
Precise estimation of the satellite pseudorange biases is of great significance in ionosphere sensing, positioning, and timing services. This study presents the GF function model to estimate the BDS pseudorange biases with high temporal resolution, including the pseudorange OSB and DCB. The feasibility, affecting factors, and necessity of the BDS pseudorange biases estimation with high temporal resolution are demonstrated. The possible BDS pseudorange biases of the B1I, B3I, B2I, B1C, and B2a signals are estimated using the global BDS observation data from the IGS MGEX network. The main conclusions can be drawn as follows.
(1) The stability of the estimated satellite pseudorange biases with high temporal resolution is at the level of sub-nanosecond and BDS-3 pseudorange bias stability is slightly better than those of BDS-2. Compared with the CAS product, the accuracy of the estimated pseudorange biases is 0.  Fig. 13 Time series of the BDS-3 C45 and C46 satellite pseudorange biases from BRD, CAS, and MGP products biases can be estimated with the product of high temporal resolution with certain reliability.
(2) Large jump and variation of the satellite pseudorange biases in the tens of nanoseconds for the BDS-2 and BDS-3 are observed in years 2021 and 2022 in the case of the two types of satellite flex power. BDS-2 redistributes the power of S6I and BDS-3 operates to simultaneously increasing the transmit power on individual signal. We stress that it's necessary to estimate the BDS pseudorange biases of high temporal resolution in the case of the satellite flex power and the low-resolution pseudorange bias products in the agencies cannot reflect the true quantity under the circumstance.