Modelling the 2012 Lahar in a Sector of Jamapa Gorge (Pico de Orizaba Volcano, Mexico) Using RAMMS and Tree-Ring Evidence

A good understanding of the frequency and magnitude of lahars is essential for the assessment of torrential hazards in volcanic terrains. In many instances, however, data on past events is scarce or incomplete, such that the evaluation of possible future risks and/or the planning of adequate countermeasures can only be done with rather limited certainty. In this paper, we present a multiidisciplinary approach based on botanical field evidence and the numerical modelling of a post-eruptive lahar that occurred in 2012 on the northern slope of the Pico de Orizaba volcano, Mexico, with the aim of reconstructing the magnitude of the event. To this end, we used the debris-flow module of the rapid mass movement simulation tool RAMMS on a highly resolved digital terrain model obtained with an unmanned aerial vehicle. The modelling was calibrated with scars found in 19 Pinus hartwegii trees that served as paleo stage indicators (PSI) of lahar magnitude in a sector of Jamapa Gorge. Using this combined assessment and calibration of RAMMS, we obtain a peak discharge of 78 m3 s−1 for the 2012 lahar event which was likely triggered by torrential rainfall during hurricane “Ernesto”. Results also show that the deviation between the modelled lahar stage (depth) and the height of PSI in trees was up to ±0.43 m. We conclude that the combination of PSI and models can be successfully used on (subtropical) volcanoes to assess the frequency, and even more so to calibrate the magnitude of lahars. The added value of the approach is particularly obvious in catchments with very scarce or no hydrological data at all and could thus also be employed for the dating and modelling of older lahars. As such, the approach and the results obtained can be used directly to support disaster risk reduction strategies at Pico de Orizaba volcano, but also in other volcanic regions.


Introduction
Lahars are torrential mass movements on volcanic slopes and have been documented in various volcanic environments around the world. Being mixtures of water and sediment, they can mobilize large amounts of weakly consolidated debris, rocks, and wood. These processes are know for their high flow velocities, and for their potential to cause major destruction and loss of lives [1][2][3]. Lahars are commonly triggered by heavy rainfall, condensate gas coming from volcanoes, excessive pore water, or glacial and snow meltdown, but the process has also been reported to occur as a result of river erosion or dam overflow [4,5].
On active volcanoes, so-called syn-eruptive lahars may be formed after heavy downpours onto unconsolidated volcanic sediment or by the sudden melting of snow and ice during volcanic eruptions. Similar processes can occur in the absence of volcanic activity; in that case, they are referred to as intra-eruptive or post-eruptive lahars [6,7]. This last type of lahars can occur on active, dormant or extinct volcanic complexes and can persist long after the last volcanic activity [8][9][10]. Intra-eruptive and post-eruptive lahars are frequently triggered by enhanced hurricane activity in subtropical environments or by cyclonic events (i.e., torrential rainfall) [8,11], but seismic activity has been shown to trigger lahars as well [4,12]. Different approaches have been used in the past to understand the occurrence and triggering of lahars. Past work included rainfall-runout analyses [8,11,[13][14][15], real-time monitoring of lahars with rainfall, seismic (geophones) and observational data [16,17], but also analyses of disturbances recorded in tree-ring series [18][19][20][21].
Lahar magnitude can be approximated by lahar discharge, i.e., the volume of water and sediment moving down a gorge per unit of time. It depends on the initial volume of water and sediment, the erosion and incorporation of secondary, exotic debris by lahars as they move downstream (bulking) [5], the incorporation of water from tributaries, and sediment deposition (debulking). Lahar discharges have been reported in the literature, and one of the largest contemporary lahars occurred in 2007 at Nevado del Huila (Colombia) [22]. Since they are widely applied to calibrate numerical models used in lahar hazard evaluation, accurate volume estimations are necessary to strengthen their results.
The magnitude of lahar events has been assessed frequently with models using different bidimensional hydraulic codes and input data; among the more commonly used modeling tools are LAHARZ [23], TITAN2D [24], FLO-2D [25] or RAMMS [26]. Thanks to these models, the magnitude of several lahars has been calculated, including, among others, the very large 2007 lahar at Nevado del Huila with an estimated maximum peak discharge of~13,000 m 3 s −1 [22]. Other studies investigated the process dynamics of the lahar caused by the 2010 Merapi eruption [27], but also looked at hypothetical lahars on the base scenario of historical lahar deposits from an event that occurred in 1877 on Cotopaxi volcano [28]. Additionally, Nocentini et al. [29] have used FLO-2D to perform debris-flow modelling in volcanic terrains of Italy.
On the most active volcanoes of central Mexico, several hydraulic models have been used in the past to study lahars, mostly with the aim of obtaining better insights into sedimentology and/or rheological parameters (discharge, velocity, lahar depth). The LAHARZ code has been used on Popocatépetl [30], Pico de Orizaba [31], Nevado de Toluca [32], San Martin at Tuxtlas [33], Tacaná [34] and Ceboruco [35] volcanoes. Subsequently, the FLO-2D code has been applied to analyze Popocatépetl's 2001 lahar [36] and in the Nexpayantla Gorge located on the northwestern slope of Popocatépetl volcano [37] to model the 2010 lahar. At Colima volcano, lahar peak discharges have been estimated for recent events as well using the FLO-2D code [16,38].
In general, all of the above-mentioned studies have contributed considerably to delimit the reach and lateral extension of lahar events and to understand the hazards and risks induced by lahars. The value and accuracy of modelling approaches will depend substantially on the availability of (post-event) field data or observations for the calibration and validation of models; however, such data are normally absent [39][40][41]. In particular, paleo-deposits (lahar terraces or deposits) are only rarely present on Mexican volcanoes due to omnipresent erosion in the high-gradient torrents. Furthermore, historical records of volcanic or geomorphic events are not generally available, especially in areas with limited population densities. In this context, the reconstruction of past events with tree-ring series can contribute significantly to an improved knowledge of the frequency and/or magnitude of past volcanic or geomorphic processes [19][20][21]42]. Hence, scarred trees have proven very useful for the dating and frequency analysis of past geomorphic events [43] and for the determination of paleostage indicators (PSI; [39]) as the level of scars on trees is considered to indicate the minimum flood level during a given event [44]. Several studies have used scar heights to calibrate two-dimensional (2D) hydraulic models and to estimate the flow discharge (magnitude) of paleofloods [45][46][47]. This methodology, originally used for torrential and river floods, has been adapted on volcanic terrains to estimate the discharge of a lahar on Canary Island, Spain, in 1997 [48].
In this paper, we aim at estimating the magnitude of a lahar that occurred in 2012 on the slopes of Pico de Orizaba volcano (Mexico) in a sector of Jamapa Gorge by combining botanical evidence (i.e., scars on tree trunks) and hydraulic modelling performed on a highly-resolved topography acquired with an unmanned aerial vehicle (UAV). To this end, we used RAMMS, a dynamic numerical model that allows simulation of debris flows (including water depth, discharge and flow velocity) in mountain basins [49][50][51].

