The 2007 caldera collapse at Piton de la Fournaise: new insights from multi-temporal structure-from-motion

We produced new multi-temporal Digital Elevation Models (DEMs) of the April 2007 summit collapse at Piton de la Fournaise from previously unused aerial photographs. This dataset reveals the precise temporal evolution of the collapsed volume and caldera morphological changes during the event. It provides a unique opportunity to study caldera formation, one of the most hazardous natural phenomena, for which relatively little scientiﬁc and quantiﬁed information is available. During this rare example of observed caldera formation, the summit started to collapse four days after the onset of a high-volume eruption at an unusually low elevation (at 20:48 UTC on April 5 th ). Our new data show that during the ﬁrst 30 hours, collapse was relatively fast (840m 3 s − 1 average), and continued for at least the following 12 days, at a slower rate (46m 3 s − 1 average), which had not previously been reported. On April 19 th , the collapse reached 96 % of its ﬁnal volume, while the remaining 4 % was probably attained by May 1 st (end of lava emission at the vent). New infrared 3D mapping of the caldera ﬂoor made a year after the event demonstrates that post-collapse hydrothermal activity in the caldera is closely associated with the main ring faults active during the collapse, which are now preferential paths for ﬂuids to reach the surface.


Introduction
In volcanic systems, rapid withdrawal of large magma volumes from reservoirs can lead to the formation of calderas -broad depressions at the top or on the sides of volcano flanks. These large but recurring events pose multiple hazards to neighboring communities (e.g. caldera collapses, summit instabilities, high gas fluxes, ash plumes, ash deposits, rapid and extended lava flows in distal areas). Even for calderas of classic morphologies on basaltic shield volcanoes, scientific records of observations during their formation are rare: Kīlauea volcano, Hawai'i, in 1960 [Delaney andMcTigue 1994]; Fernandina volcano, Galápagos, in 1968 [Simkin andHoward 1970]; Miyakejima, Japan, in 2000 [Geshi et al. 2002]; Piton de la Fournaise volcano (PdF), La Réunion, in 2007[e.g. Staudacher et al. 2009]; Bárðarbunga volcano, Iceland, in 2014 [e.g. Gudmundsson et al. 2016]; and Kīlauea again in 2018 [e.g. Neal et al. 2019]. Detailing and understanding collapse mechanisms may lead to significant breakthroughs in volcano dynamics science. In particular, detailed measurement of the collapsed volumes and precise mapping of ring faults during and after the caldera formation is of paramount importance when determining the precise geometry and location of an underlying magma reservoir. The April 2007 col-Corresponding author: allan.derrien@yahoo.fr lapse of the PdF summit is part of a wide range of collapse processes due to the withdrawal of a large volume of magma from an underlying reservoir. Its relatively small volume (approximately 0.1 km 3 according to Staudacher et al. [2009]) places it at the "small" end-member of caldera collapses, with much larger examples observed at other locations (e.g. 1 to 2 km 3 at Fernandina volcano in 1968 [Acocella et al. 2015] or 1.8 km 3 at Bárðarbunga volcano in 2014 [Gudmundsson et al. 2016]).
The style of the April 2007 event has been suggested to represent a plate or "piston" mechanism (following the caldera end-members proposed in the literature [e.g. Acocella 2007;Cole et al. 2005;Lipman 1997]), where a single main rock column above the reservoir is destabilised and collapses, coeval with voluminous magma withdrawal. Such events are triggered by magma exiting a magma reservoir through intrusions (dykes and/or sills). This intrusion may reach the surface far from the source reservoir (e.g. >45 km at Bárðarbunga in 2014 [e.g. Gudmundsson et al. 2016]), with two implications: (1) the length of the intrusion implies a significant volume of magma is tapped from the reservoir during the intrusion process and (2) the location of the vent, far from the reservoir, is often also at lower elevation, favouring rapid and voluminous lava emission (the difference in altitude between the source reservoir and the vent in this case being smaller 4D mapping of caldera collapse at Piton de la Fournaise Derrien et al., 2020 than for eruptions occurring at higher elevation close to the summit). Thus, collapses at volcano summit are often coeval with the beginning of an eruption located far from the usual active center. In later stages of the caldera formation, it is also possible that magma uses the caldera faults, sometimes rooted at the reservoir roof, as preferential pathways to reach the surface [e.g. Kennedy et al. 2018]. In these cases, the collapse process may also be coeval with, or precede, smaller volume, intra-caldera eruption(s) (e.g. the case of several eruptions inside PdF's Cratère Dolomieu in 2008[Staudacher 2010).
Piton de la Fournaise is located in the southwestern Indian Ocean ( Figure 1A). This hotspot basaltic shield volcano has undergone a series of summit caldera collapses in the last three centuries [e.g. Michon et al. 2013;Peltier et al. 2012], each forming a depression subsequently filled with the products of following summit eruptions. The latest of these events ( Figure 1B, C) began on April 5 th , 2007 (20:48 UTC) and was associated with an unusually voluminous distal eruptive phase starting on April 2 nd and ending on May 1 st . More than 140 × 10 6 m 3 of lava have been emitted according to Staudacher et al. [2009] and 240 × 10 6 m 3 according to Roult et al. [2012]. These large volumes widely emptied the shallow magma reservoir located at a depth of about 1-2 km [Michon et al. 2011;Michon et al. 2009;Muro et al. 2014;Peltier et al. 2009;Staudacher et al. 2009]. This event was well documented by ground deformation ], seismicity [Fontaine et al. 2014;], extrusion rate [Coppola et al. 2009] and gas emission [Gouhier and Coppola 2011;Muro et al. 2014;Tulet and Villeneuve 2011] datasets. No technique or tool to directly measure the evolution of the collapse morphology and volume were available in situ at PdF. In this study, we present the first quantified observations of the geomorphological evolution of a volcanic caldera during its formation. Helped by recent developments in image matching software, we produced an unprecedented series of high-resolution Digital Elevation Models (DEMs) of the summit from aerial photographs taken during its collapse. These DEMs reveal the caldera volume evolution during the first fourteen days following the initiation of the collapse. With the addition of complementary aerial infrared data from the Observatoire Volcanologique du Piton de la Fournaise acquired about a year after the collapse, we bring new light to the collapse structure and dynamics.

