We present data from two separate measurement campaigns at Stromboli volcano in Italy: September 9–11, 2019, and September 23–25, 2021 (Figure 1). Thermal emission measurements were performed using a Bruker EM27 OP-FTIR spectrometer lent by the Osservatorio Etneo dell’Istituto Nazionale di Geofisica e Vulcanologia (INGV). The EM27 consists of a “Rock Solid interferometer” (corner cube mirrors on a pendulum around the beam splitter) directing radiation towards a Stirling-cooled Mercury-Cadmium-Telluride (MCT) detector using entirely reflective optics (i.e., no external telescope was added) with an optical path difference of 1.8 cm and a field-of view (FOV) of 30 mrad. Spectra were acquired in the frequency range 600—2,000 cm⁻¹ with a resolution of ∼0.5 cm⁻¹ (or a ∼0.25 cm⁻¹ spectral sampling), and averaged with five co-adds, resulting in a sampling interval of 8–10 s. Radiometric calibration was performed using an internal blackbody target with adjustable temperature. Using two blackbody spectra at two different temperatures (typically around 20°C and 40°C), we created an instrument calibration function which can be applied to convert the single beam spectrum to calibrated radiance or brightness temperature. Clear sky spectra were collected by pointing the instrument towards an area of sky outside of the plume while keeping the inclination angle as close to that of the plume measurements as possible. The duration of acquisition periods varied between 30 min and several hours, and calibration and clear sky measurements were performed every hour or so to account for the thermal drift of the instrument and to reflect changing atmospheric conditions. In datasets for which multiple calibration sets were acquired, we account for instrument drift by interpolating the blackbody spectra between the calibration sets, thereby creating a time-dependent calibration function with unique values for each measurement. We also create a similar time-dependent representation of the expected background spectrum by interpolating between clear sky spectra. This is especially important when working with longer datasets or with datasets acquired at the time of local sunrise and sunset, when the thermal profile of the atmosphere changes rapidly.

During the first field campaign in 2019, measurements were performed from the l’Osservatorio restaurant, providing a direct view of the craters (Figure 1). Wind direction was relatively constant during the 3-day period, and the plume drifted over the Sciara del Fuoco in a NE direction. We collected multiple short sets of 30–60 min of measurements, moving the instrument field-of-view to intersect the plume at various distances from the vent. The activity during this campaign was typical for Stromboli, with explosive events occurring every 5 min on average. Events at the crater were logged manually in a notebook, with the observer entering the time, nature of the event and crater of origin as far as it was possible to observe. We also organized a series of flights with a Remotely Piloted Aircraft System (RPAS) developed and built at the University of Manchester, with the aim of sampling ash from the drifting volcanic plume during the measurement period (Figure 2). The RPAS was a quadcopter carrying a hoover-like sampling mechanism. Once the RPAS had reached the desired location within the plume, the sampling mechanism was triggered remotely, turning on a high mass-flow ducted fan to direct airflow and ash into a small vacuum cleaner bag connected to the system. We were able to collect a single sample at 17:02 local time on 11 September 2019 (analysis results in Section 3.1). The particle size distribution (PSD) of the sample collected during the RPAS flight was measured using a Retch Camsizer X2 particle size analyser at the University of Leeds. Further, the sample was mounted onto SEM stubs using clear epoxy, polished and carbon-coated for analysis with an Electron Microprobe Analyser (EMPA) at the University of Oxford, to determine the chemical composition of individual phases.

In 2021, measurements were performed from the roof of the Pedra Residence hotel in Stromboli village. The wind direction was approximately to the East, and the line of sight of the instrument intersected the plume as it passed over the summit of the island. We collected longer datasets (3-4 h) with the specific purpose of evaluating the usefulness of the method as a monitoring tool. A GoPro camera was used to capture time-lapse imagery from the observer’s vantage point, and the eruptive events were logged by reviewing the footage. Although the method allows for a continuous record of events over long periods, weaker events and those without an associated ash-rich plume were often not detectable in the footage because we did not have a direct view of the craters. Therefore, we could not assign a specific vent to any given event either. During this second campaign, we also collected simultaneous spectra with a UV spectrometer to validate the infrared SO₂ retrieval. Light was collected using a collimated telescope (diameter = 25.4 mm, f = 100 mm), connected to an Ocean Optics USB 2000+ spectrometer via an optical fibre. The spectrometer was controlled and powered by a connected laptop via USB cable. SO₂ slant column densities were retrieved using the iFit method (Esse et al., 2020). The wavelength window for the retrieval was 310—320 nm and the fit included absorption cross-sections for SO₂ at 295 K (Rufus et al., 2003) and O₃ at 223 K (Serdyuchenko et al., 2014), a Ring spectrum, a wavelength shift and stretch and a linear intensity offset. The instrument line shape (ILS) was characterised using a super-Gaussian with the shape parameters also fitted for each spectrum to account for changes with time (e.g., due to changes in temperature).