Volcanic Setting and Recorded Lahar Activity
The Jamapa Gorge is located on the northern slopes of Pico de Orizaba volcano (5640 m above sea level (a.s.l.)), also known as Citlaltépetl ("Mountain of the Star") in the Náhuatl language. Pico de Orizaba is an active, yet quiescent Quaternary stratovolcano of andesic-dacitic composition located in the eastern sector of the Trans-Mexican Volcanic Belt ( Figure 1). The volcanic complex is characterized by a long history of eruptive episodes with repeat growth and collapse of domes and other volcanic structures. According to Carrasco-Núñez [52,53], the first volcanic paleo-structure of Pico de Orizaba can be dated to~0. 65 Ma BP (= Before Present) and is referred to as "Torrecillas". This volcanic structure was destroyed between~0.29 and~0.21 Ma by a large collapse associated to debris avalanche deposits on the northeastern slope of the volcano. The second paleo-structure was dated to~0.21 Ma BP and is known as "El Espolón de Oro". The cone belonging to this episode was destroyed by a large collapse and debris avalanche associated with hydrothermal alteration of the volcanic edifice dated to~20,000 years (or 20 kyr; the abbreviation will be used hereafter). The current cone (Pico de Orizaba), is~20 kyr in age. Two Plinian eruptions have been dated to the late Pleistocene and the Holocene, between~13 kyr and~8.5 kyr ago. These events covered the flanks of the volcano with thick pyroclastic flow and tephra fall deposits known as "Citlaltépetl" ignimbrite [54][55][56]. Other volcanic events during the Holocene were recorded between 7.0-6.2 and 3.4 kyr BP [57]. Block and ash flows and lahar deposits were produced at the summit crater around 4.1-4.0 kyr BP [54] in relation to a large explosive event in the western sector of the volcano. On the southern slope of Pico de Orizaba, two lava flows were dated to between 3 and 1 kyr BP [58]. Furthermore, historical effusive and explosive volcanic events were reported in 1537, 1545, 1566 and 1613 AD [59].
Holocene and historical lahar deposits can be found in different gorges and rivers around the Pico de Orizaba volcano. For example, Carrasco-Núñez et al. [52] dated a large lahar (~16.5 kyr BP) associated with the collapse of the Espolón de Oro cone; the resulting lahar deposits were emplaced in the Tetelzingo River where they still have thicknesses of 12-20 m on average and a maximum of 100 m. Other lahar deposits were identified inside Huitzilapan Gorge and dated to~5260 yr BP [12]. In the Huitzilapan river (north-northeastern slope of Pico de Orizaba), one can observe an important lahar deposit from an event that occurred in 1920 ("Huitzilapan lahar"). The 1920 lahar was a typical compounded event, as an earthquake (magnitude 6.5) triggered several landslides from slopes saturated with water after 10 days of torrential rainfall [12]. Iverson et al. [59] reported high water marks 40 to 60 m above the channel for the 1920 event, which allowed for estimation of its volume at 4.4 × 10 7 m 3 [31].
More recently, on 5 June 2003, a catastrophic debris flow occurred on the southern slope of Pico de Orizaba after a torrential rainfall, resulting in a peak discharge of 350 m 3 s −1 [60]. On the north-northeastern slopes of Pico de Orizaba volcano (i.e., Jamapa and Seca Gorges), intra-eruptive lahars occurred in 2012, again triggered by heavy rainfalls that can be linked to the passage of Hurricane Ernesto over the Atlantic Ocean [61]. Holocene and historical lahar deposits can be found in different gorges and rivers around the Pico de Orizaba volcano. For example, Carrasco-Núñez et al. [52] dated a large lahar (~16.5 kyr BP) associated with the collapse of the Espolón de Oro cone; the resulting lahar deposits were emplaced in the Tetelzingo River where they still have thicknesses of 12-20 m on average and a maximum of 100 m. Other lahar deposits were identified inside Huitzilapan Gorge and dated to ~5260 yr BP [12]. In the Huitzilapan river (north-northeastern slope of Pico de Orizaba), one can observe an important lahar deposit from an event that occurred in 1920 ("Huitzilapan lahar"). The 1920 lahar was a typical compounded event, as an earthquake (magnitude 6.5) triggered several landslides from slopes saturated with water after 10 days of torrential rainfall [12]. Iverson et al. [59] reported high water marks 40 to 60 m above the channel for the 1920 event, which allowed for estimation of its volume at 4.4 × 10 7 m 3 [31].
More recently, on 5 June 2003, a catastrophic debris flow occurred on the southern slope of Pico de Orizaba after a torrential rainfall, resulting in a peak discharge of 350 m 3 s −1 [60]. On the northnortheastern slopes of Pico de Orizaba volcano (i.e., Jamapa and Seca Gorges), intra-eruptive lahars occurred in 2012, again triggered by heavy rainfalls that can be linked to the passage of Hurricane Ernesto over the Atlantic Ocean [61].

