© The Authors 2024

1 Introduction

Gravity (or buoyancy) waves (GWs) originating from the balance of gravity and buoyancy forces are ubiquitous in all convec-tively stable atmospheres. They have been extensively studied in the terrestrial atmosphere, where their important role in the dynamics and vertical coupling of atmospheric layers has been recognized (e.g., see reviews by Fritts & Alexander 2003; Yiğit & Medvedev 2015). A GW-like signature in the atmosphere of Mars was first detected in the entry temperature profile of Viking 2 (Seiff & Kirk 1976). Since then, GWs have been observed on Venus (Young et al. 1987), Jupiter (Young et al. 1997), Titan (Hinson & Tyler 1983), Saturn (Brown et al. 2022), and other planets (see a recent review on GWs in planetary atmospheres by Medvedev & Yiğit 2019). Recently, various GW-related phenomena in the Martian atmosphere and their impact on the whole atmosphere of the Mars system have been reviewed by Yiğit (2023).

In situ accelerometer measurements in the thermosphere of Mars performed during spacecraft aerobraking have demonstrated large amplitudes of GW-induced density disturbances and the associated wave drag (Keating et al. 1998; Creasey et al. 2006a; Fritts et al. 2006). The omnipresence of GWs in the lower Martian atmosphere and the first characterization of the GW field were revealed using remote sensing techniques (Hinson et al. 1999; Creasey et al. 2006b; Wright 2012; Nakagawa et al. 2020; Heavens et al. 2020). They showed a varying spatio-temporal structure of the GW field and evidence of a limitation of amplitude growth with height, which is an indication of wave-momentum transfer to a larger scale flow (Ando et al. 2012).

Gravity waves exist at all atmospheric heights. In situ measurements using the Neutral Gas and Ion Mass Spectrometer (NGIMS) on board the Mars Atmosphere and Evolution Mission (MAVEN) orbiter delivered a large body of GW statistics in the upper thermosphere (Yiğit et al. 2015; England et al. 2017; Terada et al. 2017; Leelavathi et al. 2020; Rao et al. 2023). In particular, the observations found an enhancement of GW activity during the global dust storm (GDS) that occurred in June 2018 (Martian year 34; MY34) (Leelavathi et al. 2020; Yiğit et al. 2021a). Conversely, observations with the Mars Climate Sounder instrument on board the Mars Reconnaissance Orbiter have shown a reduction of such activity in the lower atmosphere (Heavens et al. 2020).

Another existing controversy concerns the inverse relation between the amplitudes of GWs and the background temperature in the upper thermosphere. A number of studies have attributed it to convective instability that limits wave amplitudes and causes the so-called saturation (England et al. 2017; Terada et al. 2017; Vals et al. 2019), while Yiğit et al. (2021b) argued that the inverse relation occurs because a colder background air reduces the scale height H, thus facilitating the exponential growth of amplitude, which is proportional to 1/2H. Molecular diffusion and thermal conduction also exponentially grow with height in response to density decrease, and they eventually exceed all other damping mechanisms in the thermosphere and thereby significantly limit wave growth.

Gravity waves are generated in the lower atmosphere by a variety of mechanisms that vertically displace air parcels (e.g., flow over topography, convection, weather instabilities). While propagating, they are partially filtered out by the background mean wind. Amplitudes of the surviving harmonics grow with height. Ultimately, the harmonics reach altitudes where they are dissipated due to a combination of nonlinear interactions, molecular diffusion, and thermal conduction, and they deposit their momentum and energy into the ambient flow (Yiğit et al. 2008). This GW-mean flow interaction produces acceleration or deceleration of the large-scale circulation, which is often called “GW drag.” Its dynamical importance in the middle and upper atmosphere of Mars has been demonstrated with general circulation models (GCMs), where the effects of small-scale GWs are either parameterized (Medvedev et al. 2011a,b; Yiğit et al. 2018; Gilli et al. 2020; Roeten et al. 2022) or explicitly resolved (Kuroda et al. 2015, 2016, 2019). Such modeling studies must be validated with observations, and many of their utilized parameters need to be constrained. Therefore, an observational characterization of the GW field and its spatio-temporal variation at all heights is of great importance.

