Determination of time-and height-resolved volcanic ash emissions for quantitative ash dispersion modeling : the 2010 Eyjafjallajökull eruption

Introduction Conclusions References


Introduction
Volcanic gas and aerosol emissions influence climate (Robock, 2000), pose hazards to aviation (Casadevall, 1994;Prata and Tupper, 2009) and health (Horwell and Baxter, 2006), and iron supplied by ash fallout may enhance ocean productivity and lead Introduction

Conclusions References
Tables Figures

Back Close
Full to a drawdown of atmospheric carbon dioxide (Duggen et al., 2010;Langmann et al., 2010).These and other impacts (Durant et al., 2010) depend critically on the total mass of eruption products and the altitude at which they are effectively released into the atmosphere, neither of which is well known.Although models can calculate the long-range ash dispersion with considerable accuracy (Witham et al., 2007), robust estimates of eruption source parameters are required to obtain reliable results.Volcanoes exhibit a broad range of eruptive styles and variability (Woods, 1995), thus making theoretical attempts at predicting source parameters challenging.To date, only indirect or subjective methods are available to estimate ash emission rates (Mastin et al., 2009).For example, if eruption column heights are known, e.g., from weather radar measurements (Lacasse et al., 2004), empirical relationships may be used to estimate the mass flux of tephra (Sparks et al., 1997;Mastin et al., 2009).However, these relationships are subjective (Tupper et al., 2009), the vertical profile of ash emissions remains poorly quantified, and only a small but highly variable fraction of the tephra is fine grained enough (< 30 µm) for long-range atmospheric transport.
In this paper, we present an objective method to determine the volcanic ash emission rate, which is based on inverse modeling constrained by satellite measurement data.The eruption of the Eyjafjallajökull volcano (19.61 • W, 63.63 • N, 1666 m a.s.l.) in the year 2010 represents an ideal test case for our method, as this eruption released large amounts of ash while a wealth of measurement data have been collected during the event, which are available for model evaluation.After 18 years of intermittent seismic unrest (Dahm and Brandsdóttir, 1997), an effusive eruption of basalt on the eastern flank of Eyjafjallajökull occurred from 20 March to 12 April 2010, followed by an explosive eruption under the Eyjafjallajökull glacier on 14 April 2010.The interaction of magma and ice augmented explosive activity and generated large proportions of fine ash that were emitted into the atmosphere (Gudmundsson et al., 2010).Intense tephra discharge continued for several days and the prevailing meteorological conditions resulted in ash transport directly towards Europe, where air traffic was grounded for several days.The eruption strength increased again in May, leading to further air Introduction

Conclusions References
Tables Figures

Back Close
Full space closures.The Eyjafjallajökull eruption demonstrated how susceptible aviation is to volcanic eruptions in Iceland, as already suggested beforehand (Sveinbjörnsson, 2001).

Methods
The inverse method to determine volcanic ash emissions merges a priori information on the ash emissions, satellite observation data and simulations with a dispersion model to derive improved a posteriori ash emissions.In this section, we describe the various data sets and methods used for the inverse modeling.

Satellite data
Observational constraints on volcanic ash emissions were provided through infrared satellite retrievals of total column airborne ash loadings (Prata, 1989;Prata and Grant, 2001;Clarisse et al., 2010a) using data from two different satellite instruments, the geosynchronous Meteosat Second Generation (MSG) Spin-stabilised Enhanced Visible and Infrared Imager (SEVIRI) and the polar-orbiting MetOp Infrared Atmospheric Sounding Interferometer (IASI).These two instruments combine high temporal coverage (SEVIRI has a sampling time of 15 min but hourly averages were used here), with an enhanced sensitivity to ash (IASI has over 1000 spectral channels which can be used for ash detection).

