Tephra Mass Eruption Rate From Ground-Based X-Band and L-Band Microwave Radars During the November 23, 2013, Etna Paroxysm

The morning of November 23, 2013, a lava fountain formed from the New South-East Crater (NSEC) of Mt. Etna (Italy), one of the most active volcanoes in Europe. The explosive activity was observed from two ground-based radars, the X-band polarimetric scanning and the L-band Doppler fixed-pointing, as well as from a thermal-infrared camera. Taking advantage of the capability of the microwave radars to probe the volcanic plume and extending the volcanic ash radar retrieval (VARR) methodology, we estimate the mass eruption rate (MER) using three main techniques, namely surface-flux approach (SFA), mass continuity-based approach (MCA), and top-plume approach (TPA), as well as provide a quantitative evaluation of their uncertainty. Estimated exit velocities are between 160 and 230 m/s in the paroxysmal phase. The intercomparison between the SFA, MCA, and TPA methods, in terms of retrieved MER, shows a fairly good consistency with values up to <inline-formula> <tex-math notation="LaTeX">$2.4\times 10^{6}$ </tex-math></inline-formula> kg/s. The estimated total erupted mass (TEM) is <inline-formula> <tex-math notation="LaTeX">$3.8\times 10^{9}$ </tex-math></inline-formula>, <inline-formula> <tex-math notation="LaTeX">$3.9\times 10^{9}$ </tex-math></inline-formula>, and <inline-formula> <tex-math notation="LaTeX">$4.7\times 10^{9}$ </tex-math></inline-formula> kg for SFA with L-band, X-band, and thermal-infrared camera, respectively. Estimated TEM is between <inline-formula> <tex-math notation="LaTeX">$1.7\times 10^{9}$ </tex-math></inline-formula> kg and <inline-formula> <tex-math notation="LaTeX">$4.3\times 10^{9}$ </tex-math></inline-formula> for TPA methods and <inline-formula> <tex-math notation="LaTeX">$3.9\times 10^{9}$ </tex-math></inline-formula> kg for the MCA technique. The SFA, MCA, and TPA results for TEM are in fairly good agreement with independent evaluations derived from ground collection of tephra deposit and estimated to be between <inline-formula> <tex-math notation="LaTeX">$1.3\,\,\pm \,\,1.1\times 10^{9}$ </tex-math></inline-formula> and <inline-formula> <tex-math notation="LaTeX">$5.7\times 10^{9}$ </tex-math></inline-formula> kg. This article shows that complementary strategies of ground-based remote sensing systems can provide an accurate real-time monitoring of a volcanic explosive activity.

sedimentation which might cause hardship and damages in areas surrounding volcanoes, including threat to aviation [9]. In particular, real-time monitoring of ash-rich plumes is also crucial for initializing volcanic ash transport and dispersion models [13], [32]. Tephra dispersal from an explosive eruption is a function of multiple factors, including mass eruption rate (MER), degree of magma fragmentation, vent geometry, plume height, particle size distribution, and wind velocity [43].
Mt. Etna, located on the east coast of Sicily (Italy), is one of the most active volcanoes in Europe. The most distinctive phenomena associated with the activity of Etna are represented by the volcanic plumes, sometimes characterized by a significant tephra discharge rate [1], [3], [12]. Volcanic plumes at Etna mostly consist of sustained jets of fluid lava, propelled into the atmosphere from summit craters or lateral vents and driven by expanding gases, which commonly occur at basaltic volcanoes [7]. The fountain gains its momentum by the expansion of gas bubbles that exsolve from the magma as pressure falls while it is rising in the conduit. Height, duration, and erupted volumes of Etna volcanic plumes can greatly vary, with strong lava fountains reaching a height of several hundreds of meters. On November 23, 2013, an intense explosive eruption formed from the NSEC and lasted about an hour. This eruption has been widely analyzed in previous works, focusing on the eruptive processes and tephra volumes [3], integration of observational data [10], tephra fallout characterization [1], plume dynamics [36], and total grain-size distribution retrieval [38]. In this respect, a few instrumentbased estimates are available for the time series of the MER, that is, the amount of material erupted per unit time, a key parameter for evaluating hazard assessment and for ash plume dispersion model initialization [4], [33], [34], [45].
Near-real time MER monitoring and estimation can be provided by several techniques: 1) fixed-pointing Doppler microwave radar [16], [20], [24]; 2) optical imaging in clearair conditions [46], [48]; 3) infrasound sensor network [39]; and 4) electrical probing [6]. All techniques are affected by significant uncertainties, with considerable variations between different estimates of MER [19]. Using a centimeter and millimeter wavelength, ground-based microwave radars represent an important tool for detecting and estimating nearsource tephra MER and concentration being their wavelength comparable or larger than size of lapilli and coarse ash particles as well as less affected by two-way extinction with respect to optical sensing [16], [30], [50]. 0196-2892 © 2019 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission.
See https://www.ieee.org/publications/rights/index.html for more information. Weather radar scanning systems can be exploited to monitor a volcanic plume, measuring the reflectivity due to small lapilli and coarse ash at a fairly high spatial resolution (less than a few hundred meters) and every few minutes [22], [26]- [28]. Weather radars can provide data for estimating the plume tephra volume, total mass and height, using the volcanic ash radar retrieval (VARR) for single-polarization and dualpolarization systems at S-, C-, and X-band [30], [31], [33]. Doppler fixed-pointing radars at L-band has the antenna boresight typically oriented toward the volcano summit craters and are able to follow the plume column dynamics in near-real time providing both tephra power returns and Doppler velocities mainly due to lapilli and bombs [16], [17].
The aim of this article is to analyze the November 23, 2013, Etna volcanic plume in order to: 1) extend the applicability of the VARR methodology to L-band Doppler radar for the quantification of MER; 2) retrieve the incandescent region height and exit velocity from available polarimetric X-band radar data, Doppler L-band and thermal infrared camera; and 3) formulate MER retrieval techniques in a unified way with their own uncertainty and estimate the MER time series and total erupted mass (TEM). To reach this aim, we explore the application of three different strategies for the calculation of erupted mass based on near-surface flux, plume height, and mass continuity, respectively, and we compare the associated results with those obtained from deposit-based techniques.
This article is organized as follows. Section II briefly describes the available instruments as well as the Etna eruption. Section III discusses the proposed methodology to derive the tephra MER and erupted mass. The radars and camera data set is presented in Section IV together with the discussion of the results. Conclusions are presented in Section V.