The MIR spectrometer of the Atmospheric Chemistry Suite (ACS) experiment on board the Trace Gas Orbiter (TGO; Korablev et al. 2018) allows for this gap in the knowledge of atmospheric variability to be filled by measuring vertical profiles of density and temperature between 20 and 160–180 km in the Martian atmosphere (Belyaev et al. 2022). The algorithm of retrieving GW profiles and their characteristics, together with the first results of its application, have been presented in detail in the work by Starichenko et al. (2021). That database covered observations for 1.5 Martian years (MYs), from the middle of MY 34 to the end of MY35. The same volume of the ACS-MIR/TGO data was recently analyzed by Ji et al. (2023) applying another algorithm of extracting wave components.

In this paper, we present the results on GW activity extending over two MYs: the second half of MY 34, the whole MY 35, and the first half of MY 36. In addition to the ACS-MIR spanning altitudes from 20 to 180 km, we analyzed the profiles measured by another ACS channel, namely, near-IR (NIR; Fedorova et al. 2023). Although the channel covers altitudes only up to ~ 100 km, the number of profiles (several thousand) adds significantly to the overall statistics.

The structure of this paper is as follows. In Sect. 2, observations and temperature retrieval procedures are outlined. The methods of deriving the GW characteristics are described in Sect. 3, while Sect. 4 presents the data coverage. Latitude-altitude distributions for four Martian seasons are given in Sect. 5, and the impact of the GDS on the GW activity is discussed in Sect. 6. Conclusions are presented in Sect. 7.

2 Observations and temperature retrievals

The ACS is a part of TGO, which represents the ESA-Roscosmos ExoMars 2016 collaborative mission. The instrument consists of three infrared channels (Korablev et al. 2018): NIR (0.73–1.6 µm), middle-IR (MIR; 2.3–1.2 µm), and thermal-IR (TIRVIM; 1.7–17 µm). In this work, we use the data obtained from the MIR and NIR spectrometric channels operating in the solar occultation mode since April 2018. The ACS-MIR is a cross-dispersion echelle spectrometer that allows for retrieving vertical temperature and density profiles in the strong 2.7 µm CO₂ absorption band that covers a wide altitude range between 20 and 180 km (Belyaev et al. 2021, 2022). The ACS-NIR is an echelle spectrometer combined with an acousto-optic tunable filter, and it measures the atmospheric structure in the 1.43 µm and 1.57 µm CO₂ bands at altitudes from 10 to 100 km (Fedorova et al. 2020, 2023). The ACS MIR and NIR channels both possess a high-resolution power exceeding ~25 000, a signal-to-noise ratio of more than 1000, and sound the atmosphere with the vertical resolution of ~1 km. During simultaneous occultation measurements, the lines of sight (LOSs) of both instruments targeted identical tangent points. This provides a confidential cross-validation between the retrieved atmospheric profiles. The altitude of the tangent points was determined as the closest distance between the instrument’s LOS and the areoid surface of the planet. The atmospheric transmission spectrum at each tangential altitude was obtained as a ratio of the solar spectrum transmitted through the atmosphere to the reference solar spectrum measured at an altitude where the absorption at the given CO₂ band is negligible, that is, 200 km for the MIR case and 130 km for NIR.

