Assessing Terrestrial Water Storage Variations in Southern Spain Using Rainfall Estimates and GRACE Data

: This paper investigates the relationship between rainfall, groundwater and Gravity Recovery and Climate Experiment (GRACE) data to generate regional-scale estimates of terrestrial water storage variations in the Andaluc í a region of southern Spain. These estimates can provide information on groundwater depletion (caused by periods of low rainfall or droughts) and ground-water recovery. The spatial distribution of groundwater bodies in southern Spain is complex and current in situ groundwater monitoring methods are deﬁcient, particularly in terms of obtaining representative samples and in implementing and maintaining groundwater monitoring networks. The alternative approach proposed here is to investigate the relationship between precipitation time series and changes in the terrestrial water storage estimated from GRACE observations. The results were validated against the estimated ﬂuctuation in regional groundwater. The maximum correlation between the mean groundwater level and the GRACE observations is 0.69 and this occurs at a lag of one month because the variation in gravity is immediate, but rainfall water requires around one month to travel across the vadose zone before it reaches the groundwater table. Using graphical methods of accumulated deviations from the mean, we show that, in general, groundwater storage follows the smooth, multi-year trends of terrestrial water storage but with less short-term trends; the same is true of rainfall, for which the local trends are more pronounced. There is hysteresis-like behaviour in the variations in terrestrial water storage and in the variations of groundwater. In practical terms, this study shows that, despite the abnormal dryness of the Iberian Peninsula during the 2004–2010 drought, the depleted groundwater storage in Andaluc í a recovered almost to its pre-drought level by 2016. In addition, groundwater storage and terrestrial water storage show very similar trends but with a delay in the groundwater trend.


Introduction
The Gravity Recovery and Climate Experiment (GRACE) and the continuing GRACE-Follow-on (GRACE-FO) are space-based missions [1] designed to measure changes in the Earth's gravity field (in the form of the geoid) that are directly related to variations in surface mass [2].International research centres provide estimates of the temporal variation of the Earth's gravity field derived from GRACE observations in the form of spherical harmonic coefficients (Groupe de Recherche en Geodesie Spatiale, Geo Forchungs Zentrum (Potsdam), Centre for Space Research at the University of Texas Austin) or global mascons (NASA/Jet Propulsion Laboratory).Several websites [3,4] provide variations in the gravity field interpreted as a change in geoid height, equivalent water thickness, and viscoelastic or elastic deformation.In continental areas, variations in the gravity field, after atmospheric corrections, are mainly related to redistributions in continental water storages.These spatial and temporal variations in the surface mass signal are largely the sum of the changes in groundwater, water in the vadose zone, soil moisture and surface water stored in dams, lakes, and the snowpack [5].Previous studies have shown that the GRACE observations can be used to monitor water mass redistributions at the global scale [6][7][8], the continental scale [9,10], the regional scale [11][12][13] and large aquifer scale [14,15].For large areas in the mid-latitudes, monthly variations in the snowpack and in water in dams and lakes are small relative to the total terrestrial water storage (TWS).In these cases, the TWS can be regarded as the variation in the water in the subsurface including soil water storage (SWS), groundwater storage (GWS) and water in the vadose zone that can be considered to be in transit from the SWS to the GWS in the groundwater recharge process.Despite the importance of SWS in agriculture and evapotranspiration, GWS is the primary contributor to groundwater resources.GWS can be considered as a deposit from which water can be extracted for irrigation, industrial purposes and potable purposes.Groundwater delivers a worldwide average of 42% of the water required for irrigation [16].The natural discharge of groundwater gives rise to springs, and it is the origin of rivers that support important aquatic ecosystems [17].The widespread decline in groundwater levels over recent years [18] is generally associated with higher extraction rates and longer periods of severe drought [19,20].Reliable forecasting of the availability of water resources requires local and regional information on aquifer properties together with the impact of long dry periods on the water resources of stressed aquifer systems.The over-extraction of groundwater can cause aquifer compaction with accompanying land subsidence and consequent adverse impacts on infrastructure [21].Despite the importance of groundwater resources, the general quality and extent of monitoring networks makes it very difficult to quantify changes in the level of groundwater at regional scales [22].There are significant difficulties in processing in situ data such as the extent to which samples are representative of each aquifer; the quality and extent of the maintenance of monitoring networks; the representativeness of piezometers in relation to the spatial variability within a single aquifer; and the cost of sampling in frequency and in space [23].The use of the GRACE data could circumvent these issues and provide a direct measurement of regional groundwater variation, i.e., depletions and recoveries [20,24], whilst rendering the process more efficient.Van Loon et al. [25] used GRACE observations to study the impact of drought on spatially variable groundwater patterns in Europe.Other studies have used GRACE observations to investigate the manner in which worldwide drought is changing due to climate change [26].
The purpose of the work presented here is to assess the possibility of using GRACE data as a surrogate for groundwater variability over time (depletion and recovery of groundwater resources) at a regional scale and in a region with high hydrogeological complexity.There are many scientific papers using GRACE data in relation to the hydrology of large river basins or large aquifers.The case study presented here is in an administrative area (the Andalucía autonomous region in southern Spain with an area of 87,000 km 2 ) with high hydrogeological complexity.A total of 62% of the surface can be considered as permeable and able to store groundwater (i.e., aquifers) and 38% of the surface can be considered as impervious and without aquifers.Furthermore, the 62% of the permeable surface includes 189 aquifers with complex spatial dispositions.In this context, there is a significant difficulty, if not an impossibility in practical terms, in using in situ measurements to monitor the regional groundwater variations for the entire Andalucía region.The purpose of this paper is to analyse the possibility of monitoring the groundwater evolution (depletion and recoveries) by using the gravity data provided by GRACE.