II. CASE STUDY AND DATA
The Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo (INGV-OE) is equipped with a large set of instruments installed for real-time monitoring of the Etna eruptive activities [41]. Various sensors are included at the Etna site such as several seismic and acoustic sensors, two optical lidars, microwave radars, and thermal infrared and visible cameras installed at different times during the INGV-OE instrument site expansion [41], [48]. The Dipartimento della Protezione Civile (DPC) of Italy operates an X-band mobile weather radar installed at the Catania airport [36], [50].

A. Mt. Etna eruption on November 23, 2013
The November 23, 2013, episode represents one of the most explosive events among all the lava fountains that have occurred at Etna from 2011 to 2013 [1], [41]. The eruption started with a Strombolian activity in the afternoon of the day before. The activity moved from Strombolian to volcanic plume after 9:30 UTC and formed an eruption column that grew suddenly and reached a height of 11 km above sea level (asl) during the most intense phase [10]. The thick volcanic plume and cloud spread tephra particles up to several hundreds of kilometers from the summit craters toward the northeast direction [3]. The volcanic plume portion, characterized by the hot jet made of coarse pyroclasts [1], was clearly distinguishable from the higher eruptive column governed by buoyancy. Abundant fallout of bombs and coarse lapilli occurred on the lower north-east flanks of the volcano, whereas fine lapilli dispersed all along the Ionian coast of Sicily. The fallout of clasts of several tens of centimeters caused severe damage to buildings, solar panels, and cars, and injured a few hikers at about 5 km from the crater [1].
This volcanic plume was clearly observed from the X-band polarimetric weather radar and L-band fixed-pointing Doppler radar as well as from a visible and a thermal-infrared camera (TIC) of the INGV-OE surveillance system. For this case study, the map in Fig. 1 shows the location of the remote sensing instruments used in this article, the main characteristics of which are briefly listed here.
1) The X-band microwave weather radar (MWR) is a dual-polarization scanning radar, included in the Italian weather radar network [50]. The dual use of X-band MWR for the detection of both meteorological and volcanic clouds is possible thanks to a combination of several factors: wavelengths of about 3.1 cm (frequency of 9.6 GHz), transmitted peak power of 50 kW, half-power beamwidth of 1.3 • , and permittivity factor of ash particles (equal to 0.39 with respect to 0.93 of water particles) [29], [35]. The X-band MWR performs a 3-D scan of the surrounding scene as a function of range, azimuth, and elevation with five antenna rounds per minutes. The X-band MWR acquisitions consist of data volumes having an area of about 160 × 160 km 2 wide and 20 km height. The data volume cross sections are sampled along 12 elevation angles plus a vertical one, as shown in Fig. 1, and released every 10 min at a distance of about 32 km from the NSEC. The MWR volume, scanned near the NSEC, has a range resolution of 200 m and about 700 × 700 m 2 transverse spatial resolution. 2) The ground-based L-band Doppler radar (VOLDORAD-2B or VDR hereinafter), operating at a wavelength of 23.5 cm (frequency of 1.274 GHz) was designed by the Observatoire de Physique du Globe de Clermont-Ferrand (OPGC) for the monitoring of explosive activity [16], [17]. It can be deployed near an eruptive vent to measure in real-time the eruptive velocities and backscattered power at high rate up to 20 s −1 . The VDR signal wavelength allows to sound dense lava jets and ash-laden plumes as well as to avoid attenuation by hydrometeors because of cloudy, foggy, rainy, or snowy conditions. Owing to its modularity and limited weight (about 70 kg) the ensemble is easily transportable and, thus, can be used for short-term scientific campaigns, as well for long-term monitoring. At Etna, VOLDORAD-2B is jointly operated by INGV-OE and OPGC, sounding (fixed pointing) 13 volumes (about 28 × 10 8 m 3 ) right above the summit craters every 0.23 s. Its capacity to provide first-order MER in real-time, as well as TEM, onset and end of volcanic plumes, and eruptive crater has recently been shown [17]. Processed velocity and echo power data are stored in an open access database made available from the OPGC website and EPOS platform (see the Acknowledgement later on).
3) The TIC is located at ∼15 km southward from the craters and belongs to the INGV-OE network of video-surveillance system. TIC provides a time series of 640 × 480 pixel images with a spatial resolution of few meters considering the distance between the TIC and the NSEC crater [3], [7]. The height, width, and area of the volcanic plume can be detected by properly selecting the saturated portion of the measured brightness imagery and adopting the procedure described in [3] and [15].

