Alternative Approach for Tsunami Early Warning Indicated by Gravity Wave Effects on Ionosphere

Abstract

The rapid displacement of the ocean floor during large ocean earthquakes or volcanic eruptions causes the propagation of tsunami waves on the surface of the ocean, and consequently internal gravity waves (IGWs) in the atmosphere. IGWs pierce through the troposphere and into the ionospheric layer. In addition to transferring energy to the ionosphere, they cause significant variations in ionospheric parameters, so they have considerable effects on the propagation of radio waves through this dispersive medium. In this study, double-frequency measurements of the Global Positioning System (GPS) and ionosonde data were used to determine the ionospheric disturbances and irregularities in response to the tsunami induced by the 2011 Tohoku earthquake. The critical frequency of the F2 layer (foF₂) data obtained from the ionosonde data also showed clear disturbances that were consistent with the GPS observations. IGWs and tsunami waves have similar propagation properties, and IGWs were detected about 25 min faster than tsunami waves in GPS ground stations at the United States west coast, located about 7900 km away from the tsunami’s epicenter. As IGWs have a high vertical propagation velocity, and propagate obliquely into the atmosphere, IGWs can also be used for tsunami early warning. To further investigate the spatial variation in ionospheric electron density (IED), ionospheric profiles from FORMOSAT-3/COSMIC (F3/C) satellites were investigated for both reference and observation periods. During the tsunami, the reduction in IED started from 200 km and continued up to 272 km altitude. The minimum observed reduction was 2.68 × 10⁵ el/cm³, which has happened at 222 km altitude. The IED increased up to 767 km altitude continuously, such that the maximum increase was 3.77 × 10⁵ el/cm³ at 355 km altitude.

1. Introduction

Based on the electric charge, the Earth’s atmosphere is divided into two major layers, the troposphere and the ionosphere. The troposphere, in which the atmospheric components are electrically neutral, is the lower atmosphere ranging from the surface to about 60 km [1]. The ionosphere is the upper part of the Earth’s atmosphere, which is extended from approximately 60 km to more than 1000 km. Solar radiation produces free electrons and ions in this region that affect the propagation of electromagnetic waves [2]. Studying the coupling between these two layers has been an interesting topic for atmosphere and space weather research for several decades [3]. A tsunami displaces the atmosphere as it propagates across the open ocean, the atmosphere responds to this excitation by propagating gravity waves obliquely upward, these waves are referred to as internal gravity waves (IGWs), and grow nearly exponentially with height as they proceed into the rarefied regions of the upper atmosphere [4]. IGWs are buoyancy oscillations that can propagate horizontally and vertically, their propagation is under the gravity force of the Earth [5]. On the one hand, the ionospheric disturbances can decrease the performance of precise positioning and navigation [6]. On the other hand, the ionospheric disturbance signals induced by some tropospheric events, such as earthquakes, tsunamis, and lightning, can be used to study these events [7,8,9].

A tsunami is generated when a large oceanic earthquake or volcanic eruption causes rapid displacement of the ocean floor. Studies about the propagation of IGWs during the 1970s have suggested that the ionosphere is sensitive to IGWs by the forcing effect of a tsunami on the surrounding atmosphere [10]. Peltier and Hines inferred that IGWs might be detectable and used for tsunami warning system purposes [4]. After the Sumatra and Indian tsunami in 2004, researchers pay significant attention to observing tsunamis by ionospheric sounding. Using the high-frequency Doppler sounding network in Taiwan, Liu et al. observed ionospheric disturbances triggered by the 2004 Indian Ocean tsunami [11]. Occhipinti et al. used total electron content (TEC) data measured by the Jason-1 and Topex/Poseidon satellite altimeters to detect the signature of IGWs on the ionosphere [12]. However, both Doppler sounding networks and altimetry satellites have a low spatial and temporal resolution, which makes it difficult to study the propagation characteristics of ionospheric disturbances in detail. Artru et al. used TEC data, observed from the Geonet network in Japan, to extract the ionospheric disturbances induced by the 23 June 2001 8.2 Mw earthquake in Peru for the first time [8]. Due to the high spatial and temporal resolution, the GPS TEC data have been widely used in ionospheric monitoring and tsunami detection [13,14,15,16].

