Quantification of drainable water storage volumes on landmasses and in river networks based on GRACE and river runoff using a cascaded storage approach – first application on the Amazon

The knowledge of water storage volumes in catchments and in river networks leading to river discharge is essential for the description of river ecology, the prediction of floods and specifically for a sustainable management of water resources in the context of climate change. Measurements of mass variations by the GRACE gravity satellite or by ground-based observations of river or groundwater level variations do not permit the determination of the respective storage volumes, which could be considerably bigger than the mass variations themselves.

For fully humid tropical conditions like the Amazon the relationship between GRACE and river discharge is linear with a phase shift. This permits the hydraulic time constant to be determined and thus the total drainable storage directly from observed runoff can be quantified, if the phase shift can be interpreted as the river time lag. As a time lag can be described by a storage cascade, a lumped conceptual model with cascaded storages for the catchment and river network is set up here with individual hydraulic time constants and mathematically solved by piecewise analytical solutions.

Tests of the scheme with synthetic recharge time series show that a parameter optimization either versus mass anomalies or runoff reproduces the time constants for both the catchment and the river network τ_C and τ_R in a unique way, and this then permits an individual quantification of the respective storage volumes. The application to the full Amazon basin leads to a very good fitting performance for total mass, river runoff and their phasing (Nash–Sutcliffe for signals 0.96, for monthly residuals 0.72). The calculated river network mass highly correlates (0.96 for signals, 0.76 for monthly residuals) with the observed flood area from GIEMS and corresponds to observed flood volumes.

The fitting performance versus GRACE permits river runoff and drainable storage volumes to be determined from recharge and GRACE exclusively, i.e. even for ungauged catchments. An adjustment of the hydraulic time constants (τ_C, τ_R) on a training period facilitates a simple determination of drainable storage volumes for other times directly from measured river discharge and/or GRACE and thus a closure of data gaps without the necessity of further model runs.

1 Introduction

In the context of water resources management and climate change there is an ongoing discussion on how to assess available water resources, i.e. the storage volumes which can be used for water supply in a dynamic way beyond the limitations of sustainable extraction rates. The maximum average extraction rate for a sustainable use of water resources is limited by the long-term recharge of a catchment (Sophocleous 1997; Bredehoeft, 1997); however, this rate-based definition of groundwater stress only allows an assessment of water resources with respect to long-term sustainability and does not permit short-term management in order to satisfy specific water demands. Thus the knowledge of water resources involved in the water cycle contributing to river discharge, such as parts of the groundwater or surface water system, is essential.

Very little attention has so far been given to the quantification of the storage volumes of renewable water resources participating in the dynamic water cycle driven by precipitation P, actual evapotranspiration ET_a and river runoff R. The reason for this is seen in the problem that observations of time-variant groundwater or river levels only permit the estimation of volume changes but not absolute storage volumes, which could be considerably bigger.

Natural systems consist of many different storage components like canopy, snow/ice, surface, soil, unsaturated/saturated underground, drainage system etc. Direct measurements of storage volumes from water or pressure levels are problematic as they are based on assumptions and approximations. They are based on point measurements and quite rare on large spatial scales compared to the heterogeneity scale of the respective compartments. This leads to large interpolation errors. In addition, the storage coefficients for porous media describing the relationship between the measurable groundwater heads or capillary pressure on the one hand, and storage volume or absolute soil saturation on the other hand, are insufficiently known on large scales. Remote sensing data have been limited to near-surface water storage (open water bodies, soil) up until now and are thus of limited benefit for the quantification of water storage with respect to accuracy and coverage due to methodological constraints (Schlesinger, 2007).

In contrast to discharge-less basins and/or arid areas, which are nearly exclusively driven by precipitation and evapotranspiration, the storage dynamics of catchments draining into a river system allows the hydraulically coupled storage compartments to be addressed via their contributions to river discharge. These comprise groundwater, surface water, the river network and temporarily inundated areas. All storages draining into the river system by gravity are referred to as “drainable” storage here. So, aquifers or parts of them not draining into the river system without an energy input are not considered here.

River runoff $R\left(t\right)=Q\left(t\right)/A$ (corresponding to river discharge Q(t) from the related catchment area A) is driven by the storage height or mass density $M_{\mathrm{Storage}}=V_{\mathrm{tot}}/A$ of all superposed hydraulically coupled storages and is determined by their runoff–storage (R–S) relationship. For time periods with no recharge or losses of water (as by ET_a), i.e. with no processes affecting M_Storage other than river discharge Q, a linear runoff–storage (R–S) relationship $R\left(M\right)=M/\mathit{\tau }$ leads to an exponential decrease in river discharge or streamflow Q(t) depending on the related hydraulic time constant τ:

$$\begin{array}{}\text{(1)}& Q\left(t\right)=Q\left(t_{\mathrm{0}}\right)\cdot e^{-\cdot \frac{t-t_{\mathrm{0}}.}{\mathit{\tau }}}\end{array}$$

For this case the corresponding total drainable storage in terms of mass density M_Storage at any given time t₀ can be determined by an infinite temporal integration over river discharge Q(t) from the corresponding catchment area A starting at time t₀:

$$\begin{array}{}\text{(2)}& \begin{array}{rl}{M}_{\mathrm{storage}}\left(t_{\mathrm{0}}\right)& ={\displaystyle \frac{V_{\mathrm{tot}}\left(t_{\mathrm{0}}\right)}{A}}={\displaystyle \frac{\mathrm{1}}{A}}\cdot \underset{t_{\mathrm{0}}}{\overset{\mathrm{\infty }}{\int }}Q\left(t\right)\mathrm{d}t={\displaystyle \frac{Q\left(t_{\mathrm{0}}\right)}{A}}\cdot \underset{t_{\mathrm{0}}}{\overset{\mathrm{\infty }}{\int }}\\ & e^{-\cdot \frac{t-t_{\mathrm{0}}}{\mathit{\tau }}}\mathrm{d}t=\mathit{\tau }\cdot {\displaystyle \frac{Q\left(t_{\mathrm{0}}\right)}{A}}=\mathit{\tau }\cdot R\left(t_{\mathrm{0}}\right)\end{array}\end{array}$$

Contributions from several storage compartments (with individual time constants) superpose, if they drain in parallel and if there is no feedback from the river system. For this case, there is a wide range of time series analysis methods (Tallaksen, 1995), which allow the flow components to be separated into fast, medium or slow and the corresponding surface, interflow or groundwater flow contributions according to their individual time constants. Thus, measurements of the different time constants allow the drainable storage of the respective storage compartment and the corresponding mean drainable storage to be determined from mean runoff R or recharge N:

$$\begin{array}{}\text{(3)}& {\stackrel{\mathrm{‾}}{M}}^X=\stackrel{\mathrm{‾}}{R}\cdot {\mathit{\tau }}_X=\stackrel{\mathrm{‾}}{N}\cdot {\mathit{\tau }}_X.\end{array}$$

On global scales the absolute storage volume of the drainable storages can be determined from runoff time series directly, if there are distinct and long enough periods of negligible or even negative recharge (actual evapotranspiration ET_a > precipitation) as it occurs in seasonally dry regions (Niger, Mekong, some Amazon sub-catchments etc.). From the purely exponential decrease in river discharge the time constant can be determined directly from a curve fit, as shown in Fig. 1b for Amazon sub-catchments. If the dry period is long enough the sequence of different time constants taken from the discharge curve even permits a discrimination between the fast response by overland flow and the slow response by the groundwater system.

Figure 1 (a) R–S diagram with counter clockwise hysteresis for mean monthly observed runoff R_o versus GRACE dM for fully humid catchments including a phase adaption for the Amazon upstream from Obidos (Riegger and Tourian, 2014). (b) R–S diagram with clockwise hysteresis for mean monthly observed runoff R_o versus GRACE anomaly dM for seasonally dry catchments in the Amazon basin (Riegger and Tourian, 2014). (c) Mean monthly runoff R and recharge N for fully humid catchments in the Amazon basin (log scale for R). (d) Mean monthly runoff R and recharge N for seasonally dry catchments in the Amazon basin (log scale for R) including exponential fittings for runoff R_sim.

Catchments with permanent input, i.e. no periods of negligible recharge, however, do not show an exponential behaviour for discharge. For these cases the hydraulic time constant cannot be taken from discharge dynamics directly, but has to be estimated by hydrological models. These are intended to describe the large number of storages distributed over the catchment by the assumed processes and calibrate the involved parameters by their respective superposed flows versus the observed river discharge. The main difficulties in verifying large or global-scale hydrological models or land surface models (GHMs or LSMs) consist of the quantification of local individual storage volumes and related flows by local ground-based measurements. Thus, even though distributed hydrological models very much support an understanding of processes in the water cycle, the limitation of the calibration versus river discharge exclusively introduces an ambiguity in the impact of contributing processes and the related storages and flows.

