Next Article in Journal
Diagnostic Multidisciplinary Investigations for Cultural Heritage at Etna Volcano: A Case Study from the 1669 Eruption in the Mother Church at the Old Settlement of Misterbianco
Next Article in Special Issue
Volcanic Processes Monitoring and Hazard Assessment Using Integration of Remote Sensing and Ground-Based Techniques
Previous Article in Journal
Potential of Ultra-High-Resolution UAV Images with Centimeter GNSS Positioning for Plant Scale Crop Monitoring
Previous Article in Special Issue
Changes in the Eruptive Style of Stromboli Volcano before the 2019 Paroxysmal Phase Discovered through SOM Clustering of Seismo-Acoustic Features Compared with Camera Images and GBInSAR Data
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Comparison between Automated and Manual Detection of Lava Fountains from Fixed Monitoring Thermal Cameras at Etna Volcano, Italy

1
Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo—Sezione di Catania, Piazza Roma 2, 95125 Catania, Italy
2
Dipartimento di Ingegneria Elettrica, Elettronica e Informatica, Università degli Studi di Catania, Viale A. Doria 6, 95122 Catania, Italy
*
Author to whom correspondence should be addressed.
Remote Sens. 2022, 14(10), 2392; https://doi.org/10.3390/rs14102392
Submission received: 23 March 2022 / Revised: 12 May 2022 / Accepted: 13 May 2022 / Published: 16 May 2022

Abstract

:
The Etna volcano is renowned worldwide for its extraordinary lava fountains that rise several kilometers above the vent and feed eruptive columns, then drift hundreds of kilometers away from the source. The Italian Istituto Nazionale di Geofisica e Vulcanologia-Osservatorio Etneo (INGV-OE) is responsible for the monitoring of Mt. Etna, and for this reason, has deployed a network of visible and thermal cameras around the volcano. From these cameras, INGV-OE keeps a keen eye, and is able to observe the eruptive activity, promptly advising the civil protection and aviation authorities of any changes, as well as quantifying the spread of lava flows and the extent of pyroclastic and ash plumes by using a careful analysis of the videos recorded by the monitoring cameras. However, most of the work involves analysis carried out by hand, which is necessarily approximate and time-consuming, thus limiting the usefulness of these results for a prompt hazard assessment. In addition, the start of lava fountains is often a gradual process, increasing in strength from Strombolian activity, to intermediate explosive activity, and eventually leading to sustained lava fountains. The thresholds between these different fields (Strombolian, Intermediate, and lava fountains) are not clear cut, and are often very difficult to distinguish by a manual analysis of the images. In this paper, we presented an automated routine that, when applied to thermal images and with good weather conditions, allowed us to detect (1) the starting and ending time of each lava fountain, (2) the area occupied by hot pyroclasts, (3) the elevation reached by the lava fountains over time, and (4) eventually, to calculate in real-time the erupted volume of pyroclasts, giving results close to the manual analysis but more focused on the sustained portion of the lava fountain, which is also the most dangerous. This routine can also be applied to other active volcanoes, allowing a prompt and uniform definition of the timing of the lava fountain eruptive activity, as well as the magnitude and intensity of the event.

Graphical Abstract

1. Introduction

New data and interpretations have emerged of the geodynamics of the eastern Sicily point to Etna as a volcano, undergoing an evolutionary phase where a future increase in highly energetic explosive activity is possible [1]. As a matter of fact, the last three decades of Etna’s activity were characterized by frequent highly explosive eruptions, here called paroxysms [2,3]. Paroxysms at Etna are characterized by lava fountaining lasting 1–2 h, reaching the height of 1–3 km above the crater, and generating conspicuous and lengthy ash plumes that can drift hundreds of kilometers from the vents [4,5], often accompanied by short-lasting lava overflows from the crater rim [6,7]. During the last few decades, Etna volcano underwent several eruptions characterized by lava flows mainly from the summit vents, alternating with short-lasting but powerful explosive episodes [2,3,8,9,10]. In particular, between 2011 and 2015, Etna produced more than 50 such eruptions [3,10,11,12], releasing a cumulative erupted volume of a similar order to a major flank eruption [2], which was normally ten times greater than summit activity [13]. Given that explosive paroxysms can have a major impact on aviation [14], on road and traffic conditions, and also on the villages on the slope of the volcano [15,16,17], it is of paramount importance for a volcano observatory such as the Istituto Nazionale di Geofisica e Vulcanologia-Osservatorio Etneo (INGV-OE) to be able to raise an early warning as soon as possible, and then to advise the civil protection authorities of its possible impact on human activities [18,19,20,21]. The first and most important parameter to be detected as soon as possible is the timing of start and end of any impending activity, and this information needs to be completed with the extent of the ash plume, lava flows and lava fountains, and with the volume erupted [3,6,8,12,22,23]. The volume erupted during a lava fountain (LF) episode quantifies the magnitude of the event, whereas the eruption rate determines its intensity [24]. The information gathered from the monitoring system is then used to inform the civil protection of the magnitude and intensity of the event, and also in models routinely used for the prediction of the extent and distribution of the eruption products, which, at Etna, often comprise both lava flows [25,26] and pyroclastics [27,28,29]. An automated procedure to map the lava flows from the images of the thermal monitoring cameras was recently developed [6,7], whereas the LF detection and the estimation of the erupted volumes were normally carried out by manual analysis of the images [3,8,30]. The pyroclastic volume estimated by thermal images was compared to the total erupted volume estimated by strain and with the lava flow volume erupted during each episode as estimated by satellite [30], and an error of ~20% was estimated for the calculated fluid volume, comprising gas plus pyroclastics [31].
A sustained LF normally gives rise to an eruptive column comprising three main portions: a lower and innermost zone called gas-thrust region with the highest velocity at the exit of the erupting vent; an intermediate zone, where convective movements of the hot mixture of gas and tephra allow ingestion of the cold surrounding air, thereby slowing down the spreading hot jet, and an uppermost zone reaching the buoyancy zone and consisting of the laterally spreading umbrella region [32,33]. In addition, the eruptive columns at Etna can be distinguished into weak plumes if bent in the wind direction due to wind speeds greater than 10 m s−1, or strong to intermediate plumes when rising vertically above the vent or slightly bent in the wind direction for wind speeds lower than 10 m s−1 [3].
LF heights were determined following several different methods. One of the earliest, applied at Kilauea, was carried out by digitizing film from 8 mm time-lapse movie cameras deployed on the ground, and using a few theodolite measurements as calibration points for the film [34]. At Etna, a C++ code was developed in order to explore the INGV-OE thermal image library for image processing [35], applying appropriate thresholds, and converting the color images into a binary black and white image over which the maximum vertical extent of the LF can be easily retrieved. A less automated but sometimes more precise system involves the visual and manual analysis of each image [3,8], given that the thermal images can be affected by low gas, weather or ash clouds [36,37,38,39,40] that may reduce the automated measurement of the LF height, or the LF jet may be inclined [3,8,30,35]. Calibrated images of visible cameras can also be used to estimate the vertical extent of proximal ash plumes associated with the lava fountains [12,22].
One of the most difficult challenges in volcanology is to determine when an eruption is over, especially when it includes multiple episodes and long pauses [41], although sometimes, a gradual decline of the mass eruption rate may anticipate the end of the eruption [42]. During the lava fountain activity at Etna, the start and ending time, as well as any early warning alarm, is given on the basis of the volcanic tremor and infrasound [18,43,44,45,46]. However, a volcanic tremor does not allow us to calculate the volume erupted [3,8], and does not provide information about the extension of the lava fountains and ash plume [5,22]. Another useful device is the borehole strainmeter, which allows calculating the total erupted volume [47,48], comprising both lava flows and pyroclastics [30]. However, in order to assess the impact on the population, the amount of solely the pyroclastics component erupted during a lava fountain event needs to be established, because this affects the stability of roofs, the cleaning up of roads and motorways, the impact on the nearby Catania international airport, and the health effects on the local population [49,50]. Automated routines for volcanic activity detection and characterization have been recently developed [23,45,51,52], and they will probably resolve most of the issues related with early warning alarms. In this paper, we presented a new automated routine that, when applied to the images recorded by the thermal monitoring cameras, allowed us to calculate (1) lava fountain height, (2) area of the lava fountain jet, and (3) volume of the erupted tephra, using the formula applied by Calvari et al. [3,8]. We compared these results with those obtained by a manual analysis, discussing limits and advantages, and future possible improvements.

2. Methods