The emission spectrum retrieval algorithm (developed in Python and available at https://www.github.com/jfsmekens/plumeIR) follows the basic principles laid out by Love et al. (1998), where the simplified radiative transfer expression comprises three layers (Figure 3): 1) a lower layer between the observer and the plume, 2) a plume layer, and 3) an upper layer encompassing the atmosphere in the line-of-sight behind the plume. Given a plume height (h) and a vertical plume thickness (z _v ), the plume layer is defined between the height h 1 = h − z v / 2 and h 2 = h + z v / 2 . The forward model computes the radiance difference between a clear sky and a plume measurement, which offers the advantage of isolating the spectral features of the volcanic species. The aim is to provide a tool which can be used in the field to quantify target species in real time as the measurements are taken (i.e., with a target processing time of <10 s for each individual spectrum, comparable to the acquisition time between consecutive spectra). The algorithm was developed around an iterative non-linear least-squares fitting method. High-resolution radiative transfer calculations (.04 cm⁻¹ spectral sampling) over the entire vertical atmospheric profile are performed in a pre-processing stage to generate reference spectra used in the forward model. In contrast to recent efforts undertaken at Popocatépetl (Stremme et al., 2012; Taquet et al., 2017), where the effects of particulate species are simulated using a polynomial, we use a broad retrieval window (700–1,300 cm⁻¹) and model the radiation and attenuation due to particulate species from Mie theory, allowing us to retrieve SCDs for those species as well. The radiance difference method presents the advantage of emphasizing the spectral features due to the presence of a volcanic layer, thereby minimising the influence of atmospheric variables such as temperature and water vapour volume mixing ratios at higher altitude in the vertical profile.

The ash sample collected during the first measurement campaign in 2019 was analysed for particle size and chemical composition. This section details the results of these analyses, which were used to a) guide the selection of the most appropriate ash refractive indices in the retrieval, and b) compare particle size information between the sample and the values retrieved by the algorithm. The sample’s mass-equivalent PSD (Figure 6) is bimodal, with a population of coarser particles (200–400 μm) and a separate population of finer particles (70–100 μm), seemingly lacking very fine particles. This PSD is generally similar to that of fallout samples collected at the summit of Stromboli in 2015 (Freret-Lorgeril et al., 2019). However, this representation can be misleading, as the Camsizer measurement describes the volume fraction of particles based on their equivalent diameter (X _area : diameter of the area equivalent circle of each particle projection). In these terms, the larger particles take an outsized representation, as the very fine fraction (<10 μm) represents a negligible volume within the overall sample, despite representing the vast majority of the particles. In order to emphasize this effect, we also show the PSD in terms of fractional number density (N-density) of particles by size—calculated based on the measured volume fraction and assuming spherical particles - which is heavily shifted towards smaller particle size bins. The resulting N-density PSD is also bimodal with a large peak around 5–10 μm, and a second one around 65 μm. We argue that the relatively small number of very large particles has a limited effect on the overall optical properties of the ash within the FOV of the OP-FTIR instrument, especially for dilute plumes, and that the comparatively large number of finer particles dictates the shape of the overall extinction profile. Therefore, we believe the N-density PSD is more directly comparable with the sizes retrieved in this work, as the strongest effects on the overall shape of the particle extinction coefficients from Mie scattering largely depend on the number density of particles with diameters close to that of the measured wavelength (i.e., 8–12 μm, see Figure 5).

The chemical composition of several crystalline phases, as well as the glass matrix as determined via EMPA analysis are presented in Table. 2. The glass composition represents a basaltic trachyandesite, and the main crystalline phases are olivine (Ol), clinopyroxene (Cpx) and plagioclase (Plag). Using the backscattered images from the EPMA analysis, we estimated the relative surface area fractions for each phase (Glass: 60%; Ol: 5%; Cpx: 10%, Plag: 25%) and used these relative fractions to calculate an approximate bulk composition for the sample. This bulk composition (a trachybasalt with ∼51 wt% SiO₂) was used to parameterize the refractive index of the ash during the retrieval (see sect. 2).

