Performance of Ground-Based Global Navigation Satellite System Precipitable Water Vapor Retrieval in Beijing with the BeiDou B2b Service

: The accurate measurement of water vapor is essential for research about and the applications of meteorology, climatology, and hydrology. Based on the BeiDou PPP-B2b service, real-time precipitable water vapor (PWV) can be retrieved with the precise point positioning (PPP) software (XTW-PPP version 0.0). The experiment was conducted in Beijing in January 2023. Three solutions were designed with PPP using the BeiDou system only, the GPS system only, and the BeiDou-GPS combined solution. Real-time PWVs for the three solutions were validated with the ERA5 reanalysis data. Between the PWV values from the single BeiDou and ERA5, there was a bias of 0.7 mm and an RMSE of 1.8 mm. For the GPS case, the bias was 0.73 mm and the RMSE was 1.97 mm. The biases were less than 1 mm and RMSEs were less than 2 mm. Both the BeiDou and the GPS processing performed very well. But little improvement was found for the BeiDou-GPS combined solution, compared with the BeiDou system-only and the GPS system-only solution. This may be due to the poor handling of two different kinds of errors for the GPS and the BeiDou systems in our PPP software. A better PWV estimation with the two systems is to estimate PWV with a single system at the first step and then obtain the optimization by Bayesian model averaging.


Introduction
Water vapor is an important variable constituent of the global atmosphere.Accurate measurement of water vapor is essential for meteorology, hydrology, and climatology for weather forecasting, disaster warning, and climate prediction [1][2][3].Many observation techniques have been employed for the measurement of water vapor.The traditional radiosonde is expensive, and normally makes measurement twice a day with a low temporal resolution.Microwave radiometers and satellite remote sensing often have poor accuracy in cloudy or rainy weather.However, the emergence of GPS-based (Global Positioning System, GPS) PWV retrieval has gained increasing attention from scientists due to its all-weather capability, high accuracy, and cost-effectiveness [4,5].Many meteorological agencies have used the technique in ground GNSS meteorological network solutions in the US, Japan, Germany, and China, to name a few [6][7][8].High-quality water vapor observation from the GNSS (the global navigation satellite system, GNSS) has been validated to have an accuracy of about 2 mm.The GNSS observations can match those of the radiosonde or be even better [9,10].The scanning radiometer VISSR-2 (stretched visible and infrared spin scan radiometer) on-board FY-2 (the Fengyun satellite series 2) can also provide high spatiotemporal resolution atmospheric PWV products.But it only works on clear skies and the root mean square error (RMSE) is high, reaching over 4 mm compared with radiosonde [11].
Double-difference (DD) network solution and precise point positioning (PPP) are two widely used approaches to process ground-based GNSS measurements.The DD solution is a basic method used to calculate precise positions or satellite orbits and other parameters from the global or regional GNSS network [12].On the contrary, the PPP solution uses fixed satellite orbits and clocks and assumes that the mathematical models are consistent with those applied in the processing of the global reference network used to derive the satellite products [13].In the DD algorithm, receiver clock biases are eliminated, and the integer ambiguity resolution is easy to define.While in the PPP solution, data screening for outliers is required and the cycle slip detection and editing is more challenging.Thus, the DD performs better in terms of the system's stability and has been more widely used than the PPP in global and regional GNSS PWV analysis in recent years.
The DD network solution requires that all stations transfer their GNSS raw observations to the center for centralized processing.It is challenging work to obtain a high tempo-resolution of PWV with a short time latency.With the high reliability and availability of satellite orbit and clock products in recent years, the PPP technique much progress has been in terms of its accuracy and performance.It demonstrates the advantages of computational efficiency and timeliness over the DD network solution in data processing.With timely and high tempo-resolution PWV, bouts of severe weather in summer could be detected and targeted.
The BeiDou Navigation Satellite System (BDS-3) provides five public service signals of B1I, B1C, B2a, B2b, and B3I in three frequency bands of B1 (central frequency 1575.42MHz), B2 (1176.45MHz), and B3 (1268.52MHz).Since 2020, China's BeiDou has begun to provide an initial real-time precise point positioning (PPP) service via B2b signals for the Asia-Pacific region.The service uses three GEO satellites (C59-C61) to broadcast corrections to the orbit and clock offset for the BeiDou and GPS navigation satellites [14].

