Hydrological Modelling for Water Resource Management in a Semi-Arid Mountainous Region Using the Soil and Water Assessment Tool: A Case Study in Northern Afghanistan

To address the issues of water shortages and the loss of agricultural products at harvest in northern Afghanistan, the Soil and Water Assessment Tool (SWAT) was applied for agricultural water resource management by simulating surface runoff in the Balkhab River basin (BRB) on a monthly basis from 2013 to 2018. Elevation, slope, land cover data, soil maps, and climate data such as temperature, precipitation, relative humidity, wind speed, and solar radiation were used as inputs in the SWAT modelling. During the dry season from July to September, the water resources downstream were basically attributed to baseflow from groundwater. In the calibration, the groundwater baseflow was estimated by analyzing station-recorded discharges for 190 springs. With the estimated baseflow, the SWAT results were markedly improved, with R2 values of 0.70, 0.86, 0.67, and 0.80, Nash-Sutcliff efficiency (NSE) values of 0.52, 0.83, 0.40, and 0.57, and percent bias (PBIAS) values of 23.4, −8.5, 23.4, and 17.5 in the four different subbasins. In the validation, the statistics also indicated satisfactory results. The output of this study can be used in agricultural water resource management with irrigation practices and further in the assessment of climate change effects on the water resources in the BRB.


Introduction
Water resource management has had an inevitable link to economic development and human activities throughout history. Every continent has experienced water scarcity, and global water use has increased two times faster than the population increase in the last century [1]. Water scarcity, especially in arid and semiarid regions, highlights the importance of water resource management worldwide. Mekonnen andand Hoekstra (2016) [2] reported that 0.50 billion people faced severe water scarcity for the entire year from 1996 to 2005, while 3.97 billion people face severe water scarcity for at least one month each year. Approximately 80% of the population in Central Asia is experiencing difficulty associated with water stress, and approximately 50% of the population suffers from water shortages [3]. Studies report an increase in drought severity in Afghanistan and its neighboring countries in Central Asia [4][5][6][7]. Water resource management during severe drought, water scarcity, and shortages is essential to reduce the vulnerability of local people.
Afghanistan, with an arid to semiarid climate, receives most precipitation in the form of snow during winter and a smaller amount in spring as rainfall [8][9][10]. According model from 2003 to 2010, and model calibration and validation were conducted at a single gauging station. In Hajihosseini et al. (2016) [22], the SWAT model was utilized in the Helmand River basin from 1969 to 1979, and the model output was used to assess the 1973 treaty of the transboundary water of the Helmand River in southern Afghanistan. The Salma Dam located on the Harirod River is a newly constructed multipurpose dam in western Afghanistan. The SWAT model was used to estimate the surface runoff and sediment downstream from the Salma Dam [23]. The snowmelt runoff model (SRM) was used to estimate the daily river discharge in the Salang River basin and upper Kabul River basin by utilizing the MODIS 7-day snow cover product and metrological data from 2009 to 2011 [24,25]. Few studies have been conducted in the Afghanistan and no studies using the hydrological in the northern Afghanistan has been conducted to date. On the other hand, previous studies performed calibrations using single station discharge data for the calibration and validation of models while the current study used multiple stations, from upstream to downstream, for model verification in the study region. The result of the current study will be used for irrigation water management and as a solution for the current challenges.
The main objective for this study is to develop the hydrological model of the BRB for the sustainable water resource management practices in irrigation water management. The results of current study will be used for the evaluation of current irrigation water shortages and challenges, as well as for testing the possible policies and solutions for the proper distribution of water resources in future studies. In this study, a watershed model of the Balkhab River basin (BRB) in northern Afghanistan was developed using the SWAT, and model calibration was performed using multiple site calibration. In the mountainous watershed with arid and semiarid climates, lateral flows and springs are the main sources of surface water during the non-rainynon-rainy season. In this study, the authors used analysis of the historical discharge data at the different stations to address the issue of 190 spring discharges to the streamflow where no record of their discharge is available. In this study, the hydrological process of surface runoff within the study area was simulated using SWAT. The results of the current study can provide the core component for future work on irrigation practices and water resource management in watersheds by testing different possible scenarios. The results can also be used for studies on the impact of climate change on the water resources in northern Afghanistan.

