Seasonal and Storm Event‐Based Dynamics of Dissolved Organic Carbon (DOC) Concentration in a Mediterranean Headwater Catchment

This study investigates the spatial and temporal dynamics of Dissolved Organic carbon (DOC) concentration in a Mediterranean headwater catchment (Turbolo River catchment, southern Italy) equipped with two multi‐parameter sondes providing more than two‐year (May 2019–November 2021) continuous high‐frequency measurements of several DOC‐related parameters. The sondes were installed in two nested sections, a quasi‐pristine upstream sub‐catchment and a downstream outlet with anthropogenic water quality disturbances. DOC estimates were achieved by correcting the fluorescent dissolved organic matter—fDOM—values through an original procedure not requiring extensive laboratory measurements. Then, DOC dynamics at the seasonal and storm event scales were analyzed. At the seasonal scale, results confirmed the climate control on DOC production, with increasing background concentrations in hot and dry summer months. The hydrological regulation proved crucial for DOC mobilization and export, with the top 10th percentile of discharge associated with up to 79% of the total DOC yield. The analysis at the storm scale using flushing and hysteresis indices highlighted substantial differences between the two catchments. In the steeper upstream catchment, the limited capability of preserving hydraulic connection over time with DOC sources determined the prevalence of transport as the limiting factor to DOC export. In the downstream catchment, transport‐ and source‐limited processes were observed almost equally. The correlation between the hysteretic behavior and antecedent precipitation was not linear since the process reverted to transport‐limited for high accumulated rainfall values. Exploiting high‐resolution measurements, the study provided insights into DOC export dynamics in nested headwater catchments at multiple time scales.

be seen as "active pipelines" contributing to negative net ecosystem production (Cole et al., 2007), needs to be deeply investigated to improve the understanding of the global carbon cycle.
Hydrological factors are known to contribute to regulating the DOC balance at the reach scale (Bertuzzo et al., 2017;Parr et al., 2019).Interannual, intra-annual (seasonal) and event-based hydrological variability, particularly in headwater streams (Butman & Raymond, 2011;Rovelli et al., 2018), affects stream-hillslope organic matter exchanges and river network connectivity, leading to significant space and time variations in sources and processes regulating DOC dynamics.The impact of this interaction reflects on broader spatial scales so that, recently, regional approaches have been undertaken to evaluate the relationship between streamflow and DOC export regimes (Morison et al., 2022) or combine streamflow and DOC observations to validate catchment classification (Giesbrecht et al., 2022).
At different timescales, different processes emerge.At the seasonal scale, Viza et al. (2022) found that the intermittent flow regime of a Mediterranean river basin contributes to reducing organic matter decomposition rates.More generally, the effects of droughts on DOC transport have being extensively investigated (e.g., Ahmadi et al., 2019;Humbert et al., 2015;Mehring et al., 2013;Wu et al., 2022).Available studies have highlighted the inhibition of DOC release during low flow conditions owing to reduced network connectivity but higher DOC concentrations after droughts.
Several studies reported DOC concentration increase in the last decades (e.g., Monteith et al., 2007;Roulet & Moore, 2006;Wu et al., 2022).This trend is connected to rising temperatures that favor the DOC release (Bengtson & Bengtsson, 2007;Chen et al., 2021;Freeman et al., 2001;Zhong et al., 2020), an instance which establishes positive feedback with climate change (since DOC is eventually converted into CO 2 , a major greenhouse gas).Other causes, also linked to global warming, can concur to the observed increase of DOC export through the hydrological response, such as changes in land management, pH and sulfate, atmospheric CO 2 increase, acidic deposition decrease and runoff changes (Johnston et al., 2022;Ren et al., 2016;Worrall & Burt, 2007).
At short timescales of a few days or hours, storm events dominate DOC mobilization and transport (Parr et al., 2019).Precipitation activates direct wet deposition and indirect dry deposition deriving from vegetation canopy and stem (Song et al., 2021) and soil erosion (Chaplot & Mutema, 2021).This contribution to total DOC export is further emphasized if the wet event occurs at the end of prolonged dry periods (Blaurock et al., 2021).Fazekas et al. (2020) highlighted that anomalous events lasting overall less than 20 days in a year could define the annual behavior of the relationships between streamflow and organic matter concentration.
For a specific basin, the concentration-discharge (C-Q) relationship is a signature of the interactions between biogeochemical and hydrological processes, which in their turn depend on climatic, geological, and topographical features.C-Q relationships can reveal much of the DOC mobilization dynamics at different timescales (Chorover et al., 2017;Rose et al., 2018).At the seasonal or annual scale, null to low concentration variability in response to discharge fluctuations is called chemostasis (Basu et al., 2010;Godsey et al., 2009), indicating a homogeneous spatial distribution of DOC in the analyzed catchment.On the contrary, chemodynamic behavior identifies stronger dependence of solute concentration on streamflow (Fazekas et al., 2020;Musolff et al., 2015).This behavior is characterized by decreasing concentration with discharge if the DOC source is limited or, on the opposite, by increasing concentration with discharge if the limiting factor is the transport capacity.Under chemodynamic conditions, representing synchronous observations of DOC concentrations and flow rates in a C-Q plane yields hysteretic loops.At the event scale, the hysteretic loops' shape and direction help identify the main transport mechanisms.For example, if the DOC source is close and well-connected to the stream, clockwise hysteretic loops can be identified.On the contrary, counterclockwise loops prevail if it is far and connected by pathways with slow transport velocities (Williams, 1989).
The response of DOC dynamics is strictly connected to spatial features of heterogeneous ecosystems.Several studies showed that not only land cover type and land use (Aitkenhead-Peterson et al., 2007;Fovet et al., 2018;Seybold et al., 2019;Vaughan et al., 2017) but also local topography and geomorphic features (Weiler & McDonnell, 2006) significantly affect DOC mobilization and transport, influencing the hillslope-channel hydraulic connectivity (Botter et al., 2021).Therefore, the response in time of DOC dynamics in specific sections of a catchment is modulated by local properties of the upstream areas.Within the same catchment, significant differences can arise, which cannot be fully captured by a single downstream monitoring section that integrates heterogeneous upstream biogeochemical signals.It is a typical problem of scale (Lowe et al., 2006;McGuire 10.1029/2022WR034397 3 of 19 et al., 2014), which also affects streamwater chemistry and needs to be addressed with innovative theoretical concepts and technical approaches, including intensive spatially distributed monitoring campaigns in nested sections of the same catchment (Blaurock et al., 2021;McGuire et al., 2014).DOC dynamics monitoring across different spatial and temporal scales is possible thanks to the advancements in optical aquatic sensors technology.Through in situ continuous high-frequency measurements, such sensors catch rapidly changing concentrations during storm events and trends over more extended (seasonal to interannual) periods (Pellerin & Bergamaschi, 2014), supporting the development of accurate dynamic models (e.g., Jones et al., 2014) and, in general, providing great potential for a better understanding of aquatic ecosystems functioning (Snyder et al., 2018).Indeed, optical aquatic sensors do not measure DOC directly but rather the fraction of DOM that fluoresces (fluorescent dissolved organic matter, fDOM).fDOM data can be corrected by accounting for some physical properties of the water (e.g., Downing et al., 2012;Snyder et al., 2018;Watras et al., 2011) and related to DOC using laboratory measurements needed to calibrate the transfer function.Many studies exploit optical sensor properties integrated into multi-parameter sondes to highlight several features of coupled DOC-streamflow dynamics at different timescales.For example, Saraceno et al. (2009) analyzed a 4-week period including a short-duration storm event.Vaughan et al. (2017) and Fovet et al. (2018) focused on hysteresis in C-Q curves across many storms in catchments with different land use.Mistick and Johnson (2020)  This paper contributes to the ongoing effort to improve understanding of the related dynamics of streamflow and DOC concentration spatial variability across different timescales.Our investigation focused on a Mediterranean headwater catchment (Turbolo River, southern Italy) characterized by dry and hot summer climate enhancing network intermittency.The catchment was equipped with two multi-parameter sondes at two outlets, an upstream section closing a quasi-pristine sub-catchment and a downstream section receiving water from the upstream sub-catchment but closing a wider catchment area affected moderately by human activities (agriculture and villages).More than two-year (May 2019 to November 2021) continuous high-frequency measurements of several chemical-physical parameters were recorded, including DOC-related parameters like f DOM, streamwater temperature, and turbidity.On-site measurements were complemented by several samples collected during January-April 2021, aimed at characterizing the catchment and calibrating the f DOM-DOC transfer function.Furthermore, hydrometeorological observations were continuously performed, including discharge at the analyzed sections.
The study addresses the interrelated dynamics of DOC concentration, river discharge, and other hydrometeorological variables across multiple timescales in two nested sections of a Mediterranean headwater catchment characterized by different land uses.This general purpose was fulfilled through two specific objectives, which were addressed by exploiting a novel, simple procedure for the correction of recorded f DOM values that does not rely on extensive laboratory measurements: (a) the assessment of the seasonal variability of DOC background values related to several hydrometeorological parameters; (b) the evaluation of the DOC concentration-discharge relationships at the storm event timescale, considering season-and site-dependence, aimed at uncovering the main mobilization and transport mechanisms.For both the timescales considered in this study (storm event and seasonal), the difference in DOC response of the two nested cross-sections was analyzed to infer the dependence of DOC dynamics on scale properties and other landscape features.

