Articles | Volume 19, issue 3
Research article
08 Mar 2023
Research article |  | 08 Mar 2023

Fluvio-deltaic record of increased sediment transport during the Middle Eocene Climatic Optimum (MECO), Southern Pyrenees, Spain

Sabí Peris Cabré, Luis Valero, Jorge E. Spangenberg, Andreu Vinyoles, Jean Verité, Thierry Adatte, Maxime Tremblin, Stephen Watkins, Nikhil Sharma, Miguel Garcés, Cai Puigdefàbregas, and Sébastien Castelltort

The early Cenozoic marine sedimentary record is punctuated by several brief episodes (<200 kyr) of abrupt global warming, called hyperthermals, that have disturbed ocean life and water physicochemistry. Moreover, recent studies of fluvial–deltaic systems, for instance at the Palaeocene–Eocene Thermal Maximum, revealed that these hyperthermals also impacted the hydrologic cycle, triggering an increase in erosion and sediment transport at the Earth's surface. Contrary to the early Cenozoic hyperthermals, the Middle Eocene Climatic Optimum (MECO), lasting from 40.5 to 40.0 Ma, constitutes an event of gradual warming that left a highly variable carbon isotope signature and for which little data exist about its impact on Earth surface systems. In the South Pyrenean foreland basin (SPFB), an episode of prominent deltaic progradation (Belsué–Atarés and Escanilla formations) in the middle Bartonian has been usually associated with increased Pyrenean tectonic activity, but recent magnetostratigraphic data suggest a possible coincidence between the progradation and the MECO warming period. To test this hypothesis, we measured the stable-isotope composition of carbonates (δ13Ccarb and δ18Ocarb) and organic matter (δ13Corg) of 257 samples in two sections of SPFB fluvial–deltaic successions covering the different phases of the MECO and already dated with magnetostratigraphy. We find a negative shift in δ18Ocarb and an unclear signal in δ13Ccarb around the transition from magnetic chron C18r to chron C17r (middle Bartonian). These results allow, by correlation with reference sections in the Atlantic and Tethys, the MECO to be identified and its coincident relationship with the Belsué–Atarès fluvial–deltaic progradation to be documented. Despite its long duration and a more gradual temperature rise, the MECO in the South Pyrenean foreland basin may have led, like lower Cenozoic hyperthermals, to an increase in erosion and sediment transport that is manifested in the sedimentary record. The new data support the hypothesis of a more important hydrological response to the MECO than previously thought in mid-latitude environments, including those around the Tethys.

1 Introduction

The Middle Eocene Climatic Optimum (MECO) is a global warming event that occurred during the Bartonian (ca. 40.5–40.0 Ma) and briefly reversed the longer-term cooling trend of the middle to upper Eocene (Fig. 1; Arimoto et al., 2020; Bijl et al., 2010; Bohaty and Zachos, 2003; Bohaty et al., 2009; Bosboom et al., 2014; Galazzo et al., 2014; Henehan et al., 2020; Edgar et al., 2010, 2020; Giorgioni et al., 2019; Jovane et al., 2007; Mulch et al., 2015; Sluijs et al., 2013; Sppoforth et al., 2010; Pälike et al., 2012; van der Boon et al., 2021). Marine bulk and benthic oxygen isotope compositions (δ18O values) show a negative excursion of −1.5 ‰ over the event, which was interpreted as a gradual global warming of 3 to 6 C (Bohaty et al., 2009). In contrast, the evolution of carbon isotope composition (δ13C values), unlike earlier hyperthermals of the Cenozoic (e.g. Palaeocene–Eocene Thermal Maximum – PETM; Eocene Thermal Maximum – ETM 2), differs from site to site, showing opposite patterns between hemispheres and displaying a carbon isotope excursion (CIE) in some but not all marine records (Bohaty et al., 2009; Henehan et al., 2020; Westerhold and Röhl, 2013). This CIE suggests a rise in atmospheric partial pressure of carbon dioxide (pCO2) during the warming peak (Henehan et al., 2020; Bijl et al., 2010), and numerous potential CO2 sources have been proposed – among them, a prolonged pulse of metamorphic decarbonisation possibly linked with the Himalayan collision at that time (Bijl et al., 2010; Bohaty et al., 2009; Bouilhol et al., 2013; Sternai et al., 2020), an increase of volcanism (van der Boon et al., 2021), or lower continental weatherability (van der Ploeg et al., 2018). However, the pCO2 record remains ambiguous and difficult to link in a straightforward way to a rapid injection of exogenous carbon during the MECO (e.g. Henehan et al., 2020). In addition, regardless of the CO2 sources involved, the MECO coincides with a very long (2.4 Myr) eccentricity cycle, suggesting a possible orbital trigger (Westerhold and Röhl, 2013). Therefore, considering the unresolved MECO driving mechanism(s) and how the Earth system responded to this carbon cycle perturbation, the MECO poses a significant challenge to understanding carbon cycle variations on timescales of several hundreds of thousands of years (Sluijs et al., 2013; Henehan et al., 2020; Sternai et al., 2020). Addressing this challenge requires extensive documentation of the MECO in a range of environments and geodynamic contexts, as well as documentation of its effect on Earth surface dynamics.

Figure 1(a) Cenozoic δ18O values compilation from the Pacific Ocean, compiled in Miller et al. (2020). The continuous blue bar represents permanent ice sheet presence, and the discontinuous blue bar represents ephemeral ice sheet, modified from Zachos et al. (2001). Main climate events: EECO – Early Eocene Climatic Optimum; PETM – Palaeocene–Eocene Thermal Maximum; MECO – Middle Eocene Climatic Optimum; E–O – Eocene–Oligocene transition. (b) Carbonate δ18O values of site 702 from Bohaty et al. (2009). General climatic context of the Middle Eocene Climatic Optimum. The MECO “event”, in yellow from ca. 40.5 to 40.0 Ma in the inset, is considered the last hyperthermal of the Eocene, immediately preceding the shift to genuine Antarctic glaciation and the ice-house world of the Oligocene.

Current data converge towards the view that, during the MECO, surface and deep oceanic waters experienced a gradual and uniform warming between 3 to 6 C for all latitudes (Arimoto et al., 2020; Bijl et al., 2010; Bohaty et al., 2009; Rivero-Cuesta et al., 2019). Moreover, deep-sea carbonates were affected either by significantly reduced carbonate accumulation rates or even by dissolution, suggesting broad acidification of sea-bottom waters, involving an estimated ca. 1 km shoaling of the carbonate compensation depth (CCD; Henehan et al., 2020; Cornaggia et al. 2020; Pälike et al., 2012; Arimoto et al., 2020). However, while the temperature increase in the oceans has been inferred in multiple sites, the MECO environmental perturbation affected differently the fauna communities (Arimoto et al., 2020; Cramwinckel et al., 2019). In some locations, the warmer conditions reduced nutrient availability, decreasing the benthic productivity (Arimoto et al., 2020; Bijl et al., 2010, Galazzo et al., 2014; Moebius et al., 2015; Cramwinckel et al., 2019). In contrast, the Southern Ocean (Moebius et al., 2014) or the Neo-Tethys Ocean (Galazzo et al., 2013) recorded increased productivity during the MECO.

In contrast to the oceanic realm, the expression of the MECO in non-marine records remains scarce and variable. Mulch et al. (2015) suggested a boost of precipitation in the North American plateau derived from low δ18O values, while Bosboom et al. (2014) documented a shift towards arid conditions in the Tarim basin with a reduction in fern palynomorphs. Such a drying trend in central Asia is opposite to the Neo-Tethys Ocean dynamic, where a greater burial of organic matter (OM) immediately following the MECO may have been caused by increased nutrient runoff due to an enhanced hydrological cycle during the warm period (Galazzo et al., 2014; Giorgioni et al., 2019; Spofforth et al., 2010). These studies raise the question of the response of weathering, erosion, and sediment transport in terrestrial systems to global warming, as has also been posed for other hyperthermals recently (e.g. Chen et al., 2018; Foreman et al., 2012; Foreman and Straub, 2017; Honegger et al., 2020). This issue highlights the need for further documentation of the clastic sedimentary successions that temporally cover single and long-term climate crises (i.e. Early Eocene Climatic Optimum, Palaeocene–Eocene Thermal Maximum, etc.; Fig. 1).

In this work, we aim to understand the effects of the MECO on surface systems by exploring the interface between ocean and continent. The shallow marine settings, very sensitive to sea level changes and sediment supply, potentially provide a unique perspective of the hydrological response to climate change in the continental domain and the geochemical and isotopic evolution in the marine domain. We focus on two separate deltaic successions in the southern (Belsué locality, B) and northern (Yebra de Basa locality, YB) margins of the Jaca basin in the South Pyrenean foreland basin (SPFB; Fig. 1). The successions are characterised by excellent exposure and have already been dated by high-resolution magnetostratigraphy. Both sections reveal progradation of deltaic and fluvial systems coeval with the magnetic reversal occurring at chrons C18r and C18n.2n, near or at the zenith of MECO warmth (Edgar et al., 2010, 2020; Garcés et al., 2014; Vinyoles et al., 2021). We generated new high-resolution profiles of δ13Ccarb, δ18Ocarb, major and trace elements, clay mineralogy, and Rock-Eval parameters across the chron C18r–C18n.2n reversal to identify geochemical changes associated with the MECO onset and its recovery. We tested the possible causative links between progradation and the MECO perturbation. Finally, we discuss the sedimentary evolution of both sections to understand landscape response during the MECO and to explore the significance of its identification in the SPFB and the impact of climate shifts in source-to-sink systems as recorded at the continent–ocean interface.

2 Geological setting

The Pyrenees are a nearly E–W-trending mountain belt formed by the collision of the Iberian and European plates from the Late Cretaceous (Santonian) to the Early Miocene (Muñoz, 1992; Roure et al., 1989; Teixell, 1998; Vergés et al., 2002). The south Pyrenean zone is composed of an imbricate system of synorogenic thrusted cover sheets propagating southwards, detached above the Triassic evaporites (Labaume and Teixell, 2018; Lagabrielle et al., 2010; Mochales et al., 2012; Pueyo et al., 2002; Teixell et al., 2016, 2018). Among them, the emplacement of the South-Central Unit (SCU) by the early Eocene resulted in the partitioning of the South Pyrenean basin (SPFB) and the development of an E–W elongated deep basin draining west towards the Atlantic Ocean (Mochales et al., 2012; Muñoz et al., 2018; Puigdefàbregas, 1975; Puigdefàbregas and Souquet, 1986; Seguret, 1972). Due to the westward propagation of deformation during the middle Eocene and the differential velocity of the thrust sheets, oblique thrusted anticlines developed at the southwestern termination of the SCU (Muñoz et al., 2013). These thrusts caused the fragmentation and piggy-back transport of a wider foreland region, which included, from east to west, the Tremp–Graus basin, the Ainsa basin, and eventually the Jaca basins (Muñoz et al., 2018; Fig. 2).