Since 2013 GRACE observations of the time-variable gravity field provide monthly distributions of mass density on large spatial scales > ∼200 000 km² (Tapley et al., 2004). However, as the water storage in different compartments (snow, ice, vegetation, soil, surface water, ground water etc.) superposes with all other terrestrial (geophysical) masses, only the time-variant part of the GRACE signal can be used to quantify the terrestrial water storage (TWS) anomalies (monthly mass signals minus long-term average), but not the related absolute storage volumes. Nevertheless, this for the first time permits a direct comparison of measured TWS and observed river runoff R_o. Surprisingly some GHMs showed a considerable phase shift between measured mass anomalies by GRACE and river discharge as well as between calculated and measured runoff and an underestimation of mass signal amplitudes (Güntner et al., 2007; Chen et al., 2007; Schmidt et al., 2008; Werth et al., 2009; Werth and Güntner, 2010) even though they comprise a large number of storages like soil, surface water, groundwater etc. This is emphasized by Scanlon et al. (2019), who for tropical basins recognize the main cause of the discrepancies in insufficient storage capacity and lack of surface water inundation.

The direct comparison of GRACE anomalies and river runoff on large spatial and monthly timescales by Riegger and Tourian (2014) revealed that measured runoff–storage (R–S) diagrams show hysteresis curves of distinct form and extent (Fig. 1a, b), which are characteristic for different climatic conditions (like fully humid, seasonally dry or boreal) and can be explained by considering recharge and runoff properties (Fig. 1c, d).

Thus for example, catchments in fully humid conditions (like the full Amazon basin upstream from Obidos (295 in Fig. 1a) and some of its catchments like upstream Manacapuru (501 in Fig. 1c)) with a permanent input, i.e. only positive recharge (Fig. 1c), show a counterclockwise hysteresis (Fig. 1a). If this can be fully described by a positive phase shift, river runoff and storage behave like a linear time-invariant (LTI) system (Riegger and Tourian, 2014), i.e. the R–S relationship is linear, if the phase shift is adapted as shown in Fig. 1a. For this case the hysteresis can be purely assigned to a time lag. Once the phase shift is adapted the slope in the R–S diagram corresponds to the hydraulic time constant via $\mathit{\tau }=M/R$. The time constant and the reasonable assumption of a proportional R–S relationship (no runoff for empty storage) then facilitates the quantification of the drainable storage (Eq. 3), i.e. the volume related to the hydraulically coupled storage compartment, which drains by gravity.

In contrast, catchments with distinct periods of zero or negative recharge (like Niger, Mekong or Rio Branco (504), Rio Jurua (506) in the Amazon basin; Fig. 1b) show a clockwise hysteresis in the R–S diagram and a form which is determined by an increase in mass and runoff during wet periods, a decrease in mass and runoff with different slopes corresponding to different time constants and a possible mass loss without a related runoff (by negative recharge (Fig. 1d) by evapotranspiraton from the soil zone) during dry periods. This type of hysteresis is determined by storage changes not connected with river discharge (uncoupled) and cannot be explained by a time lag as it is not causal.

The consequence from the above discussion is that the determination of the hydraulic time constant and thus the drainable storage is only possible for catchments for which the hysteresis is fully explained by a positive phase shift; i.e. uncoupled storages are either negligible or can be separated from GRACE mass by other means (as shown below for boreal regions).

Based on this method, Tourian et al. (2018) apply an adaption of the phase shift using a Hilbert transform in order to determine the hydraulic time constants and the total drainable water storage for the sub-catchments of the Amazon basin without a consideration of the form of the R–S hysteresis. To be sure, this leads to reasonable results for the sub-catchments with permanent input (Fig. 1a, c) for which the time-dependent uncoupled storage is negligible. However, for Rio Branco (504) or Rio Jurua (506) this condition is not fulfilled as the hysteresis is determined by mass changes in the uncoupled storage and by runoff with different time constants (Fig. 1b, d). For these catchments the exclusive adjustment of the phase shift leads to negative time lags – which are not physical – and as a consequence to misleading time constants and thus to considerable errors in the determination of drainable storage volumes.

The accurate description of the R–S hysteresis of a catchment and its river network is the prerequisite for an accurate description of the system dynamics and the related storage volumes on the land masses (canopy, soil, overland flow, saturated/unsaturated underground) and in the river network.

Recent developments in river routing schemes of global hydrologic models with a hydrodynamic modelling of the flow in the river network system have successfully dealt with the description of phase shifts generated by the time lag in the river network (Paiva et al., 2013; Luo et al., 2017; Siqueira et al., 2018). Getirana et al. (2017a) emphasize the importance of integrating an adequate river routing schemes not only for an improved phase agreement with observed river discharge but also for an appropriate fit of the total mass amplitude to GRACE by the inclusion of the corresponding river network storage. Yet a hydrodynamic modelling of a complete river network system for the determination of the river network time lag and storage means a huge modelling effort (Getirana et al., 2017b).

A far more simple approach is presented by Riegger and Tourian (2014), describing the system by macroscopic variables, summarizing all coupled storage compartments on landmasses and in the river network, and analogously all uncoupled storage compartments in one respective single storage by their effect on the R–S relationship. The intention of such a “top-down” or lumped approach is to integrate the catchment-scale water balance and describe the system by large-scale variables and parameters, which are directly measurable or adjustable. For this purpose recharge based on moisture flux divergence or catchment water balance using GRACE can be used, which is quite accurate, yet limited to global scales (see below). Thus, opposite to distributed hydrological models which are based on spatially or temporally distributed data (for hydrometeorological input, local storage conditions in vegetation, soil and underground) and a detailed description of internal processes – which cannot be verified locally at present – this “top-down” approach uses measured catchment-scale input, storages and runoff. Where necessary and possible catchment-scale parameters are used to separate coupled and uncoupled storages (like MODIS snow coverage for boreal regions; Riegger and Tourian, 2014). In addition the time lag between storage and river discharge need not be explicitly described by an excessive routing scheme. Instead the related phase shift can be adapted by mathematical methods. This leads to a description of the system behaviour with high accuracy (Nash–Sutcliffe efficiency of 0.97 for the whole Amazon basin) by an adaption of only two parameters, the hydraulic timescale and the phase shift, even though the physical cause of the phase shift is not addressed explicitly.

A disadvantage of the above approaches (Riegger and Tourian, 2014; Tourian et al., 2018) is that it does not permit the individual drainable storage volumes on landmasses and in the river network to be quantified separately, but only the total drainable volume of the catchment. The information contained in the phase shift or time lag is not used for a quantification of the river network storage volume. Yet, as observations of inundated areas in river networks such as those from the GIEMS “Global Inundation Extent from Multi-Satellites” project Prigent et al. (2007); Papa et al. (2008); Papa et al. (2013) and hydrodynamic models of the river network (Paiva et al., 2013, Getirana et al., 2017b; Siqueira et al., 2018) indicate a considerable contribution of river network storage corresponding to a non negligible time lag, the river network storage must be considered in the integration of the total catchment water balance. As a sequence of storages (cascaded storages) leads to a time lag (i.e. a phase shift; Nash, 1957) and storages draining in parallel (as for overland and groundwater flow) just lead to a superposition (with no time lag), a storage cascade is considered as an appropriate description to account for a time lag.

This paper explores the accuracy and uniqueness of a lumped, top-down approach called a “cascaded storage” approach, which is based on the integration of given recharge in the water balance and utilizes a cascade of a catchment storage and a river network storage for a simple description of the observed time lag and the individual storage volumes. This permits a description of the system with a minimum number of macroscopic observation data and an adaption of only two parameters, the hydraulic time constants of the catchment and the river network. These time constants then could be used for nowcasts or even forecasts (within the time lag) of river discharge and/or drainable storage volumes directly from measurements without the need for further modelling.

The paper is structured as follows: Sect. 2 presents the mathematical framework of piecewise analytical solutions of the water balance equation for a cascade of catchment and river network storages. It also contains the description of observables, which permit the comparison of calculated and measured values. The “single storage” approach is handled as the specific case for a negligible river network time constant. In Sect. 3, the properties of the cascaded storage approach and its impact on the performance of the parameter optimization are described for synthetic recharge data and compared to the single storage approach. Based on the cascaded storage approach a fully data-driven approach is presented which permits a simplified determination of the drainable storage volumes directly from measurements without the need for further model runs. In Sect. 4 the approach is applied to data from the Amazon basin and evaluated versus measurements of GRACE mass, river runoff and flood area from GIEMS. The results are compared to GHM/LSM studies. In Sect. 5 the approach and its performance and limitation is discussed. Possible future investigations in order to overcome some of its limitations are sketched. Conclusions are drawn in Sect. 6.

2 Mathematical framework

