Evaluation of IWV from the numerical weather prediction WRF model with PPP GNSS processing for Bulgaria

Global Navigation Satellite Systems (GNSS) meteorology :::::::::: Meteorology : is an established operational service providing hourly updated GNSS tropospheric products to the National Meteorologic Services (NMS) in Europe. In the last decade through the ground-based GNSS network densification and new processing strategies like Precise Point Positioning (PPP), : it has become possible to obtain sub-hourly tropospheric products for monitoring severe weather events. In this work one year (January December 2013) of sub-hourly GNSS tropospheric products (Zenith Total Delay) are computed using the PPP strat5 egy for seven stations in Bulgaria. In order to take advantage of the sub-hourly GNSS data to derive Integrated Water Vapour (IWV) surface pressure and temperature with similar temporal resolution is required. As the surface observations are on 3 hourly basis the first step is to compare the surface pressure and temperature from ::: the numerical weather prediction model Weather Forecasting and Research :::::::: Research ::: and :::::::::: Forecasting : (WRF) with observations at three synoptic stations in Bulgaria. The mean difference between the two data-sets for 1) surface pressure is less than 0.5 hPa and the correlation is over 0.989, 10 2) temperature the largest mean difference is 1.1◦ C and the correlation coefficient is over 0.957 and 3) IWV mean difference is in range ::::::: between : 0.1-1.1 mm. The evaluation of WRF on annual bases shows IWV underestimation between 0.5 and 1.5 mm at five stations and overestimation at Varna and Rozhen. Varna and Rozhen have also much smaller correlation 0.9 and 0.76. The study of the monthly IWV variation shows that at those locations the GNSS IWV has unexpected drop in April and March respectively. The reason for this drop is likely problems with station raw data. At the remaining 5 stations a very good 15 agreement between GNSS and WRF is observed with high correlation during the cold part of 2013 i.e. March, October and December (0.95) and low correlation during the warm part of 2013 i.e. April to August (below 0.9). The diurnal cycle of the WRF model shows a dry bias in the range of 0.5-1.5 mm. Between 00 and 01 UTC the GNSS IWV tends to be underestimate


Introduction
The atmospheric water vapour is a key element of the hydrological cycle and participates in precipitation formation, energy transfer and atmospheric stability.Water vapour has a relatively short lifetime in the atmosphere, from one week to ten days and its complex life cycle includes vertical and horizontal transport, mixing, condensation, precipitation and evaporation.Due to its high temporal and spatial variability atmospheric water vapour is very demanding to observe.
An established method for monitoring water vapour is the radiosonde.In Europe, the radiosonde network consists of 93 stations operated by the National Meteorological Services (NMS) under the EUMETNET-EUCOS project (euc, 2016).The radiosonde provide high vertical resolution data but due to its high cost is operated only one or two times per day at 00 UTC, at 12 UTC or 00 and 12 UTC, respectively.To monitor the high temporal and spatial water vapour variability a new method was developed in the early 1990s using the Global Positioning System (GPS) signal delay.The method was called "GPS meteorology" but with the development of other GNSS, for example Glonass and Galileo, was renamed to "GNSS meteorology".A memorandum of understanding between E-GVAP (the EUMETSAT GNSS water vapour program) and EUPOS (the European Position determination System), which opens opportunities to use GNSS data for Bulgaria and South East Europe is in action since 2012.
Multi technique ::::::::::::: Multi-technique comparisons (Ning et al., 2012;Buehler et al., 2012;Van Malderen et al., 2014) have demonstrated that GNSS meteorology derived Integrated Water Vapour (IWV) has a root mean square error in the range of 0.4-0.6 mm.A number of studies compare IWV from GNSS and Numerical Weather Prediction (NWP) models in Europe.A recent study by Keernik et al. (2014) found that the HIRLAM NWP model underestimates the IWV by 59 % for values below 12 mm, and overestimates by 6-10 % for values over 25 mm.A study of the COSMO model diurnal IWV cycle over Germany (Tomassini et al., 2002), reports a systematic IWV underestimation larger than 1 mm in the model analysis between 06 and 18 UTC.For Switzerland, Guerova et al. (2003) report a good agreement between model analysis and GNSS in winter but in summer, a significant underestimation of IWV was found in the model, which is well correlated with significant overestimation of light precipitation.For both Germany and Switzerland a systematic underestimation of the diurnal IWV cycle between 6 and 21 UTC in both the model analysis and forecast is reported in Guerova and Tomassini (September 2003).::: For ::::::: Poland, A recent development in GNSS processing is use of the Precise Point Positioning (PPP) strategy (Zumberge et al., 1997).
The WRF model output is integrated into the SUADA.by integration over the model levels the WRF-IWV is obtained as below: where ρ w is density of liquid water, n is the number of model levels.