SEVIRI
SEVIRI observes the earth disk over a total field of view of 70 • in 12 channels from the visible to the infrared.The SEVIRI temporal sampling time is 15 min with a spatial resolution of ca. 10 km 2 at the below-satellite location which increases to ca. 100 km 2 near the edges of the scan.Only infrared channels were used in the analyses for a subregion of the SEVIRI disk covering the geographical region 30 Full to 70 • N. The retrieval relies on measurements made at 10.9 µm and 12.0 µm, that are corrected for the effects of absorption and emission by atmospheric water vapour (Yu et al., 2001) and then inverted to determine optical depths and effective particle radius using established methods (Prata, 1989;Prata and Grant, 2001;Wen and Rose, 1994).A Mie scattering program (Evans, 1988) and a discrete ordinates model (DOM) (Stamnes and Swanson, 1981) were used to estimate the top-of-atmosphere (TOA) brightness temperatures at 10.9 and 12.0 µm, assuming a plane-parallel cloud of andesitic spherical particles with uniform cloud temperature overlying a background of uniform surface temperature.Refractive indices of andesite (Pollack et al., 1973) as a function of wavelength were interpolated, convolved with SEVIRI response functions, and input to the Mie scattering code using a set of assumed modified-γ size distributions with different mean particle radii and standard deviations.The scattering code outputs a set of extinction coefficients, asymmetry parameters and single-scatter albedos for each particle size and wavelength.These values are then used in the DOM code to determine the TOA brightness temperatures.The cloud top and background surface temperatures are also required as input to the DOM code.Initial estimates of these are determined from the data by finding clear and opaque pixels and assigning these to the surface (T s ) and cloud top temperatures (T c ), respectively.As there is some error in this assignment, an ensemble of DOM calculations is performed for values of T c and T s that are ±10 K different for each initial temperature.Finally, each DOM calculation is performed in equal optical depth steps of 0.02 starting at 0 (clear field of view) and ending at 7.98 (opaque field of view).The result of these calculations is a very large table of simulated TOA brightness temperatures for the two SEVIRI channels that can be interpolated based on the water vapour corrected SEVIRI data.The retrieval provides the best fit effective particle radius (in steps of 0.02 µm) and infrared optical depth (in steps of 0.02).The mass loading (g m −2 ) can be calculated from, Introduction

Conclusions References
Tables Figures

Back Close
Full where ρ is the density of the ash, r e is the retrieved effective particle radius, τ λ is the optical depth and Q ext (λ) is the extinction efficiency, all functions of the wavelength, λ, leading to errors of 40-60% in estimated mass loading (Wen and Rose, 1994).The retrievals have a preferential sensitivity to ash with particle diameters from 2 to 32 µm.
The composition of the pixel (e.g.ash, meteorological cloud, clear) is identified prior to performing the ash retrieval.The main test requires that the water vapour corrected brightness temperature difference (BTD, which is a function of the atmospheric pathlength and so depends on the scan position) between the 10.9 and 12.0 µm SEVIRI channels must be < −0.8 K.In this case the pixel is assumed to contain volcanic ash (Prata, 1989).Other conditions applied include an opacity test where if the sensed brightness temperature at 12.0 µm is less than 230 K, the cloud is assumed to be optically thick and retrievals are not made.
For the inverse modeling, the 15 min pixel-by-pixel mass loading retrievals were binned into 0.25 • × 0.25 • grid cells and time-averaged to provide 240 longitude × 120 latitude values every hour.A parallax correction (Vicente et al., 2002) was applied to all ash-affected pixels assuming that the ash clouds were at 6 km height.This simplification results in a small error in geolocation, but is an improvement on using the data without a parallax correction.

IASI
IASI is a sunsynchronous polar orbiting infrared sounder (Clerbaux et al., 2009).With its high spectral resolution and low radiometric noise, it has proven very useful in monitoring a host of trace gases.Relatively little attention has been given to the sounding of aerosols with IASI and there are no current or planned operational products pertaining to aerosols.Infrared sounding of aerosols has, however, a number of distinct advantages, such as the availability of night time data and the high sensitivity to aerosol Introduction

Conclusions References
Tables Figures

Back Close
Full morphology.Due to its complexity, a sophisticated method for retrieving radius and mass loadings (Clarisse et al., 2010a), based on optimal estimation which iteratively fits synthetic spectra to an observed spectrum by varying radius and mass loading is not suitable for the treatment of large amounts of data.Another difficulty is that reliable results can only be obtained if the wavenumber-resolved refractive index of the observed aerosols is known.For the Eyjafjallajökull eruption, none of the few published refractive index data for the erupted ash matched the high resolution observations (for broadband sensors like SEVIRI this is not so much of an issue).
Recently, a sensitive method was presented for the detection of volcanic ash based on correlation coefficients (Clarisse et al., 2010b).We applied this method here for two months of IASI data for the region 60 • W to 50 • E and 45 • N to 90 • N. Mass loadings of sulfur dioxide in the upper troposphere and lower stratosphere were retrieved in parallel.Coincident measurements of sulfur dioxide provided good supporting evidence for the many puffs of ash detected.The ash detection algorithm does not provide mass loadings.As a measure for the total ash column, the BTD between the IASI channels at 1231.5 cm −1 and 1160 cm −1 (Clarisse et al., 2010b) was calculated (which is close to zero for clear or cloudy scenes and positive for ash scenes).Forward calculations based on basaltic ash optical properties show that this difference is approximately linearly proportional to the ash mass loading for all but the highest ash concentrations.The exact conversion factor depends mostly upon plume altitude and size distribution and was calibrated to match the SEVIRI retrievals of mass loadings from coincident measurements.

