Model modifications in Schumann resonance intensity caused by a localized ionosphere disturbance over the earthquake epicenter

This paper is a further extension of our latest observations and modeling by Hayakawa et al. (2005a), in which we discovered the anomalous behavior of Schumann resonance observed in Japan, in possible association with the Chi-chi earthquake in Taiwan. Schumann resonance intensity changes associated with a localized decrease in the lower ionospheric height over the earthquake epicenter are modeled. The knee model of the vertical conductivity profile of the ionosphere describes the regular Earth-ionosphere cavity, and the modified knee model is introduced for the disturbance. The localized ionosphere modification is of a Gaussian radial dependence; it has a 1-Mm radius, and the decrease reaches 20 km in the lower ionosphere height over the epicenter of the earthquake (Taiwan). The diffraction problem in the Earth-ionosphere cavity with a localized disturbance is resolved by using the Stratton-Chu integral equation. This solution is constructed for the case of natural resonance oscillations driven by independent random sources distributed worldwide. The data of the Optical Transient Detector (OTD) are used to introduce the source distribution. A pronounced increase in the intensity of the Schumann resonance is obtained around the fourth mode frequency (up to 20%) when thunderstorms are concentrated in Central America. The worldwide distribution of lightning strokes blurs and slightly reduces the effect (15% increase in intensity) for the observer in Japan and the localized nonuniformity positioned over Taiwan. A clear qualitative similarity is obtained in relation to the experimental data, indicating that records collected in Japan may be explained by the impact of a localized decrease in the lower ionosphere over the epicenter of the earthquake. It is admitted that the assumed conductivity decrease could only be caused by a severe change in the ionization in the middle atmosphere. It is not in the scope of this paper to discuss the possible mechanism, but rather Correspondence to: M. Hayakawa (hayakawa@whistler.ee.uec.ac.jp) to show that a closer and quantitative agreement with the experiment yields information about the form and size of the ionospheric modification and about the distribution of global thunderstorm activity during measurements.

main shock.Secondly, the goniometric estimation of the arrival angle of the anomalous signal coincides with the Taiwan azimuth (the unresolved dual direction indicates toward South America).Also, the pulsed signals, such as the Qbursts, were simultaneously observed with the "carrier" frequency around the peak frequency of the fourth Schumann resonance mode.The anomalies were similar for the second earthquake, the Chia-yi.Most likely due to a smaller magnitude, the anomaly appeared one day before and lasted one day after the main shock, with a smaller amplitude enhancement at the fourth Schumann resonance mode than in the case of the Chi-chi earthquake.Yet, the other characteristics were nearly the same, including the goniometric direction finding, frequency shift, etc.Although the emphasis of the present study is on the experimental aspects, a possible generation mechanism for this anomaly is discussed in terms of the ELF radio wave scattered by a conducting disturbance, which is likely to take place in the middle atmosphere over Taiwan.Model computations show that the South American thunderstorms (Amazon basin) play the leading role in maintaining radio signals, leading to the anomaly in the Schumann resonance.The purpose of this paper is to elaborate on the propagation model in our previous paper (Hayakawa et al., 2005a).

