An Integrated Hydrogeological and Remote Sensing Modeling Approach to Evaluate the Climate Change and Over-Irrigation Impact on Groundwater Depletion in North Jordan

Groundwater is an essential source of freshwater for the severely water-deprived country of Jordan. Thus, understanding abstracted amounts of the ground water is essential for sustaining this resource. More than 50% of the annual abstracted groundwater is used for irrigation, indicating the need to correctly estimate this usage and minimize the uncertainties. Thus, this study aims to assess actual groundwater abstraction in North Jordan by implementing a new geographic information system–remote sensing multi-method approach. The results are promising and show that the announced GW abstraction volumes provided by the Water Authority of Jordan and the Ministry of Water and Irrigation are not particularly accurate, especially for areas where extensive agricultural practices are present. This fact may explain the discrepancies between official reports on declining groundwater levels and the results of several studies that assessed North Jordan groundwater abstraction. This study showed that irrigation is a crucial factor that affects this discrepancy and that the crop water requirements are generally higher than what is reported by the governmental sector. The results showed a 31 MCM discrepancy between the official recorded groundwater abstraction volumes and the actual (crop-based) water consumption. Increasing the irrigation efficiency by 20% will aid in sustaining the groundwater in the study area, maintaining the current groundwater wells and save the farmers money


Introduction
Groundwater (GW) is an essential source of water, especially in arid and semi-arid regions (Subyani and Al-Ahmadi, 2009;Awadh et al., 2020).However, in Jordan, the second most water-poor country worldwide (Ibrahim and Al-Mashakbeh, 2016), GW represents the primary source of freshwater (Alraggad et al., 2016).GW accounts for more than 65% of Jordan's total 900 MCM/year water supply (MWI, 2013).Nonetheless, over the past few decades, Jordan's water resources in general have been overexploited to cover the needs of its rapid population growth as a result of the normal population growth rates and the refugee crises (Al-Karablieh and Salman, 2016).In addition, increased climate-induced evaporation and reduced rainfall rates have altered the GW recharge rates in Jordan as a whole (Alraggad et al., 2016).In Jordan, there are twelve main GW basins, of which two are shared with Syria in the north and Saudi Arabia in the east (Ramirez et al., 2008;Al-Harahshah et al., 2020).In Northern Jordan, particularly, the Mafraq/Azraq GW basins are primarily exploited for irrigation water to support the agricultural activities in the area (Ramirez et al., 2008).The GW basins in Jordan are generally overexploited, with the extraction rates exceeding the recharge rates.This overexploitation is mainly caused by illegal agricultural GW wells constructed without governmental permits, which makes it difficult to calculate the exact amount of water being extracted.Thus, the data provided by the Jordan Ministry of Water and Irrigation (MWI), primarily obtained from legal (permitted) GW wells, may not reflect the actual GW situation in the area.For example, in 2005, water consumption in Jordan reached 1,042 MCM/year compared with the MWI-reported annual average of 990 MCM/year.Out of the 1,042 MCM of water consumed, 679 MCM were used for agricultural activities (Ramirez et al., 2008).Accordingly, this study aims to calculate the actual GW abstraction volume in Northern Jordan and to reduce the uncertainties caused by illegal GW well pumping and the underestimation caused by legal GW well meter malfunctions through the application of a new geographic information system-remote sensing (GIS-RS) multi-method approach.
The study area (Fig. 1) is located in Northern Jordan and spreads over two main governorates, Mafraq and Zarqa.The topographic variability in the study area ranges from highlands reaching up to ca. 1,000 m asl near the Jordan-Syria borders, at close proximity to Jebel Al-Druz, to desert lowlands of ca.350 m asl.The annual average precipitation over the study area ranges from <50 mm/year in the eastern low-relief deserts to >300 mm/year over the western highlands (MWI, 2018).The study area comprises two primary GW aquifers: the upper Basalt (BSCPX) aquifer and the lower Amman Wadi Es-Sir (B2/A7) aquifer (Fig. 1b).The BSCPX aquifer extends from Jebel Al-Druz in Syria into Eastern Jordan and northern Saudi Arabia (Al-Harahshah et al., 2020;Odat, 2015).The B2/A7 aquifer primarily outcrops over the highlands, where higher precipitation is recorded, while the BSCPX aquifer covers the central and eastern parts of the study area.This study area was chosen because of the presence of extensive agricultural practices where the GW abstraction rates are usually underestimated, making it an important site to investigate the uncertainties of GW abstraction measurements.