Methodology and Assumptions for the Hydrological Cycle, Rainfall Data and GRACE Observations
The ANU GRACE visualisation tool portal [3], based on the spherical harmonics field estimated by the Groupe de Recherche en Geodesie Spatiale (GRGS), was used to estimate the equivalent water height or terrestrial water storage (TWS) time series at specific locations.The TWS time series were estimated using the integrated water change equation [3,5]: where ∆TWS(θ, λ, t) represents the variation in terrestrial water storage at the geographical location of co-latitude (θ) and longitude (λ) at time (t); R is the radius of the Earth; l, m are, respectively, the degree and order of the spherical harmonic model; P lm (cos(θ)) is the normalised associated Legendre functions at co-latitude (θ); ρ e and ρ w are, respectively, the densities of the Earth and water; k l is the elastic load Love numbers of degree l; and ∆C lm and ∆S lm are the spherical harmonic coefficient anomalies from the GRGS GRACE gravity fields at time (t).This equation is a function of the spherical harmonic coefficients and the maximum harmonic degree (L) of the spherical harmonic expansion.This parameter is important in gravity field recovery and, depending on the geographical location, an increase or decrease in its value impacts not only the stochastic noise properties of the resulting time series but may also result in the inclusion of other gravity effects from mass concentrations that are outside the spatial region that the measure represents, a phenomenon known as spatial gravity leakage [27,28].A rough estimate of the spatial resolution (SR) is given by the relationship SR = π/L [29], which implies that the GRACE measurement at a given location refers to an average value over a square area of side-length 680 km for degree 30 and of 400 km for degree 50.Using geostatistical terminology [30], this area is termed the spatial support of the GRACE data.In the same manner, the time support of the GRACE data is a month because the data are monthly accumulations of gravity anomalies.The uncertainties of the computed TWS values were calculated using the method of propagating the variances by the following equation [3] that uses the methodology of propagation of errors [31]: where σ 2 TWS , σ 2 ∆C(t) , σ 2 ∆S(t) are the variances of TWS, ∆C(t) and ∆S(t) , respectively.Usually, for each ∆TWS(θ, λ, t) there is an associated uncertainty that is evaluated by σ TWS (t), the standard error or standard deviation from Equation (2), which in general is of the order of 1 to 2 cm [3].A recent study of the uncertainty evaluation of GRACE data can be found in Boergens et al. [32].The uncertainty of the data used in the case study presented in this paper is discussed below.
One advantage of the GRGS spherical harmonic models is that they are generated in a way that does not require any filtering (i.e., spatial smoothing) or decorrelation (i.e., destriping) steps [33][34][35] before using the data [3].Some of the GRACE products use various techniques to filter the raw gravity observations recorded by satellites.We compared the results obtained for the CNES-GRGS data with those of the Centre for Space Research (CSR) at the University of Texas Austin [36] and retrieved from the EGSIEM website [4] and the results were similar.
A TWS time series is a function of the change in the Earth's gravity field due to changes in mass load, in general, water.In other words, it is the variation of the sum of terrestrial water storage that can be decomposed into surface water storage (SWS), soil moisture storage (SMS) and groundwater storage (GWS) after rainfalls in the wet season (Figure 1A).However, in the dry season, the variation in TWS is mainly proportional to GWS (Figure 1B).
Water input to the system comes from precipitation (P in Figure 1A), which may be in the form of liquid water (rainfall) or solid-state water (snow water equivalent, SWE); the latter can be seasonally accumulated in snowpacks in mountainous areas.The SWE and water stored in dams (D), lakes, wetlands and other repositories are the main components of SWS.SMS is the water stored in soil and, after accounting for evapotranspiration (ET), if the maximum storage value is exceeded, the excess water percolates and travels across the vadose zone until it reaches the groundwater zone (i.e., the recharge process).Water from rainfall enters the soil up to an infiltration value and the excess runs over the landscape (runoff water, R), entering drainage channels and forming streams that appear after a rainfall event.In the wet season, there is a transit of water across the vadose zone (VZ) from the soil zone (SZ) to the saturated zone or proper groundwater zone (SGZ).In the dry season (Figure 1B), the main contribution to the loss of TWS is from the decreasing GWS.This suggests a hysteresis-like behaviour of the GRACE signal in which the GWS lags behind the TWS and there is an obvious larger variability in the humid periods than in the dry periods.
Following [37] and our description of the water storage, the hydrological cycle and water budget can be written as: where P is the total precipitation, R is runoff, ET is evapotranspiration, ∆SWS is the variation in SWS, ∆SMS is the variation in SMS and ∆GWS is the variation in GWS.After rainfall in the wet season, it is assumed that ∆SWS, ∆SMS and ∆GWS increase.In the dry season (Figure 1B) there is a decrease in ∆SWS, ∆SMS and ∆GWS, with the latter being the most significant.Then, the water flowing in rivers in the dry season comes from aquifer discharge from the GWS.Thus, and therefore There is also the gravity or GRACE data and, as the GRGS gravity field variation (∆GF) was corrected for changes in atmospheric and ocean variations [3], it is directly proportional to the variation in TWS: Thus, during the wet season: whereas during the dry season the TWS is mainly groundwater and thus: The same reasoning can be extended from the yearly dry and wet seasons to a continuous multi-annual period of dry years with below-average rainfall (drought conditions) and to a multi-annual period of wet years (humid conditions with above average rainfall).Thus, in the variation in the gravity field, the effect lags behind its cause with more variability in the humid periods and less variability in the drought periods.During drought conditions, there is a decrease in groundwater storage that may be recovered during humid periods.However, it is possible that the rainfall in the humid periods may not be sufficient to recover to pre-drought levels, leading to aquifer depletion.and travels across the vadose zone until it reaches the groundwater zone (i.e., the recharge process).Water from rainfall enters the soil up to an infiltration value and the excess runs over the landscape (runoff water, R), entering drainage channels and forming streams that appear after a rainfall event.In the wet season, there is a transit of water across the vadose zone (VZ) from the soil zone (SZ) to the saturated zone or proper groundwater zone (SGZ).In the dry season (Figure 1B), the main contribution to the loss of TWS is from the decreasing GWS.This suggests a hysteresis-like behaviour of the GRACE signal in which the GWS lags behind the TWS and there is an obvious larger variability in the humid periods than in the dry periods.Following [37] and our description of the water storage, the hydrological cycle and water budget can be written as: where P is the total precipitation, R is runoff, ET is evapotranspiration, ∆SWS is the variation in SWS, ∆SMS is the variation in SMS and ∆GWS is the variation in GWS.After rainfall in the wet season, it is assumed that ∆SWS, ∆SMS and ∆GWS increase.In the dry season (Figure 1B) there is a decrease in ∆SWS, ∆SMS and ∆GWS, with the latter being the most significant.Then, the water flowing in rivers in the dry season comes from aquifer discharge from the GWS.Thus, and therefore There is also the gravity or GRACE data and, as the GRGS gravity field variation (∆) was corrected for changes in atmospheric and ocean variations [3], it is directly proportional to the variation in TWS: Thus, during the wet season: A useful tool in this study was the cumulative departure from the mean statistic [38] used to enhance the correlation in the short-term trends of rainfall data and TWS.Previous work, such as [38], studied the validity of the cumulative rainfall departure (CRD) and concluded that the CRD is a useful indicator of short-term rainfall trends and is an important parameter for the remainder of this paper.Given a set of n rainfall measurements {R 1 , R 2 , . . . ,R n } with mean value: the CRD is defined as: To ensure a meaningful comparison of the two datasets in the following case study, the same concept of cumulative departures from the mean statistic (CDMS) was applied to the TWS measurements and to the groundwater heads in the validation stage: where TWS is the mean value of the TWS measurements TWS i .The CDMS can be applied to the time series of any other variable.

