Study of the TEC data obtained from the DORIS stations in relation to seismic activity

Ionospheric data obtained from the DORIS system are used in this paper. The DORIS system is composed of several ground-based beacons which emit at two frequencies (400 MHz and 2 GHz) and of receivers on board several satellites (currently SPOT2, SPOT4, SPOT5, Topex-Poseidon, Jason1 and Envisat). Thanks to the density of its network coverage (∼50 stations), DORIS provides information on the ionosphere. The TEC (Total Electron Content) parameter which is the electron density integrated over the vertical could be obtained from DORIS measurements. In a first step, the paper describes the way to obtain the TEC data from the DORIS ionospheric measurements, and comparisons of the results are done with the IRI2001 model. In a second step, TEC values are used to search for correlation between ionospheric perturbations and seismic activity. Earthquakes of magnitude larger than 5 are chosen close to the ground-based DORIS stations. Among other results, the statistics show that, during the night time and at geomagnetic latitude close to the equator (<10°), TEC amplitude fluctuates at the time of the earthquakes as it is expected, but also 2 days and 5 days before. Mailing address: Dr. Michel Parrot, Laboratoire de Physique et Chimie de l’Environnement (LPCE), CNRS, Université d’Orléans, 3A Avenue de la Recherche Scientifique, 45071 Orléans Cedex 2, France; e-mail: mparrot@cnrs-orleans.fr


Introduction
The aim of this paper is to perform a statistical analysis of the variations of an ionospheric parameter called TEC (Total Electron Content) in relation to seismic activity in order to reveal a seismo-electromagnetic effect.This statistical approach was chosen due to the wide variability of the ionosphere which is mainly under the influence of solar activity.
The seismo-electromagnetic effect is the electric and magnetic perturbation caused by the natural geophysical activities (earthquake, volcanic eruption).Its existence has been underlined since the 50s and this phenomenon is very interesting because for certain cases, it occurs before the earthquakes and can be considered a precursor (see for example monographs by Hayakawa, 1999;and Hayakawa and Molchanov, 2002).
During these 50 years, scientists have formed numerous laws which are based on the mathematic descriptions of the earthquakes.Fenoglio et al. (1994) proposed a model of ruptured isolated reservoirs resulting in an electrokinetic generation of transient magnetic field.Another mechanism of ULF (Ultra Low Frequency) electromagnetic field generation based on the crack opening (micro fracturing) was studied by Molchanov and Hayakawa (1995): the creation and relaxation of electric charges at the walls of opening cracks in the earthquake hypocenter was supposed to be a possible reason of electromagnetic influence prior to the earthquakes.Surkov (1999) suggested a theory which is based on the same assumption of crack opening, but the electromagnetic perturbation is generated by MHD (Magneto Hydro Dynamics) effect from the propagating seismic wave.ULF electromagnetic emissions have also been observed prior to the occurrence of earthquakes in U.S.A., Russia and Japan (see the recent review by Hayakawa and Hattori, 2004).Molchanov et al. (2002) explain these phenomena by the following two steps: first, generation of a seismo-electromagnetic field; second, penetration of the electromagnetic fields through a dissipative ground medium to reach the ionosphere.
The above mechanisms are connected with different radiation sources existing in the lithosphere.Sorokin and Fedorovich (1982) reported that these waves propagate within thin layer of lower ionosphere along the Earth surface at low and middle latitudes with a small attenuation.Further studies have shown that there is an enhancement of DC (Direct Current) electric field which accompanies the generation of these periodic inhomogeneities (Sorokin and Yaschenko, 1999).Sorokin et al. (2002) supposed a new mechanism of seismogenic geomagnetic pulsations which is connected with the radiation source in the ionosphere; it is based on the excitation of gyrotropic waves in the presence of periodic horizontal inhomogeneities of electric conductivity in the low ionosphere.
At the same time, Tronin (2002) presented another view of the litho-atmosphere coupling in seismic processes because NOAA/AVHRR satellite thermal images indicated the presence of positive thermal anomalies that are associated with the large linear structures and fault systems of the Earth's crust.Phenomena like gas release, moisture change, and temperature increase of the soil, could be the source of thermal anomalies related with seismic activity.
In his review on models, Pulinets (2002) studied all the following phenomena: the DC vertical atmospheric electric field, ULF emissions, VLF (Very Low Frequency) emissions, etc., in relation to atmospheric plasma which can be produced by the radon and various metal aerosols emanating from the Earth's crust in seismoactive zones.The ionizing radiation and hydration process generates earlier a nearground electric field which could lead to effects of coronal discharge and triggers the electromagnetic emission.The atmospheric plasma described here contains metallic aerosols, and then, is a typical dusty plasma with characteristic plasma instabilities which can also give rise to electromagnetic emissions at different frequencies.Pulinets concludes that seismo-ionospheric effects observed at all heights of the atmosphere-thermosphere-ionosphere-magnetosphere chain is a kind of electromagnetic coupling due to large-scale anomalous electric field appearing over the seismoactive area several days before strong earthquakes.He considers this process like plasma-chemistry processes and ion-molecular reactions over the ground surface in seismo-active zones, regardless of the physical mechanism of the large scale electric field, it produces no effects within the ionosphere.
Concerning the ionospheric perturbations, Pulinets (1998) reported several seismo-ionospheric variations.These irregularities occurring in the ionosphere could be related to electrodynamic and plasma-chemical processes.The seismic activity can induce positive or negative variations in electron density, variations of h pF2, electron temperature, ion and neutral composition, and formation of sporadic E and spread-F layers.The following characteristics have been underlined: -While magnetic storms last 8 to 48 h, seismo-ionospheric variations take only 3 to 4 h, they appear every day during a period up to 5 days prior to earthquakes at the same local time.
-The seismo-ionospheric variations of the electron density could be negative or positive; at a given height the sign depends on the local time.
-The seismo-ionospheric variations result in a re-distribution of electron concentration with height.
-The modified region of the ionosphere is «tied» to the position of the anticipated earthquake though the greatest changes do not necessarily coincide with the epicenter projection.
-Ion and neutral composition changes may accompany the electron concentration variations connected with the seismic activity.
An alternative mechanism which can produce such perturbations is related to the atmospheric gravity waves which can be triggered by gas releases or thermal anomalies (Molchanov et al., 2001;Shvets et al., 2004).
At much higher altitudes above the F-layers, Afonin et al. (1999) found a correlation with the general distribution of seismic activity and ion density variations measured by satellites.
The paper is organized as follows: Section 2 will briefly summarise some recent papers related to the TEC and ionospheric parameters variations with the seismic activity.The data which will be used here will be described in Section 3. In particular the method to obtain the TEC from the DORIS system will be developed and the results will be compared with an ionospheric model.Section 4 presents the results of the statistical study about the variation of the TEC with the seismic activity and Section 5 will give conclusions.