Figure 2Geological and stratigraphic context of the study area (a) Synthetic geologic map of the Pyrenees and location of the study area in the Jaca basin (red square). Black lines represent the main tectonic structures. Modified from Teixell (1996) and Bosch et al. (2016). (b) Detailed geologic map of the Jaca basin and the study area (Belsué and Yebra de Basa). Modified from Puigdefàbregas (1975) and Boya (2018). (c) Chronostratigraphy of the Ainsa and Jaca basins, showing the westward progradation of all depositional systems. The names of the main lithostratigraphic units are represented in black. The studied sections (in white) are included in the work of Vinyoles et al. (2021), while the location of the geochemical analyses carried out in this paper is highlighted in red. The figure is modified from Vinyoles et al. (2021) and the age model GPTS 2012 from Gradstein et al. (2012).

The Tremp–Jaca basin (TJB) preserves an outstandingly exposed, complete source-to-sink system from during MECO times. In the middle Eocene, the alluvial–fluvial system of Sis–Escanilla flowed down the Tremp–Graus and Ainsa basins, draining the eroded sediments from the uplifting northeastern Pyrenees (Beamud et al., 2003; Roigé et al., 2016; Coll et al., 2020; Puigdefàbregas, 1975; Puigdefàbregas and Souquet, 1986). Sediments were transported westward into the Jaca basin, forming a mixed delta-carbonate ramp with tidal influence (Puigdefàbregas, 1975; Castelltort et al., 2003). Two main deltaic systems developed during the lower Bartonian in the southern (Belsué–Atarés) and northern (Sabiñánigo) margin of the Jaca basin that primarily fed the distal Hecho turbidites in the western sector (Mutti, 1977; Remacha and Fernández, 2003). Shallow marine environments are mainly dominated by marly facies, which, thanks to their high carbonate content and relatively deep depositional environment, are suitable for geochemical studies (Wendler, 2013). High-resolution magnetostratigraphy in the Pyrenean region correlates with different sections along the entire source-to-sink system (Vinyoles et al., 2021). We selected two lower-Bartonian sections, Belsué (BS) and Yebra de Basa (YB), because they present excellent exposition and are provided with magnetostratigraphic dating, aiming to unravel the geochemical history of the two main deltaic systems coeval to the MECO.

Figure 3Detailed geological maps from Belsué (a) and Yebra de Basa (b), modified respectively from Puigdefàbregas (1975) and Boya (2018). Red lines show the location of the studied sites, and the darker grey lines in Belsué represent the depositional sequences defined by Millán et al. (1994). Small stratigraphic logs on the bottom right represent large-scale synthetic logs with the formation name, which are modified from Vinyoles et al. (2021) for Yebra de Basa and from Castelltort et al. (2003) for Belsué.

The 130 m thick BS section is located within the External Sierras, east of the “Pico del Águila” anticline (42.30 N, 0.37 W; Fig. 3). This section has been extensively studied as a perfect case study for the tectonic-sedimentation relationship (Fig. 2; Puigdefàbregas, 1975; Puigdefàbregas and Souquet, 1986; Millán et al., 1994, 2000; Castelltort et al., 2003; Huyghe et al., 2012; Garcés et al., 2014). The lower boundary corresponds to an encrusted and ferruginous surface on top of a shallow marine bioclastic limestone (Guara Fm.), locally overlain by sandy marls rich in glauconite (Millán et al., 1994). It has been interpreted as a drowning unconformity of the Guara carbonate platform (Puigdefàbregas, 1975; Silva-Casal et al., 2019), close to the Lutetian–Bartonian boundary (Rodriguez-Pintó et al., 2013). This major unconformity led to the syntectonic deposition of the Arguís marls and the Belsué sandstones, while the Gabardiella and Pico del Águila anticlines were growing. Different authors studied the influence on the stratigraphy of local tectonic movement in the Belsué-Arguís region (Lafont, 1994; Castelltort et al., 2003), concluding that local tectonics modify the stacking pattern and position of its genetic units along with different depositional environments. The entire section covers the lower-Bartonian interval (Garcés et al., 2014) up to a maximum flooding surface MFS-2 by Muñoz et al. (1994) (Fig. 3), which corresponds to the deepest palaeobathymetry in the BS area (ca. 150 m, Sztràkos and Castelltort, 2001).

The 800 m YB section is located between the Basa anticline and the Santa Orosia syncline (42.49 N, 0.28 W; Fig. 3) and comprises one of the best outcropping sections for the Sabiñánigo sandstone (Puigdefàbregas, 1975; Lafont, 1994; Boya, 2018). It is composed of three subsections covering the lower-Bartonian interval, between the upper part of the chron C18r and the chron C18n.1r (Vinyoles et al., 2021). The first subsection is located west (Yb-C), within the Larrés marls, whose top is correlated with the base of the Vinyoles et al. (2021) magnetostratigraphic profile. From here, the next two subsections were built following the Vinyoles et al. (2021) profile, which comprises the two deltaic levels of the Sabiñánigo sandstone and most of the overlaying Pamplona marls.

3 Materials and methods

3.1 Field and sampling

We performed a complete stratigraphic study and sampling of the lower-Bartonian sections from BS and YB. The stratigraphic thickness of these sections was measured using the Jacobs staff in the field and geometric calculations when direct measurements were impossible. New mapping based on fieldwork and orthophotos was performed to correlate the different subsections. The Belsué section in the southern margin is composed of two subsections (Fig. 3); the lower one (Belsué east section, BS-E) was sampled at a resolution of 0.5–1.0 m and the upper one (Belsué west section, BS-O) every 3–6 m. In the northern margin, the higher average sedimentation rate (SR) in Yebra de Basa (>80 cm kyr−1; Vinyoles et al., 2021) motivated a high-resolution sampling of 1–3 m at the middle part of the section (Yebra de Basa high-resolution section, YB-HR) and every 9–15 m in the other intervals (lower Yebra de Basa section, YB-C; upper Yebra de Basa section, YB-sup).

A total of 101 samples in BS and 157 samples in YB were collected. Each was composed of ca. 200 g of fine-grained and fresh rock below the weathering depth to avoid alteration and grain size bias (e.g. Lupker et al., 2011). The samples were mostly marls, corresponding to rocks rich in carbonate, OM, and clays. These samples were prepared for mineralogical and geochemical analyses in the laboratories of the University of Geneva and the University of Lausanne. The sample surface was cleaned with deionised water, and the weathered material was removed and then dried at 45 C for 2–3 d. The dried samples were crushed with a hydraulic press and powdered using an agate mill. The exposure conditions were usually ideal for sampling in both sections. However, there were difficulties in four intervals, resulting in gaps in the data. The problems were due to a dominance in sandy facies at the outcrop, corresponding to moments of maximum deltaic progradation or to poor exposure because of the fine-grained nature of the marls (e.g. Quaternary cover).

3.2 Clay mineralogy

The clay mineral assemblages of 24 representative samples per section were determined by X-ray diffractometry. The used system was a Thermo Scientific ARL X-TRA diffractometer at the Institute of Earth Science of the University of Lausanne (ISTE-UNIL), following the methods described by Klug and Alexander (1974), Kübler (1983, 1987), and Adatte et al. (1996). Ground chips were mixed with deionised water (pH 7–8) and agitated. The carbonate fraction was removed by treatment with 10 % HCl at room temperature and then for 20 min or more until all carbonate was dissolved. The insoluble residue was disaggregated (ultrasonication, 3 min), washed, and centrifuged (eight times) until a neutral suspension was obtained (pH 7–8). Different grain size fractions (<2 to 16 µm) were separated by the time-settling method based on Stokes' law. The selected fraction was then pipetted onto a glass plate and air-dried at room temperature. X-ray diffraction (XRD) analyses of oriented clay samples were made after air drying at room temperature under ethylene-glycol-solvated conditions. The intensities of XRD peaks (2θ; Moore and Reynolds, 1997) characteristic of each clay mineral (e.g. chlorite, mica, kaolinite) were used for a semi-quantitative estimation of the relative percent of clay minerals present in two size fractions (<2 and 2–16 µm).

3.3 Major- and trace-element compositions

