A resilient adjustment method to weigh pseudorange observation in precise point positioning

The accurate weighting of pseudorange observations is important to improve the convergence time and positioning quality of Precise Point Positioning (PPP). Currently, the weight of a pseudorange observation is mainly determined with empirical stochastic models. However, in a complex environment, due to the inability to adapt for the dynamic changes of the user environment, the empirical stochastic models usually cannot reflect the real error level of pseudorange observations. To address this problem, a resilient adjustment method to weigh pseudorange observations is proposed, which constructs the real-time estimation and inflation model for the variance of pseudorange multipath error and measurement noise to replace the empirical stochastic model to determine the weights of pseudorange observations. A set of static and dynamic Global Positioning System (GPS) test data are used to verify the effectiveness of the proposed method. The test results indicate that the proposed real-time estimation model can provide a better representation of the pseudorange accuracy, and the positioning performance of PPP using the real-time estimation model is better than that with the empirical stochastic model. Compared with the optimal empirical stochastic model, the positioning accuracy of PPP with the real-time estimation model is improved by at least 20%, and the convergence time is reduced by at least 50%.


Introduction
In a complex environment, due to the influence of various error sources the accuracy of observation data is much related to the observed satellites, especially pseudorange observations. In order to improve the convergence time, positioning reliability and accuracy of Precise Point Positioning (PPP), it is necessary to properly weigh each observation. A stochastic model is used to describe the statistical characteristics of observation errors, which is the key to the determination of the weights of observations (Dovis et al., 2015). Therefore, accurate stochastic model is important for improving the positioning performance.
In order to determine the weights of observations, researchers have proposed a variety of empirical stochastic models, which are usually a function of the satellite elevation angle, signal-to-noise ratio, carrier-to-noise ratio, or multipath factor, to estimate the measurement variance (Parvazi et al., 2020;Zhang et al., 2019aZhang et al., , 2019b. The coefficients of an empirical stochastic model are usually fixed. However, in complex environments, measurement error is related to many factors such as the signal power, the types of receivers and observations, the design of the antenna, the satellite constellation, and the satellite elevation angle . Because those factors vary significantly during the process, it is difficult to accurately evaluate the statistical characteristics of the measurement error in real-time by using an empirical

Open Access
Satellite Navigation https://satellite-navigation.springeropen.com/ *Correspondence: zhangqieqie@nuaa.edu.cn stochastic model. Therefore, the empirical stochastic model used usually cannot reflect the real error level of observations in complex environments.
To deal with the time-varying and unpredictable observation errors in complex environments, many scholars have carried out extensive research on resilient stochastic models and put forward various resilient adjustment methods (Yang, 2019;Zhang et al., 2020), including Multiple Model Adaptive Estimation (MMAE), Variance Component Estimation (VCE), and covariance matching methods. The MMAE method (Shen et al., 2022) adopts a set of parallel filters with different observation covariance matrices to estimate the state parameters, and then fuses the state estimation according to the a posterior probability of each filter. This method can obtain the optimal state estimation theoretically, but needs to set many parallel filters, which will lead to serious memory resource consumption and affect the real-time performance. The VCE methods include Helmert Variance Component Estimation (HVCE) (Gao et al., 2016), Minimum Norm Quadratic Unbiased Estimator (MINQUE) (Zhang et al., 2019a(Zhang et al., , 2019b, the Best Invariant Quadratic Unbiased Estimator (BIQUE) (Zhang et al., 2019a(Zhang et al., , 2019b, and the Least-square Variance Component Estimator (LS-VCE) (Parvazi et al., 2020). They generally use the residuals to estimate the (co)variance of observations, which requires a lot of redundant observations to obtain reliable (co) variances. Therefore, the VCE methods are mainly used for post-processing analysis or to determine the scale factors of different types or groups of observations (Parvazi et al., 2020;Zhang et al., 2019aZhang et al., , 2019b. The covariance matching methods include the Residual-based Adaptive Estimation (RAE) (Lee et al., 2020) and the Innovationbased Adaptive Estimation (IAE) , which calculate the measurement covariance matrix based on the numerical relationship between the measurement error and measurement residual or innovation, respectively. Although the covariance matching methods have a strict mathematical basis, they require that the measurement residual or innovation vectors are of the same dimension at all epochs, which limits their application in the field of satellite navigation and positioning where the number of observations changes with time (Yang & Gao, 2006).
Although a resilient stochastic model can adjust the weights of observations adaptively, the resilient adjustment methods mentioned above are not suitable for PPP. At present, the empirical stochastic models are mainly adopted to determine the weights of observations in PPP. Since the measurement error is the most intuitive manifestation of the measurement accuracy, using the measurement error sequence extracted from the observations to analyze the measurement accuracy should best reflect the true error level of the observations (Zheng & Guo, 2016). In view of this, a resilient adjustment method for pseudorange observation weighting is proposed in this paper. In this method, the real-time estimation model is constructed, which is used to estimate the pseudorange variance based on pseudorange multipath error and measurement noise sequences, to replace the empirical stochastic model for determining the weights of pseudorange observations. This article is arranged as follows. The observation function model and the empirical stochastic model of PPP are first introduced, and the advantages and disadvantages of the empirical stochastic model are analyzed. Then, the resilient adjustment method and model to establish the pseudorange observation variance are introduced in detail. Afterwards, the performance of the proposed model is evaluated and discussed using the Global Positioning System (GPS) observation data in a static and a dynamic tests. Finally, summary and conclusion are given.

Observation model of precise point positioning
There are two types of observations in PPP, namely pseudorange and carrier phase. Considering the main error items, the pseudorange and carrier phase measurements can be modeled as (Pan, 2018) where i represents the frequency index, c denotes the speed of light in vacuum, P i and L i are the pseudorange and carrier phase measurements obtained with a receiver, respectively, ρ is the geometric range between the satellite and the receiver antenna phase centers at emission and reception time, t s is the satellite clock offsets, T r is the tropospheric delay, I i is the ionospheric delay, d mp and represent the carrier phase and pseudorange multipath errors, respectively, ε and ς are the pseudorange and carrier phase measurement errors, respectively, t r and B i are the receiver clock offset and the carrier phase bias, which are unknown parameters.
The error sources affecting the accuracy of observations include satellite orbit and clock errors, atmospheric transmission delay, multipath errors, and measurement noise. Among them, the satellite orbit and clock errors and tropospheric delay can be corrected by models (Guo et al., 2017;Ma et al., 2021;Stpniak et al., 2015). Considering that the ionospheric delay is inversely proportional to the square of the signal frequency, the ionospheric delay can be eliminated by the Ionospheric-Free (IF) observation model as follows (Pan et al., 2021). (1) where P IF and L IF are IF pseudorange and carrier P IF phase observations, which are a linear combination of the pseudorange and carrier phase observations on frequencies f j and f i , respectively. However, multipath error and measurement noise are difficult to model and usually cannot be eliminated (Seepersad & Bisnath, 2015;Lau, 2017;Zhang et al., 2018). They are the main factors that reduce the accuracy of observations, especially in a complex environment such as urban.

Empirical stochastic model
The classical method to weigh the observations is to use an empirical stochastic model to estimate the observation variance. Assuming that each error source obeys a Gaussian normal distribution with zero mean, and considering the contributions of all error sources, the observation variance for each satellite in view can be obtained as (Dovis et al., 2015) where σ i is the standard deviation of the measurement to the ith satellite. The details of the other symbols are as follows: • σ clk&eph is the standard deviation of the satellite clock and orbit errors, which is usually referred to User Range Accuracy (URA) (Bijjahalli & Sabatini, 2021). With broadcast ephemeris, the typical value of URA is 2.4 m, while with precise orbit and clock products, the URA is generally several centimeters or better (Guo et al., 2017). • σ trop is the standard deviation of troposphere correction model. Because the troposphere delay correction error is mainly caused by the hydrostatic delay σ TVE correction model error, the standard deviation of troposphere delay is modeled by (Dovis et al., 2015) where σ TVE is the standard deviation of the zenith troposphere delay, which is set to an empirical value 0.12 m, M dry (β EL ) is the dry mapping function, and β EL is the satellite elevation angle. • σ ion denotes the standard deviation of the ionosphere correction model. For a dual-frequency or multi-frequency receiver, the ionospheric delay can be eliminated by dual-frequency IF linear combination, so the standard deviation of ionospheric delay can be set to σ ion = 0. (3) • σ noise is the receiver noise standard deviation.
Because receiver measurement noise has random characteristics, it is difficult to model and eliminate it, and the noise characteristics can only be determined by statistical analysis. Researchers have proposed many stochastic models for estimating receiver noise standard deviation. Those stochastic models are generally a function of satellite elevation angle or Signal to Noise Ratio (SNR) (Parvazi et al., 2020). Among them, the sine stochastic function model based on satellite elevation angle is widely used (Parvazi et al., 2020) where σ 0 is the nominal value of measurement noise standard deviation. For example, σ L 0 = 0.003 m for GPS P code pseudorange, and σ P 0 = 0.3 m for carrier phase observable. • σ MP is the standard deviation of multipath uncertainty. The multipath error is usually treated as a noise, and an empiric stochastic function model based on the satellite elevation angle similar to the receiver measurement noise is used to estimate the standard deviation of the multipath error (Seepersad & Bisnath, 2015).
The above error model is based on the empirical model, which is simple in its implementation. However, in a complex environment, the empirical model can not reflect the real error level due to the mismatching of error model and the user environment (Dovis et al., 2015). In order to improve the reliability of the estimated Global Navigation Satellite System (GNSS) observation accuracy and the performance of the navigation system in dynamic and complex environments, we propose a resilient adjustment method to estimate the pseudorange variance by constructing a real-time estimation model.

Resilient adjustment method to weigh pseudorange observations
Because multipath error and measurement noise have random characteristics, especially multipath error is much dependent of the user environment, it is usually difficult to accurately estimate the variance using empirical stochastic models, which ultimately reduces the positioning performance of GNSS. Compared with carrier phase observations, pseudorange observations are more severely affected by multipath error and measurement noise. Accurate construction of pseudorange variance is the key to improving the convergence time and positioning quality of precise point positioning (Seepersad & Bisnath, 2015). Therefore, this paper will mainly focus on the research on resilient adjustment methods to weigh pseudorange observations. Taking the pseudorange multipath error and measurement noise together as a source of error, Eq. (4) can be rewritten as where σ mp&noise,i is estimated in real time based on pseudorange multipath error and measurement noise sequences.

Pseudorange multipath error and measurement noise estimation model
In order to extract the multipath error and measurement noise from pseudorange observations, a Code-Minus-Carrier (CMC) observable is adopted (Blanco & Haag, 2011) The CMC observations constructed by the above formula also include two times of ionospheric delay and carrier phase multipath error and measurement noise, which in general are much smaller than the pseudorange multith error and measurement noise d i mp ≪ d i MP and ς ≪ ε , so d i mp and ς can be neglected. There are two ways to handle the ionospheric divergence. For a single frequency receiver, the ionospheric term can be eliminated by fitting and removing a first or second order curve from the CMC observable (Braasch, 2017). For a dual-frequency receiver, the ionospheric delay can be estimated through the dual-frequency carrier phase linear combination (Seepersad & Bisnath, 2015): Subtracting Eq. (9) from Eq. (8) and ignoring carrier phase multipath error mp and measurement noise ς , the dual-frequen.y ionospheric delay corrected CMC observable on frequency f i can be obtained as where d i MP represents the pseudorange multipath error and measurement noise, and B i is the carrier phase ambiguities.
In addition to pseudorange multipath error and measurement noise, the CMC observation model in (10) contains a constant component associated with phase ambiguities. Assuming that no cycles slip occurs in the carrier phase observations and the mean of multipath error and measurement noise is zero, the pseudorange multipath error and measurement noise can be estimated by (Seepersad & Bisnath, 2015;Zheng & Guo, 2016) where s i,k CMC is the average of s i,t CMC from the time t = k 0 0 of the first epoch or the cycle slip occurs to current time t = k, which is the estimation of the carrier phase ambiguities B i and removed to eliminate the constant component.

Real-time estimation model for pseudorange variance
Assuming that the multipath error and measurement noise obey the Gaussian normal distribution with zero mean, the pseudorange variance at the current epoch can be estimated based on the previous multipath error and measurement noise sequences (Wang et al., 2013) where N is the number of epochs.

Variance inflation processing
In general, the variance estimated by (12) is very conservative, which cannot be used directly as the pseudorange variance. So, it needs to be enlarged to reflect the pseudorange measurement error more accurately. Considering that the multipath error and measurement noise are related to the satellite elevation angle, the following model based on satellite elevation angle is adopted where σ 2 mp&noise,k is the enlarged variance of pseudorange multipath error and measurement noise, and s is the amplification coefficient, which can be set to 2.0-3.0 (Yang & Gao, 2006).

Variance robust estimation
Since the multipath error is very susceptible to the surrounding environment, in dynamic situations, the changes in the surrounding environment sometimes cause sudden changes in multipath errors. In order to reduce the influence of the multipath error mutation on the σ mp&noise , it needs to be detected and processed.
Taking the mutation of the multipath error as an abnormal disturbance, the following hypothesis test is constructed (Zhang et al., 2019a, b) where r is test statistic, H 0 is the null hypothesis that there is no multipath and H 1 is the altervative hypothesis that there is a multipath abnormal disturbance. The multipath abnormal disturbance is detected by comparing the magnitude relationship between the test statistic r and the test threshold T. Under normal circumstances, the measurement variances at adjacent epochs should be similar. So, in order to simplify the processing, the multipath error and measurement noise can be simply regarded as a Gaussian normal distribution with a mean of zero and a variance of σ 2 mp&noise , and the test statistics can be formulated as where a is the significance level. The detection threshold T is determined by the significance level a and can be obtained from the standard normal distribution table. According to the above detection method, if multipath abnormal disturbance is detected, the Huber weight function model is used to adjust the variance of pseudorange multipath error and measurement noise as where σ 2 mp&noise,k is the updated variance of multipath errors and noise, and T 0 is the threshold, which can be set to T 0 = 1.0 (Yang & Gao, 2006).
In addition, when a cycle slip occurs in the carrier phase observations or the signal is blocked, the constant component associated with phase ambiguities need to be recalculated. At this time, the d MP is inaccurate in general. In order to avoid the influence of inaccurate d MP on the variance estimation σ 2 mp&noise , the estimation of d MP at the time of the cycle slip is not involved in the calculation. The processing of real-time estimation of pseudorange variance is shown in Fig. 1. Figure 2 illustrates an example of the results of pseudorange variance estimation of σ 2 mp&noise , σ 2 mp&noise with inflation processing, and σ 2 mp&noise with inflation processing and robust estimation. The results show that the pseudorange variance estimation obtained by the proposed real-time estimation model is more consistent with the contour of pseudorange error. Since the proposed real-time estimation model takes into account the influence of the surrounding environment on measurement noise and multipath errors, it can obtain a more accurate estimation of the pseudorange variance and provide better estimation of the pseudorange error level than the empirical stochastic model.

Experiment and analysis
In this section, two tests were conducted using a set of static data and a set of ship-borne dynamic data to evaluate the performance of the real-time estimation model. In each test, we first analyzed and compared the estimation accuracy of the real-time estimation model and the empirical stochastic model. To this end, the following three weighting schemes for pseudorange observations are designed: Scheme 1: The empirical model is adopted, the parameter σ 0 in (6) is set to a general value, σ 0 = 0.3 .
(the general empirical model) Scheme 2: The empirical model is adopted, the parameter σ 0 in (6) is set to the optimal value, which is obtained through post-processing analysis. (the optimal empirical model) Scheme 3: The real-time estimation m.del is adopted, and the amplification coefficient is set to s = 3.0 . (the real-time estimation model) The performance of the real-time estimation model was analyzed from the perspective of the positioning domain by comparing the positioning accuracy and convergence time of kinematic PPP. The positioning accuracy is described by the Root Mean Square (RMS) of positioning errors, which is calculated by n where x and x ref are the state estimate and reference, and n is the number of epochs. The convergence time is measured as the time until the positioning error in the three directions of East (E), North (N) and Up (U) directions are all less than 0.5 m, which continues for 1 h. The detailed processing strategies of PPP are listed in Table 1.

Tests with static data
The static test data was collected from 00:00:00 to 08:23:16 on Jan 07, 2019 (GPS time) by the C200-AT high-precise GNSS receiver produced by BDStar Navigation company at a sampling frequency of 1 Hz. The receiver antenna was placed on the top of a building, and its surrounding environment is shown in Fig. 3. The coordinates of the station is (− 2 610 957.4 m, 4 232 178.5 m, 3 980 666.8 m) in Earth-Centered Earth-Fixed (ECEF) coordinate system, which is calculated by Canadian Spatial Reference System (CSRS) PPP (Tétreault et al., 2005). The sky plot of visible GPS satellites is shown in Fig. 4, and the number of visible satellites in the test period is 15. It should be noted that due to the long connecting cable (about 50 m) between the receiver and antenna, and the high-rise buildings around the antenna, and the serious multipath errors, the quality of the static observation data is not very good. For this test, the parameter σ 0 of the empirical model in Scheme 2 is set to the optimal value of 1.5 m. Figure 5 shows the RMS of pseudorange measurement residuals and the standard deviation estimates at each Relativity effect Corrected by model (Pan, 2018) epoch. Figure 6 illustrates the pseudorange measurement residuals and standard deviation for each observed satellite. In figures, the blue curve represents the pseudorange measurement residuals. The red, green, and black curves represent the pseudorange Standard Deviations (SD) estimated by the general empirical model in Scheme 1, the optimal empirical model in Scheme 2, and the realtime estimation model in Scheme 3, respectively. The statistical results for each satellite are listed in Table 2. It can be seen from Figs. 5, 6, and Table 2 that the standard deviation of pseudorange estimated by the general empirical model in Scheme 1 is very conservative, and the average standard deviation is about 1.19 m, which cannot reflect the true error level of the pseudorange observations. The average pseudorange standard deviation calculated by the optimal empirical model in Scheme 2 and the real-time estimation model in Scheme 3 are 6.69 m and 8.01 m, respectively. The standard deviation of pseudorange obtained by these two models are very close to each other, and they can obtain more accurate evaluation results. The statistical results in Table 2 show that the average pseudorange standard deviation estimated by the optimal empirical model and the real-time estimation model is about 2 to 3 times the RMS of pseudorange measurement residuals. In addition, it can be seen from the graph of the pseudorange measurement residuals for each satellite in Fig. 6 that the pseudorange measurement residuals are correlated with satellite elevation angles, and the measurement residuals decrease as the satellite elevation angles increase.  As can be seen from Figs. 7, 8, 9, the convergence time of PPP using the optimal empirical model in scheme 2 and the real-time estimation model in scheme 3 are better than that of the general empirical model in scheme 1, and the positioning results of PPP using the optimal empirical model and the real-time estimation model are more accurate and stable. It can be also seen from the statistical results in Table 3 that the positioning accuracy of PPP using the general empirical model in scheme 1 is the worst, and the convergence can not be reached in the whole process. Compared with the general empirical model in scheme 1, the positioning performance of PPP using the optimal empirical model in scheme 2 and the real-time estimation model in scheme 3 are significantly improved.

Evaluation of Positioning performance
The statistical results in Table 3 show that compared with the general empirical mode in scheme 1, the RMSs of PPP using the optimal empirical model in scheme 2 are respectively, and the convergence time is reduced to 1 964 s. By comparing the positioning results of PPP with scheme 2 and scheme 3, it is found that the positioning performance of PPP with scheme 3 is slightly better than that with Scheme 2. Compared with the optimal empirical model, the real-time estimation model reduces the RMSs of PPP by 0.078, 0.038, and 0.032 m in E, N, U directions respectively, and the convergence time by 2123 s.
In summary, both the optimal empirical model and the real-time estimation model can obtain more accurate pseudorange variance estimates, and therefore can achieve better positioning performance than the general empirical model. Compared with optimal empirical model, the real-time estimation model can adaptively adjust the pseudorange standard deviation estimates according to the changes in the measurement error level at each epoch, so the real-time estimation model can

Tests with ship-borne dynamic data
To further assess the performance of the real-time estimation model in dynamic environment, a ship-borne test was carried out on January 1, 2019, which lasted for about 3.5 h. The test data was collected by the NovAtel SPAN PwrPak7D-E1 receiver, and the sampling frequency is 1 Hz. Figure 10 shows the reference trajectory, which is calculated by the NovAtel Inertial Explore 8.80 post-processing software through the carrier phase Differential Global Navigation Satellite System (DGNSS) positioning mode, with centimeter-level positioning accuracy. Figure 11 shows the GPS satellite sky plot, the number of visible satellites during the test period is 15. In this test, the parameter σ 0 in Scheme 2 is set to the optimal value of 2.5 m. Figure 12 shows the RMS of pseudorange measurement residuals and the standard deviation estimates at each epoch. Figure 13 illustrates the pseudorange measurement residuals and standard deviation for each   Table 4.

Evaluation of variance estimation
It can be seen from the statistical results in Table 4 that the RMS of pseudorange measurement residuals is about 3.63 m, while the average pseudorange standard deviation estimated by the general empirical model in Scheme 1 is about 1.26 m, which is obviously, conservative and cannot reflect the actual error level of pseudorange observations. The average pseudorange standard deviations estimated by the optimal empirical model in scheme 2 and the real-time estimation model in Scheme 3 are 8.84 m and 8.96 m, respectively. They are equivalent, about 2-3 times the RMS of pseudorange measurement residuals. It can be seen that the real-time estimation model and the optimal empirical model have the same estimation accuracy and can better reflect the pseudorange error level.

Evaluation of Positioning performance
The convergence behavior of the positioning for the ship-borne test data in E, N, U directions are illustrated in Figs. 14, 15, 16. The statistical results of positioning errors and convergence time are listed in Table 5.
It can be seen from the position errors graph in Figs. 14, 15, 16 and the statistical results of positioning performance in Table 5   In addition, compared with the model in scheme 2, the positioning accuracy of PPP with model in scheme 3 is further improved by 0.062, 0.208 and 0.146 m in E, N, U directions directions respectively, and the convergence time is reduced by 1099 s. The positioning performance of PPP with model in scheme 3 is slightly better than that with the model in scheme 2. This is because the real-time estimation model can dynamically adjust the pseudorange standard deviation according to the degree of the influence of surrounding environment on the measurement error, so it can obtain a more accurate estimation than the optimal empirical model.

Conclusions
In this paper, a resilient adjustment method to weigh pseudorange observations is proposed, which constructs a real-time estimation and inflation model for the pseudorange variance of multipath error and measurement noise to deal with the unpredictable and timevarying of pseudorange observation error in a dynamic and complex environment. Since the influence of the surrounding environment on measurement noise and multipath error is considered, the standard deviation of the pseudorange obtained by the real-time estimation model can provide a better representation of the pseudorange error level, and the positioning performance of PPP using the real-time estimation model is better than the empirical stochastic model. The experimental results illustrated that compared with the optimal empirical stochastic model, the positioning accuracy of PPP using the realtime estimation model is improved by at least 20%, and the convergence time is reduced by at least 50%. However, it should be noted that the proposed resilient adjustment method is only suitable for dual-frequency pseudorange observations, and cannot adjust the weight of the carrier phase observations, also the influence of abnormal carrier phase observations cannot be handled. In addition, since ambiguity convergence takes some time, the pseudorange variance estimated by the proposed real-time estimation model may be inaccurate during this period. Carrier phase is another important observation in GNSS. Setting appropriate weights for carrier phase observations is of great significance to improve the positioning reliability and accuracy of PPP. Therefore, how to realize the resilient adjustment of the weight of the carrier phase observation will be one of the further works.