Gravity wave propagation in the realistic atmosphere based on a three-dimensional transfer function model

In order to study the filter effect of the background winds on the propagation of gravity waves, a threedimensional transfer function model is developed on the basis of the complex dispersion relation of internal gravity waves in a stratified dissipative atmosphere with background winds. Our model has successfully represented the main results of the ray tracing method, e.g. the trend of the gravity waves to travel in the anti-windward direction. Furthermore, some interesting characteristics are manifest as follows: (1) The method provides the distribution characteristic of whole wave fields which propagate in the way of the distorted concentric circles at the same altitude under the control of the winds. (2) Through analyzing the frequency and wave number response curve of the transfer function, we find that the gravity waves in a wave band of about 15–30 min periods and of about 200–400 km horizontal wave lengths are most likely to propagate to the 300-km ionospheric height. Furthermore, there is an obvious frequency deviation for gravity waves propagating with winds in the frequency domain. The maximum power of the transfer function with background winds is smaller than that without background winds. (3) The atmospheric winds may act as a directional filter that will permit gravity wave packets propagating against the winds to reach the ionospheric height with minimum energy loss.


Introduction
Many researchers have found that severe weather phenomena in the troposphere, such as cold wave (Shen and Zi, 1986;Šauli and Boška, 2001), thunderstorm (Pierce and Coroniti, 1966), typhoon (Shen, 1982;Huang et al., 1985;Xiao et al., 2006) and vortexes (Wan et al., 1998), can excite internal gravity waves which propagate to the ionospheric heights and Correspondence to: W. Wan (wanw@mail.iggcas.ac.cn) bring forth the observable traveling ionospheric disturbances (TIDs).Row (1967) brought forward the pulse propagation theory of gravity waves in an unbounded dissipationless, planar, nonrotating, gravitationally stratified isothermal neutral atmosphere.Row's theory exhibits the major features of the observation at ionospheric heights after the occurrence of a ground impulse source, such as nuclear explosion.Bertin et al. (1975) adopted the reverse ray-tracing method to study the excitation and propagation of gravity waves and proposed that the production mechanism for these waves may involve ageostrophic perturbations of the neutral wind in a jet stream bordering the polar front.Using the transfer function method without background winds, Mayr et al. (1984Mayr et al. ( , 1990) ) simulated the excitation of gravity waves by a source at the polar thermosphere.Wan et al. (1998) statistically analyzed the wave parameters and the exciting sources of the TIDs in Central China based on the observation of an HF Doppler array.They found that the observed TIDs are excited by the tropospheric vortexes, which, in turn, are related to the topography of the Tibet plateau.Mendillo et al. (2002) suggested that the disturbances from the lower atmosphere may cause 10-30% variations in the F2-layer electron density.Rishbeth (2006) reviewed the F-region links with the lower atmosphere and emphasized the role of gravity waves.He pointed out that some of the well-known F2-layer anomalies might have their origins in the lower atmosphere, and it has been suggested that phenomena, such as storms, earthquakes and volcanic explosions, may produce F2-layer signatures.Moreover, Liu et al. (2006) investigated the ionospheric GPS total electron content (TEC) disturbances triggered by the 26 December 2004 Indian Ocean tsunami.They found that the TIDs, which were triggered by the tsunami waves, have periods of 10-20 min and horizontal wavelengths of 120-240 km in the ionosphere.
Conventionally, the researchers adopt three methods to study the propagation characteristics of gravity waves.The first method is the ray-tracing technique (Cowling et al., 1971;Hung et al., 1978;Wan et al., 1998;Yuan et al., 1997).It computes the propagation characteristics of wave energy, Published by Copernicus Publications on behalf of the European Geosciences Union.
L. Sun et al.: Gravity waves propagation in the realistic atmosphere assuming the nondissipative atmosphere.This method is easy to use and its results have clear physical meanings.Furthermore, the method can also deal with the nonplane stratified atmosphere and the time-varying medium.But compared with the realistic atmosphere, the ray theory model is known for its excessive briefness so that it can only be applied to high frequency proximateness.Moreover, the method does not include the phase and cannot calculate the real wave field.Besides, the method cannot provide the distribution characteristics of disturbance in the time-space domain.Recently, an amended ray-tracing technique which can overcome the deficiency of this method, to some extent, is put forward by Ding et al. (2003) through considering the heat conductivity coefficient and viscosity coefficient in the complex dispersion relation of gravity waves.
The second method is the full wave solution method (Francis, 1973;Lindzen et al., 1976;Liang et al., 1998).This method is very complexity.The third one is the transfer function method.Mayr et al. (1984Mayr et al. ( , 1990) investigated the characteristics of gravity waves induced by a source in the polar region.But they did not considered the background winds in their transfer function model, and they also did not study free waves excited by a source in the troposphere and did not investigate the characteristics of medium scale gravity waves induced by the local source, such as the weather process, an earthquake and nuclear detonation.
In this paper, we develop a numerical transfer function model based on the linear theory of gravity waves by considering background winds, to study the excitation and propagation of free gravity waves generated in the troposphere and the response of the ionosphere in a realistic atmosphere.The source on the Earth's surface is regarded as an input signal and the TIDs induced by the exciting source are regarded as an output signal.We consider the atmosphere between the Earth's surface and the ionosphere as a whole response system.We go into detail about the principle of the transfer function model in Sect.2, and the results are presented in Sect.3. Finally, a summary is made in Sect. 4.

The dispersion relation of internal gravity waves
In our study the movement of the atmosphere satisfies the Navier-Stokes equations (Francis, 1973;Ding et al., 2003): (1) where ρ, p, v, T are density, pressure, particle velocity and temperature for the atmosphere.Q is the quantity of heat.λ=6.71×10 −7 T 0.71 0 , µ=3.34×10 −7 T 0.71 0 (DalGarno and Smith, 1962), γ and R 0 are thermal conductivity coefficient, viscosity coefficient, specific heat ratio and universal gas constant, respectively.S is the gas viscosity tensor.g is the gravitational acceleration.The variables in Eq. ( 1) are composed of two parts: the ambient quantities and the first-rank perturbations, while higher rank perturbations are ignored.Hence, the variables can be expressed as where v 0 is neglected by Francis (1973).Under these hypotheses, Eq. ( 1) is subsequently linearized and we apply the Fourier transfer both in time and space.Then we obtain where where H =R 0 T /g is the neutral scale height.The matrix Eq. ( 3) differs from that derived by Francis (1973) in that background zero rank wind v 0 is included here.k h is the real horizontal wave number: k 2 h =k 2 x +k 2 y .k x and k y are the horizontal wave numbers in the x and y directions, respectively.It has a nontrivial solution, so long as the determinant of the coefficient matrix in Eq. ( 3) is zero.Then we can obtain the complex dispersion relation: where For a given real frequency ω and horizontal wave number k h , Eq. ( 5) has three pairs of roots for k z : one pair for gravity waves, the other two for dissipative waves.Each pair of k z corresponds to the upgoing and downgoing waves.In our calculation, only the root for gravity waves is considered.

Numerical transfer function model
According to the Eq. ( 5), the complex dispersion relation of internal gravity waves can be expressed as k z =k z (ω , k h , z) in the moving coordinate system.However, the dispersion relation can also be described as k z =k z (ω, k x , k y , z) in the  The Fourier integral of arbitrary physical quantity V (x, y, z, t) in the static coordinate system, such as temperature, wind speed, is e i(ωt−k x x−k y y−k z z) dk x dk y dk z dω.
Then we integrate k z in Eq. ( 7): The excitation source is at the Earth's surface and the response is at 300-km altitude.The atmosphere is vertically divided into 60 layers with the interval of 5 km.The background temperature and pressure in each layer are set to be constant.Through applying the Fourier transform to temporal and spatial variation of V (x, y, z, t) at the ith layer, according to Eq. ( 8), we obtain   When Eq. ( 10) is divided by Eq. ( 9), we obtain the transfer function at the ith layer So, from the Earth's surface to 300-km altitude the expression formula of the transfer function is In our paper, S(t, x, y, z 1 ) is the exciting source in the troposphere.Through using standard three-dimensional Fourier transform techniques, we obtain the expression formula of the exciting source in the frequency domain.
According to the relationship between the input signal and the output signal, we can get the response at any layer up to 300-km altitude in the frequency domain, Then we apply the inverse Fourier transform to Eq. ( 14).This method yields the temporal and spatial distribution of V (x, y, z, t) at any layer from the Earth's surface to 300-km altitude.
3 The results and discussions

The background parameters
Figure 1 displays the vertical profiles of the atmospheric winds and the temperature adopted in our model.The background wind profiles are calculated by HWM93 wind model (Hedin et al., 1991) and the temperature profile is obtained from the NRLMSISE-00 atmospheric model (Picone et al., 2002) at 06:00 UT on 22 June 2000 at Wuhan   number are most likely to reach the 300-km altitude with less energy attenuation.The amplitude of the transfer function is close to zero in the long-wave-number low-frequency region and in the short-wave-number high-frequency region.This indicates that the gravity waves cannot propagate from the Earth's surface to the 300-km altitude when the horizontal wave length is very long or short, especially when the wave frequency is relatively high or low.The shape of the shadow regions are different for gravity waves propagating along σ =135 • , σ =−45 • and σ =45 • .Through comparing the grey level of the shadow regions with the maximum amplitude in each plot, we find that the amplitude of σ =135 • is the strongest and σ =−45 • is the weakest (the value is 1.0 at z=0).All of the regions with the maximum amplitude basically correspond to wave periods of about 15-30 min and horizontal wave lengths of about 200-400 km.In other words, for the internal gravity waves, the atmosphere behaves like a band-pass filter which offers the easy passage to gravity waves with periods of about 15-30 min and the horizontal wavelengths of about 200-400 km.Furthermore, the gravity waves propagating along σ =135 • are the easiest to reach ionosphere heights, however, σ =−45 • is the most difficult.

The impulse response of the transfer function in the spatial-temporal domain
When severe weather phenomena occur, gravity waves may be triggered in the troposphere and propagate upward in the middle and upper atmosphere.Gravity waves can even reach the ionosphere (Pierce and Coroniti, 1966;Šauli and Boška, 2001).Considering a unit impulse as the excitation origin of the gravity waves, we calculate the impulse response of the atmosphere and present the two-dimensional amplitude distribution at 300-km altitude in Fig. 3.The plus denotes the horizontal position of origin.The positive values in the x axis and y axis are eastward and northward, respectively.The wave amplitude is the relative amplitude which refers to the value of the amplitude at 300-km altitude divided by the maximum amplitude at the lowest layer.In Fig. 1 we have described that the average background winds blow from southeast to northwest at 300-km altitude.Compared with the propagation of gravity waves, the southeast direction of the source is the direction against winds and the northwest is along winds.The northeast and the southwest are approximately perpendicular to the background the winds.Gravity wave disturbance can propagate away from the source toward the periphery (Fig. 3).Mayr et al. (1990) had suggested that the propagation shape of the gravity waves looks like a cone in the spatial domain and propagates as concentric circles at the same altitude under the propagating condition without the background winds.However, taking into account the background winds, we find that the distorted concentric circles develop with amplitude maxima near the source in the two-dimensional environment.This can be due to the effect of background winds which change the horizontal propagation distance of the gravity waves in different directions.The space interval between the adjacent wave crests, the horizontal propagation distance and the distribution of the wave field are different for waves propagating along different directions.For the waves propagating along winds, the horizontal propagation distance of waves are horizontally prolonged, which consequently causes the increase in wave attenuation accumulated on the propagating paths.However, for the waves propagating against winds the propagation distance is shortened, so as to decrease the energy loss.Therefore, the attenuated energy against the winds is smaller than that along winds.Those propagation characteristics of gravity waves can also be found in the low and middle atmosphere, such as at 80-km (Fig. 4), 50 km (Fig. 5) and 20-km altitude (Fig. 6).Unfortunately, we have no observational data to compare with the model output.
The time series of wave field are shown in Fig. 7 at 300km altitude.The horizontal observation point is fixed at (a) x=500 km and y=−500 km, at (b) x=−500 km and y=500 km and at (c) x=500 km and y=500 km, respectively.At the same altitude the amplitude of gravity waves increases rapidly and then decreases gradually with time.The time interval between the adjacent wave crests gradually increases and changes from 15 min to 30 min, which means the wave train includes harmonic waves with a period of 15-30 min.The decrease in wave amplitudes is caused by the increase in viscosity and thermal conduction of the atmosphere, as well as background winds.The attenuation of the wave amplitude against the winds (Fig. 7a) spends a much longer time than the other observation points (Figs.7b and c) and the   amplitude is greater than the others.The results show that the wave fields against the winds can propagate much higher than those along the winds.
Some propagation characteristics of gravity waves mentioned above are similar to those results obtained by the raytracing method (Cowling et al., 1971;Yeh, 1972;Ding et al., 2003).But, as compared with ray-tracing method which can only be used to study the monochromatic wave, the transfer function method is able not only to investigate the propagation of the wave packets but also to provide the distribution characteristic of the whole wave fields.

The filter characteristics of the transfer function
As stated above, the transfer function method can obtain the distribution of the whole wave fields in the spatial-temporal domain.In this part we will discuss the response curve of the transfer function in the frequency domain (Fig. 8), wave number domain (Fig. 9) and direction domain (Fig. 10), respectively.
The response curve of the transfer function in the frequency domain is shown in Fig. 8 at 300-km altitude.The solid line and the dotted line, respectively, denote the atmosphere with background winds and without background winds.As illustrated in the figure, the maximum power of the transfer function with background winds (including the filter effect of the background winds in different directions)     is smaller than that without background winds.This implies that the background winds filter some frequency components of waves propagating upward.The domain with the maximum power mainly focuses on the frequency domain of about 2-5 cycles/hour.It describes the characteristic of the band-pass filter for the gravity waves propagating from the Earth's surface to 300-km ionospheric altitude.Furthermore, the center frequency corresponding to the maximum power has obvious deviation under the control of the background winds.Figure 9 describes the filter characteristics of the transfer function in the wave number domain at 300-km altitude.The maximum power domain corresponds to the horizontal wave number of about 2.5-6 cycles/1000 km.It also indicates that the other components of gravity waves have been filtered by the viscosity and thermal conduction of the atmosphere, as well as the background winds when the gravity waves propagate upward to higher altitudes.Figure 10 depicts the filter characteristics of the transfer function in the direction domain at 300-km altitudes.The thick solid line and the thick dotted line, respectively, denote the distribution of the power in the atmosphere with background winds and   without background winds.The arrowhead shows the direction of the average background winds.The wave power is symmetrically distributed in all the directions without background winds.However, the dominating power concentrates on the southeast direction with background winds.Comparing with Fig. 3, we find that the gravity wave power mostly converges on the anti-windward direction.The main reason is that the upward propagating gravity waves will be affected by the background winds.
The filtering effects of the transfer function is also notable in Fig. 11 at x=500 km and y=−500 km (solid line), x=−500 km and y=500 km (dashed line) and x=500 km and y=500 km (dotted line) at 300-km altitude, respectively.For the fixed observation points, the center frequency along winds or against the winds is different from that perpendicular to the winds; moreover, the maximum power against the winds is stronger than that along the winds.Figure 12 shows the power distribution of the transfer function at the time of 60 min, 120 min and 180 min at 300 km altitudes (gravity wave excited on the surface when T=0 min).The power concentrates on the concentric circles at different times.The maximum power corresponds to the wave number of about 5 cycles/1000 km at 60 min.Then the wave number corresponding to the maximum power gradually decreases with time, which indicates that the wave length of gravity waves becomes longer.Before 120 min the maximum power is stronger than that after 120 min.The computed results also show the filter effect of the transfer function.Through analyzing the response curve of the transfer function from Fig. 8 to Fig. 12, we can conclude that the gravity waves travel easily in the anti-windward direction.This is because the propagation distance of gravity wave packets along the winds is longer than that against the winds, which leads to more energy loss when gravity waves are propagating along the winds than against the winds at the same altitude.The results show that the atmospheric winds may act as a directional filter that will permit gravity wave packets propagating against the winds to reach the ionospheric height with minimum energy loss and minimum travel time.Therefore, with directional filtering by the winds, the action of the atmospheric processes seems to favor those gravity wave packets that take the minimum time to reach the ionosphere.

Conclusions
Under the linear propagation theory of gravity waves, the transfer function model is built to study the vertical propagation of gravity waves.Our numerical results not only confirm the main conclusions obtained by ray-tracing method, e.g. the trend of gravity waves to travel in the anti-windward direction, but also yield some new propagation characteristics of gravity waves: (1) The method provides the distribution characteristic of the whole wave fields and the time series of the propagation of gravity waves.The impulse response of the transfer function in the spatial-temporal domain exhibits the propagation of the gravity wave packets in the way of the distorted concentric circles which is controlled by the background winds.
(2) The frequency and wave number response curve of the transfer function show the filter characters, where we can pick out the range of frequency and the wave number of the gravity waves (about 15-30 min periods and about 200-400 km horizontal wave lengths) which will reach ionospheric heights.The other frequency components of the waves have been filtered by the viscosity and thermal con-  duction of the atmosphere, as well as the background winds.Moreover, the center frequency of the transfer function has an obvious shift with the winds and the maximum power with the background winds is smaller than that without background winds at the same altitude.
(3) The direction response curve of the transfer function indicates that the atmospheric winds may act as a directional filter that will permit gravity wave packets propagating against the wind to reach the ionospheric heights with minimum energy loss and minimum travel time.Therefore, with directional filtering by the winds, the action of the atmospheric processes seems to favor those gravity wave packets that take the minimum time to reach the ionosphere.The maximum power of the fixed observation point in the antiwindward direction is stronger than that in the windward direction.In the development of the paper we consider the simulation of a real origin of the gravity waves, e.g.typhoon or cold wave, and its propagation process in the ionosphere.Through the comparison between the observational and modeled results, we hope to make a further study on the coupling mechanism between the lower atmosphere and ionosphere.Furthermore, the nonlinear interaction between the gravity wave and the atmosphere plays an important role in the gravity waves propagation (Fritts and Alexander, 2003), which will also be discussed in the next work.

17Fig. 1 .
Fig.1.The altitude profile of atmospheric wind speed (a, b) and temperature (c) respectively obtained from HWM93 wind model and NRLMSISE-00 atmospheric model at 0600 UT on June 22, 2000 at Wuhan (30°N, 114°E).The eastward wind and the northward wind are denoted by v0x and v0y respectively.The dotted lines in figure a, b show the average wind which blows from the southeast to the northwest.

Fig. 1 .
Fig. 1.The altitude profile of atmospheric wind speed (a), (b) and temperature (c), respectively, obtained from the HWM93 wind model and NRLMSISE-00 atmospheric model at the 06:00 UT on 22 June 2000 at Wuhan (30 • N, 114 • E).The eastward wind and the northward wind are denoted by v 0x and v 0y respectively.The dotted lines in panels a, b show the average wind which blows from the southeast to the northwest.

18Fig. 2 .
Fig.2.Contours of transfer function amplitude for gravity waves propagating at 300 km altitudes in the frequency-wavenumber domain.Fig.2a, b and c respectively describe the amplitude distribution of transfer function at the azimuth of σ=135° (against wind), σ=45° (the perpendicular direction of background winds) and σ=-45° (along wind).σ is the azimuth angle clockwise from y direction.

Fig. 3 .
Fig.3.The two-dimensional spatial distribution of amplitude for gravity waves at 300 km altitudes, which was triggered by the unit impulse source at t=0 min.The figures present the distribution of wave fields at time (a) t=60 min, (b) 120 min and (c) 180 min.The plus denotes the position of the unit impulse.X axis and y axis describe the eastward and the northward, respectively.The amplitude is the relative amplitude compared with that in the lowest layer.

Fig. 2 .
Fig. 2. Contours of transfer function amplitude for gravity waves propagating at 300-km altitudes in the frequency-wave-number domain.Panels (a), (b) and (c), respectively, describe the amplitude distribution of transfer function at the azimuth of σ =135 • (against wind), σ =45 • (the perpendicular direction of background winds) and σ =-45 • (along wind).σ is the azimuth angle clockwise from the y direction.

Fig. 2 .
Fig.2.Contours of transfer function amplitude for gravity waves propagating at 300 km altitudes in the frequency-wavenumber domain.Fig.2a, b and c respectively describe the amplitude distribution of transfer function at the azimuth of σ=135° (against wind), σ=45° (the perpendicular direction of background winds) and σ=-45° (along wind).σ is the azimuth angle clockwise from y direction.

Fig. 3 .
Fig.3.The two-dimensional spatial distribution of amplitude for gravity waves at 300 km altitudes, which was triggered by the unit impulse source at t=0 min.The figures present the distribution of wave fields at time (a) t=60 min, (b) 120 min and (c) 180 min.The plus denotes the position of the unit impulse.X axis and y axis describe the eastward and the northward, respectively.The amplitude is the relative amplitude compared with that in the lowest layer.

Fig. 3 .Fig. 4 .Fig. 4 .Fig. 4 .Fig. 5 . 19 Fig. 4 .Fig. 6 .
Fig. 3.The two-dimensional spatial distribution of amplitude for gravity waves at 300-km altitudes, which was triggered by the unit impulse source at t=0 min.The panels present the distribution of wave fields at time (a) t=60 min, (b) 120 min and (b) 180 min.The plus denotes the position of the unit impulse.The x-axis and y axis describe the eastward and the northward, respectively.The amplitude is the relative amplitude compared with that in the lowest layer. 20

Fig. 7 .
Fig.7.Time series of gravity waves amplitudes at 300 km altitudes.The obs point is horizontally fixed at (a) x=500km and y=-500km, at (b) x=-500 y=500km and at (c) x=500km and y=500km.The relative amplitude is enlarg times.

Fig. 7 .
Fig. 7. Time series of gravity wave amplitudes at 300-km altitudes.The observation point is horizontally fixed at (a) x=500 km and y=−500 km, at (b) x=−500 km and y=500 km and at (c) x=500 km and y=500 km.The relative amplitude is enlarged 1000 times.

Fig. 8 .
Fig.8.The filter characteristics of the atmosphere in the frequency domain at 300 km altitudes.The solid line and the dotted line respectively denote the power distribution of the atmosphere with background winds and without background winds in the frequency domain.

Fig. 9 .
Fig.9.The filter characteristics of the atmosphere in the wave number domain at 300 km altitudes.The solid line and the dotted line respectively denote the power distribution of the atmosphere with background winds and without background winds in the wave number domain.

Fig. 8 .
Fig. 8.The filter characteristics of the atmosphere in the frequency domain at 300-km altitudes.The solid line and the dotted line, respectively, denote the power distribution of the atmosphere with background winds and without background winds in the frequency domain.

Fig. 9 .
Fig.9.The filter characteristics of the atmosphere in the wave number domain at 300 km altitudes.The solid line and the dotted line respectively denote the power distribution of the atmosphere with background winds and without background winds in the wave number domain.

Fig. 9 .
Fig.9.The filter characteristics of the atmosphere in the wave number domain at 300-km altitude.The solid line and the dotted line, respectively, denote the power distribution of the atmosphere with background winds and without background winds in the wave number domain.

Fig. 10 .
Fig.10.The filter of the atmosphere in the direction doma altitudes.The thick solid line and the thick dashed line respectively deno distribution of the atmosphere with background winds and without back in the direction domain.The arrowhead shows the direction of backgroun radius shows the power value.

Fig. 10 .
Fig. 10.The filter characteristics of the atmosphere in the direction domain at 300-km altitude.The thick solid line and the thick dashed line, respectively, denote the power distribution of the atmosphere with background winds and without background winds in the direction domain.The arrowhead shows the direction of the background winds.The radius shows the power value.
22distribution of the atmosphere with background winds and without background winds in the direction domain.The arrowhead shows the direction of background winds.The radius shows the power value.

Fig. 11 .
Fig.11.At 300 km altitudes the power-frequency characteristics of the fixed observation point at x=500 km and y=-500 km (solid line), x=500 km and y=500 km (dotted line) and x=-500 km and y=500 km (dashed line).

Fig. 11 .
Fig. 11.At 300-km altitude the power-frequency characteristics of the fixed observation point at x=500 km and y=−500 km (solid line), x=500 km and y=500 km (dotted line) and x=−500 km and y=500 km (dashed line).

Fig. 12 .
Fig.12.The power distribution of the atmosphere in the wave number domain at 300 km altitude at (a) 60 min, (b) 120 min and (c) 180 min.

Fig. 12 .
Fig. 12.The power distribution of the atmosphere in the wave number domain at 300-km altitude at (a) 60 min, (b) 120 min and (c) 180 min.