Variations of the ionospheric electron density during the Bhuj seismic event

Abstract. Ionospheric perturbations by natural geophysical activity, such as volcanic eruptions and earthquakes, have been studied since the great Alaskan earthquake in 1964. Measurements made from the ground show a variation of the critical frequency of the ionosphere layers before and after the shock. In this paper, we present an experimental investigation of the electron density variations around the time of the Bhuj earthquake in Gujarat, India. Several experiments have been used to survey the ionosphere. Measurements of fluctuations in the integrated electron density or TEC (Total Electron Content) between three satellites (TOPEX-POSEIDON, SPOT2, SPOT4) and the ground have been done using the DORIS beacons. TEC has been also evaluated from a ground-based station using GPS satellites, and finally, ionospheric data from a classical ionospheric sounder located close to the earthquake epicenter are utilized. Anomalous electron density variations are detected both in day and night times before the quake. The generation mechanism of these perturbations is explained by a modification of the electric field in the global electric circuit induced during the earthquake preparation. Key words. Ionosphere (ionospheric disturbances) – Radio Science (ionospheric physics) – History of geophysics (seismology)


Introduction
Since the great Alaskan earthquake in 1964 (Leonard and Barnes, 1965), many evidences of ionospheric perturbations after strong earthquakes have been reported (see, for example, the review by Blanc, 1985).The perturbations occur just after the shock and are due to acoustic waves which are amplified through the atmosphere because of the decreasing atmospheric density with increasing height.
Correspondence to: M. Parrot (mparrot@cnrs-orleans.fr)Furthermore, it has been shown that electric and magnetic modifications could occur between a few hours and a few days before earthquakes in the seismic zone.Ionospheric perturbations have also been observed a few days before above the seismic zone.However, it must be noted that not all earthquakes produce such phenomena.Experiments close to some earthquakes have not recorded perturbations.This is a common factor of all precursors due to the complex nature of earthquakes.Some examples of these ionospheric perturbations can be found in Sobolev and Husamiddinov (1985), Fatkullin et al. (1989) and in the monographs edited by Hayakawa and Fujinawa (1994), and Hayakawa and Molchanov (2002).Such ionospheric anomalies are better detected during the night when the ionosphere is calm (see the review by Parrot et al., 1993).Increases as well as decreases of the critical frequencies are observed in the D-, E-and F-regions before earthquakes.Although this phenomena is not well understood, it could be related to a redistribution of the electric charges at the surface of the Earth and then in the Earth's atmospheric system.Other hypotheses concerning the mechanism of these perturbations are given in Parrot et al. (1993), Molchanov (1993), Parrot (1995), and Pulinets et al. (2003).Ionospheric effects of aerosols and metallic ions emitted in the atmosphere above seismic regions are considered in Pulinets et al. (1994Pulinets et al. ( , 1997Pulinets et al. ( , and 1999)).Pulinets et al. (2000) calculated the electric field generation on the Earth's surface and in the atmosphere.It was shown that the additional flux of metallic aerosols leads to the strengthening of an anomalous electric field.This electric field (in local region inside atmospheric layers) was evaluated.It was found that the night-time field penetration efficiency is larger than during the daytime and that it depends upon the size of the vertical field E z localization layer.Pulinets et al. (2002) introduces the concept of atmospheric plasma appearing in the near ground layer of the atmosphere under the action of radon radioactivity in seismically active regions.Such plasma could trigger or perturb different kinds of EM emissions at various frequencies observed experimentally before strong earthquakes.The paper focuses on the seismogenic variations on the higher altitude-upper ionosphere and magnetosphere-small-and large-scale irregularities within F-layers of ionosphere and magnetospheric duct formation.Pulinets and Legen'ka (2002) revealed a strong modification of the equatorial anomaly structure prior to strong earthquakes by the analysis of topside ionograms obtained on board Intercosmos-19 satellite.Equatorial anomaly (EA) may increase (deepen) or attenuate (be filled up), depending upon the local time.The most important parameter indicating that EA is modified prior to earthquake is theB2u parameter.B2u is the half thickness of the topside part of the Ne(h) profile of the F2-layer.The EA modification results on the basis of satellite data are confirmed by the data of groundbased observations.Pulinets and Legen'ka (2003) have shown that seismic ionospheric disturbances are strongly time dependent before the beginning of the main shock.Seismic ionosphere disturbances are generated weekly, several days before the first shock, but at that moment the displaced region is not located above the epicentre, but rather displaced from it.These disturbances can be transferred along magnetic field lines into the conjugate regions in the opposite hemisphere.Three earthquakes of different latitudes (High latitude, Alaska; Mid latitude, Central Italy; and Low latitude, New Zealand) have been considered.Moreover, data of vertical sounding of the ionosphere performed by the Alouette-1 and Intercosmos-19 satellites have been taken, along with the data of vertical sounding of the ionosphere by ground-based instruments.
Variations of the critical frequencies of the ionospheric layers were mainly observed with ground-based ionospheric sounders, but measurements of Total Electron Content (TEC) by satellites can also be used (Evans, 1977).The TEC gives the sum of the electron density between the altitude of the satellite and the ground.Therefore, this parameter is mainly related to the density in the F-layer which is larger than in the others.Calais and Minster (1995) have used this method to detect perturbations due to an earthquake.A statistical study with the TEC data from TOPEX-POSEIDON has been performed by Zaslavski et al. (1998) to check correlations between ionospheric perturbations and seismic activity.This work is similar to the work of Liu et al. (2002), who used a method to derive the ionospheric TEC from data recorded by a network of the Global Positioning Systems (GPS) in Taiwan.Based on the network data the latitude time-TEC (LTT) plots are constructed and showed that 1-4 days prior to three earthquake onsets, TEC values decrease and anomaly crests move towards the equator.Further, the simultaneously deduced overhead TEC and the foF2 observed at Chung-Li during the three earthquakes are compared and examined.The recent publication of Liu et al. (2004), where a statistical analysis of TEC anomalies has been performed before strong earthquakes (M>6) at Taiwan, claims a leading time of 5 days before the seismic shock for the ionospheric precursors.
All these concerns about the ionospheric anomalies have been compiled in a review paper by Pulinets et al. (2003).
Thus the aim of this paper is to study the ionosphere and associated parameters recorded by various experiments above the Gujarat active seismic zone around the time of the Bhuj earthquake.A workshop report about this event can be found in Singh and Ouzounov (2003).The experimental equipments and the data selection are described in Sect.2, whereas the results are presented in Sect.3. Discussions and conclusions are given in Sect. 4.

