typhoons and thunderstorms (Heale et al., 2019;Suzuki et al., 2013;Vadas et al., 2009Vadas et al., , 2012;;Yue et al., 2009). The all-sky airglow imager network can cover a wider spatial range and track the spatiotemporal evolution of GWs in detail, making up for the defect that the field of view of a single airglow imager can be blocked (J. Y. Xu et al., 2015).

The ever-growing availability of Global Navigation Satellite System (GNSS) measurements provides total electron content (TEC) information of the ionosphere over the globe, and great progress has been made in many ionospheric research areas (Komjathy et al., 2016;Perwitasari et al., 2022;Tang et al., 2018;Zakharenkova et al., 2016;Zhang et al., 2017Zhang et al., , 2019Zhang et al., , 2021)). Zakharov and Kunitsyn (2012) investigated the ionospheric disturbances during six hurricane events by using 92 GNSS stations in the Gulf of Mexico, proved the feasibility of GNSS for detecting hurricane-caused ionospheric disturbances. Lay et al. (2015) detected acoustic waves and GWs associated with mid-latitude thunderstorms in ionosphere altitudes (250-350 km) by using GPS-TEC measurements. Recently, beyond analysis of individual sites, two-dimensional disturbance maps of the ionosphere were used to describe the regional TIDs caused by extreme tropospheric weathers such as tornadoes, thunderstorms, and hurricanes/typhoons. Nishioka et al. (2013) detected a clear concentric wave in the ionosphere during an EF5 tornado by using dense TEC observations over the North America and confirmed that the concentric wave was caused by GWs. Based on the GNSS-TEC data from the Contiguous United States (CONUS), the characteristics (Azeem et al., 2015) and propagation parameters (Azeem & Barlage, 2018) of TIDs caused by thunderstorm-induced atmospheric gravity waves (AGWs) were investigated, and the consistency of the TID characteristics with the underlying AGWs excited by the point sources were verified (Azeem, 2021). By using ground-based TEC derived from GNSS networks in Japan and Taiwan, Chou, Lin, Yue, Chang, et al. (2017) captured two obvious TIDs with concentric and northwest-southeast wavefronts, respectively. S. Xu et al. (2019) used comprehensive satellite data to provide a detailed analysis of the propagation of GWs associated with hurricane Matthew from the troposphere to the stratosphere, then to mesosphere, and finally to the ionosphere. Song et al. (2019) observed two medium-scale traveling ionospheric disturbances (MSTIDs) during the two landfalls of Typhoon Chan-hom by using dense GPS-TEC data of China, and the horizontal phase velocities, periods, and horizontal wavelengths of the MSTIDs were calculated. Vadas and Azeem (2021) investigated TIDs caused by secondary GWs in deep convection and estimated the propagation parameters. In these studies, only horizontal propagation parameters of GWs were analyzed based on observations. The vertical propagation of the ionospheric disturbances between the bottomside and the topside ionosphere, however, has not be well established due to the lack of observation. Simultaneous information of the horizontal and vertical propagations of a TID is fundamental to understanding the GW-TID relationship and ion-neutral coupling processes.

Three-dimensional computerized ionospheric tomography (3DCIT) technology can provide critical 3D information of ionospheric electron density variations and has now become a powerful new tool for ionospheric research (Zheng et al., 2016(Zheng et al., , 2018)). Ssessanga et al. (2015) inferred the vertical and horizontal structure of electron density during nighttime MSTIDs by using 3DCIT technology. Mai and Kiang (2009) employed the 3DCIT technology to reconstruct the electron density disturbances during the Sumatra Tsunami. The latter study showed that while tsunami-induced GWs propagated upward from the sea surface, GWs-induced concentric wave-like ionospheric disturbances were identified in the vertical direction. Kong et al. (2018) analyzed the co-seismic ionospheric disturbance of the Nepal earthquake in 2015 by using the 3DCIT technology, and found that the height of the disturbances was determined to be between ∼150 and ∼300 km. Zhai et al. (2021) reconstructed the ionospheric electron density during the Japan Tohoku earthquake in 2011, and provided the detailed description of seismo-traveling ionospheric disturbances (STIDs) at different altitudes. These studies demonstrated that 3DCIT technology was an effective method for studying the propagation of GWs in the ionosphere.

