Academia.eduAcademia.edu

Advanced GNSS Processing Techniques (Working Group 1)

2019, Advanced GNSS Tropospheric Products for Monitoring Severe Weather Events and Climate

Over the last decade, near real-time analysis of GPS data has become a well-established atmospheric observing tool, primarily coordinated by the EIG EUMETNET GPS Water Vapour Programme (E-GVAP) in Europe. In the near future, four operational GNSS will be available for commercial and scientific applications with atmospheric science benefiting from new signals from up to 60 satellites observed at any one place and time, however, many challenges remain regarding their optimal combined utilization. Besides raw data streaming, recent availability of precise real-time orbit and clock corrections enable wide utilization of autonomous Precise Point Positioning (PPP), which is particularly efficient for highrate, real-time and multi-GNSS analyses. New GNSS constellation signals, products and processing methods suggest the development of advanced GNSS tropospheric products, in support of weather

Chapter 3 Advanced GNSS Processing Techniques (Working Group 1) J. Douša, G. Dick, Y. Altiner, F. Alshawaf, J. Bosy, H. Brenot, E. Brockmann, R. Brožková, Z. Deng, W. Ding, K. Eben, M. Eliaš, R. Fernandes, A. Ganas, A. Geiger, G. Guerova, T. Hadaś, C. Hill, P. Hordyniec, F. Hurter, J. Jones, M. Kačmařík, K. Kaźmierski, J. Kaplon, P. Krč, D. Landskron, X. Li, C. Lu, J. P. Martins, G. Möller, L. Morel, G. Ófeigsson, R. Pacione, C. Pikridas, E. Pottiaux, J. Resler, W. Rohm, A. Sá, J. Sammer, T. Simeonov, W. Söhne, A. Stoycheva, A. Stürze, Sz. Rozsa, F. N. Teferle, S. Thorsteinsson, P. Václavovic, H. Valentim, B. Van Schaeybroeck, P. Viterbo, K. Wilgan, L. Yang, L. Zhao, N. Zinas, and F. Zus In the following sections material is republished with kind permission: 3.2.1, 3.3.2–3.3.7, 3.4.3, 3.4.6, 3.4.7. 3.4.8, 3.5.3–3.5.6, 3.5.8, 3.5.9, 3.6.2 and 3.6.6. J. Douša (*) Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] G. Dick GFZ German Research Centre for Geosciences, Helmholtz Centre Potsdam, Potsdam, Germany e-mail: [email protected] Y. Altiner · W. Söhne · A. Stürze BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected]; [email protected]; andrea.stuerze@bkg. bund.de F. Alshawaf · Z. Deng · X. Li · C. Lu · F. Zus GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected]; [email protected]; [email protected]; [email protected]; zusfl[email protected] J. Bosy · T. Hadaś · P. Hordyniec · K. Kaźmierski · J. Kaplon · K. Wilgan Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected]; [email protected]; [email protected]. pl; [email protected]; [email protected]; [email protected] H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] E. Brockmann Swiss Federal Office of Topography swisstopo, Wabern, Switzerland e-mail: [email protected] © Springer Nature Switzerland AG 2020 J. Jones et al. (eds.), Advanced GNSS Tropospheric Products for Monitoring Severe Weather Events and Climate, https://doi.org/10.1007/978-3-030-13901-8_3 33 34 J. Douša et al. R. Brožková Czech Hydrometeorological Institute, Prague, Czech Republic e-mail: [email protected] W. Ding · F. N. Teferle University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected]; [email protected] K. Eben · P. Krč · J. Resler Czech Institute of Computer Science, Academy of Sciences, Praha, Czech Republic e-mail: [email protected]; [email protected]; [email protected] M. Eliaš · P. Václavovic · L. Zhao Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected]; [email protected]; [email protected] R. Fernandes · H. Valentim University of Beira Interior, Covilhã, Portugal e-mail: [email protected]; [email protected] A. Ganas National Observatory of Athens, Athens, Greece e-mail: [email protected] A. Geiger · F. Hurter ETH Zurich, Zürich, Switzerland e-mail: [email protected]; [email protected] G. Guerova Physics Faculty, Department of Meteorology and Geophysics, Sofia University “St. Kliment Ohridski”, Sofia, Bulgaria e-mail: [email protected]fia.bg C. Hill · L. Yang University of Nottingham, Nottingham, UK e-mail: [email protected]; [email protected] J. Jones Met Office, Exeter, UK e-mail: jonathan.jones@metoffice.gov.uk M. Kačmařík Institute of Geoinformatics, VŠB Technical University of Ostrava, Ostrava, Czech Republic e-mail: [email protected] D. Landskron · G. Möller · J. Sammer Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected]; [email protected]; julia. [email protected] J. P. Martins · P. Viterbo Instituto Português do Mar e da Atmosfera, Lisbon, Portugal e-mail: [email protected]; [email protected] L. Morel École Supérieure des Géomètres et Topographes, Le Mans, France e-mail: [email protected] 3 Advanced GNSS Processing Techniques (Working Group 1) 35 G. Ófeigsson · S. Thorsteinsson The Icelandic Meteorological Institute, Reykjavík, Iceland e-mail: [email protected]; [email protected] R. Pacione e-GEOS/Centro di Geodesia Spaziale-Agenzia Spaziale Italiana, Matera, MT, Italy e-mail: [email protected] C. Pikridas Aristotle University of Thessaloniki, Thessaloniki, Greece e-mail: [email protected] E. Pottiaux Royal Observatory of Belgium, Brussels, Belgium e-mail: [email protected] W. Rohm Institute of Geodesy and Geoinformatics, Wrocław University of Environmental and Life Sciences, Wroclaw, Poland e-mail: [email protected] A. Sá Polytechnic Institute of Guarda, Guarda, Portugal e-mail: [email protected] T. Simeonov Sofia University “St. Kliment Ohridski”, Sofia, Bulgaria e-mail: [email protected] A. Stoycheva National Institute of Meteorology and Hydrology, Sofia, Bulgaria e-mail: [email protected] Sz. Rozsa Budapest University of Technology and Economics, Budapest, Hungary e-mail: [email protected] B. Van Schaeybroeck Royal Meteorological Institute of Belgium, Uccle, Belgium e-mail: [email protected] N. Zinas Tekmon Geomatics, Ioánnina, Greece e-mail: [email protected] Abstract Over the last decade, near real-time analysis of GPS data has become a well-established atmospheric observing tool, primarily coordinated by the EIG EUMETNET GPS Water Vapour Programme (E-GVAP) in Europe. In the near future, four operational GNSS will be available for commercial and scientific applications with atmospheric science benefiting from new signals from up to 60 satellites observed at any one place and time, however, many challenges remain regarding their optimal combined utilization. Besides raw data streaming, recent availability of precise real-time orbit and clock corrections enable wide utilization of autonomous Precise Point Positioning (PPP), which is particularly efficient for highrate, real-time and multi-GNSS analyses. New GNSS constellation signals, products and processing methods suggest the development of advanced GNSS tropospheric products, in support of weather 36 J. Douša et al. numerical prediction and nowcasting will be substantially improved. Such examples are: ultra-fast and high-resolution tropospheric products available in real-time or on a sub-hourly basis, parameters monitoring tropospheric anisotropy above the station (such as horizontal gradients and tropospheric slant path delays), and indicators of severe weather such as extreme convection. Development of advanced GNSS tropospheric products within COST Action ES1206 benefited from two dedicated campaigns prepared for a collaborative effort: (1) the benchmark campaign and (2) the real-time demonstration campaign. The former served for estimating and assessing horizontal tropospheric gradients and tropospheric slant delays, estimated from GNSS, Water Vapour Radiometers and Numerical Weather Model (NWM) ray-tracing. The second campaign developed new software and strategies for realtime, multi-GNSS, high-rate tropospheric solutions including the assessment of pre-operational solutions. The impact of selected processing strategies and precise models were assessed during a long-term GNSS reprocessing campaign aimed at providing homogeneous tropospheric products for climate research. Using information from modern NWM forecasting systems, a variety of tropospheric correction models for real-time kinematic GNSS positioning were developed and assessed. Finally, a transfer of knowledge such as support for establishing new GNSS Analysis Centres and inclusion of new networks into E-GVAP were completed. 3.1 Introduction J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] G. Dick GFZ German Research Centre for Geosciences, Helmholtz Centre Potsdam, Potsdam, Germany e-mail: [email protected] The GNSS4SWEC Working Group 1 (WG1) focused on development and utilization of advanced GNSS processing techniques for the purpose of both estimation and exploitation of tropospheric parameters within geodetic, meteorological and climate applications. The main goals of the WG1 are summarized within four defined domains: 1. Coordinating the development of advanced GNSS tropospheric products in support of weather forecasting, namely such as real-time parameter estimation, troposphere asymmetry monitoring and modelling, developing severe weather indicators, assessing the impact of hydrometeors, advantage of multi-GNSS data processing. 2. Exploiting NWM data in GNSS precise positioning and real-time kinematic applications, namely supported by mapping functions or mapping factors, a priori separation of hydrostatic contributions, tropospheric horizontal gradients, tropospheric parameter scaling and conversions factors. 3. Reprocessing of GNSS data and assessing precise models for the purpose of a longterm consistent tropospheric product provision in support of a climatology research. 3 Advanced GNSS Processing Techniques (Working Group 1) 37 4. Supporting the transfer of knowledge, tools and data exchange for extending existing products and establishment of new ACs or inclusions of new networks. WG1 coordination was split into ten different sub-tasks, with many overlaps. For this reason, the structure of this chapter does not correspond to the WG1 sub-tasks, but reflects mainly goals specified above. An important role in new development was achieved by a collection of common data sets and the design of specific campaigns, which are described in Sect. 3.2. Development of advanced tropospheric products specified in the first goal are introduced in Sects. 3.3 and 3.4. Utilization of NWM-based products in precise GNSS analyses, i.e. corresponding to the second goal, is covered by Sect. 3.5. Various GNSS reprocessing and processing model assessments are then summarized in Sect. 3.6. Finally, the transfer of knowledge, the establishments of new ACs and integration of new networks is completed in Sect. 3.7. The definition of the SINEX_TRO V2 format, which was elaborated in a close collaboration with WG3, is included in appendix D. 3.2 Campaigns for Development of Advanced Tropospheric Products The main goal of this section is to introduce two campaigns suitable for developing and evaluating advanced tropospheric products. The first is the so-called GNSS4SWEC WG1 Benchmark campaign which is considered as a cornerstone that helped to accomplish various WG1 objectives within the COST Action. The second is the Real-time demonstration campaign, which helped to develop, optimize and evaluate new real-time GNSS software and products. 3.2.1 Benchmark Campaign – Common Data Set for New Product Development and Validation1 M. Kačmařík Institute of Geoinformatics, VŠB Technical University of Ostrava, Ostrava, Czech Republic e-mail: [email protected] J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] G. Dick GFZ German Research Centre for Geosciences, Helmholtz Centre Potsdam, Potsdam, Germany e-mail: [email protected] 1 Parts from this section were previously published in Douša et al. 2016 38 J. Douša et al. F. Zus GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: zusfl[email protected] R. Brožková Czech Hydrometeorological Institute, Prague, Czech Republic e-mail: [email protected] H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] A. Stoycheva National Institute of Meteorology and Hydrology, Sofia, Bulgaria e-mail: [email protected] G. Möller Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] J. Kaplon Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] Only basic information is presented here about the Benchmark campaign. For further information please see Douša et al. (2016). 3.2.1.1 Motivation An idea to create a well-prepared and extensive common data set which would enable an effective collaboration within the WG1 itself, but also with other working groups, arosed during the first WG1 meeting in Valencia, October, 2013. The natural motivation for the campaign was to support the WG1 main goals, namely the development of advanced GNSS tropospheric products in support of weather forecasting and vice versa exploitation of numerical weather data in precise GNSS positioning. The Benchmark campaign planning started with inventory of requirements based on a wide discussion within the Action members. The following requests for the Benchmark data set were summarized as follows: • Period covering a month at least to enable NWM and GNSS processing initializations and to cover different weather conditions – quiet and variable, optimally including a severe weather event. • Availability of a dense network of GNSS reference stations in Europe with a limited scale, but including flat and mountainous areas. • Availability of meteorological data from independent sources. 3 Advanced GNSS Processing Techniques (Working Group 1) 3.2.1.2 39 Description of Selected Spatial and Temporal Domain Finally, an area in central Europe was selected covering Germany, the Czech Republic, Poland and Austria. Due to its size, the whole area was firstly divided into a ‘core’ domain where below mentioned severe weather events took place and an ‘extended’ domain which surrounded the ‘core’ one. In the second step both domains were geographically divided into several clusters to allow a reasonable GNSS data handling, see Fig. 3.1. From the time perspective 2 months in June 2013 (May and June) were chosen. The weather conditions in the selected area during May 2013 were mostly quiet. On the contrary an extreme precipitation event lasting from May 31 to June 2 led to devastating flooding on Danube, Elbe and Vltava rivers. Since the event was only partly forecasted by NWM it was a suitable candidate for a GNSS meteorology benchmark campaign as it could provide additional observations for meteorological community. Significant precipitation periods hitting areas of smaller extent occurred also from June 9 to June 11 and from June 23 to June 26. 3.2.1.3 Description of Collected Data Set The Benchmark data set contains following data: GNSS observations and auxiliary products, E-GVAP operational GNSS products, synoptic meteorological observations, NWM fields, radiosonde observations, WVR observations, meteorological radar images. Collected data from individual sources are briefly described in paragraphs below. The data set is stored on an ftp server at Geodetic Observatory Pecný Fig. 3.1 Benchmark core (yellow area) and extended domains depicted together with nine clusters for GNSS stations (coloured points). The size of the points indicates height of GNSS reference stations above the WGS-84 ellipsoid 40 J. Douša et al. (GOP) and available to members of the COST ES1206 Action for research purposes. To obtain more information and access to the data set please send an email to michal. [email protected] or [email protected]. GNSS Observations observations from 430 GNSS reference stations in RINEX format with 30s sampling interval were collected in total from which 247 sites belonged to the ‘core’ domain. An average distance between two stations was about 50–70 km. From the total number of GNSS sites, 4 observed GPS, GLONASS and GALILEO satellites, 356 observed GPS and GLONASS satellites, and remaining 70 stations were equipped with GPS recievers only. Station metadata files were completed and checked carefully. A qualitative and quantitative control and a standard positioning were performed using G-Nut/Anubis software (Václavovic and Douša 2016). Besides multiple correction of metadata, 15 sites had to be rejected from the data set because of the data quality issues. E-GVAP Operational GNSS Products operational near real-time tropospheric solutions provided by 14 analysis centres for all GNSS reference stations in Europe were collected for the campaign. These products contributed routinely to the E-GVAP (http://egvap.dmi.dk) and are stored in the COST-716 format with a temporal resolution of ZTD estimates from 5 to 60 min. Synoptic Meteorological Data meteorological measurements of at least atmospheric air pressure, air temperature and relative humidity from 610 synoptic stations were collected. Original data were provided in various formats which were, additionally, converted into a single unified plain text format with a sampling interval ranging from 10 to 60 min. NWM Data and Products NWM 3D data fields from the Czech Hydrometeorological Institute’s (CHMI) local area model ALADIN-CZ were extracted in GRIB format. The horizontal resolution of the model is 4.7  4.7 km, outputs are provided at 87 model levels with a 6-h interval of analysis run (00:00, 06:00, 12:00 and 18:00 UTC) and 1-h interval for forecast range. The model didn’t assimilate any GNSS tropospheric products for delivered fields and provided necessary parameters for derivation of hydrostatic and non-hydrostatic GNSS signal path delays as well as for calculating the effect on signal due to so called hydrometeors (e.g. ice and liquid water). Radiosonde Data radio soundings from two different sources were collected providing profiles with full and reduced resolutions. Measurements with high resolution were available from two sites in the Czech Republic – Prague-Libuš and Prostějov, both provided by the CHMI. Altogether 278 files were obtained for the period of Benchmark campaign. Radiosonde data with reduced vertical resolution from 19 European stations were provided by E-GVAP based on the EUMETNET – EUREF MoU (Pottiaux et al. 2009). Water Vapour Radiometer Data observations from two Water Vapour Radiometers (WVR) situated in Germany were collected. The first one was operated by German Research Centre for Geosciences (GFZ) in Potsdam (POTS) 30 km south-westward 3 Advanced GNSS Processing Techniques (Working Group 1) 41 from Berlin and provided measurements of IWV and liquid water in two modes: slant (GPS satellite tracking) and zenith. The second was operated by Deutscher Wetterdienst (DWD) at the Lindenberg meteorological observatory (LDBG) located approximately 100 km eastward from Berlin and provided the same parameters as the first mentioned WVR but only in the zenith direction. Meteorological Radar Data raster images of combined observations from two C-band Doppler meteorological radars located in the Czech Republic and operated by the CHMI were collected. They represent a maximum reflectivity fields with side projections in horizontal resolution of 1  1 km and 30-min time interval. The area effectively covered by those two radars includes the territory of the Czech Republic and areas of approximately 100 km outside the Czech state boundary. Data Acknowledgement the members of COST Action ES 1206 thank all the institutions that provided data for the campaign. GNSS data from the Austrian network EPOSA were provided by Österreichische Bundesbahnen Infrastruktur AG; GNSS data from SAPOS network in Germany by Zentrale Stelle SAPOS in Hannover; GNSS data from several networks in the Czech Republic – (1) CZEPOS by the Czech Land Survey Office, (2) Trimble VRS Now® by GEOTRONICS Praha, s.r.o. and (3) GEONAS and VESOG stations thanks to the project CzechGeo (LM2010008) operated by the Institute of Rock Structure and Mechanics of the Academy of Sciences of the Czech Republic and the Research Institute of Geodesy, Topography and Cartography, respectively; GNSS data from Polish ASG-EUPOS network by the Head Office of Geodesy and Cartography in Poland; Synoptic data by Zentralanstalt für Meteorologie und Geodynamik, ZAMG (Austria), Deutscher Wetterdienst, DWD (Germany), Czech Hydrometeorological Institute, CHMI (the Czech Republic) and Polish Institute of Meteorology and Water Management (Poland). Finally, special thanks come to people assisting in collecting all the data and metadata, namely: Dr. Jan Řezníček (GNSS/CZEPOS), Dr. Uwe FeldmannWestendorff and Dr. Markus Ramatschi (GNSS/SAPOS), Dr. Petr Novák (RADAR/ CHMI), Dr. Martin Motl (RAOBS/CHMI), Dr. Anna Valeriánová (SYNOP/CHMI), Dr. Pavla Skřivánková (CHMI), Dr. Roland Potthast (SYNOP/DWD), Dr. Jürgen Güldner (WVR-Lindenberg/DWD) and Dr. Stefan Heise (WVR-Potsdam/GFZ). 3.2.1.4 Reference Products and Their Initial Validation After the data preparation, quality checking and cleaning, GNSS and NWM reference (tropospheric) products were generated and consequently initially validated. Results of these steps provided a first insight into variations of parameters and atmospheric conditions during the Benchmark campaign and were helpful for more detailed planning of following Benchmark-related activities. Following paragraphs provide basic information about mentioned reference products and their validation. 42 J. Douša et al. GNSS Reference Tropospheric Products two institutions delivered their GNSS tropospheric products for all GNSS reference stations within the Benchmark data set. The first reference one was generated at GOP using the BSW52 (Dach et al. 2015) with the network processing approach using double-differenced GNSS observations. The strategy for daily solutions was consistent with the GOP contribution to the EUREF Repro2 campaign (Douša et al. 2017). The second reference tropospheric product was delivered by GFZ using the GFZ EPOS software (Gendt et al. 2004; Ge et al. 2006) based on undifferenced GNSS observations and PPP approach. GOP solution used CODE precise orbit and clock products while GFZ solution was based on GFZ own precise products. The models were in case of both GOP and GFZ solutions compliant with the IERS Conventions (Petit and Luzum 2010). NWM Derived Tropospheric Products analogously to the situation with GNSS reference products also NWM derived reference tropospheric products were independently generated by GOP and GFZ. GOP used their G-Nut/Shu software (Douša and Eliaš 2014) and provided three products. Two were based on global numerical weather models – the European Centre for Medium-Range Weather Forecasts (ECMWF) ERA-Interim (Dee et al. 2011) and the Global Forecast System (GFS) of the National Centres for Environmental Prediction (NCEP) available at http:// www.ftp.ncep.noaa.gov/data/nccf/com/gfs/prod/, and one was based on regional model ALADIN-CZ which was provided within the Benchmark dataset. The following tropospheric and meteorological parameters were calculated for all stations: zenith hydrostatic and wet delays, air pressure, partial water vapour pressure, mean temperature, temperature lapse rate, water vapour pressure and zenith wet delay exponential decay rates. GFZ used their direct numerical simulation (DNS) tool (Zus et al. 2014) in order to derive the following parameters from two global NWMs (ERA-Interim and GFS): zenith hydrostatic and zenith wet delays, horizontal (1st, 2nd order) tropospheric gradients and coefficients of hydrostatic and wet mapping functions. Validation of Reference Tropospheric Products Values of ZTD, ZWD and horizontal tropospheric gradients derived from described GNSS and NWM reference products were compared to study the quality and mutual agreement. Table 3.1 summarizes comparison results of GNSS ZTDs with those derived from NWMs. Mean statistics over all 430 sites demonstrated that both GNSS reference products based on a completely different software and strategy performed very similarly when compared to NWM products. Also, both software used for derivation of GNSS related tropospheric parameters from NWM fields agreed very well. In case of individual NWM the high-resolution ALADIN-CZ model outperformed both global reanalysis models in the Benchmark domain and period mainly in terms of standard deviation. For the NCEP’s GFS products a negative mean bias of about 5 mm was observed compared to all other solutions. This bias stems from the bias in ZWD values and its possible explanation is the low vertical resolution of GFS model which resulted in larger interpolation errors. A comparable bias for GFS model was reported by Urquhart et al. (2011). NWM source (software) ERA (Shu) ERA (Shu) ERA (DNS) ERA (DNS) GFS (DNS) GFS (DNS) ALADIN (Shu) ALADIN (Shu) Grid resolution 1 deg 1 deg 1 deg 1 deg 1 deg 1 deg 4.7 km 4.7 km Analysis [hour] 6 6 6 6 6 6 6 6 Forecast [hour] 0 0 0 0 3 3 0,1,2,3,4,5 0,1,2,3,4,5 GNSS source (software) GOP (BSW52) GFZ (EPOS-8) GOP (BSW52) GFZ (EPOS-8) GOP (BSW52) GFZ (EPOS-8) GOP (BSW52) GFZ (EPOS-8) Pairs # 224 224 224 224 224 223 1343 1343 Excl # 2 3 3 3 7 7 20 22 Bias [mm] +0.0 +0.3 0.4 0.1 4.9 4.5 +0.8 +0.6 Sdev [mm] 9.6 9.7 9.4 9.6 11.0 10.9 7.6 7.3 RMS [mm] 10.0 10.0 9.8 9.8 12.0 11.8 7.8 7.5 3 Advanced GNSS Processing Techniques (Working Group 1) Table 3.1 Comparison of zenith total delays from NWM and GNSS (mean values, 430 sites) 43 44 J. Douša et al. Results of individual GNSS and NWM reference products were studied also using maps showing differences between ZTDs from GNSS GOP product and individual NWM products. Generally, a good homogeneity was observed from the statistical results. Exceptions existed mainly in relation to the orography which triggered larger differences between models particularly in mountain areas with complex terrain where regional model ALADIN-CZ performed better than both global models. 3.2.1.5 List of Benchmark Campaign Participants and Users • G. Möller (TU Wien). In the first instance, we will use RINEX data to compute STDs using different processing strategies. In addition, the GNSS GOP products are used as input for GNSS tomography and assimilation studies. • E. Pottiaux (ROB). The idea is to use the Benchmark campaign to assess the different tropospheric products (RT, sub-hourly, hourly and post-processing – ZTD, gradients and possibly slant delays from the BSW52) and refine my processing strategies. By taking part myself in the Benchmark processing, I also would like to stimulate the interaction with WG2 developments. • J. Douša et al. (GOP). Generation of the reference GNSS tropospheric products. Development of the software and the strategy for optimal provision of real-time products and asymmetry monitoring, slant delay retrievals. Assessment of NWM-derived tropospheric products and developing the combination of GNSS and NWM data. • G. Dick et al. (GFZ). Generation of the reference GNSS tropospheric products. • M. Kačmařík (TU Ostrava). Inter-comparison of tropospheric slant delays form GNSS, NWM and WVR at dual-station. • W. Rohm et al. (WUELS). We would like to join the Benchmark campaign by providing Slant Delays for selected stations (same as Michal’s) based on the GNSS observations and ray-tracing. • H. Brenot (BIRA). I would like to look at hydrometeors delays from ALADIN outputs and to show their impact for the severe flood event of June. • K. Eben (ICS ASCR). Assimilation of NRT ZTDs into the WRF mesoscale model. • T. Hadaś et al. (WUELS). The goal is to use the stored real-time products for simulated real-time troposphere monitoring (ZTD estimates using original software), to optimize the methodology and algorithms, to compare results with other real-time AC. • L. Morel (Le CNAM). To process some stations of the benchmark campaign and to deliver CNAM results (GAMIT processing). • S. Nahmani (IGN). I want to verify some of my results on these specific stations. • P. Gołaszewski (UWM). My research is focused on real time and post-processed ZTD/ZWD estimation. Using this data will allow me to present the results on the COST workshop in Potsdam, in September this year. I am interested in using 3 Advanced GNSS Processing Techniques (Working Group 1) • • • • 45 pseudo real-time demonstration (for ZTD estimation) and observation data for post-processing. S. de Haan (KNMI). My plan is to eventually assimilate the slant observations in Harmonie. But first, I am going to compare the observation with my model equivalent. Arpaci (UBIMET). We will carry out case studies to examine the impact of the assimilation on the forecast quality. Especially the catastrophic flooding from June 2013, which affected big parts of central Europe seems to be an ideal evaluation case study scenario. We will compute WRF model runs with ZTD assimilation and compare them with runs without data assimilation. UBIMET will carry out eye to eye verifications, Upper Air verifications (using sounding and aircraft data) and surface verifications using VERA analysis data. D. Kwasniak (UWM). I want to use them for my research about GPS positioning using a new positioning method called MAFA method. Results of this research I want to present on 17th Czech-Polish Workshop. Y. Altiner et al. (BKG). Want to obtain the access for a processing the Benchmark GNSS data in real-time simulated mode. 3.2.2 GNSS Real-Time PPP Demonstration Campaign P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] F. N. Teferle University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected] Providing new real-time or ultra-fast tropospheric products, such as ZTD, GRD, STD, IWV maps or other derived products estimated using data from GNSS permanent networks, is interesting for numerical and non-numerical weather nowcasting and severe weather event monitoring (Guerova et al. 2016a, b). The Precise Point Positioning (PPP) processing strategy plays a key role in the production of real-time tropospheric parameters because of its high processing efficiency, and the sensitivity to the absolute value of the tropospheric delay. It enables to exploit optimally data from all available GNSS multi-constellations, and facilitates the production of all interesting GNSS parameters such as ZTDs, GRDs or STDs. 46 J. Douša et al. Most importantly, the PPP is supported with the global orbit and clock products provided by the RTS (Caissy et al. 2012) of IGS (Dow et al. 2009). During the period 2015–2017, the COST Action ES1206 WG1 played an initiative role in the coordination of the development and the evaluation of GNSS realtime tropospheric products thanks to the design of the GNSS4SWEC Real-time Demonstration Campaign which is briefly introduced in this subsection. 3.2.2.1 Campaign Design and Contribution Specifications A list of common stations was selected from the E-GVAP supersites, EPN and IGS sites, which were limited by a maximum of 50 in total, regionally and globally distributed: • E-GVAP super-sites (5): BRST, GOPE, ONSA, YEBE, ZIM2 • EPN sites (10): CASC, HERT, HOFN, MALL, MATE, NICO, PDEL, POTS, REYK, WTZR • IGS sites (17): ADIS, ALBH, ALGO, ALIC, AUCK, DUBO, LHAZ, NKLG, NRMD, OHI3, POVE, THTI, ULAB, UNSA, WIND, YAR3, YELL Station metadata are introduced using the RINEX skeleton files which are available from EPN CB (http://www.epncb.oma.be/stations/log/skl) and IGS CB (http://igscb.jpl.nasa.gov/igscb/station/general/skel). The strategy for tropospheric estimation is generally free when optimized with respect to the individual software capabilities and following the state-of-the-art models, in particular IERS conventions, antenna phase offsets and variations, a priori tropospheric model, and mapping functions, and others. Use of GNSS systems GPS and GPS + GLONASS is recommended if supported by the software and real-time data streams, other systems are optional, however, the constellation has to be properly defined in the header of the COST-716 file using the PCD flag. The GPS and GPS + GLONASS solutions are provided in different files and different names and individual solution (including variants) has to be specified by the fourth character in the processing centre name. The character ‘G’ will be used for GPS-only results, and the character ‘R’ will be used for GPS + GLONASS. In order to support GPS + GLONASS, the IGS03 real-time product is mandatory for the utilization to guarantee a consistency of mandatory products and enable to compare the results. Thus even GPS-only solution should use the IGS03 stream, for others solutions the precise ephemeris source is optional. Parameters must be estimated with a 5-min resolution (parameter sampling rate). If the higher sampling rate is used in the processing delivered product files should be reduced to 5 min. Parameters to be estimated: ZTD (mandatory) with the product sampling rate of 5 min (processing sampling rate can be higher), horizontal tropospheric gradients (optionally) and coordinates (mandatory) estimated as static parameters. The contributors submit files with troposphere parameters every hour to the ftp-server at the Geodetic Observatory Pecny after the registering of the product. The product file is 3 Advanced GNSS Processing Techniques (Working Group 1) 47 converted to the latest COST-716 format with the file name following the COST-716 conventions, i.e. using “demo” product status. Analysis centre providing more product lines should be uploaded as separate COST-716 files when using a specific analysis centre acronym (i.e. consisting of a unique analysis centre and a product line). 3.2.2.2 Contribution and Monitoring From April 2015 till the end of the COST Action, eight agencies succeeded to start the real-time processing and provide partial contribution at least to the GNSS4SWEC Real-Time Demonstration Campaign. Real-time/ultra-fast solutions were provided using six different software and using various flavours of processing options (Table 3.2). Truly real-time solutions using the operational processing engine was provided by seven contributors. For the purpose of a feedback to such product providers, a dedicated web service for an easy monitoring and comparison of individual contributions Fig. 3.2. The access has been made available also to a wide community (in particular interested from the GNSS4SWEC WG2) at http://www.pecny.cz/COST/RT-TROPO to enable visualising site-specific time series of recently estimated ZTD and gradient parameters from real-time solutions. For the comparison purpose, also near real-time regional and global solutions from the GOP analysis centre operationally contributing to the EIG EUMETNET GNSS Water Vapour Programme, E-GVAP (http:// egvap.dmi.d) were included in the real-time demonstration monitoring. Table 3.2 Contributions to GNSS4SWEC Real-Time Demonstration campaign AC GOP Running agency Geodetic Observatory Pecný, RIGTC TUW Technical University Vienna ROB Royal Observatory of Belgium ASI ULX Agenzia Spaziale Italiana/Centro di Geodesia Spaziale, Matera University of Luxembourg Software G-Nut/ Tefnut TUW software G-Nut/ Tefnut GipsyOasis BNC TUO Technical University of Ostrava RTKLib 15.6. 2015 5.11.2015 BKG Bundesamt für Kartographie und Geodäsie GFZ German Research Centre for Geosciencies BNC 1.3.2016 EPOSRT 16.2.1017 GFZ Start 9.4. 2015 15.4. 2015 23.4. 2015 5.5. 2015 Update Realtime Realtime Realtime Hourly Realtime Realtime Realtime Realtime Solutions GPS, GLO, gradients GPS GPS, GLO, gradients GPS, gradients GPS GPS GPS, GLO GPS, GLO 48 J. Douša et al. Fig. 3.2 Web service for the monitoring of the GNSS4SWEC Real-Time Demonstration campaign 3.2.2.3 Link to 4.3.7 IAG Working Group The activity within the GNSS4SWEC project plays also a key role in the IAG Working Group 4.3.7 ‘Real-Time troposphere monitoring’ which has been established with the following objectives for the period of 2015–2019: • Stimulate the development of software that enable routine production of realtime/ultra-fast tropospheric products. • Develop optimal strategies suitable for numerical or non-numerical weather nowcasting applications, and severe weather event monitoring. • Demonstrate a reliable high-temporal resolution real-time/ultra-fast production, assess applied method, software and precise real-time orbit and clock products. • Evaluate real-time/ultra-fast tropospheric parameters and their potential for applications in meteorology. • Setting up a link to the users, review product format and requirements. 3.3 Tropospheric Asymmetry Monitoring and Advantage of Multi-GNSS This section focuses on the detection of meteorological heterogeneities surrounding ground-based GNSS stations, in particular the disturbance of GNSS signal through the neutral atmosphere as retrieved by geodetic software with gradient and residual contributions to Slant Total Delay (STD). The STD of the neutral atmosphere, measured by GNSS technique, is the result of the adjustment of two components – 3 Advanced GNSS Processing Techniques (Working Group 1) 49 isotropic and anisotropic. The ZTD of the neutral atmosphere represents the isotropic contribution above a GNSS site. To adjust the anisotropic contribution, the concept of horizontal gradients has been introduced in GNSS software. Additionally, this section included selected results of the multi-GNSS processing, which is expected to foster and improve tropospheric parameters, in particularly a possibility to monitor anisotrophy by estimating tropospheric horizontal gradients or retrieving STDs from carrier-phase post-fit residuals for more satellites in view. 3.3.1 Concept of Tropospheric Gradients H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] The accuracy of slant delay measurements of neutral atmosphere (Latm) and the number of visible satellites are critical for identifying the exact location of smallscale asymmetric tropospheric structures (Fig. 3.3). STD ¼ Latm ðє; αÞ ¼ 10 6 Z S ! N P ds ð3:1Þ Latm is commonly called slant tropospheric delay (for the elevation є and the azimuth α). Its formulation shows a dependency on the atmospheric refractivity (N ) ! and the position of source of signal ( P ). The adjustment of the slant path delay of GNSS signal (Latm) through a blob of water vapour, can be well retrieved taking into account the 1st order anisotropic character of the neutral atmosphere (model presented Fig. 3.3; see Gradinarsky 2002). Fig. 3.3 Illustration of tropospheric heterogeneity affecting path travel (S) of GNSS signal. P⃗ is the vector position along the satellite direction (S) at the elevation (є), and ρ is its projection on the horizontal plan 50 J. Douša et al. Fig. 3.4 Example of the distribution of the atmospheric refractivity (N) above a GNSS site for which a model of a flat tilted atmosphere is considered to retrieve delay gradient. The tilted line on the top, corresponds to the isoline for a constant refractivity N0 The atmospheric anisotropy reflects the fact that the speed of microwave signals emitted from satellites and recorded by ground-based receivers differs according to the azimuthal direction through the neutral atmosphere (and the asymmetrical distribution of components surroundings a GNSS station). An inclined plane model of troposphere (Fig. 3.4) schematised by linear thickness and density variations is considered to define horizontal gradients during the adjustments of tropospheric parameters (Davis et al. 1993; Gradinarsky 2002). The correction provided by gradients possesses its own mapping function, mfaz (Chen and Herring 1997). The expression of azimuthal anisotropic contribution (Laz) to the reconstruction of slant delays depends on the satellite direction (elevation and azimuth). According to Davis et al. (1993), the first order Taylor expansion (when ρ tends to 0⃗) of the atmospheric refractivity (N ) can be formulated !  ! !  ∂N ρ ; z  ! N ρ ; z ¼ N o ðzÞ þ ¼ N o þ ξ ðzÞ: ρ  !  ! ! ∂ρ ρ !0 ð3:2Þ This expresses the 1st order anisotropy around a site. The assumption of a straight line propagation is considered, with no time dependency of N (Gradinarsky 2002). ! N0 is the isotropic contribution to N with cylindrical symmetry, and ξ is the horizontal refractivity gradient. Both depend on the altitude z (function of the elevation є). Using expression Eq. 3.2 injected in Eq. 3.1, the formulation of Latm is obtained: Latm ðє; αÞ ¼ 10 6 Z S N o ðzÞ ds þ 10 6 Z ¼ Lsym ðєÞ þ Laz ðє; αÞ S ! ! ξ ðzÞ: ρ ds ð3:3Þ Lsym is the symmetric contribution to ZTD mapped in direction of satellite, and Laz is the asymmetric contribution to the total delay (with azimuthal dependency). ! Vector position ρ in the horizontal plan can be expressed using the orthonormal ! !  base uNS ; uEW of the North-South (NS) and East-West (EW) directions: 3 Advanced GNSS Processing Techniques (Working Group 1) ! ρ 51   ! ! ¼ z: cot ðєÞ: cos ðαÞ uNS þ sin ðαÞ uEW ð3:4Þ A way to link the differential of length along the path of GNSS signal (ds) and the differential of length along the vertical (dz) is to use mapping function (Niell 1996; Böhm et al. 2006a, b) which can depend on the elevation (m0(є)) ds ¼ m0 ðєÞ : dz ð3:5Þ The refractivity gradient can be expressed in the base ! ! ! uNS ; ! ! uEW  : ξ ðzÞ ¼ ξNS ðzÞ uNS þξEW ðzÞ uEW ð3:6Þ Considering Eqs. 3.3, 3.5 and 3.6, Laz can be formulated Laz ðє; αÞ ¼ 10 6 : m0 ðєÞ: cot ðєÞ:ðGNS : cos ðαÞ þ GEW : sin ðαÞÞ ð3:7Þ !  (G , GEW) are the components of the horizontal delay gradients G in the base  !NS ! uNS ; uEW expressed as ! ! ! G¼ GNS uNS þGEW uEW ¼  GNS GEW  ð3:8Þ The analytic formulation of (GNS, GEW) is following  GNS GEW  0Z 1 1 z:ξ ð z Þ dz NS B C C Z0 ¼B @ 1 A z:ξEW ðzÞ dz ð3:9Þ 0 The unit of delay gradient is the same as for the delay – the unit of length [m]. The ! delay gradient G can be defined as a correction of phase residual projections depending on the elevation and azimuth angles of visible satellites. An interpretation of the horizontal gradient is for example that, a gradient component of 1 mm at the zenith will show a higher correction for elevation of 45 and 25 (correction of about 2 mm and 4 mm, respectively). Following the implementation of tropospheric delay (ZTD), the horizontal delay gradient is the second tropospheric parameter implemented in the least-squares adjustment proceeded by geodetic software in the analysis of GNSS data (Davis et al. 1993; MacMillan 1995; Alber et al. 1997; Chen and Herring 1997; Bar-Sever and Kroger 1998). Initially, horizontal gradients were introduced into the calculations in order to improve positioning solutions, showing a 15% improvement for the horizontal repeatability (Bar-Sever and Kroger 1998) and a 25% improvement 52 J. Douša et al. between wet delays from GNSS and WVR. Afterwards, the potential of delay gradient for GNSS meteorology has been investigated. Walpersdorf et al. (2001) have showed how the use of GNSS gradient can describe the approach of a front towards Marseille in the south-east of France in 1998. Iwabuchi et al. (2003) also showed that the temporal and spatial variations of GNSS gradients matched well with the moisture field determined by ZTD and with the meteorological condition in summer 1996 over the Japan Islands (in particular during the passage of a weather front). The interest of using GNSS gradient for monitoring severe weather is still investigated, especially to improve our understanding of meteorological situation. For example, during the flash-flood event of September 2002 (southeastern France), three phases have been identified by Delrieu et al. (2005) and confirmed by GNSS data (Brenot et al. 2006). The maximum daily precipitation reached 691 mm. Figure 3.5 shows the path of a MCS over CHRN station during the Phase III (from 01:00 to 18:00 UTC on 9 September 2002). It can be seen that horizontal gradient can clearly be used to monitor tropospheric structure (MCS or blob of water vapour). Just after this rainfall event (in Autumn 2002), a dense network of 20 GPS stations has been installed by the Mediterranean Hydrometeorological Observatory Cévennes-Vivarais (OHM-CV) to proceed monitoring of convective systems and Cévenol effect, eastwards of Mont Aigoual (AIGO station located at 30 km of this dense network; see Fig. 3.5). Using different configurations in calculations (settings, geometry, constrains), the precision of NS and EW delay gradients components (zenith direction) has been evaluated to 0.35–0.7 mm and 0.2–0.5 mm respectively (Brenot et al. 2014a, b). The lack of GPS satellites at low elevation in the North for a Fig. 3.5 (left) Time-series of delay gradients over Château-Renard (CHRN) GPS station during the flash-flood event of September 2002 (over Bouche-du-Rhône, between Montpellier and Marseille, France). Purple double-arrow shows the path of a quasi-stationary Mesoscale Convective System (MCS) over CHRN station. (right) GNSS delay gradients are superposed over radar reflectivity localising the MCS close to MTPL, CHRN and VERC stations. Dash circles show representativity areas of gradients (for a cutoff angle of 10 ). 3 Advanced GNSS Processing Techniques (Working Group 1) 53 network located at a latitude of 45 N, explains why the precision of the NS component is less good than the EW one. 3.3.2 Global Validity and Behaviour of Tropospheric Gradients Estimated by GPS2 L. Morel École Supérieure des Géomètres et Topographes, Le Mans, France e-mail: [email protected] Estimation of tropospheric gradients in GNSS data processing is a well-known technique to improve positioning. Today, they are routinely estimated by several global and regional GNSS analysis centres but they are still not yet used for operational meteorology. We have studied the physical meaning of tropospheric gradients estimated from GPS observations recorded by several permanent stations located all around the world. In a first study with several stations on Corsica island, we estimated ZTD and tropospheric gradients using two software: GAMIT/GLOBK (GAMIT version 10.5) and GIPSY-OASIS II version 6.3 in order to analyse the differences in the tropospheric results (ZWD and gradients) coming from the processing strategy (double-differences for GAMIT/Globk versus zero-difference for GIPSY-OASIS). That study allowed to confirm a strong correlation between the two software for ZWD estimation (98%) and a good correlation for gradient estimation (70%). No direct correlation with elevation or geographical location has been noticed but the gradients were oriented inward land (Fig. 3.6), in opposite direction from tropospheric humidity field processed by ERA – Interim and with a direction relatively stable along the year (Morel et al. 2014). In a following study with 14 stations all around the world, selected due to their proximity of the relief, we also observed that gradient directions were stable over the time and pointed toward the relief for most of the stations selected. Correlation coefficients were processed between gradients (yearly mean values (Ge, Gn) as vector component) and direction of the steep slopes (obtained by analysing Digital Elevation Model at 20 km, 40 km and 60 km around the station), see Table 3.3. These results gave us a first step for a physical meaning to gradients when stations are close to high mountains. We can notice 10 stations with a correlation coefficient > 0.4 (60 km) and 2 stations without any correlation (BOGT and CHWK) but surrounded by mountains. Now, we are going to continue the study with more stations and years and quantify multipath effect. 2 Parts from this section were previously published in Morel et al. 2014 54 J. Douša et al. Fig. 3.6 Monthly mean gradient magnitude and direction (blue arrows estimated by GAMIT and red arrows estimated by GIPSY-OASIS) for all stations in Corsica Island over the year 2011. Gradient vectors are drawn considering monthly mean values (Ge, Gn) as vector components Table 3.3 Correlation coefficient between gradient and relief around 14 stations (red: negative correlation, orange: weak correlation and green: strong correlation) 3 Advanced GNSS Processing Techniques (Working Group 1) 3.3.3 55 Monitoring of Severe Weather from Wet Gradients, Residuals and Slants3 H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] To retrieve Slant Wet delays (SWD) in direction of GNSS satellites (for elevation є and azimuth α), two contributions are commonly considered (see Eq. 3.10): the wet isotropic contributions (Lsym ) with spherical symmetry and derived from Zenith Wet wet Delays (ZWD), and the anisotropic contributions (Lasym ) with azimuthal asymmetry and derived from horizontal wet gradients and residuals. wet wet SWD ¼ Lsym ðєÞ þ Lasym ðє; αÞ ð3:10Þ wet wet Lsym ðєÞ ¼ ZWDadjusted  mf sym ðєÞ ð3:11Þ with The GNSS technique retrieves Zenith Tropospheric Delay (ZTD) of the neutral atmosphere using an a priori ZHD (ZHDapriori) and adjusting a ZWD (ZWDadjusted ¼ ZTD – ZHDapriori). ZHDapriori is generally obtained using the formula of Saastamoinen (1972); see also Davis et al. (1985) and Elgered et al. wet depends on the elevation (є) of each satellite and using a wet mapping (1991). Lsym function, e.g. GMF (Boehm et al. 2006a, b). wet Concerning the wet anisotropic contribution (Lasym ), the 1st and the 2nd order wet contributions can be considered by using respectively the wet gradients (Laz ) and the one-way post-fit residuals (Lres), as formulated in Eq. 3.12. Generally, the contribution of residuals is not considered because it can be highly affected by multipath and artefacts in calculations. This study tries to show you the interest of using residuals, as retrieved by GAMIT software (Herring et al. 2010). The 1st order wet wet wet contribution to Laz is estimated with gradient wet components (GNS , GEW ) that is connected to the azimuth (α), and with the use of gradient mapping function  wet  mf az ¼ 1=ð sin єtanє þ CÞ which depends on satellite’s elevation and on a constant C (Chen and Herring 1997). wet wet Lasym ðє; αÞ ¼ Laz ðє; αÞ þ Lres ðє; αÞ with 3 Parts from this section were previously published in Brenot et al. (2013) ð3:12Þ 56 J. Douša et al.  wet  wet wet wet Laz ðє; αÞ ¼ mf az ðє; CÞ: GNS : cos ðαÞ þ GEW : sin ðαÞ The gradient components (GNS, GEW) retrieved by GNSS technique are total. This means there is no distinction between wet and hydrostatic gradients (Chen and Herring 1997; Flores et al. 2000). In geodetic software it is commonly considered that C ¼ 0.0032 (Herring 1992), but for the estimation of the asymmetric wet wet delay (Laz ), C ¼ 0.0031 can be used (Chen and Herring 1997). The wet gradient ! wet G is expressed by the difference of the hydrostatic to the total component. ! wet G ¼  wet GNS wet GEW  ¼  GNS GEW  Ghydrostatic NS Ghydrostatic EW ! ð3:13Þ To obtain the hydrostatic gradient components, a characterisation of the surface pressure field around each GNSS station is required. In that case, the hydrostatic gradient can be established by fitting a plane through the pressure measurements (Champollion et al. 2004; Brenot et al. 2014a, b). From the pressure field near a GNSS site, the spatial variations of the hydrostatic delay per unit of distance (km) in the north-south (Z hydrostatic ) and east-west (Z hydrostatic ) directions can be calculated. EW NS Generally, for a case study, surface pressure measurements around all GNSS stations are not available. Outputs from numerical weather model can be considered. Assuming an exponential law in the hydrostatic refractivity and considering the scale height of the gradients in the hydrostatic delays set to H ¼ 13 km (as suggested by Chen and Herring 1997), the spatial variations of the hydrostatic delay can be converted in hydrostatic gradients (Elósegui et al. 1999; Ruffini et al. 1999; Flores et al. 2000) to obtain wet gradient components (Eq. 3.14).  Gwet NS Gwet EW  ¼  GNS GEW  Z hydrostatic NS H: Z hydrostatic EW ! ð3:14Þ The 2nd order asymmetric contribution from residuals (Lres), can be considered if a station is weakly affected by multipath or these avoided by Multipath Stacking Methods (MPS), as introduced by Elósegui et al. (1995) and Shoji et al. (2004). By using MPS, low elevation measurements can be improved by identifying and avoiding non-tropospheric signature in data. Estimates of slant delays can also be improved by using Phase Centre Variation model (PCV) for the antenna (Shoji et al. 2004). This aspect is not treated in this study, as well as the questioning of artefacts in calculations. This study investigates the contribution of residuals (Lres) to wet delays without any correction of multipath and PCV model. Brenot et al. (2013) have studied in detail the rainfall event of 28–29 June 2005. This paper shows the critical role of GNSS horizontal gradients of the water vapour content to detect small scale structures of the troposphere (i.e. convective cells), and presents a strategy to identify typical water vapour configurations (dry/wet dipole in 3 Advanced GNSS Processing Techniques (Working Group 1) 57 Fig. 3.7 (a) Imaging of the 2D field of ZWD with a classic interpolation (stations are plotted using ! black circles); (b) improvement of this field by GNSS gradients. Wet gradients Gwet are plotted using grey arrows at each GNSS site. BUGG, ERPE, GERA and BRUS stations are plotted. Locations of 9 major Belgium cities (red circles) and meteorological radars (yellow triangles) are also plotted on these 2D maps time and space) and obtain preliminary signs of the initiation of deep convection. The complementary objective of this work is to investigate, step by step, the use of wet gradient to improve 2D field of ZWD and visualise water vapour blobs, and the use of SWD delays to monitor small scale tropospheric structure of convective cells. Total GNSS gradients has been used by Brenot et al. (2013, 2014a, b) to improve the spatial resolution of the 2D field of ZTD (using Hermite interpolation or pseudoobservation defined by gradients). The strategy of Brenot et al. (2013) with additional pseudo-observations has been transferred to the improvement of the 2D field ! of ZWD with wet gradients Gwet (comparison between classic interpolation and improved ZWD field is shown in Fig. 3.7). The grey arrows in Fig. 3.7b, represent ! wet G (expressed in the zenith direction) with amplitudes of about 0.01 m and more for stations GERA and BUGG at 12:00 UTC on 29 June 2005. A relevant way to visualise, in time, isotropic and anisotropic contributions to wet delays, is proposed in Fig. 3.8. Such a graph acts by superimposing horizontal wet gradient to ZWD. Four stations have been selected (ERPE, BUGG, GERA and BRUS), as being several time close or overflown by convective cells. Let’s focus our attention on the time window (10:00 to 14:00 UTC) when these four stations measure high anisotropic contributions (wet gradient higher than 0.01 m). The morning of the 29th June 2005, the wet delay is gradually decreasing for the four stations (as shown in Fig. 3.8), then a sudden increase is observed (starting at 11:00 UTC for ERPE and GERA, and 12:00 UTC for BRUS and BUGG). Figure 3.9 shows the spatial distribution of wet delay all over Belgium at 12:15, 12:30, 13:00 and 13:15 UTC. Low ZWD is observed by BUGG station at 12:15 in Fig. 3.9a. At this moment, the initiation of convection is taking place between BUGG and BRUS (dry/wet dipole). At 12:30 UTC (in Fig. 3.9b), a strong contrast of moistening is still 58 J. Douša et al. ! Fig. 3.8 Visualisation of isotropic (ZWD) and anisotropic (wet gradient, Gwet ) contributions to delays on 29 June 2005 for 4 stations (a) ERPE, (b) BUGG, (c) GERA, and (d) BRUS observed between BRUS and BUGG with significant rainfall. Wet gradients of BUGG and BRUS stations point to area with high precipitation. At 13:00 UTC, the flux of water vapour from north to south is separated by a drier area on the north and east side of BRUS, where the initiation of deep convection is taking place (on the east side of BRUS). The strength of the rapid flux of moistening is shown by a strong increase of the field of ZWD close to BUGG and ERPE and in the north-west of Belgium (see Fig. 3.9c). Wet gradient amplitudes higher than 0.015 m have been observed for several stations during this rainfall event (and especially during this time window), notably at NAMR station (south-east of BRUS in Fig. 3.9d); see Brenot et al. (2013). Looking at the improved 2D field in Fig. 3.9, dry/wet contrast of ZWD field is a good indicator of preliminary signs of deep convection and heavy precipitation. An investigation of the interest of SWD for monitoring small-scale structures (sub-kilometric size) for these four stations is presented in Fig. 3.10 (time window from 10:00 to 14:00 UTC). Using skyplots, the three contributions to SWD (isotrowet pic contribution by ZWD and Lsym , and anisotropic contributions by wet gradients wet and residuals, respectively Laz and Lres) are shown for GNSS signals from two satellites recorded by four stations. Small wet gradients and residuals are observed between 10:00 and 12:00 UTC for these four stations (except a wet structure at 11:00 on the west side of ERPE and 3 Advanced GNSS Processing Techniques (Working Group 1) ! 59 Fig. 3.9 (left) 2D fields of ZWD improved by, Gwet at (a) 12:15, (b) 12:30, (c) 13:00 and (d) 13:15 UTC on 29 June 2005; (right) Radar precipitation at (a) 12:15, (b) 12:30, (c) 13:00 and (d) 13:15 UTC on 29 June 2005. 60 Fig. 3.10 Skyplots for 4 GNSS stations of couples SWD/satellites trajectories on 29 June 2005, (a) ERPE, (b) BUGG, (c) GERA, and (d) BRUS J. Douša et al. 3 Advanced GNSS Processing Techniques (Working Group 1) 61 BUGG, 15 min later, specifically where the initiation of a convective system took place). Between 12:00 and 14:00, several tropospheric structures are identified by the four stations, showing a base of anisotropic contribution provided by the wet gradient and a precise time-space detection by the residuals. Even similarity and correlation with radar precipitation is not obvious, neither straightforward. It can be noticed that at 12:15, a blob of water vapour is detected on the west side of ERPE. For BUGG station, the negative residuals of PRN27 satellite at 12:15, decrease the anisotropy seen by the wet gradients (showing that the highest density of humidity of the wet structure is not located at 60 of elevation). No significant structure is seen by Lres for GERA station. However, a tropospheric structure is detected by residuals of PRN27 BRUS station (in the south – southeast direction for elevations between 40 and 50 ). At 12:30, the wet gradients of BRUS indicate the north, and Lres of PRN27 is negative for elevation between 50 and 70 , showing a drier structure (convective system has moved eastward). A wet structure is seen by Lres of PRN27 for ERPE station (located north-west of BRUS). Between 13:00 and 13:30, PRN27 residuals of ERPE and BRUS stations indicate wet structure in the east direction for an elevation of 70 . BUGG residuals show a wet structure in the south-east. There is a clear identification of sub-kilometric meteorological structure by one-way post fit residuals. Even no correction of multipath and PCV model has been applied (Elósegui et al. 1995; Shoji et al. 2004), a way to justify that the structures detected by residuals have properly a tropospheric origin is to compare Lres on 29 June with the one measured the day before, for which the trajectories of satellites are very similar. On 28 June, the tropospheric activity was moderate. Figure 3.11 shows clearly low residuals on 28 June, justifying that the structures detected on 29 June are due to tropospheric activity (no multipath or artefact effect). If a GNSS station is a good candidate (confirmation of the tropospheric origin of residuals), the 2nd order asymmetric contribution from residuals (Lres), can be considered in meteorological applications (assimilation in numerical weather forecasts or imaging for nowcasting). This study shows a good potential for residuals and SWD for detecting small scale tropospheric structures affecting signal propagation between GNSS satellites and stations. 3.3.4 Indicator of Tropospheric Activity Based on the Disruption of GNSS Signals4 H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] 4 Parts from this section were previously published in Brenot and Warnant (2008) 62 Fig. 3.11 Skyplots for 4 GNSS stations of one-way post fit residuals (Lres) and satellites trajectories on 28–29 June 2005, (a) ERPE, (b) BUGG, (c) GERA, and (d) BRUS J. Douša et al. 3 Advanced GNSS Processing Techniques (Working Group 1) 63 The aim of this subsection is to find a new indicator of small-scale tropospheric activity. Different candidates of tropospheric effects indicators can be considered according to GNSS carrier phase measurements (King et al. 1985; Dong and Bock 1989; Blewitt 1989; Leick 1989; Teunissen et al. 1998). However, ZTD and gradients are not the best candidate being the results of a time and space average. STDs show a good potential to detect small-scale meteorological structure (see Sect. 3.3.5). Nevertheless, double differences (L1, L2) of the ionosphere-free combination of GNSS phase observations can be used as an additional detection of the presence of small-scale structures in the troposphere. Small-scale structures induce disturbances on phase measurements. This study considers stations for which the positions and ij the geometric distancesDAB are precisely known between couples of satellites and ij ij and the ambiguity N AB ground-based receivers. The tropospheric perturbation T AB , IF remains the only unknown parameters in the double difference of phase of the ij ionosphere-free (IF) combination ϕAB , IF (for a simplified mathematical model of phase measurements (Seeber 2003; Leick 2004; Brenot and Warnant 2008). A new ij observable of phase ΦAB , IF can be estimated: ij ij ΦAB , IF ¼ ϕAB, IF f 1 ij f ij ij D ¼ 1 T AB þ N AB , IF c AB c ð3:15Þ c is the speed of electromagnetic waves, ( f1, f2) carrier frequencies of (L1 or L2). ij The ambiguity term (N AB , IF ) has the following expression:  i ij N AB , IF ¼ N A, IF N Bi , IF   N Aj, IF N Bj, IF  ð3:16Þ Ambiguities ( N Ai , IF , N Bi , IF , N Aj, IF and N Bj, IF ) are defined using ionosphere-free ij combination. The phase ambiguity term (N AB , IF ) is a real number with a constant value. ij Figure 3.12 shows an example of tropospheric perturbation T AB and ambiguity ij ij N AB, IF presented by the phase observable ΦAB, IF (called IF Double Difference and expressed in cycles) on 29 June 2005 for BRUS-GILL baseline (4 km) and the couple of satellites (27-08). Without the presence of the troposphere, a constant real ij value of the IF DD should be observed according to the ambiguityN AB , IF (which can be a real number). The error induced by the troposphere on the IF Double Difference ij observable ΦAB , IF time-series is clearly shown in Fig. 3.12 between 12:00 and 13:00 UTC. According to radar imaging important precipitations (higher than 100 mm/h) took place over and north-east of OLLN station at 12:30 UTC the 29th June 2005 (location of strong tropospheric activity). A high content of water vapour and the existence of hydrometeors induces a strong perturbation of atmospheric refractivity (Brenot et al. 2006). Perturbation of refractivity can clearly explain sudden variability of tropospheric error T Ai measured by station A for a signal emitted by a satellite i. The following expression presents 64 J. Douša et al. Fig. 3.12 IF Double Difference of BRUS-GILL baseline the 29th June 2005 (Day Of Year 180) the relation of tropospheric error T Ai (i. e. generally called STD) with neutral atmosphere refractivity (N ): T Ai ¼ 10 Z 6 ð3:17Þ N ds ds is a differential distance according to path travel of signal between satellite ij i and station A. The tropospheric error T AB induces the perturbation of the phase ij observable ΦAB, IF , as defined by Eq. 3.15, representing the double difference of phase of the ionosphere-free combination with a correction of the geometric distances. This tropospheric error has the following expression:  ij T AB ¼ T Ai T Bi   T Aj T Bj  ð3:18Þ Sudden perturbations of tropospheric errors (T Ai ,T Bi ,T Aj and T Bj ) by a small-scale ij ij structures induce direct perturbations of T AB and ΦAB , IF . Considering two epochs of measurements (epoch t0 and epoch t0 + Δt, for example Δt ¼ 5 min), Fig. 3.13 illustrates the direct impact of the occurrence of a small-scale tropospheric structure ij ij on phase measurements (observables ΦAB , IF and T AB ). The resolution of the ambiguities is required. For a selected Day Of Year (DOY), reference satellites are chosen to form double differences (DD) and maximise the time periods. The atmospheric scans by these couples of satellites are sufficient to represent the tropospheric activity. Considering NAMR-OLLN baseline DOY 180 of 2005 (couple of satellites 10–21) and BRUS-BERT baseline DOY 365 of 3 Advanced GNSS Processing Techniques (Working Group 1) 65 Fig. 3.13 Perturbation of TABij induced by a small-scale tropospheric structure for two epochs of measurements (epoch t ¼ t0 and epoch t ¼ t0 + Δt) Fig. 3.14 IF Double Difference of NAMR-OLLN baseline the 29th of June 2005 event (DOY 180) on the left, BRUS-BERT baseline on the right (no meteorological event on 31 December 2005, DOY 365). The fits of DD time-series with polynomial functions of the 3rd order are shown 2006 (couple of satellites 16–19), Fig. 3.14 shows IF Double Differences time-series ij (observable of phase ΦAB , IF of Eq, 3.15) for these two baselines (called IF DD plotted with crosses). The impact on DD depends on the elevation of considered satellites. This is a specificity of the tropospheric activity. In order to display only the influence of small-scale structures on DD time-series, fits of IF DD time-series have been assessed using polynomial functions of the 3rd order (dashed line Fig. 3.14) and 66 J. Douša et al. Fig. 3.15 1 h radar precipitation accumulation starting at (a) 16:00 UTC, (b) 17:00 UTC, (c) 18:00 UTC, and (d) 19:00 UTC on 2006/10/01 Fig. 3.16 Radar precipitation on 29 June 2005 at 14:30 UTC (left), and on 31 December 2006 at 13:15 UTC (right) the biases between IF DD and the respective fits, called IF DD Residuals, can be obtained. The estimation of the bias to the fit removes elevation effects. Then small-scale structures are clearly identified for NAMR-OLLN baseline on the 29th of June 2005 (in Fig. 3.14 between 14:00 and 15:00 UTC). Figure 3.15 shows a time-series of IF DD Index. To obtain this Index of the tropospheric activity, absolute values of IF DD Residuals (in cycles) is converted into centimetres (multiplying by the wave length: 19.029 cm). According to radar imaging of rain rate (in Fig. 3.16 on the left), the tropospheric small-scale activity around Namur (NAMR and OLLN stations) during DOY 180 of 2005 can be easily observed between 14:00 and 15:00 UTC. Note that a strong tropospheric activity was also observed between 12:00 and 13:00 UTC this day (see Figs. 3.14 and 3.15; see also Brenot et al. 2013). No tropospheric activity took place around Brussels (station BRUS) DOY 365 of 2006 at 13:15 UTC (see radar imaging in Fig. 3.16 on the right). 3 Advanced GNSS Processing Techniques (Working Group 1) 67 Fig. 3.17 Imaging of maximal IF DD Index detected (a) 29 June 2005 at 14:30 UTC; (b) 31 December 2006 at 13:15 UTC Considering all the couples of satellites for a selected baseline and all the available phase measurements, the daily tropospheric activity (superposition of all the IF DD Index of the selected satellites-stations couples) can be shown (see Brenot and Warnant 2008). Considering every baseline of the Belgian network, IF DD Index imaging are shown Fig. 3.17. In this imaging, geometric segments (each one corresponding to a baseline) are affected by the maximum IF DD Index estimated at a given moment (at 14:30 UTC on 2005/06/29 and 13:15 UTC on 2006/12/31, for the two examples presented) according to all the couples of satellites considered in our system. Note that the rainfall cell present over OLLN station does not appear with ZTD imaging due to the time and space average. Horizontal delay gradient points a direction where the local anisotropy is maximal. However gradient represents a time and space average which punctually (at 12:30 UTC) do not show exactly the location of small-scale structure (in the north-east direction of OLLN station), as seen by IF DD Index. The IF DD Index imaging (Fig. 3.17a) is clearly sensitive to sudden perturbation of tropospheric activity. That means sensitive to the occurrence of tropospheric small-scale structures which locally affect transmission of signals from GNSS satellites at a given epoch and not for a time and space average measurements. IF DD Index shows strong perturbations of GNSS signal propagation induced by the troposphere around OLLN station between 14:30 and 15:00 UTC (Brenot and Warnant 2008). The presence of water vapour and hydrometeors above OLLN and on the north-east side of this station, affects Double-Difference observations for OLLN-NAMR baseline the 29th June of 2005 (DOY 180). The deep convection process and the thermally driven turbulent mixing that moves air parcels from the lower to the upper atmosphere, shows a vertical extension up to 14 km close to 68 J. Douša et al. BRUS, MECH and BUGG stations, with occurrence of heavy rain and, for some area, hail stones at 12:20 UTC (Brenot et al. 2013). Using a dense network of GNSS stations (e.g. the Belgian dense network with baselines from 5 to 30 km), a relevant monitoring of tropospheric structure can be established with IF DD Index. As an example of severe weather, the month of September and the start of autumn 2006 was exceptionally hot and dry. During the last days of September, sea breeze was finally bringing humidity in the warm low layer. At higher altitude, a strong dynamic was taking place with strong jets maintaining powerful forcing able to generate the sturdiest thunderstorms. The differences of wind direction from low to high levels, so called wind shear, in association with strong flux of water vapour and moistening of low level, led to a critical situation on 1st October 2006. The tropospheric activity took place during all this day, initiating locally deep convection and generating specific conditions for the establishment of supercells that can stay active during few hours. One of these supercells has created a convenient meteorological situation for the formation of a tornado (see photography in Fig. 3.18). This supercell (associated with heavy rainfall from a cumulonimbus cloud, as seen on radar imaging on the south of Brussels in Fig. 3.19), has generated this tornado close to BRUS and LEEU stations. Taking the northeast direction after its creation and avoiding densely built area, the tornado nevertheless reached farms located along its trajectory. Several buildings have been seriously damaged. Figure 3.18b presents the probability of hail (maximum values) during the 1st October 2006 event (DOY 274). The passage of the supercell on the south side of Brussels can be observed. Several other supercells are also shown (wind direction oriented from south-west to north-east). Operational hail detection products are derived from the height of the freezing level and from 45 dBZ echotop values provided by single-polarization C-band weather radar (Delobbe and Holleman 2006). The supercell close to Brussels (BRUS station) and Sint-Pieters-Leeuw (LEEU station) can clearly be observed by daily IF DD Index (Fig. 3.19) applied for a baseline of 8 km. 100 90 80 70 60 50 40 30 20 10 0 Fig. 3.18 Photography of the tornado of Petit-Roeulx-lez-Braine (source: C. De Keyser) at 16:00 UTC on 2006/10/01 (left); daily probability of hail (in %) for this day, with positions and names of 6 GNSS stations; courtesy of Laurent Delobbe (right) 3 Advanced GNSS Processing Techniques (Working Group 1) 69 Fig. 3.19 Radar imaging on the 1st of October 2006 at 16:25 UTC (left); daily IF DD Index for BRUS-LEEU baseline (right) Note that the tropospheric daily IF DD Index detects strong activity DOY 274 for baseline ANTW-KALL. On the other hand, a quiet tropospheric activity is observed for the baseline VOER-TONG. The interest of looking at radar imaging of hail probability is that the possible production of hail requires an important vertical extension and a consequent amount of water vapour, linked with the existence of hydrometeors. The IF DD Index at 16:25 UTC in Fig. 3.19 shows high level of activity with Index more than 10 cm. In Fig. 3.20, 1 h radar precipitation accumulation shows that supercells have traveled from south-west to north-east. Baselines BUGG-NIKL, BREC-HERE, DIES-MOL0, LEEU-NIVL, BRUS-NIVL, BERT-OLLN and BREE-MAAS are plotted on these radar imaging. For these baselines, IF DD Index tropospheric activities are presented in Fig. 3.21. According to radar precipitation accumulation imaging (Fig. 3.20), the description of Fig. 3.21 (left) is the following: at about 16:00 UTC, tropospheric activity for baseline BUGG-NIKL (up to 7 cm) was occurring close to BUGG station (cell C of Fig. 3.20a). This same cell C was located then close to HERE station at about 17:00 UTC and induced an IF DD Index up to 7 cm for BREC-HERE baseline. Around 18:00 UTC, the cell B was approaching over BUGG station (in Fig. 3.20b) and inducing IF DD Index up to 10 cm for BUGG-NIKL baseline. The passage of the cell B can clearly be observed with IF DD Index presented in Fig. 3.21 between 18:00 and 19:30 UTC. Around 19:00 UTC, the cell B (in Fig. 3.20c) was above HERE station (IF DD Index up to 11 cm for BREC-HERE baseline), and at 19:20 UTC (in Fig. 3.20d) above MOL0 station (IF DD Index up to 8 cm for DIES-MOL0 baseline). According to radar precipitation accumulation imaging (in Fig. 3.20), the description of Fig. 3.21 (right) is the following: at 16:20 UTC, tropospheric activity has taken place around LEEU and NIVL stations (IF DD Index up to 8 cm for baseline LEEU-NIVL). Between 16:30 and 17:00 UTC the supercell A has moved from south-west to east of Brussels (successively IF DD Index of 6 cm for BRUS-NIVL and BERT-OLLN baselines). Around 18:40 UTC, the cell A was close to MAAS station and induced IF DD Index up to 7 cm for BREE-MAAS baseline. The passage of the supercell A is clearly shown from 16:00 to 19:00 UTC on IF DD Index. 70 J. Douša et al. Fig. 3.20 1 h radar precipitation accumulation starting at (a) 16:00 UTC, (b) 17:00 UTC, (c) 18:00 UTC, and (d) 19:00 UTC on 2006/10/01 The contribution of hydrometeors in association with water vapour bubble to strong IF DD Index of tropospheric activity is indisputable. The blobs of water vapour surround rainfall cells with a high vertical extension (up to 10 km for cells B and C and up to 11 for cell A, as estimated by the maximum radar reflectivity from the weather radar of Wideumont). The implementation of this index started in the frame of the GALOCAD/ESA project (2006–2008). The aim was to find relevant tropospheric and ionospheric indicators to warn the impact on NRT positioning solutions, i.e. effect in RTK-architecture for a dense network of stations (Brenot and Warnant 2008; Warnant et al. 2008, Wautelet et al. 2008; Brenot et al. 2014a, b). To summarise, this study presents GNSS indicators of meteorological activity that allow the detection of small-scale structures in the neutral atmosphere. The scope is to present a new NRT index of meteorological activity based on doubledifference of the ionosphere-free combination (so called IF DD index). Contrary to 3 Advanced GNSS Processing Techniques (Working Group 1) 71 Fig. 3.21 Daily IF DD index (2006/10/01) for BUGG-NIKL, BREC-HERE and DIES-MOL0 baselines (left), and for LEEU-NIVL, BRUS-NIVL, BERT-OLLN and BREE-MAAS baselines (right) ZTD imaging and horizontal delay gradients measurements from geodetic software (result of a mean time and space solution), the IF DD Index imaging is clearly sensitive to sudden disturbances of tropospheric activity. That means sensitive to the occurrence of tropospheric small-scale structures which locally affect couples of satellites emitted signals considered in NRT applications (i.e. GNSS meteorology or positioning) at a given epoch. The use of the IF DD Index can be planed operationally in NRT meteorological or geodetic system using dense networks, being useful for forecasters and nowcasting. The next step of this work can be to improve the time and space imaging of the IF DD Index using multi-GNSS satellites, in collaboration with forecasters. The contribution of hydrometeors to IF DD Index need to be investigated studying correlation with radar reflectivity. The use of the geometry-free combination (GF) can also be used to define an indicator of ionospheric activity (GF DD Index). This work based on DD difference can also be transferred to of L1, L2 and IF combination to obtain IF Index of the tropospheric activity between a satellite and a ground-based station. As a first investigation the flash-flood event in the Gard region, on 8–9 September 2002, has been tested (Brenot et al. 2006), showing in Fig. 3.22 strong tropospheric activity when convective cells were located close to CHRN station. 3.3.5 Validation of Slant Tropospheric Delays5 M. Kačmařík Institute of Geoinformatics, VŠB Technical University of Ostrava, Ostrava, Czech Republic e-mail: [email protected] 5 Parts from this section were previously published in Kačmařík et al. (2017) 72 J. Douša et al. Fig. 3.22 Time-series of IF combination (GPS satellite PRN06) and its fit for CHRN station on 8–9 September 2002 (left); IF Index for all the signals emitted from satellites and recorded by CHRN on 8–9 September 2002 (right) J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] G. Dick GFZ German Research Centre for Geosciences, Helmholtz Centre Potsdam, Potsdam, Germany e-mail: [email protected] F. Zus GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: zusfl[email protected] H. Brenot Royal Belgian Institute for Space Aeronomy, Uccle, Belgium e-mail: [email protected] G. Möller Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] E. Pottiaux Royal Observatory of Belgium, Brussels, Belgium e-mail: [email protected] P. Hordyniec Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] 3 Advanced GNSS Processing Techniques (Working Group 1) 73 J. Kaplon Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] L. Morel École Supérieure des Géomètres et Topographes, Le Mans, France e-mail: [email protected] Tropospheric STD represents the total delay that undergoes the GNSS radio-signal due to the neutral atmosphere along the path from a satellite to a ground receiver antenna. It is evident that STDs can provide much more information about the distribution of water vapor in the troposphere than classical ZTDs. With a continuous development of NWM forecasting and nowcasting tools a demand for highquality humidity observations with high spatial and high temporal resolutions is growing at side of a meteorological community. On the other hand, despite a set of studies which dealt with this topic there is still not a uniform consensus on how to reconstruct STDs from GNSS processing results. Therefore, we decided to realize an extensive inter-technique validation of STDs using data from Benchmark data set and try to answer at least some of the opened questions. A summary of obtained results is given here, however, we refer the reader to the publication of Kačmařík et al. (2017) for a much more detailed presentation. 3.3.5.1 Description of STD Validation Study From the complete Benchmark data set, we selected a subset of 10 GNSS reference stations situated at six different locations (Table 3.4). It also includes collocated (dual) GNSS stations playing an important role in the validation since they track GNSS satellites with the same azimuth and elevation angles, so that they should deliver the same or very similar tropospheric parameters used for STD reconstructions. Seven institutions delivered their STD solutions for this validation study, namely Ecole Supérieure des Géomètres et Topographes (ESGT CNAM), Geodetic Observatory Pecný (GOP, RIGTC), Helmholtz Centre Potsdam – German Research Centre for Geosciences (GFZ), Royal Observatory of Belgium (ROB), VŠB-Technical University of Ostrava (TUO), Vienna University of Technology (TUW), and Wroclaw University of Environmental and Life Sciences (WUELS). Principal information about individual solutions are given in Table 3.5. In total, we validated eleven solutions computed with five different GNSS processing software. Considering all available GNSS solutions, only GOP used a stochastic modelling approach to estimate all parameters. Additionally, GOP provided two solutions: (1) GOP_F using Kalman filter (forward filter only), i.e. capable of providing ZTD, tropospheric gradients and STDs in real time, and (2) GOP_S applying the backward smoothing algorithm (Václavovic and Douša 2015) on top of the Kalman filter. The latter improves the quality of all estimated parameters during 74 Table 3.4 Characteristics of 10 GNSS reference stations Name GOPE KIBG LDB0 LDB2 POTM POTS SAAL WTZR WTZS WTZZ Latitude [ ] 49.914 47.449 52.210 52.209 52.379 52.379 47.426 49.144 49.145 49.144 Longitude [ ] 14.786 12.309 14.118 14.121 13.066 13.066 12.832 12.879 12.895 12.879 Height [m] 593 877 160 160 145 144 796 666 663 666 Network IGS, EPN Dual station IGS, EPN LDB2 LDB0 POTS POTM IGS, EPN IGS IGS WTZS, WTZZ WTZR, WTZZ WTZR, WTZS Receiver TPS NET-G3 TPS GB-1000 JAVAD TRE_G2T JPS LEGACY JAVAD TRE_G3TH JAVAD TRE_G3TH DELTA TPS GB-1000 LEICA GRX1200 + GNSS SEPT POLARX2 JAVAD TRE_G3TH DELTA Antenna TPSCR.G3 TPSH TPSCR3_GGD CONE JAV_GRANT-G3T NONE LEIAR25.R4 LEIT JAV_GRANT-G3T NONE JAV_RINGANT_G3T NONE TPSCR3_GGD CONE LEIAR25.R3 LEIT LEIAR25.R3 LEIT LEIAR25.R3 LEIT J. Douša et al. Solution Name CNAM GFZ GOP_F GOP_S ROB_G ROB_V TUO_R TUO_G TUW_3 TUW_7 WUE Institution ESGT CNAM GFZ Potsdam GO Pecný GO Pecný ROB ROB TU Ostrava TU Ostrava TU Vienna TU Vienna WUELS Strategy DD PPP PPP PPP DD DD DD DD PPP PPP PPP Software GAMIT EPOS 8 G-Nut/Tefnut G-Nut/Tefnut BSW52 BSW52 BSW52 BSW52 NAPEOS NAPEOS BSW52 GNSS GPS GPS GPS GPS GPS + GLO GPS + GLO GPS + GLO GPS GPS + GLO GPS + GLO GPS Elev. cut-off 3 7 7 7 3 3 3 3 3 7 3 Mapping function VMF1 GMF GMF GMF GMF VMF1 VMF1 VMF1 GMF GMF VMF1 ZTD/gradients interval 1 h/1 h 15 min/1 h 2.5 min/2.5 min 2.5 min/2.5 min 15 min/1 h 15 min/1 h 1 h/3 h 1 h/3 h 30 min/1 h 30 min/1 h 2.5 min/1 h Post-fit residuals NO YES YES YES YES YES NO NO YES YES YES 3 Advanced GNSS Processing Techniques (Working Group 1) Table 3.5 Information about GNSS-based STD solutions used in the validation 75 76 J. Douša et al. the batch processing interval and eliminates effects such as the PPP convergence or re-convergence. Some institutions delivered also two STD solutions which differ in a single processing setting. The aim was to evaluate their impact on STDs: a) TUO_G and TUO_R exploit GPS-only and GPS + GLONASS observations respectively, b) TUW_3 and TUW_7 apply an elevation cut-off angle of 3 and 7 degrees respectively, and c) ROB_G and ROB_V use the GMF and VMF1 mapping functions respectively. Additionally, ROB solutions are the only ones based on the processing of double-difference (DD) observations and providing ZD carrier-phase post-fit residuals converted from the original DD residuals using the technique described in Alber et al. (2000). For an independent validation of STDs from GNSS processing we used STDs derived from NWM via ray-tracing and from observations of Water Vapor Radiometer (WVR). In case of NWM derived STDs, four institutions delivered their solutions based on three different NWMs: ALADIN-CZ (4.7 km resolution, limited-area hydrostatic model, operational analysis in 6-h interval with forecasts for 0, 1, 2, 3, 5, 6 h, http://www.umr-cnrm.fr/aladin/), ERA-Interim (1 horizontal resolution, 6-h reanalysis), and NCEP-GFS (1 horizontal resolution, 6-h operational analysis, https://www.ncdc.noaa.gov/data-access/model-data/model-datasets/ global-forcast-system-gfs). None of these NWM assimilates data from ground GNSS stations. For more information about the models, see Douša et al. (2016) and specifically Trojáková (2016) for the ALADIN-CZ model and Dee et al. (2011) for the ERA-Interim reanalysis. First, STD solutions using the ERA-Interim and NCEP-GFS models were delivered by GFZ Potsdam using acronym ERA/GFZ and GFS/GFZ, respectively. Two STD solutions were then delivered for the ALADINCZ model: (a) ALA/BIRA generated at Royal Belgian Institute for Space Aeronomy (BIRA), and (b) ALA/WUELS delivered by Wroclaw University of Environmental and Life Sciences. For a description of these solutions we refer to Kačmařík et al. (2017). In case of WVR we used the instrument operated at GFZ Potsdam. The instrument is situated on the same roof as the GNSS reference stations POTM and POTS. The WVR is switching between ‘zenith mode’ when it is measuring IWV and ‘slant mode’ when it is tracking GPS satellites using an in-built GPS receiver. In the latter case, Slant Integrated Water Vapour (SIWV) values are delivered for the direction of satellites. Our study focuses on the comparison of STDs, not SIWV. It was thus necessary to convert the WVR SIWV into STDs. The used conversion of WVR SIWVs to STDs aimed at minimum distorting the accuracy of original WVR observations and is described in Kačmařík et al. (2017). 3.3.5.2 Introduction to STD Estimation from GNSS Observations The tropospheric STD is being reconstructed from tropospheric parameters valid in zenith direction which are being estimated during GNSS observation processing. It can be expressed by Eq. 3.19, where ZHD and ZWD represent Zenith Hydrostatic 3 Advanced GNSS Processing Techniques (Working Group 1) 77 and Wet delay, respectively, G horizontal tropospheric gradient and ele is the elevation angle and azi is the azimuth angle of observation. STDðele; aziÞ ¼ ZHD ∙ mf h ðeleÞ þ ZWD ∙ mf w ðeleÞ þ Gðele; aziÞ þ RES MPT ð3:19Þ The elevation angle dependency of STD is described by the mapping functions, separately for the hydrostatic (mfh) and the wet (mfw) components. Additionally, post-fit residuals (RES) may contain un-modelled tropospheric effects not covered by the estimated tropospheric parameters. Obviously, residuals contain also other un-modelled effects such as multipath (MPT), however, potentially including also errors from antenna phase centre variations or systematics in satellite clocks. For eliminating such systematic effects, cleaning of post-fit residuals was applied by generating elevation/azimuth-dependent correction maps as described by Shoji et al. (2004). We thus computed mean values of post-fit residuals in 1  1 degree bins using the whole benchmark period for each solution and station. Computed means were then subtracted from the original post-fit residuals to generate solutions using cleaned residuals. Therefore, whenever zero-differenced (ZD) post-fit residuals were available for any solution delivered to the validation, three variants of the solution were used: (1) solution without residuals (nonRES), (2) solution with raw residuals (rawRES), and (3) solution with cleaned residuals (clnRES). Example maps obtained with gradient estimation, polar 1  1 degree bins for multipath determination and 2 sigma outlier rejection threshold are presented on Fig. 3.23 for LDB0. Kapłon et al. (2017) later realized a set of tests to evaluate the impact of strategy of STDs calculation on STD differences obtained from GNSS and Multipath map for station: LDB0 0 40 20 –50 –150 –100 –50 0 50 100 150 Map of std. dev. of multipath: LDB0 50 80 40 60 30 40 20 20 10 0 –150 –100 –50 0 50 Azimuth [deg] 100 Fig. 3.23 Example equal area maps of multipath effect for LDB0 station 150 StdDev. [mm] Elevation angle [deg] 60 Multipath [mm] 50 80 78 J. Douša et al. raytracing through the GFS NWM model including testing the impact of method of multipath effect calculation: polar degree bins (1  1, 2  2 and 5  5 degree) or equal area bins (Huisman et al. 2009) of 1, 2 or 5 degree height, level of outliers reduction (2 or 3 sigma) on STDs. Summarizing results of these tests, all 24 variants of multipath maps provided very similar results with a slight edge for variants using the smallest 1  1 degree bins. 3.3.5.3 Methodology of STD Comparisons Since NWM outputs are restricted to the time resolution of their predictions (typically 1, 3 or 6 h) and since WVR is able to track only one satellite at one moment, all three sources provide different numbers of STDs per day. Therefore, three different comparisons were realized: (1) results for GNSS versus GNSS comparisons, (2) results for GNSS versus NWM comparisons, and (3) results for GNSS versus WVR comparisons. All the presented results were obtained over the whole benchmark period of 56 days. No outlier detection and removal procedure was applied during the statistics computation within the study. Two variants of the comparisons are presented: ‘ZENITH’ and ‘SLANT’. ‘ZENITH’ stands for original STDs mapped back to zenith direction using 1/sin(ele) formula. Such mapping aimed at normalizing STD differences for their evaluation in a single unit. The ‘SLANT’ type of comparison denotes an evaluation of STDs at their actual elevation angles. To be more specific, STDs were grouped into individual elevation bins of 5 degrees, i.e. for example all STDs with an elevation angle between 10 and 15 degrees were evaluated as a single unit. The cut-off angle of 7 degrees was used in all GNSS versus GNSS and GNSS versus NWM comparisons. In GNSS versus WVR comparisons 15 degrees cut-off was applied to exclude problematic WVR observations from low elevation angles. 3.3.5.4 GNSS Versus GNSS: Evaluation of All GNSS Solutions Versus the Reference GNSS Solution Individual GNSS solutions were first compared to the GFZ solution in the zenith direction (ZENITH). Figure 3.24 shows all the solutions using STDs calculated from the estimated ZTD and horizontal gradient parameters, i.e. without adding post-fit residuals. Adding raw or clean residuals, applied consistently to both compared and reference solutions, provided very similar graphs (not displayed). Colours in the Figure indicate the processing software used in individual solutions. Medians of all solutions (dotted lines in each bin) are displayed for each station in order to highlight differences among the stations. These were observed mainly as biases ranging from 3.6 mm to 0.6 mm. The better agreement between GOP and GFZ solutions could be attributed to a similar strategy of both solutions compared to others. It is particularly visible for LDB0 and POTM stations where median values over all solutions differ by 2.3 mm and 3.6 mm, respectively. The reason for the 3 Advanced GNSS Processing Techniques (Working Group 1) 79 Fig. 3.24 Comparison of individual GNSS STD solutions against GFZ solution, all without using residuals (nonRES) and projected in the zenith direction: bias (left) and standard deviation (right). The median value of all solutions at each station is represented by the dotted blue line in each bin divergent behaviour at the two stations has not been identified although site metadata were cross-checked carefully. A significant difference can also be noticed for TUW_3 and TUW_7 at the station KIBG where these solutions used individual antenna calibration files while all others solution used type mean calibration (Schmid et al. 2016). Plots with standard deviations show agreements within 3–5 mm among all the stations and all solutions. The only exception is the GOP_F solution representing a simulated real-time analysis applying only a Kalman filter (not backward smoothing) and providing results by a factor of 2 worse compared to the others in terms of precision. 3.3.5.5 GNSS Versus GNSS: Evaluation in the Slant Direction Figure 3.25 provides an evaluation of the STDs at their original elevation angles for the station POTS. Four individual panels show bias (top left), normalized bias (NBIAS, top right), standard deviation (bottom left), and normalized standard deviation (NSDEV, bottom right). Normalized bias and normalized standard deviation were computed to see the dependence of relative errors in STDs at different elevations. For its computation, absolute differences of STDs from two solutions were divided by the STD values from the reference solution. We found that the agreement among individual solutions compared to the GFZ STDs is rather stable above the elevation angle of 30 degrees. Corresponding biases of individual elevation bins are within 4 mm and standard deviations are slowly increasing up to 10 mm at 30 degrees. With elevation angles decreasing below 30 degrees the biases slightly increase for some solutions. Normalized standard deviation remains almost constant over all elevation angles indicating a very consistent relative performance of STDs among all the solutions. A similar behaviour is present at all stations although the absolute values can be higher for some stations or solutions, namely GOP_F for LDB0 and WTZZ with standard deviations reaching up to 72 mm. 80 J. Douša et al. Fig. 3.25 Comparison of individual GNSS STD solutions against GFZ STD solution at station POTS, in slant directions 3.3.5.6 GNSS Versus NWM: Summary of Results A summary of the GNSS versus NWM validation is presented in Table 3.6. For each reference station a median of bias and a median of standard deviation in the zenith direction between all GNSS solutions and a particular NWM-based solution are given. If we consider ALA/BIRA and ERA/GFZ only, without the two mountainous stations KIBG and SAAL, absolute biases between NWM and GNSS solutions stay mostly below 3 mm, which represents a very good agreement between these independent sources used for retrieving slant delays. Standard deviations generally range from 8 mm to 12 mm in the zenith projection, with an exception of ALA/WUELS showing lower precision by a factor of 2.5. 3.3.5.7 GNSS Versus WVR: Summary of Results A bias of about 5.5 mm in the zenith direction was found between WVR and GNSS solutions at station POTS while the bias at station POTM was around 10 mm. The difference between stations POTM and POTS are probably related to issues with GNSS data processing at POTM. The bias between POTS and WVR roughly corresponds to 1 kg/m2 of IWV, what can be addressed as the achievable accuracy 3 Advanced GNSS Processing Techniques (Working Group 1) 81 Table 3.6 Medians of bias and standard deviation values of differences between all GNSS solutions and a particular NWM-based solution at each reference station, expressed in the zenith direction Station GOPE KIBG LDB0 LDB2 POTM POTS SAAL WTZR WTZS WTZZ Bias (mm) ERA/ ALA/ GFZ BIRA 0.3 3.3 19.3 4.9 2.0 0.7 1.6 0.9 3.4 6.3 1.7 1.4 19.4 7.8 4.8 1.5 3.5 0.9 2.1 0.9 GFS/ GFZ 8.6 9.6 5.5 6.1 12.5 7.6 11.7 4.9 4.2 6.0 ALA/ WUELS 11.5 22.5 10.6 15.1 18.9 12.5 24.3 10.2 10.8 11.6 Standard deviation (mm) GFS/ ERA/ ALA/ GFZ GFZ BIRA 8.3 10.3 7.1 11.6 17.8 11.0 9.9 10.3 8.5 9.1 10.1 8.6 8.0 10.6 9.4 7.7 10.3 9.2 12.7 17.9 11.8 11.0 11.8 8.5 11.4 12.3 8.7 11.3 12.0 8.9 ALA/ WUELS 22.4 26.7 26.2 25.4 26.2 25.8 22.9 23.1 23.7 23.7 of any technique, however, WVR accuracy is more dependent on a proper instrument calibration. Values of standard deviation around 12 mm in the zenith direction, were higher than those observed in GNSS versus GNSS comparisons and slightly higher than from GNSS versus NWM comparisons. 3.3.5.8 Results at Collocated Stations For the GNSS versus NWM and GNSS versus WVR comparisons at individual stations slightly higher values of standard deviations were always found for GNSS solutions applying raw or cleaned residuals in contrast to versions of solutions without any residuals (Kačmařík et al. 2017). However, since two erroneous techniques were always confronted to each other without knowing the true reference, these results do not tell anything about potential of post-fit residuals. For these reasons, we assessed all GNSS solutions at the collocated (dual) stations because for them we are able to provide troposphere-free differences of STDs to evaluate noise of GNSS STD retrievals. Dual stations were available in the benchmark campaign at three different locations in Germany. The first two sites collocate twin GNSS reference stations (LDB0 + LDB2 and POTM+POTS), the third location collocate three individual reference stations (WTZR+WTZS+WTZZ). During normal weather conditions, the tropospheric variation is reasonably smooth, meaning it can be well represented by GNSS STDs reconstructed only from ZTDs and horizontal gradients. However, during high temporal or spatial variabilities in the troposphere, post-fit residuals certainly contain tropospheric signals which are not modelled. If they surpass the observation noise and other residual errors from GNSS models, cleaned residuals should be considered in the GNSS STD model as described in Eq. 3.19. In order to initially address optimal STD modelling under different weather conditions within the benchmark, we tried to 82 J. Douša et al. identify days with a high variability in the troposphere. Daily standard deviations of cleaned post-fit residuals were computed individually for each day of the benchmark, for every station and GNSS solution for 1-degree elevation angle bins. We studied their daily variations considering the GNSS model applied. If cleaned postfit residuals consist of the noise of observations only, the variation in time should be negligible. However, the days showing significantly higher values, correlated at collocated stations, indicated highly variable tropospheric conditions. Three such days were identified at LDB0, LDB2, POTM and POTS stations (May 31, June 20, June 23) and 2 days at WTZR and WTZS stations (June 19, June 20). They all very well correspond to the days initiating heavy precipitations in the domain, Douša et al. (2016). Typical differences between raw and clean residuals are displayed in Fig. 3.26 for all elevations during the normal day (June 19, DOY 170) and the following day with high variability in the troposphere (DOY 171, June 20) for POTM and POTS stations using GFZ solution. Obviously, the variability of clean residuals (black dots) and their 2-sigma envelops are higher by a factor of two for the day of year 171 compared to 170. The variability is clearly visible over all elevations, but the increase is slightly higher at low elevations. The plot demonstrates the different quality of GNSS observations, particularly related to a multipath effect displayed by 2-sigma envelop (green curves). The multipath level is much lower for station POTS which is using a choke ring antenna compared to station POTM which is not using a choke ring. The similar situation was found for stations LDB0 and LDB2. Variability of 2-sigma envelopes of clean residuals (red curves) indicates a higher sensitivity of clean residuals to the weather conditions compared to station selection and observation quality, thus suggesting a significant Fig. 3.26 Comparison of individual GNSS STD solutions against GFZ STD solution at station POTS, in slant directions. Elevation-dependent variability of clean residuals (black dots) and their 2-sigma envelops (red curves) are showed for June 19 (DOY 171) and June 20 (DOY 170) and stations POTS and POTM. Additionally, plots display 2-sigma envelopes for raw residuals (blue curves) and multipath (green curves) 3 Advanced GNSS Processing Techniques (Working Group 1) 83 contribution from the troposphere to the cleaned residuals. In the same context, raw residuals show much higher sensitivity to the observation quality compared to different weather conditions, which is particularly true in case of LDB0 and LDB2 stations. In a next step elevation-dependent differences of STDs from all three versions (without residuals, with raw residuals, with clean residuals) were analysed for days with high and low variability of residuals. We noticed following: (a) STD differences are more or less similar for both days, i.e. no significantly different between days with normal and high variations in the troposphere. It suggests that increased residuals contain strong contributions from the tropospheric effect that could not have been assimilated into ZTDs and tropospheric horizontal gradients (b) STD differences using raw residuals were always the largest ones and they varied with the elevation angle, (c) relative performance of differences from STDs with clean residuals and without residuals for different days remained similar. Uncertainties of the simplified STDs at low elevations surpassed additional uncertainties due to applying clean residuals. According to the magnitude of clean residuals at low elevations (Fig. 3.26), the small uncertainties from calculated differences indicated the presence of tropospheric signals in the residuals at low elevations, roughly below 30 degrees. It seemed to be almost independent from the weather conditions and is supposed to represent mainly unmodelled horizontal asymmetry in the troposphere. Figure 3.27 displays results for comparisons of individual collocated stations in slant directions calculated from all days of the benchmark. The same statistics and plots (not displayed) were prepared also for days identified with ‘severe’ weather conditions, but only minor differences were observed. Strong variations are observed mainly in normalized biases over all elevation angles for the solutions using raw post-fit residuals (rawRES) regardless weather conditions. These are clearly related to local effects such as multipath or modelling instrumented related effects (phase centre offsets and variations) and disappear after using the cleaned residuals (clnRES). The standard deviations and normalized standard deviations at all stations are clearly the lowest for variants without using post-fit residuals (nonRES), slightly higher using cleaned residuals, and significantly higher when using raw residuals, i.e. corresponding to above performed inter-technique validations. 3.3.5.9 Future Work Three institutions (GFZ, GOP, ROB) delivered GNSS STD solutions not only for the ten GNSS reference stations but for the whole GNSS network within the Benchmark data set. Our future study will therefore focus on (a) larger comparison within the network of stations, (b) an evaluation of azimuthal dependency of post-fit 84 J. Douša et al. Fig. 3.27 Comparison of GNSS STDs at dual stations computed over whole benchmark period from individual GNSS solutions in the slant direction for dual stations from left to right: LDB0LDB2, POTM-POTS, WTZR-WTZS. Statistical parameters from top to bottom: bias, normalized bias, standard deviation, normalized standard deviation residuals under severe weather conditions and (c) an evaluation of GNSS STDs estimated from real-time and post-processed solutions using a stochastic approach. 3.3.6 Information Content in Post-fit Residuals, PPP vs DD Approach6 G. Möller Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] Based on real GNSS measurements, the tropospheric signal in post-fit residuals is difficult to assess since it is superimposed by a series of other unmodelled effects like 6 Parts from this section were previously published in Möller (2017) 3 Advanced GNSS Processing Techniques (Working Group 1) 85 observation noise, multipath, satellite clock or orbit errors. Nevertheless, observation stacking methods as applied in (Kačmařík et al. 2017; Möller 2017) allow for the reduction of common parts like multipath or clock errors, but only when longer time periods or larger GNSS networks are processed. Hence, within the COST action an initiative was carried out which addresses the general tropospheric signal in GNSS post-fit residuals, with focus on differenced data processing. Theoretically, the precise point positioning (PPP) and the double difference approach (DD) are equivalent with respect to redundancy and with respect to the estimates, in case a correct stochastic model is introduced. Practically, the DD approach has some advantages in data processing since the satellite and receiver clock errors and therewith the hardware biases cancel out, which allows the fixing of integer ambiguities. Further, also the pre-processing is less critical since the receiver clock error has to be known only with μs-accuracy. Unfortunately the greatest strength of double-difference processing, the elimination of common effects, is also a shortcoming at the same time in small networks (< 500 km). In such networks, tropospheric parameters cannot be estimated in an absolute sense but rather with respect to a reference station. Therefore, reference values (station coordinates and ZTD) have to be introduced, at least for one station, and constrained to their given values. Then the tropospheric parameters can be estimated like in PPP processing, except for the reference station. In order to analyse satellite or station specific effects in double-difference residuals (DDR), the residuals have to be converted into zero-difference residuals (ZDR), also known as pseudoZDR since certain conditions have to be applied for the reconstruction. (Alber et al. 2000) suggested a two-step approach in which the DDR vector is converted into a pseudo-ZDR vector, assuming zero-mean conditions. In order to analyse the applicability of this approach and in general of the tropospheric signal in DDR, two sets of dual-frequency GPS observations were simulated for 12 stations in Austria. Both sets differ only with respect to the applied troposphere model. While no troposphere model was applied to the first set of observations, ZTDs and East-West gradients were simulated for the second set. The observations of all 12 stations were processed in PPP and double-difference approach. If both, ZTD and gradients are estimated, the simulated STD could be recovered with sub-mm accuracy and the post-fit residuals became negligible. If only the ZTD is estimated, it is expected that an anisotropic delay remains in the post-fit residuals. It turned out that in case of PPP methods the anisotropic delay, except for a small offset which was absorbed by the ambiguity parameter, could be recovered from the PPP post-fit residuals but unfortunately not from the DDR. If only a single baseline is processed, the DDR and also the reconstructed pseudo-ZDR are almost zero since the anisotropic effects were differenced out in data processing. In best case, the resulting STD bias and standard deviation was 1 mm +/ 37 mm. This was obtained by fixing the ZTD and by taking all possible baselines between the 12 stations into account for the reconstruction of pseudo-ZDR using the method proposed by (Alber et al. 2000). However, a comparable result (0 mm +/ 38 mm) is obtained if no residuals are added to the isotropic STD. 86 J. Douša et al. In practice, the ZTD is not known and therewith cannot be fixed to its given value. Thus an additional solution was created whereby ZTD but no gradient parameters were estimated. This results in an increase of bias and standard deviation (84 mm +/ 105 mm). This example underscores the importance to estimate gradient parameters in addition to ZTDs. An unmodelled east-west gradient of 2 mm introduced a ZTD error of 35 mm +/ 13 mm. For more details the reader is referred to (Möller 2017). It becomes obvious that the applied reconstruction method proposed by (Alber et al. 2000) is less suited for the reconstruction of pseudo-ZDR in small networks. The reconstructed values are mostly too small. In addition, jumps appear in the time series every time a satellite rises or sets. The magnitude of the jumps can be reduced by downweighting of low elevation satellites; however, the reconstruction process cannot be significantly improved therewith. In consequence, for analysis of satellite or station specific effects in post-fit residuals we recommend undifferenced GNSS data processing strategies, especially in small GNSS networks. 3.3.7 Tropospheric Parameters from Multi-GNSS Analysis7 P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] Nowadays, multi-GNSS offers new satellites and signals which are expected to strengthen all estimated parameters, in particular the ZTD and horizontal linear tropospheric gradients, or to densify slant tropospheric delays for monitoring the troposphere asymmetry at individual GNSS sites. Currently, data from the US NAVSTAR Global Positioning System (GPS) and the Russian GLONASS constellation are commonly used to produce different products within scientific services. Essential models for these two systems has been already established, and their mutual combination provides better precision then the processing from any standalone system. Besides others, limitations for the use of GNSS data from other global systems, the European Galileo and the Chinese BeiDou, persist mainly in (1) incompleteness of the constellations, (2) lack of precise models and calibrations for new signals, receiver and satellite instrumentations, and (3) lack of precise orbit and clock products supporting the ultra-fast processing mode. The situation will 7 Parts from this section were previously published in Douša et al. (2018a). 3 Advanced GNSS Processing Techniques (Working Group 1) 87 change soon as both global systems will become operational in next years and the IGS Multi-GNSS Experiment (MGEX, http://mgex.igs.org, Montenbruck et al. 2017) is continuously filling the gaps in data, metadata, models, formats, standards and products for an optimal exploitation of all global satellite constellations and their regional augmentations. 3.3.7.1 Evaluation of Results from Collocated GNSS Stations The impact of using multi-constellation data on the tropospheric parameter estimation can be optimally assessed using closely collocated GNSS stations, e.g. within few meters. Although different instrumentation-specific effects, such as phase centre modelling and the quality of a receiver tracking, can affect analyses at both stations, the station should principally observe the same tropospheric delays. For the purpose of our evaluation, we selected two IGS station pairs, ZIM2-ZIMJ (Zimmerwald, Switzerland) and MAT1-MATE (Matera, Italy), all collecting data from GPS, GLONASS and Galileo systems within 10 m and 2.5 m in horizontal and vertical distances, respectively. We used the GFZ MGEX (GBM) product (Deng et al. 2017) and CODE MGEX (COM) product (Prange et al. 2017) as multi-GNSS reference solutions for the comparison. We assessed not only the impact of using more constellations but also the impact of different strategies for the parameter estimation. The first approach is the Kalman filter usable mainly for real-time, and the second strategy is the backward smoothing designed for improving the precision of parameters in post-processing (Václavovic and Douša 2015). Solutions improved by introducing multi-constellations and the backward smoothing are demonstrated in the time series of ZTD and horizontal gradient differences obtained from the two collocated stations, ZIM2 and ZIMJ, Fig. 3.28. Results of the single /multi-constellations are visualized by different colours: (1) standalone GPS in red, (2) GPS + GLO in green, and (3) GPS + GLO + GAL in blue. A positive effect is visible for all parameters, and is similar using both the Kalman filter and the backward smoothing, i.e. for both real-time and postprocessing strategy. Scatters of multi-constellation solutions are smaller compared to the standalone GPS solution. More significant effect is visible for the smoothed gradient parameters. Theoretically, zero differences are expected for the collocated stations with the same antenna height. However, a vertical difference between ZIM2 and ZIMJ is about 2 m which can cause about 0.5 mm difference in ZTDs when considering the pressure decreases approximately by 11.3 Pa/m near the geoid and the 100 Pa difference in the atmospheric pressure causes a 2.27 mm difference in ZHD (Saastamoinen 1972). As we observe a ZTD difference about 2–3 mm, it can be still attributed to remaining station-specific systematic errors, e.g. such as phase centre offset and variation models. Numerical statistics (biases and standard deviations) characterize an impact of single- and multi-constellation solutions on the estimated ZTDs and horizontal linear gradients when using the Kalman filter, Table 3.7. It should be noted, that GLONASS (R) and Galileo (E) observations were down-weighted by a factor of 88 Fig. 3.28 Time series of ZTD (top), north gradient (middle) and east gradient (bottom) differences at Zimmerwald dual-station when using the Kalman filter (left) and the backward smoothing (right) J. Douša et al. 3 Advanced GNSS Processing Techniques (Working Group 1) 89 Table 3.7 Statistics (BIAS  SDEV) for Kalman filter using GPS (G), GLONASS (R) and Galileo (E) Station pair ZIM2ZIMJ ZIM2ZIMJ ZIM2ZIMJ MAT1MATE MAT1MATE MAT1MATE GNSS G BIAS  SDEV ZTD [mm] +2.8  1.4 BIAS  SDEV N-GRD [mm] +0.08  0.17 GR +2.4  1.3 +0.02  0.14 0.02  0.12 GRE +2.0  1.3 +0.03  0.14 0.04  0.13 0.5  2.4 0.03  0.18 +0.18  0.25 GR +0.1  2.3 +0.01  0.15 +0.14  0.22 GRE +0.1  2.2 +0.00  0.15 +0.13  0.21 G BIAS  SDEV E-GRD [mm] 0.02  0.14 2 with respect to GPS (G) to reflect lower accuracy of precise products and models. A positive effect of multi-constellation is visible at all parameters, and particularly in terms of the standard deviation, while the impact of GLONASS is more significant compared to Galileo. It is expected due to a lower number of operational Galileo satellites as well as longer support of GLONASS with precise models and products in the scientific community. As already discussed, ZIM2-ZIMJ differences indicate a bias of about 2–3 mm in ZTD which has been decreased partly in multi-GNSS solutions. The improvements in all parameters reached 15–30% in terms of RMSE at both dual-stations. Table 3.8 then shows the impact of the backward smoothing on all solutions using single- or multi-constellation data. All the above mentioned characteristics are similar to the Kalman filter, and the backward smoothing then improved mainly standard deviations (by about 25%). 3.3.7.2 Carrier-Phase Post-fit Residuals and Slant Delays Figure 3.29 shows the carrier-phase post-fit residuals when using the Kalman filter PPP (left) and the backward smoothing PPP (right) for multi-GNSS solutions supported with the COM (top) and GBM (bottom) MGEX products. The carrierphase residuals are useful indicators of an overall performance of the solution including the quality of input products and models. Showing plots for the ZIM2 station only, below discussed characteristics are common to other stations too. First, we observe a common elevation-dependent pattern of characteristics of post-fit residuals when using elevation-dependent observation weighting. Second, the backward smoothing does not change the distribution of the carrier-phase post-fit residuals significantly. The main effect of the backward smoothing is thus understood mainly as improved accuracy of the estimated parameters. The tropospheric slant delays reconstructed from the model parameters and post-fit residuals will thus 90 J. Douša et al. Table 3.8 Statistics (BIAS  SDEV) for backward smoothing using GPS (G), GLONASS (R), Galileo (E) Station pair ZIM2ZIMJ ZIM2ZIMJ ZIM2ZIMJ MAT1MATE MAT1MATE MAT1MATE GNSS G BIAS  SDEV ZTD [mm] +2.7  1.1 BIAS  SDEV N-GRD [mm] +0.11  0.12 GR +2.3  1.0 +0.06  0.11 0.02  0.09 GRE +1.9  1.0 +0.07  0.12 0.04  0.09 1.3  1.6 0.04  0.15 +0.22  0.19 GR +0.6  1.4 +0.00  0.12 +0.16  0.17 GRE +0.5  1.4 0.01  0.11 +0.16  0.16 G BIAS  SDEV E-GRD [mm] 0.02  0.10 Fig. 3.29 Carrier-phase post-fit residuals from the Kalman filter (left) and the backward smoothing (right) benefit primarily from the improvement of the parameters. Third, the GPS residuals (black) are the smallest and compact compared to other systems indicating actual quality of precise models and products. Galileo shows the largest residuals, however, we had to substitute various precise models, in particular station antenna phase centre offsets and variations by using the values from GPS models. Due to the same reason, we may notice systematic changes in the elevation-dependent redistribution of Galileo residuals (red) after applying the backward smoothing. Fourth, we can notice about twice larger post-fit residuals GLONASS (green) when using the GBM product compared to the COM product. As the characteristics are common to all the stations, it indicates a lower quality of GLONASS orbits and clocks from the 3 Advanced GNSS Processing Techniques (Working Group 1) 91 GBM product or some inconsistent models used for the product generation and in the PPP software. 3.3.8 Multi-GNSS Solutions and Products Z. Deng GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] GPS PWV is considered to have observation noise of about 1~2 mm, but GPS PWV sometimes shows larger noise and jumps for some stations and/or at certain times, suggesting that a lower number of GPS satellites and poor line-of-sight condition limits the quality of ZTD estimates in such stations and at such times. Many GPS networks are now being upgraded to multi-GNSS observation networks, and this upgrade is expected to be beneficial for GNSS tropospheric monitoring. The IGS network is being upgraded to be capable to observe multi-GNSS (GPS, GLONASS, Galileo, BeiDou and QZSS) signals. GFZ has started to provide multi-GNSS orbit and clock for all the constellations (GBM) since middle 2014 (Fig. 3.30). To validate ZTD results from a global GPS-only and multi-GNSS analysis we processed multi-GNSS station data from the global ground tracking MGEX-network (Fig. 3.31) spanning a 3-month time period from December 2014 to February 2015. Color indicates maximum number of satellites available in addition to GPS per observation epoch. Multi-GNSS and GPS-only global network solutions were generated to study the impact of including additional GNSS on estimated ZTD. ZTD difference time series Fig. 3.30 Number of satellites per GNSS constellation included in the global multi-GNSS observation data processing 92 J. Douša et al. Fig. 3.31 Site selection for multi-GNSS processing derived from pre-defined sites associated to GFZ Rapid routine solution were computed for sites which provided multi-GNSS observation data. Dots depict sites suitable for comparing GPS-only versus multi-GNSS (GREC) results. Figure 3.32 shows the mean bias of the resulting ZTD differences which varies in the range of 1.5 mm. No specific latitude or longitude dependency could be identified. Moreover, there is no obvious correlation between the magnitude of the bias and the number of additional satellites used to estimate the ZTDs (compare to Fig. 3.31). Figure 3.33 shows associated standard deviations derived from the ZTD difference time series with values below 3.5 mm (Deng et al. 2015). ZTDs estimated with multi-GNSS processing are more stable than those based on GPS only. Sudden jumps observed in GPS-only ZTD are significantly reduced with the multi-GNSS processing. Because the number of satellites in multi-GNSS solution is more than twice that of only GPS observation, the noise due to rising and setting satellites is mitigated thus reducing the size of sudden jumps in ZTD. 3.4 PPP and Ultra-Fast GNSS Tropospheric Products A majority of E-GVAP ACs till now uses a double-difference observation processing in the network solution. This strategy eliminates clock errors at GNSS receiver and satellite while public products were not available in near real-time (NRT). The situation has changed in 2012, when the IGS introduced the Real-Time Service (RTS, http://rts.igs.org) providing GPS and GLONASS orbit and clock corrections by combining contributions from several IGS real-time analysis centres 3 Advanced GNSS Processing Techniques (Working Group 1) 93 Fig. 3.32 Mean bias for zenith total delay estimates between GPS-only and multi-GNSS (GRCE) solution. Differences are shown only for those sites actually providing multi-GNSS observation. Left and top subfigures show mean bias distribution w.r.t. latitude and longitude, respectively Fig. 3.33 Shows associated standard deviations (StDev) derived from the ZTD difference time series with values below 3.5 mm. In contradiction to the biases, the STDs reveal a small latitude dependency with larger magnitudes for sites below the equator 94 J. Douša et al. (Caissy et al. 2012). The IGS RTS aims at supporting real-time (RT) analyses with the Precise Point Positioning (PPP) method (Zumberge et al. 1997). The PPP is based on original observations or their linear combination without differencing between receivers or satellites. Though German Research Centre for Geosciences (GFZ) has provided a NRT PPP ZTD product (Dick et al. 2001; Gendt et al. 2004) since 2001, it was possible only thanks to their two-step processing approach consisting of (1) a global NRT solution for determining consistent satellite clock and orbit products and (2) a distributed PPP processing for ZTD estimated for each station individually. With the availability of global real-time data flow, software and standards specified for precise product dissemination, the PPP is becoming more popular for the troposphere monitoring. Compared to the traditional approach in E-GVAP dominated by the doubledifference network processing, the PPP offers several advantages: (a) an easy production in real-time or NRT fashion, (b) flexible use of central or distributed processing scheme including a receiver built-in solution, (c) an estimation of tropospheric parameters in the absolute sense with a high spatio-temporal resolution, and (d) an optimal support of all satellite constellations and new signals including multiple frequencies; all profiting from a highly efficient and autonomous processing approach. The price for mentioned advantages is however paid by several disadvantages. Compared to the strategy using double differences, all observation models need to be carefully applied to reach the best accuracy. In addition, integer ambiguity resolution is possible only if precise observation phase biases are available, thus often non-integer-fixed ambiguities are usually estimated. 3.4.1 Real-Time Data and Product Dissemination Y. Altiner BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] W. Söhne BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] A. Stürze BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] PPP is a method for high accuracy positioning using observations from a single GNSS receiver, suited for both, real-time and post-processing implementations. Traditionally, PPP is an idea of a post-processing technique for efficient evaluation of GNSS data from large scale networks. But it also enables a real-time positioning for stable (static) and movable (kinematic) objects with an accuracy of centimeter and sub-decimeter level, respectively. The fundamental advantage of PPP is that the 3 Advanced GNSS Processing Techniques (Working Group 1) 95 number of simultaneously observed stations within a global network can be significantly increased without decrease in accuracy of station coordinates. This benefit and the short duration of data processing has given the PPP method popularity to be used it efficiently also in the field of climate research, in particular to support weather forecasting techniques for medium or large scale areas. Using PPP, estimation of coordinates for positioning and determination of ZTD is possible every second (real-time PPP). However, for PPP in real-time some parameters, such as precise satellite coordinates (orbit), earth orientation parameters, and satellite clock corrections are needed from an external source, e.g., available online from the IGS RTS. The IGS RTS is generating and providing the variables of the SSR related to the orbits and clocks of GNSS satellites, an indispensable essential for real-time PPP (http://www.igs.org/rts). The RTS products created by several analysis centres contain GNSS satellite orbit and clock corrections to the broadcast ephemeris. Orbit corrections are provided as along-track, cross-track and radial offsets to the Broadcast Ephemeris in the ECF reference frame (Earth-centred and Earth-fixed). RTS corrected orbits are expressed within the ITRF implemented during the realtime GNSS observations. Clock corrections are expressed as offsets to the Broadcast Ephemeris satellite. Hereby, attention should be paid that the reference point of the satellite clocks is selected in accordance with the reference point of the satellite orbits. The IGS RTS is providing three combination solutions, IGS01, IGS02 and IGS03. While IGS01 is generated on the basis of epoch-wise combination, the Kalman filter technique is exploited for producing IGS02 and IGS03. Two different agencies are responsible for these RT products: European Space Agency (ESA) provides IGS01 and BKG provides IGS02 and IGS03. All mentioned streams include orbits/clock corrections to GPS satellites, and only IGS03 supports also GLONASS constellation. While IGS01 and IGS02 are combined from up to eight individual solutions, IGS03 has only four individual contributors. The RTS correction streams are formatted with respect to the RTCM (Radio Technical Commission for Maritime Services) standard for SSR and are transmitted using the NTRIP protocol (Networked Transport of RTCM via Internet Protocol). NTRIP was developed in co-operation between the Informatikzentrum Dortmund in Germany and BKG (http://igs.bkg.bund.de/ntrip) and initiated as an industrial standard since 2004 (Weber et al. 2005). Afterwards, NTRIP was standardized by the Special Committee 104 “DGNSS” of RTCM. The communication between the major components of the NTRIP, i.e. the server, the caster and the client are handled through HTTP ports. It is to mention that the major software components of the NTRIP are developed under “GNU General Public License”. NTRIP allows disseminating hundreds of data and product streams simultaneously for a few thousand users when applying the modified Internet Radio Broadcasting Software. It is also to note that a GNSS stream typically needs not more than 5 kbit/s bandwidth. The currently used version 2 of NTRIP, downward compatible to version 1, was completed in 2009. BKG supports the distribution of the new technology by providing the so-called Professional NTRIP Caster. This tool has been developed in 96 J. Douša et al. cooperation of BKG with Alberding Company in Germany for administration, configuration, and implementation the piece of software running using the LINUX operating system and is widely used within the IGS and EUREF (http://www.epncb. oma.be/). Meanwhile, almost every new GNSS receiver is coming with the NTRIP option. 3.4.2 BKG Real-Time Analysis Development and Contribution Y. Altiner BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] W. Söhne BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] A. Stürze BKG, Federal Agency for Cartography and Geodesy, Frankfurt, Germany e-mail: [email protected] Since April 2016, BKG is contributing to the COST real-time demonstration campaign for troposphere estimation providing solutions processed with BNC, BKG’s own software tool. The development of BNC started in 2005 by BKG in collaboration with different partners, e.g. with Technical University Prague, Czech Republic, and Alberding Company, Germany. BNC is a software for simultaneously retrieving, decoding, converting and processing or analyzing real-time GNSS data streams applying the NTRIP standard. BNC has been developed within the framework of EUREF and the IGS. Although BNC is primarily intended for real-time GNSS applications, it may also be run offline by transmission of data from an external file to simulate real-time observation conditions or for post-processing implementations. BNC provides also different modes of data evaluation like “graphics or interactive” mode to illustrate the processing state and results and “no windows” mode as well. A major module of the BNC is the option “Real-Time PPP” for positioning in real-time according to the SSR model which was first provided in 2010. To meet requirements of the PPP in real-time using the state variables of the SSR model, within this version additional messages are provided to the user, among others, satellite orbit and clock corrections for GPS as well as GPS and GLONASS combination, and ionospheric corrections as well as biases for code and phase data. The PPP module of BKG allows users a real-time positioning worldwide at sub-decimeter-level using code and phase data in ionosphere-free solutions using P3 or L3 linear combinations in static or kinematic mode within an observation time of 10 min. In 2011, the “Real-Time PPP” module was expanded by incorporating the 3 Advanced GNSS Processing Techniques (Working Group 1) 97 PPP implementation for post-processing to work offline including data from external files. In 2014, BNC comes into being able to process data of multiple stations in a single BNC job (one job for all stations). Within this version (BNC 2.12), the troposphere parameters estimated by the data processing is provided by SINEXTRO format (v0.01) to allow the usage of the TRO-file as a priori or a posteriori model for GNSS applications (Weber et al. 2016). Within the RT demonstration campaign, ZTD values were evaluated in 5-min intervals for 22 mount-points at the beginning. As of October 2017, all 32 mountpoints are implemented. As agreed, orbit and clock corrections from the IGS03 product are used. Streaming of broadcast ephemeris is taking place through the realtime access to the broadcast ephemeris stream of BKG “RTCM3EPH” for GPS and GLONASS. It should be noted that the IGS03 correction stream is a combined stream of four individual solutions providing GPS and GLONASS and is done by BKG with BNC. The evaluation of ZTD implemented by BNC is using the Kalman filtering method. One important parameter to set in the configuration file is the white noise (signal) for 1 h which can be according to the weather conditions. The standard value for the variation of the white noise is 36 mm/h1/2. Other effects or parameters influencing the accuracy of positioning so far considered by BNC are listed in Table 3.9. Processing of data within the demonstration campaign takes place using an elevation cut-off angle of 7 for observation usage and considering the float ambiguity resolution. To study the impact assessment of GNSS processing on tropospheric products, two different products are created by BKG. A GPS-only and a combined GPS plus GLONASS product is computed and submitted separately as GPS and GPS + GLO products on hourly basis to the central analysis centre of the project (http://www. pecny.cz/COST/RT-TROPO/). Each hourly solution contains 12 ZTD-values and each ZTD-value represents a time span of 5 min within the relevant hourly solution. In total, 10,234 files were submitted on hourly basis from March 14, 2016 to July 31, 2017 including solutions using the real-time GPS measurements and state Table 3.9 Important effects so far considered by BNC within data processing Effect Earth’s tide Earth’s rotation (Movement of pol) Phase-wind up Ionosphere Troposphere Multipath (phase shift of the signal) Atmospheric and hydrological loading Ocean tides Offset and phase centre variations of satellite and ground antennas (PCV) Considered Yes Yes Yes First order terms eliminated using L3 Determined No No No Yes 98 J. Douša et al. variables from the IGS03 product stream. For GPS + GLO, the quantity of the submitted hourly files reached to a total number of 8988 between March 18, 2016 and July 31, 2017. Data transmission continues also beyond the date 31.7.2017 (Fig. 3.34). The internal precision in terms of agreement of both BKG solutions – GPS-only versus GPS + GLONASS – is shown in Figs. 3.35, 3.36 and 3.37. The example covering a 15-day period in July 2017 shows a good overall agreement for the ZTD Fig. 3.34 The quantity of ZTD processed by PPP using orbit and clock corrections from IGS03 product for GPS and GPS + GLO observations. The solutions were created in 5-min intervals and combined to a single file on hourly basis to be submitted to the central analysis centre of the project (http://www.pecny.cz/COST/RT-TROPO/). Each submitted hourly solution contains 12 ZTD-values and each ZTD-value represents a time span of 5 min within the relevant hourly solution Fig. 3.35 Time series of ZTD estimates from BKG’s real-time solutions (GPS-only and GPS + GLONASS) as taken from the uploaded COST format files for stations THTI (left) and WTZR (right) 3 Advanced GNSS Processing Techniques (Working Group 1) 99 Fig. 3.36 Time series of differences of ZTD estimates between BKG’s real-time solutions “GPSonly” minus “GPS + GLONASS” for stations THTI (left) and WTZR (right). The plots are showing increased differences during periods of re-initialization of at least one of both solutions Fig. 3.37 Mean of differences between time series of ZTD estimates from BKG’s real-time solutions GPS-only minus GPS + GLONASS for 15 days in July 2017. Only three of 21 stations show a bias of almost 2 mm ZTD. The large standard deviation for station NRMD cannot be explained by frequent re-initialization; the time series show larger portions of disagreement time series for Tahiti (Fig. 3.35 left) and Wettzell. A closer look into the differences, however, shows some portion of disagreement. This can be explained by re-initialization of at least one of both solutions. Figure 3.37 shows the statistics for 21 stations (station ADIS was only occasionally available in the GPS + GLO solution). For the majority of the stations the mean bias is well below 2 mm ZTD. Regarding the submitted ZTDs of BKG within the demonstration campaign, a statistical study for station WTZR was conducted in relation to the EUREF weekly 100 J. Douša et al. combined solutions (Table 3.10). To do this, the EUREF weekly combined ZTD solutions were considered as target values and BKG real-time GPS and GPS + GLO solutions as measured values, respectively. The measured values from BKG were subtracted from the EUREF target values. In total, 8451 ZTD values from each solution were included to this statistical study. The correlation between the differences of GPS and GNSS solutions in relation to the EUREF combined solutions amounts to 85%. This suggests good coincidence between both solutions (GPS and GPS + GLO). The average differences for GPS and GPS + GLO are on the order of 11.4 mm and 8.6 mm, respectively. Contrary to the average, the standard deviation of GPS is smaller than the standard deviation of GPS + GLO (11.4 mm and 12.9 mm, respectively). To illustrate the relation of the GPS and GPS + GLO solutions to the EUREF combined solution the total 8451 ZTD values were reduced to 253 through choosing a random value for each day at 15:30 as illustrated in Fig. 3.38. It is important here to Table 3.10 Statistical aspects between ZTDs determined by BKG for station WTZR using GPS and GPS + GLO real-time corrections and the combined EUREF weekly solution. In total, 8451 ZTD values from each solution were included in the study. Averages and standard deviations, shown in bold, determined relative to the EUREF combined ZTD solution (EUREF minus BKG-GPS and EUREF minus BKG-GNSS) Solution for station WTZR BKG-GPS BKG-GPS + GLO EUREF Average in mm 11,8 8,1 Std. dev. in mm 11,4 12,9 Average of std. dev. in mm 0,7 0,7 0,5 Fig. 3.38 Relation of the GPS and GPS + GLO solutions to the EUREF combined solution for station WTZR. 8451 ZTD values from each solution were reduced to 253 through choosing of a random value for each day at 15:30. The results suggest good coincidence between the GPS and GPS + GLO solutions using the real-time orbit and clock corrections 3 Advanced GNSS Processing Techniques (Working Group 1) 101 notice that the EUREF weekly combinations include the results of several ACs which are derived using the final orbits of the IGS. Each AC does not use the same software for data processing, and network geometries of ACs are also different. Some studies conducted by BKG between the results of the real-time, near real-time und post-processing applications also suggest the good quality of real-time PPP with respect to the determination of ZTD (Altiner et al. 2009, 2010 and 2011). 3.4.3 Assessment of IGS RTS Orbits and Clock Corrections and GOP Real-Time Tropospheric Products8 L. Zhao Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] The accuracy of the RT ZTD calculated with the PPP method strongly depends on the quality of RT GNSS orbit and clock corrections (Douša and Václavovic 2014; Hadaś and Bosy 2015). We evaluated publicly available global RT products and we summarized our ZTD contributions to the RT Demonstration Campaign initiated in 2015 by this COST Action. We also studied the impact of IGS RTS (Caissy et al. 2012) on the simulated RT ZTD estimates within the GNSS4SWEC Benchmark Campaign (Douša et al. 2016). 3.4.3.1 Assessment of Real-Time Orbit and Clock Corrections We investigated the performance of four real-time products (Douša et al. 2018a) having been collected and archived at Geodetic Observatory Pecny (GOP) using the BNC Software (Weber et al. 2016) since 2013: IGS01, IGS02 and IGS03 the official IGS RTS combined products, and CNS91 (also known as CLK91) as an individual solution provided by CNES RT AC (Laurichesse 2011). Two different strategies and 8 Parts from this section were previously published in Douša et al. (2018a). 102 J. Douša et al. Fig. 3.39 Monthly statistics of availability of satellite corrections from IGS01 RT stream software are used for combining the IGS RTS (Hadaś and Bosy 2015). All mentioned streams include orbits/clock corrections to GPS satellites, and only IGS03 supports also GLONASS constellation. Navigation data from MGEX (Montenbruck et al. 2017) needs to be used together with the RT corrections to recover the precise satellite orbit and clocks. First, the availability and completeness of RT corrections were checked and, second, satellite orbit and clocks were compared to IGS final orbit and clocks, both during the period of 2013–2017. Figure 3.39 depicts monthly completeness of RT corrections for all GPS satellites from the IGS01 combined products. The others IGS products and the CNES product are generally showing similar performance. From the comparison, we can classify problems into three groups: (1) temporal unavailability period of some satellites, e.g. G03, G04, (2) source-specific unavailability, e.g. G01 for CNS91, and (3) satellite-specific incompleteness. The first group is usually caused by the loss of observations due to the upgrade of satellite, such as replacing the old Block IIA satellites with the new Block IIF satellites or a maintenance identified by satellite unhealthy status. The second and third groups of gaps are caused by data unavailability from a global network and the processing strategy including outlier detection in the product generation. The availability of the corrections is significantly lower for some months (June 2015, December 2016) compared to others which was caused by the internet connection failures at GOP when receiving the streams. The source-specific loss of data at IGS02 and CNS91 streams are visible in June 2015 and, these are mainly due to the inconsistent navigation message Issue of Date (IOD) available from the MGEX broadcast and those referred by RT corrections. It can be thus recommended to use consistent RT navigation data and precise correction streams optimally guaranteed by the same provider. In general, the availability of RT corrections is well over 90% for most satellites which agrees with findings in (Hadaś and Bosy 2015). It indicates that the RT corrections were provided continuously for use in troposphere monitoring, however, problems can be expected in a kinematic positioning which is more sensitive to the product incompleteness. 3 Advanced GNSS Processing Techniques (Working Group 1) 103 Table 3.11 RT clocks and orbits components compared to the IGS final products IGS01 IGS02 IGS03 CNS91 RMSE [mm] Radial 1.84 2.35 2.41 2.68 Along 2.83 3.71 3.82 3.07 Cross 2.38 3.04 3.10 2.47 3D 4.34 5.63 5.70 5.01 Clock 5.72 10.08 10.28 11.16 SDEV [mm] Clock 2.95 3.52 3.03 2.29 Apart from the availability of the corrections, the precision is critical for the user performance. The orbits are compared in 5-min intervals for three components: radial, along-track and cross-track while the clock comparison is based on the second order difference method. The IGS08 and the IGS14 model is used to correct satellite PCOs prior and after January 29, 2017, respectively, corresponding to the adoption of the IGS14 reference frame (Rebischung et al. 2016). The clock datum is estimated by calculating average over all satellites clocks at each epoch. The datum inconsistencies are then eliminated through single-differences between individual satellite clocks and the clock datum. The single-differences from the real-time clocks are compared to those from the IGS final product. The root-mean-square error (RMSE) of RT orbits and clocks are calculated for each day while outliers are removed using a fixed threshold. Although there is a strong correlation between clocks and radial orbit component, we haven’t corrected this dependency. Table 3.11 gives summary statistics for all products over all days. The orbit difference in radial component shows the smallest RMSE for all products, whereas the along-track and cross-track components reached slightly larger values. The IGS01 orbit shows the best agreement with respect to the IGS final orbits. Largest differences are observed for the orbits from IGS03, which might be attributed to a different outlier detection method applied when including GLONASS satellites. Time evolution of the orbit comparison for each product and specific component is shown in Fig. 3.40. Coordinate differences greater than 30 cm are plotted at the top horizontal lines of each graph. Orbits from the IGS01 stream are less affected by the outliers compared to IGS02 and IGS03 products as indicated by outliers mainly during March 2015. The switch from the IGS08 to the IGS14 PCO model (January 28, 2017) can be observed in statistics of the radial component. It seems that the CNS91 product used the new IGS14 model as of March 9, 2017, while official IGS solutions are difficult to recognize due to most likely asynchronous switches by different contributing providers. Otherwise, the orbit accuracy for all products shows an overall good consistency over the period. Table 3.11 also summarizes RMSE and standard deviation (SDEV) of the realtime clock corrections. The former represents the accuracy relevant for the processing of code pseudoranges while the latter characterizes the precision important for the carrier-phase processing. It can be also interpreted from the PPP point of view combining both observation types as follows – the former have a positive impact on the PPP convergence time while the latter enable more precise positioning within already converged solution (Ye et al. 2018). Obviously, this is the case of 104 J. Douša et al. Fig. 3.40 Daily RMS of real-time orbits with respect to IGS final orbit Fig. 3.41 Daily statistics SDEV (top) and RMSE (bottom) of real-time clocks IGS01 and CNS91 products when the first is more accurate, but the second more precise for the PPP application. The IGS02 and IGS03 products performs slightly worse in terms of both RMSE and SDEV. Figure 3.41 finally shows time series of the clock SDEV and RMSE statistics. The former (top plot) indicates a comparable high quality over the period for IGS01 3 Advanced GNSS Processing Techniques (Working Group 1) 105 and CNS91, while more outliers are observed for IGS02 and IGS03 including the problematic period in 2015 identified in the orbit availability evaluation. The clock RMSE from IGS01 is the lowest and the most stable compared to the others during the period while, the RMSE of CNS91 clocks was more accurate during 2015 when compared to the other years. 3.4.3.2 Impact of IGS RTS Products on ZTD Estimates The impact of the IGS RTS products on PPP ZTD estimates was assessed by exploiting the GNSS4SWEC Benchmark campaign with 400 GNSS stations in central Europe during the period of May–June, 2013. The ZTD was calculated using the G-Nut/Tefnut software in the post-processing mode when supported with two precise products: (1) IGS final orbits and clocks, and (2) IGS RTS orbit and clock corrections. Two reference solutions were provided for the benchmark using different software and processing strategies (Douša et al. 2016). GOP used the BSW and the double-difference processing (DD) and GFZ used the EPOS software and the PPP method. Statistics from the comparison of both testing solutions with respect to both reference products are given in Table 3.12. Generally, the results indicate a good agreement, however, the impact of the IGS RTS products (IGS01) on ZTDs is clearly visible in two aspects: a) a common systematic error of 2.4–2.8 mm, and b) a lower precision of 13–17%. Interestingly, a better agreement in terms of SDEV is reached between 10% and 20% when using two PPP solutions (G-Nut/ Tefnut vs EPOS software) compared to the processing strategies (DD vs PPP). The results also showed that input products and the processing strategy might result in a similar impact on the ZTD estimates which can reach up to 20% in terms of accuracy. Finally, it should be noted that the PPP ZTD estimation used a stochastic model and an epoch-wise filtering method in the G-Nut/Tefnut software (Václavovic et al. 2013), while a deterministic model with the least-squares batch adjustment used in the EPOS software. Table 3.12 Summary statistics from the comparison of PPP ZTD results using two inputs (IGS01 RT vs. IGS final) w.r.t. EUREF reprocessing G-Nut/Tefnut PPP Input precise products IGS final (SP3 files) IGS01 RT corrections IGS final (SP3 files) IGS01 RT corrections ZTD reference product GOP final (BSW52/ DD) GOP final (BSW52/ DD) GFZ final (EPOS/ PPP) GFZ final (EPOS/ PPP) Bias [mm] +0.9 STD [mm] 5.1 RMS [mm] 5.2 +2.4 5.8 6.4 +0.4 4.1 4.2 +2.8 4.9 5.7 106 3.4.3.3 J. Douša et al. Long-Term Quality of Operational RT ZTD Production The RT ZTD from the demonstration campaign is evaluated for 18 European stations during the initial year of the GNSS4SWEC Real-time Demonstration campaign. Two GOP solutions using the IGS03 product are compared with respect to the EUREF 2nd reprocessing combined tropospheric product (Pacione et al. 2017): (1) GOPR – standalone GPS solution, and (2) GOPQ – GPS + GLONASS solution. In Table 3.13, we can observe a systematic error in ZTD of about 2 mm in the longterm evaluation, similar as observed in the simulated real-time processing in the benchmark campaign, see the previous subsections. Although GLONASS observations are down-weighted by a factor of 2 in our solution in order to reflect the lower quality of GLONASS precise products, a small positive impact on the ZTD is observed in terms of mean bias (10%) and mean SDEV (7%), both calculated over 18 stations. Figure 3.42 shows the comparison of the GOPR solution with respect to the EUREF combined product during the first year of the RT demonstration campaign. Monthly mean ZTD biases, standard deviations and their 1-sigma scatter calculated over all 18 stations indicate a long-term stability of the operational real-time production with a small seasonal effect in SDEV due to a less accurate troposphere modelling during the summer period (Douša and Václavovic 2016). 3.4.4 Real-Time Product Development and Evaluation at ROB E. Pottiaux Royal Observatory of Belgium, Brussels, Belgium e-mail: [email protected] In the framework of COST Action ES1206 (GNSS4SWEC), the Royal Observatory of Belgium (ROB) collaborated with the Geodetic Observatory Pecny (GOP) to use their Real-Time Precise-Point-Positioning (RT-PPP) software G-Nut/Tefnut for the real-time monitoring of the troposphere, to help support nowcasting of severe weather in Belgium. On June 21, 2014, ROB started to use G-Nut/Tefnut to produce real-time tropospheric products (ZTD and horizontal gradients) using 4 different processing strategies (Table 3.14), with a particular focus on Belgium. Table 3.13 Summary statistics over 18 stations from routine RT product using GPS and GPS + GLO data Solution description GOPQ – GPS + GLO GOPR – GPS BIAS [mm] mean  sdev 1.8  2.9 SDEV [mm] mean  sdev 6.7  1.2 RMSE [mm] mean  sdev 7.5  2.5 2.0  2.8 7.2  1.0 7.9  1.5 3 Advanced GNSS Processing Techniques (Working Group 1) 107 Fig. 3.42 Monthly summary biases and standard deviations of real-time ZTDs over 18 stations Table 3.14 Setup of the different ROB’s RT-PPP Processing Common parameters to all solutions Parameter Setup Coordinates Static estimation Tropospheric parameters ZTD + horizontal gradients Tropospheric model Saastamoinen + GMF + Chen and Herring Cut-off angle 3 Time resolution 10 s Latency 100 s Ocean tide loading Coef. FES2004 Antenna model IGS08 Antex file Solution naming & differences Solution name GNSS observations CLK + ORB product ROBA GPS-only IGS02 (GPS only) ROBB GPS IGS03 (GPS + GLONASS) ROBC GPS + GLONASS IGS03 (GPS + GLONASS) ROBD GPS + GLONASS CNS91 (GPS + GLONASS) Based on these developments, ROB also participated in the RT-PPP demonstration campaign, with the main goals to extensively develop, test, and validate realtime processing methods and tropospheric products that can help supporting nowcasting and forecasting of severe weather and foster the link to WG2 activities. ROB processes thus real-time GNSS observations from the 32 GNSS sites requested to participate in the demonstration campaign, along with those from 153 additional GNSS sites located worldwide (Fig. 3.43, left), including the complete Belgian dense network (Fig. 3.43, right). In total, 185 GNSS Stations are included in ROB’s RT-PPP processing with G-Nut/Tefnut (Fig. 3.44). These 185 stations are equipped with 73 different combinations of GNSS receivers and antennas (26 receiver types, 40 antenna type, Fig. 3.45), allowing thereby to study and assess the performances of this RT-PPP processing w.r.t. the equipment, and in fine to finetune accordingly the processing strategy. All RT-PPP products are formatted in both 108 J. Douša et al. Fig. 3.43 (Left) Location of the GNSS stations included in the RT-PPP operated by ROB using the G-Nut/Tefnut software from GOP. (Right): The location of the GNSS stations of the Belgian Dense network included in this processing Fig. 3.44 Number of GNSS stations included in the RT-PPP Processing operated by ROB Fig. 3.45 Number of stations equipped with a specific type of GNSS receiver (left) or antenna (right) and included in ROB’s operational RT-PPP processing campaign 3 Advanced GNSS Processing Techniques (Working Group 1) 109 in the COST-716 v2.2a and in the SINEX_TRO v2.00 format to ease exchange and validation. Upload of these products (only the stations contributing to the WG1 RT demonstration campaign) is done every hour to a central hub at GOP, and can be visualized at http://www.pecny.cz/COST/RT-TROPO. 3.4.4.1 Monitoring and Validation Developing and maintaining such new RT-PPP processing systems and products requires regular assessment, adaptation, and fine-tuning cycles. Therefore, ROB developed its own monitoring/validation system consisting of MySQL databases, statistical assessment programs, and a web-based user-interface to monitor continuously his RT-PPP processing system (graphs, reports. . .). The system is capable of monitoring/carrying out: • The campaign setup and its evolution • Inconsistency checks (configuration, equipment, models. . .) • The performances of the products at all GNSS stations included in the RT-PPP processing (the monitoring system developed by GOP can only monitor 17% of them), at various time scale and epochs (biases, precision, geographical dependency. . .) • Alarm systems in case of problem (dataflow, processing. . .) To validate the performances of the products at all stations and all time scales (e.g. a very rapid monitoring require the use of e.g. NRT products), the monitoring system uses various reference products listed in Table 3.15, from post-processed to NRT, some computed in PPP, some in Double-Difference (DD) approach etc. In all cases, the validation starts with a screening of the RT-PPP results to reject convergence period. This screening is using the RMS of the coordinates, ZTD and horizontal gradients, as well as the GDOP values, and the number of satellite measurements used to compute the tropospheric products at each single epoch. The advantages of this approach is that it can be implemented as a “real-time Table 3.15 List of reference products chosen to validate the RT-PPP products operationally computed by ROB Solution IGS Software BSW50 ROB PPP BSW52 ROB DD BSW52 ROB NRT E-GVAP BSW52 Orb. & Clk. IGS final CODE final CODE final IGS ultrarapid MultiGNSS GPS only Time Res. 5 min Latency 3–4 weeks NB. Sta. 30 TRO Est. ZTD + GRD GPS + GLO 5 min 3 weeks All ZTD + GRD GPS + GLO 1h 3 weeks All ZTD + GRD GPS only 15 min 1h All ZTD only 110 J. Douša et al. filtering” at the production level to screen out potential performance degradations before using the product in actual meteorological applications. In most cases, this filtering works well and had a low percentage of rejected values (<1%). Having various reference products is justified by their respective advantages and drawbacks, namely IGS for a standard reference product (but having only 30 common stations), consistency for PPP approaches (but probably less precise/accurate than the double-difference approach), accuracy/precision for final products (but large latency of the monitoring), and short-latency monitoring in the case of NRT solutions. The first long-term validation has been carried out with this system over the period July 1st 2014–March 31st 2015 (9 months). Globally and at the long-term level, the 4 solutions listed in Table 3.14 agrees very well with the IGS Final products (22 stations considered altogether). The linear regressions (computed over ~600,000 samples) between RT-PPP and the IGS products have a correlation coefficient above 0.99, a slope of ~0.99, and almost a zero intercept (0.02–0.03 mm). The bias is however station dependent, consistent at the mm level for all products, and ranges from 5 to 16 mm (Fig. 3.46). The standard deviation is also station dependent and ranges typically from 5 to 10 mm (Fig. 3.47), with a maximal value Fig. 3.46 Bias observed between the ZTDs from each ROB’s RT-PPP product and the IGS final troposphere product over the period July 1st 2014–March 31st 2015 Fig. 3.47 Standard deviation observed between the ZTDs from each ROB’s RT-PPP product and the IGS final troposphere product over the period July 1st 2014–March 31st 2015 3 Advanced GNSS Processing Techniques (Working Group 1) 111 Fig. 3.48 Bias observed between the ZTDs from each ROB’s RT-PPP product and the ROB final PPP troposphere product over the period July 1st 2014–March 31st 2015. Stations are ordered by increasing longitude Fig. 3.49 Standard deviation observed between the ZTDs from each ROB’s RT-PPP product and the ROB final PPP troposphere product over the period July 1st 2014–March 31st 2015. Stations are ordered by increasing longitude observed of 26 mm (ALIC). No clear geographical dependency (longitude, latitude, and altitude) of the bias or the standard deviation could be observed in this comparison. The comparison between the 4 RT-PPP solutions and the ROB’s PPP postprocessing reference solution showed very similar results (Figs. 3.48 and 3.49), with typical biases ranging from 5 to 16 mm, typical standard deviations ranging from 4 to 7 mm (i.e. slightly better than when compared to the IGS final product), and correlation coefficient above 0.99. One can also remark that the RT-PPP solutions performs very well and homogenously in Belgium: all stations in the middle of the graphs (i.e. from VEUR to VITH) have lower biases ranging from 1 to 4 mm. We can also note the slightly better global performance obtained with the CNS91 orbit and clock product. This lower standard deviation can be due to a faster convergence time for the solution using CNS91 in GPS + GLONASS mode, but this still needs to be confirmed. Finally, processing jointly GPS + GLONASS observations (ROBC, IGS03-GLO) provides consistent solutions as processing GPS-only observations (ROBB, IGS03), in general at the 1–2 mm level. 112 J. Douša et al. In conclusion, ROB started in June 2014 to compute operationally RT-PPP products as a demonstration campaign using the G-Nut/Tefnut software from GOP. This campaign has run now for more than 3 years. Results obtained so far are very promising but leave space for further studies and improvements. This RT-PPP production and their related developments will continue in the framework of the IAG WG 4.3.7 ‘Real-time GNSS tropospheric products’, with the aim on a more longer term to be operationally provided to E-GVAP. As the next natural step in the developments, ROB participates in a case study that aims to simulate the RT-PPP processing in off-line mode. The goals of this case study are explained in the next paragraph. 3.4.4.2 Assessment of RT Products in Simulated RT Analysis of Benchmark Data As we said, the RT-PPP demonstration campaign provided promising results but also leaved space for specific studies and improvements. Such kind of studies, like studying the optimal processing settings according to e.g. weather type (severe versus normal condition), studying the influence of the station equipment, or optimizing the convergence period etc. requires to be able to replay the processing several times on the same dataset by changing solely one processing parameter at a time. In other words: simulating the real-time processing in an offline mode for finetunings and assessments. One perfect candidate for this is the WG1 Benchmark campaign, and at the end of the COST Action, we started the real-time benchmark campaign with the following main objectives: • Processing (a core group of) stations of the Benchmark campaign in simulated real-time PPP mode with various real-time orbits and clocks products to produce a final assessment of the capability of each real-time orbit and clock product. • Develop optimal strategies for the estimation of ZTDs and tropospheric gradients at high (time) resolution (e.g. 5 min), and assess them to reference solutions and/or ZTD/GRD from NWM. • Investigate the capability of dynamically constrain ZTDs and GRDs according to the weather conditions (calm, moderate, turbulent, severe, etc.). This is very important for natural hazard warnings. • Produce a dataset of simulated real-time tropospheric products (ZTD, horizontal gradients, SPD) using the final fine-tuning of all previous steps that can be re-used to assess the capability/performances of these products in real studies/applications (e.g. a nowcasting case of severe weather). • Production of IWV in almost real-time for non-numerical nowcasting applications. • Standardize methods and format, and provide guidelines towards operational production of real-time tropospheric products (link to E-GVAP). In that context, ROB processed the complete WG1 benchmark dataset in simulated real-time offline mode with various processing options and various orbit and 3 Advanced GNSS Processing Techniques (Working Group 1) 113 Orbit and Clocks Products ZTD Constraints GRD Constraints IGS01 From 0.5 (tight) to 5.0 (loose) From 0.05 (tight) to 0.4 (loose) IGS02 By Step of 0.5 By step of 0.05 IGS03 (GPS-only) Fig. 3.50 The different processing configuration tested by ROB in the context of the real-time PPP benchmark campaign (simulated real-time PPP offline processing). It includes 4 orbit and clock products configurations (no offline version of the CNS91 product could be found), 10 ZTD constraints setup, and 8 horizontal gradients constraints setup clock products, focusing on providing the final assessment of the orbit and clock products, on testing different ZTD & GRD (dynamical/optimal) constrains, on studying the dependency of the performance w.r.t. to the equipment (see text of the RT demonstration campaign), and on convergence periods. Similarly as for the regular assessment of the real-time demonstration campaign, we also reprocessed the benchmark campaign with the BSW52 to produce final reference tropospheric products both in PPP and in double-difference approach. In total, 320 flavors of the RT-PPP products are available for inter-comparison (Fig. 3.50) and assessments. These results have been produced towards the end of the COST Action but will be analyzed in the context of the IAG WG 4.3.7. 3.4.5 GFZ Real-Time Product Development and Assessment in RT Analysis C. Lu GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] X. Li GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] The multi-constellation Global Navigation Satellite Systems (GNSS) offers promising potential for the retrieval of real-time (RT) atmospheric parameters to support time-critical meteorological applications, such as nowcasting or regional short-term forecasts. In this study, we processed GNSS data from the globally distributed MultiGNSS Experiment (MGEX) network of about 30 ground stations by using the precise point positioning (PPP) technique for retrieving RT multi-GNSS tropospheric delays. RT satellite orbit and clock product streams from the International GNSS Service (IGS) were used. Meanwhile, we assessed the quality of clock and orbit products provided by different IGS RTS ACs, called CLK01, CLK81, CLK92, GFZC2, and GFZD2, respectively. Using the RT orbit and clock products, the 114 J. Douša et al. performances of the RT ZTD retrieved from single-system as well as from multiGNSS combined observations are evaluated by comparing with the U.S. Naval Observatory (USNO) final troposphere products. With the addition of multi-GNSS observations, RT ZTD estimates with higher accuracy and enhanced reliability can be obtained compared to the single-system solution. Comparing with the GPS-only solution, the improvements in the initialization time of ZTD estimates are about 5.8% and 8.1% with the dual-system and the four-system combinations, respectively. The RT ZTD estimates retrieved with the GFZC2 products outperform those derived from the other IGS RTS products, Fig. 3.51. In the GFZC2 solution, the accuracy of about 5.05 mm for the RT estimated ZTD can be achieved with fixing station coordinates. The results also confirm that the accuracy improvement (about 22.2%) can be achieved for the real-time estimated ZTDs by using multi-GNSS observables compared to the GPS-only solution. In the multi-GNSS solution, the accuracy of real-time retrieved ZTDs can be improved by a factor of up to 2.7 in the fixing coordinate mode comparing with that in the kinematic mode, Fig. 3.52. 3.4.6 Contribution to RT Demonstration Campaign from ULX9 W. Ding University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected] F. N. Teferle University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected] After initial RT solutions based on BNC were abandoned, ULX modified the PPP-Wizard developed by CNES to provide a real-time solution to the RT-Demonstration Campaign. These solutions included GPS, GLO and GAL observations and employed the real-time products from CNES CLK93, including satellite orbit, clock and code/phase biases. The latter allowed PPP ambiguity resolution for GPS using a zero-difference ambiguity resolution approach (see Table 3.16). The modifications of PPP-Wizard included: • • • • • • • 9 Apply Antenna Reference Point (ARP) correction from igs08.atx Apply receiver PCO + PCV correction from igs08.atx Solid earth tide + ocean tide loading correction (FES2004) ZTD (GPT and Saastamoinen) + ZWD (modeled as random walk process) Troposphere Mapping Function (GMF) Elevation dependent weighting strategy (Q ¼ 1/cos(zen)**2) Of particular interest was the impact of multi-GNSS on initialization times. Parts from this section were previously published in Ding et al. (2017). 3 Advanced GNSS Processing Techniques (Working Group 1) 115 50 40 30 20 10 0 CLK01 Along ( cm ) 1 2 3 5 7 6 8 50 40 30 20 10 0 1 2 C02 C03 4 3 5 6 C06 CLK92 GFZC2 GFZD2 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 GPS PRN 7 8 C08 C09 9 10 50 40 30 20 10 0 C05 CLK81 C07 12 13 14 15 GLONASS PRN C11 C10 C12 C13 17 16 E11 C14 18 E12 19 E19 20 E22 21 22 E24 23 E26 24 E30 BDS+GALILEO PRN 50 40 30 20 10 0 Cross ( cm ) 1 2 3 5 7 6 8 50 40 30 20 10 0 1 2 C02 C03 4 3 5 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 GPS PRN 6 7 8 C08 C09 9 10 50 40 30 20 10 0 C05 C06 C07 12 13 14 15 GLONASS PRN C11 C10 C12 C13 17 16 E11 C14 18 E12 19 E19 20 E22 21 22 E24 23 E26 24 E30 BDS+GALILEO PRN 20 15 10 5 0 Radial(cm) 1 2 3 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 GPS 20 15 10 5 0 1 2 3 4 5 6 7 8 9 10 12 13 14 15 16 17 18 19 E12 E19 20 21 22 23 24 GLONASS 50 40 30 20 10 0 C02 C03 C05 C06 C07 C08 C09 C10 C11 C12 C13 C14 E11 E22 E24 E26 E30 BDS+GALILEO Fig. 3.51 RMS values of the differences between IGS RTS orbits (CLK01, CLK81, CLK92, GFZC2, and GFZD2) and GFZ final orbits for the four systems (i.e., GPS, GLONASS, Galileo, and BDS) in the along (top), cross (middle) and radial (bottom) components, respectively 116 J. Douša et al. CLK01 CLK81 CLK92 GFZC2 GFZD2 GPS RT ZTD (mm) 2600 2500 2400 2300 10min 20min 30min 40min 50min 60min 10min 20min 30min 40min 50min 60min 10min 20min 30min 40min 50min 60min 10min 20min 30min 40min 50min 60min 10min 20min 30min 40min 50min 60min 10min 20min 30min 40min 50min 60min G/R RT ZTD (mm) 2600 2500 2400 2300 G/R/E/C RT ZTD (mm) 2600 2500 2400 2300 Fig. 3.52 RT ZTD estimates at station ONS1 derived from the GPS-only (“GPS”, top), the combined GPS/GLONASS (“G/R”, middle), and the combined GPS/GLONASS/Galileo/BDS (“G/R/E/C”, bottom) solutions in fixing coordinate (left panels) and kinematic processing (right panels) modes by employing different IGS RTS over the first 2 h of DOY 090, 2017 Table 3.16 List of data processing modes investigated at ULX Modes RFLT GFLT GFIX MFLT MFIX Details Float PPP solution based on GLONASS-only observations Float PPP solution based on GPS-only observations Fixed PPP solution based on GPS-only observations Float PPP solution based on GPS/GLONASS observations Fixed PPP solution based on GPS/GLONASS/Galileo observations Reproduced from Ding et al. (2017) 3 Advanced GNSS Processing Techniques (Working Group 1) 117 Fig. 3.53 Average initialization time in all data processing modes. (Reproduced from Ding et al. 2017) From Fig. 3.53 it can be found that the initialization time required by the RFLT solution is still the longest. It exceeds 30 min for two stations. Compared with that, the initialization time for the GFLT solution is shorter for all stations. The average value is 613 s (approx. 10.2 min). By applying ambiguity resolution, the initialization time becomes shorter for most of the stations, and is 583.6 s (approx. 9.7 min) on average. The initialization process can also be accelerated by utilizing GNSS observations, for which it can be achieved on average in 533 s (approx. 8.9 min). Again, this suggests that the effect of the observation geometry is larger than that of ambiguity resolution in accelerating the initialization process, especially considering that an initialization time is required to achieve the first ambiguity resolution. When applying all three GNSS and ambiguity resolution in the MFIX solution, the initialization process is finished on average in 508.3 s (approx. 8.5 min), and there are only small differences between different stations, which leads to the highest consistency in the solutions and reveals the benefit of GNSS observation and ambiguity resolution for severe weather event monitoring (for more details see Ding et al. 2017). Furthermore, the solutions from the five different RT data processing modes can be compared to the benchmark troposphere products from CODE and USNO, both being GPS + GLONASS solutions (Table 3.17). The statistics with respect to the two types of benchmark products are similar, which further validates the reliability and consistency of the reference products. The RMS of the RFLT solution for nearly all stations is smaller than 15 mm, and the average value is about 11.16 and 13.98 mm with respect to the final troposphere products from CODE and USNO. Compared with that, the RMS of the GFLT solution is better. The RMS of all stations is better than 12 mm except MOBS, and is about 9 mm on average. The worse performance of the RFLT solution may 118 J. Douša et al. Table 3.17 Mean accuracy of processing modes with respect to final troposphere products from CODE and USNO (reproduced from Ding et al. 2017) RFLT GFLT GFIX MFLT MFIX CODE Mean(mm) 0.82 0.83 2.09 0.47 1.48 STD(mm) 11.26 6.32 5.65 6.41 5.96 RMS(mm) 11.61 7.05 6.37 6.87 6.42 USNO Mean(mm) 0.61 0.59 2.03 0.41 1.52 STD(mm) 13.67 8.27 7.45 8.27 7.69 RMS(mm) 13.98 8.95 8.17 8.69 8.14 come from two points: (1) the number of GLONASS satellites is less than for GPS; (2) the accuracy of satellite products for GLONASS is worse than for GPS (Dach and Jean 2015). However, considering the accuracy requirements (10–15 mm) in updating NWP models, the RT troposphere estimates based on GPS or GLONASS only observations can both fulfill the requirements (De Haan 2006). Applying ambiguity resolution, the GFIX solution is further improved up to 0.8 mm on average compared to the GFLT solution. However, the mean bias becomes slightly bigger. Combining the observations of two systems, the MFLT solution is only 0.18 mm and 0.26 mm improved on average with respect to CODE and USNO products, which reveals that the accuracy is not greatly improved by incorporating GLONASS observations. In addition, the accuracy even becomes a little worse for some stations, which may be correlated with the weighting strategy between two systems and needs further research in the future. At last, the mean RMS of the MFIX solution is 6.42 mm and 8.14 mm with respect to CODE and UNSO products, respectively. It is the best solution among all the data processing modes, which again reveals the effect in utilizing GNSS observations and ambiguity resolution. A further comparison was performed with respect to radiosonde observations (Fig. 3.54). To summarize the accuracy of the 13 stations, we sort the results according to the distance between the GNSS station and the nearby radiosonde launch site. Since the mean bias of RT ZTD are monitored and will be corrected in the assimilation procedure, we will only calculate the standard deviation (StDev) of all stations (Bennitt and Jupp 2012). Based on the results, the accuracy of the RFLT solution is the worst, of which the StDev is especially larger and exceeds 15 mm in several stations. Among the other solutions, the StDevs are all smaller than 15 mm except for ABMF and JFNG. On average, the StDevs of the two single system solutions are 14.6 mm and 9.1 mm each, which again reveals that they can fulfill the requirements in monitoring severe weather events. However, compared with the GFLT solutions, we notice that the accuracies of the GFIX and MFLT solutions become a little lower for many stations. Since there is only one radiosonde observation in each day, this might be a consequence of the instability of GPS phase bias information and the satellite orbit/clock products for GLONASS. Additional details can be found in Ding et al. 2017. 3 Advanced GNSS Processing Techniques (Working Group 1) 119 Fig. 3.54 Standard deviation of RT ZTD errors with respect to the radiosonde observations in all data processing modes (reproduced from Ding et al. 2017) 3.4.7 New Adaptable Strategy for RT and NRT Troposphere Monitoring10 P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] In order to optimize the accuracy and the timeliness of tropospheric parameters with possibly prioritizing the product latency in real-time (RT) or its accuracy in near realtime (NRT), we have developed a new adaptable strategy in the G-Nut/Tefnut software (Douša et al. 2018a) exploiting the Precise Point Positioning, PPP (Zumberge et al. 1997). We aimed to support various users from a single processing PPP engine combining RT and NRT analysis modes. All the parameters, ZTD, GRD and STD, are estimated consistently for both product lines when supported by realtime orbit and clock products. The strategy based on data from an individual station only is also ready to exploit optimally all available multi-GNSS observations when supported by precise products. 10 Parts from this section were previously published in Douša et al. (2018a). 120 3.4.7.1 J. Douša et al. Epoch-Wise Filtering vs. Batch Processing Nowadays, a common procedure of the NRT analysis within the EIG EUMETNET GNSS Water Vapour Programme (E-GVAP, http://egvap.dmi.dk) is based on the least-squares adjustment (LSQ) analysis and a piece-wise linear function of the modelling of tropospheric parameters within the processing interval. By using the BSW52 (Dach et al. 2015), Geodetic Observatory Pecný (GOP) applies the same strategy in all long-term contributions to E-GVAP (Douša 2001a, b, Douša and Bennitt 2013). The estimated tropospheric parameters for such case are displayed in Fig. 3.55 by black dots and the piece-wise deterministic model by dash lines connecting the parameters within each hourly product update; not necessarily connected at update boundaries. According to the E-GVAP conventions, the last product value is shifted by 1 min (HR:59) in order to avoid a duplicate value with the next hour product update. So far, the tropospheric horizontal linear gradients and slant delays were not provided by GOP because of two reasons: (1) not being yet assimilated into numerical weather models (NWM), and (2) significantly increase the number of parameters in the network solution and, consequently, the computation time. From this point of view, the epoch-wise processing (e.g. Kalman filter) is an optimal strategy as it, contrary to the LSQ, estimates recurrently all unknown parameters in every epoch. One of the consequences is that only previous observations contribute to a current estimate and thus the solution needs a certain time to converge. However, an accuracy of parameters estimated during the initial convergence, or any later re-convergence, can be improved only by using the backward smoothing algorithm (Václavovic and Douša 2015) additionally using both past and following observations for improving the precision at every epoch. As a consequence, the backward smoothing was designed for the post-processing solutions and it can substitute the LSQ in a number of applications. In the first step of the filter, parameters are predicted via adding particular amount of noise to diagonal elements of the variance-covariance matrix belonging to dynamic parameters. When new observations are available the state vector with its variance-covariance matrix are Fig. 3.55 Various strategies of the troposphere modelling: (1) piece-wise model (black dots connected with dash lines) and (2) stochastic modelling by real-time Kalman filter (white points), hourly backward smoothing (blue points), and 45-min postponed hourly backward smoothing (red points) 3 Advanced GNSS Processing Techniques (Working Group 1) 121 updated. Results from the prediction as well as from the update needs to be stored for future smoothing. Since the standard Kalman filter and also smoother can suffer from numerical instabilities due to the round-off error we have implemented alternative forms of the algorithms exploiting Cholesky and Singular Value Decompositions of the variance-covariance matrix instead of the original one. 3.4.7.2 RT and NRT Combined Processing Supported by Observations from Files or Streams In E-GVAP, the NRT products are updated on hourly basis and the delivery requirement to the E-GVAP server at UK Met Office is 45 min after the last observation used in the analysis. The NRT LSQ batch processing, initiated every hour when obtaining a majority of data files, is thus a relevant solution for this purpose when using hourly data files. Although the backward smoothing approach is the most beneficial for the post-processing solutions, it can be effectively used in NRT applications too. We have thus used it for the new adaptable strategy combining a continuously running forward filter with a regularly triggered backward smoothing filter. The former is aimed for the estimating epoch-wise tropospheric parameters in real-time indicated by white points in Fig. 3.55. The backward smoothing is started periodically at a pre-defined time stamp as indicated at HR:00 in the Figure. It uses the initial state vector from the same epoch provided by the realtime filter for recalculating the past parameters from the Kalman filter. Obviously, such recalculation is able to refine significantly older parameters, but cannot improve parameters at the initial epochs of the backward smoothing. The length of the smoothing period can be set flexibly to reflect an actual user preference for a higher accuracy or a shorter latency of the product. In such way, the standard E-GVAP NRT tropospheric product can also be provided on hourly basis as shown by blue points in the Figure. Initially, the new adaptable strategy was designed to use the precise orbit and clock corrections disseminated through RT streams, it can however exploit observations coming from both real-time (streams) and near real-time (hourly or sub-hourly files). While the former is necessarily used in a simultaneous RT and NRT product generation, the latter is applicable for NRT only. In any case, both data flows can be mixed and analysed for each station independently. Additional advantage of the NRT analysis utilizing the backward smoothing and RT observations may profit from the 45-min requirement in E-GVAP for the product delivery in NRT. By starting the backward smoothing shortly before the delivery request, indicated in Fig. 3.55 by the red arrow starting at 12:40, new observations can further improve the accuracy of NRT product, in particular last parameters of the NRT product and still reduce systematic errors typical for data interval boundaries. 122 3.4.7.3 J. Douša et al. Estimating High-Resolution ZTDs and Horizontal Gradients Estimating high-resolution parameters with the batch LSQ can be more difficult due to increasing size of normal equations, in particular for the network processing, and consequently extending the processing time due to the inversion of the large matrix of the solution. On the contrary, the described combination of the Kalman filter and smother seems to be more convenient for monitoring dynamical processes when high-resolution parameters are required. Since only observations from past epochs are used in the Kalman filter processing, achieved parameters cannot react on their fast change, and batch processing can reach better precision in this case. However, when backward smoothing is applied, following observations improve the state vector estimation when preserving high resolution from the previous forward filter. Such situation is presented in Fig. 3.56 showing the ZTDs from different solutions during the fast change in the troposphere, indicated with a sudden decrease of ZTD by 6–7 cm during 2.5 h at POTS station on May 10, 2013. Obviously, the postprocessing 15-min ZTDs from the GFZ solution using the EPOS software (Gendt et al. 2004) implementing the LSQ processing (gray points) and the daily smoothed 30-s stochastic ZTDs from our software (black dots) are in a very good agreement. Note that ZTDs from our solution is resampled to 5 min in the Figure. Due to the use of past observations only, the real-time Kalman filter (red points) shows a delay in the change of estimated parameter. The random walk for ZTD was set to 5 mm/h1/2. The Figure reveals a characteristic behaviour for the real-time processing when using past observations in stochastic parameter estimation. We can observe significant improvements mainly in reducing the systematic behaviour for the 1-h backward smoothing (green points). The main improvement can be reached within 15–30 min, as indicated for the ZTD from the backward smoothing at 12:00, 13:00 and 14:00. Interestingly, the ZTD from the 2-h smoothing already shows a very good agreement with the post-processing LSQ solution. It shows a potential improvement discussed in the previous subsections and postponing the NRT smoothing by at least 30 min if RT observations are available. Fig. 3.56 Real-time (red), near real-time (green, blue) and reference (black) ZTD estimates during a fast change in the troposphere 3 Advanced GNSS Processing Techniques (Working Group 1) 3.4.7.4 123 Assessment of New Method Compared to the Existing E-GVAP Processing The new strategy has been initially developed and assessed using the GNSS4SWEC Benchmark campaign (Douša et al. 2016) and firstly compared to the GOP NRT tropospheric solution contributing operationally to E-GVAP. Thought the new strategy can provide RT and NRT products in high temporal resolution, we compared only the ZTD as a product of HH:00 and HH:59 time stamps in every hour representing the standard NRT E-GVAP product, see Fig. 3.55. Table 3.18 summarizes results of three strategies and six ZTD solutions using 13 EUREF stations selected from the benchmark campaign and exploiting the EUREF combined ZTD product as a reference for all comparisons (Pacione et al. 2017). The table shows summary statistics indicating similar improvements in terms of the standard deviation over all ZTDs estimated at HH:00 in NRT independently on applied products (IGS RTS vs. IGS final products), processing strategies and software (G-Nut/Tefnut PPP vs. BSW52 DD). Compared to the Kalman filter ZTD estimated at the last epoch (HH:59), the backward smoothing running on hourly basis showed the improvement of 20% and 24% for the IGS RTS and the IGS final product, respectively. The E-GVAP/GOP product demonstrates a similar improvement (24%) comparing ZTD from HH:00 against HH:59, which corresponds to our previous results (Douša and Souček 2005). The new adaptable PPP solution using the IGS final orbits reached the same accuracy as the E-GVAP/GOP product using the IGS ultra-rapid orbits and NRT DD network solution from the BSW. On the other hand, the use of IGS RTS products instead of IGS final products in PPP indicates a degradation of 18% in ZTD SDEV and a 2.5 mm bias. It should be finally noted, that the E-GVAP/GOP solution and the reference EUREF solution are based on a similar processing strategy and the software, while the new strategy is significantly different. Table 3.18 Summary statistics of three processing strategies and six ZTD solutions compared to EUREF combined tropospheric product Solution RT PPP (HR:59) NRT PPP (HR:00) PP PPP (HR:59) PP PPP (HR:00) NRT DD (HR:59) NRT DD (HR:00) Software G-Nut/ Tefnut G-Nut/ Tefnut G-Nut/ Tefnut G-Nut/ Tefnut BSW52 BSW52 Strategy description Kalman filter, simulated realtime solution Hourly backward smoothing in real-time Kalman filter in offline processing, IGS final Hourly backward smoothing with IGS final Last ZTD of hourly PW linear LSQ First ZTD of hourly PW linear LSQ Latency <5 min ~ 60 min < 5 min ~ 60 min ~ 90 min ~ 30 min Mean BIAS 2.4 mm Mean SDEV 5.7 mm 2.5 mm 4.6 mm 0.1 mm 4.7 mm 0.2 mm 3.6 mm 0.4 mm 4.9 mm 0.2 mm 3.7 mm 124 J. Douša et al. Fig. 3.57 ZTD standard deviations and biases for 13 EUREF stations and six processing strategies compared to the reference EUREF product Figure 3.57 shows standard deviations and biases individually for all stations comparing the first ZTDs (HH:00) and the last ZTDs (HH:59). The statistics of ZTDs from the E-GVAP/GOP solution (PPP:DD-ultra) are plotted in red and pink for HH:00 and HH:59, respectively. The results from the new strategy using the PPP with IGS final orbit and clock products (PPP:IGS-final) are shown in dark and light blue and, using IGS RTS (PPP:IGS03) in black and grey. Standard deviations for all the stations show a similar improvement in the ZTD SDEV over all the strategies, software and precise products. However, there is no significant impact of the strategy on systematic errors, and we can observe only a common positive bias attributed to the use of the IGS RTS products. 3.4.8 Optimum Stochastic Modeling for GNSS Tropospheric Delay Estimation in Real-Time11 T. Hadaś Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] K. Kaźmierski Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] P. Hordyniec Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] 11 Parts from this section were previously published in Hadaś et al. (2017). 3 Advanced GNSS Processing Techniques (Working Group 1) 125 J. Bosy Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] F. N. Teferle University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected] It is commonly accepted by the GNSS community to constrain epoch-wise ZWD estimates, usually by estimating ZWD as a random walk parameter. Wrocław University of Environmental and Life Sciences (WUELS) have shown that the optimum ZWD constraints in real-time GNSS processing, modelled as a randomwalk process, should be time and location specific. A typical approach is to perform an initial empirical testing in order to obtain the effective constraining. As an alternative, WUELS proposed to take benefit from numerical weather prediction models to define optimum random walk process noise (RWPN). Two different strategies were proposed and validated. In the first approach an archived ZTD time series can be used to calculate a grid of yearly means of the difference of ZWD between two consecutive epochs divided by the root square of the time lapsed, which can be considered as a random walk process noise. Using archived VMF-G grids, we obtained RWPN global grids for hydrostatic and wet parameter (Fig. 3.58). We noticed that grids are nearly identical year by year, with differences below 1 mm/√h for hydrostatic and wet grids. This means that a single RWPN grid can be implemented in a software as a look up table to define the optimum wet RWPN value for any station located worldwide. It was shown that RWPN values from grids are similar to those obtain with empirical testing (see Fig. 3.58). Alternatively, a short-term weather forecast can be used to perform ray-tracing in order to obtain forecast of ZTD and then to calculate RWPN dynamically in realtime. This approach was validated using forecast from GFS4 model. In this case Fig. 3.58 Hydrostatic (left) and wet (right) yearly mean RWPN grid for 2015 126 J. Douša et al. Fig. 3.59 Comparison of wet RWPN, ZTD time series, standard deviations of real-time ZTD residuals superior results were obtained, by means on the accuracy and precision of estimated tropospheric delay (see Fig. 3.59). The advantage of this approach is that the wet RWPN is regularly adjusted to the current tropospheric conditions. Its value remains low, when ZTD is stable over time, and rises when a rapid change of ZTD is expected. More details can be found in Hadaś et al. (2017). 3.5 Exploiting NWM-Based Products for Precise Positioning Numerical Weather Model (NWM) data and/or derived climatologies are increasingly used in precise geodetic applications. A prominent example is the Vienna Mapping Function, VMF (Boehm et al. 2006a), and the Global Mapping Function, GMF, (Boehm et al. 2006b). In recent years, the tropospheric models became more sophisticated including NWM tropospheric (1st/2nd order) gradient estimates, NWM tropospheric mapping factors, advantages of utilization of high-resolution NWM models and NWM predictions, and others. The use of NWM-derived parameters plays more significant role in both a long-term reprocessing and real-time kinematic positioning. The following section provides an overview of achievements in this respect within the COST Action. 3 Advanced GNSS Processing Techniques (Working Group 1) 3.5.1 127 Tropospheric Parameters from Numerical Weather Models F. Zus GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: zusfl[email protected] The tropospheric delay T is approximated as T ðe; aÞ ¼ mf h ðah; bh; ch; eÞ ∙ Z H þ mf w ðaw; bw; cw; eÞ ∙ Z W þ mf G ðC; eÞ½N cos ðaÞ þ E sin ðaފ ð3:20Þ where e and a denote the elevation and azimuth angle of the station satellite link, mfh and mfw denote the hydrostatic and non-hydrostatic Mapping Function (MF), mfg denotes the gradient MF, Zh and Zw are the zenith hydrostatic and non-hydrostatic delay, and N and E denote the (first-order) gradient components. The elevation angle dependency of the hydrostatic (non-hydrostatic) MF is based on the continued fraction form proposed by Marini (1972) and normalized by Herring (1992) to yield the unity at zenith mða; b; c; eÞ ¼ ð1 þ a=ð1 þ b=ð1 þ cÞÞÞ=ð sin ðeÞ þ a=ð sin ðeÞ þ b=ð sin ðeÞ þ cÞÞÞ ð3:21Þ The elevation angle dependency of the gradient MF is based on the form proposed by Chen and Herring (1997): mgðC; eÞ ¼ 1=ð sin ðeÞ tan ðeÞ þ C Þ ð3:22Þ where C ¼ 0.003. Therefore, provided that the tropospheric parameters ah, bh, ch, aw, bw, cw, zh, zw, N and E are known, the tropospheric delay can be assembled for any station satellite link. The tropospheric parameters are determined from ray-traced tropospheric delays. The required pressure, temperature and humidity fields are taken from a Numerical Weather Model (NWM). At GFZ Potsdam the algorithm proposed by Zus et al. (2014) is used to compute mapping factors and slant factors, i.e. the ratios of slant and zenith delays. Note that mapping factors are computed under the assumption of a spherically layered troposphere. From the mapping factors and slant factors the tropospheric parameters are estimated by least-squares fitting. This is done separately for the mapping function coefficients and the gradient components. 128 3.5.1.1 J. Douša et al. Mapping Function Coefficients For each station, 10 hydrostatic (non-hydrostatic) mapping factors are computed for elevation angles of 3, 5, 7, 10, 15, 20, 30, 50, 70, 90 degree, and the hydrostatic (non-hydrostatic) MF coefficients are determined by least-squares fitting (Zus et al. 2015a). 3.5.1.2 Gradient Components At first, 120 slant factors and corresponding mapping factors are computed at elevation angles 3, 5, 7, 10, 15, 20, 30, 50, 70, 90 degree and the azimuth angles 0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330 degree. Second, zenith delays are applied to obtain azimuth-dependent and azimuth-independent slant total delays. Third, the differences between azimuth-dependent and azimuth-independent slant total delays are computed. Finally, the gradient components are determined by leastsquares fitting (Zus et al. 2015b). 3.5.1.3 Tropospheric Model Errors The pressure, temperature and humidity fields are taken from Global Forecast System (GFS) of the National Centers for Environmental Prediction (NCEP). The NCEP’s GFS analyses are available every 6 h (00:00, 06:00, 12:00, 18:00 UTC) with a horizontal resolution of 0.5 degree on 31 pressure levels. Tropospheric parameters are derived for a grid with a resolution of 0.5 degree. For any gridpoint on Earh’s surface we examine how well the modeled tropospheric delays, i.e., the tropospheric delays assembled from the tropospheric parameters, match the true (ray-traced) tropospheric delays. This is done by calculating the elevation angle dependent postfit-residual. sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi   2 P  S e; a j T e; a j RðeÞ ¼ n ð3:23Þ where T denotes the assembled tropospheric delay and S denotes the ray-traced tropospheric delay. The elevation angles are chosen to be 3, 5, 7, 10, 15, 20, 30, 50, 70, 90 degree and the azimuth angles are chosen to be 0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330 degree. As an example we consider a single epoch (first of May 2017, 12:00 UTC). At first, we examine azimuth independent tropospheric delays. Figure 3.60 shows the residuals for an elevation angle of 7 . The residuals are well below 1 mm for any point on Earth’s surface. For comparison Figure 3.61 shows the residual for an elevation angle of 7 if we replace the NWM based MF by the GMF (Boehm et al. 2006a, b). Clearly, as the GMF is based on a climatology it cannot capture the short- 3 Advanced GNSS Processing Techniques (Working Group 1) 129 MF rmse [mm] elevation 7 deg 2 80 1.8 60 1.6 latitude [deg] 40 1.4 20 1.2 1 0 0.8 –20 0.6 –40 0.4 –60 0.2 –80 50 100 150 200 250 longitude [deg] 300 350 0 Fig. 3.60 We examine azimuth independent tropospheric delays (1st of May 2017, 12:00 UTC). The scatter plot shows the tropospheric residuals for an elevation angle of 7 GMF rmse [mm] elevation 7 deg 30 80 60 25 40 latitude [deg] 20 20 15 0 –20 10 –40 5 –60 –80 50 100 150 200 250 300 350 0 longitude [deg] Fig. 3.61 We examine azimuth independent tropospheric delays (1st of May 2017, 12:00 UTC). The NWM based MF is replaced by the GMF. The scatter plot shows the tropospheric residuals for an elevation angle of 7 130 J. Douša et al. 50 MF GMF average rmse per station [mm] 45 40 35 30 25 20 15 10 5 0 3 5 7 10 elevation angle [deg] 15 Fig. 3.62 We examine azimuth independent tropospheric delays (1st of May 2017, 12:00 UTC). The bar plot shows the average residual per grid point as a function of the elevation angle. Different colors show different options; NWM based MF and GMF term variability of the troposphere and hence the residuals are larger. Figure 3.62 shows the average residual per grid point as a function of the elevation angle. The residuals for the NWM based MF are well below 1 mm for any elevation angle emphasizing the fact that the three term continued fraction form of the mapping function works with an exquisite level of precision. Hence, if the troposphere is indeed spherically layered, the tropospheric model (its functional form) can be regarded error free. Next, we examine azimuth dependent tropospheric delays. Figure 3.63 shows the residuals for an elevation angle of 7 when no gradients are applied and Fig. 3.64 shows the residuals for an elevation angle of 7 when gradients are applied. As to expect the residuals are significantly reduced when gradients are applied. The residuals can be reduced by adding higher-order gradients to the tropospheric delay model. In essence, the tropospheric delay is approximated as T ðe; aÞ ¼ mf h ðah; bh; ch; eÞ ∙ Z H þ mf w ðaw; bw; cw; eÞ ∙ Z W þ mf G ðC; eÞ½GN cos ðaÞ þ GE sin ðaފ þ mf G ðC; eÞ½FcosðaÞ cos ðaÞ þ GcosðaÞ sin ðaÞ þ HsinðaÞ sin ðaފ ð3:24Þ where F, G and H denote higher-order gradient components. Figure 3.65 shows the residuals for an elevation angle of 7 when first- and higher-order gradients are applied. Finally, Fig. 3.66 shows the average residual per grid point as a function of the elevation angle. The application of higher-order 3 Advanced GNSS Processing Techniques (Working Group 1) 131 00 rmse [mm] elevation 7 deg 50 80 45 60 40 latitude [deg] 40 35 20 30 0 25 –20 20 15 –40 10 –60 5 –80 50 100 150 200 longitude [deg] 250 300 350 0 Fig. 3.63 We examine azimuth dependent tropospheric delays (1st of May 2017, 12:00 UTC). The scatter plot shows the tropospheric residuals for an elevation angle of 7  when no gradients are applied 01 rmse [mm] elevation 7 deg 50 80 45 60 40 latitude [deg] 40 35 20 30 0 25 –20 20 15 –40 10 –60 5 –80 50 100 150 200 longitude [deg] 250 300 350 0 Fig. 3.64 We examine azimuth dependent tropospheric delays (1st of May 2017, 12:00 UTC). The scatter plot shows the tropospheric residuals for an elevation angle of 7  when first-order gradients are applied 132 J. Douša et al. 01+02 rmse [mm] elevation 7 deg 50 80 45 60 40 latitude [deg] 40 35 20 30 0 25 –20 20 15 –40 10 –60 5 –80 50 100 150 200 250 longitude [deg] 300 0 350 Fig. 3.65 We examine azimuth dependent tropospheric delays (1st of May 2017, 12:00 UTC). The scatter plot shows the tropospheric residuals for an elevation angle of 7  when first- and higherorder gradients are applied. 90 O1+O2 O1 O0 average rmse per station [mm] 80 70 60 50 40 30 20 10 0 3 5 7 10 elevation angle [deg] 15 Fig. 3.66 We examine azimuth independent tropospheric delays (1st of May 2017, 12:00 UTC). The bar plot shows the average residual per grid point as a function of the elevation angle. Different colors show different options; no gradients applied, first-order gradients applied and first- and higher-order gradients applied. 3 Advanced GNSS Processing Techniques (Working Group 1) 133 gradients reduces the residual in particular for elevation angles below 10 . We recall that the underlying NWM has a horizontal resolution of 0.5 degree. Small scale tropospheric features cannot be represented by the low horizontal resolution. Larger residuals are expected when the horizontal resolution of the underlying NWM increases. 3.5.2 The Impact of Global and Regional Climatology on the Performance of Tropospheric Blind Models G. Möller Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] J. Sammer Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] For tropospheric delay modelling usually the concept of mapping functions is used: ΔLðεÞ ¼ ΔLhz ∙ mf h ðεÞ þ ΔLwz ∙ mf w ðεÞ ð3:25Þ It describes the total slant delay ΔL at elevation angle ε as as the sum of a hydrostatic and a wet component. Each component can be expressed as the product of a zenith delay and the corresponding mapping function. Empirical troposphere models like GPT2w (Böhm et al. 2015) can provide this information for any user position and epoch. In particular, GPT2w is based on global 1  1 gridded values of tropospheric state parameters like surface pressure for modelling the zenith hydrostatic delay or water vapour pressure, weighted mean temperature and water vapour decrease factor for modelling of the zenith wet delay. In addition, mapping coefficients are provided separately for the hydrostatic and the wet mapping function, to further reduce the mapping error, especially below 15 degrees elevation angle (Möller et al. 2014). Comparison with time series of the IGS reveal that GPT2w allows for modelling of the zenith tropospheric delay with a bias of less than 1 mm and a RMS of about 3.6 cm, see (Böhm et al. 2015). However, to further improve the model performance on regional level, the global 1  1 grid was replaced by a regional 0.2  0.3 grid. Analogous to the global grid, on the regional climatological grid tropospheric parameters are provided as mean, annual and semi-annual coefficients. The regional coefficients were derived from 3 years of ALARO model data, as provided by the Central Institute for Meteorology and Geodynamics (ZAMG), Austria on 18 pressure levels with a temporal resolution of 3 h. Figure 3.67 shows time series of pressure (p), temperature (T), temperature lapse rate (dT), mean temperature (Tm), water vapour pressure (e) and water vapour 134 J. Douša et al. Fig. 3.67 Meteorological parameters as obtained from GPT2w (blue) and GPT2w ALARO (red) at GNSS station HART, Austria. Analysed period: May–Dec 2013 Fig. 3.68 Histogram of ZTD residuals (GNSS minus model) as obtained at 45 GNSS in Austria. (Left) GNSS minus GPT2w, (right) GNSS minus GPT2w ALARO. Analysed period: May– Dec 2013 decrease factor (Lambda), as obtained from both grids, exemplary for GNSS station HART, Austria. The regional model (GPT2w ALARO) shows in general a more distinct seasonal signal than the global GPT2w. Especially in temperature lapse rate and water vapour decrease factor, significant differences can be observed. In order to evaluate the impact on the delay modelling based on both grids, ZTD time series were derived for about 45 GNSS sites in Austria and neighbouring countries. Figure 3.68 shows the ZTD residuals for GPT2w and GPT2w ALARO with respect to ZTDs derived from dual-frequency GNSS observations at GNSS site. The regional grid helps to further reduce the ZTD bias at GNSS sites, which are located north and south of the main Alpine ridge. Nevertheless, averaged over all 3 Advanced GNSS Processing Techniques (Working Group 1) 135 GNSS station and over the period May to December 2013, the mean bias (GPT2w: 2 mm and GPT2w ALARO: 3 mm) and standard deviation (GPT2w: 29 mm and GPT2w ALARO: 30 mm) of GPT2w ALARO is very comparable to the global GPT2w. In consequence, it is concluded that a higher spatial resolution does not lead consequently to a better performance of the tropospheric blind model, even not in the Alpine area. In order to further increase the performance of tropospheric blind models in future other strategies have to be discovered. 3.5.3 Refined Discrete and Empirical Troposphere Mapping Functions VMF3 and GPT312 D. Landskron Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] J. Boehm Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] The Vienna Mapping Functions 3 (VMF3) is a refinement of VMF1 with the aim of even higher precision. It eliminates shortcomings in the empirical coefficients b and c and is not only tuned for the specific elevation angle of 3 , but for the whole elevation range through least-squares adjustments. The new mapping function coefficients were determined on the basis of ray-traced delays of the ray-tracer RADIATE (Hofmeister and Böhm 2017). Comparing modeled slant delays of VMF3 and VMF1 with the underlying ray-traced delays proves the high quality of VMF3, in particular at low elevation angles. In consequence, when requiring highest precision, VMF3 is to be preferable to VMF1. For more details, the reader is referred to Landskron and Böhm (2017). Figure 3.69 shows the empirical (blind) troposphere model Global Pressure and Temperature 3, GPT3 (Landskron and Böhm 2017) as a refinement of the model Global Pressure and Temperature 2 wet, GPT2w (Böhm et al. 2015), with re-calculated mapping function coefficients and empirical horizontal gradients, available in a horizontal resolution of 5  5 and 1  1 . The meteorological quantities remain unchanged. The empirical mapping factors from GPT3 are averaged from the VMF3 applying information from 2D ray-tracing through numerical weather models of the ECMWF. GPT3 is full consistent with VMF3 and can be used for GNSS as well as VLBI analysis. It is unique in such a way as it provides the entire information which is required in order to model a priori the troposphere delay dependent on elevation angle and azimuth. 12 Parts from this section were previously published in Landskron and Böhm (2017). 136 J. Douša et al. Fig. 3.69 Hydrostatic north gradient (left) and hydrostatic east gradient (right) from GPT3. While the north gradient shows systematic features for the northern and southern hemisphere owing to the atmospheric bulge, the east gradient is rather affected by the land-ocean distribution 3.5.4 The Impact of NWM Forecast Length on ZTDs13 J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] M. Eliaš Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] P. Krč Czech Institute of Computer Science, Academy of Sciences, Praha, Czech Republic e-mail: [email protected] K. Eben Czech Institute of Computer Science, Academy of Sciences, Praha, Czech Republic e-mail: [email protected] J. Resler Czech Institute of Computer Science, Academy of Sciences, Praha, Czech Republic e-mail: [email protected] To study the impact of the GNSS tropospheric corrections on the NWP prediction length, we calculated ZTDs for stations in the Benchmark campaign using the G-Nut/Shu software and the Weather Research and Forecasting (WRF) models operated routinely by the Institute of Computer Science, Academy of Sciences, 13 Parts from this section were previously published in Douša et al. (2015a, b). 3 Advanced GNSS Processing Techniques (Working Group 1) 137 Czech Republic (ICS ASCR). The two same WRF-ICS analyses have routinely contributed to the Real-time Demonstration campaign since July 2015 and can be characterized as follows: (a) two regional domains 9  9 km (EU9) and 3  3 km (CZ3), (b) uniform horizontal grid represented using the Lambert Conformal Conic projection (LCC), (c) grid unstaggered dimensions (west-east ¼ 418 grid points; south-north ¼ 302 grid points), (d) 38 vertical levels with the top level at 50 hPa, (e) four forecasts per day: 00:00, 06:00, 12:00 and 18:00 UTC, (f) a 1-h temporal resolution in each forecast, and (g) a 13-h length of the forecast. The data used for the WRF-ICS analysis computation are collected within a 6-h window surrounding the synoptic time T. The analysis then emerges on the web site of NCEP about 3 h 25 min after the synoptic time. The global model runs the forecasts for increasing time horizons up to 14 days. The resulted forecasts are successively uploaded until about T + 5 h. The mesoscale WRF model starts its simulation run at the synoptic time T, using previous analysis (T-6 h) and the time window from T-6 h to T. The so called grid nudging is performed at this level. In this way the mesoscale model reaches a state in time T, being a downscaled analysis. Grid nudging ensures that the mesoscale model doesn’t diverge too far from the global analysis in time T. The mesoscale simulation starts at approximately T + 3 h 35 min real time and during the grid nudging phase and the subsequent spinup phase, i.e. while simulating hours (T, T + 6), the mesoscale model simulation catches up the real time, until about T + 4 h 30 min real time. About T + 7 h real time when depending of computing resources, the forecast horizon T + 78 h is produced. Based on the WRF-ICS different forecasting intervals counted from the time T: 0-6 h (spinup), 6-12 h (forecast) and 12–18 h (forecast), we calculated ZTDs for all Benchmark stations during May/June 2013 and compared then with reference GNSS ZTD parameters. Figure 3.70 shows the comparisons of the NWP-based ZTDs with Fig. 3.70 ZTD biases (top) and standard deviations (bottom) calculated for 420 stations of the GNSS4SWEC Benchmark campaign over 2 weeks in May–June 2013 NWP models compared to GOP GNSS reference solution. From left to right three prediction intervals are shown: 0–6 h (left), 6–12 h (middle) and 12–18 (right) 138 J. Douša et al. Table 3.19 Summary statistics for the prediction length NWP domain D01/EUR D01/EUR D01/EUR D02/CZ D02/CZ D02/CZ Forecast window 0–6 h 6–12 h 12–18 h 0–6 h 6–12 h 12–18 h Bias [mm] 1.50 0.89 0.51 +2.05 +2.46 +3.20 Sdev [mm] 9.93 10.95 12.91 9.14 10.33 12.83 RMS [mm] 10.42 11.55 13.48 9.84 11.90 13.50 respect to GNSS reference solution in geographical plots. Biases (top) and standard deviations (bottom) are shown for the three prediction windows (left to right). Table 3.19 then provides a summary statistics for ZTDs calculated from two WRF domains. We focused on assessing zenith tropospheric delays potentially usable as external corrections for GNSS real-time applications such as positioning and navigation. Statistics for 14 days and 420 stations of the GNSS4SWEC Benchmark campaign [6] summarize ZTDs calculated from the spinup (0–6 h) and forecast intervals (6–12 h and 12–18 h). The summary shows that the quality of ZTD within the prediction windows up to 18-h resulted in RMS of 9.5–13.5 mm. It demonstrated a slow degradation of ZTDs from NWM approximately at a rate of 1–2%/h, which can be characterized by 1–2 mm/h. A degradation of mean biases was not observed. ZTDs usable in real-time GNSS applications correspond to the prediction of 6–12 h in the standard WRF operation. 3.5.5 Dual-Layer Tropospheric Correction Model Combining Data from GNSS and NWM14 J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] M. Eliaš Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] 14 Parts from this section were previously published in Douša et al. 2018b. 3 Advanced GNSS Processing Techniques (Working Group 1) 139 K. Eben Czech Institute of Computer Science, Academy of Sciences, Praha, Czech Republic e-mail: [email protected] Various tropospheric models for GNSS real-time positioning have been developed at GOP recently taking advantages of recent enhancements in the troposphere modelling (1) analytical ZWD calculation based on the concept of Askne and Nordius (1987) when combining exponential decay parameters from the water vapor pressure vertical profile and ZWD profile (Douša and Eliaš 2014), (2) more precise vertical approximation for the ZWD parameter using new parameter for modelling of ZWD exponential decay expressed either for a dependency on pressure or altitude (Douša and Eliaš 2014), (3) more simple and accurate tropospheric model based on user parameters related to the hydrostatic and non-hydrostatic zenith delays (Douša et al. 2015a, b), and (4) a flexible parameterization of new and legacy ZWD modelling approaches in various user modes, (Douša et al. 2015a, b). Based on the abovementioned enhancements in the modelling of tropospheric corrections, GOP has developed a new concept of a dual-layer tropospheric correction model for GNSS precise real-time positioning applications which optimally benefit from the synergy between NWM and GNSS data. The idea behind is to combine and predict optimally hydrostatic and wet components of the total zenith path delay in support of real-time GNSS positioning applications. In a simple form it resembles the approach of an assimilation of GNSS ZTD into the numerical weather forecast, however, targeting the GNSS user parameters, i.e. tropospheric path delays of an electromagnetic signal up to 15GHz frequencies. The new concept aimed at enhancing the original GOP augmentation model introduced in Douša et al. (2015a, b) by combining NWM data with ZTDs estimated from GNSS permanent stations in regional networks. The first layer is represented with the background NWM-driven model available anytime for the region of interest when derived purely from a NWM forecast provided from several hours up to 1–2 days. The first-layer parameters are predicted from NWM data fields and provides the background ZHD, ZWD together with auxiliary parameters for the parameter vertical scaling. The second layer improves mainly the ZWD, or more precisely, optimized corrections complementary to the NWM-derived hydrostatic corrections of the tropospheric effect on the at radio-frequency electromagnetic signal. The description of the combination method is given in the paper by Douša et al. (2017) and includes four test cases for optimal method of ZWD combination from GNSS and NWM. Two NWM models were used: (1) global ERA-Interim reanalysis (Dee et al. 2011) provided by ECMWF (we used 6-h updates with 1  1 horizontal resolution), and (2) high-resolution regional NWM prediction (European model with 9  9 km horizontal uniform resolution, hourly updated) provided from an operational short-time prediction from the mesoscale WRF model operated by the Institute of Computer Science (ICS), Academy of Sciences of the Czech Republic. The input GNSS ZTDs stems from the analysis of GNSS data of the GNSS4SWEC WG1 Benchmark dataset (Douša et al. 2016) provided by GOP as the GNSS reference solution. 140 J. Douša et al. Fig. 3.71 Mean statistics of a single-layer model (top) and dual-layer (bottom) ZTDs with respect to GNSS ZTDs using a closed-loop test case and ERA Interim (left) and WRF-ICS (right) models during May, 2013 The concept was implemented in the G-Nut/Shu software for the testing and assessing different variants of ZWD weighting method. The closed-loop results demonstrated that the new combination concept provides highly accurate and stable results despite of several involved approximations, see Fig. 3.71. The GNSS data were able to improve significantly the non-hydrostatic part of the NWM-based tropospheric model, however it is able to correct also possible errors in the hydrostatic component coming from the underlying NWM data. The most significant improvement of using the second model layer was reached in term of ZTD standard deviations, in total over 43%, when assessed using products from independent GNSS stations. An improvement has also been observed in term of systematic errors, but these were almost negligible in the total statistics. The most important result of the combination of GNSS and NWM data was found in term of stability in the session-to-session performance of the dual-layer model ZTD statistics. The scatters in session-to-session mean standard deviations were reduced from 8.7 to 4.4 mm and from 9.6 to 3.0 mm for the WRF-ICS model when using GNSS ZTD and ZTD together with horizontal gradients, respectively, see Table 3.20. The GNSS contribution to the dual-layer augmentation model consists particularly in a local phenomena where the model resolution plays a key role for an optimal assimilation of GNSS ZTDs into the model. We demonstrated a fast decrease of the statistics when simulating a higher horizontal resolution for the ERA-Interim model and obtaining comparable or better results to those resulting from the high-resolution WRF-ICS model. Comparing several ZWD weighing methods in the combination, optimal results were achieved by replacing the original NWM ZWDs with GNSS ZWDs if these were interpolated robustly to the model grid points. Using tropospheric horizontal gradients estimated for each GNSS station, we developed a method to calculate so-called pseudo-ZTDs to densify GNSS ZTD field for its optimal contribution to the combination with NWM ZWDs stemming from the WRF-ICS high-resolution model, example is given in Fig. 3.72. The combination of GNSS ZTDs and pseudo-ZTDs with NWM showed up to 70% improvement 3 Advanced GNSS Processing Techniques (Working Group 1) 141 Table 3.20 ZTD mean statistics of GNSS and ERA-Interim data weighting within the ZWD combination and improvements with respect to the background NWM model (last line) Data reduction: Variant First layer (ERA) First layer (WRF) Second layer (ERA) Second layer (WRF) Data weighting [mm] σ GNSS ¼ 10.0; σ NWM ¼ 1 σ GNSS ¼ 10.0; σ NWM ¼ 1 σ GNSS ¼ 1; σ NWM ¼ 10.0 σ GNSS ¼ 1; σ NWM ¼ 10.0 None bias/sdev [mm] +2.6 / 9.3 R3 (33%) bias/sdev [mm] +2.6 / 9.3 R2 (50%) bias/sdev [mm] +2.7 / 9.4 0.8 / 9.6 0.7 / 9.5 0.8 / 9.6 +2.8 / 5.3 +2.7 / 5.2 +3.0 / 5.5 0.2 / 2.4 0.2 / 2.9 0.0 / 3.0 Fig. 3.72 Differences of ZTDs (points with circles) and pseudo-ZTDs at distances of 15 km and 25 km (small points) with respect to the reference GNSS ZTDs are showed for the model first layer (top panels) and for the model second layer (bottom panels) and ERA-Interim (left panels) and WRF-ICS (right panels), on May 31, 12:00 UTC, 2013 for standard deviations compared to the NWM model, and about 35% improvement compared to the initial GNSS ZTDs and NWM data combined model. As expected in this test case, we haven’t observed any improvement in case of the ERA-Interim model. As the developed method of NWM and GNSS data combination demonstrated a high accuracy and stability over time, it is very promising for implementation of a 142 J. Douša et al. service for tropospheric corrections for real-time precise positioning. Currently, we are implementing the operational prototype combining the existing GOP services for tropospheric parameters estimated (near) real-time (Douša and Václavovic 2014, 2016) provided in support of the numerical weather forecasting, and the existing demonstration prototype for the GOP model background layer utilizing the WRF-ICS operational weather forecast for Europe. More details can be found in Douša et al. 2018a, b. 3.5.6 Tropospheric Refractivity and Zenith Path Delays from Least-Squares Collocation of Meteorological and GNSS Data15 K. Wilgan Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] W. Rohm Institute of Geodesy and Geoinformatics, Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] J. Bosy Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] F. Hurter ETH Zurich, Zürich, Switzerland e-mail: [email protected] A. Geiger ETH Zurich, Zürich, Switzerland e-mail: [email protected] This subsection summarizes results of the troposphere model tested in two countries with different orography: Switzerland (mountainous) and Poland (mostly flat). The troposphere model is based on the least-squares collocation technique, where each observation is divided into a deterministic part, a correlated stochastic part (signal) and an uncorrelated stochastic part (noise). The selected parameters from different data sources are estimated simultaneously in the least-squares sense taking into account the two kinds of errors. The advantage of this method is a relatively easy implementation of additional data sources. Using the computed model coefficients it is possible to reconstruct the value of considered parameter at any time and place. 15 Parts from this section were previously published in Wilgan et al. (2017a) 3 Advanced GNSS Processing Techniques (Working Group 1) 143 Calculations were made using the software COMEDIE (Collocation of Meteorological Data for Estimation and Interpretation of tropospheric path delays), developed at ETH Zürich. For Switzerland, the profiles of total refractivity were calculated from three data sets: (1) the total refractivity calculated from ground-based meteorological measurements, (2) the total refractivity calculated from meteorological measurements and ZTD from GNSS stations and (3) the total refractivity calculated from meteorological measurements, ZTD GNSS and horizontal gradients of ZTD from GNSS stations. The data set (2) exhibits the best agreement with the reference RS data. The data set (1) based on the ground-based meteorological data also gives a good information about the tropospheric state, but only up to the height of the highest station. Above that height, it is necessary to use an additional data source, such as the ZTD GNSS, to provide an information about tropospheric parameters in the vertical profiles. Unfortunately, adding the horizontal gradients of ZTD in data set (3) did not improve the troposphere model. In the countries located mainly on lowlands such as Poland, the height distribution of ground-based meteorological stations is too flat to reconstruct the refractivity profiles with the collocation technique. Thus, the troposphere model is built based on NWP data of 10 km  10 km spatial resolution, 34 height levels and time resolution of 1 h. The total refractivity profiles and ZTD values were calculated from 4 data sets: (1) NWP WRF model, (2) WRF model integrated with ZTD GNSS, (3) ZTD GNSS only, and (4) WRF model, ZTD GNSS and ground-based meteorological measurements. The obtained total refractivity profiles were compared with the reference radiosonde (RS) data and obtained ZTD with the reference GNSS data from post-processing. For total refractivity, the best agreement with reference RS data was achieved from the WRF model integrated with GNSS data (data set 2). Including the ground-based meteorological data (data set 4) did not improve the model. Data set (1) exhibited similar agreement with RS measurements as data set (2), only in the upper layers of the model there was a displacement of 2 ppm. The data set (3) showed much worse accuracy with the discrepancies at lower altitudes even at the level of 30 ppm. Such large differences are a result of an attempt to reconstruct a whole profile of refractivity from a single ZTD value. In the next step, the model of ZTD from COMEDIE was compared with the reference near-real time GNSS data. The comparisons were made for 9 days in May 2014 that included a severe weather event. The data sets of the highest agreement with the reference GNSS data are (2) and (3) with the average biases of 3.7 mm and 3.8 mm and standard deviation of 16.7 mm and 17.2 mm, respectively. The collocation based only on the WRF model (data set 1) overestimates the ZTD values. The reason for such behavior is that the humidity values provided by the WRF model after the rainfall are too high, which directly affects the ZTD values. To sum up, the integration of NWP and GNSS data is essential to provide accurate models for both total refractivity and ZTDs. More details can be found in Wilgan et al. 2017a, b. 144 3.5.7 J. Douša et al. Improving Precise Point Positioning with Numerical Weather Models C. Lu GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] Precise positioning with the current Chinese BeiDou Navigation Satellite System is proven to be of comparable accuracy to the Global Positioning System (GPS), which is at centimeter level for the horizontal components and sub-decimeter level for the vertical component. But the BeiDou precise point positioning (PPP) shows its limitation in requiring a relatively long convergence time. In this study, we develop a numerical weather model (NWM) augmented PPP processing algorithm to improve BeiDou precise positioning. Tropospheric delay parameters, i.e., zenith delays, mapping functions, and horizontal delay gradients, derived from shortrange forecasts from the Global Forecast System (GFS) of the National Centers for Environmental Prediction (NCEP) are applied into BeiDou real-time PPP. Observational data from stations that are capable of tracking the BeiDou constellation from the IGS Multi-GNSS Experiments (MGEX) network are processed, with the introduced NWM augmented PPP and the standard PPP processing. The accuracy of tropospheric delays derived from NCEP is assessed against with the IGS final tropospheric delay products, Fig. 3.73. The positioning results show that an improvement of convergence time up to 60.0% and 66.7% for the east and vertical components, respectively, can be achieved with the NWM augmented PPP solution compared to the standard PPP solutions, while only slight improvement of the solution convergence can be found for the north component Fig. 3.74. A positioning accuracy of 5.7 cm and 5.9 cm for the east component is achieved with the standard PPP that estimates gradients and the one that estimates no gradients, respectively, in comparison to 3.5 cm of the NWM augmented PPP, showing an improvement of 38.6% and 40.1%. Compared to the accuracy of 3.7 cm and 4.1 cm for the north component derived from the two standard PPP solutions, the one of the NWM augmented PPP solution is improved to 2.0 cm, by about 45.9% and 51.2%. The Fig. 3.73 (a) The time series of NCEP and IGS ZTD at station BRST for September, 2015 – day of year 244–272. The NCEP ZTD are shown in red, and IGS ZTD in blue. (b) Distribution of ZTD differences between NCEP and IGS 3 Advanced GNSS Processing Techniques (Working Group 1) 145 Fig. 3.74 The BeiDou RT-PPP solutions at station GMSD (Japan, 30.56  N, 131.02 E) on September 1, 2015 (DOY 244 of 2015). The NWM augmented PPP solution is shown in red, the standard PPP solution that estimates gradients in blue, and the standard PPP solution that estimates no gradients in green positioning accuracy for the up component improves from 11.4 cm and 13.2 cm with the two standard PPP solutions to 8.0 cm with the NWM augmented PPP solution, an improvement of 29.8% and 39.4%, respectively. 3.5.8 Using External Tropospheric Corrections to Improve GNSS Positioning of Hot-Air Balloon16 P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] 16 Parts from this section were previously published in Vaclavovic et al. (2017). 146 J. Douša et al. J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] M. Eliaš Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] Tropospheric delay of a GNSS signal has to be considered in any precise positioning solution. When a receiver remain in a static position, which is the case of all permanent stations, coordinates are significantly constrained, and the correlation between station height and ZTD is not critical. Moreover, tropospheric parameters can be also constrained during stable weather condition. However, when a receiver is in movement, particularly in vertical direction, the height-ZTD mutual correlation can significantly decrease an accuracy of position results. It is particularly true for Precise Point Positioning because any biases are neither eliminated nor reduced. It should be noted, that known precise troposphere can improve also static receiver positioning, namely in terms of reducing convergence time due to diminishing state vector dimension and better observation model after cold start. Consequently, introducing ZTD from an external source can play important role in making the solution more robust and accurate. Sources for the tropospheric parameters can be models produced on the bases of meteorological data, tropospheric products from GNSS analyses, or their optimal combination. Geodetic Observatory Pecny has demonstrated how the correlation between height and ZTD critically influences position estimates of kinematic receiver. For this purpose, an experiment with a hot air balloon carrying a GNSS receiver together with several meteorological sensors has been arranged (Fig. 3.75). The receiver with other necessary equipment were mounted on the balloon basket and flied more than 2000 meters above the earth surface (Václavovic et al. 2017). Figure 3.76 demonstrates a vertical profile of the flight divided into five flight phases: (1) Initial, (2) Ascent, (3) Descent, (4) Landing, and (5) Finish. In order to obtain stable solution without initial convergence typical for the PPP method, the receiver remained in its static position during the initial and finish flight phases. Collected multi-GNSS (GPS, GLONASS, Galileo, BeiDou) high rate observations with meteorological measurements (atmospheric pressure, relative humidity, air temperature) were processed by the G-Nut/Geb software using the PPP method, which requires modeling all errors with the highest accuracy. This is accomplished using corrections from external precise models or products, such as satellite orbits and clocks, satellite attitude models, atmospheric delays, receiver and satellite antenna phase centre offsets and variations, relativistic and phase windup effects. Experimental campaign observations were processed utilizing real-time orbits and clocks of GPS and GLONASS satellites provided by the Centre National d’Etudes Spatiales, CNES (Laurichesse et al. 2013). The products are disseminated in a real-time stream, and 3 Advanced GNSS Processing Techniques (Working Group 1) 147 Fig. 3.75 Mounting GNSS antenna and meteorological sensors on the balloon basket (left) and the balloon during ascending (right) Fig. 3.76 Experimental vertical profile and definition of flight phases they had been stored in daily files using the BKG NTRIP Client (BNC) (Weber et al. 2016). We applied two principal strategies for ZTD modeling: (1) estimating along with all other unknowns in the adjustment, (2) introducing from an external tropospheric model. The tropospheric model used in this study was the Weather Research and Forecasting (WRF) model provided by the Institute of Computer Science, Academy of Science of the Czech Republic. A methodology for deriving ZTD from such WRF model has been developed at GOP and published in Douša and Eliaš 2014. Estimating ZTD simultaneously with the rover height resulted in a strong dependency of kinematic solution on the ZTD random walk noise setting. The effect was attributed to a strong mutual correlation described by the correlation coefficient reaching up to 0.9 for a very loose ZTD constraining spectral density 500 mm/ h1/2 of the stochastic process random walk (red line in Fig. 3.77). Depending on ZTD random walk setting, the rover height discrepancies reached up to 30–50 cm. The study suggests that ZTD should be tightly constrained in a vertically kinematic solution to stabilize the results; however the constraining influences the accuracy of estimated height in case of strong rover dynamics. The solution using not only GPS satellites but also GLONASS became more stable due to a better satellite 148 J. Douša et al. Fig. 3.77 Dependence of the correlation coefficient between ZTD and height on random walk settings geometry and a better decorrelation of estimated parameters. To achieve the best accuracy, a careful offline analysis of an optimal random walk setting for estimated ZTD should be done and, optionally, dynamically adapted. Findings of such analyses are not available in real-time, however, when a precise external tropospheric model is applied, the offline analyses is not necessary. The precise tropospheric model improved mainly solutions under the conditions of poor satellite geometry which was simulated by reducing available GPS constellation. Therefore, a significantly smaller dependence of positioning precision on the ZTD constraining was observed when using multi-constellation observation compared to single-constellation. Kinematic processing is usually difficult in real-time because receiver movement cannot be sufficiently predicted; therefore, large noise has to be introduced in the Kalman filter prediction. However, when post-processing is applicable the backward smoothing can be applied (Václavovic and Douša 2015). It the experiment, such improvement was described by the factor of two. The combination of multi-GNSS observations and the backward smoothing algorithm reached the best agreement between different solution variants using different troposphere constraining. As an optimal setting for ZTD constraining is not able to recommend, precise kinematic positioning may benefit from external tropospheric corrections. Current accuracy of NWM forecasts already provides corrections at the centimeter level, which is similar or even better when compared to values estimated from GNSS data, in particular when a standalone GNSS is used. The NWM-driven PPP solution of our vertical experiment resulted in 9–12 cm and 5–6 cm uncertainties in the rover altitude using the Kalman filter and the backward smoothing, respectively. Compared to the standard PPP, it indicates better performance by a factor of 1–2 depending on the availability of GNSS constellations, the troposphere constraining and the processing strategy used. To conclude, if observation conditions are difficult, such as with high GDOP values, external corrections from the augmented tropospheric model can significantly improve the robustness and the accuracy of the rover height in precise positioning of prevailing vertical dynamics. More details can be found in Václavovic et al. 2017. 3 Advanced GNSS Processing Techniques (Working Group 1) 3.5.9 149 Real-Time PPP Augmented with High-Resolution NWM Model Data17 K. Wilgan Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] T. Hadaś Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] P. Hordyniec Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] J. Bosy Wrocław University of Environmental and Life Sciences, Wrocław, Poland e-mail: [email protected] We proposed a high-resolution model of troposphere corrections for Poland based on numerical weather prediction (NWP) model and GNSS data from postprocessing. The chosen NWP model is Weather Research and Forecasting (WRF). In current configuration, the model outputs are given in a form of a dense horizontal grid (4 km  4 km) with 47 height levels. Another factor that can impact the positioning accuracy and convergence time is the choice of mapping functions used to reduce the zenith delay to the slant delay. The most commonly used VMF are based on the ECMWF model with spatial resolution of about 40 km and temporal resolution of 6 h. In this article, the mapping functions based on a WRF model with higher spatial and temporal resolutions were proposed. The tropospheric corrections model was applied into real-time PPP software GNSS-WARP (Wroclaw Algorithms for Real-time Positioning) developed at Institute of Geodesy and Geoinformatics, Wroclaw University of Environmental and Life Sciences. The study was conducted for 14 Polish EPN stations during three periods, with different troposphere conditions. The performance tests were conducted in six GNSS data processing variants, including two commonly used variants using a priori ZTD and mapping functions from UNB3m and VMF1-FC models, one with a priori ZTD and mapping functions calculated directly from WRF model and three variants using the aforementioned mapping functions but with ZTD model based on GNSS and WRF data used as a priori troposphere and to constrain tropospheric estimates. The application of a highresolution WRF/GNSS-based ZTD model and mapping functions results in the best agreement with the official EPN coordinates. Three types of coordinates: static, kinematic and reinitialized kinematic were estimated. In both, static and kinematic 17 Parts from this section were previously published in Wilgan et al. (2017b). 150 J. Douša et al. Fig. 3.78 Mean 3D biases and 3D StDev of static coordinate residuals (estimated – EPN official) averaged from 14 Polish EPN stations mode the application of high-resolution WRF/GNSS-based model resulted in an average reduction of 3D bias by 20 and 10 mm respectively, but an increase of 3D standard deviations by 1.5 and 4 mm respectively (Figs. 3.78 and 3.79). This approach also shortens the convergence time, e.g. for a 10 cm convergence level, by 13% for the horizontal components and by 20% for the vertical component (Fig. 3.80). More details can be found in Wilgan et al. 2017a, b. 3.5.10 Validation and Implementation of Direct Tropospheric Delay Estimation for Precise Real-Time Positioning L. Yang University of Nottingham, Nottingham, UK e-mail: [email protected] 3 Advanced GNSS Processing Techniques (Working Group 1) 151 Fig. 3.79 Mean 3D biases and 3D StDev of kinematic coordinate residuals (estimated – EPN official) averaged from 14 Polish EPN stations C. Hill University of Nottingham, Nottingham, UK e-mail: [email protected] J. Jones Met Office, Exeter, UK e-mail: jonathan.jones@metoffice.gov.uk Unmitigated tropospheric delay remains one of the major error sources in PPP. Due to the lack of real information, conventional empirical models (for both ZTD and mapping functions) have limitations in terms of positioning accuracy and convergence time, especially during active tropospheric conditions. Aiming to overcome these limitations and improve PPP performance, we investigated the feasibility of integrating the real atmospheric condition from NWM into PPP, as an external 152 J. Douša et al. Fig. 3.80 Convergence time for different levels of convergence for horizontal (up) and vertical (bottom) components and a percentage of converged solutions correction service. The main challenge is to capture and abstract the most effective information from NWP data, and transmit to the user in a bandwidth-economic way. The fundamental concept proposed is a multiple-tiered data transmission structure. The tiers will provide quality information in different levels of details. The content of the tiers is shown in Fig. 3.81. In the whole service area, the atmospheric conditions are monitored in grid cells. The higher tiers are only produced and transmitted where/when they could bring additional accuracy compared to the lower tiers. An indication type of flag map is produced in real-time to inform users of the number of active tiers available at the users’ grid cell. This tiered transmission scheme could effectively reduce the transmission bandwidth consumed, while keeping a certain guarantee of accuracy. In this project, the state of the art Unified Model (UM) produced by the UK Met Office (UKMO) is selected for the underpinning NWP. To investigate the 3 Advanced GNSS Processing Techniques (Working Group 1) 153 Fig. 3.81 Tiered troposphere delay correction structure performance of the proposed correction service, statistical analysis is carried out using data over 12 months in 2014 with a 6-h interval at 108 CORS stations in the UK. 3.5.10.1 Impact of the Multiple Tiers In Fig. 3.82, the RMS errors of different ZTD estimation methods are compared, (using DD solution as the reference). Three different empirical ZTD estimations (Tier 0, which does require any external information), i.e. the MOPS model, the GPT2 model and the ESA Blind (ESA-B) mode model, show very similar performance. While clear improvements can be observed for the higher tier solutions. Comparing to the ESA-B solution, there is a 27.3% improvement for the ESA site mode (Tier 1, which provides surface meteorological parameters), 45.2% for the ESA augmented mode (Tier 2, which provides vertical lapse rates of meteorological parameters) and 73.9% for the zenith ray tracing solution (which provides the ultimate on-site ZTD correction). Clearly the more real information provided, the better the ZTD estimation accuracy achieved. Figure 3.83 provides the slant tropospheric delay (STD) error comparison between the solutions of Tier 0 (ESA blind mode), Tier 1 (ESA site mode), Tier 2 (ESA augmented mode), Tier 3 (ZTD and mapping function correction) and Tier 4 (gradient correction). The Tier 3 solutions are further split into the ZTD only correction and the ZTD plus mapping function correction. Ray tracing is used as the reference. The mean error and standard deviation of each solution are plotted in individual figures in the left, and the RMS errors of all solutions are plotted together against elevation in the right. It can be seen that, from the blind mode to the site mode and then the augmented mode, both mean error and standard deviation are improved in general, which results in the improved RMS errors in turn. The improvement from blind mode to site mode 154 J. Douša et al. Fig. 3.82 RMS of ZTD error from different estimation methods Fig. 3.83 TSD RMS error against elevation, for different tiers is larger than the improvement from site mode to augmented mode indicating that, although the real vertical rates of change are beneficial, the real surface parameters play a more important role in the STD estimation. The ZTD correction could effectively reduce the mean error offset at higher elevations, as well as the uncertainty in the lower elevation area. This is because if the ZTD is not estimated properly, the error could be magnified by the mapping function. Even with the ZTD correction, the RMS error at 10 elevation is still higher than 10 cm. Both mapping function and gradient correction could reduce the error in the lower elevation, while the former has a stronger impact. Therefore, it is recommended that the mapping function shall be corrected earlier than the gradient. It can also be seen that at the 10 elevation, the improvement from the mapping function correction 3 Advanced GNSS Processing Techniques (Working Group 1) 155 will be larger than the improvement from the ZTD correction, and could reduce the RMS error to centimetre level. 3.5.10.2 Flag Map One of the key features in our proposed tiered troposphere delay correction scheme is the flag map. This map is updated in real time, and provides users with local information, that to which correction tier the user should use to guarantee a certain STD estimation accuracy. The size of this flag map is quite small, thus is suitable for broadcast in a transmission medium with limited bandwidth, such as a communication satellite. Figure 3.84 shows an example of the tiered structure flag map, which is consist of 60 by 90 grid squares (~17 km  17 km) for the UK region. Figure 3.85 shows the total occurrence probability of each tier, from the statistic of 12 months in the UK CORS stations. Each subfigure presents different elevation settings, ranging from 10 to 70 . In each subfigure, the x axis is various thresholds that could be pre-set as the accuracy targets, and the y axis gives the chance (in percentage) of each correction tier being required. At high elevations, a simple empirical estimation approach would be able to meet a lenient accuracy target. With an increasingly strict threshold, real information is required to replace empirical values to suppress the modelling error. When the elevation gets lower, the STD modelling assumptions based on the empirical mapping function and balanced local atmosphere profile will gradually become less valid, and the error they bring cannot be ignored. Fig. 3.84 Example of the flag map for the tiered structure (06th Oct 2014 0600 GMT) 156 J. Douša et al. Fig. 3.85 The chance of each tier to be applied, against different threshold 3.5.10.3 PPP Improvement Figure 3.86 uses four storm event cases (at each row) to indicate the PPP performance improvement from the proposed correction scheme, in which the red lines indicate the traditional PPP solution using empirical tropospheric modelling and the green lines indicate PPP with external tropospheric correction. In the left column, the convergence times for 10 cm (dash line) and 5 cm (solid line) are shown. And in the right column, the positioning standard deviation in 3D (solid line) and height component (dash line), at the end epoch after the 1-h PPP solution, are shown. The positive impact of the proposed external correction can be observed from both aspects. 3.5.10.4 Bandwidth Bandwidth is a key issue in for the proposed troposphere correction scheme, since if the bandwidth required is too high to be transmitted, the whole scheme is not practical. A customized transmission has been designed in this project, which only transmits the necessary information to where it is required and hence effectively reduce the total required bandwidth and also help in multiplexing the transmission channel. 3 Advanced GNSS Processing Techniques (Working Group 1) 157 Fig. 3.86 Comparison of PPP convergence time and standard deviation after 1 h Table 3.21 Summary of transmission design details for each tier, for UK coverage Tier Bandwidth (kbps) 1 4.74 2 2.60 3 0~5.21 4 0~2.60 5 397~531 In Table 3.21 the required transmission bandwidth is estimated. It can be seen that, for the UK coverage, the total required bandwidth for tier 1–4 will be 7.34~15.15 kbps, and 397~529 kbps for Tier 5. Meanwhile, the size of each flag map message will be 110 Bytes, which could be easily broadcasted. Tier 1–4 and the flag map message could be easily fit into communication satellite bandwidth nowadays. Tier 5 information is better to be transmitted via 4/5G network and through a two-way communication once required. 3.5.10.5 Summary • Real time ray tracing could provide an accurate STD estimation, and thus improve PPP results. • Although the required troposphere condition data for the real time 3D full ray tracing is huge, a tiered correction structure and flag maps mechanism has been 158 J. Douša et al. designed to reduce data transmission volume effectively. The data are only transmitted for the time and place when they are required. • Statistical results show that only very occasionally is the full ray tracing correction required, and for most of the time, the mapping function correction could allow the user to satisfy STD estimation accuracy. • The total bandwidth required for the UK coverage is 7.34~15.15 kbps for Tier 1–4, and 397~529 kbps for Tier 5. The flag map message size is 110B each. Apart from Tier 5, other information could be easily transmitted via communication satellite. • The current drawback of the whole proposed solution is the NWP data resolution, especially in time (currently 6-h interval). 3.6 GNSS Data Reprocessing for Climate One of the important tasks of the WG1 of the COST Action was the homogeneous reprocessing of GNSS data for climate applications in support of WG3. Main focus of WG1 was to study the impact of precise models and processing strategy for providing a long-term tropospheric parameters with an optimal accuracy, a high reliability and a consistent quality. This section gives an overview of the reprocessing activities during the COST Action period, while WG3-related activities are described in Chap. 5, including reprocessing product validation, combination, quality control and homogenization. This section introduces contributions to the GNSS second reprocessing performed within the IAG sub-commission for the European Reference Frames (EUREF) and within the IGS TIGA campaign. Although these were primarily designed for contribution to geodetic applications, within the COST Action a significant effort was additionally dedicated to an optimal estimation of tropospheric parameters on regional and global scales. A special activity supporting directly the climate community has been then started by establishing the GRUAN Central Processing Centre at GFZ for a long-term ground-based GNSS precipitable water estimation. 3.6.1 EUREF Repro2 Contribution of Swisstopo E. Brockmann Swiss Federal Office of Topography swisstopo, Wabern, Switzerland e-mail: [email protected] 3 Advanced GNSS Processing Techniques (Working Group 1) 159 Fig. 3.87 Key parameters of swisstopo repro2 The Swisstopo finished a reprocessing in 2014 covering a homogeneous processing of a time span starting beginning 1996 till almost the end of 2014. Some key parameters are given in Fig. 3.87. The number of stations increased from 20 to 170. The stations are mostly located in central Europe. Some boundary stations were included to enable a better decorrelation between troposphere parameters with station height. The number of satellites increased from 25 in 1996 to 65 in 2014. Since 2004, the number of satellites increased due to the improved GLONASS constellation. The reprocessing was performed with the BSW52 using most up-to-date models in a homogeneous way. Orbit and earth rotation parameters were used from the CODE repro2 products. The important processing options are given in Table 3.22. Figure 3.88 shows that daily solutions are even more stable in view of the Helmert transformation parameters than the previous weekly solutions. Overlapping 3-day solutions are calculated to optimize the ZTD estimates at midnight (see Fig. 3.89). Moreover, variations on the processing options were carried out in order to find the best possible modelling options and in order to do some sensitivity studies: • GMF/GPT mapping, atmosphere loading (ATL) + IGS08 group antenna PCV (submitted as LP0 to EUREF: Oct. 23. 2014) • VMF mapping, IGS08 antenna model but with individual antenna PCV where available, atmospheric tidal loading ATL, non-tidal atmospheric pressure loading ALP (submitted to EUREF as LP1: March 20, 2015). • Same as above, but GPS-only (idea borne on Thessaloniki workshop 2014, results generated in June 2015 and presented at Wroclaw COST meeting end of September 2015) 160 J. Douša et al. Table 3.22 Basic processing options used for swisstopo repro2 Software Satellite systems Elevation cutoff angle Observation weighting Antenna Troposphere Troposphere gradients Tides Conventions Ocean tides Gravity field Ionosphere Reference frame Network Time span Orbits/EOP BSW52 (+) GPS + GLO (ab 2004) 3 COSZ elevation-dependent weighting I08 absolute antenna model (group values) GMF and DRY GMF mapping for the a priori values and while estimating hourly ZPD parameters using WET GMF Chen Herring for tropospheric gradient estimation Atmospheric tidal loading applied IERS2010 FES2004 EGM08 CODE 2-h resolution; including higher order terms IGb08 Max. 180 stations DOY 007, 1996 till DOY 207, 2014 CODE reprocessing series 2011 (till DOY 106, 2011) and CODE reprocessing series 2013 (till DOY 362, 2013), CODE operational series in 2014 • The impact of additional GLONASS observations on the long-term (comparison of the before mentioned last two repro2 series) was especially analyzed in this project. The impact of the additional GLONASS observations is negligible when analyzing ZTD trends. Further conclusions can be drawn from the comparisons between different ZTD time series: • Formal errors of ZTD estimates are smaller with more GLONASS observations (max. 10%) included • Influence of GMF versus VMF: no significant rate, standard deviation 0.5–2.5 mm ZTD (109 sites with long time series) • Difference to CODE global Repro: no significant rate, standard deviation 2.0–5.0 mm ZTD (39 sites with long time series); • Influence additional GLONASS observations: no significant rate, standard deviation: 0.4–1.5 mm (111 sites) (MODA exception). This is In fact a little higher (also GPS-only time counted), but strongly dependant on the used a priori ZTD constraints, there is only a statistical effect but no significant bias. Comparisons with e.g. radio sondes have much bigger differences (see Fig. 3.90). In case of PAYE we see a standard deviation of 8–9 mm ZTD. The estimated rate of 5 mm/10 years. is also dependent on changes at the radio sonde station. In 2009, a 3 Advanced GNSS Processing Techniques (Working Group 1) 161 Fig. 3.88 Stability of the solutions expressed in Helmert Parameters between each individual solution with respect to the combined solution (operational older weekly solution top, swisstopo repro2 daily solution bottom) new humidity sensor was used. Trend estimates on the raw ZTD time series are quite sensitive – especially the various antenna changes at the GNSS stations generate quite significant jumps in station heights as well as in the ZTDs. 162 J. Douša et al. Fig. 3.89 3-day solutions to optimize the ZTD estimates at midnight Zenith Total Delay Difference for station PAYE in [mm] PAYE: R2H3 minus RS, N=7377 linear fit, std=8.3 mm 40 20 0 –20 Repro2, VMF, 3D: 8.3 mm std –40 2004 2006 2008 Radiosonde: new humidity sensor 31/03/2015 15:00 40 2010 YEAR 2012 2014 2016 PAYE: R2G1 minus RS, N=7377 linear fit, std=8.9 mm 20 0 –20 Repro2, GMF, 1D: 8.9 mm std –40 2004 2006 2008 2010 YEAR 2012 2014 2016 31/03/2015 14:50 Fig. 3.90 Swisstopo Repro2 (VMF 3 days upper diagram, GMF 1 day lower diagram) compared with radio sonde derived ZTD estimates for station PAYE 3 Advanced GNSS Processing Techniques (Working Group 1) 3.6.2 163 EUREF Repro2 Assessment of GOP Processing Variants18 J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] P. Václavovic Geodetic Observatory Pecný, Research Institute of Geodesy, Topography and Cartography, Zdiby, Czech Republic e-mail: [email protected] The results of the 2nd reprocessing of all data available from 1996 to 2014 from all stations (>300) of the European GNSS permanent network are introduced as performed at the Geodetic Observatory Pecný (GOP) (Douša et al. 2017). The reanalysis was completed during the 2nd EUREF reprocessing to support mainly the realization of a new European terrestrial reference system. The BSW52 (Dach et al. 2015) was used together with CODE precise products (Dach et al. 2014). Within the COST ES1206 project, a new set of GNSS tropospheric parameter time series was provided for applications in climate research. To achieve this goal, we improved our strategy for combining tropospheric parameters over three consecutive days that guarantees a continuity of all estimated tropospheric parameters, zenith tropospheric delays (ZTD) and tropospheric horizontal linear gradients, at all mid-nights and during transitions of GPS weeks. Basic characteristics of the 2nd reprocessing are provided in Table 3.23. Within the reprocessing, we performed additionally seven solution variants in order to study an optimal troposphere modelling, Table 3.24. We assessed variants in terms of the coordinate repeatability by using internal evalautions of coordinate repatability. Then, we compared ZTDs and tropospheric horizontal gradients with independent values obtained from ERA-Interim numerical weather reanalysis (Dee et al. 2011). Generally, the results of the GOP Repro2 yielded improvements of approximately 50% and 25% for the horizontal and vertical component repeatability, respectively, when compared to the results of the GOP Repro1 solution. Vertical repeatability was reduced from 4.14 mm to 3.73 mm when using the VMF1 mapping function (Boehm et al. 2006a, b), a priori ZHD, and non-tidal atmospheric loading corrections from actual weather data. Increasing the elevation cut-off angle from 3 to 7 /10 increased RMS errors of residuals from the coordinates’ repeatability. These findings were confirmed also by the independent assessment of tropospheric parameters using NWM reanalysis data. The differences of particular solutions were statistically analysed to demonstrate the impact of the different modelling on ZTDs, see 18 Parts from this section were previously published in Douša et al. (2017). 164 J. Douša et al. Table 3.23 Characteristics of GOP reprocessing models Analysis options Products Observations Reference frame Antenna model Troposphere Ionosphere Loading effects Description CODE precise orbit and earth rotation parameters from the IGS 2nd reprocessing. Dual-frequency code and phase GPS observations from L1 and L2 carriers. Elevation cut-off angle 3 , elevation-dependent weighting 1/cos2 (zenith), double-difference observations and with 3-min sampling rate. IGb08 realization, core stations set as fiducial after a consistency checking. Coordinates estimated using a minimum constraint. GOP: IGS08_1832 model (receiver and satellite phase Centre offsets and variations). A priori zenith hydrostatic delay/mapping function: GPT/GMFh (GO0) and VMF1/VMF1 h (GO1-GO6). Estimated ZWD corrections every hour using VMF1w mapping function; 5 m and 1 m for absolute and relative constraints, respectively. Estimated horizontal NS and EW tropospheric gradients every 6 h (GO0-GO5) or 24 h (GO6) without a priori tropospheric gradients and constraints. Eliminated using ionosphere-free linear combination (GO0-GO6). Applying higher-order effects estimated using CODE global ionosphere product (GO5). Atmospheric tidal loading and hydrology loading not applied. Ocean tidal loading FES2004 used. Non-tidal atmospheric loading introduced in advanced variants from the model from TU-Vienna (GO4-GO6). Table 3.24 GOP solution variants for the assessment of selected models and settings Solution ID GO0 GO1 GO2 GO3 GO4 GO5 GO6 Specific settings and differences GMF and 3 cut-off VMF1 and 3 cut-off ¼GO1; 7 cut-off ¼GO1; 10 cut-off ¼GO1; atmospheric loading ¼GO4; higher-order ionosphere ¼GO4; 24-h gradients Remarks and rationales Legacy solution for Repro1 New candidate for Repro2 Impact of elevation cut-off angle Impact of elevation cut-off angle Non-tidal atmospheric loading applied Higher-order ionosphere effect not applied Stacking tropospheric gradients to 24-h sampling Table 3.25. Although the VMF1 mapping function outperformed GMF in term of standard deviation for the elevation cut-off angle of 3 (GO0 vs GO1), the impact is smaller than for the others effects, i.e. due to elevation cut-off angle (GO2 vs GO1, GO3 vs GO1) and non-tidal atmospheric modelling (GO4 vs GO1). The results get worse for both ZTD and tropospheric gradients when raising the elevation cut-off angle from 3 to 7 (GO2) or 10 (GO3). No significant impact of modelling of highorder ionospheric effects (GO5) was observed; the effect is systematic for regional network and was eliminated by applying fiducial stations. Figure 3.91 shows a significant improvement for gradients (negligible impact on ZTD) when combining 3 Advanced GNSS Processing Techniques (Working Group 1) 165 Table 3.25 Median, minimum and maximum values of total ZTD biases and standard deviation over all stations Compared solution variants GO1 vs GO0 GO2 vs GO1 GO3 vs GO1 GO4 vs GO1 GO5 vs GO4 GO6 vs GO4 ZTD bias median [mm] 0.36 +0.03 +0.03 +0.05 0.02 0.02 ZTD bias min [mm] 1.52 0.81 2.22 3.29 0.31 0.23 ZTD bias max [mm] +0.70 +1.66 +2.66 +5.55 +0.07 +0.16 ZTD sdev median [mm] 2.01 0.66 1.10 1.37 0.07 1.24 ZTD sdev min [mm] 0.69 0.15 0.31 0.68 0.04 0.76 ZTD sdev max [mm] 3.82 1.29 2.04 4.72 0.30 2.46 Fig. 3.91 Monthly means of bias and standard deviation of tropospheric horizontal north (N-GRD) and east (E-GRD) gradients compared to those obtained by ERA-Interim. Error bars indicate standard errors of mean values over all compared stations plotted from the zero y-axis to emphasise seasonal variations and trends. Error bars are displayed for north gradients only, however, being representative for the east gradients too gradient parameters from the original 6-h into 24-h time resolution without applying constrains. The 6-h gradients were able to absorb some remaining errors in GNSS model. We also studied a temporal and spatial variation of ZTD differences from seven solution variants, for details see Douša et al. 2017. Figure 3.92 shows latitudinal 166 J. Douša et al. Fig. 3.92 Dependence of ZTD biases (blue) and standard deviations (red) from inter-comparisons of GOP 2nd reprocessing solution variants on station latitude. Note different y-range for the GO5 vs. GO4 comparison dependence of the ZTD differences, while dependences on height, time and geographical location are not showed here. Based on the reanalysis variants, we recommend using low-elevation observations together with the VMF1 mapping function, precise a priori ZHD values and consistent model of non-tidal atmospheric loading. The daily piecewise linear function model showed better stability for estimated gradients while did not indicate a worse repeatability of coordinates estimates. For saving the time, we could recommend such (unconstrained) approach for an optimal modelling of the firstorder tropospheric asymmetry. The main difficulties faced during the 2nd reprocessing was the quality of the GNSS historical data containing a large variety of problems. To provide high-accuracy GNSS tropospheric products, the elimination of problematic data was more critical compared to estimating coordinates on a daily only. The 3rd EUREF reprocessing should thus exploit more complex data quality control in order to optimize the reprocessing and the quality of tropospheric products. 3 Advanced GNSS Processing Techniques (Working Group 1) 3.6.3 167 CORDEX.be Reprocessing E. Pottiaux Royal Observatory of Belgium, Brussels, Belgium e-mail: [email protected] B. Van Schaeybroeck Royal Meteorological Institute of Belgium, Uccle, Belgium e-mail: [email protected] The CORDEX.be stands for “COordinated Regional Climate Downscaling EXperiment and beyond” (http://cordex.meteo.be/). This national project brings together the Belgian climate and impact modelling research groups into one research network as the first step towards the realization of climate services. The philosophy of CORDEX.be is inspired by the World Climate Research Programme (WCRP) project CORDEX (“COordinated Regional Climate Downscaling Experiment”, http://www.cordex.org/) project, but – as the “.be” in the acronym indicates – the project aims to go beyond for Belgium. CORDEX.be has 4 main targets (see http://cordex.meteo.be/meteo/view/en/ 29026726-Targets.html): (1) to contribute to the international CORDEX project, (2) to go beyond by running 4 High-Resolution Limited Area Models (ALARO-0, MAR, and two flavours of COSMO-CLM, all with a resolution of 4 km, but with different forcing strategies), (3) to go beyond by running 4 Local-Impact Models (Wave and storm surge, urban, crop Isoprene or vegetation emission models), and (4) to infer the climate uncertainties to the Belgian level. In addition, a specific task is dedicated to the validation of the high-resolution climate simulations using GNSSderived products. The main objective of this task is to go beyond the standard verification procedure of climate simulations. The traditional manner is to compare the results from climate runs with long-term ground-based surface meteorological observations. Instead, here we will implement a verification based on products estimated from continuously operating GNSS (Global Navigation Satellite Systems, such as the American GPS) stations. Therefore, a careful reprocessing of the historical observations at about 320 world-wide GNSS stations (Fig. 3.93, left) has been carried out with a focus on the period 2000–2010 (i.e. the period requested by the CORDEX.be partners for the assessment of the climate models), allowing both the validation of the standard CORDEX runs and the high-resolution runs. Unfortunately, from these 320, only 20 GNSS stations (Fig. 3.93, right) which were operating almost continuously and providing high-quality observations throughout the requested validation period (2000–2010) are located in the restricted high-resolution climate model domain, but many other European GNSS sites can still be use to validate the standard 12.5 and 50 km resolution European CORDEX runs. First results of this GNSS-based validation of the hourly ZTDs derived from the high-resolution model output is given in Fig. 3.94 (only for two climate models, MAR and COSMO-CLM, but with different forcing scenarios for MAR: NCEP, 168 J. Douša et al. Fig. 3.93 Left: All GNSS stations included in the reprocessing. Right: GNSS stations located in the high-resolution climate domains of CORDEX.be COSMO-CLM KUL (ERA-I) MAR (ERA-I forcing) Model / observed variability ZTD Anomaly correlation of hourly ZTD MAR (NCEP forcing) MAR (ERA20C forcing) 90% 80% 70% 60% 50% 40% 30% 20% Feb Apr Jun MONTH Aug Oct Dec 110% 105% 100% 95% 5 10 15 20 Hour of the day Fig. 3.94 Anomaly correlation (left) and variability ratio (right) of hourly ZTD values between different model simulations and the GNSS-derived observations for the different months of the year. The models considered are the H-Res simulations of COSMO-CLM (driven by ERA-Interim), and MAR (driven by NCEP-NCAR-v1, ERA-Interim and ERA20C). Averages are taken over seven Belgian GNSS stations ERA20C and ERA-Interim). It is seen that all model-based ZTDs correlate pretty well with the GNSS-based ZTDs. The annual cycle (Fig. 3.94, left) from all models shows a better agreement during the winter months while more pronounced departure is visible during the summer months. Very similar yearly patterns in the anomaly correlation is visible for the two climate models and the various forcing scenarios represented. For the daily cycle (Fig. 3.94, right), also very similar daily patterns in the “model over observed variability of the ZTD” is visible in all cases, with a peak drop in the mid-late afternoon. The current working hypothesis is that the higher the water vapour variability, the higher the departure between model and observation is visible. In terms of forcing the climate model runs: as expected for the 3 Advanced GNSS Processing Techniques (Working Group 1) 169 MAR results, results with ERA-Interim forcing are better than forcing with NCEPNCAR-v1 reanalysis which in turn are better ERA20C forcing. The better scores for MAR with respect to COSMO-CLM may be caused by contrasting coupling strategies. More specifically, whereas MAR is directly forced at the boundaries by reanalysis, COSMO-CLM uses an additional nesting to obtain the highest spatial resolution of 2.8 km. This validation work is still ongoing and more detailed assessments will be undertaken in the future. 3.6.4 GRUAN Reprocessing G. Dick GFZ German Research Centre for Geosciences, Helmholtz Centre Potsdam, Potsdam, Germany e-mail: [email protected] F. Alshawaf GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] The Global Climate Observing System (GCOS) Reference Upper Air Network (GRUAN) of the World Meteorological Organization (WMO) is an international reference observing network, designed to meet requirements of climate research. Upper air observations within GRUAN will provide long-term high-quality climate records. A GNSS receiver is part of the GRUAN station equipment with highest priority for deriving PWV. Due to its long-term experience in GNSS data processing, GFZ was selected by WMO as a Central GRUAN GNSS Data Processing Centre. GFZ operates also a number of GNSS stations on GRUAN sites, see Fig. 3.95. Fig. 3.95 GRUAN network with GNSS stations on GRUAN sites operated by GFZ (in red) 170 J. Douša et al. Fig. 3.96 GNSS-derived PWV results from reprocessed data for 2001–2016 compared with radiosonde measurements at GRUAN site Lindenberg (Germany) The GRUAN GNSS Processing Centre at GFZ with EPOS software covers a fully automated processing chain starting with collecting of raw GNSS observations and resulting with climate relevant validated PWV products, which are available online. PWV uncertainty estimation (Ning et al. 2016), comparisons with radiosonde (RS) measurements (Figs. 3.96 and 3.97) as well as PWV trend estimation (Fig. 3.98) are essential parts of the GNSS data analysis at GFZ for climatological applications. 3.6.5 GFZ TIGA Reprocessing Z. Deng GFZ German Research Centre for Geosciences, Potsdam, Germany e-mail: [email protected] Being a modern geodetic measuring method, GNSS (Global Navigation Satellite Systems) has reached an important role in geosciences. Within the scope of the Tide Gauge Benchmark Monitoring Working Group (TIGA-WG) of the IGS, GFZ analyses and reprocesses GNSS data of stations near tide gauges (Deng et al. 2014). This allows us to monitor tide gauges for vertical land deformations, e.g. due to postglacial uplift. TIGA also contributes to the calibration of satellite altimeters and the unification of height systems (Hunegnaw et al. 2017). In conjunction with circa 400 global IGS stations, GFZ processes data of almost 500 GNSS stations near tide gauges between 1994 and today (Fig. 3.99). In the first TIGA combination there are contributions from 3 international ACs. The GFZ TIGA solution shows the best accuracy among the three submitted solutions (Hunegnaw et al. 2017). 3 Advanced GNSS Processing Techniques (Working Group 1) 171 Fig. 3.97 GNSS-derived PWV results from reprocessed data for 2011–2016 compared seasonally with radiosonde measurements at GRUAN site Ny-Alesund (Norway). Seasonal differences between GNSS PWV and RS can be explained by different behaviour of RS during warmer and colder months Fig. 3.98 Example of PWV trend estimation from reprocessed GNSS data for 2000–2016 at GRUAN site Lindenberg. Trend value is 0.31 mm/decade, sigma of the trend is 0.075 mm/decade 172 J. Douša et al. -150˚ -100˚ -50˚ 0˚ 50˚ 100˚ 150˚ 50˚ 50˚ 0˚ 0˚ -50˚ -50˚ -150˚ -100˚ -50˚ 0˚ 50˚ 100˚ 150˚ Fig. 3.99 GNSS stations processed at GFZ in the framework of TIGA Since the GFZ TIGA processing complied with the accords of the 2nd IGS reprocessing campaign, our reprocessed solutions contributed to the determination of the ITRF2014 (Rebischung et al. 2016). In that context the GFZ solution contained the most stations among all of the nine solutions. Because it is the only solution that contains all IGS stations, the IGS Analysis Centre Working Group decided at the 2017 IGS Workshop in Paris, to routinely deliver the GFZ TIGA solution to the IGS in order to ensure, that all IGS stations are contained in the combined weekly IGS solutions. GFZ TIGA products are available via FTP servers of GFZ and CDDIS. The GFZ FTP server provides daily and weekly files for coordinates, orbits and Earth rotation parameters (ftp://ftp.gfz-potsdam.de/pub/transfer/kg_igs/igstiga/solutions/). In addition to coordinate and orbit products, troposphere parameters are provided, which can be used for climate studies (ftp://ftp.gfz-potsdam.de/GNSS/products/tiga_ repro2_tro/). 3.6.6 ULX TIGA Reprocessing19 F. N. Teferle University of Luxembourg, Luxembourg, Luxembourg e-mail: [email protected] ULX, as one of the IGS Tide Gauge Benchmark Monitoring (TIGA) Working Group AC has carried out a second reprocessing campaign in line with IGS. Using the latest 19 Parts from this section were previously published in Hunegnaw et al. (2015). 3 Advanced GNSS Processing Techniques (Working Group 1) 173 Fig. 3.100 120 selected global stations from the reprocessed TIGA solution. The selection is based on the time length of the ZTD time series and their quality. We have only selected those sites having a minimum length of 6 years available bias models and methodology the different IGS ACs re-analyzed the full history of GPS data collected by the global tracking network from 1995 to 2015. The consortium of the British Isles continuous GNSS Facility (BIGF) and the University of Luxembourg TIGA Analysis Centre (BLT) completed a new global solution using <750 GPS stations. Figure 3.100 shows a map of 120 stations. As it can be seen, the stations are globally distributed and the timeseries varies from 6 to 21 years in length. The re-processing follows a double difference network strategy using the BSW52 (Dach et al. 2015), incorporates recent bias model developments, the latest IERS 2010 conventions (Petit and Luzum 2010) and IGS recommendations. Further details are detailed in (Hunegnaw et al. 2015). The selected station network included all IGb08 core stations (Rebischung et al. 2012) and more or less the complete archive of TIGA, which encompasses a large number of GPS stations at or near the global network of tide gauges. The GPS data was re-processed using the CODE final precise orbits and Earth orientation parameters. We employed the IGS08 satellites and receiver antenna phase centre models and adopted an elevation cut-off angle of 3 (Dach et al. 2016). In our solution we make use of the VMF1 (Boehm et al. 2006a, b) that allows to describe the atmosphere with the finest detail, leading to the highest precision in the derived tropospheric parameters. In BSW52, the ZHD is parameterized as a piece-wise function variation of the delay using a piece-wise linear interpolation between temporal nodes. Observations of atmospheric pressure at the GPS station offer high precision for the ZHD estimates and minimize station height errors (Tregoning and Herring 2006). However, many of the TIGA and IGS stations do not possess integrated meteorological sensors. Thus, ZHD in units of meters was a priori obtained reliably from surface pressure data from the gridded output of the ECMWF NWP model and is provided 174 J. Douša et al. by VMF1 using the modified Saastamoinen model, which assumes that the atmosphere is in hydrostatic equilibrium (Davis et al. 1985). We estimate the ZTD parameters in an interval of 1 h with a loose constraint of 5 m. In addition, horizontal gradients in the North-South and East-West directions are estimated in a 24-h interval with the same 5 m loose relative constraint. In this manner more than two decades of ZTD time series along with station positions are available from our re-processing. Figure 3.100 shows a selection of 120 global stations for which we have carried out the further analysis described in this study. However, as the station positions are affected by on average two discontinuities per station per decade, the ZTD time series need to be homogenized before being useful for further application. The results for these 120 stations in terms of a statistical analysis of the periodic signals and stochastic porperties of the related ZWD time series can be found in Klos et al. (2018). 3.7 New Analysis Centres, Networks and Solutions In this section, the new ACs are mainly presented, however, including also new networks and strategies as well as the shared system for facilitating a collaboration between existing and new ACs. Some of the new ACs are located in Eastern and South-Eastern Europe (Bulgaria, Greece, Turkey, Romania), i.e. large areas where no operational products existed so far, others have been established in countries like Austria, Iceland and Portugal. All of them has developed or gained the expertise thanks to the COST Action. Creating of new GNSS ACs thus fulfilled one of the very important goals of the Action: to increase the observing network and to facilitate the transfer of knowledge for establishment of new GNSS ACs. Now, almost all of the new ACs contribute to E-GVAP – the EUMETNET EIG GNSS Water Vapour Programme (http://egvap.dmi.dk). 3.7.1 Trop-NET System for Collaborative Ground-Based GNSS Meteorology J. Douša Geodetic Observatory Pecný, RIGTC, Ondřejov, Czech Republic e-mail: [email protected] The Trop-NET system has been developed at the Geodetic Observatory Pecny (GOP), Czech Republic, in support of the GNSS4SWEC transfer of knowledge. The goal was to facilitate establishment of new analysis centres for near real-time troposphere monitoring in support of numerical weather prediction within the EUMETNET EIG GNSS Water Vapour Programme – E-GVAP (http://egvap.dmi. 3 Advanced GNSS Processing Techniques (Working Group 1) 175 dk). Three short-term scientific missions has been carried out within the transfer of knowledge, see the STSM summary in Appendices. The ground-based GNSS near real-time processing using a batch approach requires following aspects: (1) hourly data provision, (2) predicted orbit products and precise models, (3) efficient and robust procedure for fully automated operation, and (4) a continuous monitoring and product evaluations. Since 1997, GOP has been developing a flexible system for automated GNSS processing using the BSW and BPE. The system also included data/product flow and supports scientific applications for estimating various parameters in a flexible update rate for different purposes/services: • • • • • • • near real-time GPS regional troposphere monitoring (EGVAP) near real-time GPS + GLONASS regional troposphere monitoring (EGVAP) near real-time GPS global troposphere monitoring (EGVAP) ultra-rapid GPS and GLONASS orbit determination (IGS) rapid (daily), final (weekly) GPS solution for reference frame (EUREF) homogeneous re-processing of full EPN (EUREF) daily/hourly data flow at the local data centre (EUREF) All these applications are based on a common module library continuously extended for a higher flexibility and robustness with a near real-time GNSS analysis for the troposphere monitoring as the earliest application (Douša 2001a, b). Due to initial limits in hardware and software, frequent instabilities of data flow, low quality of 24-h predicted orbit products, the system designed had to be highly efficient, fully self-supporting and maximally robust. During the decade, enhancements were done particularly by extending the library for ultra-rapid orbit determination (Douša 2004a, b, 2010), GLONASS processing (Douša 2012), global NRT troposphere solution (Douša and Bennitt 2013) and long-term coordinate and troposphere re-analysis (Douša et al. 2017). The above mentioned experience leaded to the idea of sharing the library for a collaborative use within the GNSS4SWEC project. For this purpose, the Trop-NET package has been completed for easier dissemination, configuration and maintenance in different environments. Main goals within the COST ES1206 were: (a) facilitate the establishment of new analysis centres, (b) improve the product coverage and its homogeneity over Europe, (c) give a possibility to share future developments, and (d) enable coordinated solution updates. Currently, the TropNET pack is maintained by GOP using the Subversion repository. The system consists of several modules supporting individual settings for different user scenarios (Fig. 3.101). The distributed processing is supported by three core modules: (1) for data and product download and mirroring, (2) for GNSS data processing, and (3) for product uploading. Additionally, the system includes central components currently maintained by GOP such as software distribution and update system, NRT product monitoring system, long-term product evaluation system and information systems. Additional modules are considered for future development, for example conversion of ZTD into IWV, animation plots or local monitoring. 176 J. Douša et al. Fig. 3.101 Trop-Net modules – central (up) and disseminated (bottom) 3.7.1.1 Strategy for NRT Troposphere Monitoring Several aspects are considered generally as important for developing the near realtime GNSS troposphere estimates: (a) high efficiency and low latency of GNSS processing, (b) precise station coordinates fully consistent with troposphere estimates, and (c) robust system operated with a minimum manual interventions. For this purpose, the Trop-NET system implemented three processing levels with intermediate solutions combined into a final solution, however, still efficient in near realtime fashion (Douša 2004a, b), Fig. 3.102: 1. Processing of small network clusters using a short-term data batch (yellow and green). 2. Stacking of clusters in spatial domain into a single session network solution (grey). 3. Stacking of network solutions in long-term solutions using previous solutions (blue). Originally, the Trop-NET system applied the 1-h session because of a limited computer power in early 2000. In order to support global solution and for a reliable integer ambiguity resolution at long baselines, the processing batch has been extended to the 4-h session (Douša and Bennitt 2013). Hence, the level of processing the original GNSS data remains redundant for 3 h as visible in Fig. 3.103. The session network and sub-network solutions are temporarily saved in the form of normal equations as intermediate products archived for next 30 days at least. Various numbers of intermediate solutions are combined when estimating different parameters – tropospheric path delays, receiver coordinates or resolving integer 3 Advanced GNSS Processing Techniques (Working Group 1) 177 Fig. 3.102 Trop-NET processing in network clusters, spatial and temporal stacking solutions Fig. 3.103 Trop-NET processing redundancy phase ambiguities. Initially, all parameters are always estimated step-by-step using the 4-h session when original GNSS observations are processed and precise satellite positions introduced as known. First coordinates and tropospheric parameters are then estimated from previous solutions of 1–2 days only before these are introduced for resolving integer phase ambiguities. In later steps, the coordinates are estimated from the time span of 28 days along with introducing integer phase ambiguities and still estimating unresolved ambiguities as float values. Such coordinates are estimated within each individual hourly NRT solution and tied to the actual reference datum. Such approach keeps the system free from any external process and provides additional advantages: (a) all models being implicitly consistent for coordinate and troposphere estimates, (b) coordinates are automatically updated within the system 178 J. Douša et al. without need of an external information about their changes in time (implicitly supports solutions in tectonically active areas such as Greece, Turkey and Iceland), (c) any new station is configured once only with an implicit initialization of station coordinates which is usually done from two past days. The tropospheric parameters are finally estimated using the 12-h session when combining normal equations in temporal domain. The ambiguities are estimated along with tropospheric parameters in the final step because it guarantees more stable tropospheric product compared to the ambiguities-fixed solution when using a short data session only. Different parallelization strategies are used in various processing steps. The ‘predefinition’ regional clusters are used whenever necessarily applying a full correlation model for the sub-network solution and such solutions are usually used also for storing solution normal equations. On the other hand, adaptable strategy for parallel processing uses optimal groups of sites or baselines, often suitable for autonomous site or independent baseline processing. The examples are RINEX conversions, pseudorange smoothing, receiver clock synchronization, ambiguity resolution etc. The processing of regional clusters based on pre-defined configurations may still be automatically adapted, e.g. merging too small clusters. The adaptation of clusters uses several scenarios for generating optimal groups: (a) sorting station- or baselinespecific observation files, (b) using actually available stations, or (c) following clusters from any previous step of the processing. The processing system finally provides automatic warning and error messages either via e-mail or via SMS indicating a temporary solution problem. A warning message often informs about a temporary exclusion of station due to the incompatibility of the file header and the station metadata which always requires a station manual reconfiguration. An error indicating a solution crash often represents a temporary lack of data/products which is possibly within upcoming hours, or if caused by a system-specific reason might need a manual intervention. The status of processing solutions and monitoring indicators are archived along with the products. 3.7.2 Sofia University GNSS Analysis Centre (SUGAC): First Processing Campaign 3.7.2.1 Motivation T. Simeonov Sofia University “St. Kliment Ohridski”, Sofia, Bulgaria G. Guerova Physics Faculty, Department of Meteorology and Geophysics, Sofia University “St. Kliment Ohridski”, Sofia, Bulgaria e-mail: [email protected]fia.bg In Europe GNSS meteorology is a well-established field in both research and operation, however, large regional differences were acknowledged in GNSS4SWEC MoU, namely “while the production, exploitation and evaluation of operational 3 Advanced GNSS Processing Techniques (Working Group 1) 179 GNSS tropospheric products for NWP is well established in the Northern and Western Europe, it is still an emerging R&D field in Eastern and South-Eastern Europe”. In 2014, with the signature of national agreement between Sofia University and Bulgarian BuliPOS GNSS network, the Sofia University GNSS Analysis Centre (SUGAC, http://suada.phys.uni-sofia.bg/) was established to fill the gap of production of GNSS tropospheric products for Bulgaria and South-East Europe (Simeonov et al. 2015). The GNSS4SWEC supported STSM for knowledge transfer and processing of 1 year of tropospheric products from Bulgaria. 3.7.2.2 Main results The first SUGAC processing campaign took place during the STSM of Tzvetan Simeonov to University of Luxembourg (for details see Chap. 7). 7 stations from the BULiPOS network (red pointers in Fig. 3.104) were processed using the NAvigation Package for Earth Observation Satellites (NAPEOS) software version 3.3.1. NAPEOS is developed and maintained by the European Space Operations Centre of the European Space Agency (ESA). First SUGAC processing campaign was performed using the Global Mapping Function and 10 elevation angle cutoff. The RINEX files were processed using the PPP strategy and IGS satellite orbits and clocks. ZTD was computed every 300 s (5 min) for one year – 2013. The ZTD data was used to: (1) estimate IWV and (2) evaluate the numerical weather prediction Fig. 3.104 Map of Bulgaria with marked (red dots) ground based stations of the BuliPOS GNSS network used in first SUGAC processing campaign 180 J. Douša et al. Fig. 3.105 Monthly mean IWV from GNSS (blue markers) and WRF (red markers) for: (a) Burgas, (b) Shumen, (c) Stara Zagora, (d) Montana, (e) Varna and (f) Rozhen in 2013 (NWP) model for Bulgaria (Simeonov et al. 2016). In order to derive IWV with sub-hourly temporal resolution the surface pressure and temperature from the Weather Research and Forecast (WRF) NWP model were used. WRF model simulations were initialized at 00:00 UTC and computed on horizontal mesh of 9 km with 44 vertical levels over Bulgaria. Separately the IWV from WRF is computed by integrating the vertical profile of the water vapour density. The comparisons of monthly mean GNSS and WRF IWV at stations Burgas, Shumen, Stara Zagora, Montana, Varna and Rozhen are presented in Fig. 3.105. At all stations, with exception of Rozhen, the monthly mean IWV minimum is 10 [kg/m2] in December 2013 and the maximum is up to 25 [kg/m2] in June 2013. For station Burgas (Fig. 3.105a) good agreement between the monthly mean IWV from GNSS and WRF is seen with correlation coefficient between 0.96 and 0.84 and Root Mean Square Error (RMSE) between 1.8 and 2.8 [kg/m2] (Fig. 3.106). The maximum and minimum correlation is seen in winter and autumn, and spring and summer, respectively. Between stations Shumen (Fig. 3.105b) and Stara Zagora (Fig. 3.105c) similarities in the IWV can be observed. The two stations also share low RMSE of 2.3 and 2.5 [kg/m2] respectively. For Shumen the lowest correlation is observed in April and it remains low during the spring months. For Stara Zagora the correlation coefficient stays low in summer with minimum from April till August. Montana (Fig. 3.105d) is in Northwest Bulgaria where the influence of the Balkan mountains is significant and the interaction with synoptic flows plays a major role for the IWV distribution. The lowest GNSS and WRF IWV values are seen for December 12 [kg/m2] and the highest for June with 27 [kg/m2]. For Varna (Fig. 3.105e) of interest is the difference between GNSS and WRF, which is seen in January–April 3 Advanced GNSS Processing Techniques (Working Group 1) 181 Fig. 3.106 IWV RMSE and correlation between model and GNSS datasets for the Bulipos Network for 2013. Colors indicate: blue – first quarter of the year (JFM), green – second quarter (AMJ), red – third quarter (JAS), orange – fourth quarter (OND) and May–December period (marked with gray circle). From January to April the IWV in the WRF is lower than the GNSS and from May to December it is the opposite. Similar GNSS IWV jump between April and May is seen at Rozhen (gray circle on Fig. 3.105f). The reason for IWV jump at station Varna is identified to be change of antenna type in the RINEX file. However, for station Rozhen the reason for the IWV change needs further investigation. A summary of the correlation and RMSE for each season is presented in Fig. 3.106. For the cold part of the year (January to March and October to December) the correlation between the model and observation is highest and the RMSE is below 2.6 [kg/m2]. The warm part of the year quarter 2 and 3 (April–September) is with RMSE over 2.6 [kg/m2] and correlation below 0.9. This is expected and is related to the increased atmospheric dynamics and summer time instability and convection, which are a well-known weakness of the NWP models. A detailed investigation of WRF model performance in summer is given in Chap. 4 of this report (Slavchev and Guerova). 3.7.2.3 Future Work The first SUGAC processing campaign was a first step in building expertise in Bulgaria with processing GNSS for remote-sensing the troposphere. The work will continue by developing a pilot transnational severe weather service exploiting GNSS 182 J. Douša et al. tropospheric products to enhance the safety, the quality of life and environmental protection in the Balkan-Mediterranean region (Bulgaria, Cyprus and Greece). 3.7.3 TU Wien Near Real-Time GNSS Analysis Centre in Austria (TUW AC): First Processing Results G. Möller Department of Geodesy and Geoinformation, TU Wien, Wien, Austria e-mail: [email protected] Since March 2017 TU Wien provides near real-time ZTDs for selected GNSS reference sites in Austria and neighbouring countries, see Fig. 3.107. The processing is based on dual-frequency GPS and GLONASS observations, which are provided on a routine basis from the national reference network provider EPOSA (www.eposa.at) in hourly batches. The processing is carried out at TU Wien using the BSW52 double-difference processing strategy. The routines for processing were established within the national research project GNSS-MET Austria in the years 2009–2010 (see Karabatic et al. 2011) and were further refined during the framework of the COST action. For reliable ambiguity resolution, the observations available within the last 8 h are processed altogether. Therefore, it can be guaranteed that at least 65% (long term average) of the ambiguities can be fixed to their integer values. The accuracy of the tropospheric estimates is evaluated regularly at selected IGS sites against the final IGS tropospheric estimates. Figure 3.108 shows the results of the comparison, exemplary for station GRAZ over the first 31 days in 2018. Therefore, from each near real-time solution only the estimates of the last hour were considered. Except of a few outliers, the differences in Fig. 3.107 GNSS station distribution. (Blue) GNSS stations of the Austrian reference network EPOSA, (black) IGS and EUREF stations 3 Advanced GNSS Processing Techniques (Working Group 1) 183 Fig. 3.108 Comparison of TUW near real-time ZTDs with IGS final ZTDs at GNSS site Graz, Austria. Analysed period: First 31 days in 2018 ZTD vary between +/ 1 cm with a mean bias of 1 mm and a standard deviation of of about +/ 4 mm. A similar result is obtained for other IGS sites. 3.7.4 New Operational Solutions from ROB in Support to Global NWP Models and Rapid-Update Numerical Nowcasting E. Pottiaux Royal Observatory of Belgium, Brussels, Belgium e-mail: [email protected] Fostered by the successful developments done by ROB in the context of this COST Action ES1206, the existing near real-time regional troposphere monitoring operated continuously by ROB in the framework of E-GVAP was upgraded to the latest BSW52, using now GPS + GLONASS observations, and processing them in a double-difference batch approach with the latest modelling techniques. In that 184 J. Douša et al. process, a special attention was given to develop a highly flexible and robust processing chain using the BSW52 and the BPE (similar philosophy as the system developed by GOP and described in Sect. 3.6.1), allowing thereby using the same core processing system for various specific applications. This core systems also aims to minimize the manual intervention. As a consequence, we could extend our support to the meteorological community by developing two new troposphere monitoring systems: 1. A 15-min updated regional troposphere monitoring to support nowcasting applications, and 2. A near real-time global troposphere monitoring to support global NWP models. Similarly to the legacy ROB solution to E-GVAP, these two new contributions uses the BSW52, GPS + GLONASS observations, the latest modelling techniques, and are now fully operationally provided to all E-GVAP partners. These new monitoring systems are shortly describes below. 3.7.4.1 New Sub-Hourly GPS + GLONASS Troposphere Monitoring The sub-hourly operational solution operated by ROB includes about 235 GNSS stations providing real-time observations throughout several NTRIP broadcaster servers (Fig. 3.109). Its main objective is to enable rapid-update cycle NWP data assimilation and non-numerical nowcasting applications in the BENELUX + UK regions. This solution has an update cycle of 15 min, and is uploaded to E-GVAP (named ROBQ) with a latency of max 10–15 min after the last observation included in the processing. As expected due to the strong requirement of latency, the precision of the solution is slightly less good than for a standard regional near real-time solution, but it remains within the requirements imposed for such applications (Offiler et al. 2010). Further developments and improvements are still expected to improve this support. Fig. 3.109 Left: Location of the GNSS stations currently included in ROB’s 15-min updated operational contribution to E-GVAP. Right: ZTD time series (and its formal error below) for the Belgian EPN station located in Denterghem from this solution (ROBQ, Status: 3 October 2017) 3 Advanced GNSS Processing Techniques (Working Group 1) 3.7.4.2 185 New Near Real-Time GPS + GLONASS Global Troposphere Monitoring The global operational solutions operated by ROB includes about 315 GNSS stations (named ROBG, Fig. 3.110). Its main objective is to support data assimilation in global NWP models such as those from the U.K Met Office and Météo France. It also potentially supports meteorological agencies outside Europe but collaborating with E-GVAP (e.g. Environment Canada) to access our products. The precision of the solution is similar to the one of the standard regional near real-time solution. With ROB’s global solution, E-GVAP has now (solely) 3 global solutions (ROB, GOP and U.K. Met Office). This allows redundancy and a combination process (ASIC solution by ASI/e-Geos) at common GNSS sites, but ROBG also improves the global coverage by processing sites that are not (yet) processed by the two other ACs. Further works on this solution include performance tuning, and adding more sites in specific area such as in Antarctica to further improve the spatial coverage for global NWP models. 3.7.5 New Methods to User GNSS Vapor Estimates for Meteorology (NUVEM) R. Fernandes University of Beira Interior, Covilhã, Portugal e-mail: [email protected] H. Valentim University of Beira Interior, Covilhã, Portugal e-mail: [email protected] Fig. 3.110 Left: Location of the GNSS stations currently included in ROB’s hourly-updated global operational contribution to E-GVAP (ROBG, Status: 3 October 2017). Right: ZTD time series (and its formal error below) for the EPN station NYA1 located in Ny-Alesund, Norway 186 J. Douša et al. P. Viterbo Instituto Português do Mar e da Atmosfera, Lisbon, Portugal e-mail: [email protected] J. P. Martins Instituto Português do Mar e da Atmosfera, Lisbon, Portugal e-mail: [email protected] A. Sá Polytechnic Institute of Guarda, Guarda, Portugal e-mail: [email protected] J. Jones Met Office, Exeter, UK e-mail: jonathan.jones@metoffice.gov.uk The goal was to include GNSS PWV estimates in weather forecast of Portugal, especially in the decision process of warning dissemination of severe weather situations. For that purpose, Instituto Português do Mar e da Atmosfera (IPMA) and Space and Earth Analysis Laboratory (SEGAL) set up a scheme that provides PWV based on GNSS estimates. Scripts were developed to automatically retrieve the raw GNSS from ftp servers to the server at SEGAL each hour. NUVEM is using 146 stations from 6 GNSS networks over Portugal and Spain, see Fig. 3.111. Fig. 3.111 Stations used to estimate the PWV solutions in the framework of NUVEM project 3 Advanced GNSS Processing Techniques (Working Group 1) 3.7.5.1 187 GNSS Data Providers Portugal: RENEP (1); SERVIR (2) 1. http://www.dgterritorio.pt/cartografia_e_geodesia/geodesia/redes_geodesicas/ renep/ 2. http://www.igeoe.pt/servir/servir.asp Spain: IGN (3); Castilla (4); Extremadura (5); Andalucia (6) 1. 2. 3. 4. http://www.fomento.es http://gnss.itacyl.es http://194.224.247.162:8080/WebExtremadura/ http://www.juntadeandalucia.es/obraspublicasytransportes/ redandaluzadeposicionamiento/rap/ The product (GNSS-PWV maps) is available every 5 min with a 2 h delay and can be checked through a dedicated website (http://nuvem.di.ubi.pt/) that was created for the NUVEM project, where all relevant information, including operational results were/are being published. The ZTDs are estimated at SEGAL, which collects data provided by GNSS networks and information about the GNSS satellite orbits provided by the NASA Jet Propulsion Laboratory as inputs to the GIPSY-OASIS software that is responsible for calculating the ZTD. The conversion of ZTD to PWV is also performed at SEGAL and requires additional information about the pressure and temperature in the vicinity of the GNSS stations. In the developed scheme, these variables are extracted and interpolated for each of the stations by IPMA from the forecasts provided by the ECMWF and are sent to the SEGAL as soon as forecasts are available, with up to about 12 h in advance of their use by the GIPSY. The PWV estimates are then sent to the IPMA, where they are archived and made available to the Operational Centre of Time Forecasting in the form of maps every 15 min, although the available information allows maps every 5 min. In the moment, the products currently available for free by JPL only allow 2 h delay estimates, limiting their use in the context of nowcasting. Although the project is over, SEGAL and IPMA still maintain the operation and there are plans to improve the products (better outlier detection, quality flagging of the retrievals, reducing delay, increasing availability, etc. as well as comparison to other data sources). 3.7.6 Near Real-Time GNSS Processing at ASI/CGS, Italy R. Pacione e-GEOS/Centro di Geodesia Spaziale-Agenzia Spaziale Italiana, Matera, MT, Italy e-mail: [email protected] 188 J. Douša et al. ASI/CGS has been processing Near-Real Time data for E-GVAP since its beginning. During the years of the COST Action ES1206, the existing Near Real-Time network continuously analysed by ASI/CGS in the framework of E-GVAP was upgraded by adding as many GNSS stations as possible in order to homogenize the coverage of troposphere products over Italy. GNSS data belonging to the following regional GNSS networks: Veneto, Liguria, Piemonte, Friuli Venezia Giulia, Trentino, Umbria, Puglia, Calabria, Lazio, Abruzzo, Campania and the NetGeo commercial GNSS network for the Sardinia Island were added to the core network based on EPN and ASI stations. As of today, about 250 stations in the Central Mediterranean part of Europe are recognized by the NRT processing system. In the E-GVAP framework, ASI/CGS is participating as Analysis Centre and acts as Combination Centre delivering four tropospheric solutions: 1. Near Real Time ZTD (Operational, labelled ASI_): every hour, 150 ZTD estimates with a 1h450 latency for an European network of more than 250 sites (Fig. 3.112, left); 2. Near Real Time Combined ZTD (Operational, labelled ASIC): every hour, the 150 ZTD estimates from the contributing E-GVAP Analysis Centres are combined and made available to the project, using a combination scheme outlined in Pacione et al. 2011. On hourly basis about 550 stations on a global scale are combined (Fig. 3.112, right); 3. Near Real Time ZTD (Test, labelled ASIR): the aim of this solution is to evaluate IGS RT products in hourly PPP for NWP application; 4. Sub-hourly ZTD (Test, labelled ASIS): the aim of this solution is to test RT GNSS observation and products in sub-hourly PPP for now-casting application. For ASI_, ASIR and ASIS solutions, GIPSY-OASIS II software (Webb and Zumberge 1997) is used for data reduction. In particular, for ASI_ the standard Fig. 3.112 Left: ASI E-GVAP operational network, Right: ASI E-GVAP Combined Network 3 Advanced GNSS Processing Techniques (Working Group 1) 189 technique of network adjustment is used fixing the IGS Ultra Rapid orbits. A 4-h sliding window approach for data handling is applied with a sampling rate of 5 min and an elevation cut-off angle for the data of 100. The ZWD is estimated every 5 min with a stochastic model (random walk) and a constraint of 20 mm/h1/2. The station coordinates are kept fixed to values provided by combining 1 month of daily postprocessed solutions and are updated every 30 days taking into account the tectonic movements of the area as reported in Pacione and Vespe 2008. 3.7.7 New Analysis Centre (AUTh) and National Observatory of Athens (NOA) C. Pikridas Aristotle University of Thessaloniki, Thessaloniki, Greece e-mail: [email protected] N. Zinas Tekmon Geomatics, Ioánnina, Greece e-mail: [email protected] A. Ganas National Observatory of Athens, Athens, Greece e-mail: [email protected] 3.7.7.1 New Analysis Centre (AUTh) In the frame of a Short Term Scientific Mission on October 2014, a new analysis centre (AC) for near real-time GNSS tropospheric monitoring in Greece was established at the Department of Surveying Engineering of the Aristotle University of Thessaloniki (AUTh). Since then the AUTh Analysis Centre contributes to the EGVAP hourly ZTDs from many permanent GNSS stations in Greece (Fig. 3.113) using the Trop-NET Engine. The AC provides a unique contribution of tropospheric products to the meteorological community for the E-GVAP project that cover the whole of Greece. During the STSM the GOP’s (Geodetic Observatory Pecný) TropNET engine was installed. AUTh AC operates the BSW52 and handles GNSS data from its own and collaborated networks. Additionally, 36 GNSS stations (Fig. 3.114) from IGS and EUREF are included in the network processing scheme for datum definition and consistent absolute tropospheric estimation. The near real-time (NRT) processing engine includes the following three modules: 1. flexible data, metadata, precise product and model downloading or mirroring, 190 Fig. 3.113 GNSS stations in Greece contributing to E-GVAP Fig. 3.114 IGS and EUREF GNSS stations J. Douša et al. 3 Advanced GNSS Processing Techniques (Working Group 1) 191 2. module for GNSS processing based on BSW52 and the BPE, and 3. tropospheric product filtering module for converting in the COST-716 format version 2.2. Currently the tropospheric products are uploaded to the GOP data centre, Met-Office and AUTh ftp archive. In order to have a continuous quality monitoring of the estimated results (like coordinates) and product evaluation the AUTh AC research team (GNSS_QC) developed various scripts for automatic plots of ZTD and coordinates values for each GNSS station (Fig. 3.115). Finally, in collaboration with two other COST participating countries, Bulgaria and Cyprus, the AUTh Research team received funding under the frame of the European Territorial Cooperation Programme “Interreg V-B Balkan-Mediterranean 2014–2020” for the project BeRTISS (Balkan-Meditteranean Real Time Severe weather Service). The main objective of BeRTISS is to develop and establish a pilot transnational severe weather service by exploiting Global Navigation Satellite Systems (GNSS) tropospheric products to enhance the quality and safety of life in the BalkanMediterranean region. This monitoring service will provide continuous and uninterruptible information for nowcasting, forecasting and early warning for Fig. 3.115 Coordinate residuals for GNSS reference station AUT1 192 J. Douša et al. PWV using the GNSS derived tropospheric products and WRF (Weather Research and Forecasting) model that will be tangible and visible to the public through a dedicated web-platform. In detail, the aims of the project are: (1) Integration of networks of GNSS stations located in the three countries in a unified system, (2) Collection, processing and analysis of GNSS observations and tropospheric products, (3) Calculation of the meteorological parameter IWV/PWV for more accurate short-term prediction of severe weather events and (4) Creation of a dedicated website to provide in near real-time the National Meteorological Services and the public with PWV data and warnings of severe weather events. BeRTISS comprises the continuation of the EU-COST Action “GNSS4SWEC” in particular with respect to the expansion of GNSS tropospheric products in one of the Europe’s most remote region and vulnerable to climate change. 3.7.7.2 National Observatory of Athens (NOA) National Observatory of Athens (NOA) operates NOANET the nationwide geodetic network of 22 CORS stations (www.gein.noa.gr/gps.html). NOANET daily 30-s data are distributed via the GSAC web service http://194.177.194.238:8080/ noanetgsac/ as part of the ongoing project EPOS-IP https://epos-ip.org/. In total the NOA GSAC distributes data from 62 CORS stations in SE Europe. In addition, NOA (1) conducted several GPS field campaigns (re-measuring the position of benchmarks) of the CRL/Lefkada/Messinia network in scheduled missions (2) continued installation and maintenance of CORS stations, in the CRL area, in Rhodes (station KATC owned by UNAVCO) and in Messinia (new station ANIK) and (3) conducted installation and maintenance of six (6) continuous GNSS stations in the Ionian Sea area after the 17 November 2015 Earthquake in Lefkada (DRAN, EXAN, ASSO, FISK, KIPO, VLSM). NOA continues to apply space geodesy techniques (SAR interferometry and GNSS) as an important tool for mapping regional surface deformations due to tectonic motions and large earthquakes (e.g. Ganas et al. 2015; 2016; 2017; ongoing research by Athanassios Ganas, Panagiotis Elias, Panagiotis Argyrakis and Alexandra Moshou). In addition, NOA is engaged in research activities within the CRL project (http://crlab.eu/) with emphasis on the effect of the troposphere (ongoing research by Nikos Roukounakis and Panagiotis Elias). The troposphere introduces a path delay in the radio signal, which, in the case of GPS, can be partially removed with the use of specialized mapping functions. Moreover, tropospheric stratification and short wavelength spatial turbulences produce an additive noise to the ground deformation calculated by the (multitemporal) INSAR methodology. The objective is to further correct the vertical component in GPS measurements with the use of a high resolution meteorological model (WRF), producing a 3D tomography of the troposphere. Thus, the knowledge of the tropospheric parameters along the propagation medium can be used to estimate and minimize the effect of this noise, so that the remaining signal represents the deformation mostly due to tectonic or other geophysical processes (Roukounakis et al. 2015). The data-contributing stations 3 Advanced GNSS Processing Techniques (Working Group 1) 193 belong to the Corinth Rift Laboratory and NOANET networks, which monitor the seismicity of the region on a permanent basis. Results are compared with tropospheric delays derived from WRF re-analysis. 3.7.8 New Analysis Centre in Iceland, Icelandic Meteorological Office (IMO) S. Thorsteinsson The Icelandic Meteorological Institute, Reykjavík, Iceland e-mail: [email protected] B. G. Ófeigsson The Icelandic Meteorological Institute, Reykjavík, Iceland e-mail: [email protected] Icelandic Meteorological Office (IMO) operates most monitoring networks of natural hazards in Iceland, see Fig. 3.116. Its operations range from Meteorological monitoring, Hydrological monitoring to volcano and seismic hazard monitoring. As Fig. 3.116 Continuous GNSS network in Iceland (ISGPS). The diamonds mark the location of a GNSS station 194 J. Douša et al. a part of this monitoring effort IMO operates a continuous GNSS network mostly used for volcano monitoring. With the support of Geodetic Obseratory Pecný (GOP), IMO established an analysis centre for near real-time regional troposphere monitoring, using GOP’s Trop-NET system in March 2016. Since then the majority of Icelandic continuous GNSS stations already operated by IMO are routinely processed. The analyses is now a part of IMO’s continuous operational systems where it has become one of the analysis centres providing GNSS ZTD data to E-GVAP. IMO has started and plans to continue to do assimilation impact studies with the GNSS ZTD data as well as from 4 GNSS sites in Greenland gotten from E-GVAP in HARMONIE on the 2.5 km IGB grid domain. New decision regarding the common operational system between IMO and the Danish Meteorological Institute (DMI) is to extend the domain to cover the whole Greenland and Iceland and its surrounding islands, termed IGB domain. The HARMONIE tools and the Icelandic processing GNSS ZTD centre that we have developed in COST ES1206 to monitor convective clouds and severe weather conditions will become useful for IMO. 3.7.9 New Analysis Centre in Hungary (BUTE) Sz. Rozsa Budapest University of Technology and Economics, Budapest, Hungary e-mail: [email protected] The new Analysis Centre (SGO1) was set up in 2014. It processes the Hungarian active GNSS network (37 stations with the mean distance of 60 kms) including some permanent stations from the neighbouring countries (19) as well as some EUREF station as fiducial stations. The near realtime processing is done using BSW52 and the estimated ZTD values are automatically transmitted to the E-GVAP programme. References Alber, C., Ware, R. H., Rocken, C., & Solheim, F. S. (1997). GPS surveying with 1 mm precision using corrections for atmospheric slant path delays. Geophysical Research Letters, 24, 1859–1862. Alber, C., Ware, R., Rocken, C., & Braun, J. (2000). Obtaining single path phase delays from GPS double differences. Geophysical Research Letters, 27, 2661–2664. https://doi.org/10.1029/ 2000gl011525. Altiner, Y., Söhne, W., & Weber, G. (2009). Echtzeitnahe Schätzung troposphärischer LaufzeitVerzögerung und Koordinatenüberwachung, AVN 2/2009, pp. 42–47. Altiner, Y., Mervart, L., Söhne, W., & Weber, G. (2010). Real-time PPP results from global orbit and clock corrections. https://igsac-cnes.cls.fr/documents/egu10/EGU2010-11969_Altiner.pdf Altiner, Y., Söhne, W., & Weber, G. (2011). Von NRT zu RT. http://www.uni-stuttgart.de/gi/ research/Geodaetische_Woche/2011/Session4/S4-18-Altiner.pdf 3 Advanced GNSS Processing Techniques (Working Group 1) 195 Askne, J., & Nordius, H. (1987). Estimation of tropospheric delay for microwaves from surface weather data. Radio Science, 22(3), 379–386. Bar-Sever, Y. E., & Kroger, P. M. (1998). Estimating horizontal gradients of tropospheric path delay with a single GPS receiver. Journal of Geophysical Research, 103, 5019–5035. Bennitt, G., & Jupp, A. (2012). Operational assimilation of GPS zenith total delay observations into the UK met Office numerical weather prediction models. Monthly Weather Review, 140(8), 2706–2719. https://doi.org/10.1175/MWR-D-11-00156.1. Blewitt, G. (1989). Carrier phase ambiguity resolution for the global positioning system applied to geodetic baselines up to 2000 km. Journal of Geophysical Research, 94, 10,187–10,203. https:// doi.org/10.1029/JB094iB08p10187. Boehm, J., Niell, A., Tregoning, P., & Schuh, H. (2006a). Global mapping function (GMF): A new empirical mapping function based on numerical weather model data. Geophysical Research Letters, 33, L07304. https://doi.org/10.1029/2005GL025546. Boehm, J., Werl, B., & Schuh, H. (2006b). Troposphere mapping functions for GPS and very long baseline interferometry from European Centre for Medium-Range Weather Forecasts operational analysis data. Journal of Geophysical Research, 111, B02406. https://doi.org/10.1029/ 2005JB003629. Böhm, J., Möller, G., Schindelegger, M., Pain, G., & Weber, R. (2015). Development of an improved empirical model for slant delays in the troposphere (GPT2w). GPS Solutions, Springer, 19(3), 433–441. Brenot, H., & Warnant, R. (2008). Characterization of the tropospheric small-scale activity, Technical Report ESA, WP250, GALOCAD project. https://orbi.uliege.be/jspui/bitstream/ 2268/84849/1/GALOCAD-WP250-Report-OUT250-1_v2.pdf Brenot, H., Ducrocq, V., Walpersdorf, A., Champollion, C., & Caumont, O. (2006). GPS zenith delay sensitivity evaluated from high-resolution numerical weather prediction simulations of the 8–9 September 2002 flash flood over southeastern France. Journal of Geophysical Research, 111, D15105. https://doi.org/10.1029/2004JD005726. Brenot, H., Neméghaire, J., Delobbe, L., Clerbaux, N., De Meutter, P., Deckmyn, A., Delcloo, A., Frappez, L., & Van Roozendael, M. (2013). Preliminary signs of the initiation of deep convection by GNSS. Atmospheric Chemistry and Physics, 13, 5425–5449. https://doi.org/10. 5194/acp-13-5425-2013 (licensed under CC BY 3.0, https://creativecommons.org/licenses/by/ 3.0/). Brenot, H., Wautelet, G., Warnant, R., Nemeghaire, J., & Van Roozendael, M. (2014a) GNSS meteorology and impact on NRT position. Proceedings of the European Navigation Conference (ENC), Rotterdam, The Netherlands. Brenot, H., Walpersdorf, A., Reverdy, M., van Baelen, J., Ducrocq, V., Champollion, C., Masson, F., Doerflinger, E., Collard, P., & Giroux, P. (2014b). A GPS network for tropospheric tomography in the framework of the Mediterranean hydrometeorological observatory Cévennes-Vivarais (southeastern France). Atmospheric Measurement Techniques, 7, 553–578. https://doi.org/10.5194/amt-7-553-2014. Caissy, M., Agrotis, L., Weber, G., Hernandez-Pajares, M., & Hugentobler, U. (2012). INNOVATION-coming soon-the international GNSS real-time service. GPS World, 23, 52–58. Champollion, C., Masson, F., Van Baelen, J., Walpersdorf, A., Chéry, J., & Doerflinger, E. (2004). GPS monitoring of the tropospheric water vapour distribution and variation during the September 9, 2002, torrential precipitation episode in the Cévennes (southern France). Journal of Geophysical Research, 109, D24. Chen, G., & Herring, T. A. (1997). Effects of Atmospheric Azimuthal Asymetry on the Analysis of Space Geodetic Data. Geophysical Research Letters, 102(20), 489–420. 502. Dach, R., & Jean, Y. (2015). International GNSS service, technical report 2014. Pasadena: IGS Central Bureau. Dach, R., Schaer, S., Lutz, S., Baumann, C., Bock, H., Orliac, E., Prange, L., Thaller, D., Mervart, L., Jäggi, A., Beutler, G., Brockmann, E., Ineichen, D., Wiget, A., Weber, G., Habrich, H., Söhne, W., & Ihde, J. (2014). Steigenberger, P., and Hugentobler, U: CODE IGS Analysis 196 J. Douša et al. Center Technical Report 2013, Dach, R., & Jean, Y. (Eds.), IGS 2013 Technical Report, pp. 21–34. Dach, R., Lutz, S., Walser, P., & Fridez, P. (2015). Bernese GNSS Software Version 5.2, University of Bern, Bern Open Publishing. https://doi.org/10.7892/boris.72297. Dach, R., Stefan, S., Arnold, D., Orliac, E., Prange, L., Suśnik, A., Villiger, A., & Jäggi, A. (2016) CODE final product series for the IGS, Published by Astronomical Institute, University of Bern. https://doi.org/10.7892/boris.75876. Davis, J. L., Herring, T. A., Shapiro, I. I., Rogers, A. E. E., & Elgered, G. (1985). Geodesy by interferometry: Effects of atmospheric modeling errors on estimates of baseline length. Radio Science, 20, 1593–1607. Davis, J. L., Elgered, G., Niell, A. E., & Kuehn, C. E. (1993). Groundbased measurements of gradients in the “wet” radio refractivity of air. Radio Science, 28, 1003–1018. De Haan, S. (2006). National/regional operational procedures of GPS water vapour networks and agreed international procedures, Rep WMO/TD-No, 1340:20, KNMI, Netherlands. Dee, D. P., et al. (2011). The ERA-Interim reanalysis: Configuration and performance of the data assimilation system. Quarterly Journal of the Royal Meteorological Society, 137(656), 553–597. Delobbe, L., & Holleman, I. (2006). Uncertainties in radar echo top heights used for hail detection. Meteorological Applications, 13, 361–374. Delrieu, G., Nicol, J., Yates, E., Kirstetter, P.-E., Creutin, J.-D., Anquetin, S., Obled, C., Saulnier, G.-M., Ducrocq, V., Gaume, E., Payrastre, O., Andrieu, H., Ayral, P.-A., Bouvier, C., Neppel, L., Livet, M., Lang, M., Parent du-Châtelet, J., Walpersdorf, A., & Wobrock, W. (2005). The Catastrophic Flash-Flood Event of 8–9 September 2002 in the Gard Region, France: a First Case Study for the Cévennes-Vivarais Mediterranean Hydrometeorological Observatory, Journal of Hydrometeorology, 6, 34–51. https://doi.org/10.1175/JHM400.1. Deng, Z., Schöne, T., & Gendt, G. (2014). Status of the TIGA Tide Gauge Data Reprocessing at GFZ. International Association of Geodesy Symposia, October 2014. Deng, Z., Uhlemann, M., Fritsche, M., Dick, G., & Wickert, J. (2015). Troposphere parameters derived from multi-GNSS data processing at GFZ. Wien, EGU 2015. Deng, Z., Nischan, T., & Bradke, M. (2017). Multi-GNSS Rapid Orbit-, Clock- & EOP-Product Series. GFZ Data Services. https://doi.org/10.5880/GFZ.1.1.2017.002. Dick, G., Gendt, G., & Reigber, C. (2001). First experience with near real-time water vapor estimation in a German GPS network. Journal of Atmospheric and Solar – Terrestrial Physics, 63, 1295–1304. https://doi.org/10.1016/S1364-6826(00)00248-0. Ding, W., Teferle, F. N., Kazmierski, K., Laurichesse, D., & Yuan, Y. (2017). An evaluation of real-time troposphere estimation based on GNSS Precise Point Positioning. Journal of Geophysical Research Atmosphere, 122(5), 2779–2790. https://doi.org/10.1002/2016JD025727. Dong, D.-N., & Bock, Y. (1989). Global positioning system network analysis with phase ambiguity resolution applied to crustal deformation studies in California. Journal of Geophysical Research, 94, 3949–3966. https://doi.org/10.1029/JB094iB04p03949. Douša, J. (2001a). Towards an operational near-real time precipitable water vapor estimation. Physics andChemistry of the Earth, Part. A, 26, 189–194. https://doi.org/10.1016/S1464-1895 (01)00045-X. Douša, J. (2001b). The impact of ultra-rapid orbits on precipitable water vapor estimation using ground GPS network. Physics and Chemistry of the Earth, 26(6–8), 393–398. Douša, J. (2004a). Evaluation of tropospheric parameters estimated in various routine analysis. Physics and Chemistry of the Earth, 29(2–3), 167–175. Douša, J. (2004b). Precise orbits for ground-based GPS meteorology: Processing strategy and quality assessment of the orbits determined at geodetic observatory pecný. Journal of the Meteorological Society of Japan, 82, 371–380. Douša, J. (2010). Precise near real-time GNSS analyses at geodetic observatory pecný – Precise orbit determination and water vapour monitoring. Acta Geodynamics et Geomaterialia, 7(157), 1–11. 3 Advanced GNSS Processing Techniques (Working Group 1) 197 Douša, J. (2012). Developments of the GLONASS ultra-rapid orbit determination at geodetic observatory pecný. In S. Kenyon, M. C. Pacino, & U. Marti (Eds.), Geodesy of planet earth (IAG Symposia Series) (Vol. 136, pp. 1029–1036). Dordrecht: Springer. Douša, J., & Bennitt, G. V. (2013). Estimation and evaluation of hourly updated global GPS zenith Total delays over ten months. GPS Solutions, 17, 453–464. https://doi.org/10.1007/s10291-0120291-7. Douša, J., & Eliaš, M. (2014). An improved model for calculating tropospheric wet delay. Geophysical Research Letters, 41, 4389–4397. Douša, J., & Souček, P. (2005). The results of near real-time COST-716 GPS campaign from geodetic observatory Pecný, In J. Sledzinski (Eds.), Proceedings of the EGU G9 symposium (pp. 139–150). Reports on Geodesy, Warsaw University of Technology, Institute of Geodesy and Geodetic Astronomy, Warsaw, Poland. Douša, J., & Václavovic, P. (2014). Real-time zenith tropospheric delays in support of numerical weather prediction applications. Advances in Space Research, 53(9), 1347–1358. https://doi. org/10.1016/j.asr.2014.02.021. Douša, J., & Václavovic, P. (2016). Evaluation of ground-based GNSS tropospheric products at Geodetic Observatory Pecný. Proceedings of the IAG Symposia Series, Springer, Rizos Ch. & Willis P. (Eds.). 143:759–766. 10.1007/1345_2015_157. Douša, J., Elias, M., Veerman H., van Leeuwen, S. S., Zelle, H., de Haan, S., Martellucci, A., & Perez, R. O. (2015a). High accuracy tropospheric delay determination based on improved modelling and high resolution numerical weather model. Proceedings of the ION GNSS 2015, Institute of Navigation, Tampa, Florida, USA, September 14–18, pp. 3734–3744. Douša, J., Václavovic, P., Krč, P., Eliaš, M., Eben, E., & Resler, J. (2015b). NWM forecast monitoring with near real-time GNSS products, In Proceedings of the 5th scientific Galileo Colloquium, Braunschweig, Germany, October 27–29. http://old.esaconferencebureau.com/ docs/default-source/15a08_session1/044-dousa.pdf?sfvrsn¼2 Douša, J., Dick, G., Kačmařík, M., Brožková, R., Zus, F., Brenot, H., Stoycheva, A., Möller, G., & Kaplon, J. (2016). Benchmark campaign and case study episode in Central Europe for development and assessment of advanced GNSS tropospheric models and products. Atmospheric Measurement Techniques, 9, 2989–3008. https://doi.org/10.5194/amt-9-2989-2016. (licensed under CC BY 3.0, https://creativecommons.org/licenses/by/3.0/). Douša, J., Václavovic, P., & Eliaš, M. (2017). Tropospheric products of the second European GNSS reprocessing (1996–2014). Atmospheric Measurement Techniques, 10, 1–19. https://doi.org/10. 5194/amt-10-1-2017. (licensed under CC BY 3.0, https://creativecommons.org/licenses/by/3.0/). Douša, J., Václavovic, P., Zhao, L., & Kačmařík, M. (2018a). New adaptable all-in-one strategy for estimating advanced tropospheric parameters and using real-time orbits and clocks. Remote Sensing, 10, 232. https://doi.org/10.3390/rs10020232. (licensed under CC BY 4.0, https:// creativecommons.org/licenses/by/4.0/). Douša, J., Eliaš, M., Václavovic, P., Eben, K., & Krc, P. (2018b). A two-stage tropospheric correction combining data from GNSS and numerical weather model. GPS Solutions, 22, 77. https://doi.org/10.1007/s10291-018-0742-x. Dow, J. M., Neilan, R. E., & Rizos, C. (2009). The international GNSS service in a changing landscape of global navigation satellite systems. Journal of Geodesy, 83, 191–198. https://doi. org/10.1007/s00190-008-0300-3. Elgered, G., Davis, J. L., Herring, T. A., & Shapiro, I. I. (1991). Geodesy by radio interferometry: Water vapour radiometry for estimation of the wet delay. Journal of Geophysical Research, 96, 6541–6555. Elósegui, P., Davis, J. L., Jaldehag, R. T. K., Niell, A. E., & Shapiro, I. I. (1995). Geodesy using the global positioning system: The effects of signal scattering on estimates of site position. Journal of Geophysical Research, 100, 99219934. Elósegui, P., Davis, J. L., Gradinarsky, L. P., Elgered, G., Johansson, J. M., Tahmoush, D. A., & Ruis, A. (1999). Sensing atmospheric structure using small-scale space geodetic networks. Geophysical Research Letters, 26, 2445–2448. 198 J. Douša et al. Flores, A., Ruffini, G., & Rius, A. (2000). 4D tropospheric tomography using GPS slant wet delays. Annales de Geophysique, 18, 223–234. https://doi.org/10.1007/s00585-000-0223-7. Ge, M., Gendt, G., Dick, G., Zhang, F. P., & Rothacher, M. (2006). A new data processing strategy for huge GNSS global networks. Journal of Geodesy, 80(4), 199–203. Gendt, G., Dick, G., Reigber, C., Tomassini, M., Liu, Y., & Ramatschi, M. (2004). Near real time GPS water vapor monitoring for numerical weather prediction in Germany. Journal of the Meteorological Society of Japan, 82, 361–370. Gradinarsky, L. P. (2002). Sensing atmospheric water vapor using radio waves, Ph.D. thesis, School of Electrical Engineering, Chalmers University of Technology, Göteborg, Sweden. Guerova, G., Jones, J., Douša, J., Dick, G., de Haan, S., Pottiaux, E., Bock, O., Pacione, R., Elgered, G., Vedel, H., & Bender, M. (2016a). Review of the state of the art and future prospects of the ground-based GNSS meteorology in Europe. Atmospheric Measurement Techniques, 9, 5385–5406. https://doi.org/10.5194/amt-9-5385-2016. Guerova, G., Simeonov, Tzv.,& Vassileva, K. (2016b). Comparison of GNSS tropospheric products obtained by two processing strategies. In Proceedings of international SYMPOSIUM on modern technologies, education and professional practice in geodesy and related fields. Sofia, Bulgaria, 3–4/11/2016, 8 p. ISSN 2367-6051. Available at: http://suada.phys.uni-sofia.bg/ wordpress/wp-content/uploads/2017/07/Guerova_etal_Sofia-Symposium_2016.pdf Hadaś, T., & Bosy, J. (2015). IGS RTS precise orbits and clocks verification and quality degradation over time. GPS Solutions, 19, 93–105. https://doi.org/10.1007/s10291-014-0369-5. Hadaś, T., Teferle, F. N., Kaźmierski, K., Hordyniec, P., & Bosy, J. (2017). Optimum stochastic modeling for GNSS tropospheric delay estimation in real-time. GPS Solutions, 21(3), 1069–1081. Berlin – Heidelberg. Herring, T. A. (1992). Modeling atmospheric delays in te analysis of space geodetic data. In Proceedings of the symposium on refraction of transatmospheric signals in geodesy. Netherlands Geodetic Commission, Publications on Geodesy, Delft, the Netherlands. Herring, T. A., King, R. W., & McClusky, S. C. (2010). Documentation for the GAMIT GPS analysis software, version 10.4, Technical report. Cambridge: Massachusetts Institute of Technology. Hofmeister, A., & Böhm, J. (2017). Application of ray-traced tropospheric slant delays to geodetic VLBI analysis. The Journal of Geodesy, 91(8), 945–964. Huisman L., van der Marel H., & Teunissen P. (2009). CORS local-site finger-printing using undifferenced least squares GNSS phase residuals. In International global navigation satellite systems society IGNSS symposium 2009, Holiday Inn, Surfers Paradise, Qld, Australia, December 1–3 Hunegnaw, A., Teferle, F. N., Bingley, R., & Hansen, D. (2015). Status of TIGA activities at the British Isles continuous GNSS facility and the University of Luxembourg, C. Rizos & P. Willis (Eds.), International Association of Geodesy Symposia, 143:617–623. https://doi.org/10.1007/ 1345-2015-77. Cham: Springer. Hunegnaw, A., Teferle, F. N., Abraha, K. E., Santamaría-Gómez, A., Gravelle, M., Wöppelman, G., Schöne, T., Deng, Z., Bingley, R., Hansen, D., Sanchez, L., & Moore, M. (2017). On the Scientific Applications of IGS Products: An Assessment of the Reprocessed TIGA Solutions and Combined Products. IGS Workshop 2017, Paris. Iwabuchi, T., Miyazaki, S., Heki, K., Naito, I., & Hatanaka, Y. (2003). An impact of estimating tropospheric delay gradients on tropospheric delay estimations in the summer using the Japanese nationwide GPS array. Journal of Geophysical Research, 108, 4315. https://doi.org/10. 1029/2002JD002214. Kačmařík, M., Douša, J., Dick, G., Zus, F., Brenot, H., Möller, G., Pottiaux, E., Kapłon, J., Hordyniec, P., Václavovic, P., & Morel, L. (2017). Inter-technique validation of tropospheric slant total delays. Atmospheric Measurement Techniques, 10, 2183–2208. https://doi.org/10. 5194/amt-10-2183-2017. (licensed under CC BY 3.0, https://creativecommons.org/licenses/by/ 3.0/). 3 Advanced GNSS Processing Techniques (Working Group 1) 199 Kapłon, J., Hordyniec, P., & Rohm, W. (2017). Analysis of systematic effects in slant total delay estimation with PPP. Presentation G06-1-06 at IAG-IASPEI 2017 Joint Scientific Assembly of the International Association of Geodesy (IAG) and International Association of Seismology and Physics of the Earth’s Interior (IASPEI), Kobe, Japan, August 1. Karabatic, A., Weber, R., & Haiden, T. (2011). Near real-time estimation of tropospheric water vapour content from ground based GNSS data and its potential contribution to weather now-casting in Austria. Advances in Space Research, 47(10), 1691–1703. https://doi.org/10. 1016/j.asr.2010.10.028. King, R. W., Masters, E. G., Rizos, C., Stolz, A., & Collins, J. (1985). Surveying with GPS, Monograph 9, School of Surveying, University of New South Wales, Kensington, Australia. Klos, A., Hunegnaw, A., Teferle, F. N., Abraha, K. E., Ahmed, F., & Bogusz, J. (2018). Statistical significance of trends in Zenith Wet Delay from re-processed GPS solutions. GPS Solutions, 22, 51. https://doi.org/10.1007/s10291-018-0717-y. Landskron, D., & Böhm, J. (2017). VMF3/GPT3: Refined discrete and empirical troposphere mapping functions. The Journal of Geodesy. https://doi.org/10.1007/s00190-017-1066-2. (licensed under CC BY 4.0, https://creativecommons.org/licenses/by/4.0/). Laurichesse, D. (2011). The CNES real-time PPP with undifferenced integer ambiguity resolution demonstrator. In Proceedings of the ION GNSS 2011, Portland, Oregon. Laurichesse, D., Cerri, L., Berthias, J. P., & Mercier, F. (2013). Real time precise GPS constellation and clocks estimation by means of a Kalman filter. In Proceedings ION GNSS 2013 (pp. 1155–1163), Institute of Navigation, Nashville, Tennessee, USA, September 16–20. Leick, A. (1989). GPS satellite surveying. New York: Wiley-Interscience. Leick, A. (2004). GPS satellite surveying (3rd ed., 435p). New York: Wiley. MacMillan, D. S. (1995). Atmospheric gradients from very long baseline interferometry observations. Geophysical Research Letters, 22(9), 97–102. https://doi.org/10.1029/95GL00887. Marini, J. W. (1972). Correction of satellite tracking data for an arbitrary tropospheric profile. Radio Science, 7(2), 223–231. Möller, G. (2017) Reconstruction of 3D wet refractivity fields in the lower atmosphere along bended GNSS signal paths, Dissertation, TU Wien, Department of Geodesy and Geoinformation. http://repositum.tuwien.ac.at/obvutwhs/content/titleinfo/2268559 Möller, G., Weber, R., & Böhm, J. (2014). Improved troposphere blind models based on numerical weather data. NAVIGATION: Journal of the Institute of Navigation, 61(3), 203–211. Montenbruck, O., Steigenberger, P., Prange, L., Deng, Z., Zhao, Q., Perosanz, F., Romero, I., Noll, C., Stürze, A., Weber, G., Schmid, R., MacLeod, K., & Schaer, S. (2017). The Multi-GNSS Experiment (MGEX) of the International GNSS Service (IGS) – Achievements, prospects and challenges. Advances in Space Research, 59, 1671–1697. https://doi.org/10.1016/j.asr.2017.01. 011. Morel, L., Pottiaux, E., Durand, F., Fund, F., Boniface, K., de Oliveira, P. S., & Van Baelen, J. (2014). Validity and behaviour of tropospheric gradients estimated by GPS in Corsica. Advances in Space Research. https://doi.org/10.1016/j.asr.2014.10.004. Niell, A. E. (1996). Global mapping functions for the atmosphere delay at radio wavelengths. Journal of Geophysical Research, 101, 3227–3246. https://doi.org/10.1029/95JB03048. Ning, T., Wang, J., Elgered, G., Dick, G., Wickert, J., Bradke, M., & Sommer, M. (2016). The uncertainty of the atmospheric integrated water vapour estimated from GNSS observations. Atmospheric Measurement Techniques, 9, 79–92. https://doi.org/10.5194/amt-9-79-2016. Offiler, D. et al. (2010). EIG EUMETNET GNSS Water Vapour Programme (E-GVAP-II), Product Requirements Document, Version 1.0 – 21 December 2010. Pacione, R., & Vespe, F. (2008). Comparative studies for the assessment of the quality of near-real time GPS-derived atmospheric parameters. Journal of Atmospheric and Oceanic Technology, 25, 701–714. Pacione, R., Pace, B., de Haan, S., Vedel, H., Lanotte, R., & Vespe, F. (2011). Combination methods of tropospheric time series. Advances in Space Research, 47, 323–335. https://doi.org/ 10.1016/j.asr.2010.07.021. 200 J. Douša et al. Pacione, R., Araszkiewicz, A., Brockmann, E., & Douša, J. (2017). EPN-Repro2: A reference GNSS tropospheric data set over Europe. Atmospheric Measurement Techniques, 10, 1689–1705. https://doi.org/10.5194/amt-10-1689-2017. Petit, G., Luzum, B. (red.), IERS Conventions. (2010). IERS Technical Note No. 36, Frankfurt am Main: Verlag des Bundesamts für Kartographie und Geodäsie, ISBN 3-89888-989-6 Pottiaux, E., Brockmann, E., Söhne, W., & Bruyninx, C. (2009). The EUREF – EUMETNET collaboration: First experience and potential benefits. Bulletin of Geodesy and Geomatics, 3, 269–288. Prange, L., Orliac, E., Dach, R., Arnold, D., Beutler, G., Schaer, S., & Jäggi, A. (2017). CODE’s five-system orbit and clock solution—The challenges of multi-GNSS data analysis. Journal of Geodesy, 91, 345–360. https://doi.org/10.1007/s00190-016-0968-8. Rebischung, P., Griffiths, J., Ray, J., Schmid, R., Collilieux, X., & Garayt, B. (2012). IGS08: The IGS realization of ITRF2008. GPS Solutions, 16, 483494. https://doi.org/10.1007/s10291-0110248-2. Rebischung, P., Altamimi, Z., Ray, J., & Garayt, B. (2016). The IGS contribution to ITRF2014. Journal of Geodesy, 90, 611–630. https://doi.org/10.1007/s00190-016-0897-6. Roukounakis, N., et al. (2015). Improvement of the vertical component of GPS and INSAR measurements in the western Corinth Gulf (Greece), by the use of high-resolution meteorological modeling of the lower troposphere: The PaTrop Experiment, 12th European Conference on Applications of Meteorology (ECAM), EMS2015-544;11 September 2015, Sofia, Bulgaria. Ruffini, G., Kruse, L. P., Rius, A., Bürki, B., & Cucurull, L. (1999). Estimation of tropospheric zenith delay and gradients over the Madrid area using GPS and WVR data. Geophysical Research Letters, 26(4), 447–450. https://doi.org/10.1029/1998GL900238. Saastamoinen, J. (1972). Atmospheric correction for the troposphere and stratosphere in radio ranging of satellites. The Geophysical Monograph Series, 15, 247–251. Schmid, R., Dach, R., Collilieux, X., Jäggi, A., Schmitz, M., & Dilssner, F. (2016). Absolute IGS antenna phase center model igs08.atx: status and potential improvements. Journal of Geodesy, 90(4), 343–364. Seeber, G. (2003). Satellite geodesy (2nd ed., 589p). New York: de Gruyter. Shoji, Y., Nakamura, H., Iwabuchi, T., Aonashi, K., Seko, H., Mishima, K., Itagaki, A., Ichikawa, R., & Ohtani, R. (2004). Tsukuba GPS dense net campaign observation: Improvement in GPS analysis of slant path delay by stacking one-way post phase residuals. Journal of the Meteorological Society of Japan, 82(1B), 301–314. https://doi.org/10.2151/jmsj.2004.301. Simeonov Tzv., Sidorov, D., Teferle, N, Guerova, G., Egova, E., Vassileva, K., Milev, I., & Milev, G. (2015). Sofia University GNSS Analysis Centre, FIG working week, 17–21/05/2015, Sofia, Bulgaria. https://www.Figurenet/resources/proceedings/fig_proceedings/fig2015/papers/ts08g/ TS08G_simeonov_sidorov_et_al_7677.pdf Simeonov, T., Sidorov, D., Teferle, F. N., Milev, G., & Guerova, G. (2016). Evaluation of IWV from the numerical weather prediction WRF model with PPP GNSS processing for Bulgaria. Atmospheric Measurement Techniques Discussions. https://doi.org/10.5194/amt-2016-152. https://www.atmos-meas-tech-discuss.net/amt-2016-152/amt-2016-152.pdf. Teunissen, P., Kleusberg, A., Bock, Y., Beutler, G., Weber, R., Langley, R. B., van der Marel, H., Blewitt, G., Gload, C., & Colombo, O. L. (1998). GPS for geodesy (2nd ed.). Berlin: Springer. Tregoning, P., & Herring, T. A. (2006). Impact of a priori zenith hydrostatic delay errors on GPS estimates of station heights and zenith total delays. Geophysical Research Letters, 33(L23303). https://doi.org/10.1029/2006GL027706. Trojáková, A. (2016). The NWP activities at CHMI, Joint 26th ALADIN Workshop & HIRLAM All Staff Meeting 2016, Lisbon, Portugal, April 4–8. Urquhart, L., Santos, M., Nievinski, F., & Boehm, J. (2011). Generation and Assessment of VMF1Type Grids using North-American Numerical Weather Models, XXV IUGG General Assembly, Melbourne, Australia, June 28th – July 7th. Available at http://unb-vmf1.gge.unb.ca/Publica tions.html 3 Advanced GNSS Processing Techniques (Working Group 1) 201 Václavovic, P., & Douša, J. (2015). Backward smoothing for precise GNSS applications. Advances in Space Research, 56(8), 627–1634. https://doi.org/10.1016/j.asr.2015.07.020. Václavovic, P., & Douša, J. (2016). G-Nut/Anubis – Open-source tool for multi-GNSS data monitoring. In IAG Symposium. Berlin/Heidelberg: Springer. ISBN-13: 978-3-319-30895-1, 143, 775–782. Václavovic, P., Douša, J., & Gyori, G. (2013). G-Nut software library – State of development and first results. Acta Geodynamica et Geomaterialia, 10(4), 431–436. https://doi.org/10.13168/ AGG.2013.0042. Václavovic, P., Douša, J., Elias, M., & Kostelecky, J. (2017). Using external tropospheric corrections to improve GNSS positioning of hot-air balloon. GPS Solutions, 21(4), 1479–1489. https:// doi.org/10.1007/s10291-017-0628-3. Walpersdorf, A., Calais, E., Haase, J., Eymard, L., Desbois, M., & Vedel, H. (2001). Atmospheric Gradients Estimated by GPS Compared to a High Resolution Numerical Weather Prediction (NWP) Model. Physics and Chemistry of the Earth, 26, 147–152. Warnant, R., Wautelet, G., Spits, J., & Lejeune, S. (2008) Characterization of the ionospheric small-scale activity. Technical Report ESA, WP220, GALOCAD project. Wautelet, G., Lejeune, S., Stankov, S., Brenot, H., & Warnant, R. (2008) Effects of small-scale atmospheric activity on precise positioning. Technical Report ESA, WP230, GALOCAD project. Webb, F. H., & Zumberge, J. F. (1997). An Introduction to GIPSY/OASIS II. JPL D-11088. Weber, G., Dettmering, D., Gebhard H., & Kalafus, R. (2005) Networked Transport of RTCM via Internet Protocol (Ntrip) – IP Streaming for Real-Time GNSS Applications. In Proceedings ION-GNSS 18th international technical meeting of the satellite division (pp. 2243–2247). Weber, G., Mervart, L., Stürze, A., Rülke, A., & Stöcker, D. (2016). BKG Ntrip Client (BNC) Version 2.12, Mitteilungen des Bundesamtes für Kartographie und Geodäsie, Band 49, ISBN 978-3-86482-083-0. Wilgan, K., Hurter, F., Geiger, A., Rohm, W., & Bosy, J. (2017a). Tropospheric refractivity and zenith path delays from least-squares collocation of meteorological and GNSS data. The Journal of Geodesy, 91(2), 117–134. (licensed under CC BY 4.0, https://creativecommons.org/licenses/ by/4.0/). Wilgan, K., Hadaś, T., Hordyniec, T., & Bosy, J. (2017b). Real-time PPP augmented with highresolution NWM model data. GPS Solutions, 23(3), 1341–1353. (licensed under CC BY 4.0, https://creativecommons.org/licenses/by/4.0/). Ye, S., Zhao, L., Song, J., Chen, D., & Jiang, W. (2018). Analysis of estimated satellite clock biases and their effects on precise point positioning. GPS Solutions, 22, 16. https://doi.org/10.1007/ s10291-017-0680-z. Zumberge, J. F., Heflin, M. B., Jefferson, D. C., Watkins, M. M., & Webb, F. H. (1997). Precise point positioning for the efficient and robust analysis of GPS data from large networks. Journal of Geophysical Research, 102, 5005–5017. https://doi.org/10.1029/96JB03860. Zus, F., Dick, G., Heise, S., Douša, J., & Wickert, J. (2014). The rapid and precise computation of GPS slant total delays and mapping factors utilizing a numerical weather model. Radio Science, 49(3), 207–216. Zus, F., Dick, G., Douša, J., & Wickert, J. (2015a). Systematic errors of mapping functions which are based on the VMF1 concept. GPS Solutions, 19, 277–286. Zus, F., Dick, G., Heise, S., & Wickert, J. (2015b). A forward operator and its adjoint for GPS slant total delays. Radio Science, 50, 393–405.