According to the US Geological Survey, an earthquake (Mw = 9) occurred with the epicenter at 38.32°N, 142.37°E; it occurred on 11 March 2011 at 5:46:23 UT near Tohoku, Japan, and then induced a tsunami. In this paper, we employ various geodetic techniques to detect ionospheric disturbances and ionospheric irregularities induced by the tsunami on the west coast of the United States of America. According to the travel time map in Figure 1, the tsunami reached the west coast of the USA about 10 h after the earthquake. Various researchers have studied the arrival time of IGWs in the 2011 Tohoku tsunami, but their approach was mainly post-process. In this paper we investigated a near real-time method to detect IGWs, therefore using this procedure can lead to faster tsunami warnings compared to DARTs and tide gauges. Moreover, ionospheric perturbations caused by tsunamis were also studied by different satellite geodetic techniques.

2. Materials

In this study, 56 GPS stations were selected from the UNAVCO network, (http://www.unavco.org/data/gps-gnss.html, accessed on 20 April 2020). The sampling interval for all stations was 15 s or 30 s. To identify the peculiar signatures of the ionosphere during the tsunami, data of one day before and one day after the tsunami day were analyzed. Data from ionosonde, the deep-ocean assessment and reporting of tsunamis (DART), and tide gauges were used to validate the obtained results. Figure 2 shows the spatial distribution of the stations and the types of instruments used in the present analysis.

Solar and geomagnetic activities are the dominating factors that control the behavior of the ionosphere. Geomagnetic storms and solar activity mask the effects of tropospheric events in the ionosphere, and responses of the F region to tsunamis can be sought only under particular geomagnetic and solar conditions [18]. The index F10.7 was 129.5, 121.5, and 119.2 solar flux units (SFU) on 10, 11, and 12 March, respectively. According to the F10.7 index variation, solar radiation was in the steady state. The geomagnetic activity can be interpreted using several parameters. One of these parameters is the Kp index, which varies from 0 to 9 (values equal and greater than 5 represent a geomagnetic storm) [19]. Another parameter is the Ap index, which varies from 0 to 400 (values greater than 50 represent a geomagnetic storm) [20]. Figure 3 shows the geomagnetic activity level. The solar and geomagnetic indices were downloaded from the archives of Goddard Space Flight Center (https://omniweb.gsfc.nasa.gov/form/dx1.html, accessed on 20 April 2020).

Referring to Figure 3, the Kp index was below 5 units on March 11, from 6 to 18 UT, which indicates quiet geomagnetic conditions during the major part of the tsunami event.

3. Methods

3.1. TEC Measurements

Using dual-frequency GPS observations ($f_1$= 1575.42 MHz, $f_2$ = 1227.60 MHz), the ionospheric delay can be calculated. Then smoothed code pseudorange measurements for $f_1$ and $f_2$ are estimated using the following Equation (1) [21]:

$$\{\begin{array}{c}{\tilde{P}}_{1,i}^k={\overline{P}}_{1,i}^k+L_{1,i}^k-{\overline{L}}_{1,i}^k+2\frac{f_2^2}{f_1^2-f_2^2}(L_{1,i}^k-L_{2,i}^k-({\overline{L}}_{1,i}^k-{\overline{L}}_{2,i}^k))\\ {\tilde{P}}_{2,i}^k={\overline{P}}_{2,i}^k+L_{2,i}^k-{\overline{L}}_{2,i}^k+2\frac{f_1^2}{f_1^2-f_2^2}(L_{1,i}^k-L_{2,i}^k-({\overline{L}}_{1,i}^k-{\overline{L}}_{2,i}^k))\end{array}$$

(1)

where k and i denote the pseudo-random noise (PRN) of a given satellite and receiver, respectively; ${\tilde{P}}_{1,i}^k$ and ${\tilde{P}}_{2,i}^k$ are the smoothed dual-frequency code measurements at epoch t; $L_{1,i}^k$ and $L_{2,i}^k$ are the phase measurements observed by receiver i as seen from satellite k at epoch t; ${\overline{L}}_{1,i}^k$ and ${\overline{L}}_{2,i}^k$ are the mean phase measurements at cycle slip free time interval; ${\overline{P}}_{1,i}^k$ and ${\overline{P}}_{2,i}^k$ are the mean code measurements. Finally, slant TEC (STEC) value is calculated using the following equation:

$$STE{C}_i^k(t)=\frac{f_1^2{f}_2^2}{40.31(f_2^2-f_1^2)}({\tilde{P}}_{i,4}^k(t)-DC{B}^k-DC{B}_i+{\epsilon }_L)$$

(2)

In Equation (2) $STE{C}_i^k(t)$ is slant TEC between satellite k and receiver i, ${\tilde{P}}_{i,4}^k(t)$ is the smoothed code measurement, or the so-called geometry-free linear combination. $DC{B}^k$ and $DC{B}_i$ denote the satellite and receiver differential code biases, respectively. In this study, $DC{B}_i$ value is calculated using the Li et al. method [22]. In this method, the global ionospheric map (GIM) TEC values are interpolated to the same GPS footprint locations and times. The $DC{B}_i$ of each station is estimated by least squares adjustment using the observations of all satellites. To remove multipath effects, only the measurements with elevation angles higher than 20° are considered. STEC along the GPS line of sight (LOS) is converted into the vertical TEC (VTEC) using the mapping function as follows:

$$VTEC=STEC\cos ({\sin }^{-1}(\frac{R_e\sin z}{R_e+H_{ion}}))$$

(3)

where $R_e$ represents the mean radius of the Earth, $z$ is the zenith distance of the LOS from the receiver to GPS satellites, and $H_{ion}$ is the altitude of the ionospheric thin-layer shell, which is set to 450 km in this study.

3.2. TEC Disturbances

Hernández-Pajares et al., Ref. [23], used a first-order numerical difference method to extract ionospheric disturbances. Their method uses Equation (4), as follows:

$$\mathsf{\Delta }TEC(t)=\mathsf{\Delta }s(t)=s(t)-0.5(s(t-\tau )+s(t+\tau ))$$

(4)

where $\mathsf{\Delta }s(t)$ is the first-order numerical difference of TEC series, t and $\tau$ denote the observation epoch and time step, respectively. This method can effectively capture the TEC variations for measurements with elevation angle higher than 40°. In this paper, a second-order numerical difference method is employed to eliminate the TEC trend because the satellites at low elevation can better highlight the presence of oscillating perturbations [12]. The second-order numerical difference is calculated by Equation (5), as follows:

$${\mathsf{\Delta }}^{2}TEC(t)={\mathsf{\Delta }}^{2}S(t)=\mathsf{\Delta }S(t)-\frac{1}{2}(\mathsf{\Delta }S(t-\tau )+\mathsf{\Delta }S(t+\tau ))$$

(5)

The ${\mathsf{\Delta }}^{2}S(t)$ values and primary IGW signals have the same period ($T$), and the ratio amplitude ($R$) is calculated using Equation (6), as follows:

$$R=4{\sin }^4\left(\pi ·\tau /T\right)$$

(6)

The method is very effective when $R$ is larger than 1, whereas the sensitivity will degrade when $R$ is smaller than 0.4 [24]. If equals 300 s, the sensitive period is 6 to 26 min with the ratio of amplitude bigger than 0.4. The relation between $R$ and $T$ is plotted in Figure 4 for $\tau$ equal to 300 s.

3.3. Ionospheric Irregularities

The rate of TEC index (ROTI) defined as the standard deviation of the rate of TEC (ROT) is used to indicate the presence of ionospheric irregularities. ROTI is calculated with a sliding window for each five minute interval using Equation (7) [25,26].

$$ROTI=\sqrt{\langle RO{T}^2\rangle -{\langle ROT\rangle }^2}$$

(7)

ROT (in TECU/minute) is calculated as follows:

$$ROT(t)=\frac{L_4(t)-L_4(t-1)}{\mathsf{\Delta }t\times {10}^{16}\times 40.3\times (\frac{1}{f_1^2}-\frac{1}{f_2^2})}$$

(8)

where $L_4(t)$ is the geometry-free phase combination at epoch t; $f_1$ and $f_1$ are frequency; and $\mathsf{\Delta }t$ is the time difference between the epochs (in minutes).

4. Results

4.1. IGWs

The location of the ionospheric pierce point (IPP) is calculated for all of the satellites. Figure 5 shows the IPP paths over the west American coast, between 15:00 UT to 20:00 UT, corresponding to each station. Due to limited space and for simplification, the figure only shows one IPP path for each station. This figure shows the spatial distribution of GPS stations with different colors. Each station is related to only one IPP path, which has the same color (e.g., station P001 with red color corresponding with PRN 21).

The IGWs originated by the ocean tsunamis propagate in a conical shape through the atmosphere, with their vertex located at the epicenter of the tsunami. Therefore, the IGWs velocity has both vertical and horizontal components. As already known, by the increase in the altitude, the atmospheric density reduces exponentially, thus the amplitude of the IGWs increase as they propagate to higher altitudes [27]. One way to distinguish the signals associated with IGWs induced by the tsunami is to search for disturbances that are propagating at ∼200–300 m/s in an outward direction from the tsunami’s source [28]. The vertical propagation velocities of IGWs are 40–50 m/s [8]. Figure 6 shows the ${\mathsf{\Delta }}^{2}TEC$ time series of PRN18, PRN22, and PRN29 observed from different stations from 14:30 UT to 21:00 UT on 11 March. The ${\mathsf{\Delta }}^{2}TEC$ time series of PRN 18, 22, and 29 are shown as an example to show the ionospheric disturbances from the beginning to the end of the propagation of IGWs. The ionospheric disturbances appeared over the west American coast at about 15:10 UT, then faded at about 20:00 UT. The horizontal speed of the ionospheric disturbances was between 200 and 300 m/s. Figure 7 shows the locations of the IPPs of PRN18, PRN22, and PRN29 observed from all the GPS stations from 15:10 UT to 20:00 UT on 11 March. As it can be seen, within the land area, PRN 18 is further west compared to PRN 22 and 29, so it has detected IGW signatures earlier. Within the ocean area, PRN 22 was further west and should have been able to detect IGWs earlier, but due to the oblique propagation direction of IGWs, the signatures were not detectable there.

Figure 8 shows the ${\mathsf{\Delta }}^{2}TEC$ time series and time-frequency analysis of PRN 18, PRN 22, and PRN 29 observed from the GPS station PABH. The GPS station PABH is northwest of America and is selected as a sample. At the beginning of the propagation of IGWs, the dominant frequency and amplitude of the ${\mathsf{\Delta }}^{2}TEC$ time series was small (the dominant frequency of the ${\mathsf{\Delta }}^{2}TEC$ time series was between 0.75 and 1.5 mHz in Figure 8f), and after a while, the dominant frequency and amplitude of the ${\mathsf{\Delta }}^{2}TEC$ time series increased (the dominant frequency of the ${\mathsf{\Delta }}^{2}TEC$ time series was between 1.5 and 2 mHz in Figure 8b,d). At the end of the propagation of IGWs, the dominant frequency and amplitude of the ${\mathsf{\Delta }}^{2}TEC$ time series decreased again (the dominant frequency of the ${\mathsf{\Delta }}^{2}TEC$ time series was between 0.75 and 1.5 mHz in Figure 8d).

Figure 9 and Figure 10 show the ${\mathsf{\Delta }}^{2}TEC$ time series of PRN18, PRN22, and PRN29 observed from the same stations from 14:30 UT to 21:00 UT on 10 March (a day before the tsunami) and 12 March (a day after the tsunami). There are no obvious disturbances in the ${\mathsf{\Delta }}^{2}TEC$ time series on 10 and 12 March. The average values of ${\mathsf{\Delta }}^{2}TEC$ between 15:10 UT and 20:00 UT (time of the propagation of IGWs to the ionosphere on 11 March) on the three days of March 10, 11, and 12 of 2011 were 0.01, 0.04, and 0.01 TECU, respectively.

Figure 11 shows the locations of the DART and tide gauge stations and the GPS station CABL. Also, the location of the IPPs for PRN15 and PRN16 observed from station CABL are plotted in the figure. The IPPs of PRN15 and PRN16 were over the ocean and land, respectively. The PRN15, PRN16, and CABL stations were selected as a sample.

Figure 12 shows the ${\mathsf{\Delta }}^{2}TEC$ time series, time-frequency analysis, and ROTI time series of PRN15 and PRN16 observed from the GPS station CABL. The ionospheric disturbances and ionospheric irregularities are observed above the land and sea, but the ionospheric disturbances and ROTI values of PRN15 were more than PRN16. The ${\mathsf{\Delta }}^{2}TEC$ time series of PRN15 and PRN16 observed from the GPS station CABL have similar frequency characteristics. This phenomenon illustrates that the ocean tsunami generates IGWs that propagate obliquely in the ionosphere [15]. The results of the time-frequency analysis indicate that the average dominant frequency of ionospheric disturbances was 1.2 ± 0.3 mHz with a period of 14.42 min on 11 March.

4.1.1. Horizontal Speed Phase of IGWs

The ionospheric disturbances are observed from 15:10 UT to 20:00 UT on 11 March. According to [29], we interpolated the ionospheric disturbances in meridional and zonal wavenumbers at 30-s intervals between 15:10 UT and 20:00 UT because the sampling interval of some of the GPS stations was 30 sec. The wavenumbers were calculated in meridional and zonal wavenumbers using fast Fourier transformation (FFT) analyses of ionospheric disturbances. The horizontal phase speed of the ionospheric disturbances was calculated at 30-s intervals between 15:10 UT and 20:00 UT, using Equations (9) and (10).

$${\omega }_r=\frac{2\pi }{{\tau }_r}$$

(9)

where ${\tau }_r$ is the period of ionospheric disturbances and ${\omega }_r$ is the ground-based frequency.

$$c_H=\frac{1}{\sqrt{{\left(k/{\omega }_r\right)}^2+{\left(l/{\omega }_r\right)}^2}}=\frac{1}{\sqrt{{\left(1/c_x\right)}^2+{\left(1/c_y\right)}^2}}$$

(10)

where k and l are the zonal and meridional wave numbers, respectively, and $c_H$, $c_x$, and $c_y$ are the horizontal phase speed, zonal phase speed, and meridional phase speed of IGWs. The results indicate that the horizontal phase speed of the IGWs was 231.31 ± 44 m/s.

4.1.2. Ionosphere Response to IGWs

In this study, we calculated ROTI to investigate ionospheric irregularities in the IGWs induced by the tsunami. Figure 13 shows the average ROTI values for all of the satellites observed from all of the GPS stations used on March 10, 11, and 12. The ROTI values (red curve) suddenly increased after 15:10 UT and again decreased back to normal values at about 20:00 UT on 11 March compared to the ROTI values on 10 and 12 March (green and blue curves). The average values of ROTI between 15:10 UT and 20:00 UT (the time of propagation of the IGWs to the ionosphere on 11 March) on the three days were 0.25 ± 0.01, 0.27 ± 0.008, and 0.25 ± 0.01 TECU/min, respectively.

4.2. Tsunami Waves Activity

We compared the ${\mathsf{\Delta }}^{2}TEC$ time series with sea level time series measured by all of the active tide gauges and DART stations on the west American coast on 11 March, to detect the correlation between the propagation characteristics of the IGWs and tsunami waves. The sampling intervals of all of the stations were 60 s. The DART station data is available at (www.ndbc.noaa.gov/dart.shtml, accessed on 20 July 2020) and the tide gauge station data is available at (www.ioc-sealevelmonitoring.org/list.php, accessed on 20 July 2020). We used the same second-order numerical difference method (Equation (5)) to process the sea level time series. The arrival times of the tsunami waves to the different stations are listed in

Table 1. The time arrival of tsunami waves to different stations.

Station	Tsunami Time Arrival (UT)	Latitude (° N)	Longitude (° W)
Alameda	16:36	37.77	−122.3
Arena	15:34	38.91	−123.72
Astoria	16:24	46.21	−123.77
Cherry Point	17:08	48.86	−122.76
Jolla	16:22	32.87	−117.26
San Fran	16:15	37.81	−122.47
South Beach	15:42	44.63	−124.04
West Port	15:39	46.91	−124.11
46404	14:33	45.85	−128.78
46407	14:39	42.71	−128.83
46411	14:59	39.34	−127.07

. The DARTs stations detected the tsunami waves earlier than the tide gauges, due to the location of the DART stations, which were 300 km from the coast, while the tide gauges were near the beach. The first time, the tide gauges at Arena and West Port observed tsunami waves at 15:34 UT and 15:39 UT, respectively. The locations and names of the tide gauge and DART stations are shown in Figure 11. Figure 14 shows the time-frequency analysis and time series of the Arena and West Port station. The red dots in Figure 14a,c indicate the time of arrival of the tsunami waves at the stations. The results of the time-frequency analysis indicate that the average dominant frequency of the tsunami waves was 1.86 ± 1 mHz with a period of 10.3 min on 11 March.

The tsunami propagation velocity was computed using the equation of shallow water [8], as follows:

$$V_{tsunami}=\sqrt{gH}$$

(11)

where g is the gravity with a value of 9.8 m/s², and H is the ocean depths in meters. The ocean depth value was obtained from the DART data. The tsunami propagation velocity varies from 163.8 to 204.3 m/s.

Table 2. Tsunami propagation velocity of range.

Station	Tsunami Propagation Velocity (m/s)
46404	163.8
46407	178.9
46411	204.3

shows the variation range of the tsunami propagation velocity at each DART station. Finally, the average velocity of the tsunami waves is 182.25 ± 18 m/s.

4.3. Validation of Results

4.3.1. Using Ionosonde Data

The foF₂ data observed by the Point Arguella (PA836) ionosonde station were used for validation of the ionospheric disturbances. The data for this station are available at (ftp://ftp.ngdc.noaa.gov/ionosonde, accessed on 20 June 2020). The foF₂ can be converted to electron density using Equation (12).

$$f[Hz]=8.98\sqrt{N_e}(m^{-3})$$

(12)

The diurnal variation in the foF₂ parameter from the 10 to 12 March is shown in Figure 15 with a red color. The range of variation in foF₂ under the quiet period (10 and 12 March) was calculated using Equation (13).

$$\begin{array}{l}aver=average(fo{f}_2(10),fo{f}_2(12))\\ std=stdev(fo{f}_2(10),fo{f}_2(12))\\ upper=aver+std\\ lower=aver-std\end{array}$$

(13)

Figure 15 shows that the foF₂ values increased suddenly after 15:00 UT and then decreased back to the normal value at 20:00 UT on 11 March, compared to the quiet days (10 and 12 March). The foF2 values were not between the lower and upper limit from 15:00 UT to 20:00 UT on 11 March (Figure 15b), while on the 10th and 12th of March (Figure 15a,c) the values were completely between the lower and upper limit. Figure 16 shows the difference between the upper value and the diurnal variation in the foF2 parameter in each epoch from 14:00 UT to 22:00 UT on 11 March 2011. The ionosonde station did not have data from 17:00 to 18:15 UT on 11 March (Figure 15b). Only the PA836 station had data in the Western United States on 11 March 2011. The abnormal onset and end times of the GPS observation were 15:10 UT and 20:00 UT. Therefore, the results of ionosonde are consistent with the GPS analysis. It is also worth noting that there were no obvious disturbances on the day before (10 March) and the day after (12 March), using the ionosonde data.

4.3.2. Using Electron Density Disturbances

One of the techniques to obtain the ionospheric profiles of electron density is GPS radio occultation (RO). We used ionospheric electron density (IED) profiles obtained from F3/C RO measurements to detect the vertical ionospheric disturbances. The data are provided by CDDAC (COSMIC Data Analysis and Archive Center, http://cdaac-www.cosmic.ucar.edu, accessed on 20 August 2020). Figure 17 shows the location of the ionospheric occultation tangent points of two RO profiles and the set of four IPP tracks of GPS PRN 18 from the southern stations P001, P014, P473, and P796.

Figure 18a shows the vertical profile of electron density (red curve) measured from F3/C RO at 18:36 UT on 11 March 2011 (reference period). This is a setting occultation between the F3/C satellite (C-006) and PRN-31, starting with the occultation (tangent) point altitude of 62 km and subsequently ascending to 767 km. The location of the occultation point upon starting is 121.35°W, 28.59°N and at the ending is 106.88°W, 28.43°N. This means that the occultation point not only ascends vertically from the bottom to top, but also sweeps horizontally a large distance (~1414 km) from west to east. The tangent point projections of this profile are shown in Figure 17 with a red color. PRN31 and C006 measured another RO electron density profile at 107.13°W, 30.98°N at ~18:09 UT on 12 March 2011 (observation period). The tangent point projections of this RO profile are shown in Figure 17 with a green color. Figure 18a shows the vertical profile of electron density with a green color. Another RO electron density profile was measured at 19:11 UT on 14 March by PRN24 and C006. The tangent point projections of this profile are shown in Figure 17 with a blue color. No RO electron density profile was measured near the reference period on 10 and 13 March. Figure 18a shows the vertical profile of electron density with a blue color. As Figure 18a shows, there is little difference between the amount of electron density on March 12 and March 14, indicating no especial ionospheric event. However, there is a significant difference between the amount of electron density on day 11 and the other two days (12 and 14 March). In this study, the difference and the percentage of electron density variations between the reference and observation periods have been calculated using Equations (14) and (15). The recorded IED profiles are interpolated to a height interval of 1 km.

$$difference=N{e}_{reference}-N{e}_{observation}$$

(14)

$$Percent\_Variation=\frac{N{e}_{reference}-N{e}_{observation}}{N{e}_{observation}}\times 100$$

(15)

According to Figure 18b, the reduction in IED started from 200 km and continued up to 272 km altitude, and the maximum reduction was 2.68 × 10⁵ el/cm³ (~27%), which happened at the 222 km altitude. The IED increased up to 767 km altitude continuously, such that the maximum increase was 3.77 × 10⁵ el/cm³ (~64%) at 355 km altitude. The electron density peak was 12.22 × 10⁵ el/cm³ at 291 km altitude on 11 March, and the value was 12.18 × 10⁵ el/cm³ at 258 km altitude. Generally, it can be inferred that the tsunami induced an increase in the IED peak value and altitude.

Figure 19 shows the ${\mathsf{\Delta }}^{2}TEC$ values of PRN 18 observed from the mentioned stations on 10, 11, and 12 March 2011. There are no obvious disturbances in the ${\mathsf{\Delta }}^{2}TEC$ time series of PRN 18 observed from the four GPS stations on 10 and 12 March. Thus, it is clear that the IGWs propagated during the tsunami, reached the F2 layer peak altitude and caused electron density disturbances.

5. Discussion

In this study, the ionospheric disturbances induced by IGWs during the tsunami on 11 March 2011 on the west American coast were detected using space-based and ground-based data.

Sea level disturbances and ionospheric disturbances have similar propagation characteristics, such as speed, frequency, and arrival time. Ionosphere disturbances and ionospheric irregularities were observed in the time series of satellites that passed over the land and ocean in the region between the tsunami epicenter and the observation stations. The ionospheric disturbances and ROTI values over the land were higher than those of over the ocean. This phenomenon is the result of the IGWs propagation in the oblique direction. The ROTI, ${\mathsf{\Delta }}^{2}TEC$, and ionosonde variations were not significant because we studied the effect of IGWs 7900 km away from the epicenter. As expected, water level disturbances are greater near the epicenter than those far from the epicenter; therefore, the ROTI, TEC, and ionosonde variations would have been significant if ionospheric disturbances were detected using the GPS stations of Japan. Near the epicenter, IGWs propagate due to variations in the crust of the earth and sea level. Far from the epicenter, IGWs propagate due to variations in the sea level [28,30]. The relative contribution of the tsunami-driven IGWs, as compared to IGWs generated by the earthquake itself, deserves further research [28].

The average of ionospheric disturbances between 15:10 UT and 20:00 UT on 11 March was four times larger than the average amount on 10 and 12 March (between 15:10 UT and 20:00 UT). The average ROTI value between 15:10 UT and 20:00 UT on 11 March was 1.1 times larger than that of 10 and 12 March (between 15:10 UT and 20:00 UT).

The reduction in IED started from 200 km and continued up to 272 km altitude, and the maximum reduction was 2.68 × 10⁵ el/cm³ (~27%), which happened at 222 km altitude. The IED increased up to 767 km altitude continuously, such that the maximum increase was 3.77 × 10⁵ el/cm³ (~64%) at 355 km altitude. Generally, the tsunami induced an increase in the IED peak of value and altitude. Occhipinti et al. showed that about one hour after the occurrence of the tsunami, the main part of the IGWs energy reaches the altitude of 300 km. At this height, the electron density value becomes significant. This effect in the upward velocity adjusts the tsunami waveform as it propagates from the surface of the ocean surface to high altitudes. The ionosphere responds instantaneously to the IGW forcing, and produces a passing wave that disappears with the diffusion and chemical loss after some time. On the contrary, ion production and loss plays a crucial role; the signatures of the IGWs in the plasma are maximized in the direction of the magnetic field [31]. In addition, Occhipinti et al. demonstrated that when a tsunami occurs in the Northern Hemisphere, the perturbed electron density would not exceed 10% in both the E and F regions; however, when it travels south, the electron density could reach to up 80% [31].

The ionospheric disturbances appeared at 15:10 UT and faded at about 20:00 UT. The results of the ionosonde and GPS observations were in good agreement. The first time the tide gauges detected tsunami waves was at 15:34 UT, about 25 min after ionospheric perturbations appeared in the ionosphere. So far, various researchers have studied the arrival time of IGWs and tsunami waves, and have concluded that IGWs can be used to warn of tsunami waves [15,27]. We used the near real-time method (second-order numerical difference method) to detrend VTEC and water level measurements. As IGWs propagate obliquely in the atmosphere, with both horizontal and vertical propagation velocity components, and the vertical propagation velocity of the IGWs is high, the IGWs in the ionosphere can be detected earlier than the tsunami waves in the ocean. The vertical propagation velocity of the IGWs is about 50 m/s, which increases approximately with height due to temperature variations [4,32]. The Kp index was more than five during the tsunami. Liu et al. [33] also suggested that the magnetic storm had little effect on the ionospheric disturbances on 11 March 2011.

6. Conclusions

In this study we investigated the effects of IGWs developed by the tsunami on the bottom-side ionosphere. As the IGWs originate from the ocean, the tsunamis propagate in a conical shape through the atmosphere, and their velocity has both vertical and horizontal components. The IGWs horizontal velocity component and the tsunami’s velocity are more or less the same, but the vertical component of the IGWs velocity increases rapidly with the increase in the altitude, therefore the effect of the IGWs propagation through the bottom-side ionosphere, which is manifested as ionospheric irregularities, can be detected by the GPS ground stations in the region of the tsunami earlier than the tsunami itself.

In this investigation we studied the 2011 Tohoku tsunami and its effects on the American west coast, which is located 7900 km away from the tsunami’s epicenter. The earthquake occurred at 5:46 UT, and was detected at 15:34 UT on the American west coast for the first time. However, the first footprints of IGWs, observed at one of the GPS ground stations located in the same region, were detected at 15:10 UT, which is about 25 min earlier than the tsunami’s first detection.

This proves the fact that the effect of IGWs originated from the tsunamis on the bottom-side ionosphere can be used as an alternative approach for a tsunami’s early warning.