# Unraveling pyroclastic density current dynamics with multiparameter geophysical sensing

Francesco Biagioli, Jean-Philippe Métaxian, Eléonore Stutzmann, Pascal Bernard, Giorgio Lacanna, Anne Mangeney, Giacomo Ulivieri, Gilda Risica, Alister Trabattoni, Marc Peruzzetto, et al.

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Francesco Biagioli, Jean-Philippe Métaxian, Eléonore Stutzmann, Pascal Bernard, Giorgio Lacanna, et al. Unraveling pyroclastic density current dynamics with multiparameter geophysical sensing. Communications Earth & Environment, 2026, 7 (1), pp.111. 10.1038/s43247-025-03091-6. hal-05495939

HAL Id: hal-05495939

https://brgm.hal.science/hal-05495939v1

Submitted on 5 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_3.png" /></figure>

# Unraveling pyroclastic density current dynamics with multiparameter geophysical sensing

<figure><img src="image_4.png" /></figure>

Check for updates

Francesco Biagioli<sup>1,2✉</sup>, Jean-Philippe Métaxian<sup>2,3</sup>, Eléonore Stutzmann<sup>2</sup>, Pascal Bernard<sup>2</sup>, Giorgio Lacanna<sup>1</sup>, Anne Mangeney<sup>2</sup>, Giacomo Ulivieri<sup>4</sup>, Gilda Risica<sup>5</sup>, Alister Trabattoni<sup>6</sup>, Marc Peruzzetto<sup>7</sup>, Lorenzo Innocenti<sup>1</sup> & Maurizio Ripepe<sup>1</sup>

Detecting, locating, and characterising the dynamics of destabilised volcanic material is critical for assessing the extreme hazards posed by volcanic mass flows, such as pyroclastic density currents. Geophysical measurements of these events may offer information otherwise hardly observable at close range. Here we investigate pyroclastic density current dynamics using a multiparameter approach that combines seismic, distributed acoustic sensing, and infrasound data with thermal and visible imagery, supported by numerical simulations. We focus on two events at Stromboli volcano, Italy, that occurred in October and December 2022. By comparing visible imagery with seismic energy and applying array processing techniques, we identify different flow volumes (~23.5 ± 9.5 × 10³ m³ and ~80 ± 9 × 10³ m³, respectively) and velocities (33-42 m/s and 54-59 m/s, respectively). Simulations reveal that reproducing these velocities requires volume-dependent empirical friction angles (~27° and 21°), consistent with dry granular flow behaviour and friction weakening. These findings offer new insights into the use of distributed acoustic sensing for volcanic monitoring and underscore the value of integrating multiparameter data with modeling to better understand complex volcanic processes.

Flank instability is common at volcanoes, affecting small to large portions of their edifices<fnref n="1" />. Remobilised volcanic material can trigger gravitational flows, such as rockfalls and landslides, rapidly reshaping volcanic morphologies and posing various hazards, including the eventual generation of tsunamis<fnref n="1" />,². Failure of crater walls or tephra may trigger pyroclastic density currents (PDCs) mainly composed of previously emplaced and hot volcanic material³⁻⁵. Similar to other types of PDCs caused by the collapse of eruptive columns and lava domes, or by lateral blasts, these events are fast-moving, high-temperature, gravity-driven flows of volcanic gas and particles³,⁶,⁷. PDCs represent the primary cause of fatalities from volcanoes at close range as they may propagate up to tens of kilometres from eruptive centres and exert dynamic pressures of up to 10² kPa⁶,⁸⁻¹⁰. The role of gas or other interstitial fluids, alongside factors such as pyroclast granulometry, may influence their flow dynamics, modifying friction, run-out distance, and dynamic pressure⁶,¹¹,¹². However, quantifying these effects at the field scale is still an open question with significant implications for PDC hazard assessment.

Detecting, locating, and estimating the velocity and size of these events are also crucial aspects of volcanic hazard assessment. Multiparameter geophysical measurements of PDC-generated signals (e.g., seismic and infrasonic waves) ensure real-time detection and provide information about flow dynamics, whose internal processes are hardly observable at close distances¹³. Seismic detection and characterisation of small-volume (<10⁵ m³) mass movements such as PDCs are increasingly adopted using high-frequency signals (>1 Hz)¹⁴⁻¹⁶. Amplitude-based methods applied to networks of seismometers around volcanoes enabled PDC tracking at several volcanic sites¹³,¹⁷,¹⁸. Despite their widespread use for locating seismo-volcanic sources, such as explosions and tremors¹⁹⁻²¹, array-based techniques (such as seismic beamforming²²) using small-aperture seismic arrays are rarely applied to monitor PDCs and gravitational flows¹⁴. Small-aperture seismic arrays are typically deployed for limited time periods, often as part of temporary field experiments²⁰,²¹,²³. Therefore, the analysis of PDCs is constrained by their low occurrence frequency. Additionally, other seismo-volcanic signals, such as tremor or explosions, may obscure the detection

<footnote n="1">¹Department of Earth Science, University of Florence, Firenze, Italy. ²Université Paris Cité, Institut de physique du globe de Paris, CNRS, Paris, France. ³CNRS, IRD, Université Gustave Eiffel, Université Grenoble Alpes, University of Savoy Mont Blanc, Chambéry, France. ⁴GECO Srl, Firenze, Italy. ⁵Istituto Nazionale di Geofisica e Vulcanologia, Sezione di Pisa, Pisa, Italy. ⁶Université Côte d'Azur, IRD, CNRS, Observatoire de la Côte d'Azur, Géoazur, Valbonne, France. ⁷BRGM, Orléans, France. ✉e-mail: francesco.biagioli@unifi.it</footnote>
and tracking of gravitational flows. Heterogeneous volcanic structures also
attenuate and scatter seismic waves, reducing coherence and distorting the
wavefield, which hinders analysis and interpretation. Although small-
aperture arrays deployed near craters could improve signal-to-noise ratios
and mitigate propagation effects, their proximity also increases vulnerability
to PDC damage.

Alongside classical seismic methodologies, the distributed acoustic sensing (DAS) technology may offer a powerful tool for long-term, high-resolution volcanic surveys. DAS turns fibre-optic cables into arrays of thousands of single-component strain sensors via optical interferometry of Rayleigh back-scattered light from laser pulses sent through the fibres²⁴. Encouraging results from recent applications²⁵⁻²⁸ suggest that integrating DAS with conventional seismo-acoustic instruments could enhance the detection and comprehension of volcanic phenomena²⁶. Nevertheless, its actual use for long-term volcanic monitoring is still under evaluation.

Infrasonic arrays also enable detection and location of PDCs¹⁷,²⁹ at volcanoes. Short-range infrasonic measurements are highly effective for studying volcanic activity, owing to low atmospheric scattering and dissipation and more predictable propagation and site effects¹⁴,²⁹. Like seismic networks, infrasonic arrays are not affected by cloud or ash coverage and do not need direct visual access to the events²⁹. These features make them especially useful when ash clouds from PDCs or eruptive activity obscure key visual information, such as the dynamics of the flow front or the amount of material released.

Although not always available or suitable for analysis, high-resolution
camera imagery provides crucial data for real-time monitoring and inves-
tigation of PDCs²,³⁰,³¹. Real-time images may provide unique information to
constrain the size of an event, helping to calibrate empirical laws between
seismic energy and PDC volume. Therefore, combining seismic and
infrasonic observations with camera imagery may help discriminate

between different events (e.g., explosions, rockfalls, and PDCs) based on
their unique seismic and acoustic signatures and, when possible, direct
visual confirmation. Additionally, the integration of multiparameter geo-
physical analyses of PDC signals enables quantitative estimations of their
velocity and volume, offering a complementary approach in cases where
direct visual observations are not possible.

Finally, numerical modelling of gravitational flows and PDCs serves as a powerful tool for characterising the underlying physical processes that are essential for accurate hazard assessment. These models rely on estimates of the released volumes and rheology, i.e., the friction law and its governing parameters³², and they take into account topographic effects on the flow dynamics. Constraining the rheology of gravitational flows, particularly the basal friction, remains a key challenge for accurately modelling these complex grain-fluid flows. This difficulty is mainly due to the scarcity of field data on flow dynamics. Yet, such rheological insight is essential for assessing the tsunamigenic potential of PDCs and other gravitational flows entering water bodies³,³³.

We show here how combining seismic, DAS, and infrasonic mea-
surements with video recordings and numerical simulations makes it pos-
sible to detect and characterise PDC dynamics and frictional behaviour. At
Stromboli volcano, Italy (Fig. 1, Supplementary Note 1), PDCs are related to
explosive activity² or to the collapse of hot tephra deposits⁴,³⁴. PDCs are
generally confined within the Sciara del Fuoco depression (SdF)³¹,³⁴⁻³⁸. SdF is
an unstable, ~1 km wide sector-collapse scar with a mean slope of ~35°
(Fig. 1), extending to ~2 km below sea level. PDCs (and other gravitational
flows) reaching the sea can trigger tsunami waves, endangering the island
and neighbouring coastlines²,³⁹. When not channelled in the SdF, PDCs can
affect inhabited centres, claiming victims³⁴,³⁵. In autumn 2022, intense
volcanic activity triggered lava overflows and flank instability, leading to the
failure of the northeast (NE) crater walls and generating PDCs⁴

<figure><img src="image_5.png" /><figcaption>**Fig. 1** | **Multiparameter monitoring network.** a Map of Stromboli Volcano showing the monitoring network of the University of Florence, Italy (UniFi), consisting of broadband seismo-acoustic stations (white triangles), infrasonic pressure sensors (white diamonds), tsunami gauges (PdC and PIB, blue circles), and the LBZ (visible) and ROC (thermal) cameras. The yellow line represents the 4 km long DAS cable exhibiting four sub-arrays (i.e., DNT, DCT, DST, DSG). The DAS interrogator is hosted in the UniFi Geophysical Observatory (OSV). Black triangles are seismic</figcaption></figure>

sensors (nodes) deployed between 19 September and 20 October 2022 in two small-aperture arrays (NNT, NSG). Blue and red stars, respectively, mark the location of PDC1 and PDC2 flow fronts obtained from LBZ georeferenced frames every 2 s (see Methods and Fig. 2). **b** DNT and NNT array geometries. **c** DCT array geometry. **d** DSG and NSG array geometries. **e** EAR array geometry. The Digital Earth Model is computed from images taken in 2014 and courtesy of the Italian Civil Protection.
(Supplementary Fig. 1, Supplementary Movies 1 and 2). We focus on two PDCs (that we refer to as PDC1 and PDC2), which occurred on 9 October and 4 December 2022, respectively. These were the most intense events we recorded (Fig. 2a, b).