Figure 7 illustrates typical results of individual fits for selected spectra. Over the broad fitting window, the radiance difference spectra are dominated by a broadband contour representing the water vapour continuum. This contour can be positive or negative in absolute value, representing either an increase or a decrease in total water vapour column between the plume and clear sky measurements. The dashed line in the top panel of each figure emphasises this effect and shows the expected difference before considering the contributions from the plume. It is obtained by computing the radiance difference using best-fit parameters for water vapour scaling and observer temperature only and omitting all plume species. The shape of this baseline spectrum at the start of measurements depends on the respective positions of the line-of-sight between clear sky and plume and should approach zero for a clear sky taken at the plume location. Note that this radiance difference already contains recognisable spectral features for H₂O (continuum and narrow absorption/emission lines), O₃ and CO₂ resulting from the overall change in total atmosphere transmission when simply moving the line-of-sight of the instrument. In Figure 7A, showing a spectrum acquired in 2019, O₃ (1,000–1,080 cm⁻¹) and CO₂ (925–1,000 cm⁻¹) appear as emission lines (due to the reduced transmission associated with a H₂O scaling factor <1), and the baseline minimum value goes from −1 to −6 mW/(m²⋅sr⋅cm⁻¹) over the course of ∼30 min. This value is expected to gradually change over the course of measurements, as the background atmospheric conditions evolve. In contrast, Figure 7B shows a spectrum acquired in 2021, where the baseline is a positive radiance difference (H₂O scaling factor >1) and O₃ and CO₂ appear as absorption features. The fitting window was deliberately chosen to include areas on either side where the atmosphere is virtually opaque (<730 cm⁻¹ and >1,270 cm⁻¹). The radiance at those wavelengths represents the blackbody emission at the temperature of the most proximal layer, and we can use the value of the radiance difference to retrieve the temperature of the proximal layer (layer 1).

The bottom panel in each figure shows the contributions of each plume species (including the ad hoc corrections for atmospheric gases) superimposed over the baseline radiance difference. Each individual contribution is computed using the forward model and zeroing all quantities within the plume except for the species of interest. Note that because the background and foreground atmosphere are always an inherent part of the forward model, water vapour absorption lines still appear in the contribution from each individual species. The most instantly recognisable feature is the emission associated with SO₂. In particulate-rich spectra, a strong O₃ absorption is visible, related to the severe decrease in transmission introduced by the heavy particle burden and often requiring an ad hoc correction.

Figure 8 shows the time series of SO₂ SCDs retrieved from our FTIR measurements for one dataset acquired during our second measurement campaign in 2021. We also present SO₂ SCDs measured using a co-located UV spectrometer. Disagreement in both absolute values of the SCD and timing of individual peaks are to be expected. They may arise as a consequence of differences in a series of factors between the two methods, such as: 1) the size of the respective FOVs (.03 rad for the IR, .1 rad for the UV); 2) the alignment of the telescopes; 3) the integration time (0.2 s in the UV, 8.9 s in the IR); and 4) different random and systematic error sources between the two methods. The UV time series was smoothed using a kernel of length equal to the IR integration time and resampled to match the temporal x-axis in the IR time series. Moreover, we determined the optimal lag using a cross-correlation method and shifted the UV time series accordingly. We would expect this lag to vary over the course of the measurements, as it is tied to the plume velocity, and plumes associated with crater explosions will travel at greater velocities. However, we found that an overall lag of ∼71 s resulted in a good match between the timing of the main events. Retrieved SCDs between the IR and UV measurements generally agree within a factor of ∼2 (scatter plot in Figure 8A). The main discrepancies occur towards the beginning and end of datasets with the IR method retrieving systematically lower SCDs than the UV method at the beginning of the dataset in Figure 8A, and systematically higher SCDs towards the end of the dataset.

Figure 8 further documents explosive events, as recorded in GoPro time-lapse footage (orange stars). Here we only report ash explosions, defined as type 1 events following the terminology in Patrick et al. (2007), and mark them at the time of their first appearance above the crater rim from our vantage point. Type 2 events (essentially bubble bursts associated with ballistics) were much more difficult to identify from the 2021 vantage point without a direct view of the craters (Figure 1). It is also impossible to assign a specific crater of origin for individual events for the 2021 time series. All documented events are followed by a marked increase in retrieved SO₂, SA and ash contents, after a delay of 1-2 min representing the time needed for the plume to move into the FOV of the instruments.

