Hydrodynamics of karst aquifers in metamorphic carbonate rocks: results from spring monitoring in the Apuan Alps (Tuscany, Italy)

During the last 40 years, extensive research has characterized the hydrogeology of many karst aquifers, as they are important water resources. Despite that, a systematic investigation on metamorphic karst aquifers is still lacking. The present study investigates the functioning of marble karst aquifers by means of spring hydrological monitoring, coupled with storm-hydrograph, thermograph, and chemograph (HTC) analysis and lag time analysis. Renara and Equi springs (Apuan Alps, Italy) were selected for this investigation. These springs drain catchments that have different degrees of structural complexity. Piston flow is the common hydrodynamic response of Renara spring to infiltration. Strong dilution effects were observed only during the heaviest rainfall events. Prolonged dry conditions after a sustained recharge phase showed the delayed arrival of infiltration water about a month later. Equi spring has a more complex behaviour due to its wider and more hydrologically heterogeneous catchment but the comparison of HTC graphs during the winter dry phase helped to recognize the differential contributions of proximal and distal sectors. Both springs show a rapid discharge increase in response to impulse infiltrative events. Conversely, water temperature and specific electrical conductivity increase only slightly during floods, indicating limited chemical and thermal exchanges between the rock and the water stored in these aquifers. The hydrodynamic behaviour of these karst springs suggests that the Apuan metamorphic aquifers are characterized by the predominance of conduit porosity over fissure and matrix porosities. This is explained by a reduced interstitial porosity and fracturing of the metamorphic carbonate rocks.


Introduction
Karst aquifers are very important water resources for public, agricultural, and industrial uses in many countries worldwide but, due to their peculiar features, their management requires a deep knowledge of the hydrogeological dynamics Stevanović 2015;Goldscheider et al. 2020). Usually, karst aquifers host most of the water in fissures and small secondary voids, whereas flow is concentrated along solution-enlarged conduits that are characterized by a high hydraulic conductivity (e.g., Mangin 1994;White 2002). The metamorphic processes that transform limestone and dolostone into marbles lead to a decrease in the porosity and permeability of the rock mass. Even in lowgrade conditions, metamorphism causes the partial obliteration of sedimentary and diagenetic structures (i.e., bedding planes) that are important flow pathways and inception horizons for karstification (Klimchouk and Ford 2000;Lauritzen 2001;Filipponi et al. 2009). For these reasons, groundwater drainage systems in metamorphic carbonate aquifers appear to have peculiar structures with less infiltration and storage than the nonmetamorphosed carbonate ones (Andreo et al. 1997;Despain and Stock 2005;Antonellini et al. 2019).
Continuous monitoring is key to studying the hydrodynamics of karst springs because these features are affected by significant and rapid changes to their physical and chemical characteristics (e.g., Bonacci 1993;Mangin 1994;Goldscheider and Drew 2014). The most important monitored parameters are water discharge/level (Q/wl), temperature (T w ) and electrical conductivity (EC). Comparison of the time series of these parameters gives insights into the aquifer architecture that the hydrographs alone cannot provide (Birk et al. 2004;Mudarra et al. 2014;Vigna and Banzato 2015;Nannoni et al. 2020). Water temperature (T w ) can be considered as a nonconservative tracer due to heat exchange between water and rock (Genthon et al. 2005). The relationships between groundwater temperature, surface water temperature, geothermal flux, and climate are complex. Brookfield et al. (2016) investigated how precipitation intensity and temperature influence the groundwater temperature pattern of a merokarst in Kansas, USA. Their results demonstrated that the recharge-event temperature strongly affects the aquifer. Luhmann et al. (2011) monitored water level and temperature at 25 springs in Minnesota (USA) and identified four different thermal patterns at seasonal-and event-scale that reflect the heat exchange dynamics and permit evaluation of the storage compartments and residence times in the karst aquifers. Electrical conductivity (EC) is a nonconservative tracer that is helpful in understanding the mixing processes that occur in karst aquifers; in fact, its variation is often related to changes in localized recharge (Birk et al. 2004;Calligaris et al. 2019) and carbonate rock dissolution. Moreover, spring EC is a reliable integrated parameter to obtain a more realistic simulation of karst aquifer functioning (Birk et al. 2006;Chang et al. 2021). Hydrological monitoring is usually coupled with precipitation data to evaluate the aquifer hydrodynamic response to different recharge conditions. The karst spring hydrological and hydrochemical behaviour represents the output signal of the interaction between the infiltrating precipitation, which can be viewed as the input signal, and the water stored in the aquifer. Moreover, rainfall amount and intensity, and water supply preceding a storm event (i.e., drought or wet conditions), play a pivotal role in the spring response to infiltration (Barberá and Andreo 2012). Lags between discharge, water temperature, and EC responses to precipitation have been used to estimate the volume of karst conduit systems and to investigate karst aquifer architecture and recharge dynamics (Liñán Baena et al. 2008;Cholet et al. 2015;Adji and Bahtiar 2016;Sánchez et al. 2017).
Presently, there is a huge amount of literature about karst system functioning, but a systematic investigation on the hydrodynamics of aquifers developed in metamorphic carbonate rocks is still lacking. Apart from a few studies about the hydrodynamics of some marble aquifers, such as those developed in the Norwegian subarctic region or those found in northern Greece and Turkey (Özler 2001;Lauritzen 2013, 2017;Petalas 2017;Petalas and Moutsopoulos 2019), no general definitions/features, either qualitative or quantitative, are given for these unusual karst systems.
Due to their geological-structural setting and the high degree of hydrogeological knowledge, the Apuan Alps karst region (central Italy) can be considered a representative area to study the hydrology of worldwide metamorphosed karst aquifers (Doveri et al. 2017). In this report, the results of the monitoring of Renara and Equi springs, two of the major Apuan karst outflows, are presented. These springs were selected because (1) they drain aquifers that are representative of the functioning of a marble karst system, and (2) their hydrogeological catchments have been well-constrained by previous studies. The main aim of this study is to define the particular characteristics of karst aquifers developed in metamorphic carbonate rocks, where groundwater flow occurs mainly in well-channelled networks with a negligible contribution from matrix and fracture porosities.