The INGV-OE thermal camera network for monitoring Etna volcano includes four fixed, continuously operating thermal cameras located on the flanks of the volcano at different distances from and looking towards the summit (Figure 1, Table 1). The images recorded by these cameras are transmitted to the INGV-OE Operative Room and displayed in real-time to allow continuous monitoring of the volcano. From this perspective, operators have the task of recognizing the type of event as early as possible. Thermal cameras are remote-sensing ground-based fixed devices that have significantly improved INGV-OE’s observational capabilities. They allowed us to monitor the summit area continuously, and to identify and locate eruptive events. Only poor weather conditions, and especially thick clouds of water vapor, gas, and/or volcanic ash [37,38,39,40] may affect the visibility of the cameras, and consequently, the reliability of the acquired images, by partially or totally hiding what is happening in the monitored area. When, by manual examination of images, we recognized that visibility was limited to a few frames or interfering clouds were low and only partially obscured the lava fountain (LF), a linear interpolation was carried out on the data [3]. The manual analysis of the thermal images followed what was described by Calvari et al. [3,8].
A sequence of LF episodes occurred at Mt. Etna between 2020 and 2022, and the list of these events is given in Table 2. We have analyzed the images of the thermal cameras that recorded the event (Figure 1, Table 2), and chose for each episode the one offering the best view and the entire vertical extension of the LF, as a function of the LF size, wind direction, and consequent ash plume fallout. In cases of rotating ash plumes, we used the integration of images from more than one camera. The manual analysis of the camera images allowed retrieving the starting and ending time of each episode, the duration expressed in minutes and seconds, the maximum height of the lava fountain and its average value, as well as the erupted volume of pyroclastics and the time-averaged discharge rate (TADR, [53]). The volume of pyroclastics was calculated by following the method by Calvari et al. [3,8], based on the measurements of the LF heights at time lapses of 1 min, and considering a constant vent radius of 15 m and a pyroclastic ratio of 0.18% of the total erupted fluids comprising gas plus pyroclastics. All these results, obtained from the manual analysis of the thermal camera images, are reported in Table 2.
The dataset consisted of 66 episodes of LF recorded between 13 December 2020 and 21 February 2022 at Mt Etna by using the four thermal cameras of the INGV-OE monitoring network listed in Table 1, and whose position is shown in Figure 1. Additional technical details on the thermal sensors can be found in Calvari et al. [3,29]. Original data were provided as .avi format files, each containing about 5 min of recorded volcanic activity. Information concerning the name of the camera and the starting time of each file were embedded in its name. For instance, the filename EMOT 20210319–082500.avi, refers to a file recorded by the camera named EMOT, starting on 19 March 2021 at time 08:25:00. All times are UTC. Other information, such as the duration, the frame-rate, and other useful video properties, such as width, height, bits per pixel, and the video format were embedded in the file object. Files were pre-processed in order to crop the color bar and the information about the acquisition time and camera name, which are normally embedded in the frame.
In order to detect the presence of LF, each frame was converted from the original RGB format to grayscale and further binarized by adopting a threshold luminance value, specified as a value in the range [0, 1]. In this way, the hottest objects, such as newly erupted or cooling down products, will be represented in the binarized image as white areas, while all the others will be represented in black (Figure 2). Thus, in the absence of hot objects, such as new ejected volcanic matter, hot spots or cooling lava, the binarized images will result completely black. However, as mentioned above, it should be noted that hot objects may not be detected due to the presence of a thick cloud cover. Clearly, the choice of the threshold parameter played a crucial role. In fact, although on one hand, it allowed filtering unwanted information, due, for example, to warming of the monitored area by the Sun, on the other hand, it could also remove useful information, especially in the phase of emergence of an LF episode. Unfortunately, there were no optimal criteria for the choice of this parameter which was therefore made for each camera, adopting the traditional trial and error method. Each binarized image was processed, and detected objects, represented in white color in Figure 2c, were measured in order to obtain:
the areas A i ,   i = 1 , N , A i being the area in pixels of the i-th object and N the number of recognized objects;
the coordinates x i , y i of the centroid, of the i_th object, x i and y i being the horizontal and vertical coordinates, respectively.
The overall process, starting from the original image (Figure 2a), through the grayscale (Figure 2b), the binarized (Figure 2c), and finally, to the labeled image (Figure 2d) is reported, as an example, in Figure 2. In particular, the centroids of detected hot objects are shown by the red asterisks in Figure 2d.
The presence of multiple objects, even from an individual LF episode, was due to the fact that the volume occupied by an LF does not have a uniform temperature, as can be seen from Figure 2a where the original colors of the thermal image indicate different temperatures in a scale starting from blue (0 °C) to white (>100 °C). However, for practical reasons, as different hot objects belong to the same individual LF, it is a good practice to consider all detected objects as a single one. In our case, this was obtained by summing up the areas of all detected objects, and calculating the coordinates of a single centroid by a weighted average of the coordinates of the individual centroids, as expressed in Equation (1).
x = i = 1 N w i x i ; y = i = 1 N w i y i ,   w i = A i i = 1 N A i
where Ai is the area of the i-th object, N is the number of detected objects, and w i is the normalized area of the i-th object. It is straightforward to say that, with considerable approximation, due to the fact that the LF occurs in 3-D volume, while images refer to a 2-D area, the estimated area A is, in some way, related to the volume of hot matter, while the y-coordinate of its centroid may be related with the mean elevation. Of course, A and its centroids’ coordinates (x,y), originally expressed in pixel units, can be converted into geographical units by appropriate conversion constants, depending on the position of the considered camera with respect to the monitored area. The graphical representation of the area A (Figure 3a) and its mean altitude (Figure 3b) for an LF episode occurring on Mt. Etna on 24 February 2021, is shown in Figure 3.
Figure 3 can be interpreted as follows: during 24 February 2021, a mass of previously erupted material was cooling down around the crater area. This can be deduced from the fact that the y-coordinate of the centroids (Figure 3b) fluctuates around zero during the time interval from 00:00 to about 19:25. In Figure 3b, for convenience, the zero value of the mean altitude has been arbitrarily set to the average value of the y-coordinate, measured throughout the recording period. The lack of signal in the mean altitude (in the reported example, for instance, around 14:20 and 14:40) is due to the thick cloud cover which prevented viewing of the cooling mass. After 19:15, an LF appeared, as can be seen from both the area signal (Figure 3a) and the increasing mean altitude of the center of mass (Figure 3b). The LF continued until about 20:00, after which, the erupted material began a cooling process. It should be observed that while the value of the area slowly decreased, the value of the mean altitude of the hot objects almost instantaneously decreased when the LF ended. As can be seen from Figure 3, an LF can therefore be recognized by the characteristic bell-shaped distribution of both the area and mean altitude time series. This suggests a criterion for identifying the time mark to be associated with the beginning and the end of the paroxysmal phenomenon, by using a Change Point Detection (CPD) algorithm, which will be explained in Section 3.2.

3. Results

3.1. Manual Estimation of the Eruptive Activity

Of the 66 LF episodes that occurred between 13 December 2020 and 21 February 2022, only one could not be detected because of poor weather conditions: the episode of 27 May 2021, when there was poor visibility from all the monitoring cameras. In all the other cases, visibility was more than 80% of the duration of the episode. Thus, for the short lapses of time when there was no visibility, the LF heights were linearly interpolated. The total volume of erupted pyroclastics during the 65 episodes was ~65 × 106 m3. The average duration of the 65 episodes listed in Table 2 was of 120 min (minimum 6 min, maximum 803 min), or 7171 s (minimum 360 s, maximum 48,180 s); the average erupted volume of pyroclastics was 0.99 × 106 m3, with 41 × 103 m3 as minimum value and 4.6 × 106 m3 as maximum value; the TADR, calculated for the only pyroclastic portion of the episodes and for the whole duration of each paroxysmal event, was 145 m3 s−1 on average, spanning from a minimum of 59 m3 s−1 and a maximum of 245 m3 s−1; maximum LF heights above the vents were 1815 m on average, spanning from a peak value of 5714 m to a minimum of 333 m, and the average heights of the LFs were between 115 m and 778 m, with a peak value of 2160. All values are listed in Table 2.
From a volcanological point of view, it is worth noting that if we exclude the two outliers of 23 and 31 March 2021, which emitted more than 4 × 106 m3 of tephra, the volume of pyroclastics which erupted during the LF episodes occurring between mid-December 2020 and February 2022 increased with time. This can clearly be seen from Figure 4, where we report the distribution of the erupted volumes with time and its linear trend. In addition, the time-averaged discharge rate (TADR) and its trend increased with time (Figure 5). These observations, i.e., of the increase of erupted volume of pyroclastics and rate of eruption with time, suggested that the sequence of paroxysmal events was not yet over [31].

3.2. Change Point Detection