Our dataset consists of DAS data measured along 4 km of geolocated optical fibre (Supplementary Note 2, Supplementary Fig. 2), two temporary small-aperture arrays of nodal seismometers, and the University of Florence (UniFi) permanent monitoring network (see Methods) which includes broadband seismometers, infrasonic pressure sensors, and thermal and visible cameras (Fig. 1). We analyse the multiparameter dataset to detect and track the PDCs and quantify their flow volumes and velocities. We also aim to test the efficacy of DAS for volcanic monitoring using the nodal seismometers to validate the DAS signals and as a reference to interpret DAS array processing results. First, we calibrate our DAS strain-rate measurements, finding good phase agreement with node data despite showing lower amplitudes (Supplementary Note 3, Supplementary Fig. 3). We then assess the effectiveness of seismic, DAS, and infrasonic measurements in detecting and tracking the PDCs, and compare the derived flow velocities with estimates from video recordings. Additionally, continuous analysis of seismic and DAS data allows for the detection of other gravitational flows along the SdF. DAS array processing also reveals intense activity along the SdF preceding the destabilisation of a large mass (PDC2). Finally, we estimate the flow volumes by combining thermal imagery with seismic data from broadband stations and DAS. Numerical modelling of these PDCs

demonstrates that the Coulomb friction angles to be used in the model are inversely related to the flow volume and similar to those of dry granular flows, suggesting comparable dynamics.

## Results

### Multiparameter PDC tracking

We characterise the dynamics of the two major PDCs that occurred in 2022 by combining camera images and multiparameter geophysical measurements. First, we constrain the timing, velocity, and location of both PDCs using georeferenced images (see Methods) recorded by the LBZ visible camera (Fig. 1). We then apply array processing (see Methods) to seismic (nodes and DAS) and infrasonic arrays to retrieve the back-azimuth of seismo-acoustic waves generated by volcanic activity and the PDCs. This approach allows for characterising the activity before, during, and after the two PDCs, as well as tracking their downslope propagation. The nodal and DAS arrays are sensitive to ground-coupled seismic sources, whereas the infrasonic array detects pressure waves propagating in the atmosphere. Therefore, each technique provides complementary insights into distinct aspects of the volcanic activity, featuring multiple yet simultaneous phenomena.

We focus on three representative arrays for each method: NNT (nodes), DCT (DAS), and EAR (infrasound; Fig. 1b, c, e; Supplementary Table S1). These arrays are able to track sources migrating along the SdF, offering an azimuthal shift of ~30–40° from craters to coastline. In

<figure><img src="image_6.png" /><figcaption>**Fig. 2** | Tracking the PDCs along the flanks with visible and thermal imagery. **a, b** respectively show two snapshots recorded by the LBZ camera during PDC1 and PDC2. The empty blue and red stars, respectively, indicate the approximate initial location of the masses that generated PDC1 and PDC2. ts₁ and ts₂ mark the approximate start times of PDC1 and PDC2, respectively. The blue and red stars show the position of the pyroclastic front in the visible images at 2-second intervals. tf₁ (07:23:08 UTC) and tf₂ (15:18:12 UTC) in **a, b**, respectively, indicate the moments when PDC1 and PDC2 first appear in the video frames. te₁ (07:23:30 UTC) and te₂ (15:18:24 UTC) indicate the arrival times at sea of the two PDCs. **c** Exhibits distance</figcaption></figure>

vs. time plot for PDC1 (blue) and PDC2 (red) with the corresponding average front velocities (41.6 and 59.2 m/s, respectively). These velocities are well reproduced by numerical simulations (see the Results section). The blue and red lines represent the distance along the SdF calculated at every time step from average front velocities and different friction angles δ. **d–g** Four thermal images from the ROC camera on 4 December 2022 showing the onset of PDC2. We observe the detachment and downslope movement of the mobilised mass (white dashed profile) that breaks away from the NE crater (NEC, labelled in **d**) and descends the SdF, triggering PDC2.
<table><caption>**Table 1 | Summary of main results**</caption><thead><tr><th></th><th>Monitoring cameras (continuous recording)</th><th>Nodal seismometers (19/09/2022- 20/10/2022)</th><th>Broadband Seismometers (continuous recording)</th><th>DAS (19/09/2022- 30/01/2023)</th><th>Infrasound (continuous recording)</th><th>Granular flow numerical model (SHALTOP)</th></tr></thead><tbody><tr><td>PDC1 (2022/10/09)</td><td>Detected</td><td>Detected</td><td>Detected</td><td>Detected</td><td>Not detected</td><td>-</td></tr><tr><td>Flow velocity [m/s]</td><td>41.6</td><td>33.0 ± 2.3</td><td>-</td><td>24.3 ± 2.5</td><td>-</td><td>-</td></tr><tr><td>Flow volume [m³]</td><td>-</td><td>-</td><td>17.6 ± 3.6 × 10³</td><td>27.6 ± 5.4 × 10³</td><td>-</td><td>-</td></tr><tr><td>Friction angle</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>~27°</td></tr><tr><td>Signal envelope</td><td>-</td><td>-</td><td>1.1 × 10⁻⁶ [m²/s]</td><td>3.7 × 10⁻¹⁰ [1/s]</td><td>-</td><td>-</td></tr><tr><td>PDC2 (2022/12/04)</td><td>Detected</td><td>No data</td><td>Detected</td><td>Detected</td><td>Detected</td><td>-</td></tr><tr><td>Flow velocity [m/s]</td><td>59.2</td><td>-</td><td>-</td><td>56.6 ± 5.2</td><td>54.1 ± 10.9</td><td>-</td></tr><tr><td>Flow volume [m³]</td><td>80.0 ± 9.0 × 10³</td><td>-</td><td>80.0 ± 9.0 × 10³ (from camera images)</td><td>80.0 ± 9.0 × 10³ (from camera images)</td><td>-</td><td>-</td></tr><tr><td>Friction angle</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>~21°</td></tr><tr><td>Signal envelope</td><td>-</td><td>-</td><td>4.5 × 10⁻⁶ [m²/s]</td><td>9.6 × 10⁻¹⁰ [1/s]</td><td>-</td><td>-</td></tr><tr><td>Small-scale gravitational flows after PDC1</td><td>Detected</td><td>Detected</td><td>-</td><td>Detected</td><td>Not detected</td><td>-</td></tr><tr><td>Precursory phase before PDC2</td><td>Detected</td><td>No data</td><td>Detected</td><td>Detected</td><td>Not detected</td><td>-</td></tr><tr><td>Steam explosions at sea after PDC2</td><td>Detected</td><td>No data</td><td>Not detected</td><td>Not detected</td><td>Detected</td><td>-</td></tr></tbody></table>

Supplementary Note 4 and Supplementary Fig. 4, we assess the accuracy of the three arrays by analysing an explosion that occurred at 07:05 UTC on 9 October 2022 (15 minutes before PDC1; Supplementary Fig. 3). The EAR infrasonic array provides the most precise result, locating the acoustic source of the explosion at the active craters. The NNT and DCT arrays place the event within the crater area, ~80–350 m (~5–20°) southeast of the exploding vent. A similar source location is obtained for both arrays during a 20-minute window of volcanic tremor between 07:00 and 07:20 UTC (Supplementary Fig. 5), confirming the validity of the results and the stability of the seismic source during ordinary activity. DAS and infrasonic arrays operated continuously during both PDCs. Due to battery depletion, nodes shut down 11 days after PDC1 (20 October 2022), highlighting a key limitation in the use of these instruments. When both operate simultaneously, we use nodal array results to interpret and validate those from the DAS array, as both detect seismic waves. Below, we describe the two PDCs and their multiparameter tracking, and present results on smaller-scale gravitational flows linked to flank instability before and after PDC1 and PDC2. The results are listed in Table 1.

**PDC1 on 9 October 2022.** The first event, PDC1 (Fig. 2a), occurred on 9 October 2022 at $t_s_1$ ~07:22:50 UTC, resulting from the collapse of part of the NE craters' outer rim⁴ (Supplementary Fig. 1, Supplementary Movie 1). The collapse was likely caused by a rise in magma level within the conduit, which increased the magmatic pressure and exerted additional thrust on the craters' rim, eventually leading to its failure⁴,³⁸. Georeferenced frames from the LBZ camera show PDC1 travelling 905 m down the SdF slope in 22 s (from $t_f_1$ = 07:23:08 UTC to $t_e_1$ = 07:23:30 UTC in Fig. 2a), with a front velocity of $u_f$ = 41.6 m/s (Fig. 2a, c, Supplementary Fig. 1d, Supplementary Movie 1). PDC1 was then followed by a lava overflow starting from the NE craters and propagating along the SdF⁴ (Supplementary Fig. 1c, d, Supplementary Movie 1).

Figure 3a–c show 6-minute-long waveforms, spectrograms, and backazimuths from seismic ground velocity (nodes), strain rate (DAS), and acoustic pressure (infrasound) data during PDC1. Before the event, both seismic ground velocity and strain-rate traces predominantly feature volcanic tremor (Fig. 3a, b; Supplementary Note 1). The associated backazimuths point towards the crater area: 237–252° for the nodal array and 248–260° for the DAS array (Supplementary Fig. 5). In both cases, the backazimuths locate the sources at the active craters or up to ~350 m (≤20°)

south of them (Supplementary Note 4). The uncertainties in back-azimuth estimates remain below ~2° (see Methods and Supplementary Note 4).