Models of regular and disturbed ionosphere and relevant propagation constant
When describing radio wave propagation in a uniform Earthionosphere cavity, we apply the most advanced model: the "knee" model of the vertical conductivity profile of the atmosphere (Mushtak and Williams, 2002).It is a further development of the Greifinger and Greifinger (1978) idea of the ELF propagation depending on two characteristic altitudes in the ionosphere: the "electric" and "magnetic" heights.As was demonstrated by Jones andKnott (1999, 2004), the con-ductivity at intermediate altitudes (so-called "knee" region) is also important, and the knee model incorporates the properties of this region into consideration.One may find the details in Mushtak and Williams (2002).We show a vertical conductivity profile in Fig. 1, to outline how the characteristic heights and ELF propagation parameter are introduced.The conductivity is shown in the standard fashion: the argument, the height over the ground surface, is plotted on the ordinate in km, and the conductivity is shown along the abscissa on logarithmic scale.When obtaining the ELF propagation constant, one operates with two complex parameters, h E and h M , which are found for a given frequency f by using the following equations (Mushtak and Williams, 2002): Here h kn is the "knee" altitude corresponding to the knee frequency f kn .The knee altitude and frequency satisfy the condition: σ kn =2πf kn ε 0 , where ε 0 is the dielectric constant of vacuum, ε 0 =8.859×10 −12 F/m.We show this particular height in Fig. 1 with the asterisks.
The following values were used in our computations: f kn =10 Hz for regular and disturbed profiles, parameter h kn =55 km in the regular ionosphere and h kn =35 km in the disturbed ionosphere.Equations ( 1) and (2) physically mean that at the given frequency f the conductivity and the displacement currents are equal at the "electric" altitude h E .Parameters z a and z b are the local "electric" scale heights of the exponential conductivity profile below and above the h kn , respectively: z a =2.9 km and z b =8.3 km in our computations.The "magnetic" scale height z m is a function of frequency introduced as (Kirillov, 1993): , where z * m =4.0 km, b * m =20.0 km and f m =8.0 Hz.The upper reference altitude h * m =96.5 km was accepted based on the following equation (Kirillov, 1993;Kirillov et al., 1997;Mushtak and Williams, 2002): Here, k=ω/c is the free space wave number, ω is the wave angular frequency, ω 0 is the electron plasma frequency, ν e is the electron collision frequency, z e and z ν are the local scale heights of the electron density and electron collision frequency profiles.
We apply in the computations the following model parameters: in the regular ionosphere f kn =10 Hz, h kn =55 km, z a =2.9 km z b =8.3 km, f * m =8 Hz, h * m =96.5 km, z * m =4 km, and h kn =35 km in the disturbed ionosphere.Profiles in Fig. 1 present relevant regular and disturbed profiles.The position is shown (black dots) of the characteristic heights h E and h M at the initial 4-Hz frequency used in the Schumann resonance studies.Arrows in Fig. 1 indicate how the characteristic values are connected, and the asterisks mark the knee altitudes of regular and disturbed profiles.It is obvious from the above equations that an increase in the signal frequency f will shift the lower altitude h E upward, while the characteristic altitude h M goes down with the frequency (compare with Nickolaenko and Hayakawa, 2002).
A frequency dependent ELF radio wave propagation parameter ν(f ) is obtained from the following expression: where a is the Earth's radius in km.A propagation constant ν(f ) is used in computations of the field, see below.We use the above listed parameters in the model when computing the resonance signals for the Nakatsugawa observatory (35.45 • N and 137.3 • E).Equations (1-4) are similar to the formulas applied by Greifinger and Greifinger (1978), Nickolaenko andRabinowicz (1982, 1987), Sentman (1996), Fullekrug (2000), Nickolaenko and Hayakawa (2002), except for a slight distinction in the procedure of obtaining the characteristic heights.However, the authors of the knee model (Mushtak and Williams, 2002) insist that such details result in the most realistic predictions for the Schumann resonance curves, at least for the case of the uniform distribution of lightning activity over the globe.Therefore, we have chosen this, the most recent and advanced model, for our computations.
We analyze the impact of an ionospheric disturbance associated with pre-seismic and seismic activity; it is introduced in the simplest possible way, as in Hayakawa et al. (2005a).We suppose that modifications are driven by a "terrestrial" process and affect only the lowest part of the conductivity profile: they do not penetrate into the lower part of the region E and higher layers of the ionosphere.Thus, the lower portion of the regular conductivity profile moves downward as a whole by 20 km, stating from the heights below 80 km.The resulting profile is marked "Disturbed" in Fig. 1.The model reduces the real part of the "electric" height, as well as the "knee-height" by 20 km and leaves the rest of the model parameters unchanged.We are fully aware that this is a very severe modification of the middle atmospheric conductivity earthquake-related process could cause.At present, we do not know what kind of earthquake-related process could cause such a severe modification.A discussion of particular mechanisms and the nature of such a modification is beyond the scope of this paper.It is our goal to show that the observed changes in the Schumann resonance intensities many be explained on the basis of such an assumed profile.Diffraction of the ELF radio wave at such a disturbance is able to provide modifications similar to those observed experimentally.

