A Spatial Analysis to Define Data Requirements for Hydrological and Water Quality Models in Data-Limited Regions

The objective of the present work is a spatial analysis aimed at supporting hydrological and water quality model applications in the Canale d’Aiedda basin (Puglia, Italy), a data-limited area. The basin is part of the sensitive environmental area of Taranto that requires remediation of the soil, subsoil, surface water, and groundwater. A monitoring plan was defined to record the streamflow and water quality parameters needed for calibrating and validating models, and a database archived in a GIS environment was built, which includes climatic data, soil hydraulic parameters, groundwater data, surface water quality parameters, point-source parameters, and information on agricultural practices. Based on a one-year monitoring of activities, the average annual loads of N-NO3 and P-PO4 delivered to the Mar Piccolo amounted to about 42 t year−1, and 2 t year−1, respectively. Knowledge uncertainty in monthly load estimation was found to be up to 25% for N-NO3 and 40% for P-PO4. The contributions of point sources in terms of N-NO3 and P-PO4 were estimated at 45% and 77%, respectively. This study defines a procedure for supporting modelling activities at the basin scale for data-limited regions.


Introduction
Hydrological and water quality models have been used largely in environmental studies for identifying sources of pollutants and developing scenarios in order to test environmental measures to reduce pollution or mitigate climate change impacts on water resources [1][2][3][4][5][6].Several ecohydrological models have been developed and applied around the world for simulating hydrology and water quality, including the Soil and Water Assessment Tool (SWAT) [7] and the Annualized Agricultural Non-Point Source (AnnAGNPS) models [8].Models require that a large number of environmental parameters, spatial data, time series, and measurements be defined for their set up, calibration, and validation.A large amount of data may be required [9,10], especially for complex models, constituting a limit on model applicability in data-limited regions [11].Moreover, model results are influenced strongly by the datasets and procedures used in the model's application.Further, a model's quality depends on the use of appropriate methodologies (set up, calibration, validation) and on the possibility of verifying intermediate results [12].After defining the aims of the study and building the conceptual model, a geodatabase must be created.Data may be derived from different sources (e.g., municipalities, environmental agencies, the River Basin Authority, the Civil Protection Service, and so forth) and may not be uniform in terms of time scale and spatial resolution [13].In addition, time series frequently include several gaps.Hence, it is necessary to collect and accurately cross-check available data, which should be validated and integrated with field campaigns and farmer/citizen interviews to clarify any inconsistencies in the data or to update datasets.
Model calibration and validation procedures require streamflow and water quality data [14,15].Hence, in data-limited regions, it may be necessary to install a number of gauging stations for streamflow measurements and to define a monitoring plan for water quality samples.The creation of a geodatabase, as well as surveys and field activities, is time consuming and may require considerable economic resources [16].For these reasons, it is important to define a protocol before starting modelling activities, taking into account the aim of the study, the peculiarities of the area, as well as the time and economic resources available.A clear procedure (plan) helps modellers in making assumptions and choices consistent with the objectives of the modelling project, in reducing bias, and, at the same time, allowing other modellers to repeat or continue the project [17].Although all of these preparatory activities are fundamental for the application of a model, they are generally neglected in publications.
The general objective of the present work was to develop a methodology to support modelling activities conducted at the basin scale.With this aim in mind, a spatial analysis of the Canale d'Aiedda basin in the province of Taranto (Puglia Region, Italy), which is recognized a sensitive environmental area due to industrial and agricultural pressures, was conducted.The specific objectives of the work were to: (i) realize a geodatabase that included a large number of distributed and point data from different sources, previous studies, and from new, specifically-implemented environmental monitoring activities, and (ii) estimate the monthly nutrient loads and their knowledge uncertainty for two river sections based on observations taken over the course of a year (nutrient and streamflow), which is required for calibrating and validating models.The results, on the one hand, increase the knowledge of the area and its actual environmental status; on the other hand, they constitute a methodology for supporting modelling activities in areas with limited data.

Study Area
The Canale d'Aiedda basin is located in the Puglia Region in Southern Italy (Figure 1).The climate, land cover, and management practices are representative of the Mediterranean region.The basin drains an area of approximately 408 km 2 , and the Canale d'Aiedda stream flows into the second inlet of the Mar Piccolo (literally "Small Sea"), which is an inner, semi-enclosed sea basin (surface area of 20.72 km 2 and a coastline perimeter of 28 km).The land is used predominantly for agriculture, particularly in the form of fruit tree species (e.g., vineyards, olive groves, and citrus groves) and vegetables.There is also a significant presence of simple arable crops (e.g., durum wheat, soft wheat, barley, and oat), which are distributed uniformly throughout the basin.The presence of breeding Water 2019, 11, 267 3 of 22 farms is low, and it is mostly localized in the north of the basin.The area is characterized by flat territory, a hilly part with a maximum elevation of 520 m above sea level (a.s.l.) and the presence of morphologically depressed areas, which are subject to flooding during particularly intense rainfall events [18].This is the sub-basin located to the south of the Mar Piccolo called the Salina Grande, which, through the Canale Maestro (partially-buried river channel), enters into the Canale d'Aiedda.The river network is not well defined in the northern part of the basin since this area is karstic.