TEC variation with the seismic activity
TEC (Total Electron Content) is the integrated electron density which is measured between a satellite and a ground beacon or between two satellites.It is an important parameter which can be used to survey the ionosphere.
In recent years, many scientists have carried out investigations on perturbation of the Earth's ionosphere by earthquakes using GPS data.For example, Calais and Minster (1998) experimentally measured an ionospheric effect after the M=6.7 Northridge earthquake of 17 January 1994.Ionospheric oscillations have been also registered by Ducic et al. (2003) after the Denali Park Alaska Earthquake of 03 November 2002 (M=7.9).Concerning the pre-seismic effects, Liu et al. (2000Liu et al. ( , 2002) ) examined the TEC variations during three earthquakes: Rei-Li on 17 July 1998, Chi-Chi on 20 September 1999, and Chia-Yi on 22 October 1999.It can be seen that significant VTEC (Vertical) decreases three days before the Rei-Li earthquake: 4 days and 3 days before the Chi-Chi earthquake; 3 days and 1 day before the Chia-Yi earthquake.Liu et al. (2000) analyzed the variation of foF 2 recorded by the Chung-Li ionosonde and found a decrease of electron density 1, 3, and 4 days before the Chi-Chi earthquake, which is correlated with the TEC perturbation.They suggested that these perturbations caused by the vertical electric field are generated over the seismoactive zone.
The global variation of the ionosphere can also be studied by ionospheric sounders and recently Ondoh and Hayakawa (2002) showed anomalous sporadic ionizations in the E layer before two Japanese earthquakes with magnitude larger than 7 using a network of ionospheric sounders.They presented a model where these perturbations are attributed to seismo electrostatic discharges between the lower ionosphere, a tropospheric cloud (formed by terrestrial ions carried with upward air currents), and the ground.Liu et al. (2004) performed a statistical analysis with 20 earthquakes of magnitude ≥6.0 in the Taiwan area from September 1999 to December 2000.It shows that pre-earthquake ionospheric anomalies occur in the local time interval 18:00-22:00 during 5 days before earthquakes.