Impact of a localized disturbance on the field
Consider the scattering of ELF radio waves by a localized ionospheric disturbance.We use the known solution of the problem, see Nickolaenko (1984Nickolaenko ( , 1994) ) and Nickolaenko and Hayakawa (2002).The vertical electric field component is treated, since its orientation is independent of the mutual position of the source, observer, and the disturbance.The field at an observatory is the sum of the direct E 1 and scattered E 2 waves (see Fig. 2): (5) The direct (primary) wave propagates from the source to the observatory through the uniform cavity, so that it is found from the following equation: Here I ds(ω) is the current moment of the source, ν is the propagation constant given by Eq. ( 4) for the undisturbed ionosphere, P ν (x) is the Legendre function, and θ H is the source-observer angular distance.
The normalized field disturbance is found from: where where P 1 ν (cos θ ) is the associated Legendre function, γ is the observer-disturbance angular distance, γ i is the distance from the disturbance to the particular source (i=1, 2, 3).
The geometrical parameter M connects the propagation paths, and it is found from the relation M= ∂γ i ∂γ = sin θ i cos γ i •cos α i −sin γ i cos θ i sin γ (Nickolaenko and Hayakawa, 2002).
We suppose that seismic activity locally modifies the lower ionosphere height and thus changes the propagation constant under the disturbance (e.g.Hayakawa et al., 1996;Hayakawa and Molchanov, 2002).The modification is introduced via a complex cosine function, which is found in the regular case from the following relation: The variation δC 2 ν is a measure of the deviation of the disturbed ν D (f ) from the regular ν(f ) value: Both the regular ν(f ) and the disturbed ν D (f ) propagation constants are computed from Eq. ( 4) with characteristic electric and magnetic heights relevant to undisturbed (55 km) and reduced (35 km) knee altitudes.Parameter δC 2 ν is assumed to have a symmetric "Gaussian" radial dependence: where C 2 ν is the maximum and central disturbance, β is the angular distance form the center of ionospheric modification to the point of integration in Eq. ( 7), and d is the characteristic size of the modification measured in radians.We apply in our computations d=π /20, which is equivalent to 1 Mm (1 Mm=1000 km).
One obtains the following formula for the normalized field disturbance after integrating Eq. ( 7) in the case of a compact nonuniformity: In what follows, we show that a localized disturbance over Taiwan provides a noticeable effect on the field: the impact has much in common with the experimental results.Some additional factors must also be addressed.