Experiment Description
In order to test the performance of real-time GNSS PWV retrieval based on the BeiDou B2b service, a temporary BeiDou GNSS meteorological station, shown in Figure 1, was set up at the New Technology Base of Chinese Academy of Sciences, Dengzhuang South Road, Haidian District, Beijing on 10 January 2023.The station was located at 116.16 • E in longitude, 40.04 • N latitude, and 94.78 m above sea level in a very good observational environment with little line-of-sight obstruction.Configured with a BeiDou GNSS receiver AirPro, a choke ring antenna HX-CGX611A, a set of automatic weather station WTS400, and a set of BeiDou data transmission terminals, the observation lasted for 17 days.

The B2b Real-Time Precise Ephemeris of the BEIDOU and GPS Satellites
Messages regarding the correction of the satellite orbit, clock offset, differential code bias (DCB), and user range accuracy index (URAI) are broadcasted in the PPP-B2b signal.The orbit corrections and clock offset are used to correct the broadcast ephemeris [31].
The orbit correction of a satellite is written as Remote Sens. 2024, 16, 2902 3 of 10 where δρ is the position correction, and the subscripts r, t, and a denote the radial, the tangential, and the axial directions, respectively.The superscript T is the transpose of the vector.The precise satellite position X S C can be corrected through the broadcast ephemeris X S brdc with the following equation: where e r =

The B2b Real-Time Precise Ephemeris of the BEIDOU and GPS Satellites
Messages regarding the correction of the satellite orbit, clock offset, differential code bias (DCB), and user range accuracy index (URAI) are broadcasted in the PPP-B2b signal.The orbit corrections and clock offset are used to correct the broadcast ephemeris [31].
The orbit correction of a satellite is written as where δρ is the position correction, and the subscripts r, t, and a denote the radial, the tangential, and the axial directions, respectively.The superscript T is the transpose of the vector.The precise satellite position    can be corrected through the broadcast ephemeris    with the following equation： where ⃗⃗ is the vector of the satellite's position, and   ⃗⃗⃗ is the vector of the satellite's speed.The satellite clock is corrected by this equation: where  c  is the corrected clock offset,  brdc  is the clock offset of broadcast ephemeris,   is the clock offset correction, and c is the speed of light in a vacuum.
B2b correction message can be matched through two steps.In the first step, the IODN (Issue of Data for Navigation) is used to select the corresponding broadcast ephemeris.For BDS-3 and GPS, the broadcast ephemeris here are navigation messages CNAV1 and LNAV1, respectively.In the second step, the orbit and the clock corrections are matched by IOD CORR (Issue of Data for Correction).
Since the B2b correction messages can be delayed by propagation, the validity of the correction must be ensured.The nominal validity delays for the orbit correction, the clock offset correction, the DCB, and the URAI are 96 s, 12 s, 86,400 s, and 96 s, respectively.The satellite clock is corrected by this equation:

Real-Time PPP for ZTD with BeiDou B2b Ephemeris
where dt s c is the corrected clock offset, dt s brdc is the clock offset of broadcast ephemeris, δτ s is the clock offset correction, and c is the speed of light in a vacuum.
B2b correction message can be matched through two steps.In the first step, the IODN (Issue of Data for Navigation) is used to select the corresponding broadcast ephemeris.For BDS-3 and GPS, the broadcast ephemeris here are navigation messages CNAV1 and LNAV1, respectively.In the second step, the orbit and the clock corrections are matched by IOD CORR (Issue of Data for Correction).
Since the B2b correction messages can be delayed by propagation, the validity of the correction must be ensured.The nominal validity delays for the orbit correction, the clock offset correction, the DCB, and the URAI are 96 s, 12 s, 86,400 s, and 96 s, respectively.

Real-Time PPP for ZTD with BeiDou B2b Ephemeris
The PPP model is a dual frequency ionosphere-free linear combination (IFLC) model, in which the pseudorange and carrier observations can be expressed as where DCB s ,j and DCB s ,k can be directly used to correct the IFLC pseudorange observation because the corrections of the DCB broadcast by PPP-B2b contain the pseudorange hardware delays in the j/k frequency relative to the B3I frequency.Unlike BDS, the broadcast ephemeris clock offset of the GPS is estimated by the IFLC observations of L1 and L2, so the users do not need to correct DCB when using L1 and L2 dual frequency for positioning.
The real-time precise orbit and clock offset are brought into Equation ( 6) and then linearized to obtain the following observation equation: Then, the Kalman filter method can be used for parameter estimation, and Zenith Total Delay (ZTD) can be acquired.Table 1 shows the strategy deployed in the PPP software.
Then, ZWD can be obtained by subtracting ZHD from the ZTD.It has been proved that the ZWD is related to PWV as follows [4], Π is a mapping function that can be calculated using the following equation [4]: where T m is the weighted mean temperature, , P w is the pressure of the water vapor, and T is the Kelvin temperature; T m is often calculated empirically with the surface temperature Ts.T m = a + b•T s [4].

