The hydrological response of baseflow in fractured mountain areas

The hydrological response of baseflow in fractured mountain areas A. Millares, M. J. Polo, and M. A. Losada Environmental Fluid Dynamic Group, Andalusian Centre for Environmental Studies (CEAMA), Avda. del Mediterráneo s/n, 18006, Granada, Spain Hydrology and Hydraulics in Agriculture Research Group, University of Córdoba, Campus de Rabanales, Edif. Leonardo da Vinci, 14071-Cordoba, Spain Received: 19 March 2009 – Accepted: 28 March 2009 – Published: 22 April 2009 Correspondence to: A. Millares (mivalag@ugr.es) Published by Copernicus Publications on behalf of the European Geosciences Union.


Introduction
The analysis of subterranean responses and its quantification constitute, within the study of the hydrological cycle and all its processes, one of the commonest and, at the same time, most complex problems.Despite the heterogeneity which exists in the composition of baseflow (Hewlet, 1961;Tallaksen, 1995), it can be defined as the portion of surface flow originating from subterranean storage (Hall, 1968), where most of the input comes from saturated soil, and in most cases, from shallow, generally unconfined, aquifers (Wittenberg and Sivapalan, 1999).Their "hidden" nature and the problems in obtaining data for water storage in the subsoil present an added difficulty to their study when compared with surface techniques, where it is possible to obtain information either through direct methods (gauging sites, weather stations, etc.) or indirect methods (video image, satellite images, etc.).The problem is even more complex when it comes to identifying subterranean flow in mountainous areas and in basin headwaters, where the characteristics of the often fractured materials are very different to the standard issues concerning porous material applied in conventional hydrogeology.In particular, the flow that develops through fractures is very difficult to identify using direct methods, resulting in widely variable answers depending on the extent and depth of the fracture.It is an even more complex process in semi-arid areas, where a number of hydrological processes occur together: the presence of snow, high solar radiation, intense evaporation, the use of water for irrigation, and so on, all of which influence the final baseflow reading.
The only way to identify these elements and their relationship with other hydrological processes is by using indirect study methods.Traditionally, this type of study starts with the analysis of recession curves, on the assumption that, once the hydrogeological recharge process is over, the study of flow consists of measuring the relationship between storage and discharge, which allow us to identify the aquifer (Tallaksen, 1995).Generally, depending on the quality of data available, the hydrological system studied and on the prevailing climate, the approaches taken by these studies have differed widely in time perspectives.On the one hand, the recession study can focus on those long, continuous periods of time where the input to the system is zero or minimal, such as the dry season, evaluating the storage-discharge relationship continually, in the sure knowledge that the flows originate exclusively from the underground system (Tallaksen, 1995;Wittenberg, 1999;Chapman, 2002).On the other hand, however, in basins which are monitored infrequently, with short periods of flow and where a sufficiently lengthy dry period does not occur, recession fragments can be used, taking into account processes with a shorter time scale, such as surface runoff, episodes of snowmelt or evapotranspiration, in order to verify that their influence is not significant and that these factors do not "mask" the source of subterranean flow (Lamb and Beven, 1997;Fenicia et al., 2006).To use this methodology, it is essential to establish the relationships between selected fragments by using what is known as the Master Recession Curve (MRC), which summarizes the average characteristics of the underground system (Tallaksen, 1995;Hammond, 2006).In addition to this, when different procedures in the relationships between fragments are taken into account to generate the MRC, it is possible to identify differences in the response of these sub-surface contributions.
The main aim of this work is to propose a methodology based on MRC that can differentiate between types of response in mountain environments with fractured material.To do this, to the standard analysis of recession fragments is added a) a pre-hydrological analysis that selects the fragments representing the subterranean contributions, and b) the use of downward concatenation of these fragments to identify the faster response time scale which occurs in these materials.To achieve this, we propose the use of automatic algorithms for the selection and extraction of recession fragments, as well as concatenation algorithms for the generation of MRC, both in upward and downward form, based on the principles of automatic concatenation of recessions (Lamb and Beven, 1997), and using as support the information provided by a distributed physical model.In parallel, and in order to validate and contrast the results obtained, the recesses corresponding to the dry season will also be studied.
This methodology is applicable in watersheds that suffer from the problems stated above and where the characterization of the hydrological response in the inter-event periods is important.
This article is organized as follows: first, we will explain the methodology used to achieve the proposed objectives, along with its theoretical basis and underlying assumptions.Next, this methodology will be applied to the River Guadalfeo (in southeast Spain) as an example of high altitude headwaters, with the presence of snow every year and several melting cycles caused by the weather variations in the Mediterranean environment to which it belongs.