Geomorphological and geological setting
The Apuan Alps of central Italy have an alpine-like landscape that differentiates them from any other mountain ridge of the northern Apennines. Metamorphism enhanced the different responses to weathering and denudation processes of the lithologies, producing a rough landscape where the highest and steepest reliefs consist mainly of carbonate formations, whereas valleys and subdued areas are often made up of siliciclastic formations. Dissolution is an important geomorphic process on the carbonate outcrops. Several caves have been surveyed in this area, including some of the deepest and longest in Europe (Piccini 1998).
The local complex geological setting is characterized by the occurrence of several overlapped tectonic units, belonging to the Tuscan and Ligurian sedimentary domains. The lowermost units (Apuane and Massa units) were affected by green-schist facies metamorphism dated between 27 and 10 Ma (Kligfield et al. 1986) and were overthrusted by nonmetamorphic units and then exhumed during a late Miocene to Pliocene extensive tectonic phase (Carmignani and Kligfield 1990). The Apuane Unit consists of a metasedimentary sequence deposited since Upper Permian age on a metamorphic basement made up of Ordovician phyllites and Early Permian porphyric metavolcanics (Conti et al. 1993). The sedimentary sequence begins with conglomerates ("Verrucano Formation") on which a carbonate shelf developed from Carnic to Sinemurian. The carbonate platform sequence is formed from bottom to top by dolostone, dolomitic limestone and limestone. Since Pliensbachian, the deposition became pelagic with siliceous limestones followed by radiolarites (Tithonian). Afterwards, the sedimentation progressively shifted to shales and calcarenites followed by silico-clastic turbidites of upper Oligocene age.
The Massa Unit has the same basement of the Apuane Unit, but the sedimentary sequence displays mainly Triassic pericontinental facies. The Tuscan and Ligurian Nappes overlay both the metamorphic units; the former has a stratigraphic sequence like that of the Apuane Unit, with a higher thickness of hemipelagic and terrigenous formations (turbidites).
Large-scale NE verging isoclinal folds characterize the structural setting of the metamorphic complex. A late-Miocene-to-Pliocene extensional tectonic phase caused the uplift of the massif along NW-SE and NE-SW striking faults (Carmignani and Kligfield 1990). The exposed metamorphic complex is not significantly affected by large postorogenetic faulting, but a uniformly oriented system of minor faults and joints is recognized (Ottria and Molli 2000).