In order to investigate the impact of a non negligible river water storage on the time lag in the river system, the water balance of the total system comprising both the catchment and river network storage has to be considered. A conceptual model corresponding to a Nash cascade (Nash, 1957), called a cascaded storage approach here, is set up with individual time constants for the different storages and with the following properties:

• Surface water and shallow groundwater storages on the land mass which are draining into the river network and are being fed by recharge are summarized to a so-called “catchment” storage M^C with time constant τ_C. Overland and groundwater flow from the land masses are summarized to a “Catchment” runoff R^R.

• River runoff (river discharge per catchment area), which addresses hydraulically the flow in the river channel network including inundated areas, is determined by its hydraulic time constant τ_R. The respective river network storage M^R is assumed to be instantaneously distributed within the river network system. Internal routing effects, which might lead to an additional delay in streamflow response, are not considered.

• Any possible hydraulic feedback from the river to the catchment system is assumed to be negligible.

• Temporal variations of uncoupled storage compartments like soil or open water bodies are considered as negligible.

These conditions are chosen for the sake of conceptual and mathematical simplicity. It has to be emphasized here that for a general applicability on a global coverage several coupled storages with different time constants and different uncoupled storage compartments with their respective time dependency have to be considered, of course. For fully tropical climatic conditions with permanent recharge, however (as for the full Amazon basin), variations in the soil water storage are negligible and the different dynamics of overland and groundwater flow cannot be distinguished. Thus, applications of this first approach are limited to catchments for which the hysteresis can be fully described by a time lag; i.e. no impacts of other coupled or uncoupled storages exist.

The abbreviations used throughout the paper are described in Appendix Table A1).

The total system behaviour is described by two balance equations, one for catchment storage (Eq. 4) and one for river storage (Eq. 6).

Catchment storage:

$$\begin{array}{}\text{(4)}& {\displaystyle \frac{\partial }{\partial t}}{M}^{\mathrm{C}}\left(t\right)=N\left(t\right)-R^{\mathrm{C}}\left(t\right)=N\left(t\right)-{\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{C}}}}\cdot M^{\mathrm{C}}\left(t\right)\end{array}$$

with

$$\begin{array}{}\text{(5)}& R^{\mathrm{C}}\left(t\right)={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{C}}}}\cdot M^{\mathrm{C}}\left(t\right).\end{array}$$

River storage:

$$\begin{array}{}\text{(6)}& {\displaystyle \frac{\partial }{\partial t}}{M}^{\mathrm{R}}\left(t\right)=R^{\mathrm{C}}\left(t\right)-R^{\mathrm{R}}\left(t\right)=R^{\mathrm{C}}\left(t\right)-{\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{R}}}}\cdot M^{\mathrm{R}}\left(t\right)\end{array}$$

with

$$\begin{array}{}\text{(7)}& R^{\mathrm{R}}\left(t\right)={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{R}}}}\cdot M^{\mathrm{R}}\left(t\right),\end{array}$$

with a proportional R–S relationship for hydraulically coupled storages. N denotes the recharge as input, R^C the catchment runoff from the catchment storage M^C, which cannot be measured directly on large spatial scales, and R^R the river runoff from the river network storage M^R, which can be measured at discharge gauging stations.

The water balance equation, Eq. (4), for the catchment is generally solved by the following:

$$\begin{array}{}\text{(8)}& M^{\mathrm{C}}(t-t_{\mathrm{0}})=M^{\mathrm{C}}\left(t_{\mathrm{0}}\right)\cdot e^{-\frac{t-t_{\mathrm{0}}}{{\mathit{\tau }}_{\mathrm{C}}}}+\underset{t_{\mathrm{0}}}{\overset{t}{\int }}N\left(w\right)\cdot e^{\frac{w-t}{{\mathit{\tau }}_{\mathrm{C}}}}\cdot \mathrm{d}w,\end{array}$$

where M^C(t₀) is the initial condition and N(t) the time-dependent recharge.

For recharge N(t) being given with a certain temporal resolution in time units or by periods of piecewise constant values and arbitrary length (stress periods) the recharge time series can be described as follows:

$$\begin{array}{}\text{(9)}& \begin{array}{rl}& N\left(t\right)=\sum _{i=\mathrm{0}}^{n-\mathrm{1}}{N}_{i+\mathrm{1}}\cdot {\mathit{\gamma }}_{i+\mathrm{1}}\left(t\right)\mathrm{with}{\mathit{\gamma }}_{i+\mathrm{1}}\left(t\right)=\left\{{}_{\mathrm{0}}^{\mathrm{1}}\right.\\ & \mathrm{for}\left\{{}_{t\notin \left[t_{,}{t}_{\mathrm{1}+\mathrm{1}}\right]}^{t\in \left[t_{i,}{t}_{i+\mathrm{1}}\right]}\mathrm{for}\mathrm{each}\mathrm{interval}\left[t_{i,}{t}_{i+\mathrm{1}}\right]\right.\end{array}\end{array}$$

For calculation convenience Eq. (8) can be solved successively for each stress period using the values at the end of the last period as the starting value, which leads to the piecewise analytical solution for catchment mass for a time $t\in \left[t_{i,}{t}_{i+\mathrm{1}}\right]$ in stress period _i+1:

$$\begin{array}{}\text{(10)}& M_{i+\mathrm{1}}^{\mathrm{C}}(t-t_i)=M_i^{\mathrm{C}}\left(t_i\right)\cdot e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{C}}}}+N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{C}}\cdot \left(\mathrm{1}-e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{C}}}}\right).\end{array}$$

The respective catchment runoff R^C based on Eq. (5) and M^C from Eq. (10) is used as input for the river network water balance, Eq. (6), and leads to the general solution for the river network storage M^R:

$$\begin{array}{}\text{(11)}& M^{\mathrm{R}}(t-t_{\mathrm{0}})=M^{\mathrm{R}}\left(t_{\mathrm{0}}\right)\cdot e^{-\frac{t-t_{\mathrm{0}}}{{\mathit{\tau }}_{\mathrm{R}}}}+\underset{t_{\mathrm{0}}}{\overset{t}{\int }}{R}^{\mathrm{C}}\left(u\right)\cdot e^{\frac{u-t}{{\mathit{\tau }}_{\mathrm{R}}}}\cdot \mathrm{d}u,\end{array}$$

and the iterative solutions for time $t\in \left[t_{i,}{t}_{i+\mathrm{1}}\right]$ in stress period _i+1:

$$\begin{array}{}\text{(12)}& \begin{array}{rl}& M_{i+\mathrm{1}}^{\mathrm{R}}(t-t_i)=M_i^{\mathrm{R}}\left(t_i\right)\cdot e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{R}}}}+N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{R}}\cdot \left(\mathrm{1}-e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{R}}}}\right)\\ & +\left[M_i^{\mathrm{C}}\left(t_i\right)-N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{C}}\right]\cdot {\displaystyle \frac{{\mathit{\tau }}_{\mathrm{R}}}{{\mathit{\tau }}_{\mathrm{C}}-{\mathit{\tau }}_{\mathrm{R}}}}\cdot \left(e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{C}}}}-e^{-\frac{t-t_i}{{\mathit{\tau }}_{\mathrm{R}}}}\right).\end{array}\end{array}$$

The total mass M^T is then given by the following:

$$\begin{array}{}\text{(13)}& M_i^{\mathrm{T}}=M_i^{\mathrm{C}}+M_i^{\mathrm{R}}.\end{array}$$

The mixed term in Eq. (12) and thus the total mass are commutative in (τ_C, τ_R) and show a singularity at τ_C=τ_R with an asymptotic value. For τ_R > τ_C solutions also exist with analogous values in total mass M^T for M^R > M^C.

It has to be emphasized here that the piecewise analytical solutions for time periods of constant recharge provide a mathematical solution for an arbitrary temporal resolution without numerical limitations. Finite difference solutions are limited by stability criteria ($t_{i+\mathrm{1}}-t_i$) < τ and accuracy criteria ($t_{i+\mathrm{1}}-t_i$) < τ∕10 for the smallest τ. Analytical solutions facilitate an exact calculation of the response of the river network during the time interval of constant recharge (though the time constant of the river network could be much shorter than the time interval or the time constant of the catchment). Thus the very high temporal discretization, which otherwise would be needed using a finite difference scheme, is avoided.

The observables related to measurements by GRACE and discharge from gauging stations are the total mass anomaly dM^T and the river runoff R^R. GRACE observations with acceptable error are still limited to monthly resolution. Discharge as well as some of the meteorological variables like precipitation, evapotranspiration or moisture flux divergence are often measured in daily values, and some of the products are measured in monthly values. For an optimal adaption to the monthly resolution of GRACE products, the approach presented here is based on monthly values but could also be applied to daily data without problems.