Recession analysis
To carry out a combined analysis of recession fragments with durations varying between three days and several weeks, the relationship existing between them needs to be assessed so that they can reproduce the combined response of the subterranean flow.For this they must be concatenated into a Master Recession Curve, which summarizes the characteristics of the hydrogeological system.This graph shows the average contributions made by the fragments, and they can be analyzed in upward or downward form (Fig. 1), with a view to identifying the time differences in the response from the subterranean systems.
In the upward pattern, the resulting Master Recession Curve primarily reflects the portion of the recess which corresponds to a slow flow and greater inertia in the system.Arranging the fragments in a descending pattern, the resulting curve would allow us to evaluate those recessions which belong to the subterranean discharge and take place considerably faster.This response may be due to the mobilization of the subterranean storage in the pre-event period (Chapman and Maxwell, 1996) and seems to contradict the traditional theories about the behaviour pattern of slow flowing subterranean waters (Wittenberg, 1999).In fractured mountain systems, the response is also connected with the input of second-order fractures and more superficial sources whose relation is different to primary fractures situated at a greater depth.
To obtain the selection of recession fragments from the flow record series and the subsequent correction of the influence of hydrological processes, supported by the use of a distributed hydrological model, calibrated and validated for the study area (Aguilar, 2006;Herrero, 2007;Millares, 2008), the following steps were taken: 1. Segments of continuous flow recess were selected from the gauge records at each study area.
2. In each case, fragments shorter than 4 days were ignored, following Hammond (2006) criteria to select recession fragments only when the average for 1-day flow falls continuously for a period of at least 4 days.
3. Fragments resulting from significant snowmelt events or occurring in periods of intense evaporation were ignored.To identify such processes influence, the study areas were modelled after Herrero et al. (2009) and Aguilar and Polo (2006), and water flows due to both snowmelt and evaporation were quantified on a daily basis along the study areas.A threshold in terms of mean values was assumed for each sub-basin studied; values equal or greater than 1 mm/day for snowmelt flow and fragments where the total amount of evapotranspiration exceeded the 10% of the total discharge (following the assumptions of Lamb and Beven (1997)) were disregarded for the analysis.
4. Fragments affected by precipitation excess were identified, and the intervals of direct runoff were excluded, by using the empirical relationship (Linsley et al., 1958), where N is the number of days between the storm crest and the end of runoff and A is the area of the basin studied in km 2 .
Figure 2a summarizes the method described.At each point studied, the MRC values obtained by upward and downward methodologies from the final set of fragments were compared.For this study an algorithm based on the works described by Lamb and Beven (1997), adding downward procedures, is proposed.This algorithm (Millares, 2008), implemented in Matlab, is illustrated schematically in Fig. 2b.
Once the selection of the fragments is made, the first step is to sort these fragments, ordered in both the upward and downward direction.The concatenation procedure considers recessions in pairs, relating the flow by scanning one fragment and establishing the relationship existing when the flow is as close as possible to the flow of the next recession.Once the relation is made, the recession scanned is moved step equal to in order to obtain a continuous curve in time.
In both cases, the MRC obtained was adjusted, following linear and non-linear models of subterranean behaviour, as shown in the following section.
In addition, the recession study has been validated through the analysis of long and continuous periods of time where the inputs to the system are zero or minimal, as in the dry season.The results have been used to validate the study of recession fragments in the upward form, since the evolution of the flow during these periods is related to deep storage and slow response of the subterranean system.To do this, the period in which the inputs could be considered zero or negligible was fixed, considering the hydrological variables influencing recharge in the study area.Those years when occasional episodes of summer rain or processes of snowmelt occurred were evaluated, and those dry season recesses which were rather unrepresentative of the subterranean storage were excluded from the analysis.