During the 2019 campaign, we performed measurements from the l’Osservatorio restaurant, with a direct view of the crater area (Figure 1B). We were able to move the FOV of the instrument to intersect the plume directly above the vent. A time series of the retrieved quantities for a dataset collected in this configuration is shown in Figure 9. During this period, we captured the baseline passive degassing plume, as well as plumes resulting from explosive events. Events were documented through visual observation and logged manually. As we had a direct view of the crater area, the explosive plumes can be separated into two classes with relation to the event types defined by Patrick et al. (2007). Type 1 events (bubble bursts + ballistics) produce plumes with low ash content and elevated concentrations of SO₂ (green stars). Type 2 events, on the other hand, are ash-rich with a visibly dark plume (orange stars). The FOV of the instrument was positioned directly above the NE crater area. As a result, this dataset mainly contains events originating from the NE and central vents (which were directly behind the NE vents along the line of sight). Some emissions from the SW crater area may have entered our FOV but, given the wind direction during our measurements, this was challenging to determine visually. Type 1 events originated from both NE and C areas, and type 2 events mainly originated from the NE craters. Absolute retrieved SCDs for SO₂, ash and sulphate aerosols above the vents are relatively high (peak SO₂: 6 × 10¹⁸ molec⋅cm⁻², peak ash: 1.75 g⋅m⁻², peak SA: 1.2 g⋅m⁻²). Peaks in SO₂ within the time series generally occur shortly after the onset of each event. Retrieved SA particle sizes are shown in the bottom plot in Figures 9E, F. The probability density function, where each measurement of particle size is weighed according to the associated SCD, exhibits a mode around .8 μm, representative of the typical size retrieved during low level degassing (mean d _eff value of 1.79 ± 3.08 μm). There are also several peaks at higher diameters, representing d _eff values retrieved during the explosive events (i.e., type 2 ca 18:56, 18:58 and 19:02 and type 1 ca 18:36 and 18:43, Figure 9E).

In a scatter plot between SO₂ and SA over the entire dataset (Figure 9B), the measurements with lower amounts of particulates (light shades) define a ratio line. We calculate SO₂/SA mass ratios using a robust linear regression (Siegel, 1982) which emphasises the overall trend rather than more extreme data points, in order to minimise the effects of those measurements with high amounts of particulates during the explosive events. This is done so we can quantify the sulphur partitioning in the gas plume during passive degassing periods. We observe a ratio of ∼53 by mass in this dataset. Values from similar datasets range between 40 and 80.

Retrieved ash particle sizes are shown in the middle plot in Figures 9C, D. The retrieved size for ash particles increases during explosive events. The probability density function exhibits a clear mode at ∼10 μm during non-explosive phases (range between 5 and 15 μm). During the most intense ash-rich type 2 events (e.g., ca 18:30–18:33 and ca 18:42 in Figure 10), ash d _eff reaches much higher values, up to 200 μm, the upper bound for ash size set in the algorithm. During type 2 events of lower intensity (e.g., three events in quick succession ca. 18:56, 18:58 and 19:00), particle size also increases during each event, but to lower values of 10–30 μm.

Figure 10 shows a time series of measurements acquired from the same vantage point and on the same day (11 September 2019) as those measurements shown in Figure 9, but where the instrument FOV intersected the plume ∼0.8 km downwind (Figure 1B). The measured SCDs for all target species (peak SO₂: 2 × 10¹⁸ molec⋅cm⁻², peak ash: .3 g⋅m⁻², peak SA: .1 g⋅m⁻²) are lower than in the near-vent measurements, representing the dilution of the plume and loss processes (e.g., SO₂ oxidation and settling of PM) as it travels. Wind speed on that day was relatively low, and this distance corresponds to a plume age of approximately 3–5 min. Recorded type 1 and type 2 events at the vents appear disassociated with peaks in the SO₂ time series. However, we identify three individual events within the time series: one event with elevated SO₂ and SA, with a slow onset at ca. 15:07, assumed to be a type 2 event), and two events with sharper onsets at ca. 15:14 and 15:21 showing elevated SO₂, SA and ash contents (assumed to be type 2 events). The mean SO₂/SA mass ratio for this particular dataset is ∼17, and values for downwind datasets range between 10 and 30.

