the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 12 Jun 2025
Soil physico-chemical indicators for ecosystem services: a focus on water regulation
Abstract. This study investigates the intricate relationship between soil properties and water-related processes, with a focus on their collective impact on ecosystem service provision, particularly water regulation. Conducted in three diverse regions Marchfeld (Austria), Bologna (North Italy) and Rmel (Tunisia), the research aims to identify key soil properties that influence water infiltration (INF), groundwater recharge (GWR), and crop water stress indexes (CWSI). Key soil characteristics such as saturated hydraulic conductivity (KS ), available water content (AWC), bulk density (BD), saturated water content (θs ), organic matter (OM), clay content and soil depth were analyzed for their role in regulating water movement and the overall hydrological balance. Pairwise correlation and multiple linear regression analyses were used to assess the interactions among soil water balance processes and soil properties. The results reveal significant variations between regions in terms of the factors that control infiltration, groundwater recharge, and CWSI. For example, in Marchfeld infiltration showed a strong positive correlation with BD (r = 0.74, p < 0.001), while CWSI had the most significant negative correlation with soil depth (r = -0.35, p < 0.001). Futhermore, multiple linear regression models were developed to assess the relevance of the different soil properties and of their interactions on the components of the soil water balance. As an example, in Marchfeld, the model for infiltration (r = 0.79, p < 0.001) was highly predictive, incorporating Clay, OM and soil depth. These results emphasize the critical role of key soil properties KS , AWC, BD, OM, clay content, θs and soil depth in controlling soil water processes. The study highlights the value of using these properties in predictive models to inform water management practices to optimize crop performance and soil conservation in different agricultural settings.
- Preprint
(1424 KB) - Metadata XML
- BibTeX
- EndNote
Status: open (until 27 Oct 2025)
-
RC1: 'Comment on egusphere-2025-1927', Anonymous Referee #1, 28 Jul 2025
reply
-
AC1: 'Reply on RC1', Marialaura Bancheri, 29 Aug 2025
reply
Referee #1, 28 Jul 2025
General comments
The authors studied the link between soil physical and chemical parameters and indicators of hydrological processes and ecosystem services. They performed an extensive simulation study with the hydrological FLOWS model for 3 contrasting regions regarding soil variability and climate. They showed that the most relevant soil properties for soil functions and ecosystem services depend on the environmental conditions in a region, which is an important innovation.
The paper reads fluently, the methodological setup is sound, and the results are well presented in Tables. However, I still have questions and remarks, as described below.
A. We would like to sincerely thank the Reviewer for their careful reading of our manuscript and for their constructive comments and suggestions. Their feedback has been very helpful in improving the clarity, quality, and overall rigor of the work.
Particular comments
Q. Figure 2: For each site, a different soil legend is used. This makes comparison hard. Can you create a figure with a single uniform soil legend?
A. Although we agree with the suggestion of the reviewer, the adoption of a single legend for the three sites would result in a loss of detail in the spatial distribution of the different soil types presented in the three areas. This is due to the fact the three soil maps have a different reference scale: 1:25,000 for Marchfeld (detailed soil map, where different anthropic soils units are identified and mapped), 1:50,000 for Bologna (semi-detailed soil map, where mapping units are pedolandscapes with distinct soil series) and 1:100,000 for Rmel (reconnaissance soil map, where mapping units are in most cases associations or complexes of contrasting soil groups). This results in different levels of taxonomic details for the mapping units of the three sites, which in practice hamper the adoption of a single legend without reducing the detail of the information and the soil variability depicted in the maps.
Q. Figure 3: In the Figure of FLOWS, you show the groundwater. Does FLOWS consider capillary rise?
A. Yes, FLOWS consider capillary rise of both water and dissolved solutes. FLOWS is based on the numerical solution of the Richards equation and water fluxes are calculated according to the hydraulic gradient. In the case of a shallow water table, thus close the bottom of the simulated soil domain and interacting with the soil profile, upward water fluxes (capillary rise) will develop anytime the hydraulic head at the water table will be higher than that at the simulation node just above the water table.
Q. Also related to this Figure: “lateral drainage” is drawn in the unsaturated zone. As the majority of lateral drainage occurs in the saturated zone, better to depict “lateral drainage” there.
A. There appears to have been a misunderstanding regarding the use of the term drainage. In the FLOWS model, “drainage” does not refer to lateral drainage processes naturally occurring anytime there is a lateral hydraulic gradient. Actually, FLOWS is a 1D model and would not allow for calculating these lateral fluxes. Instead, in FLOWS this term has been used to identify fluxes to artificial drains, calculated according to the Hooghoudt’s theory, which are included in the Richards equation as a sink term.
For clarity, in the figure the term drainage was changed to artificial drainage, to distinguish it from the vertical drainage and to avoid any confusion with lateral Darcian flow.
Q. Line 130: The FLOWS model plays an important role in your analysis. Can you mention studies that showed the reliability of FLOWS?
A. The model has been used for several years in several papers related to water and solute balance, even in complex soil profiles, both at the theoretical (see for example, Coppola et al., 2009; 2012; 2015a) and applicative levels (see for example Coppola et al., 2015b; 2019a, 2019b; Hassan et al., 2024; Wang et al., 2013; Porru et al., 2024).
Coppola, A., M. Abdallah, G. Dragonetti, P. Zdruli, and N. Lamaddalena. 2019a. “Monitoring and Modelling the Hydrological Behaviour of a Reclaimed Wadi Basin in Egypt.” Ecohydrology 12, no. 4: e2084. https://doi. org/ 10. 1002/ eco. 2084.
Coppola, A., G. Dragonetti, A. Sengouga, et al. 2019b. “Identifying Optimal Irrigation Water Needs at District Scale by Using a Physically Based Agro- Hydrological Model.” Water 11, no. 4: 841. https:// doi. org/10. 3390/ w1104 0841.
Hassan, S.B.M., Dragonetti, G., Comegna, A., Lamaddalena N., Coppola A, 2024. Analyzing the role of soil and vegetation spatial variability in modelling hydrological processes for irrigation optimization at large scale. Irrig Sci (2024). https://doi.org/10.1007/s00271-023-00882-7
Coppola, A., Basile A., Comegna A., Lamaddalena N., 2009. Monte Carlo analysis of field water flow comparing uni- and bimodal effective hydraulic parameters for structured soil, Journal of Contaminant Hydrology (2008), doi:10.1016/j.jconhyd.2008.09.007.
Coppola, H. H. Gerke, A. Comegna, A. Basile, V. Comegna, 2012. Dual-permeability model for flow in shrinking soil with dominant horizontal deformation. Water Resources Research, Vol. 48, W08527, doi:10.1029/2011WR011376.
Coppola A., A. Comegna, G. Dragonetti, H. H. Gerke, A. Basile, 2015a. Simulated Preferential Water Flow and Solute Transport in Shrinking Soils. Vadose Zone J. doi:10.2136/vzj2015.02.0021
Wang X., G. Quan, Y. Pan, R. H.,Y. Zhang, A. Tedeschi, A. Basile, A. Comegna, A. Coppola, R. de Mascellis, 2013. Comparison of hydraulic behaviour of unvegetated and vegetation-stabilized sand dunes in arid desert ecosystems. ECOHYDROLOGY (ISSN:1936-0584). DOI: 10.1002/eco.1265;
Coppola A., N. Chaali, G. Dragonetti, N. Lamaddalena and A. Comegna, 2015b. Root uptake under non-uniform root-zone salinity. Ecohydrol. (2014). Wiley. DOI: 10.1002/eco.1594
Porru M.C., Hassan S. B. M., Abdelmaqsoud M., Vacca A., Da Pelo S., Coppola A., 2024. Using index and physically-based models to evaluate the intrinsic groundwater vulnerability to non-point source pollutants in an agricultural area in Sardinia (Italy). Front. Water. 6:1399170. doi: 10.3389/frwa.2024.1399170
In addition, two manuscripts specifically addressing an inter-code comparison with the HYDRUS model are in preparation and will be submitted shortly to Computational Geosciences. Their provisional titles are:
FLOWS: A model to simulate water and solutes fluxes in agricultural and environmental systems. 1. Theory and applications with inter-code benchmarking
FLOWS: A model to simulate water and solutes fluxes in agricultural and environmental systems. 2. Modules for specific applications in agro-hydrological modeling
Q. Line 147: At which soil depth did you apply the free drainage option?
A. The free drainage boundary condition was applied at the bottom of each soil profile, with the depth corresponding to the specific profile considered in the simulations. We will clarify this point in the revised manuscript.
Q. Line 159: You define here the potential soil evaporation, which in general is much larger than the actual soil evaporation, an important part of your water balance. How did you calculate the actual soil evaporation?
Q. Line 165: Is potential root uptake also assumed to be proportional to root density?
A. The two questions are somehow related. FLOWS uses the reference evapotranspiration and a crop coefficient to calculate the potential evapotranspiration, ETp. This, in turn, is separated in the potential evaporation, Ep, and potential transpiration, Tp, components by using the Beer’s law, which includes an attenuation factor for radiation penetration through the vegetation.
As for actual evaporation, Ea, FLOWS calculates the maximum upward flux at the soil surface, Emax, that the soil conditions may allow, based on the hydraulic gradient between the top boundary and the first node in the soil profile, as well as on the hydraulic conductivity in that numerical compartment. Accordingly, the code uses the following condition for calculating Ea: if Ep<=Emax, Ea=Ep; if Ep>Emax, Ea=Emax.
As for actual transpiration, Ta, FLOWS uses a macroscopic root uptake (see for example, Feddes, Kowalik, and Zaradny 1978; Feddes and Raats 2004). This approach distributes the potential transpiration, Tp, along the different simulation (numerical) nodes within the rooted depth (which is an input for the model) according to a root density distribution function, g(z). The user may opt for different root distributions, such as uniform, triangular, logistic, etc. The fraction of Tp in a node multiplied by the node compartment depth, delta_z, provides the potential sink, Sp, at that node. In other words, the integral of Sp over the rooted depth corresponds to Tp.
To obtain Ta, FLOWS firstly calculates an actual sink term, Sa, by multiplying the Sp by an uptake reduction function, α(h,hos). The sink reduction factor accounts for the possible water stress (related to soil water pressure head, h) and osmotic stress (related to the osmotic potential, hos) experienced by roots at different simulation nodes within the root zone. As the reduction uptake factor is related to the local (the specific numerical node) conditions within the root zone, there may be nodes where roots are drawing water at the potential conditions (α=1 and Sa=Sp) and stressed nodes (α<1 and Sa<Sp). The integral of Sa over the rooted depth corresponds to Ta, which may thus be <= to Tp, depending on the stress conditions along the root zone.
Accordingly, details on how the model calculates actual evaporation and transpiration will be provided in the revised text.
Q. Lines 174-183: You select for each site a median year based on the amount of rainfall. Did you also check whether the rainfall distribution over the year was representative of the site's climate?
A. In addition to selecting the median rainfall year for each site, we also examined the intra-annual rainfall distribution to ensure that the selected year was representative of the typical seasonal rainfall pattern of each site.
Q. Precipitation input: At which interval did you provide input of rainfall fluxes? In order to simulate surface runoff, daily rainfall intensities are generally not suitable.
A. The model simulations were carried out at a daily time step, with precipitation provided as daily amounts. We acknowledge that this temporal resolution does not allow for the explicit simulation of surface runoff processes driven by short-duration, high-intensity rainfall events. However, the purpose of our work was to simulate the soil–water balance at the daily scale, focusing on infiltration, evapotranspiration, and soil moisture dynamics, rather than a single event-based runoff generation. Moreover, precipitation data at shorter time steps are rarely available at the large spatial scale of our three case studies, which further justifies the use of daily inputs.
Q. Lines 184 and 185: The data of soil texture and organic matter: to which soil layer do they belong? Did you assume these data to be uniform over the entire root zone depth (max 80 cm)?
A. The data on texture, bulk density, and organic matter content were available for each horizon of the investigated soil profiles. For each horizon, the hydraulic parameters were then derived by applying the PTFs. Finally, for the purpose of the statistical analysis, we computed (i) the weighted mean of these properties over the horizon depths, and (ii) the arithmetic mean across all soil profiles within each site, in order to obtain overall spatial statistics.
Q. Line 187: The crop coefficient is defined depending on the stage of crop development. How did you derive these stages?
A. The crop development stages were determined by integrating the FAO56 guidelines with expert knowledge, literature, and site-specific field observations. This approach ensured a robust and site-appropriate parameterization of crop development for each case study. For example, based on the typical crop calendar and the growth stages described in FAO56, flowering occurs earlier at RMEL than at Marchfeld and Bologna; while sowing takes place earlier at Marchfeld and Bologna compared to the Rmel site. To further strengthen transparency, we will provide the representative development stages for each site in the supplementary material and revise the main text accordingly.
Q. Table 2: The water balances of Bologna and Rmel do more or less close. However, the fields in Marchfeld, on average gained 94 mm of water. This is a substantial amount compared to other water balance components and a bit strange for a median year. Please clarify this in the text.
A. From Table 2, the higher average water storage (91.8 mm) observed in Marchfeld soils compared to RMEL (14.4 mm) and Bologna (-3.3 mm) can be attributed to their greater soil depth (Table 1) and the resulting right-skewed distribution of soil profile depths, which allows more infiltrated water to remain within the profile. Moreover, the use of a single climatic year prevents proper closure of the water balance as would be expected when extending the simulations over longer periods.
Q. Line 274: Minimum and maximum values are not depicted in Table 2.
A. We have now updated the table to include the minimum and maximum values for each variable, providing a more complete picture of the variability across sites.
Q. Lines 417 – 433: Important that you mention the ambiguity on the role of the saturated conductivity. You might add here that the simulated amount of runoff not only depends on the saturated conductivity but also on the time intervals with which rainfall is provided. In general, daily intervals underestimate the runoff amounts.
A. Thank you for your comment. We will add the following text to clarify this aspect:
“This analysis was conducted at a daily scale; however, for certain terms of the water balance, such as potential runoff production, a finer temporal resolution would be more appropriate. This is particularly relevant in the context of increasing high-intensity, short-duration rainfall events associated with climate change (Prein et al., 2017).”
Prein, A. F., Rasmussen, R. M., Ikeda, K., Liu, C., Clark, M. P., & Holland, G. J. (2017). The future intensification of hourly precipitation extremes. Nature climate change, 7(1), 48-52.
4 Discussion: Mention important limitations of your study:
- Q. Your hydrological results are based on median weather years. Commonly, in these years, hydrology functions well. The problems might arise in extreme dry or wet weather years. A thorough analysis of soil functions and ecosystem services should include more extreme weather years.
A. We fully agree that assessing soil functions and ecosystem services under extreme weather conditions is essential. While the current analysis is based on median weather years to provide a baseline characterization, we are already working on extending the assessment to include both dry and wet extremes. This will allow us to better capture the variability and potential vulnerabilities of soil functions under climate variability and change.
Accordingly, the following text will be added to the manuscript in the discussion section at line 455:
“It should be noted that, in assessing the complex interplay between meteorological variables and soil properties, this study was conducted using a median year, which provides baseline information on climate–soil interactions. Clearly, extending such analyses over longer time spans, under more extreme climatic conditions, and at finer (e.g., hourly) temporal resolution would provide additional and complementary insights.”
- Q: You focus on the durum wheat crop. In the case of other crops or multiple crops per year, the results will be different.
A: Certainly, the results would vary with different crops or multiple cropping systems. We focused on durum wheat because it is a widely cultivated rainfed crop in the three study areas, and its growth cycle is well-documented. This allowed us to rely on robust phenological information and expert knowledge to accurately define its development stages. Moreover, the absence of irrigation simplifies the analysis by reducing variability due to water management, making it a suitable reference case for our assessment framework.
- Q: Soil layering, which is explicitly mentioned in the introduction, is not included in the analysis.
A: Soil layering was in fact considered in the analysis, but we acknowledge that this was not sufficiently explained in the manuscript. As also noted in the previous comment, we will revise the text to clearly indicate how soil layering was incorporated into the model setup and analysis.
Detailed comments
Q: Lines 110 – 114: Current text is wrongly phrased and repetitive. Define the climate of site Rmel more concisely and clearly.
A: Corrected
Q: Line 134: “Genucthen” should be “Genuchten”.
A: Corrected
Q: Line 138: The symbol “h” was already defined in line 134. Define each variable once, upon first usage in the text.
A: You’re generally right. We are going to revise and delete the repetition of the variables when it is possible, according to the specific contexts.
Q: Lines 189 and 190: Include references.
A: Corrected
Q: Line 191: Include units.
A: Corrected
Q: Line 134: “groupf” should be “group”.
A: Corrected
Q: Table 2: In the caption, clarify CWSI.
A: Corrected
Q: Table 3: In the caption, clarify SD.
A: Corrected
Citation: https://doi.org/10.5194/egusphere-2025-1927-AC1
-
AC1: 'Reply on RC1', Marialaura Bancheri, 29 Aug 2025
reply
-
RC2: 'Comment on egusphere-2025-1927', Anonymous Referee #2, 27 Aug 2025
reply
General remarks:
The authors present a study on the use of the FLOWS model for identifying the correlation between soil indicators and ecosystem services. The manuscript is properly written, describing how the study is performed adequately. However, there are a number of issues that before this study can be published.
- At lines 39 – 43 you state that the current literature is in a simplified manner showing integration of soil indicators and their overall effect on soil functions, soil processes and ecosystem services. This is a rather simplified summary of the literature review, since there are many studies actually showing an integrative, holistic approach of soil indicators (e.g. Drobnik et al., 2018, Schoenholtz et al., 2000, Dominati et al. ,2010, Dominati et al., 2016, Constantini et al., 2016, to mention a few). In fact, you are not extending your analysis to a holistic approach. You use a limited set of physical soil indicators (soil hydrological parameters, partly even taken from PTFs form soil texture data) to explain their effects on a limited set of ecosystem services (groundwater recharge, infiltration, and the calculated rest products of these). What is your contribution to your statement that current studies show a simplified manner in integrating soil indicators and ecosystem services? Clarify this in your introduction and results/discussion accordingly.
- On the methodology part: you have used the FLOWS model to calculate your identified ecosystem services (RUN, INF, GWR, CWSI). Basically, this study is a sensitivity analysis of the model, without a proper description of the model itself. You indicate that you use a ‘representative median year’ for running the model (line 170), but not what the calculation timesteps of your runs are. Coppola et al.(2019) used rain events for calibrating/validating the model, with the, among others, Ks as calibration factor. You state at line 353 that Ks has a moderate to weak correlation for infiltration (and even more pronounced at line 418). But in what way is this due to the configuration of your model-design? How does the timestep affect this, or the number of profile layers in the model setup? Is there a top-layer and a more compact bottom layer in your profile? In addition, the Ks you use is not the actual Ks, since macropore flow is not taken account for by the MvG equations. There is a limited number of publications on the use of the FLOWS model, so more information is essential to understand its concept. Moreover, you do not calibrate or validate your model runs, so it remains unclear if the magnitude of the results is realistic. Explain your model configuration, validation and perform a sensitivity analysis.
Other remarks:
Line 68: you are not evaluating chemical soil parameters, only physical. OM is often considered both chemical and physical, but still leaves you with only one (partly) chemical.
Line 75: refer to the original paper of the model (Coppola et al., 2019), although in that paper the model is also not properly described. This makes understanding the model difficult.
Line 111 – 112: Both rainfall and precipitation are mentioned, use the same phrasing.
Line 134: Typo in Genuchten
Line 199 (and elsewhere): you state at places that the runoff is calculated. But the model does not calculate runoff. Actually, at line 200 you mention the correct terminology: ‘potential runoff’. Basically, this is the outcome of the equation:
Potential runoff = Precipitation – Interception – infiltration.
The remainder is ponding, potentially available for runoff (but also for evaporation or infiltration in the next timestep of calculations). Change this term accordingly.
Line 207: Typo in ‘group’
In Results: Section 3.1 should go to ‘methods’, it describes the study area and input data.
Table 1: You present the MvG parameters, as a result of the PTFs given by Wôsten et al., 2001. I assume you are aware of the very low R2 of most of those values, making your study even more difficult without having the option for validating your model results. Describe this in your discussion, how does this fit affect your results.
Line 289: The discrepancy is not due to runoff, since runoff (better: ponding) is the final component of the calculations. So, your statement that runoff removes substantial portion of the precipitation is not true, since infiltration is calculated first.
Citation: https://doi.org/10.5194/egusphere-2025-1927-RC2 -
AC2: 'Reply on RC2', Marialaura Bancheri, 10 Sep 2025
reply
Referee #2, 27 Aug 2025
General remarks:
The authors present a study on the use of the FLOWS model for identifying the correlation between soil indicators and ecosystem services. The manuscript is properly written, describing how the study is performed adequately. However, there are a number of issues that before this study can be published.
A. The authors sincerely thank the thoughtful and constructive feedback on our manuscript. They appreciate the recognition of the overall quality of the work and the adequacy of our study description. The authors will carefully address the issues raised by the Reviewer to further improve the clarity, rigor, and contribution of our study before resubmission.
Q. At lines 39 – 43 you state that the current literature is in a simplified manner showing integration of soil indicators and their overall effect on soil functions, soil processes and ecosystem services. This is a rather simplified summary of the literature review, since there are many studies actually showing an integrative, holistic approach of soil indicators (e.g. Drobnik et al., 2018, Schoenholtz et al., 2000, Dominati et al. ,2010, Dominati et al., 2016, Constantini et al., 2016, to mention a few). In fact, you are not extending your analysis to a holistic approach. You use a limited set of physical soil indicators (soil hydrological parameters, partly even taken from PTFs form soil texture data) to explain their effects on a limited set of ecosystem services (groundwater recharge, infiltration, and the calculated rest products of these). What is your contribution to your statement that current studies show a simplified manner in integrating soil indicators and ecosystem services? Clarify this in your introduction and results/discussion accordingly.
A. At lines 39-43 we wanted to remark that the inherent complexity of soil functions resulting in the processes underpinning the delivery of soil-based ecosystem services hamper to some extent the efficacy of simplified indicators which, by definition, are static being based on a limited set of soil properties. If from the one hand holistic approaches of soil indicators have allowed the definition of theoretical frameworks for the integration of soil-based ecosystem services into policy and land planning, and actually being quite successful at doing so, from the other there is the risk to overlook the actual complexity of field dynamics and their temporal variability. Our contribution was not intended to extend the analysis of modelling results in the three case study areas to a more holistic approach, but rather to provide insight into the relationships between the outcomes of a dynamic process-based model and a set of soil properties frequently used as proxies for the assessment of regulation ecosystem services, focusing on water regulation in different climatic conditions.
We have clarified the aim of our contribution by modifying the text in the Introduction and Results/Discussion sections, as suggested by the Reviewer.
Q. On the methodology part: you have used the FLOWS model to calculate your identified ecosystem services (RUN, INF, GWR, CWSI). Basically, this study is a sensitivity analysis of the model, without a proper description of the model itself. You indicate that you use a ‘representative median year’ for running the model (line 170), but not what the calculation timesteps of your runs are. Coppola et al. (2019) used rain events for calibrating/validating the model, with the, among others, Ks as calibration factor. You state at line 353 that Ks has a moderate to weak correlation for infiltration (and even more pronounced at line 418). But in what way is this due to the configuration of your model-design? How does the timestep affect this, or the number of profile layers in the model setup? Is there a top-layer and a more compact bottom layer in your profile? In addition, the Ks you use is not the actual Ks, since macropore flow is not taken account for by the MvG equations. There is a limited number of publications on the use of the FLOWS model, so more information is essential to understand its concept. Moreover, you do not calibrate or validate your model runs, so it remains unclear if the magnitude of the results is realistic. Explain your model configuration, validation and perform a sensitivity analysis.
A: This comment from Reviewer 2 covers several topics. For clarity, we have structured our reply into two parts: the first addresses the overarching questions, while the second focuses on specific points.
A. Overarching questions:
On the comment regarding the ‘representative median year’ used for running the model (line 170), as already replied to Reviewer 1, in addition to selecting the median rainfall year for each site, we also examined the intra-annual rainfall distribution to ensure that the selected year was representative of the typical seasonal rainfall pattern of each site. As described in Section 2.3, the representative median climatic year was calculated for each of the three sites. The selected representative median year corresponds to the year within the full observation period whose annual precipitation was closest to the median of the long-term dataset at each site. This approach ensures that the chosen year reflects typical precipitation and evapotranspiration conditions, making it representative of the site’s median hydro-climatic regime.
One of the reviewer’s requests concerns additional details on the model. While it is not feasible to include all equations and processes implemented in the model within the manuscript, we will clarify several aspects that were identified as unclear in Section 2.2 and provide a link to the model manual for further technical details (Coppola et al., 2025). In addition, as already noted in our reply to Reviewer 1, two manuscripts including inter-code comparisons with HYDRUS-1D are currently in preparation and will soon be submitted.
Coppola, A., Hassan, S. B. M., & Pili, F. (2025). FLOWS - FLOw of Water and Solutes in agricultural and environmental systems (WEB Version 1.0). Zenodo. https://doi.org/10.5281/zenodo.17076523
The second general comment concerns what we consider a fundamental rationale for this work, which may not have been adequately and convincingly conveyed in the initial manuscript. This study is not conceived as an exercise in applied modelling in the strict sense. Rather, its primary purpose is to emphasize—particularly for the community of soil scientists focused on soil health and ecosystem services, but less familiar with soil hydrology—that relationships between indicators and hydrological processes are inherently complex and indirect. Such relationships are influenced by multiple interacting factors (e.g., soil layering, upper and lower boundary conditions) that cannot be captured by a simple one-to-one indicator–process assessment.
Finally, we do not consider this work reducible to a mere sensitivity analysis of model parameters, mainly because our evaluation is not based on perturbations of direct model inputs (such as van Genuchten parameters or LAI). Instead, it relies on widely used indicators of water regulation ES, which may coincide with direct model parameters (e.g., Ksat) or may not (e.g., AWC).
- Specific comments: (on the sentence Ks shows only moderate to weak correlation with infiltration)
- Time steps used in the simulation runs: Information on the daily time step used in the simulations is already reported in Section 2.2.2 Boundary and initial conditions for water flow. To improve clarity, we will restate this point when presenting the simulation results. In addition, as also discussed in our response to Reviewer 1, we acknowledge that for specific terms such as runoff generation, a finer (e.g., hourly) time step would be more appropriate, particularly in the context of high-intensity, short-duration rainfall events.
- On the role of soil layering: Each soil profile was discretized into multiple horizons based on available site data. This stratification is a key factor influencing infiltration, as water redistribution across horizons reduces the direct influence of topsoil Ks on overall infiltration. The soil layering represents a clear example of the soil complexity highlighted in the manuscript and the need of a modelling approach. For instance, a sand horizon overlying clay produces completely different infiltration and redistribution dynamics compared to the reverse arrangement, with outcomes that cannot be predicted a priori. For this reason, a physically-based modelling approach can adequately capture such interactions and provide reasonable estimates of water balance components. Conversely, identical soil properties arranged in a different layering sequence, or under a different climate, may lead to markedly different results.
Regarding the Reviewer’s question on the presence of a top layer overlying a more compact bottom layer, given the large number of soil profiles considered, we observed a variety of situations. Not all soils follow the classical sequence of increasing compaction with depth. Instead, below the ploughed Ap horizon, several horizons with different degrees of compaction are present, resulting in a range of stratifications across the profiles.
Finally, since this information on the number of soil horizons was missing from the original text, we have now added the following sentence in ‘Section 2.3 - Input dataset for model simulations’: “The number of horizons reported in the soil maps of the three areas ranges from a minimum of 1 to a maximum of 5. On average, Marchfeld, Bologna and Rmel profiles include 3 horizons.
- Ks influenced by model design/configuration: Certainly, the correlation between Ks and the different processes (i.e., infiltration) is influenced by the model design/configuration. For instance, a bucket-type model with a homogeneous soil profile would yield markedly different results compared to a physically based approach, as used in our study, where processes are simulated with state-of-the-art numerical methods. In particular, the physically-based approach with the Richards' equation explicitly accounts for real driving force, soil layering, boundary conditions, flux-to-runoff, and the dynamic redistribution of water within the profile, thereby reducing the apparent direct control of topsoil Ks alone and providing more realistic insights into infiltration and water balance dynamics.
- Ks influenced by macropores: The reviewer is certainly correct that the Ks used in our simulations does not represent the field-saturated conductivity (Ksf ). This limitation is well recognized by the authors, and for this reason we deliberately denoted it simply as Ks , in line with common practice to distinguish between the two. Furthermore, considering that our analysis involved simulations across 315 soil profiles, it is practically impossible—at least to our knowledge—to obtain consistent information on macropore presence and their parameterization at such a scale. Therefore, our focus was placed on matrix conductivity, while acknowledging that macropore flow is an important process at the local scale but remains beyond the scope and feasibility of this study.
- Implications for the realism of the magnitude of results: Explain your model configuration, validation and perform a sensitivity analysis. We appreciate the Reviewer’s concern regarding the validity of the model and the request for further information. To address this, we note that the full manual—containing all equations and details of the model design—together with the scientific papers in which the model has already been published and validated, provides a solid basis for transparency and reproducibility. In addition, Chapter 5 of the manual (Coppola et al., 2025) evaluates the reliability of FLOWS by comparing its results with those obtained from HYDRUS-1D simulations, further supporting confidence in the model’s performance.
With respect to the request for a sensitivity analysis of the model parameters, we fully agree that this is a valuable exercise. Such an analysis is indeed in preparation and will further strengthen the model’s applicability.
Other remarks:
Q: Line 68: you are not evaluating chemical soil parameters, only physical. OM is often considered both chemical and physical, but still leaves you with only one (partly) chemical.
A: We agree that our study mainly focuses on soil physical parameters (e.g., texture, bulk density), while organic matter (OM) can be considered both a chemical and a physical property. To avoid ambiguity, we revised the text accordingly and now refer to “soil parameters” instead of “chemical and physical parameters.”
Q: Line 75: refer to the original paper of the model (Coppola et al., 2019), although in that paper the model is also not properly described. This makes understanding the model difficult.
A: See our answer in the “Overarching questions” section
Q: Line 111 – 112: Both rainfall and precipitation are mentioned, use the same phrasing.
A: We will use the term rainfall throughout this work.
Q: Line 134: Typo in Genuchten
A: Thanks, we will correct the Genuchten in the text
Q: Line 199 (and elsewhere): you state at places that the runoff is calculated. But the model does not calculate runoff. Actually, at line 200 you mention the correct terminology: ‘potential runoff’. Basically, this is the outcome of the equation:
Potential runoff = Precipitation – Interception – infiltration.
The remainder is ponding, potentially available for runoff (but also for evaporation or infiltration in the next timestep of calculations). Change this term accordingly.
A. The model does not compute runoff as a residual of the water balance; rather, it explicitly simulates both Hortonian and Dunnian runoff generation processes. We will clarify this point in the model description. To prevent misunderstanding, we have replaced the term “runoff” with “flux-to-runoff,” which more accurately represents the model’s approach. FLOWS is a one-dimensional model and does not formally simulate routing runoff. Furthermore, the model allows the user to specify a maximum acceptable ponding depth (e.g., 10 cm): any flux exceeding this threshold is considered flux-to-runoff. After the rainfall ceases, the ponded water gradually infiltrates into the soil profile.
Q: Line 207: Typo in ‘group’
A: Thanks, we will correct the group in the text
Q: In Results: Section 3.1 should go to ‘methods’, it describes the study area and input data.
A: In Subsection 3.1, as indicated in the introduction to the results, we aimed to present the soil data used in the simulations, emphasizing the study area’s diverse pedological characteristics. For this reason, we considered it appropriate to include this subsection in the results. At the same time, we recognize that the table itself is more closely related to the description of the data and procedures, and therefore we have moved the table to Section 2 (Materials and Methods), while keeping the interpretative text in the results. This way, both clarity and consistency with scientific conventions are maintained.
Q: Table 1: You present the MvG parameters, as a result of the PTFs given by Wôsten et al., 2001. I assume you are aware of the very low R2 of most of those values, making your study even more difficult without having the option for validating your model results. Describe this in your discussion, how does this fit affect your results.
A: The reviewer is certainly correct in pointing out the limitations of relying on PTFs. This remains a fundamental challenge, especially when addressing problems at a regional rather than a local scale. As we have already written in our response to Reviewer 1, we will strengthen the discussion section to more explicitly acknowledge and elaborate on the limitations associated with the use of PTFs, so that the reader is fully aware of their implications in the present context.
Q: Line 289: The discrepancy is not due to runoff, since runoff (better: ponding) is the final component of the calculations. So, your statement that runoff removes substantial portion of the precipitation is not true, since infiltration is calculated first.
A: This question appears to stem from a misunderstanding of how the model calculates runoff (now referred to as flux-to-runoff). We kindly refer the Reviewer to our previous response on this issue, where the procedure is explained in detail.
Citation: https://doi.org/10.5194/egusphere-2025-1927-AC2
-
RC3: 'Comment on egusphere-2025-1927', Anonymous Referee #3, 20 Sep 2025
reply
The manuscript presents and elaborate and well structured study in a nicely written and structured way, for which I congratulate the authors. However, I do think that there several key issues that would need to be carefully considered before this submission can be published. Some of the comments call for potentially substantial amount of additional work and re-analysis of the data.
- On the description of the modeling case: I somewhat lack a proper description of the modeling process, parameterization, etc.: a clear delineation and discussion of which of the reasonably important processes are (or are not) modelled by the model. Was the full functionality of the model used or did data limitation or some other reason reduce what processes were actually considered? For instance, preferential flow is a key process to consider given the aim of the paper. I also generally lack clarity on the parameterization of the model, how many profiles have been modeled, as well as any mention of validation of the model for the purpose or location(s).
- On the matter of scale: To me it feels that there is a lack of connection established between the scale of modeling and the scale at which conclusions are made: The model is apparently a 1D simulation model that has been run on individual soil profiles. However, “spatial differences” have been referred to (around L255), while there has not been any information provided on how the upscaling to 100s of km2 has been handled. Each area is obviously heterogeneous in terms of soil type and depth. I think this should be discussed, along with the limitations of the chosen approach.
- On using only the average year in the simulations: I find working with just a single “average” year to be a weakness of the study. Most of the events of concern – be that runoff, erosion, pollution, drought, etc. – are usually the products of either extreme events or extreme seasons/years. This is unfortunately not considered in the study, and I feel that the given advice is not based on the best available information. Working with a distribution of weather conditions (years) at each site would also allow the deduction of correlations between various factors in a much deeper way, allowing the study’s conclusions to be much better supported. Given that the data are available, the workflow would just need to be extended, and the results analysed at a higher dimensionality. I think the best value for the reader would be if the analysis was rerun on all available years and the correlation analysis was deeper, accounting for interactions and distributions that look into e.g. precipitation distribution (same soil under wetter and drier conditions) and soil depth as well, which is crucial when capacity-driven metrics are of concern.
- On the diversity of soil depths and their handling: Soil depth is obviously one of the keys to the findings, whereas it has not really been addressed as such. It has mentions, but it is a key factor that I believe requires much more attention. The study mentions e.g. drainage and ground water recharge. At what depth are they defined in this study? How are these factors handled when the soil is shallow?
- On some of the obvious uncertainties that are not being acknowledged and handled: Let me pick an example. I think the conclusions made regarding Ksat vs. infiltration at L417-429 are not thought through sufficiently and therefore are a bit careless. Essentially the issue is a set of limitations that originate from possibly several sources. In this study Ksat has been estimated by a generic PTF, which smoothens differences and is generally error prone. The other variable, infiltration is a result of modeling, which is prone to its own set of uncertainties and errors. Neither of them are measured here. Ksat and infiltration are not the culprits, but the lack of our ability to account for them is. Such sets of uncertainties should be acknowledged and discussed, and the findings placed in context, and the conclusions handled accordingly.
- On advancing our understanding: The authors cite that past studies generally did not go to sufficient depth (L40-41), however, as a reader I do not really see how this study advances the state-of-the-art all that much. Much of the conclusions are rather common sense and reinforce relatively basic, textbook type knowledge that is relatable to soil depth and the “average” soil’s physical and hydrological characteristics. To refer to some examples, several very logical connections can be made in Table 2 already: infiltration is strongly influenced by Ksat (and likely soil thickness), runoff is strongly connectable with prevailing precipitation intensity and soil texture type, drainage by precipitation characteristics and Ksat. If there is no rain, there is no infiltration or runoff. If the soil is more course textured, there is good infiltration and no runoff. Examples are in e.g. lines 285-286, 291-318 and 336-337.
The quantitative confirmation by modeling is appreciated of course, but I feel that the reader does not receive sufficiently impactful new information that is of interest to THEM.
The authors should generally recognize and discuss the study’s limitations that may root from the model itself, the choices made in parameterization, the scale problems, simplifications (depth, climate “average”), and also data sources. I think the value of the reported study could be enhanced a lot by both addressing several of these with improved analysis, and discussing the rest.
Some detailed comments:
L40-41: This statement invites that the authors make a very clear statement later on how their study improves on this.
L43-48: Here I recommend consulting e.g. https://doi.org/10.1016/j.earscirev.2021.103873 and https://www.nature.com/articles/s41586-018-0463-x and citing them if they are found relevant and fitting.
L151-152: Please provide more on this – how the rule/decision was made which condition to apply. Also: was there any warm-up period considered for the model?
L155-156: How were those parameters obtained?
L167: How were the h1-h4 thresholds handled? Crop specific? Fixed? Please provide info.
L170: The choice of median year largely dismisses the significance of infiltration capacity in really wet years. Please acknowledge this appropriately.
L185-186: Please use the original 1999 Geoderma reference for the HYPRES PTFs.
L207-210: At least some of the outputs and inputs/indicators are rather closely correlated. Could you please justify their handling? Example: infiltration vs. Ksat. Also: define AWC.
L231-234: Please clarify the distinction between statistical and multivariate “relationships”. Latter is also a statistical method.
L244: Is soil depth accounted for in the modeling approach? Please describe. Table 1 suggests that there have been huge differences in soil depth _within_ each study area. Soil depth of course matters a lot when it comes to the soil’s buffering capacity in e.g. infiltration capacity, total pore volume in the profile, total available water for the crop, etc.
L264-265: These differences in coarse and fine textured soil Ks values are barely different. Usually the variability of local measurements would not even make such values statistically different. I recommend double checking and also noting this in the manuscript. PTFs tend to smoothen within the available data range, and reduce such differences normally seen in the field.
L272: “all” soil profiles? Please clarify what this really means. I also note that averaging over a huge diversity of soils within each region may simply mask the hot-spots, that matter a lot in reality.
L336-337: re: stress correlations: I think these were logically expected. Bologna has the deepest soils on average, as well as the most clayey ones, and the highest amount of precipitation. Increased clay just adds to retention, and it is likely that the “h1-h2” (wet) stress limit kicks in more often, so the reduction of clay has a positive effect. In Rmel the conditions are exactly the opposite: dry region, highest sand content, shallowest soils “on average”. Then the “h3-h4” (drought) stress limit is likely often met, so any clay present will help retain water and make it available to plants. To me this feels rather simple “on average”.
L336-339: Why not run all years and then plot and correlate precipitation, clay content and CWSI? A cloud of points and deeper analysis of the potentially available data would lend much more support to the claims made.
L340: This is against classic textbook knowledge, although it has been reported by others before. The authors could attempt to go deeper in explaining it. A starting point can be https://doi.org/10.2136/sssaj2004.0055 or https://doi.org/10.2136/vzj2016.03.0021, but there are a few other references available since.
L385-388: I think the main point is not to increase some correlation coefficients, but to help understand how variable interactions help explain why certain correlations between variables and indicators are seen and how they interact, especially in “non-average” cases. I suggest trying to focus on what can be generalized for the benefit of the independent reader, who do not associate directly with the cases of Bologna or Marchfield or Rmel.
L392: What T values are referred to? Introduce the statistic used and show the data that you are referring to.
L399-401: Without depth considered, basically the same thing has been stated already on L336-337.
L400-402: On L335-336 it has been stated that in Bologna moderate correlation was found and higher clay content has been associated with increased stress. Here it is explained that only depth (i.e. not clay) showed correlation with stress. So which one is it? Running different analyses should not be the explanation. Then one of them is the incorrect method.
L416: The variability within each of the 3 areas has been greatly masked, whereas that is what would inform the reader how to judge their own cases. I don’t really see deep enough support for that in this study.
L417-429 (and also to some extent L430-439): Ksat is being singled out here, so I come to its defense to some extent. The authors should deepen their arguments at least where prior knowledge exists. Ksat of course is known to show great spatial variability that is very difficult to capture by PTFs that use only basic physical soil properties that have very limited relationship with soil structure. The estimates by PTFs are known to be smoothened (often a lot) compared to the original data support, and be often weak descriptors of the original data, no matter what is claimed about them. I ask the authors to look up the R2 value of e.g. the used HYPRES European PTFs – they will be surprised how low it is (around 0.2, as I remember). The mean Ksat values presented on L265-266 show such sign of smoothing. The problem is by far not only with Ksat and infiltration not being strong-enough correlated, but that in this study a PTF was used – and a rather generic one. Essentially, this study did not have Ksat to work with. It had some estimations of Ksat to work with – i.e. not necessarily the correct data. Second: I feel that the model’s imperfections, simplifications, simulated processes have not been thoroughly accounted for and discussed. The authors try to draw conclusions between an estimated variable and an indicator whose simulation is dependent on a multitude of imperfectly modeled processes, under “average” climatic conditions, and smeared over varied soil depths/types. I think all such uncertainties should be acknowledged, and at least some that are feasible to handle in a more detailed way should be imroved. Relevant conclusions should of course be modified (softened?) accordingly.
L440-455: I feel that the authors could have done a deeper analysis to support this section with more results, and make it more concrete. Rainfall distribution could have been characterized, the frequency of intense (defined as…?) rain events quantified, etc.
456: This is exactly what is needed – basing the study on appropriate data. I don’t want to say that a study based on PTF estimates is illegitimate, but the conclusions drawn should account for the source of errors that originate from PTF use, model use, spatial and temporal simplifications and averaging, etc.
Citation: https://doi.org/10.5194/egusphere-2025-1927-RC3
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 675 | 67 | 22 | 764 | 17 | 29 |
- HTML: 675
- PDF: 67
- XML: 22
- Total: 764
- BibTeX: 17
- EndNote: 29
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
General comments
The authors studied the link between soil physical and chemical parameters and indicators of hydrological processes and ecosystem services. They performed an extensive simulation study with the hydrological FLOWS model for 3 contrasting regions regarding soil variability and climate. They showed that the most relevant soil properties for soil functions and ecosystem services depend on the environmental conditions in a region, which is an important innovation.
The paper reads fluently, the methodological setup is sound, and the results are well presented in Tables. However, I still have questions and remarks, as described below.
Particular comments
Figure 2: For each site, a different soil legend is used. This makes comparison hard. Can you create a figure with a single uniform soil legend?
Figure 3: In the Figure of FLOWS, you show the groundwater. Does FLOWS consider capillary rise?
Also related to this Figure: “lateral drainage” is drawn in the unsaturated zone. As the majority of lateral drainage occurs in the saturated zone, better to depict “lateral drainage” there.
Line 130: The FLOWS model plays an important role in your analysis. Can you mention studies that showed the reliability of FLOWS?
Line 147: At which soil depth did you apply the free drainage option?
Line 159: You define here the potential soil evaporation, which in general is much larger than the actual soil evaporation, an important part of your water balance. How did you calculate the actual soil evaporation?
Line 165: Is potential root uptake also assumed to be proportional to root density?
Lines 174-183: You select for each site a median year based on the amount of rainfall. Did you also check whether the rainfall distribution over the year was representative of the site's climate?
Precipitation input: At which interval did you provide input of rainfall fluxes? In order to simulate surface runoff, daily rainfall intensities are generally not suitable.
Lines 184 and 185: The data of soil texture and organic matter: to which soil layer do they belong? Did you assume these data to be uniform over the entire root zone depth (max 80 cm)?
Line 187: The crop coefficient is defined depending on the stage of crop development. How did you derive these stages?
Table 2: The water balances of Bologna and Rmel do more or less close. However, the fields in Marchfeld, on average gained 94 mm of water. This is a substantial amount compared to other water balance components and a bit strange for a median year. Please clarify this in the text.
Line 274: Minimum and maximum values are not depicted in Table 2.
Lines 417 – 433: Important that you mention the ambiguity on the role of the saturated conductivity. You might add here that the simulated amount of runoff not only depends on the saturated conductivity but also on the time intervals with which rainfall is provided. In general, daily intervals underestimate the runoff amounts.
4 Discussion: Mention important limitations of your study:
Detailed comments
Lines 110 – 114: Current text is wrongly phrased and repetitive. Define the climate of site Rmel more concisely and clearly.
Line 134: “Genucthen” should be “Genuchten”.
Line 138: The symbol “h” was already defined in line 134. Define each variable once, upon first usage in the text.
Lines 189 and 190: Include references.
Line 191: Include units.
Line 134: “groupf” should be “group”.
Table 2: In the caption, clarify CWSI.
Table 3: In the caption, clarify SD.