Advanced Hydro-Geospatial Modeling for Groundwater Potential Zoning in the Baro River Watershed: Geographic Information System and Remote Sensing Approach

1. Introduction

Groundwater is the most important and dynamic natural freshwater resource to sustain ecosystems and to serve as a significant buffer to the variability of surface water due to climate change [1,2,3]. Population growth, industrialization, and episodic droughts have led to higher demands for groundwater in various worldwide areas [4,5]. In sub-Saharan Africa, countries with little surface water infrastructure, groundwater is the major source of residential water supply consumption and irrigation demands [6,7,8]. Furthermore, the non-uniform distribution of such a resource is controlled by complex geomorphic, anthropogenic, and climatic interactions [9,10].

Ethiopia is the “Water Tower of Africa,” with an abundance of surface water but experiencing localized water scarcity [11,12]. Groundwater is found in various settings across the country, such as high-elevation volcanic plateaus and low-lying sedimentary basins [13,14,15]. However, recent studies suggest that to be climate resilient, basin authorities such as the Upper Blue Nile basin, Omo-Gibe basin, and Awash basin focus more on shifting from relying on surface water to groundwater-based projects [16,17,18]. Even though the Baro River watershed is a large-scale watershed, its subsurface potential is mostly unreported and understudied despite these consequences [19].

The Baro River watershed covers an estimated area of 23,000 km² in the southwest Ethiopian highlands, and is of great hydrological importance. It is a large tributary of the Nile system that has a substantial annual rainfall and flows into the alluvial plains of the Gambella region from broken volcanic mountains [1,20]. While the surface water has tremendous potential, the groundwater prospects are not all that clear. The application of conventional hydrogeological methods such as exploratory drilling and geophysical sounding (VES) is challenging, costly, and time-consuming in vast and massive terrains [21,22].

Hydro-geospatial modeling has developed an effective and cost-saving solution to address these challenges [23,24]. To model the groundwater potential zones (GWPZ), the latest version of Geographic Information Systems (GIS) and Remote Sensing (RS) can be integrated [2,3,8,9,17,18,19,23]. The excellent ability of GIS to handle high-resolution information and its powerful raster processing module, especially through open source software such as GRASS [10,25], makes it an effective tool for large-scale research. On the other hand, modules such as r.stream.extract and r.mapcalc can be used to derive precise secondary components needed for infiltration modeling, like drainage density and Topographic Wetness Index (TWI) [10,26].

The differentiation of GWPZ needs to take into account several issues, including the balancing of various conditioning factors [3,10,27,28]. The Analytic Hierarchy Process (AHP) is a subcategory of Multi-Criteria Decision Analysis (MCDA) [10,27,29,30], which is the best method for this. The 9 important parameters that were shown to be relevant in the hydrological aspects were lithology, slope, lineament density, land use/land cover (LULC), curvature, rainfall, geology, soil, topographic wetness index (TWI), and drainage density [2,11,14,24]. Primary and secondary porosity and the rate of runoff vs. infiltration controlled by slope and TWI were perceived using the lineament density and lithology [15,31].

On the terrestrial water balance, the main hydraulic input is rainfall, while the water availability in soil moisture and in deep infiltrations is strongly linked with the actual evapotranspiration [32]. High near surface evaporation loss in arid and semiarid environments [33,34], low recharge potential due to high evapotranspiration, and shallow water table conditions that do not allow for groundwater recharge and may result in overestimation of recharge if not included [32,35], The water balance between runoff and evapotranspiration and between evapotranspiration and precipitation depends on properties of the land surface, which are coupled to the climatic forcing by atmospheric demand and precipitation [36]. Hence, incorporating the climate outputs with land surface characteristics is required to minimize the error in water balance estimation for groundwater potential and recharge studies.

The objective of this study was to identify the groundwater potential zones (GWPZ) of the Baro River watershed through a hydro-geospatial modeling approach. The synthesis of the environmental variables can be done using Geographic Information Systems (GIS) and Remote Sensing (RS) to model the groundwater potential zones (GWPZ). GIS has high scalability, so it can be used in high-scale studies, and its raster processing capabilities with open-source software, like GRASS, are also powerful. For instance, some of the modules in the stream, including r.stream.extract and r.mapcalc, enable accurate calculation of secondary components required for infiltration modeling, such as drainage density and Topographic Wetness Index (TWI).