The Bhuj earthquake
The Bhuj earthquake shook the Indian province of Gujarat on the morning of 26 January 2001 at 03:16 GMT (08:46 LT).This devastating earthquake had a magnitude M s =7.6, and a depth of approximately 24 km.The exact location of the epicenter was 23.4 • N, 70.32 • E, approximately at 20 km from Bhuj (population 150 000).Since this earthquake occurred in a heavily industrialised region of India, it caused a lot of casualties.Numerous aftershocks exceeding M s =5 have been reported and one aftershock had a magnitude M s =5.9.It took place on 28 January 2001 at 06:32 IST.Up until 31 January 2001 hundred aftershocks with magnitudes ranging between 2.8 and 5.9 have occurred.The more important ones are displayed in the bottom panel of Fig. 4.

The TEC data from TOPEX-POSEIDON
We used the data of the TOPEX/POSEIDON satellite which is dedicated to space oceanography and has two altimetric systems for measuring the altitude of the satellite above the sea surface.Measurement of fluctuations in the integrated electron density or TEC (Total Electron Content) between the satellite and the ground is used to remove the ionospheric influences in order to obtain a better estimation of the altitude (Monaldo, 1993).Typical examples of TEC variations are shown, for example, in Robinson andBeard (1995), andZlotnicki (1994).
The TOPEX-POSEIDON satellite carries a dual frequency altimeter operating at frequencies of 13.6 GHz (Ku-band) and 5.2 GHz (C-band).The data from the dual frequency altimeter are obtained with a good time resolution (1 s) but unfortunately they are not valid above land.So, in this study we have used the data from the DORIS system (Fleury et al., 1991;Escudier et al., 1993) which is valid everywhere.However, the amplitude resolution is low because the DORIS system uses a set of ground-based stations (∼50) located at the surface of the Earth to determine the vertical sub-satellite values.
Our data are from the AVISO CD-Roms (AVISO, 1992).We used the ionospheric correction for the altitude measurement which is expressed in millimetres.This quantity is directly proportional to TEC.The time resolution of the data is one point by second.At Ku-band, this parameter can reach up to −250 mm during high solar activity.It must be noticed that this negative value, in fact, corresponds to a positive  1).increase of electron density.In the AVISO CD-Roms, expected ionospheric corrections from the Bent model are also given.The altitude of the TOPEX-POSEIDON satellite is 1340 km and it covers the geographic range from 66 • N to 66 • S. A particularity of this satellite is that its ground tracks are always identical, but it is not a Sun-synchronous satellite.Therefore, a study of parameter variations versus the local time can be performed at a given place.The satellite returns every 10 days at the same location, and at the equator, the distance between 2 ground tracks is ∼300 km, whereas the distance between 2 consecutive orbits is ∼3000 km.The data have been studied during one month around the earthquake time.Figure 1 presents the ground tracks of the satellite TOPEX-POSEIDON around the Bhuj location.In the middle, the star indicates the position of the epicentre.The satellite regularly returns on the same tracks which have been labelled A, B, C, and D. Tracks A and C are upgoing tracks, whereas B and D are downgoing tracks.The days and average hours in UT (LT=UT+5.5h) of each passes are given in the Table 1.It must be noted that tracks A and C correspond to daytime (the ionospheric correction is higher because the TEC is much important) whereas tracks B and D correspond to nighttime.Two months of data centred on the time of the earthquake have been considered.

