**2. Materials and Methods**

### *2.1. UV Camera System*

In the framework of the ERC-funded project BRIDGE (www.bridge.unipa.it), we designed a multi-instrument UV-absorption spectroscopy system for robust SO2 flux measurements (Figure S1). The system is composed of (i) an instrument module and (ii) an acquisition/processing module. The instrument module is equipped with two JAI CM-140GE-UV cameras sensible to UV-radiation, and one Ocean-Optics USB2000+ Spectrometer coupled to a telescope of rectangular, vertically-oriented Field Of View (FOV ≈ 0.3◦ × 14◦), and is spatially filtered to match the ≈12◦ vertical width. Two different band-pass optical filters with Full Width at Half Maximum (FWHM) of 10 nm, and a central wavelength of 310 and 330 nm, respectively, are applied in front of the cameras to enhance differential UV absorption in the SO2 bandwidth [45,46]. In addition, 520 × 676 pixel images are acquired at 10-bit resolution with a frame rate of 0.5 Hz.

To obtain a quantitative measure of SO2 column density within the volcanic cloud, we calculate the proportionality ratio between absorbance and SO2 concentration in a defined region of the image pointed by the Ocean-Optic USB2000+ Spectrometer [47]. The use of the UV spectrometer allows us to quantitatively measure the full UV spectrum, and then fit the theoretical SO2 absorption cross-section [48] with the differential absorption between two consecutive spectra (acquired every 5 seconds). The Ocean-Optic USB2000+ Spectrometer in use has on-board a Sony ILX511B Linear Silicon CCD Array Detector at 2048 pixels, with a wavelength response of 200–1100 nm, a dynamic range of 8.5 × 108, and a SNR of 250:1 at full signal. Calibrated SO2 column densities over the entire images are then obtained by integrating images achieved by the UV camera with information achieved by the spectrometer. The instrument module is powered with a 12 V power supply, and requires 15 W in a fully operational mode. A Fujitsu RX100 Workstation, connected to the instrument module, automatically acquires synchronous data from the instrument module, and processes data without the need of the operator. To do this, we designed algorithms to control acquisition and processing parameters, such as the automatic tuning of camera's and spectrometer's exposures, and automatic evaluation of optimal viewing conditions (see Section 2.3.2). The computer internal time drift is controlled by a specifically designed application that reads the time-stamp from an NMEA (National Marine Electronics Association) standard message coming from a GPS antenna. The instrument module communicated with the acquisition/processing module via wired or wireless TCP/IP connection.

This UV camera system was installed at the Montagnola site (on Etna Volcano, Figure 1), and designed to stream real-time SO2 flux results using a Wi-Fi data link. The objective is to capture SO2 emissions associated with diverse volcanic processes and dynamics, including quiescent (passive) degassing, explosive eruptions (strombolian activity/lava fountaining), and effusive eruptions [49]. Montagnola is located at ~3 km distance from the active summit vents and grants perfect views of the southern sector of the summit crater area (Figure 1).

#### *2.2. Seismic and Thermal Data*

We compared our SO2 fluxes with other independent geophysical parameters such as tremor and thermal radiance primary to have a benchmark in the calibration of automatic SO2 flux calculation algorithm. Last but not least, we combined these data to characterize volcanic activity.

We used seismic data recorded by the ETN station [42,50], located at Lapide Malerba, at 5 km from the summit area. ETN is equipped with a broad-band seismometer (Guralp CMG-40T, with a sensitivity of 800 V/(m/s) and eigen period of 30 s). The link between SO2 flux and volcanic tremor at Mt. Etna [22,51,52] suggests that the tremor is generated by the same degassing dynamics. We then calculated volcanic tremor amplitude from raw traces recorded at ETN station, by averaging within a 1-minute length window the maximum RMS amplitude taken within a 1-s window.