Data Requirements for Modelling Hydrology and Water Quality
Hydrological and water quality models require a large amount of spatial data, time series, and measurements to be defined for their set up, calibration, and validation.
Before starting the modelling activities, it is necessary to define the aims of the study, which also guide the model choice.Physical (topography, river network, land use, soil), meteorological (rainfall, temperature, solar radiation), non-point-(agricultural practices) and point-source pollution, and hydrological (streamflow, sediment and nutrients loads) data must be collected on the basis of the defined objectives, basin characteristics, model, and economic resources.Thus, a preliminary study of and a survey campaign regarding the basin constitute a guide for developing the conceptual model, building a specific tool or selecting an existing model, optimizing the costs associated with the required data collection, and obtaining reliable results.
In the Canale d'Aiedda basin, the aims of the future model application are to: (i) identify sources of nutrient pollution within the basin, and (ii) quantify the nutrient loads delivered to the Mar Piccolo by the river network over a long period.To accomplish these aims, the SWAT model [7] was used.There were two main reasons for selecting the SWAT model-it is well documented, and it has already been successfully applied in the cases of a number of Mediterranean basins [5,6].The current research focuses on the data requirements for setting up, calibrating, and validating the SWAT model.However, the geodatabase and the monitoring activities developed in this study can be applied to The Mar Piccolo has lagoon features, and it is divided into two inlets; the first one has an area of 8.1 km 2 and a maximum depth of 13 m, and the second one has an area of 12.6 km 2 and a maximum depth of 8 m.It is connected with the open sea, i.e., the Mar Grande (literally "Big Sea"), in the Mar Ionio (i.e., Ionian Sea) through two channels.At the bottom of this marine basin, there are submarine springs (referred to locally as "citri") which discharge groundwater from the deep limestone aquifer located in the study area [19].For this reason, the Mar Piccolo is the most important area for mussel farming in Italy [20].Furthermore, it is characterized by the presence of fisheries as well as an abandoned Italian Navy shipyard and its port [21].The surrounding area represents a sensitive environmental area due mainly to the presence of industries, including the largest Italian steel company (with a total area of 15.45 km 2 ), a refinery (2.75 km 2 ), and the tourist port [22,23].The environmental impact in terms of health and ecological risks as well as cultural and environmental heritage is very high, and the Ministry of the Environment, via National Law no.426/98, declared Taranto a "Contaminated Site of National Interest (SIN)" (i.e., "Sito di Interesse Nazionale").More specifically, the Taranto SIN consists of the SIN Mare (i.e., Sea SIN, 70 km 2 ) and SIN Terra (i.e., Land SIN, 43 km 2 ).The interventions of environmental requalification and the preliminary studies are performed within the "area of environmental crisis", which partially includes the Canale d'Aiedda basin.
The basin includes medium-sized urban areas and urban wastewater treatment plants (WWTPs), which involve three stages of sewage treatments (i.e., primary, secondary, and tertiary treatments).Treated effluents are discharged into the river system.The channel is heavily modified (concrete river banks).The stream flows into the "Regional Nature Reserve Swamp la Vela" wetland, which is a protected area (World Wide Fund, WWF oasis, 1.16 Km 2 ) included in the Site of Community Importance (SCI) "Mar Piccolo" (IT9130004).The wetland is characterized by typical Mediterranean natural vegetation (i.e., Mediterranean maquis), consisting mostly of shrubs and small trees (e.g., myrtle, mastic, and holm oak) and able to offer shelter and food to many resident species (e.g., grey herons and finch) and migratory birds (e.g., flamingos and curlew).The rich flora of the oasis is of the halophile type and is suitable for living in salty environments.Due to the urban settlements and agricultural activities, the wetland is currently strongly threatened [24].

Data Requirements for Modelling Hydrology and Water Quality
Hydrological and water quality models require a large amount of spatial data, time series, and measurements to be defined for their set up, calibration, and validation.
Before starting the modelling activities, it is necessary to define the aims of the study, which also guide the model choice.Physical (topography, river network, land use, soil), meteorological (rainfall, temperature, solar radiation), non-point-(agricultural practices) and point-source pollution, and hydrological (streamflow, sediment and nutrients loads) data must be collected on the basis of the defined objectives, basin characteristics, model, and economic resources.Thus, a preliminary study of and a survey campaign regarding the basin constitute a guide for developing the conceptual model, building a specific tool or selecting an existing model, optimizing the costs associated with the required data collection, and obtaining reliable results.
In the Canale d'Aiedda basin, the aims of the future model application are to: (i) identify sources of nutrient pollution within the basin, and (ii) quantify the nutrient loads delivered to the Mar Piccolo by the river network over a long period.To accomplish these aims, the SWAT model [7] was used.There were two main reasons for selecting the SWAT model-it is well documented, and it has already been successfully applied in the cases of a number of Mediterranean basins [5,6].The current research focuses on the data requirements for setting up, calibrating, and validating the SWAT model.However, the geodatabase and the monitoring activities developed in this study can be applied to the most commonly used hydrological and water quality models (i.e., AnnAGNPS).Model application results will be analyzed in a further research paper.
In order to simulate the physical processes associated with water, sediment, and nutrients, the SWAT model divides a watershed into sub-basins, which are then further subdivided into hydrologic response units (HRUs) that consist of areas with homogenous slopes, land uses, management practices, and soil characteristics [25].
There are two separate phases in the SWAT model for simulating hydrology and water quality.One is the landscape phase, which simulates the amount of water, sediment, and nutrient loading from each sub-basin sent to the main stream, and the other is the in-stream phase, which describes the water, sediment, and nutrient movement through the stream system to the outlet.The landscape phase requires a digital elevation model, a land use map, information regarding agricultural practices, a soil map and soil hydrological data, as well as river network and weather data.The calibration and validation of the model require that streamflow and water quality parameters be measured.
The model generates several output files and results are aggregated by different spatial (basin, sub-basin, reach, HRU) and temporal scales (daily, monthly, yearly).Thus, it allows the critical HRUs to be identified in terms of soil and nutrients losses, the nutrient loads delivered to the outlet by the river network over a long period to be quantified, and the contribution of point sources to be evaluated.

Spatial Analysis for Building a Geodatabase
A spatial analysis was performed to build the geodatabase needed for the model set up.Existing official data on topography, river network structure, land use, soil, rainfall, temperature, solar radiation, agricultural practices, and point-source pollution (e.g., discharges of urban wastewater treatment plants, WWTPs) were collected and analyzed.It should be noted that these data were derived from different sources.Indeed, in recent years, several studies have been conducted in the Taranto province by the Regional Environmental Agency, municipalities, the River Basin Authority, and the Civil Protection Service.Unfortunately, the results and data from these studies were not archived in a unique database.Moreover, they are not uniform in terms of time scale and spatial resolution.
The sources of all the data collected and used for building the geodatabase are provided in Appendix A. If replicated data utilizing different spatial resolutions were available from different sources, those with the highest resolution were selected for and implemented in the geodatabase.In addition, meteorological data time series frequently include several gaps.Thus, data from nearby meteorological stations were used to fill in any missing data regarding precipitation, temperature, and/or solar radiation.
Furthermore, the reliability of collected data concerning land use, soil type, and agricultural practices was improved by means of interviews with farmers and agricultural advisors, published studies, and field inspections.
In particular, the land use map provided by the Puglia Region was refined by means of the National Agricultural Census results [26].All of the agronomic data needed by the hydrological models (such as type, timing and amount of fertilizers used for each crop, annual crop yields, tillage operations, irrigation supply, and grazing) were identified through the analysis of the interviews (Appendix B).Both farmers and agricultural advisors were selected precisely in order to gain information covering the whole basin [27,28].Nutrients from animal manure and grazing used to fertilize fields were quantified by multiplying the live weight of each animal type by its animal-specific total nitrogen (TN) excretion rate [29].Livestock densities were available from the National Agricultural Census on a municipal scale [26].A distinction between indoor and outdoor farming was made.TN loss during manure handling and storage (27.5%) was considered in the case of manure produced by indoor farming [30].
Two databases providing soil data at the regional [31] and European levels [32] were analyzed and cross-checked.In order to calculate the soil hydraulic properties needed by the hydrological model (e.g., hydraulic conductivity, available water capacity, and bulk density), the software Soil Water Characteristics developed by the United States Department of Agriculture (USDA) Agricultural Research Service was used [33].The effects of the karstic geology on surface runoff and percolation were examined through the analysis of published studies concerning the study area [19].
Field inspections in collaboration with the Reclamation Consortium of Stornara and Tara were carried out in order to obtain an actual river network representation and understand the hydraulic operation of the tributaries (i.e., Salina Grande).Therefore, detailed land use and soil maps and data on soil hydrological parameters, river network conformation, agricultural practices, climate, and point sources were gathered.