Geomorphic Features on Jamapa Gorge
The Jamapa Gorge drains the northern flank of Pico de Orizaba volcano; its headwaters are at ~5000 m a.s.l. near the terminus of Jamapa glacier ( Figure 1). In that area, glacial landforms from the late Pleistocene and Holocene have been covered by Holocene pyroclastic deposits and lava flows. The glaciers of the Little Ice Age (LIA) covered the top of the Jamapa Gorge drainage basin and built terminal moraines ridges around 4700 m a.s.l. [62]. In addition, several paraglacial and postglacial Satellite imagery illustrating the volcano as well as the gorges mentioned. The Jamapa Gorge, where this study was realized, is located on the northern slope of the volcano. The study site is indicated with a red box.

Geomorphic Features on Jamapa Gorge
The Jamapa Gorge drains the northern flank of Pico de Orizaba volcano; its headwaters are at 5000 m a.s.l. near the terminus of Jamapa glacier ( Figure 1). In that area, glacial landforms from the late Pleistocene and Holocene have been covered by Holocene pyroclastic deposits and lava flows. The glaciers of the Little Ice Age (LIA) covered the top of the Jamapa Gorge drainage basin and built terminal moraines ridges around 4700 m a.s.l. [62]. In addition, several paraglacial and postglacial geomorphic features can be found in the area, such as debris flows or rockfall deposits and sectors characterized by fluvial erosion [63].
The study site investigated in this paper is located on the central sector of Jamapa Gorge, between 3535 and 3620 m a.s.l. (Figures 1 and 2A). Several old laharic terraces can be recognized with thicknesses of 1-4 m ( Figure 2B). Buried, tilted and scared trees are common in the laharic deposits of the low terraces ( Figure 2C). Superficial debris-flow deposits with block sizes of~3 m in diameter are embedded in a sandy matrix. Recent laharic deposits (from debris flows and/or hyperconcentrated flows) have filled the main channel with~2 m of sediments in this sector of Jamapa. Some of these lahars have reached the road connecting the federal states of Puebla and Veracruz ( Figure 2D). between 3535 and 3620 m a.s.l. (Figures 1 and 2A). Several old laharic terraces can be recognized with thicknesses of 1-4 m ( Figure 2B). Buried, tilted and scared trees are common in the laharic deposits of the low terraces ( Figure 2C). Superficial debris-flow deposits with block sizes of ~3 m in diameter are embedded in a sandy matrix. Recent laharic deposits (from debris flows and/or hyperconcentrated flows) have filled the main channel with ~2 m of sediments in this sector of Jamapa. Some of these lahars have reached the road connecting the federal states of Puebla and Veracruz ( Figure 2D).