A change point represents a transition between different states in a process that generates the time series. Change point detection (CPD) can be defined as the problem of choosing between two alternatives: no change or indeed, a change occurred. CPD algorithms are traditionally classified as online or offline [54]. Offline algorithms consider the whole data set at once and try to recognize where the change occurred. Thus, the aim in this case, is to identify all the sequence change points in batch mode. In contrast, online, or in real-time, algorithms run concurrently with the process they are monitoring, processing each data point as it becomes available, with the goal of detecting a change point as soon as possible after it occurs, ideally before the next data point arrives. In practice, no CPD algorithm operates in perfect real-time because it must wait for new data before determining if a change point occurred. However, different online algorithms require different amounts of new data before a change point can be detected. Based on this observation, an online algorithm, which needs at least ε samples in the new batch of data to be able to find a change, is usually denoted as ε-real time. Therefore, offline algorithms can be viewed as ∞-real time whereas the best online algorithm is 1-real time, because for every data point, it can predict whether or not a change point occurs before the new data point. Smaller ε values may lead to stronger, prompter CPD algorithms. To find a change point in a time series, a global optimization approach can be used with the following basic algorithm:
  • Choose a point and divide the signal into two sections.
  • Compute an empirical estimate of the desired statistical property for each section.
  • At each point within a section, measure how much the property deviates from the empirical estimate, and at the end, add the deviation for all points.
  • Add the deviations section-to-section to find the total residual error.
  • Vary the location of the division point until the total residual error attains a minimum.
As noted above, the search for a change point k can be formulated as an optimization problem where the cost function J(k) to minimize it can be written, in the general case as:
J k = i = 1 k 1 Δ x i ; χ x 1 , , x k 1 + i = k N Δ x i ; χ x k , , x N
where x 1 ,   x 2 ,   x N is the time series, χ is the chosen statistic, and ∆ is the deviation measurement. In particular, when χ is the mean, the cost function assumes the following form:
J k = i = 1 k 1 ( x i x 1 k 1 ) 2 + i = k N ( x i x k N ) 2
where the symbol ⟨·⟩ indicates the mean operator.
Another aspect to be considered, when formulating the optimization problem, is that signals of practical interest have more than one change point. Furthermore, the number of change points K is often not known a priori. To handle these features, the cost function can be generalized as:
J k = r K 1 i = k r k r + 1 1 Δ x i ; χ x k , , x r + 1 + β K
where k0 and kK are, respectively, the indexes of the first and the last sample of the signal. In the expression (4), the term βK is a penalty term, linearly increasing with the number of change points K, which avoids the problem of overfitting [55]. Here, β represents a positive coefficient that weights the number of searched change points. Indeed, in an extreme case (i.e., β = 0), J(K) reaches the minimum value (i.e., 0) when every point becomes a change point (i.e., K = N).
The algorithm described above for a univariate time series, can easily be extended to the case of a multivariate time series, which was the case, for instance, of a data set recorded by a GPS network [54]. In this case, the cost function was evaluated, of course, over the whole set of available time series. The software considered in this work was implemented in MATLAB based on the CPD algorithms described in [55,56]. The package can help the user at various levels. The lowest level is to consider the software to obtain the time series of area A and mean altitude MA of detected hot objects from images of volcanic activity, following the algorithm described in the previous section. Subsequently, the user, based on a visual inspection of these time series, can indicate presumed times for the starting and ending times of the LF, and request the software to calculate other quantities of interest such as the volumes emitted. Another possibility is to leave the software to search for the transition times from Strombolian to paroxysmal activity using the CPD algorithm, possibly selecting the statistic to be used to perform the detection (i.e., abrupt changes in the mean, in the variance or in the slope). While the latter possibility is preferable when the time series shows clear LF episodes, the former is more suitable when the automated detection of the transition from Strombolian to paroxysmal activity and vice versa can be more problematic due to external noise sources (e.g., poor visibility or interference by other kinds of hot objects).

3.3. Timing the Lava Fountains Occurring at Etna during 2020–2022

The main advantages of performing a computer process analysis of LF images are the following:
the user can quickly analyze the content of the image files recorded over days, an operation which, carried out manually, requires a considerable amount of time;
the user can speed up the computation of key quantities such as height and duration of the LF, which are necessary for the calculation of the volumes of erupted material);
it is possible to implement algorithms for automatically timing the transition from Strombolian to paroxysmal activity, which is otherwise left to human judgment, gaining in uniformity and repeatability;
in case of lack of visibility, since it is necessary to proceed with interpolation of the data, a rather difficult operation to perform manually, the user can resort to automated interpolation techniques (e.g., linear interpolation, nearest, etc.).
However, it should not be overlooked that quantitative measurements of the LF parameters, starting from the images, require overcoming several non-negligible difficulties. First of all, the aforementioned lack of visibility can make the reliability of the measurements poor. In fact, interpolation techniques can help in cases of limited amounts of data, but obviously cannot replace them when significant amounts are missing.
Furthermore, since the LFs are recognized as hot objects, they can be confused with hot objects of other kinds. Hot areas are very often formed due to the sunlight reflection of both the ground (Figure 6a) and the clouds (Figure 6c). At other times, the Sun itself was included in the images as it travels its natural orbit (Figure 6d,e). Other kinds of hot objects, which could be confused with LFs, were cooling lava flows (Figure 6b). Moreover, different hot objects can combine their effects with those of the LFs (Figure 6e,f). Some of these effects could be avoided by using special cameras, but this is not always possible. Since the current state of development of the software does not allow a reliable distinction between the noise and the LF signal, the user, for the purpose of determining the start and end times, can limit the search space, so as to exclude particularly noisy periods.
Normally, the area signal is smoother than the mean altitude signal and therefore, usually when determining the start and end time of an LF using the CPD algorithm, it is preferable to consider the time series of the areas. However, in some cases, it may be useful to consider the series of mean altitude, as described in the following example. For instance, consider the area signal shown in Figure 7a, which refers to an LF occurring on Etna on 20 February 2021. As can be seen, due to the presence of a cooling lava flow, the area signal slowly decreases (thermal hysteresis), making it difficult to accurately perform the automated detection of the ending time of this LF episode. However, this shortcoming can be overcome by performing the CPD timing on the mean altitude signal, which is not affected by the thermal hysteresis (see Figure 7b).
Here, the term ‘timing’ will be used to indicate the estimation of the starting and ending time of an LF episode. In particular, for this LF, while performing the CPD timing from the area signal (Figure 8a), the end of the LF is estimated to be tend = 06:38, whereas considering the mean altitude signal, the end of the LF is estimated tend = 00:56, which is closer to that manually estimated (tend = 01:15) and reported in Table 2. This difference was caused by the greater curvature that can be seen in Figure 8a, which was determined by the cooling lava flow area.

3.4. Timing the Lava Fountains by a Gaussian Function-Based Approach

