Volcanic forcing for climate modeling: a new microphysics-based data set covering years 1600–present

. As the understanding and representation of the impacts of volcanic eruptions on climate have improved in the last decades, uncertainties in the stratospheric aerosol forcing from large eruptions are now linked not only to visible optical depth estimates on a global scale but also to details on the size, latitude and altitude distributions of the stratospheric aerosols. Based on our understanding of these uncertainties


Introduction
For climate simulations covering the last centuries, the volcanic forcing is a crucial external influence that needs to be included as precisely as possible.It is now commonly accepted that explosive volcanism is, together with the solar variability, one of the most important pre-industrial forcings explaining a large part of the decadal variance in the past centuries (Hegerl et al., 2003).Moreover, on the centennial timescale its strong effect on oceanic heat content has been recently stressed (Gleckler et al., 2006;Stenchikov et al., 2009;Gregory et al., 2010).However, as stated in the latest IPCC report (Forster et al., 2007), the level of scientific understanding of the volcanic influence is still low, one reason being the large uncertainties in the volcanic aerosol data before the satellite period.
Volcanic stratospheric aerosols impact climate through several mechanisms related to the increased scattering of solar radiation and increased absorption of solar and terrestrial infrared radiation.These changes modify the radiative budgets in the stratosphere and troposphere and lead to important temperature changes which then impact ozone and potentially stratospheric water vapor.Dynamical responses to F. Arfeuille et al.: Volcanic forcing for climate modeling (1600-present) the final net radiative forcings further modify the climate response and can lead to regional anomalies of temperature (e.g., winter warming; Robock and Mao, 1992;Stenchikov et al., 2006) and precipitation (Robock and Liu, 1994;Joseph and Zeng, 2011).
To describe the diverse mechanisms involved as accurately as possible, general circulation models (GCMs) and chemistry-climate models (CCMs) need information about the mass and size distributions of the aerosols for every latitude and altitude.Through the Mie theory (Mie, 1908), relevant optical parameters can then be calculated: extinction, single scattering albedo and asymmetry factor.In addition, the surface area density (SAD) of the aerosols can be used to calculate heterogeneous reactions affecting ozone.
This information can be satisfactorily derived from satellite data for the most recent period, particularly for the Pinatubo 1991 eruption, despite saturation issues (Arfeuille et al., 2013); lack of extensive extinction measurements before 1979 prevents a direct inference of the volcanic aerosols optical parameters.For the 1883-1978 period, volcanic aerosol data sets are based on isolated solar optical depth data (Sato et al., 1993) and before 1883 ice core aerosol mass estimates are used (Sato et al., 1993;Ammann et al., 2003;Crowley et al., 2008;Gao et al., 2008).However, this limited information is not sufficient to derive the necessary optical parameters by itself and has to be supplemented by parameterizations and assumptions on the evolution of the size and spatial distribution of the aerosols with time.The most common approaches used in volcanic aerosol data sets and GCM modeling studies to create and implement the volcanic forcing before the satellite era are summarized in Sect. 2.
In this study we present a new volcanic stratospheric aerosol data set using the latest ice core data from Gao et al. (2008) and covering the 1600-present period.The approach is based on simulations made with the state-of-the-art twodimensional AER aerosol model developed at Atmospheric and Environmental Research Incorporation, Lexington, MA, USA (Weisenstein et al., 2007).It allows the size distributions to be directly derived from the microphysical processes described in the model for each individual eruption, while the spatial distribution is also dependent on the boundary conditions (latitude and altitude of SO 2 injection, date of eruption) and climatological winds.
This data set differs in its conception from previous volcanic forcing reconstructions in several ways which are described in Sect.3. Technical details on the boundary conditions used for each of the 26 eruptions present in the data set are given in Sect. 4 and a description of the AER model is presented in Sect. 5. Main results of the reconstruction are then shown and compared to previous studies in Sect.6.

Previous data sets and approaches
As summarized in the IPCC AR4 report, recent climate simulations implement the volcanic forcing in different ways, from the direct use of a global or latitudinally dependent radiative forcing (Crowley et al., 2003;Goosse and Bradley, 2005;Stendel et al., 2006), to prescribed aerosol optical depth (AOD) changes with interactive calculations of the long-wave and shortwave radiation budgets (Tett et al., 2007).
For climate simulations covering a period prior to 1883, the start of the using optical depth ground measurements in the Sato et al. (1993) GISS data set, direct parameterizations are often used to convert aerosol loadings to optical depths (i.e., by dividing the aerosol mass by 150 Tg as proposed by Stothers, 1984) or to convert optical depth to radiative forcing by multiplying it by −20, as in Wigley et al. (2005).An additional factor can then be used for the conversion in order to take into account the non-linearity of the shortwave forcing-aerosol mass relationship due to the decreasing scattering efficiency as particles grow through enhanced coagulation (Hegerl et al., 2011).
The Crowley et al. (2008) volcanic aerosol data set provides effective radii in addition to aerosol optical depth (AOD) data at 550 nm.For the models using this data set, one approach consists in distributing the AOD at chosen vertical levels to estimate extinction coefficients (km −1 ).The extinctions, single scattering albedos and asymmetry factors for all the wavelength bands of the model are then calculated assuming a log-normal size distribution of the aerosols.The effective radii from the data set and a constant distribution width complete the needed information.This methodology was, for instance, used in Jungclaus et al. (2010), with an assumed distribution width of 1.8 and the AOD distributed over the 20-86 hPa pressure levels.The calculated spectral extinctions are very sensitive to the choice of distribution width, which can vary from a value of 1.8 in background quiescent conditions to 1.2 for the volcanic large mode after a strong eruption, and should be chosen with care.It can be noted that using a distribution width of either 1.8 or 1.2 and the surface area density data and effective radius data calculated from SAGE II (Thomason et al., 1997;SPARC, 2006) in order to derive the aerosol size distribution leads to too large extinctions in the infrared (Arfeuille et al., 2013).Recently, some studies focused on simulation past volcanic eruptions using climate models incorporating aerosol microphysical processes (Laki eruption: Oman et al., 2006;Schmidt et al., 2010: Pinatubo eruption: Aquila et al., 2012).
Until every GCM/CCM includes coupled stratospheric aerosol microphysical modules, there is a need for an external data set providing aerosol temporal and spatial evolutions as well as size distributions (Timmreck, 2012).The data set should contain spectrally resolved optical properties, surface area densities and mean radii for chemical interactions.The forcing data sets most used in recent modeling studies and covering the pre-satellite period are now described in Sects.2.1-2.4.