The mass values used in the calculations here are assigned to the interval boundaries while the values for monthly recharge and measured runoff are constant over the interval and temporally assigned to the centre of the interval. Thus, for a comparison of the calculated mass and runoff values versus the observed monthly values of GRACE and discharge, the calculated values have to be averaged over the interval. As the dynamics follow an exponential behaviour the mean values cannot be taken from arithmetic averages at the interval boundaries but instead from an integral average over the interval.

The mean storage mass for M_X is given for each interval $\left[t_{i,}{t}_{i+\mathrm{1}}\right]$ by the following:

$$\begin{array}{}\text{(14)}& {\stackrel{\mathrm{‾}}{M}}_{i+\mathrm{1}}^X={\displaystyle \frac{\mathrm{1}}{t_{i+\mathrm{1}}-t_i}}\underset{t_i}{\overset{t_{i+\mathrm{1}}}{\int }}{M}_{i+\mathrm{1}}^X(t-t_i)\cdot \mathrm{d}t,\end{array}$$

leading to mean runoff

$$\begin{array}{}\text{(15)}& {\stackrel{\mathrm{‾}}{R}}^X\left(t\right)={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_X}}\cdot {\stackrel{\mathrm{‾}}{M}}^X\left(t\right),\end{array}$$

i.e. mean catchment mass and runoff:

$$\begin{array}{}\text{(16)}& \begin{array}{rl}{\stackrel{\mathrm{‾}}{M}}_{i+\mathrm{1}}^{\mathrm{C}}& =(M_i^{\mathrm{C}}-N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{C}})\cdot {\displaystyle \frac{{\mathit{\tau }}_{\mathrm{C}}}{(t_{i+\mathrm{1}}-t_i)}}\left(\mathrm{1}-e^{-\frac{t_{i+\mathrm{1}}-t_i}{{\mathit{\tau }}_{\mathrm{C}}}}\right)\\ & +N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{C}}\end{array}\end{array}$$

and

$$\begin{array}{}\text{(17)}& {\stackrel{\mathrm{‾}}{R}}_i^{\mathrm{C}}={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{C}}}}\cdot {\stackrel{\mathrm{‾}}{M}}_i^{\mathrm{C}}\end{array}$$

and mean river mass and runoff:

$$\begin{array}{}\text{(18)}& \begin{array}{rl}{\stackrel{\mathrm{‾}}{M}}_{i+\mathrm{1}}^{\mathrm{R}}& =(M_i^{\mathrm{R}}-N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{R}})\cdot {\displaystyle \frac{{\mathit{\tau }}_{\mathrm{R}}}{(t_{i+\mathrm{1}}-t_i)}}\cdot \left(\mathrm{1}-e^{-\frac{t_{i+\mathrm{1}}-t_i}{{\mathit{\tau }}_{\mathrm{R}}}}\right)\\ & +N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{R}}+{\displaystyle \frac{\left[M_i^{\mathrm{C}}-N_{i+\mathrm{1}}\cdot {\mathit{\tau }}_{\mathrm{C}}\right]}{(t_{i+\mathrm{1}}-t_i)}}\cdot {\displaystyle \frac{{\mathit{\tau }}_{\mathrm{R}}}{{\mathit{\tau }}_{\mathrm{C}}-{\mathit{\tau }}_{\mathrm{R}}}}\\ & \cdot \left({\mathit{\tau }}_{\mathrm{R}}\cdot e^{-\frac{t_{i+\mathrm{1}}-t_i}{{\mathit{\tau }}_{\mathrm{R}}}}-{\mathit{\tau }}_{\mathrm{C}}\cdot e^{-\frac{t_{i+\mathrm{1}}-t_i}{{\mathit{\tau }}_{\mathrm{C}}}}+\left({\mathit{\tau }}_{\mathrm{C}}-{\mathit{\tau }}_{\mathrm{R}}\right)\right).\end{array}\end{array}$$

The observables, which allow a comparison to measured data, are as follows.

• Average river runoff:

  $$\begin{array}{}\text{(19)}& {\stackrel{\mathrm{‾}}{R}}_i^{\mathrm{R}}={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{R}}}}\cdot {\stackrel{\mathrm{‾}}{M}}_i^{\mathrm{R}},\end{array}$$

  corresponding to measured monthly runoff.

  – Average total mass:

  $$\begin{array}{}\text{(20)}& {\stackrel{\mathrm{‾}}{M}}_i^{\mathrm{T}}={\stackrel{\mathrm{‾}}{M}}_i^{\mathrm{C}}+{\stackrel{\mathrm{‾}}{M}}_i^{\mathrm{R}},\end{array}$$

  corresponding to monthly GRACE data.

The equations Eqs. (10)–(20) are self-consistent, i.e. the corresponding balance equations are fulfilled with the following:

$$\begin{array}{}\text{(21)}& {\displaystyle \frac{M_{i+\mathrm{1}}^{\mathrm{T}}(t-t_i)-M_i^{\mathrm{T}}}{(t-t_i)}}+{\stackrel{\mathrm{‾}}{R}}_{i+\mathrm{1}}^{\mathrm{R}}(t-t_i)=N_{i+\mathrm{1}}.\end{array}$$

