Impact of land surface model schemes in snow-dominated arid and semiarid watersheds using the WRF-hydro modeling systems

1.   Introduction

Water resources management has a direct link to economic development and human activities. The water scarcity issue has been exacerbated significantly over the past century, demanding urgent attention and a sustainable solution. The surge in population is the key driver of water scarcity, as the increase in population corresponds to a higher demand for water ^[1]. According to FAO (2020) ^[1], the annual per capita freshwater availability declined by more than 20% in the past two decades, with particular severity in Western Asia and Northern Africa. In these regions, the per capita freshwater availability barely reaches the annual average of 1000 m³, conventionally considered the threshold for severe water scarcity ^[1]. Arid and semiarid climate zones face heightened vulnerability to drought and water scarcity. Past studies have declared a noticeable increase in drought severity in Afghanistan, Central Asia, and Iran ^[2,3,4,5]. Considering the water crisis in the region, a comprehensive approach involving extensive scientific studies on water resources is required to address the issue.

Computer-based hydrological models are a simplified representation of real-world systems consisting of a series of contemporary equations and logical sets of operation ^[6]. Modeling is a common tool in many scientific fields in general and in water sciences in particular. Hydrological models have diverse applications, such as modeling existing catchments with available input and output data (operational flood prediction, water resources management, or extension of data array for flood design of water resources assessment), coupled hydrology and meteorology (global climate models), coupled hydrology and geochemistry (nutrients and acid rain), ungauged basins runoff estimation, and the prediction of the impact of changes (land use and land cover) ^[6]. The best model is the one that can represent results closest to reality with the fewest input variables and complexity ^[7]. Hydrological models are classified into (1) empirical models or metric models, (2) conceptual models, and (3) physics-based models based on process description ^[7]. As per Singh (2018) ^[8], physics-based models can overcome many defects of the other two types of models due to the physical interpretation of the process included in the model parameters.

Most weather and climate models have adopted a one-dimensional (vertical) approach, resulting in an oversimplification of the hydrological process by ignoring lateral water movements and the subsequent re-infiltration and exfiltration processes, which could lead to errors in the representation of hydrological processes ^[9]. The hydrological enhancement of the Weather Research and Forecast (WRF) model, named WRF-Hydro ^[10], was designed to improve the simulation of land surface hydrology and energy states and fluxes at a fine spatial resolution (typically 1 km or less) with two modes: the coupled and standalone version. The WRF-Hydro model has a wide range of applications in water-related studies. Lee et al. (2022) ^[11] adopted the WRF-Hydro model to demonstrate the characteristics of recent droughts occurring between 2008 and 2015 in South Korea; the standardized soil moisture index and standardized streamflow index were estimated using WRF-Hydro simulations to evaluate the agricultural and hydrological droughts. Another study applied the coupled WRF-Hydro model to the flood early-warning system in the Ouémé river basin in Benin, West Africa, from 2008 to 2010 ^[12]. The WRF-Hydro model may also be used in water resources planning: A study performed on the Tono dam in West Africa utilized fully coupled WRF-Hydro to simulate streamflow from 1999 to 2003 at a 5-km horizontal resolution in the inner domain and used the output of the model as the input for a water balance model to simulate dam levels ^[13]. The literature describes the WRF-Hydro model as capable of simulating the hydrological processes in small-to-large basins. However, we consider the WRF-Hydro model to be a moderately computationally intensive modeling system when compared with other physics-based models, such as the Soil and Water Assessment Tool (SWAT) or the Hydrological Simulation Program-Fortran (HSPF).

Nowadays, hydrological studies in arid and semiarid mountainous watersheds have become an important research topic (e.g., ^[14,15]). Generally, modeling the snowmelt effect in large and mountainous basins is challenging due to high spatial and temporal variability in the model parameters and a theoretical simplification of the snowmelt–runoff mechanism within hydrological models ^[15]. Thus, the main objectives of this study are to (1) assess the performance of Noah multi-parameterization (Noah-MP) land-surface physics in simulating snowmelt processes at high altitudes with snow domination of the streamflow hydrograph and (2) simulate the hydrological process over a large watershed using fine spatial and temporal resolution.