Linear and non-linear recession adjustments
The relationship between storage and discharge can be expressed by the expression (1) where S is the storage system (L 3 ), Q the discharged flow (L 3 /T ), and k, whose dimensions are (L 3−3b T b ), and b, which is dimensionless, are the empirical parameters.When b=1 the system is linear and expression (Eq. 1) becomes where k is termed the storage constant and is close to average time of residence.By combining expression (Eq.2) with the continuity equation, we reach the well-known equation (Eq. 3) (e.g.Bako and Hunt, 1988), whose parameters are obtained by linear adjustment of its logarithmic expression (3) Non-linear behaviour of the system b =1 reflects coupled interaction or a lack of synchronization in the hydrological processes involved, especially under arid conditions (Hamad, 1993).By combining equation (Eq. 1) with the continuity equation, the expression for recession flows proposed by Coutagne in 1948 is obtained, In this case, the parameters were adjusted following Wittenberg (1994), who proposes the use of the squared minimums method, systematically varying the value of parameter b, and calculating parameter k for each iteration ensuring that the volume obtained during the period of time in question is equal to that of the recession curve studied, thus fulfilling the relationship (Eq.5): Where Q i is the discharge for the time period i in a recession observed.

Study area
The river Guadalfeo rises on the southern slope of the Sierra Nevada mountain range in southern Spain (see Fig. 3).In the Sierra Nevada, it receives water from tributaries -on the south-eastern slopes from the rivers Trevélez, Poquéira and Cádiar, and on the south-western slopes from the rivers Ízbor and Lanjarón.The high altitude guarantees the presence of snow from November to June every year, and this decreases the flow of direct runoff during rainfall events and levels them out with intermediate melting cycles which, in turn, refill the aquifers and maintain the surface water levels during the dry season, although in dry years this can run out.
The average annual rainfall in the basin for the period 1950-2007 is 600 mm, with minimum and maximum values of 500 and 1000 mm respectively (Millares et al., 2008).At an altitude of over 2000 meters, over 75% of the annual rainfall occurs in the form of snow (Herrero, 2007).
The study area consists of the drainage area in the metamorphic materials of the Sierra Nevada range (Fig. 3).There are numerous published studies which analyze the nature of this particular aquifer which exhibits a significant volume of water flow throughout different configurations of fractures (Gisbert, 1997;Alwany, 1997;Adarve et al., 1997;Castillo et al., 2002).Alwany (1997) remarked the importance of considering snowmelt in the recharge of this aquifer, defining three flow types which depend on the settling of these fractures; deep flow, intermediate flow and shallow flow.Castillo et al. (1999) analyzed indirectly these flows from different sources and springs, distinguishing different responses related to the geological material and the nature of the fractures.In general, we can distinguish between: 1) a fast subsurface flow linked to the fringe of the deposits, with more superficial alteration and fracturing, common in shallow thicknesses, and 2) a deep, slow flow linked to the circulation through large fractures, tectonic accidents and interstratified layers of marble and quartzite (Castillo et al., 1999).This difference in underground flows depending on the depth is evident in the river-aquifer interactions and is reflected in the flow analysis, as can be seen below.In general, it can be stated that the flow type associated with these materials is a subsurface flow through fractured and topographically controlled areas.This flow is facilitated by the existence of deposits and highly permeable strips of glacial and periglacial alterations, the successive cycles of accumulation and melting of snow throughout the year, the artificial recharge of the system through numerous man-made irrigation channels and, to a lesser extent, by the type of vegetation cover in these mountain areas that favours infiltration (Castillo and Fideli, 2002).This leads to the existence of abundant sources; some provide starting points for a shallow subsurface flow, with a quick response, which goes on for a few days after the end of the summer thaw, while others have greater volumes of flow, have slower responses, never run out and maintain the surface flow of these mountain rivers during the summer.
The method described in the previous sections has been applied to three sub-basins belonging to the same geological formation, all influenced by snow dynamics, at whose closing points there is a daily series of flow (source: Andalusian Mediterranean Basin), lasting 38, 10 and 17 years respectively, which allows us to identify the subterranean contributions.Figure 4 shows the time evolution for these records, for each series, and their relative orders of magnitude, reflecting, among other factors, their differences in surface area (66, 74 and 33 km 2 , respectively).