A priori emissions
We compiled mean eruption column heights from six-hourly (VAAC, 2010) and daily (Jakobsdóttir et al., 2010)  3-hour interval to estimate the mass flux corresponding to observed plume heights.This calculated mass flux was then vertically distributed according to model predicted magma densities and plume radius, yielding a time-height gridded inventory with 328 3-hourly intervals for the period 14 April to 24 May, and 19 layers of 650 m vertical resolution.We assumed that 10% of the erupted mass was fine ash in the size range of 2.8-28 µm to which satellite measurements are sensitive, obtaining a total fine-ash emission of 11.4 Tg for the 41-day period considered.When put into the dispersion model, this leads to simulated vertically integrated atmospheric ash loadings that are in "reasonable" overall agreement with available satellite data.Notice that our a priori emission estimate should be more accurate than common operational methods, which use statistical relationships between plume heights and total ash emission (Sparks et al., 1997;Mastin et al., 2009).These methods do not consider actual atmospheric conditions and do not model the vertical distribution of the ash.
For determining the emission uncertainties, the time-height emissions were first smoothed by applying a tri-cubic weight function.Subsequently, the emission uncertainties were set to 100% of the highest weighted emission within a window of 1300 m and 6 h half-widths, respectively.Furthermore, minimum values were prescribed in order to avoid zero uncertainties for grid cells without ash emissions.Since we cannot objectively determine the true values of the emission uncertainties, they were chosen such as to allow the inversion to substantially change the emissions, while still being guided by the a priori estimate.

Model simulations
To simulate the dispersion of volcanic ash, we used the Lagrangian particle dispersion model FLEXPART (Stohl et al., 1998(Stohl et al., , 2005)).The simulations accounted for gravitational particle settling (Naeslund and Thaning, 1991) and wet and dry deposition (Stohl et al., 2005), but ignored ash aggregation processes.FLEXPART was driven with 3hourly meteorological data from two different sources, namely the European Centre for Medium-Range Weather Forecasts (ECMWF) analyses with 0.18 • × 0.18 • resolution Introduction

Conclusions References
Tables Figures

Back Close
Full and 91 model levels, and the National Centers for Environmental Prediction (NCEP) Global Forecast System (GFS) analyses with 0.5 • ×0.5 • resolution and 26 levels.Comparisons between simulation results using these two alternative data sets were used to quantify model uncertainty (see below).
To improve the a priori emissions by the inversion algorithm, it was necessary to run the dispersion model for each one of the 6232 (328 times × 19 layers) emission grid cells.Each one of these 6232 scenario simulations explored the sensitivity of total atmospheric ash columns in our modeling domain ( 30• W to 30 • E and 40 • N to 70 • N) to the ash emissions in a single emission grid box.The simulations extended over six days and carried 360 000 particles in nine different size bins (4,6,8,10,12,14,16,18, 25 µm diameter; particle density 3000 kg m −3 ).Per simulation, a unit mass of ash weighted with the assumed size distribution (see below) and multiplied by the sensitivity of the satellite retrievals was distributed over all particles.The model output consisted of total atmospheric columns of ash for comparison with satellite observations at a spatial resolution of 0.25 • × 0.25 • .No vertically resolved ash concentrations were determined, to keep the model output at a manageable size.
To obtain vertically resolved ash concentrations, we also ran single model simulations based on the gridded a priori and a posteriori emission fluxes.These simulations produced output at high vertical (250 m) resolution and for a global domain and extended over the full period of the eruption.For these simulations, the ash size distribution was extended beyond the range measured by satellite, using a total of 25 size classes from 0.25-250 µm and releasing 24 million particles in total.Results from these simulations were used for model validation and to quantify the ash concentrations over Europe.

Ash particle size distribution
To define the ash particle size distribution, we forced the dispersion model, using the a priori source term, to reproduce the measured deposit size distribution (see Fig. 1) at a location close to the volcano.This was achieved by adjusting the emitted size Introduction

Conclusions References
Tables Figures

Back Close
Full , where M is mass and D p particle diameter.We obtained a distribution with a primary mode at 10 µm diameter (geometric standard deviation 1.3) and a secondary mode at 180 µm (geometric standard deviation 0.35), and used further distributions with the primary mode shifted towards smaller (7 µm) and larger (13 µm) particle size modes for sensitivity analyses (Fig. 1).Notice that with the chosen initial reference model size distribution (thick blue line in Fig. 1), the size distribution of the deposited ash (black asterisks in Fig. 1) can be matched nearly exactly by the model.The initial size distributions are broadly consistent with measurements of airborne ash in the 18 May 1980 Mount St. Helens eruption (Hobbs et al., 1982).Modeled size distributions with maximum modal diameter > 13 µm produced results in poorer agreement with satellite observations at long distances from the volcano due to increased gravitational settling and correspondingly reduced ash atmospheric residence time.Alternatively, modeled particle size distributions with a maximum modal diameter below 7 µm are inconsistent with downwind in situ measurements of ash particle size distributions in the Eyjafjallajökull ash cloud, which show larger modal diameters (Schumann et al., 2010).