For this purpose, the Amu River Basin (ARB), originating in northeastern Afghanistan, was selected as our study region. The ARB is the largest watershed in the country in terms of flow generation, with an extremely complex topography in the Hindu-Kush mountains, which are located in the southern part. The country has an arid to semiarid climate, receiving sparse precipitation in form of snow in the winter and rain in autumn and spring ^[16,17]. Afghanistan is extremely agriculture-dependent, with approximately 80% of the population directly or indirectly engaged in agricultural activities. The country is under development; a few decades of war and conflict destroyed the nationwide irrigation infrastructure. Based on the Ministry of Agriculture, Irrigation and Livestock (MAIL) Management Information System (MIS), the arable land includes 3.2×10⁶ ha of irrigation-fed land and 2.8×10⁶ ha of rain-fed land; 4.8×10⁶ ha of potential arable land remains uncultivated ^[18]. Afghanistan water resources are limited, but if proper management practices are adopted, there is potential for irrigating the existing farmland as well as new land development ^[19,20].

Despite its significance, only a few studies have focused on watershed modeling in Afghanistan. A few scientific works have attempted to perform hydrological modeling of Afghanistan watersheds by utilizing the SWAT to simulate the hydrological process on a monthly time scale ^{[16,21,22,23,24]}. Hydrological models with coarse spatial and temporal resolution are useful and can provide a general overview of the seasonal water-availability variation, supporting the formulation of strategies toward sustainable water resources management. So far, all hydrological modeling attempts performed in Afghanistan ^{[16,21,24,25]} have been based on a monthly time scale. However, higher temporal resolution models (daily and sub-daily time scales) are required for several purposes such as flood prediction, evaluation of catchment management consequences, reservoir management, freshwater ecology, and input provision for social, economic, and ecological models ^[26].

In this study, we assess the performance of the Noah-MP model in the simulation of snow processes over the ARB in Afghanistan. The sensitivity of the snow process to the parameterization scheme within Noah-MP land surface physics was explored by calibrating the schemes and comparing the output with the observation discharges. After optimizing the Noah-MP physics, WRF-Hydro model parameters were calibrated for the ARB to reconstruct the streamflow over three sub-watersheds in the ARB at a 3-km spatial and daily temporal resolution. This study contributes to the understanding of the community on the importance of the land surface model (LSM) of WRF-Hydro in the simulation of snow-related processes and its impact on the simulation of streamflow timing and magnitude. Dechmi et al. (2012) ^[27] emphasized the importance of continuously monitoring water quality and quantity for a better understanding of hydrological dynamics in intensely irrigated watersheds. However, the collection of long-term data is time-consuming and costly ^[28]. On the other hand, it is challenging to generalize the results of site-specific experiments to a regional level in complex watersheds with mixed land use and soil type ^[28]. In this context, the use of hydrological models could be helpful for water resources management. The findings of this study show the importance of the Noah-MP physics parameters for the simulation of the snowmelt runoff process within the WRF-Hydro model in mountainous watersheds with significant snowmelt contributions. In this study, the snowmelt runoff timing of peak flow simulation was significantly improved after calibration of the Noah-MP physics parameters.

2.   Study area and material

2.1.   Study area

Afghanistan is a mountainous country with an arid to semiarid climate. The country has a population of 34.3 million people, 24.2 million of whom live in rural areas ^[29]. The agricultural sector, industry, and services are the main sectors generating 33.7%, 16.1%, and 45.0%, respectively, of the gross domestic product (GDP). The study area is situated in the ARB, in the northeastern region of the country. The ARB is subdivided into five subbasins: Kokcha, Kunduz, Khanabad, Panj, and Ab-i-Rustaq. The study region for the current study comprises the three westernmost basins of ARB: the Kokcha (22,367.7 km²), Khanabad (11,993.5 km²), and Kunduz (28,023 km²) watersheds (Figure 1). Out of 14 stations, five were selected based on the two criteria defined in section 2.3 to eliminate outliers and possible human errors. The study area is located in a complex topographical environment with significant differences in elevation, ranging from 308 to 6847 m.a.s.l. ARB's arable land area is around 450,000 ha and the agriculture sector consumes the largest share of surface water in the basin ^[30]. Moreover, it is pertinent to mention that irrigation infrastructure is not well developed, and rivers remain largely in their natural state.