Sato et al. (1993)
The Sato et al. (1993) optical depth data set (and updates) cover the 1850-present period.It provides optical depths (in four 5 km bands from 15 to 35 km altitude) at 550 nm with varying latitude resolution over time.It is based on several sources.In the earlier time period (1850-1882), optical depths are derived from Mitchell's index (Mitchell, 1970) based on the volume of ejecta.The AOD is then spread along every latitude, except for the Askja (Iceland) 1875 eruption (30 • N-90 • N).In 1883-1959 ground-based measurements of atmospheric extinctions are used.Most of the stations then are in the mid-northern latitudes.Southern Hemisphere AOD are estimated based on the latitudes of large known eruptions.Hence for a tropical eruption the Southern Hemisphere AOD is taken to be the same as measured in the Northern Hemisphere.It can be noted that the Krakatau 1883 eruption was originally described using a single measurement station in Montpellier (44 • N), but is since a 2002 update derived from the latitudinal distribution of the Pinatubo aerosol cloud.For 1960-1978, Southern Hemisphere stations could also be used as well as lunar eclipse data.From 1979 on, satellite measurements are used.However, it can be noted that a new version of the SAGE II satellite data is now available (SPARC, 2006) which differs substantially from the one used in Sato et al. (1993) (Arfeuille et al., 2013).

Ammann et al. (2003)
The Ammann et al. (2003) volcanic data set  provides AOD at 550 nm with individual evolutions of the volcanic aerosols for each eruption.Latitudinal transport is parameterized.Transport from the tropics to the midlatitudes is assumed to take place in the lower stratosphere in the winter hemisphere only.The temporal evolution is constrained by specific e-folding decays (e.g., 12 months in the tropics) and a linear buildup.This simple model of the evolution of the aerosols cannot describe the details of the transport and sedimentation of the aerosols but compared to Sato et al. (1993) an advantage is the absence of changes in data sources which can introduce a factitious decadal variability (Ammann et al., 2003).The peak optical depths are scaled from the peak aerosol loadings (from Stothers, 1996;Hofman and Rosen, 1983;Stenchikov et al., 1998) assuming a fixed 75 wt% aerosol composition and a fixed aerosol size distribution (with an effective radius of 0.42 µm).Crowley et al. (2008) (this data set is further described in Crowley and Unterman, 2012) provides monthly aerosol optical depth at 550 nm and effective radius for four latitudinal bands for the 850-2000 period.Optical depths are de-rived using a calibration from Antarctic sulfate records of the Pinatubo and Hudson 1991 eruptions against satellite based optical depth from Sato et al. (1993).In total 22 ice cores are used.The obtained time series was validated for the Krakatau (1883), Santa Maria (1902) and Katmai (1912) eruptions against Stothers (1996).Effective radius growth and decay are derived from Pinatubo 1991 satellite observations and scaled using microphysical calculations from Pinto et al. (1989) for eruptions larger than Pinatubo.