Study Area
Afghanistan, with a semiarid to arid climate, has been divided into five major basins: The Northern, Kabul, Helmand, Harirod Murghab, and Amu Darya river basins. The study area is situated in the Northern River basin (Figure 1b), which is further divided into four river basins, namely, the Balkhab, Khulm, Sari Pul, and Shirin Tagab river basins. This study focused on the BRB with a minimum elevation of 202 m, maximum elevation of 4616 m, and an area of 28,835.2 km 2 ( Figure 1c). The Balkhab River starts from a series of six lakes in the high mountains in the south and flows through the north, and 190 springs join the stream along this river. The largest portion of water in the river is used for the irrigation of 296,774 ha of agricultural land through 3 small regulator dams and 101 traditional and engineering irrigation canals along the river. The average annual precipitation in Afghanistan is 270 mm and it ranges from 1270 mm in the high mountains to 75 mm in the southwestern area [26]. The BRB climate is semiarid with an average annual precipitation of 247 mm. Notably, there is no precipitation during the summer season. The water source for the rest of the year is the snowmelt water coming from the snow accumulated in the highland area during the winter and groundwater flow. In the northern water basin, glaciers originated in the Hindu Kush mountains in northwestern Afghanistan, and there are no glaciers located in the BRB region [27]. Snowmelt and underground water significantly contribute to the stream flow in the late spring, summer, and early autumn seasons, when precipitation is almost zero. Five hydrometeorological stations ( Figure 1c) are located in the BRB, recording the minimum and maximum temperatures, relative humidity, precipitation, and stream discharge with a daily temporal resolution and monthly average. In the current study, the daily metrological data were used for simulations of surface runoff on a monthly basis. The annual average temperature varies from 7.56 • C in the highlands of the south to 18.93 • C in the flat lands of the north based on the data from 2010 to 2018 provided by the Ministry of Energy and Water, Islamic Republic of Afghanistan. The BRB dominated land cover is dominated by range land and grasses, encompassing 62.25% of the total area. water basin, glaciers originated in the Hindu Kush mountains in northwestern Afghanistan, and there are no glaciers located in the BRB region [27]. Snowmelt and underground water significantly contribute to the stream flow in the late spring, summer, and early autumn seasons, when precipitation is almost zero. Five hydrometeorological stations ( Figure 1c) are located in the BRB, recording the minimum and maximum temperatures, relative humidity, precipitation, and stream discharge with a daily temporal resolution and monthly average. In the current study, the daily metrological data were used for simulations of surface runoff on a monthly basis. The annual average temperature varies from 7.56 °C in the highlands of the south to 18.93 °C in the flat lands of the north based on the data from 2010 to 2018 provided by the Ministry of Energy and Water, Islamic Republic of Afghanistan. The BRB dominated land cover is dominated by range land and grasses, encompassing 62.25% of the total area.

Dataset
The SWAT inputs are temperature, precipitation, relative humidity, wind speed, solar radiation, a digital elevation model (DEM), land cover map, soil-type map, slope map, and water input and outflow from the river.
The recorded daily maximum and minimum temperatures, relative humidity, and precipitation from five hydrometeorological stations from 2010 to 2018 were used as the weather input data for the SWAT model setup on a monthly time scale (Figure 1c). The wind speed and solar radiation were retrieved from the NASA Prediction of Worldwide Energy Resources (RWER) at the same geographical location as those five stations since