B. Sensor Data Processing
MWR, VDR, and TIC can be processed to derive useful near-source variables. VDR can provide an estimate of the exit velocity v ex , whereas the incandescent region height can be retrieved from MWR and TIC. From both VDR and MWR an estimate of tephra concentration, mean diameter, and MER can also be derived using the VARR algorithm.

1) MWR Data Processing:
The polarimetric MWR is capable of measuring not only the X-band copolar horizontally polarized reflectivity factor Z hh (hereinafter called reflectiv-ity), but also other polarimetric moments such as the differential reflectivity Z dr , differential phase shift K dp , and the copolar correlation coefficient r hv (elsewhere also indicated by ρ hv or ρ co ) [29], [50]. Fig. 2 shows the vertical profiles of X-band Z hh , Z dr , K dp , and r hv along the line connecting the radar antenna with horizontal maximum expansion of the plume (see [36] and [50] for similar plots). It is interesting to note the contrasting trend of X-band Z hh and r hv in the areas immediately above the crater extending vertically for a few kilometers. The X-band Z hh reaches a maximum altitude of 11-km asl, decreasing horizontally more severely after about 20 km from the summit craters, probably due to a faster fallout of large particles, a region identified by values of Z hh ≥ 50 dBZ. The area with low X-band reflectivities (Z hh < 30 dBZ) is associated with outermost edges of the plume suggesting the presence of coarse particles prone to fallout in agreement with the tephra sampling [1]. Since r hv measures the consistency of copolar signal power and phase for each received pulse pair, r hv confirms the Z hh trend, revealing a fairly apparent vertical separation between the incandescent saturated region and the convection eruptive region just above [50]. The Z dr signatures are relatively low, oscillating at around 0.02-0.005 dB, meaning that tephra particles are detected as spherical on average (the material injected from the crater is still very fluid and is characterized by some degree of anisotropy, as noted in [32]). The K dp signature shows an increase in a region, which is slightly displaced with respect to the column above the crater. Positive values of K dp typically indicate slightly horizontal orientation for oblate volcanic particles. The behavior of K dp increment could be due to the presence of falling large lapilli and bombs with a ballistic trajectory. Fig. 3 shows the maximum values of both X-band MWR Z hh and r hv along the column closest to the NSEC at about 32 km using the fourth elevation angle (see Figs. 1 and 2). This plot can be interpreted by looking at the vertical profile of Fig. 2. In particular, low values of r hv suggest non-spherical shapes and tumbling of volcanic particles [29] so that in Fig. 3 the region with r hv < 0.95 can be divided into two regions where different physical processes probably occur [36]: 1) the region above the NSEC, where the ascending gas and particles form the eruption column that is progressively inclining and thickening as it propagates downwind and 2) the region, aside the NSEC at horizontal distances less than about 5 km from the crater, which can be reasonably associated with the fallout of irregular large lapilli and bombs. The region, having r hv > 0.95, extends over the entire remaining plume, detected by the X-band MWR, thus including both lateral cloud advection by wind and fallout of tephra particles [36].
From Fig. 3 it emerges that a combined thresholding on both X-band MWR Z hh and r hv can be used to detect the incandescent region height H IR . In this case study, by doing an iterative analysis aimed at finding a relatively stable estimate, we have empirically set the combined condition Z hh ≥ 50 dBZ and r hv < 0.95 to retrieve H IR . From C-band polarimetric observations of the 2012 Mount Tongariro, for a tephra plume associated to a more viscous magma than Etna, a transition at r hv around 0.9 has been found corresponding to a stronger  (a), but for the copolar cross-correlation coefficient r hv (adim). (c) Same as (a), but for the differential reflectivity Z dr (dB). (d) Same as in (a), but for the differential phase shift K dp (adim). decrease of reflectivity upward of about 1 km above the vent. This limit reflects the upward transition to the buoyancydominated convective column that rose at about 5 km high.
Indeed, we should have a larger set of explosive eruptions, observed by a polarimetric radar, but so far the paroxysm on November 23, 2013, remains one of the best case studies at Etna where the H IR signature is quite clean in both r hv and Z hh signal. This approach may be probably improved by including the other polarimetric features in a tree-logic approach, but from the case study of November 23, 2013, it seems that the improvement is relatively negligible. The relatively good agreement of the proposed radar-based H IR retrieval algorithm with the estimates from thermal infrared camera supports the current approach. Note that, due to MWR finite antenna beam width and the distance of 32 km from the summit craters, the spatial cross-resolution along the eruption column is about 700 m (see Fig. 1).
2) VDR Data Processing: The fixed-pointing L-band VDR is measuring both the radial velocities v r and the received backscattered power P R X derived from L-band VDR [16]. From the observation geometry we can convert v r into exit velocity v ex normal to the surface of summit craters (i.e., v ex = 3.89 v r ) [17], whereas from the specifications of the L-band VDR and the radar constant, the backscattered power P R X can be transformed into the L-band horizontally polarized reflectivity factor Z hh . Fig. 3 shows the VDR Z hh corresponding to its fourth range bin, smoothed about every 17 s and down-sampled every 10 min in order to reconcile the time sampling with that of MWR. The time trend of Z hh is related to the VDR range gate closest to the NSEC. The VDR radar reflectivity factor is higher by 15 dBZ with respect to MWR during climax at 10:00 UTC. This may be due to sampling location mismatches, VDR being measuring right above the crater, effects of Mie scattering regime and impacts of particle non-sphericity, affecting the incandescent region above the crater dominated by bombs and lapilli in the proximal fallout region.
3) TIC Data Processing: The TIC measurements can be processed to extract the incandescent region height H IR from the recorded thermal-infrared brightness temperature imagery over the eruptive time interval [3]. Most techniques are based on imposing a proper threshold to the vertical spatial gradients and/or to edge-contour detection filters [15]. By selecting the TIC frames at time intervals of 1 min, it is possible to derive the incandescent region height H IR in each image.