The procedures for retrieving vertical density and temperature profiles from the transmission spectra of CO₂ are extensively described in the papers of Belyaev et al. (2021, 2022) and Fedorova et al. (2020, 2023) for the MIR and NIR channels, correspondingly. The retrieval algorithm simultaneously derives the CO₂ number density and the temperature at each target altitude while assuming hydrostatic equilibrium for pressure. The scheme is multi-iterative, with the well-known Tikhonov regu-larization of each profile. This procedure allows one to avoid possible artificial fluctuations between adjoining points in the altitude grid. Thus, our retrieved GW structures have vertical wavelengths longer than 4-5 km. Points with estimated 1-σ errors exceeding 20 K for MIR and 10 K for NIR have been removed from consideration based on expected wave amplitudes near the top of the vertical domains captured by the channels. Validation between simultaneously measured MIR and NIR profiles has demonstrated a good coincidence below 100 km, with a dispersion of 5–10 K in more than 90% of occultations (Belyaev et al. 2022). When retrieving the GW parameters (see Sect. 3), the data were checked for reasonably smooth background temperature profiles and for physically reasonable values of the wave potential energy (<1000 J kg⁻¹). Anomalous values could be found either at inflection points of temperature profiles, such as around the mesopause, or near the top of the domain (80– 100 km) for the NIR data. Such cases account for about 10–15% of all occultations and were excluded from our consideration.

In most of the simultaneous observations, both the MIR and NIR individual profiles (Figs. 1, 2a, d) along with the retrieved GW parameters (Figs. 2b, c) closely match each other. Nevertheless, in some cases, the evaluated GW amplitudes and potential energy somewhat differ between the two channels (see Figs. 2d, e). The reason for this lies in the different altitude domains used for retrieving the NIR and MIR background temperatures. Consequently, the greatest discrepancy occurs near the upper end (80–100 km) of the NIR profiles (Figs. 2d, e). Overall, statistics of occultations at the MIR 2.7 µm CO₂ band are about 10 times less frequent than those at the NIR CO₂ spectra. Thus, in the analyses to be presented, we complement the MIR profiles with those from NIR whenever measurements with MIR are not available.

Fig. 1 Vertical profiles of temperature derived from simultaneous ACS-MIR (red line) and ACS-NIR (blue line) occultations. The shaded area denotes the uncertainty of the measurements. Examples are from orbits a) 4425n1 (20 Nov. 2018, MY34, Ls=291.7°, Lat=39.2ºS) and b) 8946n1 (25 Nov. 2019, MY35, Ls=111.7°, Lat=60.9°N).

Fig. 2 Wave characteristics retrieved for two representative measurements shown in Fig. 1: orbits 4425n2 (upper row) and 8946n1 (lower row). The red and blue colors distinguish the MIR and NIR data correspondingly. In the left column (a, d), the solid lines represent wave-induced distributions of temperature T', dashed lines are for the wave packet envelopes |T'|, and shades indicate the errors. The middle column (b, e) compares the potential energy profiles, while the right column (c, f) presents the estimated GW drag.

3 Gravity wave characteristics

Separation of the observed temperature profile T(z) into the background $\overline{T}$ and the GW-induced disturbance $T^{\prime }=T-\overline{T}$ is an ambiguous procedure because no unique partition exists. The wave component depends greatly on the definition of the mean temperature $\overline{T}$. The routine used in this work was described in detail and extensively tested in the paper by Starichenko et al. (2021). It was recently applied to the retrieval of GWs in the thermosphere of Saturn (Brown et al. 2022). The mean vertical profile is determined by fitting cubic polynomials within sliding windows of 60 km in width, effectively limiting the consideration to relatively short-scale GW harmonics with vertical wavelengths smaller than 30 km. Tides, technically GWs as well, usually have longer vertical scales. Those with scales shorter than 30 km are attributed by the algorithm to GWs, thus potentially leading to some overestimation of the wave activity. The windows are first shifted from the bottom up and then downward with 7-km steps. Then, all the overlapping polynomials are averaged, and the final profile is smoothed over by a moving average procedure. The uppermost and lowest 4 km of each profile have to be dropped due to a spurious behavior of the fitted polynomials at the edges, which cannot otherwise be averaged.

After the mean and wave components for each profile are derived, the Brunt-Väisälä frequency, wave amplitude, wave potential energy, vertical flux of horizontal momentum, and GW drag can be determined. The Brunt-Väisälä frequency characterizes the convective stability of the atmosphere: $$N^2=\frac{g}{T}\left(\frac{\text{d}T}{ \text{d}z}+\frac{g}{c_p}\right),$$(1)