Study Area
The study area is the upper Turbolo creek catchment (Figure 1), closed at the Fitterizzi gauge (183 m a.s.l., referred to as FIT from hereinafter), in southern Italy, a drought-prone area (Mendicino & Versace, 2007), for which an increasing drying scenario is projected (Senatore, Fuoco, et al., 2022).The catchment area is approximately 7 km 2 , with an elevation of up to 1,005 m a.s.l.The Turbolo creek originates from the Calabrian Coastal Range, dominated by strongly altered and fractured crystalline-metamorphic rocks that entail widespread slope instability and have overall high permeability.The geology allows ample groundwater recharge and storage that 10.1029/2022WR034397 4 of 19 sustains almost perennial flow at FIT. Steep slopes characterize the catchment on the metamorphic rocks in the west.In the eastern part, slopes are less steep but affected by water erosion processes, inducing shallow landslides and soil creep (Senatore et al., 2021).
The channel network consists of two main forks: the southern one (red-contoured on the map) is the San Nicola creek, 2 km 2 wide, in whose closing section (231 m a.s.l.) a gauging station was installed (referred to as SN from hereinafter).San Nicola is a quasi-pristine sub-catchment on which only an abandoned settlement is located.The average elevation of the catchment closed at the San Nicola gauge is 600 m a.s.l.(elevation range of 774 m, from 1,005 to 231 m a.s.l.).The northern fork reflects some more relevant anthropogenic effects, with more agricultural areas (mainly non-irrigated arable land and olive groves) and the village of Mongrassano with an adjoining sewage treatment plant sized for about 500 equivalent inhabitants.The average elevation of the catchment (closed at the Fitterizzi gauge) is 491 m a.s.l.(elevation range of 822 m, from 1,005 to 183 m a.s.l.).

Continuous-Time Monitoring and Manual Sampling of Water Quality
Several chemical-physical parameters were recorded continuously with an hourly time step starting from May 2019 by two multiparametric probes in the San Nicola and Fitterizzi sites.Monitoring is ongoing, but the results presented here relate to observations performed until 25 November 2021.Two YSI EXO2™ multi-parameter water quality sondes with 7 sensor ports, including a central wiper port, were used.Among the many parameters monitored (namely, ammonia NH 3 , specific conductance SpCond, Fluorescent Dissolved Organic Matter f DOM, nitrate ions  NO3− , dissolved oxygen DO, ammonium ions  NH 4 + , redox potential ORP, salinity Sal, total dissolved solids TDS, pH, temperature T, and turbidity), Table S1 in Supporting Information S1 provided in the supporting information shows details about the availability in time for both SN and FIT of f DOM, T, and turbidity, needed for DOC calculation.Over 11,500 water samples were gathered in each site, even though only 3,760 data were acquired simultaneously at both sections, during August-November 2019 and January-March 2021.Gaps are due to sensor malfunctioning (data out of scale or not recorded data).
Continuous observations were supported by discrete monitoring carried out in January-April 2021 when 59 samples were collected on-field in both outlets and analyzed in the laboratory according to Italian (APAT-IRSA) standards to achieve reference DOC values, used to correct measured f DOM values.Besides DOC measurements, several other physico-chemical parameters were measured to characterize the catchment in a larger picture and validate further continuous observations.Finally, meteo-hydrological observations carried out at both sites were available for the whole period.Specifically, a weather station and a water stage gauge were installed in FIT, managed by the Regional Agency for the Protection of the Environment (ARPACal), while a water stage gauge operated in SN until March 2021.