4) Tephra Concentration:
Starting from Z hh derived from both X-band and L-band radars, we can apply the VARR methodology, considering the ad hoc physical-electromagnetic  Each curve is sampled every 10 min, as imposed by the MWR operational schedule. Data near to the NSEC are relative to the fourth elevation angle of MWR and to the fourth range-bin of VDR. (a) VARR-based retrievals of maximum and minimum of tephra mass concentration C t (g/m 3 ) and mean sphere-equivalent diameter D n (mm), derived from the X-band radar (MWR). (b) Same as in (a), but derived from the L-band radar (VDR). model of non-spherical ash particles in order to derive the mean sphere-equivalent diameter D n and the tephra mass concentration C t (also denoted as C a in [26]: here we prefer the notation C t dealing with near-source pyroclats). The latter is defined as [26] (see also Table I for the list of main symbols) where D is the sphere-equivalent diameter (mm), ρ t is the volcanic particle specific density (kg/m 3 ), N t is the particle size number distribution (PSD, in m −3 mm −1 ), typically characterized by three parameters (i.e., mean diameter D n , shape parameter μ, and number concentration N n ) [26]. The volcanic particle size distribution is parametrized using field and combined data [38]. Equation (1) holds if ρ t is constant and introduces the airborne-particle volumetric fraction f N providing, as a function of PSD parameters, the fraction of tephra particles per unit volume or, more generally, the degree of rarefaction of the ejected material. From (1), percentage values of f N for tephra are usually less than 0.01% [29]. The extension of VARR to L-band is quite straightforward as the backscattering model is valid for both Rayleigh and Mie regimes and the considered particle sizes range from 64 μm up to 131.072 mm. 5) VARR Processing: Fig. 4 shows the VARR-based maximum and minimum retrievals of C t and D n , obtained from X-band MWR and L-band VDR data in the range gate nearest to the NSEC (see Fig. 3). The mass concentration retrievals can reach values of 18 g/m 3 for L-band VDR and 7 g/m 3 for X-band MWR, whereas mean-diameter estimates show sizes from about 5 to 12 mm for L-band VDR and from 0.1 to 4 mm for X-band MWR. Interestingly, the modal diameter of proximal lapilli sampled immediately near the cone for the July 2011 Etna paroxysm was between 11 and 16 mm [7]. Consistent with the radar wavelength, L-band VDR is mainly sensitive to lapilli and bombs, whereas X-band MWR response is also influenced by smaller coarse particles [29]. It is worth noting that MWR peaks (around 10:00 UTC) are slightly anticipated with respect to VDR ones (around 10:10 UTC). Assuming a particle density ρ t of 2700 kg/m 3 , airborneparticle volumetric fraction f N is typically less than 10 −7 .