The seismic velocity and strain-rate traces of PDC1 are ~2 minutes long, spindle-shaped, and emerging from the volcanic tremor (Fig. 3a, b). Corresponding spectrograms exhibit energies mainly between 2 and 10 Hz (Fig. 3a, b, Supplementary Fig. 6). From a seismic analysis perspective, the nodal array appears more accurate and reliable than DAS during PDC1: it captures a clearer crater-to-sea back-azimuth migration, with timing and flow velocity more consistent with visual observations. During the flow, the nodal array reveals a marked ~40° shift in back-azimuth, indicating a migration of the source towards the shoreline (290–295°; Fig. 3a). By considering 10 different trajectories from the craters to the sea (Supplementary Fig. 7a), the backazimuth variation during PDC1 translates into an average flow velocity of 33.0 ± 3.3 m/s, which is ~10 m/s lower than the velocity derived from visual imagery (Fig. 2c). The nodal array also suggests that PDC1 reaches the sea ~17 seconds after $t_e_1$. This delay may be due to: (i) the 16-second causal window used in the array analysis, where each back-azimuth reflects the preceding signal (see Methods); (ii) the spatial extent of the seismic source⁴⁰; and (iii) ongoing material collapse along the SdF for over a minute after PDC1 ceases (Supplementary Movie 1). The DAS array also captures a seaward back-azimuth shift during PDC1. However, this variation starts ~20 s after $t_e_1$, during a phase of ongoing gravitational activity along the SdF following the main pulse of PDC1 (Supplementary Movie 1). Such a shift corresponds to an average flow speed of 24.3 ± 2.5 m/s (Supplementary Fig. 7b), consistent with slower residual movements rather than the main flow front.

Before PDC1, infrasonic detections indicate sources in the NE sector of the crater area (266–272°), consistent with visible explosive activity at the summit craters and with the array results obtained from nodes and DAS (Fig. 3c, Supplementary Movie 1). However, PDC1 does not appear in the infrasonic record, which is instead dominated by pressure transients occurring every 1–2 s, generated by vigorous activity at the NE craters (Fig. 3c, Supplementary Fig. 8). Few detections, between 07:23:40 and 07:23:44 UTC (~10 s after $t_e_1$), show back-azimuths pointing offshore (~320°), potentially indicating PDC1's impact with the sea. However, during the flow, we do not observe clear evidence of a source propagating along the SdF toward the sea. This may result from (i) an inefficient acoustic source, likely due to the relatively small scale of PDC1⁴¹, or (ii) intense activity at the NE craters, which dominates the acoustic wavefield and hinders detection of weaker signals (Supplementary Fig. 8).
<figure><img src="image_7.png" /></figure>

**Fig. 3** | PDC back-azimuth tracking with nodal, DAS and infrasonic arrays. The left and right columns respectively show 6-minute-long traces, spectrograms and corresponding back-azimuths of PDC1 and PDC2. **a**, **d** Show the seismic velocity (vertical component) recorded by the T3 node (**a**) and the ROC seismometer (**d**) in black. **b**, **e** is the DAS strain-rate trace at the offset 900 m. The black vertical lines in (**b**, **e**) indicate the scale of the individual strain-rate traces. The four vertical red lines in (**b**, **e**) numbered from 1 to 4 indicate the four segments of DCT, whose geometry is displayed in (**h**). **c**, **f** Show the acoustic pressure recorded by the MIC5 infrasonic sensor in black. Back-azimuths are calculated using the NNT, DCT, and EAR arrays. Nodal arrays are excluded from the analysis of PDC2 due to battery depletion in October 2022. The green and pink vertical dashed lines indicate the onsets

($t_1$ = 07:22:50 and $t_2$ = 15:17:50 UTC, in green) and arrival at sea ($t_e$ = 07:23:30 and $t_e$ = 15:18:24 UTC, in pink) of PDC1 and PDC2, respectively, as observed by the thermal and visible cameras (Supplementary Fig. 1, Supplementary Movies 1 and 2). The black dotted line in (**a**–**c**) indicates the end ($T_{e1}$) of activity along the SdF linked with PDC1. The black horizontal lines at the bottom right in (**a**–**c**) indicate the length of the time windows for array analysis (16, 20 and 5 s, respectively). Back-azimuths obtained with NNT and DCT array in (**a**, **b**, **e**) are smoothed with a 2 s-long window. For clarity, back-azimuths in (**e**) are shown between 210 and 310°, for a comprehensive view of the results see Supplementary Fig. 10. Back-azimuths obtained with EAR array in (**c**, **f**) are filtered for time residuals <0.05 s. **g**–**i** Show the geometries of the three arrays.

**PDC2 on 4 December 2022.** The second event, PDC2 (Fig. 2b), occurred on 4 December 2022 at ~15:18 UTC during a lava overflow and a period of intense flank instability along the SdF that began around 14:25 UTC, nearly an hour earlier (Fig. 1e, Supplementary Movies 2 and 3). Thermal imagery reveals that PDC2 originated around $t_{s2}$ ~15:17:50 UTC from the detachment, sliding, and fragmentation of a ~48 m long, ~35 m high sector of the NE crater walls (Fig. 2d, Supplementary Fig. 1, Supplementary Movie 2). As with PDC1, PDC2 was likely triggered by a magma-level rise that heightened explosive activity and pressure on the NE crater rim¹. However, the collapse generating PDC2 appears substantially larger than PDC1 (Supplementary Movie 2). The hot, dislodged

mass of PDC2 rapidly descended along the SdF, covering the final 719 m in 12 s (from $t_f_2$ = 15:18:12 UTC in Fig. 2b) with a front velocity $u_f$ = 59.2 m/s (Fig. 1b, c). PDC2 entered the sea at $t_e_2$ = 15:18:24 UTC (Fig. 2b), generating a tsunami wave of 1.2–1.4 m recorded at two gauges ~1.3 km away⁴² (Supplementary Fig. 9). The interaction between hot materials and seawater also triggered hydrothermal explosions near the coastline, as observed from camera images (Supplementary Movie 2). After PDC2, flank instability persisted for ~2 hr and was accompanied by a lava flow.

Seismic ground velocity and strain-rate waveforms generated by PDC2 are similar to those of PDC1 while exhibiting higher amplitudes (Fig. 3d, e,
**Fig. 4** | Seismic and DAS detection of granular flows. **a** Seismic velocity trace (node T3, Z component) during the same period. **b** Back-azimuth results obtained from NNT array. **c** DAS strain-rate traces (offset 900 m) between 07:00 and 08:00 UTC on 9 October 2022. **d** Back-azimuth results obtained with the DCT array. Back-azimuth values in (**b**, **d**) are smoothed with a 3 s sliding window. **e–g** Snapshots of the LBZ visible camera during the three events (Evs. 1–3). Event 1 corresponds to PDC1, while Events 2 and 3 are rockfalls along the SdF (Supplementary Movie 4).

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

Supplementary Fig. 6). During PDC2, both DAS and infrasound data reveal crater-to-sea back-azimuth shifts of ~30°, consistent with the main flow and indicating its downslope migration (Fig. 3e, f). The back-azimuth shift observed by the DAS array (Fig. 3e) indicates an average PDC2 velocity of ~56.6 ± 5.2 m/s (Supplementary Fig. 7c), again consistent with visual estimates (Fig. 2c). However, DAS results show great variability, with multiple back-azimuth fluctuations before and after PDC2 (Fig. 3e). These variations likely reflect repeated gravitational collapses along the SdF, as confirmed by video footage (Supplementary Movies 2 and 3).

Additionally, between ~15:18:30 and 15:18:57 UTC, strain rate (and seismic ground velocity) amplitudes increase by a factor of three compared to those during the main phase of PDC2 and exhibit lower frequencies, mainly below 2 Hz (Fig. 3d, e). During this interval, back-azimuth estimates show high uncertainties (up to 30°) and shift towards the southern flank of the volcano (190°; Fig. 3e, Supplementary Fig. 10). Apparent seismic velocities, typically stable around ~0.8 km/s, rise sharply to ~10 km/s, with associated errors reaching ~3 km/s (Supplementary Fig. 10). Such an abrupt velocity increase is absent from DAS or nodal array data during PDC1 or earlier phases. Given the dominant frequencies of 3.5–4 Hz (see Methods and Supplementary Note 4; Supplementary Fig. 4g), the resulting wavelengths (> 2 km) likely exceed the resolution capacity of our array processing technique and array aperture during this activity phase (see Methods). We interpret these unusually high apparent velocities as indicative of seismic waves arriving at steep incidence angles, possibly originating from deeper regions within the volcano²¹. This deeper source may be due to the detachment of PDC2 and the subsequent decompression of the volcanic edifice due to the mass released and the overflowing lava. While compelling, a thorough investigation of this phenomenon lies beyond the scope of this study and would require further research, as DAS array analysis alone is insufficient.

Infrasonic signals of PDC2 are emergent and have amplitudes up to ~5 Pa (at MC5 sensor) at frequencies <8 Hz (Fig. 3f). The infrasonic array precisely identifies the volcanic explosion at 15:16:50 UTC at the NE craters, approximately one minute before PDC2. The explosion is then followed by a marked ~30° back-azimuth shift starting at ~15:18:05 UTC. This variation tracks the propagation of PDC2 along the SdF (Fig. 3f), corresponding to an

average flow speed of 54.1 ± 10.9 m/s (Supplementary Fig. 7d), in agreement with visual observations (Fig. 2c). After PDC2, the infrasonic record is characterised by signals having higher acoustic pressure (10 Pa), and frequencies mainly between 0.2–1 Hz (Fig. 3f). During this phase, the acoustic source remains offshore until ~15:27:30 UTC (Fig. 3f), likely reflecting prolonged hydrothermal explosions where hot materials composing PDC2 and the lava front meet seawater (Supplementary Movies 2 and 3).

**Other episodes of flank instability.** Recurrent gravitational flows occurred along the SdF before and after the two PDCs, as confirmed by monitoring images (Supplementary Movies 2–4). Judging from visual observations and seismic records, these events involved smaller volumes than PDC1 and PDC2. Extending the array analysis of seismic and DAS signals beyond the main PDC windows (as in Fig. 3) reveals these additional flows.

Figure 4 shows back-azimuths from nodal and DAS array analysis between 07:00 and 08:00 UTC on 9 October 2022. For most of the hour, both arrays (NNT and DCT) point to a stable seismic source close to the summit craters, likely due to volcanic tremor and explosions. However, during three distinct events, we observe back-azimuths shifting and reaching the coastline: these three phases coincide with PDC1 (event 1 in Fig. 4b, c), and two smaller flows around 07:38 and 07:45 UTC (events 2 and 3). Camera images confirm these events as rockfalls that reach the sea (Fig. 4e–g, Supplementary Movie 4). Nodal and DAS back-azimuths are generally consistent. However, nodes offer higher precision in locating tremor and explosions, better differentiation of source types, and reduced variability in back-azimuth measurements (Fig. 4). Nevertheless, this agreement reinforces confidence in DAS outcomes during periods without node coverage, such as during PDC2. During the same period in Fig. 4, the EAR infrasonic array detects no such activity along the SdF. Infrasonic back-azimuths consistently indicate a persistent acoustic source at the NE craters (Supplementary Fig. 11).

