ZWD is defined as the integral of the wet refractivity along the vertical profile above the station.

In Eq. (1), N_w is the wet refractivity, e is the water vapour pressure in hectopascals (hPa), T is the temperature in Kelvin, $k_{\mathrm{2}}^{\prime }$ is 22.1 K hPa⁻¹, and k₃ is 373 900 K² hPa⁻¹ (Bevis et al., 1994). It can be seen from Eq. (1) that ZWD changes with height, vapour pressure, and temperature. The ZWD will decrease with increasing height due to the shorter integral length. With the profiles of water vapour pressure and temperature, one can obtain the accurate ZWD by the ray-tracing method. However, in practical applications (e.g. aircraft navigation and positioning, wide area augmentation), we usually use empirical models for ZWD corrections due to the unavailability of meteorological data profiles. Therefore, it is necessary to develop an empirical ZWD model with high precision. The temperature roughly decreases linearly with increasing height in the troposphere, while the change in water vapour is more variable, so the water vapour is the main determinant of vertical variation of ZWD. In the following content, we used the meteorological data profile of ERA-Interim pressure levels provided by ECMWF to analyse the vertical variation characteristics of ZWD and explore a suitable fitting function capable of describing the changes in ZWD with respect to height.

Figure 2ZWD vertical gradients profile (a) and linear fit with height below 2 km (b) at a grid point (0^∘ N, 0^∘ E) at 12:00 UTC on 1 January 2010.

Download

ERA-Interim can provide data at 00:00, 06:00, 12:00, and 18:00 UTC daily with a spatial resolution of not more than $\mathrm{0.75}{}^{\circ }\times \mathrm{0.75}{}^{\circ }$ and 37 pressure levels (Dee et al., 2011). The highest level data come from a height of approximately 50 km, covering almost the entire troposphere and stratosphere. We used the temperature, the geopotential height, and the specific humidity provided by the ERA-Interim pressure level data, and the discretised form of Eq. (1), to calculate the ZWD for each level height (Böhm and Schuh, 2013).

In Eq. (2), q is the specific humidity in grammes per gramme (g g⁻¹), P is the pressure in hPa, $k_{\mathrm{2}}^{\prime }$ and k₃ are empirical constants same as Eq. (1), and h is the geopotential height in metres. From Eq. (2), we can see that the ZWD at a specific level height is the sum of the ZWD portions in all layers above the specific level height. Figure 1 shows the water vapour pressure and ZWD profiles at a grid point (0^∘ N, 0^∘ E) at 12:00 UTC on 1 January 2010. From Fig. 1, it can be seen that the downward trend in the water vapour pressure varies significantly with height, and the decrease factor is different across different height intervals. The changes in ZWD with respect to height are similar to that of the water vapour pressure with respect to height: the decay is fastest up to a few kilometres height and slows down with increasing height; the ZWD values are close to zero after 10 km. Zhao et al. (2014) showed that about 50 % of the water vapour content is concentrated within 1.5 km of the surface and less than 10 % of the water vapour content remains above 5 km, leading to different ZWD decay rates within different height intervals. These results are basically consistent with our experiment results. Further, the derivative of the ZWD with respect to height (i.e. ZWD vertical gradient) is analysed to better understand the characteristic of the ZWD vertical distributions. Figure 2a shows the variation of ZWD vertical gradients with respect to height at the same grid point to Fig. 1. From Fig. 2a, it can be seen that the trends in ZWD vertical gradients at different height intervals are clearly different. Specifically, the linear fit of the ZWD gradients with height below 2 km shows a great agreement with an R² value of 0.99 (Fig. 2b). Thus we can come to a conclusion: ZWD gradients roughly change linearly below 2 km. From 2 to 5 km, and 5 to 10 km, the ZWD gradients vary non-linearly.

Figure 3ZWD gradients profiles at grid points in different latitude bands (12:00 UTC, 1 January 2010).

Download

Figure 4ZWD gradients profiles at grid points in different latitude bands (12:00 UTC, 1 July 2010).

Download