where g is the acceleration of gravity and c_p is the specific heat capacity at constant pressure. If N² approaches zero (or the temperature gradient approaches the dry adiabatic lapse rate), the stability decreases. When N² drops below zero, the atmosphere becomes convectively unstable and no longer supports GW propagation. Thus, GW harmonics experience strong dissipation and/or breaking in the regions of small or negative N².

Since GW harmonics usually propagate in wave packets, the observed instantaneous peaks and troughs do not fully characterize the wave amplitude. The wave activity is better represented by the envelope for temperature disturbances $\left|T^{\prime }\right|=\sqrt{\overline{T^{\prime 2}}}$. It is calculated by performing the Fourier decomposition in each 60-km sliding window and summation of contributions from all harmonics. The other useful measure of the wave field is the potential energy (per unit mass) E_p: $$E_{\text{p}}=\frac{1}{2}{\left(\frac{g}{N}\right)}^2{\left(\frac{\left|T^{\prime }\right|}{\overline{T}}\right)}^2.$$(2)

While |T'| and E_p describe spatio-temporal distributions of the wave field itself, the vertical flux of horizontal momentum (per unit mass) quantifies the momentum transport by propagating harmonics. This momentum flux is defined as $F=\left(F_x,F_y,0\right)=\left(\overline{u^{\prime }{w}^{\prime }},\overline{v^{\prime }{w}^{\prime }},0\right)$, where u', v', and w' are the components of wave-induced disturbances of velocity in the zonal, meridional, and vertical directions, correspondingly. The directional part of the flux cannot be inferred from a single vertical profile; however, the absolute momentum flux $F=\sqrt{F_x^2+F_y^2}$. (e.g., Ern et al. 2004, Sect. 4) can be estimated: $$F={\displaystyle \sum _{k_h,m}\frac{1}{2}}\frac{k_h}{m}{\left(\frac{g}{N}\right)}^2{\left(\frac{\left|T_{k,m}^{\prime }\right|}{\overline{T}}\right)}^2.$$(3)

The variables k_h and m in Eq. (3) are the horizontal and vertical wavenumbers, and |T'_k,m| is the amplitude of the corresponding harmonic. Summation overall k_h and m is performed. While amplitudes and vertical wavenumbers of particular harmonics are determined by Fourier analysis, the horizontal wavenum-ber cannot be derived from a single vertical profile. Instead, it serves as a scaling factor whose value has to be assigned to provide a meaningful dimension for F. The densest atmospheric footprint at a target point in occultation experiments is ~400– 500 km horizontally; therefore, harmonics with shorter wavelengths remain unresolved. In our calculations, we assigned the uniform horizontal wavelength ${\lambda }_{\text{h}}^{*}=2\pi /k_{\text{h}}^{*}=300 \text{km}$, which is consistent with the shortest possible resolved horizontal scale, and it also agrees with λ_h=100–500 km, which is commonly used in GW parameterizations implemented in numerical GCMs (e.g., Yiğit et al. 2018). The momentum lost by a given breaking and/or dissipating harmonic is transferred to the background flow, thus producing its acceleration or deceleration, or imposing the so-called GW drag $$a_{\text{h}}=\frac{1}{\overline{\rho }}\frac{\text{d}\overline{\rho }F}{ \text{d}z}.$$(4)

In Eq. (4), p denotes the mean density; the subscript h indicates that the acceleration occurs in the horizontal direction. Since the precise direction of F is not known, only absolute values of a_h can be determined from the observations. In the lower parts of the profiles, the amplitudes and the momentum flux F are small. This can lead to large differences of small values in the calculations of a_h, which can result in negative values. In order to avoid this non-physical behavior, we applied to F the iterative procedure described in the paper of Brown et al. (2022, Supporting Information S1).