Dataset
The SWAT inputs are temperature, precipitation, relative humidity, wind speed, solar radiation, a digital elevation model (DEM), land cover map, soil-type map, slope map, and water input and outflow from the river.
The recorded daily maximum and minimum temperatures, relative humidity, and precipitation from five hydrometeorological stations from 2010 to 2018 were used as the weather input data for the SWAT model setup on a monthly time scale (Figure 1c). The wind speed and solar radiation were retrieved from the NASA Prediction of Worldwide Energy Resources (RWER) at the same geographical location as those five stations since the stations did not collect these data. Monthly recorded discharge is available for all five stations, and four stations located on the mainstream were used for calibration and validation purposes. A summary of the climate data and river discharge at the five hydrometeorological stations is presented in Figure A1a-e in the Appendix A for the ground observations and data retrieved from global sources.
The ALOS PALSAR Radiometric terrain corrected (RT1) images ( Figure 1c) were downloaded from the Alaska Satellite Facility Distributed Archive Data Center (ASF DAAC) website with a spatial resolution of 12.5 m and used as the DEM in the SWAT watershed model. RT1 products are generated from high-resolution and mid-resolution DEMs [28]. The slope map for the BRB was created using this DEM in the SWAT model (Figure 2c). In the Figure 2a, the acronyms are defined as follows: water (WATR), residential-medium density (URMD), agricultural land-close-grown (AGRC), range-grasses (RNGE), forestdeciduous (FRSD), barren (BARR), agricultural land-generic (AGRR), and barren sand (SWRN). In Figure 2b the acronyms are based on the FAO/UNESCO soil type name designations and are defined as follows: Lithosols and Xerosols, strongly dissected to mountainous (I-X-c-3512); Lithosols, Cambisols and Rankers, strongly dissected to mountainous (I-B-U-2c-3503); Calcic Yermosols, rolling to hilly (Xk4-2b-3303); Haplic Yermosols, Calcaric Regosols, Calcic Yermosols, and Orthic Solonchaks, level to hilly (Yh23-2ab-3581); and Cambic Arenosols, Calcaric Regosols, Calcic Yermosols, and Solonchaks, level to gentle undulation (Qc47-1a-3542). the stations did not collect these data. Monthly recorded discharge is available for all five stations, and four stations located on the mainstream were used for calibration and validation purposes. A summary of the climate data and river discharge at the five hydrometeorological stations is presented in Figure A1a-e in the Appendix A for the ground observations and data retrieved from global sources. The ALOS PALSAR Radiometric terrain corrected (RT1) images ( Figure 1c) were downloaded from the Alaska Satellite Facility Distributed Archive Data Center (ASF DAAC) website with a spatial resolution of 12.5 m and used as the DEM in the SWAT watershed model. RT1 products are generated from high-resolution and mid-resolution DEMs [28]. The slope map for the BRB was created using this DEM in the SWAT model ( Figure 2c). In the Figure 2a, the acronyms are defined as follows: water (WATR), residential-medium density (URMD), agricultural land-close-grown (AGRC), range-grasses (RNGE), forest-deciduous (FRSD), barren (BARR), agricultural land-generic (AGRR), and barren sand (SWRN). In Figure 2b the acronyms are based on the FAO/UNESCO soil type name designations and are defined as follows: Lithosols and Xerosols, strongly dissected to mountainous (I-X-c-3512); Lithosols, Cambisols and Rankers, strongly dissected to mountainous (I-B-U-2c-3503); Calcic Yermosols, rolling to hilly (Xk4-2b-3303); Haplic Yermosols, Calcaric Regosols, Calcic Yermosols, and Orthic Solonchaks, level to hilly (Yh23-2ab-3581); and Cambic Arenosols, Calcaric Regosols, Calcic Yermosols, and Solonchaks, level to gentle undulation (Qc47-1a-3542).  The Land Cover Atlas of the Islamic Republic of Afghanistan for 2010 was prepared through the cooperation of the FAO and Ministry of Agriculture, Irrigation, and Livestock (MAIL), Islamic Republic of Afghanistan, under the "Strengthening Agricultural Economic, Market Information and Statistics Services", using medium resolution satellite imagery from SPOT-4 with a 10 m spatial resolution and the Global Land Survey (GLS) Landsat Thematic Mapper, high-resolution satellite imagery and aerial photographs as well as ancillary data [29]. The land cover map for the BRB is shown in Figure 2a. There is a lack of reliable soil-type maps in many parts of the world, for which two commonly used global soil-type maps, the FAO/UNESCO soil map of the world and the Harmonized World Soil Database (HWSD_v121), are adopted [30]. The FAO/UNESCO soil map of the world (Figure 2b) was used as the soil-type map in this watershed modelling. The FAO/UNESCO soil map of the world was prepared using topographic map series of the American Geographical Society of New York at a scale of 1:5,000,000 with a topsoil layer of 30 cm and a subsoil layer of 70 cm [30].

Point Source Data
The BRB is located in a mountainous region with 190 springs along the steam, where no records of their exact location and water yield are available. Historical discharge during the non-rainy season was used to overcome the groundwater data scarcity issue. We assumed that the precipitation and snowmelt contributions to the stream flow during July, August, and September is negligible. The monthly average precipitation at the five stations in the BRB is shown in Figure 3, where the precipitation is almost zero during July, August, and September. In the same bar plot, the snow cover area percentage was plotted as per Hussainzada, (2020) [31]. In this study, monthly snow cover percentages for the BRB were retrieved from the spatiotemporally combined MODIS Aqua and Terra daily snow cover products. The recorded monthly average discharge of the four stations Rabat-i-Bala, Pul-i-Baraq, Doshqadam, and Nazdik-i-Nayak stations in subbasins 7, 10, 28, and 30 ( Figure 1c) respectively, is shown in Figure 3 for these three dry months. The dotted lines in Figure 3 show the average flow from 2010 to 2018. The Balkhab river starts from a series of six lakes with a perennial discharge of 4.22 m 3 /s, and for subbasins 7, 10, and 28, the groundwater contribution is equal to the difference in the amount of discharge between the upper station and lower station during the non-rainy months. The groundwater contribution in the subbasins is replicated in Table 1. well as ancillary data [29]. The land cover map for the BRB is shown in Figure 2a. There is a lack of reliable soil-type maps in many parts of the world, for which two commonly used global soil-type maps, the FAO/UNESCO soil map of the world and the Harmonized World Soil Database (HWSD_v121), are adopted [30]. The FAO/UNESCO soil map of the world (Figure 2b) was used as the soil-type map in this watershed modelling. The FAO/UNESCO soil map of the world was prepared using topographic map series of the American Geographical Society of New York at a scale of 1:5,000,000 with a topsoil layer of 30 cm and a subsoil layer of 70 cm [30].