Figure 3 shows the ZWD vertical gradients with respect to height at grid points in different latitude bands. Figure 4 shows the similar ZWD vertical gradients as Fig. 3 but for a different season. The variations are similar to those in Fig. 2a, which show trend changes at about 2 and 5 km. It is worth noting that the ZWD gradients at low latitudes are much larger and water vapour is more variable than at high latitudes, resulting from the fact that the water vapour at low latitude is more variable. In addition, the ZWD gradient trends in the Southern Hemisphere are significant. In contrast, the ZWD gradients in the Northern Hemisphere are slightly complicated with respect to height: the reason for this may be that the Southern Hemisphere is mostly oceanic while the Northern Hemisphere has many sea coasts. The terrain complexity in the Northern Hemisphere contributes to the disturbances in the ZWD gradient in specific areas. According to the vertical variation characteristics of ZWD, we divided the troposphere into three height intervals (below 2 km, 2 to 5 km, and 5 to 10 km) and assumed 10 km to be the empirical tropopause beyond which the ZWD is assumed to be zero. For ZWD fitting with respect to height, TropGrid2 and GPT2w use exponential functions, while some scholars have also used a polynomial to describe the tropospheric delay with respect to height (Song et al., 2011). We used both polynomial and exponential functions to fit the variation trend of the ZWD with respect to height in the three selected intervals, respectively. The results showed that the quadratic polynomial used under 2 km and the exponential functions between 2 and 5 km and between 5 and 10 km gave the best fits. The combination of the quadratic polynomial and exponential functions for different height intervals are termed piecewise height functions. Table 1 summarises the global fitting statistics of different fit functions, demonstrating the superiority of piecewise height functions to the single polynomial function and single exponential function used for the whole troposphere.

Figure 5Decadal time series and cycle fitting results of function coefficients z₁ (a), z₂ (b), and z₃ (c) at a grid point (0^∘ N, 0^∘ E) from 2001 to 2010.

Download

Table 1Fitting RMSE of piecewise height functions, single quadratic polynomial function, and single exponential function (unit: mm).

Download Print Version | Download XLSX

From the above analysis of ZWD vertical variation and fitting, the piecewise height functions of the proposed HZWD model are:

In Eq. (3), B is the latitude in degrees; L is the longitude in degrees; H is the height in metres; function coefficients z₁, z₂, and z₃ can be regarded as the ZWD at the height of 0, 2, and 5 km, respectively. We used the monthly mean profiles of ERA-Interim pressure levels from 2001 to 2010 with a horizontal resolution of $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$ for ZWD modelling. The ZWD profiles calculated for each grid point are fitted by Eq. (3) to obtain the time series of the corresponding function coefficients: z₁, a₁, a₂, z₂, β₂, z₃, and β₃. It is worth noting that the ERA-Interim-derived ZWD data indicate that the averaged ZWD values at the three height intervals (i.e. below 2, 2 to 5 km, and 5 to 10 km) are 0.126, 0.0489, and 0.0111 m, respectively. Jin et al. (2007) found that the tropospheric delay has notable seasonal variations, mainly on annual and semi-annual cycles. Song et al. (2011) and Zhao et al. (2014) considered the temporal features of function coefficients in their troposphere models. We used the 10-year time series of those coefficients obtained to analyse their temporal variations. Figure 5 shows the time series and cycle fitting results of the function coefficients z₁, z₂, and z₃ at grid point (0^∘ N, 0^∘ E). Figure 5 shows that the time series of the function coefficients z₁, z₂, and z₃ have a significant characteristic annual cycle, and the semi-annual cycle is small but nevertheless evident.

Therefore, taking the annual and semi-annual cycles into consideration, we used Eq. (4) to fit the function coefficients derived from Eq. (3) to temporal parameters for each grid point (Böhm et al., 2015).