In the absence of comprehensive extensive drill data, the work used a multi-step theoretical validation approach to ensure consistency of GWPZ. A statistical sensitivity analysis is given to test the stability of the model and the effect of every environmental theme layer on the final results by carefully omitting some distinct theme layers as provided by [5,28,37]. Second, to verify regional alignment and provide a solid scientific foundation for groundwater management without the need for conventional well-depth measurements, the results were finally compared with earlier research by [19], geophysical survey in the Baro - Akobo- Sobat Basin.

2. Materials and Methods

2.1. Baro River Watershed Regional Settings

2.1.1. Physiography and Hypsometry

The Baro River catchment is a big tributary of the Sobat-White Nile system, and it’s located in southwest Ethiopia. The area is approximately more than 23,000 km² in a river basin (Figure 1) [19,38,39]. The regional hydraulic potential is related to the hypsometric gradient. In groundwater recharge modeling, the topographic slope (β) is determined by the slope of the digital elevation model (DEM) in both the x and y directions [19].

The steeper slopes of the highlands generate very short surface runoff time and also have a high gravitational potential, and the lowland areas offer residence time (Tr), which is used for effective deep percolation into the saturated zones. To differentiate between low-velocity and high-velocity runoff, high-resolution data is needed to identify the discontinuity features in the topography [15].

2.1.2. Hydrology and Drainage Network

The drainage system is influenced by tributaries of the Baro River watershed, including the Birbir, Geba, and Sor rivers [1,38,39,40]. After displaying a discrete dendritic system in the upper regime of the watershed, the drainage systems in the alluvial plains (lowlands of the watershed) changed to meandering forms, which are indicative of consistency resistance in the volcanic bedrock [2,4,41,42,43]. Parameters used to assess the hydrological efficiency of a watershed are: Drainage density (D_d) is calculated by stream length/total area (equation 1) [10,27].

$${\mathbf{D}}_{\mathbf{d}}={\displaystyle \frac{{\sum }_{\mathbf{i}=1}^{\mathbf{n}}{\mathbf{L}}_{\mathbf{i}}}{\mathbf{A}}}$$

(1)

where: L is stream length and A is area (23,000 Km²). The prioritized area (low drainage density) in the approach to modeling indicates high surface permeability and infiltration potential [2]. Contrary to this, the dissected highlands of the Baro River watershed have high drainage density values, caused by the high runoff coefficient and impermeable surface of dense drainage networks [14].

2.1.3. Hydrogeological Framework

The geology of the watershed of the Baro River defines the three groundwater units of occurrence in the watershed [38,39,40].

Basement rocks (Precambrian age): The origin of the Precambrian basement rocks with insignificant primary porosity is represented by the localized outcrops; the groundwater is constrained to the localized fracture zones and worn mantles [38,39,40].

Tertiary volcanic: mainly rhyolites and basalts in the highlands of the watershed, lacking substantial primary inter-granular porosity. Groundwater storage in the subsurface is controlled by secondary porosity (faults, joints, and vesicles). Accordingly, in the vertical recharge setting, these fractures at the surface appear as lineaments [4,21,23].

Quaternary alluvial are unconsolidated sediments with a considerable primary porosity that are present in Gambella's lowlands with high hydraulic permeability and significant sedimentary thickness, making them high-productive aquifers zones [38,39].

2.2. Data Type and Collection

The reliability of hydrogeospatial modeling for the assessment of Groundwater Potential Zones (GWPZ) in data-scarce areas such as the Baro watershed necessitates the integration of high resolutions and several sources of remote sensing products [44,45]. The study watershed's (23,000 km²) topography, pedologic, spectral, and climatic variations were recorded using a multisensory inventory. To improve computational accuracy and consistency, all layers were projected using EPSG: 32637-WGS84/UTM ZONE 37 and harmonized at 30 m resolutions.

2.2.1. Topographic Data

Topography is an important element of the hydrological cycle, and it regulates the flow between the rapid runoff and infiltration processes. The terrain data was high-resolution (30 m) data that was provided by the OpenTopograpty webpage (USGS, 2016). The GRASS algorithm (r.fill.dir) was used to preprocess the raw Digital Elevation Model (DEM) data in QGIS by removing false spikes and filling sinks to ensure a hydrologically reasonable surface [46].

A projected 30-meter DEM was used to extract several GWPZ factors. Including:

Slope gradient (β): The velocity of surface runoff directly is related to slope. Lower slopes will enable the water to remain on the surface in longer periods of time, an aspect that will promote the process of infiltration [47].