2.3
The TEC data from the satellites SPOT2 and SPOT4 SPOT2 and SPOT4 are two satellites devoted to the Earth's observation.The altitude of these satellites is ∼810 km and their orbit inclination is 98 • .They have both a DORIS system similar to the one on board TOPEX-POSEIDON which was presented in the above section.However, the presentation of the ionospheric data is different in such a sense that  where a is a constant (40.22/c in m 2 /(electron.s)),c is the light velocity, dt is the counting period, CET 1 ver and CET 2 ver are the vertical TEC at the beginning and at the end of the counting period, respectively, and F is the emission frequency.The geometric factor k of correction is given by where r is the position radius of the subionospheric point (r=R T +h m +50 km), R T is the Earth radius, h m is the subionospheric altitude (altitude of the maximum of electron density) , and E is the elevation angle of the satellite relatively to the station.As a slant TEC estimation is done, the subionospheric point is the point at an altitude equal to 350 km which belongs to the line between the satellite and the ground-based DORIS station.It gives the best position of the TEC estimation because 350 km is the average altitude of the maximum of the electron density.Therefore, the parameter Iono given by Eq. ( 1) is more related to a TEC time variation than to the TEC parameter itself.
Figure 2 represents the positions of these subionospheric points for each passes of the satellite SPOT2 close to the DORIS stations which are around the epicenter of the Bhuj earthquake on 24 January 2001.Using Eq. (1) with the method described in Fleury et al. (1991), the TEC parameters have been estimated for each satellite passes.This method only assumes a latitudinal variation of TEC.

2.4
The TEC data from the Bangalore ground-based station We use data from a ground-based station located at Bangalore (Lat.13.021, Long.77.57).This GPS station, named IISC, belongs to the NASA Flinn global network for continuous GPS analysis.Bangalore is quite far from the epicenter (1400 km), but it falls in the area of this earthquake preparation which depends on the magnitude (Dobrovolsky et al., 1979).