Gao et al. (2008)
The Gao et al. (2008) study which analyzed 54 polar ice cores led to the release of two distinct data sets: an aerosol injection forcing which describes the mass of aerosols released for each volcanic eruptions for both hemispheres (for the 500-2000 period), and a data set providing latitude/altitude distributions of aerosol concentrations.Values from the first data set are derived from the ice core measurements via a calibration factor.
The hemispheric calibration factors for Greenland and Antarctica cores, linking the sulfate depositions measured in the ice for each polar region (in kg km −2 ) and the total mass of stratospheric aerosol formed after the eruption (in kg).The hemispheric calibration factors cannot be calculated from the same method due to limitations in the measurements in the recent period.In Gao et al. (2008), Northern Hemisphere (NH) calibration factors for tropical eruptions are derived from nuclear bomb tests, while Southern Hemisphere (SH) factors are calculated using Mt.Pinatubo eruption data.In the Gao et al. (2008) data set, the calibration factor to convert the sulfate deposition in NH and SH to the corresponding hemispheric stratospheric loading is set to a mean value of 1×10 9 km 2 for all tropical eruptions regardless of the NH/SH partitioning (which, however, should modify the hemispheric calibration parameters Gao, 2008).Hence, ice cores hemispheric estimates can be challenged and one example of the uncertainty arising from the differences in partitioning is the discrepancy for the Agung eruption.For this eruption, the ice core calibration leads to an almost even hemispheric partitioning of the aerosols after calibration while observations indicate a much larger loading in the SH (Stothers, 2001).This discrepancy could also be attributed potentially to a concomitant high latitude eruption (Surtsey, Iceland, November 1963, Gao et al., 2007).Another example of uncertainties in hemispheric aerosol values is the difference for the estimation of the Tambora (1815) aerosol hemispheric asymmetry between Gao et al. (2008) (more aerosols in the NH) and Crowley et al. (2008) (more aerosols in the SH).
The new ice core data set of Sigl et al. (2013) also exhibits larger values of sulfate deposition in the Southern Hemisphere after the Tambora eruption.Moreover, it provides sulfate deposition values from the Pinatubo eruption from both poles which indicates stronger deposition in Greenland than in Antarctica.This would lead to a larger calibration factor for the Northern Hemisphere compared to the one for the Southern Hemisphere and hence accentuate the asymmetry towards the south for the Tambora aerosol cloud (with more than 75 % of the Tambora aerosols in the Southern Hemisphere).
As noted in Gao (2008), uncertainties in hemispheric values due to variations in hemispheric partitioning of tropical eruptions partly compensate in Crowley et al. (2008) for the calculation of the global injection masses as an underestimation of SH loading would correspond to an overestimation of the NH loading.This is because the 1 × 10 9 km 2 calibration factor used in the study is the mean of calibration factors calculated when taking into account potential asymmetry in the hemispheric distribution (from 1/3 to 2/3 of the global mass in a specific hemisphere).
The total mass injected (but not the partitioning into hemispheres) is used in our study to infer the SO 2 emission inputs for the AER-2-D model and will be described in Sect.4.1.
The second data set provided is a full latitude/altitude distribution of the aerosol concentration, varying monthly.The altitude-latitude evolution is derived from a parameterization of the transport between tropics, extratropics and poles (for instance transport from tropics to extratropics is set to a 50 % mass ratio per month).No inter-hemispheric transport is allowed and separate injections are made in the two hemisphere in order to match the ice cores mass estimates.A linear build up of the total aerosol mass during the first four months after the eruption is prescribed while an exponential decrease of the stratospheric aerosol mass is assumed, with a global e-folding time of 12 months.The e-folding time in the tropics is 36 months (little sedimentation), 12 months in the extratropics (tropopause folding), three months for the winter polar region (stronger subsidence) and six months for the summer polar region.Gao et al. (2008) use vertical distributions from the Pinatubo 1991 eruption based on the assumption that the vertical distributions of volcanic aerosols have a similar shape (Pinatubo, Agung, El Chichón and Fuego eruptions observations) and peak only at slightly different heights depending on altitudes of initial injections.The variations in initial sulfur altitude distribution from one eruption to another does not strongly affect the aerosol global surface radiative forcing later on but may affect the location of the stratospheric warming and potentially the transport of the aerosols in the stratosphere (see Sect. 6.2).

Specificities of our data set
In regards to the latest studies on the influence of volcanic eruptions on the stratosphere and climate and the current implementation of volcanic aerosol forcings in GCMs and CCMs, volcanic aerosol data sets should describe as precisely as possible: 1.The global average AOD peaks in the shortwave: in our study this is mostly dependent on the estimated global SO 2 mass release, derived from the state-ofthe-art ice core reconstruction of Gao et al. (2008).
Ice core sulfate estimates present uncertainties which can be important due to variations in sulfate depositions and volcanic signal extraction methods.The Gao et al. (2008) methodology made use of large number of ice cores in order to reduce these uncertainties (Gao et al., 2007).As seen in Fig. 1b, the relationship between aerosol mass and AOD peak is stable, as the maximum AOD does not depend to a large extent on spatial and temporal variations in individual volcanic eruptions.It should be noted that important differences in the magnitude of many eruptions are present in existing data sets.For instance, the Krakatau eruption of 1883, previously thought to have released more SO 2 than Pinatubo, is now downgraded to 2/3 of Pinatubo release by Gao et al. (2008).The Cosiguina 1835 eruption while being present in most ice core reconstructions (Zielinski, 1995;Gao et al., 2008) is absent of the Stothers (2007) data set based on lunar eclipses.An aerosol optical depth background value is also present in our data set in conformity to observations during volcanically quiescent periods.
2. The altitude-latitude distribution: in this study the spatial distribution is modeled through the AER-2-D microphysical scheme and depends on nucleation, condensation, coagulation, sedimentation and climatological transport.Initial distribution and subsequent transport depends on the specific latitudinal position of the volcanoes and dates of eruption.
Spatial distribution evolves with time as aerosol size distribution also evolves due to microphysical processes; sedimentation rates depend on calculated particle size distributions.The 1.2 km vertical resolution of the aerosol model AER-2-D leads to a good accuracy in the vertical profiles of the extinction coefficients, while data sets providing total AOD as Crowley et al. (2008) or 5 km altitude resolution AOD as Sato et al. (1993) cannot resolve the evolution of the aerosol altitude and hence the location of the stratospheric warming for instance.
Differences in the initial altitude distribution of the volcanic SO 2 is also accounted for in this study, using information from the one-dimensional Plumeria volcanic plume model (Mastin et al., 2009) when data is available.
3. The size distribution and the absorption in the infrared: it is important to characterize the mass scattering efficiency and long-wave forcing of the aerosols as well as their sedimentation speed.Using an aerosol microphysical model allows us to describe the size distributions at each location and time for every specific eruption.The AER-2-D model has a 40 size bins resolution.The absorption / scattering ratio and its dependence on the eruption size (Timmreck et al., 2009) is hence taken into account, preventing an overestimation of the shortwave insolation reduction and underestimation of the long-wave absorption following large eruptions.Figure 1a shows the variations in effective radii for each of the 26 eruptions as simulated by the AER model.Comparison of the extinction coefficients in the visible and infrared from the AER-2-D model with observations were done for the Pinatubo eruption and presented in a separate article (Arfeuille et al., 2012, AER_7 simulation).The extinctions coefficients in the visible and infrared are in reasonable agreement with the observations from SAGE II (1.02 µm), HALOE (5.26 µm), and ISAMS (12.1 µm), except at the equator in the first months after the eruption where there is an overestimation of the extinction coefficients.
The final released 1600-present volcanic aerosol data set consists in monthly, three-dimensional fields of extinction coefficients, single scattering albedos and asymmetry factors for different wavelength bands calculated upon request.Surface area densities for heterogeneous chemistry computations are also made available.