Thermal remote sensing offers a great opportunity to follow volcanic unrest from ground and space to characterize volcanic activity in near-real time [53] and to estimate near vent large pyroclastic products and lava flow discharged during eruptions. We used satellite data from the MIROVA system [43,54] and ground-based thermal cameras [44] to constrain onset, duration, and intensity through the time of eruptive events occurring at Etna during 2016. In particular, MIROVA uses the data provided by the MODIS sensor, which acquires four images per day (two daytime and two nighttime) with a spatial resolution of 1 km in the infrared bands [43]. The heat flux retrieved from MIROVA data, called Volcanic Radiative Power (VRP), is a combined measurement of the area and the integrated temperature of the hot emitters (hot vents, lava flows, etc.) at the time of each image acquisition. These data are provided without correction for cloud cover and satellite view geometry [43]. These factors may introduce noise in the dataset, which has been proven not to affect the general trend associated with the activity of Mt. Etna [43]. A mask 5 × 5 km around the volcano summit is used to filter out thermal anomalies due to wildfires.

#### *2.3. The Algorithm for Automatic Processing of UV Camera Data*

Manual processing of UV camera results [33] has the advantage that images are checked and validated manually by an operator. However, manual procedures are time-consuming, especially when dealing with huge data flows from permanent monitoring stations. To overcome this limit, we designed a MATLAB-based algorithm to automatically process images, and thus obtain SO2 fluxes in near real-time. Using a "low-cost" PC-workstation (with 8 GB RAM and Xeon E3 type CPU), we can successfully process data at ~5× speed (i.e., a 5 minutes-long record is processed in ~1 min), which allows nearly real-time monitoring. The automatic routine has been calibrated using the results of the manual procedure and includes: (1) automatic determination of background absorbance levels, (2) automatic determination of image goodness, using image quality indexes, (3) estimation of gas plume speed and its distribution across the crater area, and (4) calculation of SO2 flux distribution throughout the summit crater area. These are described in more detail below.

### 2.3.1. Image Processing and SO2 Column Densities

Relative UV absorption by volcanic plume SO2 is quantified by applying the Beer–Lambert law. Sets of synchronous images, taken by the two co-aligned cameras using different filters, are combined to obtain single absorbance images. This method, known as the "double filter method" [12,45], implies the use of two cameras with different filters, with one centered at 310 nm and the other centered on 330 nm, where UV radiation is/is not absorbed, respectively. The use of the two filters method allows compensating for aerosol attenuation/backscattering, to avoid any temporal mismatch associated with filter change while using a single camera, and to maintain the sampling rate at up to 0.5 Hz (two synchronous images every 2 s taken by two cameras [45]).

In our automatic routine, once a new raw image is acquired, a first quality check is made by the system, in order to keep the best exposure times to compensate any subtle changes in sunlight intensity. An automatic real-time tuning of exposure time is also applied to the UV spectrometer data, in order to obtain the best measurement dynamics within the UV bandwidth.

Synchronous images from the two cameras are then real-time corrected for vignetting effects associated with filters and optics, and normalized for relative exposure times. Residual intensities are then combined to obtain an un-calibrated absorbance image using the Lambert-Beer Law equation.

$$A = -\log\_{10} \frac{I^{310}}{I^{330}} - A\_0 \tag{1}$$