The ionospheric data from the Ahemadabad groundbased sounder station
We used the data of Ahemadabad ionospheric station (23.02 • N, 72.6 • E) located at 220 km from Bhuj.Data has been collected from the SPIDR (Space Physics Interactive Data Resource) NOAA Satellite and Information Services (http://spidr.ngdc.noaa.gov/).The maximum time resolution of these data is 1 h.Two parameters which characterise the F 2 ionospheric layer have been taken into account: the characteristic frequency foF2 and the virtual height h p F 2 .

The results
Figure 3a shows the variations of the TOPEX-POSEIDON ionospheric data recorded along the tracks A and C (see Table 1) as a function of the latitude.It corresponds to data recorded in the daytime but not at the same hours, because there is a time shift due to satellite passes.It is shown that the ionospheric correction increases for days which are close to 26 January, the earthquake's day, and decreases as far as we move away from this date.On 24 January, around 09:00 UT (14:30 LT), the ionospheric correction is much more important.Although all these data are recorded in local daytime (see Table 1), this tendency could be attributed to a normal time variation because TEC values are usually maximum just after noon in local time.Then a comparison has been done values of the ionospheric corrections as function of the latitude (see also Table 1), (b) differences between the ionospheric corrections and the Bent model as a function of the latitude.
with a model, and the corresponding differences with the Bent ionospheric model are plotted in Fig. 3b.Looking for differences with an absolute value larger than 20 mm, it is again observed that largest variations occur for days close to the earthquake's day, except on 13 February.Although it is a date far from the earthquake's day, the data recorded on 13 February depict a large difference with the ionospheric model.This can be explained by the magnetic activity which reaches during this day its maximum value over the time interval of two months, as it is shown in the last panel of Fig. 10.The results from SPOT2 and SPOT4 are given in Fig. 4. We have only considered the data from the KITB station which is closest to the epicentre (see Fig. 2).The top panel of Fig. 4 shows the variation of the ionospheric correction given by Eq. ( 1) for the two satellites passes during an interval of two months around the time of the Bhuj earthquake.The passes occurred between 11:00 and 13:00 LT, and for each pass the data have been averaged for latitudes between 25 • and 33 • , in order to obtain a single measurement per day.The horizontal full line indicates the mean value over the two months and the dashed lines indicate the variance.It is observed that an important TEC time variation occurred on 24 January.It corresponds to a SPOT2 pass which is in red, in Fig. 2.
Around the time of the Bhuj earthquake, another important TEC variation has been recorded on 29 January.This last anomalous variation can be attributed to aftershock activity, as it can be observed on the bottom panel of Fig. 4. Other important TEC variations have been recorded at the end of February which can be imputed to our long time interval of observations.Two factors can be involved: i) a time shift of the observations which took place around 11:00 LT at the beginning of January and around 13:00 LT at the end of February (this last local time corresponds, on average, to a maximum of TEC), and ii) a seasonal variation (Codrescu et al., 1999).The middle panel of Fig. 4 represents the TEC parameters obtained by the method mentioned in Sect.2.3.As in the top panel, the TEC values correspond to registrations between 11:00 and 13:00 LT, and they are averaged for latitudes between 25 • and 33 • .A maximum value is observed on 20 January, six days before the earthquake.Explanations concerning other values exceeding the variance at the end of February are similar to the ones given for the top panel.
Figure 5 represents the TEC data recorded at Bangalore station from 16 January until 30 January.In this time interval, two anomalous variations appear on 21 and 25 January in the afternoon -evening sector.It was shown (Pulinets, 1998) that a few days before the strong earthquakes the ionosphere variability increases.The variability of the ionosphere can be estimated by different ways (Forbes et al., 2000;Rishbeth and Mendillo, 2001) but the most appropriate and accepted is the value of the critical frequency deviation from the monthly median: where foF2 is the critical frequency, and f oF 2 is the monthly median.Statistical analysis shows that a 25% threshold is an indicator of the "normal" state of the ionosphere.Figure 7, for critical frequency deviations before the time of Bhuj earthquake, demonstrates the increase of the variability starting from 21 January (5 days before the seismic shock).
In Liu et al. (2004), a statistical process was developed to reveal possible ionospheric precursors: the interquartile range IQR is calculated to construct the upper bound X+I QR and lower bound X−I QR for the data estimation where X is the 15-day running mean.Here instead of X we use the monthly median.Liu et al. (2004) estimated IQR as 1.34σ , where σ is the standard deviation.To make the condition more rigorous we took in our processing a 2σ threshold.Figure 8 shows the results of the processing for 23 January.The black line with points represents the real measurements of the critical frequency, the green line represent the monthly median, and the dashed red and blue lines represents the upper and lower bounds of f oF 2+2σ and f oF 2−2σ , respectively.One can clearly see from the figure that the value at 10:00 UT is out of the 2σ upper bound, and that at 17:00 and 18:00 UT the value is out of the 2σ lower bound.They both show anomalous deviations of the critical frequency.To demonstrate this anomalous behavior we bring the data for the variations of the critical frequency at 10:00 UT in Fig. 8.For the month of January 2001 the value of the critical frequency at 10:00 UT is within the range between 11.0 and 11.5 MHz, and only a few days before the seismic shock one can see the sharp drop and then gradual increase of the critical frequency.This behavior is completely different in comparison with the whole month's behavior.As the epicenter was close to the station, it must be noted that a large gap is observed just after the quake.In February, fewer variations are observed that can be attributed to the aftershock activity (see the bottom panel of Fig. 4).
All critical ionospheric frequencies foF2 recorded by the Ahemadabad station during January and February 2001 are shown in Fig. 9.It represents the hourly variations of foF2 as a function of the days, which are color-coded according to the scale on the right.The missing data have been forced to zero values.An interesting information is contained in the station records: at 20:00 UT on 25 January, there is no foF2 and h p F 2 data due to the presence of a spread-F-layer.Otherwise, nothing is unusual to report from Fig. 9, except an increase in the 16:00-18:00 UT time interval (21:30-23:30 LT) which stops 3 days before the  quake.This anomalous increase has an important effect on the monthly average value, as can be seen on Fig. 7. Particular attention has been given to the nighttime and we have checked the variation of h p F 2 during the January and February months of several years.The results are presented in Fig. 10.The top panel shows the averaged values of h p F 2 which take into account the measurements at 22:00, 23:00 and 24:00 UT (03:30, 04:30, 05:30 LT).This parameter is plotted as a function of the days for January and February 1998.The thick dashed line represents the averaged values of this parameter over the two months, whereas the two thin dashed lines are related to the standard deviation σ .The second panel represents the variation of the magnetic activity index A p as a function of the same days.The two middle panels concerns the same parameters for the year 1999 whereas the two bottom panels are related to 2001, the year of the event (no ionospheric data are available for the year 2000).It is interesting to compare the variations of h p F 2 with A p for these three years.In 1998, during the two months January and February, the A p index does not exceed 40 and h p F 2 remains between +σ and −σ .In 1999, the A p index is equal to 40 and more on 13 January and 19 February and this induces at the same time a large increase in h p F 2 .In 2001, A p is much less than 40 during the 2 months and then a large variation of h p F 2 is not expected.However, three days before the quake we observe an increase larger than 2σ .Another increase on 31 January could be related to aftershocks.