Retrieved sizes for the SA downwind are 5.21 ± 4.27 μm during low level degassing, and 2.19 ± 0.36 μm on average during explosive events (type 1 and type 2). Sizes retrieved during a type 1 event with a high SA SCD (ca 15:07–15:10 in Figure 10) show elevated values, reaching the upper bound for SA size set in the algorithm (10 μm). Retrieved sizes for ash in the downwind dataset are 6.66 ± 0.57 μm outside of ash-rich events, and 7.68 ± 0.38 during ash-rich events. Both values are in good agreement with the smaller mode of the particle size distribution found in our collected sample (Figure 6).

The principal source of uncertainty in the method comes from the estimation of plume temperature in the FTIR forward model. This is illustrated by systematic disagreements between the UV and IR retrievals of SO₂ SCDs towards the beginning and end of datasets (see Figure 8A). As the measurement period evolves, the actual temperature at plume height increases or decreases (depending on time of day), while the assumed plume temperature in our model remains the same (extracted from the reference atmospheric sounding, which is taken only in 12 h intervals, and at a location >100 km away from the measurement location). If the assumed temperature is colder than the plume actually is, the amount of SO₂ retrieved by our method will be an overestimate of the actual SO₂ in the plume. Conversely, if the assumed temperature is warmer than the actual temperature, our measurements will be an underestimate. UV spectroscopy is not affected by this assumption and thus differences may arise between the retrieved quantities from the FTIR and the UV measurements. A further complication is that both plume height and atmospheric temperature at plume height may vary over the course of the measurements, leading to systematic errors in all quantities retrieved by the method. As the radiance difference is calculated using a clear sky spectrum acquired at the beginning or end of the measurements, exacerbating the differences in measured radiance which result from changing atmospheric conditions in spectra collected further away from the calibration period, using more frequent clear sky measurements will yield more consistent results. In addition, plume temperature might be higher due to the presence of the volcanic plume, i.e., the actual temperature in layer 2 may differ from the one recorded in the atmospheric profile due to the hot volcanic gas and particles.

To quantify the error linked with the plume layer temperature, we introduced a fixed temperature difference (Δt _plume ) between the temperature in the volcanic layer at the time of measurement and the clean air layer temperature in the clear sky measurement and performed the retrieval on a subset of 1,000 measurements from 25 September 2021, for various values of Δt _plume (Figure 11). The temperature of the plume at Stromboli is frequently recorded using a Multi-component Gas Analyser System (Multi-GAS) at the summit (Aiuppa et al., 2009; Aiuppa et al., 2010). Although temperature measurements are not reported in published studies, they are recorded by the instrument and typically fall within 1–3 degrees of the ambient temperature when the plume drifts over the summit (INGV personal communication). Therefore, we vary Δt _plume between values of −5 K and + 5 K. In addition, we ran the retrieval on the same data subset leaving the temperature of the volcanic layer as a free parameter. This latter approach yielded unrealistic results, with Δt _plume ranging from −30 K to + 30 K within the same dataset, leading to overestimates of up to a factor of four in the retrieved SO₂ SCD in spectra with the highest amounts of particulates. However, these results allow us to quantify the expected errors related to plume temperature. When plotting the ratio between the SO₂ SCD retrieved with a free t _plume and the SO₂ SCD retrieved with a fixed t _plume against the associated Δt _plume , they outline a very clear relationship, which can be approximated using an exponential decay function. We estimate that uncertainties in assumed plume temperature will lead to errors in the retrieved SO₂ SCDs of a factor of ∼3.5% per K. When considering a realistic range of assumed plume temperatures, this translates to overestimates of up to 19% (Δt _plume = −5 K) and underestimates of up 16% (Δt _plume = + 5 K). Measurements with high SCDs of particulates depart from this overall relationship and introduce larger overestimates. Similar values are found when examining other retrieved quantities (summarised in Table 3). These uncertainties are much larger than the typical fitting error of ∼2.2% and represent the largest limitation in the method. Local measurements of temperature at plume height and accurate visual determination of the plume height will help mitigate the potential for error.