Case Study: The Andalucía Region of Southern Spain
The study area is in the Andalucía region of southern Spain as shown in Figure 2, which has an area of 87,000 km 2 .This region is known for its high rainfall variability and includes the low-rainfall areas in the east and high-rainfall areas in the southwest.In Andalucía, dry periods can adversely impact non-irrigated crops (e.g., olive trees), especially in years when abnormally low rainfalls induce severe droughts, for example, the exceptional drought in south-eastern Spain from October 2013 to July 2014 [39].Eight geographical locations, shown in Figure 2 as large purple dots, were defined on a rectangular grid of 180 km × 90 km covering this region.The geographical co-ordinates of the grid points are shown in Table 1.Because the GRACE measurement at each point represents the mean value of a very large area (for example a square of 400 km × 400 km for maximum degree L = 50 in the spherical harmonics), the measurements overlap and it is irrelevant whether the grid was arbitrarily chosen, as the same results would be obtained with similar alternative grids.Figure 2 also shows a grid (small purple dots) of 5 km × 5 km covering the region of Andalucía that was used to obtain the rainfall estimates as explained below.The monthly TWS time series, for maximum degrees of L = 30 and L = 50, estimated at location 7 from the GRACE observations are shown in Figure 3A. Figure 3 shows 152 values, of which 144 are original GRACE data and the 8 missing monthly values were obtained by linear interpolation from their two nearest neighbours.The reason for the missing values is the failure of the satellite instruments.The completion of the GRACE time series was not strictly necessary for the purposes of this paper, but it was performed for completeness and for calculating the CDMS.It should be noted that these time series represent the average monthly accumulated variations in the gravity field for square areas of side lengths 680 km and 400 km for maximum spherical harmonic degrees L = 30 and L = 50, respectively.In particular, for the period of the anomalies below the average (dry periods denoted by (a) in Figure 3A), the peaks of the valleys are more negative, whereas during the period of anomalies above the average (humid periods, denoted by (b) in Figure 3A), the peaks of the annual maximum are higher.Nevertheless, the 0.90 coefficient of correlation between the two time series is high.Figure 3B shows, for the same maximum degree L = 50, the TWS time series for locations 6 and 7 in Figure 2. The two time-series are very similar (coefficient of correlation of 0.91) because there is a significant overlap of the spatial support of the measurements of two 400 km × 400 km areas, the centres of which are separated by 180 km.This is the conventional 400 km scale [40], which is adequate for studies of the water storage of large structures such as river basins [41], large aquifers [20], and areas with many aquifers and aquicludes as shown below.As the maximum degree of the spherical harmonics (L) increases, the spatial resolution increases (i.e., the support or area decreases) and the variability increases.As can be seen in Figure 3C for location 7a, the uncertainty of the GRACE TWS estimates also increases as the maximum number L of spherical harmonics in Equation (1) increases.In Figure 3C, the uncertainty is evaluated as the square root of the error variance of Equation (2), that is, the standard error or standard deviation of the error.The mean of the errors are 10 mm, 25 mm and 150 mm for the maximum degree of the spherical harmonics of L = 30, L = 50 and L = 80, respectively.While the standard error is acceptable for L = 30 and for L = 50, the standard error for L = 80 is so large that it renders the GRACE data almost useless and, thus, is not used in the work presented here.Figure 3C also shows the variability of the standard errors around their mean values and sometimes there is a spike in the time series, that is, a particularly high uncertainty value for a given month.The magnitude of the standard error is not correlated with the magnitude of the GRACE TWS data.The origin of these spikes may be due to the cumulative uncertainty in different GRACE instruments in that specific months possibly generating computational errors.Figure 3D shows a new perspective of Figure 3B obtained by considering the uncertainty of the measurements.The latter figure shows, for the same maximum degree L = 50, the TWS time series for locations 6 and 7 (as in Figure 3B) plus the standard error of the TWS for the location.This figure shows that the TWS of location 6 is within the interval of the TWS ± 1 standard error of location 7.In addition, the uncertainty for locations 6 and 7 is almost the same if the same maximum degree L is considered for both locations.Rainfall estimates on a 5 km × 5 km grid for the Andalucía test site were obtained from the Agencia Estatal de Meteorología (Spanish Meteorological Office (www.ign.es)accessed on 30 March 2022).Temperature estimates on a 10 km × 10 km grid were retrieved from the Spain02 climate database [42].The mean rainfall values over the test area from 2002 to 2015, the period in which the original GRACE data were recorded, were obtained by averaging the rainfall estimates on the 5 km × 5 km grid.The mean surface temperatures for the same time period were obtained by averaging the mean temperature estimates from the 10 km × 10 km grid.Rainfall estimates on a 5 km × 5 km grid for the Andalucía test site were obtained from the Agencia Estatal de Meteorología (Spanish Meteorological Office (www.ign.es)accessed on 30 March 2022).Temperature estimates on a 10 km × 10 km grid were retrieved from the Spain02 climate database [42].The mean rainfall values over the test area from 2002 to 2015, the period in which the original GRACE data were recorded, were obtained by averaging the rainfall estimates on the 5 km × 5 km grid.The mean surface temperatures for the same time period were obtained by averaging the mean temperature estimates from the 10 km × 10 km grid.