Discussion and conclusions
Table 2 sums up the ionospheric observations done before the Bhuj earthquake.The following anomalous ionospheric variations have been observed: -The TEC from TOPEX indicates a decrease on 21 January around 10:00 UT (daytime).-The TEC from SPOT2 and SPOT4 indicates an increase on 20 January (daytime) and a perturbation on 24 January around 06:30 UT (this is confirmed by the sounder because the data corresponding to this time are not available).
-The TEC from the Bangalore station shows similar afternoon -evening perturbations (decrease) on 21 and 25 January.
-The ionospheric sounder close to the epicenter indicates an increase of the foF2 variability during the night of 21 January.An anomalous behaviour of foF2 at 10:00 UT is shown several days before the quake.An increase in foF2 is observed several days before the quake during evening time (21:30-23:30 LT).An increase of h p F 2 is observed 3 days before the quake during nighttime (03:30, 04:30, 05:30 LT).On 25 January around 20:00 UT, the ground-based sounder detects a spread-F-layer at this local nighttime.
From the observations, we can conclude that the ionosphere behavior changes a few days before the seismic shock with increased variability starting on 19 January, which can be interpreted as anomalous pre-seismic ionospheric variations.
It is common knowledge that the ionosphere is mainly under the control of the solar activity and that it can be affected by TIDs (Traveling Ionospheric Disturbances).However, we have shown anomalous behaviours prior to the Bhuj earthquake relative to a large time interval around the event time.All these anomalous perturbations depend on the time and the latitude of the observations.This is in agreement with previous events reported by Pulinets andLegen'ka (2002, 2003) and Liu et al. (2001Liu et al. ( , 2002)).
A hypothesis to explain these ionospheric perturbations is related to the action of the electric field.Due to the stress of the rocks, electric charges could appear at the Earth's surface and change currents in the atmosphere -ionosphere system.The effects of these electric fields on the ionosphere were modelled in Pulinets et al. (2000).It is noteworthy that, in relation to the possible mechanisms described in Sect. 1, Virk et al. (2003) have observed a radon anomaly on 22 January 2001.

