Analyzing the Impact of Ungauged Hill Torrents on the Riverine Floods of the River Indus: A Case Study of Koh E Suleiman Mountains in the DG Khan and Rajanpur Districts of Pakistan

Abstract

Floods are one of the most destructive natural hazards in Pakistan, causing significant damage. During monsoons, when westerly winds and concentrated rainfall occur in rivers’ catchments, floods become unmanageable. Given the limited resources of Pakistan, there has been minimal effort to quantify the amount of rainfall and runoff generated by ungauged catchments. In this study, ten hill torrents in Koh e Suleiman (District Rajanpur and DG Khan), an area affected by flash flooding in 2022 due to extreme precipitation events, were investigated. The Hydrologic Engineering Centre’s Hydrologic Modeling System (HEC-HMS), a semi-distributed event-based hydrological model, was used to delineate streams and quantify runoff. Statistical analysis of the rainfall trends was performed using the non-parametric Gumbel extreme value analysis type I distribution, the Mann–Kendall test, and Sen’s slope. The results of the study show that the total inflow to the river Indus is 0.5, 0.6, 0.7, and 0.8 MAF for 25, 50, 100, and 200 years of return period rainfall, respectively. This study presents appropriate storage options with a retention potential of 0.14, 1.14, and 1.13 MAF based on an analysis of the hydrology of these hill torrents to enhance the spate irrigation potential as flood control in the future.

1. Introduction

Several mountainous regions of Pakistan are vulnerable to flash floods, which are considered catastrophic torrents [1]. During the monsoon of 2022, Pakistan experienced its worst flood in the past ten years due to extreme rains. According to Pakistan’s National Disaster Management Authority, approximately 1100 people were killed, 33 million were affected, and 1 million homes were destroyed or damaged by the floods. The worst flooding occurred along the Indus River in Punjab, Khyber Pakhtunkhwa, Balochistan, and Sindh. Approximately 150 bridges and 3500 km (2200 miles) of roads have been destroyed across the country, according to Relief Web. In addition, 2 million acres of crops and orchards, as well as more than 107,000 animals, have been destroyed. In the DG Khan District, 342 villages were damaged, 80 union councils were flooded, and 699,502 people were directly affected. Figure 1 shows the before and after flooding situation in Rajanpur as reported in the official reports, in which hill torrents affected close to 100,000 people and inundated 309,000 acres of agricultural land.

A total of 13 hill torrents have been identified in the vicinity of DG Khan and are known as Kaura, Sanghar, Vehova, Sorilund, Vidor, SakhiSarwar, Mithawan, Kaha, Chadhar, Pitok, SoriShumali, Zangi, and SoriJanubi, which could act as a conduit for floodwater from the nearby catchment [2]. The torrents, as mentioned above, enter the Indus River from the right bank of the Chashma River, the DG Khan canal, and the Kachhi canal [3]. Hill torrents from the Koh-e-Suleiman Range enter the Indus between the Taunsa Barrage in Punjab and the Guddu Barrage in Sindh. Watersheds feeding these hill torrents are ungauged, making it impossible to provide reliable information on their contribution to flood water in the Indus River [4]. Thus, the hill torrents emanating from the Koh-e-Suleman range have caused havoc for flood management officials in Sindh and Punjab.

Climate change has contributed to the increase in flood frequency and magnitude [5,6]. One of the impacts of climate change in Pakistan is concentrated rainfall that contributes to floods in the catchments of rivers during the monsoon season. These monsoon currents and wind conditions can intensify floods to an intolerable extent [4]. There is a need for a floodwater management plan to reduce the impacts of floods (hill torrents) and prevent them from recurring. It may consist of structures that can withstand large quantities of water and reduce the impact of hill torrents, especially during the monsoon season. Various models and methodologies have been used to quantify the runoff produced by ungauged catchments, which may result in flash floods caused by those torrents. Prior studies have focused more on quantifying rainfall runoff than developing mitigation strategies for flash floods. Due to Pakistan’s limited resources, it has contributed minimally to the flash flood routing and management research. In the wake of the 2022 flash flood in Pakistan, it is imperative that ungauged watershed runoff measures are quantified and that remedial measures are taken to divert flood water to the water-scarce areas. Flood water can also be conserved to be used during dry periods. This study aims to quantify hill torrents in DG Khan and Rajanpur. Furthermore, it examines remedial measures since those torrents caused significant economic damage and human casualties during the floods of 2010 and 2022. Moreover, a policy discussion raised questions as an entry point to start discussion to amend customary water laws/rules, how agroecological potential can be tapped avoiding colonial legacy, as well as what crop choices, and types of scale issues of governance paid attention.

Literature Review: Modelling Ungauged Catchments