III. ESTIMATING MASS ERUPTION RATE
The main goal of this article is to provide an estimate of MER using MWR, VDR, and TIC showing how their retrievals can be processed using a unified approach. Indeed, the capability of ground-based radars to estimate the timedependent MER, here also indicated by the symbol Q M (t), is still an open issue [31], [33], [50].
From a methodological point of view, the time-dependent MER Q M (t) can be related to the plume top height H TP and to incandescent region height H IR . Note that the incandescent region height may be similar or higher than the gas-thrust region height (depending on the eruption style where the fragmented particle momentum is coupled or not with tephra plume jet [43]). In this respect, we can apply three approaches in order to estimate Q M (t): 1) the surface flux approach (SFA), physically relating Q M (t) to the eruptive exit velocity v ex , tephra density and the geometry of the crater; 2) top plume approach (TPA), using semiempirical parametric models considering the top plume altitude and environmental parameters; the mass continuity approach (MCA), using the mass continuity equation within the erupted plume above the crater.
A. Methods 1) SFA: The tephra MER Q M (t) (kg/s) (sometimes also indicated byṀ(t)) through the crater at each instant time t can be written as where the crater has a surface S v (m 2 ) in (x,y) coordinates, ρ x (kg/m 3 ) is the density of the eruptive mixture, and v ex (m/s) is the vertical exit velocity normal to the crater surface. On the right-hand side of (2) we assume that both ρ x and v ex are constant averaged values within S v , a reasonable assumption if the crater is relatively small. The objective of the SFA approach, in order to retrieve the time series of Q M (t), is to provide an estimate of ρ x , v ex , and S v from remote sensing instruments MWR, VDR, and TIC at each instant time t, as will be discussed in the following text.
In order to estimate ρ x in (2), we can consider that Etna volcanic plumes are typically characterized by gas fractional content f g between 2% and 3% [8], a gas density ρ g between 0.10 and 0.20 kg/m 3 , and a magma density ρ m between 2500 and 3000 kg/m 3 [46]. These values result in a fragmented magma-gas mixture density ρ x at the vent given by where f m is the volumetric fraction of magma, holding on the right-hand side of (3) a linear mixing with f m =1− f g [48].
In order to estimate S v in (2), we can adopt the usual assumption of a cylindrical conduit with a circular crater with an area S v = πr 2 v , r v being the radius [39]. Experimental evidences suggest that the radius can be set to 13.5 m with an uncertainty of ∼10% [7], [8], a value which reinforces the approximation of constant ρ x . This radius estimate is confirmed by the inspection of the available thermal-infrared imagery evaluating the average size of the detected vertical column. An estimate of S v is then about 572.5 ± 57 m 2 .
In order to estimate v ex in (2), we should distinguish between the three sensors: L-band VDR provides directly v ex as a normal projection of the measured Doppler radial velocity v r , whereas both X-band MWR and TIC data can provide an estimate of the incandescent region height H IR (see Section II-B). Indeed, X-band MWR is a Doppler radar and velocity profiles can also be estimated from ad hoc data processing [36]. However, as pointed out in [36], radar estimate of the updraft velocity cannot be considered as an exit velocity but rather its proxy. This is the reason why in this article we explore the use of the MWR-based estimate of the incandescent region height. The latter can be used to retrieve v ex based on the ballistic equation, also known as Torricelli's equation [40], [45]. This equation, also deducible from energy conservation, can provide an estimate of H IR associated to the vertically directed outflow velocity v ex of a pyroclastic constant flow from the volcano crater, and vice versa, at each where g (m/s 2 ) is the Earth's gravity acceleration and the atmospheric density variation and drag effects are considered to be negligible within the incandescent region. This expression is a valid approximation within the incandescent region, when most pyroclasts are sufficiently large to be considered uniformly accelerated projectiles not entering the upper convection region of the plume [45]. This is typically the case for ballistic bombs. Note that, due to the nonlinear relation present in (4), for a H IR of 2000 m, an uncertainty of 20% (400 m) reflects into an uncertainty of 10% on v ex (19.7 m/s). The use of (4) together with (3) and (2) allows the SFA-based estimate of tephra MER from both X-band MWR and TIC data. The overall fractional uncertainty ε Q of the approximate MER Q M (t) in (2) can be estimated using the first-order error propagation theory for independent (maximum) errors, obtaining the following expression: where δρ x , δ rv, and δ H IR are the (maximum) errors on density at the crater, circular radius, and incandescent region height that are causing the overall error δ Q M on the MER. If from previous considerations we assume that δρ x = 0.15ρx, 2) MCA: The mass continuity equation can be applied to estimate the MER [31], [33]. The MER can be decomposed into two terms, the first one Q Mdi f (t) related to the timevariation of the tephra mass M t between two adjacent instants within the erupted volume, and the second one Q Madv (t) related to the advection of tephra mass with a vector velocity v through the erupted volume. In formula [33] where M t is the tephra mass within the volume enclosed within the closed plume surface S detected by the weather radar scan, C t is the tephra mass concentration, r is the range vector, n S is the unit vector normal to the surface S, and v is the tephra velocity field. The advection term in (6) can be estimated either by the Doppler field, even though radial velocities should be transformed into vector velocity with arbitrary assumptions, or by a space-time cross-correlation technique providing only a single displacement vector for the detected tephra volume [31], [33]. The contribution to the MER of the advection term Q Madv (t) is generally much less than the other term related to time-derivative of the tephra mass.
To estimate the mass eruption Q M in (6), the input data is the time series of Q Mdi f (t) and Q Madv (t) provided by scanning weather radar such as X-band MWR. As already discussed for SFA uncertainty, the overall fractional uncertainty ε (6) can be estimated by where δ Q M is the MER (maximum) error. If we assume that δ Q Mdi f = 0.20Q Mdi f (due to error in tephra volume and time undersampling) and δ Q Madv = 0.10 Q Madv (due to errors in velocity field estimation), the relative percentage error ε (MC A) Q = 22.3%. 3) TPA: Numerical 1-D models, theoretical simplified models and field-based empirical relationships relate the top-plume height H TP to the instantaneous tephra MER Q M (t) through a generalized power law [43], [45], [51], [52], here expressed in a compact form as where the coefficients a 0 , a 1 , b, and c are properly set either by experimental fitting or model analyses.
In particular, in thi article we have considered the following estimators of MER.
1) The empirical relationship, proposed by [32] and here named TPA-MA09, where in (8) a 1 = 0 and c = 0, whereas a 0 = 3.29 kg/s/m b and b = 4.15. Note that TPA-MA09 is indeed proposed in terms of volumetric eruption rate, expressed in m 3 /s, and here has been converted into MER by assuming a magma density ρ m = 2500 kg/m 3 as prescribed in [32]. 2) The analytical relationship, proposed by [14] and here named TPA-DB12, where b = 4, c = 3 and the coefficients a 0 (kg/s/mb) and a 1 (kg/s/m c ) are dependent on the air/plume density and temperature, specific heat capacity, mean buoyancy frequency, radial/wind entrainment, and mean wind velocity across the plume height [14]. In this case study, the considered value of the mean cross-wind along the plume maximum vertical extension is about 20 m/s, using the local model weather forecasts [41]. To estimate the MER Q M in (8), the input data are the time series of top-plume height H TP (t) which can be provided by visible cameras, satellite data, and a scanning weather radar [10]. As already discussed for SFA uncertainty, the overall percentage uncertainty ε Q−TPA of Q M (t) in (8) can be estimated by where δ Q M is the MER (maximum) error. If for MA09 we assume that δa 0 = 0.20 a 0 and δ H T P = 0.20H TP , the relative percentage error ε Q (TPA) = (δ H TP /H TP ) = 83%.