Results of the thus derived wave characteristics are presented in Fig. 2 for two representative occultations from Fig. 1. For the first example (orbit 4425nl), not only do the MIR and NIR temperature profiles coincide but the wave amplitudes and phases do as well (Fig. 2a). As a result, the envelopes of the wave packet determined from the MIR and NIR data agree well up to ~90 km (Figs. 2b, c). However, for the second profile (orbit 8946nl), the amplitude and estimated potential energy disagree between MIR and NIR above 50 km (Figs. 2d, e), although short vertical-scale features in the temperature profile are resolved well. The altitude distribution of the estimated wave drag agrees well up to ~80 km for both orbits, peaking with 60 m s⁻¹ sol⁻¹ around 80 km (Figs. 2c, f). Above this height, the temperature uncertainties increase significantly for the NIR data. We again emphasize that the drag and momentum flux values represent an upper limit estimate since the scaling factor is proportional to k_h = 2π/λ_h, and the unknown λ_h is equal or longer than the footprint.

4 Description of the data coverage

The analyzed NIR and MIR measurements were taken over two MYs, between the second half of MY34 and the first half of MY36 (May 2018–February 2022). For that period, the MIR 2.7 µm CO₂ band statistics encompasses ~350 occultations in each hemisphere, while the NIR observations were performed approximately 10 times more frequently. The seasonal-latitudinal coverage of individual orbits is presented in Fig. 3 (upper row) as a function of the solar longitude L_s for the northern (left column) and southern (right column) hemispheres separately. Due to the solar occultation mode, the observations were performed either during sunrises or sunsets over morning or evening twilight. However, since the local time (LT) of the solar terminator varies with orbit, latitude, and season, it may reach midday or midnight when closer to the polar regions (Figs. 3a, b).

To analyze the seasonal variability of the GW parameters, we grouped the individual vertical profiles into bins of 3° of L_s and 1 km of altitude. Thus, each bin represents an average of one to seven measured values. Since errors grow with height and profiles extend to different altitudes, the contribution of individual profiles in their top 20 km was weighted by the coefficient ranging from one to zero. The distributions for the GW drag and wave potential energy are plotted in the middle and lower panels of Fig. 3, respectively. The upper altitude spread depends on season and latitude, ranging from 10–20 km to 140–150 km at aphelion and from 20–30 km to 160–170 km at perihelion, as can be clearly seen in the southern hemisphere (Figs. 3d, f). An increased wave activity of up to 300–400 J kg⁻¹ was observed in the winter hemispheres at the mesospheric and thermospheric altitudes. In the summer hemispheres, the peaks of wave activity lie higher, with GW drag reaching maxima near or above the mesopause.

Fig. 3 Seasonal distribution of the GW activity. Upper row: coverage of the ACS measurements as a function of the solar longitude Ls and latitude. The colors indicate the local time of observations. Middle and lower rows: seasonal-altitude distribution for the GW drag and potential energy (per unit mass), correspondingly. The left and right columns present the data for the northern and southern hemispheres, respectively. The gray area denotes the period of GDS. Red dashed lines separate the Martian years MY34, MY35, and MY36.

5 Latitudinal distributions

We next turn to a more detailed examination and consider the altitude-latitude distributions of the GW characteristics. For that, we collected data in 3^° latitude bins and organized the results into four seasons centered on L_s = 0, 90°, 180°, and 270°. They represent two equinoctial seasons (L_s = 0° and 180°) and two solstitial seasons around the aphelion (L_s = 90^°) and perihelion (L_s = 270°). In order to eliminate the influence of the major dust storm of MY34 that occurred between L_s = 188° and 250°, we excluded these measurements. The differences in GW activity introduced by the GDS are explicitly considered in the next section. The cross sections of the GW potential energy are plotted in Fig. 4, where one can immediately notice that the wave activity is stronger in the first half of the MY. The maxima are located at low latitudes in the upper mesosphere and lower thermosphere during the equinoctial L_s = 0^° season (panel a) and shifted to the southern (winter) hemisphere during the solstitial L_s = 90^° season (panel b). A similar pattern is present in the second half of the year, although with a clearly smaller magnitude, especially during the northern winter solstice. We note the symmetry of the E_p distribution with respect to the equator during the equinox.