Hydrogeology
The hydrogeological setting of the Apuan Alps is strongly affected by the occurrence of wide outcrops of carbonate rocks (about 140 km 2 ; Fig. 1). The main Apuan metamorphosed karst systems consist of marble and metadolostone that display high infiltration rates and the development of Fig. 1 a Index map, the red square highlights the location of part b. b Simplified hydrogeological map of the Apuan Alps (see Piccini et al. 1999 for further details about the regional hydrogeology of the Apuan Alps). The red squares show the two investigated systems: 1 Equi, 2 Renara extensive underground conduit networks (Piccini 2002;Menichini et al. 2016;Doveri et al. 2017). Marbles and metadolostone aquifers have a low storage capacity due to low matrix porosity and to reduced fracturing at depth. For this reason, in some conditions, the epikarst can have a significant role in recharging spring baseflow during prolonged droughts (Nannoni and Piccini 2022).
Siliceous limestones (Metacalcari Selciferi and Metacalcari Selciferi a Entrochi formations) are locally well karstified, although they usually display a lower infiltration index (Piccini 2002). All the upper hemipelagic and terrigenous formations (metapelite and metaarenites) have a very low permeability. Palaeozoic basement and Cretaceous-Oligocene metapelite/metaarenites usually act as lateral limits of the carbonate aquifers. Due to the complexity of deep karst drainage, the hydrogeological and hydrographic catchments differ greatly, and subterranean flow paths do not often follow the hydraulic gradient. Therefore, it is difficult to determine the limits of the former ones if no tracer tests are available.
Despite the steep reliefs, infiltration rate usually exceeds more than 50% of rainfall. Below an elevation that ranges from 600 to 200 m above sea level (asl), groundwater flows in phreatic conditions. At the regional scale, the main karst aquifer is made up of metadolostone, dolomitic marble, marble and, locally, by cherty metalimestones; this series is limited by impermeable rocks of the Palaeozoic basement and by low-permeability rocks, mainly schists, phyllites and jasper.
Most of the karst springs are located on the seaward side, the major spring among them being Forno Spring, whose average flow rate is about 1.6 m 3 /s. The second highest average flow rate (about 0.9 m 3 /s) was measured at Pollaccia Spring, on the inner side of the ridge. Another spring with a total significant mean discharge (about 0.8 m 3 /s) is located in the north-western part of the Apuan Alps, close to the village of Equi. At least a dozen other springs have flow rates ranging between 0.4 and 0.1 m 3 /s and many others have discharge higher than 10 L/s (Piccini 2002;Doveri et al. 2017).

Investigated karst systems
The two studied springs drain water from catchments characterized by different extensions and complexity (Fig. 2). Renara spring was selected due to the relatively simple geometry of its hydrogeological structure, and the small extension of the feeding system. Equi spring has a wider catchment and a higher structural complexity. It is known that both systems have well-defined hydrogeological catchments thanks to several dye-tests (De Grande et al. 1997;Malcapi 2011;Roncioni 2011;Poggetti et al. 2015).
Renara is located at an elevation of about 295 m asl, about 6 km NE from the city of Massa. The water flows out along the Renara creek channel (Fig. 2a), just upstream of the contact between the metadolostones and the Paleozoic basement. This spring is the main outflow of a well-defined karst system that extends from M. Altissimo to the north (Fig. 2b). Before the outlet, the water flows through Buca di Renara Cave, which is located a few hundred meters upstream; it has been explored by divers for a length of about 700 m and a depth of 100 m below the water level (Fig. 2a). The spring has an average discharge of about 0.2 m 3 /s (Piccini et al. 1999), and the catchment area has an extension of about 5.5 km 2 . The hydrogeological structure consists of an isoclinal syncline, and it is delimited by the impermeable basement rocks on the eastern and western sides (Fig. 2b). The M. Pelato-M. Altissimo hydro-structure comprises two major hydrogeological systems, from NNW to SSE (Fig. 2b): the Renara catchment, and the Polla dell'Altissimo. The groundwater divides between these catchments were defined by means of combined structural, geochemical, and hydrological investigations (Menichini 2012;Doveri et al. 2013). A tracer test showed a long transit time across the Renara catchment: the first detection occurred 26 days after dye injection, with 304 and 2,200 m of vertical and horizontal distances, respectively, between injection and detection points (Menichini 2012). Another tracer test in the Abisso dei Fulmini cave helped to define the groundwater divide between Renara and Polla dell'Altissimo catchments and to recognize the fast circulation in the latter system-less than 14 h for the first dye detection (Menichini 2012).
Equi spring is located along a minor tributary of Lucido Creek. The main upper outflow is named Barrila and it emerges from a small cave hosted in dolomitic marble, near the contact with the cherty metalimestone, whereas another source discharges from Buca di Equi Cave (Fig. 2c). The Barrila outflow represents a seasonal overflow, whereas the Buca di Equi spring has a more regular regime. Some thermal springs are located along Lucido Creek, a few hundred meters upstream from the village of Equi Terme. Thermal waters (mean temperature 23 °C) are related to a deep circulation system that is controlled by a major tectonic structure, the Equi Terme Fault (Molli et al. 2015). Recent studies estimated a total average flow of about 0.8 m 3 /s for the cold-water sources, whereas the total discharge can exceed 35 m 3 /s during floods (Poggetti et al. 2015). The catchment area extends over about 20 km 2 and is now welldefined thanks to several dye tests performed in deep caves of the northern Apuan Alps (De Grande et al. 1997;Malcapi 2011;Roncioni 2011). This catchment is characterized by a complex structural setting consisting of three majororder isoclinal folds, from west to east, the Vinca anticline, the Orto di Donna syncline, and the M. Pisanino anticline (Fig. 2d). The contact between the Paleozoic basement and the metamorphic carbonate sequence limits the catchment on the west side, whereas the hinge line of the M. Pisanino anticline represents the easternmost limit. The occurrence of low-permeability rocks along the Orto di Donna syncline complicates the hydrogeologic structure, and the hydrologic connection between the Pizzo d'Uccello-M. Grondilice and the M. Pisanino-M. Cavallo areas was recognized only thanks to the dye tests. Short transit times were estimated with dye tests: 7 days passed from dye injection at Buca Nuova (SSE slope of M. Pizzo d'Uccello) and first detection at the spring, with vertical and horizontal distances of 925 and 4,500 m, respectively (Roncioni 2011). Flow velocities are high even in the central (Orto di Donna) sector of the catchment with 887 and 6,000 m of vertical and horizontal distances, respectively, travelled in 4 days from the injection point in the Abisso Panné cave to the spring. Relatively longer transit times were found for the distal sectors (M. Pisanino-M. Cavallo ridge): 7,000 m of horizontal distance and an elevation difference of 1,190 m were covered in 23 days from the injection point in the Abisso Perestroika cave to the spring but this test was performed in low-flow conditions. Flow velocities are indeed higher during high flow, as demonstrated by another test in the southeasternmost sector (Malcapi 2011).

