Stratospheric water vapor from the Hunga Tonga - Hunga Ha’apai volcanic eruption deduced from COSMIC-2 radio occultation

The eruption of the Hunga Tonga – Hunga Ha’apai (HTHH) volcano on January 15, 2022 injected large amounts of water vapor directly into the stratosphere. While normal background levels of stratospheric water vapor are not detectable in radio occultation (RO) measurements, eﬀects of the HTHH eruption are clearly observed as anomalous refractivity proﬁles from COSMIC-2, suggesting the possibility of detecting the HTHH water vapor signal. To separate temperature vs. water vapor eﬀects on refractivity, we use co-located temperature observations from the Microwave Limb Sounder (MLS) to constrain a simpliﬁed water vapor retrieval. Our results show enhancements of water vapor up to ˜1500-2300 ppmv in the stratosphere (˜29-33 km) in the days following the HTHH eruption, with propagating patterns that follow the dispersing volcanic plume. The stratospheric


Introduction
The January 2020 explosive volcanic eruption of Hunga Tonga -Hunga Ha'apai (~20 o S, 175 o W; hereafter HTHH) injected large amounts of water vapor (H2O) directly into the stratosphere.Measurements from operational radiosonde balloons over Australia in the first several days after the eruption show local H2O mixing ratios above 1000 ppmv over altitudes ~25-30 km (the top of the balloon measurements; Sellitto et al, 2022;Vömel et al, 2022), compared to typical stratospheric background H2O values of 5 ppmv.The isolated balloon measurements showed that the H2O enhancements occurred in relatively thin ~1-2 km layers.Measurements from the Microwave Limb Sounder (MLS) instrument on the NASA Aura satellite show H2O perturbations from HTHH that are unprecedented in the satellite data record, in terms of both altitude and magnitude (Millán et al, 2022).The MLS data show enhanced H2O at altitudes up to 53 km immediately after the eruption (January 15), and maxima between ~25-35 km during the following few days .
Maximum H2O values in these MLS retrievals, which represent averages over ~3 km, was ~200-400 ppmv.However, the standard MLS retrievals are not designed for the anomalously high stratospheric H2O values from HTHH, and many of these early retrievals did not pass the MLS quality screening criteria (Millán et al, 2022).The HTHH plume traveled westward and dispersed in the stratosphere, and MLS data show that local H2O maxima decreased to ~50 ppmv by early February, and to ~10-20 ppmv by late March.The MLS retrievals are much better characterized for the lower H2O values after late January.Anomalous high H2O from HTHH persists in the stratosphere through boreal summer 2022 and has spread over much of the globe (Xu et al., 2022;Legras et al., 2022;Khaykin et al., 2022).
The objective of this paper is to explore the retrieval of the extreme stratospheric H2O amounts from HTHH in the first week after the eruption using COSMIC-2 (C2) GNSS radio occultation (RO) data.RO measures the bending angle of radio waves propagating through the atmosphere, which is closely related to atmospheric refractivity (N), which in turn is dependent on temperature and moisture.Under normal conditions, H2O makes virtually no contribution to N in the stratosphere and thus is not detectable by RO measurements.
However, effects from HTHH are clearly observed as anomalous stratospheric N profiles from C2 that follow the HTHH plume during the first week after the eruption (shown below).These anomalous N profiles can potentially be associated with both temperature and H2O effects from the HTHH eruption, and we separate these influences by using independent temperatures from nearby MLS measurements.Our results include comparisons with the limited radiosonde H2O measurements over Australia (up to ~30 km) to help evaluate the H2O profiles derived from C2 measurements.

