Assimilation of zenith total delays in the AROME France convective scale model: a recent assessment

The impact of assimilating GPS zenith total delays (ZTD) in the convective scale model AROME is assessed over a 1-month period in summer 2013. The experimental set-up is similar to the current operational usage at MétéoFrance where the observing system has been expanded in July 2013 in a three-dimensional variational (3D-Var) data assimilation scheme with a 3-hour cycling. Three experiments are performed. In a baseline experiment the GPS ZTD provided through the E-GVAP programme are withdrawn from the observing system (NOGPS). In a second experiment, GPS ZTD from E-GVAP are included in the observing system, representing the operational configuration at Météo-France (EGVAP). The last experiment is similar to EGVAP but new ZTD observations processed by the University of Luxembourg are also assimilated on top of all other observations (UL01). In the first stage, it has been verified through a systematic comparison with model counterparts that the quality of ZTD data processed by the University of Luxembourg is similar to the one provided by other analysis centres from the E-GVAP programme. After a number of quality controls, it has been possible to assimilate around 90 additional observations on top of around 600 stations from E-GVAP every 3 hours. Despite the small fraction of observations assimilated in AROME that ZTD represent (B2%), it is shown that they systematically improve the atmospheric humidity short-range forecasts by a comparison with other observing systems informative about water vapour (radiosoundings, satellite radiances, surface networks) even though it is by small amounts. When examining objective precipitation scores over France, the improvement brought by the UL01 stations on top of E-GVAP is systematic for all daily precipitation thresholds. Examination of several case studies reveals the ability of the ZTD observations to modify the intensity and location of precipitating areas in accordance with previous studies. The addition of ZTD from UL01 is also found to be beneficial, by improving rainfall patterns. Planned improvements to the AROME forecasting and assimilation systems with higher horizontal resolution and hourly cycling of 3D-Var assimilation will be of benefit to ZTD observations.


Introduction
Data assimilation systems aim at providing accurate initial conditions for numerical weather prediction (NWP) models. During the last decade, measurements from polar orbiting satellites have become the dominant observing system of data assimilation for global NWP models, particularly with the emergence of hyperspectral infrared sounders having thousands of channels such as the Interferometer Atmospheric Sounding Instrument (IASI) on board METOP-A and B (Hilton et al., 2012). In the meantime, there has been significant progress in the design of data assimilation systems for limited area mesoscale NWP models (Park and Zupanski, 2003;Sun, 2005). An accurate prediction of mesoscale phenomena such as severe convective storms, squall lines, wind gusts, sea-breezes or fog events requires information on the atmosphere at small temporal and spatial scales, particularly regarding the moisture field, that polar orbiting satellites can hardly provide. Therefore, in addition to conventional (radiosoundings, surface stations, aircraft data) and satellite observations, dedicated observing systems have been explored for operational purposes. One can cite reflectivity and Doppler winds from ground-based weather radar networks and Zenith Total Delays (ZTD) from groundbased GPS receivers (De Pondeca and Zou, 2001;Cucurull et al., 2004;Vedel and Huang, 2004;De Haan, 2013).
At Me´te´o-France, a three-dimensional variational (3D-Var) data assimilation system with a 6-hourly cycle was developed by Fischer et al. (2005) for the initialisation of the limited area NWP model ALADIN (Bubnova et al., 1995) with a 9.5 km horizontal resolution. A convection-permitting version of ALADIN was developed under the AROME project, leading to the AROME model (Seity et al., 2011) with a 2.5 km horizontal resolution that became operational in December 2008. The ALADIN 3D-Var assimilation system was adapted to the smaller scales of AROME with a 3-hour Rapid Update Cycle (RUC) . The first version of the 3D-Var AROME took advantage of a number of feasibility studies regarding the usage of dedicated observations suitable to the mesoscale with high spatial resolution and frequent temporal availabilities. This has been the case for SEVIRI infrared radiances for the geostationary satellite METEOSAT (Montmerle et al., 2007), radar Doppler winds from the French network ARAMIS (Montmerle and Faccani, 2009) and ZTD from ground-based GPS stations (Yan et al., 2009b). Since then, the number of data assimilated has significantly increased due to two main factors: the use of satellite at higher spatial density (the spatial thinning of the data has been reduced from 125 to 80 km) and the use of radar reflectivities through a Bayesian inversion (Caumont et al., 2010;Wattrelot et al., 2014). The number of ZTD observations has also been increased with revised quality controls.
Many studies have been undertaken to demonstrate the interest of GPS ZTD observations for the meteorological community, starting from comparisons with Integrated Water Vapour (IWV) from radiosoundings to impact and validation studies in NWP models (Yang et al., 1999;Cucurull et al., 2000;Vedel and Huang, 2004;Macpherson et al., 2008;Boniface et al., 2009;Bennitt and Jupp, 2012). They have shown the important value of such data available at high temporal frequency providing accurate information on precipitable water and surface pressure in all weather conditions, thanks to the use of the L-band microwave frequency. In Europe, the joint effort of the meteorological and geodetic communities has led to the EUMETNET programme E-GVAP, 1 providing in near real-time (timeliness around 90 minutes) about 1000 ZTD measurements every 15 minutes over a large part of the European continent. However, when feasibility studies were undertaken at Me´te´o-France both in the global model ARPEGE (Poli et al., 2007) and in the limited area models ALADIN (Yan et al., 2009a) and AROME (Yan et al., 2009b), the number of observation types was much lower (including ZTD data).
The goal of the present paper is two-fold. First we want to provide a recent assessment of the impact of GPS ZTD data in the current version of the operational 3D-Var data assimilation system for AROME. Indeed, the last previous similar study was undertaken by Yan et al. (2009b) in 2008 with a pre-operational version of the AROME 3D-Var and, as stated above, many components of the observing system have changed since then. Moreover, AROME has been improved regarding its vertical resolution (60 levels instead of 43) and the horizontal domain has been extended. Such an impact is examined by withdrawing GPS ZTD observations from the current observing system of the AROME 3D-Var over a 1-month period in summer 2013. Secondly, we want to examine the impact of increasing the GPS network through the addition of a new processing centre from the University of Luxembourg (UL01). Indeed, the current GPS network available through E-GVAP is expanding and the interest of new partners in the E-GVAP programme providing additional surface stations and alternative ZTD solutions can be assessed. Moreover, the systematic comparison of ZTD solutions from an analysis centre with model counterparts provides a way to quantify the quality of ZTD values that could lead to improved data processing.
The paper is organised as follows. Section 2 summarises the main features of the numerical model AROME. In Section 3, the 3D-Var assimilation system is described with a special emphasis on the pre-processing of GPS ZTD data. This section also provides information on the additional dataset processed by the University of Luxembourg. The experimental set-up of the data assimilation experiments is explained in Section 4 together with the results in terms of analysis and forecast impacts. Finally in Section 5, the main conclusions of the study are drawn together with a number of perspectives for an improved usage of GPS ZTD observations for mesoscale NWP models.