Methods
Renara and Equi springs were both equipped with a multiparametric probe (model Eijkelkamp CTD-Diver) measuring water pressure, temperature (T w ), and electrical conductivity (EC), and with a barometric probe (model Eijkelkamp Baro-Diver) for atmospheric pressure and air temperature. The pressure resolution was 2 mm with an accuracy of 5 mm. The two pressure sensors measured the water load on the immersed probe and the atmospheric pressure, respectively. The difference between the pressure measured by the diver and the atmospheric pressure provides the height of the water level above the immersed probe (wl). The characteristics of Renara spring do not allow flow rate measurements; therefore, only the water level height (wl) was considered. At Equi spring, three discharge (Q) measures were performed by means of a propeller flow meter, and a discharge curve was obtained. These measures were taken in low to medium flow conditions, because during floods (Q > 5 m 3 /s) the water flux is too irregular and turbulent to be measured with certainty. However, the discharge values in high flow conditions are reliable data for the definition of the storm hydrograph shape and magnitude. The water temperature sensor has a resolution of 0.01 °C and an accuracy of ±0.1 °C. The electrical conductivity is measured with a resolution of 0.1%, and accuracy of 1% of reading (about 2-3 μS/ cm). Conductivity data are converted to specific conductivity (EC sp ) at the reference temperature of 25 °C. A sampling frequency of 15 min was set for all the probes. Precipitation data (P) were collected by two meteorological stations of

Storm hydro-thermo-chemo-graphs and lag measurements
The components of the storm hydrographs of the two springs were analysed, with a focus on the rising limb and its relationship with the precipitation pattern (e.g., Kresic and Stevanovic 2009). The hydrograph parameters considered in this study are described in the electronic supplementary material (ESM1). Two types of storm hydrographs can be broadly distinguished: simple and complex (multiple)-the former consists of single peaks, whereas the latter comprises consecutive peaks because of multiple precipitation events. Rainfall rate and the duration of dry conditions (i.e., without significant precipitation) before each infiltration event were considered. The comparison between wl/Q, T w , and EC sp was carried out by means of hydro-thermo-chemo-graphs (henceforth abbreviated to HTC-graphs), and the lag times were measured between P and wl/Q peaks and between wl/Q and T w /EC sp peaks (Fig. S1a of ESM1). Figure 3 shows the HTC-graphs-water level wl for Renara (Fig. 3a), discharge Q for Equi (Fig. 3b)-of the two springs compared to precipitation data (P). The data are reported with a 30-min time step and concern only the period from autumn to late spring, because from the beginning of summer to early autumn both springs were affected by dry conditions that led to frequent emersions of the probe. The statistics of the hydrological and meteorological monitoring are reported in Table 1.