Field Data
The calibration and validation of hydrological and water quality models require that the model performance be evaluated for both hydrology (streamflow) and water quality.The appropriate performance criteria should be selected in advance [34,35]; consequently, measurements of streamflow and water quality are needed, possibly covering several years.
In the study area, historical data on water quality were not available.Hence, a monitoring scheme was defined while keeping in mind the final aim of the monitoring, i.e., characterizing the basin in terms of hydrology and water quality and creating a basic dataset for model calibration and validation (streamflow and water quality).Several aspects were considered before selecting the locations for and the number of gauging stations.The study budget, environmental characteristics of the area, and accessibility, as well as the visibility of the sites in order to avoid acts of vandalism, were the most important factors influencing the choices made [36].
When the time series for the streamflow is short, as in this case study, an alternative that can be utilized is the "split-in space" calibration strategy, which is based on the measurements recorded at a number of gauging stations [9].Thus, at the beginning of the monitoring activity, it was decided to monitor streamflow and water quality at two gauging stations (A, B) located on the main course of the stream and one (C) located on a tributary, taking into account that the sum of the nutrient loads measured in sections A and B provided nearly the entire load delivered to the Mar Piccolo by the Canale d'Aiedda (Figure 1).
However, after the first monitoring campaign at these three stations, new river sections were selected for monitoring the water quality in order to integrate data provided by local authorities (e.g., discharges from springs and WWTPs) and verify whether or not there were any further point sources beyond those provided by local authorities [37] discharging into the river network.The treated wastewater discharged by the WWTPs of Monteiasi (about 2 km upstream of B) and San Giorgio Ionico (about 3 km upstream of A) was sampled (W2, W3), as was the spring water (S3).In addition, samples were also taken from the river section upstream of the two WWTPs (O1, O2), the basin mouth section (O4), and a river section (O3), which was 250 m upstream of O4 and not subject to Mar Piccolo saltwater intrusion.Thus, surface water samples were defined on a fortnightly or monthly time scale for 10 sections (A, B, C, W2, W3, S3, O1, O2, O3, O4-Figure 1).Point sources upstream from the river sections O1 and O2 were not sampled since O1 and O2 were found to be dry for most of the year.Moreover, it was not possible to take water samples in the Salina Grande river network (partially buried) since its open channels were found to be dry during all of the field inspections except for the open channel that received water from a spring (S3).
The concentrations of total nitrogen (TN), nitrate (N-NO 3 ), nitrite (N-NO 2 ), ammonia (N-NH 4 ), total phosphorus (TP), mineral phosphorus (P-PO 4 ), and total suspended solids (TSS) were determined.The analyses were carried out in the laboratories of the Agricultural and Environmental Science Department (DISAAT) of the University of Bari Aldo Moro and the Water Research Institute-National Research Council (IRSA-CNR, Bari) according to the Agency for Environmental Protection and Technical Service's (APAT-IRSA/CNR) standard analytical method [38].
For streamflow measurements, a data logger (MDS Dipper-PT, © 2019 Seba Hydrometrie, Kaufbeuren, Germany) was installed and programmed following the operation guide (Figure 2).The sensor measured water level and temperature.It was not possible to employ an event-based sampling strategy by means of automatic sampler [16] due to limited economical resources.
The water level was registered with a sub-hourly time scale whenever it was lower than 0.30 m.However, the time step was reduced to 15 min during flood events (water level greater than 0.30 m).
Using these water level measurements, streamflow was estimated using the Gauckler-Strickler equation [39].A survey of the concrete river sections (A, B, and C) was carried out, and water velocity measurements were taken (taking into consideration different flow conditions) via a flow module (ISCO 750 Area-Velocity Flow Module, © 2018 Teledyne ISCO, Lincoln, NE, USA) in order to calibrate the equation's parameters (i.e., roughness coefficients and channel slope).When the water level was lower than 0.04 m and it was not possible to use the flow module, float measurements were performed [40].The water level was registered with a sub-hourly time scale whenever it was lower than 0.30 m.However, the time step was reduced to 15 min during flood events (water level greater than 0.30 m).
Using these water level measurements, streamflow was estimated using the Gauckler-Strickler equation [39].A survey of the concrete river sections (A, B, and C) was carried out, and water velocity measurements were taken (taking into consideration different flow conditions) via a flow module (ISCO 750 Area-Velocity Flow Module, © 2018 Teledyne ISCO, Lincoln, NE, USA) in order to calibrate the equation's parameters (i.e., roughness coefficients and channel slope).When the water level was lower than 0.04 m and it was not possible to use the flow module, float measurements were performed [40].

Riverine Nutrient Load Export
Nutrient load measurements are required for model calibration and validation.The calculation of nutrient loads (L) passing through a river section in a time interval is not easy when streamflow measurements are continuous (daily) and nutrient concentrations are discrete (a few measurements per month), as in this case.However, several methodologies have been developed [41][42][43][44] that require a continuous streamflow time series and discrete values of nutrient concentrations.In this study, the "Loads Tool" [45] was used to estimate monthly and annual load with different methods.Previous studies [42,44,46,47] report criteria and methods for estimating loads.Based on the characteristic of the dataset (fortnightly or monthly time scale), three methods were selected to calculate loads.Method 1 uses the inter-sample mean concentration (Equation ( 1)), Method 2 calculates the intersample mean concentration using mean flow (Equation ( 2)), and Method 3 requires a linear interpolation of the concentration (Equation (3)).Following the suggestions reported by Quilbé et al. [42], regression methods were not applied due to the low correlation between nutrient concentration and streamflow (R 2 < 0.3).
The inter-sample mean concentration technique uses averaged concentration values to estimate concentrations on non-sampled days and streamflow measured daily.
where Cj is the jth sample concentration, Qj is the jth streamflow, and n is the number of concentration measurements.