To minimize the destruction caused by flash floods, it is necessary to quantify them to make early warnings, timely preparations, and appropriate adaptation strategies. Research studies examined various methods for quantifying flash floods in ungauged catchments. For example, Sene (2013) has explored multiple conceptual, data-driven, physically based, and probabilistic flood forecasting models [7]. A study conducted by Isabelle Braud et al. (2010) examined a flood using two distributed hydrological models: CVN and MARINE. The CVN model is found to have a greater range of uncertainty than the MARINE model [8]. There is a higher sensitivity of the CVN model to Manning’s roughness coefficient than the MARINE model. Neither model evaluated the relevance of calculating the runoff coefficient during post-maximum discharges. Mishra et al. (2008) developed an empirically based hydrological model for paddy agricultural watersheds with limited hydrometeorological data [9]. There was a deviation in the results of flood peaks in the range of 9–33%, and a deviation in the runoff coefficient was 4–11%. A lack of data prevented the model from efficiently comparing observed and modeled results, as Nash–Sutcliffe coefficient values (0.10–0.55) and correlation coefficient values (0.45–0.66) were calculated.

Camarasa-Belmonte (2016) analyzed flash flood events in five Mediterranean ephemeral streams in Spain to better understand the semi-arid fluvial system [10]. The percent accumulation curves of rain and flow showed strong similarities at the beginning of the flood when rainfall intensities were higher. On the other hand, higher intensities at the end showed dissimilarities in both curves. Additionally, high rainfall intensity shortens the response time of the basin, while high amounts result in flood peaks. Black box models perform better in high-intensity events, while distributed or semi-distributed models perform better in low-intensity events. Adamovic et al. (2016) incorporated a simple dynamical system approach into the distributed hydrological model and named it SIMPLEFLOOD [11]. The study produced satisfactory results over the entire period. In the wet years, the model simulation performed well; however, in the dry years, the simulation performed poorly. Rozalis et al. (2010) used an uncalibrated hydrological model to simulate the watershed of the Mediterranean Sea, covering an area of 27 km² [12].

The study examined limited data usage for rainfall–runoff modeling, land-use-change impacts on runoff, and the impact of rainfall distribution on flash floods. Based on the study’s results, the model performed well in predicting peak flow discharges, but depended on the storm type. The importance of curve number (CN), rainfall amount, and rainfall intensity has been noted to be crucial in simulating runoff production since these factors affect the magnitude of a runoff flow. With the use of a distributed hydrological model, Zoccatelli et al. (2010) have investigated the dependency of rainfall variability on flash flood modeling [13]. An examination of three extreme flash floods that occurred in Romania between 2005 and 2007 is presented in this study. According to the study, rainfall’s spatial variability significantly affects the flash floods prediction since the Nash–Sutcliffe coefficients are less than 0.8 in two cases and 0.6 in one. Using geology and rainfall variability to determine the impact of flash floods, Zanon et al. (2010) studied the flash flood event in Western Slovenia that occurred on 18 September 2007 [14]. The study found errors in modeled flood peaks in rainfall volume. A low runoff coefficient is found due to the low soil moisture in the initial conditions, ranging from 0.17 to 0.24. Yasin and Nabi (2014) analyzed the impact of the Mithawan hill torrent in the DG Khan area by using the semi-distributed hydrological model HEC-HMS to assess the degree of damage caused by this hill torrent to the Kachhi canal [15]. Research has demonstrated that the HEC-HMS model successfully determines the flow peaks within acceptable limits.

An evaluation of rainfall runoff patterns in the upper Baitarani River Basin in east India was carried out by Verma et al. (2010) using the HEC-HMS and WEPP hydrological models [16]. Based on the study’s results, both hydrological models simulated lower stream flows during the validation period. During the simulation range of different ef-ficiency coefficients NSE, R² was found to be between 0.63–0.83 and 0.73–0.84, respec-tively. Based on the results, the HEC-HMS model performed better than the WEPP model in simulating daily stream flows for the upper Baitarani River basin. Using the Refh rainfall–runoff model, Joo et al. (2014) developed a comparison between the Refh rainfall–runoff model and HEC-HMS model in two catchments of Korea (Bukil and Jeungpyeong) located within the Guem River basin [17]. The model lumped characteristics allow it to perform well only in small catchments. The semi-distributed nature of HEC-HMS makes its performance more reliable.

Yan et al. (2015) have examined the flow routing of two rivers, (1) the Yuan River and (2) the Danube River using a Generalized Nash Model (GNM). This model uses Laplace transformation and mathematical induction [18]. It was concluded from the study that GNM is more accurate at predicting flows than the traditional IUH model. By improving the current information, it is concluded that IUH’s performance in forecasting short lead times can be enhanced. Paiva et al. (2011) conducted a hydrological and hydrodynamic modeling study in a region of the Purus River Basin with limited data and used the IPH-IV and MGB-IPH models [19]. The study concluded that little uncertainties and errors occur due to vegetation and cross-section geometry limitations of DEMs.

In another study, Chatterjee et al. (2014) estimated the runoff volume and peak discharge in the India—Damodar watershed using HEC-HMS [20]. According to the results, both runoff volume and peak discharge are affected by the impervious area and infiltration rate. A study by Urias et al. (2007) determined the probability of precipitation and flood return in Juarez, Mexico [21]. The statistical precipitation distribution is determined using the Hazen plotting position method. Additionally, it calculates the annual precipitation, which is then arranged in ascending order. Then, each event is ranked, which leads to the calculation of precipitation probability and return period.