Inversion algorithm
In previous studies (Eckhardt et al., 2008;Kristiansen et al., 2010), we developed an inversion algorithm to calculate the vertical distribution of sulfur dioxide emission rates for instantaneous volcanic eruptions.The algorithm is based on the theoretical work of Seibert (2000) and was used in modified form also for determining the spatial distribution of greenhouse gas emissions (Stohl et al., 2009).A full description of the algorithm was given previously (Eckhardt et al., 2008;Stohl et al., 2009) and is not repeated here.A method to derive a posteriori uncertainties by propagating the uncertainties in the a priori emissions, the observations and the model calculations through the inversion algorithm was developed by Seibert et al. (2011).Sensitivity tests and Introduction

Conclusions References
Tables Figures

Back Close
Full evaluation with measurement data were also reported in these previous studies.In this paper, we use this algorithm for the first time to yield volcanic ash emission rates.Furthermore, it is the first time we derive emissions both as a function of altitude and time for a six-week long eruption.
To establish the sensitivity of atmospheric ash loadings to spatially and temporally resolved emissions, 6232 different emission scenarios were established.For each scenario, a unit amount of ash was emitted in one of 19 vertical layers stacked up to 12.3 km altitude and during one of 328 3-h time intervals between 14 April and 24 May 2010.For each scenario, FLEXPART was run to evaluate the atmospheric ash total column loadings.The model results for all scenarios were matched with about 2.3 million 0.25 • × 0.25 • gridded satellite observations, which were available hourly for the geosynchronous platform and twice daily for the polar orbiter.The matched data set was then fed into the inversion algorithm which optimally merges satellite observations, a priori emissions and simulated sensitivities (Eckhardt et al., 2008;Kristiansen et al., 2010).The result of the inversion is vertically and temporally resolved a posteriori emissions on the original emission grid, obtained as a linear combination of all scenario source terms, which optimize the agreement with both the a priori emissions and the satellite observations when considering the uncertainties of both data sets.An important improvement over our previous work (Eckhardt et al., 2008;Kristiansen et al., 2010;Seibert, 2000) is that the model error is approximated for every individual grid cell as the difference between a posteriori model results based on GFS and ECMWF meteorological input data, respectively.These model results were obtained from simulations over the full study period using, by means of an iterative loop, the a posteriori emissions.While a larger model ensemble would be preferable to characterize the model error, this is an improvement over assuming an arbitrary constant model error.The model results were also used to identify regions where ash older than six days contributed more than 20% of the total ash.Data from such locations were not used for the inversion because the sensitivity calculations extended only over six days and older ash would contaminate the inversion.All other satellite observations Introduction

Conclusions References
Tables Figures

Back Close
Full detecting ash were used for the inversion but 75% of the observed zero concentration values were removed, applying a random data thinning scheme which was weighted towards keeping cases providing a strong emission constraint.In total, 2.3 million observation cases were used for the inversion.While the inversion method formally propagates stochastic errors in the input data to the calculated a posteriori emissions, the overall uncertainty of the emissions is driven mainly not by stochastic errors but by the 40-60% errors of the satellite ash retrievals, which are partly systematic.Therefore, we report all ash emission and concentration errors as 50%.In particular, low-altitude ash clouds with loadings < 0.5 g m −2 are often below the detection limits of the two satellite instruments used and, thus, emissions during episodes with less intense eruptions may be biased low.On the other hand, errors for thick high-altitude ash plumes may be smaller than 50%.

A posteriori emissions
From the inversion, we obtain a total fine ash (diameter 2.8-28 µm) emission of 8.3 ± 4.2 Tg, about 73% of the a priori estimate.Extrapolating the emissions to the size range of 0.25-250 µm using the size distribution shown in Fig. 1 yields a total ash emission of 11.9 ± 5.9 Tg but we must keep in mind that the mass outside the 2.8-28 µm range is not actually constrained by satellite data.On average, the a posteriori emissions are shifted towards higher altitudes (Fig. 2a) and the temporal evolution of the emissions is considerably different from the a priori values (Fig. 2d).For example, high-intensity emissions between 16-18 April are reduced by more than a factor of four, whereas emission intensity on 12-13 and 15 May is substantially increased, due to satellite observations that are inconsistent with the a priori emissions.The results for the GFSand ECMWF-based inversions are very similar, both with respect to the temporal and vertical distribution.Generally, the a posteriori emissions are more variable than the Introduction

Conclusions References
Tables Figures