Riverine Nutrient Load Export
Nutrient load measurements are required for model calibration and validation.The calculation of nutrient loads (L) passing through a river section in a time interval is not easy when streamflow measurements are continuous (daily) and nutrient concentrations are discrete (a few measurements per month), as in this case.However, several methodologies have been developed [41][42][43][44] that require a continuous streamflow time series and discrete values of nutrient concentrations.In this study, the "Loads Tool" [45] was used to estimate monthly and annual load with different methods.Previous studies [42,44,46,47] report criteria and methods for estimating loads.Based on the characteristic of the dataset (fortnightly or monthly time scale), three methods were selected to calculate loads.Method 1 uses the inter-sample mean concentration (Equation (1)), Method 2 calculates the inter-sample mean concentration using mean flow (Equation ( 2)), and Method 3 requires a linear interpolation of the concentration (Equation ( 3)).Following the suggestions reported by Quilbé et al. [42], regression methods were not applied due to the low correlation between nutrient concentration and streamflow (R 2 < 0.3).
The inter-sample mean concentration technique uses averaged concentration values to estimate concentrations on non-sampled days and streamflow measured daily.
where C j is the jth sample concentration, Q j is the jth streamflow, and n is the number of concentration measurements.
The calculation of the inter-sample mean concentration using the mean flow method assumes that the average daily concentration on non-sampled days is calculated by a simple average of the concentrations from the last sample data point and the next sample data point.The flow is assumed to be the average flow up to the next sampled concentration value.
Water 2019, 11, 267 where Q <j+1 is the average flow to the end of the j + 1 period.The linear interpolation concentration technique assumes that the concentration on non-sampled days is determined by interpolating linearly between fortnightly or monthly sampled concentrations.
where Q j is the inter-sample mean flow evaluated between periods j and j + 1.
It is expected that monthly loads will differ from one method to another.Indeed, the load estimate is affected by an overall level of uncertainty due to different sources, such as the estimation technique employed (defined as knowledge uncertainty), the variability of sampled measurements (stochastic uncertainty), and unrepresentative measurements or sampling (measurement uncertainty).Hence, we assume here that the mean value acted as a measure of nutrient monthly load and that only the knowledge uncertainty was estimated, which was expressed as the coefficient of variation (ratio of standard deviation to the mean).

Monitoring Streamflow and Water Quality
The Canale d'Aiedda stream has a flow regime characterized by typical semi-arid features with a seasonal pattern of a low-flow period or drought from May to October and a wet season from November to April.Intense storms occurring in summer often lead to intense flash floods [4,48,49].
The mean daily flow measured at the A and B gauging stations from August 2017 to July 2018 is shown in Figure 3. From August 2017 to July 2018, the total amounts of precipitation recorded at the San Giorgio Ionico and Grottaglie gauging stations were 530.8 mm and 544.8 mm, respectively, corresponding to 90.6 % and 89.8%, respectively, of the average values for these stations (calculated from 2000 to 2013).Concerning the water quality (Figure 3) during the dry period, TN and N-NO3 concentrations were found to be sensibly lower than those registered during the wet period (winter and spring) since, due to the rainfall events, an increase in the transport of nutrients to the stream channels occurred.These results confirm that the Canale d'Aiedda has the behavior typical of a Mediterranean intermittent stream [4,50].In summer, nutrient loads are mostly due to point sources, as rainfall is limited to a few sporadic events, as Figure 4 shows.Due to tourism, the volume of the effluent is higher than that discharged in winter.In the wet season, the contribution of diffuse sources can be relevant depending on flood events.In section A of the river, TP and P-PO4 concentrations did not present a pattern linked to the hydrological regime.In section B, in the dry season, the concentrations were higher than those recorded in the wet period, for which it seems that a dilution effect contributes to a reduction in concentrations.The correlation coefficient between concentration and streamflow was 0.22 for N-NO3 and 0.06 for P-PO4 in section A. The concentration of N-NO3 and streamflow were not correlated in section B; however, after excluding an outlier, the correlation was 0.58 (Figure 5).Meanwhile, the P-PO4 concentration was negatively correlated with streamflow (−0.2897).The corresponding coefficient of determination (R 2 ), as reported in Figure 5, indicates that most of the variability in the concentrations of both N-NO3 and P-PO4 cannot be explained by streamflow but instead depends on many other factors, such as point source discharges, rainfall intensity and its spatial distribution, vegetation cover, and agricultural practices.The Canale d'Aiedda river network shows a natural intermittent character due to the climate and the geological features of the region.However, the natural flow regime was recorded only at gauging station C, for which no point sources were identified upstream.Conversely, a continuous flow was registered at both gauging stations A and B due to the presence of inlets (WWTP discharges); therefore, the hydrological regime of these river segments was altered heavily.The hydrographs recorded at the three gauging stations show a short lag time (time between peak rainfall and peak discharge) and a very rapidly rising stage.Flow regimes reflect the precipitation patterns.Flood duration is typically only a few hours, during which streamflow increases and decreases rapidly.This information is relevant for setting hydrological parameters in the model project.The river sections O1 and O2 (upstream of the WWTPs) were found to be dry for most of the year.
Concerning the water quality (Figure 3) during the dry period, TN and N-NO 3 concentrations were found to be sensibly lower than those registered during the wet period (winter and spring) since, due to the rainfall events, an increase in the transport of nutrients to the stream channels occurred.These results confirm that the Canale d'Aiedda has the behavior typical of a Mediterranean intermittent stream [4,50].In summer, nutrient loads are mostly due to point sources, as rainfall is limited to a few sporadic events, as Figure 4 shows.Due to tourism, the volume of the effluent is higher than that discharged in winter.In the wet season, the contribution of diffuse sources can be relevant depending on flood events.In section A of the river, TP and P-PO 4 concentrations did not present a pattern linked to the hydrological regime.In section B, in the dry season, the concentrations were higher than those recorded in the wet period, for which it seems that a dilution effect contributes to a reduction in concentrations.The correlation coefficient between concentration and streamflow was 0.22 for N-NO 3 and 0.06 for P-PO 4 in section A. The concentration of N-NO 3 and streamflow were not correlated in section B; however, after excluding an outlier, the correlation was 0.58 (Figure 5).Meanwhile, the P-PO 4 concentration was negatively correlated with streamflow (−0.2897).The corresponding coefficient of determination (R 2 ), as reported in Figure 5, indicates that most of the variability in the concentrations of both N-NO 3 and P-PO 4 cannot be explained by streamflow but instead depends on many other factors, such as point source discharges, rainfall intensity and its spatial distribution, vegetation cover, and agricultural practices.
Measured N-NO 3 and P-PO 4 concentrations in ten river sections along the Canale d'Aiedda river network are summarized in Figure 6.From August 2017 to July 2018, 17 campaigns were planned and implemented, and a large number of samples were analyzed (n = 127).Relatively high concentrations of P-PO 4 (greater than 0.1 mg L −1 ) were found in several sections.Given that phosphorus concentrations greater than 0.1-0.2mg L −1 are generally perceived to be sufficiently high as to result in freshwater eutrophication [51], the high values observed, especially at A, B, and O3, were of particular concern.In these river sections, the N-NO 3 concentrations were found to be relatively high as well.However, more than 20% of rivers with concentrations exceeding 3.6 mg N L −1 are found in many EU countries [51].Measured N-NO3 and P-PO4 concentrations in ten river sections along the Canale d'Aiedda river network are summarized in Figure 6.From August 2017 to July 2018, 17 campaigns were planned and implemented, and a large number of samples were analyzed (n = 127).Relatively high concentrations of P-PO4 (greater than 0.1 mg L −1 ) were found in several sections.Given that phosphorus concentrations greater than 0.1-0.2mg L −1 are generally perceived to be sufficiently high as to result in freshwater eutrophication [51], the high values observed, especially at A, B, and O3, were of particular concern.In these river sections, the N-NO3 concentrations were found to be relatively high as well.However, more than 20% of rivers with concentrations exceeding 3.6 mg N L −1 are found in many EU countries [51].
Nutrient concentrations measured at the treated wastewater sites (W2, W3) and the downstream river sections (A, B, O3, O4) were found to be highly variable over time (Figure 6).Measured concentrations in the O4 river section (basin mouth) were lower than those recorded at O3, which was located 250 m upstream, due to the effect of sea intrusion (dilution effect).Differences in N-NO3 and P-PO4 concentrations were found between the river sections located upstream and downstream of the WWTPs (W2 and O1, W3 and O2, respectively).In section W2, as Figure 6 shows, the median value of the P-PO4 concentrations was higher than that measured at O1, as well as the variability recorded over time.Indeed, the interquartile range (75th-25th quartile) and 95th percentile were much higher at W2 than those measured at O1.The variability in nutrient concentrations registered at the WWTP outlets (W2, W3), which depends on the treatment efficiency at the time of sampling, influences the water quality of the downstream river sections.Hence, at the gauge stations A and B, nutrient concentrations showed a higher variability (high 95th percentile and high interquartile range) than those recorded for the sampling sites C, O1, and O2, which are not fed directly by WWTPs.
This different degree of variability in concentrations among the river sections is much more evident for P-PO4.As Figure 6b shows, the concentrations at C, O1, and O2 exhibited low variability.However, N-NO3 concentrations showed a certain degree of variability in the river sections C, O1, and O2 due to agriculture being the most relevant source of nitrates and a number of other factors (e.g., agricultural practices, climate, vegetation cover, hydrological regime) influencing its export and riverine transport [50].Thus, in this basin, TP and P-PO4 are influenced mainly by point sources, and agriculture is a secondary source of phosphorus.
The high variability of nutrient concentrations related to the quality of WWTP discharge is an important piece of information that must be considered since it could introduce a source of uncertainty in the modelling results.Hence, in order to perform satisfactory model calibration and validation, detailed data concerning point sources discharge (flow and load) are needed.
The streamflow in the river section S3 (spring) is sustained by deep groundwater [19] and, consequently, the surface water quality is influenced by the groundwater regime.P-PO4 concentrations here were found to be the lowest recorded in the watershed.Meanwhile, the median N-NO3 concentration (5.25 mg L −1 ) is indicative of groundwater pollution.Data concerning nitrate contamination of groundwater published for the Apulia Region confirmed the high level of nitrates [52].Given the results obtained from the monitoring activities, the presence of other point sources in addition to those provided by local authorities can be excluded.