The 2012 Hydrometeorological Event and Rainfall Triggering
On 6 August 2012, tropical storm "Ernesto" attained a hurricane strength at 1200 UTC (= Coordinated Universal Time) in the Caribbean, east of Chetumal, Mexico. The hurricane then turned westward and continued to strengthen as it approached the coast of the southern Yucatan Peninsula. The hurricane made landfall at Cayo Norte in the Banco Chinchorro Islands of Mexico around 0100 UTC on 8 August [64]. On 9 August, it caused torrential rainfalls in the range of 112 to 300 mm 24 h −1 in the states of Oaxaca, Veracruz, Hidalgo and Puebla. On 10 August, rainfall decreased from 126 to 180 mm 24 h −1 in Puebla, Guerrero, Chiapas, and Veracruz [65]. Based on records from the Huatusco meteorological station (located on the NE slopes of Pico de Orizaba, 1186 m a.s.l.) the maximum rainfall was recorded on 9 August 2012 with 133 mm, and reached 153 mm in two days (9-10 August 2012), 159 mm in three days (8-10 August 2012), and 212 mm in five days (9-13 August 2012) [66]. This hydrometeorological event triggered regional landslide and debris-flow events on the northern slopes of the volcano, i.e., in the Jamapa and Tliapa Gorges where they caused massive structural damage [61]. Field recognition showed several deposits left by the geomorphic events, namely in the form of terraces, debris and block deposits, scarps, and scars on several trees (Figure 2A-C). . On 9 August, it caused torrential rainfalls in the range of 112 to 300 mm 24 h −1 in the states of Oaxaca, Veracruz, Hidalgo and Puebla. On 10 August, rainfall decreased from 126 to 180 mm 24 h −1 in Puebla, Guerrero, Chiapas, and Veracruz [65]. Based on records from the Huatusco meteorological station (located on the NE slopes of Pico de Orizaba, 1186 m a.s.l.) the maximum rainfall was recorded on 9 August 2012 with 133 mm, and reached 153 mm in two days (9-10 August 2012), 159 mm in three days (8-10 August 2012), and 212 mm in five days (9-13 August 2012) [66]. This hydrometeorological event triggered regional landslide and debris-flow events on the northern slopes of the volcano, i.e., in the Jamapa and Tliapa Gorges where they caused massive structural damage [61]. Field recognition showed several deposits left by the geomorphic events, namely in the form of terraces, debris and block deposits, scarps, and scars on several trees (Figure 2A-C).

Digital Terrain Model
A UAV was used to construct a digital superficial model (DSM) and a digital terrain model (DTM) based on the structure-from-motion (SfM) photogrammetry of the landforms with high-resolution topography [67]. For the modelling of the 2012 lahar discharge in the study reach within Jamapa Gorge, we constructed a DTM based on 306 digital images taken with a DJI Phantom 4 and a Mavic Pro Platinum ( Figure 3A). Flight height was 70 m on average above the terrain surface with an image resolution of 4000 × 3000 dpi. Side overlap between images was 75% and front overlap was 40%. The processing and DTM extraction were performed using Pix4DMapper software (Pix4D SA, Switzerland) with a cell size of 3.6 cm ( Figure 3B).