Injection data
Table 1 shows the sulfate aerosol mass estimates from Gao et al. (2008) for Northern and Southern Hemispheres for each eruption recorded.
From this data, the Laki 1783 eruption was not used because according to previous studies, the SO 2 release was either mainly tropospheric (Highwood and Stevenson, 2003) or distributed in the upper troposphere/lower stratosphere region with, however, 85 % of the formed aerosols being removed in the summer/fall 1783 (Thordarson et al., 2003).Due to difficulties to model accurately such an eruption with our approach, we chose not to implement it in the final data set.Further work could be done in the future to add this particular eruption, for instance following the work of Oman et al. (2006) and Schmidt et al. (2010).The eruptions in 1925 and 1943 were removed as they do not appear in any available reconstruction and no known large eruption were recorded around these periods.The Pinatubo eruption is not derived from ice core measurements as no ice core data is available from Greenland and was added in Gao et al. (2008) based on satellite data.In the same manner, we added the El Chichón eruption of 1982 and the 1976 eruption attribution was changed to the large tropical eruption of Mt.Fuego of October 1974 as dates and aerosol amounts correspond well.

Calculation of SO 2 amount
First, the aerosol mass is converted to a SO 2 injection mass (needed by the AER-2-D model) for each volcanic eruption: where m is the mass of the compounds, M their molecular weight and wt% the percentage in mass of H 2 SO 4 in the aerosols.Sulfuric aerosol total amount (NH + SH) is used to calculate SO 2 initial release assuming a 100 % conversion efficiency and a H 2 SO 4 average composition of 75 wt%.This is used for consistency with the average composition assumed for the ice core parameterization in Gao et al. (2008).For the optical properties calculations, the H 2 SO 4 wt% is set to match high temperatures after a strong volcanic eruption (Pinatubo 1991 eruption simulation) in order to minimize the absolute error arising from not using the actual composition (which depends on the temperature).

Attribution to specific eruptions
Some tentative attributions of sulfate peaks in ice cores to specific historical eruptions were made in this study when large historical eruptions (data from the Global Volcanism Program of the Smithsonian Institution) have a consistent location and timing of eruption compared to the loadings given by Gao et al. (2008).This is the case for the Parker (January 1641, Philippines), Gamkonora (May 1673, Indonesia), Serua (April 1693, Indonesia), Katla (October 1755, Iceland), Makian (September 1760, Indonesia) and Babuyan Claro (1831, Philippines) eruptions.The Sigl et al. (2013) study with new high-resolution ice core data finds bi-polar signals in the deposition following 1693 and 1831 which correspond well to our attribution to the tropical eruptions of Serua and Babuyan Claro, respectively (while these 2 eruptions are only recorded in one hemisphere in Gao et al., 2008).The results from Sigl et al. (2013) would indicate then higher total emissions than prescribed here for these two eruptions.

Timing and latitude of injection
The latitudinal distribution of the initial SO 2 injection is given by the location of the volcanoes.For tropical eruptions, it was observed that a few weeks after the eruption, SO 2 spread quickly inside the tropical pipe region.Injection in the AER-2-D model was made in the 15 • S-0 • and 0 As in Gao et al. (2008), the location of an eruption is set to tropical if unknown but recorded in both hemispheres.This is the case only for the Unknown 1809 eruption.In the same way as in Gao et al. (2008), the timing of an eruption is set to April if the eruption is unknown, except for the eruption date of the 1809 eruption which was set to February (Cole-Dai et al., 2009).As only high latitude eruptions in our data set have unknown dates of eruptions, no tropical eruption was artificially set to an April eruption, hence no bias towards a specific hemispheric partitioning is introduced in the data set.