B. TEM
In summary, the instantaneous MER Q M (t) is obtained. 1) For SFA, from (2) after deriving v ex (t) directly from L-band VDR and through H IR in (4) where Q V (m 3 /s) is the DRE eruption rate, given by the MER divided by the magma density ρ m , and the event duration. Note that M I (t) expresses the mass erupted till a given time t starting from the eruption onset, whereas M T is the total mass erupted during the whole event. The latter M T is converted into V T through magma density ρ m . The total erupted volume is introduced to allow a straightforward comparison with results available in the literature on the same case study.
As previously discussed, all the methods SFA, MCA, and TPA are affected by uncertainties and so is the TEM M T . If we assume that the Q M error fraction ε f for each method is time invariant, it can be easily shown that the erupted mass fractional error is given by where t i are the sampling time instant and t is the time step. Note that the time series of the MER Q M (t) and volumetric eruption rate Q V (t) is sometimes summarized by the timeaveraged MERQ M and volumetric eruption rateQ V [1], [3]. In our notation they are defined as where t is the event duration.Q M andQ V are, indeed, the time-average of M T and V T , respectively, that is the TEM and volume divided by the event duration t.

IV. RESULTS
With the SFA, MCA, and TPA methodologies defined in the previous section, we can apply them to estimate the MER and its derived parameters in (10) and (12) for the case study of the November 23, 2013, Etna volcanic plume.
SFA-based techniques are dependent on the estimate of the tephra exit velocity at the crater, derived from each available sensor. Fig. 5 shows estimates of the exit velocity v ex directly derived from L-band VDR, as well as those derived from X-band MWR and TIC data using H IR and (4). Estimated exit velocity v ex shows a behavior similar to H IR with a maximum at 10:00 UTC with values of around 240 m/s from L-band VDR, 230 m/s from TIC data, and 235 m/s from X-band MWR. Fig. 5 also shows estimates of H IR directly derived from of X-band MWR and TIC data as well as those derived from L-band VDR using v ex and (4). The maximum H IR is reached at 10:00 UTC, which for the L-band VDR is around 2600 m above the crater level (acl), for the TIC data around 2500-m acl and for the X-band MWR around 2550-m acl. These H IR estimates are consistent in terms of both values and trends.
On the left panel (a) of Fig. 6 the time series of the retrieved MER Q M (t) are shown, sampled every 10 min (the lowest temporal sampling due to the X-band radar), using SFA, TPA, and MCA methods from TIC data, X-band MWR and L-band VDR. The TPA method, based on the DB12 parametric model, uses an average wind velocity of 20 m/s averaging the vertical wind profile inside the eruption column.
MER estimates from all sensors are in fairly good agreement in the paroxysmal time step from 09:50 till 10:20 UTC with maximum values between 1.0 × 10 6 and 2.3 × 10 6 kg/s. At the beginning of the eruption L-band VDR and TIC data tend to provide MER estimates higher than X-band MWR ones, the latter probably being affected by the uncertainties in H IR discrimination due to a relatively poor cross-correlation coefficient ρ hv signal as well as to the low transverse resolution of the radar beam. At the climax of the eruptive episode all retrievals methods are in fairly good agreement, the TIC-based MER being the largest. The impact of wind velocity on DB12 model tend to provide an MER which is slightly higher (about 25% more) than the SFA-based ones. On the other hand, the MCA-based approach is very close to the  DB12 estimates during the paroxysm climax, probably due to a similar approach based on the erupted plume features such as the plume top height (for TPA) or airborne plume mass (for MCA). As expected for bent-over plumes, the MA09 strategy tends to underestimate the MER with respect to the other methods due to strong wind advection; such an effect is taken into account into DB12 TPA-based method.
MER estimates can be used to evaluate the accumulated ejected mass during the eruption temporal evaluation. The right panel (b) of Fig. 6 also shows the time-integrated erupted mass M I (t) in (10) using the same methods as shown in Fig. 6(a), i.e., SFA and MCA from X-band MWR, L-band VDR, and TIC data as well as TPA from DB12 and MA09. In this plot, the TEM M T is represented by value at the last time step at 10:40 UTC, as deducible from (9). Where MER starts to decrease around 10:10 UTC, M I (t) values tend to saturate. Indeed, from (10), the last value at 10:40 UTC of M I (t) provide the erupted mass retrieved values between 3.6 × 10 9 kg/s from VDR up to 4.7 ×10 9 kg/s from TIC data, except for MA09 (1.7 × 10 9 kg/s) affected by more uncertainty (estimated to be a factor of about 50 at 95% [20]). Note that at the beginning of the volcanic eruption around 9:30 UTC, VDR seems to have already detected some extra mass rate, which could be added to values reported here [20].
The uncertainty of each MER estimation technique, introduced in (5), (7), and (9), can suggest the confidence interval of the obtained results. Fig. 7 shows the trend of the MER estimated value and the respective uncertainty for all methods (SFA, TPA, and MCA) at each sampling time. Between 09:50 and 10:10 UTC, that is, the interval of the largest eruptive activity, the uncertainties (error bars) are, in general, bigger. Estimates using the SFA method show a more significant departure, even if their uncertainty is considered, within the ending tail of the eruption, a feature probably related to the uncertainty in the use of the Torricelli equation (to estimate of incandescent region height or exit velocity, as shown in Fig. 5) and the summit crater which may even change within the eruption itself. The two DB12-based and MA06-based TPA methods show an expected discrepancy which is not accounted for their relative uncertainties. MCA-based values have uncertainties comparable to DB12-based ones. Finally, the average and standard deviation of all-method retrieved MERs show an overall fairly consistent increasing and decreasing trend with a paroxysm MER standard variability between 1.5 × 10 6 and 2.4 × 10 6 kg/s. In most cases the only way to validate MER estimates is to compare the TEM, derived from MER, with available ground deposits. The estimates of the TEM M T (kg) is shown in Table II for the November 23, 2013, Etna volcanic plume. The results refer to erupted mass derived from the literature data (see rows a, b, c1, and c2) together with TPA-MA09 and TPA-DB12 estimates (see rows d and e), SFA retrievals from TIC data, X-band MWR, and L-band VDR (see rows f , g, and h) and X-band MWR MCA (see row i ), respectively. Uncertainties, estimated by (5), (7), and (9), are also indicated. All retrievals show the same order of magnitude around few kilotons in agreement with the erupted mass derived from satellite data (3 × 10 9 kg in b and 5.7 × 10 9 kg in c 1 ) as well as from wind-driven TPA-DB12 (4.3 × 10 9 kg in d). The SFA estimates are interestingly very similar among them with estimates 3.8 × 10 9 kg (from L-band VDR) and 4.7 × 10 9 kg (from TIC data). These MER retrievals are within the uncertainty of the November 23, 2013, Etna eruption field values (see rows a and c 1 ) providing erupted mass values between 1.3 ± 1.1 × 10 9 kg and 5.7 × 10 9 kg, obtained by estimating the fallout deposit using the Weibull distribution [1], [38]. Moreover, they are higher than TPA estimates    [1], [3], AND [20] from MA09 with a erupted mass of 1.7 × 10 9 kg, but in fairly good agreement with DB12 and MCA ones showing values around 4.3 (with a mean wind velocity of 20 m/s) and 3.9 × 10 9 kg, respectively. Independent estimates of time-average mass and volumetric rates can also be used for comparison [1], [3]. Using magma density ρ m = 2700 kg/m 3 [8] and t = 4200 s [70 min, from Fig. 6(a)], Table III shows erupted mass and volumetric eruption rates, computed as defined in (10) and (12), from SFA methods using TIC data, X-band MWR, and L-band VDR, from the X-band MWR MCA method and derived from the literature. As expected from the discussion on erupted mass estimates, the three SFA retrievals are in good agreement with time-average MERs between 1.6 × 10 6 and 0.9 × 10 6 kg/s, time-average volumetric discharge rates between 627 and 335 m 3 /s and DRE volumes between 1.9 × 10 6 and 1.4 × 10 6 m 3 . The lowest values are shown in MCA approach but similar to that reported in [1]. By assuming ρ m = 2700 kg/m 3 (accepting an uncertainty of 10%) and t = 3000 s (50 min, disregarding the first and last 10 min), SFA-based timeaverage volume discharge rates are between 627 and 335 m 3 /s, the latter value close to 360 m 3 /s as reported in [3] using an estimate based on a TIC for the same event. The MCA-based time-averaged volumetric eruption rate is characterized by a value between 348 and 526 m 3 /s. The total erupted volumes, derived from SFA methods, are in fairly good agreement with those provided in [3] and in [20] of 1.6-1.7 × 10 6 m 3 .

