Sediment yield model implementation based on check dam infill stratigraphy in a semiarid Mediterranean catchment

Soil loss and sediment transport in Mediterranean areas are driven by complex non-linear processes which have been only partially understood. Distributed models can be very helpful tools for understanding the catchment-scale phenomena which lead to soil erosion and sediment transport. In this study, a modelling approach is proposed to reproduce and evaluate erosion and sediment yield processes in a Mediterranean catchment (Rambla del Poyo, Valencia, Spain). Due to the lack of sediment transport records for model calibration and validation, a detailed description of the alluvial stratigraphy infilling a check dam that drains a 12.9 km2 sub-catchment was used as indirect information of sediment yield data. These dam infill sediments showed evidences of at least 15 depositional events (floods) over the time period 1990–2009. The TETIS model, a distributed conceptual hydrological and sediment model, was coupled to the Sediment Trap Efficiency for Small Ponds (STEP) model for reproducing reservoir retention, and it was calibrated and validated using the sedimentation volume estimated for the depositional units associated with discrete runoff events. The results show relatively low net erosion rates compared to other Mediterranean catchments (0.136 Mg ha −1 yr−1), probably due to the extensive outcrops of limestone bedrock, thin soils and rather homogeneous vegetation cover. The simulated sediment production and transport rates offer model satisfactory results, further supported by in-site palaeohydrological evidences and spatial validation using additional check dams, showing the great potential of the presented data assimilation methodology for the quantitative analysis of sediment dynamics in ungauged Mediterranean basins.


Introduction
Modelling sediment yield is a complex task due to the nonlinearity of natural processes intervening at slope and basin scale (Schumm and Lichty, 1965;Coulthard et al., 1998;Roering et al., 1999).Recent computing advances, together with a better understanding of hydrodynamic processes involved in the surface runoff, sediment production and sediment transport, have stimulated the development of physically based and distributed parameter models (e.g.WEPP, EUROSEM and LISEM).The reliability of such sediment yield models depends on a robust calibration and/or validation process that, at ungauged catchments, as it is the case of most of small basins around the world, may limit a broad use of such models.Different authors have used the sediment volume accumulated in lakes and reservoirs as an indirect validation method for modelling sediment yield at regional scale (Van Rompaey et al., 2003;Grauso et al., 2008).Reservoir sediment volumes have been used since the 1950s as an estimate of the catchment mean sediment yield for comparison with the results of empirical equations.Some examples are Geiger (1957), Ackermann and Corinth (1962), Rohel (1962), Farnham et al. (1966), Callander and Duder (1979), Jolly (1982), Le Roux and Roos (1982), Duck and McManus (1993), Avendaño Salas et al. (1995, 1997) and Verstraeten et al. (2003).
Recently, sediment volumes stored in water-retention dams have also been used for distributed mathematical model validation, as shown in De Vente et al. (2005, 2008) and Alatorre et al. (2010).In De Vente et al. (2005), two Published by Copernicus Publications on behalf of the European Geosciences Union.