Digital Terrain Model
A UAV was used to construct a digital superficial model (DSM) and a digital terrain model (DTM) based on the structure-from-motion (SfM) photogrammetry of the landforms with highresolution topography [67]. For the modelling of the 2012 lahar discharge in the study reach within Jamapa Gorge, we constructed a DTM based on 306 digital images taken with a DJI Phantom 4 and a Mavic Pro Platinum ( Figure 3A). Flight height was 70 m on average above the terrain surface with an image resolution of 4000 × 3000 dpi. Side overlap between images was 75% and front overlap was 40%. The processing and DTM extraction were performed using Pix4DMapper software (Pix4D SA, Switzerland) with a cell size of 3.6 cm ( Figure 3B).

Dendrogeomorphic Analysis and PSI
Within the study sector in Jamapa Gorge, we identified and sampled all trees with visible scars inflicted by the 2012 lahar between ~3535 to 3620 m a.s.l., at a point located ~7 km away from the Pico de Orizaba crater (Figure 1). To this end, we used increment borers and an electric chainsaw to obtain cores or wedges to properly date the injuries [43]. A total of 47 Pinus hartwegii trees were sampled, of which 38 were disturbed trees (with scars) and 9 were undisturbed trees for the construction of a reference chronology ( Figure 3B). In addition to their dating, scars were used as paleostage indicators (PSI) by measuring the height of impacts with respect to the channel bed [40]. In the case that a tree

Dendrogeomorphic Analysis and PSI
Within the study sector in Jamapa Gorge, we identified and sampled all trees with visible scars inflicted by the 2012 lahar between~3535 to 3620 m a.s.l., at a point located~7 km away from the Pico de Orizaba crater (Figure 1). To this end, we used increment borers and an electric chainsaw to obtain cores or wedges to properly date the injuries [43]. A total of 47 Pinus hartwegii trees were sampled, of which 38 were disturbed trees (with scars) and 9 were undisturbed trees for the construction of a reference chronology ( Figure 3B). In addition to their dating, scars were used as paleostage indicators (PSI) by measuring the height of impacts with respect to the channel bed [40]. In the case that a tree had scars located at different heights, different samples were taken to distinguish between different lahar events. We considered the uppermost point of the highest scar as the peak of lahar discharge ( Figure 4) [45,46]. In addition, for each tree, we also recorded the geomorphic and geographic position with a Global Positioning System (GPS) device. Additional measures were taken in the field to delineate the distance, direction and inclination of trees with respect to topographic features (such as big blocks, old terraces, or volcanic slopes) using a compass, tape measure and inclinometer.
After fieldwork, samples were prepared in the laboratory according to the standard procedures described in Bräker [68] and Stoffel and Bolslchweiler [69]. After the dating of samples, ring widths were measured by means of a microscope and a sliding Velmex [70] measurement plate connected to a computer using the time series analysis software TSAPWin TM (Time Series and Analysis Program for Windows) [71].
In addition to the dating of scars to the year [45,46], we also aimed to identify the seasonal timing of lahars based on the position of growth anomalies (i.e., injuries and callus tissue) within the tree ring according to Stoffel et al. [72]. The growing season was adapted to central Mexico conifers [20].
Water 2020, 12, 333 7 of 15 had scars located at different heights, different samples were taken to distinguish between different lahar events. We considered the uppermost point of the highest scar as the peak of lahar discharge ( Figure 4) [45,46]. In addition, for each tree, we also recorded the geomorphic and geographic position with a Global Positioning System (GPS) device. Additional measures were taken in the field to delineate the distance, direction and inclination of trees with respect to topographic features (such as big blocks, old terraces, or volcanic slopes) using a compass, tape measure and inclinometer. After fieldwork, samples were prepared in the laboratory according to the standard procedures described in Bräker [68] and Stoffel and Bolslchweiler [69]. After the dating of samples, ring widths were measured by means of a microscope and a sliding Velmex [70] measurement plate connected to a computer using the time series analysis software TSAPWin TM (Time Series and Analysis Program for Windows) [71].
In addition to the dating of scars to the year [45,46], we also aimed to identify the seasonal timing of lahars based on the position of growth anomalies (i.e., injuries and callus tissue) within the tree ring according to Stoffel et al. [72]. The growing season was adapted to central Mexico conifers [20].