Correction of Measured f DOM Values
It has been shown that measurements of raw f DOM signal by field-deployable sensors need to be corrected by instrument-specific characteristics and are strongly affected by environmental conditions and site-specific characteristics.pH is one of the environmental factors that could affect fluorescent emission by dissolved organic matter (f DOM signal).However, in the typical pH range of natural aquatic systems (5 < pH < 9), it has a minimal impact on the f DOM measurements (Spencer et al., 2007).On the other hand, the f DOM signal decreases as a function of temperature due to the increase of non-radiative deactivation pathways.Thus, the raw f DOM data must be corrected to account for this effect.Watras et al. (2011) have shown that the fluorescent intensity measured by different fluorometers decreases linearly with the temperature in a relatively wide T range (0-40°C).Notably, the above authors showed that the ratio between the slope and the intercept of the lines describing the dependence of the f DOM value versus T, obtained at different concentrations of fluorescence organic matter, is almost constant within the experimental uncertainty.Therefore, they proposed the following compensation equation for the temperature effect: (1) In the above equation,   and   are measured and reference temperatures in °C, and ρ is a specific temperature attenuation coefficient (°C −1 ) equal to −0.01°C −1 (Exo User Manual, 2020).Equation 1shows that the compensation for T is restricted to ±5% in the 20-30°C range.Larger deviations occur for temperatures that differ more than 5°C from the reference temperature.We measured a range of temperatures between 6 and 26°C in the Turbolo basin waters.Therefore, Equation 1 corresponds to compensations of the raw f DOM value ranging from about −16% at 6°C to +1% at 26°C when the reference T is 25°C.
The other two most important phenomena that significantly attenuate the fluorescence signal of field-deployable sensors are the inner filtering effect and the suspended particles.The absorption of the excitation wavelengths directly causes the first one.The filtering effect is due to the organic matter dissolved in the water standing in the path between the excitation source and the fluorescence detector of the instrument (Carstea, 2012;Henderson et al., 2009) and indirectly by absorption of the fluorescence light emitted.This effect is significant only at concentrations of DOC larger than approximately 4 mg/L and exponentially increases with DOC following a Beer-Lambert law (Downing et al., 2012, and references therein).
Attenuation due to suspended particles is caused by light scattering and absorption (Saraceno et al., 2017, and references therein).This effect depends on the chemical nature and granulometric distribution of the suspended particles.It is relevant already at a concentration of suspended particles, measured in terms of turbidity, as low as 40 NFU.Attenuation of the f DOM signal increases exponentially with the sample turbidity following a Beer-Lambert law (Downing et al., 2012;Saraceno et al., 2017, and references therein).
To use f DOM as a proxy for the concentration of DOC, a linear equation could be set up, which reads: where f DOM corr is the corrected f DOM value, measured under ideal conditions, A and B are, respectively, the slope and intercept of the linear function, and DOC is the concentration of dissolved organic matter measured in the laboratory.Ideal conditions for this equation can be reasonably defined based on previous works as DOC < ∼4 mg/L and turbidity <∼40 FNU, that is, when the attenuation effects of the fluorescence due to inner filtering and suspended particles are both almost negligible.
If Equation 2 is valid, then the following proportionality relation holds: The   DOM  DOC ratio contains all the deviations of the measured f DOM corrected for T ( f DOM T ) from its ideal value ( f DOM corr ), caused by the above effects.
In our study campaign on the Turbolo basin, we detected min, max, and medium DOC values on 59 discrete samples collected under different meteorological conditions, respectively of 1.18, 3.88, and 2.2 (±0.8) mg/L.According to Downing et al. (2012), the inner filtering should have a negligible effect on the fluorescence signal at these DOC values.Thus, most of the f DOM deviations should be attributed to the attenuation effects caused by the suspended particles.These effects can be evaluated by plotting the above ratio as a function of the sample turbidity, which is a proxy for the concentration of suspended particles.Figure 2 shows that the f DOM signal non-linearly decreases with the turbidity.By inspection of the dispersion plot, the decrease follows an exponential trend, as expected from the Beer-Lambert behavior.Therefore, we modeled the above trend with an exponential function: where C is a pre-exponential parameter (ppb QSU/(mg L −1 )) and k is the particle suspended attenuation factor.The fitting returned R 2 = 0.59, and the following values for the two fitting parameters were achieved: C = 11.1 ± 0.5 ppb QSU/(mg L −1 ) and k = 0.004 ± 0.001 FNU −1 .Even though the attenuation factor is site-specific because it depends on the chemical composition and size distribution of the suspended particles, the value of the attenuation factor k fitted by our procedure is similar to those reported by others in previous works (Downing et al., 2012;Lee et al., 2015;Snyder et al., 2018), which were obtained using complex procedures requiring extensive laboratory measurements.(5) The corrected f DOM becomes, in this way, a proxy for the concentration of DOC in our natural aquatic system.In order to calculate the corresponding DOC value, the coefficients A and B of Equation 2 must be determined.
Figure 3 shows the raw f DOM data (corrected only for T) and those corrected applying Equation 5 on a data set of 59 matched samples.The raw f DOM data are randomly distributed in the whole DOC concentration range.In contrast, a good linear correlation was obtained between the measured DOC and the corresponding corrected f DOM values, with a Pearson coefficient r equal to 0.76 (R 2 = 0.58), which is relatively high considering that values between 0.3 and 0.8, were reported previously, in dependence on the aquatic environment analyzed (Snyder et al., 2018).
The linear fitting provides the following parameters, which were used to calculate the DOC values from the measured f DOM in situ: with c = 0.8 ± 0.2 mgC/L and m = 0.054 ± 0.007 mgC/L ppb QSU −1 .
For turbidity values much higher than those used to estimate the factor k, Equation 5 could not be applied.Therefore, for continuous measurements with turbidity values higher than 600 FNU, DOC values were extrapolated using multiple linear regressions of discharge and accumulated precipitation from 1 up to 12 previous hours, reaching R 2 values of 0.7 compared to observations.