Altitude of injection
In the tropical reservoir region where most of the SO 2 mass from large tropical eruptions is released, even a small altitude difference can theoretically affect the timing of the transport to the extratropics.The altitude of injection is taken from the Plumeria model (Mastin et al., 2009) when input parameters as mass of ejecta and duration of eruptions or eruption rate estimates are available.Otherwise, the same initial SO 2 vertical distribution as Pinatubo is used.Results summarized in Table 2 indicate for instance that the neutral buoyancy height of the SO 2 release by the Agung 1963 eruption was significantly lower than that of the Pinatubo eruption.The Plumeria model is useful in order to get estimates of the neutral buoyancy height of eruption cloud but results are subject to caution as dynamical parameters used to infer the mass flux of the eruptions are uncertain.A limitation is the lack of aerosol heating impact on the SO 2 lofting in our approach (Aquila et al., 2012).However, initial SO 2 heights given here represent rather the altitude of the eruption cloud after the large initial lofting occurring in the first weeks (Pinatubo: 24-25 km), when the volcanic cloud has not yet circled the globe.Our calculated heights are useful to take into account variations between eruptions which are due to differences in the eruption dynamics, independently of the subsequent aerosol induced lofting.Since the power of an explosive eruption is not directly linked to the amount of SO 2 released, and hence to the amount of induced lofting, these two effects do not necessarily add up.An additional difficulty in estimating the SO 2 altitude distribution arises from the fact that some eruptions involve distinct plinian phases as well as a co-ignimbrite phase (Herzog and Graf, 2010).This was the case for Mt.Tambora 1815 eruption.Ejecta mass released for each of the phases are then hard to infer and an additional difficulty is that the specific SO 2 contents of each phase depend on the sulfur partitioning in the subvolcanic magma reservoir.

AER model
The two-dimensional sulfate aerosol model developed at AER (Atmospheric and Environmental Research Incorporation, Lexington, MA, USA) was used to model 26 volcanic eruptions in the 1600-present period.For validation of the model against observations see Weisenstein et al. (1997), SPARC (2006) and Weisenstein et al. (2007).The model includes the following species: SO 2 , OCS, DMS, H 2 S, and CS 2 .For the volcanic eruptions the sulfur is put into the stratosphere in the form of SO 2 .The model uses precalculated OH concentrations and photolysis rates, derived from a model calculation with standard stratospheric chemistry (Weisenstein et al., 1997).Reaction rates are according to Sander et al. (2000).Sulfuric acid aerosols (H 2 SO 4 /H 2 O) microphysics is modeled in a global domain which extends from the surface to 60 km with 1.2 km vertical resolution, though thermodynamics ensures that sulfuric acid aerosols evaporate above about 35-40 km.The horizontal resolution is 9.5 • .The model accounts for significant recycling of gaseous H 2 SO 4 into SO x in the upper stratosphere via photolysis.The sulfate aerosols are treated as liquid binary solution droplets.The size distribution (with 40 size bins spanning the range 0.4 nm-3.2 µm by volume doubling) of aerosol particles is defined by the microphysical processes of nucleation, condensation, evaporation, coagulation, and sedimentation.The shape of the size distribution is unconstrained.Aerosols are removed below the tropopause by washout and surface deposition.Transport parameters are based on Fleming et al. (1999) and were calculated from observed O 3 , water vapor, zonal wind, and temperature for climatological conditions (transient for Pinatubo eruption).The lack of interaction of the sulfate aerosol transport with induced radiative perturbations and the lack of information on past QBO phases can limit the accuracy of the aerosol transport in the model.The planetary and gravity waves drag is included.

Size distributions
The influence of the various volcanic SO 2 injection masses on the effective radii of the aerosols and subsequently on the evolution of the AOD is discussed here.Figure 1a shows the increase in effective radius with increasing eruption size due to intensification in the coagulation process from small to large eruptions because of increased aerosol concentrations in the stratosphere.We can note that Northern Hemisphere extra-tropical eruptions tend to exhibit smaller effective radii www.clim-past.net/10/359/2014/Koyagushi et al. (1993); Pallister (1992); Baker (1996); Newhall (1996); Self et al. (2004); Stix and Kobayashi (2008).
than tropical eruptions for a given SO 2 erupted mass.This could be due to the upwelling in the tropics which can reduce the removal rate of large aerosols.The peak AOD (global average, Fig. 1b) ratio between two eruptions follows very closely the relation: where M represents the mass of stratospheric aerosols and AOD the peak of the optical depth at 550 nm.This idealized relationship is a consequence of the total aerosol mass being proportional to the volume of the aerosols while the optical depth at 550 nm is proportional to the total surface area.
The actual mean equation linking our calculated global mean peak AOD at 550 nm to the aerosol masses given in Gao et al. (2008) is The AOD values for each individual eruption can differ from the mean relationship between aerosol mass and AOD by 20-30 % due to differences in the latitude and altitude of the injections.The AOD at 550 nm of the 26 volcanic eruptions simulated for the 1600-present period have large variations in terms of altitude/latitude distribution and time evolution (see Sects.5.2 and 5.3).
Following Toohey et al. (2011), the cumulated optical depths for the two years after each eruption are computed (sum of 24 monthly averages).Figure 2a shows the relation: The diminished power law compared to Eq. ( 3) shows the effect of the decreased scattering efficiency and faster sedimentation of the large aerosols formed as SO 2 emissions get larger.Using Eq. ( 4) to calculate the cumulative AOD for 17 and 700 Mt SO 2 eruptions as simulated in Toohey et al. (2011) using MAECHAM5-HAM, we obtain optical depths of 2.76 and 18.63, respectively, which is in good agreement with the results from Toohey et al. (2011) for the 17 Mt eruption and lower than their calculation for the 700 Mt case (global mean 2 yr cumulative AOD of approximately 20-24 depending on season of eruptions in Toohey et al., 2011).A tendency of winter eruptions to lead to smaller southern extratropical mean AOD than spring and summer eruptions is seen in Fig. 2c.In the same way, the winter hemisphere produce larger AOD in the Northern Hemisphere extratropics compared to other seasons (Fig. 2b).This agrees well with the results of Toohey et al. (2011).