For the single storage approach the above piecewise analytical solutions of the cascaded storage approach (Eqs. 8–21) are used for τ_R≪τ_C (here ${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$ months). For this case the river network mass is negligible compared to the catchment mass.

3 Properties and optimization performance

For the evaluation of the parameter optimization performance of the cascaded storage approach an example with synthetic recharge as input is investigated. This permits the quantification of the uniqueness and accuracy of the parameter estimation undisturbed by noise. It also facilitates the discrimination of errors in the calculation scheme itself and impacts arising from undescribed processes when compared to real-world data. For an application to GRACE measurements the main question is if and why the time constants τ_C and τ_R can be determined independently by an optimization versus anomalies in total mass and/or river runoff. Thus, in order to understand the optimization results with respect to uniqueness the general properties of the approach are presented and discussed first. For the synthetic case a recharge time series of sinusoidal form with a period of 12 arbitrary time units and length units with an amplitude and mean value of 1 is used as the driving force, and the calculation is run until equilibrium is reached. The example in Fig. 2 shows the effect of a non negligible river network time constant τ_R=2.5 time units for a catchment time constant τ_C=3 time units, which leads to an increase in total mass $M^{\mathrm{T}}\left(t\right)=M^{\mathrm{C}}\left(t\right)+M^{\mathrm{R}}\left(t\right)$ with respect to the average level and signal amplitude and to a phase shift between total mass M^T and river mass M^R, i.e. the corresponding river runoff R^R.

Figure 2 Time series of recharge N, catchment mass M^C, river network mass M^R and total masses M^T for the synthetic case at equilibrium.

In order to describe the general behaviour of the mass and runoff time series and their dependence on τ_C and τ_R, their properties are summarized here in the form of statistical values for the synthetic case with the sinusoidal recharge in equilibrium. This helps our understanding of why unique values for the time constants are achieved in the parameter optimization process. The values of time constants τ_C and τ_R used for the statistical description cover a wide range from 0.1 to 100 time units and are combined independently.

3.1 Catchment and river mass

Based on the mean mass values, Eqs. (14), (16), (18), of each stress period the long-term averages for the storage compartments are given by the following:

$$\begin{array}{}\text{(22a)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}^{\mathrm{C}}=\stackrel{\mathrm{‾}}{N}\cdot {\mathit{\tau }}_{\mathrm{C}},\text{(22b)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}^{\mathrm{R}}=\stackrel{\mathrm{‾}}{N}\cdot {\mathit{\tau }}_{\mathrm{R}},\text{(22c)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}^{\mathrm{T}}=\stackrel{\mathrm{‾}}{N}\cdot ({\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}).\end{array}$$

For τ_R≪τ_C (here ${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) the river network mass is negligible and the solution corresponds to a single storage approach. For a non negligible river network storage the given average values for total mass M^T mean that the effective “total” time constant is given by the sum of the catchment and river time constants ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$, which means that the total mass M^T observed by GRACE is bigger than the mass M^C calculated for the catchments alone. However, Eq. (22c) cannot be used for the determination of ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ from GRACE measurements directly as GRACE only provides mass anomalies.

The relative signal amplitudes (standard deviations normalized with those of the respective input) of both the catchment mass M^C or river mass M^R show the same functional form ${\mathit{\sigma }}_{\mathrm{MC}}/{\mathit{\sigma }}_{\mathrm{N}}\sim {\mathit{\sigma }}_{\mathrm{MR}}/{\mathit{\sigma }}_{\mathrm{RC}}=$ SD(M^C) ∕ N for the respective time constants τ_C or τ_R (Fig. 3, ${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) with a monotonous increase to an asymptotic value ${\mathit{\sigma }}_{\mathrm{MC}}/{\mathit{\sigma }}_{\mathrm{N}}\sim {\mathit{\sigma }}_{\mathrm{MR}}/{\mathit{\sigma }}_{\mathrm{RC}}=\mathrm{2}$, which is reached at about one full period of the input. The superposition of the signal amplitudes for the observable total mass $M^{\mathrm{T}}\left(t\right)=M_{\mathrm{C}}^{\mathrm{C}}\left(t\right)+M^{\mathrm{R}}\left(t\right)$ leads to a complex behaviour for σ_MT∕σ_N(τ_C, τ_R) (Fig. 3), if the river time constant τ^R is not negligible (${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) and especially if it gets close to τ_C.

Figure 3 Signal amplitudes of total mass normalized by recharge: σ_MT∕σ_N versus total mass time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ for different river time constants τ_R.

3.2 Catchment and river runoff

The calculated long-term averages of the runoff contributions R^C and R^R correspond to the ones of the water balance equations, Eqs. (4), (6), given by the mean recharge and thus are not dependent on the time constants.

$$\begin{array}{}\text{(23)}& {\stackrel{\mathrm{‾}}{R}}^{\mathrm{R}}\left(t\right)={\stackrel{\mathrm{‾}}{R}}^{\mathrm{C}}\left(t\right)=\stackrel{\mathrm{‾}}{N}\end{array}$$

Thus, an observed long-term average of runoff does not permit the determination of the time constant and hence the storage volume (Eq. 22).

The relative signal amplitudes of both catchment and river runoff (normalized with the respective input σ_RC∕σ_N and σ_RR∕σ_RC) show the same functional form corresponding to a single storage approach (Fig. 4, ${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) and decrease monotonously with the respective time constants τ_C and τ_R to an asymptotic zero. However, the signal amplitude of the observable river runoff ${\mathit{\sigma }}_{\mathrm{RR}}/{\mathit{\sigma }}_{\mathrm{N}}({\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$), normalized with recharge N, shows a deviations for different τ_C and τ_R with the same τ_C+τ_R (Fig. 4).

Figure 4 Signal amplitudes in standard deviations for river runoff normalized by recharge: σ_RR∕σ_N versus total mass time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ for combinations in (τ_C, τ_R).

Both observables, total mass and river runoff, show a non unique behaviour with respect to combinations in (τ_C, τ_R) for the same ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ and considerable deviations from the single storage approach (${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$). Measurements of the signal amplitudes thus only provide coarse estimates of the total time constant τ_T, yet do not permit distinction between τ_R and τ_C and between catchment and river network storage.

However, so far, only the signal amplitudes are examined, but not the specific properties of the time series, i.e. the dynamic response to input signals in form and phase. The convolution in the solution of the balance equation, Eqs. (8) and (11), leads to a different phasing with respect to the input N(t), which can be utilized for a separation of the respective time constants.

3.3 Phasing

For the synthetic example with a sinusoidal recharge time series N(t) as input the phasing ω of the different response signals is determined by the fit of a sinusoidal function (Fig. 5). This facilitates the easy determination of the phasing and thus the relative phase shift Δω between the signals. Masses and the related runoffs are in phase for the same storage compartments (Eq. 15). For a negligible river network time constant (${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) river runoff R^R is in phase with the catchment storage M^C.

Figure 5 Phasing of river network mass with respect to recharge time series displayed versus τ_C for different τ_R.

The functional form of the phasing ${\mathit{\omega }}_M^{\mathrm{C}}$ for the catchment mass M^C or the corresponding runoff R^C relative to recharge N(t) (Fig. 5) can be empirically described by the monotonous function:

$$\begin{array}{}\text{(24)}& {\mathit{\omega }}_{\mathrm{MC}}\left({\mathit{\tau }}_{\mathrm{C}}\right)={\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}}{\mathit{\lambda }}}\right),\end{array}$$

with the empirical parameters ${\mathit{\omega }}_{max}=\mathrm{2.8}$ and λ=2.7 and an error ε < ∼2 % relative to the maximum.

As the catchment runoff R^C with the phasing ${\mathit{\omega }}_{\mathrm{M}}^{\mathrm{C}}$ serves as input into the river system, the phasing of the river system with respect to catchment runoff R^C, which has the same functional form as Eq. (24), is added on top of it (Fig. 5). The resulting phasing of the river network storage or river runoff is thus given by a superposition in the following form:

$$\begin{array}{}\text{(25)}& {\mathit{\omega }}_{\mathrm{RR}}({\mathit{\tau }}_{\mathrm{C}},{\mathit{\tau }}_{\mathrm{R}})={\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}}{\mathit{\lambda }}}\right)+{\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}}\right)\end{array}$$

for any combination (τ_C, τ_R) and with the same empirical parameters as in Eq. (24).

As total mass $M^{\mathrm{T}}\left(t\right)=M^{\mathrm{C}}\left(t\right)+M^{\mathrm{R}}\left(t\right)$ is the superposition of the signals with the respective amplitudes and phasing, the phasing of total mass M^T(t) is situated between catchment and the river system mass according to τ_R. This means that for non negligible river network mass (τ_R > 0) a phase shift between total mass (GRACE) and observed river discharge and also between total mass and modelled catchment mass must occur. The phasing of total mass M^T(t) for all combinations (τ_C, τ_R) (Fig. 6) shows the same functional form as ω_MC and ω_MR (Eqs. 24, 25) if displayed versus the total time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$.

Figure 6 Phasing of total mass versus total time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$.

It can be approximated by the fitting function M^T fit:

$$\begin{array}{}\text{(26)}& {\mathit{\omega }}_{\mathrm{MT}}({\mathit{\tau }}_{\mathrm{C}},{\mathit{\tau }}_{\mathrm{R}})={\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}}\right),\end{array}$$

with the empirical parameters ω_max=2.95 and λ_T=3.2.

The phase shift between GRACE total mass and river runoff is thus given by the following:

$$\begin{array}{}\text{(27)}& \begin{array}{rl}\mathrm{\Delta }\mathit{\omega }({\mathit{\tau }}_{\mathrm{C}},{\mathit{\tau }}_{\mathrm{R}})& ={\mathit{\omega }}_{\mathrm{RR}}-{\mathit{\omega }}_{\mathrm{MT}}={\mathit{\omega }}_{max}(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}}{\mathit{\lambda }}})\\ & +{\mathit{\omega }}_{max}(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}})-{\mathit{\omega }}_{max}(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}}).\end{array}\end{array}$$

The empirical phase shift Δω from Eq. (27) corresponds to the one determined by a phase adaption Δω_adapt (Eqs. 38 and 39) of total mass and runoff within <∼5 % (see Supplement). This in principle allows for a determination of τ_C and τ_R separately from the adapted phase shift Δω_adapt and the total mass time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ according to Eq. (27). However, errors introduced by the linear interpolation used for the adaption of the phase shift lead to a much lower accuracy than the parameter estimation via the time series.

3.4 Parameter estimation

The analytical solutions for synthetic recharge time series permit the evaluation of the uniqueness and accuracy of the parameter optimization for given observables independent from limitations in the accuracy of numerical schemes and independent from noise in real-world data sets. For given combinations (τ_C, τ_R) the analytical solutions are used as synthetic measurements and are fitted with the same algorithm in order to retrieve the fit parameters (τ_C, τ_R).

As the total mass M^T (Eq. 20), and the phasing, (Eqs. 25–27), are commutative in (τ_C, τ_R), either of the data ranges τ_R < τ_C or τ_R > τ_C has to be used for a unique optimization. This is realized via an additional constraint in the optimization. For the discussion here the condition τ_R < τ_C is used, which hydrologically reflects the more frequent situations that the inundation volume is smaller than the catchment storage, but the results can also be applied to τ_R > τ_C, which might be the case in flat areas with a dense river network (such as the Amazon), which typically leads to temporarily inundated areas.

As absolute signal values are not relevant for the determination of the time constant from runoff or not available for GRACE data, the optimization versus the respective time series is based on signal amplitudes and the phasing. Thus, for a unique determination of (τ_C, τ_R) the following conditions have to be fulfilled.

• a.

  Optimization versus runoff:

  $$\begin{array}{}\text{(28)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{RR}}/{\mathit{\sigma }}_{\mathrm{N}}({\widehat{\mathit{\tau }}}_{\mathrm{C}},{\widehat{\mathit{\tau }}}_{\mathrm{R}})={\mathit{\sigma }}_{\mathrm{RR}}/{\mathit{\sigma }}_{\mathrm{N}}({\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}})\text{(29)}& {\displaystyle }{\mathit{\omega }}_{\mathrm{RR}}({\widehat{\mathit{\tau }}}_{\mathrm{C}},{\widehat{\mathit{\tau }}}_{\mathrm{R}})={\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}}{\mathit{\lambda }}}\right)+{\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}}\right).\end{array}$$