Event Selection and Indices Calculation
A set of 29 focus events for each sub-catchment were identified during the study period following the approach proposed by Ladson et al. (2013).At first, the baseflow was separated from quickflow using the following basic filter equations: Where q f (i), q(i), and q b (i) are the quickflow, the streamflow, and the baseflow response at the ith sampling time (hourly), and α the filter parameter.This iterative method must be run multiple times (called passes) forward and backward.Thereafter, the Baseflow Index, defined as the ratio between the baseflow and the streamflow volume, was used to identify the events.The events selection based on this method was performed through the "hydroEvents" package in Software R using α = 0.925, while the appropriate number of passes to separate the baseflow for hourly data was chosen equal to 9, as suggested by Ladson et al. (2013).The method adopted not always detected events in both catchments during the same storm.Furthermore, some minor events, though automatically detected, were discarded due to the very low flow associated, especially in the small San Nicola creek.
The hysteresis index (HI) and the flushing index (FI) were calculated to evaluate the dynamics of DOC concentration in the analyzed basin, which is mobilized and transported by storm events.
The HI indicates a clockwise or counterclockwise behavior in the concentration-discharge (C-Q) relationship (Lloyd et al., 2016;Vaughan et al., 2017).For each event, the HI index was calculated starting from the normalized values of discharges and DOC concentrations: where Q i and C i are the discharge and the DOC concentration at the ith time step, Q min and Q max are the maximum and minimum discharge values, respectively, and C min and C max are the maximum and minimum DOC concentrations of the storm event.The density of the measured points on the C-Q plane was subsequently increased by linear interpolation considering 2% intervals between two adjacent measurements (i.e., 49 new points between two adjacent measurements).For the same intervals (called j), the hysteresis index HI j was calculated as follows: where C j,rising and C j,falling are the DOC concentrations in the rising and falling limb, respectively.The final hysteresis HI (ranging from −1 to +1) of each storm event was obtained by averaging all HI j values.Negative values indicate a counterclockwise hysteresis, while positive HI values indicate a clockwise hysteresis.Counterclockwise behavior in the concentration-discharge (C-Q) graphs means that DOC concentration is higher in the falling limb of the hydrograph than in the rising.This behavior generally occurs when the primary DOC sources are relatively far and hydraulically disconnected from the active river network at the beginning of the storm or when DOC transport is lower than water flux into the channel.Negative values of HI might also imply that the export process is transport-limited, that is, the process ceases because of the reduced transport capacity when the water drains toward the watercourse.On the other hand, positive HI values mean higher concentrations in the rising limb, the proximity of the major DOC sources to the active network, and source-limited process, that is, despite a still sustained water flow from the catchment to the river network, DOC concentration reduces in time.Therefore, the HI index allows quantifying event hysteresis dynamics, even with complex patterns (Williams, 1989) that are not easily interpretable.
The FI evaluates the increase of concentration, that is, flushing effect (positive values), or the decrease of concentration, that is, diluting (negative values) effect of DOC concentration on the rising limb (Butturini et al., 2008;Vaughan et al., 2017).The FI index is defined as: where C Qpeak, norm and C initial, norm are the normalized DOC concentrations at the peak of discharge and the beginning of the storm, respectively.

Results
In the following, the main results achieved in the two nested catchments were analyzed at different timescales.
First, the seasonal variability of the background (i.e., related to non-rainy days) DOC concentration was evaluated using data at the daily and monthly time steps.Then, the relationship between discharge and DOC export and between different DOC concentrations in the nested catchments was assessed by considering hourly data.Finally, the importance of storm events and related high flows on DOC dynamics was analyzed.

Seasonal Variability of Background DOC Concentration
The mean DOC concentration across the monitoring campaign (2019-2021) for FIT and SN were 1.7 ± 0.3 and 2.1 ± 0.5 mg l −1 , respectively.These are relatively low values in agreement with typical DOC concentrations in freshwater (≤5 mg l −1 ) (Stumm & Morgan, 1996).Concentrations between the two sites might not be directly comparable because the two recording periods do not entirely overlap.However, reducing the statistics only to the 91 non-rainy days with overlapping observations, DOC concentrations were higher in the SN section (1.9 ± 0.4 and 1.8 ± 0.3 mg l −1 for SN and FIT, respectively).Slightly higher background values for the upstream monitoring section can be due to the enhanced biomass production and decomposition and the steeper topography of the upstream forested catchment.
Box and whisker plots of the DOC showing seasonal trends for the background values (days in the absence of rain) are given in Figure 4 for FIT and SN sites.Table S2 in Supporting Information S1 adds descriptive statistics for the entire monitoring campaign.The analysis was conducted on average daily values, and the set of data corresponding to each site was also divided into four categories corresponding to calendar seasons-spring (Sp), summer (S), autumn (A), and winter (W).The two sites showed partially contrasting behavior.For FIT, which is the downstream site, the highest average DOC was recorded in autumn, with a descending order of A > S > W > Sp.In SN (the upstream site) the highest average DOC was observed in summer (S > A > W > Sp) instead of autumn, probably due to the greater availability of organic material in the surrounding pristine area that during the summer is converted into DOC by photochemical processes.The one-way ANOVA statistical analysis showed that the four data sets associated with each season were significantly different (p < 0.001 both in FIT and SN).To find potential correlations between DOC and key hydrometeorological parameters recorded in FIT, the data were further processed by multivariate data statistical analysis.Figure 5 shows the principal component analysis, including DOC and water temperature, air temperature, solar radiation, and discharge.For both sites, PC1 and PC2 altogether explained more than 80% of the total variance, and seasonal clusters could be distinguished with a clear temporal trajectory.Both PCA plots confirm the univariate statistical analysis presented before, indicating higher DOC concentrations were found mainly in summer and, in addition, provide us with further interesting insights.Indeed, while the seasonality observed for both sites reflects a multiparametric response rather than that of a single variable, the influential variables on the above patterns and their correlation differ between the sites.DOC and discharge are negatively correlated in SN while they lie in the same hyperplane in the PC model in FIT.It may effectively reflect differences in the limitation on DOC transport and release between the two sites, as already pointed out.In contrast, solar radiation is inversely correlated to DOC in FIT but not in SN.Interestingly, temperature and DOC seem to affect the multivariate pattern in both sites in the same way, thus highlighting the role of temperature as an alternative parameter to DOC for intra-annual variability of water quality.In summary, the DOC controlling factors do not have the same influence on the two sites; however, the temperature seems to contribute similarly to DOC seasonality.The PCA models strongly support the seasonality, but the differences highlighted above lead to different temporal trajectories in the two sites, with the DOC decreasing in a counterclockwise direction from summer to spring in FIT whilst clockwise in SN.The concentration effect the catchment undergoes during summer is one of the most important factors controlling the above patterns.However, discharge, solar radiation, and temperature may strongly affect the quantity and quality of the DOC through the activation of site-specific DOC sources and the enhancement of several photocatalytic degradation processes of dissolved organic matter (Stumm & Morgan, 1996).
Figure 6 provides more insights into the relationship between the DOC, discharge and air temperature at the monthly time scale.Specifically, Figures 6a and 6c show the seasonal hysteretic cycles in the Q-C plane for Fitterizzi and San Nicola.Data refer to the years 2021 and 2020, respectively, for which an almost complete series of data across the entire year were available.They confirm the expected pattern for the two sites under investigation, even though an anticipated DOC peak occurred in San Nicola, where the maximum DOC concentration was observed in late summer/early autumn.Figures 6b and 6d show the corresponding hysteretic loops for the monthly mean atmospheric temperature, highlighting counterclockwise patterns to discharge for both locations.
The peak temperature lagged about 1 month behind DOC in SN, while there was a larger inertia in the broader and lower elevation FIT catchment.