2.2.   Dataset

The minimum input data for the WRF-Hydro model include air temperature, precipitation, surface wind speed, specific humidity, shortwave radiation, longwave radiation, and surface air pressure. For these variables, the current study utilized Global Land Data Assimilation System Version 2.1 gridded data with 3-h temporal resolution as the atmospheric boundary conditions for the model. After data quality control (please refer to section 2.3), the observed discharge data from five gauging stations were shortlisted for model accuracy assessment. The MODIS-based land use data was used as input to the model. The MODIS Terra+Aqua Combined Land Cover product incorporates a 12-month input to identify two classes by the International Geosphere Biosphere Program (IGBP) based on 2001 data ^[31]. The 5-min FAO soil texture default dataset of WRF with 16 soil categories was used for the soil type ^[32]. The digital elevation model (DEM) was used for the generation of the WRF-Hydro routing grid inputs. Table 1 summarizes the input data for the simulation and accuracy assessment purposes.

2.3.   Observed data quality control

The velocity-area method is used to measure streamflow in the ARB. The Ministry of Energy and Water (MEW) uses a current meter for measuring the velocity and wading or cableways for measuring the water depth. Errors in the observed data are unavoidable and are classified as observational errors or rough errors ^[33]. The observational errors are divided into random and systematic errors. Random errors refer to misreading caused by human error, while systematic errors are the difference in the average of all observations from the true values ^[33]. On the other hand, rough errors can be caused by problematic observations, manual errors, or corrupted data during transmission. The main purpose of quality control is to remove errors from the data. In this study, after inspection of the data, two measures were considered for error removal; suspect records were deleted. We did not attempt to fill missing data or suspect data that was removed based on the following two criteria:

(1) Daily observations with the exact same value for more than 10 consecutive days were flagged as suspect data ^[34].

(2) To detect outliers, the Z-score approach was adopted ^[35]. The Z score was estimated based on Eq. (1), and values greater than 3 and smaller than −3 were considered outliers.

where ${x}_{i}$ is the observation, $\stackrel{-}{x}$ is the mean value of the observations, and $\sigma$ is the standard deviation. After removing the outliers and records with the exact same value for 10 continuous days, data with less than 10% missing values in the total observation length were considered reliable. Table 2 presents a summary of the quality control stations.

3.   Methodology

3.1.   WRF-Hydro standalone model configuration

WRF-Hydro is an open-access hydrometeorological modeling system originally developed to couple Weather Research and Forecast (WRF) atmospheric models for simulating surface and subsurface lateral water movement and shallow aquifer processes ^[10]. WRF-Hydro can also be applied in the uncoupled mode with external forcing data in offline mode or in a standalone version. In this study, the WRF-Hydro standalone mode was adopted to simulate the hydrological processes in the study region. WRF-Hydro is a multi-process, fully distributed, and multi-scale three-dimensional land surface hydrological simulation system depicting the surface, subsurface, channel, and groundwater lateral redistribution ^[9]. These interfaces enable a better representation of the relationships between water and energy fluxes at the atmospheric-terrestrial level. Moreover, for such a complex basin with high terrain, the WRF-Hydro represents a realistic simulation that considers the thermal processes and complete dynamics of the watershed ^[36].

WRF-Hydro is coupled with two land surface models: Noah and Noah-MP. The purpose of Noah-MP is to improve upon some of the limitations of its successor ^[37]. In this configuration, the Noah-MP LSM represents the energy flux and water from upstream to downstream, including surface runoff, snowmelt and accumulation, evapotranspiration, aquifer recharge, and soil water drainage and storage ^[38]. The model configuration is summarized in Table 3. In the physics option, surface overland flow routing, subsurface routing, and channel routing modules were activated. The baseflow bucket model was also activated with the exponential option. A detailed description of the WRF-Hydro modules and groundwater bucket model is available in Gochis et al., (2020) ^[39].

