https://orcid.org/0000-0002-3379-9213
SUMMARY
Earth observation satellites yield a wealth of data for scientific, operational and commercial exploitation. However, the redistribution of mass in the system Earth is not yet part of the standard inventory of Earth Observation (EO) data products to date. It is derived from the Gravity Recovery and Climate Experiment (GRACE) mission and its Follow-On mission (GRACE-FO). Among many other applications, mass redistribution provides fundamental insights into the global water cycle. Changes in continental water storage impact the regional water budget and can, in extreme cases, result in floods and droughts that often claim a high toll on infrastructure, economy and human lives. The initiative for a European Gravity Service for Improved Emergency Management (EGSIEM) established three different prototype services to promote the unique value of mass redistribution products for Earth Observation in general and for early-warning systems in particular. The first prototype service is a scientific combination service to derive improved mass redistribution products from the combined knowledge of the European GRACE analysis centres. Second, the timeliness and reliability of such products is a primary concern for any early-warning system and therefore EGSIEM established a prototype for a near real-time service that provides dedicated gravity field information with a maximum latency of 5 d. Third, EGSIEM established a prototype of a hydrological/early warning service that derives wetness indices as indicators of hydrological extremes and assessed their potential for timely scheduling of high-resolution optical/radar satellites for follow-up observations in case of evolving hydrological extreme events.
1 INTRODUCTION
Until the end of the last century, global monitoring of temporal changes in the distribution of environmental mass had been a challenge. Before the advent of satellite missions dedicated to global gravity field recovery such as the CHAllenging Minisatellite Payload (CHAMP; cf. Reigber et al.1998), the Gravity Recovery and Climate Experiment (GRACE; cf. Tapley et al.2004), and the Gravity field and steady-state Ocean Circulation Explorer (GOCE; cf. Drinkwater et al.2006), geodetic satellites tracked by Satellite Laser Ranging (SLR) and altimeter and other satellites tracked by Doppler Orbitography and Radiopositioning Integrated by Satellite (DORIS) made the biggest contribution to the determination of Earth’s gravity field models (e.g. Cheng et al.1997; Lemoine et al.1997; Bianco et al.1998). However, since the launch of the dedicated GRACE mission in 2002, intersatellite ranging has been established as the state-of-the-art technique to globally observe mass variations in the system Earth from space, as evidenced by several thousand scientific publications in the meantime. Although GRACE enabled spectacular science results, for example, summarized by Wouters et al. (2014) and Tapley et al. (2019), numerous important questions regarding the changes and dynamic processes in the continental hydrology, cryosphere, ocean, atmosphere and solid Earth remain unresolved (Pail et al.2015). Sustained and improved observation systems such as the recently launched GRACE Follow-On (GRACE-FO; Flechtner et al.2013) and next generation gravity missions (e.g. Cesare et al.2013) are needed to extend the available time-series and to provide another leap in accuracy and spatial resolution. This will be a prerequisite to, for example, better separate human-induced changes from natural climate changes. However, a full exploitation of the potential of future gravity missions will only be possible if this may also be achieved for the less precise data of the current GRACE mission and if, in parallel, limiting error sources such as accelerometer errors and temporal aliasing errors may be further reduced (Loomis et al.2012; Flechtner et al.2016).
More than 17 yr after the launch of the GRACE satellites, the exploitation of the GRACE data is still ongoing. A growing number of Analysis Centres (ACs) inside and outside the GRACE Science Data System (SDS) is still challenged with the Level-1B data processing. To date each new release of monthly gravity fields, provided to the scientific user community as Level-2 products in a (usually) spherical harmonic (SH) representation, still represents a substantial and significant improvement with respect to previous releases. Official SDS solutions are provided by the Center for Space Research (CSR; Bettadpur 2012) and the German Research Center for Geosciences (GFZ; Dahle et al.2012). Additional solutions are computed by the Jet Propulsion Laboratory (JPL) for validation (Watkins & Yuan 2012). Apart from these official monthly GRACE gravity field solutions, a growing number of ACs are providing additional solutions, for example, the Centre National d’Etudes Spatiales (CNES; Bruinsma et al.2010), the Technical University of Graz (TUG, solutions labelled ITSG Mayer-Gürr et al.2010) and the Astronomical Institute of the University of Bern (AIUB; Meyer et al.2016). Unfortunately, however, the release of different ACs are not in every case fully comparable. This is mainly observed in terms of noise (Jean et al.2018), but sometimes differences may also result in terms of signal due to different processing strategies. A well-known pitfall was the GFZ RL05 solution (Dahle et al.2014) as illustrated in Fig. 1 for Greenland mass change estimates. At the end of the shown time period substantial differences of a few hundred Gigatons have been observed with respect to the other solutions due to an unintended (hidden) regularization, which was caused by separating the inherently linked orbit determination of the two GRACE satellites from the process of recovering the gravity field (Meyer et al.2015). As a lesson learned one may conclude that a certain agreement on common processing standards will help to ensure the compatibility of the different solutions in view of a proper geophysical interpretation.
Trends of Greenland mass change estimated from different releases of SDS solutions (RL05) based on different processing strategies (left-hand panel) and with respect to the CSR solution (right-hand panel). No filtering and no further corrections (GIA, leakage) are applied.
On the other hand, the example of the GFZ RL05 solution also nicely shows that the availability of several competing solutions is very beneficial to detect outliers. The current situation is, however, sub-optimal from the user perspective to fully exploit the information content offered by the data of the past GRACE and the recently launched GRACE-FO mission. Today the end user is essentially left with the difficult choice of deciding (1) which solutions from which AC to use and (2) what actions to perform for converting SH coefficients (Level-2 products) into gridded mass values (Level-3 products) appropriate for their study when not selecting one of the three official SDS solutions.1
In contrast to other Earth Observation (EO) data, satellite-based measurements of gravity represent total water storage variations, that is, variations of all surface and subsurface water storage compartments. As such the past GRACE and the recently launched GRACE-FO mission provide unique information on the wetness state of a river basin with regard to its actual flood generation potential or its susceptibility to a drought. Reager & Famiglietti (2009) estimated flood potential at the regional scale by means of determining repeated maxima in water storage anomalies, which suggest an effective storage capacity in a region, beyond which additional precipitation must be met by increases in runoff or evaporation. Thomas et al. (2014) presented a quantitative approach for measuring hydrological drought occurrence and severity based on GRACE data by calculating the magnitude of the deviation of regional, monthly total water storage anomalies from the time-series’ monthly climatology. Humphrey et al. (2016) surveyed key features of temporal variability in the GRACE record by decomposing gridded times-series of monthly equivalent water height (EWH) into linear trends, interannual, seasonal and intra-annual components, with an additional focus on extreme dry anomalies and their relation to documented drought events.
Today it takes approximately two months from the time that the Level-1B data is collected onboard the GRACE satellites to the time when scientists can access and examine the Level-2 or Level-3 products provided to the scientific community. The temporal sampling of the gravity field solutions is at best 10 d but most reliably restricted to one month when not using regularization techniques (Bruinsma et al.2010). Both of these time constraints currently limit the potential of using the results from satellite gravimetry in time-critical monitoring applications (Pail et al.2015). This applies in particular for early-warning and forecasting systems of extreme hydrological events. Flood forecast models need, for example, near real-time (NRT) information to estimate the probable development of the event in terms of flood stage or river discharge with typical lead times of a few days for larger river basins (see overview by Emerton et al.2016). Also the usefulness of high-resolution follow-up observations such as optical and radar EO data for emergency management is strongly influenced by the time-span from alert reception, satellite programming, satellite acquisition and data reception (Voigt et al.2016).
Given these limitations, the GRACE ACs in Europe have joined forces and established the European Gravity Service for Improved Emergency Management (EGSIEM) project, an initiative supported by the Horizon 2020 Framework Program for Research and Innovation of the European Commission and the Swiss State Secretariat for Education, Research and Innovation (SERI). In this article the EGSIEM concept is presented, the current implementation status is reported and future perspectives are outlined. The article is structured as follows: Section 2 introduces the general concept of the EGSIEM initiative. Section 3 presents the complementary data prepared within EGSIEM, whereas Section 3.2 is dedicated to the GRACE Level-1B data processing. Sections 4–6 introduce the EGSIEM prototype services and present the current achievements. Section 7 closes the article with a summary and outlook.
2 CONCEPT OF EGSIEM AND ITS PROTOTYPE SERVICES
The main driver of the EGSIEM initiative was that gravity-based observations of the redistribution of water and ice masses provide valuable and complementary information with respect to more traditional EO data as provided by, for example, Copernicus, the European Programme for the establishment of a European capacity for Earth Observation.2 We, therefore, defined three key objectives: (1) demonstrating that superior gravity products can be derived by combining solutions from different ACs, (2) providing NRT gravity products and (3) assessing the potential of these NRT products for flood and drought monitoring and forecasting. These objectives are achieved by setting up three dedicated prototype services, namely a scientific combination service, a NRT service and a hydrological service. Fig. 2 illustrates the EGSIEM service structure, the interaction between the individual services and the input data used. The services and the achieved results are described in more detail in Sections 4, 5 and 6, respectively.
General concept of the EGSIEM: Satellite data from Altimetry, Gravity, GNSS, SLR and Copernicus missions were used to create three prototype services.
2.1 Scientific combination service
The main input data to the first two EGSIEM prototype services are GRACE Level-1B data. For combination purposes, Level-1B data needs to be consistently processed by several ACs in the post-processing mode to derive monthly gravity field solutions. These are then combined by the Scientific Combination Service to derive high quality scientific products with a time resolution of one month and latencies of up to 60 d. Since all EO systems must refer to one and the same reference frame, the consistent use of the reference frame is an essential step in the EGSIEM concept. To the extent possible a generic formulation of processing standards has been developed and consequently been applied for the GRACE processing at all EGSIEM ACs (EGSIEM 2015). It is expected that the measures pave the way for a long awaited standardization in the gravity field community and to significantly increase the quality, robustness and reliability of products derived from satellite gravimetry.
A dedicated GPS reprocessing campaign was conducted for a consistent integration of data from the global tracking network of the International GNSS Service (IGS; cf. Dow et al.2009) as described in Section 3.1 and to provide one single set of up-to-date GPS products to all EGSIEM ACs (see Section 3.1). Subsequently, ACs may then apply their preferred, individual approach to process the GRACE input data. The approaches differ by the used observables, the parametrizations, the stochastic noise models applied and the relative weighting of the different observables. All ACs are required to provide free solutions, that is, solutions that do not depend on an a priori gravity model. This is of crucial importance to avoid biases in the combined solutions. If all solutions are unbiased, the strengths and weaknesses as well as the different background models of the different solutions are expected to average out in the combination and the combined solution will be statistically better and more robust. Details about the processing strategies and the specific background models used at the individual EGSIEM ACs are documented in EGSIEM (2015).
2.2 Near real-time service
The main input data are again GRACE Level-1B data for the recovery of the Earth’s gravity field but this time individual ACs derive in their NRT mode rapid gravity products for flood and drought alerting with latencies down to 5 d. The increased temporal resolution from typically 1 month to 1 d is achieved by using sophisticated regularization techniques to compensate for the loss of spatial resolution.
GFZ converts the measured intersatellite tracking data by means of spacecraft velocities into daily differences of the dynamic forcing acting on the twin satellites that are then reduced by forces stemming from geophysical background models and non-conservative origin. The reduced data are then projected by means of functions of geopotential gradient differences to Earth gravity potential variations at ground level. The TU Graz approach directly relates the measured kinematic orbit positions and intersatellite tracking data to the unknown gravity field parameters. After reduction of all modelled forces, the daily gravity variations are recovered using the methodology of Kurtenbach et al. (2012). Further details are provided in Section 5.
The provision of the corresponding NRT gravity field information within 5 d is a unique product that has the potential to translate into added value for warning and forecasting the onset of natural hazards.
2.3 Hydrological service
The main input data of the Hydrological Service are the rapid gravity products from the NRT service to quantify water storage anomalies and to derive indicators of hydrological extreme events such as floods and droughts, whereas the post-processed combined monthly solutions are used for comparison purposes. Their value is especially assessed in view of supporting operational satellite-based flood information services such as those within the framework of the DLR’s (Deutsches Zentrum für Luft-und Raumfahrt) Center for Satellite Based Crisis Information (ZKI3). A large variety of other EO data such as medium- to high-resolution Synthetic Aperture Radar (SAR) and optical satellite data, for example, from Envisat ASAR, TerraSAR-X, MODIS, Sentinel-1 and the Laboratoire d’Etudes en Géophysique et Océanographie Spatiales (LEGOS) hydroweb data base from altimetry as well as global flood data bases such as the Dartmouth Flood Observatory (Brakenridge 2016) are analysed to compile suitable historical and large-scale flood situations to assess the value of the gravity-based flood and drought indicators.
3 REPROCESSING OF INPUT DATA
3.1 GNSS reprocessing
Besides ultraprecise K-Band intersatellite ranging, the GRACE satellites, as well as all other dedicated gravity missions, are equipped with onboard dual-frequency Global Positioning System (GPS) receivers to allow for precise GPS high-low satellite-to-satellite tracking (GPS hl-SST). In case of the GRACE and the recently launched GRACE-FO mission the analysis of GPS data is thus an inherent part of the Level-1B data processing to infer monthly gravity field solutions (e.g. Beutler et al.2010). For this purpose high-quality information on GPS satellite orbits and GPS satellite clock corrections is necessary. State-of-the-art products are routinely provided by the International GNSS Service (IGS; cf. Dow et al.2009). The operational IGS products are, however, continuously improved by taking the latest developments into account, for example, updates in the underlying reference frame (Altamimi et al.2011, 2016) or the conventions recommended by the International Earth Rotation and Reference Systems Service (IERS; cf. Petit & Luzum 2010). As a consequence the operational products are inevitably inhomogeneous over time. Reprocessing campaigns are thus needed to serve applications where long and homogenously processed GPS product series are crucial (e.g. Steigenberger et al.2006; Fritsche et al.2014).
Because previous IGS reprocessings did not usually provide the necessary high-rate (5 s) GPS satellite clock information, as it is needed for a most precise determination of the orbits of low Earth orbiting (LEO) satellites (Bock et al.2009), a dedicated reprocessing effort was initiated within the frame of EGSIEM (Sušnik et al., in preparation). This allowed us to take into account the latest improvements in GNSS orbit modelling (Arnold et al.2015), which in turn allowed for a further improvement in kinematic LEO orbit determination and the subsequent recovery of the Earth’s gravity field from kinematic positions. For the first time 5 s high-rate satellite clock corrections are not only provided for the GPS but also for the GLONASS satellites. This reprocessing campaign thus not only assisted spaceborne GPS but also terrestrial multi-GNSS precise point position (PPP) applications (Zumberge et al.1997). GLONASS satellite clock corrections are provided after the year 2008 (30 s) when the full GLONASS constellation was complete, and with 5 s sampling from 2010 onwards. For a detailed description of the EGSIEM reprocessing campaign we refer to Sušnik et al. (in preparation).
3.2 GRACE reprocessing
Each associated AC was requested to provide monthly gravity fields of 2006-2007 that agreed with the EGSIEM standards concerning reference frame, Earth orientation, satellite geometry/antenna reference points, relativistic effects and third bodies perturbing the satellites’ motion. The processing had to be based on common GPS orbits and clock corrections (see Section 3.1). All gravity field contributions were requested to be free solutions (i.e. not regularized) complete to degree and order 90. The normal equations of the monthly GRACE gravity fields were generated in the Solution (Software/technique) INdependent EXchange (SINEX) format.
This reprocessing was used by the individual analysis centres to revise their processing strategy in general. Major changes compared to the pre-EGSIEM gravity field time-series were
AIUB: revision of observation screening, low-pass filtering of geometric K-band correction,
GFZ: 3 hr accelerometer scale factors and biases in three directions, down-weighting of GPS phase observations from 0.7 to 1.0 cm mean error, ocean tide model update to FES2014 (Carrere et al.2016),
CNES: free and complete solution instead of single-value-decomposition (SVD) approach, extension to degree and order 90, revision of relative weighting GPS/KRR, reduction of resolution of GPS-contribution to degree and order 40.
TUG: introduction of observation screening, use of star camera and accelerometer sensor fusion attitude data and geometric K-band correction, accelerometer parametrization (e.g. full scale matrix; see Klinger et al.2016), orbit integration based on elliptic reference orbit instead of linear motion, increase of correlation length of empirical covariances from 1 to 3 hr, improved constraining of co-estimated daily variations based on geophysical models.
In Fig. 3 we compare the EGSIEM releases of the individual gravity fields with their respective predecessors for the 2 yr demonstration period of 2006-07. It illustrates global grids of the RMS of EWH anomalies. Anomalies are defined by subtraction of deterministic models of secular and seasonal changes per grid cell from monthly grids of EWH variations. While over the continents strong non-seasonal signals of mainly hydological origin are visible, the ocean areas are very quiet, with few exceptions like in the South Atlantic. We therefore concentrate on the oceans to assess the noise content of the gravity field time-series. A significant reduction of noise is visible in the case of TUG and GFZ, a minor reduction of noise in the case of AIUB. For CNES no corresponding comparison can be shown since no free solution comparable to the EGSIEM release exists.
RMS of anomalies of the filtered predecessor releases (left-hand column) and the EGSIEM releases (right-hand column) of the monthly gravity field solutions computed at the different ACs.
4 SCIENTIFIC COMBINATION SERVICE
The goals of the scientific combination service as described in Section 2 are achieved by careful standardization, screening and combination of the individual AC solutions described in Section 3.2. Relative weights are derived on solution level by variance component estimation. The final combination is performed on normal equation level to produce optimally combined solutions. Eventually external validations are investigated to independently assess the quality of the derived solutions.
4.1 Combination approach
To ensure the quality of the combined products, the individual contributions are compared in terms of signal and noise content. The signal content is assessed by analysing the amplitude of seasonal variations of EWH in three selected river basins and by ice mass trends in Greenland and West Antarctica (e.g. Meyer et al.2016). An attempt to characterize the noise levels is made by the evaluation of anomalies, that is, the monthly variability after subtraction of a deterministic model of secular and seasonal variations. Anomalies are studied either spatially in regions of little variability, that is, over the oceans, or spectrally for a range of spherical harmonic coefficients with a strong noise contribution, for example, between about degrees 50–90. The purpose of the comparison is twofold. Time-series of monthly gravity fields with reduced signal content, for example, due to regularization, are rejected in order to avoid a bias in the combined solution. Furthermore, individual monthly gravity fields with increased noise level, for example, due to bad observational data, are also screened out in order to not deterioriate the combined solution.
4.1.1 Level-2 products
After individual screening all ACs use a subset of the same GRACE K-Band range-rate observations. No additional information is therefore introduced by the combination. But the noise levels of the individual ACs’ monthly gravity fields are today still dominated by background model errors (Flechtner et al.2016) and individual analysis noise (Jean et al.2018), but not due to observation errors (see baseline accuracy in Kim 2000). These errors are specific to the individual ACs and therefore reduced in the combination, which in turn increases the signal to noise ratio. In the frame of EGSIEM combinations are performed on both the solution level and the normal equation (NEQ) level. As opposed to the findings of Sakumura et al. (2014), who could not demonstrate weighting to be beneficial when combining the three official SDS solutions, a weighted average was found to be crucial to account for the very different noise levels of the individual gravity field solutions. As described in detail by Jean et al. (2018), relative weights are iteratively derived per solution for each month by comparing the individual solutions to a weighted average by using the mathematical framework of variance component estimation (e.g. Koch 2007).
The combination on NEQ-level is superior in taking into account the correlations between the estimated corrections to the gravity field model and satellite- and arc-specific instrument and orbit parameters. But due to very diverse noise models and formal error characteristics, a standard weighting scheme based on variance factors (e.g. Koch 2007) does not lead to satisfactory results. To reduce the effect of the different noise models the individual NEQs are scaled until equal contribution of the individual NEQs to pairwise combinations is achieved and consequently the noise-based relative weights derived on solution level are applied as described in detail by Meyer et al. (2018).
Within the frame of the EGSIEM project monthly combined solutions were produced for the demonstration period 2006 January to 2007 December. The individual as well as the combined solutions for 2006 January are shown in Fig. 4 in terms of the median (2006-07) of degree amplitudes of anomalies of spherical harmonic coefficients, expressed in EWH. In order to limit the comparison to the physically meaningful part of the spectrum, Fig. 4 (bottom) only considers the orders 0–29. Fig. 4 confirms that the quality of the four individual contributions of the EGSIEM ACs is quite diverse, with the ITSG contribution being dominate in terms of low noise. This can be attributed to the different approaches to absorb noise, either by pseudo-stochastic accelerations (AIUB), K-band and accelerometer instrument parameters (GFZ and GRGS) or a sophisticated empirical noise model and observation-dependent weighting (as in case of ITSG). In Fig. 5. the monthly RMS of anomalies over the oceans are shown to assess the noise of the individual and combined solutions in the spatial domain.
Median (2006-07) of degree amplitudes of coefficient wise anomalies taking into account all orders (top) or only orders up to 29 (bottom).
Weighted RMS over the oceans of filtered monthly solutions of the EGSIEM ACs and of the combined monthly solutions.
The combined solutions outperform the individual contributions in the part of the spectrum with a strong noise component, that is, beyond degree 31 (corresponding to the second resonant order 31). Compared to the combination on solution level, the combination on NEQ-level turns out to be more robust, for example, in the vicinity of resonant orders, because accuracy information and correlations are taken correctly into account. A more in-depth analysis may be found in Meyer et al. (2018).
The combined products may be visualized by the dedicated EGSIEM plotter and are distributed via the EGSIEM webpage4 and the International Center for Global Earth Models (ICGEM), either in spherical harmonic coefficients (Level-2 products) or as post-processed grids of EWHs (Level-3 products). The latter are generated individually for hydrological and oceanographic applications.
4.1.2 Level-3 products
To derive Level-3 (Table 1.) from Level-2 products, as a first step the full (non-tidal) signal content is reconstructed by adding the monthly means of the atmosphere and ocean de-alising products back to the monthly fields. These are combined from the monthly means provided by the individual ACs using the same relative weights as derived for the monthly fields by VCE on solution level. As a next step, spherical harmonic degree 1 coefficients derived by Satellite Laser Ranging (SLR) are added to transform from a centre of mass to the centre of figure frame. C20 is not replaced, because no clear advantage of SLR derived values over the combined GRACE solution could be found.
Definition of L3-products.
| L3-product | Constituents |
|---|---|
| GRACE non-tidal | GRACE + monthly means of atmosphere and ocean de-aliasing products |
| Hydrology | GRACE non-tidal − atmosphere − ocean − GIA |
| Oceanography | GRACE non-tidal − atmosphere − ocean − GIA − hydrology + OBP |
| L3-product | Constituents |
|---|---|
| GRACE non-tidal | GRACE + monthly means of atmosphere and ocean de-aliasing products |
| Hydrology | GRACE non-tidal − atmosphere − ocean − GIA |
| Oceanography | GRACE non-tidal − atmosphere − ocean − GIA − hydrology + OBP |
Definition of L3-products.
| L3-product | Constituents |
|---|---|
| GRACE non-tidal | GRACE + monthly means of atmosphere and ocean de-aliasing products |
| Hydrology | GRACE non-tidal − atmosphere − ocean − GIA |
| Oceanography | GRACE non-tidal − atmosphere − ocean − GIA − hydrology + OBP |
| L3-product | Constituents |
|---|---|
| GRACE non-tidal | GRACE + monthly means of atmosphere and ocean de-aliasing products |
| Hydrology | GRACE non-tidal − atmosphere − ocean − GIA |
| Oceanography | GRACE non-tidal − atmosphere − ocean − GIA − hydrology + OBP |
The gravity signal induced by the last glaciation about 20 000 yr ago overlaps signals of other sources such as hydrology, and thus needs to be effectively removed with a model of glacial isostatic adjustment (GIA). Within EGSIEM a new GIA model has been developed based on a series of regional ice history models (Tarasov et al.2012; Briggs et al.2014; Lambeck et al.2014; Lecavalier et al.2014; Hughes et al.2015; Nordman et al.2015; Root et al.2015) which were combined temporally and spatially to form a global ice model called LM17.3 (Steffen et al.2017), and the laterally homogeneous VM5a Earth model (Argus et al.2014). The ice model has a spatial resolution of 0.5 × 0.5 degrees and consists of 52 time steps of different length related to distinct past changes. The GIA model filtered with DDK5 (Kusche et al.2009), called EGSIEM GIA LM17.3 DDK5, can be downloaded from the EGSIEM plotter. It is planned to provide an update with a slightly improved ice model and a laterally heterogeneous (3D) Earth model soon.
For ocean-related applications the monthly mean of ocean bottom pressure (OBP) is added back and a hydrological model is subtracted to avoid signal leakage from the continents (the specific models applied are listed in Table 2). Finally the spherical harmonic coefficients are transformed to global grids with one degree resolution, applying specific filters for either hydrological or oceanographic applications. The filter matrices are derived in analogy to the DDK filters (Kusche et al.2009), but instead of a simulated, average GRACE noise model, the monthly calibrated errors of EGSIEM-ITSG are used. All monthly filter matrices are available to the users via the EGSIEM webpage.
Models applied for the definition of Level-3.
| Degree 1 SHC | AIUB-SLRa (Sośnica et al.2015) |
| Atmosphere | AOD1B-RL05: GAA |
| Ocean | AOD1B-RL05: GAB |
| OBP | AOD1B-RL05: GAD |
| GIA | LM17.3 (Steffen et al.2017) |
| Hydrosphere | WGHM, evaluated at the epochs of the monthly GRACE gravity fields (Döll et al.2003) |
Models applied for the definition of Level-3.
| Degree 1 SHC | AIUB-SLRa (Sośnica et al.2015) |
| Atmosphere | AOD1B-RL05: GAA |
| Ocean | AOD1B-RL05: GAB |
| OBP | AOD1B-RL05: GAD |
| GIA | LM17.3 (Steffen et al.2017) |
| Hydrosphere | WGHM, evaluated at the epochs of the monthly GRACE gravity fields (Döll et al.2003) |
4.2 External validation of gravity field solutions
Ensuring data quality is an essential part of EGSIEM. This is achieved by comparing gravity-field derived station displacements with vertical displacement time-series derived from GNSS according to Davis et al. (2004) and van Dam et al. (2007). Besides the individual solutions derived by the four individual ACs in Section 3.2 and the EGSIEM combined solution from Section 4.1, we additionally include also the three SDS gravity field solutions. Note that we have validated all the EGSIEM Level-2 and Level-3 gravity products within the EGSIEM project, while we present here only results of validating Level-2 combined gravity solutions. For results of validating other EGSIEM gravity products, readers are referred to Chen et al. (2018).
The validation procedure described here closely follows Chen et al. (2018). This includes post-processing of the GRACE data, processing of three different GNSS data sets used for the validation, as well as the computation of performance evaluation indicators, which are degree and cumulative degree WRMS reduction measures. For detailed information regarding these aspects we refer to Chen et al. (2018). Here we only mention some slight differences with respect to Chen et al. (2018). For the GRACE data post-processing we add back SLR-derived degree 1 coefficients from Sośnica et al. (2015) and we restore the AOD1B atmospheric and oceanic dealiasing products from the EGSIEM combined dealiasing products (see Section 4.1). Regarding the three GNSS station time-series, we use only those station time-series that do not suffer from any data gap within the 2-yr period of our EGSIEM gravity field solutions. As such, we have 178 station coordinates from the EGSIEM GNSS reprocessing products (see Section 3.1), (2) coordinate solutions of 388 stations from the ITRF2014 realization (Rebischung et al.2016) and (3) 377 station time-series of the JPL GNSS products (Bock et al.2012).
Fig. 6 shows the mean degree WRMS reductions and mean cumulative degree WRMS reductions of the eight different gravity field solutions in comparison to the ITRF2014 displacements, which demonstrates similar characteristics as in Chen et al. (2018) in terms of degree contributions. In the context of this article, we are primarily interested in the relative performance of the various solutions. Both the mean degree WRMS and especially the mean cumulative degree WRMS reduction indicate that the combined EGSIEM solutions generally provide the strongest WRMS reductions. Also two of the individual AC solutions, AIUB and ITSG, outperform the three SDS solutions in the low degrees, which can be assigned to the consequent implementation of the processing standards to ensure a higher consistency with the IGS processing standards. Nevertheless, it has to be noted that the ACs solutions from GFZ and GRGS currently do not provide good results at certain low degrees, in particular the imperfect performance at degree 2. Via the degree WRMS reduction analysis, we are able to track the degree 2 problem of the EGSIEM-GRGS solution down to be caused by C21 and S21 terms. Further investigation will be needed to identify the reason of this abnormal behaviour. It is interesting to note that the combined solutions were able to overcome the degree 2 problems from the EGSIEM-GRGS solution and perform best among all EGSIEM solutions, which confirms once more the hypothesis that the combination of solutions may yield an improvement in the gravity field solutions.
Mean degree WRMS reductions (top) and mean accumulative degree WRMS reductions (bottom) of the eight gravity solutions with respect to the ITRF2014 GNSS residuals on the full signal level at 388 GNSS stations. As no significant contributions come from degrees beyond 30, we truncate the plots at degree 40. EGSIEM-COMB indicates the EGSIEM combined solutions at the NEQ level. EGSIEM-AIUB, EGSIEM-GFZ, EGSIEM-GRGS and EGSIEM-ITSG represent gravity solutions from the four individual ACs.
Fig. 7 shows in the spatial domain the WRMS of the reduction at the 388 GNSS stations of the two time-series, that is, one gravity field solution and the ITRF2014 station displacements. Visually, yellow to red colours dominate the spatial patterns indicating strong agreements between GNSS and GRACE. Up to |$72.13\, {\rm per\, cent}$| of WRMS reduction is observed using the EGSIEM-COMB solution at CUIB located in Cuiabá, Brazil and similar patterns are also seen for other gravity solutions. The negative WRMS reductions (light blue to black dots) mostly appear in the GNSS stations located in islands or along the coast which are possibly affected by non-tidal oceanic effects (Tesmer et al.2011) or even spurious long-period signals due to unmodelled short-periodic displacements (Penna et al.2007).
Spatial plots of the WRMS reductions from all the gravity solutions compared to the ITRF2014 GNSS residuals on the full signal level at 388 GNSS stations.
To be more specific, a summary of statistics is displayed in Table 3 which shows the mean and median WRMS reductions on the full signal level as well as the annual signal level from eight gravity solutions with respect to the three GNSS products. With respect to the ITRF2014 and EGSIEM-reprocessed GNSS time-series, the EGSIEM combined solution provides the best statistic numbers in terms of both the full signal and the annual signal level. For instance, up to 26.65 per cent mean WRMS reduction at the full signal level is observed for the EGSIEM-COMB solution comparing to the EGSIEM-reprocessed GNSS time-series. In comparison to the JPL GNSS time-series, the EGSIEM-COMB also demonstrates the top performance as the EGSIEM-ITSG and CSR RL05 gravity solutions.
Mean WRMS reduction on the full signal level and median WRMS reductions on the annual signal level between eight GRACE gravity products and three GNSS solutions. Best performances are indicated in bold letters.
| ITRF2014 | EGSIEM-reprocessed | JPL | ||||
|---|---|---|---|---|---|---|
| Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | |
| EGSIEM-COMB | 25.31 | 72.72 | 26.65 | 72.35 | 23.59 | 67.72 |
| EGSIEM-AIUB | 24.50 | 70.64 | 26.22 | 70.73 | 22.81 | 66.33 |
| EGSIEM-GFZ | 22.17 | 68.28 | 23.08 | 67.44 | 20.27 | 64.65 |
| EGSIEM-GRGS | 16.95 | 63.64 | 22.07 | 66.64 | 18.41 | 62.36 |
| EGSIEM-ITSG | 24.78 | 71.63 | 26.42 | 69.25 | 23.83 | 68.89 |
| GFZ RL05a | 22.61 | 70.20 | 24.96 | 72.01 | 22.16 | 65.87 |
| CSR RL05 | 23.78 | 70.20 | 25.94 | 70.31 | 23.61 | 67.75 |
| JPL RL05.1 | 22.56 | 69.16 | 24.41 | 65.94 | 21.86 | 64.14 |
| ITRF2014 | EGSIEM-reprocessed | JPL | ||||
|---|---|---|---|---|---|---|
| Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | |
| EGSIEM-COMB | 25.31 | 72.72 | 26.65 | 72.35 | 23.59 | 67.72 |
| EGSIEM-AIUB | 24.50 | 70.64 | 26.22 | 70.73 | 22.81 | 66.33 |
| EGSIEM-GFZ | 22.17 | 68.28 | 23.08 | 67.44 | 20.27 | 64.65 |
| EGSIEM-GRGS | 16.95 | 63.64 | 22.07 | 66.64 | 18.41 | 62.36 |
| EGSIEM-ITSG | 24.78 | 71.63 | 26.42 | 69.25 | 23.83 | 68.89 |
| GFZ RL05a | 22.61 | 70.20 | 24.96 | 72.01 | 22.16 | 65.87 |
| CSR RL05 | 23.78 | 70.20 | 25.94 | 70.31 | 23.61 | 67.75 |
| JPL RL05.1 | 22.56 | 69.16 | 24.41 | 65.94 | 21.86 | 64.14 |
Mean WRMS reduction on the full signal level and median WRMS reductions on the annual signal level between eight GRACE gravity products and three GNSS solutions. Best performances are indicated in bold letters.
| ITRF2014 | EGSIEM-reprocessed | JPL | ||||
|---|---|---|---|---|---|---|
| Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | |
| EGSIEM-COMB | 25.31 | 72.72 | 26.65 | 72.35 | 23.59 | 67.72 |
| EGSIEM-AIUB | 24.50 | 70.64 | 26.22 | 70.73 | 22.81 | 66.33 |
| EGSIEM-GFZ | 22.17 | 68.28 | 23.08 | 67.44 | 20.27 | 64.65 |
| EGSIEM-GRGS | 16.95 | 63.64 | 22.07 | 66.64 | 18.41 | 62.36 |
| EGSIEM-ITSG | 24.78 | 71.63 | 26.42 | 69.25 | 23.83 | 68.89 |
| GFZ RL05a | 22.61 | 70.20 | 24.96 | 72.01 | 22.16 | 65.87 |
| CSR RL05 | 23.78 | 70.20 | 25.94 | 70.31 | 23.61 | 67.75 |
| JPL RL05.1 | 22.56 | 69.16 | 24.41 | 65.94 | 21.86 | 64.14 |
| ITRF2014 | EGSIEM-reprocessed | JPL | ||||
|---|---|---|---|---|---|---|
| Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | Full (per cent) | Annual (per cent) | |
| EGSIEM-COMB | 25.31 | 72.72 | 26.65 | 72.35 | 23.59 | 67.72 |
| EGSIEM-AIUB | 24.50 | 70.64 | 26.22 | 70.73 | 22.81 | 66.33 |
| EGSIEM-GFZ | 22.17 | 68.28 | 23.08 | 67.44 | 20.27 | 64.65 |
| EGSIEM-GRGS | 16.95 | 63.64 | 22.07 | 66.64 | 18.41 | 62.36 |
| EGSIEM-ITSG | 24.78 | 71.63 | 26.42 | 69.25 | 23.83 | 68.89 |
| GFZ RL05a | 22.61 | 70.20 | 24.96 | 72.01 | 22.16 | 65.87 |
| CSR RL05 | 23.78 | 70.20 | 25.94 | 70.31 | 23.61 | 67.75 |
| JPL RL05.1 | 22.56 | 69.16 | 24.41 | 65.94 | 21.86 | 64.14 |
5 NRT-SERVICE
The goal of the NRT-Service was to provide daily gravity field solutions with a latency of maximum 5 d. For this purpose, two independent approaches were developed at the German Research Centre for Geosciences (GFZ) and TU Graz (TUG) based on a Kalman filter approach, first introduced by Kurtenbach et al. (2012).
Common to both approaches is the input data consisting of rapid GNSS products provided by the centre of orbit determination (CODE; Dach et al.2009, 2017), JPL’s GRACE Level-1B quick-look data and the GFZ Release 6 GRACE Level-1B atmosphere and ocean de-aliasing product (AOD1B RL06; Dobslaw et al.2017) that are both provided within the GRACE Science Data System. These data sets feature a latency of 17, 24 and 10 hr, respectively, which meets the above requirements to produce NRT gravity field solutions. Since data pre-processing, outlier detection and gravity field computation in both approaches differ significantly, two suitably independent solutions are produced on a daily basis.
GFZ converts the measured intersatellite range-rates and range-accelerations by means of spacecraft velocities derived from precise orbit determination (POD) into differences of the dynamic forcing acting on the twin satellites. These derived dynamic observations can be reduced by forces stemming from geophysical background models, representing the lithosphere density contrast (static gravity field), third body attractions, ocean tides, atmospheric and oceanic non-tidal mass variations (AOD1B) and forces of non-conservative (mainly drag and solar radiation pressure) origin. The reduced data is then projected by means of functions of geopotential gradient differences to Earth gravity potential variations at ground level and expressed in equivalent water layer thickness. The spherical mean Earth surface is discretized into 2 × 2 arc degree equal area surface tiles that are post-processed to regular 2x2 arc degree grids and can be further decomposed in spherical harmonic coefficients for the generation of additional (e.g. Level-3) products. Daily updates are computed in a Kalman filter approach, where external data sets provide the constraints for a prediction in time, driven by the expected average temporal variation of daily water storage in the hydrological system derived from WGHM (Döll et al.2003) and atmospheric and oceanic mass variations derived from AOD1B. Secular variations, for example, of the cryosphere, are contained in the background gravity field model (EIGEN-6C; Shako et al.2014) and are therefore reduced from observations beforehand. This results in statistically smooth and stationary signal. Multistep noise estimates for the observed dynamic forcing completes the measurement update in the Kalman filter. The main idea behind this is the whitening of the observation residuals after the surface signal computation step, which is achieved for the respective autocorrelations if they reduce to Dirac-like pulses. A detailed description of the methodology can be found in Gruber & Gouweleeuw (2019).
TUG relies on intersatellite range rate and kinematic orbit positions that are related to the unknown gravity field parameters by means of variational equations. The gravity field is parametrized by a spherical harmonic series up to degree and order 40. The processing of the GRACE Level-1B data closely follows the strategy used for TUG’s monthly solutions (Klinger et al.2016; Mayer-Gürr et al.2016; Ellmer et al.2017). A key component of this processing chain is the determination of weighting between the individual observation groups by analysing the post-fit residuals of an unconstrained monthly gravity field solution. Using this weighting, unconstrained daily normal equations, which serve as input for the Kalman filter, are computed. The signal content of these daily normal equations is the residual gravity field left in the GRACE observations, after reducing all background models. Consequently, the major constituents left are sub-seasonal continental hydrosphere and cryosphere, as well as errors in the dealiasing model. This is reflected in the Kalman filter process dynamics, which is composed of empirical auto- and cross-covariance estimates of geophysical model output of exactly these signal sources. In the case at hand, the hydrosphere is represented by WGHM while cryosphere and dealiasing errors are represented by the ‘I’ component and dealiasing model error estimates contained in the ESA Earth System Model (Dobslaw et al.2015). The models were only evaluated up to the beginning of 2002 to avoid any overlap with the GRACE time-series. Furthermore, an isotropic noise component of approximately 7 mm EWH which corresponds to 5 per cent of the average global signal was added to account for unmodelled effects. A detailed description of the methodology may be found in Kvas & Mayer-Gürr (in preparation).
The performance of the daily gravity solutions has been tested using historical hydrological extreme events (Gouweleeuw et al.2018) and relative intercomparisons between both NRT time-series based on global RMS differences (2.53 cm EWH for 2002–2015 and 2.43 cm EWH for 2004–2010) and absolute comparisons with GNSS-derived vertical displacements (2.67 mm for 2002–2015 and 2.51 mm for 2004–2010 in case of TUG; 2.67 and 2.57 mm, respectively, for GFZ). As a threshold the 2002–2015 EWH RMS difference between monthly ITSG-Grace2016 and CSR RL05 solutions is 1.9–2.8 cm depending on the applied filter. The comparison of CSR RL05 with GNSS-derived vertical displacements results in 2.49 mm for 2002–2015 and 2.32 mm for 2004–2010 (Flechtner et al.2017), which demonstrates the excellent performance of both daily time-series.
An operational test run of the NRT service was foreseen in the final phase of the EGSIEM project between 2017 April and September. However, the deteriorating health of GRACE-B posed a serious challenge in GRACE data processing and is the limiting factor in the quality of the gravity field solutions. Following a battery cell failure on 2016 October 25, the accelerometer aboard GRACE-B was turned off. Its data stream has subsequently been replaced with a ‘transplant’ data product, where the accelerometer measurements of GRACE-A are shifted in time and rotated to substitute the missing measurements on GRACE-B. Furthermore, to shed load on the GRACE batteries following the cell failure, intersatellite K-Band ranging (KBR) data were only collected in orbit segments where the satellites were fully exposed to the sun. Since 2011 the satellites were subject to regular position switching maneuvers to conserve fuel on GRACE-B. These events are correlated with the β angle between the orbit plane and the Earth–Sun direction and have a periodicity of 161 d. When β is smaller than approximately 15°, the intersatellite link was switched off and the KBR data are completely missing.
Following this procedure, intersatellite observations between the GRACE spacecraft were once more available since 2017 March 17. To alleviate the accelerometer transplant process, the pitch angle of both spacecraft relative to the line of sight was increased from zero to approximately one degree (Himanshu Save, Gerhard Kruizinga, personal communication 2017) on March 30 (Fig. 8, left-hand panel). The resulting increased KBR antenna phase centre correction (PCC) magnitude (Fig. 8, right-hand panel) causes errors in spacecraft attitude more prominently propagating into the gravity field solutions and manifest in east-west striping patterns.
Pointing angle (pitch component) of GRACE-A relative to the line of sight and antenna phase centre corrections (PCC, expressed as range-rate) during the operational NRT service run.
KBR data collection was steadily increased from sunlight-only orbital segments (approximately 40 per cent reduction of observation count) to full revolutions. On May 2 the accelerometer on GRACE-B was switched on again for 22 d at the expense of the amount of KBR data collected during this period. Nevertheless during this time span the nominal set of GRACE science data was available.
The operational phase of the NRT-Service started on April 1 with both GFZ and TUG producing daily GRACE solutions in an automated manner until 2017 June 29, when the instruments have been switched-off again. In fact these were the last gravity measurements of GRACE and the mission was officially declared ended on 2017 October 27.5
The problems with KBR antenna PCCs described above impacted the quality of the TUG and GFZ models in a different manner. While GFZ relies fully on JPL’s PCCs provided within the Level-1B data records, TUG calculates their own correction using a smoothed attitude product based on a combination of star camera observations with angular accelerations (Klinger et al.2014). This combination acts as a low pass filter and reduces high frequency noise in the derived PCC’s compared to the Level-1B data set. The discrepancies between the solutions become obvious when looking at daily total water storage maps (Fig. 9) and derived flood indices (Fig. 10) within the operational test run (e.g. 2017 June 5) and a daily solution at ‘normal’ satellite conditions (e.g. 2008 June 5). In terms of difference RMS over the continents, TUG and GFZ solutions differ approximately 2.8 cm during the NRT service run, which is an increase of about 30 per cent compared to the same three months time span in 2008. Nevertheless, from the historical data it can be concluded that both, the gravity maps (9, top) and the flood indices (Fig. 10, top), look very similar for 2008 and again demonstrate the quality of the TUG and GFZ solutions. Operational latency of the NRT gravity field solutions at TUG is on average below 1 d from the last data epoch collected to the upload of the final solution. GFZ solutions needed about 1 additional day, but in general the latency is also well below the requirement of 5 d.
Total water storage anomaly derived from TUG (left-hand panel) and GFZ (right-hand panel) for 2008 June 5 (top) and 2017 June 5 (bottom). A three month mean around the evaluation epoch is subtraced from both time-series.
Wetness indices for TUG (left-hand panel) and GFZ (right-hand panel) NRT gravity products for 2008 June 5 (top) and 2017 June 5 (bottom).
6 HYDROLOGICAL SERVICE
Fig. 10 shows the above defined WI derived from the TUG and GFZ NRT gravity field solutions presented in Section 5. Similar to the water storage maps in Fig. 9 (top), the geographical patterns of the WI look similar for the two NRT solutions derived from historical data. Fig. 10 (top) shows exceptionally high (Central South America, Siberia) and low (Central Australia, Turkey, Central Asia) values of the wetness index. In contrast, for deficient satellite conditions during the operational NRT test run, the global WI patterns derived from the two solutions considerably differ (Fig. 10, bottom) as also the daily water storage maps differ from each other (Fig. 9, bottom).
6.1 Test case Danube basin
Fig. 11 shows a retrospective analysis of daily global gravity solutions for the years 2002–2015 for the Danube basin. Both the gravity solutions from TUG and GFZ and the WIs derived from these solutions indicate increased values in particular for the flood events in the Danube basin in 2002, 2006 and 2010, and less pronounced for the flood events in 2005, 2009, 2013 and 2014. Particularly relevant with respect to early flood warning is the build-up of basin-wide water storage of several weeks duration prior to the larger flood events of 2006 and 2010, which were triggered by a combination of (early season) snowmelt (resulting from unusual early high temperature peaks in 2006) and intense (2006) or excessive (2010) rainfall (e.g. ICPDR 2007, 2012).
Basin-average daily water storage anomalies (TWSA), gravity-derived WIs and river discharge at Ceatal Izmail gauging station at the basin outlet of the Danube basin (2002–2015) for (a) TUG and (b) GFZ, (c) for the combined WI. The vertical lines indicate the days of peak discharges for flood events.
Fig. 12 illustrates the flood and WI dynamics in more detail for the 2006 and 2010 events. River flow at the Ceatal Izmail station at the outlet of the Danube Basin peaks on 2006 April 26. Elevated WI values of about 2, corresponding to the 95 per cent percentile of WI values in the Danube basin, were recorded from early March onwards, possibly indicating a level close to basin storage saturation. A first WI peak over the presumed threshold of 2 on March 14, marked as a (green) vertical line in Fig. 12 (left-hand panel) provides a lead-time to river peak flow at the basin outlet of 43 d (or over 6 weeks). This lead time exceeds the traveltime of flood waves from the upstream to the downstream parts of the basin which is in the order of 3–4 weeks. Fig. 12 (right-hand panel) illustrates the setting for the 2010 Danube flood. River discharge at the basin outlet peaks on July 6. WI values clearly exceeding the indicative threshold of 2 are recorded from May 30 onwards. This leads to a flood warning lead-time prior to peak flow at the basin outlet of 37 d (or over 5 weeks).
Basin-average gravity-derived WIs for the Danube basin for the 2006 (left-hand panel) and the 2010 (right-hand panel) flood together with discharge at the basin outlet, gauging station Ceatal Izmail. The vertical dashed lines indicate time of peak flow (cyan) and early flood warning (green) based on WI exceeding a threshold value of 2.
A closer look at the figures discussed above reveals different dynamics of the two WIs derived from the individual daily gravity solutions. While the ITSG-derived WI dominates (i.e. generates the more extreme WI values) for the 2006 flood (Fig. 12, left-hand panel), the GFZ-derived WI does so for the 2010 flood (Fig. 12, right-hand panel). This tentatively points to a strength of the combined WI, which selects the most extreme absolute value. Interestingly, in the longer term perspective (2002–2015), the similarity of the two WIs is striking, despite differences in dynamics and/or noise levels of the individual daily gravity solutions they are derived from (Fig. 11c). Here, apart from the larger floods in 2006 and 2010 (exceeding the indicative WI threshold of 2), smaller and shorter duration floods from 2013 and 2014 are detected with shorter lead times, possibly indicating a different flood triggering mechanism with less water storage accumulation. A low of the wetness index reflects dry conditions during the 2003 and 2015 European heat waves. A gradual increase of the WI from 2003 to 2006 can be an indication of the persistence of low water storage as a consequence of the drought event, and a reason why the floods in 2004 and 2005 are not captured by the gravity-based index.
Fig. 13 shows maps of total water storage anomalies and of the WIs for the Danube basin on 2006 April 26 when flow peaks at Ceatal Izmail station at the basin outlet. Increased gravity TWSA values (dark blue) in both gravity solutions are reflected in elevated WI values, in particular in the lower parts of the river basin, with a spatially focused pattern of high WI by the GFZ gravity solution along the downstream reaches of the Danube where major flooding occurred. For the example of the Achleiten gauging station located at the German–Austrian border, the Upper Danube sub-basin size of just under 80 000 km2 is at the limit of what is a physically detectable gravity signal by GRACE. In fact, lead-times prior to peak flow in the Upper Danube are critically reduced (2006) or even negative (2010). Thus, the potential of early flood warning by a detectable WI signal of elevated storage is markedly reduced for smaller upstream sub-basins.
Daily gravity-based total water storage anomalies (TWSA) and gravity-based Wetness Indices (WI) for the Danube basin for 2006 April 26 based on TUG gravity solution (a, b), GFZ gravity solution (c,d) as well as the combined WI (e). The red cross indicates the location of the gauging station Ceatal Izmail.
6.2 Potential for early warning applications
Early warning by the wetness index presented here is expected to improve the programming and the efficient use of high resolution satellites that are used for disaster management activities such those carried out at DLR’s Center for Satellite-Based Crisis Information (ZKI), being one of the main value adders of the International Charter (Voigt et al.2007; Martinis et al.2017). The International Charter ‘Space and Major Disasters’, which is an agreement of space agencies and satellite operators worldwide with the aim of providing a uniform system for the quick acquisition and delivery of EO data in case of disasters, is a major international activity in this respect (Voigt et al.2016). Generally, the Charter mechanism is activated through authorized users. In some cases in the past, when user requests came in relatively late, satellite tasking could not be put into effect until the flood peak had already passed the area of interest. Here, for the example of the Danube flood in 2006, we assess the value of early proactive satellite tasking based on external information from gravity-based wetness indicators before the actual activation of the Charter. As an example, Fig. 14 shows a map of flooded areas in the delta region of the river Danube in eastern Romania which was derived from MODIS (Moderate Resolution Imaging Spectroradiometer) data acquired on 2006 April 26. This corresponds to the time of the flood peak at Ceatal Izmail station at the outlet of the Danube basin (see Fig. 11). MODIS is a medium resolution optical satellite sensor with a spatial resolution of up to 250 m used for large-scale environmental monitoring (Masuoka et al. 1998). MODIS acquires images of the Earth’s surface continuously every 1–2 d, provided that there are no clouds. Higher resolution satellites (particularly SAR satellites), which are more suitable for flood delineation even for cloudy conditions, need to be programmed prior to the flood peak in order to acquire up-to-date high resolution images of the flood affected regions to be used for disaster management activities. For the 2006 flood, the International Charter was activated for Romania on 2006 April 18 (Call No. 121, Table 4). The area affected in the Danube delta (see Fig. 14) was one of the two areas of interest (AOI) of the Charter call for which a number of rapid mapping products were produced based on the Charter satellite data, shown in Fig. 15, apart from the routinely available MODIS data acquired on 2006 April 26.
Map of flooded areas in the Danube delta based on MODIS data acquired on 2006 April 26 (source: https://disasterscharter.org/web/guest/-/flooding-in-romania).
Possible overpasses and acquisitions of ENVISAT ASAR (wide swath mode with 150 m spatial resolution), Radarsat-1 (standard mode with 30 m spatial resolution), SPOT-4 (HRVIR sensor with 10 m spatial resolution) and SPOT-5 satellite (HRG sensor with 10 m spatial resolution) for the Danube delta region shown in Fig. 14. Acquisitions that have been realized during Charter Call 121 are marked with dark green symbols.
Satellite data acquisitions for the Danube delta in Romania during Charter Call 121.
| Satellite | Sensor | Spatial Resolution (m) | Date |
|---|---|---|---|
| SPOT-4 | Optical | 10 | 23.04.2006 |
| ERS-2 | SAR | 30 | 23.04.2006 |
| Radarsat-1 | SAR | 30 | 25.04.2006 |
| ENVISAT-A SAR | SAR | 30 | 28.04.2006 |
| Satellite | Sensor | Spatial Resolution (m) | Date |
|---|---|---|---|
| SPOT-4 | Optical | 10 | 23.04.2006 |
| ERS-2 | SAR | 30 | 23.04.2006 |
| Radarsat-1 | SAR | 30 | 25.04.2006 |
| ENVISAT-A SAR | SAR | 30 | 28.04.2006 |
Satellite data acquisitions for the Danube delta in Romania during Charter Call 121.
| Satellite | Sensor | Spatial Resolution (m) | Date |
|---|---|---|---|
| SPOT-4 | Optical | 10 | 23.04.2006 |
| ERS-2 | SAR | 30 | 23.04.2006 |
| Radarsat-1 | SAR | 30 | 25.04.2006 |
| ENVISAT-A SAR | SAR | 30 | 28.04.2006 |
| Satellite | Sensor | Spatial Resolution (m) | Date |
|---|---|---|---|
| SPOT-4 | Optical | 10 | 23.04.2006 |
| ERS-2 | SAR | 30 | 23.04.2006 |
| Radarsat-1 | SAR | 30 | 25.04.2006 |
| ENVISAT-A SAR | SAR | 30 | 28.04.2006 |
As an example, Fig. 15 shows four of the Charter satellites, for which acquisitions would have been possible for the Danube delta area (Fig. 14) based on the available overpasses during April 2006. This analysis has been carried out with the ‘SaVoir Charter - Swath Acquisition Planner’ V4.5.4.0 (©Taitus Software).
Taking into account the lead time of the gravity-based wetness index of 43 d prior to the flood peak at Ceatal Izmail station at the outlet of the Danube basin (see Fig. 12), a large number of additional satellite overpasses (35 for ENVISAT-ASAR, 32 for Radarsat-1, 26 for SPOT-4 and 27 for SPOT-5) would have been possible to acquire data for flood monitoring and emergency mapping if the wetness indices were available at that time. For the Danube flood 2006 the Charter was activated four times (Czech Republic/Slovakia: April 1; Austria: April 7; Hungary: April 14; Romania: April 18). The WI Early flood warning time stamp of the Upper Danube was 2006 March 13. Regarding the time when the Charter was first activated for the Danube flood on April 1 for the Czech Republic and Slovakia being the point in time when it was considered as a major disaster, an early flood warning lead time of 18 d can be determined in this case.
Requirements expressed by the users of satellite rapid mapping products focus on timely and high frequency flood monitoring from the onset of a flood event with special focus on mapping the flood extent at peak level until water levels have receded to near normal stages. Based on the example of the Danube flood 2006 it could be demonstrated that gravity-based early-warning indicators for total water storage anomalies can improve the programming and the efficient use of EO satellites for rapid flood mapping tremendously. For the case of widespread flooding along the Danube in southern Romania after a dam break close to the village of Bistret on 2006 April 17 (not shown in detail here), proactive early satellite tasking would have offered the possibility to have satellite data available for disaster response teams a few hours after the dam break instead of a few days.
7 SUMMARY
The European GRACE ACs have joined forces to establish the European Gravity Service for Improved Emergency Management (EGSIEM) initiative. During the European Commission backed funding phase of 2015–2017, three prototype services were set up: (1) a scientific combination service, (2) an NRT service, and (3) a hydrological/early warning service. With the scientific combination service EGSIEM aimed to pave the way to a long-awaited standardization in the gravity field community. For a demonstration period of 2 yr normal equations of monthly GRACE gravity fields in terms of spherical harmonics were generated in the SINEX format by four ACs (AIUB, TUG, GFZ, CNES) by adopting different approaches but agreed upon processing standards. The resulting SINEX files were generated without regularizations and combined on the normal equation level to derive one consolidated monthly gravity field solution with significantly increased quality, robustness and reliability. With the end of the EGSIEM project in 2017 December, the scientific combination service is currently being extended to also include non-European ACs, to extend the time-series of combined monthly gravity fields covering the entire GRACE time period, and to be prepared for the upcoming release of GRACE-FO data. The service will be continued as the International Combination Service for Time-variable Gravity Fields (COST-G). COST-G will be officially launched as the product centre for time-variable gravity fields of IAG’s International Gravity Field Service (IGFS) in 2019 July at the 27th General Assembly of the International Union of Geodesy and Geophysics (IUGG).
Within the NRT service, EGSIEM generated for the first time daily gravity field solutions in NRT using two different approaches at the ACs TUG and GFZ. The performance of the daily gravity solutions was tested using historical hydrological extreme events and relative intercomparisons between the two generated NRT time-series. An excellent agreement of the two solutions was found for these historical flood events, which was the proof-of-concept of this service. The operational test run took place from 2017 April 1 until June 29 but the deteriorating health of GRACE-B posed a serious challenge to the data processing. In terms of difference RMS over the continents, the two solutions differed approximately 2.8 cm during the operational test run, which is an increase of about 30 per cent compared to a similar three month time span in the offline test for 2008. Thus, the instrument performance was found to be the limiting factor in the quality of the gravity field solutions but we found no limiting factor in the procedure and workflow of the NRT service. The test can therefore be considered a success and the service is prepared for the soon to be released GRACE-FO data.
Within the hydrological service a wetness index was developed indicating departure of the GRACE-derived total water storage anomaly from a mean seasonal cycle (after removing the long-term trend). A retrospective analysis of both daily gravity field solutions showed increased values of the wetness index for the flood events in the Danube basin in 2002, 2006 and 2010. Of particular relevance with respect to early flood warning is the build-up of basin-wide water storage over several weeks prior to the larger flood events of 2006 and 2010. Elevated wetness indices of about two were registered with a lead-time of more than 6 and 5 weeks before the river peak flow at the basin outlet for the 2006 and 2010 floods, respectively, significantly exceeding the traveltime of flood waves from the upstream to the downstream parts of the basin. The gravity-based wetness index of the hydrological service has been incorporated into existing mapping and alerting systems for testing purposes. In particular, it has been included as a data layer into the European Commission’s Copernicus Global Flood Awareness System (GloFAS) platform. At DLR/ZKI, an interactive web client has been developed to visualize the wetness index together with other DLR/ZKI data sources such as with results from DLR’s operational Sentinel-1 (Twele et al.2016) and TerraSAR-X Flood Services (Martinis et al.2015). Being incorporated into the operational ZKI workflow the daily global gravity-based wetness index may serve as one of the triggers for an early and improved satellite tasking that offers high resolution (e.g. TerraSAR-X) acquisition planning for supporting disaster response and disaster management.
ACKNOWLEDGEMENTS
This research was supported by the European Union‘s Horizon 2020 research and innovation program under the grant agreement No. 637010 and the Swiss State Secretariat for Education, Research and Innovation. All views expressed are those of the authors and not of the Agency. Discharge station data are kindly provided by the Global Runoff Data Centre, 56068 Koblenz, Germany, and by the National Institute of Hydrology and Water Management, Bucharest, Romania.
Footnotes
The US Tellus website http://grace.jpl.nasa.gov/data/monthly-mass-grids is only dedicated to monthly gravity field models of the GRACE SDS ACs