Hydrological Controls on DOC Export in the Nested Catchments
Most of the export through streamflow occurred relatively quickly, that is, during high flow events.Discharges above Q10 (i.e., the flow equaled or exceeded only 10% of the time) were responsible for 79% of the total yield in FIT (which is equal to 15.3 • 10 3 kg in 11,567 hr, i.e., 1.32 kg hr −1 ) and 69% in SN (4.4 • 10 3 kg in 11,732 hr, i.e., 0.37 kg hr −1 ). Figure 7 shows the steep slopes of the flow duration curves (hourly data, normalized to the maximum recorded hourly streamflow in both sites), emphasizing the high flow variability and the significant impact of quickflow.The corresponding normalized accumulated DOC load curves are represented in the same graph.These curves were created considering the DOC load values associated with each hourly streamflow value.The loads were accumulated starting from that corresponding to the lowest streamflow value, that is, from the lowest value of the flow duration curve.Therefore, to the lowest streamflow, the corresponding DOC load was associated; to the second lowest streamflow, the sum of the DOC loads of the lowest and second lowest streamflow was associated, and so on.The accumulated DOC curves were normalized to the total accumulated values at the end to allow comparison between the two sites.This representation makes it possible to directly display the percentage of total DOC export as a function of dimensionless duration.The normalized accumulated DOC load curves are also convex, even though they are less than the corresponding duration curves.The lower convexity of the SN DOC load curve compared to FIT DOC highlights the relatively more consistent DOC contribution in this more forested catchment with moderate to high discharges.
The behavior of DOC concentrations in FIT and SN was evaluated by comparing the 3,760 hourly data acquired simultaneously at both sections.Figure 8 shows that DOC concentration was generally higher in FIT (61% of the time) than in SN.This result, which overturns the indications obtained in the previous background DOC analysis showing higher non-rainy days concentrations in SN, highlights the importance of the rain/discharge events determining the highest concentration values.The total yield measured in the overlapping measurements period was equal to 3.32 and 0.61 kg hr −1 , respectively, for FIT and SN, leading to a ratio between the two total yields of 5.4, which is higher than the ratio between the two catchment areas, approximately equal to 3.5.The FIT catchment is characterized, overall, by a flatter topography.It is influenced by the significant contribution of the upstream northern fork, having different (less forested, more agricultural) land use and more relevant water erosion processes.Indeed, for the overlapping period, DOC concentrations were more positively correlated to the discharge observed in FIT than SN (correlation coefficient r equal to 0.65 and 0.34, respectively).This result is consistent with the behavior of the accumulated DOC load curves in Figure 7, confirming the more substantial impact of high flows in FIT on DOC yield.