Figure 2 depicts the overall flow process of the model setup, calibration, and validation. Based on the available historical discharge data, the model output was set to daily discharge. The model requires an initial input, routing functions, and meteorological forcing data. The following steps were executed during data preprocessing:

(1) Defining the mode domain: Data required to execute the Noah-MP model were defined in the netCDF file generated using the GEOGRID component of the WRF preprocessing system (WPS). The GEOGRID component of WPS automated the procedures for defining the georeferencing and space and attributing most land surface data for the execution of Noah-MP ^[39]. The horizontal resolution for the Noah-MP model was set to 3 km.

(2) Initial surface conditions: Initial surface conditions such as soil temperature, soil moisture, and snow states were generated in netCDF (wrfinput_d0x.nc) with the aid of an R script by the National Centre of Atmospheric Research (NCAR).

(3) Hydrologic routing input: The WRF-Hydro routing functions are executed on the sub-grid of the LSM grid ^[40]. A suite of Python-based utilities, the WRF-Hydro GIS preprocessing toolkit ^[40], is available online in the NCAR website. The inputs for this step are the GEOGRID file, Digital Elevation model (DEM), and geographical coordinates for the generation of point data for the observation station sites. A detailed description of the procedures and output is available in Sampson & Gochis, (2018) ^[40]. The hydrologic routing grid was set to 250 m.

(4) Land surface model and lateral routing parameters file: Parameters for the Noah-MP land surface model are general, vegetation, and soil parameters associated with land cover and soil-type classes. These are available in the model directory in.txt format with a TBL extension ^[39]. As part of the National Water Model implementation effort, these land surface model parameters can be presented in a single netCDF file (soil_properties.nc) with the SPATIAL_SOIL option on during model compilation. The soil property files were generated using the R script provided by NCAR. The parameters for the lateral routing component of WRF-Hydro are defined in the text files as the LSM parameters or a netCDF file hydro2dtble.nc.

(5) Meteorological forcing data: Simulation of land–atmosphere exchanges and terrestrial hydrological processes in the standalone version of WRF-Hydro requires meteorological forcing data ^[39], obtained using the Earth System Modeling Framework (ESMF) with scripts using NCAR Command Language (NCL).

The period from January 1, 2010, to December 31, 2013, was considered as the spin-up period for model warmup and was excluded from further analysis. The period from January 1, 2014, to December 31, 2016, was selected as the model calibration period, and the period from January 1, 2017, to December 31, 2019, was selected as the model validation period.

3.2.   Calibration and validation

3.2.1.   Noah-MP physics parameters sensitivity analysis

The default LSM for WRF-Hydro version 5 is the Noah-MP, which is an augmented Noah-LSM model with multiple parameterization options ^[38]. The Noah-MP model offers an expanded range of user control process schemes, unlike its predecessor Noah ^[9]. The calibration of the Noah-MP model focuses on improving and matching the simulation discharge in the snow melting season. A reasonable representation of the snowmelt and accumulation within land surface hydrological models is key to a successful setup of the model over mountainous basins with significant snowmelt contributions. The Noah-MP model parameters were subjected to sensitivity analysis, and the best setup to represent the snowmelt timing and runoff peak was chosen. The recommended setup for the Noah-MP and the modified option for the current study model are presented in Table 4.

3.2.2.   WRF-Hydro model parameterization

The manual calibration approach was adopted for the model parameter calibration. The model was compared with the observed discharge at the gauging stations. The model calibration process attempts to improve the simulation output in terms of discharge magnitude and timing. Sensitive parameters were selected based on a literature review ^[9,41,42,43] and an accuracy assessment of the default model run. After sensitivity analysis, five of the parameters were adjusted to the new optimized value; the default value for the three parameters found to be the best were kept as default values, as summarized in Table 5. The performance of the hydrological model was evaluated using four statistical indicators as shown in Eq 2–7, namely, the correlation coefficient (R), coefficient of determination (R²), Nash-Sutcliffe efficiency (NSE), and Kling-Gupta efficiency (KGE) ^[44].