Topographic Wetness Index (TWI): After sinks were filled and the DEM's flow directions were determined, the topographic wetness index was computed to determine the local topography-based locations where water prefers to accumulate (Equation 2).

$$\mathbf{T}\mathbf{W}\mathbf{I}=\mathbf{l}\mathbf{n}\left({\displaystyle \frac{\mathsf{\alpha }}{\mathbf{t}\mathbf{a}\mathbf{n}\mathsf{\beta }}}\right)$$

(2)

where: α is the local upslope contributing area /contour length, and β is the slope angle.

Curvature: Flow across the surface is accelerated and decelerated by using the r.slope.aspect GRASS module from the DEM curvature calculation.

Lineament density (L_d): The lineament density, which represents structural weakness (fractures/faults) used for principal conduits for groundwater transport, was successfully calculated using line density from the GRASS module [48]. Computed using equation 3 below:

$${\mathbf{L}}_{\mathbf{d}}={\displaystyle \frac{\sum {\mathbf{L}}_{\mathbf{l}\mathbf{i}\mathbf{n}}}{\mathbf{A}}}$$

(3)

where: $\sum {\mathbf{L}}_{\mathbf{l}\mathbf{i}\mathbf{n}}$ is the total length of identified lineaments (cumulative sum of length of all mapped)

2.2.2. Lithology and Land Use/Land Cover Data

To determine the Groundwater Protection Zone (GWPZ) for the study area, surface properties, including land use and land cover (LULC), were integrated with subsurface geological properties through modeling. This integration utilized vector-based geographical data alongside high-resolution land cover images. The combination of high-resolution land cover imagery and vector-based geological data enabled the modeling of both subsurface and surface properties in the research area, which were essential for defining the GWPZ.

Geology: Structural and lithological information from Ethiopia-Geoportal was utilized to assess the hydrogeological potential of the watershed. This assessment involved defining primary permeability and aquifer storage capacity, thereby characterizing the region's water resource potential effectively [6].

LULC: The land use and land cover products from LULC, ESRI Sentinel-2 products, were available to provide data on the land use and land cover, which affects the surface hydraulic roughness and the tendency to evapotranspiration [49].

2.2.3. Soil and Climate Data

The recharge input of the system was determined based on the distribution of rainfall across the region by CHIPRS (chc.ucsb.edu/data/chirps) and the physical properties of the soil by FAO.

Soil: FAO had a shape file containing soil data. Soil texture and drainage class parameters were used to calculate the surface interface's infiltration capacity (fao.org/soils-portal).

Rainfall is the primary source of groundwater recharge, as determined by Climate Hazards Group Infrared Precipitation with Station data (CHIRPS). The mean annual rainfall was computed using data from 12 gridded meteorological stations (Gambella, Bure, Uka, Yibdo, Dembidolo, Mettu, Shebel, Guliso, Simbo, etc.) in the delineated study area. Inverse Distance Weighting (IDW) interpolation techniques were used to compute aerial precipitation [50].

To account for the impact of evapotranspiration, the precipitation input from the macroclimatic model was combined with the characteristics of land use/land cover in data-scarce areas, using Multi-Criteria Decision Analysis (MCDA) [51,52]. High-resolution surface hydraulic properties and localized transpiration demands from the LULC layer from ESRI Sentinel-2. High-density vegetation and wetlands have high actual evapotranspiration (AET) and reduce water infiltration to deep aquifers [53]. When the weightage of the land use land cover parameter with high evapotranspiration is reduced, the model becomes a significant balance between gross precipitation input and atmospheric loss to the system [51,52].

Table 1. Summarized data type and sources.

Data category/type	Specific layer	Primary sources	Access URL/Portal
Land cover	LULC	ESRI Sentinel-2	livingatlas.arcgis.com
Elevation	DEM	OpenTopography/ALOS PALSAR/	opentopography.org
Geological	Lithology, Lineament density	EthiGeoportal	ethiogeoportal.gov.et
Pedological	Soil texture, type	FAO	fao.org/soils-portal
Climatic	Mean annual rainfall	CHIRPS	chc.ucsb.edu/data/chirps

Table 1. Summarized data type and sources.