In Eq. (4), c₀ is the annual mean, c₁ and b₁ are the annual cycle parameters, c₂ and b₂ are the semi-annual cycle parameters, and doy is the day of the year. The fittings show that the annual means and annual and semi-annual amplitudes of z₁, z₂, and z₃ are distinct. For instance, the cycle fitting results at a grid (0^∘ N, 0^∘ E) (Fig. 5) indicate that the temporal parameters (i.e. c₀, c₁, b₁, c₂, and b₂) of z₁ are 0.2911, 0.0237, 0.0312, −0.0006, and −0.0227 m, respectively; the temporal parameters of z₂ are 0.1215, 0.0118, 0.0203, 0.0004, and −0.0146 m, respectively; the temporal parameters of z₃ are 0.0255, 0.00031, 0.0070, −0.0019, and −0.0044 m, respectively. It should be noted that the fitting results of coefficients a₂, β₂, and β₃ show that all their annual means and annual and semi-annual amplitudes are small. However, below 2 km, the lack of cycle terms in a₂ would cause centimetre-level errors in the ZWD estimates, so these terms have been retained. For β₂ and β₃, ZWD itself is small at heights above 2 km, so the annual mean suffices for a desirable ZWD estimate. The experiment reveals that the loss of accuracy due to the lack of annual and semi-annual terms in β₂ and β₃ for the ZWD estimates is less than 0.1 mm. Therefore, only the annual means are retained for these two coefficients.

Figure 6 shows the global distributions of annual means of model coefficients z₁, z₂, and z₃. From Fig. 6 we can see that the extremum of ZWD annual means at 0 m height occur near the Equator and the maximum exceeds 0.36 m. The ZWD annual means decrease with increasing latitude. The distributions of ZWD annual means at 2 and 5 km heights are similar to that at 0 m, but the areas with the large values near the Equator decrease in extent and the ZWD distributions tend to be uniform, indicating that the water vapour content near the Equator is greater than that in other regions, and the ZWD value is also larger in low-altitude regions. As the height increases, the difference in water vapour content or ZWD between the Equator and other areas begins to decrease, but remains significant. Overall, there are some differences in the ZWD distribution at different heights, and it is necessary to model the spatio-temporal variations of ZWD at different heights.

Figure 6Global distributions of annual means of HZWD model coefficients z₁ (a), z₂ (b), and z₃ (c).

Download

After the fitting processes involving Eqs. (3) and (4), the global ZWD model HZWD, using piecewise height functions, is established. The spatial resolution of the HZWD model is $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$. Each grid point contains seven primary coefficients: z₁, a₁, a₂, z₂, β₂, z₃, and β₂. Among these coefficients, z₁, z₂, z₃, a₁, and a₂ are further expressed by Eq. (4) with five temporal parameters, respectively. Therefore, there are 27 parameters for each grid point and total 68 094 parameters for the HZWD model. As a comparison, the GPT2w model has a number of 77 760 parameters, which is 14 % more than that of the HZWD model. It is worth noting that the UNB3m model only has 50 parameters due to its coarse spatio-temporal resolution. When the HZWD model is applied, the four grid points surrounding the station are determined according to the horizontal position (latitude and longitude) of the station, and then the model coefficients of the corresponding height intervals at the four selected points are calculated according to Eq. (4). The ZWD of the four grid points are extrapolated to the station height by using Eq. (3), and finally the ZWD at the station location is obtained by using bilinear interpolation. The HZWD model only needs time, latitude, longitude, and height as input parameters. It can calculate ZWD without meteorological data and can provide wet delay correction products for navigation and positioning at different heights.

To test the precision of the HZWD model and analyse the model correction performance compared to other troposphere models, we used the ERA-Interim pressure level data and radiosonde data from the year 2015 as external data sources and compared the results with the commonly used models UNB3m and GPT2w. The parameters used for the validation are bias and root mean square error (RMSE), expressed as follows.

In Eqs. (5) and (6), ZWD${}_i^M$ is the value estimated by the HZWD model developed in this study and ZWD${}_i^{\mathrm{0}}$ is the reference value.

For the UNB3m model, the ZWD at mean sea level (m.s.l.) is first calculated, then a vertical correction is applied to transform the ZWD to the target height. The formulae are as follows (Leandro et al., 2006):

where e₀, T₀, and ZWD₀ are the water vapour pressure, temperature, and ZWD at MSL, respectively; R_d is the specific gas constant for dry air (278.054 J kg⁻¹ K⁻¹); γ and λ are the temperature lapse rate and water vapour decrease factor, respectively; g_m is the gravity acceleration at the mass centre of the vertical column of the atmosphere and can be computed with geodetic latitude φ and height h by

and T_m is the weighted mean temperature computed by