Setting Inputs for Hydrological and Water Quality Models
The geodatabase derived from the spatial analysis was built by including all of the input data required by the SWAT model, but it is also able to cover the input requirements of other models (i.e., AnnAGNPS).Table 1 summarizes the procedures adopted for developing the geodatabase.Figure 7a shows the reclassified land use map, and Figure 7b shows the soil map to be used as model input with the exception of those areas (e.g., outcropping limestone) where percolation is prevalent with respect to the superficial runoff.Thus, the basin (408 km 2 ) was clipped appropriately to include only the area contributing to surface runoff, which resulted in a total area of 309 km 2 .Moreover, the basin was divided into sub-basins based on the dominant hydrological processes.With this aim, the drainage area of the Salina Grande (Basin S, 86.9 km 2 ) was separated from the remaining part of the basin (Basin N, 222.1 km 2 ) due to the peculiarities of this sub-basin, as characterized by the presence of buried channels.For this reason, two different hydrological simulations were carried out for these two basins.

Input Description Format Digital Terrain Model
The Digital Terrain Model (DTM) provided by the Puglia Region was used [53].grid

River network
The stream features provided by the Hydro-Geomorphological and Regional Technical maps [53] were cross-checked and integrated with information obtained by means of field inspections.The actual river network representation was obtained with particular reference to the interaction between the Salina Grande and the Canale d'Aiedda streams.Nutrient concentrations measured at the treated wastewater sites (W2, W3) and the downstream river sections (A, B, O3, O4) were found to be highly variable over time (Figure 6).Measured concentrations in the O4 river section (basin mouth) were lower than those recorded at O3, which was located 250 m upstream, due to the effect of sea intrusion (dilution effect).Differences in N-NO 3 and P-PO 4 concentrations were found between the river sections located upstream and downstream of the WWTPs (W2 and O1, W3 and O2, respectively).In section W2, as Figure 6 shows, the median value of the P-PO 4 concentrations was higher than that measured at O1, as well as the variability recorded over time.Indeed, the interquartile range (75th-25th quartile) and 95th percentile were much higher at W2 than those measured at O1.The variability in nutrient concentrations registered at the WWTP outlets (W2, W3), which depends on the treatment efficiency at the time of sampling, influences the water quality of the downstream river sections.Hence, at the gauge stations A and B, nutrient concentrations showed a higher variability (high 95th percentile and high interquartile range) than those recorded for the sampling sites C, O1, and O2, which are not fed directly by WWTPs.
This different degree of variability in concentrations among the river sections is much more evident for P-PO 4 .As Figure 6b shows, the concentrations at C, O1, and O2 exhibited low variability.However, N-NO 3 concentrations showed a certain degree of variability in the river sections C, O1, and O2 due to agriculture being the most relevant source of nitrates and a number of other factors (e.g., agricultural practices, climate, vegetation cover, hydrological regime) influencing its export and riverine transport [50].Thus, in this basin, TP and P-PO 4 are influenced mainly by point sources, and agriculture is a secondary source of phosphorus.
The high variability of nutrient concentrations related to the quality of WWTP discharge is an important piece of information that must be considered since it could introduce a source of uncertainty in the modelling results.Hence, in order to perform satisfactory model calibration and validation, detailed data concerning point sources discharge (flow and load) are needed.
The streamflow in the river section S3 (spring) is sustained by deep groundwater [19] and, consequently, the surface water quality is influenced by the groundwater regime.P-PO 4 concentrations here were found to be the lowest recorded in the watershed.Meanwhile, the median N-NO 3 concentration (5.25 mg L −1 ) is indicative of groundwater pollution.Data concerning nitrate contamination of groundwater published for the Apulia Region confirmed the high level of nitrates [52].
Given the results obtained from the monitoring activities, the presence of other point sources in addition to those provided by local authorities can be excluded.