a) C2 GNSS-RO observations
By measuring the phase delay of radio waves from GNSS satellites as they slow and bend in Earth's atmosphere, profiles of bending angles and refractivity can be obtained (e.g.Anthes et al, 2008).At microwave frequencies in the troposphere and stratosphere, the refractivity varies due to contributions from the dry air and water vapor.Specifically, the atmospheric refractivity N can be related to atmospheric pressure p, temperature T, and water vapor partial pressure e (Smith and Weintraub 1953): In this study, we obtain level-2 C2 GNSS-RO profiles for January 2022 processed by the COSMIC Data Analysis and Archive Center (CDAAC) at the University Corporation for Atmospheric Research (UCAR).C2 is providing ~6000 profiles/day between ~40°N-40°S.
We use the 'atmPrf' product, which provides refractivity and dry temperature (retrieved under the assumption of dry air) from near the surface up to ~60 km; data above 40 km are strongly influenced by climatology, and we focus on altitudes 25-35 km in this study.The profiles are quality controlled at CDACC by assigning 'bad' flags based on several metrics, including deviation from the climatology by a specific threshold.We note that a few C2 measurements intersect the early HTHH plume on January 15 and indicate large N anomalies, but these are flagged as 'bad' retrievals.These profiles are briefly discussed below and shown in the Supplement.The effective vertical resolution of RO soundings varies from ~200 m in the lower troposphere to ~1-2 km in the upper stratosphere (more details below) while the horizontal footprint (horizontal scale represented by a single observation) is ~200 km (e.g.Anthes et al, 2008).
At CDAAC, the neutral atmosphere profiles are retrieved in the geometric optics approximation above 20 km.In this approximation, the resolution is physically limited by the 1st Fresnel zone  [Kursinski et al., 1997] which depends on altitude and is about 1.4 km at 30 km.At CDAAC, data smoothing is performed using the Savitzky-Golay filter; details can be found in [Zeng et al., 2019].The half-width of the filter response function used for the standard processing is 1.1 ~ 1.5 km at 30 km.Specifically for this study, we investigated the possibility of modifying the filtering in order to better resolve those  structures in the stratosphere induced by HTHH eruption.By modeling and testing we found that reducing the half-width of the response function by one-half, to 0.75 km at 30 km, results in improvement of the vertical resolution at 30 km without substantial increase of the effect of observational noise (mainly from the ionosphere), while further reduction increases the noise and creates artifacts in the profiles.Thus, C2 RO data used in this study were processed with twice increased vertical resolution compared to the standard retrieval.
b) MLS temperatures MLS temperatures are retrieved using measured limb emission from atmospheric oxygen (O2), and these data are not strongly influenced by the HTHH eruption.We use MLS temperature profiles based on the v4.2 retrieval (Livesey et al, 2022).MLS provides global sampling with approximately 3500 measurements per day.Temperature retrievals are provided on standard pressure levels; data extend from 261 hPa (~9 km) to above 1 hPa (~48 km, with an effective vertical resolution of 3-4 km over the main region of interest here (~30 km).

c) Radiosondes
We include comparisons with radiosonde measurements of stratospheric H2O during the first few days after the HTHH eruption.These data are discussed in detail in Vömel et al (2022), based on measurements from the operational upper air network using the Vaisala RS41 radiosonde.These data can detect the large stratospheric H2O from HTHH, and provide measurements up to ~ 30 km (depending on the balloon burst altitude of individual flights) with vertical resolution of ~100 m.

a) RO sensitivity to stratospheric H2O
N decreases strongly with altitude in the troposphere and stratosphere, and the H2O term in (1) contributes little to N above the middle troposphere, where H2O < 100 ppmv (e.g.Johnston et al, 2022).Sensitivity tests based on (1) can determine the magnitude of GNSS RO refractivity variations expected from isolated H2O and temperature anomalies in the stratosphere.The sensitivity of N to an isolated large stratospheric H2O perturbation (1000 ppmv anomaly at 30 km) is illustrated in Fig. 1, showing a positive N anomaly of ~3.5% compared to an unperturbed background.While this is a relatively small fractional N anomaly, it is not small compared to the observed N standard deviation in the stratosphere measured by C2, with typical values of 1% background levels (as shown in Fig. 1; these were calculated from the pre-volcanic C2 measurements in the vicinity of HTHH).Hence, an N anomaly due to a 1000 ppmv increase in H2O is comparable to a (rare) positive 3-sigma C2 noise event.Figure 1 also shows the N anomaly associated with a local -8 K temperature perturbation at 30 km, which is similar (~3.5%) to the effect of 1000 ppmv H2O.These calculations suggest that C2 RO measurements may be able to detect a large stratospheric H2O perturbation on the order of 1000 ppmv, but that similar N signals can arise from large negative temperature anomalies.