GNSS processing strategy and tropospheric products
Archived in SUADA are GNSS tropospheric products like Zenith Total Delay (ZTD over 12 000 000 individual observations) and derivatives like IWV (over 55 000) from five GNSS processing strategies and 37 stations in Bulgaria/Southeast Europe for the period 1997-2013.The temporal resolution of the GNSS data is from 5 minutes to 6 hours.
In this work we use GNSS tropospheric products from the BULgarian intelegent the GMF (Global Mapping Function) (Boehm et al., 2006) and 10 • elevation cut-off angle.The data were processed using the PPP strategy employing IGS satellite orbits and clocks.The computed ZTDs are with a temporal resolution of 300 s (5 min).
The ZTD data is archived in the GNSS SUBSCRIPTNBIN :: IN : SUADA table (Fig. 1).To derive GNSS-IWV the WRF model surface pressure (p s , [hPa]) and temperature (t s , [K]) are used in Eq. 2 (Davis et al., 1985) and Eq.3,4 (Bevis et al., 1992) (2) where ] are constants derived first by Thayer (1974) and is the height and θ is the latitude variation of the gravitational acceleration.
The pressure at the GNSS station altitude is calculated using the model pressure at the nearest model grid point.The pressure difference between the GNSS station altitude and the nearest NWP model grid point is calculated using the polytropic barometric formula Sissenwine et al. (1962): where P g is the pressure at the GNSS station altitude, P m is the pressure at meteorological station altitude, T [K] is the temperature in meteorological station, L = 6.5 K/km is tropospheric lapse rate, H m [km] is the altitude of the meteorological station, H g [km] is the altitude of the GNSS station, g 0 = 9.806 m s 2 is the gravitational acceleration, M 0 = 28.9644g mol is the molar mass of air and R = 8.31432 N m (molK) is the universal gas constant.

Surface observations
Archived in SUADA are also surface observations of: 1) pressure, 2) 2 m temperature and 3) precipitation (PP).The measurements are from the surface observation network (SYNOP) of the National Institute of Meteorology and Hydrology (NIMH) in Bulgaria and are collected manually every 3 hours (00, 03, 06, 09, 12, 15, 18 and 21 UTC).The data is available from the OGIMET weather information server (ogi, 2016).The surface pressure and temperature are used for derivation of IWV from the GNSS tropospheric products as described in Sect.2.2.Using surface observations the IWV is derived every 3 hours and is referred to as IWV* in Table 1.