2. Materials and Methods

2.1. Study Area

The study focus on the Koh e Suleiman region of Dera Ghazi Khan and Rajanpur in Punjab, Pakistan. D.G. Khan is located at latitude 30.0489° N and longitude 70.6455° E, at an average altitude of 124 m above mean sea level, as shown in Figure 2.

It extends from the southern part of the Hindu Kush Mountain system in FATA (Southern Federally Administered Tribal Areas) and Afghanistan. In addition, some of its parts are located in Khyber Pakhtunkhwa and the southwest of Punjab, while most of its parts originate from northern Baluchistan. Arid to semi-arid conditions prevail in Dera Ghazi Khan. There is an erratic precipitation pattern in this region’s hill torrent areas. The region usually experiences little or no rainfall following a heavy precipitation pattern.

2.2. Types of Data

2.2.1. Meteorological Data

Pakistan’s Meteorological Department (PMD) provides precipitation throughout the country. This study utilized 30 years’ worth of daily rainfall data from 8 PMD rainfall gauge stations from 1989 to 2008 (see

Table 1. Inventory of Rainfall Stations and Data.

S. No.	Station	Period of Record	Years of Record	Time Scale of Data
1	Barkhan	1989–2018	30	Daily
2	DG Khan	2003–2018	16	Daily
3	Multan	1989–2018	30	Daily
4	Khanpur	1989–2018	30	Daily
5	Rahimyar Khan	2003–2018	16	Daily
6	Jacobabad	1989–2018	30	Daily
7	DI Khan	1989–2018	30	Daily
8	Zhob	1989–2018	30	Daily

). A rainfall storm frequency analysis was performed to calculate the rainfall and runoff return periods.

2.2.2. Soil Data

The main governing parameter in rainfall–runoff modeling is soil texture. Soil textures were acquired through an open-source online database https://soilgrids.org/ (accessed on 8 August 2022). The database has a grid cell spatial resolution of 250 m for soil texture. Geotiff files containing clay and sand content at a depth of 15 cm are downloaded and utilized in the generation of Curve No.

2.2.3. Topographic Data