Analysis and Discussion of the Results
The areal monthly rainfall estimates for a 200 km × 200 km area centred on location 7 (red square in Figure 1) and for the whole of Andalucía (small purple dots in Figure 1) are shown in Figure 4A.These estimates show that, for large areas, whilst the variability in the rainfall time series is similar, the magnitudes are different.In Figure 4B, these rainfall estimates are compared with the TWS time series from GRACE for maximum degree L = 50, from which it can be concluded that, in general, the GRACE signal for maximum degree L = 50 captures the bimodality of rainfall together with other small variations that can be attributed to noise or to the filtering technique applied to the raw gravity field [36].They may also be the result of spatial leakage, a major limitation on the quantitative interpretation of satellite gravity measurements at regional scales.The latter phenomenon has been well documented over the past decade [40,43].Thus, in assessing the variation in TWS, the GRACE signal for maximum degree L = 30 (Figure 3A) is used when smoother variation is required and the signal for maximum degree L = 50 is used when the bimodality of each annual cycle is required.In addition, it is evident that for several peaks in Figure 4B there is a delay of around one to two months between the TWS time series and the rainfall data.This delay has been noted in various studies, e.g., [20,44], and is caused by the GRACE data registering the gravity anomalies introduced by the mass of water supplied by rainfall stored on the surface or underground.In other words, the TWS time series is derived from rainfall and if the rainfall data were a perfect sinusoid, the TWS signal (a cosine obtained by a 90 • horizontal shift) would have a delay of 90 • , or a quarter of the three-month wavelength.Figure 5 shows the TWS time series estimated with GRACE data at maximum degree L = 50, together with the mean surface temperature.The two curves are in opposite phase.The correspondence between minimum temperature and maximum TWS, and the reverse, is a consequence of rainy seasons that occur in winter and spring and dry seasons that occur during the summer and early autumn.
A spectral analysis using the Lomb-Scargle periodogram [45] was used to check the cyclical signal recorded in the TWS and the rainfall data.The Lomb-Scargle periodogram is specially designed for the spectral analysis of time series with uneven spacing of the observations such as the TWS time series where there are missing data.This avoids the artefacts introduced in the estimated spectral functions when the time series are completed by interpolation to yield even time series.The estimated periodogram is shown in Figure 6A in which the main cyclicities are highlighted.The annual hydrological cycle is evident, with the most prominent peak for a one-year period.The wet season in the study area is bimodal [46], meaning that it appears to have an almost half-year cycle with peaks at 5.4 months due to the variability of the bimodality changing from year to year. Figure 6A also shows longer cycles of 3.2 and 11.9 years.Some studies have discussed the relationship between the 3.2-year cycle and the North Atlantic Oscillation phenomenon [47].The 11.9-year cycle may be related to the North Atlantic Oscillation cycle in the six to tenyear range [48] or with the 11-year sunspot cycle.The time series are too short to draw any conclusive conclusions.Comparable peaks can be seen in the GRACE data and the rainfall data.The relationship between the TWS time series and rainfall is clear from the cross-spectral analysis as shown in [45].The squared coherence plots in Figure 6B show high-value peaks (the maximum value of the squared coherence is 1) at frequencies that correspond to periods of 3.6 years, 1 year and 7.8 months.The statistical confidence levels of these coherence values shown in Figure 6C were evaluated using the permutation test [45].The estimated phases are shown in Figure 6D in which the phases of 3.6 years, 1 year and 7.8 months as periodicities have values of 55 • , 27 • and 24 • , respectively, which imply delays between TWS and rainfall of 6.6 months, 1 month and 0.5 months, respectively.While the 1 month or 1.5-month delays are visually evident in Figure 4B, the delay of almost 7 months implies a delay in the trends of the signals as discussed below.
The general and local trends for the CRD and the accumulated TWS time series for maximum degree L = 50 and location 7 are shown in Figure 7.The results show a possible short increasing trend between (a) and (b) followed by a significant negative slope between (b) and (c) and a significant positive slope between (c) and (e).The variations in the accumulated TWS time series are smoother than those of the rainfall data.The TWS time series has more inertia, i.e., less variability and a smoother response, than the large oscillations in rainfall, although both variables show the same trend.Rainfall, in particular, includes the runoff water mass and evapotranspiration that are not included in the TWS variations.Assuming the variations in the gravity field recorded in Andalucía are mostly related to GWS levels, the minimum CRD and the accumulated TWS in (c) around the end of 2010 must have resulted in minimum GWS in the area, which is verified below.Figure 7 clearly shows the significant drought in Andalucía that started in 2004 and ended in 2010 (Junta de Andalucía, Medio Ambiente, http://www.juntadeandalucia.es/medioambiente/(accessed on 30 March 2022)).It also shows that, by 2016, the rainfall since 2010 allowed the depleted GWS to recover to a level almost comparable to that of the pre-drought period.According to [49], the southern half of the Iberian Peninsula received less than 40% of its usual precipitation by June 2005, with the consequent depletion of large aquifers in the peninsula [50].
series is derived from rainfall and if the rainfall data were a perfect sinusoid, the TWS signal (a cosine obtained by a 90° horizontal shift) would have a delay of 90°, or a quarter of the three-month wavelength.Figure 5 shows the TWS time series estimated with GRACE data at maximum degree L = 50, together with the mean surface temperature.The two curves are in opposite phase.The correspondence between minimum temperature and maximum TWS, and the reverse, is a consequence of rainy seasons that occur in winter and spring and dry seasons that occur during the summer and early autumn.A spectral analysis using the Lomb-Scargle periodogram [45] was used to check the cyclical signal recorded in the TWS and the rainfall data.The Lomb-Scargle periodogram is specially designed for the spectral analysis of time series with uneven spacing of the observations such as the TWS time series where there are missing data.This avoids the artefacts introduced in the estimated spectral functions when the time series are completed by interpolation to yield even time series.The estimated periodogram is shown in Figure 6A in which the main cyclicities are highlighted.The annual hydrological cycle is evident, with the most prominent peak for a one-year period.The wet season in the study area is bimodal [46], meaning that it appears to have an almost half-year cycle with peaks at 5.4 months due to the variability of the bimodality changing from year to year. Figure 6A also shows longer cycles of 3.2 and 11.9 years.Some studies have discussed the relationship between the 3.2-year cycle and the North Atlantic Oscillation phenomenon [47].The 11.9-year cycle may be related to the North Atlantic Oscillation cycle in the six to tenyear range [48] or with the 11-year sunspot cycle.The time series are too short to draw any conclusive conclusions.Comparable peaks can be seen in the GRACE data and the rainfall data.The relationship between the TWS time series and rainfall is clear from the cross-spectral analysis as shown in [45].The squared coherence plots in Figure 6B show high-value peaks (the maximum value of the squared coherence is 1) at frequencies that correspond to periods of 3.6 years, 1 year and 7.8 months.The statistical confidence levels of these coherence values shown in Figure 6C were evaluated using the permutation test [45].The estimated phases are shown in Figure 6D in which the phases of 3.6 years, 1 year and 7.8 months as periodicities have values of 55°, 27° and 24°, respectively, which imply delays between TWS and rainfall of 6.6 months, 1 month and 0.5 months, respectively.While the 1 month or 1.5-month delays are visually evident in Figure 4B, the delay of almost 7 months implies a delay in the trends of the signals as discussed below.The general and local trends for the CRD and the accumulated TWS time series for maximum degree L = 50 and location 7 are shown in Figure 7.The results show a possible short increasing trend between (a) and (b) followed by a significant negative slope between (b) and (c) and a significant positive slope between (c) and (e).The variations in the accumulated TWS time series are smoother than those of the rainfall data.The TWS time series has more inertia, i.e., less variability and a smoother response, than the large oscillations in rainfall, although both variables show the same trend.Rainfall, in particular, includes the runoff water mass and evapotranspiration that are not included in the TWS variations.Assuming the variations in the gravity field recorded in Andalucía are mostly related to GWS levels, the minimum CRD and the accumulated TWS in (c) around the end of 2010 must have resulted in minimum GWS in the area, which is verified below.Figure 7 clearly shows the significant drought in Andalucía that started in 2004 and ended in 2010 (Junta de Andalucía, Medio Ambiente, http://www.juntadeandalucia.es/medioambiente/(accessed on 30 March 2022)).It also shows that, by 2016, the rainfall since 2010 allowed the depleted GWS to recover to a level almost comparable to that of the pre-drought period.According to [49], the southern half of the Iberian Peninsula received less than 40% of its usual precipitation by June 2005, with the consequent depletion of large aquifers in the peninsula [50].d) is a significant inflexion in rainfall that does not appear as strongly in the TWS accumulated time series; the variability in TWS is less than that of rainfall because rainfall includes TWS and the water that was lost by runoff and evapotranspiration.
The value of the TWS provided by GRACE can only be fully appreciated by considering the difficulty in estimating the state of the groundwater resources at the regional level of the Andalucía region (Figure 1). Figure 8 shows the distribution of groundwater bodies.Each groundwater body may comprise a single aquifer or several aquifers.In Figure 8, it is also possible to see large zones of low permeability (aquicludes).There is thus a problem in monitoring groundwater levels across a region with a high number of  d) is a significant inflexion in rainfall that does not appear as strongly in the TWS accumulated time series; the variability in TWS is less than that of rainfall because rainfall includes TWS and the water that was lost by runoff and evapotranspiration.
The value of the TWS provided by GRACE can only be fully appreciated by considering the difficulty in estimating the state of the groundwater resources at the regional level of the Andalucía region (Figure 1). Figure 8 shows the distribution of groundwater bodies.Each groundwater body may comprise a single aquifer or several aquifers.In Figure 8, it is also possible to see large zones of low permeability (aquicludes).There is thus a problem in monitoring groundwater levels across a region with a high number of groundwater bodies and aquicludes.An expert selection of eleven piezometers, which are generally considered representative of the natural behaviour of the corresponding aquifers and are not affected by pumping, are shown as solid black dots in Figure 8.For illustrative purposes, Figure 9 shows the evolution of water heads for four of these piezometers.The groundwater variability was estimated as the difference between each water-head time series and its mean followed by averaging the eleven time series.In addition, the time series was rescaled as an anomaly with respect to the period 2004-2009, which is the period to which GRACE data are usually referred.The time series of the anomalies are shown in Figure 10.The rainfall and TWS (GRACE, maximum degree L = 30) time series are referred to the same baseline to calculate the anomalies.The accumulated differences with respect to the mean are shown in Figure 11 together with the accumulated differences with respect to the mean for the groundwater time series in Figure 10.This Figure shows that the groundwater time series register the multi-year trends, and their behaviour is much smoother with respect to the rainfall and TWS time series, which also register the same trends together with other short-range trends as discussed in relation to Figure 7.The groundwater levels begin to decrease in 2005, and the minimum is reached in 2010, followed by recovery.mélange, metamorphic rocks) that are not proper aquifers but aquicludes and aquitards.The figure shows the spatial complexity of the distribution, which makes detailed monitoring of groundwater levels very difficult.In addition, most of the aquifers are subjected to pumping for irrigation and potable water and only a small number of piezometers can be regarded as being unaffected and measuring the real evolution of water levels.An expert selection of 11 representative piezometers (black dots) were used in this study; they are irregularly distributed and monitor only a small number of groundwater bodies.11.This is our best estimate of the evolution of GW in Andalucía.First, all anomalies for each piezometer were calculated with respect to each mea and these were then averaged.Next, the anomalies were referred to the mean of the 2004-2009 p riod that is covered by the GRACE data.11.This is our best estimate of the evolution of GWS in Andalucía.First, all anomalies for each piezometer were calculated with respect to each mean and these were then averaged.Next, the anomalies were referred to the mean of the 2004-2009 period that is covered by the GRACE data.
Figure 10.Mean groundwater head anomaly calculated as the average of the groundwater head anomaly of the 11 piezometers shown in Figure 11.This is our best estimate of the evolution of GWS in Andalucía.First, all anomalies for each piezometer were calculated with respect to each mean and these were then averaged.Next, the anomalies were referred to the mean of the 2004-2009 period that is covered by the GRACE data.The groundwater time series and the TWS data and rainfall data are compared further in Figure 12. Figure 12A shows good agreement between the evolution of the TWS anomalies and the groundwater anomalies with a delay in the peaks.The delay is because Figure 11.Accumulated differences from the mean for rainfall (orange line), TWS data (black line) for location 7 and groundwater (blue line), for anomalies with respect to the mean for the period 2004-2009.The groundwater curve shows the general trends described in Figure 7, with a decrease (depletion) over a five-year period and a partial recovery in the following six years.
The groundwater time series and the TWS data and rainfall data are compared further in Figure 12. Figure 12A shows good agreement between the evolution of the TWS anomalies and the groundwater anomalies with a delay in the peaks.The delay is because the TWS anomaly is due to an excess of water mass that includes the water in the soil and in the vadose zone prior to the recharge; this reaches the water table at the end of the recharge, thus affecting the groundwater level.The figure includes the coefficient of correlation and the values of the correlation function of different lags, where lag 1 is equal to one month and lag 2 is equal to two months.The coefficient of correlation is the correlation function for lag 0. The maximum correlation is for lag 1 and is equal to 0.69.Thus, there is an average time delay of one month between the gravity signal and the groundwater level because of the previously explained traveling time of the part of the rainfall water that, infiltrated into the soil, must travel across the vadose zone until it reaches the groundwater table.The effects of this traveling water on gravity are instantaneous, but they can only be detected by in situ monitoring when they arrive and increase the groundwater level.the TWS anomaly is due to an excess of water mass that includes the water in the soil and in the vadose zone prior to the recharge; this reaches the water table at the end of the recharge, thus affecting the groundwater level.The figure includes the coefficient of correlation and the values of the correlation function of different lags, where lag 1 is equal to one month and lag 2 is equal to two months.The coefficient of correlation is the correlation function for lag 0. The maximum correlation is for lag 1 and is equal to 0.69.Thus, there is an average time delay of one month between the gravity signal and the groundwater level because of the previously explained traveling time of the part of the rainfall water that, infiltrated into the soil, must travel across the vadose zone until it reaches the groundwater table.The effects of this traveling water on gravity are instantaneous, but they can only be detected by in situ monitoring when they arrive and increase the groundwater level.Figure 12B compares the rainfall anomalies and the groundwater anomalies and the same delay, for the same reasons, is evident between the peaks.There must be a positive Figure 12. (A) TWS (GRACE maximum degree L = 30) time series and groundwater head anomalies in which there is a good correspondence between both signals over the annual cycle and the multiyear trends.(B) Time series of the rainfall and groundwater head anomalies in which there is a deficit of rain (drought) and a negative anomaly of the groundwater head that may imply depletion of groundwater resources.(C) TWS (GRACE maximum degree L = 30) time series and rainfall anomalies in which there is a good correspondence between both signals but with higher variability in rainfall.The very large spikes in rainfall do not correspond to the acute spikes in TWS because there is a proportional loss of water mass by runoff.
Figure 12B compares the rainfall anomalies and the groundwater anomalies and the same delay, for the same reasons, is evident between the peaks.There must be a positive rainfall anomaly (humid period) if the groundwater levels recover, whereas there must be a negative rainfall anomaly (drought) if the groundwater levels decrease.Figure 12C compares the TWS anomalies and the rainfall anomalies from which similar conclusions can be drawn to those in the discussion of Figure 4B.
Finally, the long-term trend between TWS and GWS (groundwater heads) was assessed by filtering the annual trend using a moving average filter with a half-bandwidth of 12 months.The result of the filtering is shown in Figure 13, which clearly shows that the GWS trend is a delayed version of the TWS.Thus, the TWS can be used to monitor the GWS of the Andalucía region.This is a significant outcome for a region that has a complex distribution of aquifers and aquicludes and very limited monitoring.
rainfall anomaly (humid period) if the groundwater levels recover, whereas there must be a negative rainfall anomaly (drought) if the groundwater levels decrease.Figure 12C compares the TWS anomalies and the rainfall anomalies from which similar conclusions can be drawn to those in the discussion of Figure 4B.
Finally, the long-term trend between TWS and GWS (groundwater heads) was assessed by filtering the annual trend using a moving average filter with a half-bandwidth of 12 months.The result of the filtering is shown in Figure 13, which clearly shows that the GWS trend is a delayed version of the TWS.Thus, the TWS can be used to monitor the GWS of the Andalucía region.This is a significant outcome for a region that has a complex distribution of aquifers and aquicludes and very limited monitoring.