where *A* is the absorbance, *I* <sup>310</sup> and *I* <sup>330</sup> are pixel intensities associated with cameras mounting the 310 or 330 nm filter, while *A*<sup>0</sup> is the absorbance level associated with a clear background sky sub-area of the image (assumed to be unaffected by SO2 absorption, *I*<sup>0</sup> <sup>310</sup> and *I*<sup>0</sup> <sup>330</sup> of Figure 2) and calculated as −log10(*I* <sup>310</sup>/*I* 330), following Kern [55].

In the automatic processing module, this background sky sub-area is automatically selected for each image by monitoring a distal sky horizontal section with respect to the vent position, and selecting the sector with the lowest absorbance intensity.

Residual absorbance is then converted into SO2 column density integrating data from the co-located ultraviolet spectrometer, which is pointing to a known sub-area within the camera field of view. This procedure yields, in real-time, the proportionality ratio between absorbance and SO2 column densities using the method described in McGonigle [47].

**Figure 1.** (**a**) The Montagnola site where the UV camera system was installed (position of ETN seismic station is also indicated). Sites from where pictures (**c**), (**d**), and (**e**) have been taken are shown. Black inset identifies the zoomed area of Figure 1b. (**b**) Etna summit area (redrawn from Reference [56]) shows the active summit craters (BN: Bocca Nuova; VOR: Voragine; NEC: North-East Crater; SEC: South-East Crater; NSEC: New South East Crater), the 7 August degassing vent, and the graben-like structure discussed in the text. The site where picture f has been taken is also shown. (**c**) Photo showing the Montagnola site and the summit area. (**d**) Thermal snapshot from INGV-OE monitoring camera capturing the 18 May lava fountaining episode. (**e**) Vigorous degassing from the vent opened on 7 August. (**f**) Picture taken from Reference [56] showing BN crater collapse occurred on 10 October, and marking the end of enhanced degassing activity (see text).

**Figure 2.** Calculation of the SO2 column densities using the dual camera method with 330 nm (**a**) and 310 nm (**b**) optical filters. Absorbance image (**c**) is obtained using the Lambert-Beer equation after image normalization with respect to background absorbance intensities (black circles in (**a**) and (**b**), automatically located). Black circles in (**a**) and (**b**) represent the sky areas where absorbance is assumed to be SO2 free (minimum absorbance level). Black rectangles (*Ssky Sgrn*) in (**a**) represent the areas used for calculation of the visibility index. The red circle in (**c**) shows the FOV of co-located UV-scanning spectrometer used to convert un-calibrated absorbance intensities into SO2 column densities.

#### 2.3.2. Automatic Determination of the Optimal Viewing Condition

The optimal plume viewing conditions, and the presence of a clear sky, are required for reliable SO2 density measurements. However, weather conditions are extremely variable on Etna's summit, and often prevent optimal SO2 observation. To minimize uncertainties due to poor weather conditions, we set-up a sub-routine for real-time calculation of two visibility indexes. The first visibility index (Fog index) is calculated as the unsigned ratio between the mean pixel intensity associated with the camera FOV's portions capturing sky and ground, respectively. Tests we conducted on real and synthetic images show that the higher this ratio is, the better the visibility condition is. SO2 measurements are then selected by setting a threshold on the visibility index, and discarding measurements below the threshold (e.g., biased by a poor visibility condition).

Detecting a "sky" signal well above the "ground" signal (Figure 2) is a required but not sufficient condition for reliable SO2 measurements. This is especially true in the presence of a highly condensed plume, where the SO2 absorbance signal can be masked. Thus, a second automatic procedure was developed and run in real-time, which allows us to select only images with a clear SO2 signal above atmospheric noise.

This latter procedure is based on the principle of combining absorbance and 310 nm images associated with the plume. A well detected and measurable SO2 signal requires that lower intensities in the 310 nm image are measured in the plume relative to its surroundings (because SO2 is absorbing solar radiation), and that higher intensities are consistently obtained in the absorbance image (see the Lambert Beer equation). This condition is only verified if the SO2 signal is high enough to emerge above the atmospheric noise. To discriminate this condition in real-time, we defined a correlation index (Figure 3) as the correlation coefficient between absorbance and 310-nm pixel intensities over a cross-section intersecting the plume (Figures 3 and 4). The correlation coefficient is defined as *C(i,j)*/*((C(i,i)\*C(j,j))ˆ(1*/*2)*, where *C* is the covariance matrix, *<sup>i</sup>* and *<sup>j</sup>* are pixel intensities over the cross-section of absorbance and 310-nm images, respectively. In such plots, the closer the correlation coefficient is to the value −1, the more absorbance can be related to gas. Images that do not satisfy this condition (e.g., that have a Correlation Index < −0.5) are disregarded by the automatic computation (Figure 4).

**Figure 3.** The image quality calculation method using the correlation coefficient between the 310 nm filter image (**a**) and the corresponding absorbance image (**b**). An intensity profile associated with a section (dashed line) crossing the volcanic plume, for both the 310-nm filter and absorbance images, is obtained (**c**). If these profiles are negatively correlated with a high correlation coefficient (**d**), then the gas is visible within the plume.

**Figure 4.** Example of output of the quality indexes sub-routine. The visibility (fog) and correlation indexes fluctuate through time as visibility conditions change (see snapshots on top of the figure). Gas is visible only when the fog index is greater than 4 and the correlation index is less than −0.5.

#### 2.3.3. Plume Velocity Field

A robust plume velocity field is mandatory for reliable SO2 flux measurements. Errors in plume velocity have shown to contribute to 40% or more of the overall error in the determined fluxes [13,57].

The UV camera approach offers the unique opportunity to track the gas while dispersing right after atmospheric emission, which minimizes errors in plume speed determination of yet more established DOAS and COSPEC methods. These methods indirectly infer plume speed from either on-site measurement of wind velocity or from meteorological models [11,15].

The UV camera approach allows us to derive the velocity profile over the summit craters by applying an optical flow algorithm that tracks gas fronts in consecutive frames [58].

Optical flow consists ofthe apparent motion pattern of image objects between two consecutive frames, caused by the movement of either the object or the camera, and is valid under the assumptions that the pixel intensities of an object do not change significantly between consecutive frames, and that the neighboring pixels have similar motion.

If a pixel, with intensity *I*(*x*,*y*,*t*), where (*x*,*y*) are the pixel coordinates and *t* is the time in first frame, moves by distance (*dx*,*dy*) in the next frame taken after time *dt*, it can be assumed that:

$$I(\mathbf{x}, y, t) = I(\mathbf{x} + d\mathbf{x}, y + dy, t + dt) \tag{2}$$

Then, from the Taylor series approximation of the right-hand side, removing common terms and dividing by dt, one gets the following equations.

$$
\rho \frac{\partial I}{\partial x} u + \frac{\partial I}{\partial y} v + \frac{\partial I}{\partial t} = 0 \tag{3}
$$

$$u = \frac{d\mathbf{x}}{dt}; v = \frac{dy}{dt} \tag{4}$$

where <sup>∂</sup>*<sup>I</sup>* <sup>∂</sup>*<sup>x</sup>* and <sup>∂</sup>*<sup>I</sup>* <sup>∂</sup>*<sup>y</sup>* are the image gradients, <sup>∂</sup>*<sup>I</sup>* <sup>∂</sup>*<sup>t</sup>* is the gradient along time, and u and v are horizontal and vertical velocities that are unknown. Lucas & Kanade [59] provide a method (LK) to derive these unknown velocities, by solving the basic optical flow Equations (3) and (4) for all the pixels using the least squares criterion, and by combining information from nearby pixels. We then applied the LK algorithm, included within the Open-CV toolbox, to our dataset [60]. We tested the performance of this method by applying it to artificial images with known particle velocities. The method has successfully determined velocity field with an error of <5%.

Absorbance images, obtained using the two filters method, contain gas-rich and ash-free portions of the plume with a higher absorbance relative to the background, and/or to ash-rich or particle-rich plume portions. We exploit this feature to track only gas moving fronts in consecutive frames by filtering them from other moving features such as lapilli and ash and by applying the LK method to the absorbance images rather than to the raw images directly acquired by our dual camera system. Velocities are then calculated by selecting the best features to track within the image, and which correspond to the areas with the highest pixel intensities (i.e., high SO2 column densities) and high spatial coherence in consecutive frames (taken every 2 s).