Model description
The AROME spectral mesoscale NWP model is derived from the global model ARPEGE (Courtier et al., 1991) and is operational at Me´te´o-France since December 2008. It has been adapted to the limited area geometry and the hydrostatic hypothesis has been relaxed in the adiabatic prognostic equations (Seity et al., 2011). The package of physical parameterisation schemes is derived from the Meso-NH research model (Lafore et al., 1998). Stratiform and deep convective clouds are explicitly represented by a set of six prognostic equations for the water species (vapour, cloud water, rain water, ice crystals, snow and graupel) with a bulk mixed-phase microphysical scheme (Pinty and Jabouille, 1998) that describes the conversions between species. For shallow convection, a sub-grid parametrisation scheme is considered using the Eddy Diffusity Mass Flux (EDMF) concept to describe the vertical transport by thermals in the boundary layer (Pergaud et al., 2009). The turbulence scheme is based on a prognostic Turbulent Kinetic Energy (TKE) equation combined on a diagnostic mixing length (Cuxart et al., 2000). The radiative transfer scheme is derived from the ECMWF model: the shortwave spectrum is described with six spectral bands (Fouquart and Bonnel, 1980) and the longwave spectrum is computed by the Rapid Radiative Transfer Model (RRTM) of Mlawer et al. (1997). The surface fluxes are described using the externalised module SURFEX (Masson et al., 2013) where the ISBA scheme (Noilhan and Mahfouf, 1996) is used for soil and vegetation transfers and the Town Energy Budget (TEB) scheme (Masson, 2000) is used over urban areas.
The AROME domain covers a significant part of Western Europe (Fig. 1). The horizontal resolution is 2.5 km on a Lambert projection with its centre at (46.48N, 2.28E) with 750 and 720 physical grid points in the eastÁwest and northÁ south directions respectively. The domain is vertically divided into 60 layers using a hybrid pressure terrain-following coordinate system. The centre of the uppermost layer is located at 1 hPa. The height of the lowest layer centre is about 10 m above the ground. Between the surface and this first model level, a set of six levels has been added to resolve explicitly the vertical turbulent transfers in the surface boundary layer with a TKE scheme for temperature, specific humidity, and wind components (Masson and Seity, 2009). The use of a semi-implicit semi-Lagrangian temporal integration scheme allows a model time-step of 60 seconds. Lateral boundary conditions are provided by hourly ARPEGE Fig. 1. Map of the domain of the AROME model with the GPS stations selected from assimilation on the 18 July 2013 at 03 UTC. The black squares correspond to the stations provided through the E-GVAP network (EGVAP). The grey circles correspond to the new stations processed by the University of Luxembourg (UL01). Note that except for stations located in Belgium, the UL01 stations are also processed by other analysis centres.
forecasts. This configuration of AROME was adopted for operational applications in April 2010 and is used for the experiments described hereafter.