The topographic data was obtained from the website of the Japan Aerospace Exploration Agency (http://www.eorc.jaxa.jp, accessed on 8 August 2022). Advance Land Observing Satellite (ALOS) DEM is more suitable because of its higher spatial resolution (30 m × 30 m) and capability to cover high-altitude steep mountain regions. The ALOS DEM provides the basis to delineate the catchment and sub-watersheds.

2.2.4. Land Use/Land Cover Data

Land use and land cover data are used to determine which types of surfaces will produce more runoff and which will produce less. Landsat 8 imagery from the USGS website https://earthexplorer.usgs.gov/ (accessed on 8 August 2022) has been used to classify the study area’s land use and land cover. The classified land cover from the image contains the following classes:

• Open Water
• Barren Land
• Cultivated Crops
• Grassland
• Pasture
• Subshrubs

2.3. Data Map Preparation

To delineate the watershed of the study area, spatial data are prepared using Arc GIS 10.5.1 software. ArcHydro and the HEC-GeoHMS toolbar in ArcGIS are used to preprocess data for watershed delineation. In order to calculate different parameters, such as shape area, shape length, and other parameters, the downloaded DEM is projected. Since the study area lies within UTM zone 42, the DEM is projected in UTMzone42N using the ArcGIS project tool. Using the Extract by mask tool, the study area of interest (AOI) is clipped on the projected DEM. Every Digital Elevation Model (DEM) has a few sink issues that require adjustment following the neighboring grid. For this purpose, the tool “Fill sinks” is used, which adjusts values in the DEM. Using the “Fdr” tool in HEC-GeoHMS, the flow direction was calculated to understand the location of drainage patterns. Using flow direction results, “Fac” calculates the flow accumulation grid.

• 60,000 for large watersheds
• 5000 for small watersheds

The tool “StrLnk” links/segments the stream generated in the above process. Cat generates a grid of catchments draining into each other based on stream segmentation and flow direction. The already prepared catchment grid can be converted into polygons using the catchment tool. Next, a drainage line processing tool was used to create a drainage pattern to define the watershed’s boundaries further. In the end, small catchments are prepared using adjoint catchment processing tools.

2.3.1. Hec-Geo-HMS Characteristics

It is necessary to define the characteristics of the watersheds after they are generated, including the length of the river, the slope of the river, the slope of the basin, the longest flow path, the centroid of the basin, the centroid elevation, and the centroidal longest flow path. The soil maps are prepared using downloaded soil grids. Clay and sand content layers are converted to soil texture using the open-source software QGIS. The shape files are imported into ArcMap and clipped to the study area. The following textures characterize the soil type:

• Clay Loam
• Loam
• Sandy Clay Loam
• Sandy Loam

Further, the texture is assigned a hydrological soil group A, B, C, D. According to [22], assigning group A low runoff rate to D High runoff rate, as shown in

Table 2. Soil texture having different soil groups.

Object ID	Soil Type	Soil Code
1	Clay Loam	D
2	Loam	B
3	Sandy Clay Loam	C
4	Sandy Loam	A

.

Landsat 8 imagery was used to produce the land use and land cover maps. The raster data was then clipped to the study area using the Extract by Mask tool. Using ArcMap’s raster to polygon tool, raster data is converted to vector data. In hydrology, the runoff curve number is an empirical parameter used to calculate the amount of direct runoff and infiltration following a precipitation event. CN uses land cover values to predict the amount of runoff being produced. Paved areas are assigned a value of 0, while open areas are assigned a value of 100. The CN grid is constructed first by intersecting soil type, land use, and land cover data. In addition, a CNLOOKUP table containing soil type information is used. This

Table 3. Curve Number LOOKUP Table (Source: [22]).

Object Id	LUVALUES	Land Cover	A	B	C	D
1	1	Open Water	100	100	100	100
2	2	Barren Land	77	86	91	94
3	3	Cultivated Crops	67	77	83	87
4	4	Grassland	49	69	79	84
5	5	Pasture	49	69	79	84
6	6	Subshrubs	63	77	85	88

contains the object ID, land use values (LUVALUE), land cover data, and hydrological soil parameters. A CN grid is generated after the three-soil type, land cover, and CNLOOKUP tables are converted within HEC-GeoHMS in the utility of the CN grid. Following the above processes, input files are generated for exporting these datasets to HEC-HMS.

2.3.2. Hydrological Modeling

A hydrological model is a dynamic process that calculates all parameters associated with the water cycle, from evaporation to runoff. In this study, HEC-HMS is used. Many researchers have used HEC-HMS to quantify rainfall runoff due to its event-based nature. Through the use of soil type and land cover information, the loss method determines how much water will infiltrate after rainfall has occurred. Runoff modeling requires information on how much is infiltrated and lost. In HEC-HMS, eleven different methods are used to calculate losses. This study uses the soil conservation service curve developed by the United States Department of Agriculture (USDA). The formulas for different computing parameters for runoff using the SCS-CN method are provided below:

$$Q=\frac{{\left(P-Ia\right)}^2}{\left(P-Ia\right)+S}$$

(1)

$$Ia=0.2 S$$

(2)

$$S=\frac{1000}{\mathrm{CN}}-10$$

(3)

where “Q”, “P”, “Ia”, “S”, and “CN” are Discharge (IN), Rainfall (IN), Initial Abstraction (IN), Potential maximum retention after runoff begins (IN), and Curve Number.

The lag time is calculated using the transform method. A lag time is the amount of time between when maximum rainfall occurs and when peak discharge occurs. In this study, the SCS unit hydrograph method was used. For each Hill torrent and its sub-basins, lag time has been computed using the CN grid, longest flow paths, and basin slopes, as shown in

Table 4. Curve Number Values using Loss Method and Lag Time Computed.

Hill Torrent	Sub Basin	CN	Lag Time (Min)
Kaha	W970	79.56	579.46
	W910	84	461.11
	W1040	87.3	436.48
Sanghar	W830	86.2	316.27
	W770	82.62	436.44
Vidore	W110	88.6	154.32
	W120	90	93.411
	W130	87.2	153.50
Kaura	W140	87	128.6
	W120	91.4	85.70
Chachar	W1560	91.64	268.77
	W1530	88.16	184.58
Pitok	W300	89.27	100.96
	W410	90.17	99.417
Vehova	W410	84.5	422.49
	W450	86.4	240.11
SakhiSarwar	W230	90.8	64.47
	W300	91	57.73
RakhiMounh	W110	83.57	86.58
	W120	83.1	104.3
	W130	89	43.76
Mithawan	W380	83.74	102.94
	W430	90.26	47.07
	W440	87	80.335

.

2.4. Statistical Analysis

2.4.1. Frequency Analysis

The Barkhan rain gauge was selected for frequency analysis since it was the only gauge representative of all hill torrents spatially. Frequency analysis is performed to determine the return periods and probability of the rainfall events. Gumble Extreme Value Analysis Type 1 Distribution has been used for frequency analysis of the Barkhan rain gauge for 56 years (1963–2018). There have been many statistical analysis methods used for storm rainfall frequency analysis, but generally, Gumble Extreme Value Analysis Type 1 Distribution is selected based on the previous study used for its best fit to storm rainfall frequency data in the area [15]. However, in this analysis, a cross comparison is not drawn between different distribution to determine the best fit for study area rainfall data.

2.4.2. Rainfall Trend Analysis

Nonparametric Mann–Kendall and Sen’s slope methods are used to determine positive and negative trends in precipitation data with their statistical significance [23,24]. Three scenarios are analyzed, i.e.,

• Annual Maximum, Minimum, and Mean variation;
• Seasonal, Pre-Monsoon, Post Monsoon, and Monsoon variation;
• Monsoon Monthly variation.

2.5. Sensitivity Analysis

Various parameters used to generate runoff should be analyzed for their influence on runoff generation. In sensitivity analysis, lower and upper bounds are used to determine model sensitivity. Therefore, changes are made to one parameter while keeping the others unchanged. The parameters that HEC-HMS uses in simulating rainfall runoff generation are Curved No. (CN), Initial abstraction (Ia), impervious percent, and lag time. Analysis was carried out on 2, 5, and 10% changes in CN and Initial abstraction. Impervious percent and lag time has no significant impact because the catchment has a small portion of built-up area. Ahmad et al. (2021) [25] study the Al-Adhaim River catchment in Northern Iraq using HEC-HMS, and noticed that the main parameters which affect runoff quantities were the curve number and initial abstraction.

2.6. Storage Availability Analysis

The storage availability is assessed to identify the potential for spate irrigation for the local community residing adjacent to these torrents, as they divert some rainwater for agriculture making kaccha bunds. Based on basin area, three large basin area hill torrents have been chosen for storage availability, i.e., Kaha, Sanghar, and Vehova. DEM contours are used to analyze storage at all three hill torrents.

3. Results

3.1. Watershed Characteristics

Depending on the characteristics of the watershed, the quantity of water to be drained varies from one watershed to another. Several factors contribute to the characteristics of a watershed, including its topography, land use/cover, climate, and soil type.

3.1.1. Topography

Topography refers to features on the surface of the earth, such as mountains, rivers, valleys, and built-up areas. A watershed’s topography includes the basin’s area, shape, and slope. A Digital Elevation Model is used for topography, which has elevation values ranging between 75–3023 m, as illustrated in Figure 3a.

3.1.2. Basin Area

In a watershed, the basin area is the total area draining towards a common outlet. In order to model rainfall runoff, it is necessary to delineate the basin area of watersheds accurately. Using Digital Elevation Models (DEMs) and ArcMap software, basin areas for each of the ten major hill torrents are calculated as shown in

Table 5. Hill torrents showing their respective areas in Sq Miles. (Mi²) and four groups of soil texture and their respective areas.

S. No	Hill Torrent (Watershed)	Area (Mi2)
1	Kaha	2122
2	Sanghar	1848
3	Vidore	291
4	Vehova	1011
5	Mithawan	93
6	Pitok	90
7	RakhiMounh	40
8	Chachar	298
9	Sakhisarwar	41
10	Kaura	197
Soil Texture	Soil Type	Area (Mi2)
Loam	B	3168
Clay Loam	D	2033
Sandy Clay Loam	C	826
Sandy Loam	A	1.83

.

3.1.3. Basin Slope

The slope of a watershed basin provides a better understanding of water movement. The slope of the basin affects the time of concentration, which directly impacts the runoff volume. Runoff takes less time to reach the outlet on higher or steeper slopes. The slopes of the basins of each hill torrent in this study were computed using DEM, which has values ranging from 0–377 m, as shown in Figure 3b.

3.1.4. Land Use/Land Cover (LULC)

The characteristics of the watershed influence hydrologic responses. In rain-fall-runoff modeling, land cover and soil type are key determinants of a healthy and unhealthy watershed. A watershed’s land cover contributes to assigning curve numbers, facilitating runoff calculation after rainfall. Using Landsat 8 imagery, six land use and land cover classes are prepared for each hill torrent’s catchment: open water, barren land, cultivated crops, grassland, pastures, and subshrubs. Furthermore, LULC classification areas were calculated and analyzed for each hill torrent separately, which enabled a better understanding of the response of each portion of the land cover class. Based on the results of the analysis of individual hill torrents in the Koh e Suleiman mountainous range, it was determined that the majority of the land cover area consists of barren land with a range of 47–98%, grassland with a range of 2–50%, crops ranging from 13% to 17%, pasture ranging from 17% to 17%, and shrubs comprising 1–5%, as shown in Figure 3c. Land cover within the study area indicates a higher value of the Curve Number, which suggests that rainfall will result in a robust runoff pattern.

3.1.5. Soil type

The soil is divided into four types: A, B, C, and D, each with its runoff potential. As a result of the study, soil data of an average 15 cm depth was used to determine the presence of the following soil textures in the study area, as shown in

Table 6. Annual rainfall showing maximum, minimum, and mean trends.

Time Series	1963–1988	1989–2013	1963–2013
Annual Maximum	0.93	0.0001	0.20
Annual Minimum	−2.73 **	−0.05	−0.84
Annual Mean	−1.85 +	−0.91	−1.39
Pre-Monsoon	2.83 **	0.37	1.28
Post Monsoon	1.64	−1.48	−0.32
Monsoon	0.12	−0.12	−0.01

** Significance level _99%. + Significance level _90%. Bold = Significant Negative Trend.

and Figure 3d.

3.2. Statistical Analysis

3.2.1. Annual Maximum, Minimum, and mean Rainfall Trend Analysis

The Mann–Kendall non-parametric method was used to analyze rainfall trends. The rainfall data is divided into three-time series: annual maximum, annual minimum, and annual average. According to trend analysis, maximum time series data show an increasing trend, while minimum time series data show a decreasing trend. The annual mean rainfall also shows a downward trend. Over the period 1963–2013, rainfall patterns in the study area decreased with an increase in extreme events, as shown in Table 6.

3.2.2. Seasonal Rainfall Trend Analysis

Rainfall is divided into three seasons: pre-monsoon, post-monsoon, and monsoon season. Most rainfall occurs during the monsoon season. The pre-monsoon rainfall has been increasing historically from 1970–1993. There is a slight shift in the rainfall pattern toward the pre-monsoon season, as illustrated in Figure 4d. The rainfall during the monsoon season is relatively stable in comparison to the pre-and post-monsoon seasons. An analysis of the monsoon seasons shows that the pattern of rainfall is significantly shifting towards the months of June, August, and September, although this trend is not statistically significant.

3.2.3. Frequency Analysis

Extreme rainfall events are the cause of flood events. Life loss and economic damages are their adverse consequences. To better understand this phenomenon, extreme frequency of occurrence and probabilities are calculated. The return period of storm rainfall helps to understand future rainfall extremes better and to predict/simulate the runoff produced by those events. Using the Gumble Extreme Value Analysis Type 1 Distribution, the frequency of rainfall at Barkhan rain station has been analyzed, and precipitation has been calculated over four return periods.

3.3. Simulation Using HEC-HMS

3.3.1. Peak Flow Simulation

HEC-HMS simulates extreme events while taking into account the characteristics of the watershed. Hydrological simulations in HEC-HMS can be conducted using different methods of loss and transformation of rainfall. This study uses SCS-CN as a loss and transforms method.

Table 7. Hill torrents showing their peak discharge in Cusecs values at the given return periods.

Hill Torrents Name	Return Period (Years)				ADB Study (2007)	Maximum (Cusec)
	25 (116) *	50 (143)	100 (151)	200 (170)
Kaha	136,311.8	167,637.4	197,866.8	232,213.0	…	96,000 (2010)
Sanghar	159,435.6	193,478.6	225,901.5	262,329.1	122,730	229,000 (2010)
Kaura	22,255.7	26,153.1	29,816.8	33,894.7	45,160	128,500 (2010)
Vehova	67,526.1	78,848.9	89,620.2	101,728.6	85,520	110,500 (2010)
Vidore	32,955.9	38,726.4	44,150.1	50,186.2	…	…
Mithawan	10,287.9	12,139.3	13,879.9	15,817.0	…	…
SakhiSarwar	12,456.7	14,540.9	16,501.0	18,684.2	…	…
Chachar	24,926.1	29,956.5	34,742.7	40,119.1	…	…
RakhiMounh	11,114.2	13,194.9	15,151.9	17,329.8	…	…
Pitok	26,415.3	30,995.3	35,301.8	40,096.5	…	5000 (2010)

* In parenthesis rainfall value at the relative return periods.

provides peak discharges at the outlets of each hill torrent based on the calculated rainfall to the adjacent return period. Peak discharges are compared with grey literature of the area and historical maximum discharge data collection of the irrigation department. For that purpose, the available report used prepared by Asian Development Bank (ADB) in 2017 [26] through consultant team. However, the report does not mention the modeling parameter details and used different years of data for hydrological modeling. The ungauged nature of these hill torrents and data scarcity is one of the limitations for the study area.

3.3.2. Volume Simulation

Local communities in DG Khan and Rajanpur divert runoff water during low flows into temporary Bunds; however, the volume of water required to construct a permanent storage structure still needs to be assessed. In terms of storage potentials on these hill torrents, the volume of runoff is computed at each return period as shown in

Table 8. Hill torrents showing the volume of runoff in Acre-Ft at the given return periods.

Hill Torrents Name	Return Period (Years)
	25 (116)	50 (143)	100 (151)	200 (170)
Kaha	94,355.4	118,423.3	141,984.8	169,072.3
Sanghar	121,706.1	152,614.4	182,792.0	217,404.6
Kaura	27,760.2	33,686.0	40,367.5	45,767.0
Vehova	99,975.5	111,204.8	121,887.2	133,895.6
Vidore	39,664.9	48,161.9	56,298.0	65,485.1
Mithawan	12,257.0	14,983.5	17,603.4	20,569.7
SakhiSarwar	16,391.0	19,652.9	22,758.6	26,250.3
Chachar	19,371.0	23,728.6	27,927.0	32,691.9
RakhiMounh	12,773.4	15,791.4	18,706.6	22,020.8
Pitok	31,108.8	37,573.5	43,750.7	50,714.6
Total (ACRE-FT)	475,363.3	575,820.3	674,075.8	78,3871.9
Total (MAF)	0.475363	0.57582	0.674076	0.783872

.

3.3.3. Sensitivity Analysis

The sensitivity of peak runoff was examined by changing two parameters: CN and Initial abstraction. The percentage change increased and decreased by 2, 5, and 10% used for sensitivity analysis.

Curve Number Sensitivity

The land cover of the area determines the curve number. The more impervious the area, the higher the curve number. Initially, peak discharges were calculated using the values assigned by CN-Grid to each sub-basin of the watershed in HEC-GeoHMS. Afterward, increasing and decreasing changes are made while maintaining other parameters, and discharges are simulated. Changes in curve numbers have a significant impact on discharges. It was found that the upper and lower bounds of 10% change in CN were 7.8% change and 13% change in simulated discharges, respectively.

Initial Abstraction

Initially, rainwater is lost due to infiltration, which is hydrologically referred to as initial abstraction. The initial abstraction depends on the basin’s CN and the soil type. In the initial abstraction, the same changes of 2, 5, and 10% were made by keeping other parameters unchanged.

Table 9. Sensitivity analysis Based on CN, and Initial Abstraction.

	Curve No	Initial Abstraction (Ia)
Increase	% Change	% Change
2%	2.081873	−0.0182
5%	4.8061197	−0.0489
10%	7.7970864	−0.0982
Decrease
2%	−2.323105	0.021581
5%	−6.124852	0.051673
10%	−13.20128	0.099394

illustrates that there was no significant impact on runoff generation.

3.4. Storage Availability Analysis

Three torrents were selected to determine the storage potential of the three hill torrents. ArcScene creates depth versus area tables and area versus volume tables. As shown in Figure 5, Kaha, Vehova, and Sanghar have estimated potentials of 0.14, 1.14, and 1.13 MAF, respectively. Numerous studies on how to best use flash flood water for irrigation and lessen the effects of drought on hill torrents came to the conclusion that building storage reservoirs, diversion structures, and cross-drainage systems is crucial [2]. Growers in hill torrent locations exploit low torrent flows for traditional agriculture by constructing little embankments known as “Gandaz” locally. Flash floods have been observed to destroy minor embankments, preventing farmers from using this priceless water. Punjab supports 80% of Pakistan’s agricultural output [27]; hence it plays a significant part in the country’s economic life. Punjab produces more than 90% of its food through irrigated land [28], utilizing a sizable portion of its traditional land and water resources. Today’s population growth and diminishing land resources necessitate the use of wastelands that can only be irrigated by floodwater [29]. According to the degree and frequency of flooding, it has been observed that cropping intensity varies greatly in the piedmont area of the hill torrents. On most of the cropped area during flood wetness, sorghum and millets are seeded, and then oil seeds are grown.

4. Policy Discussion: Socio-Economic Consequences of Water Infrastructure

Water accounting of Pakistan’s total water resources is somewhat uncertain due to limited data collection (especially for internally generated water resources from the hill torrents of Balochistan and parts of Sindh and Punjab) and the need for real-time water accounting instrumentation required [30]. National Water Policy (NWP) and National Water Conservation Strategies also discuss the potential of these untapped water resources and asserted the profiling of untapped water resources (hill torrents) for efficient management [31]. There is a need to plan the construction of water infrastructure. The importance of these infrastructures is well understood in the policy circle, and it is well recognized by Khyber Pakhtunkhwa irrigation official Mr. Zubair in a Pakistan Water Week organized by IWMI from 24–28 October 2022. Mr. Zubair pointed out “Gomal Zam Dam has covered its cost by over 10 times. It has saved the communities of Tank and DI Khan from the impacts of August 2022 floods. The dam has stored 1.1-million-acre ft water in these floods” [32]. Abdul Wahab Kakar, Director General on Farm Water Management (OFWM) Balochistan, also endorsed Mr. Zubair’s point of view and said, “In Balochistan, floods were caused by hill torrents. Such surface water could be used in the rainfed areas for productive use, and capacity building of farmers also needs to go alongside” [33].

There is no doubt about the economic and flood protection benefits, but the development of these infrastructures also reshaped the social structure of society and groundwater hydrology contour. The reshaping of land and water resources create new beneficiaries and losers. This is clearly a case for politics of scale.

The politics of scale in water resources refers to the ways in which different levels of government, organizations, and individuals interact and make decisions about water management [34,35]. This can include issues such as allocation of water resources, development of water infrastructure, and implementation of regulations and policies. The politics of scale can also involve conflicts and negotiations between different stakeholders with different interests, such as agricultural and urban users, as well as local and regional governments. It also deals with how different levels of government and other actors work together to address complex water management issues that often extend beyond political boundaries [34]. In the case of hill torrent management in Punjab, Pakistan, the politics of scale is particularly relevant due to the complex and multi-faceted nature of the issue. In Punjab, the provincial government is responsible for overall flood management, including the development of infrastructure and policies to mitigate the impact of hill torrents. However, local governments, such as district and city councils, also play a key role in responding to and recovering from these events. Additionally, public sector organizations, such as irrigation have a vested interest in managing hill torrents as they can affect their operations and infrastructure. Furthermore, the local communities who are most affected by hill torrents, also have a significant role to play in the management of these events. They have the knowledge of the area and their participation is crucial in the decision-making process. However, the involvement of local communities in the decision-making process is not always guaranteed and it can be a source of conflicts between different stakeholders [36]. Overall, managing hill torrents in Punjab, Pakistan requires coordination and cooperation between different levels of government, private sector organizations, and local communities. The politics of scale is therefore an important consideration in this context, as it can affect the ability of different actors to work together effectively and make decisions that are in the best interest of all stakeholders. In the case of surface storage for hill torrent management, several types of water conflicts may arise. For example, the following conflicts may arise:

• Allocation conflicts: Surface storage projects such as dams and reservoirs can affect the distribution and availability of water, leading to conflicts between different water users, such as irrigation command area distribution between upstream and downstream communities.
• Environmental conflicts: Surface storage projects can have significant impacts on the natural environment, such as altering flow patterns, flooding areas, and disrupting aquatic ecosystems. This can lead to conflicts between those who support the projects and those who are opposed to them due to environmental concerns.
• Decision-making conflicts: Local communities who are most affected by hill torrents and surface storage projects, also have a significant role to play in the management of these events. However, their participation in the decision-making process is not always guaranteed and it can be a source of conflicts between different stakeholders.

There is a need to rethink the past experiences of irrigation management during a colonial-era canal infrastructure and plan future water resource development accordingly. Building these hill torrents infrastructures and its command area development need to understand the existing customary water use practices and crop choices for its use. How to update these customary practices (Governance) for improving this precious water resource use. Flood irrigation is traditionally used to redirect the flow of a hill torrent into areas that need to be irrigated to cultivate seasonal crops. The agricultural operations are designed to withstand heavy floods and droughts, and a special irrigation system for hill torrent areas is used locally called “Kamara Irrigation.” This method disregards the frequency and volume of flows produced by storms and enforces successive water rights, dictating irrigation patterns from higher to lower riparian areas. This strategy prevented far lower riparians from receiving irrigation water during a year of low flow. People have been observed to prepare fields by building around 1.8 m high embankments to hold the water while considering the local soil type, water availability, and several other considerations. As soon as the field’s water supply stopped, crops were sowed and thrived because of the moisture the soil had retained. Other than the rain, no more watering is conceivable if it comes. Discussion on points like; what are the ways to intervene in these customary water laws/rules required, how their agroecological potential tapped but not repeating the colonial legacy of canal colonies, and what are ways these water resources provide an economical and just economic dividend to its inhabitants keeping in view the loser and winner of this new irrigation development.

Based on the previous experience of the Indus Basin irrigated canal development crop choices need to fixed—based on local climate and soil topography, so that the water logging and salinity menace is avoided. In such areas, it is important to select crops that are more drought-tolerant and require less water to grow. Some examples of drought-tolerant crops that may be more suitable for areas with limited water resources include:

• Millets: Millets such as pearl millet (Bajara) and finger millet (Ragi) are drought-tolerant crops that can be grown in areas with low rainfall and poor soil fertility.
• Pulses: Pulses such as lentils and chickpeas are drought-tolerant crops that can be grown in areas with low rainfall and poor soil fertility. They are also good source of protein.
• Oilseeds: Oilseeds such as groundnut and sunflower are drought-tolerant crops that can be grown in areas with low rainfall and poor soil fertility. They are also good source of oil and protein.
• Fodder crops: Fodder crops such as oat and barley can be grown as a source of feed for livestock and require less water than other crops.
• Agroforestry: Agroforestry systems, that integrate trees and crops on the same land, can be a good choice for command area development in these areas. The trees can provide shade, prevent soil erosion, and improve soil fertility while the crops can provide food and income.

It is important to note that the crop selection will also depend on the specific conditions of the area and the preference of the farmers. In addition, the farmers will have to be trained on different crop management practices, irrigation techniques, and post-harvest management to ensure the successful cultivation of the crops. Additionally, the water management strategy, such as rainwater harvesting, water conservation and efficient irrigation systems can also be implemented to make the best use of the available water resources.

5. Conclusions and Recommendations

This study aimed to estimate the contribution of hill torrents to the river Indus flow. The hydrological modeling analysis revealed that for return periods of 25, 50, 100, and 200 years, the total volume of flow that contributes to the river Indus is approximately 0.5, 0.6, 0.7, and 0.8 million acre-feet, respectively. A sensitivity analysis indicates that the CN and initial abstraction have a significant impact on runoff. Changes in CN by 10% with upper and lower bounds may affect runoff generation by 7.79% and −13.20%, respectively. The storage potential of three large basin hill torrents, Kaha, Vehova, and Sanghar, is 0.14, 1.14, and 1.13 MAF, respectively.

The study indicates that a substantial amount of water is available for spate irrigation development in the region and hydraulic structures need to be built to make the most use of this resource. To examine the effectiveness of hydraulic structures, a sediment transport analysis must be conducted. Furthermore, to obtain a more precise understanding of rainfall patterns, a network of rain gauges should be established.