where ${Q}_{0}$ is the measured discharge (${m}^{3}/s)$, ${Q}_{m}$ is the simulated discharge (${m}^{3}/s)$, ${\stackrel{-}{Q}}_{0}$is the average observed discharge (${m}^{3}/s)$, ${\stackrel{-}{Q}}_{m}$ is the average simulated discharge (${m}^{3}/s)$, $\beta$ is the bias ratio (dimensionless), $\gamma$ is the variability ratio (dimensionless), ${\sigma }_{s, o}$ is the standard deviation of the simulated and observed discharge, and $n$ is the number of days.

4.   Results

4.1.   Noah-MP physics scheme calibration

Snow accumulation and melting simulation are highly important in the reconstruction of the streamflow signal. The sensitivity test of the Noah-MP LSM physics options shows that the simulation of snow accumulation and melting are particularly responsive to surface layer drag coefficients (SFC) and snow/soil temperature time scheme (TEMP). Figure 3a–e illustrates four different combinations of these two physics options. The hydrographs in Figure 3 depict the significance of Noah-MP LSM model configuration in the simulation of the snow process in the WRF-Hydro model. There is a huge lag between the simulated and observed discharge in all five stations using the recommended physics option. The calibrations of the Noah-MP physics option significantly improved the model simulation by eliminating the delay in the start of melting season and peak flow matched with observations.

As shown in Figure 3, the change of the SFC scheme between option 1 (M-O) and option 2 (original Noah) slightly alters the daily snow accumulation period and early melting season. Meanwhile, the SFC scheme affects the snow ablation period. A faster snow ablation was simulated using the original Noah option in the Noah-MP model. The TEMP option 1 (semi-implicit) and option 2 (fully implicit) shows the most significant improvement in the simulation result. The TEMP semi-implicit scheme simulated the snowmelt process with an earlier start of melting and early end of the snow ablation regarding the observations in all stations. The TEMP fully-implicit scheme captured well the start and end of the melting season on a daily basis. The TEMP is the option used to solve the thermal diffusion equation in the Noah-MP model ^[45,46]. TEMP has a large influence on the simulation of snow cover and melting in the model; as a result, the discharge simulated by WRF-Hydro was highly improved after selecting the proper options for SFC and TEMP. These two options are extremely important in the simulation of the snow depth and snow-water equivalent in the modeling process. In regions with limited observations of snow depth and snow-water equivalent, recorded hydrographs can be used to depict the snow melting process.

4.2.   Calibration

After optimization of the Noah-MP model parameters and improvement of the snowmelt simulation, we calibrated the WRF-Hydro model parameters. The simulated discharges from the model calibration were assessed for accuracy against the ground truth data from five stations on the three rivers of Kokcha, Kunduz, and Khanabad. Figure 4 shows the simulated hydrographs after calibration and the observed discharge at the five stations (for the locations of the observation stations and the rivers, please refer to Figure 1). The Nazdik-i-Keshem station is located in the midstream, and the Khawjaghar station originates downstream of the Kokcha River (Figure 4a and 4b). A summary of the model accuracy assessment is presented in Table 6. Overall, the statistics reveal an acceptable result, except for the NSE value at the Nazdik-i-Keshem station, which is 0.33. According to Moriasi et al., (2015) ^[47], this result is unsatisfactory, but all the other indicators are within acceptable ranges. The hydrographs for simulated and observed discharge over the Khanabad River at the Nazdik-i-Taluqan (right branch of the Khanabad River) and Pul-i-Bangi (left branch of the Khanabad River) stations are shown in Figure 4c and 4d, respectively. The statistical indicators show good agreements between the observations and simulations at the Nazdik-i-Taluqan station located in the right branch of the Khanabad River upstream with respect to Pul-i-Bangi. The NSE and R² values for the Nazdik-i-Taluqan station are not within a satisfactory range, while the KGE is reasonable for this station. At the Pul-i-Bangi station, the NSE and KGE values are not within the acceptable range, but the R² value is. The model performance for the Gerdab station, located downstream of the Kunduz River (Figure 4d), did not significantly improve and became closer to the observation. The statistical indicators show a correlation of 0.59 and 0.35 for r and R², respectively. However, the NSE and KGE values are not significant and represent negative values. The model failed to represent the peak flow in the Gerdab station and overestimated the discharge during the peak season. The presence of missing data in the peak flow season in this station may affect the statistical indicators, especially the NSE value.