Altitudinal distributions
For each eruption, our approach provides altitude-resolved extinction coefficient values.Figure 3a and c shows the distribution of extinctions at 550 nm with altitude and time in the equatorial region for the Pinatubo (top) and Tambora eruptions (bottom).At this latitude, high extinctions values span approximately between 90 and 15 hPa.One interesting feature to point out is the faster sedimentation of the aerosols in the case of the Tambora eruption compared to Pinatubo due to a more intense coagulation forming larger particles.In Fig. 3b and d, the cumulative monthly extinctions for the first two years after the eruptions are shown.The transport to the extratropics and polar regions follows the Brewer-Dobson circulation and is influenced by the gravitational settling of the particles.The region of high extinctions in the extratropics is in the 150-40 hPa range.The tropical pipe region is visible with strong extinction gradients at its limits.
The exact influence of the initial altitude distribution on the aerosol transport is still unclear.Using different SO 2 injection altitudes does not change the aerosol altitude distribution to a large extent; however, some variations in the altitude of the aerosol cloud in the first months following the eruption.are observed.Also, while small differences of 1-2 km can theoretically influence the temporal and latitudinal evolution of the aerosols in the tropical reservoir, limitations in the model (e.g., vertical resolution, uncertainties in transport) prevent a precise testing and analysis of this influence.In addition model microphysical processes will tend to dampen altitudinal sensitivity, as sedimentation velocities are a strong function of air pressure.Sensitivity studies made for the Tambora eruption with initial injection at 23-25 km and 27-29 km tend to indicate a non-negligible impact on the aerosol distribution, with shorter time residence in the tropical region for the 27-29 km simulation and stronger extratropical transport.Increasing the altitude of injection from 24 to 28 km hence decreases the mass peak loading by 10-15 % in the tropics and increases it by 13 % and 8 % at 50 • S and 50 • N, respectively.

Latitudinal distributions
The season and latitude of volcanic eruptions highly impact the latitudinal distribution of the aerosols especially through the hemispheric partitioning of tropical eruptions.In the tropics, aerosols first spread in the tropical pipe region (Plumb, 1996) where exchange to the extratropics is limited and varies seasonally.The transport to the extratropics is made preferentially towards the winter hemisphere as the increased meridional gradient leads to increased planetary wave propagation into the stratosphere (Holton et al., 1995).Hence, in a first approxi-mation, a spring-summer tropical eruption tends to produce high aerosol concentrations in the Southern Hemisphere.Additionally, the latitudinal position of the tropical pipe in the first weeks after the eruptions is also of strong importance as it influences the early latitudinal distribution of the aerosols inside the tropics.
The latitudinal evolution of the most important eruptions are shown in Figs.4-10.The specific hemispheric partitioning of the Tambora eruption is for instance illustrated in Fig. 6, where the eruption which occurred in April 1815 produced an asymmetrical aerosol cloud in our model, with most of the mass in the Southern Hemisphere.The data set from Crowley et al. (2008) also depict this behavior for this eruption.
The AOD distribution for the largest volcanic eruptions in the last 400 yr are shown in Fig. 4 1963, El Chichón 1982) for our reconstruction and Crowley et al. (2008).From these 12 large volcanic eruptions, the El Chichón eruption, Krakatau, Unknown 1809, Serua and Huaynaputina eruptions show sizably different latitudinal evolution in the two data sets while the others eruptions show similar developments.It can also   be noted that the different ice core parameterizations in Gao et al. (2008) lead to a different hemispheric partitioning than in Crowley et al. (2008) for Santa Maria, Serua and Tambora (Table 1 and Figs . 4-6).Sigl et al. (2013) provide sulfate deposition values in both Greenland and Antarctica for the Pinatubo eruption.Results from the Sigl et al. (2013) ice core data would indicate a Huaynaputina (1600) eruption with an almost even hemispheric distribution (slightly to the south), and the Serua (1693) and Unknown (1809) eruptions with larger values in the Southern Hemisphere if we use the sulfate deposition they obtain from Pinatubo to derive new hemispheric calibration factors.

Agung 1963
A strong asymmetry was observed for the Agung March 1963 eruption by the solar extinction measurements compiled by Stothers (2001) as shown in Fig. 13.The extent of the asymmetry is only captured in our data set while Ammann et al. (2003) overestimates the Northern Hemisphere loading and Sato et al. (1993) and Crowley et al. (2008) underestimate the Southern Hemisphere value.