For the GPT2w model, the modelled meteorological parameters at the four grid points surrounding the target location are extrapolated vertically to the desired height, and then the Askne and Nordius formula (10) is used to calculate the wet delays at those base points: finally the wet delays are interpolated to the observation site in horizontal direction to get the target ZWD.

In GPT2w model, T_m is an empirical parameter modelled with seasonal components and g_m is simplified to a constant (9.80665 m s⁻²). It should be noted that the GPT2w model provides both $\mathrm{1}{}^{\circ }\times \mathrm{1}{}^{\circ }$ and $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$ resolution versions. Since the horizontal resolution of the HZWD model is $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$, we used the GPT2w model with the same resolution for validation.

4.1 Validation with ECMWF data

Modelling of the HZWD model is based on the monthly mean profiles of ERA-Interim pressure level data from 2001 to 2010, while we used the ERA-Interim pressure level data with the full time resolution of 6 h in 2015 for the model validation. This is to validate the model performance on the daily scale. Regarding the ZWD profiles calculated from these data as reference values, we calculated the global annual average bias and RMSE of the ZWD for three models (HZWD, GPT2w, and UNB3m) within the three height intervals: below 2 km, 2 to 5 km, and 5 to 10 km (Table 2).

Table 2Error statistics for the three models compared to the 2015 ECMWF data (unit: mm).

Download Print Version | Download XLSX

From Table 2, it can be seen that the HZWD model is the most accurate model across all three intervals, followed by the GPT2w model, and the UNB3m model has the worst performance. The annual average biases of the HZWD model are lower than that of the GPT2w model and the UNB3m model, except below 2 km. Compared with the RMSEs in the GPT2w model, those of the HZWD model are decreased by 1.4, 0.9, and 1.2 mm within the three height intervals, corresponding to improvements of about 6 %, 6 %, and 32 %, respectively. The improvements of the HZWD model over GPT2w model will result in precision improvements of 2.8, 1.8, and 2.4 mm respectively in height estimates in real-time aircraft positioning. The correction performance improvement from 5 to 10 km height is particularly evident. Figure 7a shows the ECMWF ZWD profile and the ZWD profiles of the three models at 12:00 UTC on 1 January 2015 at a representative grid point (0^∘ N, 20^∘ E). More clearly, Fig. 7b shows the differences between the ZWD profiles of the three models and ECMWF ZWD profile at different heights. It can be seen that HZWD is the most stable model, showing the best agreement with the ECMWF ZWD data, which is superior to both the GPT2w and the UNB3m models.

Figure 7The ZWD profiles (a) of ECMWF and the three models (HZWD, GPT2w, and UNB3m) and corresponding biases (b) at a grid point (0^∘ N, 20^∘ E) at 12:00 UTC on 1 January 2015.

Download

The variation of the troposphere has a strong correlation with latitude. To analyse the correction performances of the three models in different regions around the world, we calculated the three models' errors in different latitude bands (10^∘ intervals). Figures 8 and 9 show the correction performances at different latitudes. It can be seen from Fig. 8 that the bias of the UNB3m model is basically positive in the three height intervals, indicating that its ZWD estimates are relatively large compared to the ECMWF data. Moreover, the bias in the Southern Hemisphere is significantly larger than that in the Northern Hemisphere, indicating systematic deviations in the Southern Hemisphere. Both the GPT2w model and the HZWD model have large biases in the low latitudes. The biases of the GPT2w model are positive from 2 to 5 km and 5 to 10 km height, indicating that the ZWD is overestimated by the GPT2w model with increasing height. For the HZWD model, the bias in each latitude band is relatively small with few exceptions, resulting in a global average bias close to zero (see Table 2).

Figure 8Bias comparisons between the three models (HZWD, GPT2w, and UNB3m) in different latitude bands over the year 2015.

Download

Figure 9 shows the RMSEs of the three models. It can be seen from Fig. 9 that the precision of the HZWD model is significantly better than that of the UNB3m model across the three height intervals and all latitude bands, which is better than the GPT2w model in general. The precision of the three models declines with decreasing latitude, because the active change in water vapour in these areas limits the precision of the model. Corresponding to Fig. 8, the errors in UNB3m are asymmetric: the main reason for this is that the meteorological parameters of UNB3m are interpolated from the coarse look-up table with a latitude interval of 15^∘, and UNB3m does not consider the longitudinal variations of any meteorological elements. It should be pointed out that the UNB3m model is based on the simple symmetric assumption of the Northern and Southern Hemispheres, and its modelling data source only comes from the atmospheric data collected over North America, which leads to poor precision in the Southern Hemisphere, especially in the high latitudes thereof.