Conclusions
The hypothesis considered in this work is that GRACE observations and, in particular, the TWS time series, provide a means of investigating the regional evolution of groundwater resources based on rainfall estimates and, for the target region, most of the TWS signal is determined by variations in groundwater.The relationship between the TWS time series and rainfall data was analysed in detail.Using the accumulated rainfall differences and the accumulated TWS curves, our hypothesis is supported by the correlation between rainfall and the TWS time series in terms of temporal trends as shown graphically in Figure 7.More importantly, there is good correlation (numerical values in Figure 12A) between the mean groundwater and TWS.The cumulated differences of the TWS time series show that the GWS depletion over the long drought from 2004 to 2010 mostly been recovered by 2016.Furthermore, the estimated regional groundwater variation validates the previous findings and shows that the GRACE TWS is a smooth version of rainfall and, in turn, the groundwater is a smooth version of the TWS.TWS shows hysteresis-like behaviour between drought periods and humid periods.Our analyses of the periodogram of the rainfall and of the TWS time series show similar peaks at known frequencies representing cyclicities related to seasonal cycles and climate events (e.g., NAO phenomenon, sunspot cycles).Finally, we show that the pluri-annual GWS trend is a delayed version of Figure 13.Time series of trends for the TWS (GRACE maximum degree L = 30) and groundwater head anomalies.A moving average filter, with a half-bandwidth of 12 months, was used.The trend shows similar behaviour between TWS and GWS but with a delay in the groundwater signal.