Previous measurements have shown the PM-related optical properties of the Stromboli plume to vary with the occurrence of individual events (Sellitto et al., 2021). Our results further support these observations, with type 1 events generally associated with a peak in SO₂ SCD, and an increase in SA SCD, whereas type 2 events are associated with increased amounts of both ash and SA SCD (Figure 9). Some type 1 and type 2 events were not accompanied by particle increases and in one case no SO₂ peak was observed. This may reflect the fact that some plumes drifted above the FOV of our instrument, either because the explosion that produced them injected significant momentum into the plumes, or because shifts in wind direction or speed caused the plume to rise more vertically in some cases than others. As the dataset shown in Figure 9 was acquired in the run up to sunset (approximately 19:15 local time), it is also possible that our manual observations are incomplete or inaccurate, particularly towards the end of the dataset as it became darker.

The presence of high amounts of particulates within a plume (e.g., ash-bearing plumes after explosions or those with a high optical thickness due to large amounts of aerosols) can negatively affect UV measurements and lead to large errors in retrieved SO₂ SCDs if realistic radiative transfer is not taken into account (Kern et al., 2010; Kern et al., 2012; Kern et al., 2013; Varnam et al., 2020). Although both over- and underestimations are possible, depending on the specific conditions within a plume (i.e., magnitude of the SO₂ column amounts and nature of the scattering aerosol), underestimates are more likely in plumes with high optical thicknesses, particularly if the particulate species within the plume also exhibit strong absorbing properties, as is the case for ash-laden plumes captured close to the vent (Kern et al., 2013). In such cases, IR retrievals may offer an advantage, as they are less likely to be impacted by multiple scatterings (typical sizes for volcanic particles are often much smaller than the 8–12 μm wavelength range of the method) and by light dilution effects. Our method offers a means of quantifying the aerosol optical thickness and SO₂ SCDs from a single measurement. For optically thick plumes containing very large amounts of ash or PM, underestimates in our IR retrievals should still be expected, as the radiance collected at the instrument would originate mainly from the surface of the plume. However, for moderate particle concentrations, our method offers a way to compensate for the expected loss of radiance from within the plume layer. This is apparent at several points in the time series presented in Figure 8 and highlighted in subpanels B-C, where ash-rich plumes are generally associated with higher retrieved SO₂ SCDs in the IR, by factors ranging from 2 to 5. We propose that this represents an underestimate of the true SO₂ burden by the UV method, and that the SCDs retrieved in the IR are closer to the concentrations within the plume at the time of measurement. Alternatively, it is possible that plumes with high amounts of ash retain heat more efficiently, and as a result the assumption of a plume in thermal equilibrium with the atmosphere is no longer valid. In such a case, the retrieved SO₂ SCD with the IR method would constitute an overestimate of the real values, as warmer plumes emit more radiation. These discrepancies merit further investigation within a dedicated, systematic side-by-side experiment in which multiple scattering and light dilution effects are accounted for in both UV and IR retrievals.

Table 4 summarises sulphur speciation values previously reported and measured in proximal plumes (<1 km from source) at Stromboli and elsewhere. Sulphur partitioning between the gas and aerosol phases is generally reported as the molar ratio of SO₄ ^2- ions to SO₂ (SO₄ ^2-/SO₂) when measured with filter packs. Alternatively, when retrieved with spectroscopic techniques, the reported value is often the mass ratio of the gas phase over the aerosol phase, including the water (SO₂/SA). The conversion between the two values depends on the acidity of the aerosols themselves, which is not always reported. For a simple comparison, we assume an acidity of 65 wt% H₂SO₄ (the value used to generate our reference extinction spectra), and present both ratios. The SO₂/SA measured in our near-vent measurements is on the lower side of the range observed in proximal plumes at other volcanoes, suggesting the presence of large amounts of primary sulphate aerosols (i.e., those condensed directly from the high-temperature gas emitted from the magma) in the passive plume and during type 1 events (which typically have lower particulate loads than type 2 events). Alternatively, this could reflect the rapid formation of sulphate via oxidation of SO₂ condensing onto the primary particles in the plume or a higher water content of the SA in our measurements. This is supported by the observations that the sizes measured here (Figure 9) are rather large compared to those expected for newly formed particles (<.1 μm for particles within the “nucleation mode”) (Whitby, 1978; Mather et al., 2003b).