On 4 December 2022, the formation of PDC2 at 15:18 UTC was preceded by intense and rapidly evolving slope activity. To evaluate whether DAS alone can capture this precursory phase, we analyse strain-rate data, spectrograms, and back-azimuths from 13:50 to 15:40 UTC, ~1 hour and 20
**Fig. 5 | DAS monitoring of flank instability before PDC2.** a Cumulative strain-rate time series along offsets 890–910 between 13:50 and 15:40 UTC on 4 December 2022. The three colours represent three distinct phases of increasing volcanic activity observed from video (Supplementary Movie 3), spectrogram and back-azimuth results.

b Spectrogram of (a). c Back-azimuth results obtained with the DCT array. The grey and blue horizontal bands in c represent the back-azimuths of the craters and the shoreline, respectively. The black dashed vertical lines in (b, c) represent the two times when a change in activity is observed (i.e., ~14:25 and ~15:18 UTC). d–f Three snapshots of the LBZ camera during the three phases of activity.

<figure><img src="image_10.png" /></figure>

minutes before and 20 minutes after PDC2 (Fig. 5, Supplementary Movie 3). In Fig. 5, back-azimuths from the DCT array reveal a sustained phase of gravitational instability culminating in the detachment of the large PDC2 mass. During the initial phase (Phase 1 in Fig. 5, between 13:50–14:25 UTC), strain-rate amplitudes remain moderate, with dominant frequencies between 1–5 Hz, and stable back-azimuths (235–265°) pointing to sources near the crater area. During this interval, visible camera imagery shows no evident activity along the flanks (Fig. 5d). From ~14:25 UTC (Phase 2), around 50 minutes prior to PDC2, cameras detect repeated gravitational flows along the SdF (Fig. 5e, Supplementary Movie 3). Consequently, strain-rate amplitudes and dominant frequencies increase (up to ~15 Hz), waveform coherence declines (from ~0.8 to 0.6), and back-azimuths fluctuate along the SdF, consistent with frequent landslide activity. This gravitational activity is not detected by the EAR infrasonic array, which instead locates acoustic sources at the NE craters during this phase (Supplementary Fig. 12). Following the onset of PDC2 at 15:18 UTC (Phase 3), strain-rate amplitudes surge across all frequency bands. Back-azimuths repeatedly shift towards the coastline, indicating intense gravitational activity over the subsequent 15 minutes, as corroborated by visual observations (Fig. 5f, Supplementary Movies 2 and 3).

## PDC volumes

The two PDCs exhibit different downslope velocities and appear to involve different volumes. Flow volume is essential for numerical simulations and extremely useful for the definition of the possibly associated tsunami risk. We estimate the volume from camera images and seismic signals of PDC1 and PDC2 recorded with broadband seismometers and DAS. Thermal images show PDC2 originating from a ~48 m long, ~35 m high sector of the NE crater rim, which detaches and slides downward (Fig. 2d, Supplementary Fig. 1, Supplementary Movie 1). A paraboloid geometry with a height of ~35 m, a minor axis of 48 m, and a major axis of 90 ± 10 m (corresponding to the NE crater diameter; Supplementary Fig. 13), superimposed on the SdF topography, yields an estimated volume for PDC2 of $V_2 \approx 80 \pm 9 \times 10^3 \text{ m}^3$ (Fig. 6d; see Methods).

As there are no clear thermal or visible images showing the PDC1 mass, as in the case of PDC2, we estimate its volume $V_1$ on the assumption that the

PDC volume $V$ is proportional to the seismic energy $E_s$ released during the flow¹⁵:

$$V = \frac{3E_s}{R_{s/p} \rho_f g L (\text{tany } \cos\theta - \sin\theta)}, \quad (1)$$

where $R_{s/p} = E_s/\Delta E_p$ is the seismic (E_s) to potential (ΔE_p) energy ratio, $\rho_f$ the flow density, g the gravitational acceleration, L the slope length, $\theta$ the mean slope angle, and $\gamma$ the mean angle of the deposit¹⁵. Assuming identical slope conditions (constant angle $\theta$ and length L), similar flow materials, and a constant seismic-to-potential energy ratio $R_{s/p}^{16}$, we retrieve a linear relationship between flow volume and seismic energy¹⁵,⁴³. Accordingly, we estimate the volume of PDC1 using the known volume $V_2$ of PDC2 and the ratio of their seismic energies $R_E = E_{s2}/E_{s1}$, such that $V_1 = V_2/R_E$.

We calculate the seismic energies by integrating the envelopes of the signals generated by PDC1 and PDC2 (Supplementary Fig. 6), as recorded by broadband seismic stations and DAS (see Methods). We adopt window lengths of 44 s for PDC1 and 38 s for PDC2, reflecting their crater-to-sea propagation time (Fig. 2; see Methods). At the permanent broadband stations, the seismic energies between 1–15 Hz are $E_{s1} = 0.30 \pm 0.16 \text{ MJ}$ for PDC1 and $E_{s2} = 1.20 \pm 0.62 \text{ MJ}$ for PDC2, consistent with values reported for small rockfalls and PDCs in Montserrat¹⁶. The average energy ratio between the two events is $R_E = 4.1$, yielding a PDC1 volume of $V_1 = 19.5 \pm 2.2 \times 10^3 \text{ m}^3$. DAS strain-rate data in the same frequency band give a lower ratio of $R_E = 2.7$ (Supplementary Fig. 14), corresponding to $V_1 = 29.6 \pm 3.3 \times 10^3 \text{ m}^3$.

Since seismic energy (and thus the energy ratio and estimated volume) depends on the selected time window (see Methods), we tested an extended window that starts 10 s before the PDC onset ($t_1$ and $t_2$) and continues 10 s after the arrival at sea ($t_e$ and $t_2$). The window stops before the high-energy phase following PDC2 (Supplementary Fig. 15). This approach yields higher $R_E$ values, i.e., smaller PDC1 volumes, by ~4–5 × $10^3 \text{ m}^3$ (see Methods). We finally determine the volume of PDC1 $V_1 = 23.5 \pm 9.5 \times 10^3 \text{ m}^3$ by averaging the results from seismic stations and DAS across different window lengths, and accounting for uncertainties in PDC2 volume (Fig. 6a). This volume is ~3.4 times smaller than that of PDC2.
<figure><img src="image_11.png" /><figcaption>**Fig. 6** | Numerical simulations of the two PDCs. Snapshots of the simulations that best match the observed front velocities of PDC1 (a–c) and PDC2 (d–f). We consider friction angles $\delta_1 = 27^\circ$ and $\delta_2 = 21^\circ$, respectively. **a**, **d** show the initial location of the</figcaption></figure>

flow mass with colours representing its thickness. **b**, **c**, **e**, **f** show the propagation of PDC1 and PDC2 at different time intervals, compared with real observations. The colours in (a, b) and (e, f) indicate the flow velocity.

## Numerical simulations of the PDCs

From geophysical measurement and camera imagery of PDC1 and PDC2 we estimate their velocities (~33–42 m/s and ~54–59 m/s, respectively) and volumes (~23.5 ± 9.5 × 10³ m³ and 80 ± 9 × 10³ m³, respectively). We then use these velocities and volumes to constrain numerical flow simulations using the SHALTOP model44–47 to investigate the physical processes controlling PDC dynamics. SHALTOP simulates dry granular flows over complex topography. In its simplest form, all energy dissipative effects are controlled by a basal Coulomb friction force, with a constant friction coefficient $\mu = \tan(\delta)$ (with $\delta$ the friction angle, see Methods). Using the estimated PDC volumes, we aim to find the basal friction coefficient that allows us to reproduce both the front velocity and the duration of the two PDCs.

The simulated PDC1 fits visual observations by considering a friction angle $\delta_1 = 27^\circ$ (Fig. 6a, Supplementary Fig. 16). Such a coefficient yields a flow front velocity of ~43.2 m/s in the lower part of the SdF (Fig. 2c). This velocity is comparable to the value of 41.6 m/s obtained from the visible camera (Fig. 2c), and is broadly consistent with the average flow velocity estimated from geophysical measurements (~33 m/s; Fig. 3, Table 1). The simulated PDC1 reaches the shoreline in 41 s, consistent with visual observations of ~40 s.

The friction angle that best matches the velocity and timing of PDC2 is instead $\delta_2 = 21^\circ$ (Fig. 6d, Supplementary Fig. 16). Simulated PDC2 reaches the shoreline in 32 s (close to the observed flow duration of ~34 s) with a final front velocity of ~58.4 m/s. This velocity matches estimates from imagery (59.2 m/s; Fig. 2c) and geophysical data (~54–57 m/s; Fig. 3, Table 1).

## Discussion

In autumn 2022, several PDCs at Stromboli posed significant hazards and offered a unique case study. We combine seismic, DAS, and infrasonic signals from PDCs, with imagery from permanent monitoring cameras and numerical simulations. Multiparameter analysis using seismic, DAS, and infrasonic arrays tracks the seismo-acoustic source migration of the PDCs from the craters to the sea, which allows for the estimate propagation velocities ranging from ~33 m/s (PDC1) to ~57 m/s (PDC2), increasing with flow volumes. These velocities align with independent estimates from

camera imagery (41.6 and 59.2 m/s, respectively) and are also consistent with those of previous PDCs at Stromboli (~45–50 m/s)2,31,36 and other volcanoes13,17,18,29,30. Differences in velocity estimates (of ~5–10 m/s) between geophysical data and monitoring cameras likely reflect that the former yield average flow velocities along the slope, whereas the latter provide direct estimates of terminal front velocity.

The three arrays (nodes, DAS, and infrasonic) provide coherent yet different (and complementary) results in tracking the PDCs, depending on the flows' characteristics and their seismo-acoustic source mechanisms. Infrasonic array tracking depends on the coupling between the flow and the atmosphere30. However, the infrasonic source mechanisms of mass movements remain poorly understood, though frequency and amplitude generally seem correlated with the volume and thickness of the sliding body14,48,49. The large volume and high flow velocity (up to ~59 m/s) of PDC2 is thus associated with clear infrasonic signals showing a distinct back-azimuth shift linked with the seaward PDC propagation. PDC1, by contrast, produced weak or undetectable infrasound, likely masked by the intense crater activity (Supplementary Fig. 8, Supplementary Movie 1). This is probably due to lower flow volume and velocity. Such a pattern mirrors that of debris flows, where lower velocities and shallower flow depths yield weaker acoustic emissions41.