Conclusions
The hypothesis considered in this work is that GRACE observations and, in particular, the TWS time series, provide a means of investigating the regional evolution of groundwater resources based on rainfall estimates and, for the target region, most of the TWS signal is determined by variations in groundwater.The relationship between the TWS time series and rainfall data was analysed in detail.Using the accumulated rainfall differences and the accumulated TWS curves, our hypothesis is supported by the correlation between rainfall and the TWS time series in terms of temporal trends as shown graphically in Figure 7.More importantly, there is good correlation (numerical values in Figure 12A) between the mean groundwater and TWS.The cumulated differences of the TWS time series show that the GWS depletion over the long drought from 2004 to 2010 mostly been recovered by 2016.Furthermore, the estimated regional groundwater variation validates the previous findings and shows that the GRACE TWS is a smooth version of rainfall and, in turn, the groundwater is a smooth version of the TWS.TWS shows hysteresis-like behaviour between drought periods and humid periods.Our analyses of the periodogram of the rainfall and of the TWS time series show similar peaks at known frequencies representing cyclicities related to seasonal cycles and climate events (e.g., NAO phenomenon, sunspot cycles).Finally, we show that the pluri-annual GWS trend is a delayed version of the TWS trend, and thus the filtered GWS time series can be used as a surrogate of the state of the groundwater resource, its depletion, and its recoveries without the need for an exhaustive monitoring network.It is unlikely that such a network will ever be implemented in this region (and in many other similar regions around the world) because of the costs of design, implementation, maintenance and the routine collection of measurements at piezometer sites.