• b.

  Optimization versus mass anomalies:

  $$\begin{array}{}\text{(30)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{MT}}/{\mathit{\sigma }}_{\mathrm{N}}({\widehat{\mathit{\tau }}}_{\mathrm{C}},{\widehat{\mathit{\tau }}}_{\mathrm{R}})={\mathit{\sigma }}_{\mathrm{MT}}/{\mathit{\sigma }}_{\mathrm{N}}({\mathit{\tau }}_{\mathrm{C}},{\mathit{\tau }}_{\mathrm{R}})\text{(31)}& {\displaystyle }{\widehat{\mathit{\omega }}}_{\mathrm{MT}}({\widehat{\mathit{\tau }}}_{\mathrm{C}},{\widehat{\mathit{\tau }}}_{\mathrm{R}})={\mathit{\omega }}_{\mathrm{MT}}({\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}})={\mathit{\omega }}_{max}\left(\mathrm{1}-e^{-\frac{{\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}}{\mathit{\lambda }}}\right).\end{array}$$

With the constraints τ_R < τ_C or τ_R > τ_C there is only one (τ_C, τ_R) fulfilling the respective conditions, thus leading to unique solutions. The optimization delivers RMSE errors for the time series in the range ${\mathrm{10}}^{-\mathrm{8}}-{\mathrm{10}}^{-\mathrm{7}}$ and estimated time constants (τ_C, τ_R) with a relative error ε(τ_X)∕τ_X which does not depend on absolute values of (τ_C, τ_R) but on their ratio τ_R∕τ_C (Fig. 7).

Figure 7 Relative error ε(τ_X)∕τ_X of the time constants τ_C and τ_R for the cascaded storage approach with respect to optimizations versus total mass M^T or versus river runoff R^R.

For the synthetic case relative errors ε(τ_X)∕τ_X are very small ($\sim {\mathrm{10}}^{-\mathrm{7}}$ at ${\mathit{\tau }}_{\mathrm{R}}/{\mathit{\tau }}_{\mathrm{C}}\sim \mathrm{0}$) and show an exponential increase to a maximum of ∼1 % at τ_R∼τ_C. The error for τ_R < τ_C is analogous to τ_R > τ_C and equal for an optimization versus runoff or mass anomalies.

For catchments showing a phase shift between total mass and runoff the description of the system by a single storage approach (${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$) leads to a considerably higher relative error ε(τ_X)∕τ_X in the estimated time constant ${\mathit{\tau }}_{\mathrm{C}}\sim ({\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}})$ and thus also in drainable storage volume. It follows a power function and corresponds to ε < 10 % for τ_C < 3 and ε > 40 % for τ_C > 6. For this case the optimization versus river runoff or mass anomalies leads to different total time constants (relative difference ε > 7 % for τ_C > 5). Even though this might look like an acceptable result for τ_C < 3, there are still inevitable deviations in signal amplitudes (10 %–20 %) and phasing between the modelled and measured signals for both total mass and river runoff time series.

It can be summarized that in contrast to the single storage approach the cascaded storage approach permits the determination of both time constants (τ_C, τ_R) independently in a unique, highly accurate way for optimizations with respect to either total mass anomalies or river runoff. However, it has to be mentioned that even though the theoretical error in time constants remains below 1 % for τ_R∼τ_C, the ambiguity for τ_R < τ_C or τ_R > τ_C cannot be solved without further information on the volume of the river network.

3.5 Fully data-driven determination of drainable storage volumes

For the case that river discharge is available for a sufficient period of time the cascaded storage approach facilitates a simple determination of the drainable water storage volumes both for the catchment and for the river network directly from observations without the necessity of new model runs. The two time constants (τ_C, τ_R) adapted during a training period permit the quantification of the drainable water storage volumes M^T, M^C and M^R at other times directly from observations of GRACE mass anomalies and river discharge. With a simple numerical adaption of the phase shift Δω resulting from the time constants (τ_C, τ_R) according to Eq. (27) a quite accurate determination of the total drainable storage volume from measured river discharge exclusively (or the other way round, of runoff from GRACE mass anomalies) is possible.

These can be determined by the following calculations.

1. Long-term averages of drainable storage volumes from observed runoff R_o.

   $$\begin{array}{}\text{(32)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}_{\mathrm{sim}}^{\mathrm{C}}={\mathit{\tau }}_{\mathrm{C}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}},\text{(33)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}_{\mathrm{sim}}^{\mathrm{R}}={\mathit{\tau }}_{\mathrm{R}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}},\text{(34)}& {\displaystyle }{\stackrel{\mathrm{‾}}{M}}_{\mathrm{sim}}^{\mathrm{T}}={\mathit{\tau }}_{\mathrm{T}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}}\mathrm{with}{\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}}=\stackrel{\mathrm{‾}}{N}\mathrm{and}{\mathit{\tau }}^{\mathrm{T}}={\mathit{\tau }}^{\mathrm{C}}+{\mathit{\tau }}^{\mathrm{R}}\end{array}$$

   according to Eqs. (22a, 22b, 22c) and (23).

2. Time series of drainable storage volumes from GRACE and observed runoff R_o without the need for a phase adaption.

   $$\begin{array}{}\text{(35)}& M_{\mathrm{sim}}^{\mathrm{R}}\left(t\right)={\mathit{\tau }}_{\mathrm{R}}\cdot R_{\mathrm{o}}\left(t\right),\end{array}$$

   (NS_S 0.961, NS_R 0.576, corr_R 0.859 vs. M^R from Eq. 18),

   $$\begin{array}{}\text{(36)}& \begin{array}{rl}{M}_{\mathrm{sim}}^{\mathrm{T}}\left(t\right)& =\mathrm{d}{M}^{\mathrm{T}}\left(t\right)+{\stackrel{\mathrm{‾}}{M}}^{\mathrm{T}}=\mathrm{GRACE}\left(t\right)+{\stackrel{\mathrm{‾}}{M}}^{\mathrm{T}}\\ & =\mathrm{GRACE}\left(t\right)+{\mathit{\tau }}_{\mathrm{T}}\cdot \stackrel{\mathrm{‾}}{N}\cdot \\ & =\mathrm{GRACE}\left(t\right)+{\mathit{\tau }}_{\mathrm{T}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{0}},\end{array}\end{array}$$

   (NS_S 0.973, NS_R 0.751, corr_R 0.901 vs. M^T from Eq. 20),

   $$\begin{array}{}\text{(37)}& \begin{array}{rl}{M}_{\mathrm{sim}}^{\mathrm{C}}\left(t\right)& =M^{\mathrm{T}}\left(t\right)-M^{\mathrm{R}}\left(t\right)\\ & =\mathrm{GRACE}\left(t\right)+{\mathit{\tau }}_{\mathrm{T}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}}-{\mathit{\tau }}_{\mathrm{R}}\cdot R_{\mathrm{o}}\left(t\right),\end{array}\end{array}$$

   (NS_S 0.906, NS_R −0.065, corr_R 0.607 vs. M^C from Eq. 16).

   The simplified calculations directly based on observations lead to accurate equivalences to the fully calculated time series of total storage and the river network storage volumes M^T and M^R and to a reasonable description of the catchment storage volumes M^C.

3. Time series of total drainable storage volumes M^T, directly from observed runoff R_o or simulated river runoff $R_{\mathrm{sim}}^{\mathrm{R}}$ from GRACE with a numerical phase adaption of Δω.

   Use of the phase shift Δω_adapt adapted between GRACE and observed river runoff by a linear temporal interpolation (Riegger and Tourian, 2014) permits a simple description of river runoff directly from GRACE (Eq. 38) or of total drainable water storage M^T directly from observed runoff (Eq. 39) and corresponds to Δω from Eq. (27) within < ∼10 %. Both lead to very similar fitting performances.

   $$\begin{array}{}\text{(38)}& \begin{array}{rl}{R}_{\mathrm{sim}}^{\mathrm{R}}\left(t_i\right)& ={\displaystyle \frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{T}}}}\cdot \left[(\mathrm{1}-\mathrm{\Delta }\mathit{\omega })\cdot \mathrm{GRACE}\left(t_i\right)\right.\\ & +\mathrm{\Delta }\mathit{\omega }\cdot \mathrm{GRACE}\left(t_{i-\mathrm{1}}\right)+{\mathit{\tau }}_{\mathrm{T}}\cdot {\stackrel{\mathrm{‾}}{R}}_{\mathrm{o}}\left.\right]\end{array}\end{array}$$

   (NS_S 0.943, NS_R 0.698, corr_R 0.864 vs. measured R_o,)

   $$\begin{array}{}\text{(39)}& M_{\mathrm{sim}}^{\mathrm{T}}\left(t_i\right)={\mathit{\tau }}_{\mathrm{T}}\cdot \left[(\mathrm{1}-\mathrm{\Delta }\mathit{\omega })\cdot R_{\mathrm{o}}\left(t_i\right)+\mathrm{\Delta }\mathit{\omega }\cdot R_{\mathrm{o}}\left(t_{i+\mathrm{1}}\right)\right]\end{array}$$

   (NS_S 0.946, NS_R 0.483, corr_R 0.859 vs. GRACE)