When available, the nodal array outperforms DAS in tracking PDCs through seismic signals. Our findings demonstrate the effectiveness of small-aperture arrays in tracking fast-moving volcanic flows and highlight a valuable opportunity to enhance future monitoring strategies. Nevertheless, we must clarify that, unlike the DAS system or the permanent monitoring network, the nodal arrays lack telemetry capabilities, meaning that we can analyse their data only a posteriori. Real-time analysis is possible with seismic arrays composed of autonomous, synchronised, and telemetered stations. However, such configurations are more complex and costly to deploy and maintain. The interest of DAS is that, while offering less precision in resolving individual event dynamics, it effectively detects gravitational flows and changes in volcanic activity linked to increasing flank instability, consistent with node observations. We point out that, unlike nodes, we continuously recorded DAS measurements over 4 months (from September 2022 to January 2023) with real-time access and minimal data
loss. This demonstrates DAS reliability for long-term monitoring and its potential for real-time application. Nevertheless, further strategies to store large amounts of data and to analyse them in short time periods need to be undertaken$^{50,51}$.

DAS signals carry information about the precursory flank activity preceding the destabilisation of the large PDC2 mass. Such precursory signals are absent from the infrasonic array, which continues to locate acoustic sources at the NE craters (Supplementary Fig. 12). As with PDC1, the small size of the gravitational flows may produce high-frequency seismic radiation (e.g., from clast impacts and friction$^{14,40}$) detected by DAS but not infrasonic pressure waves$^41$. Conversely, the infrasonic array detects the steam explosions at sea following PDC2 (Fig. 3f), which are absent in the DAS array results (Fig. 3e). These findings highlight the need to combine both seismic and infrasound measurements for a comprehensive understanding of different and simultaneous volcanic phenomena.

Our DAS array technique involves analysing individual fibre segments, as done in previous studies$^{7,52}$. Obtaining reliable back-azimuth measurements requires manual tuning of the DAS segment selection. We observe slight variations in back-azimuth values depending on the specific configuration and orientation of the selected segments (Supplementary Fig. 17). Additionally, array analysis of DAS data may be influenced by local heterogeneities$^{53}$ and directivity effects. We suggest developing more robust methods capable of leveraging the entire cable$^{54}$ along with innovative strategies to reduce the computational resources required to analyze the large volumes of DAS data. De-noising techniques$^{55,56}$ could also enhance results in sections with low signal-to-noise ratios due to optical signal attenuation (Supplementary Fig. 18).

Integrating seismic and DAS data with monitoring images reveals two distinct volumes for the two PDCs: PDC1 has a smaller volume ($\sim 23.5 \pm 9.5 \times 10^3 \text{ m}^3$) compared to PDC2 ($\sim 80 \pm 9 \times 10^3 \text{ m}^3$). These values are broadly consistent with independent estimates$^{38}$ for the same PDCs, namely $6.4 \pm 1 \times 10^4 \text{ m}^3$ for PDC1 and $88.9 \pm 26.7 \times 10^3 \text{ m}^3$ for PDC2. Discrepancies in these estimates, particularly for PDC1, may originate from assuming identical flow material for both PDCs. Variations in material properties can influence the seismic energy generated by each flow. Additionally, using shorter analysis windows than the full flow duration (from onset to final deposition) may lead to overestimating the seismic energy and the volume of PDC1 with respect to PDC2. Finally, volume estimates based solely on seismic energy also carry uncertainties of up to a factor of $2-3^{15}$. These uncertainties could be reduced by refining the relationship between PDC volume and seismic energy through additional, optically documented events. Nevertheless, the estimated PDC volumes align with those of PDCs that generated tsunamis of similar amplitude in summer 2019 and spring 2021$^2$, i.e., $\sim 70-200 \times 10^3 \text{ m}^3$.

Numerical simulations of dry granular flows$^{44-47}$ relying on the estimated volumes reproduce well the velocity and timing of the two PDCs. The friction angles of the two PDCs ($\sim 27^\circ$ and $21^\circ$ for PDC1 and PDC2, respectively) are inversely related to the flow volumes and consistent with other estimates (of $\sim 20^\circ$) for PDCs of similar size at Stromboli$^2$. The values of friction angles resemble those associated with landslides, implying similar dry flow dynamics$^{32}$. They also agree with previous studies on PDCs and landslides, which show that reproducing the mobility of large landslides requires lower effective friction angles than those derived from laboratory experiments$^{11,16,32,57-60}$. The origin of the friction weakening for increasing volumes (and velocities) remains an open question, though mechanisms such as flash heating, acoustic lubrication, fluid pore pressure effects, and the influence of fine particles have been suggested$^{32,61-63}$. Despite existing uncertainties, the good agreement between modelled and observed flow velocities and timing validates the use of constant Coulomb friction angles in thin-layer numerical models (SHALTOP) and shows the potential of this approach for volcanic hazard assessment.

Our results demonstrate that integrating seismo-acoustic monitoring with visual observations provides valuable insights into PDC dynamics and associated pre- and post-collapse processes. This approach lays the foundation for (i) numerical simulations that can uniquely constrain PDC

behaviour and (ii) near-real-time estimates of gravitational flow volumes from continuous seismic or DAS monitoring. It also enables independent assessment of flank erosion rates over time. Although further calibration is required to refine flow parameters derived from geophysical observations, such data provide crucial insights into PDC behaviour, especially when imagery is unavailable (e.g., due to cloud cover) or insufficient (as in the case of PDC1 volume). Detection, tracking, and volume estimates with seismo-acoustic sensors may also offer an independent means of evaluating the tsunamigenic potential of PDCs entering water bodies, since flow volume and velocity are key factors controlling tsunami generation and wave heights$^{54-66}$. Despite their higher uncertainties, geophysical measurements provide valuable additional constraints on flow parameters. These can support further calibration of PDC frictional behaviour and potentially improve modelling of PDC propagation and resulting tsunami amplitudes. Therefore, these findings underscore the importance of integrating multi-parameter datasets with physical simulations to improve volcano monitoring and deepen understanding of event dynamics, with significant implications for hazard assessment.

## Methods

### Stromboli monitoring network and DAS experiment

This study integrates multiparameter data from the UniFi permanent monitoring network, a fibre-optic DAS array deployed on the volcanic flanks, and 22 temporary nodal seismometers (Fig. 1). The UniFi network includes broadband seismometers, infrasound sensors, tiltmeters, tsunami gauges, and thermal and visible cameras$^{67}$. Broadband seismic stations (STR, ROC, SCI) use Guralp digitisers with CMG-40T 30 s seismometers, sampling at 100 Hz. They also feature differential acoustic pressure transducers with a 0.01–100 Hz flat frequency response, a 250 Pa full-scale range, and a 100 Hz sampling rate. Two stand-alone acoustic sensors (MC1, MC5) complement the STR station to form the EAR infrasound array (Fig. 3i). The visible camera (LBZ) is a FOSCAM model with a 4 mm focal length, $70^\circ$ horizontal and $75^\circ$ diagonal fields of view, sampling every 2 s. The thermal camera (ROC) is an FLIR-A310 with a $25^\circ \times 18.\,8^\circ$ optical lens, 0.05 $^\circ C$ thermal sensitivity, and 2 Hz sampling rate.

DAS data were recorded along a 4 km fibre-optic cable buried 20–30 cm deep on the volcano’s NE flank. The cable comprises a polyethylene jacket, a loose Kevlar layer, and eight tight-buffered single-mode fibres. A FEBUS A1-R DAS interrogator at the UniFi Geophysical Observatory (OSV) operated with a gauge length $L_G = 10 \text{ m}, 2 \text{ m}$ channel spacing, and recorded strain-rate data at 50 Hz. The cable forms several sub-arrays: three triangular (DNT, DCT, and DST) and one grid-like (DSG) with four $\sim 30 \text{ m}$ squares (Fig. 1, insets). DAS data were acquired from 19 September 2022 to 30 January 2023. DAS channels are named by their distance from fibre start (OSV). We performed a geolocation of the DAS channels via tap-testing and differential GPS tracking of the cable layout (Supplementary Fig. 2, Supplementary Note 2). The calibration of DAS strain-rate wave-forms with co-located nodes shows good phase correlation between the two datasets despite showing lower amplitudes recorded by DAS$^{27}$ (Supplementary Fig. 3, Supplementary Note 3).

The 22 SmartSolo IGU-16HR three-component (3C), short-period (5 Hz) nodal seismometers were co-located with the optical cable in two arrays: NNT (nodes T1-7) and NSG (nodes G1-15), arranged in triangular and grid geometries, respectively (Fig. 1a). The nodes were used to monitor volcanic activity and calibrate DAS data (Supplementary Note 3). The nodes operated from 19 September 2022 to 20 October 2022, sampling at 250 Hz. Node data are corrected for the instrumental response. DAS and nodal array names follow a systematic convention: the first letter indicates the instrument (D for DAS, N for nodes), the second denotes the location (N for North, C for Central, S for South), and the final letter refers to the geometry (T for Triangle, G for Grid).

### Camera analysis

Visible images are recorded every 2 s (0.5 Hz of sampling rate) by the LBZ camera (Fig. 1) with a resolution of $640 \times 480$ pixels. The LBZ camera is
located at an altitude of ~190 m a.s.l. and has a view direction of 215°. We use these parameters to georeference each frame with the IMGRAFT software⁶⁸ via a high-resolution (0.5 m) DEM resampled at 5 m (Fig. 1). When visible, we select the pyroclastic front in each frame. The selected pixels are then converted into 3-D UTM coordinates (Fig. 1, Supplementary Fig. 1).

## Seismic and DAS array processing

Seismic signals generated by volcanic tremor, explosions, and PDCs have no marked P and S wave arrivals, challenging the analysis with traditional source-location methods²³. Conversely, array techniques allow for the identification of both back-azimuth and apparent velocity of seismic waves generated by these sources¹⁴,⁶⁹⁻⁷¹. To track the downslope movement of the seismic sources of PDCs, we use array processing techniques that rely on delay times between sensors, assuming a unique and dominant plane wave in each time window.