Point Source Data
The BRB is located in a mountainous region with 190 springs along the steam, where no records of their exact location and water yield are available. Historical discharge during the non-rainy season was used to overcome the groundwater data scarcity issue. We assumed that the precipitation and snowmelt contributions to the stream flow during July, August, and September is negligible. The monthly average precipitation at the five stations in the BRB is shown in Figure 3, where the precipitation is almost zero during July, August, and September. In the same bar plot, the snow cover area percentage was plotted as per Hussainzada, (2020) [31]. In this study, monthly snow cover percentages for the BRB were retrieved from the spatiotemporally combined MODIS Aqua and Terra daily snow cover products. The recorded monthly average discharge of the four stations Rabati-Bala, Pul-i-Baraq, Doshqadam, and Nazdik-i-Nayak stations in subbasins 7, 10, 28, and 30 ( Figure 1c) respectively, is shown in Figure 3 for these three dry months. The dotted lines in Figure 3 show the average flow from 2010 to 2018. The Balkhab river starts from a series of six lakes with a perennial discharge of 4.22 m 3 /s, and for subbasins 7, 10, and 28, the groundwater contribution is equal to the difference in the amount of discharge between the upper station and lower station during the non-rainy months. The groundwater contribution in the subbasins is replicated in Table 1.   Table 1. The recorded average discharge at the stations during the non-rainy months and an estimation of the groundwater contribution based on the difference between records in the non-rainy seasons.

Irrigation Canals
Information on irrigation canals is used as the seasonal irrigation practice for crop cultivation in watershed modelling. The agricultural sector in the BRB extracts water through 101 irrigation canals and 3 regulator dams in the downstream region, extracting water from the river. Three regulator dam functions distribute the water to the irrigation canals, and there is no storage behind the dam. Therefore, the dam downstream discharge is equal to the dam upstream discharge minus the water withdrawn by irrigation canals. Data on the locations of the irrigation canal intakes and their discharge are provided by the Northern River Basin Authority Report [32]. Figure 4 depicts the irrigation canal intake along the Balkhab River. Thirty-four canals out of hundred and one canal in the BRB originate from the main Balkhab River and are designated by the green color in Figure 4; these canals irrigate 258,305.5 ha of agricultural land, and other canals are small canals that feed from springs or small temporary streams mostly located upstream in the southern part of the study area. The sources of water for the other sixty-seven canals designated by the red color in Figure 4 are springs or temporary small rivers. In this study, thirty-four canals and their effect on river discharge were considered negligible.

Irrigation Canals
Information on irrigation canals is used as the seasonal irrigation practice for crop cultivation in watershed modelling. The agricultural sector in the BRB extracts water through 101 irrigation canals and 3 regulator dams in the downstream region, extracting water from the river. Three regulator dam functions distribute the water to the irrigation canals, and there is no storage behind the dam. Therefore, the dam downstream discharge is equal to the dam upstream discharge minus the water withdrawn by irrigation canals. Data on the locations of the irrigation canal intakes and their discharge are provided by the Northern River Basin Authority Report [32]. Figure 4 depicts the irrigation canal intake along the Balkhab River. Thirty-four canals out of hundred and one canal in the BRB originate from the main Balkhab River and are designated by the green color in Figure 4; these canals irrigate 258,305.5 ha of agricultural land, and other canals are small canals that feed from springs or small temporary streams mostly located upstream in the southern part of the study area. The sources of water for the other sixty-seven canals designated by the red color in Figure 4 are springs or temporary small rivers. In this study, thirty-four canals and their effect on river discharge were considered negligible.