For the representativeness of the fitting performance the fully data-driven approach (Eqs. 35–39) is compared to the respective masses and runoff from the cascaded storage approach applied to the Amazon basin (see below) and not to synthetic data. The related calculations are accessible in the Microsoft Excel workbook provided in the Supplement.

This performance means that the determination of the two time constants (τ_C, τ_R) by the cascaded storage approach during a sufficient training period facilitates a simple quantification of drainable storage volume or runoff time series directly from measured river discharge or GRACE anomalies. This provides a possibility to close data gaps in river discharge or GRACE directly from measurements with high accuracy.

4 Application to the Amazon basin

The R–S diagram of the full Amazon basin shows a hysteresis (Fig. 1a) corresponding to a phase shift, which can be interpreted as the time lag of river discharge. The Amazon basin upstream from Obidos is situated in a fully humid tropic environment with permanent, yet variable recharge and is large enough (4 704 394 km²) for low noise levels in the signals of GRACE and moisture flux divergence. With permanent recharge, flow contributions from overland flow and groundwater cannot be distinguished in the discharge curve. Also, on a spatial average over the full Amazon basin, with permanent recharge the uncoupled storages (like soil water storage, open water bodies etc.) are not time variant, i.e. there is no dry-out effect. Any contribution from time-dependent, uncoupled storages could be recognized in the R–S diagram as it would appear as a hysteresis, which does not correspond to a time lag, or by the respective deviations in the scatter plots of calculated versus measured runoff or storage volumes (see the Supplement). This is not the case.

Generally recharge from different approaches and products can serve as input to the system, such as the following:

1. $$\begin{array}{}\text{(40)}& N\left(t\right)=P\left(t\right)-\mathrm{ET}\left(t\right),\end{array}$$

   from the hydrometeorological products precipitation P and actual evapotranspiration ET_a.

2. $$\begin{array}{}\text{(41)}& N\left(t\right)=-\mathrm{\nabla }\cdot \mathit{Q},\end{array}$$

   from atmospheric data, with monthly vertically integrated moisture flux divergence viMFD.

3. $$\begin{array}{}\text{(42)}& N\left(t\right)={\displaystyle \frac{\partial }{\partial t}}M\left(t\right)+R\left(t\right),\end{array}$$

   from the terrestrial water balance with monthly temporal derivatives of GRACE measurements and measured river runoff R_o of the basin.