The technique we use to analyse node and DAS data consists of estimating time delays between all possible N(N-1)/2 pairs of N sensors composing an array. Time delays are computed on subsequent sliding time windows with a coherence-weighted linear fit of the cross-spectrum phase between the two traces in each sensor pair²³,⁷². Estimated time delays have an accuracy below half the sampling rate, making this method suitable for small-aperture arrays with respect to the wavelengths of interest²³. For a unique plane wave, time delays $\tau_{ij}$ between sensors i and j are:

$$ \tau_{ij} = s \cdot r_{ij}, \qquad (2) $$

with s the slowness vector ($s = (-||s||\sin\theta, -||s||\cos\theta)$, $\theta$ the back-azimuth and $r_{ij}$ the relative position vector between the sensors. Back-azimuth and apparent velocity for each time window are calculated via inversion of N(N-1)/2 Eqs. (2) in a least-squares approach. The technique also provides uncertainties in the estimated back-azimuth and apparent velocity. We apply the technique on the vertical (Z) seismic signals recorded at the NNT array (Fig. 1) having an aperture of ~129 m.

Array analysis of DAS strain-rate signals is not as straightforward as in the case of seismometers. DAS data exhibit a direction-dependent sensitivity and an averaging effect of the gauge length⁵⁴. Furthermore, DAS is more sensitive to slower, sub-horizontal waves⁵²,⁷³. The technique we propose here has already been applied to DAS records of volcanic tremor and explosions at Stromboli²⁷. Here, we focus on one DAS array, DCT. It exhibits more favourable location to analyse sources along the SdF, showing a craters-to-sea azimuth shift of 40°. It consists of four 30-m-long linear (±10°) segments with different orientations and an aperture of ~100 m. (Figs. 1 and 3h). Estimating time delays requires high waveform coherence. We observe a sharp coherence drop between DAS channels at distances >30m and having different orientations²⁷. Therefore, we estimate the delay times separately for each segment. Under the assumption of a dominant wave type, we then obtain back-azimuth and apparent velocity by inverting all delay times together.

We apply the technique to node and DAS signals between 1 and 5 Hz. In this frequency band, traces exhibit high energy and coherence (Fig. 3, Supplementary Fig. 19). We use causal time windows, respectively, of 16 and 20 s for nodes and DAS and sliding every 1 s. Each back-azimuth value is calculated based on the 16 or 20 s of signals preceding its associated time-stamp. We selected these windows after several manual tests, as a balance between the stability of the results and the shortest possible length. Shorter time windows yield higher resolution but produce more scattered results with too few points for robust cross-spectrum phase interpolation, while longer windows reduce resolution but enable more reliable interpolation and results²³,⁷² (Supplementary Fig. 20).

In Supplementary Note 5 and Supplementary Fig. 21, we provide the results using additional array geometries: NSG for the nodes, DST and DSG for the DAS. Nodes and DAS channels that compose all the arrays are summarised in Supplementary Table S1.

## Infrasonic array processing

We handle infrasonic data using multi-channel correlation analysis (MCCA)⁴⁹,⁷⁴,⁷⁵. For a three-sensor array and a given time window, delay times $\tau_{ij}$ between sensor pairs (i, j) are determined from the peak of the cross-correlation function. The back-azimuth $\theta$ and apparent velocity $c_a$ of an incoming plane acoustic wave crossing an array triplet (i, j, k) are then computed by solving the following system:

$$ \left\{ \begin{array}{l} \frac{L_{ij} \cos(\beta_{ij} - \theta)}{\tau_{ij}} = c_{ij}, \\ \\ \frac{L_{ik} \cos(\beta_{ik} - \theta)}{\tau_{ik}} = c_{ik}, \end{array} \right. \qquad (3) $$

where L and $\beta$ are the horizontal distance and azimuth between sensors. For coherent signals, the time residual $\tau_n$ in the sensor triplet,

$$ \tau_n = |\tau_{ij} + \tau_{jk} + \tau_{ki}|, \qquad (4) $$

approaches zero ($\tau_n \to 0$), while uncorrelated signals yield positive residuals. We classify infrasonic detections when mean residuals fall below an empirical threshold (0.05 s in our case). Each detection provides estimates of $\theta$ and $c_a$. MCCA is applied to 1-10 Hz infrasonic signals using 5 s windows sliding every 1 s. We use the EAR array (Figs. 1 and 3i), which has an aperture of ~290 m.

## Estimating seismic energies

To retrieve the volume of PDC1, we rely on the proportional relationship between flow volumes and seismic energies generated by the two flows and recorded with UniFi broadband seismometers (STR, ROC, and SCI) and DAS. We compute the seismic energies using the equation⁷⁶,

$$ E_s = \int_{t_1}^{t_2} 2\pi r \rho c u_{env}(t)^2 e^{\alpha r} dt, \qquad (5) $$

being $t_1$ and $t_2$ the starting and ending times of the seismic signal generated by the event, r the source-receiver distance, $\rho$ the ground density, h the thickness of the propagation layer for surface waves with speed c, $u_{env}(t) = \sqrt{u(t)^2 + Ht(u(t))^2}$ the envelope of the ground velocity derived from the Hilbert transform $Ht$, and $\alpha = f\pi/Qc$ a frequency-dependent damping factor with Q the quality factor¹⁵,⁷⁶.

To compute the seismic energies generated by the flows, we take the signal onset from visible cameras, corresponding to $t_s_1$ and $t_s_2$ provided in Fig. 2 for PDC1 and PDC2, respectively. The window onset aligns with the increase in seismic amplitude (Supplementary Fig. 6). Since PDC2 seismo-acoustic signals are followed by a high-intensity phase probably unrelated to the flow dynamics (Fig. 3), we consider the final time $t_2$ in Eq. (5) as 4 s after $t_1$ and $t_2$, i.e. when the flows reached the coastline (Fig. 2, Supplementary Fig. 6). We include a 4 s delay to account for the travel time of waves propagating from the shoreline to the most distant parts of the seismic and fibre-optic network (~2.3 km) at velocities as low as 600 m/s²⁷. This assumption is expected to underestimate seismic energies while roughly preserving consistent energy ratios between the two flows. However, the method¹⁵ is designed to describe the total loss of potential energy from the initial mass to the final deposit. Using a time window shorter than the full propagation of the flow up to its resting position (which, in our case, is prevented by the PDCs entering the sea) may introduce errors in the seismic energy estimate and, consequently, in the derived volumes, potentially leading to overestimation of smaller volumes.

Being $f = 5$ Hz the approximate centre of the frequency range with the highest energies (Supplementary Fig. 6), we consider $h = 50$ m as one third of the wavelength of Rayleigh waves with a speed of 700 m/s²⁷. The distance r is calculated between each station and a point halfway down the SdF, $\rho = 2000$ kg/m³, and $Q = 40$ m³/s²⁷. We consider here the vertical Z component of the broadband seismic seismometers to calculate the energy as in Levy et al.¹⁶. Moreover, since seismic waves from both PDCs travel along similar paths to the stations and the other terms in Eq. (5) remain constant, we may
estimate the volume of PDC1 relative to that of PDC2 directly from the ratio of the signal envelopes integrated over the PDC durations.

Finally, we use Eq. (5) to estimate seismic energy from DAS strain rate. DAS strain-rate measurements are spatially integrated along the fibre to obtain the horizontal seismic velocity<sup>78</sup>. The envelope ratio between the two PDCs from DAS measurements remains unaffected by choosing either strain rate or seismic velocity (Supplementary Fig. 14).

## Numerical simulations

We simulate granular flows using the SHALTOP numerical model<sup>44–47</sup>. The model describes dry granular flows over a 2D topography, b = z(x, y), and relies on the thin-layer approximation (i.e., small thickness compared to longitudinal extension) for depth-averaged equations of mass and momentum conservation<sup>16,32</sup>. SHALTOP models the evolution of (i) the flow thickness perpendicular to the topography, h(x, y, t), and (ii) the depth-averaged velocity along the topography, u(x, y, t), where (x, y) are the Cartesian coordinates and t is the time. SHALTOP has successfully reproduced both natural landslides and laboratory experiments using Coulomb-type friction laws with either a constant friction coefficient μ or a variable coefficient μ(I) dependent on the inertial number I<sup>15,16,32,46,47,79–81</sup>. In this study, we adopt the Coulomb friction law with a single empirical parameter, i.e., the friction coefficient μ = tan(δ), with δ the friction angle. This friction coefficient serves as a proxy for the effective energy dissipation during the flow<sup>32,79,82</sup>. Although more complex laws could describe these granular flows, such as the μ(I) and Φ(I) rheology (with Φ the solid volume fraction) or formulations accounting for dilatancy and pore fluid pressure effects<sup>63,83,84</sup>, we adopt a simple Coulomb friction law here for its practical suitability in hazard assessment.

The topography we use is a 5 m resolution Digital Earth Model of the SdF (Fig. 6). The flow volumes used in the simulations are those estimated from camera images and the seismic signals of the PDCs recorded by broadband seismometers and DAS, namely V₁ = 23.5 ± 9.5 × 10³ m³ and V₂ = 80 ± 9 × 10³ m³ for PDC1 and PDC2, respectively. We locate the mobilised mass of both PDCs below the NE craters at an altitude of ~730 m (Fig. 6a, d).

We first retrieve an initial friction coefficient μ using an empirical relation that links it to the flow volume V. This relation was obtained by fitting run-out distances of small to large landslides, mostly dry, across various contexts<sup>32</sup>:

$$ \mu = \tan(\delta) = V^{-0.0774}. \quad (6) $$

The estimated PDC volumes V₁ and V₂ correspond to friction angles of δ₁ = [24.1°, 25.5°] and δ₂ = [22.5°, 22.8°], respectively. For each PDC, we then test several values of δ around those provided by Eq. (6) to find the ones that best describe the observed PDC velocities and durations. We calculate the simulated flow front velocity by averaging flow velocity within a 150-metre radius of the front, considering only areas where the flow thickness is >0.1 m. The values of δ that fit our observations are ~27° and 21° for PDC1 and PDC2, respectively. These values are slightly higher and lower than those obtained directly from Eq. (6) but are within the uncertainty range of δ ± 5.7° of the estimated friction coefficients<sup>32</sup>.

## Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

## Data availability

All data used in this article can be found at the following public repository: https://doi.org/10.5281/zenodo.17610889.

## Code availability

We handled seismic and infrasonic data with the Obspy Python library<sup>85</sup>, while DAS data with the Xdas Python library<sup>86</sup> (https://doi.org/10.5281/zenodo.15913005). Numerical simulations were performed with the

SHALTOP numerical model<sup>44–47,87</sup> (https://doi.org/10.5281/zenodo.10964107). We processed the GPS data with the RTKLIB software (https://rtklib.com), and we used the Interparc Matlab library (https://www.mathworks.com/matlabcentral/fileexchange/34874-interparc, John D’Er-rico, 2012, MATLAB Central File Exchange) to interpolate the fibre layout for geolocating DAS channels. Geolocating the visible imagery was done with the IMGRAFT Matlab software<sup>88</sup>. Other codes used to generate results and figures are available upon request to the authors.

Received: 16 July 2025; Accepted: 2 December 2025; Published online: 29 January 2026

## References

1. Acocella, V. Volcano Flank Instability and Collapse. In Acocella, V. (ed.) Volcano-Tectonic Processes, 205–244 (Springer International Publishing, Cham, 2021).

2. Ripepe, M. & Lacanna, G. Volcano generated tsunami recorded in the near source. Nat. Commun. 15, 1802 (2024).

3. Risica, G., Rosi, M., Pistolesi, M., Speranza, F. & Branney, M. J. Deposit-derived block-and-ash flows: the hazard posed by perched temporary tephra accumulations on volcanoes; 2018 Fuego Disaster, Guatemala. J. Geophys. Res. Solid Earth 127, e2021JB023699 (2022).

4. Di Traglia, F. et al. Generation of deposit-derived pyroclastic density currents by repeated crater rim failures at Stromboli Volcano (Italy). Bull. Volcanol. 86, 69 (2024).

5. Alvarado, G. E. & Soto, G. J. Pyroclastic flow generated by crater-wall collapse and outpouring of the lava pool of Arenal Volcano, Costa Rica. Bull. Volcanol. 63, 557–568 (2002).

6. Branney, M. J. & Kokelaar, P. Pyroclastic Density Currents and the Sedimentation of Ignimbrites. https://pubs.geoscienceworld.org/gsl/books/book/1536/Pyroclastic-Density-Currents-and-the-Sedimentation. (Geological Society of London, 2002).

7. Lube, G., Breard, E. C. P., Esposti-Ongaro, T., Dufek, J. & Brand, B. Multiphase flow behaviour and hazard prediction of pyroclastic density currents. Nat. Rev. Earth Environ. 1, 348–365 (2020).

8. Brosch, E. et al. Destructiveness of pyroclastic surges controlled by turbulent fluctuations. Nat. Commun. 12, 7306 (2021).

9. Tanguy, J.-C., Ribière, C., Scarth, A. & Tjetjep, W. S. Victims from volcanic eruptions: a revised database. Bull. Volcanol. 60, 137–144 (1998).

10. Brown, S. K., Jenkins, S. F., Sparks, R. S. J., Odbert, H. & Auker, M. R. Volcanic fatalities database: analysis of volcanic threat with distance and victim classification. J. Appl. Volcanol. 6, 15 (2017).

11. Delannay, R., Valance, A., Mangeney, A., Roche, O. & Richard, P. Granular and particle-laden flows: from laboratory experiments to field observations. J. Phys. D: Appl. Phys. 50, 053001 (2017).

12. Roche, O., Attali, M., Mangeney, A. & Lucas, A. On the run-out distance of geophysical gravitational flows: insight from fluidized granular collapse experiments. Earth Planet. Sci. Lett. 311, 375–385 (2011).

13. Watson, L. M., Cannata, A. & Andronico, D. Tracking a pyroclastic density current with seismic signals at Mt. Etna (Italy). J. Geophys. Res. Solid Earth 128, e2022JB026114 (2023).

14. Allstadt, K. E. et al. Seismic and acoustic signatures of surficial mass movements at volcanoes. J. Volcanol. Geotherm. Res. 364, 76–106 (2018).

15. Hibert, C., Mangeney, A., Grandjean, G. & Shapiro, N. M. Slope instabilities in Dolomieu crater, Réunion Island: from seismic signals to rockfall characteristics. J. Geophys. Res. Earth Surf. 116 (2011).

16. Levy, C. et al. Friction weakening in granular flows deduced from seismic records at the Soufrière Hills Volcano, Montserrat. J. Geophys. Res. Solid Earth 120, 7536–7557 (2015).

17. Yamasato, H. Quantitative analysis of pyroclastic flows using infrasonic and seismic data at Unzen Volcano, Japan. J. Phys. Earth 45, 397–416 (1997).
18. Jolly, A. D., Thompson, G. & Norton, G. E. Locating pyroclastic flows on Soufrière Hills Volcano, Montserrat, West Indies, using amplitude signals from high dynamic range instruments. *J. Volcanol. Geotherm. Res.* **118**, 299–317 (2002).

19. Chouet, B. et al. Source and path effects in the wave fields of tremor and explosions at Stromboli Volcano, Italy. *J. Geophys. Res.: Solid Earth* **102**, 15129–15150 (1997).

20. Inza, L. A. et al. Analysis of dynamics of vulcanian activity of Ubinas volcano, using multicomponent seismic antennas. *J. Volcanol. Geotherm. Res.* **270**, 35–52 (2014).

21. Métaxian, J.-P. et al. Migration of seismic activity associated with phreatic eruption at Merapi volcano, Indonesia. *J. Volcanol. Geotherm. Res.* **396**, 106795 (2020).

22. Rost, S. & Thomas, C. Array seismology: methods and applications. *Rev. Geophys.* **40**, 2–1–2–27 (2002).

23. Métaxian, J.-P., Lesage, P. & Valette, B. Locating sources of volcanic tremor and emergent events by seismic triangulation: application to Arenal volcano, Costa Rica. *J. Geophys. Res.: Solid Earth* **107**, ECV 10–1–ECV 10–18 (2002).

24. Jousset, P. et al. Dynamic strain determination using fibre-optic cables allows imaging of seismological and structural features. *Nat. Commun.* **9**, 2509 (2018).

25. Klaasen, S., Paitz, P., Lindner, N., Dettmer, J. & Fichtner, A. Distributed acoustic sensing in volcano-glacial environments— Mount Meager, British Columbia. *J. Geophys. Res. Solid Earth* **126**, e2021JB022358 (2021).

26. Jousset, P. et al. Fibre optic distributed acoustic sensing of volcanic events. *Nat. Commun.* **13**, 1753 (2022).

27. Biagioli, F. et al. Array analysis of seismo-volcanic activity with distributed acoustic sensing. *Geophys. J. Int.* **236**, 607–620 (2024).

28. Nishimura, T. et al. Source location of volcanic earthquakes and subsurface characterization using fiber-optic cable and distributed acoustic sensing system. *Sci. Rep.* **11**, 6319 (2021).

29. Ripepe, M., De Angelis, S., Lacanna, G. & Voight, B. Observation of infrasonic and gravity waves at Soufrière Hills Volcano, Montserrat. *Geophysical Research Letters* **37** (2010).

30. Delle Donne, D. et al. Chapter 9 Thermal, acoustic and seismic signals from pyroclastic density currents and Vulcanian explosions at Soufrière Hills Volcano, Montserrat. *Geol. Soc. Lond. Mem.* **39**, 169–178 (2014).

31. Calvari, S. et al. Overflows and pyroclastic density currents in March-April 2020 at stromboli volcano detected by remote sensing and seismic monitoring data. *Remote Sens.* **12**, 3010 (2020).

32. Lucas, A., Mangeney, A. & Ampuero, J. P. Frictional velocity-weakening in landslides on Earth and on other planetary bodies. *Nat. Commun.* **5**, 3417 (2014).

33. Marbœuf, A. et al. Rheology and bathymetry effects captured by a multilayer landslide-tsunami model. *Geophys. J. Int.* **242**, ggaf232 (2025).

34. Bevilacqua, A. et al. Major explosions and paroxysms at Stromboli (Italy): a new historical catalog and temporal models of occurrence with uncertainty quantification. *Sci. Rep.* **10**, 17357 (2020).

35. Di Roberto, A., Bertagnini, A., Pompilio, M. & Bisson, M. Pyroclastic density currents at Stromboli volcano (Aeolian Islands, Italy): a case study of the 1930 eruption. *Bull. Volcanol.* **76**, 827 (2014).

36. Calvari, S. et al. Multi-parametric study of an eruptive phase comprising unrest, major explosions, crater failure, pyroclastic density currents and lava flows: Stromboli volcano, 1 December 2020–30 June 2021. *Front. Earth Sci.* **10**, https://doi.org/10.3389/feart.2022.899635 (2022).

37. Giudicepietro, F. et al. Seismic and thermal precursors of crater collapses and overflows at Stromboli volcano. *Sci. Rep.* **13**, 11115 (2023).

38. Di Traglia, F., Borselli, L., Nolesini, T. & Casagli, N. Crater-rim collapse at Stromboli volcano: understanding the mechanisms leading from the failure of hot rocks to the development of glowing avalanches. *Nat. Hazards* **115**, 2051–2068 (2023).

39. La Rocca, M. et al. Seismic signals associated with landslides and with a Tsunami at Stromboli Volcano, Italy. *Bull. Seismol. Soc. Am.* **94**, 1850–1867 (2004).

40. Arran, M. I. et al. Laboratory landquakes: insights from experiments into the high-frequency seismic signal generated by geophysical granular flows. *J. Geophys. Res. Earth Surf.* **126**, e2021JF006172 (2021).

41. Belli, G., Walter, F., McArdell, B., Gheri, D. & Marchetti, E. Infrasonic and seismic analysis of debris-flow events at Illgraben (Switzerland): relating signal features to flow parameters and to the seismo-acoustic source mechanism. *J. Geophys. Res. Earth Surf.* **127**, e2021JF006576 (2022).

42. Schindelé, F. et al. A review of tsunamis generated by volcanoes (TGV) source mechanism, modelling, monitoring and warning systems. *Pure Appl. Geophys.* **181**, 1745–1792 (2024).

43. Farin, M., Mangeney, A., de Rosny, J., Toussaint, R. & Trinh, P.-T. Relations between the characteristics of granular column collapses and resultant high-frequency seismic signals. *J. Geophys. Res. Earth Surf.* **124**, 2987–3021 (2019).

44. Bouchut, F., Mangeney-Castelnau, A., Perthame, B. & Vilotte, J.-P. A new model of Saint Venant and Savage–Hutter type for gravity driven shallow water flows. *C. R. Math. 336*, 531–536 (2003).

45. Bouchut, F. & Westdickenberg, M. Gravity driven shallow water models for arbitrary topography. *Commun. Math. Sci.* **2**, 359–389 (2004).

46. Mangeney-Castelnau, A. et al. On the use of Saint-Venant equations to simulate the spreading of a granular mass. *J. Geophys. Res. Solid Earth* **110** (2005).

47. Mangeney, A., Bouchut, F., Thomas, N., Vilotte, J. P. & Bristeau, M. O. Numerical modeling of self-channeling granular flows and of their levee-channel deposits. *J. Geophys. Res. Earth Surf.* **112**, F02017 (2007).

48. Barfucci, G. Modeling infrasonic sources related to density currents. https://flore.unifi.it/handle/2158/1170483 (2019).

49. Marchetti, E. et al. Infrasound array analysis of debris flow activity and implication for early warning. *J. Geophys. Res. Earth Surf.* **124**, 567–587 (2019).

50. Arrowsmith, S. J. et al. Big data seismology. *Rev. Geophys.* **60**, e2021RG000769 (2022).

51. Quinteros, J., Carter, J. A., Schaeffer, J., Trabant, C. & Pedersen, H. A. Exploring approaches for large data in seismology: user and data repository perspectives. *Seismol. Res. Lett.* **92**, 1531–1540 (2021).

52. van den Ende, M. P. A. & Ampuero, J.-P. Evaluating seismic beamforming capabilities of distributed acoustic sensing arrays. *Solid Earth* **12**, 915–934 (2021).

53. Capdeville, Y. & Sladen, A. DAS sensitivity to heterogeneity scales much smaller than the minimum wavelength. *Seismica Library*.mcgill.ca/article/view/1007 (2024).

54. Näsholm, S. P. et al. Array signal processing on distributed acoustic sensing data: directivity effects in slowness space. *J. Geophys. Res. Solid Earth* **127**, e2021JB023587 (2022).

55. Isken, M. P., Vasyura-Bathke, H., Dahm, T. & Heimann, S. De-noising distributed acoustic sensing data using an adaptive frequency–wavenumber filter. *Geophys. J. Int.* **231**, 944–949 (2022).

56. Chen, Y. et al. Denoising of distributed acoustic sensing seismic data using an integrated framework. *Seismol. Res. Lett.* **94**, 457–472 (2022).

57. Bayarri, M. J. et al. Using statistical and computer models to quantify volcanic hazards. *Technometrics* **51**, 402–413 (2009).

58. Yamada, M., Mangeney, A., Matsushi, Y. & Matsuzawa, T. Estimation of dynamic friction and movement history of large landslides. *Landslides* **15**, 1963–1974 (2018).

59. Borykov, T. et al. Empirical investigation of friction weakening of terrestrial and Martian landslides using discrete element models. *Landslides* **16**, 1121–1140 (2019).

60. Martin, H. A. et al. Ultrasound-induced dense granular flows: a two-time scale modelling. *J. Fluid Mech.* **1004**, A10 (2025).
61. Léopoldès, J., Jia, X., Tourin, A. & Mangeney, A. Triggering granular avalanches with ultrasound. *Phys. Rev. E* **102**, 042901 (2020).

62. Roche, O., Henry, C. D., Azzaoui, N. & Guillin, A. Long-runout pyroclastic density currents: analysis and implications. *Geology* **50**, 1172–1176 (2022).

63. Breard, E. C. P. et al. The fragmentation-induced fluidisation of pyroclastic density currents. *Nat. Commun.* **14**, 2079 (2023).

64. Darvenne, A., Viroulet, S. & Lacaze, L. Physical model of landslide-generated impulse waves: experimental investigation of the wave-granular flow coupling. *J. Geophys. Res. Oceans* **129**, e2024JC021145 (2024).

65. Løvholt, F., Pedersen, G., Harbitz, C. B., Glimsdal, S. & Kim, J. On the characteristics of landslide tsunamis. *Philos. Trans. R. Soc. A: Math. Phys. Eng. Sci.* **373**, 20140376 (2015).

66. Heinrich, P. et al. Simulation of water waves generated by a potential debris avalanche in Montserrat, Lesser Antilles. *Geophys. Res. Lett.* **25**, 3697–3700 (1998).

67. Ripepe, M., Delle Donne, D., Legrand, D., Valade, S. & Lacanna, G. Magma pressure discharge induces very long period seismicity. *Sci. Rep.* **11**, 20065 (2021).

68. Messerli, A. & Grinsted, A. Image georectification and feature tracking toolbox: ImGRAFT. *Geosci. Instrum. Methods Data Syst.* **4**, 23–34 (2015).

69. Capon, J. High-resolution frequency-wavenumber spectrum analysis. *Proc. IEEE* **57**, 1408–1418 (1969).

70. Schmidt, R. Multiple emitter location and signal parameter estimation. *IEEE Trans. Antennas Propag.* **34**, 276–280 (1986).

71. Goldstein, P. & Archuleta, R. J. Array analysis of seismic signals. *Geophys. Res. Lett.* **14**, 13–16 (1987).

72. Poupinet, G., Ellsworth, W. L. & Frechet, J. Monitoring velocity variations in the crust using earthquake doublets: an application to the Calaveras Fault, California. *J. Geophys. Res. Solid Earth* **89**, 5719–5731 (1984).

73. Trabattoni, A. et al. Microseismicity monitoring and site characterization with distributed acoustic sensing (DAS): the case of the irpinia fault system (Southern Italy). *J. Geophys. Res. Solid Earth* **127**, e2022JB024529 (2022).

74. Cansi, Y. An automatic seismic event processing for detection and location: the P.M.C.C. method. *Geophys. Res. Lett.* **22**, 1021–1024 (1995).

75. Ulivieri, G. et al. Monitoring snow avalanches in Northwestern Italian Alps using an infrasound array. *Cold Reg. Sci. Technol.* **69**, 177–183 (2011).

76. Vilajosana, I. et al. Rockfall induced seismic signals: case study in Montserrat, Catalonia. *Nat. Hazards Earth Syst. Sci.* **8**, 805–812 (2008).

77. Petrosino, S., Cusano, P., Saccorotti, G. & Del Pezzo, E. Seismic attenuation and shallow velocity structures at Stromboli Volcano, Italy. *Bull. Seismol. Soc. Am.* **92**, 1102–1116 (2002).

78. Trabattoni, A. et al. From strain to displacement: using deformation to enhance distributed acoustic sensing applications. *Geophys. J. Int.* **235**, 2372–2384 (2023).

79. Favreau, P., Mangeney, A., Lucas, A., Crosta, G. & Bouchut, F. Numerical modeling of landquakes. *Geophys. Res. Lett.* **37**, 2010GL043512 (2010).

80. Moretti, L. et al. Numerical modeling of the Mount Steller landslide flow history and of the generated long period seismic waves. *Geophys. Res. Lett.* **39** (2012).

81. Moretti, L. et al. Numerical modeling of the Mount Meager landslide constrained by its force history derived from seismic data. *J. Geophys. Res. Solid Earth* **120**, 2579–2599 (2015).

82. Mangeney, A. et al. Erosion and mobility in granular collapse over sloping beds. *J. Geophys. Res. Earth Surface* **115** (2010).

83. Bouchut, F., Fernández-Nieto, E. D., Mangeney, A. & Narbona-Reina, G. A two-phase two-layer model for fluidized granular flows with dilatancy effects. *J. Fluid Mech.* **801**, 166–221 (2016).

84. Bouchut, F., Drach, E., Fernández-Nieto, E. D., Mangeney, A. & Narbona-Reina, G. A series of two-phase models for grain-fluid flows with dilatancy. *J. Fluid Mech.* **1008**, A43 (2025).

85. Beyreuther, M. et al. ObsPy: a Python toolbox for seismology. *Seismol. Res. Lett.* **81**, 530–533 (2010).

86. Trabattoni, A. et al. Xdas: a Python framework for distributed acoustic sensing. *Seismol. Res. Lett.* **96**, 3221–3230 (2025).

87. Bouchut, F., Mangeney, A., Perruzzetto, M., Lucas, A. & Sainte-Marie, J. Shaltop2d - A SHALlow water depth-averaged numerical model for geophysical granular flows over complex TOPography https://zenodo.org/records/10964107 (2024).

## Acknowledgements

This study depends on the fundamental fieldwork done by S. Oliva and his team, who dug the trench for the fibre-optic cable on the volcano's steep slopes. We also acknowledge the help offered in the field by G. Diana, M.-P. Bouin, R. Longo, Y. Capdeville, S. Kiel, M. Masini, and J.P. Sepulveda Birke, and the fruitful discussions with D. Delle Donne and J.-B. De Chabalier. We thank C. Jestin, G. Calbris, and V. Lanticq of FEBUS Optics for their support and two anonymous reviewers for their comments, which improved the quality of the manuscript. We also thank S. Cannavò for providing the drone images used in this work. This work is financially supported by the ANR-19-CE04-0011-01 MONIDAS (Natural Hazard Monitoring using DAS) Project. This study received funding through the service activity enhancement project "Sviluppo del sistema unico (INGV-Università) di monitoraggio vulcanico e rilevamento precoce dei maremoti e delle esplosioni parossistiche di Stromboli" funded by the Dipartimento della Protezione Civile and the Istituto Nazionale di Geofisica e Vulcanologia (INGV). The study/service/product does not necessarily reflect the policy and position of the INGV and Dipartimento della Protezione Civile, Italy.

## Author contributions

J.-P.M., E.S., P.B., M.R., F.B., and A.T. conceived the experiment and the main ideas and objectives. J.-P.M., E.S., P.B., M.R., A.T., L.I., and F.B. conducted the experiment in 2022. F.B., J.-P.M., E.S., P.B., and M.R. analysed DAS and node data together with the results. F.B., G.U., M.R., and G.L. analysed infrasonic and seismic broadband data and results. F.B., G.L., M.R., and G.R. analysed visible and thermal images. F.B., A.M., and M.P. performed simulations of granular flows and interpreted the results. F.B. wrote the first draft of the manuscript and produced graphs and plots. All authors discussed the results, reviewed the manuscript, and provided critical feedback.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s43247-025-03091-6.

**Correspondence** and requests for materials should be addressed to Francesco Biagioli.

**Peer review information** *Communications Earth and Environment* thanks the anonymous reviewers for their contribution to the peer review of this work. Primary Handling Editors: Domenico Doronzo and Alireza Bahadori. A peer review file is available.

**Reprints and permissions information** is available at http://www.nature.com/reprints

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2026