V. CONCLUSION
Three different approaches have been presented and compared to determine MER from microwave radars at L-band (23.5-cm wavelength) and X-band (3.1-cm wavelength), namely, the SFA, the MCA, and the TPA. These approaches exploit the radar Doppler or polarimetric capabilities as well as fixed-pointing or scanning mode and both radar data have been processed by means of the model-based VARR methodology. We have also discussed the overall formulation and some assumptions behind both SFA and TPA methods, showing how these uncertainties can reflect into the estimate of the TEM as well as time-average discharge rates. As a reference we have taken into account the estimate of the MER from a video TIC, exploiting its capability to detect the incadescent region height in this event and applying the Torricelli equation to estimate the exit velocity. The latter, when estimated from X-band MWR, L-band VDR, and TIC data, are between 160 and 230 m/s in the paroxysmal event within a difference among the various sensor retrievals less than 25%.
The intercomparison between the SFA and TPA methods, in terms of both MER and erupted mass, shows fairly good agreement. Estimated erupted mass is between 3.7 × 10 9 and 3.8 × 10 9 kg for SFA applied to L-and X-band radar data, respectively, and between 1.7 × 10 9 and 4.7 × 10 9 kg based for TPA, slightly less than the camera-based estimates equal to 4.7 × 10 9 kg. MCA-based erupted mass estimates are comparable to SFA ones. These SFA, MCA, and TPA results for erupted mass are in good agreement with the tephra fall deposit mass estimates between 1.3 ± 1.1 × 10 9 and 5.7 × 10 9 kg. Moreover, SFA-based time-averaged MERs and DRE volume eruption rates from the three remote sensors are in agreement with other independent estimates, available in the literature.
The analysis of this case study indicates that ground-based radars can be exploited to provide a self-consistent monitoring of the time-varying activity of explosive volcanic eruptions. Polarimetric weather radars can offer the capability of 3-D scanning instruments thus providing a monitoring of the plume dynamics. By combining radar at different wavelengths (23.5 cm at L-band, 3.1 cm at X-band, and 0.9 cm at Ka-band) together with lidar monitoring at visible nearinfrared wavelengths (0.5 and 1.1 μm) to gain a sensitivity to finer particles [34], the total grain size of the tephra plume could be retrieved. Further work is needed to assess the SFA methods, using more explosive eruption cases with a set of instruments at least comparable to the one used in this event and deposit-based volumes. A more robust selfconsistent approach to the near-real-time estimates of MER should be able to remove some arbitrary assumptions in the SFA formulation (e.g., crater geometry) by exploiting different methodologies and multiple sensor data. She is a Volcanologist with the Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo, Catania. Her research focuses on the analysis of the dispersal and fallout processes of eruptive plumes during explosive eruptions; calibration, sensitivity analysis and uncertainty estimation of ash dispersal models; laboratory and field experiments; development of a multidisciplinary system for the detection and monitoring of volcanic plumes and analysis of explosive activity using different remote sensing techniques (e.g. radar, lidar, satellites).
Dr. Scollo coordinated several projects and one of them, the VAMOS SEGURO project, was selected in 2012 as a "best practice" among several European Cooperation Projects. In 2011, she obtained the Rittmann Medal for young researchers in volcanology and in 2010 the paper Scollo et al. He held an NSF Post-Doctoral Fellowship at Penn State University, State College, PA, USA. He then spent one year at the Observatoire de Physique du Globe de Clermont-Ferrand (OPGC), Université Clermont Auvergne (UCA), Clermont-Ferrand, to specialize in sounding of volcanic explosive eruptions using dedicated transportable radars. Since 2002, he has been a Physician Adjoint with the Laboratory Magmas et Volcans, OPGC, and was awarded the HDR diploma in 2017. As part of the French SNOV (CNRS-INSU), he is in charge of a unique service (VOLDORAD) dedicated to volcanological Doppler radars comprising four instruments (three UHF, one mm-wave) involved in ten campaigns, one radar being permanently monitoring Etna in collaboration with INGV-OE. He has devoted most of his research to better quantifying volcanic jets and plumes dynamics and source parameters from radar remote sensing. He was a Coordinator of the European MED-SUV program from 2013 to 2016 and led the ash plume radar project of the ClerVolc Laboratory of Excellence at UCA from 2015 to 2018. He also teaches courses in geophysics.
Dr. Donnadieu is a member of the French CNAP and SNOV national committees the IAVCEI, the scientific board on Soufrière de Guadeloupe, and the OPGC administration board. He received three PEDR and one UCA awards for scientific excellence.