Potential groundwater recharge in north-western India vs spaceborne GRACE gravity anomaly based monsoonal groundwater storage change for evaluation of groundwater potential and sustainability

Journal Pre-proof Potential groundwater recharge in north-western India vs spaceborne GRACE gravity anomaly based monsoonal groundwater storage change for evaluation of groundwater potential and sustainability R.S. Chatterjee, Pranshu Pranjal, Sujit Jally, Bipin Kumar, Vinay K. Dadhwal, S.K. Srivastav, Dheeraj Kumar PII:

S2352-801X(19)30033-5

DOI:

https://doi.org/10.1016/j.gsd.2019.100307

Reference:

GSD 100307

To appear in:

Groundwater for Sustainable Development

Received Date: 26 January 2019 Revised Date:

29 August 2019

Accepted Date: 20 November 2019

Please cite this article as: Chatterjee, R.S., Pranjal, P., Jally, S., Kumar, B., Dadhwal, V.K., Srivastav, S.K., Kumar, D., Potential groundwater recharge in north-western India vs spaceborne GRACE gravity anomaly based monsoonal groundwater storage change for evaluation of groundwater potential and sustainability, Groundwater for Sustainable Development (2019), doi: https://doi.org/10.1016/ j.gsd.2019.100307. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier B.V.

Potential groundwater recharge vs. spaceborne GRACE gravity anomaly based GWS change in North-Western India during the observation period (2005-2012)

7

6

3 2

1

(a) (b)

(d) 4 (c)

9

8

(e) 10

5

Drainage basins (a) Upper Ganga (b) Yamuna (c) Chambal (d) Luni-Ghaggar (e) Narmada Indian states 1. Delhi 2. Haryana 3. Punjab 4. Rajasthan 5. Uttar Pradesh 6. Uttarakhand 7. Himachal Pradesh 8. Gujarat 9. Madhya Pradesh 10. Maharashtra

High potential groundwater recharge and high positive groundwater storage change Low potential groundwater recharge but high positive groundwater storage change High potential groundwater recharge but low positive groundwater storage change Low potential groundwater recharge and low positive groundwater storage change High potential groundwater recharge but negative groundwater storage change Low potential groundwater recharge and negative groundwater storage change

Potential Groundwater Recharge in North-Western India vs Spaceborne GRACE Gravity Anomaly based Monsoonal Groundwater Storage Change for Evaluation of Groundwater Potential and Sustainability R.S. Chatterjee1*, Pranshu Pranjal1,2, Sujit Jally1, Bipin Kumar1, Vinay K Dadhwal1,3,S.K. Srivastav1 and Dheeraj Kumar2 1

Indian Institute of Remote Sensing (Indian Space Research Organization), 4-Kalidas Road, Dehradun-248001, Uttarakhand, INDIA

2.

Indian Institute of Technology (Indian School of Mines), Dhanbad-826004, Jharkhand, INDIA

3.

Indian Institute of Space Science and Technology, Thiruvananthapuram-695547, Kerala, INDIA

*

Corresponding Author

Tel.: 91-135-2524156; Fax: 91-135-2741987, 2748041 Email: [email protected]

Abstract

The north-western part of India has been affected by severe ground water depletion. In this study, an attempt has been made to assess potential ground water recharge of North West (NW) India and its comparison with Gravity Recovery and Climate Experiment (GRACE) gravity anomaly based monsoonal groundwater storage change. An assessment of potential groundwater recharge was done in terms of terrain suitability for groundwater recharge and average monsoonal recharge. Terrain suitability was assessed based on hydrogeology and terrain parameters which were then integrated with the spatial variability in monsoonal recharge (expressed as the surface water column available for underground infiltration). The aaverage monsoonal groundwater storage (GWS) change was derived from GRACE terrestrial water storage (TWS) change. An evaluation of 1

potential groundwater recharge vs. terrain suitability for groundwater recharge and potential groundwater recharge vs. monsoonal GRACE GWS change shows that there is high terrain suitability but low potential groundwater recharge in the north-western part of NW India essentially due to very low precipitation. In the north-western part, terrain suitability is high primarily due to the presence of sandy alluvium, deeper sandy loam soil and flat terrain. The remaining alluvial plains, lithologically covered by unconsolidated sediments with loamy-silty-sand or clayey-loamysilty-sand soils, show moderate terrain suitability. Rainfall induced surface water infiltration capacity is high in Upper Ganga drainage basin and in major parts of Narmada drainage basin but low to very low in Luni-Ghaggar drainage basin and in major parts of Yamuna and Chambal drainage basins. The average monsoonal GRACE GWS change is low in the north-western part of NW India which is in agreement with potential groundwater recharge scenario. The comparison between precipitation, terrain infiltration potential and GRACE gravity based average groundwater storage observations enables to evaluate prospective and vulnerable areas in terms of groundwater availability and sustainability in NW India.

Keywords: Potential groundwater recharge; GRACE TWS change; groundwater storage change; NW India.

1. Introduction Groundwater depletion or lowering of groundwater level is one of the major issues controlling fresh water availability (Khurana and Mahapatra, 2008; NWP, 2012). The demand for water has increased rapidly in the recent years which leads to water crisis not only in India but also

2

in many other parts of the world (Hinrichsen and Tacio, 2002; NWP, 2012). Rapid population growth, urbanization, irregularity in rainfall and groundwater-intensive agriculture are major causes of groundwater depletion (Garduño et al., 2011; Wyrwoll, 2012; UNICEF 2013; Steward, 2014). India has 16 percent of the world’s total population but only 4 percent of the world’s fresh water resources. Groundwater plays an important role for high agricultural productivity in Punjab, Haryana, parts of Rajasthan and Uttar Pradesh. As a result, groundwater exploitation has been increasing at an alarming rate leading to severe groundwater depletion in this region. In India, groundwater extraction in agricultural, domestic and industrial sectors constitutes approximately 89%, 9% and 2% respectively (Suhag, 2016; CGWB, 2019). The primary cause of groundwater extraction has been the high demand in the agricultural sector. In fact, in the last four decades, there is an increase of roughly 84% net irrigated area through groundwater (Suhag, 2016). Groundwater extraction scenario can be indirectly assessed through present groundwater development stages in different parts of the country (Figure 1a). Groundwater development (in percentage) is a ratio of the annual gross groundwater extraction for all uses to the net annual groundwater availability multiplied by 100. Considering groundwater development in the past 20 years, groundwater extraction has been assessed into four categories such as safe, semi-critical, critical and overexploited categories for groundwater development ≤ 70%, 70-90%, 90-100% and > 100% respectively. Groundwater extraction is very high in NCR Delhi, Haryana, Punjab, Rajasthan, UT of Chandigarh and some parts of Gujarat, Uttar Pradesh and Madhya Pradesh with the level of groundwater development > 100% (Suhag, 2016; CGWB, 2019). The Gravity Recovery and Climate Experiment (GRACE) satellite mission was launched by National Aeronautics and Space Administration (NASA) and the German Aerospace Center (DLR) for monitoring the Earth’s gravity field and its temporal variation. Monthly gravity variations obtained from GRACE mission 3