SO₂/SA ratios in diluted downwind datasets (range 10–30, mean of ∼17 in the dataset illustrated in Figure 10) are lower than those found near-vent (range 40–80, mean of ∼53 in the dataset illustrated in Figure 9), suggesting either an increase in sulphate aerosol mass and/or a decrease in SO₂ mass as the plume ages. Several mechanisms exist which could explain the observed decrease in SO₂/SA ratio with plume age. SO₂ depletion by oxidation into sulphuric acid is a commonly proposed mechanism (Eatough et al., 1994) and has been documented in volcanic plumes, with SO₂ loss rates in tropospheric plumes ranging between 10⁻⁷ and 10⁻³ s⁻¹ (Oppenheimer et al., 1998b; McGonigle et al., 2004; Rodríguez et al., 2008). As SO₂ is oxidized, the decrease in SO₂/SA ratio then results from both a loss of gaseous SO₂ in the plume, and from condensation of sulphate formed as SO₂ is oxidized, either increasing the size of existing particles, or nucleating new particles, and ultimately causing an increase in aerosol mass. However, even assuming a fast depletion rate at the upper bound of the range documented above (10⁻³ s⁻¹), it is not clear that this oxidation process would be fast enough to explain the changes observed here, which took place over distances of <1 km and for plume ages of only a few minutes. Alternatively, the mass of SA may also be increased by condensation of water into the aerosol, a process by which the existing sulphate aerosols lower the supersaturation threshold required to grow droplets (Mather et al., 2004a; Mather et al., 2004b), potentially ultimately acting as seeds for the formation of clouds (Pianezze et al., 2019). This process would increase the overall mass of the SA while decreasing their relative acidity. Given the fixed parameter for SA acidity in our retrieval, our model is not able to measure this change in acidity. However, hygroscopic growth is consistent with both the slight increase in retrieved SA size (during periods of passive degassing) and the decrease in SO₂/SA mass ratio in the downwind measurements and is therefore a plausible mechanism to explain our observations. Two recent studies have investigated the evolution of SO₂ and SA in the Stromboli and Etna plume (Pianezze et al., 2019; Sahyoun et al., 2019). Both studies found an increase in the number concentration of SA particles with distance from the vent, combined with an increase in the size of the measured particles. Although the mechanisms investigated by these studies take place over much larger distances and targeting older plumes, their interpretations are broadly consistent with our observations, and suggest the mechanisms of particle growth that ultimately can result in cloud condensation occur even in very young plumes such as those we measured. Further investigation of these mechanisms might elucidate volcanic impacts and interactions with cloud cover in the troposphere (e.g., Ebmeier et al., 2014).

As mentioned previously, the sensitivity of the retrieval to particle sizes above 20 μm is not well established, and it is possible that the retrieved values do not represent the true size within the plume. Nevertheless, our observations of significant increases in optical depth in the long-wave infrared spectrum and rapid shifts in particle size during type 2 events are consistent with the increase in aerosol optical depth and decrease in the Angstrom exponent observed at shorter wavelengths by Sellitto et al. (2021) and are the likely the result of an injection of coarser particles from the explosions. Sizes retrieved over the entirety of the dataset capture the full range of ash particle sizes found in the sample collected with the RPAS, with the fine fraction often found during low-level degassing and the coarser particles during type 2 events.

The size of ash retrieved in diluted plumes (either downwind or above the vent during non-explosive activity) is remarkably consistent, with values of 5–10 μm in both datasets. This value is also consistent with the fine size fraction in our collected sample. Coarser particles (65–200 μm) retrieved from near-vent measurements are not detected in the downwind measurements, even when a high SCD of ash is retrieved, suggesting that they have already been removed by sedimentation at that point. This is consistent with our visual observations during the measurement period of large “fingers” of ash raining from the drifting plumes and depositing the coarser particles on the slopes of the volcano, as documented in Freret-Lorgeril et al. (2020). It is also consistent with the observations of Sellitto et al. (2020) who found a decrease in aerosol size with increasing distance from the vent in the proximal (<20 km) plume at Etna volcano, which they interpreted as due to the sedimentation of ash (and possibly coarser SA). It should be noted that the sensitivity of the IR retrieval to particle size is most pronounced in the range of particle sizes 1–10 μm, where the size of the particles approaches the wavelength of radiation. The spectral shape for ash shows relatively little variation with increasing sizes above 20 μm (see Figure 5), and the detection of coarser particle sizes in this work should be taken as an indication of their presence rather than an accurate determination of their actual size.