Back Close
Full a priori emissions, both in time and altitude (Fig. 2b,c).The a posteriori emissions are released mainly in a few strong pulses, typically close to or even above the top of the a priori eruption column.During these pulses, little ash was emitted below 4 km above sea level.
The inversion reduced the root-mean-square (RMS) error between the gridded ash total columns in the ECMWF-data-driven simulation and the satellite data by 28% when comparing a priori and a posteriori results.Visual inspection suggests that the a posteriori ash dispersion patterns are more consistent with those observed by the satellites and a large fraction of the remaining RMS error is due to noise in the satellite data, as shown below.

Emission uncertainty reduction
The a priori uncertainties of the volcanic ash emissions are assumed to be 100% of the highest emission value in the vicinity of a space-time grid cell.These errors are reduced by the inversion due to the incorporation of observation information.The relative error reduction (Fig. 3) is typically largest where ash emissions (Fig. 2) -and, thus absolute values of a priori uncertainties -are highest.The reason for this is that the signal-to-noise ratio is larger in grid cells containing a lot of ash.Error reductions can be close to 100%, meaning that the formal methodological a posteriori errors are very small.However, real errors are likely much larger, due to probable systematic satellite retrieval errors of 40-60% and error correlation, which is not taken into account by the inversion method.Furthermore, the inversion scheme assumes errors to be uncorrelated, thus yielding a too large error reduction.Still, Fig. 3 provides valuable information as it indicates where the satellite data provide a good constraint on the emissions and where such constraints are probably weaker (e.g., during the period 6-11 May).Introduction

Conclusions References
Tables Figures

Back Close
Full

Sensitivity to ash particle size distribution
While dispersion model results are sensitive to the assumed ash particle size distribution, the impact on the estimated source term is remarkably small for the tested size distributions (Fig. 4).One reason for this is that the satellite retrievals are only sensitive to fine volcanic ash in a rather limited range of particle sizes (Fig. 1).Thus, only differences in the shape of the size distribution in this range can affect the results.Absolute differences in sedimentation velocities are relatively small for particle sizes smaller than about 10 µm.Furthermore, the inversion is guided most strongly by observations of large volcanic ash columns in the vicinity of the volcano, which further limits the effects of different sedimentation rates for the chosen size distributions.

Sensitivity to changing the satellite data set
We performed ECMWF-based inversions also for subsets of the satellite data used, namely for either SEVIRI or IASI data alone.The a posteriori total fine ash emissions when using only SEVIRI (IASI) data were 7.9 (10.4) Tg, 7% less (22% more) than the 8.5 Tg obtained for the ECMWF-based inversion when using both data sets.The inversion using IASI data only is closest to the a priori emission of 11.4 Tg.The main reason for this is that the number of gridded IASI observations is about an order of magnitude smaller than the number of SEVIRI observations, thus providing a weaker constraint on the emissions, which therefore remain closer to the a priori values.
The emission changes made by the inversion are, however, relatively consistent for the two data sets separately and for the combined data set.For instance, the high a priori emissions from 16-18 April are considerably reduced in all cases, even though the reduction is smallest when using IASI data only (Figs. 5 and 2).All inversions also lead to substantial emission increases for 12-13 May and to a general shift of the ash emissions to higher altitudes.Introduction

Conclusions References
Tables Figures

Back Close
Full

Model evaluation
Model results were evaluated against a large number of independent observational datasets, which are presented in the following.