Hydrological monitoring
At Renara spring, the highest water levels were reached during January and February 2016, with fluctuations up to 5 m of the internal sump (Fig. 3a). The highest flow oscillations were due to abundant, long-lasting precipitation, with minor snowmelt contributions. In fact, the Renara catchment is mostly located at a relatively low elevation (mean 1,010 m asl) and on the seaward side of the Apuan ridge, where the local climate does not allow significant snow accumulation. From October to the end of November 2015, the water-level rises were characterized by sharp increases and partially undisturbed recession phases. These peaks were related to intense but shortlived storm events. The whole of December was characterized by no precipitation, causing low flow conditions. Water level fluctuations were accompanied by T w and EC sp oscillations. Both parameters showed slight seasonal trends with the highest mean values occurring in December 2015, whereas lower values lasted from January to April 2016. EC sp significantly decreased from about 230-240 μS/cm to roughly 200 μS/cm after the strongest storms. Such marked shift was not observed for the T w seasonal trend. During storm events, EC sp was constant or showed a significant variation (up to 30-40 μS/cm), either decreases or increases, whereas T w always showed rapid and clear decreases when flow peaks occurred. Such behaviour is unusual, since temperature normally decreases due to the infiltration of meteoric waters, which normally results also in a marked lowering of EC sp .
The Equi spring hydrograph showed rapid, pulse-like responses to storm events (Fig. 3b). Heavy flood events occurred throughout the whole monitored period and the two highest flow peaks were observed in December 2012 (Q = 34 m 3 /s) and March 2013 (Q = 31 m 3 /s). The shape of the hydrographs depended on both the intensity of precipitation and the snowmelt contribution that can be relevant in the highest and northern sectors of the recharge area. Moreover, this catchment is in the internal side of the Apuan Alps where the climate conditions are colder and wetter than those on the seaward side. In fact, the yearly precipitation of the Orto di Donna station is much higher than that registered at the Pian della Fioba station: the former station showed a mean yearly precipitation of 2,816 mm for the 2010-2020 interval, whereas 1,923 mm is the mean value calculated for the latter in the same time interval. Accordingly, the two meteorological stations collected different precipitation amounts during the monitoring interval: 2,657 and 1,110 mm for Orto di Donna and Pian della Fioba, respectively. Two seasonal drought phases occurred: the first one happened between mid-December 2012 and mid-January 2013, when only a couple of weak precipitation events occurred. The  driest interval lasted the whole of February 2013 with very weak precipitation that did not cause any flow rate variations. Both T w and EC sp showed decreasing trends from autumn to spring, with a more pronounced decrease for EC sp . These trends are related to snow melting during spring that is absent in autumn. The small-scale, diurnal oscillation of T w are an artifact related to the influence of diurnal air temperature oscillations. However, these short-period oscillations are smaller in intensity than the T w variability observed during storm events; therefore, the thermograph can be used to evaluate the spring response to recharge events. During peak discharge, T w showed an impulsive behaviour characterized by consecutive positive-negative peaks with total oscillations over a 1.0 °C of amplitude. EC sp response to flow-rate increases was also impulsive but it was not characterized by the double peak oscillations observed for T w . The fluctuations were either positive or negative, but the positive peaks were usually smaller in amplitude than the latter ones. The strongest EC sp oscillation was observed during March 2013 with a decrease of about 100 μS/cm following a heavy multiple storm event.

HTC-graph analysis
The storm HTC-graph analysis was performed for both springs on most of the flood events that occurred during the monitoring periods. The measured parameters are listed in ESM2. Two storm HTC-graphs (either single or complex) are described here in detail for each spring, being representative of the hydrodynamic behaviour in different hydrological conditions.

Renara spring
The storm events monitored at Renara spring showed a lower variability (range = 1.5-8 h) and a lower mean value of lag times Δt lag with respect to Equi (mean = 3.52 h, see ESM2). Similar observations can be done for the starting time Δt s ; Renara spring always had quick responses to precipitation (range = 0.5-5.5 h, mean = 2.15 h). The time lags between Q C (i.e., the flow rate/water level peak) and the first T w peak ranged between -4 h (i.e., T w peak precedes the Q peak) and 5.5 h, with an average value of 0.5 h, which is higher than that observed for Equi spring. EC sp time lags varied between -3 and 13 h, with a mean value of 2.75 h. Many small to medium recharge events showed very low to no variation of EC sp . These recharge events happened during winter-e.g., 14/01/2016, 31/01/2016, 03/02/2016 (ESM2). The first graph (Fig. 4a) resulted from two consecutive storm events. The first flow increase caused the highest water-level (wl) fluctuation of the autumn recharge phase (Q C -Q A = 377.2 cm) that was observed 4 h after the precipitation pulse. The second peak had a smaller amplitude, but it happened only 2 h after the last P impulse. T w sharply decreased by 0.5 °C during the wl increase, with secondary peaks during the two wl maxima. Overall, EC sp showed a decreasing trend throughout the multiple recharge event (-12 μS/cm variation) but two small positive peaks occurred immediately before or after the two wl maxima (Δt lag(ECp1,2-Qc) = -3 h, 2.5 h).
The second storm graph (Fig. 4b) was characterised by a regular, strong water-level increase (Q C -Q A = 383.5 cm) following a short-lived but intense winter storm event (P rate = 5.63 mm/h). The response to precipitation was fast with a time lag of 2.5 h. T w decreased by 0.17 °C and the lowest value was measured half an hour before the maximum water level. T w returned to preflood values after 24 h. EC sp first slightly increased (+5 μS/cm) during the hydrograph rising limb, then it gradually decreased of about 13 μS/cm, reaching the lowest value 15 h after the peak flow. EC sp returned to preflood values 2 days after the start of the discharge increase.
The two examples suggest the occurrence of a piston flow effect of water stored inside the aquifer followed by mixing phenomena that are characterized by small EC sp decreases and marked T w decreases. The on-phase and off-phase T w and EC sp trends will be examined in the 'Discussion' section.