Numerical Simulations
RAMMS has been developed in the Swiss Federal Institute for Forest, Snow and Landscape research (WSL) [73] to simulate the runout of muddy and debris-laden flows in complex terrain. The model uses the Voellmy-Salm fluid flow continuum model [74] based on the Voellmy fluid flow law [75]. In the model, the frictional resistance is divided into two parts: a dry-Coulomb type friction, proportional to the normal stress at the flow bottom (coefficient μ) and a viscous resistance turbulent friction depending on the square of the velocity (coefficient ξ (m/s 2 )). The model yields maximum water depth, sediment deposition, flow velocity, and pressure. For further details about RAMMS and the equations used, see [76,77].
In our study, the input parameters of the model were set up as follows. The lahar density was initially fixed to 1400 g/m 3 . The μ coefficient was defined on the basis of field observations at a value

Numerical Simulations
RAMMS has been developed in the Swiss Federal Institute for Forest, Snow and Landscape research (WSL) [73] to simulate the runout of muddy and debris-laden flows in complex terrain. The model uses the Voellmy-Salm fluid flow continuum model [74] based on the Voellmy fluid flow law [75]. In the model, the frictional resistance is divided into two parts: a dry-Coulomb type friction, proportional to the normal stress at the flow bottom (coefficient µ) and a viscous resistance turbulent friction depending on the square of the velocity (coefficient ξ (m/s 2 )). The model yields maximum water depth, sediment deposition, flow velocity, and pressure. For further details about RAMMS and the equations used, see [76,77].
In our study, the input parameters of the model were set up as follows. The lahar density was initially fixed to 1400 g/m 3 . The µ coefficient was defined on the basis of field observations at a value of 0.15, representing the average slope of the deposition zones [75]. The ξ coefficient was set to 400 m/s 2 according to the nature of the flux observed in similar contexts in Mexico [36]. The lahar peak discharge was obtained by applying an iterative step-backwater procedure [44]. For each of the simulation runs, we released a triangular hydrograph at the beginning of the river reach in which we defined the maximum peak discharge to occur after 60 s and the event to stop after 850 s [36]. Then, in order to estimate the peak discharge of the 2012 lahar, the maximum water depth modelled by RAMMS was compared with the maximum scar height observed in the 19 trees analyzed (Figures 3B and 4). The magnitude of the lahar event was then defined as the peak discharge that minimized the deviation between the observed scar height and the modelled lahar stage. Once the peak discharge of the 2012 lahar was fixed, we performed a sensitivity model analysis by varying the density, µ and ξ coefficients by ±25% with respect to the initial values and evaluating the changes in the deviation between observed and modelled water heights.

Tree Ring Evidence of the 2012 Lahar Event
In this study, a total of 38 trees with ages of up to 150 years have been sampled, all showing disturbances related to the 2012 lahar event. Most disturbances were in the form of scars, stem tilting and decapitation. Among these disturbances, 19 Pinus hartwegii trees had scars that were suitable to be used as PSI. The passage of the lahar left by far the largest amount of severe damage in trees that were growing in the main channel as well as in those individuals standing on the lower terraces (13 scars at heights 22-120 cm), whereas in six cases, scarred trees were located in a secondary channel (with scar heights of 3-70 cm). Scars were recorded on the upslope side of trunks, facing the direction of the lahar flow ( Figure 4). The analysis of scar seasonality showed that the injuries were located in the early latewood, corresponding to the late summer of 2012 ( Figure 5) and likely of the passage of hurricane "Ernesto".