Major- and trace-element concentrations from 24 representative samples per section were determined by X-ray fluorescence (XRF) spectrometry using a PANalytical PW2400 spectrometer from ISTE-UNIL. The major elements were analysed on fused glass discs. First, 2.7 to 3 g of sample powder was heated in a crucible oven at 1050 C for 1 night to obtain the loss on ignition (LOI) value. Then, 1.2000±0.0005 g of the calcinated sample was mixed with 6.0000±0.0005 g lithium tetraborate (Li2B4O7) to prepare the fused glass disc using an automated glass-bead-casting machine (Pearl-X'3) at the University of Geneva.

The trace elements were analysed using pressed powder discs prepared at the University of Geneva from 3.000±0.0005 g of wax and 12.000±0.0005 g of non-calcinated sample powder. The mixture was homogenised and pressed (1 t hydraulic press, 20 s). Accuracy of the XRF analyses, assessed by analyses of standard reference materials, was 0.4 wt % for the major elements and 1 to 3 ppm for the trace elements.

3.4 Rock-Eval pyrolysis

The quality and quantity of the organic matter (OM) were determined in 237 bulk rock powders using a Rock-Eval 6 instrument at ISTE-UNIL, following the method described by Behar et al. (2001) and using the IFP 160 000 standard. Aliquots of samples were placed in an oven, heated at 300oC under an inert atmosphere, and then gradually pyrolysed up to 650 C. After the pyrolysis, the samples were transferred into another oven and gradually heated up to 850 C in the presence of air, analysing the CO2 and hydrocarbon (HC) concentration during the entire process. The calculated parameters included total organic carbon content (TOC, wt %), hydrogen index (HI in mg HC g−1 TOC), oxygen index (OI in mg CO2 g−1 TOC), and Tmax (C) according to Espitalié et al. (1985) and Behar et al. (2001).

3.5 Carbonate carbon and oxygen stable isotopes

Carbonate carbon and oxygen stable-isotope ratios (δ13Ccarb and δ18Ocarb) of whole-rock powders containing >10 wt % CaCO3 (n=237) were determined at the laboratories of the Institute of Earth Surface Dynamics of the University of Lausanne (IDYST-UNIL) using a Thermo Fisher Scientific GasBench II carbonate preparation device connected to a Delta Plus XL isotope ratio mass spectrometer. The CO2 extraction was done by reaction with phosphoric acid at 70 C. The carbon and oxygen stable-isotope ratios were reported in the delta (δ) notation as the per mil (‰) relative to the Vienna Pee Dee belemnite standard (VPDB), where δ=(Rsample-Rstandard)/Rstandard×1000 and R=13C/12C or 18O/16O. The δ13Ccarb and δ18Ocarb values were standardised relative to the international VPDB scale by calibration of the reference gases and working standards with the international reference materials NBS 18 (carbonatite, δ13C=-5.04 ‰, δ18O=-23.00 ‰) and NBS 19 (limestone, δ13C=+1.95 ‰, δ18O=-2.19 ‰). Analytical uncertainty (1σ), monitored by replicate analyses of the international calcite standard NBS 19 and the laboratory standard Carrara Marble (δ13C=+2.05 ‰, δ18O=-1.7 ‰), was not greater than ±0.05 ‰ for δ13C and ±0.1 ‰ for δ18O.

3.6 Organic carbon stable isotopes

The organic carbon stable-isotope ratios (δ13Corg values in ‰ vs. VPDB) were determined in 155 samples, which were previously decarbonated by treatment with 10 % v/v HCl, thoroughly washed with deionised water, and dried (40 C, 48 h). The δ13Corg measurements were performed at the IDYST-UNIL by elemental analysis and isotope-ratio mass spectrometry using a Carlo Erba 1108 (Fisons Instruments, Milan, Italy) elemental analyser connected to a Delta V Plus isotope-ratio mass spectrometer via a ConFlo III split interface (both of Thermo Fisher Scientific, Bremen, Germany) operated under continuous helium flow (Spangenberg and Herlec, 2006). The calibration and normalisation of the measured 13C to the VPDB scale was performed with international reference materials and UNIL in-house standards (Spangenberg and Herlec, 2006; Spangenberg, 2016). The repeatability and intermediate precision were better than 0.1 ‰ for δ13Corg.

4 Results

4.1 Stratigraphy and sedimentology

Belsué (BS) stratigraphic succession records the interfingering between prodelta (Arguís Fm.) and deltaic sediments (Belsué Fm.). The Arguís Fm. are highly bioturbated marls and silts, often rich in glauconite, with sparse bioclasts (e.g. bivalves) and oxidised OM fragments. Sandstone beds (Belsué Fm.) are interlayered within the marls, forming two major coarsening and thickening upward sequences that consist of medium sandstone beds (5–10 m thick) with sharp erosion bases, parallel stratification, undifferentiated ripples, and glauconite-rich horizons (Fig. 4). The Arguís marls are interpreted as prodelta deposits in a poorly circulated and relatively deep marine environment (Millán et al., 1994). The marls prelude deltaic mouth bars (Belsué Fm.) where the fluvial component predominates, although local effects from storms and tides are observed (Millán et al., 1994; Castelltort et al., 2003). Calculated palaeocurrents show a corrected east–southeast sediment supply source, in agreement with previous studies (Puigdefàbregas, 1975; Lafont, 1994; Millán et al., 1994; Garcés et al., 2014). Both formations are interpreted as a mixed delta-carbonate ramp system prograding westward, spanning from Bartonian to Priabonian (Castelltort et al., 2003).

Figure 4Stratigraphic logs of Belsué – East (BS-E), Belsué – West (BS-O), and Yebra de Basa, with a more complete high-resolution log (Yebra de Basa HR). Red lines represent flooding surfaces (FS), and green lines represent the maximum flooding surfaces (MFS) correlation. Facies interpretation of the sedimentary logs is represented by grey bars, being more proximal the grey bar and more distal the white. Abbreviations used are as follows: Magneto. for magnetostratigraphy, Lithostrat. for lithostratigraphy, and Seq. Strat. for sequence stratigraphy. The poor-exposure zones in Belsué are covered with a semi-transparent white rectangle with a black cross.


On the northern margin of the basin, the Yebra de Basa (YB) section starts with laminated blue marls (Larrés Fm.) that are interlayered by sparse siltstone beds and two dark levels rich in OM (56–60 m in Yebra-HR). The Larrés Fm. transitions to the Sabiñánigo sandstone (SS), composed of two thickening and coarsening upward sequences. The sandstone beds present planar and, through cross-stratification, erosion scours, sigmoidal beds, and flasser-wavy stratification. The upper boundary (ca. 205 m in Yebra HR) is marked by a sharp contrast towards a highly bioturbated and fossiliferous horizon (hard ground), leading to the deposition of laminated grey-blue marls (Pamplona Fm.), less interlayered with siltstones beds but richer in fossiliferous horizons (Turritella sp. mainly). The system is interpreted as a deltaic siliciclastic shelf prograding W–SW (Roigé, 2018), defined as a fluvial delta with local tidal rework. The deltaic system ended abruptly with a major flooding that formed a hard ground and led to the Pamplona marls deposition (Lafont, 1994; Puigdefàbregas, 1975).

As mentioned in Sect. 3.1, although the exposure conditions were excellent for sampling, we had four intervals without continuous sampling, leading to local data gaps. Three of these data gaps were due to the dominance of sandy facies. In YB, the sandy intervals correspond to the Sabiñánigo sandstone deltaic bodies located at approximately 100–120 and 180–200 m (YB-HR section; Fig. 4). In BS, the Belsué sandstone interval is placed between 55 and 60 m (Belsué – E section; Fig. 4). The fourth data gap located at 85–115 m in Belsué – E is the result of a lack of sufficient exposure within the marls and the presence of a coarse-grained sandy interval (Fig. 4).

Facies associations were described by combining the observations in the field and information from previous studies in the Jaca basin (Millán et al., 1994, 2000; Castelltort et al., 2003; Lafont, 1994; Puigdefàbregas, 1975; Boya, 2018). Using the vertical variations of facies in the studied sections, we defined the depositional sequences that record the shoreline progradation and retrogradation (P–R) cycles. Here, we used the smallest correlatable sequences, which are termed parasequences when bounded by the two shallowest facies (flooding surface, FS; Van Wagoner, 1998; Van Wagoner et al., 1990), or genetic units when bounded by the two deepest facies (maximum flooding surface, MFS; Homewood et al., 1992). The sequence stratigraphic interpretation is summarised in Fig. 4, where parasequence thickness from Belsué and Yebra de Basa varies from a few to tens of metres (5–50 m), and its stacking pattern defines two main P–R cycles in both sections.

4.2 Clay mineralogy

In both sections, Belsué and Yebra de Basa, more than 90 % of the total recorded clay mineral assemblage correspond to the sum of chlorite, chlorite–smectite (CS), mica, or illite–smectite (IS) (Fig. 5). This association of clay minerals is characteristic of dominant physical erosion (Adatte et al., 2000). Mica is the most common clay in both sections (40 %–65 %), followed by chlorite (10 %–36 %). In contrast, the percentage of kaolinite is very low (<5 %). In Belsué, both progradations show different chlorite contents, being higher in the upper part. At Yebra de Basa, we observe an increase in mica that coincides with the OM peak. The absence of smectite indicates it has been transformed into CS or IS mixed layers during diagenesis. Its percentage, 20 %–30 % in the studied sections, can be used as a burial estimation (Kübler and Jaboyedoff, 2000). Kaolinite content positively correlates with the deltaic progradation, indicating that kaolinite could be mainly transported.

Figure 5Schematic stratigraphic log of Belsué – E (a) and Yebra de Basa HR (b) with clay mineral assemblages, IS (illite–smectite), CS (chlorite–smectite), mica, chlorite, and kaolinite. Highlighted in pale grey is the OM-rich interval in Yebra de Basa, and in pale yellow are the main sandstone levels. The progradation–retrogradation cycles (P–R) are drawn with grey and white triangles. The dashed lines represent non-sampled intervals.


4.3 Geochemistry

4.3.1 Major and trace elements

The major and trace elements have been normalised to aluminium (Al) to limit the dilution effect caused by different proportions of terrigenous sediment components (van der Weijden, 2002; Fig. 6). At BS, two increasing pulses of detrital major (Si, Fe, K, and Ti) and trace elements (Mn and Sr) are concomitant with both deltaic progradations (Fig. 6). The similar trend of calcium (Ca) suggests an extra-basinal origin, likely from the eroded Mesozoic or Palaeocene carbonate platforms. Only K shows a negative trend compared to the detrital elements, likely related to clay abundance. At YB, the high TOC interval (depicted in red in Fig. 6) coincides with a relative decrease of the major Si, Ca, Ti, and K and the trace Sr, Zr, and Sn, normalised to Al. In contrast, the OM-related element increase, such as the V/Cr ratio, which is related to OM-rich and suboxic–anoxic conditions (van der Weijden et al., 2006), or the Ni/Co ratio, which is related to biogenic production (Tribovillard et al., 2006).

Figure 6Stratigraphic profiles of Belsué – E and Yebra de Basa HR with the normalised major element concentrations (Si, Fe, Mg, Ti, and K), trace element concentrations (Sr and Zn), and trace element ratios (Ni/Co, V/Cr, and U/Th). The progradation–retrogradation cycles (P–R) are drawn with grey and white triangles. They are highlighted: in pale yellow, the sandstone levels, and in pale grey, the OM-rich interval in Yebra de Basa. The dashed lines represent non-sample intervals.


Figure 7Stratigraphic logs of Belsué and Yebra de Basa including the results of total organic carbon (TOC), carbonate stable isotopes (δ18Ocarb and δ13Ccarb), and organic carbon stable isotopes (δ13Corg). The two progradation–retrogradation cycles referred to in the text are drawn next to the stratigraphy; note the scale change. Highlighted in grey is the MECO isotopic signature, and highlighted in pale red is the OM-rich interval coeval to the MECO thermal peak. The coloured broad lines correspond to the three-point moving average, whilst the central part of the Yebra de Basa section has a seven-point moving-average curve due to the high sampling resolution. Magnetostratigraphic logs from Vinyoles et al. (2021) and Garcés et al. (2014). The dashed lines represent non-sample intervals.

4.3.2 Organic matter content, type, and evaluation

The total organic carbon content (TOC) is low in both sections (average < 0.5 wt %; Fig. 7). At BS, the TOC values range from 0.06 % to 0.38 % (average 0.2±0.07 wt %). The lower one shows a decreasing trend of TOC that follows the deltaic progradation, showing that the OM concentrations decrease with increasing clastic material. The upper section also depicts this trend, with the higher and more stable TOC values (∼0.3 wt %) within marly prodelta deposits. Conversely, YB TOC values range from 0.14 % to 1.3 % (average 0.3±0.13 wt %). The most prominent feature is a dark-grey marl interval with a TOC spike of >1 wt %, which is associated with negative carbonate carbon and oxygen isotope excursions (Fig. 7). Apart from this significant excursion, most TOC values are close to 0.3 wt % and demonstrate no other significant variation along the section. In the high-resolution part of the section, there are small oscillations varying up to ±0.1 wt %.

Figure 8(a) The scatter plot of hydrogen index (HI) vs. Tmax values allows one to discriminate between the different kerogen types, (B), whereas the scatter plot of hydrogen (HI) vs. oxygen index (OI) allows one to assess the OM source and thermal maturity. Here we display the kerogen types (a) and the OM origin (b) from Yebra de Basa (orange) and Belsué (green). Samples with TOC values lower than 0.2 wt % have been excluded. Reference lines for kerogen types and maturity after Espitalié (1986).

The Tmax values and the HI (hydrogen index) / OI (oxygen index) ratios serve to classify the OM in terms of type (origin) and thermal maturity (Espitalié et al., 1985). The Tmax values of the samples range between 422 and 445 C (Fig. 8), with some exceptions reaching almost 460 C. This indicates that the character of the preserved OM is generally immature or within the oil window (Fig. 8). The HI values in YB and BS are generally <100 mg HC g−1 TOC (average 65 mg HC g−1 TOC), corresponding to OM of types III and IV, indicative of a high input of terrestrial plants (Espitalié et al., 1985, Fig. 8). Some samples in Belsué (BS-W) have HI > 150 mg HC g−1 TOC, which is probably related to a slightly different depositional condition between the sections (more distal in the W). The OI values show a similar trend to the HI values, having values below 100 mg CO2 g−1 TOC (average 82 mg CO2 g−1 TOC) but being more dispersed than HI values. In summary, the Rock-Eval parameters indicate a recycled source and/or terrestrial origin of the organic matter in both sections (YB and BS).

4.3.3 Carbonate carbon and oxygen stable isotopes

At BS, the δ13Ccarb values range from −2.6 ‰ to −0.2 ‰ (-1.5±0.4 ‰) and the δ18Ocarb values from −5.7 ‰ to −2.9 ‰ (average -4.4±0.6 ‰; Fig. 7). The δ18Ocarb values show a gradual decrease (of −1.5 ‰) during the first two deltaic progradations, coinciding with a very gradual and small positive trend of the δ13Ccarb values. After the second progradation, δ18Ocarb rapidly returns to more positive values, which are maintained until the top of the section. This steadiness is not followed by the δ13Ccarb values, which record a pronounced positive shift of +1 ‰ between 125 and 150 m height. The δ18Ocarb values show a gradually decreasing trend during the first two deltaic progradations (−1.5 ‰), coinciding with a gradual and small positive trend of the δ13Ccarb values. After the second progradation, δ18Ocarb rapidly returns to more positive values, maintained until the top of the section. This steadiness is not followed by the δ13Ccarb results that record a pronounced positive shift between 125 and 150 m height (of +1 ‰). At YB, the δ18Ocarb values range from −6.4 ‰ to −4.2 ‰ (average -4.9±0.4 ‰), and the δ13Ccarb values range from −1.8 ‰ to 0.6 ‰ (-0.8±0.4 ‰; Fig. 7). Small oscillations (±0.5 ‰) of the δ18Ocarb dominate the lower part of the section and end with a significant negative shift at 285 m. There, the δ13Ccarb values decrease by 0.8 ‰, and the δ18Ocarb values decrease by 1.3 ‰. This level is OM rich (1.0 wt %–1.5 wt % TOC) and also shows a decrease of 2 ‰ in the δ13Corg values (see below). Above the OM-rich interval, the δ13Ccarb and δ18Ocarb values return to pre-event background values, with a shift to higher values towards the base of the Sabiñánigo sandstone, where the δ13Ccarb values reach the maximum value of 0.6 ‰. A gradual decrease of the δ18Ocarb values coincides with the recurrent progradation events evidenced by the Sabiñánigo sandstone. In contrast, the δ13Ccarb follows a stable trend around −1 ‰ until the top of the section. Above the Sabiñánigo sandstone and within the Pamplona marls, the δ18Ocarb values show no variance until the top of the section.

4.3.4 Organic carbon isotopes

At BS section, the δ13Corg values range from −26.3 ‰ to −22.6 ‰ (-24.6±0.5 ‰; Fig. 6), including two groups of samples with different δ13Corg values. The first group is formed by the samples before the first siliciclastic progradation (0–50 m), having relatively low TOC content (0.1 wt %–0.3 wt %) and relatively high δ13Corg values (−25 ‰ to −24 ‰). The second group is formed by samples between the deltaic propagations (65–85 m), which have higher TOC content (0.3 wt %–0.5 wt %) and lower δ13Corg values (-26 ‰). At YB, the δ13Corg values vary between −27.2 ‰ and −23.7 ‰ (-24.9±0.8 ‰; Fig. 6), with the lowest value of −27.2 ‰ having been measured within the OM-rich level at 285 m. The δ13Corg values increase upward till the base of the Sabiñánigo sandstone, where they show a negative spike, coinciding with the positive shift of the δ13Ccarb and δ18Ocarb values. Then the δ13Corg values first return to the pre-event background values and then show a negative excursion of 2 ‰ in the last two samples of the Sabiñánigo sandstone.

Figure 9Oxygen isotope (δ18Ocarb) correlation panel for the studied sections (Belsué and Yebra de Basa) with MECO target curves from Alano (Italy, Tethys Ocean; Spofforth et al., 2010), ODPS 1051 (N Atlantic Ocean; Edgar et al., 2010), ODPS 702 (S Atlantic Ocean; Bohaty et al., 2009), and ODPS 738 (S Indic Ocean; Bohaty et al., 2009). Data from the bulk- and fine-sediments fractions. Highlighted in red is the OM-rich interval (TOC peak) in Yebra de Basa. The two progradation–retrogradation cycles referred to in the text are drawn with grey and white triangles. The data are scaled according to magnetostratigraphic tie points, between C18r–18n.2n and C18n.2n–C18n.1r chrons.

5 Discussion

5.1 MECO isotopic record

In summary, at Belsué, the oxygen isotope record shows a general trend towards more negative values from the base to the middle sandstone units. This trend peaked with a negative δ18Ocarb shift of ∼1 ‰ before the sandstone unit, just in the chron transition C19r–C18n.2n (Fig. 9). The δ18Ocarb values in Yebra de Basa show a small-scale variability consistent with local effects. Then a negative shift of ∼1.2 ‰ occurs close to the deltaic progradation (Fig. 9). The MECO zenith around the magnetic reversal C19r–C18n.2n (Bohaty et al., 2009; Edgar et al., 2010; Henehan et al., 2020) is represented by the progradation of the deltaic facies over the prodelta, i.e. the deltaic facies of Belsué Fm. in Belsué and the Sabiñánigo sandstones in Yebra de Basa (Fig. 9). No isotopic data were obtained for this interval. Nevertheless, the onset of the main thermal event, just before the sandstone occurrence, is preserved as the negative excursion. After the sandstone progradation, the δ18Ocarb values in both sections return closely to those before the excursion (Fig. 9).

The trend observed in the studied sections of the South Pyrenean basin (SPFB) is shared with most high-resolution offshore and nearshore isotopic records of the MECO (Bohaty et al., 2009; Bohaty and Zachos, 2003; Edgar et al., 2010, 2020; Spofforth et al., 2010; Jovane et al., 2007; Giorgioni et al., 2019; Galazzo et al., 2014). They all show the same trend towards more negative δ18Ocarb values, which is intensified during the MECO peak (Fig. 9). The end of the event is marked in both sections by a rapid increase of the δ18Ocarb values by ∼1 ‰, as reported in other sections worldwide (Bohaty et al., 2009; Galazzo et al., 2014; Edgar et al., 2010, 2020; Giorgioni et al., 2019; Spofforth et al., 2010).

Contrarily to the agreement between our new data and most of the available oxygen isotope records, the δ13Ccarb results do not show a clear correlation with global curves (Fig. 7). On one hand, the Belsué section records a positive δ13Ccarb excursion, with a delay with respect to the δ18Ocarb minimum. The Yebra de Basa section shows a prominent positive δ13Ccarb excursion just before the main deltaic progradation (320–340 m from YB). On the other hand, most of the oceanic geochemical records show a small carbon isotope excursion (CIE) at the MECO peak of warming (∼40 Ma; Westerhold and Röhl, 2013; Bohaty et al., 2009; Spofforth et al., 2010), but before and after the δ13Ccarb values are highly variable, showing opposite trends between hemispheres (Henehan et al., 2020; Giorgioni et al., 2019). This CIE, like in other Northern Hemisphere sites (Sppofforth et al., 2010; Giorgioni et al., 2019), is not well represented in our data. Therefore, our δ13Ccarb data seem to confirm that the MECO is not associated with an extensive input of depleted 13C in the environment, as suggested in another study (Henehan et al., 2020). However, an alternative explanation for this discrepancy would be that the sandstone progradation masked the CIE.

Instead of a global origin, however, the δ13Ccarb variations could be caused by local changes in the carbon isotope composition of the dissolved inorganic carbon (DIC), pH, rate and temperature of carbonate precipitation, and its mineralogy (e.g. Swart, 1995). Here, given the proximity of the continent in the shelf environment of the studies section, we cannot rule out that spatial and temporal variations in freshwater input of different components (i.e. riverine, estuarine, and groundwater sources) could have altered the isotopic composition of the DIC and the carbonate δ13C and δ18O values (e.g. Marshall, 1992; Saltzman and Thomas, 2012; Wendler, 2013; Läuchli et al., 2021; Castelltort et al., 2017). This freshwater input could produce carbonate depleted in 13C and 18O (see detailed discussion in Sect. 5.2). The fact that the Belsué and Yebra de Basa isotopic records preserve the MECO excursion in δ18Ocarb suggests that the δ13Ccarb values could also record a global signal. However, this is difficult to appreciate because of the small δ13Ccarb variations in the studied sections and the somewhat variable and peculiar published carbon isotope record during the MECO. Therefore, the correlations were based on the δ18Ocarb records (Fig. 9). In summary, our results record a decoupling between oxygen and carbon isotopes. The δ18Ocarb values seem to follow the global trend, while the variation of the δ13Ccarb values remains ambiguous. The particular carbon isotope record during the MECO, with variable and opposite trends between hemispheres (Henehan et al., 2020) and a brief negative carbon isotope excursion recorded just at some sites (Westerhold and Röhl, 2013; Bohaty et al., 2009; Spofforth et al., 2010), could be an explanation. However, given the location of the studied sections on a continental shelf, it is important to check the possible diagenetic influence or alteration of the primary isotopic signal.

5.2 Primary versus diagenetic signals

The carbonate primary carbon and oxygen isotope compositions may be affected by post-depositional processes, including the neoformation of authigenic and diagenetic phases. Therefore, before the palaeoenvironmental interpretation of δ13Ccarb and δ18Ocarb records from shallow marine environments, it is necessary to determine primary versus diagenetic signal components. This discrimination requires an understanding of the factors controlling the primary marine isotopic composition and an evaluation of potential diagenetic overprints on the original geochemical signatures (e.g. Marshall, 1992; Schrag et al., 1995).

Oxygen isotopes in carbonates are controlled by the temperature of formation, the δ18O value of the carbonate-precipitating fluid (δ18Ow), the mineralogy (e.g. higher δ18O in dolomite vs. calcite), and any environmental parameter (e.g. pH, salinity) affecting the rate of carbonate precipitation (Swart, 1995). The effect of diagenetic alteration is more pronounced in the case of oxygen isotopes than in carbon isotopes due to the high amount of oxygen relative to carbon present in post-depositional fluids and their variable δ18O values (e.g. Marshall, 1992; Schrag et al., 1995; Fio et al., 2010). Carbonate with low δ18O values can be produced by increasing temperature, freshwater input, and meteoric diagenesis, whereas 18O enrichment could indicate either lower temperature or evaporation (e.g. Marshall, 1992; Patterson and Walter, 1994; Schrag et al., 1995). In contrast, carbon isotopes are not thought to be directly influenced by temperature and are generally more resistant to diagenetic processes (Patterson and Walter, 1994; Schrag et al., 1995; Swart, 1995). However, δ13C values are also controlled by kinetic effects, mineralogy, and mainly by the δ13C value from the DIC (Wendler, 2013). The primary diagenetic process that affects the δ13C values of the DIC is the oxidation of the organic matter, which produces CO2 (and DIC species) depleted in 13C (low δ13C values). Therefore, the δ13C values of the DIC and derived carbonates indicate the source of carbon, including the type of degraded or oxidised organic matter (OM) of different types, original seawater carbon, and skeletal and non-skeletal carbonate sources (e.g. Swart, 1995). In proximal depositional environments, however, the δ13C values could be modified by (1) OM source, productivity, and burial rate; (2) extra-basinal carbonate input; (3) water circulation, stratification, and evaporation; and (4) terrestrial runoff and weathering (Saltzman and Thomas, 2012; Läuchli et al., 2021). Considering this, δ13C is usually used as a global correlation tool, since it can register eustatic sea level fluctuations, changes in weathering flux, or significant perturbations in the global carbon cycle (e.g. volcanic CO2 input; Wendler, 2013, and references therein).

Figure 10δ18Ocarbδ13Ccarb scatterplot of all Belsué and Yebra de Basa samples. The size of the symbols is small for samples with TOC < 0.2 wt % and big for samples with TOC > 0.2 wt %. The Pearson correlation coefficients (r) and regression lines are shown.


The degree of diagenetic alteration was assessed through three different approaches. First, we evaluated the relationship between δ13C and δ18O values (Brasier et al., 1996). Statistically, a non-significant correlation (Pearson correlation coefficient r<0.6) indicates that a diagenetic overprint of the primary isotopic signature can be excluded (e.g. Fio et al., 2010). In both sections, no statistically significant correlation (r<0.3) was found between the δ18Ocarb and δ13Ccarb values. This lack of relationship suggests that no or minor diagenetic modifications affected the primary isotopic compositions (Fig. 10). The second approach used to assess the degree of alteration uses clay mineralogy. Kübler and Jaboyedoff (2000) defined four diagenetic zones by comparing illite crystallinity with mineral assemblages and organic-matter type. The Belsué and Yebra de Basa samples have 20 %–30 % smectite within the illite–smectite (IS) mixed layers and are within the third diagenetic zone of Kübler and Jabeyedoff (2000), i.e. shallow diagenesis (ca. 60–80 C). Another diagenetic indicator is the maximum temperature (Tmax) reached during the Rock-Eval pyrolysis (S2), which marks the maturity of the OM. The Tmax values obtained in samples with relatively high OM content (TOC > 0.5 wt %; S2 > 0.2) were <440C (Fig. 8), which corresponds to the beginning of the oil window (ca. 60 C; Espitalié et al., 1985). This maturity level of the organic matter agrees with vitrinite reflectance and Raman measurements in the studied area (Labaume et al., 2016). In summary, the three approaches for assessment of the diagenetic degree – i.e. carbonate δ13C and δ18O values, illite crystallinity, and thermal maturation of the organic matter (Tmax) – suggest that the diagenetic overprint in the studied Belsué and Yebra de Basa rocks is low. The primary isotopic signal is preserved largely in both sections. It can be safely used to study palaeoenvironmental conditions and can be compared to global key isotopic curves during the MECO event.

5.3 The organic matter peak in Yebra de Basa

In the Yebra de Basa (YB) section, an increase in TOC content at the 280–290 m interval (up to 1.5 wt % TOC) is associated with a negative isotope excursion of −1.5 ‰ for δ18Ocarb, −2.0 ‰ for δ13Corg, and −0.8 ‰ for δ13Ccarb values (Figs. 7 and 9). The OM-rich interval occurs 50 m below the main Sabiñánigo sandstone progradation in Yebra de Basa. It is not coincident with the prominent increase of detrital input, marked by an increase in grain size. This boost in OM burial is also observed in the Neo-Tethys region, like in Italy (Spofforth et al., 2010) and the Crimea–Caucasus (Benyamovsky, 2012), which may have played an important role in carbon drawdown and rapid cooling after the MECO event (Bohaty et al., 2009; Henehan et al., 2020).

Several possibilities could explain the presence of an OM-rich interval before a deltaic progradation. First, a significant freshwater input in a restricted basin can lead to water stratification where anoxic conditions are favoured, increasing OM preservation independently of its source. Nevertheless, the slight increase in redox-sensitive elements (V and Mo, Fig. 6) is too limited to support the development of water stratification and the resulting suboxic–anoxic conditions (Tribovillardet al., 2006). Second, the enhanced freshwater input could have increased nutrient availability and marine productivity. We, however, reject this hypothesis because the geochemical data point to the main terrestrial components of the organic matter (low HI–OI) and little to no nutrient availability increases (low Ni concentration; Tribovillard et al., 2006). Therefore, the most probable explanation is that the OM peak could be related to a significant increase in detrital input and terrestrial OM. The presence of several westward dark-marl beds suggests that it was not a unique episode but a series of recurrent events (Boya, 2018). In addition, the terrestrial origin is also supported by the strong correlation (r>0.7) observed between the siliciclastic elements (Al, Ti, Fe) and the TOC or all the OM-related trace elements (V, Mo, Ba, and Th; Tribovillard et al., 2006). Despite this, the isotopic results do not agree with this correlation because pre-Miocene marine OM had lower δ13C than terrestrial OM (Popp et al., 1989). Thus, an alternative explanation for the negative δ13Corg and δ13Ccarb excursion may be an increased input of organic matter released from soils containing bacterial biomass with low δ13Corg values (Fio et al., 2010). This agrees with the Rock-Eval results, which point towards a terrestrial origin of this organic matter (Fig. 8).

As a result, the Sabiñánigo sandstone represents a singular deltaic event embedded in long-lasting prodelta conditions (Vinyoles et al., 2021) in which no evident organic events occur. Therefore, we interpret the occurrence of the OM-rich level just before the Sabiñánigo sandstone as the first indicator of a shift towards a setting with more fluvial conditions, this being the first evidence of the main MECO excursion in the region.

5.4 MECO response in the South Pyrenean foreland basin

The integration of available age constraints (Garcés et al., 2014; Vinyoles et al., 2021) and the new high-resolution isotopic record show that MECO's warming peak (∼40 Ma) is associated with isochronous progradation, which can be followed all along the SPFB source-to-sink system (Fig. 11; Vinyoles et al., 2021). In the Tremp–Graus basin, the Escanilla fluvial system was fed by the Sis–Gurp and Pobla alluvial systems, where a grain size increase is recorded at ca. 40 Ma (Whittaker et al., 2011). Downstream, in the time-equivalent sections in the Ainsa basin, an anomalous amalgamated Olsón sheet stands out from the landscape as a continuous and thick conglomeratic bed, interpreted as a stacking of several braided river channels (Fig. 11; Verité, 2019; Labourdette, 2011; Puigdefàbregas, 1975; Vinyoles et al., 2021). In the deltaic counterparts (Jaca basin), a significant progradation of deltaic deposits on top of slope marls is observed in the studied sections (BS and YB; Lafont, 1994; Puigdefàbregas, 1975; Vinyoles et al., 2021). Finally, in the deeper sink environments of the Jaca and Pamplona basins, the correlation with the turbiditic systems is still debated and needs further research.

Figure 11Correlation panel between Belsué, Yebra de Basa, and Olsón section with the GPTS 2016 (Ogg et al., 2016). The stratigraphic sections are modified from Garcés et al. (2014) and Vinyoles et al. (2021). The oxygen isotopic (δ18Ocarb) record from ODPS 738 and the MECO age constraints defined by yellow and blue bars are modified from Bohaty et al. (2009). The oxygen isotopic record (δ18Ocarb) from Belsué and Yebra de Basa corresponds to our results, and the dashed lines represent non-sample intervals. The sedimentation rates (SRs) from Belsué, Yebra de Basa, and Olsón are the average SR between chron C18n.2n and C17r, as from the data of Vinyoles et al. (2021).

Previous works (Puigdefàbregas, 1975; Lafont, 1994) interpreted these deltaic sequences as eustatic fluctuations of the relative sea level, which can be related to different possibilities, such as thermal expansion or glacioeustasy. Ephemeral ice sheets in Antarctica during the Middle Eocene are likely, and it seems plausible that the progressive shift towards ice-house conditions could have significant implications during the MECO (Edgar et al., 2007; Huyghe et al., 2012; Baatsen et al., 2020). However, considering the temperature increase interpreted during the MECO zenith (+4 to 6 C; Bohaty et al., 2009), we should expect a sea level rise (ice caps melting and thermal expansion) instead of the observed regression and system progradation.

Alternatively, an abrupt sediment supply increase can also explain a progradation of deltaic systems. Several studies observed that the main Palaeogene hyperthermals are often associated with an enhanced flux of terrigenous material, interpreted as a boost of the hydrological cycle and higher seasonality (Schmitz and Pujalte, 2007; Chen et al., 2018; Foreman and Straub, 2017; Pujalte et al., 2015; Schmitz and Pujalte, 2007). Although the MECO is not an abrupt event like other hyperthermals, instead being a more extended period of gradual warming (ca. 500 kyr; Bohaty et al., 2009), we also observe this progradation being focused during the warming peak (ca. 40 Ma). Accordingly, an explanation for the progradation is that the MECO prolonged warming produced an enhanced hydrological cycle that favoured sediment production and transport, thus leading to an increase in sediment supply and favouring the system progradation at the peak of the event. The nature of a greater sediment provision (Qs) should be considered to have originated upstream, for instance, being linked to enhanced sediment remobilisation (e.g. floodplain) or accelerated hillslope processes (Foreman et al., 2012).

Therefore, the coincidence in time of a basin-wide progradation in the SPFB and the MECO might implicate a link between them. Our geochemistry analyses also suggest a terrestrial origin for this OM, which points towards an increase in soil remobilisation, erosion, and transport in continental environments during the MECO event.

5.5 Global implications and correlation

The global impact of the MECO event in continental settings remains poorly documented, with only a few studies in continental environments having been performed around the globe (e.g. Bosboom et al., 2014; Mulch et al., 2015). In the North American plateau, a boost of precipitation during the MECO is derived from lower δ18Ocarb values (Mulch et al., 2015). In contrast, in the Tarim basin (China), a shift towards arid conditions has been interpreted from a reduction in fern palynomorphs (Bosboom et al., 2014). This aridification trend in central Asia differs from the documented Neo-Tethys Ocean dynamic, where marine records show an increase in organic-matter (OM) burial during the MECO peak and part of the post-MECO recovery (Spofforth et al., 2010; Giorgioni et al., 2019 Benyamovskiy, 2012). Increased sediment supply due to enhanced erosion and transport provides a mechanism for the more efficient burial of OM during this and other hyperthermals (Galy et al., 2007). If this enhanced OM burial is global or sufficiently widespread (it is absent in several sections, including Belsué in this study), it could represent an important mechanism to explain the carbonate δ13C increase that is recorded globally during the post-event recovery and the associated rapid return to pre-event conditions, maybe playing an essential role in the drawdown of atmospheric carbon (e.g. Bohaty et al., 2009; Henehan et al., 2020; Sluijs et al., 2013; Edgar et al., 2020; Giorgioni et al., 2019; Spofforth et al., 2010).

Considering the long duration of the MECO event (ca. 500 kyr; Bohaty et al., 2009), some of the most important effects in the ocean occur during its peak phase, e.g. ocean acidification (Bohaty et al., 2009; Henehan et al., 2020; Arimoto et al., 2020) or OM burial (Giorgioni et al., 2019; Spofforth et al., 2010). In the SPFB, the continental progradation also occurred at the end of the event, supported by the sedimentological and geochemical evidence that shows an increase of sediments delivered to the sea, including large amounts of organic matter of terrestrial origin. Hence, our work suggests a link between enhanced hydrological cycles and enhanced OM transport and burial, possibly accounting for the observations of enhanced OM burial around the Neo-Tethys region. This response in sediment delivery rate (OM burial in shallow and restricted basins) has been previously documented for other early Eocene hyperthermals (Chen et al., 2018; Foreman et al., 2012; Pujalte et al., 2015; Foreman and Straub, 2017; Honegger et al., 2020). Hence, despite its important differences in relation to the early Eocene hyperthermals, the MECO shares several attributes with them around the warming peak. In summary, our results point to a more intense hydrological cycle perturbing rainfall patterns in the Pyrenean region during the MECO peak and leading to increased sediment supply, expressed by a major progradation of sedimentary systems and, eventually, an increase in OM burial in the nearby oceanic basins.

6 Conclusions

In the South-Pyrenean foreland basin, an important progradation affected the entire sediment-routing system from fluvial to deltaic environments at the time of the Middle Eocene Climatic Optimum (MECO). Here we present a new high-resolution multiproxy dataset, including stable isotopes, Rock-Eval, XRF, and clay minerals, covering the different MECO phases from two well-dated key sections. The new stable-isotope records from Belsué (BS) and Yebra de Basa (YB) sections show a significant negative shift in the shallow marine sediments around the main warming peak of the MECO event, reported for the first time in the Pyrenean region. In Yebra de Basa, an organic-rich interval of terrestrial origin is found before the main deltaic progradation. It is associated with a negative excursion in oxygen and carbon isotopes. The correlation between the MECO and the basin-wide progradation, as well as the new geochemical results, presents compelling evidence for a climatic driver, suggesting an enhanced hydrological cycle in the Pyrenean region that caused a boost in sediment and carbon export. This is in agreement with previous studies from the Neo-Tethys Ocean that recorded an increase in organic-matter burial during the peak of the MECO and early post-MECO.

Although the duration of the MECO and its isotopic signature differ with respect to early Eocene hyperthermals (e.g. PETM), there are similarities around the warming peak that trigger a comparable response, including ocean acidification, organic-matter burial, or a boost in sediment supply export from land to sea. Nevertheless, further work is needed to understand the role of potential sediment supply increase from the proximal continental environments towards the deeper oceanic basins and, importantly, to quantify sediment and organic export and its relationship with carbon burial and silicate weathering.

Our results support the view that high-accommodation settings in foreland basins are important recorders of palaeoenvironmental signals, even in shallow marine environments. Although certainly noisy, the fact that climate signals are preserved in these settings provides a range of potentially expanded sections that can be an interesting complement to high-resolution but more condensed deep-sea palaeoclimatic records, particularly during high-CO2 globally warm episodes of the Earth's history when the carbonate-rich oceanic records may undergo intervals of non-deposition or dissolution.

Data availability

All the data (stable isotopes, clay minerals, organic matter, major and trace elements) can be found in the Supplement.


The supplement related to this article is available online at:

Author contributions

SPC led fieldwork, sampling, sample preparation, data interpretation, and writing. LV contributed to the fieldwork, data interpretation, discussion, and writing. JES performed stable-isotope analyses, data interpretation, discussion, and writing. TA performed XRD analyses, data interpretation, and discussion. JV and AV contributed to the fieldwork preparation, sampling, and discussion. MT, SW, and NS contributed to the discussion and writing. CP contributed to fieldwork, discussion, and writing. AV and MG helped with magnetostratigraphic data interpretation and discussion. SC supervised the project, funding, interpretation, and writing.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher' note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


The authors would like to acknowledge the Societé de Physique et d'Histoire naturelle de Genève and Equinor (grant to Castelltort) for financing part of the field missions, as well as the Swiss National Science Foundation for funding this project. We would like to thank the Faculty of Science of the Université de Genève for the Excellence Master Fellowship (grant to Peris Cabré), as well as the ELSTE Master programme between Université de Genève and Lausanne. We also acknowledge Marta Roigé and Salvador Boya for their help during field campaigns and the long scientific discussions, and Antoine de Haller from Université de Genève for his help during XRF analyses. Finally, we would like to thank the anonymous reviewer and Eric Barefoot for their critical and insightful comments, as well as the editor Denis-Didier Rousseau for his helpful comments and editorial handling.

Financial support

This research has been supported by the Swiss National Science Foundation (SNSF grant no. 200020_182017: Earth Surface Signalling Systems 2).

Review statement

This paper was edited by Denis-Didier Rousseau and reviewed by Eric Barefoot and one anonymous referee.


Adatte, T., Stinnesbeck, W., and Keller, G.: Lithostratigraphic and mineralogic correlations of near K/T boundary clastic sediments in northeastern Mexico: Implications for origin and nature of deposition, Special Paper 307: The Cretaceous-Tertiary Event and Other Catastrophes in Earth History, Geological Society of America, 211–226,, 1996. 

Adatte, T., Bolle, M. P., Kaenel, E. D., Gawenda, P., Winkler, W., and Von Salis, K.: Climatic evolution from Paleocene to earliest Eocene inferred from clay-minerals: A transect from northern Spain (Zumaya) to southern (Spain, Tunisia) and southeastern Tethys margins (Israel, Negev), inL: Vol. 122, Taylor & Francis, 7–8,, 2000. 

Arimoto, J., Nishi, H., Kuroyanagi, A., Takashima, R., Matsui, H., and Ikehara, M.: Changes in upper ocean hydrography and productivity across the Middle Eocene Climatic Optimum: Local insights and global implications from the Northwest Atlantic, Global Planet. Change, 193, 103258,, 2020. 

Baatsen, M., von der Heydt, A. S., Huber, M., Kliphuis, M. A., Bijl, P. K., Sluijs, A., and Dijkstra, H. A.: The middle to late Eocene greenhouse climate modelled using the CESM 1.0.5, Clim. Past, 16, 2573–2597,, 2020. 

Beamud, E., Garcés, M., Cabrera, L., Muñoz, J. A., and Almar, Y.: A new middle to late Eocene continental chronostratigraphy from NE Spain, Earth Planet. Sc. Lett., 216, 501–514,, 2003. 

Behar, F., Beaumont, V., and Penteado, H. D. B.: Rock-Eval 6 technology: performances and developments, Oil Gas Sci. Technol., 56, 111–134,, 2001. 

Benyamovskiy, V. N.: A high resolution Lutetian–Bartonian planktonic foraminiferal zonation in the Crimean–Caucasus region of the Northeastern Peri-Tethys, Aust. J. Earth Sci., 105, 117–128, 2012. 

Bijl, P. K., Houben, A. J., Schouten, S., Bohaty, S. M., Sluijs, A., Reichart, G. J., Sinninghe Damsté, J. S., and Brinkhuis, H.: Transient Middle Eocene atmospheric CO2 and temperature variations, Science, 330, 819–821,, 2010. 

Bohaty, S. M. and Zachos, J. C.: Significant Southern Ocean warming event in the late middle Eocene, Geology, 31, 1017–1020,, 2003. 

Bohaty, S. M., Zachos, J. C., Florindo, F., and Delaney, M. L.: Coupled greenhouse warming and deep-sea acidification in the middle Eocene, Paleoceanography, 24, PA2207,, 2009. 

Bosboom, R. E., Abels, H. A., Hoorn, C., van den Berg, B. C., Guo, Z., and Dupont-Nivet, G.: Aridification in continental Asia after the middle Eocene climatic optimum (MECO), Earth Planet. Sc. Lett., 389, 34–42,, 2014. 

Bosch, G. V., Teixell, A., Jolivet, M., Labaume, P., Stockli, D., Domenech, M., and Monie, P.: Timing of Eocene–Miocene thrust activity in the Western Axial Zone and Chaînons Béarnais (west-central Pyrenees) revealed by multi-method thermochronology, Comptes Rendus Geoscience, 348, 246–256,, 2016. 

Bouilhol, P., Jagoutz, O., and Hanchar, J. M.: Dating the India–Eurasia collision through arc magmatic records, Earth Planet. Sc. Lett., 366, 163–175,, 2013. 

Boya, S.: El Sistema deltaico de la Arenisca de Sabiñánigo y la continentalización de la cuenca de Jaca, PhD thesis, Universitat Autònoma de Barcelona, Barcelona, 207 pp., ISBN 9788449083006, 2018. 

Brasier, M. D., Shields, G., Kuleshov, V. N., and Zhegallo, E. A.: Integrated chemo-and biostratigraphic calibration of early animal evolution: Neoproterozoic–early Cambrian of southwest Mongolia, Geolog. Mag., 133, 445–485,, 1996. 

Castelltort, S., Guillocheau, F., Robin, C., Rouby, D., Nalpas, T., Lafont, F., and Eschard, R.: Fold control on the stratigraphic record: a quantified sequence stratigraphic study of the Pico del Aguila anticline in the Southwestern Pyrenees (Spain), Basin Res., 15, 527–551,, 2003. 

Castelltort, S., Honegger, L., Adatte, T., Clark, J. D., Puigdefàbregas, C., Spangenberg, J. E., Dykstra, M. L., and Fildani, A.: Detecting eustatic and tectonic signals with carbon isotopes in deep-marine strata, Eocene Ainsa Basin, Spanish Pyrenees, Geology, 45, 707–710,, 2017. 

Chen, C., Guerit, L., Foreman, B. Z., Hassenruck-Gudipati, H. J., Adatte, T., Honegger, L., Perret, M., Sluijs, A., and Castelltort, S.: Estimating regional flood discharge during Palaeocene-Eocene global warming, Sci. Rep., 8, 1–8,, 2018. 

Coll, X., Gómez-Gras, D., Roigé, M., Teixell, A., Boya, S., and Mestres, N.: Heavy-mineral provenance signatures during the infill and uplift of a foreland basin: An example from the Jaca basin (southern Pyrenees, Spain), J. Sediment. Res., 90, 1747–1769,, 2020. 

Cornaggia, F., Bernardini, S., Giorgioni, M., Silva, G. L., Nagy, A. I. M., and Jovane, L. : Abyssal oceanic circulation and acidification during the Middle Eocene Climatic Optimum (MECO), Sci. Rep., 10, 1–9,, 2020. 

Cramwinckel, M. J., Van Der Ploeg, R., Bijl, P. K., Peterse, F., Bohaty, S. M., Röhl, U., and Sluijs, A.: Harmful algae and export production collapse in the equatorial Atlantic during the zenith of Middle Eocene Climatic Optimum warmth, Geology, 47, 247–250,, 2019. 

Edgar, K. M., Wilson, P. A., Sexton, P. F., and Suganuma, Y.: No extreme bipolar glaciation during the main Eocene calcite compensation shift, Nature, 448, 908–911,, 2007. 

Edgar, K. M., Wilson, P. A., Sexton, P. F., Gibbs, S. J., Roberts, A. P., and Norris, R. D.: New biostratigraphic, magnetostratigraphic and isotopic insights into the Middle Eocene Climatic Optimum in low latitudes, Palaeogeogr. Palaeocl. Palaeoecol., 297, 670–682,, 2010. 

Edgar, K. M., Bohaty, S. M., Coxall, H. K., Bown, P. R., Batenburg, S. J., Lear, C. H., and Pearson, P. N.: New composite bio-and isotope stratigraphies spanning the Middle Eocene Climatic Optimum at tropical ODP Site 865 in the Pacific Ocean, J. Micropalaeontol., 39, 117–138,, 2020. 

Espitalié, J.: Use of Tmax as a Maturation Index for Different Types of Organic Matter. Comparison with Vitrinite Reflectance, in: Thermal Modeling in Sedimentary Basins, 44. Editions Technip, edited by: Burrus, J., Publications de l'Institut Français du Pétrole/Institut Français du Pétrole, Paris, 475–496, (last access: 1 March 2023), 1986. 

Espitalié, J., Deroo, G., and Marquis, F.: Rock-Eval pyrolysis and its applications, Revue De L'Institut Français Du Petrole, 40, 563–579,, 1985. 

Fio, K., Spangenberg, J. E., Vlahović, I., Sremac, J., Velić, I., and Mrinjek, E.: Stable isotope and trace element stratigraphy across the Permian–Triassic transition: A redefinition of the boundary in the Velebit Mountain, Croatia, Chem. Geol., 278, 38–57,, 2010. 

Foreman, B. Z., Heller, P. L., and Clementz, M. T.: Fluvial response to abrupt global warming at the Palaeocene/Eocene boundary, Nature, 491, 92–95,, 2012. 

Foreman, B. Z. and Straub, K. M.: Autogenic geomorphic processes determine the resolution and fidelity of terrestrial paleoclimate records, Sci. Adv., 3, e1700683,, 2017. 

Galazzo, F. B., Giusberti, L., Luciani, V., and Thomas, E.: Paleoenvironmental changes during the Middle Eocene Climatic Optimum (MECO) and its aftermath: The benthic foraminiferal record from the Alano section (NE Italy), Palaeogeogr. Palaeocl. Palaeoecol., 378, 22–35,, 2013. 

Galazzo, F. B., Thomas, E., Pagani, M., Warren, C., Luciani, V., and Giusberti, L.: The middle Eocene climatic optimum (MECO): A multiproxy record of paleoceanographic changes in the southeast Atlantic (ODP Site 1263, Walvis Ridge), Paleoceanography, 29, 1143–1161,, 2014. 

Galy, V., France-Lanord, C., Beyssac, O., Faure, P., Kudrass, H., and Palhol, F.: Efficient organic carbon burial in the Bengal fan sustained by the Himalayan erosional system, Nature, 450, 407–410,, 2007. 

Garcés, M., López-Blanco, M., Valero, L., Beamud, E., Pueyo-Morer, E., and Rodríguez-Pinto, A.: Testing orbital forcing in the Eocene deltaic sequences of the South-Pyrenean Foreland Basins, in: Vol. 16, EGU General Assembly 2014, 27 April–2 May 2014, Vienna, Austria, EGU2014-10681-1, 2014. 

Giorgioni, M., Jovane, L., Rego, E.S., Rodelli, D., Frontalini, F., Coccioni, R., Catanzariti, R., and Özcan, E.: Carbon cycle instability and orbital forcing during the Middle Eocene Climatic Optimum, Sci. Rep., 9, 9357,, 2019. 

Gradstein, F. M., Ogg, J. G., Schmitz, M., and Ogg, G.: The Geologic Time Scale 2012, in: Vol. 2, Elsevier, Cambridge University Press, Cambridge, ISBN 978-0-44-459431-1, 2012. 

Henehan, M. J., Edgar, K. M., Foster, G. L., Penman, D. E., Hull, P. M., Greenop, R., and Pearson, P. N.: Revisiting the Middle Eocene Climatic Optimum `Carbon Cycle Conundrum' with new estimates of atmospheric pCO2 from boron isotopes, Paleoceanogr. Paleocl., 35, e2019PA003713,, 2020. 

Homewood, P., Guillocheau, F., Eschard, R., and Cross, T. A.: Corrélations haute résolution et stratigraphie génétique: une démarche intégrée, Bulletin des Centres de Recherches Exploration-Production Elf-Aquitaine, 16, 357–381, 1992. 

Honegger, L., Adatte, T., Spangenberg, J. E., Caves Rugenstein, J. K., Poyatos, M., Puigdefàbregas, C., and Harlaux, M.: Alluvial record of an early Eocene hyperthermal, Castissent Formation, the Pyrenees, Spain, Clim. Past, 16, 227–243,, 2020. 

Huyghe, D., Castelltort, S., Mouthereau, F., Serra-Kiel, J., Filleaudeau, P. Y., Emmanuel, L., Berthier, B., and Renard, M.: Large scale facies change in the middle Eocene South-Pyrenean foreland basin: The role of tectonics and prelude to Cenozoic ice-ages, Sediment. Geol., 253, 25–46,, 2012. 

Jovane, L., Florindo, F., Coccioni, R., Dinarès-Turell, J., Marsili, A., Monechi, S., and Sprovieri, M.: The middle Eocene climatic optimum event in the Contessa Highway section, Umbrian Apennines, Italy, Geol. Soc. Am. Bull., 119, 413–427,, 2007. 

Klug, H. P. and Alexander, L. E.: X-ray diffraction procedures: for polycrystalline and amorphous materials, John Wiley & Sons Inc., p. 992, ISBN 978-0-471-49369-3, 1974. 

Kübler, B.: Dosage quantitatif des minéraux majeurs des roches sédimentaires par diffraction X, in: série ADX, Vol. 1, Cahiers Institut Géologie, Neuchâtel, Suisse, p. 12, 1983. 

Kübler, B.: Cristallinité de l'illite, méthodes normalisées de préparations, méthodes normalisées de mesures, in: série ADX, Vol. 1, Cahiers Institut Géologie, Neuchâtel, Suisse, p. 13, 1987. 

Kübler, B. and Jaboyedoff, M.: Illite crystallinity, Comptes Rendus de l'Académie des Sciences-Series IIA – Earth and Planetary Science, 331, 75–89,, 2000. 

Labaume, P. and Teixell, A.: 3D structure of subsurface thrusts in the eastern Jaca Basin, southern Pyrenees, Geol. Acta, 16, 477–498,, 2018. 

Labaume, P., Meresse, F., Jolivet, M., and Teixell, A.: Exhumation sequence of the basement thrust units in the west-central Pyrenees. Constraints from apatite fission track analysis, Geogaceta, 60, 11–14, 2006. 

Labourdette, R.: Stratigraphy and static connectivity of braided fluvial deposits of the lower Escanilla Formation, south central Pyrenees, Spain, AAPG Bull., 95, 585–617,, 2011. 

Lafont, F.: Influences relatives de la subsidence et de l'eustatisme sur la localisation et la géométrie des réservoirs d'un système deltaïque, Exemple de l'Eocène du bassin de Jaca, Pyrénées espagnoles, PhD thesis, Université Rennes, Rennes, (last access: 25 February 2023), 1994. 

Lagabrielle, Y., Labaume, P., and de Saint Blanquat, M.: Mantle exhumation, crustal denudation, and gravity tectonics during Cretaceous rifting in the Pyrenean realm (SW Europe): Insights from the geological setting of the lherzolite bodies, Tectonics, 29, TC4012,, 2010. 

Läuchli, C., Garcés, M., Beamud, E., Valero, L., Honegger, L., Adatte, T., and Castelltort, S.: Magnetostratigraphy and stable isotope stratigraphy of the middle-Eocene succession of the Ainsa basin (Spain): New age constraints and implications for sediment delivery to the deep waters, Mar. Petrol. Geol., 132, 105182,, 2021. 

Lupker, M., France-Lanord, C., Lavé, J., Bouchez, J., Galy, V., Métivier, F., and Mugnier, J. L.: A Rouse-based method to integrate the chemical composition of river sediments: Application to the Ganga basin, J. Geophys. Res.-Earth, 116, F04012,, 2011. 

Marshall, J. D.: Climatic and oceanographic isotopic signals from the carbonate rock record and their preservation, Geolog. Mag., 129, 143–160,, 1992. 

Millán, H., Aurell, M., and Meléndez, A.: Synchronous detachment folds and coeval sedimentation in the Prepyrenean External Sierras (Spain): a case study for a tectonic origin of sequences and systems tracts, Sedimentology, 41, 1001–1024,, 1994. 

Millán, H., Morer, E. P., Cardona, M. A., Aguado, A. L., Urcia, B. O., and Peña, B. M.: Actividad tectónica registrada en los depósitos terciarios del frente meridional del Pirineo central, Revista de la Sociedad Geológica de España, 13, 279–300, 2000. 

Miller, K. G., Browning, J. V., Schmelz, W. J., Kopp, R. E., Mountain, G. S., and Wright, J. D.: Cenozoic sea-level and cryospheric evolution from deep-sea geochemical and continental margin records, Sci. Adv., 6, eaaz1346,, 2020. 

Mochales, T., Barnolas, A., Pueyo, E. L., Serra-Kiel, J., Casas, A. M., Samsó, J. M., and Sanjuán, J.: Chronostratigraphy of the Boltaña anticline and the Ainsa Basin (southern Pyrenees), Bulletin, 124, 1229–1250,, 2012. 

Moebius, I., Friedrich, O., and Scher, H. D.: Changes in Southern Ocean bottom water environments associated with the Middle Eocene Climatic Optimum (MECO), Palaeogeogr. Palaeocl. Palaeoecol., 405, 16–27,, 2014. 

Moebius, I., Friedrich, O., Edgar, K. M., and Sexton, P. F.: Episodes of intensified biological productivity in the subtropical Atlantic Ocean during the termination of the Middle Eocene Climatic Optimum (MECO), Paleoceanography, 30, 1041–1058,, 2015. 

Moore, D. M. and Reynolds, R. C.: X-ray diffraction and the identification and analysis of clay minerals, Oxford University Press, Oxford, New York, p. 400, ISBN 9780195051704, 1997. 

Mulch, A., Chamberlain, C. P., Cosca, M. A., Teyssier, C., Methner, K., Hren, M. T., and Graham, S. A.: Rapid change in high-elevation precipitation patterns of western North America during the Middle Eocene Climatic Optimum (MECO), Am. J. Sci., 315, 317–336,, 2015. 

Muñoz, J. A.: Evolution of a continental collision belt: ECORS-Pyrenees crustal balanced cross-section, in: Thrust Tectonics, edited by: McClay, K. R., Springer, Dordrecht,, 1992. 

Muñoz, J. A., McClay, K., and Poblet, J.: Synchronous extension and contraction in frontal thrust sheets of the Spanish Pyrenees, Geology, 22, 921–924,<0921:SEACIF>2.3.CO;2, 1994. 

Muñoz, J. A., Beamud, E., Fernández, O., Arbués, P., Dinarès-Turell, J., and Poblet, J.: The Ainsa Fold and thrust oblique zone of the central Pyrenees: Kinematics of a curved contractional system from paleomagnetic and structural data, Tectonics, 32, 1142–1175,, 2013. 

Muñoz, J. A., Mencos, J., Roca, E., Carrera, N., Gratacós, O., Ferrer, O., and Fernández, O.: The structure of the South-Central-Pyrenean fold and thrust belt as constrained by subsurface data, Geolog. Acta, 16, 439–460,, 2018. 

Mutti, E.: Distinctive thin-bedded turbidite facies and related depositional environments in the Eocene Hecho Group (South-central Pyrenees, Spain), Sedimentology, 24, 107–131,, 1977. 

Ogg, J. G., Ogg, G., and Gradstein, F. M.: A concise geologic time scale: 2016, Elsevier, Cambridge University Press, Cambridge, p. 234, ISBN 978-0-444-63771-0, 2016. 

Pälike, H., Lyle, M. W., Nishi, H., Raffi, I., Ridgwell, A., Gamage, K., and Baldauf, J.: A Cenozoic record of the equatorial Pacific carbonate compensation depth, Nature, 488, 609–614,, 2012. 

Patterson, P. P. and Walter, L. M.: Depletion of 13C in seawater ΣCO2 on modern carbonate platforms: Significance for the carbon isotopic record of carbonates, Geology, 22, 885–888,<0885:DOCISC>2.3.CO;2, 1994. 

Popp, B. N., Takigiku, R., Hayes, J. M., Louda, J. W., and Baker, E. W.: The post-Paleozoic chronology and mechanism of 13C depletion in primary marine organic matter, Am. J. Sci., 289, 436–454,, 1989. 

Pueyo, E. L., Millán, H., and Pocovı, A.: Rotation velocity of a thrust: a paleomagnetic study in the External Sierras (Southern Pyrenees), Sediment. Geol., 146, 191–208,, 2002. 

Puigdefàbregas, C.: La sedimentación molásica en la cuenca de Jaca, Pirineos, CSIC, PhD Thesis, (last access: 25 February 2023), 1975.. 

Puigdefàbregas, C. and Souquet, P.: Tecto-sedimentary cycles and depositional sequences of the Mesozoic and Tertiary from the Pyrenees, Tectonophysics, 129, 173–203,, 1986. 

Pujalte, V., Baceta, J. I., and Schmitz, B.: A massive input of coarse-grained siliciclastics in the Pyrenean Basin during the PETM: the missing ingredient in a coeval abrupt change in hydrological regime, Clim. Past, 11, 1653–1672,, 2015. 

Remacha, E. and Fernández, L. P.: High-resolution correlation patterns in the turbidite systems of the Hecho Group (South-Central Pyrenees, Spain), Mar. Petrol. Geol., 20, 711–726,, 2003. 

Rivero-Cuesta, L., Westerhold, T., Agnini, C., Dallanave, E., Wilkens, R. H., and Alegret, L.: Paleoenvironmental changes at ODP Site 702 (South Atlantic): anatomy of the Middle Eocene Climatic Optimum, Paleoceanogr. Paleocl., 34, 2047–2066,, 2019. 

Rodriguez-Pintó, A., Pueyo, E. L., Barnolas, A., Pocoví, A., Oliva-Urcia, B., and Ramón, M. J.: Overlapped paleomagnetic vectors and fold geometry: a case study in the Balzes anticline (Southern Pyrenees), Phys. Earth Planet. Inter., 215, 43–57,, 2013. 

Roigé, M.: Procedència i evolució dels sistemes sedimentaris de la conca de Jaca (conca d'avantpaís Sudpirinenca): Interacció entre diverses àrees font en un context tectònic actiu, Tesis Doctoral, Universitat Autònoma de Barcelona, Barcelona, p. 315, ISBN 9788449079177, 2018. 

Roigé, M., Gómez-Gras, D., Remacha, E., Daza, R., and Boya, S.: Tectonic control on sediment sources in the Jaca basin (Middle and Upper Eocene of the South-Central Pyrenees), Comptes Rendus Geoscience, 348, 236–245,, 2016. 

Roure, F., Choukroune, P., Berastegui, X., Munoz, J. A., Villien, A., Matheron, P., and Deramond, J.: ECORS deep seismic data and balanced cross sections: Geometric constraints on the evolution of the Pyrenees, Tectonics, 8, 41–50,, 1989. 

Saltzman, M. R. and Thomas, E.: Carbon isotope stratigraphy, in: The geologic time scale:, edited by: Gradstein, F., Ogg, J., Schmitz, M. D., and Ogg, G., Elsevier, Oxford, UK, 207–232,, 2012. 

Schmitz, B. and Pujalte, V.: Abrupt increase in seasonal extreme precipitation at the Paleocene-Eocene boundary, Geology, 35, 215–218,, 2007. 

Schrag, D. P., DePaolo, D. J., and Richter, F. M.: Reconstructing past sea surface temperatures: Correcting for diagenesis of bulk marine carbonate, Geochim. Cosmochim. Ac., 59, 2265–2278,, 1995. 

Seguret, M.: Étude tectonique des nappes et séries décollées de la partie centrale du versant sud des Pyrénées, Pub. Ustela, Série Géologie structurale 2, PhD thesis, Montpellier, p. 155, 1972. 

Silva-Casal, R., Aurell, M., Payros, A., Pueyo, E. L., and Serra-Kiel, J.: Carbonate ramp drowning caused by flexural subsidence: the South Pyrenean middle Eocene foreland basin, Sediment. Geol., 393, 105538,, 2019. 

Sluijs, A., Zeebe, R. E., Bijl, P. K., and Bohaty, S. M.: A middle Eocene carbon cycle conundrum, Nat. Geosci., 6, 429–434,, 2013. 

Spangeberg, J. E.: Bulk C, H, O, and fatty acid C stable isotope analyses for purity assessment of vegetable oils from the southern and northern hemispheres, Rapid Commun. in Mass Spectrom., 30, 2447–2461,, 2016. 

Spangeberg, J. E. and Herlec, U.: Hydrocarbon biomarkers in the Topla-Mezica zinc-lead deposits, northern Karavanke/Drau range, Slovenia: paleoenvironment at the site of ore formation, Econ. Geol., 101, 997–1021,, 2006. 

Spofforth, D. J. A., Agnini, C., Pälike, H., Rio, D., Fornaciari, E., Giusberti, L., Luciani, V., Lanci, L., and Muttoni, G.: Organic carbon burial following the middle Eocene climatic optimum in the central western Tethys., Paleoceanography, 25, PA3210,, 2010. 

Sternai, P., Caricchi, L., Pasquero, C., Garzanti, E., Hinsbergen, D. J. J., and Castelltort, S.: Magmatic Forcing of Cenozoic Climate?, J. Geophys. Res.-Solid, 125, e2018JB01646,, 2020. 

Swart, P. K.: The geochemistry of carbonate diagenesis: The past, present and future, Sedimentology, 62, 1233–1304,, 1995. 

Sztrákos, K. and Castelltort, S.: La sédimentologie et les foraminifères bartoniens et priaboniens des coupes d'Arguis (Prépyrénées aragonaises, Espagne). Incidence sur la corrélation des biozones à la limite Bartonien/Priabonien, Revue de Micropaléontologie, 44, 233–247,, 2001. 

Teixell, A.: The Ansó transect of the southern Pyrenees: basement and cover thrust geometries, J. Geol. Soc., 153, 301–310,, 1996. 

Teixell, A.: Crustal structure and orogenic material budget in the west central Pyrenees, Tectonics, 17, 395–406,, 1998. 

Teixell, A., Labaume, P., and Lagabrielle, Y.: The crustal evolution of the west-central Pyrenees revisited: inferences from a new kinematic scenario, Comptes Rendus Geoscience, 348, 257–267,, 2016. 

Teixell, A., Labaume, P., Ayarza, P., Espurt, N., de Saint Blanquat, M., and Lagabrielle, Y.: Crustal structure and evolution of the Pyrenean-Cantabrian belt: A review and new interpretations from recent concepts and data, Tectonophysics, 724, 146–170,, 2018. 

Tribovillard, N., Algeo, T. J., Lyons, T., and Riboulleau, A.: Trace metals as paleoredox and paleoproductivity proxies: an update, Chem. Geol., 232, 12–32,, 2006. 

van der Boon, A., Kuiper, K. F., van der Ploeg, R., Cramwinckel, M. J., Honarmand, M., Sluijs, A., and Krijgsman, W.: Exploring a link between the Middle Eocene Climatic Optimum and Neotethys continental arc flare-up, Clim. Past, 17, 229–239,, 2021. 

van der Ploeg, R., Selby, D., Cramwinckel, M. J., Li, Y., Bohaty, S. M., Middlelburg, J. J., and Sluijs, A.: Middle Eocene greenhouse warming facilitated by diminished weathering feedback, Nat. Commun., 9, 2877,, 2018. 

van der Weijden, C. H.: Pitfalls of normalization of marine geochemical data using a common divisor, Mar. Geol., 184, 167–187,, 2002. 

van der Weijden, C. H., Reichart, G. J., and van Os, B. J.: Sedimentary trace element records over the last 200 kyr from within and below the northern Arabian Sea oxygen minimum zone, Mar. Geol., 231, 69–88,, 2006. 

Van Wagoner, J. C.: Sequence stratigraphy and marine to nonmarine facies architecture of foreland basin strata, Book Cliffs, Utah, USA, Reply, AAPG Bull., 82, 1607–1618, 1998. 

Van Wagoner, J. C., Mitchum, R. M., Campion, K. M., and Rahmanian, V. D.: Siliciclastic sequence stratigraphy in well logs, cores, and outcrops: concepts for high-resolution correlation of time and facies, American Association of Petroleum Geologists, Tulsa, Oklahoma, p. 55, ISBN 978-0891816577, 1990. 

Vergés, J., Fernàndez, M., and Martìnez, A.: The Pyrenean orogen: pre-, syn-, and post-collisional evolution, J. Virtual Exp., 8, 55 74,, 2002. 

Verité, J.: Enregistrement sédimentaire et climatique d'un hyperthermal en domaine continental: l'Optimum Climatique de l'Éocène Moyen dans le domaine Sud-Pyrénéen, Formation d'Escanilla, Espagne, MS Thesis, Observatoire des Sciences de l'Univers de Rennes, Rennes, France, 2019. 

Vinyoles, A., López-Blanco, M., Garcés, M., Arbués, P., Valero, L., Beamud, E., and Cabello, P.: 10 Myr evolution of sedimentation rates in a deep marine to non-marine foreland basin system: Tectonic and sedimentary controls (Eocene, Tremp–Jaca Basin, Southern Pyrenees, NE Spain), Basin Res., 33, 447–477,, 2021. 

Wendler, I.: A critical evaluation of carbon isotope stratigraphy and biostratigraphic implications for Late Cretaceous global correlation, Earth-Sci. Rev., 126, 116–146,, 2013. 

Whittaker, A. C., Duller, R. A., Springett, J., Smithells, R. A., Whitchurch, A. L., and Allen, P. A.: Decoding downstream trends in stratigraphic grain size as a function of tectonic subsidence and sediment supply, Geol. Soc. Am. Bull., 123, 1363–1382,, 2011.  

Westerhold, T. and Röhl, U.: Orbital pacing of Eocene climate during the Middle Eocene Climate Optimum and the chron C19r event: Missing link found in the tropical western Atlantic, Geochem, Geophy, Geosy,, 14, 4811–4825,, 2013. 

Zachos, J., Pagani, M., Sloan, L., Thomas, E., and Billups, K.: Trends, rhythms, and aberrations in global climate 65 Ma to present, Science, 292, 686–693,, 2001. 

Short summary
The Middle Eocene Climatic Optimum (MECO) was a global warming event that took place 40 Myr ago and lasted ca. 500 kyr, inducing physical, chemical, and biotic changes on the Earth. We use stable isotopes to identify the MECO in the Eocene deltaic deposits of the Southern Pyrenees. Our findings reveal enhanced deltaic progradation during the MECO, pointing to the important impact of global warming on fluvial sediment transport with implications for the consequences of current climate change.