Equi spring
This spring exhibited a broader Δt lag range (1.5-13 h, ESM2) than the one observed at Renara spring. Moreover, the average Δt lag has a value of 5.24 h. The spring showed also the highest mean Δt s (4.67 h) among the two springs and a range of variation comprised between 0.5 and 9.5 h. The lag times between the peak discharge and the first temperature peak (Δt lag(Tp1-Qc) ) showed a mean value of -0.35 h. The measured Δt lag(Tp1-Qc) varied between -10.5 and 13 h; therefore, this spring had the highest variability of this parameter. EC sp lag time (Δt lag(ECp1-Qc) ) ranged from -9.5 h to 30 h and its mean value was equal to 3.83 h. The latter value is higher than that of Renara, therefore Equi showed, on average, the slowest EC sp arrival time during flood events.
The first presented HTC-graph (Fig. 5a) is a late autumnal event that was made up by multiple discharge peaks as a response to three consecutive precipitation events. The first Q peak (Q C ) occurred 2.5 h after the rainfall pulse, and it was characterized by a strong Q increase-Q C -Q A = 19.9 m 3 /s (ESM2). The precipitation pattern strongly influenced the shape of the hydrograph peaks. In fact, the P rate was equal to 13.7 mm/h during the first storm, whereas it reached 3.5 and 3.2 mm/h during the second and third P impulses, respectively. Consequently, the corresponding Q increases were small but long-lasting. T w showed a sequence of positive/negative impulses with a total amplitude of 0.8 °C. T w almost returned to preflood values 2.5 days after the last discharge peak. EC sp had a trend like that of T w , with a positive/ negative sequence. The total amplitude of the oscillations was equal to 42 μS/cm. EC sp started to continuously increase after the flood event exceeding prestorm values.
The second HTC-graph (Fig. 5b) has a Q C value close to that of the previous one, but it was caused by a storm event with a lower P rate . However, the latter occurred during winter conditions (February 2013, Fig. 5b) when snow accumulation is significant in the highest parts of the catchment. The flow rate increase was slow during the rising limb (concentration time Δt c = 20 h, ESM2) and the Δt lag was higher than those measured for the previous example (7.5 h). T w and EC sp showed again in-phase patterns with positive-negative consecutive oscillations. Despite some secondary oscillations, T w and EC sp first increased by 0.4 °C and 15 μS/cm, respectively. Then, a clear decrease occurred for both parameters during the peak flow and the beginning of the falling limb, although with a small offset between their negative peaks (see also Δt lag(Tp2-Qc) = 2 h, and Δt lag(ECp2-Qc) = 6 h). Then, the conductivity started to rise slowly during the recession phase to return to previous values. The presented storm HTC-graphs show that the common behaviour of Equi spring is a small-scale short-lived piston flow effect that is followed by a significant dilution effect. The balance between these two phenomena depends on several factors such as the climate conditions (i.e., snow availability in the catchment), the storm event characteristics, and the recharge provenance (proximal vs distal sectors).