b) C2 HTHH refractivity observations and H2O retrievals
Based on the sensitivity calculations shown in Fig. 1, our analyses search for positive N anomalies in C2 data that are larger than 3-sigma.Anomalies are calculated as % differences from the non-volcanic background January average.Figure 2  Examples of the vertical profiles of N anomalies for the cases identified on January 16 (Fig. 2a) are shown in Fig. 3a.As selected, each of the profiles has a maximum N anomaly above 3%, with several local maxima above 6%, and the observed profiles peak over altitudes ~29-33 km.The vertical thickness of the N anomaly profiles is ~2 km.
We calculate a simplified estimate of H2O based on Eq. 1, incorporating C2 N observations and MLS temperatures as independent information, and using dry pressure p from the C2 atmPrf files.Specifically, for each profile we take the full C2 N profile (not anomalies) and find the nearest co-located MLS temperature profile (within 600 km and 6 hours).We then solve Eq. 1 for water vapor pressure e and convert to H2O mixing ratio.
Figure 3b shows the derived H2O profiles on Jan. 16, with maximum H2O values of ~1500-2300 ppmv over altitudes ~29-33 km; the altitude maxima for H2O approximately match the N anomaly maxima in Fig. 3a.We note that these large H2O amounts are below ice saturation with respect to the background stratospheric temperatures that increase with height (e.g.Sellitto et al., 2022;Vömel et al., 2022).Figures 3c-d show derived H2O profile results for the following days, showing a decrease in the maximum H2O amounts to ~1500-2000 ppmv on January 17 and ~1000 ppmv on January 18.The measured altitudes of the H2O maxima also systematically decrease between January 16-18.
Analyses of the MLS temperatures show that there are not systematic temperature anomalies associated with the identified C2 N anomalies for the volcanic plume events (Fig. 4), so that the stratospheric RO N anomalies are arising primarily from enhanced H2O.This statement is qualified by the lower vertical resolution of MLS temperatures (~4 km) compared to C2 retrievals (~1-2 km), in addition to imprecise co-locations of the MLS and C2 data.In contrast to the lack of temperature anomalies in the MLS data, temperatures derived from the C2 dry retrieval (atmPrf, neglecting the H2O term in Eq. 1) show negative temperature anomalies up to -20 K (Fig. 4).These are likely spurious results related to the large positive N anomalies when not accounting for H2O effects.As a note, we have also examined co-located temperatures from ERA5 reanalysis (Hersbach et al., 2020) and find small systematic temperature anomalies for the co-locations, similar to the MLS results in Fig. 4.
RO H2O sensitivity inferred from Fig. 1 together with observed C2 noise levels, we expect that C2 will lose sensitivity to stratospheric H2O amounts below ~700 ppmv (which corresponds to positive N anomalies ~2-sigma).
We have also investigated several C2 measurements that intersected the HTHH plume early after the eruption on January 15, and the associated N anomalies and retrieved H2O amounts are shown in Fig. S1.The results show large positive N anomalies (~10-25%) and extreme H2O values (~ 3,000 -6,000 ppmv) over altitudes 30-40 km, and smaller maxima over 45-50 km.However, these retrieved N profiles are mostly classified as 'bad', and the reason is not simply understood because there are numerous potential causes.While we are less confident in these retrievals, the vertical structure is reasonable and the large H2O values in Fig. S1are still below ice saturation and so may represent actual geophysical structure (as noted in Khaykin et al., 2022).