SWAT Model
The SWAT model is a physical-based semi-distributed hydrological model developed to predict the impact of water management, sediment, and agricultural chemical yield in a large watershed at a daily time step [16,33,34]. Based on the daily precipitation, air temperature, relative humidity, solar radiation, and wind speed, the SWAT simulates the daily, monthly, and yearly hydrological processes in a watershed.
The SWAT estimates the surface runoff based on two different methods: (1) the soil conservation service (SCS) runoff curve number and (2) the Green-Ampt infiltration method [34]. In this study, the SCS curve number method was used for modelling purposes. The SWAT divides the basin into subbasins that are composed of smaller components called hydrological response units (HRUs). HRUs are divided based on homogenous land use, soil type, and topography. First, the surface runoff was estimated for the HRUs and then the total runoff for the entire watershed was obtained using the hydrological model [16,35]. The model predicts the surface runoff based on the water balance equation (Equation (1)).
where SW t and SW are the respective final and initial soil water contents at time t; R, Q, ET, P, and QR are the daily amounts of precipitation, runoff, evapotranspiration, percolation, and return flow, respectively; all units are mm of H 2 O. The model estimates potential evapotranspiration (PET) based on three functions: Hargreaves, Penman-Monteith, and Priestley-Taylor. The Hargreaves method estimates PET as a function of extraterrestrial radiation and air temperature [16]. In this study, the Hargreaves method was utilized for PET calculations by estimating the PET as a function of extraterrestrial radiation and air temperature [16]. The Hargreaves temperature-based method is the most commonly used method for the estimation of PET when weather data are unavailable [36]. The Hargreaves method has been modified for SWAT by replacing the extraterrestrial radiation with RAMX (maximum possible solar radiation at the Earth's surface) and through the use of temperature exponents from 0.5 to 0.6 [16].
where T mx , T mn , and T are the daily maximum, minimum, and average air temperatures in • C, respectively, and HV is the lateral heat of vaporization (J g −1 ). Snow present in the watershed will melt when the second soil layer temperature exceeds 0 • C, and snowmelt will account for the same amount of rainfall during the day to estimate runoff and percolation [16]. The model uses the degree-day method for estimation as a linear function of air temperature by setting the snowmelt threshold [34].

Model Setup
ArcSWAT 2012 with the ArcGIS interface was used for BRB watershed model development. In the first step of the modelling, watershed delineation using ALOS PALSAR DEM with 12.5 m spatial resolution for simulating the streamflow and subbasins was conducted. In the watershed delineation step, SWAT divided the BRB into 31 subbasins based on the flow accumulation and DEM file (Figure 4). The watershed delineation does not cover all low land areas in the northern BRB because those area slopes are far away from the Balkhab River, and their surface runoff does not contribute to the stream flow of the Balkhab River. Second, the HRU definition was manipulated using the FAO/UNESCO soil map of the world, land cover map of the study site, and slope of the site using the DEM (Figure 2). The study area was divided into 1532 HRUs for all 31 subbasins based on the input dataset. Then, the weather records at five hydrometeorological stations from 2010 to 2018 were defined as the input weather data for the model. The SWAT model was run for the BRB with a three-year warm-up period from 2010 to 2012, and outputs were generated from 2013 to 2018 on a monthly basis. In the current study, the calendar year was adopted for model simulation starting in January and ending in December because the data are available from January 2010, even though the hydrological years in the northern hemisphere start in October and end in September. The model simulation output was divided into two periods for calibration and validation. The simulation results from 2013 to 2015 were used for calibration purposes, and the simulation results from 2016 to September 2018 were selected as the validation period. The major working procedure for the model setup and model calibration is summarized in the workflow diagram ( Figure 5). to 2018 were defined as the input weather data for the model. The SWAT model was run for the BRB with a three-year warm-up period from 2010 to 2012, and outputs were generated from 2013 to 2018 on a monthly basis. In the current study, the calendar year was adopted for model simulation starting in January and ending in December because the data are available from January 2010, even though the hydrological years in the northern hemisphere start in October and end in September. The model simulation output was divided into two periods for calibration and validation. The simulation results from 2013 to 2015 were used for calibration purposes, and the simulation results from 2016 to September 2018 were selected as the validation period. The major working procedure for the model setup and model calibration is summarized in the workflow diagram ( Figure 5). The BRB starts from six lakes upstream with a perennial flow from the reservoir of these lakes, while 190 springs along the river also contribute to the baseflow of the river, especially during the dry season with almost zero precipitation. On the other hand, the irrigation canals along the river withdraw a considerable amount of water during the entire year from the Balkhab River. All these water contributions and water withdrawals for thirty-four canals were included in the model as point sources to develop a reliable hydrological model for the BRB.
The model was calibrated for the period from January 2013 to December 2015 based on the monthly discharge recorded at four stations, Rabat-i-Bala, Pul-i-Baraq, Doshqadam, and Nazdik-i-Nayak, which are located in subbasins 7, 10, 28, and 30, respectively. Then, the model was validated for the period from January 2016 to September 2018 with the remaining data. The BRB starts from six lakes upstream with a perennial flow from the reservoir of these lakes, while 190 springs along the river also contribute to the baseflow of the river, especially during the dry season with almost zero precipitation. On the other hand, the irrigation canals along the river withdraw a considerable amount of water during the entire year from the Balkhab River. All these water contributions and water withdrawals for thirty-four canals were included in the model as point sources to develop a reliable hydrological model for the BRB.
The model was calibrated for the period from January 2013 to December 2015 based on the monthly discharge recorded at four stations, Rabat-i-Bala, Pul-i-Baraq, Doshqadam, and Nazdik-i-Nayak, which are located in subbasins 7, 10, 28, and 30, respectively. Then, the model was validated for the period from January 2016 to September 2018 with the remaining data.
The parameter optimization was conducted using SWAT-CUP (SUFI2) at multiple sites for the parameters listed in Table 2. The monthly average discharges at four stations in subbasins 7, 10, 28, and 30 were included in the calibration process. In the model parameterization, a careful combination of the parameters should be considered before starting the calibration. The rainfall parameters as driving parameters or snow parameters (SFTMP, SMTMP, SMFMX, SMFMN, TIMP), which introduce water into the system, should not be calibrated with other parameters [40]. Non-uniqueness drives more than one solution for the SWAT model, which means that the model can generate the same discharge signal with two completely different sets of parameters, as described in Karim C. Abbaspour et al., (2017) [40]. The non-uniqueness issue can be reduced by increasing the number of variables and using multiple outlets for calibration purposes. An increase in the number of observation stations during the calibration can lead to a more reliable model after the calibration process. However, the use of more variables (i.e., discharge, ET, sediments, or crop yield) in the calibration process can ensure a better model representation of the original site condition. In this study, we could not access datasets other than discharge, and other datasets must be included in future studies.
The variability of the temperature and precipitation with elevation in the mountainous basin is an important phenomenon that can affect the model performance. Five elevation bands were used to adjust the temperature and precipitation variation to the elevation in each subbasin. Then, the temperature and precipitation lapse rate were calibrated, and the fitted values were replaced in the model. Snow and its contribution to the streamflow in mountainous basins is the second most important process and should be calibrated before the model parameters are tuned. The calibration process was continued with optimization of the seven snow parameters in the model. SWAT-CUP was run for the snow parameters with 200 simulations, and the fitted parameters were replaced in the model.
The model calibration requires knowledge of hydrology, while a correct understanding of the site and its hydrological condition is also necessary to choose the appropriate parameters. The model parameters were chosen based on a literature review of previous studies and the output of the calibrated model for the temperature lapse rate, precipitation lapse rate, and snow parameters. The chosen parameters were subjected to a sensitivity test using the SWAT-CUP platform to determine their sensitivity and changes on the hydrograph. Then, seven more parameters in Table 2 were used for the final model calibration with 1000 simulations.
The model performance for the calibration and validation was computed separately for the Rabat-i-Bala, Pul-i-Baraq, Doshqadam, and Nazdik-i-Nayak stations in the subbasins 7, 10, 28, and 30, respectively. The simulated runoff from the SWAT with monthly observed discharge based on the data availability were used for model verification. The Nash-Sutcliffe efficiency (NSE) [41], coefficient of determination (R 2 ), r-factor [42], p-factor, and percent bias (PBIAS) [43] were used as statistical indicators to determine the model goodness-of-fit for both the calibration and validation periods and are shown in Equations (3)-(5), respectively. The p-factor (percentage of measured data bracketed by 95% prediction boundary), often referred to 95PPU is used to quantify the model uncertainty of the SWAT. The p-factor ranges from 0-1, and values closer to 1 show a higher model performance and efficiency. The r-factor is the average width of the 95PPU band divided by the standard deviation of the measured variable.
where Q i is the measured discharge (m 3 /s), P i is the simulated discharge (m 3 /s), Q is the average observed discharge (m 3 /s), P is the average simulated discharge (m 3 /s), n is the number of months, y M t i ,97.5% and y M t i ,2.5% are the upper and lower boundaries of the 95PPU; and σ obs is the standard deviation of the observed data. The model performance can be evaluated with the general criteria presented in Table 3 [44,45]. The relationship of groundwater contribution to stream flow is a poorly understood phenomenon in hydrology [47,48]. In the dry season, groundwater discharge is a principal source of water in the lowland valley [47][48][49]. Figure 6 illustrates the hydrographs of two different simulations with and without constant groundwater input to the model. As shown in Figure 6, the simulations after inclusion of the spring flow contribution, which were estimated from the historical discharge analysis in the gauging stations presented in Figure 3 and summarized in Table 1, show a significant improvement in the baseflow by eliminating underestimations. The hydrographs present a constant gap between the simulation without groundwater input to the observations and the simulation after ground-water contribution is included in the model. In this section, the simulation result before and after adding the ground water contribution using point source data explained in Section 2.2.2 is presented. tion 2.2.2 is presented. However, the model without groundwater consideration still captured the temporal variation and pattern of the observation, but the large underestimation makes the model unrealistic. In terms of the statistics for model performance, the model can be evaluated as satisfactory to very good with respect to the R 2 values, while the NSE and PBAIS values were unsatisfactory in subbasins 7, 10, and 28, where high underestimation decreased the NSE range between −0.71 and −0.15 and the PBAIS range between 56.1% and 69.7%. In subbasin 30, the model did not differ with and without groundwater input since the subbasin is located far upstream, where six lakes feed a perennial discharge to the streamflow. The model statistical summary without the groundwater contribution case is presented in the validation section.  However, the model without groundwater consideration still captured the temporal variation and pattern of the observation, but the large underestimation makes the model unrealistic. In terms of the statistics for model performance, the model can be evaluated as satisfactory to very good with respect to the R 2 values, while the NSE and PBAIS values were unsatisfactory in subbasins 7, 10, and 28, where high underestimation decreased the NSE range between −0.71 and −0.15 and the PBAIS range between 56.1% and 69.7%. In subbasin 30, the model did not differ with and without groundwater input since the subbasin is located far upstream, where six lakes feed a perennial discharge to the streamflow. The model statistical summary without the groundwater contribution case is presented in the validation section.

Calibration Result
The simulated discharges from the model calibration were compared with the ground truth at the four gauging stations located in subbasins 7, 10, 28, and 30. Figure 7 demonstrates the calibration result between the observed and simulated monthly discharges at the four gauging stations in addition to monthly precipitation records at each site. Uncalibrated hydrographs show that the peak discharges appeared earlier in March and April than the observed peak in May. The uncalibrated peaks corresponded well with the peaks of observed precipitation. Thus, this shifting can be caused by faster snowmelt or misclassifying snow as rainfall in the model. Calibration of the temperature lapse rate, precipitation lapse rate, and snow parameters can improve the resulting hydrograph by considering the snowmelt conditions in the model. mation of the peak flows.
The model uncertainty results also depicted acceptable output with p-factors of 0.89, 0.61, 0.71, and 0.56 for subbasins 7, 10, 28, and 30 and r-factors of 0.75, 1.27, 1.32, and 1.06, respectively. The P-factor and r-factor are the percentage of observations covered by the 95 PPU uncertainty range and the relative width of the 95% probability band, respectively [50]. The model performance can be evaluated satisfactorily to very well with respect to the R 2 values as per the ranges described in Table 3. The NSE indicator classified the result as unsatisfactory for subbasin 28 and satisfactory to very good for the remaining sites. In addition, the PBIAS values classified the results as from satisfactory to good for all stations.  The simulated discharges after calibration of those parameters showed good agreement with the observed discharges at the gauging stations. The model captured variations in the streamflow in all the subbasins during the steady-flow conditions in winter, summer, and autumn with an almost constant underestimation gap between the observed and simulated discharges in subbasin 7. However, the model could not reproduce the peak flows well during the spring in subbasins 7, 28, and 30 in 2014 and 2015, with underestimation of the peak flows.

Validation
The model uncertainty results also depicted acceptable output with p-factors of 0.89, 0.61, 0.71, and 0.56 for subbasins 7, 10, 28, and 30 and r-factors of 0.75, 1.27, 1.32, and 1.06, respectively. The P-factor and r-factor are the percentage of observations covered by the 95 PPU uncertainty range and the relative width of the 95% probability band, respectively [50]. The model performance can be evaluated satisfactorily to very well with respect to the R 2 values as per the ranges described in Table 3. The NSE indicator classified the result as unsatisfactory for subbasin 28 and satisfactory to very good for the remaining sites. In addition, the PBIAS values classified the results as from satisfactory to good for all stations.

Validation
Then, model validation was conducted using the available observations of the monthly average discharge from January 2016 to September 2018 in subbasins 7, 10, 28 and 30. The validation results are presented in Figure 8 for the subbasins. There are missing data for the discharge observations in subbasin 28 from October 2017 to September 2018 and in subbasin 30 from October 2016 to September 2017. In the validation process the model runs in a single simulation with the same calibrated parameters.
(R 2 ) exhibited satisfactory to good model performance for all four subbasins. The PBIAS showed that the model performance was satisfactory to very good. The NSE showed satisfactory results in the two downstream stations in subbasins 7 and 10 and relatively lower values for subbasins 28 and 30, with values of 0.47 and 0.43, respectively.
The comparison between the statistical indicators for the calibration and validation depicted a slight decrease in the validation with respect to the calibration except for subbasins 7 and 28, which had slightly higher NSE values with respect to their calibration period. These slight increases were achieved by the better simulation of the peak flows while minimizing the underestimation in the validation period.

Effects of Ground Water Contribution to the SWAT Model Result
Creating a reliable model that can present the watershed condition largely depends on the input dataset and correct definition of the hydrological process within the watershed. Proper understanding of the hydrological process and water management activity of a study site can lead to a more accurate and decisive hydrological model, thereby obtaining more reliable modelling results. In a high mountainous basin with steep terrain, the main source of water comes from snowmelt from high altitude areas and water infiltration, which will join the streamflow elsewhere in the lower downstream area in the form of lateral flow or springs. In developing countries with a lack of observations and The comparison between the statistical indicators for the calibration and validation depicted a slight decrease in the validation with respect to the calibration except for subbasins 7 and 28, which had slightly higher NSE values with respect to their calibration period. These slight increases were achieved by the better simulation of the peak flows while minimizing the underestimation in the validation period.

Effects of Ground Water Contribution to the SWAT Model Result
Creating a reliable model that can present the watershed condition largely depends on the input dataset and correct definition of the hydrological process within the watershed. Proper understanding of the hydrological process and water management activity of a study site can lead to a more accurate and decisive hydrological model, thereby obtaining more reliable modelling results. In a high mountainous basin with steep terrain, the main source of water comes from snowmelt from high altitude areas and water infiltration, which will join the streamflow elsewhere in the lower downstream area in the form of lateral flow or springs. In developing countries with a lack of observations and ground truth data, finding accurate and solid records of hydrological processes is difficult and sometimes impossible.
The SWAT divides underground storage into shallow and deep aquifers, while the water entering the deep aquifer is not accounted for in the water budget calculation since it would not get back to the water cycle [51]. In the study area, approximately 190 springs contribute to the baseflow. As described in   [51], the inclusion of the deep aquifer contribution to the baseflow can improve the model output.
Considerable amounts of precipitation occurred in the BRB during the winter and spring seasons, while the summer and autumn seasons were relatively dry. The hydrograph in Figure 3 illustrates the streamflow during the dry months in July, August, and September, where the precipitation is very low and the snow cover area in the mountains of the BRB almost disappears. In such an environment, the contribution of the deep aquifer to the baseflow is very important.
Current estimation using the recorded average baseflow during the dry season where the contribution of rainfall, snowmelt, and shallow aquifers reaches the minimum can improve the modelling result by minimizing the underestimation. Since a systematic underestimation in Figure 6 is noticeable for subbasins 7, 10, and 28, this constant underestimation of the baseflow can come from a constant source of water, which is the deep aquifer water returning downstream as the springs.

Peaks in the Hydrograph
The model predicted the baseflow variation well, but the peaks were not simulated as well as the baseflow. The underestimation of the peak in the two upstream subbasins 28 and 30 was more noticeable. In spring, the increase in rainfall and snowmelt causes a considerable increase in streamflow. The snowmelt process accelerated by the increase in air temperature and rainfall also contributed to rapid melting on rainy days. Rapid melting and rainfall together with steep slopes in the mountainous area upstream may cause flooding during early spring. On the other hand, flash floods in mountainous regions are another cause of the peak flow during rainy days at high altitudes when high-intensity rainfall occurs during a short period of time.
The SWAT model snowmelt algorithms are based on the degree-day factor. In this method, the snowmelt estimates are based on the increase in the temperature above the snowmelt threshold temperature. The contribution of rainfall to rapid snowmelt is not considered in the model; this issue caused underestimation of the peak flow during the early spring season when most floods were recorded at the BRB stations. On the other hand, the coarse weather and sparse stations could not capture the high precipitation in the mountainous parts in the upstream area and can result in the underestimation of peak flows.

Effect of Snowmelt: Comparison of the SWAT Snowmelt and Monthly MODIS Snow Cover Area
The SWAT calculates the snowmelt based on the degree-day factor method. Figure 10 presents a comparison between the snowmelt output from the SWAT model in the BRB and snow cover area (%) retrieved from the MODIS snow cover product C6. The snow cover area (%) is estimated based on the MODIS Aqua and Terra snow cover area after a spatial and temporal three step combination for reduction of the cloud obstruction. The snow cover area (%) data is adopted from the previous study in the BRB [31].

Conclusions
In the current study, the SWAT model successfully reconstructed the mon charges at the four stations along the BRB. The model statistics showed a satisfa

Conclusions
In the current study, the SWAT model successfully reconstructed the monthly discharges at the four stations along the BRB. The model statistics showed a satisfactory to very good result with respect to the coefficient of determination (R 2 ) and PBIAS, while the evaluation showed an unsatisfactory to very good result based on the NSE value. The hydrologic regime in mountains with arid and semiarid climates can change seasonally, and BRB surface runoff is dominated by three different sources of water during the year. The hydrological regime is dominated by an increase in precipitation in the late autumn and winter seasons. The surface runoff during the spring and early summer seasons is strongly affected by snowmelt. During May especially, a spring peak will appear due to rapid snowmelt in addition to the contribution of rainfall to the snow melting process. In late summer and early autumn, the river discharge is dominated by groundwater flow through the springs where the water infiltrates upstream and joins the river in the flat lands downstream. In such cases, SWAT is highly capable of capturing hydrological processes that vary temporarily and spatially. However, SWAT model calibration is challenging in mountainous catchments with low annual precipitation. In mountainous catchments, the conditions can change spatially and temporarily, and calibration for accurate parameters is necessary to apply the model, which can represent the actual site conditions.
Hydrologic models are valuable for water resource management and water-related practices in Afghanistan and other data-scarce countries. The results of this study can be useful in decision-making processes for short-and long-term water resource planning and management in the BRB. Studies can expand further to investigate the effects and the role of irrigation practices under different scenarios, not only in the BRB, but also in Afghanistan.