Hydrology 2023 , 24 Figure 1 .
Figure 1.Relationship between rainfall and terrestrial water storage (TWS).TWS comprises surface water storage (SWS) in the surface zone, soil moisture storage (SMS) in the soil zone (SZ) and groundwater storage (GWS) in the saturated groundwater zone (SGZ).Water storage in the vadose zone (VZ) is considered to be in transit from the SMS to the GWS.P is precipitation, R is runoff, I is infiltration, ET is evapotranspiration, SWE is snow water equivalent, D is water in dams and Gd is groundwater discharge.(A): Following rainfall events, TWS increases as the result of changes in SWS, SMS, GWS and the water in transit across the VZ.Changes in GWS by recharge from rainfall are delayed by the transit through the VZ.(B): In the dry season, the TWS decreases because there is a decrease in SWS, SMS and GWS.In (A,B), grey color represents the soil zone (SZ), green color represents the vadose zone (VZ) and blue color represents the saturated zone (SGZ).

Figure 1 .
Figure 1.Relationship between rainfall and terrestrial water storage (TWS).TWS comprises surface water storage (SWS) in the surface zone, soil moisture storage (SMS) in the soil zone (SZ) and groundwater storage (GWS) in the saturated groundwater zone (SGZ).Water storage in the vadose zone (VZ) is considered to be in transit from the SMS to the GWS.P is precipitation, R is runoff, I is infiltration, ET is evapotranspiration, SWE is snow water equivalent, D is water in dams and Gd is groundwater discharge.(A): Following rainfall events, TWS increases as the result of changes in SWS, SMS, GWS and the water in transit across the VZ.Changes in GWS by recharge from rainfall are delayed by the transit through the VZ.(B): In the dry season, the TWS decreases because there is a decrease in SWS, SMS and GWS.In (A,B), grey color represents the soil zone (SZ), green color represents the vadose zone (VZ) and blue color represents the saturated zone (SGZ).

Figure 2 .
Figure 2. Geographical location of the study area in Andalucía, Southern Spain.The small purple dots are the centres of a regular grid of square cells with side lengths of 5 km.The large purple dots are eight locations at which the gravity-derived TWS time series are calculated.The red square is a 200 km × 200 km square centred at point 7.The dashed blue square represents a 400 km × 400 km square centred at point 7.The dashed green square represents a 400 km × 400 km square centred at point 6.

Figure 2 .
Figure 2. Geographical location of the study area in Andalucía, Southern Spain.The small purple dots are the centres of a regular grid of square cells with side lengths of 5 km.The large purple dots are eight locations at which the gravity-derived TWS time series are calculated.The red square is a 200 km × 200 km square centred at point 7.The dashed blue square represents a 400 km × 400 km square centred at point 7.The dashed green square represents a 400 km × 400 km square centred at point 6.

Figure 2 .Figure 3 .
Figure 2. Geographical location of the study area in Andalucía, Southern Spain.The small purple dots are the centres of a regular grid of square cells with side lengths of 5 km.The large purple dots are eight locations at which the gravity-derived TWS time series are calculated.The red square is a 200 km × 200 km square centred at point 7.The dashed blue square represents a 400 km × 400 km square centred at point 7.The dashed green square represents a 400 km × 400 km square centred at point 6.

Figure 3 .
Figure 3. (A): Terrestrial water storage (TWS) for the 144 months between 2002 and 2015 for location 7 and maximum degrees for the spherical harmonics of L = 30 and L = 50.There is an evident annual cycle superimposed on the other cycles of longer periods.The high-frequency variations increase from a maximum degree of L = 30 to L = 50.(B): Terrestrial water storages (TWSs) for the same maximum degree L = 50 for two different locations 6 and 7. Note that the point location is a reference at which the support of the measurement is centred.The values are averages of squares of around 400 km side length for maximum degree L = 50 and 680 km for maximum degree L = 30.(C): Time series of the uncertainties of the TWS time series for location 7 and where the uncertainty was evaluated by the standard error (SE) or square root of the variance of the error described by Equation (2).(D): Time series of TWS for locations 6 and 7 (as in (B)) but where the TWS ± SE is shown for location 7. The TWS for location 6 is within the TWS ± SE of location 7 and thus both time series of locations 6 and 7 are indistinguishable if the uncertainties are considered.

Figure 3 .
Figure 3. (A): Terrestrial water storage (TWS) for the 144 months between 2002 and 2015 for location 7 and maximum degrees for the spherical harmonics of L = 30 and L = 50.There is an evident annual cycle superimposed on the other cycles of longer periods.The high-frequency variations increase from a maximum degree of L = 30 to L = 50.(B): Terrestrial water storages (TWSs) for the same maximum degree L = 50 for two different locations 6 and 7. Note that the point location is a reference at which the support of the measurement is centred.The values are averages of squares of around 400 km side length for maximum degree L = 50 and 680 km for maximum degree L = 30.(C): Time series of the uncertainties of the TWS time series for location 7 and where the uncertainty was evaluated by the standard error (SE) or square root of the variance of the error described by Equation (2).(D): Time series of TWS for locations 6 and 7 (as in (B)) but where the TWS ± SE is shown for location 7. The TWS for location 6 is within the TWS ± SE of location 7 and thus both time series of locations 6 and 7 are indistinguishable if the uncertainties are considered.

Figure 4 .
Figure 4. (A).Mean areal monthly rainfall around location 7 using a window of 200 km × 200 km and total rainfall for the entire Andalucía region.(B).Comparison of the mean rainfall for the 200 km × 200 km (red rectangle in Figure2) and the TWS (mm) derived from the GRACE data using maximum degree L = 50.

Figure 4 .
Figure 4. (A).Mean areal monthly rainfall around location 7 using a window of 200 km × 200 km and total rainfall for the entire Andalucía region.(B).Comparison of the mean rainfall for the 200 km × 200 km (red rectangle in Figure2) and the TWS (mm) derived from the GRACE data using maximum degree L = 50.

Figure 4 .
Figure 4. (A).Mean areal monthly rainfall around location 7 using a window of 200 km × 200 km and total rainfall for the entire Andalucía region.(B).Comparison of the mean rainfall for the 200 km × 200 km (red rectangle in Figure2) and the TWS (mm) derived from the GRACE data using maximum degree L = 50.