In some cases, assuming that a typical LF has a time distribution of area and mean altitude, which roughly has a bell shape, it might be useful to approximate the measured data by using a Gaussian function. This can be useful, for example, when the images are affected by thick clouds passing through the field of view of the camera, generating a trend such as those shown in Figure 9a, which refers to the LF episode on 13 December 2020 and was observed from the ENT camera.
For this LF episode, the manually estimated starting and ending times were 22:00 and 22:48, respectively. However, from Figure 9a, it can be seen that precisely between these two times, the recorded signal is discontinuous due to poor visibility, but nevertheless, it can be seen that the area signal shows a well-detectable peak. Fitting the area samples, it is possible to obtain the results shown in Figure 9b and thus estimate the start and end dates of the LF to be 22:02 and 22:20, respectively, which are closer to the manually assessed corresponding times.
Timing of an LF episode, after having carried out the approximation of the curve by means of a Gaussian function, is simply established by using a threshold approach: the starting time is set as the one in which the recorded data exceed, for the first time, a threshold of the normalized function height. Similarly, the end time is established as the one in which the recorded signal falls, for the first time, below the threshold. In this paper, the threshold value has been set, after a trial and error approach, to be 25% of the maximum value.
Of course, the Gaussian approach can be used for timing LF as an alternative to the CPD one, even when the visibility problems described above do not exist, as shown in the example of Figure 10.
The software package presented in the previous sections was considered to perform the timing of the 65 LF episodes in the data set reported in Table 2. Moreover, for each episode, the heights of the LF at 1 min intervals were used to calculate the total fluid erupted volume from Equations (5) and (6), which included both gas and pyroclastics [8,31]:
U = 2 gH 0.5
V = U · A v · D  
In expression (5), U is the mean fluid exit velocity at the vent, H is the LF height, and g is the gravity acceleration, while in expression (6), V is the fluid volume (gas + pyroclastic) erupted by the LF, A v is the vent section area, and D is the duration of the LF in seconds. The vent surface area was calculated by assuming a circular vent with a diameter of 30 m [8] and supposed to be constant. Moreover, we have computed the volume V 2   of pyroclastics from the total erupted volume V (gas + pyroclastics), considering 0.18% as the ratio between the volumes of magma and fluids within the eruptive column as typical for Etna’s fountains [3].
According to expression (6), the estimated volumes depend on the mean fluid exit velocity U and LF duration D, for assigned values of the vent surface area. Thus, the performances of the automated approach will depend on its readability to estimate height H and duration D of the LF episodes. Concerning the reliability of the automated estimating H, the comparison with the corresponding manual reading, for a few episodes of the data set, is shown in Figure 11.
It can be seen in Figure 11 that manual height readings normally have a greater range than those automatically estimated. Here, it should be borne in mind that, as expressed in Equation (1), the automatically measured heights are a weighted average of the centroids, while those measured manually are normally taken as the maximum height of the lava fountain jet taken along the spreading direction. Considering that the heights represent the only geometric element on which the volumes of erupted material depend, it follows that with the automated estimation, these will normally be slightly less when compared to the manual ones, but with the advantage of immediacy. It is also necessary to bear in mind that for the purposes of estimating the volume of erupted matters, it is not so much the precise values of the instantaneous heights that are relevant, but their average value, which therefore also depends on the estimated duration for each individual episode. The comparison among the mean heights for the whole LF episodes after estimating the duration D by using both the CPD and Gaussian approaches are shown in Figure 12. In this figure, the abscissa is the integer N ranging from 1 to 66, i.e., number of LF episodes in the considered dataset (Table 2). The episode 32 is lacking because of poor visibility from all the monitoring cameras (Table 2).
In more detail, Figure 12 shows that for most of the episodes, there is a good agreement between the average mean heights estimated for each episode, not only between those obtained by using the CPD and Gaussian approaches, but also between these and those manually estimated. Altitudes obtained with the manual approach are, on average, 12% lower than those obtained by using the CPD approach, and 7% lower than the Gaussian one.
As regards to the duration D for each episode, the comparison between the automated and the manual estimation is shown in Figure 13.
Figure 13 shows that the durations obtained by using both the methods for automated timing are generally in good agreement with each other, as well as lower than those obtained manually. In more detail, duration manually estimated was on average about 34% higher, with a standard deviation of 77%, than those estimated by using the CPD, and about 3% higher, with a standard deviation of 93%, than the Gaussian one. To justify the discrepancy, it is worth noting that the automated approaches generally identify the sustained part of each LF, while the manual approach is not able to clearly distinguish the threshold of the intermediate phase preceding the sustained portion of the LF [3,8].
The comparison between estimated volumes and TADR by the three methods are reported in Figure 14a,b. A good agreement between the automated and manual estimation is apparent, bearing in mind that the automated volumes are usually smaller than the manual ones, because the durations refer to the sustained phase of the LF, while the manual and automated TADR are in good agreement because this feature is computed as the ratio between volumes and duration of each LF episode. In more detail, manual estimated volumes are, on average, about 26% higher, with a standard deviation of about 77% than those estimated by using the CPD approach, and 13% lower than the Gaussian, with a standard deviation of about 107%. The TADR manually estimated is about 15% lower that the CPD, with a standard deviation of 77% and 8% lower than the Gaussian, with a standard deviation of 33%.
Figure 15 shows the differences between the values obtained by the manual and automated routines, and indicates the good agreement between heights of the lava fountain (Figure 15a) and TADR (Figure 15d), and the discrepancies between duration (Figure 15b) and calculated erupted volume of pyroclastics (Figure 15c).
Considering the results obtained from the automated routines, the total volume of erupted pyroclastics during the 65 episodes computed by using the automated CPD approach was ~34.5 × 106 m3, with a minimum of ~0.012 × 106 m3 and a maximum of ~3.25 × 106 m3. The average duration of the sustained part of LF was of ~65.7 min (minimum 1 min, maximum 509 min). The TADR, calculated only for the pyroclastic portion of the episodes and for the whole duration of each paroxysmal event, was 159 m3 s−1 on average, spanning from a minimum of 84 m3 s−1 and a maximum of 379 m3 s−1; LF average height above the vents was 745 m, spanning from a peak value of 1834 m to a minimum of 223 m.

4. Discussion

In this paper, we have presented an automated routine that might help volcanic observatories such as INGV-OE to detect (1) the starting and ending time of an LF episode, (2) the LF heights, (3) the erupted volumes, and (4) the TADR, saving time and especially providing consistency and uniform data extraction from thermal monitoring videos. This would allow a prompt understanding of the state of the volcano, and of the magnitude and intensity of each explosive paroxysm as soon as it ended, allowing a timely volcanic hazard assessment. In addition, both our automated routines, based on the CPD and Gaussian interpolation, proved to be reliable in constraining the climax phase of the paroxysm leading to a sustained eruptive column, which is the phase posing the greater hazard for its impact on aviation and the population. Conversely, the manual analysis had clear difficulties in distinguishing the threshold between the intermediate phase and the LF sustained phase [3,8,30].
However, in order to routinely use the algorithms proposed here, it is necessary to overcome some limits that we described earlier and illustrated in Figure 6. The first shortcoming arose from the detection of unwanted objects falling in the field of view of the eruption, such as the Sun or the surfaces it irradiated. In distinguishing this anomalous pattern, the Gaussian interpolation might help, which would reveal and remove any deviation from the normal trend. A more common problem, and one that is hard to handle, is the cloud interference, with clouds resulting from water droplets, ash or gases filtering or obscuring the thermal images [36,37,38,39,40]. In the cases of clouds partially obscuring the field of view, it was still possible to interpolate the missing data, provided that they represent a small percentage of the total duration of the episode, which was the procedure also carried out with manual analysis. However, when the cloud cover was too continuous and extended in many directions, such as the episode of 27 May 2021 (Table 2), there was no way to retrieve any useful data, and an estimation of the erupted volume can only be performed by considering the timing obtained from the seismicity or from borehole strainmeters [17,43,44,45,47,48], and multiplying this for the average TADR estimated for each single episode occurring during the whole period lasting from 13 December 2020 to 21 February 2022. Thus, considering for the episode # 32 of 27 May 2021, the duration of 60 min (=3600 s) estimated from the seismicity, and multiplying this time for the TADR averaged over the 65 episodes (Table 2; ~146 m3 s−1), we obtained an estimated volume of ~0.53 × 106 m3, which was in line with, and slightly below, the average of the other LF events here considered (Table 2). This brought the total erupted volume of pyroclastics or tephra, erupted between 13 December 2020 and 21 February 2022, to ~65.2 × 106 m3.
Considering the manually obtained results, from a volcanological point of view, it is worth noting that, if we excluded the two outliers of 23 and 31 March 2021, which emitted more than 4 × 106 m3 of tephra (Table 2), the volume of pyroclastics erupted during the LF episodes, which occurred between mid-December 2020 and February 2022, increased with time (Figure 4), and also, the time-averaged discharge rate (TADR) increased with time (Figure 5). Figure 4 and Figure 5 display a wide variability of values, and although this variability might hide shorter eruptive cycles, it is however clear that the trend of TADR and erupted volume increased with time. These observations, i.e., of the increase of erupted volume of pyroclastics and rate of eruption with time, would suggest that the sequence of paroxysmal events was not yet over [31], and urges reliable and automated routines to be promptly developed, tested, and applied to the analysis of the LF episodes.

5. Conclusions

In this paper, we have presented the timing of start and end for 65 of the 66 LF episodes which took place at Etna volcano between 13 December 2020 and 21 February 2022, together with their duration, maximum and average LF heights, erupted volume of pyroclastics, and TADR (Table 2), obtained by manual analysis of the monitoring thermal images recorded by the INGV-OE network. We have then presented two automated routines, based on the CPD and Gaussian interpolation, that analyzed the thermal images and provided a fast and reliable way to obtain the same parameters acquired manually, in a timely way. The results obtained with the automated and manual routines are comparable (Figure 15), thus suggesting that a complete automation of the process is feasible. However, our analysis also highlighted important shortcomings arising from the presence of unwanted hot objects comprised in the field of view of the explosive episode that may lead to false results. Moreover, the presence of ash, weather, and gas clouds caused important interference with the data analysis, and might have reached the point of a complete obscuration of the field of view, as in the case of the episode #32 of 27 May 2021 (Table 2). We have shown that the Gaussian interpolation may limit the errors caused by a partial view, but more studies are necessary before this analysis can be routinely used for monitoring purposes. The results of our study showed an increasing magnitude (erupted volume) and intensity (TADR) of the explosive events in the period here considered (see Figure 4 and Figure 5), and this issue would urge a faster and reliable analysis to be obtained as soon as possible, thus motivating the work presented here.

Author Contributions

Conceptualization, S.C. and G.N.; methodology, S.C. and G.N.; software, G.N.; validation, S.C. and G.N.; formal analysis, S.C. and G.N.; investigation, S.C. and G.N.; resources, S.C. and G.N.; data curation, S.C. and G.N.; writing—original draft preparation, S.C. and G.N.; writing—review and editing, S.C. and G.N.; visualization, S.C. and G.N.; supervision, S.C. and G.N.; project administration, S.C.; funding acquisition, S.C. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Project FIRST—ForecastIng eRuptive activity at Stromboli volcano: timing, eruptive style, size, intensity, and duration, INGV-Progetto Strategico Dipartimento Vulcani 2019 (Delibera n. 144/2020; Scientific Responsibility: S.C.).