Storm and High Flow Regulation on DOC Concentration Dynamics
The importance of storms and consequent high flows in the regulation of DOC export stands at the basis of the event-based analysis shown in the following.29 events were selected for each station, 19 representing the response produced in the two sections by the same storm.Mean and maximum discharge values were higher in FIT (average values of 0.627 and 1.517 m 3 s −1 , respectively) than in SN (average values of 0.208 and 0.487 m 3 s −1 , respectively).Consistently with results shown in Figure 8, mean and maximum DOC concentrations were higher in FIT (average values of 4.03 and 8.45 mg l −1 , respectively) than in SN (average values of 3.31 and 6.78 mg l −1 , respectively).Tables S3 and S4 in Supporting Information S1 provide further details and statistics on the selected events.
Analyzing the indices accounting for DOC concentration during storms helps understand the nature of the relevant mobilization and export processes driving DOC dynamics in stream water.Figure 9 provides a comprehensive overview of the catchment response during these events by comparing the two sites' hysteresis (HI) and flushing (FI) indices.In FIT and SN, positive FI values largely prevailed (only SN showed 5 out of 29 slightly nega tive values).A positive FI means that the DOC sources in the regions that contribute to the fast response of the catchment (root zone and riparian areas) are abundant enough to increase concentration when discharge increases (rising limb of the hydrograph).
HI results provided more contrasting information.In FIT, HI values were equally subdivided into positive and negative, with few values practically equal to zero.In SN, 8 values out of 29 were positive, but the general behavior was not very different from FIT. HI values for the two sections during the 19 simultaneous events were correlated quite well (r = 0.63).Overall, the slightly higher number of negative HI values in SN can be correlated to lower hydraulic connectivity of the upstream mountainous, steep catchment, presenting more accentuated flow spatial intermittency than FIT, as is typical for the headwater catchments.
Beyond the observed differences between the upstream and downstream sections, no general rules exhaustively explained the occurrence of clockwise/counterclockwise hysteresis during flow events.In SN, positive HI values were observed at the end of the winter, consistently with enhanced hillslope-channel hydraulic connectivity at the end of the wet season.Three positive HI events were detected from 13 February 2021 to 21 March 2021, also providing the highest loads due to the corresponding high discharges.Nevertheless, four events with HI > 0 occurred in early autumn (07 October 2019, 03 November 2019, 25 September 2020, and 15 October 2020) and one even in summer (08 August 2020).However, all these events were characterized by very low flows (the maximum peak flow overcame 0.1 m 3 s −1 only in one case), therefore providing a relatively low contribution in terms  of DOC load.The maximum discharges of the events were weakly negatively correlated to HI values (r = −0.21),while correlations with FI were stronger (r = −0.43).The latter result can be explained by considering that DOC removal might become supply-limited during high-intensity events.Indeed, the average maximum discharge of the five diluting events in SN was almost double (0.865 m 3 s −1 ) that of all events.
The positive HI values in FIT events occurred only in autumn and winter.The only exception was given by one positive HI event in spring (23 April 2021), which, however, took place immediately after the wet winter season.Furthermore, the most negative index values occurred all in late spring/summer, with generally dry conditions.HI values were positively correlated with the maximum discharges of the events (r = +0.32),while FI correlations were weaker than SN (r = −0.12).
Figure 10 shows examples of hydrographs, DOC chemographs and the corresponding C-Q relations for FIT and SN in the case of positive and negative HI values.Specifically, the 17 January 2021 event (Figures 10a-10d) concerned both catchments, with opposite HI signs in SN and FIT,respectively).This event is peculiar because the DOC concentration evolution had similar behavior and was synchronous in the two sections, occurring at the same time as the FIT discharge peak, while the discharge peak in SN was brought forward by 1 hr.In SN, despite the flow reduction, DOC concentration increased for the hour following the peak flow, contrasting the decrease in the total load.Rainfall peak intensity in the Fitterizzi rain gauge (only 4.4 mm hr −1 ) occurred in the same hour as the rainfall peak in SN.However, this rain gauge is not within the SN catchment; hence some rainfall features in this catchment, like the exact amount and timing, could have been missed.
Figures 10e and 10f show a summer event (15 July 2019) for FIT, with a negative HI value (−0.26) and significantly high DOC concentration.This event was characterized by a higher rainfall amount than the previous case (32.6 mm hr −1 2 hr before the peak flow, 8.6 mm hr −1 1 hr before and 15.8 mm hr −1 at the peak flow time).The suddenly increased hydraulic connection in the river network given by this typical summer rain shower following a dry period contributed to higher DOC concentration values in the falling limb.
Finally, an event in SN with a positive HI value (+0.27) is shown (Figures 10g and 10h).This event occurred at the end of the winter (wet) period (20 March 2021), concatenating two consecutive smaller events with clockwise evolution (Figure 10h).Interestingly, discharge was lower and DOC concentration higher in the first event, consistently with the assumption that, for positive HI values, the export processes are source-limited.In this case, rainfall intensity was low (maximum 3.8 mm hr −1 at the Fitterizzi gauge station).
The higher correlation of HI with discharges in FIT and the more pronounced seasonality suggests that its variability can be partially explained by focusing on the antecedent weather conditions, influencing soil moisture and hydraulic connectivity.On the other hand, it can be tested if similar mechanisms were activated even at SN despite the lower correlation with discharge.Figure 11 shows the HI correlation with the precipitation accumulated over different time intervals, ranging from 2 hours to 10 days before the discharge peak value.For both SN and FIT, a relative correlation peak was found for approximately 0.5 days of rainfall accumulation.After that, SN correlation decreased, while HI variability in FIT was more clearly explained by precipitation accumulated in the previous 6-7 days, up to r = 0.38 (p < 0.001).
The scatter plots of HI versus accumulated precipitation for the most highly correlated time intervals (i.e., 6 days for FIT and 0.5 days for SN; Figure 12) highlight the non-linear nature of the correlation.Scattered points were interpolated through generalized additive models (GAMs, Hastie & Tibshirani, 1986, 1990), which are smooth, nonparametric functions.Especially in FIT (Figure 12a), GAMs highlighted the non-linear relation between HI and antecedent precipitation.In general, for low accumulated precipitation values, DOC export was  transport-limited (HI < 0).Then, for higher accumulation values, meaning continuous (not necessarily intense) precipitation in the considered interval, DOC sources tended to be flushed, and the process became source-limited (HI > 0).Nevertheless, when accumulated precipitation was high enough, it could mobilize other DOC sources, and the process tended to return transport-limited.The event with the highest antecedent 6-day rainfall in FIT started on 25 January 2021 after several other events had just happened.Indeed, the C-Q graph of this event (not shown) includes more than one loop and is characterized as complex, according to Rose et al. (2018).The relation between antecedent precipitation and hysteresis described for FIT is much more roughly sketched out in the smaller SN catchment (Figure 12b), for which a smaller accumulation interval was also considered.