Data category/type	Specific layer	Primary sources	Access URL/Portal
Land cover	LULC	ESRI Sentinel-2	livingatlas.arcgis.com
Elevation	DEM	OpenTopography/ALOS PALSAR/	opentopography.org
Geological	Lithology, Lineament density	EthiGeoportal	ethiogeoportal.gov.et
Pedological	Soil texture, type	FAO	fao.org/soils-portal
Climatic	Mean annual rainfall	CHIRPS	chc.ucsb.edu/data/chirps

2.3. Geospatial Resources Analysis Support System (GRASS) GIS

The hydro-geospatial modelling in the Baro River watershed has been performed by the GRASS GIS 8.x (Geospatial Resources Analysis Support System) environment. It can be used with more flexibility than traditional GIS applications with GUIs, such as ArcGIS [46], due to its open-source modular structure, high accuracy of the raster engine even for very large watershed discretization, and good ability to process topological data structures.

Groundwater potential zones need to be modeled to achieve proper hydrological characterization of the study area, which involves applying proper processing lines or techniques to convert basic topography data into important hydrological factors or predictors [5,17,23]. The workflow started with using the r.fill.dir GRASS module to remove the sinks in order to create a hydrologically consistent surface. The groundwater recharge zones were then identified (or located) to calculate the flow accumulations and drainage networks using the r.watershed module. Additionally, the r.slope.aspect module in GRASS was used to govern surface residence time (a gentle slope offers vertical percolation), where the local evapotranspiration and soil moisture retention are influenced by aspects, in order to comprehend the infiltration potential of the study area's slope gradient (β). Using a combination of these output results, in particular upslope contributing area (α) of r.watershed and slope gradient of r.slope.aspect, one can compute the Topographic Wetness Index (TWI) in order to map or define the spatial distributions of potential zones and groundwater vulnerability throughout the study area watershed [45,46,54].

2.4. Multiple Criteria Decision Analysis (MCDA)

This interaction among climate, geology, topography, soil, and LULC governs the complexity of the groundwater occurrence in the study area and needs multiple decision-making frameworks. The Analytical Hierarchy Process (AHP) was used to weigh and incorporate several environmental factors into a single potential map.

2.4.1. Development of Factor Architecture

Based on the impact of recharge, storage, and transitions, this study identified nine environmental primary variable elements that were identified from Ethiopian highlands literature (usually utilized 6 to 7 factors only) [1,2,9,11,14,15,39,55]. Each of these was considered as a theme rating (1–5) from very low to very high potentials by the authors of QGIS (GRASS). The following features were used: rainfall (orographic effects) [21,56], lithology (geology) (define aquifer storage capacity and hydraulic conductivity) [31,56], lineament density (structural conduits, fractures, and faults) [55,57], slope (control infiltration to runoff ratio) [7], drainage density (inverse indicator of permeability) [2], soil texture (infiltration rate and soil atmosphere interface) [15], topographic wetness index (TWI-quantify topographic controls) [37], curvature (decelerate flow and concentrate flows) were used.

2.4.2. Analytic Hierarchy Process (AHP) and Consistency Verification

The AHP approach, developed by [29,58]), was used to assign the relative weight of these nine criteria by using the Pairwise Comparison Matrix (A). Each component was compared with each other, and the comparison was made on a scale of 1 to 9, in which 1 is equal importance and 9 is extreme importance [29], To ascertain the impact of the nine environmental elements, four condensed steps on AHP approaches were used.

Because of the variation among individual geological formations, the technique used in the study is a cross-section between lineament density layers and lithology layers as opposed to assuming each formation was homogeneous. In fractured zones, for instance, another volcanic (tertiary) layer was given a rank of 4 due to its secondary porosity and infiltration rate characteristics, whereas in an un-fractured zone it was given a lower rank. This method suggests that it is capable of representing localized structural deformations and has a strong representation of aquifers.

Step 1: Pairwise Comparison Matrix (A): The reciprocal matrix is compared to scales developed by [29] and can calculated using equation 4 below.

$${\mathit{\alpha }}_{\mathit{i}\mathit{j}}={\displaystyle \frac{1}{{\mathit{\alpha }}_{\mathit{j}\mathit{i}}}}$$

(4)

Table 2. Pairwise comparison matrix.