Data Availability Statement

The videos of eruptive activity used in this paper belong to the Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo – Sezione di Catania and are used for monitoring purposes. Selected videos can be made available upon request to the first author of this paper.

Acknowledgments

We would like to thank the INGV-OE scientists and technicians for the monitoring network maintenance, and especially Michele Prestifilippo for providing information essential for this work. We thank Stephan Conway for revising the English style.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Barreca, G.; Branca, S.; Corsaro, R.A.; Scarfì, L.; Cannavò, F.; Aloisi, M.; Monaco, C.; Faccenna, C. Slab detachment, mantle flow, and crustal collision in eastern Sicily (southern Italy): Implications on Mount Etna volcanism. Tectonics 2020, 39, e2020TC006188. [Google Scholar] [CrossRef]
  2. Bonaccorso, A.; Calvari, S. Major effusive eruptions and recent lava fountains: Balance between expected and erupted magma volumes at Etna volcano. Geophys. Res. Lett. 2013, 40, 6069–6073. [Google Scholar] [CrossRef]
  3. Calvari, S.; Cannavò, F.; Bonaccorso, A.; Spampinato, L.; Pellegrino, A.G. Paroxysmal Explosions, Lava Fountains and Ash Plumes at Etna Volcano: Eruptive Processes and Hazard Implications. Front. Earth Sci. 2018, 6, 107. [Google Scholar] [CrossRef]
  4. Calvari, S.; Bilotta, G.; Bonaccorso, A.; Caltabiano, T.; Cappello, A.; Corradino, C.; Del Negro, C.; Ganci, G.; Neri, M.; Pecora, E.; et al. The VEI 2 Christmas 2018 Etna Eruption: A Small but Intense Eruptive Event or the Starting Phase of a Larger One? Remote Sens. 2020, 12, 905. [Google Scholar] [CrossRef] [Green Version]
  5. Pailot-Bonnétat, S.; Harris, A.J.L.; Calvari, S.; De Michele, M.; Gurioli, L. Plume Height Time-Series Retrieval Using Shadow in Single Spatial Resolution Satellite Images. Remote Sens. 2020, 12, 3951. [Google Scholar] [CrossRef]
  6. Ganci, G.; Harris, A.J.L.; Del Negro, C.; Guéhenneux, Y.; Cappello, A.; Labazuy, P.; Calvari, S.; Gouhier, M. A year of lava fountaining at Etna: Volumes from SEVIRI. Geophys. Res. Lett. 2012, 39, L06305. [Google Scholar] [CrossRef]
  7. Ganci, G.; James, M.R.; Calvari, S.; Del Negro, C. Separating the thermal fingerprints of lava flows and simultaneous lava fountaining using ground-based thermal camera and SEVIRI measurements. Geophys. Res. Lett. 2013, 40, 5058–5063. [Google Scholar] [CrossRef]
  8. Calvari, S.; Salerno, G.G.; Spampinato, L.; Gouhier, M.; La Spina, A.; Pecora, E.; Harris, A.J.L.; Labazuy, P.; Biale, E.; Boschi, E. An unloading foam model to constrain Etna’s 11-13 January 2011 lava fountaining episode. J. Geophys. Res. 2011, 116, B11207. [Google Scholar] [CrossRef]
  9. De Beni, E.; Behncke, B.; Branca, S.; Nicolosi, I.; Carluccio, R.; Caracciolo D’Ajello, F.; Chiappini, M. The continuing story of Etna’s New Southeast Crater (2012–2014): Evolution and volume calculations based on field surveys and aerophotogrammetry. J. Volcanol. Geotherm. Res. 2015, 303, 175–186. [Google Scholar] [CrossRef]
  10. Andronico, D.; Cannata, A.; Di Grazia, G.; Ferrari, F. The 1986–2021 paroxysmal episodes at the summit craters of Mt. Etna: Insights into volcano dynamics and hazard. Earth-Sci. Rev. 2021, 220, 103686. [Google Scholar] [CrossRef]
  11. Freret-Lorgeril, V.; Donnadieu, F.; Scollo, S.; Provost, A.; Fréville, P.; Guéhenneux, Y.; Hervier, C.; Prestifilippo, M.; Coltelli, M. Mass eruption rates of tephra plumes during the 2011–2015 lava fountain paroxysms at Mt. Etna from doppler radar retrievals. Front. Earth Sci. 2018, 6, 73. [Google Scholar] [CrossRef]
  12. Corradini, S.; Guerrieri, L.; Lombardo, V.; Merucci, L.; Musacchio, M.; Prestifilippo, M.; Scollo, S.; Silvestri, M.; Spata, G.; Stelitano, D. Proximal Monitoring of the 2011–2015 Etna Lava Fountains Using MSG-SEVIRI Data. Geosciences 2018, 8, 140. [Google Scholar] [CrossRef] [Green Version]
  13. Neri, M.; Casu, F.; Acocella, V.; Solaro, G.; Pepe, S.; Berardino, P.; Sansosti, E.; Caltabiano, T.; Lundgren, P.; Lanari, R. Deformation and eruptions at Mt. Etna (Italy): A lesson from 15 years of observations. Geophys. Res. Lett. 2009, 36, L02309. [Google Scholar] [CrossRef]
  14. Guffanti, M.; Casadevall, T.J.; Budding, K. Encounters of Aircraft with Volcanic Ash Clouds: A Compilation of Known Incidents, 1953–2009; Version 1.0, Series 545, U.S. Geological Survey Data; U.S. Geological Survey: Reston, VA, USA, 2010; 12p., plus 4 appendixes including the compilation Database. Available online: http://pubs.usgs.gov/ds/545 (accessed on 22 March 2022).
  15. Small, C.; Naumann, T. The global distribution of human population and recent volcanism. Glob. Environ. Chang. Part B Environ. Hazards 2001, 3, 93–109. [Google Scholar] [CrossRef]
  16. Horwell, C.J.; Baxter, P.J. The respiratory health hazards of volcanic ash: A review for volcanic risk mitigation. Bull. Volcanol. 2006, 69, 1–24. [Google Scholar] [CrossRef]
  17. Pardini, F.; Corradini, S.; Costa, A.; Esposti Ongaro, T.; Merucci, L.; Neri, A.; Stelitano, D.; de’ Michieli Vitturi, M. Ensemble-Based Data Assimilation of Volcanic Ash Clouds from Satellite Observations: Application to the 24 December 2018 Mt. Etna Explosive Eruption. Atmosphere 2020, 11, 359. [Google Scholar] [CrossRef] [Green Version]
  18. Alparone, S.; Andronico, D.; Lodato, L.; Sgroi, T. Relationship between tremor and volcanic activity during the Southeast Crater eruption on Mount Etna in early 2000. J. Geophys. Res. 2003, 108, 2241. [Google Scholar] [CrossRef] [Green Version]
  19. Alparone, A.; Andronico, D.; Sgroi, T.; Ferrari, F.; Lodato, L.; Reitano, D. Alert system to mitigate tephra fallout hazards at Mt. Etna Volcano, Italy. Nat. Hazards 2007, 43, 333–350. [Google Scholar] [CrossRef]
  20. Cannavò, F.; Aranzulla, M.; Scollo, S.; Puglisi, G. A New GNSS-Based Approach for Volcanic Crater Location During Lava Fountains. IEEE Geosci. Remote Sens. Lett. 2019, 16, 697–701. [Google Scholar] [CrossRef]
  21. Cannavò, F.; Cannata, A.; Cassisi, C.; Di Grazia, G.; Montalto, P.; Prestifilippo, M.; Privitera, E.; Coltelli, M.; Gambino, S. A multivariate probabilistic graphical model for real-time volcano monitoring on Mount Etna. J. Geophys. Res. 2017, 122, 3480–3496. [Google Scholar] [CrossRef]
  22. Scollo, S.; Prestifilippo, M.; Pecora, E.; Corradini, S.; Merucci, L.; Spata, G.; Coltelli, M. Eruption column height estimation of the 2011–2013 Etna lava fountains. Ann. Geophys. 2014, 57, S0214. [Google Scholar] [CrossRef] [Green Version]
  23. Corradino, C.; Ganci, G.; Cappello, A.; Bilotta, G.; Calvari, S.; Del Negro, C. Recognizing Eruptions of Mount Etna through Machine Learning using Multiperspective Infrared Images. Remote Sens. 2020, 12, 970. [Google Scholar] [CrossRef] [Green Version]
  24. Walker, G.P.L. The Taupo pumice: Product of the most powerful known (Ultraplinian) eruption? J. Volcanol. Geother. Res. 1980, 8, 69–94. [Google Scholar] [CrossRef]
  25. Del Negro, C.; Cappello, A.; Ganci, G. Quantifying lava flow hazards in response to effusive eruption. Geol. Soc. Am. Bull. 2016, 128, 752–763. [Google Scholar] [CrossRef]
  26. Del Negro, C.; Cappello, A.; Bilotta, G.; Ganci, G.; Hérault, A.; Zago, V. Living at the edge of an active volcano: Risk from lava flows on Mt. Etna. Geol. Soc. Am. Bull. 2020, 132, 1615–1625. [Google Scholar] [CrossRef]
  27. Scollo, S.; Prestifilippo, M.; Spata, G.; D’Agostino, M.; Coltelli, M. Monitoring and forecasting Etna volcanic plumes. Nat. Hazards Earth Syst. Sci. 2009, 9, 1573–1585. Available online: http://www.nat-hazards-earth-syst-sci.net/9/1573/2009/ (accessed on 22 March 2022). [CrossRef] [Green Version]
  28. Scollo, S.; Prestifilippo, M.; Bonadonna, C.; Cioni, R.; Corradini, S.; Degruyter, W.; Rossi, E.; Silvestri, M.; Biale, E.; Carparelli, G.; et al. Near-Real-Time Tephra Fallout Assessment at Mt. Etna, Italy. Remote Sens. 2019, 11, 2987. [Google Scholar] [CrossRef] [Green Version]
  29. Scollo, S.; Boselli, A.; Corradini, S.; Leto, G.; Guerrieri, L.; Merucci, L.; Prestifilippo, M.; Zanmar Sanchez, R.; Sannino, A.; Stelitano, D. Multi-Sensor Analysis of a Weak and Long-Lasting Volcanic Plume Emission. Remote Sens. 2020, 12, 3866. [Google Scholar] [CrossRef]
  30. Calvari, S.; Bonaccorso, A.; Ganci, G. Anatomy of a Paroxysmal Lava Fountain at Etna Volcano: The Case of the 12 March, 2021, Episode. Remote Sens. 2021, 13, 3052. [Google Scholar] [CrossRef]
  31. Bonaccorso, A.; Calvari, S. A new approach to investigate an eruptive paroxysmal sequence using camera and strainmeter networks: Lessons from the 3–5 December2015 activity at Etna volcano. Earth Planet. Sci. Lett. 2017, 475, 231–241. [Google Scholar] [CrossRef]
  32. Sparks, R.S.J.; Bursik, M.I.; Carey, S.N.; Gilbert, J.S.; Glaze, L.S.; Sigurdsson, H.; Woods, A.W. Volcanic Plumes; John Wiley: New York, NY, USA, 1997; 574p. [Google Scholar]
  33. Bonadonna, C.; Phillips, J.C. Sedimentation from strong volcanic plumes. J. Geophys. Res. 2003, 108, 2340. [Google Scholar] [CrossRef]
  34. Orr, T.R.; Ulrich, G.E.; Heliker, C.; De Smither, L.G.; Hoffmann, J.P. The Pu‘u ‘Ō‘ō Eruption of Kīlauea Volcano, Hawai‘i—Episode 21 Through Early Episode 48, June 1984–April 1987; U.S. Geological Survey Scientific Investigations Report; United States Geological Survey: Reston, VA, USA, 2018; Volume 5109, 107p. [CrossRef]
  35. Carbone, D.; Zuccarello, L.; Messina, A.; Scollo, S.; Rymer, H. Balancing bulk gas accumulation and gas output before and during lava fountaining episodes at Mt. Etna. Sci. Rep. 2015, 5, 18049. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  36. Andò, B.; Pecora, E. An advanced video-based system for monitoring active volcanoes. Comput. Geosci. 2006, 32, 8591. [Google Scholar] [CrossRef]
  37. Ball, M.; Pinkerton, H. Factors affecting the accuracy of thermal imaging cameras in volcanology. J. Geophys. Res. 2006, 111, B11203. [Google Scholar] [CrossRef]
  38. Sawyer, G.M.; Burton, M.R. Effects of a volcanic plume on thermal imaging data. Geophys. Res. Lett. 2006, 33, L14311. [Google Scholar] [CrossRef]
  39. Spampinato, L.; Calvari, S.; Oppenheimer, C.; Boschi, E. Volcano surveillance using infrared cameras. Earth Sci. Rev. 2011, 106, 63–91. [Google Scholar] [CrossRef]
  40. Harris, A.J.L. Thermal Remote Sensing of Active Volcanoes; Cambridge University Press: Cambridge, UK, 2013; 717p. [Google Scholar]
  41. Sheldrake, T.E.; Sparks, R.S.J.; Cashman, K.V.; Wadge, G.; Aspinall, W.P. Similarities and differences in the historical records of lava dome-building volcanoes: Implications for understanding magmatic processes and eruption forecasting. Earth-Sci. Rev. 2016, 160, 240–263. [Google Scholar] [CrossRef] [Green Version]
  42. Kauahikaua, J.; Mangan, M.; Heliker, C.; Mattox, T. A quantitative look at the demise of a basaltic vent: The death of Kupaianaha, Kilauea Volcano, Hawai’i. Bull. Volcanol. 1996, 57, 641–648. [Google Scholar] [CrossRef]
  43. Cannata, A.; Montalto, P.; Patané, D. Joint analysis of infrasound and seismic signals by cross wavelet Atmospheric transform: Detection of Mt. Etna explosive activity. Nat. Hazards Earth Syst. Sci. 2013, 13, 1669–1677. [Google Scholar] [CrossRef] [Green Version]
  44. Cannata, A.; Montalto, P.; Privitera, E.; Russo, G.; Gresta, S. Tracking eruptive phenomena by infrasound: May 13, 2008 eruption at Mt. Etna. Geophys. Res. Lett. 2009, 36, L05304. [Google Scholar] [CrossRef]
  45. Nunnari, G. Clustering activity at Mt Etna based on volcanic tremor: A case study. Earth Sci. Inform. 2021, 14, 1121–1143. [Google Scholar] [CrossRef]
  46. Moschella, S.; Cannata, A.; Di Grazia, G.; Gresta, S. Insights into lava fountain eruptions at Mt. Etna by improved source location of the volcanic tremor. Ann. Geophys. 2018, 61, VO446. [Google Scholar] [CrossRef]
  47. Bonaccorso, A.; Currenti, G.; Linde, A.; Sacks, S.; Sicali, A. Advances in Understanding Intrusive, Explosive and Effusive Processes as Revealed by the Borehole Dilatometer Network on Mt. Etna Volcano. Front. Earth Sci. 2020, 7, 357. [Google Scholar] [CrossRef] [Green Version]
  48. Bonaccorso, A.; Carleo, L.; Currenti, G.; Sicali, A. Magma Migration at Shallower Levels and Lava Fountains Sequence as Revealed by Borehole Dilatometers on Etna Volcano. Front. Earth Sci. 2021, 9, 740505. [Google Scholar] [CrossRef]
  49. Barsotti, S.; Andronico, D.; Neri, A.; Del Carlo, P.; Baxter, P.J.; Aspinall, W.P.; Hincks, T. Quantitative assessment of volcanic ash hazards for health and infrastructure at Mt. Etna (Italy) by numerical simulation. J. Volcanol. Geotherm. Res. 2010, 192, 85–96. [Google Scholar] [CrossRef]
  50. Horwell, C.J.; Sargent, P.; Andronico, D.; Lo Castro, M.D.; Tomatis, M.; Hillman, S.E.; Michnowicz, S.A.K.; Fubini, B. The iron-catalysed surface reactivity and health-pertinent physical characteristics of explosive volcanic ash from Mt. Etna, Italy. J. Appl. Volcanol. 2017, 6, 12. [Google Scholar] [CrossRef] [Green Version]
  51. Corradino, C.; Amato, E.; Torrisi, F.; Calvari, S.; Del Negro, C. Classifying Major Explosions and Paroxysms at Stromboli Volcano (Italy) From Space. Remote Sens. 2021, 13, 4080. [Google Scholar] [CrossRef]
  52. Fenner, D.; Rümpker, G.; Li, W.; Chakraborty, M.; Faber, J.; Köhler, J.; Stöcker, H.; Srivastava, N. Automated Seismo-Volcanic Event Detection Applied to Stromboli (Italy). Front. Earth Sci. 2022, 10, 809037. [Google Scholar] [CrossRef]
  53. Harris, A.J.L.; Dehn, J.; Calvari, S. Lava effusion rate definition and measurement: A review. Bull. Volcanol. 2007, 70, 1–22. [Google Scholar] [CrossRef]
  54. Nunnari, G.; Cannavó, F. Online detection of offsets in GPS time series. Earth Sci. Inform. 2021, 14, 267–276. [Google Scholar] [CrossRef]
  55. Killick, R.; Fearnkead, P.; Eckley, I. Optimal detection of changepoints with linear computational cost. J. Am. Stat. Assoc. 2012, 107, 1590–1598. [Google Scholar] [CrossRef]
  56. Lavielle, M. Using penalized contrasts for the change-point problem. Signal. Process. 2015, 85, 1501–1510. [Google Scholar] [CrossRef] [Green Version]