Figure 5 provides further insight into the climatology of GWs. It presents the latitude-altitude distributions of the zonal GW drag (shaded) along with the mean zonal wind (red contours) simulated with the Martian Atmosphere Observation and Modeling (MAOAM) Martian general circulation model (MGCM; Hartogh et al. 2005; Medvedev & Hartogh 2007) for MY34 and 35 and accordingly averaged. The figure shows that the regions of large GW drag, in general, align with the areas of relatively weak zonal wind, which agrees with the theoretically expected propagation and dissipation characteristics of gravity waves. Harmonics, especially those with relatively slow ground-based horizontal phase speeds c, are substantially damped when their phase speeds approach the mean wind speed $\overline{u}$. This decrease in the intrinsic horizontal phase speed $|c-\overline{u}|$ causes absorption of a significant portion of GWs propagating along the mean wind $\overline{u}$. Harmonics having $c>\overline{u}$ or traveling in the opposite direction to the wind can avoid filtering and propagate higher, grow in amplitude, and eventually break down when wave-induced wind fluctuations |u'| approach the intrinsic phase speed $|c-\overline{u}|$. The weaker the wind $\overline{u}$, the smaller the amplitude |u'| that is required for breaking and/or saturation. This is illustrated by the enhanced momentum deposition in the regions of the weak mean wind shown in Fig. 5.

The equinoctial circulation consists of two prograde (eastward) jets centered in middle-to-high latitudes of each hemisphere and the region of weak winds at low latitudes. The inferred distribution of the estimated GW drag reflects this pattern of inter-hemispheric symmetry. Weaker mean winds in low latitudes allow for GW breaking at lower altitudes. The mean wind changes direction to retrograde (westward) above the mesopause. This causes harmonics traveling westward (c < 0) to break and/or dissipate and deposit their momentum there. In fact, the wind reversal itself is the result of the GW drag (Medvedev et al. 2011a,b, 2013).

The solstitial circulation features eastward (westward) jets in the winter (summer) hemispheres. The jets are stronger during the perihelion solstices, as seen from the comparison ofFigs. 5b and d, due to greater insolation and larger meridional temperature and pressure gradients in this season. Assuming that the GW harmonics excited in the lower atmosphere have a broad range of phase speeds c and travel in all horizontal directions, stronger background winds $\overline{u}$ filter out more waves propagating in the same direction. The remaining harmonics must acquire larger amplitudes in order to break down and, therefore, have to propagate higher. This mechanism explains why the GW drag is localized in a relatively narrow altitude range during the perihelion solstice (Fig. 5d) compared to that during the aphelion (Fig. 5b).

Fig. 4 Latitude-altitude distributions of the retrieved wave potential energy (per unit mass) Ep for four representative seasons centered around (a) Ls = 0°, (b) 90°, (c) 180°, and (d) 270°. The period of the MY34 GDS has been excluded.

6 Impact of the global dust storm

As mentioned previously, the GDS in MY34 occurred between L_s = 188° and L_s = 250°, while the dust load over the same period in MY35 was close to normal, with a minor enhancement between L_s = 235° and 250° (Olsen et al. 2021). Therefore, in order to make a comparison with the dusty MY34, only observations between L_s = 188^° and L_s = 235^° were considered to represent the “dustless” MY35. Figure 6 shows the vertical profiles of potential energy averaged over the southern high-latitude region (southward of 60°S, Fig. 6a), the low and middle latitudes (between 60^°S and 60^°N, Fig. 6b), and the northern polar region (northward of 60^°N, Fig. 6c). The number of individual profiles in each bin varies from 60 to 300 below 100 km and from 7 to 30 above 100 km. In order to determine whether the difference between mean profiles for 2 yr is statistically significant, we conducted Student’s t-test in each bin. The white areas in Fig. 6 indicate regions of confidence where the mean profiles differ with a probability greater than 95%. Essentially, these are the altitudes where standard errors of the mean (drawn in reddish and blue colors) do not overlap. The gray shaded areas show the regions where the robust evaluation of the difference between the means is not possible due to an insufficient number of observations.