Factor	R	G	L	S	D	So	Lu	T	C
Rainfall (R)	1	2	3	4	5	6	7	8	9
Geology (G)	1/2	1	2	3	4	5	6	7	8
Lineament (L)	1/3	1/2	1	2	3	4	5	6	7
Slope (S)	1/4	1/3	1/2	1	2	3	4	5	6
Drainage (D)	1/5	1/4	1/3	1/2	1	2	3	4	5
Soil (So)	1/6	1/5	1/4	1/3	1/2	1	2	3	4
LULC (Lu)	1/7	1/6	1/5	1/4	1/3	1/2	1	2	3
TWI (T)	1/8	1/7	1/6	1/5	1/4	1/3	1/2	1	2
Curvature (C)	1/9	1/8	1/7	1/6	1/5	1/4	1/3	1/2	1

Table 2. Pairwise comparison matrix.

Factor	R	G	L	S	D	So	Lu	T	C
Rainfall (R)	1	2	3	4	5	6	7	8	9
Geology (G)	1/2	1	2	3	4	5	6	7	8
Lineament (L)	1/3	1/2	1	2	3	4	5	6	7
Slope (S)	1/4	1/3	1/2	1	2	3	4	5	6
Drainage (D)	1/5	1/4	1/3	1/2	1	2	3	4	5
Soil (So)	1/6	1/5	1/4	1/3	1/2	1	2	3	4
LULC (Lu)	1/7	1/6	1/5	1/4	1/3	1/2	1	2	3
TWI (T)	1/8	1/7	1/6	1/5	1/4	1/3	1/2	1	2
Curvature (C)	1/9	1/8	1/7	1/6	1/5	1/4	1/3	1/2	1

Step 2: Normalization and weight computation (W): The matrix was normalized by dividing each cell by the total of the columns, and the relative weights (W_i) were obtained by averaging the rows of the normalized matrix (equation 5) (Table 3).

$${\overline{\mathbf{a}}}_{\mathbf{i}\mathbf{j}}={\displaystyle \frac{{\mathbf{a}}_{\mathbf{i}\mathbf{j}}}{{\sum }_{\mathbf{i}=1}^{\mathbf{n}}{\mathbf{a}}_{\mathbf{i}\mathbf{j}}}}\mathbf{a}\mathbf{n}\mathbf{d}{\mathbf{w}}_{\mathbf{i}}={\displaystyle \frac{{\sum }_{\mathbf{j}=1}^{\mathbf{n}}{\overline{\mathbf{a}}}_{\mathbf{i}\mathbf{j}}}{\mathbf{n}}}$$

(5)

Step 3: Consistency verifications (CR): To assess the judging reliability, the consistency ratio (CR) was computed. A valid model is shown by CR < 0.10 (Table 3) (equation 6).

$$\mathbf{C}\mathbf{I}={\displaystyle \frac{{\mathsf{\lambda }}_{\mathbf{m}\mathbf{a}\mathbf{x}}-\mathbf{n}}{\mathbf{n}-1}} \mathbf{a}\mathbf{n}\mathbf{d} \mathbf{C}\mathbf{R}={\displaystyle \frac{\mathbf{C}\mathbf{I}}{\mathbf{R}\mathbf{I}}}$$

(6)

Where: $n=9andRI=1.49$

2.5. Overlay Analysis and Raster Algebra

Using raster calculators in the GRASS GIS module (r.mapcal), which enable sophisticated, pixel-based map algebra, the GWPZ was synthesized by Weighted Linear Combinations (WLC) (equation 7) as the final way of modeling the hydro-geospatial of the research area. Because it is effective at preserving the statistical subtleties of the AHP weights throughout the whole study area and capable of performing floating point operations with high accuracy, this module is preferable to basic overlay tools [29,58].

$$\mathbf{G}\mathbf{W}\mathbf{P}\mathbf{Z}={\displaystyle \sum _{\mathbf{i}=1}^{\mathbf{n}}}\left({\mathbf{W}}_{\mathbf{i}}\ast {\mathbf{R}}_{\mathbf{i}}\right)$$

(7)

Where: GWPZ – the final Groundwater Potential Index for each cell, W_i – relative weight of factor i, and R_i- the reclassified rating of factor i.

Figure 2. Methodological framework of GWPZ.

Figure 2. Methodological framework of GWPZ.

3. Results and Discussion

The results of hydro-geospatial modeling provide a thorough spatial synthesis of the nine environmental parameters influencing groundwater potential and occurrence in the Baro River watershed research area. The detailed reclassified thematic layers and their subsequent integrations with the final Groundwater Prospect Zone (GWPZ) map will be shown in this part.

3.1. Thematic Factors Reclassifications and Analysis