Magnitude Reconstruction of the 2012 Lahar
We simulated the 2012 lahar in a trial-and-error approach with maximum peak discharges ranging from 25 to 500 m 3 s −1 to approximate and ultimately reduce the deviation between scar heights found in the analyzed trees and maximum lahar depths obtained in the model. Our results suggest a mean deviation ranging from −0.43 m for 25 m 3 s −1 to 3.08 m for 500 m 3 s −1 . Refining the approach further, the minimum absolute mean deviation was reached for a peak discharge of 78 m 3 s −1 , translating into a lahar volume of 33,000 m 3 based on the assumed shape of the hydrograph (see Figure 6D). Varying lahar density and the ξ coefficients by ±25% indicate that the mean deviation can change by up tõ ±0.43 m and~±0.47 m, respectively (Table 1). Figure 6A shows the relationship between the maximum water depth and scar heights, as well as flow velocity and stream power for the estimated peak discharge. Deviations between the lahar model and scar heights were between −0.68 m and 1.18 m. Relatively small deviations (±0.25 m; 50th percentile) were observed in 52% of the cases. Moderate deviations (ranging between ±0.25 and ±0.35 m; 50th-75th percentile) could be found in 26% of the cases, whereas 22% of the scars found in trees showed deviations exceeding ±0.35 m (>75th percentile). Only one tree showed deviations exceeding 1 standard deviation (σ = 0.77 m; Table 2). Noteworthy, most of the trees for which deviations were more relevant were growing at sites where they tend to be quite exposed to the flow, i.e., especially in those parts of the reach with the steepest slopes and unruffled bedrock and where specific kinetic energy tends to be higher as well. The maximum lahar velocity simulated by RAMMS ( Figure 6B) occurred in the central part of the main channel with~18 m/s. of 0.15, representing the average slope of the deposition zones [75]. The ξ coefficient was set to 400 m/s 2 according to the nature of the flux observed in similar contexts in Mexico [36]. The lahar peak discharge was obtained by applying an iterative step-backwater procedure [44]. For each of the simulation runs, we released a triangular hydrograph at the beginning of the river reach in which we defined the maximum peak discharge to occur after 60 s and the event to stop after 850 s [36]. Then, in order to estimate the peak discharge of the 2012 lahar, the maximum water depth modelled by RAMMS was compared with the maximum scar height observed in the 19 trees analyzed ( Figures 3B  and 4). The magnitude of the lahar event was then defined as the peak discharge that minimized the deviation between the observed scar height and the modelled lahar stage. Once the peak discharge of the 2012 lahar was fixed, we performed a sensitivity model analysis by varying the density, μ and ξ coefficients by ±25% with respect to the initial values and evaluating the changes in the deviation between observed and modelled water heights.

Tree Ring Evidence of the 2012 Lahar Event
In this study, a total of 38 trees with ages of up to 150 years have been sampled, all showing disturbances related to the 2012 lahar event. Most disturbances were in the form of scars, stem tilting and decapitation. Among these disturbances, 19 Pinus hartwegii trees had scars that were suitable to be used as PSI. The passage of the lahar left by far the largest amount of severe damage in trees that were growing in the main channel as well as in those individuals standing on the lower terraces (13 scars at heights 22-120 cm), whereas in six cases, scarred trees were located in a secondary channel (with scar heights of 3-70 cm). Scars were recorded on the upslope side of trunks, facing the direction of the lahar flow ( Figure 4). The analysis of scar seasonality showed that the injuries were located in the early latewood, corresponding to the late summer of 2012 ( Figure 5) and likely of the passage of hurricane "Ernesto".