Remote Sensing for Determining Crop Evapotranspiration Demand
Several irrigated crops such as grapes, olives, pomegranates, and different vegetables are being planted in the study area.To estimate the crop water requirements (CWRs) for each crop type, a Landsat 8 Operational Land Imager (OLI), terra satellite image was acquired from the United States Geological Survey (USGS), the cloud-free image was selected on July 11, 2018, to represent the study area.Different preprocessing techniques were applied to prepare the Landsat 8 image prior to mapping the different crop types (Fig. 2).Similar ground resolutions (30 m) of visible and near-infrared (NIR) bands of Landsat 8 were combined by implementing layer stacking using the PCI Geomatica (2017) RS software.Then the image was geometrically corrected by collecting ground control points (GCPs) at crossroads using a global positioning system (GPS) device, with an accuracy of ±3 m.A third-order transformation model was applied to correct the image; then the image was clipped to the borders of the study area.The satellite images available for 2018 and 2019 were directly downloaded from Earth Explorer.Later on, a correction for cloud effect was applied using the F-mask algorithm afterward (Zhu and Woodcock, 2012).The ArcMap 10.3 software was used to reclassify cloud mask layers and to check the quality of the satellite images.Normalized difference vegetation index (NDVI) values were derived from the satellite images.The corrected bands were then used to derive the NDVI values based on the following equation: where B5 is the NIR band and B4 is the red band, both representing spectral reflectance at the ground level.Field observations were done on selected spots and merged with GW field office data and observations, which were all transferred into a GIS shapefile in the form of an initial crop map.The map was overlaid with the NDVI values map to assign NDVI values for each crop type.The NDVI values were then extracted to show olive trees, fruit trees, vegetables, and alfalfa.As it was difficult to extract NDVI values for different fruit tree types, field data was used to reclassify the fruit trees into grape, stone fruit, pomegranate, and date palm.Accordingly, the main outcomes of the NDVI profiles in the study area were as follows: Alfalfa was found to show the highest NDVI values (up to 0.99).However, these values were found to be highly variable because of the impact of alfalfa harvesting, which takes place once every 1.5 to 3 months (Fig. 2).Olives showed a more stable pattern of NDVI values all year round given their evergreen pattern (Fig. 3).Vegetables showed moderate NDVI values (0.43 to 0.57) during the high season (June and July), while they fell below 0.3 during the early growing stages.Fruit trees in general showed low NDVI values of less than 0.3 during early spring, while they reached up to 0.73 during the summer season (June and August).The NDVI values and profiles were used to determine the growing seasons of different crops and to validate water abstractions for irrigation recorded by the MWI.Also, it helped in the estimation of KC and ETc.

Digitizing of Land Parcels and Map Verification
As the resolution of the satellite images is coarser than tree spacing and parcel borders, highresolution Google Earth Pro (GE) images were used to extract the exact on-ground irrigated area.The NDVI images were used as a background for the GE images, and onscreen digitizing was conducted, as shown in Fig. 4.

Calculating Net Crop Water Requirements
The daily climatic records from weather stations managed by the MWI were used to calculate the net crop water requirements (NCWRs) following the standard method of FAO56 (Allen et al., 1998) to calculate ETc as follows: where ETc is the actual crop evapotranspiration (mm), Kc is the mean monthly crop coefficient computed following the standard method of FAO56, and ETo is the grass reference evapotranspiration.The daily records of MWI weather stations were used as an input to the FAO-ETo calculator, which is a free open-source software by the Food and Agriculture Organization of the United Nations or FAO.The resulting daily ETo values were processed to calculate the monthly ETo, while the Kc values were derived from the Kc standard tables given by FAO56 (Allen et al., 1998).Other data-such as planting time, growing season length, and harvesting dates-was obtained from the field data of the MWI.

Estimation of Groundwater Abstraction for Irrigated Crops
Defining GW abstraction for irrigation purposes is usually done by the Water Authority of Jordan (WAJ) through the WAJ regional offices (ROs).Over the past thirty years, GW abstraction has been managed through direct water meter readings.In the case of no meter installed or some technical issues in the existing meter, the WAJ would make an estimation of the abstraction based on their experience in the area and the present cropping patterns.Nonetheless, many wells are not monitored as a result of the presence of hundreds of illegal wells in the study area.Based on the previous work of the MWI and upon observing the extreme decline in GW levels in the study area, we can conclude that the actual abstraction is still underestimated.The CWR represents an efficient tool to estimate agricultural abstraction.However, in the study area, most of the farmers over-irrigate their crops, which means that even the abstraction estimated by the CWR is still less than the actual abstraction, but it is certainly more accurate than the WAJ records.Accordingly, ETc values and calculations from FAO56 were used by attributing the CWR for each crop area generated in the previous sections.MODFLOW-Model-Muse modelling was implemented to construct the GW numerical model for this study (McDonald and Harbaugh, 1988;Harbaugh and McDonald, 1996;Winston, 2009), at 30 × 30 m horizontal grid block dimensions and 50 m to 700 m vertical scale layer discretization.The time unit used was days (d), and meters for the length.Thus, the head measurements are given in meters above mean sea level (amsl), the hydraulic conductivity values are in meters per day (m/d), and the flow values are given in cubic meters per day (m3/d).The stratigraphic settings of the area were derived from the BGR (2018) crosssections and, more importantly, from the isopach maps of all the geological formations included in the BGR reports (Isopach Maps 1.1 to 1.16, Volume 3, Part 2; Appendix 1, Margane and Hobler, 1994;2018).Fig. 5 shows the subsurface settings of the study area with different aquifers and aquitards (for the cross-section locations, see Fig. 7).

Conceptual Model of the Groundwater Flow System
The GW system in the study area is modeled as a three-dimensional multi-layer system.The model layers are mapped to different hydrogeological formations using equivalent parameter values.Following this approach, it is possible to separate the choice of number of model layers from the number of hydrogeological formations.Several phenomena have to be simulated in this system, including areal recharge, evaporation and evapotranspiration, extraction by wells, and relation to the surrounding aquifer systems.
Various pump tests have been performed in the study area.Available pumping tests have been analyzed with the AquiferTest software v.3.5.The drawdown in the pumping wells is compared with the analytical drawdown for relevant hydrogeological conditions: confined (Theis equation) and unconfined (Neuman equation) aquifers in which the aquifer parameters were adjusted to give the best fit.By visually comparing the fitted drawdown curves, the proper hydrogeological condition was chosen, and the hydraulic conductivity based on that curve was selected.
The storativity values here are considered unreliable as the drawdown has only been measured in the pumping wells and not in a monitoring well at some distance to the pumping well (see also the manual of Aquifer Test).Unfortunately, Aquifer Test does not work properly for pumping tests that increase the pumping rates stepwise.These pumping tests have been fitted using the drawdowns of the first period only before the pumping rates were increased.For a number of wells, no reasonable fit could be obtained.The drawdown for these pumping wells does not follow any of the analytical functions, and no hydraulic conductivity could be obtained using the pumping test data.A summary of the pumping test results is given in Table 2, which shows that the pumping tests indicate a wide range of hydraulic conductivity values from 0.02 m/day to 358 m/day.These values are initially used to set reasonable minimum and maximum values of the hydraulic conductivity in the model calibration.However, as will be shown later, it was still required by the inverse solution algorithm to increase the maximum and decrease the minimum allowed values to allow better calibration.Acknowledging the fact that a large variability is observed in the estimated hydraulic conductivity from the pumping test analyses, increasing the upper and decreasing the lower bounds, as will be discussed later, is reasonable.
The boundary conditions were set according to the geological model constructed in the previous sections and the available time series water level measurements near the boundaries.The model boundaries in Mafraq were found to be more complicated given the complex aquifer settings and flow regime variability within the BSPCX and B2/A7 aquifers.The transboundary flux coming from Syria was a big challenge given the lack of data.However, the MWI (2013) estimated a total flux of 27 MCM from Jebel Al-Druz toward the Mafraq area.This value was assigned for the northern model boundary as specified flux.A GW flow map was prepared by the BGR (2013) for the major aquifers in the study area.It shows that the GW flows from the Jebel Al-Druz basaltic area radially toward the surrounding areas.This map was updated for the year 2017 to show similar trends but with a 5-15 m lowered head (Fig. 6).Thus, the eastern and western boundaries were assigned as "no flow" boundaries (Fig. 7), where the eastern boundary represents the divide to the Azraq basin and the western boundary represents the limit of saturation of the aquifer (BGR, 1995).
Within this model, the hydrological processes are simulated based on the water balance parameters independency.Firstly, the parameters are individually estimated and then used as a combined balance to calculate the recharge (Fig. 8).In addition, the precipitation isohyets were calculated based on Thiessen polygons method (Chow, 1988) until adequate statistical certainty is achieved.Daily climate data for the period  was acquired from the MWI meteorological stations.Longterm precipitation data interpolation was conducted based on inverse distance weighing (IDW) method, alongside corrections for elevation.At the central parts of the study area, the estimated mean precipitation was 112 mm/year.Moreover, following the Penman-Monteith ( 1948) equation ( Jensenand Allen, 1990), evapotranspiration was calculated using the J2000 ETP module, which proved valid for calculating evapotranspiration for green cover (FAO, 1989).
Accordingly, the results obtained from the J2000 model showed that the variability in the amounts of recharge over the study area follows the precipitation trends, where higher recharge is recorded in the northwestern parts of the study area (Fig. 9).

Model Calibration
The model was initiated under steady state conditions for the year 1980, conditions before extensive GW extraction, followed by a transitory historical simulated period of 40 years, wherein historical GW alterations were applied.Calibration was done using both trial-and-error and automatic calibration (inverse problem) methods for the first aquifer (the B2/A7 aquifer complex) in the first model layer.The available time series have shown a high degree of noise.As a result, data was necessary.Thus, the data was averaged per year, and any single year with less than four measurements or with high variance was ignored.Following to this, zonation parameterization was used to map the spatial variability of the governing hydraulic parameters.Zones indicating uniform hydraulic parameters were chosen primarily based on surface geology before being subdivided according to the density of the pumping wells, and then available head measurements were used to discriminate among areas of different hydraulic parameters.Finally, a subset of the filtered head time series and head data (MWI) was used in the calibration given higher weight compared with the latter (Fig. 10).

Climate Statistical Downscaling Model
Historical climate records have recorded a 1 mm increase in evaporation and a 3°C increase in temperature for Jordan since 1965 (MWI 2020).Thus, to assess the climate variability and its consequent impact on the basin's local area, statistical downscaling model (SDSM) v. 4.2, a hybrid model that uses transfer function methods to generate probabilistic weather (Wilby et al., 2002), was implemented to evaluate the site-specific atmospheric conditions variability.This model allows for calculating downscaled daily local weather based on multiple linear regressions of predictor-predictant associations (Hashmi et al., 2009).We calibrated and tested the model for twelve rainfall stations and five meteorological stations.Daily mean temperatures and rainfall amounts for the period 1970-2080 were downscaled based on data from the HadCM3 (UK Met Office) climate model.The model runs this estimation in five processes: (1) predictor variables screening, (2) model calibration, (3) observed data (PAR) production, (4) scenario generation, and (5) statistical analyses.We then used MWI (2017) and the Jordan Meteorological Department (JMD, 2017) data from stations collected since their starting dates.Furthermore, missing values, double-mass curves and outliers data processing were conducted prior the final time series analysis which included trend and climate analyses (Fig. 11).The climate change scenarios were adopted as an input to the GW recharge layer in the GW flow model.

Crop Evapotranspiration Demand Using Remote Sensing
The spatial and geometric calculation showed that the total irrigated area is 106,000 dunums, as given in Table 4.The analyses showed that stone fruits are the major crop in the area, with 40,000 dunums representing 37% of the total crops.The second most dominant crop is olives, with 39,000 dunums, followed by vegetables, with around 14,000 dunums.This was also validated in the field (Fig. 12).Monthly evapotranspiration for the main crop types in the area was found to be very high compared with the region given the low relative humidity enhancing more evapotranspiration, as listed in Table 5.Based on the irrigated areas, crop type, and value of ETc, GW pumping was estimated (Table 6).Field validation showed that the irrigation efficiency is no more than 70% (Fig. 13).This percentage was used to correct the pumping rates and account for this additional amount of abstraction.Considering the 70% irrigation efficiency, the MWI records clearly and extremely underestimate the actual GW abstraction volume in the study area.A total of ca. 31 MCM is being left out of the GW abstraction calculations of the MWI, which may account for the lower GW drawdown records usually published by the MWI.Increasing the efficiency of the irrigation practices would significantly decrease the GW abstraction rates and sustain the GW for longer periods in the study area.The simulated effects of climate variability, as measured by the SDSM, show increased evaporation starting in 2008 and ending in 2010; this increase is followed by a decreasing trend toward the year 2018.There is a slight variation in the mean maximum and minimum temperatures, as shown in Fig. 32.These records are similar to the long-term rainfall and temperatures recorded by global climate models (Fig. 14).These decreasing precipitation trends, accompanied by increasing temperatures and consequent evaporation, imply the true effect of climate change on such a vulnerable GW system.The decrease in precipitation will generate less surface flow, while the increase in temperature will allow more evaporation, which will result in less GW recharge.The GW level simulations for the years 2019, 2021, and 2023 showed that the average drawdown is expected to be ca.9.12 m by 2023.The maximum drawdown of 14.21 is expected to occur along Baghdad Road (eastern parts of the study area), where the most extensive agriculture is taking place.Fig. 15 shows the simulated GW heads in the Mafraq area for the years 2019, 2021, and 2023.The GW drawdown can also be illustrated by a (E-W) cross-section showing the decrease in water levels (Fig. 16).

Discussion
GW is an essential source of fresh water for Jordan; thus, understanding the abstracted amounts of GW is essential for sustaining this resource.More than 50% of the annual abstracted GW is used for irrigation, indicating the need to correctly estimate this usage and minimize uncertainties.As mentioned earlier, the presence of illegal GW wells, associated with frequent meter malfunctions and rough estimations of GW abstraction done by the WAJ staff, are all factors affecting the estimation of the actual GW abstraction volumes and the sustainability of the GW aquifers, especially in Northern Jordan.The results of this study showed that the WAJ and MWI estimations of GW abstraction volumes in agricultural areas are underestimated, even when considering 70% irrigation efficiency.The GW level varied from 420 m in the eastern parts to about 520 m in the western parts for the year 2019.However, for the year 2021, there is a clear drawdown of about 10 m spread over the study area and particularly across the eastern parts, where the aquifer base is at about 380 m.This drawdown increases for the year 2023 simulation, considering a scenario of business as usual (BAU), and more severe drops are noted, reaching up to 14 m around areas where GW is being abstracted for both agricultural and domestic uses (Fig. 17).Looking further into the future, the BAU scenario is expected to lead to a severe drop in GW levels and aquifer "total dryness" in the eastern parts of the study area by 2050.Nonetheless, when applying a 20% enhancement of irrigation efficiency (Scenario 1 in Fig. 18), we can see a lower drawdown and better sustainability of this endangered aquifer.This irrigation enhancement will not only be beneficial in maintaining the GW aquifer but also help maintain the current GW wells and prevent the need for deepening the well and changing the pump and associated parts, hence saving the farmers financial expenses.West: Zatari area.

Conclusions
This study aimed to assess the actual GW abstraction volume in North Jordan by implementing a new GIS-RS multi-method approach.The results were promising and showed that the announced GW abstraction volumes provided by the WAJ and the MWI are not particularly accurate, especially for areas where extensive agricultural practices are present.This fact may explain the discrepancies between official reports on the decline of GW levels and the results of some studies.This study showed that irrigation is a crucial factor affecting this discrepancy and that the CWR is generally higher than what is reported by the governmental sector.The results showed a discrepancy of about 31 MCM between the official recorded GW abstraction volumes and the actual (crop-based) water consumption.Increasing the irrigation efficiency by 20% will aid in sustaining the GW in the study area, maintaining the current GW wells, and save the farmers money.

Fig. 1 .
Fig.1.(a) Location of the study area in Jordan (inset) and the digital elevation model of the study area (m asl); (b) The main hydrological units (GW aquifers) in the study area

Fig. 4 .
Fig.4.GE image for a sample farm showing the variation in agricultural activities, tree spacing, and crop coverage (in Mafraq)

Fig. 5 .
Fig.5.Subsurface settings of the study area

Fig. 7 .
Fig.7.GW model boundaries and hydrological units in the Mafraq area.

Fig. 10 .
Fig.10.Computed versus observed head values for head measurements for the study area

Fig. 11 .
Fig.11.(a) Open pan ETP, (b) Max temperature, and (c) Min temperature in the study area

Fig. 13 .
Fig.13.GW pumping based on calculated CWR as mm, assuming 70% irrigation efficiency in the study area

Fig. 15 .
Fig.15.(a) Simulated GW level in the Mafraq area in 2019.(b) Simulated GW level in the Mafraq area in 2021.(c) Simulated GW level in the Mafraq area in 2023

Fig. 16 .
Fig.16.Simulated GW levels in the study area for the years 2019, 2021, and 2023

Table 1 .
Summary of Pumping Test Results (MWI Files)

Table 2 .
Total Area of Irrigated Crops in the Study Area

Table 3 .
Monthly evapotranspiration for main crop types in the study area

Table 4 .
Estimated GW abstraction based on crop types considering 100% and 70% irrigation efficiency along with the actual amount billed by the MWI (MCM)