Figure 1. (a) Google map of Sicily, with the white rectangle showing the area magnified in (b), which is the Google map of Etna volcano showing the position of the INGV thermal monitoring cameras used in this study. The red circle indicates the position of the summit craters producing the lava fountain activity.
Figure 1. (a) Google map of Sicily, with the white rectangle showing the area magnified in (b), which is the Google map of Etna volcano showing the position of the INGV thermal monitoring cameras used in this study. The red circle indicates the position of the summit craters producing the lava fountain activity.
Remotesensing 14 02392 g001
Figure 2. (a) The original RGB thermal image from the EMOT camera, with temperatures in °C comprised between blue (0 °C) and white (>100 °C); (b) the cropped grayscale; (c) the binarized image, and (d) the binarized image with computed centroids (red asterisks), representing the center of mass of the hot objects detected by the threshold process.
Figure 2. (a) The original RGB thermal image from the EMOT camera, with temperatures in °C comprised between blue (0 °C) and white (>100 °C); (b) the cropped grayscale; (c) the binarized image, and (d) the binarized image with computed centroids (red asterisks), representing the center of mass of the hot objects detected by the threshold process.
Remotesensing 14 02392 g002
Figure 3. Area (a) and mean altitude (b) of the lava fountain episode at Mt Etna on 24 February 2021.
Figure 3. Area (a) and mean altitude (b) of the lava fountain episode at Mt Etna on 24 February 2021.
Remotesensing 14 02392 g003
Figure 4. Volume of pyroclastics (y-axis, expressed in m3) against time (x-axis, date) for the 65 explosive episodes listed in Table 2, together with a linear trend and its formula.
Figure 4. Volume of pyroclastics (y-axis, expressed in m3) against time (x-axis, date) for the 65 explosive episodes listed in Table 2, together with a linear trend and its formula.
Remotesensing 14 02392 g004
Figure 5. Time-averaged discharge rate of pyroclastics (TADR, y-axis, expressed in m3 s−1) erupted during the lava fountain episodes against time (x-axis, date). The graph displays the data for the 65 explosive episodes listed in Table 2, together with a linear trend and its formula.
Figure 5. Time-averaged discharge rate of pyroclastics (TADR, y-axis, expressed in m3 s−1) erupted during the lava fountain episodes against time (x-axis, date). The graph displays the data for the 65 explosive episodes listed in Table 2, together with a linear trend and its formula.
Remotesensing 14 02392 g005
Figure 6. (a) Effect of the Sun reflection on the ground slope (green area), EMOT camera; (b) cooling lava flow (white area), EMOT camera; (c) Effect of the Sun on the clouds (pink area), ENT camera; (d) the Sun in the camera field of view (white circle), EMCT camera; (e) Combined effects: lava fountain in presence of the Sun in the field of view, EMCT camera; (f) lava fountain and Sun reflection on the vegetation in the foreground, EBT camera. (a1f1) are the corresponding gray images, and (a2f2) are the corresponding black and white binarized images.
Figure 6. (a) Effect of the Sun reflection on the ground slope (green area), EMOT camera; (b) cooling lava flow (white area), EMOT camera; (c) Effect of the Sun on the clouds (pink area), ENT camera; (d) the Sun in the camera field of view (white circle), EMCT camera; (e) Combined effects: lava fountain in presence of the Sun in the field of view, EMCT camera; (f) lava fountain and Sun reflection on the vegetation in the foreground, EBT camera. (a1f1) are the corresponding gray images, and (a2f2) are the corresponding black and white binarized images.
Remotesensing 14 02392 g006
Figure 7. (a) Area (in pixel2) and (b) mean altitude (in pixels above the crater rim) against time of the lava fountain on 20–21 February 2021 at Mt. Etna, retrieved from the EMCT camera.
Figure 7. (a) Area (in pixel2) and (b) mean altitude (in pixels above the crater rim) against time of the lava fountain on 20–21 February 2021 at Mt. Etna, retrieved from the EMCT camera.
Remotesensing 14 02392 g007
Figure 8. (a) Area (in pixel2) and (b) mean altitude above the crater rim (in pixels, measured above the crater rim) with the timing (gray vertical lines) of start and end of the lava fountain episode. See text for explanation.
Figure 8. (a) Area (in pixel2) and (b) mean altitude above the crater rim (in pixels, measured above the crater rim) with the timing (gray vertical lines) of start and end of the lava fountain episode. See text for explanation.
Remotesensing 14 02392 g008
Figure 9. (a) Area (in pixel2) and mean altitude (in pixels, measured above the crater rim) against time of the 13 December 2020 lava fountain episode as retrieved from the ENT camera, and (b) Gaussian interpolation (red line) of the pyroclastic area against time of the 13 December 2020 lava fountain episode retrieved from the ENT camera.
Figure 9. (a) Area (in pixel2) and mean altitude (in pixels, measured above the crater rim) against time of the 13 December 2020 lava fountain episode as retrieved from the ENT camera, and (b) Gaussian interpolation (red line) of the pyroclastic area against time of the 13 December 2020 lava fountain episode retrieved from the ENT camera.
Remotesensing 14 02392 g009
Figure 10. Graphs showing the timing of the LF episode N. 11, which occurred on 22–23 February 2021 (Table 2) obtained from (a) the CPD method and from (b) the Gaussian method. The black vertical lines indicate the start and end time, and the red line in (b) is the Gaussian interpolation. Table 2 shows that the manual method indicated the LF episode took place between 21:17 on 22 February and 00:03 of 23 February 2021. In accordance with the manual method, (a) shows that the CPD automated timing indicated that the sustained phase of the LF developed between 23:26 and 23:55 on 22 February, and similarly, the Gaussian method (b) indicated a timing comprised between 23:21 and 23:55 on 22 February.
Figure 10. Graphs showing the timing of the LF episode N. 11, which occurred on 22–23 February 2021 (Table 2) obtained from (a) the CPD method and from (b) the Gaussian method. The black vertical lines indicate the start and end time, and the red line in (b) is the Gaussian interpolation. Table 2 shows that the manual method indicated the LF episode took place between 21:17 on 22 February and 00:03 of 23 February 2021. In accordance with the manual method, (a) shows that the CPD automated timing indicated that the sustained phase of the LF developed between 23:26 and 23:55 on 22 February, and similarly, the Gaussian method (b) indicated a timing comprised between 23:21 and 23:55 on 22 February.
Remotesensing 14 02392 g010
Figure 11. Time series of heights estimated by using the automated and the manual approaches for a few episodes of the data set against time, with the blue line for the automated, and the orange line for the manual detection: (a) Episode 7, 17–18 February 2021; (b) Episode 12, 28 February 2021; (c) Episode 65, 10 February 2022, and (d) Episode 66, 21 February 2022.
Figure 11. Time series of heights estimated by using the automated and the manual approaches for a few episodes of the data set against time, with the blue line for the automated, and the orange line for the manual detection: (a) Episode 7, 17–18 February 2021; (b) Episode 12, 28 February 2021; (c) Episode 65, 10 February 2022, and (d) Episode 66, 21 February 2022.
Remotesensing 14 02392 g011
Figure 12. Mean altitude for each of the 65 lava fountain episodes estimated after having established the timing of start and end for each recorded time series by using both the CPD (blue line) and the Gaussian (red line) approach. Their values were compared with the corresponding manually estimated values (yellow line). Episode 32 was lacking because of poor visibility from the cameras (see Table 2).
Figure 12. Mean altitude for each of the 65 lava fountain episodes estimated after having established the timing of start and end for each recorded time series by using both the CPD (blue line) and the Gaussian (red line) approach. Their values were compared with the corresponding manually estimated values (yellow line). Episode 32 was lacking because of poor visibility from the cameras (see Table 2).
Remotesensing 14 02392 g012
Figure 13. Estimated duration performed after timing the LF episodes by the CPD (blue line) and Gaussian (red line) approaches, in comparison with the manual (black line) ones.
Figure 13. Estimated duration performed after timing the LF episodes by the CPD (blue line) and Gaussian (red line) approaches, in comparison with the manual (black line) ones.
Remotesensing 14 02392 g013
Figure 14. (a) Comparison among volumes, and (b) comparison among TADR obtained from the manual and automated approaches.
Figure 14. (a) Comparison among volumes, and (b) comparison among TADR obtained from the manual and automated approaches.
Remotesensing 14 02392 g014
Figure 15. Comparison among automated and manual analysis of (a) heights of the lava fountain; (b) duration of each episode; (c) estimated erupted volume of pyroclastics, and (d) TADR (time-averaged discharge rate [53]).
Figure 15. Comparison among automated and manual analysis of (a) heights of the lava fountain; (b) duration of each episode; (c) estimated erupted volume of pyroclastics, and (d) TADR (time-averaged discharge rate [53]).
Remotesensing 14 02392 g015
Table 1. List of the INGV monitoring thermal cameras used in this paper and their main features. The field of view is considered at the crater rim.
Table 1. List of the INGV monitoring thermal cameras used in this paper and their main features. The field of view is considered at the crater rim.
LabelTypeLocationDistance from the Craters (km)Frame RateField of View
ENTThermal
FLIR A40M
Nicolosi, South flank
730 m a.s.l.
15.02 frames/s320 × 240 pixels
EBTThermal
FLIR A320
Bronte, NW flank
85 m a.s.l.
13.52 frames/s25° × 18.8°
EMCTThermal
FLIR A320
Mt. Cagliato, East flank
1160 m a.s.l.
8.32 frames/s320 × 240 pixels
EMOTThermal
FLIR A320
Montagnola, South flank
2600 m a.s.l.
3.01 frame/s320 × 240 pixels
Table 2. Paroxysmal explosive episodes occurring at Etna volcano between 2020 and 2022 (updated to 23 February 2022). The first column reports the episode’s progressive number. All times are expressed in UTC. LF = lava fountain. TADR = time averaged discharge rate [53]. Data from episode # 32 are lacking because there was no visibility from any of the thermal monitoring cameras.
Table 2. Paroxysmal explosive episodes occurring at Etna volcano between 2020 and 2022 (updated to 23 February 2022). The first column reports the episode’s progressive number. All times are expressed in UTC. LF = lava fountain. TADR = time averaged discharge rate [53]. Data from episode # 32 are lacking because there was no visibility from any of the thermal monitoring cameras.
Ep. #DateStarting Time
(hh:mm)
Ending Time
(hh:mm)
Duration (in Minutes and in Seconds)Label of Cam-era UsedMax LF
Height
(m above the Crater)
Mean LF
Height
(m above the Crater)
LF Volume
(×106 m3)
TADR
(m3 s−1)
113 December 202022:0022:4848 min.
2880 s
ENT5142310.2484
213–14 December 202023:5800:1113 min.
780 s
ENT4002350.0790
314 December 202001:0201:4038 min.
2280 s
ENT2861210.1359
421 December 202009:1109:5948 min.
2880 s
EBT308012960.58201
522 December 202003:0505:13128 min.
7680 s
ENT8002950.7293
618 January 202119:3821:0385 min.
5100 s
ENT10673430.52101
716 February 202116:1117:0251 min.
3060 s
EMCT15607570.46150
817–18 February 202122:3200:51139 min.
8340 s
EMCT12303580.8298
919 February 202108:1610:06110 min.
6600 s
EMCT13654920.78118
1020–21 February 202121:3201:15223 min.
13,380 s
EMCT15003861.43107
1122–23 February 202121:1700:03166 min.
9960 s
ENT
EMOT EMCT
36676861.19120
1223 February 202103:5204:5058 min.
3480 s
ENT
EMOT EMCT
9003370.3292
1324 February 202118:5621:41165 min.
9900 s
ENT18006491.37139
1428 February 202107:3108:3463 min.
3780 s
ENT360013760.70185
152 March 202111:2314:50207 min.
12,420 s
EMOT6062781.0383
164 March 202101:3004:10160 min.
9600 s
ENT6002040.7578
174 March 202107:1109:32141 min.
8460 s
ENT323312751.58186
187 March 202104:1007:01171 min.
10,260 s
EMOT
EBT
40006381.07104
199–10 March 202123:5502:46171 min.
10,260 s
ENT
EMCT
18606551.44140
2012 March 202107:3509:45130 min.
7800 s
ENT
EBT
240011491.63209
2114–15 March 202123:1201:42150 min.
9000 s
ENT13336700.5359
2217 March 202101:3004:57207 min.
12,420 s
ENT15335381.58128
2319 March 202108:1810:13115 min.
6900 s
EMOT6291710.75108
2423–24 March 202121:3308:19646 min.
38,760 s
EMCT13334564.56118
2531 March– 1 April 202119:3008:53803 min.
48,180 s
EMCT6302414.1085
2619 May 202100:5004:25215 min.
12,900 s
ENT6674821.59124
2721 May 202100:5002:44114 min.
6840 s
EMCT15336830.99145
2822 May 202120:2722:08101 min.
6060 s
ENT12006490.87143
2924 May 202120:2521:4984 min.
5040 s
ENT14678310.82162
3025 May 202118:2018:5333 min.
1980 s
ENT5333170.20102
3126 May 202110:2011:1050 min.
3000 s
ENT12676270.41137
3227 May 202112:0013:0060 min.
3600 s
EMCTPoor visibilityPoor visibilityPoor visibilityPoor visibility
3328 May 202106:3007:2757 min.
3420 s
ENT8004330.40116
3428 May 202116:0516:116 min.
360 s
ENT3332950.04112
3528 May 202119:4820:5062 min.
3720 s
ENT10006010.51138
3630 May 202104:2005:4484 min.
5040 s
ENT10005890.69137
372 June 202108:3010:46136 min.
8160 s
EBT16409241.38170
384 June 202116:1217:4088 min.
5280 s
EMCT11706650.76143
3912 June 202120:0021:46106 min.
6360 s
EMCT8104380.73115
4014 June 202121:1522:2166 min.
3960 s
EMCT8704190.45112
4116 June 202111:3712:3861 min.
3660 s
ENT17336730.52142
4217 June 202122:4023:5575 min.
4500 s
EMCT11404040.50111
4319 June 202118:4719:3548 min.
2880 s
EMCT11405720.38131
4420 June 202122:4023:4464 min.
3840 s
ENT24678920.63165
4522 June 202102:3003:4575 min.
4500 s
ENT20008480.71159
4623 June 202102:4403:1935 min.
2100 s
ENT186710350.38183
4723 June 202118:0019:1272 min.
4320 s
ENT293313010.60153
4824 June 202109:4510:2641 min.
2460 s
ENT
EMOT
17338250.39157
4925 June 202100:3801:4870 min.
4200 s
ENT23338760.64153
5025 June 202118:4019:2040 min.
2400 s
ENT11336910.36149
5126 June 202116:0016:3838 min.
2280 s
ENT16007720.35155
5227 June 202108:5309:4350 min.
3000 s
ENT16006740.43143
5328 June 202114:2515:3065 min.
3900 s
EBT239012110.75193
541–2 July 202122:5000:2797 min.
5820 s
ENT18008040.91156
554 July 202115:1516:5095 min.
5700 s
ENT14678730.94164
566 July 202122:1623:4488 min.
5280 s
EBT327016731.20227
578 July 202120:3522:1297 min.
5820 s
EBT271012421.12192
5814 July 202110:4512:30105 min.
6300 s
EBT223010971.15183
5920 July 202105:1008:11181 min.
10,860 s
EBT351015332.26208
6031 July 202119:4423:37233 min.
13,980 s
EBT383014732.83202
619 August 202102:0704:11124 min.
7440 s
EBT239012801.48199
6229 August 202116:2417:5591 min.
5460 s
EBT239013101.10202
6321 September 202107:2108:3574 min.
4440 s
ENT233312340.88199
6423 October 202108:2010:17117 min.
7020 s
ENT400018441.63232
6510 February 202218:4021:56196 min.
11,760 s
ENT571421602.88245
6621 February 202211:1112:5099 min.
5940 s
ENT405718651.39234
Average
Duration (min./s)
Average Max LF
height
(m)
Average Mean LF
height
(m)
Average
LF volume
(×106 m3)
Average
TADR
(m3 s−1)
120/7171 18157840.993144.75
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Calvari, S.; Nunnari, G. Comparison between Automated and Manual Detection of Lava Fountains from Fixed Monitoring Thermal Cameras at Etna Volcano, Italy. Remote Sens. 2022, 14, 2392. https://doi.org/10.3390/rs14102392

AMA Style

Calvari S, Nunnari G. Comparison between Automated and Manual Detection of Lava Fountains from Fixed Monitoring Thermal Cameras at Etna Volcano, Italy. Remote Sensing. 2022; 14(10):2392. https://doi.org/10.3390/rs14102392

Chicago/Turabian Style

Calvari, Sonia, and Giuseppe Nunnari. 2022. "Comparison between Automated and Manual Detection of Lava Fountains from Fixed Monitoring Thermal Cameras at Etna Volcano, Italy" Remote Sensing 14, no. 10: 2392. https://doi.org/10.3390/rs14102392

APA Style

Calvari, S., & Nunnari, G. (2022). Comparison between Automated and Manual Detection of Lava Fountains from Fixed Monitoring Thermal Cameras at Etna Volcano, Italy. Remote Sensing, 14(10), 2392. https://doi.org/10.3390/rs14102392

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop