the Creative Commons Attribution 4.0 License.

the Creative Commons Attribution 4.0 License.

# Investigation of gravity waves using measurements from a sodium temperature/wind lidar operated in multi-direction mode

### Alan Z. Liu

A narrow-band sodium lidar provides high temporal and vertical resolution observations of sodium density, atmospheric temperature, and wind that facilitate the investigation of atmospheric waves in the mesosphere and lower thermosphere (80–105 km). In order to retrieve full vector winds, such a lidar is usually configured in a multi-direction observing mode, with laser beams pointing to the zenith and several off-zenith directions. Gravity wave events were observed by such a lidar system from 06:30 to 11:00 UT on 14 January 2002 at Maui, Hawaii (20.7° N, 156.3° W). A novel method based on cross-spectrum was proposed to derive the horizontal wave information from the phase shifts among measurements in different directions. At least two wave packets were identified using this method: one with a period of ∼ 1.6 h, a horizontal wavelength of ∼ 438 km, and propagating toward the southwest; and the other one with a ∼ 3.2 h period, a ∼ 934 km horizontal wavelength, and propagating toward the northwest. The background atmosphere states were also fully measured and all intrinsic wave properties of the wave packets were derived. Dispersion and polarization relations were used to diagnose wave propagation and dissipation. It was revealed that both wave packets propagate through multiple thin evanescent layers and are possibly partially reflected but still get a good portion of energy to penetrate higher altitudes. A sensitivity study demonstrates the capability of this method in detecting medium-scale and medium-frequency gravity waves. With continuous and high-quality measurements from similar lidar systems worldwide, this method can be utilized to detect and study the characteristics of gravity waves of specific spatiotemporal scales.

- Article
(14167 KB) - Full-text XML
- BibTeX
- EndNote

Atmospheric gravity waves are generated when air parcels are perturbed vertically, and gravity or buoyancy acts as the restoring force. They can propagate vertically up to the thermosphere and horizontally over a considerable distance (up to several thousand kilometers). The most common wave sources include convection, orography, and fronts, among others (Fritts and Alexander, 2003, and references therein). The momentum and energy transported by gravity waves dramatically impact the general circulation and thermal structure of the middle and upper atmosphere. Phenomena and processes include, but are not limited to, the cold summer mesopause (Holton, 1982; Siskind et al., 2012), the quasi-biennial oscillation (QBO) in the tropical lower stratosphere (Ern et al., 2014) and the semiannual oscillation (SAO) in the upper stratosphere/lower mesosphere (Ern et al., 2015), instability and turbulent mixing in the atmosphere (Fritts, 1984; Fritts et al., 2013), as well as irregularities and traveling ionospheric disturbances (Fritts and Lund, 2011; Liu and Vadas, 2013), which are all related to gravity waves. Therefore, understanding gravity wave generation, propagation, and breaking has significant impacts on weather and climate applications.

Theoretically, gravity waves are governed by the fluid Euler equations for a set of fundamental variables including pressure *p*, density *ρ*, temperature *T*, and zonal, meridional and vertical winds (*u*, *v*, *w*). Many theoretical and observational studies of gravity waves are based on the linear wave theory, and it is commonly used to quantify the propagation and dissipation characteristics of gravity waves. The linearization of the Euler equations is implemented under different assumptions regarding wave and background atmosphere properties. Except for waves with very large horizontal scales, the effect of earth's rotation is often ignored. Zhou and Morton (2007) derived the Euler equations for a compressible atmosphere with altitude-varying background temperature and wind. Taylor (1931) and Goldstein (1931) derived the 2-D Euler equations with the Boussinesq approximation in a continuous shear flow without temperature variations. These specific Euler equations are referred to as Taylor–Goldstein equations (Nappo, 2012). Fritts and Alexander (2003) derived the Euler equations without wind shear but considered the Coriolis effect. Linearized wave solutions require that the vertical wavenumber *m* be independent of altitude. Strict independence is not likely since background temperature and wind vary with altitude. If the variations are relatively slow within the range of the vertical wavelength, the Wentzel–Kramer–Brillouin (WKB) approximation is applied. A monochromatic gravity wave can be assumed to be a traveling plane wave and represented by a sinusoidal function in the form of

of which *x*, *y*, and *z* are zonal, meridional, and vertical distances in a local Cartesian coordinate, and *k*, *l*, and *m* are zonal, meridional, and vertical wavenumbers. *ω* and *ϕ* are the observed (Eulerian) angular wave frequency and initial phase, *H*_{s} is the atmospheric density scale height, and $\mathrm{1}/\left(\mathrm{2}{H}_{\mathrm{s}}\right)$ corresponds to the wave amplitude growth rate with increasing altitude for upward propagating gravity waves without dissipation.

Many remote-sensing and in situ techniques have been developed to observe atmospheric gravity waves and their influences in past decades. There is an inherent “observation filter” (Alexander, 1998; Gardner and Taylor, 1998; Alexander et al., 2010) in any of these observation instruments such that they can only measure part of the gravity wave spectrum. Single-site ground-based techniques like lidar (Hu et al., 2002; Li et al., 2007; Lu et al., 2009; Cai et al., 2014; Chen et al., 2016), radar (Nastrom and Eaton, 2006; Fritts et al., 2010; Liu et al., 2013), and space-based limb-sounding satellites (Jiang et al., 2005; Alexander et al., 2009) are limited to providing vertical profiles and can only resolve vertical structures of the wave field. Other techniques like nadir sounding satellites (Gong et al., 2012; Alexander and Grimsdell, 2013; Hoffmann et al., 2014) and airglow imaging (Taylor, 1997; Espy et al., 2006; Li et al., 2011; Fritts et al., 2014; Cao and Liu, 2022) can only retrieve the horizontal structures over a certain area. In some cases, the unobserved horizontal or vertical information can be estimated by indirect methods based on the polarization and dispersion relationships (Hu et al., 2002; Lu et al., 2015). For reliable estimates of intrinsic wave parameters and characterization of the dissipation process, it is necessary to observe gravity waves fully in both horizontal and vertical directions, i.e., in 3-D space. Measurements from multiple complementary instruments, such as co-located lidar and the airglow imager, provide good opportunities to resolve gravity waves in both vertical and horizontal directions (Bossert et al., 2014; Lu et al., 2015; Cao et al., 2016). More space-borne limb-sounding observations were also analyzed with special techniques to resolve the gravity wave as fully as possible. Two-dimensional (2-D) limb-sounding measurements along a single satellite track, such as temperature variations from the High-Resolution Dynamics Limb Sounder (HIRDLS), were used to estimate the horizontal wavelength and momentum flux for individual wave events (Alexander et al., 2008). However, there were considerable uncertainties in the derived wave properties because the 2-D sampling measures only the apparent wavelength, which is likely longer than the actual one. This was also illustrated by combining rarely occurring co-locations of two HIRDLS profiles and one COSMIC radio occultation (RO) profile (Alexander et al., 2015), which showed that the 2-D method overestimated the occurrence of long horizontal wavelength waves relative to the 3-D method and underestimated momentum flux. Multiple methods have been proposed to resolve the intrinsic horizontal propagation and momentum flux assuming a dominant wave mode within a spatial (latitude and longitude) and temporal window large enough that it is coherent across more than three RO profiles (Wang and Alexander, 2010; Faber et al., 2013). Several studies have analyzed the horizontal wave information from closely spaced profiles in serendipitous geometries, for example, from COSMIC satellites just after launch and before they have separated into their final orbits. Alexander et al. (2018) used four spatially and temporally close RO profiles to derive intrinsic wave propagation and wavelength. These previous studies demonstrate the challenges and importance of resolving gravity waves as fully as possible.

The measurement profiles retrieved from narrow-band sodium lidar systems were typically used to observe the temporal and vertical variations of the gravity waves. This study provides a prospective solution to partly address the “observation filter” effect of these lidars by resolving often neglected horizontal wave information. We propose a novel method using cross-spectrum to retrieve the slight phase shift among the measurements from different laser beams and derive the horizontal wave information. Temperature and wind measurements on the night of 14 January 2002 from a lidar deployed in the Maui/Mesosphere and Lower Thermosphere (Maui/MALT) campaign are used to demonstrate the proposed method. The paper is organized as follows: Sect. 2 describes the instrumentation and methodology. Section 3 presents the observational results and diagnostic analysis. Section 4 presents a sensitivity study to demonstrate the capability of this method. Finally, the discussions and summary are presented in Sect. 5. Extra figures are included in the Appendix as supporting information.

A narrow-band sodium lidar measures the atmospheric temperature and wind in the mesopause region (80–105 km) based on the thermal broadening and Doppler shift of atomic spectral lines of the sodium atoms. The sodium lidar transmits pulsed laser tuned to the sodium D2a line at 589.158 nm into the sky and the fluorescence scattered photons are collected by optical telescopes. The temperature and line-of-sight (LOS) winds are retrieved based on the well-known shape of the Na atomic spectrum using a three-frequency technique (She and Yu, 1994; Krueger et al., 2015). In order to measure the full 3-D wind vectors, the laser beam is configured to point to multiple directions, at zenith and off-zenith at several cardinal directions. A sodium lidar system operated by the University of Illinois at Urbana-Champaign (UIUC) was deployed at Air Force Maui Optical Station (AMOS) at Maui, HI (20.7° N, 156.4° W) from January 2002 to June 2007 (Liu, 2023). The laser transmitter was fixed, but the laser beam was directed by multiple mirrors to point to five directions: zenith (Z), 30° off zenith to the north (N), south (S), east (E), and west (W). The laser beam has an average output power of 1 W at 30 pps (pulses per second). A steerable astronomical telescope of 3.67 m diameter was coupled with the laser beam. The return photons were collected by the telescope pointing in the same direction as the laser beam. Figure 1a shows the orientation of laser beams in five directions. With a 30° off-zenith angle, there is a ∼ 50 km separation distance between the off-zenith and zenith directions at 90 km altitude. The laser beam was directed to rotate in a Z–N–E–Z–S–W sequence and the photon integration time was 1.5 min at each direction with a 1.7 min cadence (extra ∼ 0.2 min for telescope steering). The resulting measurement intervals at zenith and any off-zenith directions were 5.1 and 10.2 min; however, some irregularities existed because of mechanical issues. In order to accommodate the data processing of filtering and spectral analysis, the raw data of all directions was interpolated into the regular 6 min interval. The targeted waves presented in this study have much longer periods and interpolation would not yield any artifacts. An example of the time series of the temperature perturbations of one target wave without interpolation is shown in Fig. A3. The spatial resolution of measurements is 500 m along the laser beam in all directions. The temperature and winds from off-zenith directions are interpolated to the same altitude grids as the Z direction with a uniform 500 m resolution. The lidar-measured temperature and wind accuracies largely depend on the sodium atom density; thus, they vary with altitudes. With the temporal and spatial resolutions in this study, the accuracies are ∼ 2 K for temperature and ∼ 4 m s^{−1} for horizontal winds near the peak sodium density altitudes (Li et al., 2012; Krueger et al., 2015). The most reliable measurements are mainly within the 85–105 km altitude. The mean sodium atom density and the altitude-varying accuracies for temperature as well as horizontal and vertical winds are given in Fig. A1.

The relationship between LOS winds (*V*_{E}, *V*_{W}, *V*_{N}, *V*_{S}, *V*_{Z}) and zonal, meridional, and vertical winds (*u*_{x}, *v*_{y}, *w*_{z}) in different directions (*x* = E, W, *y* = N, S, *z* = E, W, N, S, Z) is described by Gardner and Liu (2007, Eq. A1) as

where *α* is the off-zenith angle of laser beams. When measurements from this type of lidar were analyzed, the separations of laser beams among different directions were generally ignored. Under the assumption of the atmosphere being nearly homogenous among different laser beams and vertical winds being smaller than horizontal winds, the zonal *u*, meridional *v*, and vertical *w* winds are derived from LOS winds as follows:

Note that the derived zonal wind (*u*_{E} and *u*_{W}), meridional wind (*v*_{N} and *v*_{S}), and vertical wind (*w*_{Z}) are available in different directions and time steps, while the scalar temperature measurements are available in all five directions. The derived zonal, meridional, and vertical winds have relatively large errors. The latter wave analyses are mainly based on temperature measurements, while the winds are only shown for reference.

If a wave defined by Eq. (1) propagates through the five laser beams from a certain direction, the laser beams probe different phases of the wave, as is illustrated by Fig. 1a. There should be some phase shifts in the wave-induced perturbations of different laser beams, as shown by the simulations in Fig. 1b. In theory, phase shifts can be determined from these perturbations and further used to derive horizontal wave information. The key lies in the accurate determination of the phase shift. In this study, we utilize the cross-spectrum to retrieve the phase shifts in the wave-induced perturbations from different directions. Two time series of the same frequency *ω* but different phases *ϕ* are described by

The corresponding frequency spectra are *Y*_{1}(*ω*)=FFT[*y*_{1}(*t*)] and *Y*_{2}(*ω*)=FFT[*y*_{2}(*t*)], and the cross-spectrum between the two time series is ${I}_{\mathrm{12}}\left(\mathit{\omega}\right)={Y}_{\mathrm{1}}\left(\mathit{\omega}\right)\cdot {Y}_{\mathrm{2}}^{\star}\left(\mathit{\omega}\right)$. The argument of the complex cross-spectrum is the phase difference between two waves ($\mathrm{Arg}\left[{I}_{\mathrm{12}}\right(\mathit{\omega}\left)\right]={\mathit{\varphi}}_{\mathrm{1}}-{\mathit{\varphi}}_{\mathrm{2}}$). Here, we treat the target waves as quasi-monochromatic, which might not reflect the realistic dispersive waves existing in the form of wave packets but preserve the main characteristics. In the wave parameter estimations, there might exist ambiguity in the phase differences. Therefore, the wave horizontal wavelength must be larger than the laser beam separation that yields a small phase difference for this method to be effective. Furthermore, the validity of the derived wave speed could be used to eliminate some unphysical results. There is also ambiguity in distinguishing the wave propagation direction. In this study, measurements from the third available direction were utilized to resolve the ambiguity. Considering the uncertainties in the measurements, the phase differences could not be too small to be distinguished. Therefore, the method is immune to the aliased signals of tides as they have much larger scales and yield tiny phase shifts within 100 km. A sensitivity study was carried out to verify the effectiveness of the method in resolving waves of different scales and periods.

This type of sodium lidar usually operates at night, with a typical duration of no longer than 10 h. The resulting resolution ($\mathrm{1}/\mathrm{10}$ h^{−1}) in the spectral domain is relatively coarse. There exist possibilities of spectral leakage and the true peaks falling between two adjacent integral spectral points. We improved the identification of the spectral peaks by applying a non-linear fitting of a parabolic function on the magnitudes of three integer spectral points close to the apparent peak. With the refined method, we hope to acquire the non-integer spectral peaks closer to the actual values and estimate the wave amplitudes *A*, periods *τ* (*T* is reserved to abbreviate “temperature”), and phase shifts *ϕ* more accurately. Given the off-zenith angle *α*, the spatial separation between zenith and any off-zenith directions can be calculated as $\mathrm{\Delta}=h\cdot \mathrm{tan}\left(\mathit{\alpha}\right)$ at altitude *h*. Once the phase differences in the zonal and meridional directions (*ϕ*_{x} and *ϕ*_{y}) are determined by the cross-spectral methods, the horizontal wavenumbers in the zonal and meridional directions are derived from

Then, the full set of horizontal wave parameters of wavelength *λ*_{H}, propagation azimuth angle *θ*, and observed (ground-based) phase speed *c*_{H} can be calculated as

## 3.1 Temperature/wind perturbations and background states

On the night of 14 January 2002, the lidar was operated in five-direction mode from around 06:30 to 11:00 UT. Sodium density, temperature, and winds were continuously observed within the period as the laser was rotated in all five directions. A second-order polynomial fitting was used to remove the mean state and tidal signals. The detrended temperature measurements in different directions are shown in Fig. 2. Abundant perturbations of various periods are identified from measurements in all directions, and distinct downward phase progression is seen in the perturbations, which implies an upward wave propagation. There might exist some tide residues in the perturbations that could change the wave amplitudes, but they will not influence the estimate of the phase shift. The amplitudes of temperature perturbations reach ±10 K and there is a layered structure in the vertical direction. The wave patterns of the perturbations in different directions are similar, so they are likely the same wave packets spreading over a larger area and being captured by the laser beams in different directions. Close inspection of the temperature perturbation peaks in different directions revealed some shifts in time, resulting from the spatial separation of laser beams in different directions. The detrended perturbations of different wind components are shown in Fig. 3, in which similar wave patterns with a downward phase progression can be identified in the zonal and meridional winds, with an amplitude of up to ±20 m s^{−1}. The wave pattern is still clear in the vertical wind perturbation, with an amplitude ±2 m s^{−1}. However, the downward phase progression is less evident due to the uncertainty of vertical winds. For both temperature and winds, the measurements at the top and bottom sides are associated with larger uncertainties and should be interpreted with caution.

To identify the dominant wave modes from the temperature and wind measurements, the frequency spectra were calculated from detrended temperature and wind perturbations at all altitudes. Figure 4 shows the spectra of temperature perturbation in five directions, which all show a similar pattern. The average spectrum of all five directions is shown in the upper right corner. The identified spectral peaks at each altitude using the fitting method were also denoted on the average spectrum. There might be one to several peaks at each altitude. Overall, there exist two prominent peaks; one has a period of about 3.2 h and the other one about 1.6 h. The 3.2 h period component is persistent along with the whole altitude range, and the 1.6 h peak exists mostly below 90 km and reaches a maximum at 90 km. The spectra of wind perturbation are shown in Fig. A2. The spectral peak around 3.2 h is also dominant in the horizontal (zonal and meridional) winds. However, the 1.6 h one is less evident. This is likely due to the spectral leakage of the 3.2 h wave component with much larger magnitudes which overwhelms the 1.6 h period one. The two dominant wave components need to be separated from the mean background states, variations with longer periods, and from each other for further cross-spectral analysis. Cut-off periods/frequencies are determined for desired digital filters based on the mean spectra of temperature and wind perturbations. Two spectral peaks are quite close in the frequency domain, so we used Chebyshev type II filters with flat passband and steep transition to stopband. To filter out the 1.6 h wave component, the cut-off period of a high-pass filter is selected as 2.2 h (0.46 h^{−1}).

To fully understand the propagation condition of waves, the background atmosphere states were analyzed. Figure 5a–c show the background temperature *T*_{0}, zonal wind *u*_{0}, and meridional wind *v*_{0} retrieved by detrending as defined above. The background atmosphere states show clear modulation of tides, indicated by a slow downward phase progression in both temperature and winds. The horizontal winds are quite strong, with a magnitude of ∼ 100 m s^{−1} toward the northeast above 95 km and a magnitude of ∼ 50 m s^{−1} toward the southeast around 85 km, and there is a relatively calm layer around 90 km. Squared buoyancy frequency *N*^{2} is calculated from background temperature *T*_{0} through

where *g* is the gravity acceleration constant and *c*_{p} is the specific heat at constant pressure. Larger values of *N*^{2} indicate a more statically stable atmosphere, while values of negative *N*^{2} imply an unstable atmosphere. The squared buoyancy frequencies *N*^{2} shown in Fig. 5d reveal that the background atmosphere was layered with a convectively stable layer between 90 and 97 km, with larger *N*^{2} values about 6–8 × 10^{−4} s^{−2}. This stable layer gradually moved downward as modulated by tides. Multiple relatively unstable layers existed in between the stable layers around 87, 98, and 106 km. These layers with positive but smaller *N*^{2} can also be unfavorable for wave propagation, as the resulting buoyancy period might be longer than the wave intrinsic period.

The Richardson number *Ri* is commonly used to characterize the dynamical (shear) instability and is calculated through

where S is the vertical shear of the horizontal wind. The atmosphere is considered to be dynamically unstable when $\mathrm{0}<\mathit{Ri}<\mathrm{1}/\mathrm{4}$. A strong horizontal wind shear and a negative vertical temperature gradient make the atmosphere dynamically unstable. As shown in Fig. 5f, the atmosphere is in an overall stable status, with *Ri* approaching $\mathrm{1}/\mathrm{4}$ in a few thin layers near 87, 95, and 102 km where dynamical instability was likely to occur. Large wind shears are the main factors in these unstable areas. The layer near 90 km was relatively stable with a large *N*^{2} and a small wind shear, resulting in larger *Ri*. This layer is favorable for the propagation of the atmospheric wave or allows wave amplitudes to reach larger values. As shown by the spectral analysis, there are two isolated wave components presented at the same time. The tides dominate the background states shown in Fig. 5. However, the longer-period (3.2 h) wave effectively acts as the background for the shorter-period (1.6 h) wave. The perturbation resulting from the longer-period wave might change the temperature gradient, and wind shear leads to a different stability condition for the shorter-period wave. The background states that include the longer-period (3.2 h) wave component are given in Fig. A6. Noticeable differences could be identified in the *Ri* where the relatively unstable area largely expands, especially above 95 km, which might lead to the dissipation of the shorter-period wave toward higher altitudes. However, the stable layer around 90 km also expands to a slightly wider altitude range, which improves the propagation condition for the shorter-period wave. The corresponding background states are used to diagnose each wave component in later sections.

## 3.2 Filtered wave components and wave parameters

After applying the desired filters to the temperature and wind perturbations, the two dominant wave components are isolated. Using the improved spectral peak determination method, the best-estimated periods of the two dominant components are determined to be 3.35 h (0.2986 h^{−1}) and 1.63 h (0.6148 h^{−1}). These two wave packets are referred to as wave #1 and wave #2 in the subsequent analysis. Figure 6 shows the filtered temperature perturbations of wave #1 in all five directions, with *Ri* overlapped by contour lines for some values (0.25 and 0.5). The atmosphere is generally in a stable condition favorable for upward wave propagation, with only a few thin layers with *Ri* less than 0.25, and the layers only persist for a short time. However, the wave perturbations have larger amplitudes (about ±5 K) around these layers with a larger Richard number, making the waves appear like nodal structures. Figure 7 shows the filtered perturbations for wave #2 with one local minimum at 90 km highlighted (white crosses). The slight shift in time is visible by comparing perturbations of different directions, and the wave packet roughly moves from the northeast toward the southwest. The wave pattern and downward phase progression are visible in all directions. The perturbations reach a maximum of ±10 K at about 90 km altitude, confined within both layers above and below with smaller magnitudes of *Ri*. Another maximum can be seen at below 85 km. Even though there are wave patterns at higher altitudes above 100 km, the wave pattern is less consistent in different directions and shows up with varying periods. Multiple thin unstable layers exist in this altitude range, so upward propagation waves might undergo non-linear wave-mean flow interaction, resulting in wave dissipation.

Filtered wind perturbations are shown in Figs. A4 and A5. The signatures of both waves are evident in the horizontal winds but with less evident downward phase progression. The node structure is clear in the vertical direction with at least two maxima at different altitudes. The wave patterns are slightly different in the zonal, meridional, and vertical winds, as the relations are determined by the polarization formulas. Even though the wave pattern is clear, the phase shifts among measurements of different directions are imperceptible in horizontal wind perturbations for both wave packets. The measurement uncertainties of 5–10 m s^{−1} are too large compared with the wave amplitudes of 10–20 m s^{−1}, making the slight phase shift hard to distinguish. In later analysis, only the temperature measurements are used for the cross-spectral method to estimate wave parameters.

Using the proposed method detailed in Sect. 2, we choose the temperature measurements at four off-zenith directions to derive the phase shifts in the zonal and meridional directions, and horizontal wavelength, propagation azimuth, and observed phase speed were then calculated. In Fig. 8, the background wind and wave phase speed at different altitudes are shown in the “wind barb” manner, with the length representing speed and the staff pointing to the direction to which the wind blows and toward which the wave propagates. In Fig. 8b and c, wind barb does not exist at some altitudes because the derived wave phase speed is too large ($c-{u}_{\mathrm{0}}>\mathrm{0.5}{c}_{\mathrm{s}}$) (Zhou and Morton, 2007), invalidating the results. Here, the speed of sound *c*_{s} in the atmosphere is estimated as ${c}_{\mathrm{s}}=\sqrt{\mathit{\gamma}R{T}_{\mathrm{0}}}$, where *γ* is the ratio of specific heat, *R* is the ideal gas constant, and *T*_{0} is the background temperature. The typical speed of sound *c*_{s} is calculated to be around 280 m s^{−1} in the lidar observation altitude range. During the lidar observation period and altitude range, the background wind is mostly toward the east, with strong northeastward winds above 95 km and moderate southeastward winds around 85 km. The calculated propagation azimuth and phase speed are relatively consistent with altitude ranges within blue boxes in Fig. 8b and c, including 85–97 km for wave #1 and below 90 km for wave #2. In addition, both wave packets are estimated to propagate opposite to the background winds, with wave #1 propagating mostly toward W and wave #2 propagating toward W–SW. In ideal conditions, these ground-observed wave parameters are invariant if the wave pattern sustains and does not dissipate. The wave parameters of frequency/period, propagation azimuth, phase speed, and horizontal wavelength, as well as their uncertainties, represented by the mean and standard derivation within the most reliable altitude range (within blue boxes in Fig. 8), are listed in Table 1. Both frequency/period and phase speed are estimated in a ground-based frame. The horizontal wavelengths of wave #1 and wave #2 are estimated to be around 934 and 438 km, both with a ∼ 20 % uncertainty. The propagation azimuth angles are estimated to be 282 and 233° for the two waves, both with a 15–25° uncertainty. These wavelengths and azimuths correspond to phase shifts of −36 and 8° between measurements of E–W and N–S for wave #1, and phase shifts of −65 and −49° for wave #2. The 8° phase shift in the N–S direction for wave #1 is too small and hardly visible in Fig. 6. The observed phase speeds of the two wave packets are quite similar, around 80 m s^{−1}. Both waves propagate through the background atmosphere with varying stability and potentially undergo some dissipation and dispersion, especially at higher altitudes. Therefore, the monochromatic wave assumption is no longer satisfied there, and wave speed and azimuth determined from the phase shift show large fluctuations at these altitudes. The wave pattern shows downward phase progression in the vertical direction, and not a single complete wave cycle is identified within the altitude range. Therefore, the vertical wavelengths were roughly estimated from the phase slope to be larger than 10 km for both wave packets.

## 3.3 Wave diagnosis: dispersion and polarization relations

In linearized gravity wave theory, the dispersion relation links the vertical wavenumber to the horizontal wave parameters and background states. It is often used to diagnose gravity wave propagation, reflection, and ducting. For the acoustic–gravity waves in a compressible atmosphere, Eqs. (9) and (10) in Zhou and Morton (2007) are full descriptions of the dispersion relation. For waves with a small intrinsic horizontal phase speed ($|c-\stackrel{\mathrm{\u203e}}{u}|<\mathrm{0.5}{c}_{\mathrm{s}}$), which is valid for most observed gravity waves, the dispersion relation can be described as

where ${H}_{\mathrm{s}}=R{T}_{\mathrm{0}}/g$ is the density scale height and *γ* is the ratio of specific heat, and *c*, $\stackrel{\mathrm{\u203e}}{u}$, and *c*_{s} are the observed horizontal phase speed, background wind speed in the direction of wave propagation, and speed of sound, respectively. Horizontal wavenumber *k*_{H} is related to *k* and *l* through ${k}_{H}^{\mathrm{2}}={k}^{\mathrm{2}}+{l}^{\mathrm{2}}$. The term $c-\stackrel{\mathrm{\u203e}}{u}$ is the intrinsic horizontal phase speed, usually denoted as $\widehat{c}$. The corresponding intrinsic wave frequency is related to the observed wave frequency by $\widehat{\mathit{\omega}}=\mathit{\omega}-{k}_{H}\cdot \stackrel{\mathrm{\u203e}}{u}$. When the atmosphere is treated as incompressible where the acoustic wave is eliminated, and background temperature varies slowly within the vertical wavelength of the wave, we have *c*_{s}→∞ and $\mathrm{d}{H}_{\mathrm{s}}/\mathrm{d}z\to \mathrm{0}$.

The dispersion relation (Eq. 9) is reduced to the following form that is derived based on Taylor–Goldstein equation (Nappo, 2012, Eq. 2.34):

The coefficient of the last term in Eq. (10) is different from the original one due to a correction for a compressible atmosphere based on discussions in Zhou and Morton (2007). If the wind shear terms are further neglected, the dispersion relation (Eq. 10) is simplified to Eq. (24) in Fritts and Alexander (2003) but without the Coriolis term and is also the same as the dispersion relation derived by Hines (1960):

Through the dispersion equations, the vertical wavenumber *m* is related to the wave characteristics, including the horizontal wavenumbers *k*_{H} and phase speed *c*, and the background states, including the projected wind on wave propagation direction $\stackrel{\mathrm{\u203e}}{u}$ and wind shear, as well as the background temperature *T*_{0} and its gradient as reflected by scale height *H*_{s} and buoyancy frequency N. In the regions of the atmosphere where *m*^{2}>0, gravity waves are able to propagate freely and are characterized by corresponding *m*, *k*, *l*, and *c*. Regions of *m*^{2}<0 indicate evanescence for gravity waves, where wave amplitudes decay exponentially. When a propagating wave encounters a region where *m*^{2}<0, partial or total reflection can occur depending on the depth of the evanescent region. Gravity waves whose propagation is restricted in a region with reflective regions of evanescence, above and below, are said to be ducted. At altitudes where $c=\stackrel{\mathrm{\u203e}}{u}$, the vertical wavenumber approaches infinity, and the waveform is overturned, so wave breaking or dissipation occurs. The phenomenon is called wave critical-layer filtering, which could be partial or total. In this case study, the altitude range is limited to within 85–105 km, during which the background temperature gradient is about −5 K km^{−1} and the wind shear is as strong as 30 m s^{−1} km^{−1}. The contribution of terms related to the temperature gradient and wind shear is not insignificant. Since all wave parameters and background states are explicitly determined, we evaluated all three forms of dispersion relations to diagnose the propagation of the identified wave packets.

Figure 9a and d show the calculated *m*^{2} for the two retrieved wave packets based on the full dispersion relation of Eq. (9). The *m*^{2} shows up in layered structures for both waves, potentially creating ducts for the gravity waves. For the 3.35 h wave (#1), there are major layers of negative *m*^{2} around 85–86, 95, 98, and 103 km altitudes. There are also several thin layers of negative *m*^{2} at other altitudes lasting shorter times. When upward propagating waves encounter the layers, their amplitudes will attenuate on top of the growth due to decreased atmospheric density. For the evanescent layer around 85–86 km, the waves are supposed to dampen at a fast rate (larger value of *i*×*m*). However, as this layer is thin, with a maximum thickness of 2 km, plenty of wave energy could penetrate the layer to higher altitudes. This also applies to other thin evanescent layers where the attenuation is even less. The upward propagating waves were partially reflected and refracted at each evanescent layer, which could change the wavefront orientation. The partial refraction/reflection was also presented by the previous numerical simulations in similar scenarios (Heale et al., 2014; Cao et al., 2016). For the temperature perturbation in Fig. 6, the wave pattern shows up with a nodal structure with clear amplitude maximums around 85, 92, and 102 km, and discontinuities in the phase progression in the vertical direction. These features directly result from the partial reflection and refraction at the multiple evanescent layers where the perturbations observed by the lidar are the superposition of incident and reflected waves between the evanescent layers. For the 1.63 h wave (#2), the overall layered structure in *m*^{2} remains similar to the 3.35 h wave. However, the 3.35 h wave effectively changes wind and temperature, as well as creates differences in *m*^{2}, and leads to different propagation conditions for the 1.63 h wave. The evanescent layer around 86 km reduced the attenuation rate from 09:00 UT onwards, which largely increased the portion of the transmission of wave energy. In addition, the 3.35 h wave increased the evanescent layer thickness to 5 km around 95 km, which could limit the transmissible wave energy to further altitudes. As shown in Fig. 7, there is a visible maximum below 85 km, and another maximum was around 90 km, whose wave amplitudes largely increased after 09:00 UT. The amplitudes decreased dramatically right above 95 km; however, some wave peaks can be seen above 100 km. The average *m*^{2} profiles of the whole period are calculated for both wave packets and the corresponding vertical wavelengths are estimated and shown in Fig. 9b and d. The positive vertical wavelengths for the freely propagating waves are about 15–20 km for wave #1 and 10–15 km for wave #2, similar to the ones estimated from the downward phase progression slope. The negative wavelength corresponds to the evanescent region and is equivalent to the scale height of the wave amplitude decay. The overall decay/growth rate of the wave amplitude is described by ${e}^{(\mathrm{1}/\mathrm{5}{H}_{\mathrm{s}}-im)z}$, where decay occurs when *m* is an imaginary number in the evanescent region. The factor of 5*H*_{s}, rather than 2*H*_{s}, was chosen to mimic wave dissipation; otherwise, the wave amplitude would grow too fast. In Fig. 9c and f, the relative wave amplitudes are estimated, assuming one unit of wave amplitude at the lowest altitude. The predicted wave amplitude shows fluctuation and several maxima at different altitudes, which are the combined efforts of evanescent decay and conserved growth. The evanescent layers at different altitudes would attenuate the wave amplitude, but layers with the time-varying *m*^{2} could reduce the dampening efforts and spare more energy penetrating to higher altitudes. Bossert et al. (2014) presents case studies using lidar observations and simulations to show high-frequency (period shorter than 15 min) gravity waves propagating to higher altitudes through alternating regions of evanescence and free propagation over a few kilometers. In this study, the two medium-frequency wave packets observed by the lidar propagate through multiple thin and time-varying evanescent layers with a good portion of wave energy penetrating to higher altitudes, and partial reflection and refraction occur with the observed amplitude maxima found between these evanescent layers.

The *m*^{2} estimated by Eqs. (10) and (11) are shown in Figs. A7 and A8. The results of Eq. (10) show overall similarity with the ones of Eq. (9), but some differences exist. The evanescent layers estimated by Eq. (10) are slightly thinner, which shows the contribution of neglected temperature gradients. The *m*^{2} calculated by Eq. (11) fails to capture most of the layered structures and underestimates all evanescent layer thickness as wind shear is a significant contributing factor for both waves. As discussed above, inconsistency exists among different dispersion relations, and some simplifications fail to capture the authentic characters. In the application of dispersion relations, the full background temperature and wind measurements might not be available in all cases or are sometimes limited only to a few altitudes. Nevertheless, simplified formulas (Eqs. 10 and 11) can be best utilized to diagnose wave propagation; however, the results should be interpreted with caution.

Another important relation derived from linearized wave equations is the polarization relation that describes the relative phase differences and amplitude ratios of various wave quantities. If gravity waves do not undergo dissipation, the complex wave amplitudes of the relative temperature $\stackrel{\mathrm{\u0303}}{T}\left({T}^{\prime}/\stackrel{\mathrm{\u203e}}{T}\right)$, zonal wind $\stackrel{\mathrm{\u0303}}{u}$, meridional wind $\stackrel{\mathrm{\u0303}}{v}$, and vertical wind $\stackrel{\mathrm{\u0303}}{w}$ should satisfy the following polarization relations (Fritts and Alexander, 2003; Vadas, 2013; Lu et al., 2015):

The complex amplitudes of $\stackrel{\mathrm{\u0303}}{T}/\stackrel{\mathrm{\u0303}}{w}$, $\stackrel{\mathrm{\u0303}}{T}/\stackrel{\mathrm{\u0303}}{u}$, and $\stackrel{\mathrm{\u0303}}{u}/\stackrel{\mathrm{\u0303}}{v}$ describe the amplitude and phase relations among different wave quantities. On the one hand, the missing quantities of observed gravity waves can be estimated through these relations assuming non-dissipation. On the other hand, the discrepancies between observed and theoretical values can be used to indicate wave dissipation. It is also possible to estimate higher-order statistical quantities, such as gravity wave momentum flux $\left(\stackrel{\mathrm{\u203e}}{{u}^{\prime}{w}^{\prime}}\right)$ and heat flux $\left(\stackrel{\mathrm{\u203e}}{{w}^{\prime}{T}^{\prime}}\right)$, from these relationships with limited observations (Liu, 2009; Guo et al., 2017).

^{a} $\stackrel{\mathrm{\u0303}}{T}$ here is relative temperature perturbation ${T}^{\prime}/\stackrel{\mathrm{\u203e}}{T}$ expressed in percentage. ^{b} The mean and standard derivation of the quantities are calculated within the corresponding altitude range. ^{c} The altitude range of 85–100 km was used to calculate the mean and standard deviation.

The predicted values of $\stackrel{\mathrm{\u0303}}{T}/\stackrel{\mathrm{\u0303}}{w}$, $\stackrel{\mathrm{\u0303}}{T}/\stackrel{\mathrm{\u0303}}{u}$, and $\stackrel{\mathrm{\u0303}}{u}/\stackrel{\mathrm{\u0303}}{v}$ can be calculated from Eq. (12) with the retrieved wave and background parameters. They are complex numbers, with their absolute magnitudes representing the wave amplitude (*A*) ratio and phases representing the phase (*φ*) difference between any two quantities. Two sets of amplitude ratios and phase differences are calculated, one for positive *m*^{2} and the other for negative *m*^{2}, corresponding to the wave free-propagating and evanescent regions. Table 2 lists the predicted results for both wave packets. In the region of free propagating, the values are relatively constant; the variations are mainly because of the change in intrinsic frequency along altitudes. However, large variations exist for the phase differences in the evanescent region. In this case, the wave packets propagate through multiple evanescent layers where partial reflection occurs. The lidar observed perturbations are the superposition of incident and reflected waves, and the propagating waves also undergo dissipation. It is difficult to accurately estimate the amplitude ratios and phase differences from the observed perturbations. Another difficulty is that the horizontal winds are not retrieved at the zenith direction and have to be indirectly estimated from horizontal winds at off-zenith directions by correcting the phase shift. Therefore, we only estimate the averaged amplitude ratios from the observations. In Table 2, the observed results are estimated with larger uncertainties, and all show discrepancies with the predicted ones. However, the observed values of $A\left(\stackrel{\mathrm{\u0303}}{T}\right)/A\left(\stackrel{\mathrm{\u0303}}{u}\right)$ are closer to the predicted ones in the evanescent region. The uncertainties in the wind measurements and wave dissipation could also contribute to this. For the ratio $A\left(\stackrel{\mathrm{\u0303}}{T}\right)/A\left(\stackrel{\mathrm{\u0303}}{w}\right)$, it is much smaller than the predicted values, which implies that the observed perturbations in the vertical wind might contain large errors. The presented results reveal the complexity of gravity waves propagating in the atmosphere. The polarization relation is good for diagnosing the free propagating waves without being reflected and refracted, and dissipation is not severe. It might not be proper when complicated wave-mean flow interaction occurs, such as in this wave case.

Using the proposed cross-spectral method, we identified two gravity wave packets and retrieved all the wave parameters, as well as determined the background states. The fully retrieved information was used to validate the linear gravity wave theories with the least assumptions. Consistent results are obtained from the diagnostic analyses using the dispersion relation, partly explaining the wave observations. However, raw lidar measurement uncertainties exert difficulties in some wave parameter estimations and diagnostic analyses using the polarization relation. In the next section, we implement a sensitivity study to evaluate the general usage of this method in detecting gravity waves.

It is well known that most observation techniques are restricted by the “observation filter” effect in resolving atmospheric waves. These techniques are sensitive to certain parts of the spatial and temporal spectra of the waves. The effect also applies to the lidar and the wave extraction method presented in this study. To find a spectral range of gravity waves that is more favorable to be identified by this method, we did a sensitivity study using a forward simulation. The accuracy of the cross-spectral methods in recovering the wave parameters of amplitude, period, and phase shift is quantified.

The off-zenith angle of laser beams coupled with steerable telescopes could be adjusted. However, this configuration is no longer available, and newer lidar systems are equipped with multiple fixed telescopes pointing to certain directions; thus, the off-zenith direction and angle are fixed. The photon integration time at each direction and the laser beam rotating sequence could be changed based on applications, such as only one direction in the zonal and meridional directions. There are also some lidar systems with one master laser beam being split and shooting to multiple directions simultaneously. The typical photon integration time at each direction is about 1 min to several minutes. A configuration similar to the lidar deployed in Maui is utilized in the simulation where a 30° off-zenith angle corresponds to a separation of ∼ 100 km at 90 km altitude between two off-zenith directions (W and E, or S and N). For simplicity, the temporal resolution was selected as a uniform 6 min in all directions. Because the retrieval of wave parameters in the zonal and meridional directions (*k* and *l*) is independent, the simulations only consider measurements from two laser beams aligned in one direction. The perturbations of a traveling gravity wave observed by two laser beams are described by

of which the definitions of the terms are the same as in Eq. (4). Extra terms *δ*_{1} and *δ*_{2}, two independent (associated with different “seeds” of the pseudorandom number generator) sets of Gaussian-distributed random numbers, are introduced to mimic the uncertainties in measurements. The cross-spectral method is applied on the data (*y*_{1}, *y*_{2}) to retrieve the wave amplitude *A*, period $\mathit{\tau}=\mathrm{2}\mathit{\pi}/\mathit{\omega}$, and phase shift $\mathrm{\Delta}\mathit{\varphi}={\mathit{\varphi}}_{\mathrm{1}}-{\mathit{\varphi}}_{\mathrm{2}}$. To qualify the accuracy, the percentage errors of wave amplitude *A*, period *τ*, and phase difference Δ*ϕ* between the retrieved results and preset values are calculated.

In the simulation, the wave amplitude *A* was selected as two (with an arbitrary unit) and the random noise *δ* had a standard deviation of 0.5. Therefore, the following simulation and discussion only apply to this signal-to-noise ratio ($A/\mathit{\delta}$). In the simulation, 5 h of data (50 data points) were used. In order to recover the spectral peaks falling between two integer spectral points, the aforementioned fitting method was used to best retrieve the actual peaks. The horizontal wavelength was selected to vary from 300 to 3000 km (3–30 times the laser beam separation) and the period range was from 1 to 5 h. In Fig. 10, the percentage errors of wave amplitude (*A*), period (*τ*), and phase shift (Δ*ϕ*) are shown. The retrieved wave amplitudes have an error of 10 % from the true value for all periods and wavelengths. A slight negative bias of 5 % and a wave pattern are identified in amplitude errors. The retrieved periods have a 5 % error for periods shorter than 3 h and up to 10 %–20 % for waves with longer periods and shorter wavelengths. The phase shift has the largest error among these wave parameters, with errors up to 20 % for waves with wavelengths shorter than 2000 km and up to 60 % for the waves with longer horizontal wavelengths (2000–3000 km). The errors in phase shift are independent of the wave periods. The cross-spectral method is likely unable to retrieve the small phase shifts (≤ 18°) resulting from waves with larger scales, which include the residues or aliased signals from tides that often have planetary scales.

The proposed method relies on cross-spectral analysis and inherits the typical drawback of any spectral method. In this simulation, there is only one frequency being simulated. If spectral leakage is not an issue when the actual peak is close to the integer spectral points, the wave parameters could be determined accurately. However, with limited samples, the resolution in the spectral domain is coarse, and there is a chance of spectral leakage that the actual peaks fall between two integer spectral points. The negative bias of the retrieved amplitude and the wave pattern in the errors of the wave parameters (*A* and *τ*) are related to spectral leaks. When multiple wave components are mixed together, the spectral leakage will be more complicated and limit wave identification. With measurements from only two directions, the retrieved wavelength and phase speed are apparent ones along that direction. The actual wavelength and speed in any direction can be the composition of apparent ones in the zonal and meridional directions. As narrow-band resonance lidars usually provide measurements over an altitude range between 80 and 110 km, the wave propagating in the vertical direction can be captured. The retrieved period, wavelength, and phase speed are invariant if the waves do not undergo severe dissipation. This could provide another level of verification of the retrieved wave parameters at different altitudes.

Here, we provide an introductory assessment of the proposed wave extraction method. The data duration is critical, as more data points could help alleviate the spectral leakage and retrieve accurate wave parameters. However, the duration is limited by the operation of such type of lidar that mostly works at night unless extra filters are used to facilitate daytime operation (Chen et al., 1996). There is a trade-off between longer datasets for higher spectral resolution and the duration of wave presence with an invariant period and amplitude. The most favorable wave period of this method is limited to be shorter than the data duration, to minimize the spectral leakage, and longer than the temporal resolution to achieve enough data points to resolve the variation in one wave cycle. The phase shifts need to be large enough to be distinguished by the cross-spectral method, which makes the method unsuitable for waves with very long wavelengths. If observations have higher signal-to-noise ratios ($A/\mathit{\delta}$), the method should be able to identify smaller phase shifts and determine longer wavelengths. It is hard to give an applicable wavelength range for this method, as it is highly subject to the actual data quality and off-zenith angle. Such a sensitivity analysis can be implemented for specific lidar system configurations. In practical terms, the proposed method is adept at detecting gravity waves of medium scale and medium frequency for a lidar operated at similar configurations and with comparable measurement uncertainties.

With measurements from a single ground-based instrument, a narrow-band sodium lidar operated in multiple-direction mode, two gravity wave packets were detected on the night of 14 January 2002 at Maui, Hawaii, and were resolved in 3-D space by a cross-spectral method. The retrieved phase differences among measurements in different directions enable the retrieval of the horizontal wave information. With this method, the horizontal wavelength and phase speed were estimated with ∼ 20 % uncertainties. One wave with a horizontal wavelength of 934 km and a period of 3.35 h propagated toward a 282° azimuth at a phase speed of 78 m s^{−1}, and the other one with a horizontal wavelength of 438 km and a period of 1.63 h propagated toward a 233° azimuth at a phase speed of 76 m s^{−1}. Both waves propagated toward the nearly opposite direction of the background winds and larger wave amplitudes were found in the relatively stable regions, as indicated by the squared buoyancy frequency *N*^{2} and Richardson number *Ri*. With full sets of wave parameters and background states determined, multiple versions of dispersion relation, some with simplifications, are examined in this study. Both wave packets are found to propagate through multiple thin evanescent layers where partial reflection and refraction possibly occur around those evanescent layers. However, a good portion of wave energy penetrates to higher altitudes where waves undergo further dissipation and dispersion. The longer-period (3.35 h) wave effectively changes the background and leads to a different propagation condition for the shorter-period (1.63 h) wave. The comparisons among different versions of dispersion relations show that the effects of the background temperature gradient and wind shear are important in the linearized wave theory, and diagnostic analysis based on simplified dispersion relations should be interpreted with caution. Polarization relations are also examined among terms $\stackrel{\mathrm{\u0303}}{T}$, $\stackrel{\mathrm{\u0303}}{u}$, $\stackrel{\mathrm{\u0303}}{v}$, and $\stackrel{\mathrm{\u0303}}{w}$. However, the complexity of the wave propagation conditions and uncertainties in the measurements, especially of the winds, make it difficult to retrieve accurate results from the data.

Continuous lidar measurement profiles with proper sampling rate and duration can capture a wide variety of periods of waves. However, those profiles can only resolve the vertical variations, and horizontal information is often complemented by other observations or inferred indirectly. In this study, we propose a novel method using data from a single-site lidar configured in multiple-direction observing mode to fully resolve gravity waves in 3-D space, with both horizontal and vertical wave information retrieved directly. The sensitivity study reveals the capability of this method in detecting medium-scale and medium-frequency gravity waves. This partially makes up the spectral gaps in the “observation filter” for this type of lidar. The proposed method could provide extra opportunities for the gravity wave studies based on lidar systems with similar configurations that were deployed at other sites, either in the past or still in operation (Hu et al., 2002; Hildebrand et al., 2012; Cai et al., 2014; Ban et al., 2015; Liu et al., 2016; Li et al., 2020). Unlike the lidar presented in this study, newer lidar systems are equipped with two to four telescopes that are pointed in several directions depending on research requirements, and the off-zenith angles are often fixed. If their sampling interval and rotating sequence are properly configured, it is feasible to use this method to detect more medium-scale and medium-frequency gravity wave events. The determination of 3-D wave parameters combined with background atmosphere states would also enable backward ray tracing to identify the wave source location (Vadas et al., 2009; Krisch et al., 2017; Krasauskas et al., 2023). To further examine wave propagation, reflection, and dissipation, a numerical simulation that takes account of the complete wave parameters and background states would provide important insights into the interpretation of the observed results and those unobserved beyond the field of view.

The narrow-band sodium lidar data, including sodium density, temperature, and winds from the Maui/MALT campaign and other deployments of an upgraded lidar system operated by UIUC and Embry-Riddle Aeronautical University, can be found at http://alo.erau.edu/data/nalidar/index.php (Liu, 2022), and the data used for this study can also be downloaded at https://doi.org/10.5281/zenodo.8124900 (Liu, 2023).

AZL processed the raw lidar data and administrated the data curation. BC developed the methodology and analyzed the data. BC prepared the paper with the contribution of AZL.

The contact author has declared that neither of the authors has any competing interests.

The work by Alan Z. Liu is supported by (while serving at) the U.S. National Science Foundation (NSF). Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the NSF.

Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.

This research has been supported by the U.S. National Science Foundation (NSF) (grant nos. AGS-1110199, AGS-1115249, and AGS-1759471).

This paper was edited by Wen Yi and reviewed by three anonymous referees.

Alexander, M. J.: Interpretations of observed climatological patterns in stratospheric gravity wave variance, J. Geophys. Res., 103, 8627–8640, https://doi.org/10.1029/97JD03325, 1998. a

Alexander, M. J. and Grimsdell, A. W.: Seasonal cycle of orographic gravity wave occurrence above small islands in the Southern Hemisphere: Implications for effects on the general circulation, J. Geophys. Res.-Atmos., 118, 11589–11599, https://doi.org/10.1002/2013JD020526, 2013. a

Alexander, M. J., Gille, J., Cavanaugh, C., Coffey, M., Craig, C., Eden, T., Francis, G., Halvorson, C., Hanningan, J., Khosravi, R., Kinnison, D., Lee, H., Massie, S., Nardi, B., Barnett, J., Hepplewhite, C., Lambert, A., and Dean, V.: Global estimates of gravity wave momentum flux from High Resolution Dynamics Limb Sounder observations, J. Geophys. Res., 113, D15S18, https://doi.org/10.1029/2007JD008807, 2008. a

Alexander, M. J., Geller, M., McLandress, C., Polavarapu, S., Preusse, P., Sassi, F., Sato, K., Eckermann, S., Ern, M., Hertzog, A., Kawatani, Y., Pulido, M., Shaw, T. A., Sigmond, M., Vincent, R., and Watanabe, S.: Recent developments in gravity-wave effects in climate models and the global distribution of gravity-wave momentum flux from observations and models, Q. J. Roy. Meteorol. Soc., 136, 1103–1124, https://doi.org/10.1002/qj.637, 2010. a

Alexander, P., de la Torre, A., Schmidt, T., Llamedo, P., and Hierro, R.: Limb sounders tracking topographic gravity wave activity from the stratosphere to the ionosphere around midlatitude Andes, J. Geophys. Res.-Space, 120, 9014–9022, https://doi.org/10.1002/2015JA021409, 2015. a

Alexander, P., Schmidt, T., and de la Torre, A.: A Method to Determine Gravity Wave Net Momentum Flux, Propagation Direction, and “Real” Wavelengths: A GPS Radio Occultations Soundings Case Study, Earth and Space Science, 5, 222–230, https://doi.org/10.1002/2017EA000342, 2018. a

Alexander, S. P., Klekociuk, A. R., and Tsuda, T.: Gravity wave and orographic wave activity observed around the Antarctic and Arctic stratospheric vortices by the COSMIC GPS-RO satellite constellation, J. Geophys. Res., 114, D17103, https://doi.org/10.1029/2009JD011851, 2009. a

Ban, C., Li, T., Fang, X., Dou, X., and Xiong, J.: Sodium lidar-observed gravity wave breaking followed by an upward propagation of sporadic sodium layer over Hefei, China, J. Geophys. Res.-Space, 120, 7958–7969, https://doi.org/10.1002/2015JA021339, 2015. a

Bossert, K., Fritts, D. C., Pautet, P. D., Taylor, M. J., Williams, B. P., and Pendleton, W. R.: Investigation of a mesospheric gravity wave ducting event using coordinated sodium lidar and Mesospheric Temperature Mapper measurements at ALOMAR, Norway (69° N), J. Geophys. Res., 119, 9765–9778, https://doi.org/10.1002/2014JD021460, 2014. a, b

Cai, X., Yuan, T., Zhao, Y., Pautet, P.-D., Taylor, M. J., and Pendleton, W. R.: A coordinated investigation of the gravity wave breaking and the associated dynamical instability by a Na lidar and an Advanced Mesosphere Temperature Mapper over Logan, UT (41.7° N, 111.8° W), J. Geophys. Res.-Space, 119, 6852–6864, https://doi.org/10.1002/2014JA020131, 2014. a, b

Cao, B. and Liu, A. Z.: Statistical Characteristics of High-Frequency Gravity Waves Observed by an Airglow Imager at Andes Lidar Observatory, Earth and Space Science, 9, e2022EA002256, https://doi.org/10.1029/2022EA002256, 2022. a

Cao, B., Heale, C. J., Guo, Y., Liu, A. Z., and Snively, J. B.: Observation and modeling of gravity wave propagation through reflection and critical layers above Andes Lidar Observatory at Cerro Pachón, Chile, J. Geophys. Res.-Atmos., 121, 12737–12750, https://doi.org/10.1002/2016JD025173, 2016. a, b

Chen, C., Chu, X., Zhao, J., Roberts, B. R., Yu, Z., Fong, W., Lu, X., and Smith, J. A.: Lidar observations of persistent gravity waves with periods of 3–10 h in the Antarctic middle and upper atmosphere at McMurdo (77.83° S, 166.67° E), J. Geophys. Res.-Space, 121, 1483–1502, https://doi.org/10.1002/2015JA022127, 2016. a

Chen, H., White, M. A., Krueger, D. A., and She, C. Y.: Daytime mesopause temperature measurements with a sodium-vapor dispersive Faraday filter in a lidar receiver, Opt. Lett., 21, 1093, https://doi.org/10.1364/ol.21.001093, 1996. a

Ern, M., Ploeger, F., Preusse, P., Gille, J. C., Gray, L. J., Kalisch, S., Mlynczak, M. G., Russell, J. M., and Riese, M.: Interaction of gravity waves with the QBO: A satellite perspective, J. Geophys. Res.-Atmos., 119, 2329–2355, https://doi.org/10.1002/2013JD020731, 2014. a

Ern, M., Preusse, P., and Riese, M.: Driving of the SAO by gravity waves as observed from satellite, Ann. Geophys., 33, 483–504, https://doi.org/10.5194/angeo-33-483-2015, 2015. a

Espy, P. J., Hibbins, R. E., Swenson, G. R., Tang, J., Taylor, M. J., Riggin, D. M., and Fritts, D. C.: Regional variations of mesospheric gravity-wave momentum flux over Antarctica, Ann. Geophys., 24, 81–88, https://doi.org/10.5194/angeo-24-81-2006, 2006. a

Faber, A., Llamedo, P., Schmidt, T., de la Torre, A., and Wickert, J.: On the determination of gravity wave momentum flux from GPS radio occultation data, Atmos. Meas. Tech., 6, 3169–3180, https://doi.org/10.5194/amt-6-3169-2013, 2013. a

Fritts, D. C.: Gravity wave saturation in the middle atmosphere: A review of theory and observations, Rev. Geophys. Space Phys., 22, 275–308, 1984. a

Fritts, D. C. and Alexander, M. J.: Gravity wave dynamics and effects in the middle atmosphere, Rev. Geophys., 41, 1003, https://doi.org/10.1029/2001RG000106, 2003. a, b, c, d

Fritts, D. C. and Lund, T. S.: Gravity wave influences in the thermosphere and ionosphere: Observations and recent modeling, in: Aeronomy of the Earth's Atmosphere and Ionosphere, Springer, 109–130, https://doi.org/10.1007/978-94-007-0326-1_8, 2011. a

Fritts, D. C., Janches, D., and Hocking, W. K.: Southern Argentina Agile Meteor Radar: Initial assessment of gravity wave momentum fluxes, J. Geophys. Res., 115, D19123, https://doi.org/10.1029/2010JD013891, 2010. a

Fritts, D. C., Wang, L., and Werne, J. A.: Gravity Wave–Fine Structure Interactions. Part I: Influences of Fine Structure Form and Orientation on Flow Evolution and Instability, J. Atmos. Sci., 70, 3710–3734, https://doi.org/10.1175/JAS-D-13-055.1, 2013. a

Fritts, D. C., Pautet, P.-D., Bossert, K., Taylor, M. J., Williams, B. P., Iimura, H., Yuan, T., Mitchell, N. J., and Stober, G.: Quantifying gravity wave momentum fluxes with Mesosphere Temperature Mappers and correlative instrumentation, J. Geophys. Res.-Atmos., 119, 13583–13603, https://doi.org/10.1002/2014JD022150, 2014. a

Gardner, C. S. and Liu, A. Z.: Seasonal variations of the vertical fluxes of heat and horizontal momentum in the mesopause region at Starfire Optical Range, New Mexico, J. Geophys. Res., 112, D09113, https://doi.org/10.1029/2005JD006179, 2007. a

Gardner, C. S. and Taylor, M. J.: Observational limits for lidar, radar, and airglow imager measurements of gravity wave parameters, J. Geophys. Res., 103, 6427–6437, 1998. a

Goldstein, S.: On the Stability of Superposed Streams of Fluids of Different Densities, P. Roy. Soc. Lond. A, 132, 524–548, https://doi.org/10.1098/rspa.1931.0116, 1931. a

Gong, J., Wu, D. L., and Eckermann, S. D.: Gravity wave variances and propagation derived from AIRS radiances, Atmos. Chem. Phys., 12, 1701–1720, https://doi.org/10.5194/acp-12-1701-2012, 2012. a

Guo, Y., Liu, A. Z., and Gardner, C. S.: First Na lidar measurements of turbulence heat flux, thermal diffusivity, and energy dissipation rate in the mesopause region, Geophys. Res. Lett., 44, 5782–5790, https://doi.org/10.1002/2017GL073807, 2017GL073807, 2017. a

Heale, C. J., Snively, J. B., and Hickey, M. P.: Numerical simulation of the long-range propagation of gravity wave packets at high latitudes, J. Geophys. Res.-Atmos., 119, 11116–11134, https://doi.org/10.1002/2014JD022099, 2014. a

Hildebrand, J., Baumgarten, G., Fiedler, J., Hoppe, U.-P., Kaifler, B., Lübken, F.-J., and Williams, B. P.: Combined wind measurements by two different lidar instruments in the Arctic middle atmosphere, Atmos. Meas. Tech., 5, 2433–2445, https://doi.org/10.5194/amt-5-2433-2012, 2012. a

Hines, C. O.: Internal atmospheric gravity waves at ionospheric heights, Can. J. Phys., 38, 1441–1481, 1960. a

Hoffmann, L., Alexander, M. J., Clerbaux, C., Grimsdell, A. W., Meyer, C. I., Rößler, T., and Tournier, B.: Intercomparison of stratospheric gravity wave observations with AIRS and IASI, Atmos. Meas. Tech., 7, 4517–4537, https://doi.org/10.5194/amt-7-4517-2014, 2014. a

Holton, J. R.: The role of gravity wave induced drag and diffusion in the momentum budget of the mesosphere, J. Atmos. Sci., 39, 791–799, 1982. a

Hu, X., Liu, A. Z., Gardner, C. S., and Swenson, G. R.: Characteristics of quasi-monochromatic gravity waves observed with lidar in the mesopause region at Starfire Optical Range, NM, Geophys. Res. Lett., 29, 2169, https://doi.org/10.1029/2002GL014975, 2002. a, b, c

Jiang, J., Eckermann, S., Wu, D., Hocke, K., Wang, B., Ma, J., and Zhang, Y.: Seasonal variation of gravity wave sources from satellite observation, Adv. Space Res., 35, 1925–1932, https://doi.org/10.1016/j.asr.2005.01.099, 2005. a

Krasauskas, L., Kaifler, B., Rhode, S., Ungermann, J., Woiwode, W., and Preusse, P.: Oblique Propagation and Refraction of Gravity Waves Over the Andes Observed by GLORIA and ALIMA During the SouthTRAC Campaign, J. Geophys. Res.-Atmos., 128, e2022JD037798, https://doi.org/10.1029/2022JD037798, 2023. a

Krisch, I., Preusse, P., Ungermann, J., Dörnbrack, A., Eckermann, S. D., Ern, M., Friedl-Vallon, F., Kaufmann, M., Oelhaf, H., Rapp, M., Strube, C., and Riese, M.: First tomographic observations of gravity waves by the infrared limb imager GLORIA, Atmos. Chem. Phys., 17, 14937–14953, https://doi.org/10.5194/acp-17-14937-2017, 2017. a

Krueger, D. A., She, C.-Y., and Yuan, T.: Retrieving mesopause temperature and line-of-sight wind from full-diurnal-cycle Na lidar observations, Appl. Opt., 54, 9469–9489, https://doi.org/10.1364/AO.54.009469, 2015. a, b

Li, J., Williams, B. P., Alspach, J. H., and Collins, R. L.: Sodium Resonance Wind-Temperature Lidar at PFRR: Initial Observations and Performance, Atmosphere, 11, 98, https://doi.org/10.3390/atmos11010098, 2020. a

Li, T., She, C.-Y., Liu, H.-L., Leblanc, T., and McDermid, I. S.: Sodium lidar–observed strong inertia-gravity wave activities in the mesopause region over Fort Collins, Colorado (41° N, 105° W), J. Geophys. Res.-Atmos., 112, D22104, https://doi.org/10.1029/2007JD008681, 2007. a

Li, T., Fang, X., Liu, W., Gu, S.-Y., and Dou, X.: Narrowband sodium lidar for the measurements of mesopause region temperature and wind, Appl. Opt., 51, 5401–5411, https://doi.org/10.1364/AO.51.005401, 2012. a

Li, Z., Liu, A. Z., Lu, X., Swenson, G. R., and Franke, S. J.: Gravity wave characteristics from OH airglow imager over Maui, J. Geophys. Res., 116, D22115, https://doi.org/10.1029/2011JD015870, 2011. a

Liu, A. Z.: Estimate eddy diffusion coefficients from gravity wave vertical momentum and heat fluxes, Geophys. Res. Lett., 36, L08806, https://doi.org/10.1029/2009GL037495, 2009. a

Liu, A. Z.: Andes Lidar Observatory Dataset, Andes Lidar Observatory [data set], http://alo.erau.edu/data/nalidar/index.php, last access: 1 February 2022. a

Liu, A. Z.: Na Lidar Data acquired at Haleakala, Maui, Hawaii from 2002–2005, Zenodo [data set], https://doi.org/10.5281/zenodo.8124900, 2023. a, b

Liu, A. Z., Lu, X., and Franke, S. J.: Diurnal variation of gravity wave momentum flux and its forcing on the diurnal tide, J. Geophys. Res., 118, 1668–1678, 2013. a

Liu, A. Z., Guo, Y., Vargas, F., and Swenson, G. R.: First measurement of horizontal wind and temperature in the lower thermosphere (105–140 km) with a Na Lidar at Andes Lidar Observatory, Geophys. Res. Lett., 43, 2374–2380, https://doi.org/10.1002/2016GL068461, 2016. a

Liu, H.-L. and Vadas, S. L.: Large-scale ionospheric disturbances due to the dissipation of convectively-generated gravity waves over Brazil, J. Geophys. Res.-Space, 118, 2419–2427, https://doi.org/10.1002/jgra.50244, 2013. a

Lu, X., Liu, A. Z., Swenson, G. R., Li, T., Leblanc, T., and McDermid, I. S.: Gravity wave propagation and dissipation from the stratosphere to the lower thermosphere, J. Geophys. Res.-Atmos., 114, D11101, https://doi.org/10.1029/2008JD010112, 2009. a

Lu, X., Chen, C., Huang, W., Smith, J. A., Chu, X., Yuan, T., Pautet, P.-D., Taylor, M. J., Gong, J., and Cullens, C. Y.: A coordinated study of 1 h mesoscale gravity waves propagating from Logan to Boulder with CRRL Na Doppler lidars and temperature mapper, J. Geophys. Res.-Atmos., 120, 10006–10021, https://doi.org/10.1002/2015JD023604, 2015. a, b, c

Nappo, C.: An Introduction to Atmospheric Gravity Waves, International Geophysics, Elsevier Science, ISBN 9780123852236, 2012. a, b

Nastrom, G. D. and Eaton, F. D.: Quasi-monochromatic inertia-gravity waves in the lower stratosphere from MST radar observations, J. Geophys. Res.-Atmos., 111, D19103, https://doi.org/10.1029/2006JD007335, 2006. a

She, C. Y. and Yu, J. R.: Simultaneous three-frequncy Na lidar measurements of radial wind and temperature in the mesopause region, Geophys. Res. Lett., 21, 1771–1774, 1994. a

Siskind, D. E., Drob, D. P., Emmert, J. T., Stevens, M. H., Sheese, P. E., Llewellyn, E. J., Hervig, M. E., Niciejewski, R., and Kochenash, A. J.: Linkages between the cold summer mesopause and thermospheric zonal mean circulation, Geophys. Res. Lett., 39, L01804, https://doi.org/10.1029/2011GL050196, 2012. a

Taylor, G. I.: Effect of Variation in Density on the Stability of Superposed Streams of Fluid, P. Roy. Soc. Lond. A, 132, 499–523, https://doi.org/10.1098/rspa.1931.0115, 1931. a

Taylor, M. J.: A review of advances in imaging techniques for measuring short period gravity waves in the mesosphere and lower thermosphere, Adv. Space Res., 19, 667–676, 1997. a

Vadas, S. L.: Compressible f-plane solutions to body forces, heatings, and coolings, and application to the primary and secondary gravity waves generated by a deep convective plume, J. Geophys. Res.-Space, 118, 2377–2397, https://doi.org/10.1002/jgra.50163, 2013. a

Vadas, S. L., Yue, J., She, C.-Y., Stamus, P. A., and Liu, A. Z.: A model study of the effects of winds on concentric rings of gravity waves from a convective plume near Fort Collins on 11 May 2004, J. Geophys. Res., 114, D06103, https://doi.org/10.1029/2008JD010753, 2009. a

Wang, L. and Alexander, M. J.: Global estimates of gravity wave parameters from GPS radio occultation temperature data, J. Geophys. Res.-Atmos., 115, D21122, https://doi.org/10.1029/2010JD013860, 2010. a

Zhou, Q. and Morton, Y. T.: Gravity wave propagation in a nonisothermal atmosphere with height varying background wind, Geophys. Res. Lett., 34, L23803, https://doi.org/10.1029/2007GL031061, 2007. a, b, c, d