The data
3.1.The DORIS system DORIS (Doppler Orbitography and Radio positioning Integrated by Satellite) is a radioelectrical system for high accuracy orbit determination and station positioning.
The DORIS system is composed of the following elements: -A world-wide permanent network of transmitting stations, called «orbitography network», it is part of the ground segment of the DORIS system and is composed of 56 stations evenly distributed on the Earth's surface (see fig. 1).
-So-called «ground localization stations», whose position is unknown a priori.
-A control center that performs system monitoring, instrument programming, and data processing and archiving.

Calculation of the TEC parameter
The ionospheric data files are downloaded from the ftp server dedicated to the DORIS ionospheric products ftp://cddisa.gsfc.nasa.gov/pub/doris/products/iono/ssa where we can find out the data for several satellites during given time intervals shown in table I.
The files are organized by the satellite passes over the different beacons and contain the following information: -for the header: satellite, beacon, number of observations, max of satellite pass elevation, local time, meteorological data: pressure (mb), temperature (deg), humidity (%); -for the data: CNES julian date, second in the day of the observation (TAI), elimination criteria, count interval (2 GHz channel), count interval (400 MHz channel), tropospheric correction (2 GHz channel), tropospheric correction (400 MHz channel), ionospheric correction where i is the number of measurements, mes is the ionospheric correction in m/s (0.02 cycle → →0.3 mm/s at 2 GHz); c i is a constant related to the frequency of emission and the count interval; TEC i is the vertical TEC value at the subionospheric point for each measurement during a cross coupling; and ki is a geometric factor to convert a slant TEC to a vertical TEC for each measurement during a cross coupling r is the radius of the position of the subionospheric point relative with the Earth's centre; RE is the radius of the Earth; Ei is the elevation angle of the satellite relative with the station.It is assumed that the TEC values can be approximated by a polynomial of degree N function of the latitude of the subionospheric point.Several values of N have been tested and N=5 appeared to be enough to reproduce the TEC variations.Then the TEC parameter can be expressed as (2 GHz channel), ionospheric correction (400 MHz channel), elevation angle (°), azimuth (°), station-satellite distance (m), acquisition mode, received power level (400 MHz), received power level (2 GHz), ponderation (0 if eliminated or 1), Doppler count (2 GHz).Among these data, the following parameters related to each pass close to a station are used: date, time, count interval (2 GHz), ionospheric correction (2 GHz), elevation angle, azimuth, and station-satellite distance.
For each flying over each station, we calculate at each time the coordinates (latitude and longitude) of the satellite track and then the coordinates of the subionospheric point at which the TEC will be estimated with the parameters: elevation angle, azimuth, and station-satellite distance (see fig. 2).
The value of the maximum of ionization Hm which is the peak of electron density content along altitude is chosen to be 300 km, and the value of the offset is fixed in the processing parameters to 50 km.
One example of the results is given in fig.
3 which shows the positions of the station AMSB, the subionospheric points and the trace of the satellite SPOT2 during a visibility on January 7, 2001.
Then, TEC which represents the electrons content of a 1 m 2 vertical column between the    ) ) ) ) ) This over-determined system was solved with a least square method in order to get the coefficients a i for each satellite passes over a beacon, and then, the TEC values are obtained with eq.
(3.3).Some results are presented and compared with the IRI2001 model in fig.4a-f.The International Reference Ionosphere (IRI) is an international project which produces an empirical standard model of the ionosphere, based on all available data sources.For given location, time and date, IRI describes the electron density, electron temperature, ion temperature, and ion composition in the altitude range from about 50 km to about 2000 km, and also the electron content.The 6 panels (a) to (f) are related to different satellites and different stations.The data is shown as function of latitude for a complete flying over a beacon.It shows that the TEC estimation is in good agreement with the IRI2001 model for the different stations and satellites.However this polynomial estimation gives smoothed TEC values and the detailed variations which are in the original DORIS measurements have been lost.

The earthquakes
The data related to the earthquakes were downloaded the USGS web server: <http:// neic.usgs.gov/neis/epic/epic_global.html>.The following parameters were used: date, hour, latitude, longitude, magnitude.Events were selected from January, 1, 2001 until April 30, 2004 according to the conditions below: -The magnitude M s of the earthquakes must be larger than 5.0.This minimum value was chosen in conformity with Pulinets et al. (2004) who show that M∼5 is the magnitude threshold to detect ionospheric earthquake precursors.
-The distance between the seismic centre and the stations must be less than 400 km.This is in agreement with the empiric formula given by Dobrovolsky et al. (1979) concerning the radius of the earthquake preparation zone.
Following these constraints, only 17 DORIS stations were considered in the statistical work.Their locations are given in fig.5, and table II gives the number of earthquakes for each of them.
In total there are 303 earthquakes which meet these conditions for the 5 satellites: Jason1, Spot2, Spot4, Spot5, Topex1 and for the 17 stations.In order to observe the TEC perturbations above these stations, we take the data during an interval of time which is 8 days before the earthquakes and 3 days after the earthquakes.It must be noted that for a given earthquake and then for a given station, several satellites may pass over this station (but not at the same time).

The statistical results
To process the statistics, we use the method of superposed-epoch and the data base is organized around three parameters: -The time difference between the seismic time and the time of TEC measurements.The seismic time is chosen to 0, and then this time difference is between −8 days and +3 days.A resolution of 8 h is taken.
-The geomagnetic latitude of which values are between 0 and 60°with a resolution of 10°.
These periods were chosen to follow the main steps of the average TEC diurnal fluctuations.The data base not only contains the TEC estimated from the DORIS measurements but also the TEC from the IRI2001 model.It must be noted that this model takes into account the global magnetic activity.
Then the system is organized in 33×6 ×4 cells in which the TEC measurements are averaged and the plots in fig.6a-d display some results.Figure 6a presents the averaged measured TEC as functions of the time and of the geomagnetic latitude.The data is colour-coded according to the scale on the right.This plot corresponds to night time (20:00-6:00).Figure 6b displays with the scale of fig.6a, similar data but for the IRI2001 TEC.  Figure 8c presents the corresponding variances of the average values plotted in fig.8a.The average differences between the measured TEC and the IRI2001 TEC are plotted in fig.8d.

Conclusions
This work has been done in relation to the main scientific objectives of the DEMETER micro-satellite project which is to study the ionospheric perturbations in relation with the seismic activity (Parrot, 2002).TEC were initially calculated from the DORIS ionospheric files, and a very good agreement was found between the estimated TEC and the TEC of the IRI2001 model.A statistical study was done with these estimated TEC data to check if the measurements can be influenced by the seismic activity.The statistical results show that there are perturbations of TEC several days before the earthquakes and right after the time of earthquakes.These two points are well in agreement with previous works.Ducic et al. (2003) presented observation of TEC variations just after the time of a powerful earthquake.In this work, TEC variations are found 2, 2.5, and 5 days before earthquakes, at the end of afternoon and at night-time, for invariant latitude between 0°a nd 10°.It corresponds with the results of Liu et al. (2004) who observed variations between 18:00 and 22:00 LT within 5 days before earthquakes in Taiwan area.With another set of data, Pulinets (1998) claims that ionospheric perturbations can appear up to 5 days prior to earthquakes.

Fig. 1 .
Fig. 1.The stations of the DORIS system (from the IGN web server).

Fig.
Fig. 4a-f.Examples of comparisons between the TEC estimation and the TEC determined from the IRI2001 model for different satellites and different stations.Date and local time are indicated.
of the M DORIS measurements; is the vector of unknown polynomial coefficients; and B is a matrix of size N * 6 of which elements are

Fig. 5 .
Fig. 5. Map of the DORIS stations used in the statistics.
Figure6apresents the averaged measured TEC as functions of the time and of the geomagnetic latitude.The data is colour-coded according to the scale on the right.This plot corresponds to night time (20:00-6:00).Figure6bdisplays with the scale of fig.6a, similar data but for the IRI2001 TEC. Figure 6c presents the corresponding variances of the average values plotted in fig.6a.The average differences between the measured TEC and the IRI2001 TEC are plotted in fig.6d.The main features of the panels (a), (c) and (d) of fig.6ad are observed for geomagnetic latitudes between 0 and 10°where three anomalies are seen.The first one which occurs just after the 00-20:00).Only one anomaly is observed about 2.5 days before earthquakes for geomagnetic latitudes less than 10°.At geomagnetic latitudes larger than 50°, variations in fig.7a,cand fig.7dare due to the small numbers of events at these high invariant latitudes.Figure 8a-d contains the most important data shown in figs.6a-d and 7a-d but presented in a different way.All panels are plotted for an invariant latitude less than 10°and, for example, fig.8a presents the averaged measured TEC as functions of the time (before and after the quakes) and of the local time of the TEC.General feature of the usual diurnal TEC variation is observed on fig.8a and 8b which is similar to fig.8a but for TEC obtained from the IRI2001 model.The TEC values are maximum in the noon sector (10-15).Anomalies seen in figs.6ad and 7a-d are again detected on fig.8a in the end-of-afternoon sector (15-20) and in the night sector (20-6).

Fig
Fig.7a-d.Results of the statistics for end-of-afternoon.

Fig
Fig. 8a-d.Results of the statistics for invariant latitude 0°and 10°.

Table II .
Numbers of earthquakes for each selected DORIS station.