Results
Figure 5 shows the results obtained with the hydrological model for rain, snowmelt and evapotranspiration for each sub-basin studied.
With the criteria described, a total of 171, 31 and 61 fragments at the observation points of Cádiar, Trevélez and www.hydrol-earth-syst-sci.net/13/1261/2009/ Hydrol.Earth Syst.Sci., 13, 1261-1271, 2009 Lanjarón, respectively, were chosen and corrected.A detailed study of the recessions obtained proves the existence of two distinct sections; one segment, with a steeper gradient, that coincides with quick responses in the subterranean system and another segment, with a lesser gradient, where only slower contributions, originating from deeper areas, occur.Figure 6a and b show both types of flow, from the concatenation of the chosen fragments, characterized by the MRC produced in upward and downward form at each gauging site studied.The linear and non-linear adjustments of these MRC curves are also shown in the figures and their parameter values are listed in Table 1.
For the slow flow (see Fig. 6 and Table 1), very clear differences can be observed between the Cádiar and Trevélez sub-basins, with k linear storage values of 17 and 15 days, respectively, and the Lanjarón sub-basin, with values of around 10 days.For the quick response, this linear adjustment shows similar results in the three sub-basins, with k values of around 10-12 days.
The results fit within the estimated interval for subterranean input in rivers with little flow, i.e. between 7 and 21 days (Smakhtin, 2001).As well as being within the same order of magnitude, we have been able to identify differences in the responses with the help of the methodology proposed in this work.
In addition, a set of dry-season recessions was selected at each of the stations studied, excluding those not valid for analysis.To be exact, 17, 5 and 9 recessions were obtained, which correspond to Cádiar, Trevelez and Lanjarón, respectively, as shown in Fig. 7.The results of the individual analysis of each recession were averaged to calculate an average storage coefficient for all the recessions corresponding to one station, and together, from downward concatenation processes.
Table 2 shows an overall view of the results of the dryseason study corresponding both to the average value of the individual adjustments of each recession and to obtaining the summer CRM.
As can be seen, the average results obtained from the individual study of each recession are consistent with the values corresponding to the CRM.In addition, the values for the sub-basins of and Cádiar and Trevélez, with values between 16-22 days, are clearly very different from the Lanjarón subbasin, which has lower storage values at around 7 days.

Discussion
The results obtained show very little inertia in the subsurface response.Other authors have reached similar conclusions in the same study area.The analysis of recession curves in the dry seasons on the western slopes of Sierra Nevada, made by Fideli and Castillo in 2002, showed a very small component of inertia, which is not surprising given the unique features of this mountain aquifer.Another important finding   was the tendency of the subterranean flow to impulsive responses after the occurrence of periods of intense snowmelt, with subsequent quick depletion.These results confirm the high permeability of the materials studied.In the study mentioned, conducted at a range of altitudes between 1000 and 1500 m, the depletion rates varied between 12 and 67 days, showing great heterogeneity in the response relating to the maintenance of baseflow from the many springs of different types and location.
Based on the results obtained in the study of recession fragments, the quick response in the three sub-basins shows a notable similarity.This uniformity in the quick response reflects the very similar characteristics in the surface fracturing of the soil throughout the area of study.The linear results for the Cádiar and Trevélez sub-basins allow us to identify a quick response (11 days) and a slightly slower response (17 days), characterized by different values of upward and downward MRC.For the Lanjarón sub-basin, however, there is no such difference in response, with similar values of around 10.5 days obtained in both concatenation criteria.
The results obtained from the dry-season recessions bear a great resemblance to the study of recession fragments of upward methodology, taking into account the significant methodological differences in the selection process of recessions from the same data set.Only the Trevélez sub-basin shows differences, ranging from 15 to 22 days, taking into account the few recessions selected in the dry season at this station.
In this sense, the homogeneity between the estimated fast and slow response for the Lanjarón sub-basin, the most westerly of the three, points to a possible dominance of its geomorphology in the subterranean inputs, with a lesser amount of deep storage.
Different authors have estimated average values of recession constants for different geomorphologic regions or geological formations, in some cases establishing quantitative A. Millares et al.: The hydrological response of baseflow in fractured mountain areas relationships between the recession parameters and the characteristics of the sub-basins (Tóth, 1970;Zecharias and Brutsaert, 1988;Demuth, 1989;Tallaksen, 1989).In these studies, the main indices with influence on the storage-discharge relationship are related to geological features and climate.However some authors (Mwakalila et al., 2001) remark the influence of topographic and geomorphologic indices on the baseflow patterns, like the drainage density, defined as the ratio of the length of the main stream to the square root of the catchment area, or the Gravelius index, defined as I c =0.28×P / √ S (where P and S are the perimeter (km) and catchment area (km 2 ) respectively).These indices can be linked to the storage parameter k, calculated for a sub-basin where linearity is assumed.In Fig. 8, the values of area, Gravelius index and drainage density are plotted, and its relation to the storage coefficient k for the three basins studied can be observed.The marked differences between the Lanjarón sub-basin and the other two can be seen, and there is a clear link between the area of these sub-basins and the storage coefficient, as well as in the shape of the sub-basins and drainage density.
However, the number of sub-basins studied is too small to draw conclusions about the relationship between the geomorphology of the sub-basins and the characteristics of the recession in terms of linearity of the system, and so it would be interesting to extend the study to other sub-basins that share the same geological features in order to study this relationship more deeply.
In both sets of recessions, the linear fit provides an adequate approximation, but the results of the non-linear adjustment allow us to analyze the behaviour of the subsurface flow of each zone.
The non-linear adjustment shows values for the dimensionless parameter b of around 1.25, always greater than 1, which reflects a very low concavity of the adjusted MRC, and differs considerably from the results found by other authors, with b values of about 0.5 (Wittenberg, 1999;Chapman, 2002); this shows that the characteristics of the aquifer are very different from those of a standard unconfined aquifer in porous isotropic media, with homogeneous conductivity and river-aquifer convergence.The end causes of this unusual behaviour are many; the existence of a subsurface flow caused by fracturing could be affected by problems of turbulence similar to those caused by macropores (Wittenberg, 1999).
The negligible difference found between the linear and non-linear adjustment allows us to confirm how close the system is to displaying linear system behaviour, taking into account that the chosen fragments have been corrected for those external factors which are mainly responsible for nonlinearity in the subterranean response, such as evapotranspiration, which induces underestimation of k values (Wittenberg and Sivapalan, 1999), or snowmelt, which causes overestimation.