Discussion
The seasonal pattern observed for DOC concentration in both sites (Figures 4-6) agrees with other studies that reported freshwater DOC concentration peaking between late summer (e.g., Roig-Planasdemunt et al., 2017) and autumn (e.g., Aubert et al., 2013).In most cases, this seasonality is dependent on discharge (negative correlation), leading to clockwise hysteretic loops with respect to seasons due to lower DOC concentration in winter/ spring and higher concentration in summer/autumn periods (Aubert et al., 2013;Dawson et al., 2008;Fovet et al., 2018;Mulholland & Hill, 1997).In particular, intra-annual hysteretic loops in the Q-C plane (Figure 6) reflect seasonality in the catchment streamflow dynamics and in-stream variation in the net nutrient release and uptake rate.Therefore, they may show a site-dependent behavior (Aubert et al., 2013;Brooks et al., 1999;Dawson et al., 2001Dawson et al., , 2008;;Mulholland & Hill, 1997), as confirmed by our data set, which shows an anticipated DOC peak in the upstream SN site.Differences between DOC cycles in the two sites may be explained through the unlike contributions of the dominant hydrologic pathway and in-stream processes.On the other hand, temperature is a critical hydro-climatic parameter affecting the seasonal pattern of stream water chemistry.Specific climate conditions can even lead temperature and discharge relations to be equivalent with respect to DOC concentrations (Aubert et al., 2013).Lower autumn temperatures in SN than in FIT due to the higher mean altitude of the contributing catchment inhibit DOC production and export.
The baseflow carries a relatively low amount of DOC, primarily mobilized by individual storms.The seasonal analysis suggested that the difference between mean and median values is larger in winter for both sites (Figure 4).This result is connected to the more frequent high-discharge events observed in winter, whose effects can be seen even on the non-rainy days considered in the seasonal analysis.In agreement with previous studies, DOC peaks were observed during flood events (Blaurock et al., 2021;Rose et al., 2018;Vaughan et al., 2017), when significant DOC enhancement could be measured.As an example, Fovet et al. (2018) recorded an increase ΔC = 5 ± 4 mg C l −1 compared to a concentration of DOC < 1 mg C l −1 during baseflow, Blaurock et al. (2021) observed peaks of 10.2-18.6 and 8.5-16.9mg l −1 for two sites in comparison to a concentration of 2-3 mg l −1 during the baseflow.This case study confirms the influence of the topography on the mechanism of DOC mobilization and export during storms.Similarly to the catchment located in southeastern Germany and analyzed by Blaurock et al. (2021), the DOC was monitored in two different sub-catchments situated in an upper position with steep slopes and a lower and flatter site.While the DOC concentration was higher during the background periods in the upstream sub-catchment (i.e., SN), greater concentrations were recorded in the downstream site (i.e., FIT) during the storm events.Presumably, such events are the only ones capable of connecting to the river network rich DOC sources present in the flatter downstream area, more subject to agriculture and other human activities.The DOC average values in response to storm events confirm the results of Blaurock et al. (2021).They found average values for the 4 events reported equal to 3.88 and 1.75 mg l −1 in the lower and upper sub-catchment, respectively.They correlated this behavior to topography, highlighting that saturated soils are needed in flatter areas to allow efficient lateral water transport through DOC-rich soil layers toward the active river network.Therefore, it is plausible that DOC sources' connection to the active drainage network in FIT depends more on rain events.On the other hand, the SN discharge regime, mainly controlled by groundwater sources with lower DOC concentration (mountain springs), can reduce the effect of DOC sources' contribution activated by rain events.Moreover, like in Blaurock et al. (2021), in our study site, DOC mobilization was generally delayed in the flat lower catchment, as confirmed by the slighter decrease after the peak and hysteretic loops wider than the upper catchment (larger absolute values of HI in 13 cases out of the 19 simultaneous events).
In literature test areas, there is a prevalence of counterclockwise loops (negative values of HI) for DOC hysteresis (e.g., 51 counterclockwise compared to only 3 clearly clockwise and 46 complex events in Rose et al. (2018); 6 cases out of 8 counterclockwise in Blaurock et al. (2021)).Instead, the FI FI is mainly positive (Vaughan et al., 2017).In this study, the negative HI values were 21 over 29 and 13 over 29 for SN and FIT, respectively, while negative FI values were only 5 out of 29 for SN and were undetected in FIT.The clockwise loop is likely linked to the total catchment wetness seasonal pattern.Indeed, during or immediately after the end of the wet season, when the catchment water storage is high, the hillslope-channel hydrological connections are favored compared to other periods, as shown by previous studies in the catchment (Micieli et al., 2022;Senatore et al., 2021).Therefore, the DOC peak anticipates the discharge peak describing a clockwise hysteresis.In our study area, this effect proved more frequent in the flatter downstream FIT catchment, characterized by higher hydraulic connectivity.An example is given by the 17 January 2021 event (Figures 10a-10d), in which the time needed to reach the peak flow in FIT corresponded to the time required to mobilize the primary DOC sources in both catchments, while the peak flow in SN was attained earlier.The seasonal dependence of hysteresis is also in line with the results of Fovet et al. (2018).They showed clockwise hysteresis for 62% of events at a high flow period typical of the wet season in a Brittany, Western France catchment.
Finally, the study of the correlation between catchment wetness conditions and hysteresis direction was already addressed.Blaurock et al. (2021) found a positive correlation with the rainfall accumulated 14 days before the event started.Our analysis based on a variable accumulation time window applied in the two nested catchments confirmed such influence and highlighted different process timings depending on the catchment's size and morphology (Figure 11).The lower correlation of SN with more extended accumulation periods can be explained by its smaller extent (hence, faster response to storm events) and the relatively higher contribution of groundwater to its discharge.On the other hand, network connectivity in FIT looks more sensitive to the precipitation accumulated over a longer time interval.Such a result could be reasonably generalized since in larger, flatter catchments, complete hydraulic activation generally requires extended periods (ca.6-7 days in our case, for a 7 km 2 wide catchment), while in smaller, steeper catchments, the hysteretic processes are affected by shorter precipitation durations (less than 1 day for a 2 km 2 wide catchment).Nevertheless, precipitation characteristics (i.e., intensity and amount) are also important.We demonstrated that by expanding the analysis to increasingly high precipitation amounts, DOC export moves from transport-limited (HI < 0) to source-limited (HI > 0) and then again to transport-limited due to the connection of new DOC sources far from the stream.Figure 12 clearly shows that such a non-linear correlation can be identified in the steep upstream and, most clearly, in the flat downstream catchment, where the effect is probably amplified by the richer DOC sources due to agricultural and other human activities, which, however, are not well connected to the river network.