Setting Inputs for Hydrological and Water Quality Models
The geodatabase derived from the spatial analysis was built by including all of the input data required by the SWAT model, but it is also able to cover the input requirements of other models (i.e., AnnAGNPS).Table 1 summarizes the procedures adopted for developing the geodatabase.Figure 7a shows the reclassified land use map, and Figure 7b shows the soil map to be used as model input with the exception of those areas (e.g., outcropping limestone) where percolation is prevalent with respect to the superficial runoff.Thus, the basin (408 km 2 ) was clipped appropriately to include only the area contributing to surface runoff, which resulted in a total area of 309 km 2 .Moreover, the basin was divided into sub-basins based on the dominant hydrological processes.With this aim, the drainage area of the Salina Grande (Basin S, 86.9 km 2 ) was separated from the remaining part of the basin (Basin N, 222.1 km 2 ) due to the peculiarities of this sub-basin, as characterized by the presence of buried channels.For this reason, two different hydrological simulations were carried out for these two basins.

Digital Terrain Model
The Digital Terrain Model (DTM) provided by the Puglia Region was used [53].grid

River network
The stream features provided by the Hydro-Geomorphological and Regional Technical maps [53] were cross-checked and integrated with information obtained by means of field inspections.The actual river network representation was obtained with particular reference the interaction between the Salina Grande and the Canale d'Aiedda streams.
shp Mask Map of areas contributing to surface runoff.Areas with outcropping limestone, caves, and sinkholes (Figure 7b) were identified throughout the basin, including with information from Zuffianò et al. [19]. grid

Land use map
The land use map provided by the Puglia Region [53] was reclassified by considering the crop data provided on a municipal scale by the National Agricultural Census [26].In particular, arable land was reclassified and divided between durum wheat, set-aside land, and herbage.Thus, a detailed land cover map was obtained, and 23 soil use classes were identified throughout the basin (Figure 7a).grid

Soil map
The 20 polygonal features of the soil map provided by the Puglia Region [53] were further divided, taking into account the differences between the percentages of coarse content provided by JRC-ESDAC [32].Hence, 24 soil types were identified throughout the basin. grid

Soil database
The soil information provided by the Puglia Region for each soil type (e.g., number of soil layers, soil layer depth, texture, organic carbon content) [31] was integrated with information provided by JRC-ESDAC (e.g., rock fragment content) [32] and physical parameters computed with the software Soil Water Characteristics (e.g., bulk density of each soil layer, available water capacity, saturated hydraulic conductivity) [33].mdb

Point sources
Groundwater database [54] and urban wastewater treatment plant discharge databases [37] were used and integrated with the results of the surface water sampling campaign (points W2, W3, S3) in order to obtain information on the specific point sources (springs and WWTP discharges) present in the study area.The mean annual volume of treated sewage was provided for each plant by the plant manager (Acquedotto Pugliese, AQP), and a new database was built (Table 2).mdb

Meteorological data
Daily data (precipitation, temperature, solar radiation, wind speed and direction, humidity) were acquired from the Civil Protection Service-Puglia Region (13 meteorological stations) [55] and Assocodipuglia Consortium (5 meteorological stations) [56].The coordinates were related to climatic data, and a new database was built.mdb Agricultural practices A database was implemented including for each crop (land use) agricultural practices, such as tillage operations (type and date) as well as fertilizer and irrigation applications (i.e., timing, amount, and type).Grazing was also considered on pasture land.Interviews with farmers and agricultural advisors were also used (Appendix B).
mdb Water 2019, 11, 267 Agricultural practices A database was implemented including for each crop (land use) agricultural practices, such as tillage operations (type and date) as well as fertilizer and irrigation applications (i.e., timing, amount, and type).Grazing was also considered on pasture land.Interviews with farmers and agricultural advisors were also used (Appendix B).  Figure 8 shows the boundaries of the two basins, as well as the river network, the point sources, and the weather stations included in the geodatabase.
Water 2019, 11, x FOR PEER REVIEW 13 of 23 Figure 8 shows the boundaries of the two basins, as well as the river network, the point sources, and the weather stations included in the geodatabase.Constant values for daily volume and nutrient load were considered for each spring.Meanwhile, daily loads from WWTP discharges were considered to be not constant.For each nutrient compound, the load was computed by multiplying the daily treated water volume (a constant value provided by AQP) by the sampled concentration.The latter was measured in our sampling campaign or the value published by Arpa Puglia [37] was considered.For those days when measurements were missing, a mean monthly value was used.
Table 2 summarizes the mean daily volume of water discharged (Q) as well as the mean daily loads of nutrients (N-NO3, N-NH4, P-PO4) associated with the point sources.Constant values for daily volume and nutrient load were considered for each spring.Meanwhile, daily loads from WWTP discharges were considered to be not constant.For each nutrient compound, the load was computed by multiplying the daily treated water volume (a constant value provided by AQP) by the sampled concentration.The latter was measured in our sampling campaign or the value published by Arpa Puglia [37] was considered.For those days when measurements were missing, a mean monthly value was used.
Table 2 summarizes the mean daily volume of water discharged (Q) as well as the mean daily loads of nutrients (N-NO 3 , N-NH 4 , P-PO 4 ) associated with the point sources.When considering agricultural practices, according to the information provided by the farmers, durum wheat is fertilized mainly in December and February, while vineyards are fertilized in November and February.Fertilizer operations in olive groves take place in April and August.Manure produced by animals that are free to graze is applied directly to pastures.Conversely, manure produced by indoor farming is spread on fallow land and herbage after a 90-day storage period, which aims to reduce the bacterial load.Irrigation is provided only during the dry season (i.e., from May to September) to vineyards, olive groves, and orchards.Drip irrigation systems are supplied by groundwater pumping wells.Despite the fact that the current study gathered high-quality data due to the collaboration of agricultural advisors and farmers, it was not possible to assign actual agricultural practices to each field due to the high fragmentation of land use and the large differences in adopted management practices within the study area.Therefore, the same averaged values were assigned to specific crops throughout Basin N and Basin S whenever site-specific information was not available.