Here recharge (mm month⁻¹) is taken either from the water balance (Eq. 42), or from moisture flux divergence (Eq. 41), provided by ERA-INTERIM of ECMWF and processed by the Institute of Meteorology and Climate Research, Garmisch, Germany. For GRACE mass anomalies, data from GeoForschungsZentrum GFZ Potsdam Release 5 are used in millimetre equivalent water height. Both are handled as described in detail in Riegger and Tourian (2014). Their spatial resolution limits the application of the approach to global scales ≫200 000 km². River discharge is taken from the ORE HYBAM project (http://ore-hybam.org, last access: April 2016) and converted to runoff (mm month⁻¹) by normalization with the basin area. For a comparison of the calculated river network storage with observations from the Global Inundation Extent from Multi-Satellite (GIEMS; Prigent et al., 2001), flood area (km²) is used. As GRACE mass anomalies are most accurate for a monthly time resolution at present, the other data sets are aggregated to a monthly resolution as well. For the parameter optimization, time series of river runoff and GRACE mass anomalies are used for the time period from January 2004 until January 2009. Monthly runoff and the storage volume of the basin and river network are calculated for Amazon based on different recharge products here and optimized either versus runoff or GRACE mass anomalies. The results calculated with recharge from the terrestrial water balance optimized versus GRACE are shown in Figs. 8–10 for both (a) the monthly signal and (b) the monthly residual (monthly value minus mean monthly value) for January 2003–2009.

Figure 8 Time series of river runoff $R_{\mathrm{m}}^{\mathrm{R}}$ for the Amazon and optimization versus GRACE for the signal (a) for the residual (b).

Figure 9 Time series of total mass anomalies dM^T for the Amazon and optimization versus GRACE for the signal (a) for the residual (b).

Figure 10 Time series of river network storage $M_{\mathrm{m}}^{\mathrm{R}}$ and inundated area from GIEMS for the Amazon and optimization versus GRACE for the signal (a) for the residual (b).

The calculated river runoff R^R, total mass anomaly dM^T and river network mass M^R fit very well with the measured river runoff, GRACE and the flooded area from GIEMS both with respect to the signal and the de-seasonalized monthly residual.

The cascaded storage approach reproduces the phase shift between measured runoff R_o and total mass dM^T from GRACE. The calculated river network mass M^R of about 50 % of the total mass M^T for Amazon is proportional to observed runoff R_o without any phase shift (Fig. 11).

Figure 11 R–S relationships for observed runoff R_o versus the mass anomalies of GRACE, calculated total mass dM^T and river network mass dM^R.

Calculated hydraulic time constants, mean values and signal amplitudes for the absolute storage volumes are provided in Table 2 for the full Amazon basin upstream from Obidos. In addition the performance of optimizations either versus river runoff (column A) or versus GRACE (column B) and for different recharge products (column D, E) is displayed. This shows that the optimization versus different references leads to a very similar results while the fitting performance for the two recharge products (columns A, B and D, E) is quite different. For recharge from water balance (Eq. 42) the resulting time constants and thus the storage masses differ in a range of ∼5 % for the different references while they vary ∼10 % for recharge from moisture flux divergence.

In order to illustrate the benefits of the cascaded versus a single storage approach even in the fitting quality, results for a fixed ${\mathit{\tau }}_{\mathrm{R}}={\mathrm{10}}^{-\mathrm{3}}$, which correspond to a single storage, are shown (column C, F) for different recharge products. With the single storage approach – besides the much worse fitting performance – the resulting time constant ${\mathit{\tau }}_{\mathrm{T}}={\mathit{\tau }}_{\mathrm{C}}+{\mathit{\tau }}_{\mathrm{R}}$ is overestimated (corresponding to the investigations in Sect. 3) and the modelled signal amplitude is about 20 % less than that measured from GRACE. In addition a non negligible phase shift remains between the modelled runoff and measured discharge.

Table 1 The statistical characteristics are listed for calculated river runoff R^R, total mass M^T, basin mass M^C, river network mass M^R, observed river runoff R_o, GRACE mass anomalies and flood areas from GIEMS using the following: root mean square error (RMSE) of simulated versus measured values, Nash–Sutcliffe coefficient of the signal (NS_S) (simulated values versus long-term mean of measured values), Nash–Sutcliffe coefficient of monthly residuals (NS_R) (simulated values versus monthly mean of measured), correlation of simulated versus measured signals (corr_S), correlation of simulated versus measured monthly residuals (corr_R). Avg and SD are the long-term mean and standard deviations. The prefix “d” is used for anomalies related to the long-term mean. Results are compared for the different optimization references: runoff (R_o) or GRACE for recharge from water balance ($R+\mathrm{d}M/\mathrm{d}t$) (A, B) and for atmospheric input (−divQ) (D, E). The cascaded storage approach is compared to the single storage approach in (C, F).

The cascaded storage approach with recharge from the water balance (Eq. 42) leads to high-accuracy fits between calculated and measured river runoff and total storage mass for the signals (NS${}_{\mathrm{S}}{R}^{\mathrm{R}}-R_{\mathrm{o}}=\mathrm{0.96}$, NS_S dM^T− GRACE =0.98) and for the residuals (NS${}_R{R}^{\mathrm{R}}-R_{\mathrm{o}}=$ 0.74, NS_R dM^T− GRACE =0.74; the respective calculations are available in the Excel workbook provided in the Supplement).

The comparison of the water budgets for 14 different GHMs and LSMs (Getirana et al., 2014) for the Amazon basin permits the results of the cascaded storage approach to be sorted into those of the GHMs and LSMs (Fig. 14 of Getirana et al., 2014). With a Nash–Sutcliffe coefficient, NS_R (R= with respect to the mean seasonal cycle), of 0.74 and a correlation corr_R=0.90 compared to an NS_R of 0.58 and a correlation corr_R=0.84 for the best LSM, the cascaded storage approach outperforms the GHMs/LSMs for the Amazon basin upstream of Obidos. Even the fully data-driven approach (Eq. 39) leads to a performance comparable to the best GHM–LSM models tested by Getirana et al. (2014) with an NS_R of 0.483 and a correlation corr_R of 0.859 for simulated mass anomaly versus GRACE. The related calculations are accessible in the Excel workbook provided in the Supplement.

This is partly seen as the result of the simplicity of the lumped approach averaging out errors that emanate from the large number of different processes described by the GHMs and LSMs. However, the main reason for the better performance is seen in the quality of recharge data taken from the water balance using GRACE and river runoff, as the use of moisture flux divergence for this purpose leads to much worse performance.

The calculated river network mass M^R of the Amazon varies in the range of 40 %–65% of total mass M^T with an average of ∼50 %, corresponding to the values found by Paiva et al. (2013) and Papa et al. (2013) or ∼41 % by Getirana et al. (2017a). The correlation versus the observed flood area from GIEMS is higher for the calculated river network mass M^R (0.96 for the signal and 0.76 for the monthly residual) than for GRACE (0.92 and 0.65 respectively). The consistency of the calculated river network mass (and the corresponding observed river runoff $R^{\mathrm{R}}=M^{\mathrm{R}}{\mathit{\tau }}_{\mathrm{R}}^{-\mathrm{1}}=\mathrm{0.742}{M}^{\mathrm{R}}$) with the flood areas is seen much more clearly in the phasing (Fig. 12), which shows a clear phase shift for GRACE versus GIEMS (see also Papa et al., 2008), yet none for calculated M^R. As Getirana et al. (2017a) already emphasized, only an appropriate description of the river network storage permits a correct description of the total storage in amplitude and phasing for a comparison with GRACE.

Figure 12 GRACE mass, calculated river network mass dM^R and observed river runoff dR_o versus flood area dGIEMS, all displayed as anomalies (please consider d$M^{\mathrm{R}}=\mathrm{d}{R}^{\mathrm{R}}{\mathit{\tau }}^{\mathrm{R}}=\mathrm{1.53}$dR^R).

5 Discussion

Distributed hydrological models use a lot of detailed local information in order to address a large number of involved processes for each grid cell. In this way they provide a spatially distributed and a very detailed composition of the involved storages and flows. However, it is very difficult to discriminate the respective processes locally with the consequence that only their superposition can be compared to measurable data like river discharge. This creates a kind of ambiguity between the different contributions, thus losing some of the benefits of a detailed description. As has also been pointed out by Getirana et al. (2017a), for a comparison of the superposed storages to GRACE anomalies the river network storage changes have to be quantified as well, as only the total storage changes are measured by GRACE. This means that an appropriate description of the river network storage and the time lag is an inevitable prerequisite for an appropriate adaption of model parameters. A hydrodynamic modelling of the river network facilitates the quantification of its storage to be sure, yet it involves a real computational challenge.

The cascaded storage approach permits the quantification of the drainable storage volumes for both land masses and river network directly from GRACE and measured river discharge for gauged catchments. For ungauged catchments for which GRACE and recharge data can be used, it provides good estimates for the respective storages and also for the otherwise unknown river discharge. This is achieved by adapting only two parameters, the time constants. Neither detailed information on local vegetation, surface, unsaturated/saturated zone, etc., and related flow processes nor a hydrodynamic modelling with detailed hydraulic information of the river network on river roughness, cross section, gradient or backwater effects is needed.

At present, the cascaded storage approach is limited to climatic and physiographic conditions for which the hysteresis is completely explained by a time lag, i.e. that no impacts of uncoupled storage components are visible in the R–S diagram. For global coverage the cascaded storage approach has to be extended by an explicit integration of coupled and uncoupled storage compartments to account for other regional climatic and physiographic conditions. The uncoupled storage components then have to be quantified either by their absolute storage volume or by their relative contribution to total storage.

As Riegger and Tourian (2014) have shown for boreal catchments, this can be done by means of remote sensing and a conceptional description. Boreal catchments are temporarily dominated by snow, leading to a huge hysteresis due to a superposition of masses from fully coupled (liquid) and uncoupled (solid) storage compartments. Remote sensing of the catchment snow coverage by MODIS facilitates the separation of the coupled liquid storage (proportional to river runoff) on the uncovered areas and the uncoupled frozen part on the snow-covered areas. The coupled liquid storage determined in this way actually constitutes a LTI system; i.e. the hysteresis can be fully explained by a phase shift. This fulfils the prerequisites for the cascaded storage approach and thus permits an application to boreal catchments as well. As a consequence, the principle of the cascaded storage approach is not limited to fully humid climatic conditions. It permits an application to other climatic regions as well, provided that the coupled and uncoupled storage compartments can be separated.

The description of monsoonal regions for example, which play an important role in the global water budget, is a considerable challenge. For these regions with seasonally dry periods, high-precipitation events during the wet season lead to distinct runoff in parallel from overland and groundwater, with different time constants τ_Surface and τ_GW, and to time-dependent uncoupled storage compartments like soil or isolated open water bodies, which do not contribute to discharge. All these storages have to be addressed for an adequate description of the drainable storage volumes. The uncoupled storage compartments have to be quantified by remote sensing (soil moisture and open water body altimetry from satellites) and subtracted from the total catchment mass measured by GRACE, which of course is a major task.

6 Conclusions

The test of the cascaded storage approach with synthetic recharge data has shown that the parameter optimization versus either mass anomalies or runoff reproduces the time constants (τ_C, τ_R) for both the catchment and the river network in a unique way with high accuracy, yet with an ambiguity for τ_R < τ_C or τ_R > τ_C, and thus in the related storage volumes. This problem can only be solved by reasonable assumptions or, preferably, by additional information on the volume of river network or flood areas, which can be taken from ground-based observations or remote sensing like GIEMS flood areas and water levels from altimetry. Numerical tests have also shown that the description of a system (showing a phase shift) by a single storage approach can only address the total drainable storage and thus leads to phasing differences between the calculated and measured runoff or storage and to considerable errors in the time constant of the total system.

The application to the full Amazon basin shows that the system behaviour including the time lag can be described by a simple conceptual model with a catchment and a river network storage in sequence and an adjustment of only two parameters, the time constants. The storage amplitudes for the total drainable water storage and the time lag to runoff are described with high precision. Calculated river network volume and the observed flood area correspond to GIEMS observations and the newest model results (including river routing; Getirana et al., 2017a) and are in phase with river discharge. This independent quantification of the river network volume permits an investigation of the relationship between flood areas, flood volumes, river runoff and calculated river network with its additional information and might provide insights into river hydraulics, i.e. routing times and the mass–area–level relationships of flooded areas.

As the optimization performance is comparable for either reference, the observed river runoff or GRACE anomalies, a calculation with given recharge and an optimization versus measured GRACE data can be used to determine both the river discharge and the drainable storage volumes, even for ungauged basins. For these cases the availability of accurate recharge data determines the accuracy of runoff and storage calculations at present. However, for ungauged basins the use of moisture flux divergence still provides quite acceptable results based on remote sensing and atmospheric data exclusively.

For the case that river discharge is available for a sufficient period of time in order to adapt the two time constants (τ_C, τ_R) sufficiently, the cascaded storage approach facilitates a simple “fully data-driven” determination of the drainable water storage volumes M^T(t), M^C(t) and M^R(t) at other times directly from observations of GRACE or river discharge without the necessity of new model runs. This permits data gaps to be closed or even forecasts within the period of the time lag to be made.

As the spatial resolution of GRACE and the accuracy of moisture flux divergence is limiting the applicability of the cascaded storage approach to global-scale catchments at the moment, any improvement in the spatial or temporal resolution and accuracy of GRACE and hydrometeorological data products will tremendously increase the number of catchments which can be described by this approach in future.

Appendix A: Abbreviations in the mathematical descriptions

Abbreviation	Description	Units: general or for application
N	Recharge = (precipitation − actual evapotranspiration)	Volume area−1 time−1 (mm month−1)
MC	Storage mass catchment	Mass density in equivalent water height (mm)
τC	Time constant catchment	Time unit (month)
RC	Runoff catchment	Volume area−1 time−1 (mm month−1)
ωMC	Phasing catchment mass	Time unit (month)
MR	Storage mass river network	Mass density in equivalent water height (mm)
τR	Time constant river network	Time unit (month)
RR	Runoff river network	Volume area−1 time−1 (mm month−1)
ωRR	Phasing river runoff	Time unit (month)
MT	Storage mass total system	Mass density in equivalent water height (mm)
τT	Time constant total system	Time unit (month)
ωMT	Phasing mass total system	Time unit (month)
Ro	Observed river runoff	Volume area−1 time−1 (mm month−1)
GRACE	GRACE mass anomaly	Mass density in equivalent water height (mm)
GIEMS	Flood area	Area (km2)
Prefix “d”	Indicates signal anomalies from long-term mean (anomalies)
Suffix “m”	Indicates mean values on the intervals