Fig. 1 .
Fig. 1.Orbit tracks of the satellite TOPEX-POSEIDON around the position of the Bhuj earthquake (see also Table1).

Fig. 2 .
Fig. 2. Positions of the SPOT2 subionospheric points around three different DORIS stations during 24 January passes.

Fig. 3 .
Fig. 3. Daytime variations observed by TOPEX-POSEIDON, (a)values of the ionospheric corrections as function of the latitude (see also Table1), (b) differences between the ionospheric corrections and the Bent model as a function of the latitude.

Fig. 4 .
Fig. 4. Top panel: TEC variations in time as a function of the latitude for the two satellites SPOT2 and SPOT4 during two months around the time of the Bhuj earthquake.It concerns the data from the KITB station (see Fig. 2).Middle panel: Corresponding TEC values in TEC units (10 16 m −2 ).Bottom panel: Magnitude and occurrence of the main aftershocks.In the two first panels the horizontal full line represents the mean value, whereas the dashed lines indicate the variance, and the vertical line indicates the time of the Bhuj earthquake.

Fig. 5 .Fig. 6 .
Fig. 5. Variations of the TEC data recorded by the Bangalore station during the second half of January 2001.

Fig. 7 .
Fig. 7. foF2 variations on 23 January 2001 in a continuous line with points.In the green line, it is the monthly median, in the red dashed line it is the upper bound f oF 2+2σ , and in the blue dashed line it is the lower bound f oF 2−2σ .

Fig. 9 .
Fig. 9. Hourly variations of foF2 in UT at the Ahemadabad ionospheric station during January and February 2001.The data are color-coded according to the scale on the right.

Fig. 10 .
Fig.10.Variations of an average h p F2 value taken between 22:00, 23:00, and 24:00 UT at the Ahemadabad ionospheric station and variations of the magnetic index A p as a function of the January and February days.From the top to the bottom, the plots are done for the year 1998, 1999, and 2001 (the quake occurred on 26 January 2001).

Table 1 .
Days and hours corresponding to the TOPEX-POSEIDON passes along the tracks A, B, C, and D of Fig. 1.
slant TEC values between satellites and ground-based stations are given.The available parameter is an ionospheric correction in m/s given by

Table 2 .
Main ionospheric variations observed by various experiments before the Bhuj event.Ahemadabad f 0 F 2 variations f 0 F 2 variations f 0 F 2 variations f 0 F 2 variations f 0 F 2 variations