Conclusions
The selection and correction of recession fragments with a hydrological model allow us to use these fragments outside the dry season.This makes more initial information available in basins such as the Mediterranean basin, where, in addition, it has enabled us to understand more about a fast system response originating from a subterranean storage area where it is difficult to explain the evolution of the base flow throughout the year due to the fractured nature of the materials.
The downward concatenation of these fragments proposed in this article has proven to be a valid method for quantifying the dynamics of these quick responses in the basins studied, and has allowed us to identify common patterns associated with the surface morphology of the fractures in the materials that make up the study area.The use of traditional upward criteria to obtain the MRC has yielded coherent coefficient storage values in the three sub-basins studied.All in all, the comparison of methodologies starting from recession fragments or dry-season recessions enables us to compare and validate the results for the storage-discharge relationship from two different approaches.In the study area, the values for the Cádiar sub-basin were very similar when seen from both approaches with a value of k=17 days, with slight differences in the Lanjarón sub-basin (7>k>10 days) and the most marked differences in Trevelez (15>k>22), which is most likely to be due to the smaller amount of data available.
All in all, the methodology chosen to differentiate quick response subsurface flow through fractured materials has yielded satisfactory results and has identified this flow, with 10-day patterns following the recess produced in deeper areas of the aquifer, which exhibit slightly slower behaviour.
The linear adjustment results also allow us to adopt this model to characterize both quick and slow responses in these fractured materials.In addition, the high rates of non-linear b obtained in this work (between 1.1 and 1.3), may indicate artificial recharging processes via the network of irrigation channels.This information could be used as a first step towards quantifying these flows, which have been diverted from their natural courses, and about which we have no information whatsoever at present.
Finally, the relationship between BF I (Base Flow Index), defined as the relationship between river flow and the fraction corresponding to baseflow, and geo-morphological indices could be applied in the future, once these indices, as well as the relationship between them and the storage coefficient k have been obtained in the largest number of sub-basins.
The methodology proposed is valid and applicable in basins with fractured materials that cause a delay in generating runoff due to the initial storage which the size of the fracture and its subsurface movement produce.The results obtained in this work will be included in a full distributed hydrological model, aimed at Mediterranean basins with steep topography; in particular, they will allow us to configure the model for subterranean inputs in the form of baseflow throughout the year.

Fig. 1 .
Fig. 1.Different methodologies of concatenation depending on the relationship between the fragments.

Fig. 4 .
Fig. 4. Daily flow registered at the three observation points in the study area.

Fig. 5 .
Fig. 5. Daily results for the period 1969-2004 for rainfall, snowmelt and evapotranspiration obtained with the hydrological model used in this study.

Fig. 8 .
Fig. 8. Relation between area, Gravelius index and drainage density and the slow response storage coefficients in the three sub-basins studied.

Table 1 .
Results obtained after analysis of recession fragments.

Table 2 .
Linear and non-linear adjustment of the dry-season recession study by sub-basins.