Position of sources and interference with scattered signals
Let us discuss the general role of the geometry among the source, observer, and the disturbance, including the possible influence of the disturbance on the scattered field amplitude.Initially, we suppose that global thunderstorms are concentrated at point source places in Africa, America, or in Asia.
The assumption allows us to separately estimate their relative contributions to the spectra in the regular and disturbed Earth-ionosphere cavities.
There is an obvious reason why the South American thunderstorms must play an outstanding role in the Schumann resonance observations performed in Japan; it was addressed in detail in the Hayakawa et al. (2005a) paper: the lightning discharges are found close to the antipode of the observatory.We show in Fig. 3 the distance dependence of the amplitude of the vertical electric field at the 26-Hz frequency (the fourth SR mode).The source-observer distance is plotted on an abscissa in Mm and the field amplitude in arbitrary units is shown on the ordinate.A noticeable increase in the field amplitude is observed when the source approaches the antipode of the observatory.Since the geographical coordinates of the Nakatsugawa field site are 35.4• N, 137.5 • E, its antipode (35.4 • S, 42.5 • W) is found close to South America.Keeping in mind the width of the antipode maximum, we must expect a prominent contribution from the South American thunderstorms to the Schumann resonance records collected in Japan.
The great amplitude of the direct radio wave suggests that the amplitude of the wave scattered by a nonuniformity is also great.The interference pattern of direct and reflected waves depends on their mutual phase shift; the latter varies when we alter the positioning of the paths.To illustrate this property, consider a point source located at: 10 • N and 105 • E -Asia, 0 • N and 25 • E -Africa, and 10 • N and 60 • W -America. Figure 4 illustrates the geometry of the relevant paths, in which the observer-disturbance distance is 2.4 Mm.Particular propagation distances of direct waves are equal to 4.3, 12.0, and 14.6 Mm for Asia, Africa, and America, correspondingly.Relevant distances from the source to disturbance are 2.4, 10.7, and 16.8 Mm.The path projection parameters M i =∂γ i /∂γ are equal to 0.99, 0.68, and −0.96.Hence, the position of the American source corresponds to the backscatter of the radio wave from the irregularity, and the terms are summed in Eq. ( 8).As a result, the Q i (ν) amplitude increases, and reflections become more pronounced.An important role of wave backscatter was described by Nickolaenko (1984Nickolaenko ( , 1994)), Nickolaenko and Hayakawa (2002) and Hayakawa et al. (2005a).Thus, the role of a localized nonuniformity is enhanced when the source is found in America, owing to an "optimal" negative and close-to-unity value of the parameter M i .Physically, the effect is conditioned by the "co-linearity" of the propagation paths to the observer and the disturbance, as seen in Fig. 4.
Geometry considerations suggest that American thunderstorms play a decisive role when treating the ionospheric irregularity over Taiwan.It is easy to evaluate the frequency where the interference of direct and scattered waves reaches its maximum.It occurs when the direct and reflected paths deviate by a quarter of the wavelength.The American source (see Fig. 4) corresponds to the approximate equality λ/4=2.4Mm, so that λ=9.6 Mm.Since the wavelength of the basic Schumann resonance mode is equal to the Earth's circumference (40 Mm), the above condition corresponds to the fourth Schumann resonance mode.Thus, we may expect that a signal arriving from American thunderstorms causes pronounced interference around the frequency of 26 Hz.Similar considerations show that the combination of the direct and scattered signals from compact Asian or African sources increases the field significantly at frequencies positioned far above the Schumann resonance.
We plot in Fig. 5 the results of the computations.The left frames depict the power spectra at the Nakatsugawa observatory.Black curves refer to the regular waveguide.Spectra are shown by the curves with black circles in the disturbed cav- ity.The lower left plot illustrates that the disturbance at Taiwan provides a distinct increase in the Schumann resonance signal around the fourth and fifth modes when the source is in America.
The right frame of Fig. 5 depicts the frequency variation of the normalized disturbance defined by Eq. ( 12).It shows how strong the relative amplitude modifications might be in comparison with the regular cavity.The field enhancement may reach 150% in our model.However, such great effects occur in the vicinity of 14-17 and 32-40 Hz, where regular amplitudes are small, owing to the source distances rather than to large reflections.Moderate deviations from the regular case in these frequency bands are seen in the left frames of the figure.
It is clear that the model we use contains the desired property, and it is able to explain the experimental data even when idealistic point source models are applied.The real Schumann resonance signal is driven by many lightning strokes occurring all over the globe at a rate of 100 events per second.Thunderstorm activity moves around the planet during the day, causing variations in the source-observer-disturbance geometry.The level of thunderstorm activity also varies in time.Therefore, it is desirable to include all of these factors into the model by using the time dependent global distribution of lightning strokes.We exploit here the recently acquired space-borne optical observations by OTD.

