Tracking aerosols and SO₂ clouds from the Raikoke eruption: 3D view from satellite observations

The 21 June 2019 eruption of the Raikoke volcano (Kuril Islands, Russia; 48^∘ N, 153^∘ E) produced significant amounts of volcanic aerosols (sulfate and ash) and sulfur dioxide (SO₂) gas that penetrated into the lower stratosphere. The dispersed SO₂ and sulfate aerosols in the stratosphere were still detectable by multiple satellite sensors for many months after the eruption. For this study of SO₂ and aerosol clouds we use data obtained from two of the Ozone Mapping and Profiler Suite sensors on the Suomi National Polar-orbiting Partnership satellite: total column SO₂ from the Nadir Mapper and aerosol extinction profiles from the Limb Profiler as well as other satellite data sets. We evaluated the limb viewing geometry effect (the “arch effect”) in the retrieval of the LP standard aerosol extinction product at 674 nm. It was shown that the amount of SO₂ decreases with a characteristic period of 8–18 d and the peak of stratospheric aerosol optical depth recorded at a wavelength of 674 nm lags the initial peak of SO₂ mass by 1.5 months. Using satellite observations and a trajectory model, we examined the dynamics of an unusual atmospheric feature that was observed, a stratospheric coherent circular cloud of SO₂ and aerosol from 18 July to 22 September 2019.

1 Introduction

An eruption of the Raikoke volcano (see Fig. 1; 0.55 km summit altitude; uninhabited island at 48.29^∘ N, 153.4^∘ E) occurred on 21 June 2019 at 18:00 UT. The eruption was so strong that a cloud of ash and volcanic gases was ejected to a height of 17–19 km (Gorkavyi et al., 2020), putting a significant part of the Raikoke volcanic plume into the stratosphere and above the heights of meteorological water and ice clouds and the tropopause, which is ∼10–12 km for such northern latitudes (see Fig. 1). The Microwave Limb Sounder (MLS) satellite data for 23–27 June 2019 show that the observed parts of the dispersed SO₂ cloud had heights between 11 and 18 km, with a peak concentration at 14 km. The maximum volcanic cloud top height increased by more than 6 km within 4 d after the eruption due to aerosol–radiation interaction (Muser et al., 2020).