Magnitude Reconstruction of the 2012 Lahar
We simulated the 2012 lahar in a trial-and-error approach with maximum peak discharges ranging from 25 to 500 m 3 s −1 to approximate and ultimately reduce the deviation between scar heights found in the analyzed trees and maximum lahar depths obtained in the model. Our results suggest a mean deviation ranging from −0.43 m for 25 m 3 s −1 to 3.08 m for 500 m 3 s −1 . Refining the approach further, the minimum absolute mean deviation was reached for a peak discharge of 78 m 3 s −1 , translating into a lahar volume of 33,000 m 3 based on the assumed shape of the hydrograph (see   [36]. The dashed black line represents the assumed shape of the hydrograph used in this study.

Discussion and Conclusions
In this paper, we provide insights into a tree-ring based lahar discharge reconstruction for an event that occurred in the Jamapa Gorge (northern slopes of Pico de Orizaba volcano, Mexico) during August 2012, likely as a result of the passage of hurricane "Ernesto". To this end, we used scar heights found in 19 Pinus hartwegii trees to calibrate simulation runs realized with the two-dimensional RAMMS model. Even if the combined use of tree-ring evidence and model outputs has been used previously for peak discharge estimation of fluvial processes [40], this study likely is the first of its kind to reconstruct peak discharge of a lahar event with botanical evidence and models.
The intra-annual dating of 19 scars found in trees (Figures 3B) indicate that the lahar event must have taken place in early latewood which locally corresponds to late summer ( Figure 5); this timing is consistent with the passage of hurricane "Ernesto" on 9 August 2012 at the study site. Evidence of lahars triggered by the 2012 hurricane can be found elsewhere in central Mexico, namely on Popocatépetl [20]; Iztaccíhuatl [78] and La Malinche [41] volcanoes. The timing of the 2012 lahar at our study site is further supported by evidence from local newspapers [79]. However, the 2012 lahars in Jamapa Gorge could also have been favored by water added to the system as a result of progressive ice-melting at Jamapa glacier in recent years [61].
In Mexico, dendrogeomorphic approaches have been used to date and reconstruct lahars [18][19][20] and debris flows [21,41] in the past. By contrast, studies focusing on the magnitude of hydrogeomorphic events are still critically lacking. The approach that we present here has allowed for retrospective estimation of the peak discharge of the 2012 lahar at 78 m 3 s −1 at the level of the study reach with a ~3 km 2 catchment area. For this modelled peak discharge, the difference between the maximum height of the scars and lahar stage modelled were minimized. The reconstructed peak discharge is somehow lower than the lahars of Nevado del Huila [22], which was estimated at 13,500  [36]. The dashed black line represents the assumed shape of the hydrograph used in this study.

Discussion and Conclusions
In this paper, we provide insights into a tree-ring based lahar discharge reconstruction for an event that occurred in the Jamapa Gorge (northern slopes of Pico de Orizaba volcano, Mexico) during August 2012, likely as a result of the passage of hurricane "Ernesto". To this end, we used scar heights found in 19 Pinus hartwegii trees to calibrate simulation runs realized with the two-dimensional RAMMS model. Even if the combined use of tree-ring evidence and model outputs has been used previously for peak discharge estimation of fluvial processes [40], this study likely is the first of its kind to reconstruct peak discharge of a lahar event with botanical evidence and models.
The intra-annual dating of 19 scars found in trees ( Figure 3B) indicate that the lahar event must have taken place in early latewood which locally corresponds to late summer ( Figure 5); this timing is consistent with the passage of hurricane "Ernesto" on 9 August 2012 at the study site. Evidence of lahars triggered by the 2012 hurricane can be found elsewhere in central Mexico, namely on Popocatépetl [20]; Iztaccíhuatl [78] and La Malinche [41] volcanoes. The timing of the 2012 lahar at our study site is further supported by evidence from local newspapers [79]. However, the 2012 lahars in Jamapa Gorge could also have been favored by water added to the system as a result of progressive ice-melting at Jamapa glacier in recent years [61].
In Mexico, dendrogeomorphic approaches have been used to date and reconstruct lahars [18][19][20] and debris flows [21,41] in the past. By contrast, studies focusing on the magnitude of hydrogeomorphic events are still critically lacking. The approach that we present here has allowed for retrospective estimation of the peak discharge of the 2012 lahar at 78 m 3 s −1 at the level of the study reach with ã 3 km 2 catchment area. For this modelled peak discharge, the difference between the maximum height of the scars and lahar stage modelled were minimized. The reconstructed peak discharge is somehow lower than the lahars of Nevado del Huila [22], which was estimated at 13,500 m 3 s −1 (maximum peak discharge) with a~26 km 2 catchment area (Páez River). By contrast, the calculated peak discharge is higher than the 48 m 3 s −1 estimated for "Patrio" lahar at Colima volcano [16] (Montegrande catchment area:~6.2 km 2 ). In addition, the shape of the channel [30], lahar rheologic characteristics [36] and resolution of the digital terrain model [48] can affect the lahar scenarios.
Documented differences between observed scar heights and modeled flow heights are comparable with results obtained in past work focusing on the fluvial domain. The largest deviations (exceeding ±0.35 m or the >75th percentile of deviations) were found in trees located on cut banks and/or in straight and exposed positions along the investigated river reach at the lower end of the channel. By contrast, moderate and smaller deviations were typically found in trees located in intermediate and low-energy landforms, such as alluvial terraces located in the central and lower end of the channel. These results are in agreement with previous studies carried out in the fluvial domain [40,[80][81][82][83], suggesting that intermediate-energy landforms such as alluvial terraces could be considered appropriate geomorphic locations for tree-ring-based reconstructions of lahar magnitude. Here, instead of using a non-Newtonian flow-like model for high concentrated flows [48], we employed a single-phase model to simulate inertial flows of solid-fluid mixtures such as debris flows [84]. In this study, the calibration of RAMMS has been based exclusively on scar heights and the extension of fresh deposits observed during field recognition. Despite the fact that we performed a sensitivity analysis on the main parameters, i.e., ±25% on lahar density and ξ coefficients (µ was fixed according to field recognition), uncertainties remain regarding the total volume of the event.
The combination of tree-ring based lahar reconstruction and process modelling with RAMMS is the first of its kind and has helped considerably to estimate the magnitude of the 2012 lahar in a sector of Jamapa Gorge. As such, we conclude that the combination of 2D dimensional dynamic modeling and tree-ring analysis can be very useful to assess and calibrate the magnitude of lahars in other temperate forests growing on Mexican volcanoes, especially in areas with scarce or inexistent hydrologic data and/or for the dating and assessment of older lahars.