Figure 9RMSE comparisons between the three models (HZWD, GPT2w, and UNB3m) in different latitude bands over the year 2015.

Download

Table 3Error statistics for the three models validated by 2015 radiosonde data (unit: mm).

Download Print Version | Download XLSX

Summarising the distributions of bias and RMSE across different latitude bands, we can see that the HZWD model performs best with the ECMWF data as reference values. Compared with the models GPT2w and UNB3m, the HZWD model basically eliminates systematic error in the 5 to 10 km height interval and the correction performance is stable at all heights and regions. To investigate the model's performance over time, Fig. 10 shows the time series of RMSE for the three models at 6 h intervals throughout the year 2015 at grid point (0^∘ N, 20^∘ E). We can see that the HZWD model has the best overall performances within the three height intervals over the year 2015. We noticed the significantly large RMSE for all three models across all three height intervals around the doy 19 and doy 195 of 2015. This can be attributed to the sharp short-term ZWD variations in the Equator area. The short-term variations are hardly accounted for by all three models which only consider the seasonal variations of ZWD. Moreover, the GPT2w model has the worst performance from 5 to 10 km height, which is also identified by Figs. 8 and 9.

Figure 10RMSEs in ZWD estimates of the three models (HZWD, GPT2w, and UNB3m) compared to the ECMWF data over the year 2015 at grid point (0^∘ N, 20^∘ E).

Download

Figure 11Uncertainty of ZWD with respect to height at station 01241 (63.70^∘ N, 9.60^∘ E at 10 m) at 11:00 UTC on 1 January 2015.

Download

Figure 12Global distributions of bias for the three models (HZWD, GPT2w, and UNB3m) compared to 2015 radiosonde data.

Download

4.2 Validation with radiosonde data