Figure 1 Left: Raikoke volcano (Kuril Islands, Russia; 48^∘ N, 153^∘ E) (Google Maps, Imagery: NASA). Right: International Space Station (ISS) photo of the Raikoke eruption (21 June 2019, 22:45 UTC), which shows that the volcanic plume rose much higher than the layer of tropospheric clouds (photo: iss059e119250 22 June 2019, ISS/NASA, https://visibleearth.nasa.gov/images/145226/raikoke-erupts/145230w, last access: 29 November 2021).

The study of dispersed volcanic clouds and the transformation of their SO₂ into sulfuric acid aerosols is interesting for several reasons. Firstly, clouds of ash and sulfate aerosol are aviation hazards (Carn et al., 2009; Krotkov et al., 2014). Secondly, the lifetime of volcanic SO₂ depends on its injection height, SO₂ uptake on ash, and the concentrations of oxidants that vary with season and location (Guo et al., 2004; Krotkov et al., 2010; Beirle et al., 2014; Zhu et al., 2020). SO₂ from volcanic eruptions has a longer lifetime if it is injected into the stratosphere (Carn et al., 2016). SO₂ in the stratosphere is converted to sulfate aerosol, which interacts with the long-lived Junge layer – the naturally occurring background stratospheric aerosol layer. Furthermore, the dynamics of the Junge layer itself is not well understood. In addition, stratospheric aerosols are an important climate forcing factor, because aerosols modify both the shortwave and longwave radiation in the atmosphere and reaching the Earth's surface (Toohey et al., 2019; von Savigny et al., 2020). Thus, each case of volcanic injection of large amounts of SO₂ into the stratosphere is of great interest to climate scientists (Robock 2000; Foster et al., 2008). Our ability to study volcanic SO₂ has evolved in recent decades along with satellite remote sensing technology (Krueger 1983; Krueger et al., 2000, 2008; Bovensmann et al., 1999; Guo et al., 2004; Carn et al., 2003, 2007, 2008, 2009, 2017; Clerbaux et al., 2009; Clarisse et al., 2010, 2013; Penning de Vries et al., 2014; Sandvik et al., 2019; Hedelt et al., 2019; Theys et al., 2019; Fisher et al., 2019). The spectral data obtained from satellites have made it possible to analyze in detail transient volcanic SO₂ clouds, such as of the eruptions of Pinatubo (Bluth et al., 1992; Guo et al., 2004; Fisher et al., 2019), El Chichón (Krueger 1983; Krueger et al., 2008), and Kasatochi (Bourassa et al., 2010; Prata et al., 2010; Krotkov et al., 2010; Clarisse et al., 2011).

The 2019 eruption of Raikoke was well observed by a number of satellite sensors including the Ozone Monitoring Instrument (OMI), the Tropospheric Monitoring Instrument (TROPOMI), the Ozone Mapping and Profiler Suite (OMPS) Nadir Mapper (NM) and Limb Profiler (LP), and the Cloud–Aerosol Lidar with Orthogonal Polarization (CALIOP). Studying the nadir passive data of OMPS NM and TROPOMI, our attention was drawn to an interesting phenomenon, the presence of compact long-lived stratospheric SO₂ clouds, hereafter referred to as coherent circular clouds or CCC, each about 300 km in diameter. The first was noted on 29 June 2019 centered at 51^∘ N, 157^∘ W (near Alaska), and the second was observed starting on 18 July 2019 near Kamchatka. The first of these moved towards the North Pole while the second moved south to ∼20–30^∘ N, where it remained for more than 2 months, having made almost three complete revolutions around the Earth.

Figure 2 TROPOMI image of total column SO₂ (log scale), a week after the Raikoke eruption (image: courtesy Copernicus TROPOMI/Sentinel-5 Precursor SO₂ data, https://so2.gsfc.nasa.gov/, last access: 29 November 2021). The triangle represents the location of the Raikoke volcano. A coherent circular cloud is seen on 29 June 2019 centered at 51^∘ N, 157^∘ W. Note that this plot crosses the date line and so includes data from 29 (where the CCC is located) and 30 June (UTC).

A dense cloud of SO₂ and aerosol that formed after the eruption of the Raikoke volcano spread over a week at latitudes above 40^∘ N (see TROPOMI SO₂ map on 29 June 2019 in Fig. 2). Atmospheric currents stretched the cloud into long jets, twisted them with spirals, and even formed CCCs, as shown in Fig. 2 in the region near the coordinates 52^∘ N, 156^∘ W. To analyze in detail the evolution and dynamics of the Raikoke volcanic cloud, we examine (1) the amount and height of SO₂ and aerosol emitted during the eruption using satellite nadir and limb data; (2) the mutual evolution of SO₂ and aerosol for 100 d, starting from the moment of the eruption with limb data using an adjustment correcting for limb viewing geometry effects; and (3) the propagation of the CCCs with a case study in July 2019 using a trajectory model.

2 Data and methods

In this section, we detail the instruments and methods used to (1) examine SO₂ and aerosol distributions obtained from passive nadir sensors for SO₂ and limb and lidar for aerosol, (2) calculate zonal mean values of SO₂ and aerosol as a function of a height using nadir data for SO₂ and limb data for aerosol, and (3) employ trajectory modeling for analysis of individual clouds.

2.1 Satellite mapping of SO₂

We use SO₂ total column data obtained from two satellite spectrometers: (1) TROPOMI (see Fig. 2; Theys et al., 2017, 2019) on the Sentinel-5 Precursor satellite and (2) OMPS Nadir Mapper (NM) on the Suomi National Polar-orbiting Partnership (SNPP) satellite (Li et al., 2013, 2015; Zhang et al., 2017). Both instruments observe the Earth's backscattered radiance and solar irradiance at hyperspectral ultraviolet (UV) wavelengths with two-dimensional charge-coupled devices (CCDs), measuring in the spectral domain in one dimension and in the spatial domain (across a satellite track) in the other dimension. The satellite motion provides measurements along the satellite track. The spatial resolution of TROPOMI (3.5 km by 7 km) is much finer than that from OMPS NM measurements made in the nominal mode (50 km by 50 km), and the spectral resolution is higher for TROPOMI as well (∼0.5 nm as compared with ∼1 nm for OMPS NM). The satellite swaths are wide (2700 km for TROPOMI and 2800 km for OMPS NM), providing nearly global daily coverage. Column SO₂ in Dobson units ($\mathrm{1}\mathrm{DU}=\mathrm{2.69}\times {\mathrm{10}}^{\mathrm{16}}$ molecules cm⁻²) is retrieved from hyperspectral solar backscatter at UV wavelengths (312–390 nm). TROPOMI uses a spectral fitting algorithm based on differential optical absorption spectroscopy (DOAS) (Theys et al., 2017, 2019). OMPS uses a spectral fitting approach with a principal component analysis (PCA) scheme (Li et al., 2013, 2017). For large volcanic SO₂ signals like Raikoke, comparisons between TROPOMI and SNPP/OMPS so far (for several eruptions) show little bias, with the total SO₂ mass estimates from the two normally agreeing to within 5 %–10 % (with the exception of the very early stages of large eruptions, where the density of SO₂ and/or volcanic ash is too high to be fully accounted for in operational algorithms). For retrieval noise on a pixel-to-pixel basis, SNPP/OMPS SO₂ (for stratospheric clouds) is less than 0.1 DU. TROPOMI's noise on a pixel-by-pixel base is several times greater, but once TROPOMI pixels are averaged to OMPS footprints, the noise is reduced by ∼30 %. We also used the Ozone Monitoring Instrument (OMI) on NASA's Aura satellite as an additional data source.

2.2 Aerosol profile data

We use backscatter data from CALIOP (Fairlie et al., 2014), the OMPS Limb Profiler (LP) (Loughman et al., 2018; Chen et al., 2018), and the Stratospheric Aerosol and Gas Experiment instrument aboard the International Space Station (SAGE III/ISS). The OMPS LP V1.5 aerosol retrieval algorithm is described Sects. 2 and 3 of Chen et al. (2018). The height of an aerosol cloud can be estimated with both OMPS LP and CALIOP. We use the CALIOP lidar images for aerosol, specifically 532 nm total attenuated backscatter signal with a spatial resolution of 40 km along track and a vertical resolution of 120 m (altitude<20 km) and 360 m (altitude>20 km).

OMPS LP hyperspectral measurements capture aerosol data with a sampling of ∼1 km in altitude. The actual cross-track field of view (FOV) of each OMPS LP slit is ∼4 km (Dittman et al., 2002). OMPS LP views the atmosphere in a backward direction along the orbit track with three vertical slits, one (central) aligned with the orbit track and the other two (left and right) separated by $\pm \mathrm{4.25}{}^{\circ }$ horizontally. Here, we use aerosol data from version 1.5 (Chen et al., 2018).

Figure 3 The diagram for observations of a limb sensor, showing the relationship between the observed (h) and real (H) heights of one cloud in five positions, as well as φ – the angular displacement of the position of the cloud from the tangential point T. Ψ is a latitude. The figure is not to scale.

Limb measurements of the atmosphere that view an altitude H at the tangent point also view altitudes above H in the foreground and background of the line of sight. However, the converse situation is also true: an object at a fixed altitude H will appear to be at a lower altitude h if it is located closer or farther than the tangent point. This is shown schematically in Fig. 3, where a single cloud is located at true altitude H and five successive LP measurement events (each event generates a new profile) are illustrated. For event C, the cloud location is at the tangent point, and the apparent altitude h equals the true altitude H. For events A and B, the cloud position is closer to the LP instrument than the tangent point, and h is less than H. For events D and E, the cloud position is farther from the LP instrument than the tangent point, and again h<H.

Figure 4 Ratio of real and observed heights for a cloud 1 km thick and 226 km long, the center of which is located at different distances from the tangential point.

Plotting successive altitude profiles of LP (e.g., extinction coefficient) as function of tangential latitude when a vertically limited feature (such as a cloud) is present gives the impression of an “arch” in the data. The arch effect is observed when the length of the visible part of the cloud is less than ∼1100 km (at a cloud height of 25 km). Figure 4 shows cloud F₀G₀, 1 km thick and 226 km long, centered above tangential point T. Due to the curvature of the globe, such a cloud has an observed thickness of 2 km (see Fig. 4). If we take a cloud 226 km long and with a real thickness of N km, then the observed cloud thickness will be N+1 km. Thus, the real average height of a thin (1–2 km) cloud is underestimated by 0.5–1 km even under the most optimal observation conditions. Consider a cloud FG (see Fig. 4), the center of which is displaced from the tangential point by 273 km (or by ϕ∼2.5^∘). The real height of the FG cloud is 24–25 km, but its observed height varies from 13 to 22 km. If we consider the F₀G cloud with a length of 499 km, then its real height above the Earth's surface will be 24–25 km, and the observed height is 13–25 km. Let us take into account that the Limb Profiler assigns the latitude of the tangential point to any extended cloud. Therefore, a single cloud shown in Fig. 3 in five different observed positions, instead of one real geographic latitude, receives several observed latitudes, which creates an arch effect. Let the region F₀G₀ be a gap in a continuous cloud. Then this gap, together with the arch effect, will lead to a decrease in the maximum observed height of the cloud layer by 1 km (see Fig. 4).

Figure 5 (a) The extinction coefficient profiles (674 nm) from OMPS LP for the orbit 40 636 (central slit, same orbit as in Fig. 12, which shows the profiles for events 102, 105, and 108) for 31 August 2019. To accurately calculate the extinction coefficient it is necessary to take into account the arch effect, which is clearly visible (OMPS LP data). (b) (Left) An arch that appears during the limb observation of an aerosol cloud 2 km thick (H=23–25 km) and 2^∘ long. Each frame corresponds to a satellite orbital shift by 1^∘. Radiance from aerosol is proportional to the number of observed particles, taking into account the distance to them (the left part of the arch is slightly brighter than the right, because on frames 17–21 the cloud was closer to the satellite than on frames 23–27). The radiance units are arbitrary. (b) (Right) The red line is the profile of a cloud of 400 particles, which is observed on frame 22 (middle of the arch in Fig. 5a. The black line is the number of visible particles, summed over all frames, that is, over the entire arch.

Figure 5, showing an OMPS LP extinction profile track across the volcanic plume, illustrates the arch effect. The maximum values (yellow) appear to vary in altitude by 2–3 km due to the projection effect shown in Fig. 3. Note that this effect will be present for high, thin clouds such as cirrus clouds as well as for aerosol plumes. If we believe that these lower-altitude values do not represent a true aerosol signal, we need to apply a correction in order to accurately determine overall aerosol loading.

Table 1 The observed height h and angular displacement φ for the cloud with true height H.

The relationship between true cloud height H and apparent cloud height h in Fig. 3 is given by a simple function (DeLand and Gorkavyi, 2021):

$$\begin{array}{}\text{(1)}& \cos \mathit{\phi }={\displaystyle \frac{R+h}{R+H}},\end{array}$$

where R is the radius of the Earth and φ is the angular displacement. Table 1 gives some specific examples of the angular displacement predicted for specific combinations of H and h. This displacement is approximately equal to a shift in latitude for a sun-synchronous satellite such as SNPP ($\text{inclination}=\mathrm{98.9}{}^{\circ }$) outside the polar regions. From this table, we see that for a true cloud at H=25 km, the apparent cloud signal at h=15 km would be displaced by ∼3.2^∘ both earlier and later along the orbit. This total separation of ∼6.4^∘ is consistent with the arch width at 15 km shown in Fig. 4. Note also that because of differences in overall path length, we expect events A and B in Fig. 3 (left branch of arch in Fig. 4) to have a stronger signal than the corresponding events D and E at the same apparent altitude. We can therefore use Eq. (1) to calculate and apply a correction for determining the magnitude and position of an aerosol cloud (we believe that all parts of the arch below the real height are artifacts, so the value of the extinction coefficient for them should be equal to zero).

The arch in Fig. 5a is formed by several measurements of one CCC, which has a true height of 23–25 km (see the analysis of this CCC in Sect. 3.2). Parts of the arch (between events 102 and 108), indicated with arrows in Fig. 5a below 23 km, are artifacts of multiple registrations of the same object and should be removed when calculating the total amount of aerosol at given latitudes or heights. The low stratospheric aerosol layer at an altitude of 18–20 km decreasing with latitude to 15–16 km (blue color) is the background layer. It merges in frames 101–109 with an apparent volcanic cloud arch. Note that the CCC has been the southernmost part of the Raikoke plume since late July 2019 (see Sect. 3.3). The tropical cloud at latitude ∼8^∘ N extends up to 16 km altitude.

The reality of the discussed arch effect is confirmed by simple modeling. Arch model in Fig. 5b was obtained by direct modeling of a cloud of 400 particles located at the nodes of a uniform grid (20 particles are distributed at a length of 2^∘, and 20 rows of particles are uniformly distributed along the radius in the range of 23–25 km). The model does not use radiative transfer models. In Fig. 5b (left) only the distance between the particles and the satellite is taken into account, while Fig. 5b (right) shows the simply visible (for the limb sensor) number of particles at a given altitude h (in 1 km step). The red line corresponds to a one-time observation of a cloud located at the tangent point at an altitude of h=22–25 km (an increase in the apparent thickness of the cloud by 1 km is associated with the curvature of the Earth, which is why the cloud itself turns out to be curved – see Figs. 3 and 4). That is, it is the most realistic observation of the cloud at its optimal location. The black line shows the sum of the cloud particles observed at different times. As a result of this summation, the number of particles visible at a given observed height h (which differs from the real constant cloud height H=23–25 km) turns out to be overestimated. Therefore, when analyzing the picture in Fig. 5, we must remember that it is composed of frames received at different times, so the one cloud will be registered many times. The same effect of multiple registrations will be observed for clouds of any complexity and configuration, including a uniform aerosol layer, because any cloud can be divided into a large number of elementary pieces, similar to a simple compact cloud considered in Fig. 5b. Obviously, the arch effect for a spherically uniform aerosol layer should be fully compensated for by the 1D retrieval model. But the further the real system is from the spherical symmetry, the more difficult it will be to take into account the arch effect.

The considered model of the arch effect does not depend on the specific model of radiation transfer and the methods of retrieval of the spatial distribution of aerosol. Therefore, for specific limb sensors (OMPS-LP, SAGE III, OSIRIS), it is necessary to evaluate how accurately the available retrieval packages handle compact clouds. This is especially true for 1D retrieval methods, which assume spherical symmetry of the atmosphere and which are used to obtain aerosol extinction in the OMPS/LP. To estimate the possible retrieval uncertainty due to the arch effect, we apply a simple compensation method for one specific case of the Raikoke aerosol cloud (see next section). This compensation method assumes that the considered aerosol clouds form compact clusters or a highly heterogeneous aerosol layer and that the arch effect was not taken into account in retrieval. Thus, this example should be regarded as a maximal estimation for the arch effect. Where this assumption is not valid, our correction will be overestimated, as, for example, happened with the correction of the background aerosol value (see Fig. 6), which was observed before the Raikoke eruption.

Figure 6 The daily zonal mean (45–85^∘ N) SO₂ mass (assuming a cloud height of 13 km) and the average aerosol extinction coefficient at 675 nm (summed up over 13–18 km and divided by 6 (km) to get the average extinction coefficient). The arrow bar shows the maximum possible error in calculating the aerosol extinction due to the arch effect.

Figure 7 Time series of the logarithm of the estimated total SO₂ mass (in metric tons) within (a) the main plume (45–85^∘ N) and (b) the CCC.

The principle of our posterior algorithm for estimating the arch effect is as follows: the data in the uppermost pixel of the arch are considered true. For a given pixel, artifact extinctions are calculated for pixels in the side branches of the arch. These artifact extinctions are subtracted from the initial arch extinction. The procedure is repeated for the second-highest pixel in the already modified arch image. This continues until the side branches of the arch completely disappear. The top pixels with the corrected extinction are summed into the cloud with the corrected data. We consider this procedure of correction only as an estimate, which shows the possible significance of the arch effect.

The arch effect is a specific example of the effects of inhomogeneity along the line of sight that is an issue for all types of limb sounders. One way to account for such effects is to use a two-dimensional (2D) radiative transfer model (RTM) that is able to account for such effects along with multiple observations in a tomographic retrieval (e.g., Livesey et al., 2006; Zawada et al., 2018; Loughman et al., 2018). Instead, we have developed an a posteriori adjustment method that is effective in correcting for the arch-type effects that appear with isolated features within the lines of sight for a series of observations.

The arch effect characteristic of the limb observations should be taken into account when calculating the optical thickness of aerosol clouds, determining a characteristic cloud height, or using calibration heights of ∼45 km. Calibration heights should be free of aerosol, and they are indeed free of aerosol but can be contaminated with polar stratospheric and mesospheric clouds (PMCs) and other stratospheric aerosols. For example, PMCs located at an altitude of 80–85 km, due to the arch effect, can be projected onto the calibration heights. Such contamination of the calibration heights distorts the whole picture of aerosol distribution with height.

SAGE III/ISS makes direct measurements of aerosol extinction through the attenuation of the solar beam during sunrise and sunset. The vertical resolution of these measurements is about $\mathrm{1}/\mathrm{2}$ km. Because they are occultation measurements, SAGE observations are less frequent than OMPS-LP.

Using the color ratio (extinction at 512 nm to extinction at 1022 nm) along with the extinction at 1022 nm, clouds can be removed from the SAGE observations to produce average aerosol extinctions (Thomason and Vernier, 2013).

Figure 8 Latitudinal distribution of the aerosol extinction coefficient (×10⁴, km⁻¹) at 675 nm at three altitude levels – (a) 14.5, (b) 18.5, and (c) 20.5 km – averaged every 5^∘ of latitude. The Raikoke eruption is represented by “o”, Alberta fires are represented by “*”, Siberian fires are represented by “x”, and Ulawun eruptions are represented by “#”.

Figure 9 (a) Daily location of the maximum altitude reached by the Raikoke aerosol plume as detected by OMPS LP, colored by the day number since 1 June and plotted every other day until day 150. (b) Daily zonal mean aerosol extinction profiles at 675 nm (km⁻¹) between 20^∘ S–90^∘ N (which is close to hemispheric coverage) measured by OMPS LP from June 2019 to June 2020 and smoothed spatially using a five-point boxcar averaging. Only profiles measured above tropopause +1 km are used. The white line is the aerosol plume maximum altitude in kilometers (km). (b) (Top) Stratospheric aerosol optical depth (×10³, sAOD) at 675 nm, derived by integrating aerosol extinction profiles above the tropopause to 35 km.

2.3 The NASA ftraj trajectory model

The “ftraj” trajectory model from NASA's Goddard Space Flight Center Atmospheric Chemistry and Dynamics Laboratory uses a fourth-order Runge–Kutta integration scheme to track parcels isentropically, with optional diabatic adjustments (Schoeberl and Sparling, 1995). The model is driven with winds at $\mathrm{1}/\mathrm{4}$^∘ horizontal resolution and spaced every 6 h from the Goddard Earth Observing System (GEOS) forward-processing system produced by the NASA Global Modeling and Assimilation Office (GMAO). Each model run was initialized with 3000 parcels distributed at random horizontal and vertical positions within a cylinder of radius 150 km centered about the location of the observed cloud and stretching between the lowest and highest potential temperature (Θ) values of the observed cloud. Those Θ values were derived by interpolating the GEOS fields to the parcels' initial longitude, latitude, altitude, and time. Both the starting and ending cloud observations extended through a range of altitudes and hence Θ values.

Figure 10 The ratio of the zonal mean, daily average aerosol extinction for the indicated dates after the Raikoke eruption to the averaged aerosol extinction immediately prior to the eruption shows the evolution of the volcanic aerosol in the background aerosol layer (13–20 km) and the occurrence of the CCC near 25 km from OMPS LP data.

Figure 11 (a) SAGE aerosol extinction at 1022 nm showing enhanced extinction spreading from high northern latitudes toward the tropics. The elevated part of the plume is visible at an altitude of 23–25 km. (b) SAGE profile locations during September 2019. Red circle shows where SAGE crossed the elevated portion of the Raikoke plume.

3 Results

The sensitivity of the satellite data (OMPS LP, CALIOP) we examined is such that the aerosol cloud from the Raikoke volcano was observable for many months following the eruption. Thus, a synergistic study of Raikoke volcanic emissions using various satellite-based instruments provides a good opportunity to study the dynamics of the volcanic cloud long-term dispersion. For quantitative analysis, we will use 2D SO₂ vertical column density (VCD) maps obtained by OMPS NM, as well as vertical profiles of aerosol extinction that are obtained from the OMPS LP.

3.1 SO₂ and aerosol evolution

The retrieved Raikoke SO₂ mass increased for 2 d after the eruption and then exponentially decayed with an e-folding timescale of 8–18 d. The amount of stratospheric SO₂ decreases due to photochemical conversion to sulfuric acid through gas-phase reaction with the hydroxyl radical, OH. Sulfuric acid nucleates new particles and condenses onto pre-existing particles to form long-lived stratospheric sulfate aerosol. Figure 6 shows the change in the total SO₂ mass (kt) retrieved by the OMPS NM and the altitude-averaged (13–18 km) aerosol extinction coefficient retrieved with the OMPS LP in the latitude range 45–65^∘ N after the Raikoke eruption. The aerosol extinction data in Fig. 6 are shown before and after removing the arch effect. These estimates show the potential significance of the arch effect. As we discussed earlier, the arch effect estimation assumes high cloud heterogeneity and should be considered the maximum estimation, as indicated by the arrow bars on Fig. 6. For more accurate calculations, the arch effect should be investigated using a 2D RTM.

Figure 6 shows that retrieved SO₂ mass increases for 2 d after the eruption during rapid dispersion of the fresh opaque volcanic cloud. Following this initial dispersion, the SO₂ mass decreases exponentially due to chemical conversion to sulfate. Aerosol extinction has opposite behavior; it increases by approximately a factor of 4 as the concentration of SO₂ decreases, reaching a maximum ∼50 d after the eruption. After that, the aerosol extinction starts to decrease due to gravitational sedimentation and other processes, but very slowly (from OMPS LP and NM data). Adjustments for the arch effect do not change the overall temporal pattern.

As can be seen from Fig. 7, the total SO₂ mass decreases exponentially with a variable (8–18 d) timescale. Krotkov et al. (2010) derived a volcanic SO₂ mass decay of ∼9 d following the 8 August 2008 Kasatochi eruption (30–90^∘ N) using the Aura OMI. The overall apparent e-folding time for Raikoke is ∼8–10 d (Fig. 7, solid black lines), in good agreement with the Kasatochi e-folding time. However, we see two different regimes in Raikoke SO₂ clouds with significantly different estimated e-folding time, as discussed below.

Figure 7a shows that during the first 20 d after the eruption, the SO₂ mass decreases with a longer apparent e-folding lifetime of ∼19 d, compared to the later period (days 21 to 60), when the decrease in SO₂ mass accelerates with the apparent e-folding lifetime ∼8 d. The initial estimate of the e-folding lifetime of ∼19 d, derived from the earlier stage after the eruption, probably represents the true timescale for chemical conversion of SO₂ into sulfate aerosol on altitudes>20 km (Carn et al., 2016). Following the drop of the SO₂ column amounts, the SO₂ detection limit of the OMPS instrument becomes a more important factor, leading to faster apparent decay of the total observed mass of SO₂. As the plume spreads out to larger areas, and more and more pixels with SO₂ fall below the detection limit of the OMPS NM sensor (OMPS in the stratosphere can typically detect 0.2–0.3 DU of SO₂), the apparent decay rate becomes larger, likely reflecting the combined effects of the chemical SO₂ loss and the diminished OMPS sensitivity to weaker SO₂ signals. Figure 7b shows a similar pattern for change in SO₂ mass within the CCC. For the period between 28 and 45 d after the eruption, the e-folding lifetime is ∼19 d. This suggests that the time after the eruption does not play a big role here – only the compactness of the CCC with a high level of SO₂ is essential. As in the case of Fig. 7a, the second phase of the evolution of the CCC is characterized by a faster apparent (the e-folding lifetime ∼8 d) rate of decrease in the recorded SO₂.

3.2 Stratospheric aerosol from limb observations

Figure 8 shows the time–latitude evolution of the zonal average aerosol extinction at different heights (without adjustment for the arch effect). Most of the aerosol is transported poleward at 14.5 km (Fig. 8a), and the effects of the eruption lasted for almost a year. In addition, Figs. 8b and c also show the aerosol transport to subtropics and tropics at higher altitudes. Increased aerosol loading in the lower stratosphere can also be attributed to two pyrocumulonimbus (pyroCb) events that took place before and after the eruption, Alberta fires (18 June) and Siberian fires (2 July) (Kloss et al., 2021). OMPS LP detected both plumes in the stratosphere at 12–13 km, although it became difficult to separate them from the Raikoke plume once it spread around the Northern Hemisphere. The zonal mean aerosol extinction profile between 20 and 90^∘ N shows the vertical transport of the Raikoke plume to higher altitudes and its persistence in the lower stratosphere (Fig. 9). Figure 9a shows that the maximum altitudes of the plume are around 25 km, when the plume penetrates the tropics. The plume altitude is derived using the OMPS LP cloud algorithm, which can identify enhanced aerosol layers in the stratosphere.

Figure 10 shows the altitude–latitude profile of a daily zonal mean ratio of the aerosol extinction coefficient for 1 and 21 July and 1 and 30 September 2019 to a similar daily zonal mean for a quiet period (16–20 June 2019) before the eruption. The means were computed without taking into account the arch effect (see DeLand and Gorkavyi, 2021, and Sect. 2.2). Therefore, the lower parts of the layers shown may have an overestimated extinction coefficient. In addition to increasing the density of the background aerosol layer (the contribution of seasonal changes is possible here), an interesting CCC appears at an altitude of about 25 km that will be discussed in Sect. 3.3.

The SAGE 1022 nm extinctions for the month of September 2019 are shown in Fig. 11 along with the profile latitudes. Figure 11a shows that enhanced aerosol extinction extended from the high-latitude tropopause to the tropics during September, revealing the extent of the aerosol plume. The SAGE aerosol distribution is similar to the OMPS-LP distribution shown in Fig. 10. The SAGE III profile latitudes (Fig. 11b) show that SAGE had reasonable coverage from 60^∘ S to 60^∘ N during the month. SAGE III observed the high-altitude portion of the Raikoke plume (red circle) at 30^∘ N as the occultations were moving southward. The latitudes for this observation are indicated with the red circle in Fig. 11b. The location of this detached plume is indicated in Fig. 11a, which agrees well with the OMPS-LP observations shown in Figs. 5 and 10.

Figure 12 The aerosol formed as a result of the Raikoke eruption is clearly visible in the extinction coefficient profiles obtained by the OMPS LP aerosol extinction vertical profile months after the eruption.

Figure 12 shows three LP aerosol extinction profiles that were measured by OMPS LP on 31 August 2019, more than 2 months after the Raikoke eruption. The extinction peak at 22–25 km represents the remaining signal of the stratospheric plume from the eruption. The apparent altitude variation between these peaks also shows the impact of the LP viewing geometry on the retrieved extinction profile. Each measured profile (also termed “event”) gives an along-track separation of ∼125 km. The along-track field of view integrates over a distance of ∼180 km for each 1 km vertical sample (DeLand and Gorkavyi, 2021). The OMPS LP version 1.5 algorithm has a restriction on how much the retrieved aerosol extinction is allowed to grow per iteration at each altitude relative to the first guess, which may cause an underestimation of the retrieved aerosol where the plume is concentrated (README Document for the Suomi-NPP OMPS LP L2 AER675 Daily Product, 2019).

3.3 Coherent circular cloud (CCC)

The CCC is part of the original aerosol/SO₂ cloud resulting from the eruption of the Raikoke volcano. It formed 1 month after the eruption in the same area where the volcano is located (Fig. 13). Perhaps, the formation of the CCC is associated with a tropospheric vortex, which was observed at the same time in the area (Fig. 13). This vortex could have affected the temperature or pressure in the stratosphere, which caused the CCC. On 24 July 2019, as the CALIOP backscatter data show, the CCC had a height of 19–20 km and quickly shifted to the south (Fig. 13). On 30 July 2019, it reached a latitude of 30^∘ N and was extended across China, fell into easterly winds, and moved west at about 17 m s⁻¹ (Fig. 14). On 4 August, it was recorded over the Persian Gulf region (Fig. 15), on 5 August over Egypt, and on 8 August near the Azores.

Figure 13 (a) TROPOMI SO₂ STL (lower stratosphere) columns (upper tropospheric and stratospheric SO₂ column with the center of mass altitude 17 km) for 20–21 July (UTC) plotted using a log scale. Only the Raikoke volcano is marked (triangle); (b) OMPS SO₂ STL columns plotted using a linear scale for 24 July 2019. The Raikoke volcano is marked (triangle). Figures show the formation of CCC, which appeared near Kamchatka and shifted to latitude ∼30^∘ N and moved to the west at a speed of 1400 km d⁻¹. From OMPS NM data we estimated the initial mass of SO₂ in this cloud to be ∼20 kt; (c) CALIOP total attenuated backscatter data (${\mathrm{km}}^{-\mathrm{1}}{\mathrm{sr}}^{-\mathrm{1}}$) reveal the height of the CCC to be 18.5–20 km.

Figure 14 OMPS SO₂ STL columns plotted using a linear scale for 30 July 2019. The data show the westward drift of the CCC with easterly stratospheric winds at ∼14–15^∘ d⁻¹, or ∼1400 km d⁻¹, ∼60 km h⁻¹.

Figure 15 OMPS SO₂ STL columns plotted using a linear scale for 4 August 2019.

Figure 16 The CALIOP total attenuated backscatter data (${\mathrm{km}}^{-\mathrm{1}}{\mathrm{sr}}^{-\mathrm{1}}$) for 18 August and 1 September 2019, where the CCC is visible at the indicated latitudes and heights.

Using the data of CALIOP, OMPS LP and NM, and OMI, the CCC was tracked until the end of August, while making 1.5 orbits around the Earth. CALIOP and OMPS LP (see Fig. 16 with CALIOP data from 18 August and 1 September 2019 and Figs. 5 and 12 from 31 August 2019 – OMPS LP data) were able to track the movement of the accompanying aerosol part of this CCC up to 22 September as it made ∼3 orbits around the Earth.

Figure 17 The observed path of the CCC as inferred from OMI/OMPS/TROPOMI and CALIOP observations. The individual observations are shown by dots color-coded by month.

Figure 17 shows satellite observations of the CCC from 18 July to 22 September 2019. During this period, the height of the CCC has increased significantly. There are more points by the end of the CCC observation period, because the cloud spreads out and is more often observed by CALIOP. Figure 17 was published in December 2019 at AGU-2019 (Gorkavyi et al., 2019). On 24 September 2019 the CCC was observed by lidars in Hawaii (Chouza et al., 2020). Chouza et al. (2020) traced the trajectory of this cloud back to 17 July 2019. Although the two studies were done independently, they came up with very similar results. Chouza et al. (2020) consider this cloud as a Raikoke plume, but we prefer to call it CCC because it is a very small part of the Raikoke plume.

Figure 18 Backward trajectory modeling. Panel (a) shows several discrete time snapshots from the Goddard “ftraj” back-trajectory modeling for 3000 parcels within a 150 km radius of the observed CCC position initialized on 18 August 2019. (b) The trajectory model run was initialized with the observed cloud position on 19 September 2019 and then run backwards to 18 August 2019. Panel (c) shows the same model results as for panel (a) but using a different visualization. Here the parcels are binned and counted in longitude and time, regardless of altitude or latitude. (d) Similar to panel (c), except that the model run is initialized on 19 September 2019 and run backwards to 18 August 2019.

Figure 18 shows the results of backward trajectory modeling of the CCC for the 1-month period with August 18 through 19 September 2019. The basic idea is to use the trajectories to see if this long-lived CCC can be correctly transported using assimilated meteorological observations.

We have two hypotheses to evaluate. The first hypothesis (“the great red spot” hypothesis) is that the observed CCC is contained within a gyre (i.e., a whirlpool) that existed from around the time of the eruption on 21 June to 1 September or even later, advected around the globe by the large-scale circulation. The second hypothesis (the “dead fish” hypothesis) is that the CCC was simply advected around by the large-scale winds without being sheared apart because the environment of the summer stratosphere is fairly quiescent.

We inspected the meteorological fields from NASA's Global Modeling and Assimilation Office (GMAO) Modern-Era Retrospective analysis for Research and Applications, version 2 (MERRA-2) reanalysis (Gelaro et al., 2017) and from the Global Data Assimilation System (GDAS) by the National Centers for Environmental Prediction (NCEP) at the location and altitude of the observed volcanic CCC on 21 and 24 July; the winds show no discernible gyre. This is to be expected, as the satellite data that go into these assimilation products are unlikely to pick up features on sufficiently small scales that they would show up clearly and prominently in the winds. Therefore, it is not possible to test the great red spot hypothesis using a trajectory model with the wind products available to us.

We can examine the dead fish hypothesis. If a trajectory model run shows a cloud of parcels maintaining its integrity over a long period of time, then this hypothesis becomes plausible. There are difficulties, however. First, the winds used to drive the trajectory model have spatial and temporal resolutions that are coarse enough that trajectories cannot be modeled with perfect accuracy. One can expect that a cloud of parcels will lose many of its parcels along the way. It would be unrealistic to expect to duplicate the volcanic cloud's integrity and position as a function of time. Instead, the test is whether the simulated cloud of parcels can maintain its existence over the course of a month or more during this time period. Second, a more serious problem is that of altitude. We can expect a volcanic CCC to experience self-lofting, and the CCC indeed is observed to increase in altitude with time. However, the winds vary with altitude, and thus the parcel trajectories will depend critically on getting their altitude changes right.

A carbonaceous aerosol plume associated with wild fires in British Columbia in August 2017 reached the stratosphere a few days following the initial injection and resulting from self-lofting triggered by solar heating (Torres et al., 2020; see similar effect for Australian fires – Khaykin et al., 2020). The case of a volcanic CCC, which we are considering, is interesting in that the effect of an increase in height from 19 to 26 km is observed for a cloud presumably consisting of sulfate aerosol. In this case, the main increase in the height of the CCC occurred during its movement at approximately the same latitude. Early on (19 July) it is at 80 hPa (450 K), but by 1 September it had risen to 30 hPa (590 K). This is a diabatic self-lofting heating rate of about 3 K d⁻¹. Note that this CCC, increasing its height, simultaneously contained a detectable amount of SO₂, which is also confirmed by independent satellite data from the Atmospheric Chemistry Experiment–Fourier transform spectrometer (ACE-FTS) (Cameron et al., 2021).

The thermodynamic processes within the cloud are not well known, so rather than attempting to model the cloud's self-lofting and the overall vertical advection, we impose an approximation to the cloud's vertical motion. We determined an average daily Θ adjustment that would map the column of parcels from its initial Θ range to its final range. Back-trajectories were then computed isentropically for one day at a time. After each day, the parcels' potential temperature values were changed by their daily adjustment. In this way, the vertical motion of the CCC is guided by observations, while the horizontal motions are unconstrained. Note that this adjustment procedure implies two assumptions: that the self-lofting is linear in time and that it is linear in Θ as well. The map (see Fig. 18a) displays the parcel population count as a heat map or density display: the more intense the color at a given point, the more parcels are near that point (the resolution of the map is $\mathrm{1}{}^{\circ }\times \mathrm{1}{}^{\circ }$). The plot below the map shows a similar heat map of parcel population count as a function of longitude (x axis) and altitude (y axis); the resolution here is 1^∘ longitude by 0.5 km altitude.

Figure 18a superimposes several discrete time snapshots from the trajectory model output, starting with the observed cloud position on 18 August 2019. The modeled positions are generally close to the observed, except that the model fails to pick up the effects of the synoptic weather feature beginning around 23 July. Figure 18b shows the output from another model run, initialized with the observed CCC position on 19 September 2019 and then run backwards to 18 August 2019. Figure 18c and d show the same model results but using a Hovmöller-type plot.

The CCC appears not to be a whirlpool but rather a CCC with easterly flow and drifting (dead fish model). Neither the turbulent mixing nor gravity waves penetrate well into the summer stratosphere to mix features, so they tend to move intact if injected at the right time.

4 Conclusion

We studied the process of converting SO₂ to aerosol for the 21 June 2019 eruption of the Raikoke volcano (Kuril Islands, Russia; 48^∘ N, 153^∘ E). The peak sulfate aerosol extinction is 1.5 months behind the time of the release of SO₂ in the Raikoke eruption. We see two different regimes in Raikoke SO₂ clouds with significantly different e-folding times. Initially, the SO₂ mass decreases with a longer apparent e-folding lifetime of ∼19 d, compared to the later period when the decrease in SO₂ mass accelerates with the apparent e-folding lifetime ∼8 d. The e-folding lifetime of ∼19 d, derived from the earlier stage after the eruption, probably reflects the true conversion of SO₂ into sulfate aerosol. After that, the detection limit of the OMPS instrument becomes a more important factor, leading to faster apparent decay of the total observed SO₂.

We evaluated the influence of the arch effect on the calculation of the extinction coefficient of a finite volcanic cloud using a homogeneous spherical shell atmosphere assumption adopted for processing the data of the OMPS Limb Profiler. We have shown that this effect is significant; therefore, it should be investigated by more accurately accounting for inhomogeneous atmospheric composition along the line of sight within the framework of a 2D model of radiative transfer.

We examined the unusual coherent circular cloud (CCC) of SO₂ and aerosol, which was observed for more than 2 months (from 18 July to 22 September). The CCC was embedded in the summer easterly flow in the stratosphere and demonstrated an unprecedented structural cohesion along with a diabatic self-lofting heating rate of about 3 K d⁻¹. For 2 months, the CCC circled the globe almost three times at the latitude of 30^∘ N (during September, the CCC shifted south to 15–20^∘ N) and increased its height from 19 to 26 km.