Radiosonde and ERA5 Reanalysis Data
There is an L-band radiosonde station (#54511) in Beijing that carried out observations twice a day to obtain the meteorological profiles of temperature, pressure, humidity, and wind with a resolution of every minute.The observation errors are ±0.3 • C for temperature, ±0.5 h Pa for pressure, and ±5% for relative humidity.The radiosonde station is located at a longitude of 116.46 • E and a latitude of 39.81 • N, about 40 km away from the BeiDou GNSS meteorological station.
ERA5 is the latest reanalysis dataset from ECMWF's fifth generation [33], with finer spatial grids (79 km to 31 km), a higher temporal resolution (3-6 h to 1 h), a higher vertical resolution (60-137), new numerical weather prediction (NWP) models, and an increased assimilation data volume compared to its predecessor ERA-Interim [34].The dataset covers the period from 1950 to the present.As a three-dimensional atmospheric parameter reanalysis dataset providing a global hourly resolution, ERA5 serves as a good baseline for comparison and validation with other remote sensing observations.The precipitable water vapor can be calculated with the same integration method as for radiosonde.

The PWV Accuracy Evaluation
When it comes to assessing the accuracy of PWV, three metrics, the average deviation (bias), the root mean square error (RMSE), and the correlation coefficient (r), are often used for the validation.Their formulas are written as follows: where Y i , Y i are the values and the mean values of PWV, respectively, and n represents the sample number.

Validation Data for the GNSS PWV
Radiosondes are often used for the comparison and validation of GNSS PWV measurements [35].However, radiosonde data as a reference also have the disadvantage that the profile is not vertical at the specific site, as the sounding balloon moves due to the wind.Furthermore, it takes about 3 h to finish a profiling, which means it has a poor temporal resolution and timeliness.Although it is often used to validate the GNSS PWV traditionally, the ERA5 is more appropriate for time and space matching.Meteorological reanalysis uses advanced numerical prediction models and data assimilation systems to integrate model forecasts and historical observation data to obtain long-series historical weather data with rich variables, complete spatial coverage, and stable time uniformity.It enables meteorological agencies to obtain more detailed and timely data for weather analysis and study.The ERA5 reanalysis is currently one of the best datasets for weather and climate research, though errors may exist in some regions and some times [36,37].
Figure 2 shows the PWV values calculated using the radiosonde and ERA5 during the experiment.Between the radiosonde and ERA5, the bias of PWV is −0.27 mm, and the root mean square error is 0.41 mm.In terms of relative values, the relative bias is −5% and the relative RMSE is 9%.This complies with the relative humidity error of the radiosonde, which indicates that ERA5 has a very good analysis quality and ERA5 should be a reliable dataset for PWV validation.

Validation Data for the GNSS PWV
Radiosondes are often used for the comparison and validation of GNSS PWV measurements [35].However, radiosonde data as a reference also have the disadvantage that the profile is not vertical at the specific site, as the sounding balloon moves due to the wind.Furthermore, it takes about 3 h to finish a profiling, which means it has a poor temporal resolution and timeliness.Although it is often used to validate the GNSS PWV traditionally, the ERA5 is more appropriate for time and space matching.Meteorological reanalysis uses advanced numerical prediction models and data assimilation systems to integrate model forecasts and historical observation data to obtain long-series historical weather data with rich variables, complete spatial coverage, and stable time uniformity.It enables meteorological agencies to obtain more detailed and timely data for weather analysis and study.The ERA5 reanalysis is currently one of the best datasets for weather and climate research, though errors may exist in some regions and some times [36,37].
Figure 2 shows the PWV values calculated using the radiosonde and ERA5 during the experiment.Between the radiosonde and ERA5, the bias of PWV is −0.27 mm, and the root mean square error is 0.41 mm.In terms of relative values, the relative bias is −5% and the relative RMSE is 9%.This complies with the relative humidity error of the radiosonde, which indicates that ERA5 has a very good analysis quality and ERA5 should be a reliable dataset for PWV validation.

High Tempo-Resolutional PWV with the PPP Solution Based on B2b Service
In the China Meteorological Administration, the traditional network processing method is used by the GAMIT software (version 10.50).The zenith tropospheric delay and

High Tempo-Resolutional PWV with the PPP Solution Based on B2b Service
In the China Meteorological Administration, the traditional network processing method is used by the GAMIT software (version 10.50).The zenith tropospheric delay and water vapor content are estimated hourly after a latency of 50 min with the sliding window (12 h) method [38].
Figure 3 shows the contrast of the PWV values measured using the traditional network solution (black dot) and the PPP solution with BDS (black line) on 15 January 2023.The temporal resolution of the traditional network solution is 1 h.On the contrary, with the PPP solution based on BeiDou's B2b service, PWV can be obtained every 2 min, which is near real time.A very good agreement is shown for PWV values obtained using the two solutions, but more details about the variation in the PWV levels are revealed from the PPP resolution.

Assessment of PWVs Processing with Single BeiDou and Single GPS
Figure 4 shows that the PWVs processing with single BeiDou (PWV_C) and single GPS (PWV_G) contrast with the PWV calculation of ERA5 (PWV_ERA5).Since the ERA5 data have a 1 h resolution, they change more smoothly, and so many peaks and great values are lost.Between the PWV values from BeiDou and ERA alone, there is a bias of 0.7 mm and a RMSE of 1.8 mm.For the GPS case, the bias is 0.73 mm and RMSE is 1.97 mm.Both the BeiDou and the GPS processing performed very well in the estimation of PWV.
The temporal resolution of the traditional network solution is 1 hour.On the contrary, with the PPP solution based on BeiDou's B2b service, PWV can be obtained every 2 min, which is near real time.A very good agreement is shown for PWV values obtained using the two solutions, but more details about the variation in the PWV levels are revealed from the PPP resolution.

Assessment of PWVs Processing with Single BeiDou and Single GPS
Figure 4 shows that the PWVs processing with single BeiDou (PWV_C) and single GPS (PWV_G) contrast with the PWV calculation of ERA5 (PWV_ERA5).Since the ERA5 data have a 1 h resolution, they change more smoothly, and so many peaks and great values are lost.Between the PWV values from BeiDou and ERA alone, there is a bias of 0.7 mm and a RMSE of 1.8 mm.For the GPS case, the bias is 0.73 mm and RMSE is 1.97 mm.Both the BeiDou and the GPS processing performed very well in the estimation of PWV.The ZTDs show a very similar situation, as shown in Figure 5. Referenced with the ERA5 calculation (ZTD_ERA5), for the BeiDou processing (ZTD_C), the bias is 9.27 mm

Assessment of PWVs Processing with Single BeiDou and Single GPS
Figure 4 shows that the PWVs processing with single BeiDou (PWV_C) and single GPS (PWV_G) contrast with the PWV calculation of ERA5 (PWV_ERA5).Since the ERA5 data have a 1 h resolution, they change more smoothly, and so many peaks and great values are lost.Between the PWV values from BeiDou and ERA alone, there is a bias of 0.7 mm and a RMSE of 1.8 mm.For the GPS case, the bias is 0.73 mm and RMSE is 1.97 mm.Both the BeiDou and the GPS processing performed very well in the estimation of PWV.The ZTDs show a very similar situation, as shown in Figure 5. Referenced with the ERA5 calculation (ZTD_ERA5), for the BeiDou processing (ZTD_C), the bias is 9.27 mm

PWV from the BeiDou-GPS Combined Solution
In the PPP processing, we treated the BeiDou MEO Satellites as GPSs and obtained the BeiDou-GPS combined solution.As can be seen, the trend in Figure 6 is very similar to that in Figure 4. Referenced with the ERA5 calculation, the bias is 0.74 mm and the RMSE is 1.97 mm for the BeiDou-GPS combined processing.Little improvement was found for the assisting the BeiDou with a GPS and vice versa.This may be due to the poor handling of two different kinds of errors for the GPS and the BeiDou systems in our PPP software.

PWV from the BeiDou-GPS Combined Solution
In the PPP processing, we treated the BeiDou MEO Satellites as GPSs and obtained the BeiDou-GPS combined solution.As can be seen, the trend in Figure 6 is very similar to that in Figure 4. Referenced with the ERA5 calculation, the bias is 0.74 mm and the RMSE is 1.97 mm for the BeiDou-GPS combined processing.Little improvement was found for the assisting the BeiDou with a GPS and vice versa.This may be due to the poor handling of two different kinds of errors for the GPS and the BeiDou systems in our PPP software.

PWV from the BeiDou-GPS Combined Solution
In the PPP processing, we treated the BeiDou MEO Satellites as GPSs and obtained the BeiDou-GPS combined solution.As can be seen, the trend in Figure 6 is very similar to that in Figure 4. Referenced with the ERA5 calculation, the bias is 0.74 mm and the RMSE is 1.97 mm for the BeiDou-GPS combined processing.Little improvement was found for the assisting the BeiDou with a GPS and vice versa.This may be due to the poor handling of two different kinds of errors for the GPS and the BeiDou systems in our PPP software.We optimized the PWV estimation after the PWV values had been estimated by both the GPS and BeiDou by using Bayesian model averaging [39].
where   is the fusion of _ , _ , and _ as illustrated in Section 3.3, and   and   are their fitting errors, respectively.Table 2 shows the biases and RMSEs for the GPS-and BeiDou-only solutions, the combined solution, and the fusion of PWVs, compared with PWV_ERA5.We optimized the PWV estimation after the PWV values had been estimated by both the GPS and BeiDou by using Bayesian model averaging [39].
where PWV opt is the fusion of PWV_G, PWV_C, and PWV_G as illustrated in Section 3.3, and σ c and σ G are their fitting errors, respectively.Table 2 shows the biases and RMSEs for the GPS-and BeiDou-only solutions, the combined solution, and the fusion of PWVs, compared with PWV_ERA5.

Conclusions and Discussion
Here, three solutions were designed with PPP using the BeiDou system only, the GPS system only, and the BeiDou-GPS combined system.Real-time PWV values for the three solutions were validated with the ERA5 reanalysis data.The biases are less than 1 mm and RMSEs are less than 2 mm.The BeiDou system performed very well, matching the GPS system.But little improvement was found for the assistance of GPS to BeiDou and vice versa, which might be due to the poor handling of two different kinds of errors between the GPS and the BeiDou in our PPP software.Furthermore, there were many error losses during the software processing.In the ZTD estimation, the main errors are satellite orbit and clock errors as well as model errors.In the conversion of ZTD to PWV, there are errors in the pressure and temperature observations and the calculation of the weighted mean temperature (Tm).Further improvements will be needed in the future to reduce these errors.A better PWV estimation strategy using the two systems is to first estimate the PWV values with a single system only, and then to obtain a weighted mean estimation.
The PPP-based algorithm can obtain higher temporal PWV values than the doubledifference network solution.Referenced with ERA5, the deviation and the RMSE are relatively small.Therefore, it shows an advantage over the traditional network solution for severe weather analysis and warnings.

Figure 3 .
Figure 3. PWVs from the traditional network processing method and the PPP solution.

Figure 4 .
Figure 4. Comparison of PWV values from the BeiDou solution, the GPS solution, and the ERA5.

Figure 3 .
Figure 3. PWVs from the traditional network processing method and the PPP solution.

Figure 3 .
Figure 3. PWVs from the traditional network processing method and the PPP solution.

Figure 4 .
Figure 4. Comparison of PWV values from the BeiDou solution, the GPS solution, and the ERA5.

Figure 4 .
Figure 4. Comparison of PWV values from the BeiDou solution, the GPS solution, and the ERA5.The ZTDs show a very similar situation, as shown in Figure5.Referenced with the ERA5 calculation (ZTD_ERA5), for the BeiDou processing (ZTD_C), the bias is 9.27 mm and the RMSE is 10.82 mm.For the GPS processing (ZTD_G), the bias is 11.63 mm and the RMSE is 12.68 mm.

Figure 5 .
Figure 5. ZTDs from the BeiDou, the GPS processing and calculation of ERA5.

Figure 5 .
Figure 5. ZTDs from the BeiDou, the GPS processing and calculation of ERA5.

Figure 6 .
Figure 6.PWVs from the BeiDou-GPS combined solutions and the ERA5.

Figure 6 .
Figure 6.PWVs from the BeiDou-GPS combined solutions and the ERA5.
The superscript s represents the satellite, the subscript r represents the receiver, and j and k represent the signal frequency.dt r and dt s denote the clock offset of the receiver and the satellite.T s r is the tropospheric delay; d r,IF jk and b r,IF jk are the pseudorange and carrier hardware delays experienced by the receiver.d s IF jk are the ambiguity and the wavelength.ε P and ε L are the noises in pseudorange and carrier observations.Since BDS broadcast their ephemeris based on B3I frequency, the corrected real-time satellite clock offset should include the hardware delay of B3I signal.Letδ ts = dt s + d s ,B3I , and δt r = dt r + d r,IF jk ; then, jk are observations of the pseudorange and carrier phases, respectively.

Table 2 .
Biases and RMSEs for different PWVs referenced with EAR5 calculation.