Estimating Riverine Nutrient Loads
Monthly N-NO 3 and P-PO 4 loads were estimated in the river sections A and B (Figure 8).Their sum was nearly the entire amount of nutrients exported by Basin N and was delivered to the Mar Piccolo.Only the contribution of a 5.9 km 2 area was excluded.N-NO 3 load shows an intra-annual variability, with the highest value in April due to the rainfall events after fertilizer applications in both sections (Table 3).The temporal pattern of P-PO 4 differs from that of the nitrates (Table 4).It shows a high intra-annual variability, which seems to be related not to the streamflow regime (Figure 5) but rather to the WWTP discharges.
Annual N-NO 3 loads range between 13.044 t year −1 and 16.014 t year −1 at gauge A and between 24412 t year −1 and 30.118 t year −1 at gauge B, and the annual load of P-PO 4 was found to be in the range of 0.474 and 0.554 t year −1 at A and between 1.534 and 1.684 t year −1 at B. These values were computed by summing the monthly contributions, minimum and maximum, respectively, among the three above-mentioned averaging methods (par.2.5.) (Tables 3 and 4).Indeed, for each month, the three values can be considered to be samples from a probability distribution whose central value is the true unknown monthly load.
Based on the above-mentioned calculations, the average annual load of N-NO 3 delivered to the Mar Piccolo was 41.744 t year −1 , and the average annual load of P-PO 4 was 2.123 t year −1 .The contributions from point sources (Table 2) were 45% and 77% for N-NO 3 and P-PO 4 , respectively.When considering the monthly time scale, load estimates can be significantly different depending on the method used for computation [57].Tan et al. [41], who described the data requirements and applicability of the methods for load estimation, suggested using averaging techniques if no significant relationship between streamflow and nutrient concentration exists.However, it is important to take into account that the accuracy of the averaging methods depends on the number of samples and on the time between two successive samples [28].In this study, the number of samples (<100) and the lack of a significant relationship between concentration and streamflow (Figure 5) did not allow other methods (ratio or regression techniques) to be used or a data stratification strategy meant to improve the estimation to be applied [58].Only the three methods classified as averaging were used, which provided similar biased results on a yearly basis and generally underestimated loads if events were not adequately sampled [57].
The knowledge uncertainty associated with the monthly value (due to the methods used for determining the estimation), expressed here as the coefficient of variation, was very high both for N-NO 3 (up to 25%) and P-PO 4 (up to 40%) (Tables 3 and 4, respectively).In the model simulation, this uncertainty should be considered in the discussion of the results.
The knowledge uncertainty could be reduced by increasing the understanding of hydrological and pollutant transport processes.The results obtained here constitute a first assessment of nutrient loads, based upon which the monitoring program could be improved.
Figure 9 shows the contributions per hectare of the N-NO 3 and P-PO 4 monthly loads.
number of samples (<100) and the lack of a significant relationship between concentration and streamflow (Figure 5) did not allow other methods (ratio or regression techniques) to be used or a data stratification strategy meant to improve the estimation to be applied [58].Only the three methods classified as averaging were used, which provided similar biased results on a yearly basis and generally underestimated loads if events were not adequately sampled [57].
The knowledge uncertainty associated with the monthly value (due to the methods used for determining the estimation), expressed here as the coefficient of variation, was very high both for N-NO3 (up to 25%) and P-PO4 (up to 40%) (Tables 3 and 4, respectively).In the model simulation, this uncertainty should be considered in the discussion of the results.
The knowledge uncertainty could be reduced by increasing the understanding of hydrological and pollutant transport processes.The results obtained here constitute a first assessment of nutrient loads, based upon which the monitoring program could be improved.
Figure 9 shows the contributions per hectare of the N-NO3 and P-PO4 monthly loads.Monthly loads per unit area of N-NO 3 (kg ha −1 ), estimated for the sub-basin draining in gauge A (sub-basin A) (36.2 km 2 ), were higher than the loads computed for the area draining in gauge B (sub-basin B) (180.0 km 2 ) (Figure 9).This result was expected, as sub-basin B includes natural areas and forests (9%), and the agriculture is less intensive than that of sub-basin A. The productive lands for sub-basins A and B consist of about 76 and 78% of the total land, respectively.In both sub-basins, higher values for the N-NO 3 and P-PO 4 loads were recorded in winter and early spring.In particular, the highest monthly loads of N-NO 3 were recorded in April, with about 1.3 kg ha −1 and 0.6 kg ha −1 recorded at A and B, respectively, and was due to fertilizer applications in olive groves.In the winter, the specific monthly loads of P-PO 4 recorded at A were higher than the loads recorded at B; in the summer, except for the month of July, the specific P-PO 4 loads at B were higher than those at A. This behavior was attributable to the WWTP discharges.
The contributions per hectare of N-NO 3 and P-PO 4, computed for the area contributing to the surface runoff in both sections A and B (approximately the area of Basin N, excluded 5.9 km 2 ), were estimated to be 1.9 kg ha −1 year −1 and 0.1 kg ha −1 year −1 , respectively, which are quite low values if compared with results from other studies.The N-NO 3 load per unit area was in the range of 0-2 kg N ha −1 year −1 , as estimated by Bouraoui et al. [58] for the Apulia Region for the year 2000.However, similar studies reported a wide range of N losses, from 1.5-19 kg N ha −1 year −1 for an extensive agricultural watershed [59,60] to more than 100 kg N ha −1 year −1 for very intensive agricultural areas [61,62].Also, the P-PO 4 load per unit area was of the same order as those loads identified in studies conducted in the Apulia Region [58,63].The majority of the TP load per unit area estimates ranged between 0.05 and 0.65 kg TP ha −1 year −1 [60,63,64].

Improving Monitoring and Load Assessment
Point sources (WWTPs) and diffuse pollution from agriculture are the main anthropogenic pressures in the Canale d'Aiedda basin.Based on the first year of the monitoring activities carried out in the basin, an improvement in future monitoring activities could be made in order to save time and money.Considering the main aim of future modelling applications, which is to quantify the nutrient loads delivered to the Mar Piccolo by the Canale d'Aiedda river network, it seems necessary to monitor the nutrient concentrations in at least the river sections A and B where the streamflow is also currently recorded.Meanwhile, measurements of nutrient concentrations in section C are not required since it is located upstream of gauge B and the contribution of the drainage area upstream of section C to the nutrient load is limited to flood events.In addition, it seems necessary to characterize the WWTP discharges (W2, W3) through high-frequency measures of nutrient concentrations since the latter were found to be highly variable over time.The WWTPs in this basin constitute a relevant pressure affecting the water quality downstream of the discharges.
Concerning the sampling frequency, it seems that a few measurements per month are sufficient to estimate loads, but the uncertainty associated with these values can be significant.More accurate estimations could found by increasing the number of samplings in winter and spring.In small basins, most of the sediment and nutrient loads are delivered to the river during floods [50, 57,65], which should be monitored in order to increase the accuracy of the load estimations.However, doing so is an expensive task requiring the installation of automatic samplers with a specific program based on fixed water level changes or on fixed flow rates [16].Currently, economic resources and the choice of a good location at which to install the instrument (protected from acts of vandalism and the power of floods) are the major limiting factors to overcome in order to increase the monitoring in the study area.The possibility of using passive rising stage samplers in river sections A and B should be also evaluated.This sampling methodology allows one to take continuous sampling on the rising limb of flood events without requiring a power supply [66,67].These samples can then be collected manually several days later when access is possible.Thus, it may be a cheaper and effective alternative to grab and automated sampling.
In order to acquire data to calibrate and validate the model relative to Basin S (Figure 8), an additional gauging station should be installed that takes into account the fact that the river network is partially buried.The sampling strategy and the method used for load computations have a direct influence on model results since nutrient loads are used to calibrate and validate models.Consequently, the overall uncertainty affecting loads will influence the model results.