There are many efficiency criteria for the assessment of hydrological models, each of which associated with some limitations. The greatest disadvantage of the NSE is that the differences between the simulation and observation values are estimated as square values; as a result, larger values in a time series are overestimated, and lower values are neglected ^[48]. As shown in Legates and Davis (1997) ^[49], correlation-based measures are more sensitive to outliers than values near the mean. Each efficiency criterion reveals different information on the model's performance. On the other hand, the study regions cover a large watershed subjected to spatial and temporal variability. A single set of model parameters cannot accurately represent all hydrological processes. The quality of the observation dataset and the presence of missing data could be potential reasons for inadequate model efficiency criterion values. Hydrologists must make subjective and objective estimates of the "closeness" of the simulated behavior of the model ^[50]. In this approach, the hydrologist may formulate the subjective assessment of the model behavior (e.g., over- and underestimation, raise limb, falling limb, and baseflow) parallel to the objective assessment using mathematical formulations ^[50].

Overall, the model properly captured the baseflow and seasonal fluctuations at all stations except for the Gerdab station (Figure 4e). The model reproduced well the seasonal fluctuation in discharge by reconstructing the streamflow during the start and end of the melting season during model calibration. The most sensitive model parameter was the soil pore size distribution index (BEXP), which controls infiltration into the soil column and baseflow as well as peak flow. The BEXP controls the actual hydraulic conductivity of the soil column according to the saturated hydraulic conductivity of the soil ^[9]. A higher BEXP during the simulation allowed more water to infiltrate the soil column, reducing the surface runoff, and vice versa. The melt factor for snow depletion curve (MFSNO) controls the snowmelt characteristics and can determine a delayed or rapid melting in the model. The maximum bucket depth (Zmax) and exponent of bucket model (EXPON) are the groundwater bucket model parameters controlling the calculation of baseflow. The WRF-Hydro model uses a simple conceptual baseflow bucket model that connects the baseflow and overlaying channel in one way ^[39]. This simple bucket model is a highly conceptualized and abstracted representation of the groundwater process where the parameters and depth of the bucket do not represent any physical meaning. On the other hand, linear scaling of the "openness" of the bottom (Slope) controls the simulated interaction between the water infiltrated into the soil column and the aquifer at the bottom of the soil column.

4.3.   Validation

The results of the validation are shown in Figure 5 for all five stations in the ARB, and a summary of the statistical analysis is presented in Table 6. A comparison of the statistical data during the calibration and validation periods shows slight differences. The statistical indicators for the Nazdik-i-Keshem, Khawjaghar, and Nazdik-i-Taluqan stations show slight differences between the calibration and validation periods. The R and R² values decreased slightly during the validation period, except for the R² for the Nazdik-i-Talquan station, which increased by 0.03. The NSE values show slight improvements for Nazdik-i-Keshem (0.13), Khawjaghar (0.06), and Nazdik-i-Taluqan (0.08). The hydrographs for the Nazdik-i-Keshem, Khawjaghar, and Nazdik-i-Taluqan stations are presented in Figure 5a, 5b, and 5c, respectively. The KGE values decreased with respect to those during the calibration period. On the other hand, the statistical indicator of the Pul-i-Bangi station decreased during the validation period compared with the calibration period. For the Gerdab station, R² (0.50) is acceptable, but KGE and NSE are not within the satisfactory range. Overall, the model captured the start and end of the melting season except for the year of 2018, when a delay in the start of the melting season was visible at the Nazdik-i-Keshem and Khawjaghar stations. The model was less successful at simulating the peak flows in the hydrographs, over-representing them most of the time.