Discussion
The monitoring of these two karst systems shows complex hydrodynamic behaviour that is affected by structural and meteorological factors. The Renara spring hydrodynamics can be considered representative of the functioning of the Apuan metamorphosed karst aquifers characterized by catchments with simple hydrogeological structures.
Most of the Renara storm chemographs show a slight increase of EC sp . The lack of marked dilution effects during the storm events of the autumnal recharge phase (a slight one was observed only during the strongest event, see Fig. 4a) suggests that the infiltrated water pushes out water already present in the karst conduits; therefore, the piston flow phenomenon is the common hydrodynamic response to infiltration. Marked dilution effects were observed only during the strongest flood events. The simultaneous decrease of T w and the EC sp increase during flood peaks could be caused by the removal of colder water trapped in the depressed sections (e.g., at the bottom of the downward segments) of large phreatic conduits in slow flow conditions due to its higher density. During flood events, the entire water stored in the aquifer is pushed out and mixed, with the result of a cooling accompanied by a modest positive variation of EC sp . A similar observation was made by Drysdale et al. (2001) at Cartaro Spring (in the Frigido River basin, north-western Apuan Alps), even if in this case a dilution effect was observed a few hours after the flow peak, thus indicating a shorter transit time for the infiltration water to reach the spring. Fig. 4 Examples of storm HTC-graphs observed at Renara spring compared to precipitation at the Pian della Fioba station. Each one highlights some features of the spring hydrodynamic response to precipitation: a response to autumn recharges with characteristic T w and EC sp patterns; b storm HTC-graph observed at the end of winter The chemograph of Renara spring shows EC sp variations not directly related to storm events that can be explained by admitting that the infiltrating water propagates slowly through the phreatic conduit system and without significant mixing with the water already stored in the aquifer (Fig. 6a). The portion of the chemograph from 18 to 25 November 2015 shows a net increase in conductivity of about 25 μS/ cm following heavy rainfall ("1" in Fig. 6b). After the flood peak, the EC sp remained on values higher than the pre-event ones and decreased very slowly. A month later, a net decrease of EC sp was observed, which lasted for about 4 days, followed by an increase and then by other minor decreases until the occurrence of a new infiltrative event ("2" in Fig. 6b). These decreases, which were not related to recent rainfall and dilution phenomena, probably correspond to the arrival at the spring of the water infiltrated during the rainfall event that occurred about 1 month earlier (Fig. 6c). Such long transit times agree with the results of tracer tests performed in this system during low-flow conditions (see section 'Investigated karst systems). This behaviour implies the existence of a distinctive hierarchical feeding system consisting of regularly shaped conduits where laminar flow occurs, and the water-rock chemical exchanges are limited. Finally, it should be noted that these EC sp variations are not accompanied by significant temperature variations.
Equi spring showed complex hydrodynamic responses to infiltration characterized by in-phase EC sp /T w oscillations, usually positive-negative. Commonly, EC sp shows short-lasting increases during the rising limb of the storm hydrograph, followed by marked decreases (from -20 to -75 μS/cm variation) during the peak flow and the recession phase. Such behaviour indicates a limited "piston effect" that occurs in the first 15-20 h, due to the removal of water contained in the conduits close to the outflow, the volume of which can be roughly estimated at about one million cubic meters (assuming a mean flow rate of 12 m 3 /s for 24 h), followed by the arrival of water diluted by infiltration all over the catchment area.
The return to preflood conditions typically requires 4-5 days. As mentioned before, this system is made up of multiple subcatchments divided by several first and second order, tight isoclinal folds where the noncarbonate rocks act as lateral barriers to water flow (Fig. 7a). The differential contributions from proximal and distal areas can be observed during the winter dry phase of February 2013 (Fig. 7b). After the storm event on the 2nd of February, EC sp returned to prepeak values ("1" in Fig. 7b), only to decrease again without the input of infiltrating water or flow rate increases ("2" in Fig. 7b). EC sp did not return to prepeak values after that.
This mineralization decrease was probably related to the arrival of water that infiltrated in the more distal part of the system (i.e., the M. Cavallo-M. Pisanino sector, Fig. 7a,c). This water took less than 10 days to move through the aquifer in low-flow conditions, covering roughly 6-7 km and an elevation difference of ~1,000-1,190 m). Such short transit times could be related to the fast free-surface vadose flow in a few master conduits that quickly transfer water to the proximal saturated zone, as already suggested by the results Fig. 5 Examples of storm HTCgraphs observed at Equi spring: a spring multiple response to an autumnal recharge event; b storm HTC-graph observed at the end of winter of the tracer tests (De Grande et al. 1997;Malcapi 2011;Roncioni 2011;Poggetti et al. 2015).
Overall, the two studied aquifers show some common features such as: (1) strong structural control on the aquifer architecture, (2) vauclusian springs, and (3) fast, impulsive response to precipitation, as demonstrated by the short P-Q and EC sp /T w -Q lag times of both single and multiple storm events. The P-Q lag analysis allows the comparison of the values obtained for the Apuan aquifers with those observed in some nonmetamorphic karst aquifers. Few works dealt in detail with the single-event P-Q lags (and even less regarding the T/ EC-Q lags); however, they show that the nonmetamorphic karst aquifers have longer lags than those observed for the Apuan metamorphic karst aquifers (Fig. 8). Only the Petoyan system shows lags like those of Renara and Equi, but it is characterized by a very small catchment (3 km 2 ; Adji and Bahtiar 2016).
The very short P-Q and EC sp /T w -Q lags observed in the Equi system, regardless of its catchment dimension, can be explained if water moves in a few enlarged master conduits predominantly as free-surface flow, with small fracture storability (Liñán Baena et al. 2008;Sanchéz et al. 2017;Skoglund and Lauritzen 2017). Renara has a smaller catchment and water is transferred from a mostly vertical vadose zone to a few, preferential phreatic conduits in which the water moves with laminar flow.
The opposing oscillation of the EC sp chemograph and T w termograph during a storm event is another peculiar characteristic of the Apuan metamorphic aquifers. There are few other examples of springs that showed an increase in mineralization and a simultaneous temperature decrease during storm events. Calligaris et al. (2019) reported such patterns for Timavo Spring and the Abisso Samar underground river during a couple of storm events in February 2016. Sanchéz et al. (2017) found the same contrasting patterns at a karst spring in southern Spain during a storm event, though it was only a single case of such hydrodynamic behaviour. The recurrency of this phenomenon in the Apuan karst aquifers could be a general characteristic of metamorphosed karst aquifers with a well-hierarchized conduit system.