Analysis of observed plume top heights
Eruption column top heights were estimated from archived images of two webcams viewing the Eyjafjallajökull summit when the view was not obscured by clouds.Both webcams are made by MOBOTIX and were installed by Mila (http://www.mila.is).Webcam 1 was located at a distance of approximately 10 km and has a field of view of 8.2 × 6.1 km.The maximum visible altitude is about 4.8 km.Webcam 2 was positioned roughly 15 km from summit and has a 17.25 × 12.3 km field of view.The maximum visible altitude is about 7.8 km.When the plume was clearly visible, plume top altitudes were estimated from geometrical principles, taking into account the camera characteristics.Errors in the plume top altitudes can occur when the plume is tilting from or towards the cameras.For a few cases, plume top heights in the vicinity of the volcano were also inferred from stereographic observations made with the Multi-angle Imaging Spectro-Radiometer (MISR) (http://www-misr2.jpl.nasa.gov/EPA-Plumes/projectArea.cfm?ProjectArea=19) and thermal emissivity observations from the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER).Furthermore, a few observations were made from a research aircraft (Schumann et al., 2010).Differences between estimated plume heights from the various data sources are substantial (Fig. 6), which is partly due to strong temporal variability of eruption column heights.There are also biases between the various data sets, for example webcam 2 gives consistently higher plume tops than webcam 1 (notice, however, that webcam 1 cannot observe plume tops above 4.8 km).The a posteriori ash emissions are broadly consistent with the observed plume tops.With few exceptions, high plume tops are observed only for periods for which the inversion resulted in relatively large ash emissions, and the modeled plume tops are within the range of observed values.Introduction

Conclusions References
Tables Figures

Back Close
Full In particular, there is often very close agreement between the modeled and webcam 2 plume tops.Furthermore, in general, the a posteriori plume heights are in better agreement with the webcam observations than the a priori plume heights.One period with very different a priori and a posteriori emissions is 12-13 May, for which eruption intensity was reported to have declined slightly from previous days with plume tops at 4-5 km (Jakobsdóttir et al., 2010).Dense low clouds obscured the eruption plume from the ground on 12 May, making reports of visually observed plume tops inherently uncertain.The inversion adjusted both the emission rate and height to higher values (Fig. 2), which is confirmed by plume top heights from independent satellite data (5.5 km) and a short period of webcam observations unaffected by clouds (6-7 km).Due to the emission changes on 12-13 May, the a posteriori simulation reproduces an ash cloud observed by satellite over Great Britain on 14 and 15 May, which the a priori simulation misses nearly completely (see movies described in the appendix).

Evaluation with lidar data
For evaluation of the vertical structure of the modeled ash cloud, we used data from the CALIOP (Cloud-Aerosol Lidar with Orthogonal Polarization) lidar on the CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations) platform (Winker et al., 2007), which has a vertical resolution of 30-60 m.We analyzed total attenuated backscatter at 532 nm, which is a primary Level 1B profile product (version 3.01).This signal responds to aerosols (including volcanic ash) as well as meteorological clouds which in many cases can be distinguished.
Volcanic ash emitted between 6-10 May wrapped around an anticyclone centered south of Iceland and formed several ash bands extending throughout much of the Eastern North Atlantic (Fig. 7).The satellite observations show detailed horizontal distributions but also some important differences between the SEVIRI and IASI data sets due to retrieval problems, for instance in the vicinity of the volcano and near 37 • N where ash loadings are uncertain, partly because of clouds.CALIPSO passed over the ash Introduction

Conclusions References
Tables Figures

Back Close
Full cloud and scanned along a vertical cross-section through several ash bands, allowing an evaluation of the modeled ash vertical distribution for ash aged from a few hours to several days.Since the measured total attenuated backscatter cannot be easily converted to ash mass concentrations, the comparison with the model is qualitative and further complicated by meteorological clouds, which produce similar backscatter signals as ash clouds.Nevertheless, the fresh plume (age ≈ 4 h) observed emanating from the volcano (near 60 • N) is clearly identifiable as volcanic ash and modeled at the correct position, height and vertical extent (2-4 km).The cross section through aged (≈ 4 days old) ash from 49-57 • N shows modeled ash and measured backscatter maxima at various heights between 3 and 9 km.Not all features are in perfect agreement but the modeled plume maximum near 51 • N is very close to the observed backscatter maximum.The ash band from 46-50 • N (age ≈ 40 h) is well simulated, with the ash cloud extending from 3-6 km and tilting upwards to the south.Finally, the ash band near 38 • N below 5 km (age ≈ 3 days), with faint features tilting northward up to 9 km, is also reproduced by the model.
In the appendix, we show a comparison between the a priori and a posteriori model results for the overpass discussed above (Fig. A3) and four more comparisons with other CALIPSO overpasses.While meteorological clouds often complicate the comparisons, qualitatively we find that the a posteriori results are in better agreement with the CALIPSO data (Figs.A1-A5).We also evaluated our results against quantitative vertical ash concentration profiles obtained from lidar measurements over Europe (Figs.A6-A8) and find that the modeled a posteriori ash concentrations are similar to the observed concentrations.

Evaluation with aircraft data
The volcanic ash cloud was also probed by the DLR Falcon research aircraft carrying instrumentation for aerosol and other measurements.We use the ash mass concentrations derived from the measured particle size distribution obtained for the most realistic assumption of the ash particles' optical properties, as described in detail by 5558 Introduction

Conclusions References
Tables Figures