5.   Discussion

5.1.   Baseflow simulation

The contribution of groundwater from aquifers to streamflow is highly important for hydrological modeling, especially in the simulation of baseflow in arid and semiarid regions with limited precipitation. Several studies suggest that the baseflow bucket model poorly represents the contribution of groundwater to streamflow, and due to a lack of information, the baseflow bucket model was excluded from the model structure in their studies ^[13,42]. However, for long-term simulations of streamflow where the baseflow and low-flow processes are considerable contributors, this model is useful ^[39]. The inclusion of the bucket model in the simulation process requires proper and accurate calibration of the parameters.

The calibration of groundwater parameters should be considered not only for the calibration period but also for the spin-up period. The model starts from the initial depth of water in the groundwater bucket (Zinti), and the depth increases during the spin-up period. The spin-up period should be long enough for the depth of water in the ground bucket to reach an equilibrium state. The slow bucket discharge is a function of the depth of water in the soil column to the maximum depth of water in the soil column (Z/Zmax) and of the model parameters. When the calculated depth of water exceeds the maximum proscribed depth of the soil column, the excess water essentially produces a simulated spring, producing a high discharge value as a result. It is essential to calibrate the groundwater parameters prior to the other model parameters; since it takes a few model years for the bucket model to reach a stable state, parameter tuning should also be considered during the spin-up period.

5.2.   Noah-MP model physics schemes for the snowmelt season

Snow is a crucial component of the hydrological cycle. Snow accumulation and melting strongly impact soil moisture and surface runoff, especially in mountainous regions. Snow accumulates during the winter and contributes to streamflow as the temperature rises above the melting point. Meltwater contributes to surface flow through direct runoff on the land surface or through infiltration of the soil and recharging of groundwater. The simulations of snow accumulation and melting are strongly affected by two physic schemes: (1) the surface layer drag coefficient (SFC) and (2) the snow/soil temperature time (TEMP). The TEMP option is highly influential on the snow accumulation and melting process. Switching between the TEMP options in the Noah-MP physics option significantly changes the shape of the hydrograph as highlighted in Figure 3. The TEMP option 1 simulates the snow melting season approximately one month earlier with respect to the observation data at the site. The SFC option highly affects the hydrograph during the snow ablation period with faster or slower ablation.

6.   Conclusions

In this study, the WRF-Hydro model reasonably matched the streamflow signal timing over the ARB. The findings of this study enhanced the performance of the WRF-Hydro model for the simulation of the snowmelt process and improved the melting timing and magnitude of the simulated discharge. This study recommends a careful selection of the LSM model physics schemes in watersheds with high altitude and snowmelt contribution to surface runoff. The baseflow bucket model parameters are key to a successful model calibration. Users should consider the calibration of the groundwater parameters during the spin-up period, while the remaining parameters could be calibrated during the calibration period. The simulation of the groundwater contribution to streamflow is challenging in mountainous regions due to the high topographic relief, hydrological heterogeneity, and limited data availability ^[51]. Naturally, timing of water infiltrating into the ground can produce a large effect, entering streamflow after a considerable amount of time, which makes its simulation challenging. The performance of the groundwater bucket flow model in the WRF-Hydro is reasonably good despite the simplification of the process for the long-term simulation. The soil pore size distribution index (BEXP) was the most sensitive parameter in the WRF-Hydro model parameters. BEXP improves the simulation baseflow and peak flows by controlling the amount of infiltration into soil column.

Author contributions

WH: Conceptualization, Methodology, Software, Data curation, Writing Original draft, Visualization, Investigation. HSL: Supervision, Conceptualization, Methodology, Investigation, Writing- Review and Editing.

Acknowledgments

The first author is supported by the Japanese Government Scholarship MEXT, Japan. Authors would like to acknowledge the support of the Panj-Amu River Basin Authority, especially Mr. Payanda Naimi, the Head of the Water Resources, Panj-Amu River Basin, for providing information and data for the study area.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Conflict of interest

The authors declare no conflicts of interest.