Figure 5 .
Figure 5. Variation in the TWS using GRACE data for maximum degree L = 50 together with variation in mean temperature.The maximum temperature in summer corresponds to the minimum TWS, whereas the minimum temperature in winter corresponds to the maximum TWS due to the rainy season occurring in winter and spring.

Figure 5 .
Figure 5. Variation in the TWS using GRACE data for maximum degree L = 50 together with variation in mean temperature.The maximum temperature in summer corresponds to the minimum TWS, whereas the minimum temperature in winter corresponds to the maximum TWS due to the rainy season occurring in winter and spring.

Figure 6 .
Figure 6.(A).Periodogram for rainfall and TWS calculated from GRACE data for maximum degrees L = 30 (blue solid line) and L = 50 (red solid line).The main periodicities represent the yearly hydrologic cycle with a bimodal wet season and the climatic cycles.(B).Squared coherence spectrum between rainfall and TWS calculated from GRACE data for maximum degrees L = 50.(C).Achieved confidence level (%) of the squared coherence.(D).Phase spectrum (degree) between rainfall and TWS calculated from GRACE data for maximum degrees L = 50.

Figure 6 .
Figure 6.(A).Periodogram for rainfall and TWS calculated from GRACE data for maximum degrees L = 30 (blue solid line) and L = 50 (red solid line).The main periodicities represent the yearly hydrologic cycle with a bimodal wet season and the climatic cycles.(B).Squared coherence spectrum between rainfall and TWS calculated from GRACE data for maximum degrees L = 50.(C).Achieved confidence level (%) of the squared coherence.(D).Phase spectrum (degree) between rainfall and TWS calculated from GRACE data for maximum degrees L = 50.Hydrology 2023, 10, x FOR PEER REVIEW 16 of 24

Figure 7 .
Figure 7. Accumulated differences from the mean for rainfall (blue line) and TWS data (black line) for location 7, exhibiting three significant trends in the time series (red arrows).The first increases from the start of the time series at (a) to the maximum at (b).The second decreases from (b) to the minimum at (c).The third is from (c) to the end of the time series at (e).(d) is a significant inflexion in rainfall that does not appear as strongly in the TWS accumulated time series; the variability in TWS is less than that of rainfall because rainfall includes TWS and the water that was lost by runoff and evapotranspiration.

Figure 7 .
Figure 7. Accumulated differences from the mean for rainfall (blue line) and TWS data (black line) for location 7, exhibiting three significant trends in the time series (red arrows).The first increases from the start of the time series at (a) to the maximum at (b).The second decreases from (b) to the minimum at (c).The third is from (c) to the end of the time series at (e).(d) is a significant inflexion in rainfall that does not appear as strongly in the TWS accumulated time series; the variability in TWS is less than that of rainfall because rainfall includes TWS and the water that was lost by runoff and evapotranspiration.

Hydrology 2023 , 24 Figure 8 .
Figure8.The 189 polygons (shaded grey) are the Andalucía groundwater bodies, each comprising one or more aquifers and, in total, comprising 62% of the total surface of the region.Between the groundwater bodies there are large zones with low-permeability lithologies (e.g., flysch, detritic mélange, metamorphic rocks) that are not proper aquifers but aquicludes and aquitards.The figure shows the spatial complexity of the distribution, which makes detailed monitoring of groundwater levels very difficult.In addition, most of the aquifers are subjected to pumping for irrigation and potable water and only a small number of piezometers can be regarded as being unaffected and measuring the real evolution of water levels.An expert selection of 11 representative piezometers (black dots) were used in this study; they are irregularly distributed and monitor only a small number of groundwater bodies.

Figure 8 .
Figure8.The 189 polygons (shaded grey) are the Andalucía groundwater bodies, each comprising one or more aquifers and, in total, comprising 62% of the total surface of the region.Between the groundwater bodies there are large zones with low-permeability lithologies (e.g., flysch, detritic mélange, metamorphic rocks) that are not proper aquifers but aquicludes and aquitards.The figure shows the spatial complexity of the distribution, which makes detailed monitoring of groundwater levels very difficult.In addition, most of the aquifers are subjected to pumping for irrigation and potable water and only a small number of piezometers can be regarded as being unaffected and measuring the real evolution of water levels.An expert selection of 11 representative piezometers (black dots) were used in this study; they are irregularly distributed and monitor only a small number of groundwater bodies.

Figure 9 .
Figure 9.Time series of the water head evolution for four of the piezometers shown in Figure8for the purpose of illustrating the complex nature of the groundwater level evolution in different aquifers.

Figure 9 . 2 Figure 10 .
Figure 9.Time series of the water head evolution for four of the piezometers shown in Figure 8 for the purpose of illustrating the complex nature of the groundwater level evolution in different aquifers.Hydrology 2023, 10, x FOR PEER REVIEW 18 of 2

Figure 10 .
Figure 10.Mean groundwater head anomaly calculated as the average of the groundwater head anomaly of the 11 piezometers shown in Figure11.This is our best estimate of the evolution of GWS in Andalucía.First, all anomalies for each piezometer were calculated with respect to each mean and these were then averaged.Next, the anomalies were referred to the mean of the 2004-2009 period that is covered by the GRACE data.

Figure 11 .
Figure11.Accumulated differences from the mean for rainfall (orange line), TWS data (black line) for location 7 and groundwater (blue line), for anomalies with respect to the mean for the period 2004-2009.The groundwater curve shows the general trends described in Figure7, with a decrease (depletion) over a five-year period and a partial recovery in the following six years.

Figure 12 .
Figure 12.(A) TWS (GRACE maximum degree L = 30) time series and groundwater head anomalies in which there is a good correspondence between both signals over the annual cycle and the multiyear trends.(B) Time series of the rainfall and groundwater head anomalies in which there is a deficit of rain (drought) and a negative anomaly of the groundwater head that may imply depletion of groundwater resources.(C) TWS (GRACE maximum degree L = 30) time series and rainfall anomalies in which there is a good correspondence between both signals but with higher variability in rainfall.The very large spikes in rainfall do not correspond to the acute spikes in TWS because there is a proportional loss of water mass by runoff.

Figure 13 .
Figure 13.Time series of trends for the TWS (GRACE maximum degree L = 30) and groundwater head anomalies.A moving average filter, with a half-bandwidth of 12 months, was used.The trend shows similar behaviour between TWS and GWS but with a delay in the groundwater signal.

Table 1 .
Geographical co-ordinates of the 8 purple points represented in Figure2.

Table 1 .
Geographical co-ordinates of the 8 purple points represented in Figure2.