With scientific insights from various prior findings, this section presents the chosen analysis of the nine parameters with rank, area proportions, geographic zoning, and detailed descriptions (Figure 3 and Figure 4 and Table 6). These parameters are integrated by the GRASS GIS environment.

Rainfall: Rainfall is a key input for regional hydrology and the main source of groundwater recharge [1,2,14,20,31,43] in the Baro River watershed. The Southeast regions have substantial moisture availability, as shown by the mean annual precipitation of the area with notable orographic effects from Figure 3a and Table 6. This maximum moisture availability of zoning for deep percolations, which makes up 6.05% of the entire watershed area, was rated at 5 (extremely high). On the other hand, due to high potential evapotranspiration and restricted infiltration, which account for 13.29% of the total area, the western lowlands of the watershed have the lowest precipitation available, represented by a rank of 1. Rainfall is designated as the primary function for recharging in the AHP model, and its weight is 0.31 (Table 3). This result showed that rainfall plays a significant role in deep percolations, which is consistent with regional water balance analysis by [20], geochemical and isotope validation by [1] and recent findings by [42] for hydrological extremes in the Baro-Akobo basins and other basins in Ethiopia, such as Upper Blue Nile (Chemoga watershed) for GWPZ delineation by [2] and drought-prone area analysis by [31].

Vertical recharge from local precipitation is limited by plains in the western part of Gambella. The very high groundwater potential in this low-lying area persists despite the high potential evapotranspiration. The variations are due to the difference in storage capacity of the quaternary alluvial deposits and the hydrologically active area of high primary porosity in the east and the fractured highland escarpments as a source of groundwater flow regions [1]. Then, the geometry of the basin and the underground geological properties are important modulating factors that protect the aquifer from localized losses by evaporation from the atmosphere [19].

Geology (Lithology): The geological architecture defines the subsurface environment's principal permeability and aquifer storage capacity. Precambrian basement rocks, Tertiary volcanic, and Quaternary alluvial deposits (found in river valleys, covering 2.97% of the watershed and ranked as 5 due to their highest permeability and storage capacity due to their unconsolidated nature) make up the complex lithology of the Baro River watershed, which is mapped from Figure 3b and Table 6 (supplementary materials). Due to exceptional secondary porosity resulting from fracture networks, 49.82% of the lithology covered by weathered volcanic found in eastern plateaus was graded as 4. Due to their non-porous nature and function as a regional aquitard [19], the remaining 45.2% of the area was made up of enormous bedrocks that were primarily found in enteral highlands and were ranked 1. These results were compared with those of other studies conducted in the Baro River watershed and other river basins in Ethiopia, such as the hydrogeological characterization of the Baro River basins by [1] which revealed that quaternary alluvium is the most prominent unit; other researchers discovered that weathered volcanic and alluvial deposits in the Upper Omo Gibe basins by [11] were prioritized for GWP; and, in the Chemoga River, water lithology is the second significant factor in geospatial and AHP modeling by [2]. Lithology holds 0.25 model weights (Table 3).

Lineament density (L_d): In hard-rock settings, lineaments, structural flaws like faults and fractures, act as the main channels for groundwater movement. The top level, ranked by 5 (tectonic faults with 2.85% of the watershed), offered considerable secondary porosity [21], as shown by the geographical distributions of the lineament density shown in Figure 3c and the characterization of L_d in Table 6 (supplementary materials). 7.68% of the watershed has fracture networks that enhance permeability and high-density zones at the highland borders [23,38], which are graded as 4. Compact strata with minimal structural porosity, especially in the western plains (41.72%), were classed as low density [4]. Lineament density has model weights of 0.16 (Table 3).

Slope: One important factor that establishes the limits between surface runoff and vertical infiltrations is slope gradient. The large hypsometric gradient in the research region is seen in Table 6 below. The Gambella plains, which made up 30.22% of the basins, had a 0–2% slope and a level topography. Because it maximizes the water residence period and promotes deep percolations to saturated zones, this location received a score of 5 (very high). With a rating of 4 (high) and moderate slopes (2–8%), 43.73% of the higher pediments are covered by high infiltration with minimal erosion danger. In the eastern escarpments, more than 35% of the slope, which is categorized as steep slopes and classed as 1, produces fast surface runoff that stops infiltration into subseries. Table 3 shows the weighted topographic slope of 0.11.