To understand TID vertical propagation and wave characteristics, fine height resolutions are important. Our study attempts to provide critical vertical ionospheric disturbance information along with simultaneous horizontal information during Category 4 hurricane Matthew on October 7, 2016. In the following sections, Section 2 presented the information on hurricane Matthew, the ground-based GNSS-TEC and ionosonde data, as well as solar/geomagnetic condition. The 3DCIT technique was described in Section 3. In Section 4, the frequency range of disturbances was determined by spectrum analysis, and the spatiotemporal features of the concentric traveling ionospheric disturbances (CTIDs) were analyzed in detail by using two-dimensional TEC disturbance maps and 3DCIT reconstruction results. Finally, the main results were discussed and summarized in Section 5 and Section 6, respectively.

Hurricane Matthew formed in the Caribbean Sea at 12:00 Universal Time (UT) on September 28, 2016, and continued to move westward as a tropical storm. It was upgraded to a Category 1 hurricane (on the Saffir-Simpson hurricane wind scale) on September 29 and intensified to a Category 4 hurricane at 18:00 UT on the following day. The wind reached its peak speed at ∼145 knots (kt) at 0:00 UT on October 1. Since then, Matthew made landfall in Haiti, Cuba, Western Grand Bahama Island, and South Carolina, and caused serious economic losses and casualties. It was downgraded to a Category 1 hurricane on October 8, and then gradually weakened and turned into an extratropical low-pressure system.

