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.


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 5 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 10 participating in the dynamic water cycle driven by precipitation P, actual evapotranspiration ETa 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 yet no 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 15 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 20 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 to address the hydraulically coupled storage compartments via their contributions to river discharge. These comprise groundwater, surface water, the 25 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.
For time periods with no recharge hydraulically coupled storage components with a linear runoff storage (R-S) relationship 30 lead to an exponential decrease in river discharge or streamflow depending on the related hydraulic time constant : The corresponding total Drainable Storage in terms of mass density MStorage for any given time t0 is then given by an infinite temporal integration over river discharge Q(t) from the corresponding catchment area A -called river Runoff R(t) = Q(t)/A here -starting at time t0 : 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 to separate the flow components 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 to determine the Drainable Storage of the respective storage compartment and the corresponding mean 10 Drainable Storage : (3) from mean runoff R or recharge N.
On global scales the absolute storage volume of the Drainable Storages can be determined from runoff time series directly, if 15 there are distinct and long enough periods of negligible or even negative recharge (actual evapotranspiration > 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 subcatchments. 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. 20 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 intend 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 25 discharge. The main difficulties in verifying large or Global Scale Hydrological Models or Land Surface Models (GHMs / LSMs) consist in the quantification of local individual storage volumes and related flows by local groundbased 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. 30 Since 2013 GRACE observations of the time-variable gravity field provide monthly distributions of mass density on large spatial scales >~ 200000km 2 (Tapley et al., 2004). However, as the water storage in different compartments (snow, ice, vegetation, soil, surface-, 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 5 comparison of measured TWS and observed river runoff Ro. Surprisingly some global hydrological models (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 et al., 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 10 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 time scales by Riegger and Tourian, 2014, revealed that measured runoff-storage (R-S) diagrams show hysteresis curves of distinct form and extent (Fig1a,b), which are characteristic for different climatic conditions (like fully humid, seasonally dry or boreal) and can be explained considering recharge and runoff properties (Fig1c,d) . 15 Thus for example, catchments in fully humid conditions (like the full Amazon basin upstream Obidos (295) and some of its sub-catchments like upstream Manacapuru (501)) 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 shiftriver 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 20 the slope in the R-S diagram corresponds to the hydraulic time constant via =R/M. 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 5 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). 5 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 shure, this leads to reasonable results for the sub catchments with permanent input (Fig.1 a, c) for which the time dependent uncoupled storage is negligible. However, for Rio Branco 10 (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.

15
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 20 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). 25 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 "topdown" or lumped approach is to integrate the catchment scale water balance and describe the system by large scale variables 30 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 are quite accurate, yet limited to global scales (see below). Thus, opposite to distributed hydrological models which are based on spatially / 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 5 (Nash Sutcliffe 0.97 for full Amazon) by an adaption of only two parameters, the hydraulic time scale and the phase shift, even though the physical cause of the phase shift is not addressed explicitly.
A disadvantage of the above approaches Tourian, 2014, Tourian et al., 2018) is that it does not permit to quantify the individual Drainable storage volumes on landmasses and in the river network separately, but only the total 10 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 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 15 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 "Cascaded storage approach" based 20 on the integration of given recharge in the water balance utilizing 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 to describe 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 25 further modelling.
The paper is structured as follows: Section 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 30 for a negligible river network time constant. In section 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 of further model runs. In section 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. Conclusions are drawn in section 5and discussed with respect to global applications of GSMs/LSMs. Possible future investigations in order to overcome limitations of the approach are sketched.

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 5 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 "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 10 groundwater flow from the land masses are summarized to a "Catchment" runoff R R .
• River runoff (river discharge / 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 MR 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. 15 • 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 20 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. 25 The following abbreviations are used in the mathematical description (  Table.1: Abbreviations in the mathematical descriptions: The total system behaviour is described by two balance equations : 5 with a proportional R-S relationship for hydraulically coupled storages. N denotes the recharge as input, RC C the catchment runoff from the catchment storage M C , which cannot be measured directly on large spatial scales, and RR R the river runoff from the river network storage MR R which can be measured at discharge gauging stations. 5 The water balance equation, Eq.(4), for the catchment is generally solved by: where MC C (t0) is the initial condition and N(w) 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 10 arbitrary length (stress periods) the recharge time series can be described as : For calculation convenience Eq. (8) The respective catchment runoff RC C based on Eq. (5) and M C MC 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 MR R : and the iterative solutions for time The total mass MT T is then given by : (13) 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 MT T for MR R > MC C .
It has to be emphasized here, that the piecewise analytical solutions for time periods of constant recharge provide a 5 mathematical solution for an arbitrary temporal resolution without numerical limitations. Finite Difference solutions are limited by stability criteria (ti+1-ti)< and accuracy criteria (ti+1-ti)< 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 10 The observables related to measurements by GRACE and discharge from gauging stations are the total mass anomaly dMT T and the river runoff RR R . GRACE observations with acceptable error are still limited to monthly resolution. Discharge as well as some of the meteorological inputs like precipitation, evapotranspiration or moisture flux divergence are often measured in daily values, some of the products in monthly values. For an optimal adaption to the monthly resolution of 15 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 20 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 MX is given for each interval   i.e. mean catchment mass and runoff : and mean river mass and runoff : The equations Eq. (10) For the Single Storage approach the above piecewise analytical solutions of the Cascaded Storage approach, Eq. (8) -Eq.
(21), are used for R<<C (here R = 10 -3 months). For this case the river network mass is negligible compared to the catchment mass 3 Properties and optimization performance 15 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 20 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 one 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 MT T (t) = MC C (t)+MR R (t) with respect to the average level and signal amplitude and to a phase shift between total mass MT T and river mass MR R i.e. the corresponding river runoff RR R . 5

Catchment and river mass
Based on the mean mass values, Eq.(14), (16), (18), of each stress period the long term averages for the storage compartments are given by : For R<<C (here R = 10 -3 ) the river network mass is negligible and the solution corresponds to a Single Storage approach. 20 For a non negligible river network storage the given average values for total mass MT T mean that the effective "total" time constant is given by the sum of the catchment and river time constants T = C + R, which means that the total mass MT T observed by GRACE is bigger than the mass M C MC C calculated for the catchments alone. However, Equation (22c) cannot be used for the determination of T = C + 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 MC C or river mass MR R show the same functional form MC / N ~ MR / RC = stdev(M C ) / N for the respective time 5 constants C or R (Fig.3, R = 10 -3 ) with a monotonous increase to an asymptotic value MC / N ~ MR / RC = 2 which is reached at about one full period of the input. The superposition of the signal amplitudes for the observable total mass MT T (t) = MC C (t) +MR R (t) leads to a complex behaviour for MT / N (C, R) (Fig.3), if the river time constant  R is not negligible (R = 10 -3 ) and especially if it gets close to C.

Catchment and river runoff
The calculated long term averages of the runoff contributions RC C and R R correspond to the ones of the water balance equations, Eq.(4), (6), given by the mean recharge and thus are not dependent on the time constants.  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, R = 10 -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 RR / N (C + R), normalized with recharge N, shows a deviations for different C and R with the same C+R (Fig.4). 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, Eq.(8) and (11), leads to 15 a different phasing with respect to the input N(t), which can be utilized for a separation of the respective time constants.

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 20 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 (R = 10 -3 ) river runoff RR R is in phase with the catchment storage MC C . The functional form of the phasing MC C for the catchment mass M C MC C or the corresponding runoff RC C relative to 5 recharge N(t) (Fig.5) can be empirically described by the monotonous function : with the empirical parameters max = 2,8 and = 2.7 and an error  < ~2% relative to the maximum.
As the catchment runoff RC C with the phasing MC C serves as input into the river system, the phasing of the river system with 10 respect to to catchment runoff RC 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 form:  (25) for any combination (C, R) and with the same empirical parameters as in Eq. (24).

15
As total mass MT T (t) = M C MC C (t) + MR R (t) is the superposition of the signals with the respective amplitudes and phasing, the phasing of total mass MT 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 MT T (t) for all combinations (C, R) Fig.6 shows the same functional form as MC and MR, Eq.(24), (25) if displayed versus the total time constant T = C + 20 R.
with the empirical parameters max = 2,95 and = 3.2 The phase shift between GRACE total mass and river runoff is thus given by : 10 The empirical phase shift  from Eq.27 corresponds to the one determined by a phase adaption adapt (Eq. 38 & 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 T = C + R according to Eq. (27). However, errors 15 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.

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). 5 As the total mass MT T , Eq. (20), and the phasing, Eq. (25)(26)(27), are commutative in (C, R), either the data range 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. 10 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 15 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 10 -8 -10 -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). 25 For the synthetic case relative errors (X) /X are very small (~10 -7 at R / C ~ 0) and show an exponential increase to a 5 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 (R = 10 -3 ) leads to a considerably higher relative error (X) /X in the estimated time constant C ~ (C + R) and 10 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 (rel. Diff.  > 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. 15
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. 20

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 to quantify the Drainable water storage volumes MT T , M C MC C and MR R at other times directly from observations of 5 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 : 10 ) (NSS 0.906, NSR -0.065, corrR 0.607 vs M C MC C from Eq.16) 25 The simplified calculations directly based on observations lead to accurate equivalences to the fully calculated time series of total and the river network storage volumes MT T and MR R and to a reasonable description of the catchment storage volumes MC C . 30 3. Time series of total Drainable storage volumes MT T , directly from observed runoff Ro or simulated river runoff R R sim 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 MT T directly from observed runoff (Eq.39) and corresponds to  from Eq.27 within 5 <~10%. Both lead to very similar fitting performances.  For the representativeness of the fitting performance the fully data driven approach (Eq.35-Eq.39) is compared to the respective masses and runoff from the Cascaded storage approach applied to the Amazon basin (see below) and not to 15 synthetic data. The related calculations are accessible in the xls-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 20 directly from measurements with high accuracy.

Application to the Amazon catchmentbasin
The R-S diagram of the full Amazon basin shows a hysteresis (Fig.1b, d) corresponding to a phase shift, which can be interpreted as the time lag of river discharge. The Amazon basin upstream Obidos is situated in a fully humid tropic 25 environment with permanent, yet variable recharge and is large enough (4704394km 2 ) for low noise levels in the signals of Generally recharge from different approaches and products can serve as input to the system such as : from the hydrometeorological products precipitation P and actual evapotranspiration ETa from atmospheric data, with monthly vertically integrated moisture flux divergence viMFD 3.
from the terrestrial water balance with monthly temporal derivatives of GRACE measurements and measured river 10 runoff Ro of the catchmentbasin.
Here recharge [mm/month] is taken either from the water balance, Eq. (34)  The calculated river runoff RR R , total mass anomaly dMT T and river network mass MR 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 Ro and total mass dMT T (or GRACE respectively). The calculated river network mass MR R of about 50% of the total mass MT T for Amazon is proportional to observed runoff Ro without any phase shift (Fig.11) Table.2 for the full Amazon basin upstream 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.(34), the resulting time constants 15 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 R = 10 -3 , which correspond to a Single Storage, are shown (Column C, F) for different recharge products. With the Single Storage approach -beside the much worse fitting performance -the resulting time constant T = C + R is 20 overestimated (corresponding to the investigations in section 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.   The Cascaded Storage approach with recharge from the water balance, Eq. (34), leads to high accuracy fits between calculated and measured river runoff and total storage mass for the signals (NSS RR R -Ro=0.96, NSS dMT T -GRACE = 0.98) and for the residuals (NSR RR R -Ro = 0.74, NSR dMT T -GRACE = 0.74, the respective calculations are available in the xlsworkbook provided in the supplement).
The comparison of the water budgets for 14 different GHMs / LSMs (Getirana et al. (2014)) for the Amazon basin permits to 5 sort in the results of the Cascaded Storage Approach into those of the GHMs/LSMs (Getirana et al. (2014) , Fig. 14). With a Nash-Sutcliffe coefficient NSR (R = with respect to the mean seasonal cycle) of 0.74 and a correlation corrR = 0.90 compared to an NSR of 0.58 and a correlation corrR=0.84 for the best LSM, the Cascaded Storage approach outperforms the GHMs/LSMs for full Amazon. Even the fully data driven approach (Eq.39) with an NSR of 0.483, and a correlation corrR=0.859 for simulated mass anomaly versus GRACE is comparable to the best performances in the GHM/LSM test by 10 Getirana et al. (2014). The related calculations are accessible in the xls-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/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. 15 The calculated river network mass MR R of the Amazon varies in the range of 40-65% of total mass MT T with an average ~50%, corresponding to the values found by Paiva et al., 2013and Papa et al., 2013or ~41% by Getirana et al., 2017a. The correlation versus the observed flood area from GIEMS is higher for the calculated river network mass MR 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 20 river network mass (and the corresponding observed river runoff RR R = MR R R − = 0.742 MR 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 MR R . As already Getirana et al., 2017a 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

Conclusions and Discussion
The test of the Cascaded storage approach with synthetic recharge data has shown that the parameter optimization either 5 versus 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 better 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 10 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 catchment 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 15 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 newest model results incl. 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 20 hydraulics i.e. routing times and the mass-area-and level-relationships of flooded areas.

Formatted: Not Highlight
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 as well as 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. 5 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 MT T (t), MC C (t) and MR R (t) at other times directly from observations of GRACE or river discharge without the necessity of new model runs. This permits to close data gaps or even do forecasts within the period of the time lag. 10

Discussion
Comparing global hydrological or land surface models with the presented Cascaded storage approach it can be summarized: 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 15 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 loosing parts of the benefits of a detailed description. As it 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 20 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 means a real computational challenge.
The Cascaded storage approach permits the quantification of the Drainable storage volumes both for land masses and river 25 network directly from GRACE and river discharge for gauged catchments or from GRACE and recharge for ungauged catchments by adapting only two parameters, the time constants. No detailed information on local vegetation, surface, unsaturated / saturated zone, etc. and related flow processes nor a hydrodynamic modelling with detailed hydraulic information on river roughness, cross section, gradient or backwater effects is needed.

30
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 Formatted: Heading 1 a 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.

5
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 catchments 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 10 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. In 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.

15
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 20 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.
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 / temporal resolution 25 and accuracy of GRACE and hydrometeorological data products will tremendously increase the number of catchments which can be described by this approach in future.

Data availability
In the supplement calculations and data are provided in an EXCEL workbook for the synthetic case and for the Amazon catchment.