Conclusions
In the Apuan Alps, there are 12 karst springs with a flow rate exceeding 100 L/s. The major of them have been monitored in recent years by the Tuscan Speleological Federation and the Department of Earth Sciences of Florence, Italy. Among these, two springs were selected in this work as representative of karst systems in Apuan metamorphic rocks.
During floods, piston flow is the common hydrodynamic response of Renara spring to infiltration. Strong dilution effects were observed only during the heaviest rainfall events. Prolonged dry conditions after a sustained recharge phase permitted recognition of the delayed arrival of infiltration water (~1 month later) evidenced by a strong decrease in the electrical conductivity. This suggests that water moves with laminar flow in a distinctive hierarchical system of phreatic conduits in which the water-rock chemical exchanges are limited. The temperature decreases observed during some flood events can be explained by a thermal stratification of the waters in the deepest sectors of the phreatic zone.
Equi spring has a more complex behaviour, due to its wider catchment, which includes sectors with different hydrodynamic conditions. The T w and EC sp oscillations Fig. 7 a Schematic profile across the Equi aquifer (see Fig. 2d for its location); α, β, and γ refer to the three different cave systems in which the mentioned tracer tests were performed (the numbers in brackets are the caves' identification codes in the Tuscan Speleologi-cal Federation archive); b Q, T w and EC sp graph from 30 January 2013 to 4 March 2013 showing a period with no significant precipitation events; c Groundwater flow models during steps 1 and 2 of b, the labels in brackets refer to the catchment sectors in a and their relationship with the infiltration dynamics through the various catchment sectors are complex to interpret. However, the comparison of HTC-graphs during the winter dry phase helped to recognize the differential contributions of proximal and distal sectors: the infiltrating water does not reach chemical equilibrium with the rock because it does not have time to reside either in the former sectors (small distances) or in the latter ones, due to vadose fast flow in the master conduits.
Both springs show a rapid response of the discharge to impulse infiltrative events. This behaviour is typical of springs characterized by a drainage network in a saturated zone capable of instantaneously transmitting the hydraulicload increases. Therefore, the delay between flow peaks and precipitation impulses is essentially due to the transit time of the infiltrating waters in the vadose zone, which usually is limited to a few hours. On the other hand, the two sources show only slight increases in T w and EC sp during floods, indicating that chemical and thermal exchanges are limited between the rock and the water stored in these aquifers. The hydrodynamic behaviour of these karst springs suggests the occurrence of a saturated zone with water reserves stored almost entirely in the network of conduits. The low impedance of such a network allows for rapid dissipation of the abrupt increases of hydraulic head that occur during storm events. The percolation through the unsaturated zone does not occur in a diffuse mode, but chiefly along some preferential pathways whose distance from the springs affects the lag times. All these features seem to indicate that the Apuan metamorphic aquifers are characterized by the predominance of conduit porosity over fissure and matrix porosities. This is explained by assuming that metamorphism is responsible for a reduced interstitial porosity and fracturing of carbonate rocks. The drainage through well-structured conduit networks determines a high hydraulic "sensitivity" and limits the phenomena of mixing, homogenization, and attenuation of meteorological signals. It was possible to recognize these characteristics by means of HTC-graph analysis during the prolonged periods of low flow not influenced by significant infiltrative events. Further research is necessary to strengthen this conceptual model of metamorphic karst aquifers; in particular, flow, solute, and heat transport modelling could be applied to long time-series data to better constrain conduit network patterns.