The track data of hurricane Matthew was obtained from the National Hurricane Center (NHC) (https://www.nhc. noaa.gov). Figure 1 shows the track of Matthew from 12:00 UT on September 28, 2016, to 12:00 UT on October 9. The time interval of each point on the track is 6 hr, and the colors of points represent the wind speeds. This study focused on the time period from 4:00 to 8:00 UT on October 7, 2016, when hurricane Matthew made landfall on western Grand Bahama Island.

The ground-based GNSS-TEC data were generated by Massachusetts Institute of Technology (MIT)'s Haystack Observatory, as part of the Millstone Hill Geospace Facility program to analyze observations of dense networks of worldwide GNSS receivers. The vertical TEC as well as line-of-sight (LOS) TEC products were published at Madrigal database. Over the CONUS, LOS TEC data from about 2000 GNSS stations (black dots in Figure 1) can be used to investigate ionospheric disturbances. For this study, the cutoff elevation angle was set to 20°, sampling interval of GNSS-TEC was 30 s, and the altitude of the ionospheric single layer model (SLM) (i.e., the height of ionospheric pierce points [IPPs]) was set to 350 km. The magenta diamonds in Figure 1 are the locations of ionosondes.

Figure 2 shows the F10.7, Dst, and Kp indices during the period of hurricane Matthew, and the orange dashed box indicates the indexes on the day of this study. The daily F10.7 are stable from October 5 to October 10. Dst indexes are more than -30 nT and Kp values are less than 3 on October 7. The influence of solar and geomagnetic activity on this study can be excluded.

Each slant TEC (STEC) can be expressed as the integral of the electron density along the LOS from a satellite to a receiver. The reconstruction region was divided into voxels ("voxel" is the smallest unit on the division of 3-D ionosphere reconstruction region). In our study, each voxel had a geographical extent of 0.5° × 0.5° (longitude × latitude) and a vertical height of 10 km. The electron density of each voxel was assumed to be the same; therefore, STEC can be expressed as follows (Zhai et al., 2021):

where m and n represent the number of LOS (i.e., the number of STECs) and the total number of voxels in the reconstruction region, respectively; matrix A is the matrix of intercepts of LOS in voxels, which can be determined by the position of LOS and reconstruction region. In addition, x is the column of electron density in voxels, ε indicates the measurement noise.

In 3DCIT technique, STEC is directly used to reconstruct the vertical electron density distribution of the ionosphere. However, the magnitude of TIDs is usually much smaller than the ionospheric background TEC. To extract TIDs, the Butterworth band-pass filter is employed to de-trend the background variation and the resulting de-trended STEC (dSTEC) is then used to reconstruct the 3-D TIDs:

where 𝐴𝐴 𝐴𝐴 𝑇𝑇 𝑇𝑇𝑇𝑇𝑇𝑇 𝑛𝑛×1 represents the electron density of the TIDs in voxels, and the rest of the variables are the same as those in Equation 1. In this study, the simultaneous algebraic reconstruction technique (SART) was used to iteratively obtain the electron density of voxels (Andersen & Kak, 1984).

Six LOS sequences of TEC measurements were selected to analyze the frequency ranges of the ionospheric disturbances in this event. Figure 3 shows the IPPs trajectories and corresponding GNSS stations of these TEC sequences. For the same TEC sequences, the GNSS satellite number, IPPs trajectories and GNSS stations are marked with the same color (e.g., AL62 represents GNSS station, the black G23 arrow represents the IPPs of GPS satellite No. Twenty three tracked by AL62).

Figure 4 shows the de-trending process of the six TEC sequences mentioned above. (a)-(c) are the analysis results of the stations NCRB, ALHC, and GACR, and (d)-(f) are the analysis results of the stations AL62, TN28, and MSIN. The black lines are the original TEC sequences, the spectrograms show the spectrum analysis results of TEC sequences, and the red lines represent the dTEC sequences obtained after de-trending. dTEC is the vertical value of dSTEC.

Combining Figures 3 and4, it can be seen that, with the decrease of the latitude of IPPs, the TEC shows an increasing trend. The de-trending process was used to eliminate the impact of changes in the geographical location of the IPPs on the TEC. Meanwhile, the amplitude of the ionospheric disturbances caused by hurricane-induced GWs was small, so spectrum analysis was used to determine the frequency range of disturbances. As is shown in the middle panel of Figure 4a, amplitudes of the disturbances between ∼0.5 and 1.4 mHz enhanced significantly from 5:00 to 6:00 UT. The peak of the disturbances was ∼0.1-0.2 TECU (1 TECU = 10 16 el/m 2 ). The spectrum analysis of other stations showed similar results. Thus, a fourth-order Butterworth band-pass filter in the above-mentioned frequency range was applied to obtain the dTEC sequences (Chou, Lin, Yue, Tsai, et al., 2017;Rolland et al., 2011). After de-trending by band-pass filter, the influence of the geographical differences of IPPs on the TEC sequences was well removed and higher frequency ionospheric background noise was filtered. The dTEC sequences showed obvious disturbances from 5:00 to 6:00 UT, and the amplitude was ∼0.2 TECU (red lines in Figure 4). These disturbances were periodic with a period of 20-30 min, which may be caused by GWs. The dense GNSS network in CONUS provides a large amount of observations to investigate the temporal and spatial evolutions of TIDs. The 2-D TEC disturbance maps were obtained by using nearest neighbor interpolation after de-trending the TEC data with the band-pass filter. The temporal and spatial resolution of disturbance maps was 30 s and 0.2° × 0.2° (longitude × latitude). Figure 5 shows the 2-D TEC disturbance maps in the region of 20°-50°N in latitude and 75°-105°W in longitude. The location of the hurricane eye at 6:00 UT is marked with a red cross mark. A remarkable CTID segment was detected at ∼5:40-6:10 UT, within 1,000-1,500 km northwest of the hurricane eye. Its wave-like spatial structure centered on the hurricane eye indicated that it was likely caused by the CGWs excited by hurricane Matthew. In this event, small wave-like disturbances were captured at ∼5:30 UT, and the more pronounced positive and negative waves appeared at ∼5:40 UT. The maximum To study the traveling features of the CTIDs, we took the position of the hurricane eye at 6:00 UT as the coordinate origin and recorded the time-distance relationship of the disturbances. Figure 6a gives the ionospheric disturbances at IPPs from observations of satellite G03, and the red arrow indicates the direction of the IPPs movement. As shown in Figure 6b, the horizontal phase velocities of these CTIDs are ∼150.7-183.1 m/s, the periods are ∼20-30 min, and the horizontal wavelengths are ∼251.6-276.1 km. In comparison, the above-mentioned remarkable CTID (Figure 5) had the horizontal phase velocity of ∼153.4 m/s, average period of ∼30 min, and average horizontal wavelength of ∼276.1 km. The CTID period from the time-distance analysis is consistent with that Based on these dSTEC data, the 3-D structure of disturbances was reconstructed using the 3DCIT technique described earlier. The reconstruction region was set as 20°-50°N in latitude, 75°-105°W in longitude, and 100-1,000 km in altitude. The horizontal resolution was 0.5° × 0.5° (longitude × latitude), the altitude interval was 10 km and the LOS cutoff elevation angle was 15°.

Figure 7 gives the slices of reconstructed electron density disturbances from 100 to 400 km between 05:10 UT and 06:30 UT. The results showed obvious wave-like ionospheric disturbances, especially at altitudes from 200 to 300 km. At 5:30 UT, both positive and negative wavefronts of the disturbances were detected, but the features were dominated by positive wavefronts. At ∼5:50 UT, the positive and negative wavefronts exhibited comparable magnitudes and some disturbances even reached 400 km altitude. The magnitude of the CTIDs started to decrease at ∼6:10 UT, and the disturbances were dominated by negative wavefronts. At ∼06:30 UT, the wave-like disturbances tended to fade away in the reconstruction region. Compared with 2-D TEC disturbance maps in Figure 5, 3-D reconstruction results of disturbances provided more details at different altitudes.

Figure 8 shows the altitude slices of disturbance reconstruction results with 100 km altitude interval at 05:40, 05:50, and 06:00 UT. The disturbances were most obvious at altitudes from 200 to 300 km, the maximum magnitude of disturbances was ∼0.02 × 10 11 el/m 3 , and the maximum height of upward propagation was ∼400 km. The black dotted line in each figure shows that the spatial structure of the CTID was an inverted (upside-down) cone since the disturbances had a cone shape and propagate upward obliquely. Similar results were also shown in previous simulation (Meng et al., 2018) and observational (Chou, Lin, Yue, Chang, et al., 2017;Zhai et al., 2021) studies.

In addition to the disturbances, the ionosphere background was also reconstructed using 3DCIT for the calculation of percentage disturbances in the following (Figure 11b). Figure 9 compares the ionosphere background reconstruction results based on STEC data and the electron density profiles observed by four ionosondes in the North America at 05:30 UT. It should be noticed that the values above peak height were extrapolated using Chapman function since the ionosondes can only detect the electron density values below ionospheric peak height. The ionosphere background reconstruction results showed consistent peak heights (∼300-350 km) with all ionosonde profiles and the electron density profiles were close to ionosonde observations.

To study the vertical propagation parameters of CTIDs, the time-altitude-disturbance relationships at profile point (88°W, 34°N) were investigated. As shown in Figure 10, the red cross mark (79.7°W, 27.7°N) is the center of the CTIDs since the hurricane eye was located near 79.7°W, 27.7°N between 04:00 and 08:00 UT. The profile point is marked with a purple dot and the black line indicates the direction from the CTIDs center to the profile point.

Figure 11a gives the vertical profile of electron density disturbances at 88°W, 34°N from 4:00 to 8:00 UT. Positive and negative disturbances alternated from 4:30 to 6:30 UT and the vertical (downward) phase velocities were ∼203.7-277.8 m/s. The observed periods of the positive and negative phases from 5:00 to 6:00 UT were ∼26-28 min, and vertical wavelengths were ∼367.1-411.4 km. The disturbances were most remarkable from 150 to 400 km altitudes, but also reached 450-500 km altitudes. The same disturbance was captured earlier at higher altitudes, which was consistent with the spatial structure of the CTIDs in Figure 8. Figure 11b shows the percentage of electron density disturbances (ratio of disturbance to ionospheric reconstruction results) from 4:00 to 8:00 UT. The CGWs caused ∼5% disturbances in electron density. Since the background values of electron density were smaller at 100-250 km altitudes, the percentage of electron density disturbances was significant in this altitude range, which was comparable to the simulation results of Zhao et al. (2020). Moreover, some small-amplitude wave-like disturbances during 6:30 and 8:00 UT were more salient in Figure 11b.

As shown in Figure 12, the time-distance-dTEC relationships at 200 and 300 km altitudes in the direction from the CTIDs center to the profile point (black line in Figure 10) from 4:30 to 6:30 UT were investigated. The horizontal phase velocities and wavelengths of the disturbances were ∼145.1-178.5 m/s and ∼243.8-278.5 km at 200 km altitude, and ∼149.1-181.5 m/s and ∼250.5-283.1 km at 300 km altitude. The horizontal phase velocities of the CTIDs at 300 km altitude were slightly larger than those at 200 km, and the disturbances propagated further away from the hurricane eye. Combined with Figure 11a

Severe convective weathers conditions can excite CGWs with widely varying spatial scales and frequencies (Yue et al., 2009). Some CGWs can reach the mesosphere and lower thermosphere (MLT), and cause air density variations (Azeem et al., 2015;Horinouchi et al., 2002). Secondary GWs can be also excited in these regions due to body forcing associated with the primary GWs (Vadas & Azeem, 2021;Vadas & Becker, 2018;Vadas & , 2009, 2013). When these GWs propagate upward into the ionosphere, they cause oscillations in the ions via neutral-ion collisions, which create TIDs. Yue et al. (2009) observed CGWs near the mesopause that were associated with deep convection from local thunderstorms; however, their scales (horizontal wavelengths of ∼40-80 km at z = ∼90 km) tended to be smaller than what we observe here. This is largely due to the damping of smaller scales by molecular viscosity in the thermosphere, which filters out GWs with smaller horizontal wavelength (λ H ) as a function of altitude (see fig. 9b of Vadas, 2007). Thus, the value of λ H (∼240-280 km) obtained from 3DCIT was comparable to the GWs from individual convective plumes shown in fig. 10a of Vadas and Liu (2013). It turns out that the higher the altitude, the larger the λ H in order that the GWs survive filtering by molecular viscosity. Additionally, the overall updraft of a large number of organized convective plumes can also contribute to the formation of larger-scale CGWs. Hurricane Matthew reached category 4 on October 7, 2016, the prevalent deep convective cells in hurricane could form many GW excitation source regions. CGWs with certain scales can be excited and generate remarkable CTIDs in the ionosphere when the overall upward of organized convective plumes from the same source region hit the tropopause. Some studies detected the propagation parameters of the CTIDs in detail. The horizontal phase velocities, horizontal wavelengths, and periods of CTIDs observed during Typhoon Nepartak were ∼161-200 m/s, ∼160-270 km, and ∼15-22 min, respectively (Chou, Lin, Yue, Chang, et al., 2017). The CTIDs with horizontal phase velocities of ∼143 and ∼268 m/s were detected during typhoon Chan-hom on July 11 and 13, 2015 (Song et al., 2019). The horizontal propagation parameters of the CTIDs obtained using GNSS-TEC in this study are comparable to those in previous studies. Zhao et al. (2020) obtained a similar time-distance-disturbance relationship (Figure 6) in the simulation of tropical cyclone-induced CGWs.

In ionosphere altitudes, CGWs with larger horizontal phase speeds would enable them to propagate higher into the thermosphere before dissipation (Azeem et al., 2017;Vadas, 2007). When CGWs travel upward in the thermosphere, the amplitudes increase rapidly until they reach the dissipation height, and then decay rapidly (Vadas & Liu, 2013). S. Xu et al. (2019) analyzed the same hurricane event in this study using GW dispersion relation and propagation theory, and detected several CTIDs within 100-400 km altitudes. They also inferred that the CGWs in this event traveled up to ∼400-450 km from the convective source of the hurricane. The disturbance reconstruction results in this study are comparable to their conclusions. In addition, some discreet red "blobs" were captured in the slices of reconstructed electron density perturbations (Figure 7). Based on the abundant GNSS observations and the consistency shown in the 2-D perturbation maps, these blobs were likely not caused by discreet slant TEC satellite measurements. Instead, constructive/destructive interference between waves generated by different individual convective plumes within the hurricane may be responsible for their generation.

When the primary CGWs propagate from the excitation source to the ionosphere, the TIDs can appear as concentric rings in the ionosphere if the initial GW amplitudes of the fast, large-λ H GWs are small and if wind filtering is not very important (Vadas et al., 2009). However, most GWs have relatively small phase velocities and large amplitudes, which makes them prone to break and dissipate in the MLT (Alexander et al., 1995;Chun & Kim, 2008). Local body forces from deep convection can accelerate the neutral wind in the direction of primary GWs propagation and excite secondary GWs. Secondary GWs in the thermosphere can further induce CTIDs with partial concentric rings structure (Vadas & Azeem, 2021;Vadas et al., 2018;Vadas & Liu, 2009, 2013). Although the spatial structure and propagation parameters of CTIDs were obtained, the launch altitude (in the lower atmosphere or the thermospheric body force altitude) had a significant effect on the propagation and dissipation of GWs. Due to the lack of observations, we could not determine whether the perturbations were caused by primary or secondary waves here.

In this study, the remarkable CTID segment was detected ∼1,000-1,500 km away from the hurricane eye in the northwest direction. The neutral wind velocity vector component of GWs in the direction of the magnetic field lines is an important driving factor for TIDs. The magnetic declination in the east CONUS is ∼-10°, so the 12 of 15 northwestward propagating CGWs have larger components in the direction of magnetic field and are expected to induce CTIDs with significant amplitudes. The magnetic inclination of the reconstruction region is ∼60°-70° downward. The northwestward propagating CGWs were ∼30° and ∼10° from the magnetic field lines in the horizontal and vertical plane, respectively. Therefore, these CTIDs appeared in the magnetic meridian with appreciable TID amplitudes.

TEC observations have been often used to characterize TIDs although TEC represents only integrated electron content. Our results from this event indicate that TEC horizontal propagation speeds (∼150-183 m/s) are very close to electron density horizontal propagation speeds near 250-300 km (∼149-181 m/s). Figure 13 shows these altitudes were close to the F2 peak heights as measured by several nearby ionosondes.

In this study, the 4-D ionospheric disturbances during hurricane Matthew were investigated. The frequency range of the disturbance was determined by spectrum analysis and a Butterworth band-pass filter was used to de-trend the TEC sequences to detect TIDs. A remarkable CTID segment was captured on October 7, 2016 by high spatial-temporal resolution 2-D TEC disturbance maps. It moved northwestward at a horizontal speed of ∼153.4 m/s during ∼5:40-6:10 UT. The average period was ∼30 min and horizontal wavelength was ∼276.1 km.

The CTID segment was reconstructed by 3DCIT technique involving dSTEC data. The results showed that the wave-like disturbance was first captured at ∼300 km altitude at 5:10 UT and tended to disappear at all altitudes at 6:30 UT. The features of the disturbances were dominated by positive wavefronts at ∼5:30 UT. Then, the positive and negative wavefronts exhibited comparable magnitudes at ∼5:50 UT, and were dominated by negative wavefronts at ∼6:10 UT.

In the reconstruction results from 4:00 to 8:00 UT, the disturbances were most remarkable from ∼200 to ∼300 km altitudes, and the maximum height of upward propagation was about 450 km. The vertical (downward) phase velocities of the CTIDs were ∼203.7-277.8 m/s, the intrinsic periods were ∼26-28 min, and the vertical wavelengths were ∼367.1-411.4 km. The horizontal phase velocities at 200 and 300 km altitudes were ∼145.1-178.5 and ∼149.1-181.5 m/s, respectively. The positive and negative phases of the CTID reconstructed by 3DCIT exhibit an inverted cone structure in the vertical slices and profile results. Such CTID information could be incorporated into GW theories to understand GW propagation properties; however, this task is beyond the scope of the present research. In further, validation of the 3-D structures of CTIDs using other independent observations, such as ground-based ionospheric radars and ionosondes as well as in situ satellite data, should be carried out.