Conclusions
Hydrological and water quality models are used commonly in river basin management for quantifying water resources, identifying source areas of pollutants, and testing environmental measures designed to reduce pollution or mitigate climate change impacts on water resources.When the models are used in decision-making processes, which can have significant implications for humans and the environment, detailed inputs and a clear procedure for a correct application are necessary.Hence, before applying a model at the basin scale, it is necessary to define a modelling project and to build a geodatabase.This work is time-consuming and could require special efforts, especially in regions with limited data availability.The present study reports on a procedure for a spatial analysis that is aimed at supporting hydrological and water quality model applications through a case study of the Canale d'Aiedda basin, where existing data were collected, checked, and integrated with field campaigns and interviews of farmers, agricultural advisors, and citizens.
This study shows that existing datasets exhibit missing data and differing information.Hence, collecting, checking, and elaborating on data constitute a fundamental task when applying hydrological models.The spatial analysis developed in the GIS environment and the database was the result of merging the data with previous studies and field surveys, and new studies are a valid support for model applications.Continuous monitoring of streamflow in the Canale d'Aiedda resulted in high-quality information pertaining to this intermittent river system, compensating for the gaps in historical data and allowing modelling to be implemented using a "split-in space" calibration strategy while under time and budgetary constraints.
Basic water quality monitoring (four samples per year) is not adequate for calibrating models in small basins since streamflow and nutrient concentrations are highly variable.Monthly or fortnightly samples allow the estimation of monthly and annual loads by using averaging methods; however, the related uncertainties could be significant.The sampling strategy and the methods used for load computations have a direct influence on model results since measured nutrient loads are used to calibrate and validate models.An increase in the sampling frequency, especially in winter and spring, and a monitoring plan that includes flood events could reduce this uncertainty, even if an increase in

Water 2019 , 23 Figure 1 .
Figure 1.Study area: Canale d'Aiedda basin boundaries and location of point sources (i.e., spring and wastewater treatment plants), flow gauge stations, and sampling points.The extents of the Land and Sea contaminated Sites of National Interest (i.e., SIN Terra and SIN Mare) and of the Area of environmental crisis are also shown.

Figure 1 .
Figure 1.Study area: Canale d'Aiedda basin boundaries and location of point sources (i.e., spring and wastewater treatment plants), flow gauge stations, and sampling points.The extents of the Land and Sea contaminated Sites of National Interest (i.e., SIN Terra and SIN Mare) and of the Area of environmental crisis are also shown.
Water 2019, 11, x FOR PEER REVIEW 9 of 23

Figure 4 .
Figure 4. Monthly runoffs measured at gauges A (a) and B (b), and monthly rainfall measured at the

Figure 4
Figure 4 compares monthly rainfall recorded at the above-mentioned precipitation gauge stations with the monthly runoff measured at gauges A and B. During the study period, dry conditions were never registered at gauges A and B. The mean daily flow ranged from 0.017 to 0.181 m 3 s −1 (mean 0.045 m 3 s −1 ) at the gauge A (drainage area 36.2 km 2 ) and from 0.009 to 0.838 m 3 s −1 (mean 0.080 m 3 s −1 ) at gauge B (drainage area 180.0 km 2 ).Conversely, for gauging station C (drainage area 48.0 km 2 ), a completely dry condition was recorded in August 2017.The mean daily flow measured at gauge C ranged from 0 to 0.186 m 3 s −1 (mean 0.004 m 3 s −1 ).

Figure 4 .
Figure 4. Monthly runoffs measured at gauges A (a) and B (b), and monthly rainfall measured at the nearest gauges (i.e., San Giorgio Ionico station and Grottaglie station, respectively).

Figure 4 .
Figure 4. Monthly runoffs measured at gauges A (a) and B (b), and monthly rainfall measured at the nearest gauges (i.e., San Giorgio Ionico station and Grottaglie station, respectively).

Figure 5 .
Figure 5. Measured N-NO3 (a, c) and P-PO4 (b, d) concentrations versus measured streamflow at gauging stations A and B.

Figure 5 .
Figure 5. Measured N-NO 3 (a,c) and P-PO 4 (b,d) concentrations versus measured streamflow at gauging stations A and B. Water 2019, 11, x FOR PEER REVIEW 11 of 23

Figure 6 .
Figure 6.Box plots of N-NO3 (a) and P-PO4 (b) measured in ten river sections along the Canale d'Aiedda river network.The horizontal line within each box plot indicates the median; lower and upper boundaries of each box indicate the 25th and 75th percentiles, respectively.The lower and the upper whiskers indicate the minimum and the 95th percentile of the values for each box, respectively.

Figure 6 .
Figure 6.Box plots of N-NO 3 (a) and P-PO 4 (b) measured in ten river sections along the Canale d'Aiedda river network.The horizontal line within each box plot indicates the median; lower and upper boundaries of each box indicate the 25th and 75th percentiles, respectively.The lower and the upper whiskers indicate the minimum and the 95th percentile of the values for each box, respectively.

Figure 7 .
Figure 7. Reclassified land use map (a); map of soil texture and main karst units (b).Figure 7. Reclassified land use map (a); map of soil texture and main karst units (b).

Figure 7 .
Figure 7. Reclassified land use map (a); map of soil texture and main karst units (b).Figure 7. Reclassified land use map (a); map of soil texture and main karst units (b).

Figure 8 .
Figure 8.Initial basin (River Basin Authority), basins (Basin N, Basin S) contributing to surface runoff, point sources (springs and wastewater treatment plants), flow gauging stations, and weather stations.

Figure 8 .
Figure 8.Initial basin (River Basin Authority), basins (Basin N, Basin S) contributing to surface runoff, point sources (springs and wastewater treatment plants), flow gauging stations, and weather stations.

Figure 9 .
Figure 9. Mean monthly N-NO3 loads (kg ha −1 ) (a) and P-PO4 loads (kg ha −1 ) (b) estimated via the three different methods at gauging stations A and B. Error bars show the standard deviations.

Figure 9 .
Figure 9. Mean monthly N-NO 3 loads (kg ha −1 ) (a) and P-PO 4 loads (kg ha −1 ) (b) estimated via the three different methods at gauging stations A and B. Error bars show the standard deviations.

Table 1 .
Data included in the geodatabase.Researchers may obtain access to the geodatabase via the corresponding author.

Table 1 .
Data included in the geodatabase.Researchers may obtain access to the geodatabase via the corresponding author.

Table 2 .
Point source data.Q is the mean daily volume of water discharged into the Canale d'Aiedda.

Table 2 .
Point source data.Q is the mean daily volume of water discharged into the Canale d'Aiedda.S1, S2, W1, W2, and W3 are located in Basin N; S3, W4, and W5 are located in Basin S.

Table 3 .
Monthly N-NO 3 loads estimated by using the inter-sample mean concentration (Method 1), inter-sample mean concentration using mean flow (Method 2), and linear interpolation of concentration (Method 3) methods.CV is the coefficient of variation.

Table 4 .
Monthly P-PO 4 loads estimated by using the inter-sample mean concentration (Method 1), inter-sample mean concentration using mean flow (Method 2), and linear interpolation of concentration (Method 3) methods.CV is the coefficient of variation.