G. Bussi et al.: Sediment yield model implementation
semi-quantitative models for mean annual net erosion rates estimation are compared, and their results contrasted versus reservoir sedimentation rates.In De Vente et al. (2008), the reservoir sedimentation rates were used to compare the results of three distributed approaches for soil erosion rates and long-term sediment yield estimation: the WATEM/SEDEM model (Van Rompaey et al., 2001), the PESERA model (Kirkby et al., 2008) and the SPADS model (De Vente et al., 2008).In Alatorre et al. (2010), the WATEM/SEDEM model is calibrated using the depositional record of the Barasona reservoir (NE Spain) and then used for Ésera River catchment sediment yield modelling, providing mean annual erosion and sediment yield.
Not only large reservoirs, such as the Spanish data set mentioned above (Avendaño Salas et al., 1997), can be used for sediment yield estimation, but also smaller reservoirs like check dam reservoirs, irrigation and water supply ponds, etc., can be in a similar way, a source of information.Verstraeten and Poesen (2000) quantify the number of these structures around the world in a few million dams, with small dams being defined as retention structures with a storage capacity of between 50 and 5 × 10 6 m 3 .This large number of small reservoirs in the world is a high potential source of information for sediment yield assessment and modelling.Some examples of sediment yield studies based on sedimentation volumes in small reservoirs include McManus and Duck (1985), Van den Wall Blake (1986), Neil and Mazari (1993), Foster and Walling (1994), White et al. (1996), Romero-Díaz et al. (2007), Boix-Fayos et al. (2008), Sougnez et al. (2011) and Bellin et al. (2011).Verstraeten and Poesen (2002) calculated the error on sediment yield estimation for 21 catchments located in central Belgium using small reservoir deposits and concluded that this is a suitable methodology for medium-sized catchments (10 1 -10 4 km 2 ) and for midterm sediment yield estimations (10 0 -10 2 yr).Errors on topographical surveys, sediment dry bulk density and reservoir trap efficiency must be taken into account, although the mean accuracy of this information is comparable to other methodologies used for sediment yield estimation, such as sediment rating curves or suspended sediment sampling.
Reservoir sedimentation rates are a very helpful tool for estimating catchment sediment yield, but indeed this methodology has some weaknesses (Foster, 2006): (i) the quality of the reservoir storage capacity estimation is sometimes questionable, especially for the starting reservoir capacity (i.e. when the reservoir was built); (ii) the calculated sediment yield is averaged over an extended time (Alatorre et al., 2012); (iii) the total deposition volume does not give information about temporal patterns (event sediment production) and their variability.In the case of large reservoirs and artificial lakes, sediment coring and paleolimnological techniques, including geochronological dating (Cs-137, Pb-210), have been used for temporal characterisation of sediment rates, as in the artificial Lake Matahina in New Zeland (Phillips and Nelson, 1981) or in the Brno reservoir in Czech Republic (Nehyba et al., 2011).
In Mediterranean ephemeral streams a large number of check dams have been built to prevent or reduce sediment inputs into perennial streams during the first winter or rainy season following a wildfire (Boix-Fayos et al., 2008) or to correct local channel slope (Romero-Diaz et al., 2007).In these check dams, infill deposits record historical pulses of sediments produced during discrete flood events after their construction.The coarse texture of the deposited material prevents coring, but allows the use of fluvial palaeohydrological techniques.The detailed analysis of their alluvial stratigraphy may provide quantitative information for specific events, such as the number of events, timing, and deposited volume(s) of an individual flood or floods.Similar techniques have been used in the reconstruction of the magnitude and frequency of past floods using geological evidence (Kochel and Baker, 1982;Baker, 2008;Benito et al., 2010;Machado et al., 2011).
The aim of this paper is to present a new methodology for model calibration and/or validation in catchments with no sediment data availability, by taking advantage of the flood sediment proxy information obtained from check dam infills.An existing hydrological and sediment model, the distributed TETIS model (Francés et al., 2002(Francés et al., , 2007;;Montoya, 2008), is used to illustrate the proposed procedure.With this model we aim to reproduce the hydrological and sediment regime of a small semi-arid Mediterranean basin (Rambla del Poyo, Valencia, Spain).Downstream-measured water discharge records at a gauge station (with a catchment area of 184 km 2 ) are employed to calibrate and validate the hydrological sub-model; the sediment trapped in a small reservoir, draining a 12.9 km 2 upstream subcatchment, is used to calibrate the sediment sub-model.Time-variable trap efficiency is taken into account by coupling the STEP model (Verstraeten and Poesen, 2001) with TETIS.The sediment yield temporal validation is carried out by comparing the model output with the results of the stratigraphical description (unit thickness, minimum sediment volume and texture) of depositional sequences observed in two trenches dug across the reservoir sediment infill.The sedimentary infill shows evidence of 15 flood events that occurred following the dam construction in 1992 up to year 2009.A date is assigned to each flood unit based on layers containing anomalous, high charcoal accumulation due to well-documented historical wildfires.Finally, the model is spatially validated by estimating the sedimentation of seven check dams located in the catchment headwaters.

The TETIS model
The TETIS model is a distributed and conceptual model for water and sediment cycle simulation, whose majority of parameters have a physical meaning.The TETIS model was developed in 1995 and has also been extensively used in Spain and Latin America, both for research and consulting purposes (the most recent examples are Vélez et al., 2009;Andrés-Doménech et al., 2010;Francés et al., 2011;Salazar et al., 2013).It has proven to be a reliable hydrological model for semi-arid catchments and other climate conditions.The TETIS model is composed of a hydrological sub-model and a sediment sub-model.
The hydrological component of TETIS has already been presented in many literature studies, such as Francés et al. (2002Francés et al. ( , 2007)), Morales de la Cruz and Francés (2008) and Vélez et al. (2009).Each cell of the catchment spatial grid is conceptualised by five connected tanks, as shown in Fig. 1.The connections between tanks are represented by linear reservoirs and flow threshold schemes.
The first tank (T1) corresponds to the sum of the upper soil capillary retention and surface and vegetation interception and it is called static storage; its only exit is evapotran-spiration.The second tank (T2) reproduces the surface water, i.e. the part of precipitation which generates overland flow.The third tank (T3) corresponds to the gravitational storage of the upper soil; it generates the interflow.The fourth tank (T4) corresponds to the saturated soil or aquifer, which produces the base flow.The percolation process is modelled according to both soil saturation conditions and vertical hydraulic conductivity, and the remaining water in T3 is available to feed the interflow.The last tank (T5) represents the gully or river channel, i.e. the stream network storage within each cell.The flow-routing in the stream network is carried out by the Geomorphologic Kinematic Wave methodology, employing nine geomorphologic parameters obtained from power laws (Francés et al., 2007) and estimated by geomorphologic regional studies.The TETIS hydrological component also includes an automatic calibration module based on the SCE-UA algorithm (Duan et al., 1992(Duan et al., , 1994)).The automatic calibration of the hydrological parameters is carried out adjusting up to nine correction factors (called CFs) in order to fit the observed hydrographs.CFs globally correct the parameter maps (static storage, hydraulic conductivity of different soil layers, and surface and channel flow velocities), reducing the number of variables to be calibrated (Francés et al., 2007).It is expected that this split-structured effective parameters will help simplifying the model calibration process, which is a major advantage in poorly gauged catchments, such as many semi-arid catchments.Nevertheless, the methodology presented in this study should be applicable with any model, without the need to reformulate.
The TETIS sediment component was implemented by Montoya (2008) and it is based on the conceptualisation of the CASC2D-SED model (Johnson et al., 2000;Ogden and Heilig, 2001;Rojas, 2002).This conceptualisation is based on the balance between sediment availability and flow transport capacity.Fine sediment transport is limited by sediment availability, while coarse material transport is limited by flow transport capacity (Julien, 2010).The TETIS model divides sediment flow into three textural classes (sand, silt and clay), assigning to each of them a representative diameter and settling velocity.
The hillslope sediment erosion and transport processes are described by means of the Kilinc and Richardson equation (Kilinc and Richardson, 1973) for the total transport capacity.The sediment discharge per unit width in terms of weight is given by where Q [m 3 s −1 ] is the cell overland discharge, W [m] is the cell width, S o is the terrain slope [m m −1 ] and α a dimensional and empirical parameter (around 25 000 for sandy bare soil with the expressed units).This equation follows a standard formulation for overland flow sediment transport equations, as illustrated by Julien and Simons (1985) and Prosser and Rustomji (2000).Julien and Simons (1985) found that the exponents varies between 1.2 and 1.9 (slope exponent), and 1.4 and 2.4 (discharge exponent), being the mean values 1.66 and 2.035 respectively, i.e. the values used by the Kilinc and Richardson equation.Prosser and Rustomji (2000) recommended the ranges 1-1.8 and 0.9-1.8,suggesting the median values (1.4 for both exponents) in case of using a single value.In the TETIS model, the original exponents of the Kilinc and Richardson equation were maintained, while the α coefficient, which was found to be the most influential parameter of the equation in this study, was chosen to be calibrated, as mentioned below.
The Kilinc and Richardson equation has been modified by Julien (2010) in order to consider land use, cropping management and soil characteristics.The modified equation is the following: where γ s is sediment-specific weight [t m −3 ], and K, C and P are the USLE (Wischmeier and Smith, 1978) soil erodibility, cropping management and support practice factors [-], respectively.Notice that 0.15 is the K value for sandy soil.Q h is divided into three parts, depending on the textural composition of the transported material, and each transport capacity part is used to route downstream the corresponding soil textural class of suspended sediment.If there is residual transport capacity, the deposited sediments are mobilised; if there is still residual transport capacity, the parental material is eroded.In the downstream cell, the sediments are separated into suspended or deposited depending on settling velocity, cell residence time and flow depth.
The stream network erosion and transport processes are described by the Engelund and Hansen equation (Engelund and Hansen, 1967), where the transport capacity depends on hydraulic radius, flow velocity, friction force and grain characteristics.The maximum sediment concentration is given by where and β is a non-dimensional calibration coefficient (not existing in the original expression).So, the streamflow transport capacity for the textural class i is expressed as follows: where Q is the stream channel discharge [m 3 s −1 ].Sediments are eroded and routed downstream following the same scheme as for hillslope processes.However, stream network parental material is not considered, because in many cases the most relevant source of channel sediments are deposits left by previous floods (Piest et al., 1975), which can be introduced as initial conditions in the model.In any case, parental material can be simulated as a large sediment deposit at the beginning of the simulation.The calibration of the sediment sub-model is carried out by adjusting the α and β values.Only these two values are calibrated within the whole sediment sub-model in order to avoid overparametrization.
In order to compute the effects of small retention structures such as check dams on the sediment transit, the problem of estimating trap efficiency was taken into account.For large reservoirs, empirical formulae such as Brown (1943) or Brune (1953) can provide very likely values, especially if the reservoir trap efficiency is close to 100 % (Bangqi Hu et al., 2009).Nevertheless, the uncertainty in calculating sediment trap efficiency is higher for smaller reservoirs.Moreover, trap efficiency decreases when deposition takes place and this phenomenon it is not taken into account by empirical formulae.For this reason, theoretical models have been developed to predict trap efficiency; a review can be found in Verstraeten and Poesen (2000).In this study, TETIS was coupled with the Sediment Trap Efficiency for Small Ponds model (STEP, Verstraeten and Poesen, 2001) in order to reproduce the sedimentation dynamics of a small reservoir or check dam during a flood event.The advantages of this model are: (i) its conceptualisation is simple and parsimonious, and requires a small amount of data; (ii) it was developed for small ponds, similar to the check dam analysed in this study; (iii) it was developed for long-term trap efficiency estimation and for continuous simulation; (iv) it is compatible with TETIS, given that it only requires the inlet water and sediment discharge as input, which can be provided by this model.As in other models, STEP divides the pond into several finite volumes (originally called chambers) on the basis of equal surfaces.Water is routed from volume to volume by means of a simple mass balance equation.The original STEP model can take into account the effect of the reservoir bottom slope, although complex hydraulic conditions, such as an irregular pond shape, cannot be reproduced (Cheng, 2008).In this study, the STEP model was adapted in order to model any reservoir by introducing an elevation-volume curve.
Sediment routing is also carried out by solving a mass balance equation for each volume, considering settling velocity of each particle and residence time (Chen, 1975).Then the proportion of deposited/transmitted sediments is computed for each volume.In TETIS, the STEP model was implemented as an online sub-routine, interacting at each time step.

The study area
The study area is the Rambla del Poyo catchment, a Mediterranean ephemeral stream located 30 km west of Valencia (Spain), as showed in Fig. 2. A stream gauge, located at the basin outlet, and a rain gauge, located at the same place, provide 5 min resolution discharge and precipitation series; the available data begin in 1988.
The geology consists of dolomites and limestones in the headwaters and marls in the lower part of the catchment (Camarasa Belmonte and Segura Beltrán, 2001).The mean annual precipitation is 450 mm and the mean annual potential evapotranspiration is 1100 mm, i.e. it is clearly a catchment with semi-arid conditions.The Rambla del Poyo catchment at the stream gauge station has an area of 184 km 2 (Fig. 2).The upper part, or headwater, located at west, is formed by high slopes and reliefs up to 1080 m a.s.l.Continuous land abandonment during various decades of the 20th century and the previous intensive livestock use, along with frequent wildfires, favoured the development of a rather homogeneous and dense shrubland cover (matorral), leaving only a little portion of pine forest.This phenomenon was already described in literature, e.g. by Cerdà (1998a), Rey-Benayas et al. ( 2007) and Baeza et al. (2007).
The intermediate part of the catchment is mainly nonirrigated arable land with complex cultivation patterns and transitional zones, with presence of terraced fields.The catchment lower part is dominated by agricultural land and alternation of urban zones and peri-urban agriculture, mainly orchards and citrus (Salazar et al., 2013).The stream network is composed by three major water courses: the Barranco Grande in the north, the Barranco de la Cueva Morica in the centre and the Barranco del Gallo in the south (Fig. 2).
Soil data for estimating hydrological and sedimentological model parameters was taken from LUCDEME project (Rubio et al., 1995).53 soil profiles collected within the LUCDEME project and located in or close to the Rambla del Poyo catchment were selcted, within an area of around 1500 km 2 .The topsoil texture is silt loam for headwater soils (around 20-60-20 % of sand-silt-clay percentage respectively, following the USDA classification).The content of sand gradually increases from west to east.In the lower catchment, the mean texture is clay loam (around 30-40-30 %), although sandy loam areas can also be found.Soil texture data, organic matter content and soil salinity data were used to feed the Saxton and Rawls (2006) pedotransfer functions and to obtain available water content and saturated infiltration capacity.Percolation capacity was estimated by reclassifying the lithological map, considering permeability values taken from literature (Fig. 3).
For C and K values, we employed data from a previous study done by Antolín (1998).The corresponding maps are also shown in Fig. 3.In that study, the C value was estimated as a function of the vegetation type and the cover density, following the guidelines given by Wischmeier and Smith (1978) and Dissmeyer and Foster (1984).The K factor was estimated from soil analysis (texture, organic matter and salinity) of a data set covering the whole Valencia Region, using the equation proposed by Wischmeier and Mannering (1969), and then interpolated in space.The C and K value were not changed during our calibration process.
Within the Rambla del Poyo catchment, the highest C values (0.32) are located in the headwaters, and correspond to the less dense shrubland areas.The lowest values (0.1) are also located in the headwaters, and correspond to small areas of pine forest that have survived the historical wildfires.The intermediate and lower parts of the catchment are characterised by values around 0.2.The K values decrease from the lower floodplain towards the headwater, mainly due to the variation of the soil sand content: up to 50 % in the lower part and decreasing up to 10 % in the headwater, following Rubio et al. (1995).The main statistics of C and K factor are shown in Table 1.The P factor of the USLE was set to 1, because no support practice is implemented.All maps were resampled to a 100 × 100 m mesh.This mesh was chosen as a compromise between map accuracy and computational time.
The effect of two wildfires that occurred in 1994 and 2000 was taken into account by modifying the C USLE factor during a "window of disturbance" (as defined by Prosser and Williams, 1998) after each wildfire.As demonstrated by Campo et al. (2006) for a plot located close to the Rambla del Poyo catchment, this effect seriously increases the soil erosion, especially when a severe rainfall occurs a few days or weeks after the wildfire.The length of the windows of disturbance was chosen following Andreu et al. (2001), so that the highest susceptibility to soil erosion of burnt areas takes place at the first severe rainfall event after the wildfire, usually in autumn/winter (Shakesby, 2011).A wildfire usually increases the erosion rates due to the reduction of the infiltration rate, the increase of surface runoff and the increase of soil erodibility (Cerdà et al., 1998b), among other effects.It has been demonstrated that, for Mediterranean shrubland plots, the site recovery is fully achieved after 2-4 yr (Cerdà, 1998c) and that the most important runoff and erosion rates alteration occurs within a few months after the fire, depending on the precipitation (Cerdà, 1998b, c, Andreu et al., 2001).An example of this phenomenon was described by Cerdà and Lasanta (2005).These authors demonstrated that in a shrubland catchment under natural conditions the erosion rates are low (0.04-0.1 Mg ha −1 yr −1 ), but they can increase up to 10 times within the 2-3 yr after the fire.For the Rambla del Poyo catchment, no information about fire intensity, duration or ash production was available.In order to reproduce the effect of erosion increase during a window of disturbance, the 1994 and 2000 wildfires were modelled by increasing the C factor of the USLE (similarly to what was done by Rulli et al., 2005) to 0.9 for the extreme events that occurred in the following 2 yr (December 1995 for the 1994 wildfire and October 2000 for the 2000 wildfire), corresponding with the highest peak of erosion increase.The error introduced with this approximation will be corrected by model calibration (by adjusting α and β coefficients).
The infiltration capacity was not modified during the windows of disturbance.This was done because the hydrological sub-model was correctly calibrated and validated (as will be shown below) without taking into account the eventual infiltration capacity decrease.Furthermore, there was no evi- dence of such an effect when comparing simulated and observed water discharge series.A possible reason for the absence of infiltration decrease is the one suggested by Cerdà and Doerr (2008), who stated that, after a forest fire, the layer of ash may compensate the effect of infiltration capacity reduction.In the Rambla del Poyo catchment we observed a similar behaviour: the comparison between the model results and the observed water discharge suggests that no sudden increase in runoff took place immediately after the wildfires.

The check dams
Indirect evidence of sediment production at the study area was provided by sediments stored on a small concrete check dam (Fig. 4), built in 1992 and draining a catchment area of 12.9 km 2 (the grey area in Fig. 2).The reservoir maximum storage capacity is 3000 m 3 , and at the time of the field survey was about half of its total volume capacity.This dam was chosen for its high siltation rate and for its accessibility.
A field study was carried out to survey the dam body and to describe the infill flood stratigraphy, including the collection of sediment samples for textural analysis of the sedimentary sequences.A topographic survey was also carried out at the sedimentation area and the surrounding slopes, in order to better estimate the deposited volume.An RTK (Real Time Kinematic) differential GPS was employed.The global scale error was around 6 m, while the relative accuracy was 1.5 cm (horizontal accuracy) and 1.8 cm (vertical accuracy).
Two trenches were dug across the reservoir sedimentation infill, at 9.5 and 22 m from the dam respectively, called BG-1 and BG-2 (Fig. 5).Detailed stratigraphic panels were carried out at centimetre resolution, using a one meter side vertical grid over the trench, where stakes with reference numbers were put at regular intervals.In these panels all the depositional contacts were traced laterally, with emphasis on breaks that indicate sedimentary interruption and post-flood surface exposure.The panels allowed a better detection of lateral interfingering beds, potential erosion and depositional gaps.Correlation between these two panels was possible due to anomalous charcoal content of some reference alluvial beds.Sediment samples of each unit were collected for determination of organic matter and micro-charcoal content (as indicator of fires) and for a complete textural analysis.
Starting from the reservoir geometry, the GPS survey and the flood unit morphology, the volume of each alluvial layer deposited by an individual flood was estimated by using two different methodologies: (i) "wedge" approach: given the layer depths inside the trenches, the sedimentation length and the distance between trenches, every layer volume was calculated as if each flood unit had a pyramidal shape (such as a wedge); (ii) proportional approach: given the surface shape and the layer depths inside the trenches, each layer volume was estimated by subtracting to the actual deposits the average accumulated layer depth, considering the thickness difference of each layer in each trench and approaching it to a pyramidal shape.
A complete textural analysis was also carried out for each flood unit in BG-1, showing that sediments are mainly composed by sand, whose percentage varies between 77 and 99 %.
The dry bulk density of the infill deposits was estimated in 1.195 gr cm −3 using the approach suggested by Lane and Koelzer (1943) based on textural data and coefficients for dry reservoirs.In order to validate this value, five measurements of the dry bulk density were carried out at different depths (from 10 to 90 cm).The results range between 1.014 and 1.389 gr cm −3 , and the mean value is 1.150 gr cm −3 .It is expected that the average dry bulk density of the whole deposit should be slightly higher than the measured value, given that the total depth is around 2 m.For this reason, the dry bulk density value calculated by means of the Lane and Koelzer approach can be considered more adequate for the purposes of this study.
Other seven check dams, located in the headwater, can be identified in Fig. 2. No stratigraphical analysis was carried out on these check dams, as they will only be used for spa-tial validation (see Sect. 3.3).Their characteristics are similar to the check dam described above.Their drainage areas ranges between 2.3 and 16.6 km 2 , the storage capacity ranges between 1200 and 23 000 m 3 and the dry bulk density of their deposits, calculated with the Lane and Koelzer formula, ranges between 1.19 and 1.25 gr cm −3 .The check dams show different degrees of filling, comprised between 2 and 90 % of their total capacity.For the sake of clarity, these check dams will be referred to as "secondary check dams" while the check dam subject to the stratigraphical description will be referred to as "main check dam".

Hydrological calibration and validation
The first step of this study was the calibration and validation of the TETIS hydrological sub-model.It is proven that, in Mediterranean climate, only a few events scattered over a large time period are responsible for most of the total sediment load (see for example Gallart et al., 2005).Given that our aim is sediment yield modelling, calibration and validation processes focused on the reproduction of heavy rainfall events.Two models with different time discretization were implemented within this study: with 5 min time step and with a daily time step.The first one (the main model) was used for reproducing with a fine time step only the flood events.The second one simulated the no-flood periods, with the objective of estimating the soil moisture initial condition for each flood event.This strategy was employed in order to save computational time.
The daily scale model was automatically calibrated from October 2000 to October 2003, and used for the estimation of the initial soil moisture state for 38 flood events.Its results are shown in Fig. 6.The calibration gave a Nash-Sutcliffe index (NSE, Nash and Sutcliffe, 1970) equal to 0.82, and the validation (from 1998 to 2010 excluding the calibration period) provided a 0.72 NSE index, which can be considered a "very good" daily calibration and a "good" daily validation, following the performance classification given by Moriasi et al. (2007).
Using the initial soil moisture given by the daily model simulation, the 5 min time step model automatic calibration was also carried out.The 5 min time step model was calibrated on a single storm event in October 2000, by means of the TETIS automatic calibration algorithm and validated on 37 storm events from 1990 to 2009.The simulated hydrograph for the calibration event is shown in Fig. 7 (left).The obtained NSE index is 0.78, with a volume error of −10 %, which can be considered as a very accurate model output (a "very good" model performance, following Moriasi, 2007).The temporal flood validation also provided good results.In Fig. 7 (right) the January 1998 event is shown.For this event, the NSE index is 0.5 and the volume error was  24 %; although there is a time shift between the observed and simulated peaks, maybe due to a poor description of the rainfall spatial distribution, this validation result can be judged as "good", given that the hydrograph shapes are very similar and the peak discharge error, which is very relevant for sediment yield modelling, is relative small (15 %).For all other flood events the model also obtained acceptable performances in terms of NSE index, volume and peak discharge errors and hydrograph visual fit.No spatial validation of the hydrological sub-model was carried out due to the absence of other stream gauges located within the Rambla del Poyo catchment.
The model also reproduced satisfactorily the ephemeral behaviour of the catchment; the base flow is absent, and the channel flow is composed mainly by overland flow with a little contribution of interflow (for heavy flood events, only 0.1 % of total flow), which is in accordance with our prior catchment hydrological knowledge.The model tended to provide a good estimation of the high peak flows, while the error on small intensity events was greater, probably due to the use of spatially uniform precipitation (only one rain gauge, located at the gauge station, was available for the whole catchment).The initial soil moisture estimation by a warmup simulation period at the daily scale was proven to be suitable, although some small error can be detected, as for example in the Fig. 7 (left), where the first peak of the flood is underestimated, probably due to an underestimation of the initial soil moisture.

Alluvial infill description and volume estimation
The geometry of the dam alluvial infill can be described as a sedimentary wedge with a triangular plan view (Fig. 8).The active channel is bordering the right margin of the reservoir area, partially undercutting the slope deposits.In the uppermid reach of the reservoir the most relevant morphosedimentary feature is a lateral gravel bar (43 m in length) attached to the left valley side, with a prominent 1m high frontal scarp, indicating a progradation over the fine deposits located closer to the dam.This alluvial bar is composed by poorly sorted gravels and boulders in a matrix of sand and silt, with a lack of structure, suggesting a deposition by flash flow(s) associated with detrital heavy load and loss of energy due to slope reduction, caused by the previous dam infill.Closer to the dam wall, the alluvial infill comprises a 2.5 m-thick deposit composed by multiple layers of well sorted sands and silts with ripples, planar and cross-bedded lamination and parallel lamination.The geometry of the layers is horizontally close to the reservoir centre but increases its elevation and decreases its thickness towards the valley side.Two trenches were dug across the reservoir sedimentation infill.The first trench (BG-1: Figs. 5 and 8) was about 10 m in length, covering from the left valley side to the main channel at the right margin, by 2.5 m in depth exposing sequences of multiple fine-grained flood deposits linked to the development of an eddy flow behind of the over-elevated left part of the dam wall.The stratigraphic sequences found in BG-1 (Fig. 9) provided evidence of at least 15 individual floods units post-dating the dam, as the dam infill fine sediments overlay old slope and stream channel gravels (Fig. 8).The lower seven flood units suggest a period of relatively small floods, on the basis of the very fine and fine sand grain size and thin stratigraphic layers.The upper part of the sequence is represented by eight flood layers of medium to coarse sand within units of 20 to 60 cm in thickness, with parallel and planar cross-stratification indicating a higher energy and sediment load than the lower flood units.The flood units 3 to 10 contain a large amount of charcoal debris concentrated on distinct 1 to 2 cm-thick laminae, which were probably deposited after severe wildfires that occurred in 1994 and 2000.
The second trench (BG-2, Fig. 8) is about 8.5 m long by 2 m deep, and it was excavated across the reservoir infill 12.5 m upstream of BG-1.In the lower 1 m section, at least eight flood units were distinguished and correlated with the lower ten flood units in BG-1, with exception of units 5 and 8 that pinched out at some point between both trenches (Fig. 8).The differentiation between layers was done by recognising Fig. 10.Reconstruction of deposited volumes with indication of the flood dates (dd/mm/yyyy).The six larger modelled events were assigned to a single or multiple stratigraphic units (numbers preceding the flood date) whose numbers are indicated in Fig. 9.The last three events were assigned to the surface gravel body.all kind of discontinuity elements such as mud cracks, root marks, changes in the sedimentary structure, organic matter, non-natural materials, etc.These characteristics were also used for establishing a correlation between the two trenches.Texture and charcoal content were also used to corroborate this correlation.The upper one meter of BG-2 is composed by gravels with cross-bedding at the base and massive nonstructured gravels at the top, the latter being the frontal lee face of the lateral gravel bar described previously.The relationship of these gravels with the fine-grain deposits of BG-1 is not obvious and either the gravels are the proximal facies of the floods that deposited units 11 to 15 in BG-1, or they correspond to a later large magnitude flood whose coarse sediments are prograding over the upper five units in BG-1.However, the lack of stratigraphic breaks in the gravel unit prevents a detailed correlation with the last five events described in BG1, and for practical purposes the gravel volume was considered as a sum of the last five events described in BG1.
The volume of each layer was computed following the two methodologies presented in the Methods section, and the results are shown in Table 2.As shown in Fig. 10,  Table 2. Flood unit volumes, calculated by two techniques: (i) "wedge" approach: every layer volume was calculated as if each flood unit had a pyramidal shape (such as a wedge); (ii) proportional approach: by subtracting from the actual deposits the average accumulated layer depth, considering the thickness difference of each layer in each trench and approaching it as a pyramidal shape.The "surface gravel body" layer represents the sedimentation produced by multiple events without clear stratigraphic contacts.

Flood
Volume (i) Volume (ii) unit (m 3 ) (m 3 ) the depositional sequence shows the predominance of seven flood events, which account for the 86 % of the total deposits of the "main check dam" reservoir for the time period 1990-2009.For the sediment regime reconstruction, a date of occurrence was assigned to the flood units described in the stratigraphy considering (1) the relative stratigraphic order of layers since dam construction in 1990, (2) largest rainy events were more probable to produce the largest sediment yields, and moreover (3) sedimentary units containing 1-2 cm charcoal debris lamina were deposited during rainfalls following wildfire events.Two major wildfires have affected the "main check dam" catchment since early 1990s, dating summer 1994 and 2000, and the first floods following the fires took place in December 1995 and October 2000.Lamina with high content of charcoal debris was detected in layer 3 and 9, with decreasing concentration on the overlying beds.Hence, flood unit 3 was related to the December 1995 flood event and flood units 1 and 2 to the previous two floods (December 1992 andApril 1994, respectively).As a partial confirmation of this statement, flood unit 3 is one of the thickest layers, and following the hydrological model results, the December 1995 flood event had the second highest peak discharge of the simulated series (1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009).Flood units 4 to 7 were assigned to four consecutive minor flood events Fig. 11.Sediment temporal validation using the first calibration sediment parameter set: 13 events (out of 38 modelled flood events) were associated with the 15 detected flood units.The surface gravel body corresponds to the last three modelled events.The model was calibrated using the total sedimentation volume.
(January 1996, January 1997, January 1998 and July 1999, respectively).Flood units 8, 9 and 10 were all related to the October 2000 flood event; in fact, this flood event presented three peaks (Fig. 7 -left), and for this reason, deposited three flood units.Flood units 11 to 15 were related to the following 5 flood events (April 2001, May 2002, September 2003, November 2006and April 2007).The upper layer of undistinguished sediments (the surface gravel body) was assumed to be produced by the last three flood events of the time series (October 2007, October 2008and September 2009).

Sediment sub-model calibration and validation
TETIS sediment sub-model was first calibrated using the total sediment volume accumulated behind the "main check dam".The variability in time of the reservoir trap efficiency was taken into account by the STEP approach.The reservoir was divided into 10 finite volumes and the reservoir routing time step was 1 second; the incoming water and sediment discharges were calculated by means of TETIS model.The input sediment yield was divided into three textural classes (sand, silt and clay), and a representative reservoir settling velocity was assigned to each particle size, according to the values recommended by Julien (2010).The calibration was carried out by trial and error adjusting the values of α (for hillslopes) and β (for stream channel network) coefficients within a range of feasible values, following the modellers' expertise.The objective function used in this calibration process was the total volume error.The best value of the objective function (i.e.0 %) was provided by the parameter set α = 350 and β = 0.05.For temporal validation, the chronological assumptions made before on the dam's alluvial stratigraphy and the sediment volume estimated for each flood unit (i.e. the observed Fig. 12. Sediment temporal validation with the second calibration sediment parameter set: 13 events (out of 38 modelled flood events) were associated to the 15 detected flood units.The surface gravel body was not considered.Notice that model was calibrated using the sum of sedimentation volumes from 1 to 15. volume of each layer, Table 2) were compared with the simulated trapped volumes.As Fig. 11 shows, the results are reasonably acceptable.Nevertheless, it is clear that the model tends to overestimate the observed values, especially for the high magnitude events.Since the model was calibrated using the total volume deposited in the reservoir, the overestimation is compensated by the underestimation of the remaining sedimentation volume; i.e. the upper gravel and sand massive deposit on BG-2, with an error of −66 %.This error is probably due to an incorrect reproduction of the sediment reservoir dynamics when the reservoir filling overcomes a certain level.In this case, the conditions for fine alluvial deposition are not fulfilled anymore, and erosion and mixing processes take place.The STEP model does not take into account these phenomena, since it only considers sediment deposition.
In order to overcome this problem, the model was finally calibrated using the period 1992-2007; in this way, the simulated period does not cover the last years, when the STEP model is supposedly not correctly reproducing the reservoir dynamics.The resulting parameters were α = 268 and β = 0.05.In this case, the validation results show a better agreement between the estimated and the simulated volumes, as shown in Fig. 12.The model performance can be described as satisfactory, since the volume error for the thickest flood units ranges between −50 and 50 %, which, given the high uncertainty involved in the modelled process, can be considered a very positive result.This statement confirms that small check dams can be a very important source of information for model calibration and validation, and shows that palaeoflood techniques can help to improve model performance and to estimate sediment yield both for long and short term.The NSE index, calculated with the observed/simulated flood unit sediment volumes, is 0.872 if observed volumes are estimated by the wedge approach and 0.609 if the proportional approach is used.
In order to ensure the robustness of the model implementation, a spatial validation was also carried out.Due to the lack of other flow gauge stations than the one at the catchment outlet, the model was spatially validated by comparing the observed and simulated filling of the seven remaining check dams ("secondary check dams").The model results, shown in Fig. 13, substantially agree with the observed filling within a reasonable approximation.The only anomalous behaviour was noticed at check dam 7, where the model underestimated the observed sedimentation.This could be due to an error in the model parameterization of the sub-catchment draining to check dam 7, or to a poor characterisation of the deposit volume (because of bad accessibility and dense vegetation covering the reservoir surface).Nevertheless, the model can be considered satisfactorily validated in space.

Discussion of sediment results
The simulated texture of deposited sediments is sandy (between 87 and 100 % of sand), agreeing with the field measurements, as shown in Fig. 14.The mean texture of the "main check dam" catchment (20-60-20 % of sand-silt-clay) does not correspond with the reservoir deposit texture, much sandier, because only the coarser material is trapped into the "main check dam".Nevertheless, the simulated deposit texture of the deposit is very similar to the observed one.This observation suggests that the sediment trap efficiency submodel STEP is working properly.
The model provided an average sediment trap efficiency of 51 %, ranging from 29 to 100 %.The trap efficiency varies depending on the flood magnitude and the reservoir capacity, which changes in time due to reservoir filling.Given the mean annual simulated flow (2.05 Hm 3 ) and the reservoir storage capacity (3000 m 3 ), the Brune curves (Brune, 1953) provide a trap efficiency value ranging between 44 and 68 %, with a median of 57 %, which is reasonably close to the value provided by the model.On the other hand, the Brown equation (Brown, 1943), used for example in Bellin et al. (2011) and in Boix-Fayos et al. (2008), provides a lower value, equal to 33 %, very different from the value obtained by the model.The STEP model was able to reproduce the event-to-event trap efficiency variability, providing a time-variable value which could not be computed with simpler sediment trap efficiency equations such as the Brown or Brune approaches.
The resulting mean specific sediment yield (SSY) at the "main check dam" catchment simulated by TETIS is 0.136 Mg ha −2 yr −1 .On the other hand, from the observed reservoir sediment deposit, the following formula can be employed for calculating the area-specific sediment yield (Mg ha −1 yr −1 ): where γ is the mean bulk density (gr cm −3 ) of the deposited sediment, V is the estimated volume of trapped sediment (m 3 ), A the drainage area (ha −1 ), Y the period over which sediment has accumulated (years), and TE the trap efficiency (%).Using the sediment trap efficiency value provided by the model and the density previously computed, the SSY provided by Eq. ( 5) is 0.137 M ha −1 yr −1 .Therefore, a simple lumped approach has obtained the same result in terms of specific sediment yield at the "main check dam".However, SSY computed in this simple way provides limited information, especially in a Mediterranean catchment, due to the high interannual and interevent variability of runoff and sediment yield.In fact, the model results suggest that annual sediment yield varies between 0.0032 Mg ha −1 for the year 2004 and 1.1682 Mg ha −1 for the year 2000.The above-estimated mean SSY can be considered a low erosion rate for Mediterranean catchments (Boix-Fayos et al., 2005;González-Hidalgo et al., 2007;Romero-Díaz et al., 2007;Bellin et al., 2011;and Sougnez et al., 2011).In partic-ular, Boix-Fayos et al. (2005) found that 1 Mg ha −1 yr −1 was one of the lowest erosion rates recorded in the SE of Spain at catchment scale.Nevertheless, these studies were carried out in more erosion-prone catchments.In fact, low erosion rates in shrubland catchments with limestone geology, like Rambla del Poyo, were also observed by other authors (e.g.Kosmas et al., 1997).The main reasons for this significant difference with other Spanish areas are the land cover and the lithological origin of the soil.No degraded areas such as marl gullies or badlands are present in the Rambla del Poyo catchment, mainly because the dominant rock is the limestone with relevant outcrop areas, and the vegetation cover is rather homogeneous and denser than other Mediterranean catchment studied in the previously cited papers (mostly located in more arid areas such as the SE of Spain).As stated by Cerdà (1997), soils originated from limestone have high infiltration rates especially during the dry season, reducing the direct flow on the hillslopes and thus decreasing soil erosion.The homogeneous shrubland cover has also a positive effect on land degradation (Cerdà et al., 1998b), although it increases the risk of fire.This dynamic is typical of many Mediterranean catchments, which suffered strong land abandonment during the 1960s, inducing accelerated land degradation and the development of a shrubland cover, as is the case of the Rambla del Poyo catchment.This behaviour was mentioned in various studies, such as Cerdà et al. (1998a), Rey-Benayas et al. (2007) and Baeza et al. (2007).
The greatest flood event, in terms of peak flow and sediment yield, was the October 2000 flood, which accounted for the 40 % of the total deposited volume and the 43 % of the total sediment yield of the "main check dam" sub-catchment.The SSY for this event was 1.16 Mg ha −1 , a high sediment yield value for shrubland catchments, and the trap efficiency was 35 %.The largest four events accounted for the 80 % of the total sediment yield, and the largest eight for the 90 %.This phenomenon, which was noticed in many ephemeral streams (e.g.Gallart et al., 2005), is due to the extremely variable rainfall regime and the well-known highly non-linear relationship between water discharge and sediment yield.

Conclusions
Deposits stored in check dams are direct evidence for sediment produced at the upstream catchment since dam construction, for both the short (event scale) and long term.Mean annual sediment yield can be calculated from the total volume of sediment retained behind check dams, as done by McManus and Duck (1985), Van den Wall Blake (1986), Neil and Mazari (1993), Foster and Walling (1994), White et al. (1996), Romero-Díaz et al. (2007), Boix-Fayos et al. (2008), Sougnez et al. (2011) and Bellin et al. (2011).Other authors also used this total volume to calibrate or validate mathematical models (e.g.De Vente et al., 2005, 2008;Alatorre et al., 2010).In this paper, the model Hydrol.Earth Syst.Sci., 17, 3339-3354, 2013 www.hydrol-earth-syst-sci.net/17/3339/2013/ implementation considered not only the total volume retained in the "main check dam" for calibration and spatial validation, but also volumes associated with individual flood layers for temporal validation.Detailed alluvial stratigraphy was analysed in two parallel trenches across the dam infill.At least 15 flood layers associated to single flood or flood pulses were identify based on evidences of aerial exposure of sediment contacts (e.g.mudcracks and rootmarks).These palaeoflood sedimentary units were traced along the trenches and correlated between the two trenches on the basis of charcoal debris and dark mud content.A detailed differential GPS survey together with the thickness and geometry of individual flood layers provided a total estimated accumulation volume that were deposited over the time period 1990-2009.The TETIS sediment sub-model was calibrated using the total sediment volume accumulated in the "main check dam" infill.The variability in time of the reservoir trap efficiency was taken into account by coupling the TETIS and the STEP models.This approach represents a very innovative technique for implementing distributed mechanistic models in sediment ungauged catchments, which is the most frequent situation in practical studies.The advance in the scientific and technical knowledge provided by this paper is represented by the possibility of implementing sediment models at ungauged catchments by reproducing the proposed strategy.
The simulated results show good agreement with the estimated sediment volumes retained behind the check dams, both for the short and long term.The model provides a specific sediment yield of 0.136 Mg ha −1 yr −1 for a 12.9 km 2 sub-catchment of Rambla del Poyo, which is lower than other sediment yield rates recorded or estimated in the east of Spain, probably due to extensive limestone bedrock outcrops, thin soils and dense vegetation cover (shrublands with a little portion of pine forest).The model reproduces the hydrological ephemeral behaviour of the stream, and the Mediterranean character of the sediment yield, mainly associated to flow pulses during a limited number of storm events.Almost 90 % of total deposited volume behind the "main check dam" is due to only 8 events in 20 yr.The greatest flood event (October 2000) accounted for the 40 % of the total deposited volume and the 43 % of the total sediment yield, following model results.

Fig. 2 .
Fig. 2. Location of the catchment and check dams.The grey area is the check dam subcatchment used for sediment sub-model calibration and temporal validation.

Fig. 7 .
Fig. 7. Hydrological calibration (left) and one of the temporal validations (right) at 5 min time resolution.

Fig. 8 .
Fig. 8. Upper left: geomorphological sketch with location of the dam structure and trenches BG-1 and BG-2.Bottom left: longitudinal profile across the dam infill and distribution of the two main depositional bodies, i.e. gravel bar and fine deposits (sand and silt).Right: synthetic stratigraphical profiles from trenches BG-1 and BG-2.

Fig. 13 .
Fig. 13.Sediment sub-model spatial validation at seven check dams ("secondary check dams"), comparing the observed and simulated degree of reservoir filling.

Fig. 14 .
Fig. 14.Observed vs. simulated sediment texture of all layers found inside the analysed deposit.