General features
The AROME 3D-Var assimilation system (Seity et al., 2011) is derived from the regional ALADIN-France 3D-Var scheme (Fischer et al., 2005) inherited from the ARPEGE-IFS software (Courtier et al., 1991). The atmospheric data assimilation system is configured with a 3-hour cycle. The two components of the wind, temperature, specific humidity and surface pressure are analysed at the model resolution (2.5 km) with an incremental version of the 3D-Var system. The other prognostic variables are left to adjust themselves during the first steps of the subsequent forecasts.
The observations used stem from radiosondes, wind profilers, ship and buoys reports, aircraft reports and automatic land stations for conventional data. Regarding remotesensing observations, the infra-red and microwave radiances from the following instruments available on various satellites are assimilated: HIRS, AMSU-A, AMSU-B/MHS, SSMI/S, ATMS, AIRS, IASI, CrIS. A spatial thinning of 125 km is used for all radiances except for AMSU-A and IASI where an 80 km sampling is adopted. Radiances from these instruments are assimilated over all surfaces (but with different channel selections) except for the radiances from the microwave imager SSMI/S that are only considered over oceans. Horizontal winds derived from the imagery of geostationary and polar orbiting satellites and from the ASCAT and OSCAT scatterometers are also included. The infra-red radiances from SEVIRI on board METEOSAT-10 are assimilated (Montmerle et al., 2007) at 50 km resolution with a recent increased usage of channels sensitive to the surface over land (Guedj et al., 2011). AROME assimilates radial velocities from Doppler radars of the French network (Montmerle and Faccani, 2009). Radar reflectivities are also assimilated with a two-step approach where reflectivities are first converted into relative humidity profiles with a 1D Bayesian inversion that are then considered as pseudoobservations in the 3D-Var (Wattrelot et al., 2014). The spatial thinning of radar data is 15 km to avoid spatial correlations of observation errors. Finally the 3D-Var assimilates GPS ZTD observations from the E-GVAP network, to be described more extensively in the next section. The background error statistics of the 3D-Var with a multivariate formulation have been derived from an ensemble data assimilation technique described in Brousseau et al. (2011). AROME also includes a dedicated surface assimilation scheme based on local optimum interpolation for the prognostic variables of the land surface schemes (Giard and Bazile, 2000;Mahfouf et al., 2009). This soil analysis scheme has been recently extended from the ISBA scheme to the other surface schemes available within the externalised modelling platform SURFEX (Masson et al., 2013;Hamdi et al., 2014). As input, this assimilation scheme uses increments of an analysis of screen-level temperature and relative humidity based on a bi-dimensional optimum interpolation called CANARI, considering SYNOP observations from the surface weather network. The CANARI system is also used to perform a sea surface temperature analysis from buoy and ship reports. Figure 2 displays the number of observations assimilated in the AROME operational system that went operational on the 2nd of July 2013 for a particular day (26 January 2014) with significant rainfall events over the domain. Four main observation types dominate: radiances from the two IASI instruments (with more than 100 channels each) (27%), radial winds and pseudo-humidity profiles from radars (22%), data from surface weather stations (16%), and aircraft reports (15%). Then, satellite data (dominated by SEVIRI and ATOVS radiances) and radiosoundings represent about 17% of the observations, whereas GPS ZTD are only 1.6% of the total data. This percentage reaches 6% when considering only observations sensitive to atmospheric moisture. This small amount is partly compensated by the fact that these measurements are available every 3 hours in the AROME 3D-Var RUC in all weather conditions with a relatively good accuracy. In terms of information content on the analysis that can be measured by the Degree of Freedom for Signal (DFS) as proposed by Chapnik et al. (2005), the GPS ZTD data have a contribution between 3 and 4%. When restricting this contribution to observations sensitive to atmospheric humidity, this percentage increases to 12%. From these numbers, even though the impact of GPS ZTD measurements on AROME analyses and forecasts is expected to be positive as demonstrated in previous sensitivity studies, it cannot be huge, given the small fraction they represent in the AROME 3D-Var system, both in terms of number and DFS value. However, we expect them to play an instrumental role in specific severe weather situations. Note that a DFS value provides information on the capacity of an observing system to modify the a priori information provided by the background (short-range forecast) but it does not indicate if these changes are beneficial or detrimental.