Pinatubo 1991
The latitudinal evolution of the Pinatubo 1991 eruption is also well captured by our method.For this eruption, transient winds instead of a climatology was used.This reduces the equatorial optical depth after the eruption but does not sizably impact the hemispheric partitioning of the aerosols.Figure 14 shows a comparison of the time evolution of Mt.Pinatubo optical depth for this study and Crowley et al. (2008) compared to Sato et al. (1993) data derived from satellite observations.The hemispheric partitioning is reproduced correctly in both data sets.
All data sets tend to overestimate the aerosol optical depth in the tropics (Fig. 14 and 12) compared to AOD derived from SAGE II data (optical depth data taken from SAGE_4λ in Arfeuille et al., 2013, to avoid gap-related issues).Our data set is showing the largest overestimation until end of 1991 while Ammann et al. (2003) data shows the largest values afterwards.A reasonably well-reproduced feature for the Pinatubo eruption is the temporal evolution of the AOD in both the tropics and extratropics compared to Sato et al. (1993) (Fig. 10).Crowley et al. (2008) do not capture precisely the transport timing to the extratropics, with too long of a residence time in the tropics and no time delay in the tropical to extra-tropical AOD.Their building of the AOD in the tropical region in June-August 1991 is too slow, with a time shift of approximately one month.Comparison with SAGE II data (Fig. 12) show also a general tendency of all data sets to overestimate the AOD of the Pinatubo eruption in the extratropics.
While the recent eruptions of Agung and Pinatubo are well reproduced, a bias towards the Southern Hemisphere is seen for the El Chichón eruption in our data set, where the observed higher loading in the Northern Hemisphere is   not reproduced by the model.Erroneous latitudinal distributions arise from the lack of knowledge on the specific stratospheric dynamics.Using a coupled aerosol-GCM model, Aquila et al. (2012) noted for instance the importance of the aerosol radiative influence on the transport, which is not taken into account here, but also general difficulties to infer transport across the equator.They found that for the El Chichon eruption, the low altitude of the aerosol cloud (which they explain by the limited SO 2 release and hence the limited increased uplift) may have prevented the transport to the SH through the mid-stratosphere pathway.In general, however, they note that the cross-equator transport is strongly dependent on the season of eruption.Figure 12 shows that both our reconstruction and Ammann et al. (2003) overestimate the Southern Hemisphere AOD compared to Sato et al. (1993).Crowley et al. (2008) which follows the Sato et al. (1993) data reproduce well the hemispheric partitioning of El Chichón (Fig. 6).
Despite differences in the latitude and time distribution of the aerosols, global average AOD of our reconstruction show similar values as Crowley et al. (2008) for many past eruptions as shown in Fig. 11a and b.Some small eruptions in 17th century and 19th century are present in Crowley et al. (2008) but not in the Gao et al. (2008) ice core estimates and hence not in our data set.Conversely, some small to medium eruptions in the 18th century are not present in Crowley et al. (2008) or of a reduced magnitude.Uncertainties due to the lack of information from the vertical wind structure in the tropics and QBO phases for historical eruptions is a limitation to obtaining accurate hemispheric partitions in the past.However, ice core estimates do not always provide a robust partitioning between hemispheres.Peak values of AOD and its temporal evolution are generally well reproduced in our reconstruction when compared to Sato et al. (1993), Stothers (2001) and SAGE II although it overestimates the optical depth values in the tropical region.Tropical mean AODs are generally higher than in Sato et al. (1993) and Crowley et al. (2008), and close to values from Ammann et al. (2003).Hemisphere extratropics (c).Black line is the AOD from this study, green line is the AOD from Crowley et al. (2008), blue line is the AOD of Sato et al. (1993), red line from Ammann et al. (2003), magenta line is the AOD derived from SAGE II observations.used in climate model simulations covering the last centuries in order to assess the extent of model deficiencies.The uncertainties in the volcanic forcing are still high due notably to lack of knowledge on the aerosol spatial and size distributions after large eruptions and this can lead to strong differences in the representation of both global and regional climatic anomalies in climate models.
In this study, the AER-2-D aerosol model has been used to simulate volcanic eruptions in the 1600-present period.Results for global mean optical depths at 550 nm are in overall agreement with previous reconstructions from Crowley et al. (2008), Ammann et al. (2003) and Sato et al. (1993).Previous studies exhibit large differences in hemispheric partitioning and tropical values as well as differences in detec-tion of some historical eruptions from the different ice core data sets.The AER-2-D model simulation of the Pinatubo eruption showed an overestimation of the tropical optical depth at 550 nm relative to SAGE II satellite data.The Agung eruption hemispheric partitioning is better captured by our method than in previous reconstructions as shown by comparisons with solar extinction measurements from Stothers (2001).However, some errors in the aerosol hemispheric partitioning for past eruptions probably remain due to the coarse resolution of the AER model and discrepancies in the climatological transport compared to actual wind fields.The latter includes a lack of the QBO and of the dynamical response to the radiative perturbation of the aerosols.Hence, information about timing of eruption, latitude of the volcano  and initial SO 2 height used in our study are valuable but not sufficient to constrain fully the latitudinal evolution of the aerosol clouds.Each individual eruption reconstruction gives a potential aerosol spatial distribution partially constrained by the relevant information available but which may deviate from reality in a number of cases.Despite these limitations, this method represents a credible alternative approach to the latitudinal parameterization constrained by hemispheric calibration from ice cores, which are found to be uncertain.
Our new volcanic forcing data set provides spectrally resolved optical properties and surface area densities as a function of latitude and altitude, derived from microphysical simulations which account for transport and sedimentation of the aerosols.This can help to refine the representation of the radiative, dynamical and chemical influence of the volcanic aerosols on climate compared to using data sets based on total optical depths at 550 nm (and radii), which do not give information on the altitude distribution of the aerosols and also require additional assumptions on the aerosol size distributions (e.g., distribution width) in order to provide spectral optical properties and SAD.The aerosol absorption at infrared wavelengths, which leads to a stratospheric warming and a non-negligible long-wave forcing at the surface, is proportional to the aerosol volume and is directly described by our data set.Similarly to extinction coefficients at 525 nm (compared to SAGE II observations), comparisons with data available for the Pinatubo eruption show an overestimation of the extinction coefficients at 5.26 and 12.1 µm for the tropics for the first months after the eruption but a good agreement with HALOE and ISAMS measurements later on (Arfeuille et al., 2013).For the extratropics, results are in good agreement with the observations at both visible and infrared wavelengths.
Further development could be applied to the volcanic aerosol data set.First, the development and testing of fully interactive stratospheric aerosol-climate models could lead to more realistic aerosol transport for past eruptions as it includes atmospheric background state and response to the aerosols forcing.Also, the recent increase in the background stratospheric optical depth since 2004-2005 could be implemented, either by increasing tropospheric SO 2 emissions or, following Vernier et al. (2011) and Solomon et al. (2011), by including several small volcanic eruptions.These recent volcanic eruptions include Ruang (September 2002), Manam (January 2005), Montserrat (May 2006), Kasatochi (August 2008), Sarychev (June 2009) and Puyehue (June 2011).Haywood et al. (2010) estimate for instance that the Sarychev eruption released 1.2 Tg of SO 2 , Carn et al. (2009) estimates around 1.4 Tg SO 2 release for Kasatochi.Additionally, future scenarios of volcanic aerosols forcing can be computed using the methodology described here and following the statistical approach by Ammann and Naveau (2010).A preliminary analysis showed that the Gao et al. (2008) ice core data follows well a Pareto distribution.Extension to the period before 1600 could also be done as ice core data improve and attribution of the deposits to specific eruptions can be done.Another future implementation could also involve the description of stratospheric water vapor (Joshi and Jones, 2009), and halogens release.A modeling study of the volcanic plume dynamics and its implications for stratospheric release of halogens showed that up to 25 % of the initial halogen release may cross the tropopause (Textor, et al., 2003); however, large uncertainties arise from the percentage of halogens reaching the stratosphere, which itself depends strongly on the humidity profile above the volcano and the dynamics of the plume.