Drainage density (D_d): Drainage density is an inverse measure of permeability; a well-developed stream network efficiently removes water as runoff from a high-density concept [1,21]. The research area's dissected highlands have a very high drainage density (ranked by 1), which replicates the zone of impermeable surfaces that quickly drain water away, as seen in Figure 4e and noted in Table 6 below. On the other hand, 4.74% of the research area is situated in western lowlands with spatially low drainage density (ranked 5), making it more water-absorbing than sheds. 51.66% of the research area's enteral watershed is covered by the most extensive areas with moderate drainage density [13,59]. The hydrological efficiency of watersheds is determined by these parameters, which have model weights of 0.08.

Soil and LULC: Soil texture and land use and land cover (LULC) serve as essential surface regulators for groundwater recharge by regulating infiltration rates and hydraulic roughness throughout the research region. From below, Figure 4f and Table 6 show that the lowland regions are mostly covered by 30.20% of the area of sandy and alluvial soils (ranked 5), which permit strong hydraulic conductivity, while ranked 4 (sandy loams) facilitate vertical seepage throughout 43.8% of the pediments. On the other hand, thick clays (ranked 1) at highland peaks prevent infiltration and cause excessive surface runoff. The majority of the land covered by forest and woodland (72.73% of the area, ranked 4) uses root systems to encourage deep percolation, while 1.43% of the land covered by urban areas (ranked 1) forms impermeable surfaces that completely prevent recharge. LULC patterns interact with soil qualities [21,38].

TWI and Curvature: Surface curvature and the Topographic Wetness Index (TWI) enhance the identification of water accumulation zones and flow dynamics. The TWI rank 5, covering 15% of western floodplains, indicates high saturation potential and moisture retention. Figure 4h shows valley bottoms rated 4 and concave depressions rated 5, serving as key points for surface water convergence and recharge. These features indicate significant groundwater prospects through the integration of favorable vegetative cover, soil texture, and concave topography in floodplains and lower pediments.

3.2. Groundwater Potential Zoning (GWPZ)

In the study area the GWPZ was synthesized using the Weighted Linear Combination (WLC) module and the GRASS GIS environment (r.mapcal). The r.mapcal Grass module was employed to ensure the precision performance operations [46,60]. Figure 5, of groundwater potential distribution of the Baro River watershed below, shows the potential level. The classification was done using their capacity to replenish aquifers and generate sustainable yields through a combination of hydrological, climatic, lithological, and geomorphological factors.

Very high potential (Zone 5): these are mostly found in the northern Gambella lowlands, which have a flat slope to maximize residence time for deep percolations. They are distinguished by a large region of quaternary alluvial sediments that give high primary porosity [61].

High potential (Zone 4): Expansive regions, which are primarily found in the Eastern and Southern volcanic highlands, are areas of secondary porosity created by cracks and joints. The systems are more susceptible to seasonal rainfall variability than low-land alluvial aquifers because of the presence of broken networks in geometry, which result in minimal storage even when maximum rainfall occurs.

Moderate potential (Zone 3): The center and eastern parts of the watershed are referred to as transition areas and fall into moderate potential (Zone 3); there is a balance between highland runoff and low land deposition. The strong lineament-driven recharge from the volcanic highlands makes this area less porous than that of the Gambella plains, which has a consistent and moderate inflow yield [21].

Low to very low potential (Zones 1 and 2): The western and northwest portions of the research area are comprised of large Precambrian basement rocks, with low to very low potential (Zones 1 and 2). In contrast, with the high slope and limited permeability, most of the annual rainfall in this area will drive to create rapid surface runoff instead of aquifer recharge, as would occur in the permeable alluvial lowlands.

The hydrogeological profile analysis and groundwater distribution (Figure 5) show that over 90% of the examined region has moderate to high groundwater potential. The region's primary potential zone was made up of broken volcanic hills and alluvial sediment plains. The result will offer a solid hydro-geospatial foundation for water resource management in the basin and other data basins. This study advances large-watershed Groundwater Potential Zone (GWPZ) mapping by incorporating additional environmental factors and leveraging cutting-edge GIS tools, including GRASS modules within QGIS 3.40.15. It newly identifies high-potential fractured volcanoes in the southeastern region, building on [62]. Previous research accurately delineated high-yield zones in smaller western alluvial plains, attributing this to their flat topography and sediment deposits [63].