GPS observations and pre-processing
In this study, two GPS ZTD datasets are considered. The first one is provided through the EUMETNET E-GVAP programme. Currently a set of 16 analysis centres is contributing to E-GVAP. Each centre provides data at a temporal frequency between 5 and 60 minutes (even though the observation frequency is 15 minutes for most of the centres).
Over the AROME domain, about 700 GPS stations are selected for assimilation every 3 hours, as shown in Fig. 1 for a particular assimilation cycle. Note that measurements from a given station can be processed by several analysis centres at the same time which leads to more than one solution (i.e. ZTD value) per station. The second dataset comes from the University of Luxembourg (UL01) and is currently being monitored by E-GVAP but not distributed yet. It provides about 220 solutions every 15 minutes, with 28 new stations located in Belgium and Luxembourg and about 200 existing stations already processed by other analysis centres. Only a sub-set of these stations is displayed in Fig. 1. They correspond to the ones that have been retained for assimilation after a pre-processing step described next.
The pre-processing step allows selecting the stations and processing centres of good quality that can enter a 'white list' in order to be potentially assimilated. Any station-centre pair that is not present in the 'white list' cannot be assimilated. On the other hand, before any station-centre pair can be assimilated, several quality controls are performed. The initial 'white lists' proposed by Poli et al. (2007) and Yan et al. (2009aYan et al. ( , 2009b were set-up when only few stations having data of good quality were available over Europe. The 'white list' has been significantly extended in 2010 since the network has expanded and the quality of solutions from analysis centres has also considerably improved. Currently the operational 'white list' for AROME contains 806 station-centre pairs from E-GVAP. This list has been extended by 155 stations processed by the University of Luxembourg (UL01). These numbers have been obtained through a careful examination of the statistical behaviour of the differences between ZTD observations and the model equivalent, named hereafter background departures, computed from 3-hour AROME forecast fields of surface pressure, temperature and specific humidity using the following observation operator formulation: where p is the pressure, T the temperature, e the water vapour pressure, k 1 00.776 Pa (1 K, and k 3 03730 Pa (1 K 2 . This simple formula has been compared by Yan et al. (2008) with a more recent one proposed by Brenot et al. (2006) without significant differences in terms of relative humidity analysis increments. When the surface pressure P stat at a GPS station is between two pressure model levels P k and P k'1 , the contribution of the last layer in eq. (1) is reduced by a factor: (P k Á P stat )/(P k Á P k'1 ). Similarly, when the surface pressure at GPS station is higher than the model surface pressure P surf , the contribution of the last layer in eq. (1) is increased by a factor: (P 1 Á P stat )/(P 1 Á P surf ) where P 1 is the atmospheric pressure of the first model level above ground. The underlying assumption is that atmospheric properties are kept the same between P surf and P stat . Complementary details can be found in Poli et al. (2007). The contribution to the ZTD above model top z top at pressure P top is currently accounted for in the bias correction. In the future, we will use the following correction: where R d is the gas constant for dry air and g top the acceleration of gravity at z top . It has been proposed in Poli et al. (2007) and Yan et al. (2009a) to account explicitly for the missing part of the upper atmosphere. First a spatial thinning is done to ensure that the distance between two stations is no less than 10 km. The observation error is set to 12 mm for all E-GVAP stations. For UL01, this error has been assigned for each station based on the standard deviation of background departures and assuming a model background error of 6 mm, leading to actual values for individual stations ranging between 10 and 15 mm. The mean value of background departure statistics is used to provide a constant bias correction (accounting mostly for the difference in altitude between the model and the station and for the atmospheric contribution to ZTD above z top ) that is applied to each station-centre pair before assimilation. The time series of background departures have been computed over the 1-month assimilation period performed in this study for the UL01 ZTD observations whereas the statistics for E-GVAP ZTD observations correspond to the period of April 2013. We have checked for E-GVAP stations that the biases computed over the two periods are very similar (not shown). Station-centre pairs with a toolarge bias (larger than 30 mm) are not selected in the 'white list'. During the assimilation, for a given station, the processing centre is chosen on the basis of its availability and by comparing the background departure statistics, instead of being selected a priori in the 'white list'. For each 3-hourly analysis, when several observations are available from a given station-centre pair, the GPS observation closest to the analysis time and within '/ ( 1.5 hours around the analysis time is selected. Like other observations assimilated in the AROME 3D-Var, ZTD data are subject to a First-Guess Quality Control check that rejects data too far from the model background. This check is based on background and observation errors with a threshold value of 3.5. In practice, it means that ZTD values leading to background departures (OmB) larger than 40 mm are systematically rejected, accounting for too-large mismatches between model and station altitudes. From the initial list of GPS stations from UL01 (about 200), only 155 have been retained for assimilation based on the above criteria. In particular the set of six stations in Luxembourg has been excluded since they presented too large biases (around (30 mm). The reason for this behaviour is unclear and is the topic of on-going investigations at the processing centre. On the other hand, a set of 22 stations in Belgium has been retained since the mean bias was only 0.63 mm.

The experimental design
Three AROME assimilation experiments have been set up in order to investigate the impact of the assimilation of GPS ZTD observations with the most recent version of the 3D-Var system and also to evaluate the impact of additional ZTD data processed by the University of Luxembourg (Ahmed et al., 2014). A first experiment (NOGPS) assimilates all the observations from the operational AROME 3D-Var since 2 July 2013, but excluding all GPS ZTD observations. The second assimilation (EGVAP) is similar to NOGPS but with EGVAP ZTD observations (this is a research experiment identical to the operational configuration). The last assimilation experiment (UL01) is similar to EGVAP but with ZTD solutions provided by the University of Luxembourg on top of operational E-GVAP solutions.
The three assimilation experiments started at 0300 UTC 18 July 2013 from the AROME operational analysis and ended at 2100 UTC 20 August 2013. Short-range forecasts are run every 3 hours to provide the background of the next analysis and once a day at 0300 UTC, a 30-hour forecast run is launched. This experimental set-up is summarised in Table 1.

Impact on analyses.
A first assessment of the assimilation of GPS ZTD in the AROME 3D-Var is done by examining the ZTD distributions of differences between observations and their model equivalent, both for the background (3-hour forecast) (OmB: Observation minus Background) and for the analysis (OmA: Observation minus Analysis) over the period of interest corresponding to 240 analyses. Figure 3 displays the distributions of (OmB) and (OmA) from experiments EGVAP and UL01. Only stations processed by the University of Luxembourg are selected Table 1. Summary of the 3D-Var assimilation experiments in the AROME model undertaken with different observing systems. E-GVAP corresponds to ZTD data provided operationally through the EUMETNET programme E-GVAP. UL01 corresponds to a data set processed by the University of Luxembourg for the current study.

Experiment name
Observations assimilated in 3D-Var NOGPS All observations minus E-GVAP ZTD EGVAP All observations (operational configuration) UL01 All observations plus UL01 ZTD from UL01 experiment. The size of the sample is 147591 for EGVAP and 20572 for UL01. For EGVAP, the shapes of the distributions are close to Gaussian curves with small values for the means ((0.34 mm before and (0.09 mm after analysis) indicating that the various quality controls and bias corrections that have been done beforehand have allowed the data to satisfy the necessary criteria of the 3D-Var system to use the observations optimally (unbiased observations with Gaussian error statistics). The narrower distribution of (OmA) with respect to (OmB) (reduction of the standard deviation from 11.49 to 2.33 mm) indicates that the assimilation has brought the model state, mostly in terms of humidity field, closer to GPS ZTD observations.
A posteriori diagnostics on innovations and analysis residuals  provide values for observation and background errors that are much smaller than the prescribed ones: 4.4 and 3.7 mm, respectively. This means that the 3D-Var system could use, in future experiments, smaller errors and in particular a ratio of actual observation to actual background errors closer to 1. From experiment UL01 and considering only the GPS ZTD values processed by the University of Luxembourg, the histograms have also a Gaussian shape for both (OmB) and (OmA), but they are less smooth due to the smaller size of the data sample (20572 vs. 147591). The mean and standard deviation of the (OmB) and (OmA) distributions compare well to those obtained with the whole set of selected EGVAP stations. This agreement reveals that the selected stations processed by UL01 have the same quality as the EGVAP stations and that the 3D-Var AROME system can assimilate them efficiently. This is also a very good check that the bias correction and the observation error specifications for UL01 are fully consistent with values imposed for the EGVAP stations. Another useful aspect of the (OmB) and (OmA) statistics to examine is their temporal evolution during the period of interest. Times series of standard deviations for (OmB) and (OmA) are compared for the two datasets EGVAP and UL01 in Fig. 4. There is a good level of agreement between these two datasets both in terms of daily and synoptic variations. Despite a much smaller number of stations coming from UL01, their geographical distribution covers a large part of the AROME domain, allowing them to capture almost the same variability as the full EGVAP network. A diurnal cycle of the departures is also present on both time series. This is a feature present for all analysis centres: maximum values take place at 21 UTC whereas minimum values are noticed at 06 UTC, and is likely a signature of the diurnal cycle of water vapour in the boundary layer. A dependency of the bias correction with the diurnal cycle could be envisaged in future developments. Synoptic events show up on time series: the anticyclonic situations are characterised by lower values of (OmB) and (OmA) (around day 14 and day 24) whereas perturbed situations exhibit higher values (around days 6, 10 and 16). The previous analysis of (OmA) and (OmB) distributions and time series allows to check that the 3D-Var AROME behaves as expected regarding the use of GPS ZTD observations. This is a necessary condition, but not sufficient to demonstrate that the analysis is improved by the GPS ZTD observations. Since they modify mostly the humidity field, by examining the fit of the model short-range forecasts (3-hour) to independent observations (i.e. that have not been assimilated yet), it is possible to check if the analysis from the previous cycle has been able to bring the model predicted state closer to the truth. Figures 5, 6 and 7 show the fit of the AROME 3-hour forecasts to other observations sensitive to humidity (radiosoundings and moisture sensitive channels from satellite radiometers). The assimilation of GPS ZTD reduces the standard deviation of (OmB) for specific humidity from radiosoundings below 500 hPa (Fig. 5). The UL01 observations enhance this reduction at 1000 hPa but are slightly detrimental at 700 hPa with respect to EGVAP. Regarding the window channels of the SSMI/S instrument (Fig. 6), non-negligible improvements are noticed by the use of EGVAP for the weak water vapour absorption band at 22 GHz and also for the 85 GHz channel in vertical polarisation (that has the largest errors). The influence of UL01 is either neutral or slightly negative with respect to EGVAP, but always positive with respect to NOGPS. Since SSMI/S is only available over ocean surfaces, it indicates that the corrections brought by the GPS ZTD data over land have propagated over the ocean through error correlations and model dynamics (AROME has 40% of ocean surfaces). The (OmB) of the sounding channels from the ATMS microwave instrument in the water vapour absorption band at 183 GHz are displayed in Fig. 7. Surprisingly, the largest impact shows up for the high peaking channels (that are closer to the centre of the absorption line at 183 GHz). It means that ZTD observations have the potential to also modify the mid-tropospheric humidity (between 500 and 300 hPa), probably through the vertical correlations imposed in the background error covariance matrix of the 3D-Var system. This signal is more pronounced when the model state is projected onto the brightness temperature space of ATMS than on radiosonde data. The impact of UL01 is either neutral or slightly detrimental with respect to EGVAP, but still positive when compared to NOGPS. However, the signal coming from radiosoundings is certainly more robust given the fact that in the midtroposphere the number of data from that observing system is about three times larger than ATMS radiances (see Fig. 2). Regarding the (OmB) biases for specific humidity from radiosoundings, negative values noticed in the boundary layer (i.e. the model is too moist) in experiment NOGPS (around 0.2 g/kg) are slightly reduced by about 10% through the assimilation of GPS ZTD observations (not shown). Such reduction is consistent with previous findings by Yan et al. (2009a).
In order to examine the impact of the background error covariance matrix on the assimilation of ZTD observations for a given analysis, we compare the analyses produced after the first cycle of the three 3D-Var experiments. Indeed, since each experiment starts from the same background field, the differences in analysis indicate the contribution of the GPS ZTD from EGVAP to modifications of the humidity field when compared to the experiment NOGPS, and also the contribution of UL01 on top of EGVAP (through differences between UL01 and EGVAP). The differences in IWV between EGVAP and NOGPS and UL01 and EGVAP are shown in Fig. 8, and can be interpreted as analysis increments produced in the 3D-Var by the EGVAP and UL01 ZTD networks, respectively. The geographical distribution of the increments follows the density of each network (see Fig. 1). Isolated stations, such as the one in south Sardinia, provide information about the correlation errors chosen in the background error statistics of the 3D-Var. The statistical model has homogeneous and isotropic correlations with horizontal length scales about 50 km. This particular coastal station shows how terrestrial ZTD stations can modify humidity over oceans. The standard deviation of IWV increments for EGVAP is 0.87 mm that can be compared to the mean IWV value of 23 mm over the whole domain. The corresponding value for UL01 is 0.38 mm.    Therefore the mean corrections are small. Even though mean values are very close to zero for EGVAP and UL01, the distribution of increments for UL01 is skewed towards positive corrections. This is specific to this date: through the whole period there is no tendency for the UL01 dataset to systematically moisten or dry the model. Maximum corrections are about 4Á5 mm, whereas with EGVAP over the western part of Switzerland there is a local drying of 7 mm. With UL01, large corrections are taking place over the Wallonie region thanks to the density of the network. As seen over this particular area, when the network is dense enough, the shape of the increments is less dominated by the hypotheses of the error correlation model. In the northeastern part of France (Alsace and Lorraine regions) UL01 brings corrections around 3 mm whereas almost no corrections were provided by EGVAP. This stems from the fact that the UL01 solutions have been chosen instead of the SGN (name of analysis centres from the French National Geographical Institute) ones, on the basis of a better statistical behaviour.
In summary, the GPS ZTD have a small but positive impact on the humidity analyses and short-range forecasts of AROME. The behaviour of the UL01 stations is very consistent with the EGVAP stations giving confidence on the quality of the data selected for assimilation.

Impact on forecasts
. The impact on screen-level relative humidity RH 2m forecasts is presented in Fig. 9. The AROME forecasts are examined in terms of mean and standard deviation errors with respect to a reference assumed to be close to the truth. The reference is given by screen-level analyses based on the optimum interpolation CANARI presented in Section 3.1. Screen-level analyses are done independently from the atmospheric 3D-Var analyses and are used for diagnostic purposes and to correct soil temperatures and moisture contents. Around 1600 observations are used to perform the screen-level analyses, with rather dense surface networks except over Italy and Spain. The RH 2m mean errors have a strong diurnal cycle with values close to zero during daytime and negative ones during night-time. The nocturnal surface boundary layer is too moist but it is also too warm as revealed by T 2m biases (not shown), contributing to reduce the RH 2m bias. The standard deviation increases with the forecast range up to 15 hours (18 UTC) and then remains around 9% error until the end of the forecast (after 30 hours). The assimilation of GPS ZTD has a small positive impact on the bias by reducing it after 15 hours of forecasts. The impact on the standard deviation error is also rather small, with nonnegligible values only at 15 UTC (12-hour forecast).
The spatial distribution of these errors is examined more precisely at 15 UTC where the largest differences are noticed between the baseline NOGPS and the experiments. This can be more easily performed on gridded analyses than with actual observations at irregular locations. The mean and standard deviation 12-hour forecast errors for RH 2m at 15 UTC are displayed in Fig. 10. We first examine the reference simulation NOGPS. The most striking pattern of the mean errors is a positive bias (model too moist) over areas with significant orography (mostly the Alps but also the Pyrene´es and the Massif Central). Values are above 10% over the western part of Switzerland. The eastern part of France and Germany are characterised by moist positive biases around 5%. On the other hand, Spain and Italy exhibit a slight negative bias (model too dry) whereas over Sardinia values are around (10%. Regarding the standard deviation errors, there are only few areas with values below 5% (Southern Spain), most of the errors being between 10 and 15%. Regions with the largest errors around 20% appear over Western Germany, eastern France and around the Pyrene´es. Over oceans the biases and the standard deviations are small since the surface boundary layer is constrained by a saturated surface for which the temperature is imposed from a dedicated surface analysis.
The assimilation of GPS ZTD data has a tendency to reduce the positive biases noticed over Germany and eastern France (for this last region it has been replaced by a small negative bias). The maxima over western Switzerland is also decreased. No noticeable differences exist between EGVAP and UL01 experiments. The areas of lower biases are also regions where the standard deviation errors are significantly decreased. The band of large values (above 20%) oriented SouthÁWest/NorthÁEast from the Massif Central to Western Germany in the NOGPS experiment is not present in the experiments when ZTD data are assimilated (both EGVAP and UL01). A number of rainfall events have passed through the domain during the chosen period, most of them having an orientation SouthÁWest/NorthÁEast (discussed later), that can explain why significant moisture errors can be found in the boundary layer at these locations.
The precipitation forecasts are compared over France against precipitation analyses developed at Me´te´o-France at 0.0128 resolution at hourly time scales from a blend between radar derived rain rates for small scale precipitation (convective events) and rain gauges that describe accurately the large-scale precipitation events and that can be spatially interpolated using a classical kriging technique. This analysis called ANTILOPE can be considered as the closest representation of the truth over France. In the following, only 24-hour accumulations from 06 UTC on day 1 to 06 UTC on day 2 from AROME are evaluated against ANTILOPE. First an objective comparison with classical categorical scores derived from contingency tables is performed and summarised in Fig. 11 for the Frequency Bias Index (FBI), the Probability of Detection (POD), the False Alarm Ratio (FAR) and the Equitable Threat Score (ETS) for the following thresholds: 0.2, 1, 2, 5, 10 and 20 mm. The actual definition of these scores can be found at: http://www.cawcr.gov.au/projects/ verification. The NOGPS experiment exhibits a slight positive bias up to 15% for large precipitation amounts (above 5 mm). Both ZTD assimilation experiments reduce this bias. The bias has the smallest value with UL01 for all rates. When examining the POD, the three experiments produce very similar values; however, EGVAP appears to slightly degrade the POD above 5 mm, whereas the quality of UL01 is always as good as or better than NOGPS. The assimilation of ZTD reduces the FAR for all thresholds, but UL01 is better above 1 mm. Finally the ETS, that accounts within a single measure for the POD and FAR skills, reveals higher values with the assimilation of ZTD, but the UL01 is clearly better since the improvement with EGVAP is only noticeable below 5 mm. The small detrimental effect of ZTD assimilation for the lowest threshold (0.2 mm) on the POD is also present on the ETS. These results are consistent with previous feasibility studies and  Mean and standard deviation of screen-level relative humidity 12-hour forecast errors (15 UTC) in percentages for experiments NOGPS (top panels), EGVAP (middle panels) and UL01 (lower panels) averaged from 20 July to 20 August 2013. also reveal that increasing the density of the GPS network even by a small percentage can be beneficial to the quality of precipitation forecasts from NWP models, despite an almost neutral impact on other predicted quantities.
We examine a number of case studies that took place during the period of interest in order to highlight the sensitivity of the assimilation of GPS ZTD data on the prediction of severe convective events.
On 25 July 2013, two intense rainy bands are present in the ANTILOPE analysis that are embedded in a frontal system crossing France in a synoptic southwesterly flow (Fig. 12). A first rain band is located south of Paris crossing the 488 latitude and located to the East of the Greenwich meridian, and a second one is located in southwestern part of France (Bordeaux region) with also a small band of precipitation above (Charentes region). Precipitation maxima are above 50 mm. In the NOGPS experiment, AROME has the capacity to simulate these two rainy systems, but the northern one is located too far west from its actual position with an excessive southerly orientation. The intensity of the system is underestimated. The system over the Bordeaux region is better described, but in terms of orientation and precipitation maxima, it extends too far east. The EGVAP experiment intensifies the northern precipitation band but the wrong orientation and the extension too far west are amplified. The intensity of the southern precipitation band is globally reduced. This can be seen as an improvement for the spurious cell located on the east, but as a degradation for the intense cell located near the coast of the Atlantic Ocean. In the UL01 experiment, the northern system is displaced northwards and eastwards, therefore improving significantly its location with respect to the truth, while keeping the larger amounts already produced by the EGVAP experiment. Similarly, Fig. 11. Categorical scores (FBI: Frequency Bias Index; POD: Probability of Detection; FAR: False Alarm Ratio; ETS: Equitable Threat Score) for daily accumulated precipitation AROME forecasts from the three experiments: NOGPS, EGVAP, UL01 compared with ANTILOPE precipitation analyses over France (period of interest: 18 JulyÁ21 August 2013).
the southern system is closer to the observed location, being less elongated to the east.
On 26 July 2013, the frontal system described above has only moved slightly northwards (Fig. 13). The northern precipitating area is located above Paris (Picardie and Nord regions) with three narrow bands: two with a WestÁEast orientation and another one above them with a SouthÁ West/NorthÁEast orientation. The southern precipitating area has a NorthÁSouth orientation. In the NOGPS experiment, there is a hint of the three bands observed in the northern part of France, but the third one to the East is too wide and too intense particularly over Luxembourg. The southern precipitating system has a correct orientation but is situated too far west over the Atlantic Ocean. The EGVAP experiment reduces the intensity of the rainy cell over Luxembourg but is displaced to the South where it is not observed. For the southern precipitating system, there is a reduction of the rainy band over the ocean and an intensification of the band over the continent parallel to this one. When considering UL01, in the northern part of the domain, the too-wide precipitating area southwest of the Luxembourg is not simulated anymore, and it has been replaced by narrower rainy bands that, even though not exactly at the proper location, resemble the actual precipitation patterns in this area. In the southwestern part of France, the northern branch of the precipitating system is improved but not the southern one.
On 02 August 2013, a rainy system moved over France in a southwesterly large-scale flow leading to precipitation around 15 mm over the Lorraine region (northeastern part of France) and in the Centre region (Fig. 14). A band of more intense precipitation with values above 30 mm and a Daily accumulated precipitation in mm as analysed by ANTILOPE and simulated by the AROME model starting from different atmospheric analyses on 25 July 2013: analysis (upper left), experiment NOGPS (upper right), experiment EGVAP (lower left), experiment UL01 (lower right). more westerly orientation is noticed over the Aquitaine region. Several narrow precipitation bands are captured by the ANTILOPE analysis. In experiment NOGPS, even though the model simulates a rainy system with the correct orientation and maxima located in the southwestern part of the domain, their intensity is underestimated by a factor of 2 (around 25 mm) and the precipitation band does not extend enough westwards near the Atlantic Ocean coast. In experiment EGVAP, the area of precipitation maxima located over the Aquitaine region is strongly enhanced with maximum values close to the observed ones (40 mm). The assimilation of ZTD data from UL01 leads to an even wider precipitating area in the southwestern part of France. However, none of the simulations can describe the smallscale patterns of the precipitation displayed in the analysis.
In experiment UL01, the underestimation of observed precipitation located east of Paris, already present in experiment EGVAP is somewhat amplified with respect to experiment NOGPS.

Conclusions and perspectives
In this paper, we have performed an impact assessment of the assimilation of GPS ZTD observations within the 3D-Var system of the convective scale model AROME developed at Me´te´o-France for issuing operational mesoscale short-range forecasts over Western Europe. Experiments have been undertaken during a 1-month period in JulyÁAugust 2013, corresponding to a recent upgrade of the operational forecasting suites with additional observing systems at Me´te´o-France. Even though many impact studies have been performed in the past to examine the interest of GPS ZTD on the skill of NWP models, an original aspect concerns the design of the observing system where GPS ZTD data representB2% of the total observation assimilated in AROME. This small fraction is a consequence of the capacity of the 3D-Var to ingest many other data types, in particular radar reflectivities and radial winds, that provide significant information on moisture and dynamical fields at mesoscale. However, contrary to other remote sensing observations that are only useful when it rains (weather radars) or in clear-sky or above clouds (satellite radiances), GPS ZTD data are available in all weather conditions at high temporal sampling. First it was important to verify that the GPS ZTD still contribute in a positive manner to the analyses and to the forecasts of AROME, in particular since a recent revision of the 'white list' and of the data selection. Then, we wanted to examine the sensitivity of the forecasting system to additional GPS ZTD observations provided by the University of Luxembourg (UL01). We have shown that the static bias correction scheme allows to remove most of the bias and that the data selection leads to unbiased and Gaussian statistics that are necessary conditions for an optimal assimilation. This finding is true for both E-GVAP and UL01 ZTD datasets. This is particularly encouraging for the use of UL01 within the E-GVAP network: after removing few stations having too large biases, the quality of the remaining stations was found similar the one of current operational stations provided through E-GVAP.TheAROME short-range forecasts have been found closer to other Daily accumulated precipitation in mm as analysed by ANTILOPE and simulated by the AROME model starting from different atmospheric analyses on 2 August 2013: analysis (upper left), experiment NOGPS (upper right), experiment EGVAP (lower left), experiment UL01 (lower right).
observations types sensitive to humidity when GPS ZTD are assimilated. This is a clear signal that through the assimilation of GPS ZTD the 3-hour prediction of the humidity field by AROME is closer to the truth. The signal is mostly seen on standard deviation errors but low level humidity biases are also slightly reduced. The new set of ZTD observations provided by the University of Luxembourg also satisfies these criteria and their assimilation has been successful in AROME. Regarding the impact on objective forecast skill scores, a small positive benefit has been noticed both on screen-level relative humidity and 24-hour precipitation accumulations. The categorical scores are systematically improved when the UL01 data are assimilated on top of E-GVAP ZTD observations. When examining case studies, it has been confirmed that GPS ZTD observations affect the predicted location and intensity of rainy systems that generally improves the quality of the numerical forecasts. For such specific situations the additional ZTD data processed by the University of Luxembourg significantly modify rainfall patterns with, most of the time, a better location and intensity of precipitating cells.
Despite the overall positive impacts summarised above, the use of GPS data in mesoscale data assimilation systems can still be improved. Several areas are currently explored both at Me´te´o-France and through European collaborations. The temporal availability of most ZTD observations every 15 minutes makes the assimilation of one observation every 3 hours clearly sub-optimal. At Me´te´o-France, GPS ZTD observations are also assimilated in the global model ARPEGE that has a four-dimensional variational (4D-Var) data assimilation system allowing to have several observations within a 6-hour window. A direct consequence is that the number of GPS ZTD observations is 10 times more important in ARPEGE than in AROME, even though the spatial coverage is only slightly increased (i.e. the EGVAP network is complete but no data outside Europe are currently exchanged on the Global Telecommunication System). AROME does not have a 4D-Var system, but a new version of 3D-Var RUC with hourly cycling is currently under evaluation that would increase the number of ZTD data assimilated by a factor of 3. With this new system, AROME has a finer grid mesh (1.25 km instead of 2.5 km) leading to a better description of orography that reduces the mismatch between model and station altitudes, allowing more stations to be assimilated. The assimilation of observations that are made more frequently would require an improved timeliness delivery of the data, in particular for future nowcasting applications based on NWP models. Another activity concerns the bias correction scheme where we want to replace the current static strategy (bias computed on a past period and applied to a more recent one), with a dynamical approach that could automatically account for changes in the bias as a function of relevant predictors with coefficients estimated in the minimisation of the variational assimilation (i.e. adaptive bias correction scheme).
Finally, another area of development that is currently under study within the COST ES1206 action 2 is to assimilate in NWP models additional information provided by the geodetic community that could give access to horizontal and vertical gradients of atmospheric humidity, such as horizontal gradients of zenith delays, slant delays, vertical profiles of refractivity or humidity deduced from tomography. All of these developments will help to make better use of tropospheric delays in NWP models, with expected improvements in terms of forecast skill scores.