Back Close
Full April 2010 (Fig. 8).At that time, the ash concentrations over Europe were already much lower than from 15-17 April 2010 (see movies described in the appendix) but most of the European air space was still closed because the ash concentrations over Europe were not known at the time.However, both measured and modeled total ash concentrations along the flight track were considerably lower than the current threshold (< 0.2 mg m −3 ) for the "Normal" flying zone (European Commission, 2010).Simulations using either meteorological input data reproduced the observed ash layers near Leipzig, Stuttgart and Munich.Despite model over-and underestimates for sections of the flight, the overall bias is small and within systematic measurement uncertainties (Schumann et al., 2010).The a posteriori mass concentrations are furthermore closer to the measured values than the a priori mass concentrations.Higher ash concentrations were observed by the Falcon on other flights.In particular, extensive ash layers were probed over the North Sea on 17 May (Fig. 9) and over the North Sea and Germany on 18 May (Fig. 10).Besides near Iceland, these two flights detected the highest ash concentrations of the entire airborne campaign (Schumann et al., 2010).The model captured these ash layers and there is relatively good quantitative agreement between the a posteriori model results and the measurements.Plumes of volcanic ash were encountered in 34 events on 10 flights in different regions over Europe and the North Atlantic, allowing a statistical comparison between the model and the measurements.Schumann et al. (2010) reported mean ash concentrations for 12 ash encounters in their Table 3 and referred to the remaining 22 ash encounters in the text.One case was excluded from our analysis, namely the descent into the top part of a dense ash plume relatively close to Iceland on 2 May.The descent was abandoned for safety reasons when high ash concentrations were encountered and, thus, the aircraft spent only three minutes in that part of the plume (Schumann et al., 2010).The model produces orders of magnitude too low ash concentrations at the exact location of the ash encounter.However, this is not surprising given the sampling strategy and, therefore, this case was excluded from the comparison.Introduction

Conclusions References
Tables Figures

Back Close
Full The average ash concentration during the remaining 33 plume encounters was 64.3 µg m −3 .The ECMWF-based model sampled along the flight legs with ash encounters produces mean a priori concentrations of 29.5 µg m −3 and mean a posteriori concentrations of 50.9 µg m −3 (Fig. 11), even though, overall, the total a posteriori ash emissions are lower.The improved agreement with the observations comes from the more accurate a posteriori plume positions, which is also reflected by an increase in the Pearson correlation coefficient from an insignificant value of 0.18 to a significant (at better than the 0.05% level) value of 0.61.Notice that a small underestimation by the model is expected and does not necessarily indicate a model bias.The data for the comparison was selected by screening the entire observation data set (including gasphase measurements) for volcanic plumes.An unbiased but imperfect model would underestimate the observations in such a comparison, since slight displacements of the modeled plumes would lead to sampling lower-concentration parts of the plume in the model as compared to the observations.

The area over Europe affected by volcanic ash
Based on the a posteriori simulations, we determine the area over Europe (defined here as the longitude-latitude box 10 • W-30 • E, 36 • N-60 • N) where volcanic ash was present somewhere between the surface and 13 km altitude at concentrations exceeding thresholds set (ex post facto) by the European Commission (2010) (Fig. 12).During three episodes in April and May, volcanic ash concentrations at some altitude in the atmosphere exceeded the limits for the "Normal" flying zone in up to 14% (6-16%), 2% (1-3%) and 7% (4-11%), respectively, of total area over the European domain, which triggered partial closures of European air space.For a limit of 2 mg m −3 only two episodes with fractions of 1.5% (0.2-2.8%) and 0.9% (0.1-1.6%) occurred, while the current "no-fly" zone criterion of 4 mg m −3 was rarely exceeded.Clearly, the area affected by volcanic ash depends strongly on the choice of the threshold value.at present there are no agreed values of ash concentration which constitute a hazard to jet aircraft engines."Patches of highly concentrated ash were present over Europe (10 • W-30 • E, 36 • N-60 • N) during both April and May, and it is important for aviation to avoid them.

Conclusions
In this paper we have, for the first time, objectively determined the ash emissions of a volcanic eruption as a function of time and altitude.For this, we have used an inverse method which optimally merges a priori information on the emissions based on observed plume heights and an eruption column model, satellite observations of total atmospheric ash columns and sensitivity calculations with a dispersion model.We applied the method to the eruption of Eyjafjallajökull in April-May 2010, which caused severe problems for aviation over Europe because of its extensive ash emissions.We evaluated the model simulations using webcam observations of the eruption column, ground-based and space-based lidar observations, and aircraft observations and found that the inversion consistently improved the quality of the model simulation.
From the inversion, we obtain a total fine ash emission of 8.3 ± 4.2 Tg for particles in the size range of 2.8-28 µm diameter.Extrapolation to the size range of 0.25-250 µm yields a total ash emission of 11.9 ± 5.9 Tg.
We find that during three episodes in April and May, volcanic ash concentrations at some altitude in the atmosphere exceeded the limits for the "Normal" flying zone in up to 14% (6-16%), 2% (1-3%) and 7% (4-11%), respectively, of the European area.For a limit of 2 mg m −3 only two episodes with fractions of 1.5% (0.2-2.8%) and 0.9% (0.1-1.6%) occurred, while the current "No-Fly" zone criterion of 4 mg m −3 was rarely exceeded.
Our methodology has broad applicability.It is efficient enough for real-time application and could supply ash forecasting models (Witham et al., 2007) with an objectively derived quantitative source term.Improved forecasts would then allow more effective Introduction

Conclusions References
Tables Figures