Fig. 1 .
Fig. 1.Total mass of aerosol (H 2 SO 4 + H 2 O) in Mt versus (a) calculated peak of effective radius by the AER model for each of the 26 modeled eruptions over the 1600-2000 period (effective radius values are global mean average above 100 hPa, and (b) peak AOD at 550 nm.Black line shows the logarithmic fit for the tropical eruptions.Stars indicate spring eruptions, triangles summer eruptions, diamonds fall eruptions and crosses winter eruptions.Black color indicates tropical eruptions, while blue and green colors indicate extratropical eruptions in the Southern and Northern Hemispheres, respectively.

Fig. 2 .
Fig. 2. Total mass of aerosol (H 2 SO 4 + H 2 O) in Mt versus 2 yr cumulated AOD at 550 nm for global average (a), 30 • N-90 • N (b) and 90 • S-30 • S (c).Black line shows the logarithmic fit to the tropical eruptions.Stars indicate spring eruptions, triangles summer eruptions, diamonds fall eruptions and crosses winter eruptions.Black color indicates tropical eruptions, while blue and green colors indicate extratropical eruptions in the Southern and Northern Hemispheres respectively.

Table 1 .
Gao et al. (2008)osol data after ice core calibrations for Northern and Southern Hemispheres fromGao et al. (2008).Attribution of the aerosol peaks to specific volcanic eruptions in this study.Tg amounts correspond to the total aerosol mass (H 2 SO 4 + H 2 O) assuming a 75 wt% of H 2 SO 4 .
The latitudinal position of the tropical pipe given the time of year then affects the transport inside the tropics.When unknown, eruptions recorded only in the Northern/Southern Hemisphere are set respectively to high northern latitudes (55 • -65 • N) and Southern Hemisphere mid-latitudes (40 • S).These latitudes were chosen to best approximate common extra-tropical eruption locations (e.g., Iceland, Alaska, Chile, New Zealand, etc.) and, for instance, do not take into account potential local eruptions in Antarctica.
• -15 • N regions for Southern and Northern Hemisphere tropical eruptions, respectively.

Table 2 .
Altitude (km) of maximum injections for previous studies and Plumeria simulations, as well as the neutral buoyancy heights (NBH) in Plumeria.