Atmospheric pressure
To study atmospheric signatures due to the 2022 Hunga Tonga eruption, I use atmospheric pressure data from a dense array of weather sensors “Soratena”, composed of ~ 3,000 stations deployed by Weather News Inc. throughout Japan (global.weathernews.com). Short-term enhancements of atmospheric pressure were recorded to pass over the Japanese Islands from SE to NW ~ 7 hours after the eruption, together with smaller signals ~ 38 hours after the eruption (global.weathernews.com/news/16551). They indicate the passages of LW propagated from the volcano along the short path (Fig. 1 inset), and the same wave after traveling round the Earth. In this study, I select 168 points uniformly distributed in Japan (Fig. 1) and analyze the 1-minute sampling atmospheric pressure time series in Jan.15–17 in local time (from Jan 14, 15 UT to Jan. 17, 15 UT). These data shall be compared with the ionospheric disturbances recorded by GNSS stations.
Ionospheric TEC
Global Positioning System (GPS), the oldest GNSS, has been widely used for studying ionospheric TEC. GPS satellites employ orbital periods of a half sidereal day. Hence, they cannot stay above a reasonable elevation angle from a ground station longer than 4–5 h, and GPS-TEC data are composed of numerous short arcs. On the other hand, Quasi-Zenith Satellite System (QZSS), the Japanese satellite system for positioning, is composed of three satellites with quasi-zenith orbits (J01, J02, J03) and one geostationary orbit satellite (J07). These satellites stay longer within the view of Japanese stations (~ 8 and 24 h a day for J01-03 and J07, respectively). Especially, J07 enables us to observe TEC changes during the passages of air waves over Japan for days without disruptions.
I use QZSS raw data from GEONET. I use the whole set made of ~ 1,500 points in drawing maps of ionospheric anomalies, but use a subset made of ~ 150 stations in drawing time series of wavelet transformed TEC anomalies. Line-of-sights to J07 are stationary and intersect the layer of maximum ionization (assumed at altitude 300 km) of the ionosphere at points shifted a few hundreds of kilometers toward SSW from stations (Fig. 1). On the other hands, J01, J02, and J03 come in turn and stay around the local zenith for ~ 8 hours a day so that we can monitor ionosphere just above stations with nearly vertical line-of-sights using either one of the three satellites. We converted the L1 and L2 carriers into slant TEC (STEC). STEC with J07 is converted to vertical TEC (VTEC) by multiplying with the cosine of the incidence angle of line-of-sight with the hypothetical thin layer at 300 km altitude. Basic procedures of the GNSS-TEC study follow Heki (2021).
Geomagnetic activities during this period were rather high, and the studied period (2022 Jan. 15–18) corresponds to the recovery phase from a geomagnetic storm that occurred on January 14. The daily Kp index was in the range 2.0–3.3 during this period (omniweb.gsfc.nasa.gov).
Passages Of The Atmospheric Waves
To isolate atmospheric signals associated with the LW passages, I approximate temporal variations of the pressure using polynomials and extract the departure from these reference curves as anomalies. Figure 2 gives four snapshots of their spatial distribution when the first LW passed over Japan. The signals look uniform throughout the country. The positive pressure anomalies lasting for ~ 20 minutes traveled from SE to NW with a speed of ~ 0.3 km/s. The positive peak is as strong as ~ 1.5 hPa, larger than ~ 0.1 hPa caused by the LW observed in Japan after the 1980 St. Helens eruption (Ogawa et al., 1982). Considering the similar distance of the Hunga Tonga and the St. Helens volcanoes from Japan (~ 8,000 km), this difference would reflect the intensity contrast of the two eruptions.
Figure 3 shows the pressure anomalies as the function of time and distance from the volcano. The first passage of LW occurred 11–12 UT, Jan. 15, as seen in clear red colors along black dashed lines indicating the expected arrival times assuming 0.3 km/s propagation speed. In the close-up view (Fig. 3a), we can see numerous weaker waves follow the initial pulse. In contrast, there are no significant signals before the initial arrival of LW.
There seems to be another passage of a wave with a slower velocity at 14–15 UT (Fig. 3a), which might be the internal gravity wave (IGW). Their propagation velocity of ~ 0.225 km/s is close to the tsunami speed through an ocean as deep as 5,200 m, and this IGW may have amplified the meteo-tsunami that hit Japanese coasts 3–6 hours after the LW arrival (www.jamstec.go.jp/j/jamstec_news/tonga/column03/).
Here I define t indicating time in hour counted from the beginning of Jan. 15. The second passage of LW propagated along the long path occurred ~ 8 UT on Jan. 16., i.e. t = ~ 32, as characterized by linear features with negative slopes. Its arrival times are earlier than expected by an hour (Fig. 3b). The second passage of IGW from NW to SE is expected t = ~ 44, but it seems to have occurred somewhat earlier, around t = ~ 41 (Fig. 3d). The third passage of LW from SE to NW is clearly seen at t = 48, early on Jan. 17, about a half hour earlier than expected.
Such differences in the arrival time are not due to the sphericity of the Earth. I calculate the distance between Japan and the volcano by approximating the Earth as a sphere with the average radius (6371.012 km). The real short-path distance along the ellipsoid surface is somewhat shorter (e.g., vldb.gsi.go.jp/sokuchi/surveycalc/surveycalc/bl2stf.html), but the difference does not exceed 20 km. Atmospheric structure (tropospheric temperature and sound velocity) difference would be mainly responsible for these arrival time differences.
Passage Of The Ionospheric Waves
TEC time series Jan. 15–18
Figure 4 shows continuous VTEC over four consecutive days observed with QZSS J07 at seven stations located along the LW propagation path. In addition to diurnal changes, we can see strong positive pulses 11–12 UT on Jan. 15. Their amplitudes are comparable to the background VTEC (~ 10 TECU, and 1 TECU implies 1016 electrons/m2). This amplitude is remarkable considering that the atmospheric pressure anomalies remain in the order of 1/1000 of the background values. These initial anomalies arrive later for stations farther from the volcano reflecting the propagation of the anomalies from SE to NW.
The anomalies associated with the second LW passage around t = 32–34 is not so clear as the first one because the passage occurred in the afternoon when the background TEC rapidly decays. The anomalies appear as small undulations of the smooth decays and occur earlier at northern stations due to LW propagation from NW to SE. The anomalies associated with the third passage are recognized as weak and strong peaks around the noon Jan. 17, t = ~ 48 and t = ~ 50, respectively. The occurrences of these peaks at the two ends have a time lag of ~ 0.5 hours, corresponding to the time of LW propagating for ~ 600 km from SE to NW. Passage of the fourth wave t = ~ 72 is not clear in VTEC time series.
Figures 2 and 3 showed that atmospheric pressure signatures are quite uniform throughout Japan. On the other hand, Fig. 5 indicates that TEC signatures of the first LW passage are highly variable along the small circle with the same distance from the volcano. For example, peaks are smaller, and the shapes are simpler in northeastern stations, e.g., 0019, than those in southwest, e.g., 0099. We can identify two peaks at 0165 and 0049 while only the earlier peak is clear at 0019 and 0054.
Detection of wave passages using wavelet
Because it is not easy to recognize the wave passages directly in VTEC plots (Figs. 4 and 5), I extract disturbance signals with a certain period (~ 30 minutes in this case) from TEC time series by wavelet transformation (Heki and Ping, 2005). Figure 6e, f shows the results for J07 continuous VTEC records during Jan.15–18 at 42 selected stations. The J01-J03 data around the four passages of waves are magnified in Fig. 6a-d. Time windows of Fig. 6a-c correspond to those in Fig. 3a-c (the atmospheric pressure data not available for the fourth passage).
Figure 6 are much noisier than the barometer data (Fig. 3). Nevertheless, we can recognize components traveling from SE to NW around t = 12–14 and t = 48–52, and from NW to SE around t = 32–34 and t = 70–72. The strength of four passage signals cannot be simply compared because they depend on background VTEC, which varies by an order of magnitude during a day (gray curves in Fig. 6e,f). For example, the first two passages occurred during evening hours when background VTEC is less than 10 TECU while the third and the fourth passages during daytimes. The fifth passage should be around t = 84–86, but Fig. 6f does not show any signatures for that.
There are several distinct differences in the ionospheric data (Fig. 6) from the barometer data (Fig. 3). An important difference is that significant ionospheric anomalies (Fig. 6a) seem to start well before the arrival of the pressure anomaly (Fig. 3a) by the LW passage. This is discussed in the next chapter. In the third passage, LW signatures are clear at t = 48–49 in barometer records (Fig. 3c). However, ionospheric signals at t = 48–49 are weaker than those around t = 50–52.
Distinction from MSTID
To interpret Fig. 6, it is important to distinguish the eruption-origin wave passage signatures from those by medium scale traveling ionospheric disturbances (MSTIDs). In Fig. 6f, we recognize multiple clear stripes with negative slopes during Jan. 18 daytime (t = 74–80). Here I demonstrate that these are the passage of winter daytime MSTIDs irrelevant to the Hunga Tonga eruption.
A statistical study by Otsuka et al. (2011) showed that MSTIDs are classified into summer nighttime MSTID and winter daytime MSTID. The latter has nearly EW wavefronts and propagate southward by 0.1–0.2 km/s as IGW. In January daytime, the occurrence rate of such MSTID exceeds 40%. Their signatures might be misinterpreted as the wave passages from NW to SE caused by the Hunga Tonga eruption.
Figure 7 shows four snapshots of the positive and negative TEC anomalies around the noon, Jan. 18 (t = ~ 75). Their wavefront strikes nearly EW, which is significantly different from the LW wavefront running NE-SW (Fig. 2). Moreover, their propagation speed is ~ 0.1 km/s, much slower than LW. Hence, stripes appearing on Jan. 18 daytime (Fig. 6f) would merely show MSTID irrelevant to LW from the volcano.
In contrast, the second passage (Fig. 6b) occurred late in the afternoon (t = 32–33, corresponding 5–6 PM in local time), and it is unlikely that MSTID occurred at such local times (Otsuka et al., 2011). The fourth passage seen in Fig. 6d would not be due to MSTID, either, because it occurred early in the morning (t = 70–71, corresponding to 7–8 AM in local time) when the occurrence probability of MSTID is still low.
Relation Between The Ionospheric And The Atmospheric Waves
Two TEC peaks in the first passage and their arrival times
After the 1980 St. Helens eruption, Ogawa et al. (1982) compared the barometric and ionospheric signatures in Japan. However, they could not discuss precise time lags between them because the barometers and ionospheric SIPs are separated by > 200 km. Dense coverage of GNSS and barometers in the present case (Fig. 1) enables us to discuss exact time lags between the disturbances in the two spheres.
The VTEC anomalies propagate with the LW velocity of ~ 0.3 km/s for all the four passages (Fig. 6a-d). However, the occurrences of the strong peaks in TEC are not simultaneous with those of the Lamb pulses. Figure 8 compares the pressure and TEC changes in the same time windows at four positions, where the J07 SIPs of the GNSS stations and barometers nearly overlap. Figure 5 showed that the TEC disturbances are highly variable from place to place, and this is the case for these four stations. We recognize the existence of two peaks in Fig. 8a and Fig. 8d, the first peak comes ~ 40 minutes before the LW peaks and the second peak occurs ~ 40 minutes later than the LW. The second peak is not always clear and is almost absent in Fig. 8b and Fig. 8c. This is irrelevant to the difference in the irregularity of the Earth’s surface because their SIPs are located near the coast, that is, the LW has travelled above flat ocean surface and just arrived at the coast in all these cases.
Figure 9 shows the VTEC changes at almost the same points in ionosphere as Fig. 8. I selected the satellite J03, which was near zenith at that time, and plotted VTEC using GNSS stations close to the barometers. VTEC shown in Figs. 8 and 9 should nearly coincide with each other. However, they show subtle differences especially around sharp peaks of TEC anomalies, possibly due to complex spatial structures of the electron density anomalies. Nevertheless, existence of the two peaks in (a) and (d) and absence of the second peak in (b) and (c) are common in Figs. 8 and 9.
Origins of the two TEC peaks
Although the emergences of the first TEC peaks are earlier than the LW arrivals (Figs. 8 and 9), its propagation velocity is consistent with the LW velocity of ~ 0.3 km/s (Fig. 6). Thus, the ionospheric disturbance would be a secondary feature made by energy leakage from the atmospheric LW. This is also supported by HF Doppler observations after the 1980 St. Helens eruption. Ogawa et al. (1982) found a subtle phase lag of a few minutes between the HF Doppler records with two different frequencies (5 and 8 MHz) and considered that the phase lag reflects the upward propagation of the disturbance. Here, I assume that ionospheric disturbances are caused by upward energy leakage from LW propagating through the Earth’s troposphere.
Interpretation of the backward energy leakage, i.e., TEC changes occurring after the arrivals of pressure waves, is straightforward. The sharp peaks of atmospheric pulses (Figs. 8, 9) include various frequency components and would excite atmospheric waves over a wide spectral range, including both acoustic wave (AW) and IGW. They would reach the ionospheric F region and disturb the TEC in 8–10 minutes (AW) or 20–30 minutes (IGW). TEC peaks occurring after the LW passage in Figs. 8 and 9 would originate from such atmospheric waves.
In contrast, forward energy leakage, i.e., the TEC increases occurring ~ 40 minutes prior to the LW arrival (Figs. 8 and 9), may seem enigmatic. I propose a simple qualitative model as sketched in Fig. 10. There, I assume that the LW traveling northwestward from Tonga just reached the Pacific coast of Japan at coordinate zero (white lines in Fig. 10a-c). I then simulate how AW is excited by LW and propagate.
To make it simple, the LW is supposed to excite N-shaped AW pulse every 5 minutes (I assume the period of AW 2 minutes). The sound speed profile shown in Fig. 10e does not allow AW to propagate straight but makes the ray paths refracted downward. Figure 10a, b, and c show the ray paths of AW generated 10, 15, and 20 minutes ago, respectively. There I calculated electron density anomalies following Kundu et al. (2021). AW excited 10 minutes ago (Fig. 10a) has just reached the ionospheric F region. On the other hand, a part of the ray paths of AW excited 15 minutes ago (Fig. 10b) is refracted and extends almost horizontally. Because of higher velocity in ionosphere, its northeastern front is located beyond the LW wavefront. AW excited 20 minutes ago (Fig. 10c) has already overtaken the LW wavefront and disturbs the ionosphere located ~ 200 km ahead of it.
Figure 10d shows 8 sets of acoustic ray paths excited every 5 minutes in the past back to 40 minutes ago. The northwestern tip of wavefronts of such “old” waves goes further beyond the LW wavefront by propagating within high velocity ionosphere. Because of the geometric decay of old AW, intense ionospheric disturbances would occur only within a certain distance from the LW wavefront. In the case of the 2022 Jan. 15 Tonga eruption, the ionospheric disturbance extended ~ 700 km beyond the LW wavefront, causing the apparent delay of LW arrival by ~ 40 minutes. As seen in Fig. 5, the amplitudes and waveforms of the ionospheric disturbances are highly variable along the wavefront. The forward energy leakage would make ionospheric electron density anomalies quite complex in space. A realistic numerical simulation is beyond the scope of this paper, but it would validate this hypothesis in the future.