A radiosonde is used in a sounding technique that regularly releases balloons to collect atmospheric meteorological data at different heights: it can obtain profiles of various meteorological data with high accuracy. At present, the Integrated Global Radiosonde Archive (IGRA) website (ftp://ftp.ncdc.noaa.gov/pub/data/igra/, last access: 1 November 2018) provides free downloads of global radiosonde data. We used radiosonde data from 318 stations collected in 2015 to test the HZWD model. After data pre-processing, the data with gross errors have been removed and a total of 163 671 radiosonde data epochs remained. With the provided profiles of geopotential height, temperature, and water vapour pressure, the data forms of the radiosonde data are very similar to the ECMWF pressure level data, and thus the radiosonde ZWDs can be calculated using the same method by Eq. (2). Before the validation, we conducted an assessment of the uncertainty of ZWD derived from radiosonde data. Rózsa (2014) showed that the uncertainty of ZWD is ±1.5 mm in the case of the Vaisala RS-92 radiosondes in central and eastern Europe. However, this uncertainty is only valid for the ZWD calculated from the height of the lowest layer and is limited to Europe. Using the same uncertainties of radiosonde meteorological data given by the technical specification of the radiosonde and the algorithm proposed by Rozsa (2014), we calculated the ZWD uncertainty for all heights in all radiosonde stations. Figure 11 shows the uncertainty of ZWD with respect to the height for radiosonde station 01241 located in Orland, Norway (63.70^∘ N, 9.60^∘ E at 10 m). We can see that the uncertainty of ZWD is less than ±1.5 mm near the height of 0 m and decrease quickly with increasing height. The global mean uncertainties of ZWD in all stations in the three height intervals are ±1.3, ±0.7, and ±0.2 mm, respectively, indicating the high accuracy of ZWD derived from radiosonde data.

Figure 13Global distributions of RMSE for the three models (HZWD, GPT2w, and UNB3m) compared to 2015 radiosonde data.

Download

Figure 14RMSEs in ZWD estimates of the three models (HZWD, GPT2w, and UNB3m) for radiosonde station 01241 over the year 2015.

Download

Taking the radiosonde ZWDs as reference ZWD values, we validated the ZWDs from the models HZWD, GPT2w, and UNB3m. Table 3 shows the statistical results of the three models. It can be seen from Table 3 that the HZWD model has the best overall stability of the average bias and RMSE, indicating the best precision, and the UNB3m model is the worst. Compared with the GPT2w model, the RMSEs in HZWD in the three height intervals are reduced by 0.6, 0.9, and 1.7 mm, which equates to precision improvements of 2 %, 5 %, and 33 %, respectively. Moreover, these improvements correspond to an error reduction of 1.2, 1.8, and 3.4 mm respectively in height estimates in geodetic techniques. Taking the uncertainty of radiosonde ZWD into account, the improvement of the HZWD model over GPT2w model below 2 km seems to be insignificant. Nevertheless, we can reasonably think that the ZWD predicted by HZWD is closer to true ZWD due to its smaller RMSE. It is worth noting that the bias and RMSE of the HZWD model and the GPT2w model are both larger than those of the results from ECMWF data in Table 2. The reason is that the HZWD model and the GPT2w model are based on ECMWF data, and thus the test results with radiosonde data are slightly worse than those using ECMWF data. On the contrary, the bias of the UNB3m model decreases, and the RMSE between 2 and 5 km and between 5 and 10 km are less than those in Table 2. It may be due to the fact that most of the radiosonde stations are in the Northern Hemisphere, accounting for more than 60 % (192∕318) of the total, which has a positive impact on the test results for the UNB3m model based on North American meteorological data.

Figure 12 shows the global distributions of bias for the three models within the three height intervals, and Fig. 13 shows the global distributions of RMSE for the three models. As can be seen from Fig. 12, the three models show a poorer performance in low-latitude areas than in mid- and high-latitude areas for all height intervals, similar to the results in Sect. 4.1. Within the 5 to 10 km interval, the bias of the GPT2w model is large and positive in the equatorial region, indicating that the ZWD of the GPT2w at this height is significantly overestimated, and the global bias of the UNB3m model in this height interval is positive, also indicating an overestimate of the ZWD in the UNB3m model. The bias of the HZWD model does not show obvious regional differences with respect to height, and the overall distribution of HZWD model bias has no tendency to either the positive or negative. Figure 13 further illustrates the precision of the HZWD model. The global RMSE distributions of the HZWD model are similar to that of GPT2w model below 2 km and between 2 and 5 km, but the precision of the HZWD model is slightly better. Combining this with the bias distribution of the GPT2w model in Fig. 12, the GPT2w model also has a large RMSE near the Equator in the 5 to 10 km interval, which shows that the GPT2w model is unstable at high height in low-latitude areas. The precision of the UNB3m model is poorer than that of both the HZWD and GPT2w models. Below 2 km, the UNB3m model reaches decimetre-level precision near the Equator, and even exceeds 12 cm in some areas: the distribution of north–south heterogeneity remains obvious.

These results validate the spatial stability of the precision of the HZWD model; furthermore the temporal stability of the model precision is verified next. Figure 14 shows the results of ZWD corrections of the three models for the radiosonde station 01241 for the whole of 2015. It can be seen from Fig. 14 that the HZWD model and the GPT2w model are relatively stable throughout the year, while the correction performance of the UNB3m model in 2015 is worse than those of the HZWD and GPT2w models. The probable reason for this is that the UNB3m model only takes into account the annual variations in the metrological elements with a fixed phase, resulting in precision instability throughout the year. The improved performance arising from use of the HZWD model, compared to that arising from use of the GPT2w model, is more apparent with increasing height: this shows that modelling ZWD piecewise with height can effectively approximate the real ZWD profile and improve the precision of ZWD estimation.

The ERA-Interim data can be downloaded from ECMWF (Dee et al., 2011). The radiosonde data can be downloaded from the IGRA (ftp://ftp.ncdc.noaa.gov/pub/data/igra/, last access: 1 November 2018) (NOAA, 2018).

YY and YH contributed to the conception of the study. YY contributed significantly to the data analysis and manuscript preparation. YH performed the model validation and wrote the manuscript.

The authors declare that they have no conflict of interest.