Abstract
There is a multitude of wave-like phenomena in Earth’s ionosphere and thermosphere such as acoustic waves, gravity waves, planetary waves, tides, and Traveling Ionospheric Disturbances (TIDs) which are the ionospheric manifestation of atmospheric waves. These phenomena are often difficult to study since measurements are typically irregular in time and space due to geographic constraints for deploying ground instruments and the natural orbital motion of satellites. This frequently precludes Fourier methods such as the Fast Fourier Transform (FFT) from being used. The Lomb-Scargle Periodogram (LSP) provides FFT-like analysis when measurements are irregular. To our knowledge, all prior use of the LSP in space science has been one-dimensional. This paper uses a N-Dimensional extension of the LSP (ND LSP) to study traveling ionospheric disturbances in four dimensions on a quiescent day near solar minimum. We use an exquisite dataset consisting of 12 ionosondes over Australia on June 29, 2019. The ND LSP resolves the full 3-dimensional wave vector as well as the period for many discrete TIDs. To the degree possible, we validate our findings from ionosonde data processed with the ND LSP by using an FFT-based method on line-of-sight TEC data from the same period and find similar wavelengths and periods for the large TIDs. We show that TIDs occur preferentially near elevation and could be missed or mischaracterized if using TEC data in the thin-shell approximation.
1 Introduction
The ionosphere is a charged region of earth’s atmosphere extending from 100 km to 1,000 km in altitude. The ionosphere has sources of variability from both above, such as the sun and magnetosphere, and below such as the lower atmosphere. One class of variability is a Traveling Ionospheric Disturbance (TID) which is a change in the electron density that moves through space in a quasi-periodic way. TIDs are generally caused by waves in the neutral atmosphere. If the restoring force is gravity, the wave is called an Atmospheric Gravity Wave (AGW). If the restoring force is pressure, the wave is called an Acoustic Wave (AW). AGWs have longer periods than AWs, with the shortest AGW period given by the Brunt-Vaisala frequency (15 min at ionospheric altitudes), and the longest period AW given by the acoustic cutoff frequency ( 12 min).
The TIDs that “follow” AGWs can be broken down into Small-Scale TIDs (SSTIDs) with horizontal wavelengths shorter than 100 km, Medium-Scale TIDs (MSTIDS) with horizontal wavelengths in the range 100–300 km, and Large-Scale TIDs (LSTIDs) with wavelengths over 300 km (e.g., Harris et al., 2012). LSTIDs can be generated at high latitudes and travel equatorward and westward with periods of 20 min to multiple hours. MSTIDs have less uniform propagation directions, and have periods in the range of 10–60 min.
TIDs can be launched by a variety of mechanisms. Joule heating in Auroral storms is known to generate LSTIDs. Tropospheric phenomena such as large convective thunderstorms Azeem et al. (2018), Azeem et al. (2015), Vadas and Azeem (2021) can generate AGWs and TIDs. Orographic forcing is a common creator of TIDs Becker and Vadas (2018), Becker et al. (2022). Earthquakes Sanchez et al. (2022), tornadoes Nishioka et al. (2013), and Azeem et al. (2017), Crowley et al. (2016) have also been shown to generate TIDs. TIDs have also been observed following explosions both above and below ground Huang et al. (2019) and references therein, and rocket launches Li et al. (2022), Chou et al. (2018), Mabie and Bullett (2019), Mabie et al. (2016). Electrified MSTIDs are launched by the Perkins instability rather than neutral atmosphere processes. Electrified MSTIDs are distinct from other MSTIDs because they occur at nighttime, carry polarization electric field, are northwest-southeast aligned (in the northern hemisphere) and preferentially propagate southwestward (in the northern hemisphere) Makela and Otsuka (2012), Otsuka et al. (2021).
TIDs are studied using a variety of datatypes. For example, Martinis et al. (2011) use All Sky Imagery (ASI) to study hemispheric asymmetry of TIDs. Incoherent Scatter Radar data has been used to study TIDs in Nicolls et al. (2004), Goodwin and Perry (2022). High Frequency (HF) Doppler radars such as the TID Detector, Built In Texas (TIDDBIT) Crowley and Rodrigues (2012); Vadas and Crowley (2010) have been used to study TIDs. Additionally, the Super Dual Auroral Radar Network (SuperDARN) has also been used to study TIDs Oinats et al. (2015), Frissell et al. (2014). Ionosonde data has also been used to study TIDs Morgan et al. (1978), Tedd and Morgan (1985), Reinisch et al. (2018), Emmons et al. (2020). Radio Telescopes have been used Loi et al. (2016) often in collaboration with ionosondes Obenberger et al. (2019). By far, the most common data type with which to study TIDs is GNSS Total Electron Content (TEC) data Azeem et al. (2018), Azeem et al. (2015), Vadas and Azeem (2021), Sanchez et al. (2022), Crowley et al. (2016), Li et al. (2022), Chou et al. (2018). This data type has many advantages—it is widely distributed over North America and Europe which provides excellent horizontal resolution, it often has a very high cadence, and it is publicly available. However, TEC is an integrated measurement which limits the investigation of vertical structure without tomography.
The numerical techniques used to analyze these datatypes include: interpolating the data onto a regular spatial/temporal grid and using Fourier methods e.g., Azeem et al. (2017), manuallly curve fitting the raw data, e.g., Li et al. (2022) or using the Lomb-Scargle Periodogram (LSP) which provides Fourier-like analysis without the need to interpolate the data. This has been done for studies of longer-period phenomena like Intraseasonal Oscillations in the Mesosphere/Lower Thermosphere Gong et al. (2022), and thermospheric tides Gong et al. (2013), Gong et al. (2021). The LSP has also been used in ionospheric studies. For example, Goodwin and Perry (2022) used the LSP to study small-scale structures in high-latitude electron density from spatially irregular ISR data.
This paper builds off these prior studies and uses a multivariate (N-Dimensional) extension of the LSP (ND LSP) to study TIDs using ionosonde data. This technique resolves the full four-dimensional structure of the TIDs: 3-dimensional vector (latitudinal, longitudinal, and vertical wavelength) and period.
We compare these results to FFT analysis of interpolated TEC. This comparison is not apples to apples because it uses different data (ionosondes and GNSS TEC) and different processing (ND LSP and interpolation and FFT). Although we do not expect these two analyses to agree completely, large perturbations in the ionosphere should be reflected in both analyses.
2 Data
Our work uses both ionosonde and TEC data from June 29, 2019 in Australia. This day is quite calm with staying below 3, the sunspot number being 5, the DST index varying between −10 and +10 nT, and F10.7 being near 70 SFU. The ionosonde data was collected during a dedicated experimental campaign, and the TEC data are from the Madrigal daily line of sight TEC product cos (Author Anonymous, 2019). Each dataset is described separately below.
2.1 Ionosonde data
The ionosonde dataset contains 1,097 ionograms collected from 12 separate ionosondes on June 29, 2019 over Australia. These are mostly Lowell Digisondes (DPS-1), with the exception of Curtin VIS, Alice Springs QVIS, and Humpty Doo QVIS, which are DSTG PRIME ionosondes Harris et al. (2016). Ionosondes measure power as a function of delay and frequency. Then a “trace” of delay as a function of frequency is extracted, and scaled to true height. Trace extraction and scaling can be performed with auto-scaling software e.g., Galkin and Reinsch (2008) or hand-scaling by a human expert. Trace extraction and scaling contribute the majority of the error to ionosonde data, so we use hand-scaling for this dataset. The locations of the ionosondes and their names are shown in panel (a) of Figure 1. The size of the colored dot indicates the total number of measurements provided by that ionosonde. The times of each profile collected by each ionosonde are shown in panel (b). The nominal cadence is 15 min, but there are gaps depending on the ionosonde. For example, Alice Springs and South Hedland are very reliable, but data from Ajana is not available for the second half of the day.
FIGURE 1
For each profile, we discard measurements that are not in the frequency or height range that is directly measured by the ionosonde. For example, even in these hand-scaled profiles, there is data extending above hmF2 which is model-driven. To remove points like this, we identify all measurements with electron densities lower than the maximum electron density at lower altitudes and delete them. This removes the topside and valley region. We also remove measurements below the minimum frequency actively measured by the ionosonde. This removes most low-altitude measurements during local night when the E region disappears.
Our goal with this work is to capture deviations from a smooth background. We use a Thermosphere-Ionosphere-Electrodynamics General Circulation Model (TIE-GCM) Richmond et al. (1992); Roble et al. (1977) run for the same day performed by the Community Coordinated Modeling Center (CCMC) as a background model for its accessibility and reproducibility. This background model has a native 2.5° resolution and 20 min timestep. We interpolate the electron density from the background model to the precise latitude, longitude, altitude, and time of each measurement. This is shown in Figure 1C which shows the electron density as a function of time for both the Alice Springs ionosonde in blue and the TIE-GCM electron density interpolated to the same location and time. We use measurements near 205 km here, knowing that the actual measurement altitude varies with each profile. Note the large gap in measurements at local night, and contrast the smoothness of the TIE-GCM electron density with the variability of the ionosonde measurements. To quantify the deviations from the smooth TIE-GCM background we use the noise parameter defined in Equation 1 to follow Hughes et al. (2022):where is the ionosonde measurement, and is the smooth background, both in [#/]. Since a log axis is used for Figure 1C, is the distance between the blue and orange lines. For this ionosonde and altitude, is positive all day, but when considering all ionosondes and altitudes is well-modeled by a Gaussian distribution with and . The sign of indicates that the real measurements are on average higher than TIE-GCM by a factor of on this day and in this region.
In summary, our dataset contains 1,097 ionograms from 12 unique ionosondes. Each ionogram contains an electron density profile with many measurements of electron density. Each measurement in the region directly measured by the ionosonde is compared to the electron density from a TIE-GCM model run to compute . In total, this results in 14,039 measurements of all at different latitudes, longitudes, altitudes, and times.
2.2 TEC data
The TEC data is the line-of-sight TEC data product publicly available from the Madrigal Database cos (Author Anonymous, 2019). This global data product is derived from GPS and GLONASS observations from more than 6,000 GNSS receivers, including about 60 in Australia that are used for this study. The line-of-sight TEC values along the raypaths from each receiver to the multiple GNSS satellites in view are reported at a 30 s cadence. The data product also includes the vertical TEC for each line-of-sight TEC value, calculated using a mapping function and placing the ionospheric piercing point (IPP) at 350 km altitude. For each station in the vicinity of Australia, we store the transmitting satellite, time, vertical TEC, and location of the IPP.
The vertical TEC data are further processed in a manner similar to Azeem et al. (2017). This study detrended the TEC by subtracting a running 40 min mean to yield dTEC. This detrending would remove the long-period TIDs which we wish to study in this analysis, so we do not detrend our TEC. Portions of the TEC arc with the GNSS satellite elevation angle below are discarded. Figure 2 shows the density of the Madrigal data.
FIGURE 2
The TEC data are then binned into 15 min in time and in latitude and longitude, within a geographical boundary of latitude and longitude (the red boundary in Figure 2). This boundary is chosen to approximately match the geographical extent of ionosondes shown in Figure 1A. There are some gaps in the binned data due to the irregular distribution of IPP locations. For each time bin, the gaps in latitude and longitude bins are filled using 2D linear interpolation. No extrapolation is necessary for the chosen geographical boundary.
3 Methods
3.1 ND LSP
Many different methods have been used to study periodic structures in the ionosphere. If the data are regularly spaced, Fourier analysis using the Fast Fourier Transform (FFT) provides robust and easily interpretable results. Furthermore, FFT analysis is supported in virtually every scientific computing language making it readily available. If the data are not regularly spaced, they can be interpolated onto a regular grid, but this can affect the spectra depending on the size of the gap being interpolated through Munteanu et al. (2016). Techniques that do not require interpolation include phase folding, least squares, and Bayesian methods. One very powerful member of the least squares family of spectral analysis techniques is the Lomb-Scargle Periodogram (LSP) named for Lomb (1976) and Scargle (1982). This technique is very common in astronomy and astrophysics where measurements are typically irregular. An excellent review of this technique and its use can be found in VanderPlas (2018). The 1D LSP is a well-studied numerical technique Press et al. (1992) and is supported in many common scientific programming languages such as Matlab, Python, and IDL. However, to the author’s knowledge, there is only one multivariate (N-D) implementation Seilmayer et al. (2022), and it is in the R programming language.1 This is in contrast to Fourier analysis which has broad support for N-dimensional analysis. Their derivation is briefly repeated here. We have implemented this in python, with a broadcasting extension. The python package as well as a demo and visualization script are provided on github.
3.1.1 Derivation
Assume that we have a dependent variable that is a function of different independent variables which we group into the vector . For our case, is the electron density, which is a function of the four independent variables of time, altitude, latitude, and longitude, so . Further, assume that we have samples of our independent and dependent variables so that is of shape and is of shape [, ). Assuming that our signal is periodic, we are interested in finding the coefficients and in a spectral decomposition defined in Equation 2:where each is a size vector containing a unique set of frequencies. For our case, they are all of the form for periods and wavelengths . The coefficients , describe the amplitude and phase of at frequency . Equation 2 differs from the standard discrete Fourier decomposition in two ways; first, the frequencies are not a priori fixed; and second, the parameters are frequency dependent. These differences cause the basis functions in the decomposition of Equation 2 to be not necessarily orthogonal. Hence, the calculation of such a decomposition requires a method to choose the parameters and , preferably such that each of , , , and are pairwise orthogonal for each . Therefore, instead of computing the full decomposition, we calculate the parameter such that the two functions and are orthogonal on our grid of observations, and subsequently project our observations onto these functions for a dense grid of choices for .
We ensure orthogonality by first defining a new parameter which we note is a scalar for each . We must find the value of such that the sine and cosine term are orthogonal for our sampling, This leads to Equation 3:
Trigonometric identities are applied to render this equation only in sines and cosines of either or . Our goal is to isolate and solve for . This results in Equation 4
After some algebraic manipulation, one can separate the in Equation 5:Where this sum is performed across each of the samples. Now that we have the value of that makes the sine and cosine terms orthogonal on our sampling domain, we can find the and parameters by taking the dot product with our dependent variable This is done in Equation 6:
We can now find the amplitude and phase of at the frequency This is done in Equation 7:where the signs of must be used to determine the quadrant of the phase.
3.1.2 Frequency selection
Fourier and LSP analysis both assume that a wave is coherent over the domain. Since TIDs have wavelengths and periods that change with altitude e.g., Nicolls et al. (2004) the data are limited to altitudes above 150 km to minimize the effects of differing vertical wavelengths. TIDs change direction with time. For example, Crowley and Rodrigues (2012) showed that the prevailing TID direction is offset from the prevailing neutral wind direction which rotates through per day. If we used the full day of data simultaneously, TIDs of different azimuths might cancel out and not appear in the analysis. To minimize these types of errors, we only use data from 15 to 21 UT (1-7 LT) to focus on a particularly active time when the TIDs can be assumed to be more coherent. This time range also means that all measurements considered occur during local night. We apply the LSP to a subset of the ionosonde dataset by providing one-dimensional arrays of , , altitude, time, and the noise parameter . The variables are the North and East distance (in km) from a central point at (−, ) in the center of Australia.
In Fourier analysis, the observable frequencies are given by Nyquist criterion. For LSP analysis, the frequencies must be chosen by the user. In general, uneven sampling allows for much higher frequencies to be resolved than any “pseudo-Nyquist” limits would indicate VanderPlas (2018). In fact, Eyer and Bartholdi (1999) shows that the theoretical maximum frequency is set by the number of significant digits in the data rather than the spacing. Instead of using these pseudo-Nyquist frequencies, we choose frequencies based on the TIDs that we hope to observe, and then use a shuffling approach to determine observability and significance. The periods for each dimension are shown in Table 1. Negative values for the spatial dimensions indicate motion in the opposite direction. For example, a −300 km northward period indicates a wave traveling south with a 300 km wavelength.
TABLE 1
| Temporal period [hrs] | , 6, 3, 2.88, 2.75, 2.62, 2.50, 2.37, 2.25, 2.12, 2, 1.88, 1.75, 1.62, 1.50, 1.37, 1.25, 1.20, 1.12, 1, 0.86, 0.75, 0.67, 0.6, 0.55, 0.5 |
| Vertical Period [km] | −100, −200, −300, −400, −500, −600, −700, −800, −900, −1,000, , 1,000, 900, 800, 700, 600, 500, 400, 300, 200, 100 |
| Eastward Period [km] | −300, −600, −900, −1,200, −1,500, −1800, −2,100, −2,400, −2,700, −3,000, , 3,000, 2,700, 2,400, 2,100, 1800, 1,500, 1,200, 900, 600, 300 |
| Northward Period [km] | −300, −600, −900, −1,200, −1,500, , 1,500, 1,200, 900, 600, 300 |
Periods probed by LSP analysis.
We probe up to 1,000 km vertical wavelengths, despite only having measurements from 150 to 350 km. These long wavelengths would not be resolvable in only one dimension, as the vertical range is too short to see an entire period. However, we are able to discern these waves if they have a resolvable period in another dimension. For example, in van de Kamp et al. (2014), one can visually extrapolate where the constant phase contour would intersect the axis in their Figure 2 or 8.
To assess the significance of the amplitude, we empirically calculate the noise amplitude by shuffling the variable and re-running the ND LSP. We do this 10 times and average the two highest amplitudes for every combination of frequencies. This estimates the 95% significance threshold—amplitudes higher than this from the unshuffled data have less than a 5% probability of occurring simply due to the structure of the frequency grid and likely have physical significance. We also refer to the 95% significance threshold as the “noise amplitude”.
3.2 FFT TEC analysis
FFT analysis is applied to the binned dTEC data to identify periodic structures. As previously mentioned, the observable frequencies in Fourier analysis are set by the Nyquist criterion that depends on the sampling rate. The grid for binning dTEC was chosen to make the largest observable frequencies (equivalently, shortest observable period or wavelength) similar to those probed by the ND LSP of the ionosonde data. For the FFT analysis, the shortest period is 30 min and the shortest horizontal wavelength is 3° (approximately 330 km). Of course, the FFT analysis of the dTEC data cannot provide any information about the vertical frequencies that the ND LSP can reveal in the ionosonde data.
The amplitudes at frequencies , , and are determined by applying a 3D FFT (scipy.fft.fftn) to the binned dTEC data. The spatial frequencies in latitude and longitude in units of [1/degree] are converted to units of [1/km] using approximations and . Our procedure for determining the 95% significance threshold is similar to the ND LSP. The 3D FFT is applied to 100 shuffled versions of the binned dTEC data and the 95th percentile for each point in frequency space is the significance threshold. We note that the shuffling is applied amongst the binned dTEC data prior to the 2D linear interpolation step, leaving the gaps in the same place. This ensures that the locations and times of the gaps (which occur more often in the interior of Australia where there are fewer IPP tracks) and their impact on the FFT results are properly accounted for in the 95% significance threshold.
4 Analysis
The spectra of the ionosonde and TEC data are four and three-dimensional, respectively, which makes visualization difficult. However, higher-dimensional datasets can be visualized by taking planar slices through a central point. For example, a 3-dimensional dataset can be visualized by taking three orthogonal planar slices through a central point as . We use this approach to visualize the four-dimensional amplitude of the ionosonde data as a function of wavelength and period. While there are three orthogonal slices through a 3D dataset, there are six orthogonal planes that slice the four-dimensional amplitude array for the ionosonde data. In Figure 3, we show these six orthogonal slices for a central point corresponding to a TID moving northwest with a wavelength of km and a period of 2.4 h. This wave has the highest amplitude of any of the significant waves with an amplitude more than twice the noise amplitude at this frequency. The amplitude of this wave is 0.12, which corresponds to variations in electron density on the order of compared to the smooth background.
FIGURE 3
Each panel shows a slice through a different plane. The central point is identified by the pair of thin white lines forming a crosshair in each panel. Panel (g) shows a map with the wave at this central point overlaid as color. The arrow shows the distance the wave travels in 1 h. This central point is a clear peak in panels (a–f), although it is more prominent in some dimensions than others. For example, panel (a) shows only two peaks, the central point being much higher than the secondary one while panel (e) shows a rich sea of other local peaks corresponding to horizontal waves with the same period but other azimuths and wavelengths.
Next, we apply a two-part threshold to identify significant waves and characterize them statistically. First, the amplitude must be larger than its neighbors in all dimensions. Points on the edge of the domain (such as temporal periods of 30 min or ) are not included, as one cannot know whether they are indeed a local max. Second, the amplitude must be larger than the empirically determined noise amplitude. This classifies 110 of the 126,126 points (0.087%) as statistically significant waves.
These 110 points are classified as either horizontal or non-horizontal waves. This binary distinction is based on a gap in the measured elevation angles. An elevation angle of 0 means that the TID propagates horizontally, and an elevation angle of 90 means the TID propagates vertically. Defining and , the elevation angle is defined in Equation 8:
Horizontal waves are observable, but the next lowest elevation wave comes from the largest vertical wavelength ( km) and the shortest horizontal wavelength ( km). This gives an elevation angle of arctan (0.212) . Therefore, waves with are classified as horizontal and all others are classified as non-horizontal. The periods, wavelengths, elevations, and travel directions are shown in Figure 4.
FIGURE 4
Panel (a) of Figure 4 shows a scatter plot of the period and wavelength for the non-horizontal waves. The color and size of the dots indicate the amplitude of the waves. A typical amplitude is 0.075 which translates to a amplitude deviation in electron density. Thin black lines on both margins show normalized histograms of the wavelengths and periods measured. These histograms aid interpretation by showing whether the lack of a scatter point indicates that there was no TID with those characteristics in the ionosphere, or whether that wavelength and period were simply not probed. Panel (c) shows a normalized histogram of the elevation angle for all 126,126 points (blue) and the 110 TIDs (orange). The blue histogram shows the gaps between 0 and as well as a preference for angles near . The difference between the shape of the orange and blue distribution shows that TIDs are most common near elevation angles. There is an asymmetry in the distribution with elevation angles being much more common than elevation angles, suggesting that on this day more waves originated from above than from below. Additionally, the lowest elevation angle for TIDs is despite frequencies leading to lower elevation angles being measured which indicates a real absence rather than a lack of measurement.
There are many more non-horizontal waves than horizontal ones. The period and wavelength of the eight horizontal waves are shown in Figure 4B. The black dashed line shows a speed of 300 m/s. Points to the left of this line have speeds 300 m/s, and points to the right have speeds 300 m/s. Figure 4D shows the speed as a function of azimuth for the horizontal waves.
Our instrument geometry and cadence as well as the frequencies probed preclude observation of some waves. We have previously discussed the limitations in elevation angle, but our 15 min measurement cadence and corresponding 30 min minimum resolvable period prevents all acoustic waves and many small-scale and medium-scale TIDs from being measured. We also lack the ability to resolve short horizontal wavelength waves. Even with these limitations, we can still comment on the amount of variance that is explained by waves in the probed frequency domain. We do this by summing the amplitudes of all 110 of the significant waves, and comparing this to the sum of the amplitudes at all probed frequencies. This indicates that waves account for approximately 0.157% of the total deviation from a smooth background.
Now we compare the results of the ionosonde ND LSP analysis to the TEC analysis for the horizontal waves. Figure 5 shows a scatter plot of the waves found using the TEC data as large colored circles. The size and color of the circle indicate the amplitude of the waves in TECU. Thin black lines in the margins show the measured wavelengths and periods in the TEC data. The black dashed line shows a speed of 300 m/s, just as in Figure 4. Additionally, this plot shows red x markers for the eight significant horizontal TIDs from the LSP analysis of the ionosonde data.
FIGURE 5
Because the measured periods and wavelengths are slightly different, we do not see exact matches for the waves found using each dataset/method. However, the largest TID in both analyses is very similar, having a period and wavelength of (2.4 h, 1,073 km) in the ND LSP analysis and (2 h, 1,095 km) in the TEC FFT analysis. Because the FFT analysis does not oversample the frequencies, it is possible that the peak at 2 h is the shoulder of a peak that is truly closer to 2.4 h. There are two TIDs in the TEC analysis near (3 h, 600 km) that do not have a match in the ionosonde ND LSP analysis. However, there is another pair near (2 h, 600 km) that do match a single ionosonde ND LSP TID. There is a single TEC TID near (1 h period, 800 km) that is near an ND LSP TID with a slightly lower period and slightly longer wavelength. There are two clusters of ND LSP TIDs near (45 min, 1,200 km) and (1 h, 550 km) that do not have corresponding TEC TIDs. In total, there are three close matches (2 h, 1,100 km), (2 h, 600 km), (1 h, 800 km) of which the first is the strongest TID in either analysis. There are two clusters in the ND LSP analysis (45 min, 1,200 km), (1 h, 600 km) with no match from the TEC analysis. There is one cluster in the TEC analysis (3 h, 600 km) that has no match on the ND LSP ionosonde analysis.
Perhaps the most likely reason for these discrepancies is that we are measuring the TEC signatures of non-horizontal TIDs. A full mapping of amplitude and elevation angle to the TEC signature would require topside densities which are not measured in the ionosonde data used here. However, we can say that near horizontal TIDs with the vertical wavelength greater than the “thickness” of the ionosphere (500 km) which contributes significant TEC, will impact vertical TEC measurements. More specifically, a vertical ray has no chance of cancellation by traveling through a peak and trough. Smaller values may also affect vertical TEC measurements, but cancellation will muddle their interpretation. We use a minimum TEC elevation angle of in this work which could lead to complex aliasing of these non-horizontal TIDs and could explain why there are TIDs present in the TEC analysis absent from the ionosonde ND LSP horizontal analysis.
5 Discussion
Comparing the ND LSP results to previous observations, we see general agreement in vertical TID characteristics. Although most TID observations focus only on horizontal characteristics, ISR observations have previously been used to measure the vertical structure of TIDs. Shibata and Schlegel (1993) used EISCAT Tromso to study vertical and horizontal TID propagation. During a geomagnetically quiet period (9/7/1988), they found AGWs they observed AGWs with an average elevation angle of , consistent with the TID parameters observed in this study. Additionally, Chum et al. (2021) found a maximum elevation angle near ( in their coordinate system) for results collected over a year at solar min and solar max. The azimuthal angle of propagation was also similar to that of the largest amplitude waves found in the current study. Rice et al. (1988) examined TIDs using multiple ISRs during the WAGS campaign (10/18/1985, moderate geomagnetic activity) and observed TIDs with elevation angles 75.5–82 .
These results also show overall agreement with previous studies on TID propagation and horizontal characteristics. The observed horizontal TIDs have speeds from 100–400 m/s. These speeds are slower than average LSTID speeds Hocke and Schlegel (1996), Hajkowicz (1991), but are within the expected range of speeds. The largest observed TID in this analysis had a period of 2.4 h. Similar TIDS with a dominant mode of 2.5–3.5 h have been observed in the southern hemisphere Habarulema et al. (2013).
The LSP technique is novel in that it can examine both horizontal and vertical propagation of TIDs using ionosonde data. Most observational studies have only examined horizontal propagation. Characterizing both vertical and horizontal propagation gives a more accurate picture of TID dynamics.
6 Conclusions and future work
In this paper, we have used a N-dimensional Lomb Scargle Periodogram and demonstrated its utility on a dataset consisting of high-cadence ionograms with minimal processing and no interpolation. The resulting analysis is compared to a concurrent analysis of line of sight TEC data. The two analyses agree on the horizontally-propagating TIDs which are resolvable by both methods. The ionosonde LSP analysis additionally resolves a spectrum of propagation elevation angles. The ND LSP is a complementary technique to other spectral analysis methods that can resolve periodic structures from unstructured data. Although shown here for ionosonde data, the ND LSP could be applied to TEC data as well.
We show for this dataset that TIDs are more common at than other elevation angles. Although this population of near-vertical TIDs has been observed before, it has received much less study than its horizontal counterpart. This population can easily be mischaracterized or completely lost if using TEC data with the thin-shell approximation.
Future work has many possible avenues. First, this technique could be used with ionosonde data for active periods with known TIDs such as tsunamis, large convective storms, and geomagnetically active periods where the aurora launches TIDs equatorward. Secondly, this technique could be used with dispersed space in-situ data such as mIVM measurements from COSMIC-2 and other missions.
7 Appendix
The code for computing the ND LSP on general data is provided here https://github.com/joe-hughes26/ND-LSP along with a demo script and visualization.
Statements
Data availability statement
Publicly available datasets were analyzed in this study. Simulation results for the TIE-GCM run have been provided by the Community Coordinated Modeling Center (CCMC) at Goddard Space Flight Center through their publicly available simulation services https://ccmc.gsfc.nasa.gov. The specific run is hosted at https://ccmc.gsfc.nasa.gov/results/viewrun.php?domain=IT&runnumber=Joe_Hughes_050423_IT_1.
Author contributions
JH: Conceptualization, Investigation, Methodology, Software, Visualization, Writing–original draft. IC: Data curation, Formal Analysis, Methodology, Software, Validation, Visualization, Writing–original draft. AN: Methodology, Validation, Visualization, Writing–review and editing. RK: Software, Writing–review and editing. WW: Funding acquisition, Project administration, Resources, Supervision, Writing–review and editing. KO: Funding acquisition, Project administration, Resources, Writing–review and editing. RL: Funding acquisition, Project administration, Resources, Supervision, Writing–review and editing. JC: Funding acquisition, Project administration, Resources, Supervision, Writing–review and editing. JM: Funding acquisition, Project administration, Resources, Writing–review and editing.
Funding
The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This work was performed as part of the Impact of Ionospheric Irregularities program at AFRL/RVB under the leadership of KO on contract FA9453-19-C-0400.
Acknowledgments
The authors extend many thanks to Federico Gasperini, Jeff Steward, and Ryan Nguyen for many fruitful discussions on the ND LSP. We thank the Defense Science and Technology Group (DSTG) for providing the ionosonde data. The views expressed are those of the authors and do not reflect the official guidance or position of the United States Government, the Department of Defense or of the United States Air Force. The appearance of external hyperlinks does not constitute endorsement by the United States Department of Defense (DoD) of the linked websites, or the information, products, or services contained therein. The DoD does not exercise any editorial, security, or other control over the information you may find at these locations.
Conflict of interest
Authors JH, IC, AN, RK and WW were employed by Orion Space Solutions.
The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that no Generative AI was used in the creation of this manuscript.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
References
1
Author Anonymous (2019). Data from the cedar madrigal database.
2
AzeemI.VadasS. L.CrowleyG.MakelaJ. J. (2017). Traveling ionospheric disturbances over the United States induced by gravity waves from the 2011 tohoku tsunami and comparison with gravity wave dissipative theory. J. Geophys. Res. Space Phys.122, 3430–3447. 10.1002/2016JA023659
3
AzeemI.WalterscheidR. L.CrowleyG. (2018). Investigation of acoustic waves in the ionosphere generated by a deep convection system using distributed networks of gps receivers and numerical modeling. Geophys. Res. Lett.45, 8014–8021. 10.1029/2018GL078107
4
AzeemI.YueJ.HoffmannL.MillerS. D.StrakaW. C.CrowleyG. (2015). Multisensor profiling of a concentric gravity wave event propagating from the troposphere to the ionosphere. Geophys. Res. Lett.42, 7874–7880. 10.1002/2015GL065903
5
BeckerE.VadasS. L. (2018). Secondary gravity waves in the winter mesosphere: results from a high-resolution global circulation model. J. Geophys. Res. Atmos.123, 2605–2627. 10.1002/2017JD027460
6
BeckerE.VadasS. L.BossertK.HarveyV. L.ZülickeC.HoffmannL. (2022). A high-resolution whole-atmosphere model with resolved gravity waves and specified large-scale dynamics in the troposphere and stratosphere. J. Geophys. Res. Atmos.127, e2021JD035018. 10.1029/2021jd035018
7
ChouM.-Y.LinC. C. H.ShenM.-H.YueJ.HubaJ. D.ChenC.-H. (2018). Ionospheric disturbances triggered by spacex falcon heavy. Geophys. Res. Lett.45, 6334–6342. 10.1029/2018GL078088
8
ChumJ.PodolskaK.RuszJ.BaseJ.TedoradzeN. (2021). Statistical investigation of gravity wave characteristics in the ionosphere. Earth, Planets, Space73, 60. 10.1186/s40623-021-01379-3
9
CrowleyG.AzeemI.ReynoldsA.DulyT. M.McBrideP.WinklerC.et al (2016). Analysis of traveling ionospheric disturbances (tids) in gps tec launched by the 2011 tohoku earthquake. AGU Radio Sci.51, 507–514. 10.1002/2015RS005907
10
CrowleyG.RodriguesF. S. (2012). Characteristics of traveling ionospheric disturbances observed by the tiddbit sounder. Radio Sci.47. 10.1029/2011RS004959
11
EmmonsD. J.DaoE. V.KnipplingK. K.McNamaraL. F.NavaO. A.ObenbergerK. S.et al (2020). Estimating horizontal phase speeds of a traveling ionospheric disturbance from digisonde single site vertical ionograms. Radio Sci.55, e2020RS007089. 10.1029/2020RS007089
12
EyerL.BartholdiP. (1999). Variable stars: which nyquist frequency?Astron. Astrophys. Suppl. Ser.135, 1–3. 10.1051/aas:1999102
13
FrissellN. A.BakerJ.RuohoniemiJ. M.GerrardA. J.MillerE. S.MariniJ. P.et al (2014). Climatology of medium-scale traveling ionospheric disturbances observed by the midlatitude blackstone superdarn radar. J. Geophys. Res. Space Phys.119, 7679–7697. 10.1002/2014JA019870
14
GalkinI.ReinschB. (2008). “The new artist 5 for all digisondes,” in Ionosonde Network advisory group bulletin (Sydney, New South Wales, Austrailia).
15
GongY.LvX.ZhangS.ZhouQ.MaZ. (2021). Climatology and seasonal variation of the thermospheric tides and their response to solar activities over arecibo. J. Atmos. Solar-Terrestrial Phys.215, 105592. 10.1016/j.jastp.2021.105592
16
GongY.XueJ.MaZ.ZhangS.ZhouQ.HuangC.et al (2022). Observations of a strong intraseasonal oscillation in the mlt region during the 2015/2016 winter over mohe, China. J. Geophys. Res. Space Phys.127, e2021JA030076. 10.1029/2021JA030076
17
GongY.ZhouQ.ZhangS. (2013). Atmospheric tides in the low-latitude e and f regions and their responses to a sudden stratospheric warming event in january 2010. J. Geophys. Res. Space Phys.118, 7913–7927. 10.1002/2013JA019248
18
GoodwinL. V.PerryG. W. (2022). Resolving the high-latitude ionospheric irregularity spectra using multi-point incoherent scatter radar measurements. Radio Sci.57, e2022RS007475. 10.1029/2022RS007475
19
HabarulemaJ. B.KatamziZ. T.McKinnellL.-A. (2013). Estimating the propagation characteristics of large-scale traveling ionospheric disturbances using ground-based and satellite data. J. Geophys. Res. Space Phys.118, 7768–7782. 10.1002/2013JA018997
20
HajkowiczL. A. (1991). Auroral electrojet effect on the global occurrence pattern of large scale travelling ionospheric disturbances. Planet. Space Sci.39, 1189–1196. 10.1016/0032-0633(91)90170-F
21
HarrisT. J.CerveraM. A.MeehanD. H. (2012). Spice: a program to study small-scale disturbances in the ionosphere. J. Geophys. Res. Space Phys.117. 10.1029/2011JA017438
22
HarrisT. J.QuinnA. D.PederickL. H. (2016). The dst group ionospheric sounder replacement for jorn. Radio Sci.51, 563–572. 10.1002/2015RS005881
23
HockeK.SchlegelK. (1996). A review of atmospheric gravity waves and travelling ionospheric disturbances: 1982-1995. Ann. Geophys.14, 917–940. 10.1007/s00585-996-0917-6
24
HuangC. Y.HelmboldtJ. F.ParkJ.PedersenT. R.WillemannR. (2019). Ionospheric detection of explosive events. Rev. Geophys.57, 78–105. 10.1029/2017RG000594
25
HughesJ.ForsytheV.BlayR.AzeemI.CrowleyG.WilsonJ.et al (2022). On constructing a realistic truth model using ionosonde data for observation system simulation experiments. AGU Radio Sci.57. 10.1029/2022RS007508
26
LiX.DingF.YueX.LiJ. (2022). Short-period concentric traveling ionospheric disturbances excited by the launch of China’s long march 4b rocket detected by 1 hz gnss data. AGU Space Weather20. 10.1029/2021SW003003
27
LoiS. T.CairnsI. H.MurphyT.EricksonP. J.BellM. E.RowlinsonA.et al (2016). Density duct formation in the wake of a travelling ionospheric disturbance: murchison widefield array observations. J. Geophys. Res. Space Phys.121, 1569–1586. 10.1002/2015JA022052
28
LombN. (1976). Least-squares frequency analysis of unequally spaced data. Astrophysics Space Sci.39, 447–462. 10.1007/BF00648343
29
MabieJ.BullettT. (2019). Infrasonic wave-induced variations of ionospheric hf sounding echoes. Radio Sci.54, 876–887. 10.1029/2019RS006826
30
MabieJ.BullettT.MooreP.VieiraG. (2016). Identification of rocket-induced acoustic waves in the ionosphere. Geophys. Res. Lett.43, 11024–11029. 10.1002/2016GL070820
31
MakelaJ. J.OtsukaY. (2012). Overview of nighttime ionospheric instabilities at low- and mid-latitudes: coupling aspects resulting in structuring at the mesoscale. Space Sci. Rev.168, 419–440. 10.1007/s11214-011-9816-6
32
MartinisC.BaumgardnerJ.WrotenJ.MendilloM. (2011). All-sky imaging observations of conjugate medium-scale traveling ionospheric disturbances in the american sector. J. Geophys. Res. Space Phys.116. 10.1029/2010JA016264
33
MorganM. G.CalderónC. H. J.BallardK. A. (1978). Techniques for the study of tid’s with multi-station rapid-run ionosondes. Radio Sci.13, 729–741. 10.1029/RS013i004p00729
34
MunteanuC.NegreaC.EchminM.MursulaK. (2016). Effect of data gaps: comparison of different spectral analysis methods. Ann. Geophys.34, 437–449. 10.5194/angeo-34-437-2016
35
NicollsM. J.KelleyM. C.CosterA. J.GonzálezS. A.MakelaJ. J. (2004). Imaging the structure of a large-scale tid using isr and tec data. Geophys. Res. Lett.31. 10.1029/2004GL019797
36
NishiokaM.TsugawaT.KubotaM.IshiiM. (2013). Concentric waves and short-period oscillations observed in the ionosphere after the 2013 moore ef5 tornado. Geophys. Res. Lett.40, 5581–5586. 10.1002/2013GL057963
37
ObenbergerK. S.DaoE.DowellJ. (2019). Experimenting with frequency-and-angular sounding to characterize traveling ionospheric disturbances using the lwa-sv radio telescope and a dps4d. Radio Sci.54, 181–193. 10.1029/2018RS006690
38
OinatsA.KurkinV.NishitaniN. (2015). Statistical study of medium-scale traveling ionospheric disturbances using superdarn hokkaido ground backscatter data for 2011. Earth Planets Space67, 22. 10.1186/s40623-015-0192-4
39
OtsukaY.ShinboriA.TsugawaT.NishiokaM. (2021). Solar activity dependence of medium-scale traveling ionospheric disturbances using gps receivers in Japan. Earth, Planets Space73, 22. 10.1186/s40623-020-01353-5
40
PressW. H.TeukolskyS. A.VetterlingW. T.FlanneryB. P. (1992). Numerical recipies: the art of scientific computing. Cambridge University Press.
41
ReinischB.GalkinI.BelehakiA.PaznukhovV.HuangX.AltadillD.et al (2018). Pilot ionosonde network for identification of traveling ionospheric disturbances. Radio Sci.53, 365–378. 10.1002/2017RS006263
42
RiceD. D.HunsuckerR. D.LanzerottiL. J.CrowleyG.WilliamsP. J. S.CravenJ. D.et al (1988). An observation of atmospheric gravity wave cause and effect during the october 1985 wags campaign. Radio Sci.23, 919–930. 10.1029/RS023i006p00919
43
RichmondA. D.RidleyE. C.RobleR. G. (1992). A thermosphere/ionosphere general circulation model with coupled electrodynamics. Geophys. Res. Lett.6, 601–604. 10.1029/92gl00401
44
RobleR. G.DickinsonR. E.RidleyE. C. (1977). Seasonal and solar cycle variations of the zonal mean circulation in the thermosphere. J. Geophys. Res.82, 5493–5504. 10.1029/JA082i035p05493
45
SanchezS.KheraniE.AstafyevaE.de PaulaE. (2022). Ionospheric disturbances observed following the ridgecrest earthquake of 4 july 2019 in California, USA. Remote Sens.14, 188. 10.3390/rs14010188
46
ScargleJ. (1982). Studies in astronomical time series analysis ii. statistical aspects of spectral analysis of unevenly spaced data. Astrophysical J.263, 835–853. 10.1086/160554
47
SeilmayerM.GonzalezF. G.WondrakT. (2022). The multivariate extension of the lomb-scargle method. arXiv.
48
ShibataT.SchlegelK. (1993). Vertical structure of agw associated ionospheric fluctuations in the e- and lower f-region observed with eiscat-a case study. J. Atmos. Terr. Phys.55, 739–749. 10.1016/0021-9169(93)90017-S
49
TeddB. L.MorganM. G. (1985). Tid observations at spaced geographic locations. J. Geophys. Res. Space Phys.90, 12307–12319. 10.1029/JA090iA12p12307
50
VadasS. L.AzeemI. (2021). Concentric secondary gravity waves in the thermosphere and ionosphere over the continental United States on march 25–26, 2015 from deep convection. J. Geophys. Res. Space Phys.126. 10.1029/2020JA028275
51
VadasS. L.CrowleyG. (2010). Sources of the traveling ionospheric disturbances observed by the ionospheric tiddbit sounder near wallops island on 30 october 2007. J. Geophys. Res. Space Phys.115. 10.1029/2009JA015053
52
van de KampM.PokhotelovD.KauristieK. (2014). Tid characterised using joint effort of incoherent scatter radar and gps. Ann. Geophys.32, 1511–1532. 10.5194/angeo-32-1511-2014
53
VanderPlasJ. T. (2018). Understanding the lomb–scargle periodogram. Astrophysical J. Suppl. Ser.236, 16. 10.3847/1538-4365/aab766
Summary
Keywords
traveling ionospheric disturbances, ionosondes, GNSS TEC, spectral analysis, Lomb-Scargle Periodogram
Citation
Hughes J, Collett I, Newheart A, Kelly R, Wilson WJ, Obenberger K, Landry R, Colman J and Malins J (2024) N-dimensional Lomb Scargle Periodogram analysis of traveling ionospheric disturbances using ionosonde data. Front. Astron. Space Sci. 11:1519436. doi: 10.3389/fspas.2024.1519436
Received
29 October 2024
Accepted
28 November 2024
Published
24 December 2024
Volume
11 - 2024
Edited by
Michel Blanc, UMR5277 Institut de recherche en astrophysique et planétologie (IRAP), France
Reviewed by
Juha Vierinen, UiT The Arctic University of Norway, Norway
Yun Gong, Wuhan University, China
Updates
Copyright
© 2024 Hughes, Collett, Newheart, Kelly, Wilson, Obenberger, Landry, Colman and Malins.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Joe Hughes, joe.hughes@orionspace.com
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.