is inverted to obtain TWS change on spatial grids of few hundred kilometres. Using the monthly gravity anomaly, GRACE can sense the changes in the water storage at surface and sub-surface level including groundwater storage change. Satellite-based GRACE gravity anomaly study of Rodell et al (2009) and Tiwari et al (2009) highlighted severe ground water depletion in the northwestern part of India. Subsequently, several researchers conducted similar studies and corroborated groundwater depletion in NW India and Indo-Gangetic basin (Chen et al., 2014; Dasgupta et al., 2014). GRACE gravity observation on Groundwater storage depletion has also been successfully studied in other parts of the globe such as Mississippi River Basin, USA (Rodell et al., 2007), California Central Valley, USA (Famiglietti et al., 2011), Middle East and Arabian countries (Voss et al., 2013; Richey, 2014), Canning Basin, Australia (Richey, 2014), Central Asia (Li et al., 2018), and in different parts of China (Feng et al., 2018; Gong et al., 2018; Chen et al., 2019). To avoid the situation of severe groundwater depletion, the usage of groundwater should be controlled at par the rate of replenishment (De vries, 1984). Towards this mission, an in-depth knowledge on potential groundwater recharge is utmost important. It is indeed essential to develop a geospatial framework model for generating and organizing the data layers such as geomorphology, lithology, digital surface model (DSM), land use land cover (LU LC), soil depth and texture, and monsoonal precipitation as available from Earth observation and field-based data sources. In the present study, hydrogeology and terrain parameters were coupled with monsoonal surface water column available for underground infiltration for deriving potential groundwater recharge of the study area. It may be noted that monsoonal rainfall is the primary source of groundwater recharge in the country (Sinha and Sharma, 1988). In NW India, groundwater recharge dominantly occurs during the southwest monsoon (June - October). GRACE gravity anomaly based 1°×1° grid data of monthly terrestrial water storage (TWS) change during the monsoonal months of 2005-2012 was used to 4

derive average monsoonal groundwater storage (GWS) change. Finally, potential groundwater recharge in NW India was compared with average monsoonal GWS change for comparative evaluation of groundwater sustainability in NW India.

2. Materials and methods 2.1 Study area The study area includes the north-western part of India covering the states of National Capital Region (NCR) Delhi, Haryana, Rajasthan, Uttar Pradesh (U.P.), Uttarakhand, south-eastern parts of Punjab and Himachal Pradesh (H.P.), Gujarat, Madhya Pradesh (M.P.) and north-western part of Maharashtra. We considered five major drainage basins of NW India such as Upper Ganga, Yamuna, Chambal, Luni-Ghaggar, and Narmada drainage basins as hydrological units for assessing groundwater recharge scenario in the region (Figure 1b). In general, the NW India is divided into three physiographic zones such as Indo-Gangetic Plain (in north and east), Deccan Plateau (in south and west) and Thar Desert (in west). The region receives moderate to low rainfall and essentially during the monsoon season (June - October). During the remaining period of the year, the region remains more or less dry. The extreme western part of the study area namely the western part of Rajasthan state experiences very low rainfall (<100 mm). The major rock types in the study area include unconsolidated sediments, sedimentary rocks, igneous rocks, and combination of igneous, metamorphic and sedimentary rocks belonging to the geological groups and formations such as

5

6

Indo-Gangetic alluvium, Deccan trap and Vindhyans. Indo-Gangetic alluvium consists of unconsolidated Quaternary sediments comprising recent alluvium, older alluvium and aeolian sediments. The sediments are made of sands, silts, clays, pebbles etc. They form alluvial aquifers with high infiltration capacity, high storage and transmissivity. Deccan trap consists of alternate layers of compact and vesicular basaltic lava flows and inter-trappean clastic sedimentary beds. It has, in general, low infiltration capacity, low storage and transmissivity except for highly weathered, fractured and vesicular basalts occurring in limited areas, it has moderate infiltration, storage and transmissivity. Vindhyan rocks consist of consolidated clastic and non-clastic sedimentary an meta-sedimentary rocks such as sandstone, shale, conglomerate, limestone, slate and quartzite with primary and/or secondary porosity. The Vindhyan rocks are characterized by low to moderate infiltration capacity, low to high storage and transmissivity depending on the availability of secondary porosity in addition to primary porosity of consolidated Proterozoic sandstone and conglomerate. The study area is characterized by a wide variety of soils such as alluvial soil, desert soil, sandy soil, black soil and red soils. A variety of geomorphological processes have been operated in the study area which lead to the development of diverse landforms such as fluvial landforms, aeolian landforms, volcanic landforms and structural landforms. Due to diverse climate and physiographic conditions, diverse vegetation types including tropical dry deciduous forest, tropical moist deciduous forest, tropical semi evergreen forest, tropical thorn forest, sub-tropical pine forest and littoral and swamp vegetation occur in the region (ISFR, 2011).

2.2 Data preparation and analysis 2.2.1 Earth Observation (EO) data

7

Indian Remote Sensing (IRS) Resourcesat-2 Advanced Wide Field Sensor (AWiFS) multispectral data available in visible and near infrared (VNIR) and short wave infrared (SWIR) bands were used for false colour composite (FCC) image generation. The images with a spatial resolution of 56 m and a combined swath of 740 km were acquired during November-December, 2013. Enhanced FCC images were used for refinement and updation of the existing thematic data layers such as geomorphology, lithology, and land use land cover maps to update the boundaries between lithological units, to refine geomorphological classification scheme in terms of groundwater potential and to address LU LC dynamics in the study area.

2.2.2 Digital Elevation Model (DEM) and terrain analysis data We used Shuttle Radar Topography Mission (SRTM) 1 arc second (approximately 30 m spatial resolution at equator) DEM available in 1 degree x 1 degree tiles (USGS, 2015). The DEM represents a seamless and complete elevational surface of the Earth after applying a number of post processing steps by the data provider and by spatial interpolation of the remaining no-data regions. The DEM is available in geographic (lat/long) projection system with the WGS84 horizontal datum and the EGM96 vertical datum. The relative vertical accuracy of the DEM is reported to be 12-15 m which allows to achieve a Digital Terrain Elevation Data Level 2 (DTED-2) land surface model and precise terrain-related thematic data layers. From the DEM, specific catchment area and topographic slope were determined using ArcGIS tools. Finally, Topographic Wetness Index (TWI) map was generated using specific catchment area and topographic slope of the study area as follows: TWI = ln(Specific catchment area)/slope (Figure 2). TWI is a key element for identifying ground water recharge potential as it measures the effect of local topography on hydrological

8

processes and describes the spatial distribution of soil moisture and surface saturation (Wilson and Gallant, 2000).

2.2.3 Geomorphology Geomorphology describes Earth surface processes and resulting landforms with emphasis on their shape, size, relief characteristics, material composition, spatial distribution and evolution. Geomorphology governs surface runoff, groundwater recharge and discharge potential of the landforms as a function of their relief characteristics and material composition. The existing geomorphological map (1:2,000,000 scale) of Geological Survey of India (GSI, 2002) was used after needful updation and modification based on visual analysis of AWiFS standard FCC image. 9

The geomorphological classes were regrouped as a function of ground water recharge potential to obtain simplified applied geomorphological classes (Figure 3). The geomorphological classes in the study area were broadly combined into 3 groups: (1) hills and plateau, (2) plains and (3) valleys