Production of multi-temporal DEMs of the caldera
The complete 2007 aerial photography dataset used in this study comprises (1) seven aerial photography cam-paigns (pre-and syn-collapse on January 10 th , April 7 th [two sets], 10 th , 12 th , 17 th and 19 th ), (2) a postcollapse set of 52 pictures (A.I.G.L. ; taken April 12 th , 2008 and also used to georeference the other sets) and (3) a set of 130 post-collapse infrared images (May 5 th , 2008). Except in the case of the A.I.G.L. campaign (which used automated photography and a camera set beneath an airplane), the photographs were acquired manually by different operators (Table 1) from an ultralight airplane flying roughly between 200 and 600 meters above the caldera rim (these values are different from the flying altitude values given in Photoscan reports, which are computed from the average point cloud to camera distance). Photographs with sufficient overlap were selected for structure-from-motion (SfM) processing. SfM is a technique that enables 3dimensional digital models to be built of real objects from multiple-view photographs [e.g. Koenderink and van Doorn 1991]. We processed the 9 datasets with Agisoft Photoscan Pro (V.  (Table 1). Georeferencing followed three steps: (1) to georeference the April 2008 survey, the field coordinates of 39 ground control points (GCPs) were measured by GNSS, and used in the SfM processing.
(2) From the resulting model, positions of 25 different features were extracted. (3) The coordinates of these features were used to georeference the other surveys. All models were thus co-referenced so that the exported DEMs could be compared. Note that GCPs were not included as control measurements in the bundle adjustment, but only to scale, rotate and translate the models (i.e. there were no independent check points). Seven high-resolution (0.5 to 2 m) subsidence maps (or differential DEMs) during and after the collapse process were ultimately produced ( Figure 2).

Volumes and uncertainties
Caldera volumes at different time steps were computed from the differential DEMs, by integration of the values over the caldera surface. Due to the uncertainty in GCP locations on the photographs and to errors in the photogrammetric models, the final co-referencing of the models yielded a mean absolute position error (i.e. the mean distance, in a model, between input coordinates and virtual coordinates for GCPs) of 0.52 m to 3.01 m, mostly in the z component. Considering the 0.795 km 2 surface of the summit caldera, this implies that the caldera successive caldera volumes were estimated with a 2.4 × 10 6 m 3 maximum uncertainty.
"Acquisition Information Géo-Localisée" (former french private aerial photogrammetry company). [C] Timeline of the Piton de la Fournaise activity (eruptions and collapse) during the period of study. The dates of the aerial surveys used in this study are also reported.

Morphological evolution of the caldera from the multi-temporal DEMs
During the first 6 days of the collapse (April 5 th to 10 th ), the whole Cratère Dolomieu area subsided (Figure 2A-D). After April 10 th , some terraces to the south and to the east did not display further subsidence ( Figure 2E-G), but most of the caldera area continued to subside until at least April 19 th (e.g. 6 m 3 s −1 between April 17 th and April 19 th ). In some parts, the newly formed cliff became unstable, leading to rockfalls within the caldera walls (from 1 × 10 3 up to 809 × 10 3 m 3 ). These rockfalls occurred mainly in the northwest and west of the caldera. After April 19 th , only the northern half of the caldera ( Figure 2H) continued to subside at a very low rate (<6 m 3 s −1 ) until an unknown date. Between April 19 th , 2007 and April 12 th , 2008, the caldera volume had increased by 4.6 × 10 6 m 3 , and a few sec-tors in the northwestern cliff had detached and created deposits inside the caldera. There was a clear difference in behaviour between the northwest of the caldera (a quick collapse within the perimeter of the inner ring fault, not associated with terraces) and the south and east (a slower collapse associated with terraces and outer ring faults). These differences helped to map the surface trace of the faults at different stages of the caldera formation process (e.g. Figure 2I).

Temporal evolution of the caldera volume
The DEMs of the caldera showing the morphology at different stages of the April 2007 collapse allow us to measure a persistent subsidence several days after the onset of the event (Figure 3). After April 5 th and 6 th , the subsidence rate followed an exponential decay trend (i.e. slowing of the caldera volume increase as shown on Figure 3C). On April 7 th , 2007, at 2:40 UTC, or 30

Post-collapse thermal mapping
The thermal mapping on May 5 th , 2008, revealed further morphological clues relating to the caldera's structure. The spatial distribution of temperature inside the caldera showed a ring of thermal anomalies. This thermal ring, where temperature was between 20 and 39 • C greater than the surrounding rock average temperature (Figure 4), circled the northern half of the caldera. It formed a kidney-like shape (as seen from above) running along the cliff of the Cratère Dolomieu mostly between 2300 and 2340 m a.s.l. Furthermore, we noted three distinct secondary thermal zones, in the southwest and south of the caldera Figure 4C). This thermal belt and the associated secondary thermal zones were still observed in 2017 (e.g. Figure 4, inset of December 31 st , 2017).

Caldera volume and main caldera faults
Our new dataset enables us to quantify the caldera's final volume at 126.6 ± 2.1 × 10 6 m 3 , higher than, but in relatively good agreement with previous studies (96 × 10 6 m 3 from Aster stereo-imagery and 100-120 × 10 6 m 3 from field surveys Urai et al. [2007]). Lava emission volumes from the vent at low elevation were available at the time of the eruption, but were found not reliable enough to provide a scientifically relevant comparison with collapsed volume measurements (as they were derived from satellite radiance measurements instead of direct SfM 3D modeling; [Coppola et al. 2009]). Many causes of caldera collapse can be found in the literature. The two main models link the caldera to a failure of the structure lying above a magma reservoir as a consequence of (1) an overpressure in the reservoir leading to fractures in the roof [e.g. Gudmundsson et al. 1997] or (2) an underpressure following a large magma withdrawal [e.g. Roche et al. 2000;Staudacher et al. 2009]. In the case of the April 2007 collapse at PdF, numerous geophysical parameters suggest underpressure from magma withdrawal ]. Unusually high fountaining at the low-elevation vent and a strong increase in the eruptive tremor on April 5 th , a few hours before the onset of the collapse ], suggest unusually high magma extrusion rates prior to the collapse. Furthermore, GNSS permanent stations showed a significant deflation of the summit during the days preceding the collapse (April 2 nd to 5 th ; ]), strengthening the hypothesis of a large amount of magma being withdrawn from the shallow reservoir located below the summit. Calderas are generally considered to be delimited by a set of concentric faults (ring faults) linking the failure of the reservoir roof to the collapse structure at the surface [e.g. Acocella 2007;Gudmundsson 1988;Michon et al. 2009;Walker 1984]. These faults can be normal (inward-dipping), vertical or reverse (outward-dipping), depending on the models and location [Burchardt and Walter 2010;Gudmundsson et al. 2016;Roche et al. 2000]. They can be divided into   showed an initial ring fault with a sub-circular shape centred in the northwestern part of the caldera, associated later with outer ring faults to the south and east, which eventually define the final extent of the caldera With our new data, we were able to map the ring faults active during the collapse, and link them to the aforementioned broader collapse models. Figure 4 shows the spatial evolution of these faults and their associated terraces with a selection of time steps, deduced from the April 6 th aerial photograph ( Figure 4A) and the DEMs we produced for April 7 th and 12 th (Figure 2). On April 6 th , the IRF ( Figure 4A, B), centred in the north-western part of the caldera, was completed by a set of ORFs in the south and east, delimiting three main terraces (Et, St1 and St2). On April 7 th St1 had subsided and moved northward, giving the IRF a kidney shape. The other terraces had also moved towards the centre of the caldera ( Figure 4B) and the fault at the base of St2 had significantly propagated eastward. On April 12 th , two new terraces appeared in the southwest and all the structures had again moved towards the centre of the caldera. The final area of Cratère Dolomieu (0.788 km 2 , versus 0.729 km 2 before the col-Volcanica 3(1): 55 -65. doi: 1 .3 9 9/vol. 3. 1.5565 lapse) was defined on the April 12 th , 2008 DEM, after complete stabilization of the caldera floor ( Figure 5E).
The evolution of the caldera morphology was characterized by the development of a strong N-S and weak E-W asymmetry. The southern terraces moved sufficiently northward to cover the initial surface trace of the IRF. This caldera morphology is consistent with models proposed for collapses with high roof aspect ratios (i.e. upward propagation of vertical faults surrounding the reservoir, becoming slightly to significantly reverse (outward-dipping) close to the surface, [e.g. Holohan et al. 2011;Roche et al. 2000]), as would be expected at PdF (where the reservoir extent is thought to have a radius of about 500 m at around 2 km below the summit [e.g. Michon et al. 2009;Peltier et al. 2008]). After the end of the caldera formation process, its rim remained unstable and still produced hazardous rockfalls, especially in its northwestern quarter [details in ].

Active faults during the caldera formation and
their relationship with the post-collapse hydrothermal activity Based on our DEMs, we propose the following model for development of the shallow subsurface caldera structure ( Figure 5). First, following the propagation of the IRF to the former caldera floor, a collapse event was initiated on April 5 th (first collapse increment at 20:48 UTC, ). On April 6 th , the central structure, delimited by the reverse IRF, had subsided significantly (as proposed by previous studies based on geophysical parameters, it behaved as a coherent block [e.g. Michon et al. 2009;Staudacher et al. 2009]). This behavior is in good agreement with those observed in sandbox analogue models of this type of caldera formation (in the case of magmatic reservoirs with high roof aspect ratios [e.g. Acocella 2007;Holohan et al. 2011;Roche et al. 2000]). The ring fault network was completed in the south and east by a normal ORF [Michon et al. 2009]. On April 12 th , the normal ORF network increased in density so that the collapse extended further (in particular in the south). Debris filled the bottom of the caldera, so that the initial trace of the IRF was no longer exposed. We propose that a low-dipping reverse fault in the south, possibly as an extension of the southern IRF, was responsible for the significant northward displacement of St1 ( Figure 5C). The thermal mapping of the caldera's floor carried out on May 5 th , 2008, revealed a thermal ring closely matching the IRF in the west, north and east, and the base of St1 in the south ( Figure 4B, C). We thus propose that this thermal ring follows the IRF at depth. In March 2008, field observers reported that 70-80 % of the thermal ring displayed fumarolic activity [Staudacher 2010], so that the thermal anomaly seems linked to hydrothermal fluid circulation in the caldera. This is consistent with models detailed by Stix et al. [2003] and 4D mapping of caldera collapse at Piton de la Fournaise Derrien et al., 2020Garden et al. [2017, where flow of hydrothermal fluids preferentially follows pre-existing faults. The three other distinct secondary thermal zones ( Figure 4C) seem to be located on the southern ORFs, and close to the location of four older pit-crater events (1953, 1961, 1986 and 2002) at the summit of PdF ( Figure 5E). These secondary thermal zones are more spatially diffuse than the main thermal ring, therefore inter-granular fluid flow may happen inside intra-caldera breccia at these locations (a process also evidenced by Garden et al. [2017]). A well-developed hydrothermal system is inferred below the summit [e.g. Barde-Cabusson et al. 2012;Dumont et al. 2019;Lénat et al. 2012]. It thus seems likely that hydrothermal fluids preferentially use the April 2007 faults to reach the surface, and mainly the most damaged areas (i.e. the IRF and the location of previously collapsed structures). As suggested by previous authors [Chaput et al. 2019;Peltier et al. 2012], in the case of Piton de la Fournaise, this hydrothermal circulation would result in rock alteration and the reduction of cohesion and material strength at depth, leading to zones more prone to fail. These zones of weakness are clearly visible in resistivity models [Dumont et al. 2019;Lénat et al. 2000], and would therefore have been active over a long time and could destabilize regularly.

Conclusion
Our new multi-temporal Digital Elevation Models (DEMs) of the April 2007 summit collapse at Piton de la Fournaise allows us: 1. to assess the spatio-temporal evolution of the caldera throughout the event (both its volume increase with a 2.1 × 10 6 m 3 uncertainty and its morphology), making it possible to demonstrate that the collapse process lasted longer (at least 12 days) than previously thought; 2. to assess the precise location, geometry and role of the caldera ring faults. The surface trace of the initial ring fault was still mostly visible in 2019 in the west, south and east of the caldera, but was masked early in the collapse process by terraces and debris in the south; and 3. to identify a thermal anomaly that closely matches the caldera ring faults. Volcanic fluids preferentially used these faults to propagate from depth to the surface. Whether these fluids were directly transferred from the magmatic reservoir or just produced by very shallow hydrothermal circulation remains, however, to be determined.
Overall, very little dynamic geomorphological information is currently available for collapse processes on volcanoes, which are amongst the most spectacular and hazardous natural phenomena. This study provides multiple new insights (i.e. details of the temporal evolution of the caldera volume, precise mapping of faults active during the collapse and of thermal anomalies inside the caldera) into the formation process of the Piton de la Fournaise summit caldera, helping to further understand the relationship between faulting, hydrothermal circulation and collapse dynamics for this type of process worldwide.