Conclusions
This study presented the results of a long-term monitoring campaign to unveil space and time DOC dynamics in a Mediterranean headwater catchment, relating them to meteorological and hydrological drivers.The different DOC dynamics observed in two nested sites were linked to spatially heterogeneous catchment properties (extent, orographic features, land uses).Two multi-parameter sondes were used to achieve that aim, and high-resolution continuous timeseries of several biogeochemical parameters were obtained.
The analysis relied on an original correction method, requiring water temperature and turbidity measurements to convert the observed f DOM into DOC values.Then, analyses performed at seasonal and storm event timescales provided several insights into DOC mobilization and export processes: • At the seasonal scale, univariate and multivariate statistical analysis confirmed the climate (seasonal) control on DOC production, with background concentrations increasing in hot and dry summer months due to the combined effect of enhanced photocatalytic degradation and reduced discharge in the channels; • Comparison of DOC concentrations taken simultaneously over 91 non-rainy days led to slightly higher values in the forested upstream catchment, having steeper topography and, of course, smaller streamflow; • However, observations made clear the importance of the hydrological regulation of DOC export, significantly activated by high-flow events, with discharge above Q10 being associated with 69% of the total yield in the upstream and 79% in the downstream site; • Also, the increased hillslope-channel connectivity all over the downstream catchment triggered by hydrological processes overturned the results of the seasonal background analysis, with DOC concentration higher in the downstream site considering the 3,760 simultaneous observations at the hourly time scale; • DOC sources proved to be plentiful in the zones contributing to the catchment's fast response in both sites, being able to increase concentration during almost all the storm events.Instead, the limiting factor of DOC export processes varied by season and location.In the steep upstream catchment with accentuated spatial intermittency, generally, such processes were transport-limited, while in the downstream catchment, more source-limited processes were observed; • Therefore, the HI was more positively correlated to antecedent precipitation in the downstream catchment.
However, such correlation was not linear since new DOC sources were activated with exceptionally high accumulated rainfall values, and the process tended to be transport-limited again.
Using high-resolution measurements has proven crucial to explain DOC dynamics at multiple time and spatial scales with a quantitative approach.However, though supported by laboratory measurements, the on-site recording showed some inherent weaknesses, primarily when high discharge was associated with high turbidity values, requiring the statistical retrieval of DOC peak values.Such a drawback can be partially overcome with increased on-site discrete automatic sampling during storm events and subsequent laboratory analysis.One of the further developments of the research goes toward this direction.Furthermore, it will be necessary to focus more on the processes' scaling properties, taking advantage of both the measurements in this and other sites, to support modeling approaches and contribute to a better understanding of the global carbon cycle.
analyzed seasonaland storm-scale DOC responses in clear-cut and forested headwater streams.Blaurock et al. (2021) highlighted the dependency on topography and antecedent wetness conditions.Koenig et al. (2017), Werner et al. (2019), Shogren et al. (2021), and Fazekas et al. (2020) performed multi-year investigations of the C-Q behavior across multiple sites and timescales.

Figure 1 .
Figure 1.The Turbolo catchment closed at the Fitterizzi outlet in southern Italy (inset).Colors show the altitude from 183 m a.s.l.(Fitterizzi gauge, labeled as FIT) to 1,005 m a.s.l.Isolines are at intervals of 50 m.The settlements of Mongrassano and Cavallerizzo are pointed out with orange shapes.A dashed red line contours the sub-catchment closed at the San Nicola gauge (labeled as SN).

Figure 2 .
Figure 2. Plot of the   DOM  DOC ratio versus the concentration of suspended particles measured in terms of turbidity.The   DOM  DOC contains all the deviations of the measured f DOM corrected for the temperature ( f DOM T ) from its ideal value ( f DOM CORR ), caused by inner filtering and suspended particle effects.

Figure 3 .
Figure 3. Correlation analysis between the laboratory Dissolved Organic carbon values measured on discrete samples and the corresponding f DOM corr values (gray circles).The corresponding f DOM T values are also reported (red circles).The increases of the f DOM values of two samples are highlighted as examples.

Figure 4 .
Figure 4. Box and whisker plots of temporal variation of background (i.e., no-rain days) seasonal Dissolved Organic carbon daily values (S: Summer, A: Autumn, W: Winter, and Sp: Spring) at Fitterizzi and San Nicola stations.Boxplot shows the first quartile (lower line), the median value (bold line), and the third quartile (upper line).The lower and upper whiskers show the minimum and maximum values, respectively; the dots beyond the whiskers indicate outliers.Black triangles represent seasonal means.

Figure 5 .
Figure 5. Principal component analysis for background seasonal data sets collected in (a) Fitterizzi and (b) San Nicola.The color points show the seasons: orange for autumn, blue for winter, pink for spring, and ochre for summer.The arrows indicate the eigenvectors and the related variables.

Figure 6 .
Figure 6. Background monthly data sets: mean monthly concentrations of (a, c) Dissolved Organic carbon and (b, d) mean monthly atmospheric temperature plotted against mean monthly specific discharge for (a, b) Fitterizzi (2021) and (c, d) San Nicola (2020).

Figure 7 .
Figure 7.For both SN (in red) and FIT (in blue), flow duration curves (hourly time step, continuous lines) expressed as Q/Q max (streamflow normalized to the maximum recorded hourly streamflow in both sites) and corresponding normalized accumulated Dissolved Organic carbon (DOC) load curves DOC acc /DOC acc max (dashed lines).The latter curves were created considering the DOC load values associated with each hourly streamflow value, which were accumulated and normalized to the total accumulated value (more details in the text).

Figure 8 .
Figure 8.Comparison between simultaneous FIT and SN Dissolved Organic carbon observations.

Figure 9 .
Figure 9. Storm hysteresis index (HI) versus storm flushing index (FI) for (a) FIT and (b) SN.For each quadrant, the description of HI (clockwise/anticlockwise) and FI (flushing/diluting) with the corresponding number of events is provided.The color points show the seasons: orange for autumn, blue for winter, pink for spring, and ochre for summer.

Figure 10 .
Figure 10.Dissolved Organic carbon concentrations and hydrographs (a, c, e, and g) and corresponding DOC-Q hysteresis (b, d, f, and h) during selected events: (a, b) positive hysteresis index (HI) event in FIT; (c, d) negative HI event in SN; (e, f) negative HI event in FIT; (g, h) positive HI event in SN.The spinning arrows and the ± signs in the DOC-Q graphs highlight the HI values (clockwise/positive and anticlockwise/negative). Positive flushing index characterizes all events.

Figure 11 .
Figure11.Hysteresis index correlation with precipitation accumulated in different time intervals ranging from 2 hours to 10 days before the discharge peak value, respectively, for FIT (blue) and SN (red).Correlation values were calculated considering all the selected events for each station.