with an interface group between ‘hills and plateau’ and ‘plains’ namely ‘slopes and badlands’. The geomorphological group ‘hills and plateau’ comprises the geomorphological classes with higher elevations under hills as well as plateau. The geomorphological group ‘plains’ occupies more than 66 % of the study area. It comprises a diverse variety of geomorphological classes such as rocky terrace, rocky desert, planation surface with or without duricrust, low-lying flats, ranns (low lying salt deserts), prodelta, islands and lagoons, alluvial plains and depositional terraces, aeolian plains with sand dunes and sand sheets. There are limited number of geomorphological classes under the geomorphological group ‘valleys’ such as playas, valleys and infilled valleys.

10

2.2.4 Lithology Lithology or rock types describe hydrogeological characteristics in terms of porosity and permeability of rock materials to assess groundwater recharge potential. Loose unconsolidated sediments or alluvium generally possess high porosity and permeability and are therefore favourable for high level of groundwater recharge. Primary porosity is essentially found in clastic sedimentary rocks whereas secondary porosity may be observed in igneous, metamorphic and sedimentary rocks due to mechanical processes such as brittle deformation and fracturing and/or geochemical processes such as dissolution. Geologically younger sedimentary rocks have generally higher primary porosity than that of the same type of geologically older rocks because of lesser degree of compaction. The lithological map of Geological Survey of India on 1: 2,000,000 scale (GSI, 1998) was visually compared with AWiFS FCC image and the broad lithological formation contacts were updated. Various lithological classes were regrouped based on similar rock types and

11

geological age as a function of their favourability for groundwater recharge. The broad rock types include unconsolidated sediments or alluvium, sedimentary rocks comprising clastic and non-clastic rocks, igneous rocks comprising ultrabasic-ultramafic rocks, intermediate-basic rocks and granitic rocks, a combination of igneous and metamorphic rocks, and a combination of igneous, metamorphic and sedimentary rocks and ophiolites. Among these, unconsolidated Quaternary sediments or alluvium occupies approximately 50% of the study area and primarily contributes groundwater recharge in the study area. Finally, based on the similarity in hydrgeoological properties, broad applied lithologic classes were made from a large number of lithological classes available in the study area (Figure 4).

2.2.5 Soil texture and soil depth Soil texture and soil depth map of India prepared by National Bureau of Soil Survey and Land Use Planning, Nagpur, India on 1: 2,000,000 scale (NBSS-LUP, 2002) was used in this study. We regrouped a wide variety of soil texture classes into four broad classes in terms of groundwater recharge potential: (1) fine texture (e.g., Clay, loamy clay, silty clay, sandy clay, sandy clay loam, silty clay loam), (2) medium texture (e.g., loam, silty loam, silt, sandy loam), (3) coarse texture (e.g., loamy sand, sand), and (4) others (e.g., rocky outcrop, built-up) with very coarse texture or without soil cover (Figure 5a). Similarly, we regrouped a wide range of soil depth into four classes: (1) deep soil (>50 cm), (2) moderately deep soil (25-50 cm), (3) shallow soil (10-25 cm), and (4) others with very shallow (<10 cm) or no soil cover (Figure 5b). Finally, the regrouped soil texture and soil depth classes were combined as a function of groundwater recharge potential. In this process, we obtained 3 major soil texture classes with their subclasses of soil depth except for the fourth class named as ‘others’ class. 12

2.2.6 Land use land cover The land use land cover map was essentially adapted from Agrawal et al. (2003) and modified based on IRS Resourcesat-2 AWiFS multispectral data of November-December, 2013. In the study 13

area, 13 different types of land use land cover classes were mapped such as settlement, wood or forest land (good cover), wood or forest land (poor cover), pasture range (poor), desert shrub (poor), meadow (good), meadow (poor), agriculture, bare soil, mangrove, mud flat, salt pan, water body/ snow cover (Figure 6). Except for water body/ snow cover, the LU-LC classes in the study area possess characteristically different groundwater recharge potential.

2.2.7a Rainfall Daily and monthly precipitation data for the study area were generated by combining satellite based and ground based observations. Satellite based daily rainfall measurements with spatial resolution 0.25°x0.25° was obtained from the Tropical Rainfall Measuring Mission (TRMM), a joint space mission of NASA, USA and JAXA, Japan.

14

Ground based rainfall data at 0.25°x0.25° grid was obtained from the Indian Meteorological Department (IMD, 2014). The gridded rainfall data for the entire country was generated by IMD by interpolating ground based rainfall measurements available from 6995 rainfall stations. We used essentially IMD rainfall data for those grids where one or more rainfall stations exist. For the remaining grids, we took weighted average (correlation coefficient based weighted average between IMD and TRMM rainfall data) of TRMM and IMD gridded rainfall data where both the data are available. Finally, the combined product was passed through a 3x3 low-pass median filter to smoothen sharp edges between neighbouring grids. Daily rainfall data at 0.25°x0.25° grid was later used for run-off estimation. Daily rainfall data was also used to obtain cumulative monthly rainfall during the observation period 2005-2012.

2.2.7b Runoff In the present study, the SCS runoff curve number method (Kent, 1973) was adopted to calculate runoff from daily rainfall data at 0.25°x0.25° spatial grids. The curve number depends on the hydrologic soil group, land use/ treatment and hydrologic condition. The hydrologic soil groups, land use/ treatment and slope maps were prepared from the soil map (modified after NBSS-LUP, 2002), land use land cover map (modified after Agarwal et al., 2003) and SRTM 1 arc second DEM.

2.2.7c Evapotranspiration Evapotranspiration plays a crucial role in hydrological water balance. MOD16 global evapotranspiration product is derived from Earth land surface by using Terra/Aqua MODIS remote sensing data. It provides global evapotranspiration (ET) and potential evapotranspiration (PET) at 1 15

km x 1 km regular grids for the 109.03 million km2 global vegetated land areas at 8 days, monthly and annual intervals (Mu et al., 2013). Terrestrial ET includes water vapour fluxes from soil evaporation, wet canopy evaporation and plant transpiration at dry canopy surface. MODIS 8-day fraction of photosynthetically active radiation (FPAR) is used as vegetation cover faction to quantify how much surface net radiation is allocated between soil and vegetation. MODIS 8-day albedo and daily downward solar radiation and air temperature available from daily meteorological reanalysis data are used to calculate surface net radiation and soil heat flux. Besides, daily air temperature, vapor pressure deficit (VPD) and relative humidity data, and 8-day MODIS leaf area index (LAI) are used to estimate surface stomatal conductance, aerodynamic resistance, wet canopy, soil heat flux and other key environmental variables. MOD16 ET data products are essentially generated by the algorithm of Mu et al (2007, 2011) using MODIS land cover, FPAR/LAI data and global surface meteorology data available from the Global Modeling and Assimilation Office (GMAO v. 4.0.0). In the present study, MOD 16 monthly ET data was used after upscaling into 0.25°x0.25° grids by resizing based on pixel aggregated averaging. Finally, the volume of surface water available for underground infiltration was calculated by subtracting monthly runoff and monthly ET from the monthly rainfall data products. We used the average surface water column thickness in each grid (0.25°x0.25°) potentially available for underground infiltration during the monsoon (Figure 7) to translate terrain suitability for groundwater recharge map into potential groundwater recharge map of the NW India.

16

2.2.8 Soil moisture change In this study, soil moisture data of Global Land Data Assimilation System (GLDAS) Noah land surface