Overall, the GWPZ shows moderate to extremely high stability for 90% of the watershed. Geology-derived sustainable extractability and precipitation-derived recharge potential should be separated. Actual extraction amounts are determined by the storage capacity, which is dependent on local variations in aquifer characteristics, seasonal rainfall differences, and subsoil storage capacity. The study identifies two contrasting geographical areas, namely the Eastern Volcanic Highlands and the Gambella low-lying plains, both of which have maximum rainfall, but the Eastern Volcanic Highlands has lower sustainability potential with higher slopes and rapid drainage, while the Gambella low-lying plains have high sustainability potential with unconsolidated Quaternary alluvial sediments with high porosity, thickness, and high storage capacity [3,11].

3.3. Validation Framework

Since the well yield data were not available, the researcher used a multi-faceted approach by using consistency verification based on AHP and comparative analysis in the regional studies to evaluate the dependability of the model based on groundwater potential zoning. This approach ensures resilience and application, especially when data is limited, as is the case in Ethiopia [21].

The above table 3, which is produced from the pairwise comparison matrix, is calculated through the Consistency Ratio (CR), which is less than the threshold level (<0.10), meaning that the criterion weights for elements such as slope, geology, and rainfall have been validated. The recent studies carried out in Ethiopia have established CR ranging between 0.04 and 0.09 for groundwater prospecting by the professionals, which confirms the above-mentioned [11,15]. Studies conducted in the Gidabo watershed (75% moderate–high potential) [64] and Chemoga watershed [2] are consistent with this study, which found that the Baro River watershed's high potential zones were driven by nine elements (rainfall, geology, and slope) during AHP mapping. When applied to the Omo Gibe basin, recharges and geology were identified as key factors by the ROC test, which was performed against the transmissivity value obtained from the integrated AHP-WetSpass model [11]. This shows the effectiveness of AHP in the tectono-geomorphic highlands of Ethiopia [11,14,21,55,59,65].

An independent approach using physically based temporal simulation outputs was added to improve the validation program for the SWAT+ model. The analysis in Figure 6 shows that the simulated groundwater recharge is correlated with the base flow and that groundwater recharge peaks determine the amount and timing of base flow. The GWPZ map integration shows that there is an active recharge process in the high potential zones. Previous results are also confirmed by temporal trends, with high correlation between base flow responses and fluctuations in recharge, suggesting hydrologically active spatial areas [2,3,11,36,41,43]. This is an additional validation method that improves the efficiency and accuracy of the water resources model in the region.

4. Conclusion and Recommendation

4.1. Conclusion

To design groundwater potential possibilities in the Baro River watershed, a hydro-geospatial technique that included remote sensing, GIS, and Multi-Criteria Decision Analysis (MCDA) via Analytical Hierarchy Process (AHP) was used. Nine environmental parameters were included in this study to effectively represent the complex interplay between climatic, geomorphological, and land surface characteristics that govern groundwater occurrence and movement. Significant groundwater potential is indicated by the moderate to very high groundwater potential suitability grades found in over 90% of the examined region. Because of their level topography and loose alluvial deposits that promote infiltration and groundwater recharge, the Gambella lowlands in the northern parts of the watershed have the most potential. However, the eastern and southern parts of the watershed have substantial potential, mostly because of fractured volcanic strata where groundwater storage is regulated by secondary porosity. However, cautious basin management is required since these volcanic fractured formations are more vulnerable to precipitation variability. The GIS and AHP integration approach worked well and was appropriate for places with little data. The calculated Consistency Ratio (CR = 0.067), which shows how consistent the weighting scheme and decision-making process are, is less than the permissible threshold value of 0.10. This article has shown that despite the absence of large amounts of filed information, GIS-MCDA methods can produce scientifically complete groundwater potential maps. Altogether, this study contributes to the baseline hydrological understanding of understudied of the study watershed by providing a systematic, open-sources, geospatial workflow that can assist water resources planners in data-limited regions.

4.2. Recommendation

According to the findings of the research, a number of specific recommendations can be made. Its growth must concentrate on regions that are highly viable, such as the Southern Volcanic, and the Gambella lowlands. The use of conservation methods of water, including contour bunds and check dams, is necessary in improving the availability of groundwater, particularly in areas of high runoff. The model validation requires field verification and geophysical surveys. Basin-scale planning should include the potential of groundwater to help minimize the drought vulnerability and encourage the use of surface and ground water in an integrated manner. Having a centralized GIS database will facilitate the continuous assessment and decision-making. Further study is necessary to employ hydro-geochemical and isotopic techniques to understand better how to manage groundwater flow and sustainability, especially in fractured aquifer regions.