Precipitation efficiency
In order to study water availability Tuller (1971) proposed ::: The : Precipitation Efficiency (PE) :: is : a :::::: value, :::::::: calculated ::: for ::::: each :::::: region :::::: climate :::: and :::::: climate ::::::::: variations.:: PE :: is : expressed as percentage of the IWV that is converted and measured as precipitation.Bordi et al. (2015) proposed to use GNSS IWV to compute PE.In this work the daily PE is computed as following: where PP and IWV are daily averaged precipitation and IWV at the station.Below :::::::::: Precipitation :::::::: efficiency ::::: gives : a ::::::::: long-term The annual mean, standard deviation and correlation between the surface pressure and temperature from WRF and SYNOP are presented in Table 1.The WRF pressure and temperature are extracted with the temporal resolution of the SYNOP i.e. every 3 hours.The correlation coefficient between the two data sets for atmospheric pressure for station Lovech (LOVE) is 0.99 with the mean difference 0.5 hPa.The correlation coefficient for the temperature is 0.96.The largest differences in the two data sets are observed for December 2013.For station Varna (VARN) the correlation coefficient for the pressure is 1 and for the temperature 0.96 with a mean difference between the SYNOP and NWP-WRF of 0.2 hPa for the pressure and 0.2 • C for the temperature.For station Burgas (BURG) the correlation coefficient for the pressure is 1 and for the temperature 0.96.The mean difference for the pressure is 0.1 hPa and 0.2 • C for the temperature.The NWP-WRF surface pressure shows an agreement of 0.5 hPa or better with the SYNOP data-set.This allows to take advantage of deriving IWV with the temporal resolution of the GNSS tropospheric products.A comparison between IWV and IWV* for station Burgas is seen in Fig. 3.

WRF-GNSS IWV: annual and monthly mean
A comparison between the GNSS and WRF IWV is presented in Table 3.For Burgas, Lovech, Montana, Shumen and Stara Zagora the correlation coefficient is very high and lies between 0.95 and 0.96.The mean IWV difference is between 0.5 and 1.8 mm.The smallest mean difference is obtained for Shumen and Burgas and is a consequence of the small altitude difference between GNSS station and WRF model height (Table 2).The altitude difference for station Lovech is 107 m and there the largest mean difference of 1.8 mm is obtained.For Varna and Rozhen the correlation coefficient is 0.9 and 0.76, and the mean IWV difference is negative with -0.9 and -3.2 mm respectively.

WRF-GNSS comparison: precipitation efficiency in 2013
Figure 6 showns the monthly mean precipitation efficiency computed with GNSS and WRF for Burgas and Lovech.The very good agreement between the model and GNSS is apparent at both locations.At Burgas (Fig. 6a) the PE has minimum in August less than 1 % and maximum in May 14 %.For Lovech (Fig. 6b) the maximum PE is in May but is slightly smaller than in Burgas (12 %) and the minimum is in September of about 1 %.The PE at the two stations also shows differences with are expected as the two stations are located in different climatic regions in Bulgaria.While Burgas is in south-east Bulgaria close to the Black Sea, Lovech is in north-west Bulgaria.The atmospheric circulation in Bulgaria is dependent on the Balkan mountains in middle of the country i.e. south of the range the Mediterranean cyclones are the main source of precipitation, while the north of the range their influence is largely reduced.The annual PE in both stations is in the range of 5.5-5.8 % from GNSS and 5.9-6.0 % from WRF, which is in agreement with the range of 5-10 % for the region found in Tuller (1971).

Conclusions
In this work GNSS tropospheric products (ZTD) :::: with temporal resolution 5 min .are derived using the PPP processing for one year period.In order to take advantage of the high temporal resolution of GNSS products for derivation of IWV the surface pressure and temperature from the NWP WRF model is used.The WRF surface pressure and temperature was evaluated against surface observations from three synoptic stations in Bulgaria.The mean difference for surface pressure between the two datasets is less than 0.5 hPa and the correlation is over 0.989.For the temperature the largest mean difference is 1.1 • C and the correlation coefficient is over 0.957.The IWV computed with this two data-sets has a mean difference is in range of 0.1-1.1 mm.

Figure 2 .
Figure 2. Map of the ground based stations of the Bulipos GNSS network.The red markers show the station locations.

:
Guerova et al. (2014)rs GNSS meteorology was developed in Bulgaria within a Marie Curie funded project.As apart : a :::: part : of this project a regional database for Bulgaria and Southeast Europe the Sofia University Atmospheric Data Archive (SUADA,Guerova et al. (2014)) was developed to facilitate the use of GNSS tropospheric products (see Sect. 2.2).In this study for the first time a PPP GNSS processing strategy is applied to seven stations in Bulgaria and GNSS IWV is derived using the Weather Forecasting System ::::::: Research :::: and ::::::::: Forecasting : (WRF) NWP model.