model

(Rodell

and

Beaudoing,

(https://mirador.gsfc.nasa.gov/collections/GLDAS_NOAH025_M__001.shtml) 17

was

2007) used.

Soil

moisture is available in 4 layers: first layer 0-10 cm, second layer for 10-40 cm, third layer 40-100 cm and fourth layer 100-200 cm. By summing up the values in four layers, the volumetric soil moisture (volumetric soil water content) was obtained. Soil moisture being a state function, the monthly soil moisture change was estimated with respect to a long-term average.

2.2.9 GRACE gravity observation and TWS change The GRACE satellite mission was launched by NASA and DLR for monitoring the Earth’s dynamic gravity field in terms of monthly gravity anomaly. The gravity anomaly of the Earth surface may be inverted to derive spatio-temporal variation of mass and therefore terrestrial water storage (TWS) change around the planet. By measuring the changing distance between the two satellites and combining that with precise geocentric positions obtained from the onboard GPS instruments, the gravity anomaly over the Earth's surface is produced. GRACE level-2 monthly product consists of a series of spherical harmonic coefficients which describes the shape of the Earth's gravity field. The level-2 products are converted to level-3 products describing terrestrial water storage change in terms of the deviations from the series mean by averaging kernels (e.g., Wahr et al., 1998, 2004; Seo et al., 2006; Swenson and Wahr, 2006; Landerer and Swenson, 2012). Terrestrial Water Storage (TWS) represents a land surface state which may be compared with other hydrological states or water storage components of the land surface. It is the sum of surface water storage, groundwater storage, soil moisture, accumulated snow, and plant canopy water storage. However, GRACE observation cannot provide an absolute measurement of TWS but a relative measurement of the anomaly with respect to a long-term average. In the present study, GRACE Level-3 TWS change data products available in 1° x 1° grids were used. We downscaled the GRACE monthly TWS change from 1° x 1° grids to 0.25° x 0.25° grids using high resolution EO18

based and in-situ hydrological layers such as surface run-off, ET, soil moisture and groundwater level fluctuation data by linear regression based predictive analysis. For presentation of the results, we further downscaled GRACE monthly TWS change data from 0.25° x 0.25° grids to 10 km x 10km grids by inverse distance weighted (IDW) algorithm based statistical interpolation of grid centre values.

2.3 Spatial data integration The data layers used in this study were made available from different sources such as Resourcesat-2 AWiFS EO data, SRTM DEM, geomorphology and lithological data from Geological Survey of India, soil texture and soil depth data from NBSSLUP, existing land use land cover map prepared from SPOT Vegetation data, TRMM and IMD rainfall based surface water column thickness available for underground infiltration, GRACE gravity based TWS change data. As a result, the scale and projection systems are different those are not compatible for spatial comparison and integration. All the derived data layers were reprojected into Albers Conical Equal Area projection and subsequently converted to raster layers with 10 km x 10 km spatial grids for spatial comparison and integration.

2.3.1 Layer weight assignment For determination of layer weight, we adopted Analytic Hierarchy Process (AHP) (Saaty, 1980, 1994, 2008; Coyle, 2004; Saaty and Vargas, 2012). It is an effective tool for setting the priorities in a complex multi-variate decision making problem based on a series of pairwise comparisons. It overcomes the bias in the decision making to a great extent by ensuring consistency in the process. The first step in this process is to construct a pair-wise comparison matrix for all the elements such 19

as topography (or topographic wetness describing surface water flow accumulation per unit cell), geomorphology, lithology, soil texture and soil depth, and land use land cover those are contributing to “terrain suitability for groundwater recharge”. In pair-wise comparison, an index from 1 to 9 indicates the degree of importance: (1) equally important – 1, (2) somewhat more important (moderate) – 3, (3) much more important (strong) – 5, (4) very much more important (very strong) – 7, (5) absolutely more important (extreme) – 9 and intermediate values 2, 4, 6 and 8 when some compromise is needed. The basic principle of the rating scale is that if an element ‘M’ is absolutely more important than another element ‘N’ then ‘N’ must be absolutely less important than ‘M’ and the ratings for ‘M’ and ‘N’ would be 9 and 1/9 respectively. Thus for 5 layers, a square reciprocal matrix A = (aij) of order 5 is formed for aij = 1/ aij for i ≠ j and aii = ajj = 1 for all i and j equal to 1 to 5. The next step is the calculation of a list of the relative weights or Eigen vectors such that Aω = λω, where ω is said to be an Eigen vector of order n = 5 and λ is an Eigen value. The elements of the Eigen vector of relative importance (A) is calculated as follows: for each row, multiply all the row elements and take the nth root of the product which is then divided by the sum (sum of the nth root of the product of the row elements for all the rows). The product Aω vector with five elements gives five estimates of λmax by dividing the five elements of the product vector Aω with the corresponding element of the Eigen vector. All the five estimates of λmax are always ≥ 5, else there has been error in the calculation. Finally, the mean of these five estimates gives the value of λmax. The difference between λmax and n is an indication of the inconsistency in the judgements. Therefore, Consistency Index (CI) followed by Consistency Ratio (CR) are calculated to check inconsistency. When the value of CR is zero, the judgements are considered to be perfectly

20

consistent. In practice, when CR is < 0.1, the judgements are considered to be consistent and therefore trustworthy (Raharjo and Endah, 2006).

2.3.2 Class weight assignment Individual class weightage for each of the data layers was determined on 1-9 linear scale based on their relative significance towards groundwater recharge. For this purpose, we first arranged the classes of a layer in terms of their relative significance to groundwater recharge and regrouped into hydrogeologically similar classes. Now, based on the number of groups or categories in a data layer, the entire 1-9 weight range was linearly distributed. Further, depending on the number of classes in a group or category, the weight range was linearly allocated based on next level criteria such as process of formation, grain size, geological age or degree of fracturing etc. for lithology layer, geomorphic process, ruggedness, slope or relative age of formation etc. for geomorphology layer. In this process, the recharge potential behavior of different classes in a data layer becomes more or less linear and subjectivity in the weight scale is avoided. The procedure for assigning class weightage for individual layers is elaborated below.

2.3.3 Integration of data layers Finally, the data layers were integrated by Index Overlay method using AHP based layer weights and knowledge-based linear class weights to obtain integrated terrain suitability map of the study area for groundwater recharge (Figure 8). Terrain suitability for groundwater recharge represents the capability of the study area for potential groundwater recharge but it does not take into account the availability of surface water column available for underground infiltration as a function of rainfall and terrain condition. TRMM and IMD rainfall based monsoonal surface water column 21

available for underground infiltration was therefore integrated with terrain suitability layer to obtain potential groundwater recharge during the observation period.

3. Results and discussion 3.1 Terrain suitability for groundwater recharge Terrain suitability for groundwater recharge was prepared using five terrain parameters such as topographic wetness index (TWI), applied geomorphology, lithology, soil texture and depth, and land use land cover of the study area by AHP based layer weights and linear class weights based on relative significance (Figure 8). It is observed that the terrain suitability for groundwater recharge is high in Luni-Ghaggar drainage basin in the western part of the study area primarily due to the presence of sandy alluvium, deeper sandy loam soil and flat terrain. It is high in some patches of Upper Ganga and Yamuna drainage basins covering valleys and infilled valleys. Besides, a large part of the study area, geomorphologically occupied by alluvial plain, lithologically covered by unconsolidated sediments with loamy-silty-sand or clayey-loamy-silty-sand soil types shows moderate terrain suitability for groundwater recharge. The area dominantly falls under Upper Ganga and Yamuna drainage basins covering U.P., Haryana and NCR Delhi; northern part of LuniGhaggar drainage basin and central part of Chambal drainage basin covering western part of Haryana, southern part of Punjab and eastern part of Rajasthan states; and Sabarmati-Mahi, Narmada and Tapi alluvial plains of Narmada drainage basin covering parts of Gujarat, M.P. and Maharashtra states. The remaining area is occupied by adverse geomorphology such as hills and plateau or slopes and badlands, adverse lithology such as igneous rocks or a combination of igneous metamorphic and sedimentary rocks, adverse soil types such as rocky non-soil or thin clayey-

22

loamy-silty soil and unfavourable land covers. The area shows fair to poor terrain suitability for groundwater recharge.

3.2 Potential groundwater recharge The monsoonal surface water column available for underground infiltration was integrated with the terrain suitability for groundwater recharge layer to obtain potential groundwater recharge in the study area. For this purpose, the surface water column available for underground infiltration during the monsoon was rescaled into 1-9 linear scale similar to terrain suitability for groundwater recharge layer. The two layers were integrated as an algebraic product raised to the power ‘γ’ and (1-γ) as follows (Bonham-Carter, 1994):
 γ  γ

(2)

23

Where, ‘µ’ represents the integrated value of two layer variables ‘A’ and ‘B’ (e.g., terrain suitability for groundwater recharge and surface water column available for underground infiltration), and the value of ‘γ-operator’ varies in the range of 0.1 – 0.9 to consider the effect of extreme values in the integration product. For example, where one of the variables (say, layer ‘A’) is very low (rescaled values 1-3) and the other variable is sufficiently high (rescaled values >3-9), the effect of the variable ‘A’ becomes critical and we assume γ = 0.9-0.7 for layer ‘A’ and 1-γ = 0.1-0.3 for layer ‘B’ respectively. But, where both the variables are very low (rescaled values 1-3) or both the variables are sufficiently high (rescaled values >3-9), the effect of the variables will be essentially governed by their values and we assume γ = 1-γ = 0.5 (Bonham-Carter, 1994; Marjanovic and Caha, 2011). Finally, we obtained potential groundwater recharge map of the study area (Figure 9). It is observed that potential groundwater recharge is high in the major part of Upper Ganga drainage basin covering U.P. and Uttarakhand states excluding the northern hilly region of Uttarakhand and south-eastern part of U.P., and in parts of Narmada drainage basin covering Narmada-MahiSabarmati alluvial plains in M.P. and Gujarat. It is also found to be high in the northern part of Yamuna drainage basin covering south-western part of Uttarakhand state essentially due to high rainfall and therefore high surface water column available for underground infiltration. Potential groundwater recharge is found to be low in Luni-Ghaggar drainage basin particularly in the areas covering desert region of Rajasthan state and mud flat and salt pan covered Kutch region of Gujarat state. It is also found to be low in the hilly northern part of Upper Ganga drainage basin in Uttarakhand state, sloping rugged terrain in the north-western part of Narmada drainage basin in Rajasthan state and Chambal-Yamuna gullied badlands in Chambal and Yamuna drainage basins. In the remaining areas, potential groundwater recharge is found to be moderate. While comparing with terrain suitability for groundwater recharge, we observed a disagreement in Luni-Ghaggar 24

drainage basin, in parts of Narmada drainage basin and in the northern part of Upper Ganga drainage basin essentially due to inversely proportional surface water availability (very low surface water column thickness in areas of high terrain suitability and vice versa) for underground infiltration in these regions. In the Kutch region of Luni-Ghaggar drainage basin, both the factors are unfavourable which lead to extremely low potential groundwater recharge.

3.3 GRACE monsoonal GWS change GRACE TWS change data represent mass anomaly in centimetres of equivalent water thickness relative to the 2004-2009 time-mean baseline. Post-processing of GRACE observation was done using a destriping filter and Gaussian filter followed by truncation of the gravity signal at a spectral 25

degree and order 60. The quantitative estimates of signal attenuation and leakage error due to postprocessing filters were determined with respect to simulated monthly TWS anomalies obtained from GLDAS-Noah land surface model data (Swenson and Wahr, 2006; Landerer and Swenson, 2012). Besides, GRACE gravity variations are caused due to the changes in the mass distribution in the solid Earth and at the Earth’s surface. The major sources of gravity variation due to solid Earth processes include glacial isostatic adjustment (GIA) and earthquakes (at local scale). The effect of GIA has been incorporated by the data provider based on the contemporary GIA correction model (Peltier, 2004; Chambers et al., 2010; Wahr and Zhong, 2013). After correcting the effects of solid Earth processes, GRACE gravity variation successfully represents the water mass exchanges at the Earth surface. The GRACE TWS change represents the sum of the changes in surface water storage, groundwater storage, soil moisture, accumulated snow and plant canopy water storage. In the present study, the monthly surface water storage change was considered as the sum of monthly aggregated surface runoff and ET. The monthly soil moisture change was derived with respect to a long-term average covering the same period as GRACE TWS long-term average. It may be noted that the changes in accumulated snow and plant canopy water storage are not appreciable in the study area due to their limited spatial coverage (accumulated snow cover is confined to high altitude NE boundary region of Uttarakhand and H.P. states of the study area and vegetation cover is restricted to rocky terrain over hills, plateau, slopes and bad land regions those are not favourable for groundwater recharge). From the monthly TWS change, we subtracted monthly surface water storage change and monthly soil moisture change to obtain spaceborne GRACE gravity anomaly based monthly GWS change. Finally, GRACE monthly GWS change layers were combined for the monsoonal period (June-October) to produce average monsoonal GRACE GWS change of the study area for each year. It was observed that spaceborne GRACE gravity anomaly based average 26

monsoonal GWS change is low in Luni-Ghaggar drainage basin and north-western part of Yamuna drainage basin covering desertic arid region of Rajasthan, whole of Haryana and NCR Delhi, southern part of Punjab, western part of U.P., hilly northern part of Uttarakhand and Kutch region of Gujarat states. On the other hand, it is moderate to high in Narmada drainage basin, southern part of Chambal and Yamuna drainage basins, and eastern part of Upper Ganga drainage basin covering eastern and southern parts of Gujarat and northern part of Maharashtra states, M.P. (excluding northern part), and eastern part of U.P. states respectively.

3.4 Comparison of GRACE GWS change with in-situ groundwater level based GWS change We used groundwater level (GWL) data collected by Central Groundwater Board (CGWB) at 4307 monitoring wells scattered across the states of NW India. We prepared year-wise GWL change 27

point maps for pre-monsoon (May) - post monsoon (November) periods. Finally, each GWL change point map was interpolated by Kriging and Inverse Distance Weighted (IDW) algorithms using 90% data points. The remaining 10% data points (those are spread over the study area) were used for calculation of RMSE between observed and interpolated values. Two interpolated GWL change maps for each year were then integrated by RMSE based weighted average algorithm. Insitu groundwater storage change during the monsoon was estimated using groundwater level fluctuation method as follows: Change in groundwater storage per unit area = ∆S = h * Sy ….. (1)

(CGWB, 2014)

Where, h represents groundwater level change during the monsoon period and Sy represents specific yield of aquifer materials. Lithologically, the NW India is primarily made of unconsolidated and semi-consolidated sediments followed by sedimentary rocks, extrusive igneous rocks (basalts), intrusive igneous rocks, metamorphic rocks and their combinations (granites, gneisses, phyllites, schists and quartzites). We used a standard specific yield value of 0.12 applicable for unconsolidated and semi-consolidated materials of Indo-Gangetic alluvium (CGWB, 2014). Porosity acts as the governing factor for specific yield and it ranges 30-50% for unconsolidated sediments, 10-30% for sedimentary rocks and basalts, and 5-20% for fractured partially-weathered crystalline rocks (granites and granite gneisses) (Zhang, 2017). We used a normalization factor of 1.0 for unconsolidated and semiconsolidated sediments (specific yield: 0.12), 0.5 for sedimentary rocks and basalts (specific yield: 0.06), and 0.25 for fractured partially-altered igneous rocks such as granites, gneisses and igneous metamorphic rocks combinations (granites, gneisses, phyllites, schists and quartzites) (specific yield: 0.03) to obtain GWL fluctuation-based average monsoonal GWS change (Figure 11a). It is observed that GWL fluctuation based monsoonal GWS change is very low in Luni-Ghaggar 28

drainage basin, north-western parts of Yamuna, Chambal and Upper Ganga drainage basins, northern Kutch of Gujarat state, rugged hard rock terrain of eastern Rajasthan, alluvial plains in Punjab, Haryana, NCR Delhi and western U.P. Overall, the relative spatial distribution of in-situ GWL level fluctuation based GWS change is found to agree well with the GRACE satellite based monsoonal GWS change. We used Pearson product-moment correlation coefficient to measure the nature and strength of correlation between two data layers. Grid-based Pearson correlation coefficient values (generated on 1 degree x 1 degree spatial grids followed by grid-centre point interpolation by IDW algorithm) show good agreement between two results (Figure 11b). We found a positive correlation between GRACE GWS change and GWL fluctuation based GWS change in 11,43,050 km2 out of the total 11,70,375 km2 area which amounts to 97.7% of the study area. Besides, the correlation is moderate to high in 77.4 % area and it is low in 20.3 % area. We found a negative correlation in 27,325 km2 area covering isolated patches in Uttarakhand and U.P. states in the northern part of the study area and in Rajasthan, Gujarat, M.P. and Maharashtra states in western and southern parts of the study area. This is perhaps due to high surface runoff and base flow and large differences in the spatial resolutions between GRACE and GWL fluctuation data in hydrogeologically diverse area with contrasting GWS capacity.

29

30

3.5 Comparison of potential groundwater recharge and GRACE monsoonal GWS change We made a comparison between potential groundwater recharge essentially during the monsoon and average monsoonal GRACE monthly GWS change during the observation period. We made two potential groundwater recharge classes such as high and low potential groundwater recharge and three monsoonal GWS change classes such as high and low increase in GWS (positive GWS change) and decrease in GWS (negative GWS change) classes based on spatial frequency distribution. In this comparison, the following combinations were made:

(1) high potential

groundwater recharge and high increase in monsoonal GWS, (2) low potential groundwater recharge but high increase in monsoonal GWS, (3) high potential groundwater recharge but low increase in monsoonal GWS, (4) low potential groundwater recharge and low increase in monsoonal GWS, (5) high potential groundwater recharge but decrease in monsoonal GWS (negative GWS change in monsoon), (6) low potential groundwater recharge and decrease in monsoonal GWS (negative GWS change in monsoon) (Figure 11). The spatial coverage of the six classes are 1,82,606 km2, 3,37,236 km2, 1,04,179 km2, 4,75,244 km2, 10,534 km2 and 60,868 km2 respectively. Both the potential classes with high increase in monsoonal GWS may be considered safe and prospective in terms of sustained availability of groundwater. On the other hand, the classes with high and low potential groundwater recharge but with low increase in monsoonal GWS are poorly prospective and vulnerable to groundwater depletion. Moreover, the classes with decreasing monsoonal GWS (negative GWS change in monsoon) are at severe risk in terms of groundwater sustainability and are vulnerable to excessive groundwater depletion. It is observed that Luni-Ghaggar drainage basin covering Rajasthan, Haryana, Punjab and NW part of Gujarat states is worst affected in terms of low potential groundwater recharge condition and low increase or decrease in monsoonal GWS. The major part of Yamuna drainage basin covering Delhi, 31

part of Haryana, Rajasthan and U.P. states is badly affected in terms of low potential groundwater recharge condition and low increase in monsoonal GWS. Similarly, the western part of Narmada drainage basin covering a large part of Gujarat state is characterized by low potential groundwater recharge condition with low increase in monsoonal GWS. Besides, the western part of Upper Ganga drainage basin and northern part of Yamuna and Luni-Ghaggar drainage basins covering western U.P., part of Haryana, Punjab, Uttarakhand and H.P. are characterized by high potential high or low groundwater recharge condition but low increase in monsoonal GWS. Such areas are vulnerable to groundwater depletion and therefore at risk for groundwater sustainability in in the near future. The remaining areas such as eastern part of Upper Ganga drainage basin, southern part of Yamuna and Chambal drainage basins, eastern and central part of Narmada drainage basin covering eastern U.P., M.P., eastern part of Gujarat and northern part of Maharashtra states are characterized by either of potential groundwater recharge condition but high increase in monsoonal GWS. Such areas may be considered safe and prospective in terms of sustained availability of groundwater.

4. Conclusions Terrain suitability for groundwater recharge and potential groundwater recharge in 5 major drainage basins of NW India such as Upper Ganga, Yamuna, Chambal, Luni-Ghaggar and Narmada drainage basins were assessed based on hydrogeology and terrain parameters and surface water column available for underground infiltration as a function of monsoonal precipitation, hydrological soil group, land use land cover and topographic slope. The potential groundwater recharge is high in the major part of Upper Ganga drainage basin and in parts of Yamuna and Narmada drainage basins. It is found to be low to very low in Luni-Ghaggar drainage basin including desertic region of Rajasthan and Kutch region of Gujarat states. In the remaining areas including Chambal, Yamuna 32

and Narmada drainage basins, it ranges from low to moderate. We found a sharp disagreement between terrain suitability for groundwater recharge and potential groundwater recharge in LuniGhaggar drainage basin essentially due to very low surface water column available for underground infiltration. GRACE gravity observation provides high precision quantitative estimate of TWS change on monthly interval. In post-processing of GRACE data, a number of uncertainties and leakage errors 33

are involved. The primary sources of uncertainty in the analysis include filtering, leakage correction and glacial isostatic adjustment in addition to the uncertainties associated with geocenter motion and the Earth oblateness. The uncertainty in the global water budget was found to be close to ±0.27 mm/year during 2005-2015 at 90% confidence level (Blazquez et al., 2018). In the present study, GRACE TWS change data could be successfully used to assess monsoonal groundwater recharge in terms of satellite-based GWS change. The GRACE gravity anomaly based average monsoonal GWS change was found to be low in Luni-Ghaggar drainage basin, northern part of Yamuna drainage basin, north-western parts of Upper Ganga, Chambal and Narmada drainage basins. It is moderate to high in the remaining areas of Upper Ganga, Yamuna, Chambal and most of Narmada drainage basins. Moreover, the north-western and western parts of the study area covering LuniGhaggar, Yamuna, Chambal and Narmada drainage basins are dominantly characterized by low potential groundwater recharge condition with low increase or decrease in monsoonal GWS which suggests that the given area is vulnerable to excessive groundwater depletion and is at severe risk towards sustained groundwater. Besides, the western part of Upper Ganga drainage basin and northern part of Yamuna and Luni-Ghaggar drainage basins are characterized by high potential groundwater recharge condition but significantly low increase in monsoonal GWS which suggests that such areas are also vulnerable to significant groundwater depletion over the years. The remaining areas such as eastern part of Upper Ganga drainage basin, southern part of Yamuna and Chambal drainage basins, eastern and central parts of Narmada drainage basin may be considered safe in terms of groundwater sustainability but differ in groundwater recharge potential.

34

Acknowledgement The study was funded by Indian Space Research Organization (ISRO) under the Earth Observation Application Mission (EOAM) project scheme. GRACE monthly mass grid data available from http://grace.jpl.nasa.gov were used in the study.

References Agrawal, S., Joshi, P.K., Shukla, Y., Roy, P.S., 2003. Spot Vegetation multi-temporal data for classifying vegetation in south central Asia. Current Science 84, 1440–1448. Blazquez, A., Meyssignac, B., Lemoine, J.M., Berthier, E., Ribes, A., Cazenave, A., 2018. Exploring the uncertainty in GRACE estimates of the mass redistributions at the Earth surface: Implications for the global water and sea level budgets. Geophysical Journal International 215 (1), 415-430. Bonham-Carter, G.F., 1994. Geographic Information Systems for Geoscientists: Modelling with GIS, Pergamon, Ottawa, Ontario, Canada, 416p. CGWB, 2014. Dynamic ground water resources of India (as on 31st March 2011). Central Ground Water

Board,

Faridabad,

India,

299p.

http://www.cgwb.gov.in/Documents/Dynamic-GW-

Resources-2011.pdf (accessed 10 March 2014). CGWB, 2019. National compilation on dynamic ground water resources of India, 2017. Central Ground Water Board, Faridabad, India, 298p. Chambers, D.P., Wahr, J., Tamisiea, M.E., Nerem, R.S., 2010. Ocean mass from GRACE and glacial isostatic adjustment. Journal of Geophysical Research (Solid Earth) 115 (B11), B11415.

35

Chen, H., Zhang, W., Nie, N., Guo, Y., 2019. Long-term groundwater storage variations estimated in the Songhua River Basin by using GRACE products, land surface models, and in-situ observations. Science of the Total Environment 649, 372-387. Chen, J., Li, J., Zhang, Z., Ni, S., 2014. Long-term groundwater variations in Northwest India from satellite gravity measurements. Global Planetary Change 116, 130–138. Coyle, G., 2004. The Analytic Hierarchy Process (AHP), Practical Strategy, Open Access Material (ID: #S8RP9920). Pearson Education Limited, Harlow, U.K. Dasgupta, S., Das, I.C., Subramanian, S.K., Dadhwal, V.K., 2014. Space-based gravity data analysis for groundwater storage estimation in the Gangetic plain, India. Current Science 107, 832– 844. De Vries, J.J., 1984. Holocene depletion and active recharge of the Kalahari groundwaters - A review and an indicative model. Journal of Hydrology 70, 221–232. Famiglietti, J.S., Lo, M., Ho, S.L., Bethune, J., Anderson, K.J., Syed, T.H., Swenson, S.C., de Linage, C.R., Rodell, M., 2011. Satellites measure recent rates of groundwater depletion in California’s Central Valley, Geophysical Research Letters, 38, L03403. Feng, W., Shum, C.K., Zhong, M., Pan, Y., 2018. Groundwater Storage Changes in China from Satellite Gravity: An Overview. Remote Sensing 10 (5), 674. Garduño, H., Romani, S., Sengupta, B., Tuinhof, A., Davis, R., 2011. India groundwater governance case study. http://water.worldbank.org/sites/water.worldbank.org/files/GWGovernance India.pdf (accessed 13 April 2014). Gong, H., Pan, Y., Zheng, L., Li, X., Zhu, L., Zhang, C., Huang, Z., Li, Z., Wang, Z., Wang, H., Zhou, C., 2018. Long-term groundwater storage changes and land subsidence development in the North China Plain (1971–2015). Hydrogeology Journal 26 (5), 1417-1427. 36

GSI, 1998. Geological Map of India, 7th ed. (Scale 1: 2,000,000). Geological Survey of India, Kolkata, India. GSI, 2002. Geomorphological Map of India, 1st ed. (Scale 1: 2,000,000). Geological Survey of India, Kolkata, India. Hinrichsen, D., Tacio, H., 2002. The coming freshwater crisis is already here, in finding the source: The linkages between population and water. Environmental Change and Security Program, Woodrow

Wilson

Center.

Washington,

D.C.

http://www.wilsoncenter.org/topics/pubs/

popwawa2.pdf (accessed 10 April, 2014). IMD, 2014. Daily gridded rainfall data of India (0.25o x 0.25o), National Climate Centre, India Meteorological Department, India. ISFR, 2011. India State of Forest Report. Forest Survey of India, Ministry of Environment and Forests, Dehradun, India, 270p. Jenson, S.K., Domingue, J.O., 1988. Extracting topographic structure from digital elevation data for geographic information system analysis. Photogrammetric Engineering and Remote Sensing 54 (11), 1593–1600. Kent, K.M., 1973. A method for estimating volume and rate of runoff in small watersheds. U.S. Soil Conservation Service, SCS-TP-149, 19p. Khurana, I., Mahapatra, R., 2008. Water and sanitation perspective. Wateraid India, Macrographics Pvt. Ltd., New Delhi. http://www.wateraid.org/~/media/Publications/drought-drinking-water-crisisbundelkand-india.pdf (accessed 14 March 2014). Landerer, F.W., Swenson, S.C., 2012. Accuracy of scaled GRACE terrestrial water storage estimates. Water Resources Research 48 (4), W04531.

37

Li, X., Gao, X., Chang, Y., Mu, D., Liu, H., Sun, Z., Guo, J., 2018. Water storage variations and their relation to climate factors over Central Asia and surrounding areas over 30 years. Water Science and Technology: Water Supply 18 (5), 1564-1580. Marjanovic, M., Caha, J., 2011. Fuzzy approach to landslide susceptibility zonation, V. Snasel, J. Pokorny, K. Richta (Eds.): Dateso 2011, 181–195. Mu, Q., Heinsch, F.A., Zhao, M., Running, S.W., 2007. Development of a global evapotranspiration algorithm based on MODIS and global meteorology data. Remote Sensing of Environment 111, 519–536. Mu, Q., Zhao, M., Kimball, J.S., McDowell, N.G., Running, S.W., 2013. A Remotely Sensed Global Terrestrial Drought Severity Index. Bulletin of the American Meteorological Society, 94 (1), 83-98. Mu, Q., Zhao, M., Running, S.W., 2011. Improvements to a MODIS global terrestrial evapotranspiration algorithm. Remote Sensing of Environment 115, 1781–1800. National Water Policy (NWP), 2012. Ministry of Water Resources, Government of India. http://mowr.gov.in/writereaddata/linkimages/DraftNWP2012_English9353289094.pdf (accessed 12 April, 2014). NBSS&LUP, 2002. Soils of India (Soil Resource Maps), 1st ed. (Scale 1: 1,000,000). NBSS Publication 94, NBSS&LUP, Nagpur, India. Peltier, W.R., 2004. Global Glacial Isostasy and the Surface of the Ice-Age Earth: The ICE5G(VM2) model and GRACE. Annual Review of Earth and Planetary Sciences 32, 111-149. Raharjo, H., Endah, D., 2006. Evaluating relationship of consistency ratio and number of alternatives on rank reversal in the AHP. Quality Engineering 18, 39-46.

38

Richey, A. S., 2014. Stress and Resilience in the World’s Largest Aquifer Systems: A GRACEBased Methodology, PhD Dissertation, University of California, Irvine. Rodell, M., Beaudoing, H.K., 2007. GLDAS Noah Land Surface Model L4 Monthly 0.25 x 0.25 degree V001, NASA/GSFC/HSL (accessed 14 May 2017). Rodell, M., Chen, J., Kato, H., Famiglietti, J.S., Nigro, J., Wilson, C.R., 2007. Estimating groundwater storage changes in the Mississippi River basin (USA) using GRACE. Hydrogeology Journal 15 (1), 159-166. Rodell, M., Velicogna, I., Famiglietti, J.S., 2009. Satellite-based estimates of groundwater depletion in India. Nature 460, 999–1002. Saaty, T.L., 1980. The Analytic Hierarchy Process, 1st ed. McGraw Hill, New York. Saaty, T.L., 1994. How to make a decision: the analytic hierarchy process. Interfaces 24(6), 19–43. Saaty, T.L., 2008. Decision making with the analytic hierarchy process. International Journal of Services Sciences 1(1), 83-98. Saaty, T.L., Vargas, L., 2012. Models, methods, concepts and applications of the analytic hierarchy process, International series in operations research and management science (Series vol. 175), 2nd ed. Springer, New York. Seo, K.W., Wilson, C.R., Famiglietti, J.S., Chen, J.L., Rodell, M., 2006. Terrestrial water mass load changes from Gravity Recovery and Climate Experiment (GRACE). Water Resources Research 42(5), W05417, 15p. Sinha, B.P.C., Sharma, S.K., 1988. Natural groundwater recharge estimation methodologies in India, in: Simmers, I. (Ed.), Estimation of Natural Groundwater Recharge. D. Reidel Publishing Co., Dordrecht, pp. 301-311.

39

Steward, D.R., 2014. Groundwater depletion and irrigated agriculture in the semi-arid grasslands of the central USA, GWF discussion paper 1404, Global Water Forum, Canberra, Australia. http://www.globalwaterforum.org/wp-content/uploads/2012/04/GWF-1404-Groundwaterdepletionand-irrigated-agriculture-in-the-semi-arid-grasslands-of-the-central-USA.pdf (accessed 16 May 2014). Suhag, R., 2016. Overview of ground water in India, PRS Legislative Research, 11p. Swenson, S., Wahr, J., 2006. Post‐processing removal of correlated errors in GRACE data. Geophysical Research Letters (Hydrology and Land Surface Studies) 33 (8), L08402. Tiwari, V.M., Wahr, J., Swenson, S., 2009. Dwindling groundwater resources in northern India, from satellite gravity observations, Geophysical Research Letters 36(18), L18401, 5p. UNICEF, FAO, SaciWATERs, 2013. Water in India: Situation and prospects. UNICEF, FAO, New Delhi, India. http://www.unicef.org/india/ Final_Report.pdf (accessed 14 April 2014). USGS, 2015. Shuttle Radar Topography Mission (SRTM) 1 arc-second global elevation data. https://lta.cr.usgs.gov/SRTM1Arc (accessed 25 June 2015). Voss, K.A., Famiglietti, J.S., Lo, M., de Linage, C., Rodell, M., Swenson, S.C., 2013. Groundwater depletion in the Middle East from GRACE with implications for transboundary water management in the Tigris‐Euphrates‐Western Iran region. Water Resources Research 49 (2), 904-914. Wahr, A.G.J., Zhong, S., 2013. Computations of the viscoelastic response of a 3-D compressible Earth to surface loading: an application to Glacial Isostatic Adjustment in Antarctica and Canada. Geophysical Journal International 192 (2), 557–572. Wahr, J., Molenaar, M., Bryan, F., 1998. Time variability of the Earth’s gravity field: Hydrological and oceanic effects and their possible detection using GRACE. Journal of Geophysical Research 103(B12), 30205-30229. 40

Wahr, J., Swenson, S., Zlotnicki, V., Velicogna, I., 2004. Time-variable gravity from GRACE: First results. Geophysical Research Letters 31(11), L11501, 4p. Wilson, J.P., Gallant, J.C., 2000. Terrain Analysis: Principles and Applications. John Wiley & Sons, New York, USA, 479p. Wyrwoll, P., 2012. India’s groundwater crisis, GWF discussion paper 1228, Global Water Forum, Canberra,

Australia.

http://www.globalwaterforum.org/2012/07/30/indias-groundwater-crisis

(accessed 15 April 2014). Zhang, P., 2017. EAS 44600 Groundwater Hydrology Lecture 4: Porosity and Permeability. https://www.coursehero.com/file/6061543/EAS446lec4/ (accessed 16 August 2017).

41

Appendix-1

Derivation of Topographic Wetness Index (TWI) from Digital Elevation Model (DEM) From SRTM 30 m DEM, flow direction map was prepared by eight direction (D8) flow model (Jenson and Domingue, 1988) using the ‘Flow Direction’ tool of ArcGIS software. In this algorithm, the direction of flow which could be towards any of the eight adjacent pixels is determined by the direction of steepest descent or maximum drop in the DEM with respect to the central pixel as follows: Maximum Drop (%) = (Change in Height value/Distance) * 100. When the direction of steepest descent is found, the output cell (central pixel) is coded with the direction value. When two pixels flow towards each other, they are termed as sink pixels. So, ‘Fill Sinks’ operation needs to be implemented to remove local depressions by assigning lowest pixel height from the adjacent pixels before executing flow direction step. Subsequently, flow accumulation map was prepared from flow direction map using ‘Flow Accumulation’ tool of ArcGIS. In flow accumulation map, the number of pixels flow into each pixel is determined and all the pixels are then coded with the corresponding numbers. From the flow accumulation map, specific catchment area map was prepared using flow accumulation area for each pixel per unit length as follows: Specific Catchment area = (Flow Accumulation area)/(Size of the pixel). Finally, Topographic Wetness Index (TWI) was generated using specific catchment area and topographic slope of the study area as follows: TWI = ln(Specific catchment area)/slope.

42

Highlights • Reports of severe ground water depletion in different parts of North-Western India • Potential groundwater recharge based on terrain suitability and surface water infiltration • Comparison of potential groundwater recharge vs monsoonal groundwater storage change • Evaluation of groundwater sustainability and categorization into different classes