Model of the worldwide distribution of thunderstorm activity and its daily variations
When computing expected power spectra, we apply the Optical Transient Detector (OTD) data presenting the dynamics of lightning strokes worldwide (Christian et al., 2003;Hayakawa et al., 2005b  of spatial distribution of lightning strokes over the globe.Each map corresponds to a particular hour UT.The flash rate (in strokes per square km) is given in the cells of 2.5×2.5 • grid.The flash rate was converted into the number of flashes within individual cells, and the cells size was extended to a 10×10 • grid.Such an enlargement is insignificant for the accuracy of the Schumann resonance computations, as the cell remains much smaller than the wavelength; simultaneously, the number of elementary sources is substantially reduced, and the computations are accelerated.
The cell centers have coordinates varying along the longitude: −175 • , −165 • , . . ., +175 • (from west to east), and along the latitude: −85 • , −75 • , . . ., +85 • (from south to north).Computations are organized in the following way.We fix the UT hour and pick the particular map of the stroke distribution.In computations, every cell is taken, and the distance θ H is computed from its center to the observer and to the center of ionospheric disturbance γ i , involved in Eq. ( 8).All necessary parameters of propagation paths are also computed, for instance, the current value is found for the derivative parameter M i .

Power spectra expected in uniform and nonuniform cavities with lightning sources distributed worldwide
Afterwards, we compute the direct wave at the given frequency from Eq. ( 6).The wave E 1 travels from the elementary source of a unit current moment I ds(ω)=1 to the observer in the uniform Earth-ionosphere cavity.In parallel, the normalized disturbance B is computed from Eq. ( 12).The disturbed field E D , i.e. the field in the presence of the nonuniformity is found by The power spectrum in the uniform cavity P U i is associated with the sources within the given cell; it is directly proportional to the number of strokes N i in this cell: The same is valid for the nonuniform cavity: The complete power spectra are computed in the regular and disturbed cavities by summing up the contributions from separate cells: and The above computation procedure suggests that the lightning stroke occurrence is subject to the Poisson law.Therefore, mutual time delay of individual discharges has an exponential distribution (the Ehrlang law), so that no pulse interference takes place in the frequency domain: only the energies of elementary excitations are summed in the power spectrum (Nickolaenko and Hayakawa, 2002).As a result, we obtain the current power density at a given frequency f in both uniform (Eq.16) and disturbed (Eq.17) cavities.Time (UT) is a parameter.By varying frequency and repeating computations, we obtain the expected power spectra as functions of frequency for the fixed time.

Results and discussion
We present the results of computations in Fig. 6, where four pairs of spectra are plotted, each corresponding to a characteristic time of a day: UT=03:00, 09:00, 15:00, or 21:00 h.The three latter times correspond to a maximum in the global thunderstorm activity in Asia, Africa, and America, while UT=03:00 h corresponds to the diurnal minimum in the activity.The frequency from 5 to 38 Hz is shown along the abscissa in each frame, and the power density is shown on the ordinate in arbitrary units.Black curves correspond to the spectra in the uniform Earth-ionosphere cavity.Curves with dots show spectra when the disturbance is present over Taiwan.Computations predict that the localized depression of the lower ionosphere tends to increase the amplitude of the Schumann resonance signal, especially at the fourth and fifth modes, provided that the observer is placed in Japan.The intensity of resonance oscillations and spectral pattern vary during the day, which is caused by a change in the level and position of the global lightning activity.Spectra in Fig. 6 must be regarded as a feasible picture during the earthquake, as we do not know exactly which Ann. Geophys., 24, 567-575, 2006 www space-time distribution of thunderstorms was pertinent to a particular day and at the specific moment of observations, and the spectral pattern and its modifications depend on the source position relative to the disturbance and observatory.The character of modifications becomes more apparent when we use the normalized modification of the power spectrum P D (f ) P U (f ) for the fixed times UT shown in Fig. 7.
The ratio of power in the disturbed cavity to that of the regular duct is shown in Fig. 7 for the same set of time moments.Here, the ordinate depicts the ratio, and the signal frequency is shown on the horizontal axis.As we expected from elementary considerations, a substantial effect (reaching 15%) is observed around the fourth Schumann resonance mode, conditioned by mutual positioning of the field site and the disturbance.The worldwide distribution of lightning strokes pertinent to any hour of a day blurs the pronounced resonance enhancement exclusively over the fourth Schumann resonance mode, say, a noticeable increase is also found around the fifth mode, compare with Fig 5 .A "resonance" enhancement is observed when we place the strokes into a relatively small isolated zone, so that relevant geometry conditions are satisfied.Spacious allocation of thunderstorms causes an increase in the field intensity in a wider frequency band in comparison with the case of a single compact source.One might expect "up and down" oscillations around the undisturbed value of 1, but computations show that this is not the case for the globally distributed sources.Figures 6  and 7 demonstrate that the general effect of a localized irregularity depends on the time of observations, as it varies with the motion of the global lightning activity.
After computing the regular and disturbed spectra for every hour of the day, we combine them into the dynamic spectra of Fig. 8. Two frames of this figure present diurnal variations of the field power in the uniform cavity (the left frame) and in the cavity with a nonuniformity (the right plot).An increase in the intensity is apparent at higher Schumann resonance modes.Spectral modification is the most pronounced in the morning hours before the UT noon when Asian thunderstorms provide the major contribution to the electromagnetic signal.Simultaneously, computations show that our model suggests a moderate enhancement, lower than 20%.A closer correspondence between the model data and modifications observed experimentally might be obtained with greater modifications of the conductivity profile and specific source distribution.

Conclusion
This paper is a further extension our latest paper by Hayakawa et al. (2005a).The essential improvements of this paper are (1) the use of a more quantitative knee model for the vertical conductivity profile, and (2) taking the dynamics of global lightning activity into consideration.
To model the impact of the seismic process on ELF radio propagation, we applied the solution of the diffraction problem in the Earth-ionosphere cavity with the localized ionosphere disturbance.The solution of the Schumann resonance problem in the uniform cavity was based on the knee model of the lower ionosphere.A severe localized modification of the ionosphere was introduced as a 20 km depression of the middle atmospheric conductivity profile.The disturbance was centered over the Taiwan earthquake.We assumed that the disturbance had a Gaussian radial dependence of 1 Mm.Wave scattering from the ionosphere irregularity was treated with the help of an analytic solution of the Stratton-Chu integral equation in the Born approximation.The effect of a disturbance is sensitive to the mutual position of the observer, source, and non-uniformity.The greatest impact of a disturbance over Taiwan is observed in Japan around the fourth Schumann resonance mode when sources of electromagnetic waves were placed in Central America.
Since the global lightning activity is distributed over the globe, and it moves around the planet in time, we had to account for these features.A formal solution was constructed for the intensity of the Schumann resonance field in the case of random independent lightning strokes occurring worldwide.Numerical simulations of ELF experimental data exploited the source distributions acquired by the Optical Transient Detector (OTD).We used 24 maps of global lightning activity, each corresponding to a particular hour UT.
Computations showed a pronounced increase in the intensity of the Schumann resonance oscillations around the fourth SR mode (up to 20%) when a point source was placed in Central America.The globally distributed lightning activity blurs the effect, reduces its amplitude to 15%, and increases the third mode oscillations.Our treatment showed that experimental observations of the Schumann resonance in Japan qualitatively agree with the idea of diffraction and scattering of radio waves by a localized depression of the lower ionosphere over the epicenter of the Taiwan earthquake.Additional information is necessary for obtaining a quantitative reciprocity of the measured and model data, information about the particular form and size of the ionosphere modification and especially concerning the space-time distribution of the global thunderstorm activity during measurements.
The major characteristics of the present work have supported our previous conclusion that numerous observations of distortions in the Schumann resonance recorded in Japan prior and during the earthquakes in Taiwan may be caused by the irregularities in the lower ionosphere located over the shock epicenter.

Fig. 1 .Fig. 1 .
Fig. 1.Quiet and disturbed conductivity profiles with a scheme of how the propagation parameters are obtained.

Fig. 3 .Fig. 3 .
Fig. 3. Distance dependence of the amplitude of vertical electric field component at the frequency of 26Hz.Fig. 4. Scheme of ELF radio propagation in the Earth-ionosphere cavity with a localized disturbance placed over the earthquake epicenter.

Fig. 4 Fig. 4 .
Fig.4Scheme of ELF radio propagation in the Earth-ionosphere cavity with a localized disturbance placed over the earthquake epicenter.

Fig. 6 .Fig. 6 .
Fig.6.Regular and disturbed fields at the Nakatsugawa observatory on particular times of the day.

Fig. 8 .Fig. 8 .
Fig. 8. Dynamic spectra in the knee model and the disturbance of 1 Mm radius placed over TaiwanThe left panel refers to the quiet cavity, while the right, disturbed cavity.
574 A. P. Nickolaenko et al.: Model modifications in Schumann resonance intensity Relevant alterations in the ELF propagation constant were obtained at different frequencies and used in the solution of diffraction problem.