Summary and Discussion
C2 data demonstrate that RO measurements are sensitive to the extremely large stratospheric H2O perturbations during the first week following the HTHH eruption.This is tied to the large observed refractivity anomalies associated with the HTHH plume, combined with the low stratospheric background noise in the C2 measurements (with standard deviations of ~1% background values).Sensitivity calculations (Fig. 1) suggest that stratospheric H2O anomalies of ~1000 ppmv or larger should be detectable in C2 refractivity measurements (equivalent to 3-sigma noise levels), with magnitude similar to a localized -8 K temperature anomaly.
To separate H2O and temperature effects on N we have used independent temperature observations from the MLS satellite.This has some limitations in terms of different vertical sensitivities and imperfect co-locations between C2 and MLS measurements, but is a necessary step to separate H2O and temperature influences on N. Comparison of MLS temperature anomalies vs. derived H2O amount (Fig. 4) suggests small systematic temperature anomalies for the H2O plume during the first week after the eruption, so that the C2 N anomalies are primarily due to H2O influence.As noted above, we find similar results using ERA5 temperatures.In contrast, C2 dry temperature retrievals (atmPrf) for the anomalous HTHH profiles analyzed here show large negative temperature anomalies (up to - shows the C2 daily sampling density in the region of the HTHH eruption on January 16, 17 and 18 (grey dots) along with the locations of local N anomalies with maximum values above 3-sigma over altitudes 25-35 km (colored symbols).The stratospheric plume from HTHH moved westward in the background stratospheric easterly winds (e.g.Sellitto et al, 2022;Millán et al, 2022), and the positive N anomalies >3-sigma in Fig.2track the plume movement downstream from the January 15 eruption.The altitudes of the N anomaly maxima for these days generally range from 27 to 33 km.

Figure 1 .
Figure 1.Vertical profiles of refractivity (N) anomalies associated with a localized 1000 ppmv H2O perturbation at 30 km (blue) and a -8K temperature perturbation (red).Results are expressed in terms of percentage anomalies with respect to a standard background N profile.The solid black line shows the observed standard deviation () for C2 N measurements based on non-volcanic conditions in the vicinity of HTHH, and the dashed and dotted lines show 2and 3 values.

Figure 2 .
Figure 2. Grey dots show the C2 sounding locations in the region of the HTHH eruption on January 16, 17 and 18. Colored symbols indicate locations with positive N anomalies over altitudes 25-35 km exceeding 3-sigma background values, with color and size denoting the altitude and magnitude of the maximum N anomaly.

Figure 3 .
Figure 3. (a) Vertical profiles of refractivity (N) anomalies for cases selected on January 16(colored symbols in Fig.2a).(b-d) Vertical profiles of C2-derived H2O mixing ratios (grey lines) for the cases selected on January 16-18 (colored symbols in Figs.2a-c).The blue lines in (b-d) show radiosonde H2O measurements over Australia during these days.Note that the radiosonde data only extend upwards to the balloon-burst altitude near ~30 km; the highest measurements in each profile are indicated by the blue asterisks.

Figure 4 .
Figure 4. Red symbols show MLS temperature anomalies as a function of the C2 H2O retrievals in the HTHH plume over January 16-18 (the same profiles as identified in Fig. 3).Blue symbols show the corresponding dry temperature anomalies calculated in the C2 dry retrievals (atmPrf).Black dashed line corresponds to a slope of (-8K/1000 ppmv), consistent with the T-H2O sensitivity results shown in Fig. 1.

Figure 5 .
Figure 5. (a) Longitude vs. time diagram showing locations of H2O profiles with maximum values >700 ppmv over 25-35 km, over the domain 10-30 o S, for the first week following the HTHH eruption (noted by the 'X' on January 15).The colors and symbols denote the altitudes and magnitudes of the H2O maxima.The grey dots indicate all of the C2 measurement locations.(b,c) Time series of plume height and H2O maximum amount for C2 H2O retrievals in the HTHH plume.

Figure S1 .
Figure S1.Vertical profiles of (a) refractivity anomaly and (b) derived H2O mixing ratios for C2 observations within the HTHH plume on January 15.Most of these N profiles are associated with a 'bad' flag in the retrieval, although the structures are reasonable and may be consistent with actual geophysical information.