The difference between the years shows several systematic features introduced by the dust storm. First, the decrease of the wave activity during the GDS occurs in the high-latitude northern hemisphere below 80 km (Fig. 6c). A similar reduction during the MY34 dust storm, particular on the dayside, agrees well with the Mars Climate Sounder observations in the lower atmosphere below ~30 km (Heavens et al. 2020). Simulations with a high-resolution global circulation model also reproduced the approximately factor of two decrease of the GW potential energy in the lower and middle atmosphere (Kuroda et al. 2020). The decrease is related to the reduction of wave generation caused by convective and baroclinic stabilization of the atmosphere induced by the storm. The same simulations predicted a gradual increase of the GW activity with height such that E_p exceeds the “low dust” values near the top of the model domain at around 80 km. The presented observations very accurately validate this prediction: Wave activity during the GDS surpasses that without the storm exactly above 80 km. Further observational evidence of the enhancement of GW activity in the upper atmosphere during dust storms was provided in the work by Yiğit et al. (2021a) based on NGIMS measurements on board the MAVEN orbiter. These measurements covered the thermo-sphere between 160 and 230 km, which is above the upper limit of the MIR data presented here.

The results in Fig. 6a show an additional feature that was neither observed, modeled, nor otherwise predicted before. It is an increase of the wave activity in the south polar region during the GDS. Thus, high-resolution simulations of Kuroda et al. (2020) demonstrated a stronger reduction of E_p in the southern extratropics, and the MRO-MCS observations of Heavens et al. (2020, 2022) did not see the increase at all. It is worth noting here that the ACS observes relatively long horizontal and short vertical wavelengths, whereas the MRO-MCS, MGS-TES (Pankine et al. 2024), and ODY-THEMIS (Battalio et al. 2023) sample relatively shorter horizontal and slightly vertically broader scales. These waves can be excited by different sources that respond to dust storms differently.

Fig. 5 Same as in Fig. 4 but for the estimated momentum forcing (or GW drag) (shaded). Simulations with the MAOAM MGCM of the mean zonal wind corresponding to the same intervals of Ls are shown with contour lines. The solid and dashed lines represent eastward and westward winds, respectively. Data for the MY34 GDS period have been excluded.

Fig. 6 Altitude profiles of potential energy Ep averaged over three latitude intervals for the period of the GDS of MY34 (Ls=188°–250°, red profiles) and the corresponding low-dust period of MY35 (blue profiles). The colored regions around the profiles show the standard error of the mean. White regions indicate altitudes where the difference between the means was determined with statistical confidence (with a significance level of 95%, according to the Student’s t-test).

7 Conclusions

We presented the results of an analysis of GWs retrieved from solar occultation measurements by NIR and MIR channels of the ACS instrument on board TGO taken over two MYs (mid-MY34 to mid-MY36). The retrieved temperature profiles spanning altitudes up to 100 km (NIR) and 160 km (MIR) were separated into mean component and disturbances, which were used for characterizing the wave field. In particular, we focused on the wave activity represented by wave potential energy E_p and the dynamical impact on the mean flow in terms of the GW momentum deposition (drag). The main inferences of this study are as follows:

• 1.

  Gravity waves are present at all times and in all places in the Martian atmosphere. Within the considered data set, we found no time period or location (except for a few profiles of questionable quality) when and where GW disturbances were absent;

• 2.

  Wave activity is distributed symmetrically with respect to the equator during the equinoctial seasons, whereas the maximum is shifted to the winter hemisphere during solstices;

• 3.

  The maxima of GW drag align with the areas of weak zonal wind along the edges of seasonally varying zonal jets. This feature agrees well with the physics of GW-mean flow interactions;