Back Close
Full emergency response.For retrospective analysis, more accurate knowledge about the spatial distribution of volcanic emissions in the atmosphere would improve the quantification of their radiative (Robock, 2000) and environmental (Durant et al., 2010) impacts.For instance, improved estimates of ash deposition into the ocean would allow a better quantification of ocean fertilization (Duggen et al., 2010), which could be relevant especially for the Icelandic Sea which may be iron-limited (Nielsdóttir et al., 2009).

Movies of ash dispersion
A Movies comparing ash total columns retrieved from SEVIRI with FLEXPART a priori and a posteriori simulations are available from http://zardoz.nilu.no/~sabine/APRIL.gif and http://zardoz.nilu.no/~sabine/MAY.gif.These movies show that the a posteriori results are always in better agreement with the SEVIRI observations than the a priori observations, although sometimes differences are relatively small, showing that our a priori emissions are also of high quality.Introduction

Conclusions References
Tables Figures

Conclusions References
Tables Figures

Back Close
Full Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | reports, and three-hourly radar data (Petersen and Arason, Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | distribution.We then fitted two lognormal curves to the optimized emitted size distribution to specify our initial ash size distribution d M d logD p Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Schumann et al. (2010).The first research flight was performed over Germany on 19 Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | While the European Commission (2010) has set thresholds during the eruption, the International Civil Aviation Organization states in the year 2007 (ICAO, 2007): "Unfortunately, Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Figures A1-A5 show four additional comparisons between cross-sections of CALIPSOtotal attenuated backscatter and simulated ash concentrations based on both the a priori and a posteriori emissions.FigureA3is the same case as shown in the main part of the paper (Fig.7) but here also the a priori model results are shown, allowing to assess the model improvement by the inversion.The model simulations based on our a posteriori emissions were compared with observations from ground-based lidars in Germany (Leipzig, Munich, Hohenpeissenberg) on 16-17 April.FiguresA6-A8show the model results, which are compared in the following with published measurement data.The ash cloud mostly extends from about 2-6 km altitude, and a clearly simulated maximum occurs over Leipzig on 16 April at 15:00 UTC (Fig.A6b,c) as also observed by the Leipzig lidar (Fig.1inAnsmann et al., 2010).Especially the shape of the distinctly tilted ash cloud is in very good agreement with the lidar-observed ash cloud.Furthermore, for the ECMWF-based simulations the simulated ash concentrations (≈ 500-1200 µg m −3 ) over Leipzig and Munich (Fig.A7) agree with the mass concentrations retrieved from the lidar observations(Figs.4 and 5  in Ansmann et al., 2010), where the indicated mass concentrations were on the order of 1000 µg m −3 in the densest ash layer and short-term peak values were close to 1500 µg m −3 .Other reported(Gasteiger al., 2010)  maximum mass concentration of atmospheric dispersion models using the 1 November 2004 Grimsvötn eruption,

Fig. 1 .Fig. 2 .Fig. 4 .Fig. 5 .Fig. 6 .
Fig. 1.Ash mass size distributions.Particle size distributions measured in an ash sample collected at the ground at a distance of 60 km from the vent at 11:30 UTC on 15 April 2010 (black asterisks) (see http://www.earthice.hi.is/page/ies_EYJO2010_Grain) and by the DLR Falcon aircraft 450 km downwind of the vent at 15:00 UTC on 2 May 2010 as described by Schumann et al. (2010) (violet plusses).The blue lines are a fit through the initial size distribution required to match the ground sample distribution with the model (thick blue line) as well as shifted distributions used for sensitivity analyses (thin blue lines).The red line shows the sensitivity of the SEVIRI and IASI satellite retrievals to ash particle size.The yellow and grey boxes indicate the size range used for the inverse modeling and for other modeling, respectively.Only the ash mass in the yellow size range (ca.70% of the total mass for the fitted reference distribution) is constrained by the satellite measurements.

Fig. A1 .
Fig. A1.Comparison between satellite-retrieved and modeled ash on 8 May at 4 UTC.Top row: Ash total columns retrieved from geosynchronous satellite observations between 4-5 UTC (left panel) and a posteriori (middle panel) and a priori (right panel) modeled ash total columns.The red line shows the track of CALIPSO overflying the ash plume.Middle row: Vertical cross-sections through the a posteriori (left panel) and a priori (right panel) modeled ash concentrations.Bottom row: Vertical cross-section along the red line shown in top panels of total attenuated backscatter at a wavelength of 532 nm obtained from CALIPSO (colored) superimposed with 20 µg m −3 (black line) and 200 µg m −3 (white line) isolines of modeled volcanic ash concentrations.The left panel shows the a posteriori, the right panel shows the a priori model results.Notice that many features seen in the CALIPSO data are meteorological clouds, not ash.