• 4.

  During the MY34 GDS event, GW activity decreased in the northern hemisphere below 80 km and increased above this height in the high-latitude northern latitudes. This behavior agrees with previous observations and simulations;

• 5.

  The increase of wave activity during the MY34 GDS in the southern polar latitudes below 80 km that we report was not previously observed or predicted.

The climatology of the GW activity and drag in the middle and upper atmosphere based on two MYs of ACS observations confirms theoretical or modeling predictions. However, it also reveals new features (such as the enhancement of the wave activity in the southern polar latitudes below 80 km) that were not anticipated before. Given the dynamical importance of GWs, these features have to be accounted for in numerical models. As ACS continues to perform observations, new data will help further elucidate the spatio-temporal behavior of the GW field.

Acknowledgements

The ExoMars is a joint mission of the European Space Agency (ESA) and Roscosmos. The authors affiliated with LATMOS acknowledge funding from the Centre National d'Etudes Spatiales (CNES) and the Centre National de la Recherche Scientifique (CNRS). Data availability statement The ACS data are available from the ESA Planetary Science Archive (PSA; https://archives.esac.esa.int/psa/%23!Table%20View/ACS=instrument#!Home%20View). The temperature vertical profiles retrieved from ACS-NIR and ACS-MIR measurements are described in Fedorova et al. (2023); Belyaev et al. (2022) and are available at Fedorova (2022); Belyaev (2022), respectively. The most recent MAOAM model output can be accessed at https://mars.mipt.ru. The data supporting our figures are available at https://data.mendeley.com/datasets/drcmwsbhdb/1.

References

All Figures

	Fig. 1 Vertical profiles of temperature derived from simultaneous ACS-MIR (red line) and ACS-NIR (blue line) occultations. The shaded area denotes the uncertainty of the measurements. Examples are from orbits a) 4425n1 (20 Nov. 2018, MY34, Ls=291.7°, Lat=39.2ºS) and b) 8946n1 (25 Nov. 2019, MY35, Ls=111.7°, Lat=60.9°N).
In the text

Fig. 2 Wave characteristics retrieved for two representative measurements shown in Fig. 1: orbits 4425n2 (upper row) and 8946n1 (lower row). The red and blue colors distinguish the MIR and NIR data correspondingly. In the left column (a, d), the solid lines represent wave-induced distributions of temperature T', dashed lines are for the wave packet envelopes |T'|, and shades indicate the errors. The middle column (b, e) compares the potential energy profiles, while the right column (c, f) presents the estimated GW drag.

In the text

Fig. 3 Seasonal distribution of the GW activity. Upper row: coverage of the ACS measurements as a function of the solar longitude Ls and latitude. The colors indicate the local time of observations. Middle and lower rows: seasonal-altitude distribution for the GW drag and potential energy (per unit mass), correspondingly. The left and right columns present the data for the northern and southern hemispheres, respectively. The gray area denotes the period of GDS. Red dashed lines separate the Martian years MY34, MY35, and MY36.

In the text

	Fig. 4 Latitude-altitude distributions of the retrieved wave potential energy (per unit mass) Ep for four representative seasons centered around (a) Ls = 0°, (b) 90°, (c) 180°, and (d) 270°. The period of the MY34 GDS has been excluded.
In the text

	Fig. 5 Same as in Fig. 4 but for the estimated momentum forcing (or GW drag) (shaded). Simulations with the MAOAM MGCM of the mean zonal wind corresponding to the same intervals of Ls are shown with contour lines. The solid and dashed lines represent eastward and westward winds, respectively. Data for the MY34 GDS period have been excluded.
In the text

Fig. 6 Altitude profiles of potential energy Ep averaged over three latitude intervals for the period of the GDS of MY34 (Ls=188°–250°, red profiles) and the corresponding low-dust period of MY35 (blue profiles). The colored regions around the profiles show the standard error of the mean. White regions indicate altitudes where the difference between the means was determined with statistical confidence (with a significance level of 95%, according to the Student’s t-test).

In the text