Detection of ionospheric signatures from GPS-derived total electron content maps

Abstract The processing of measurement data from satellite constellations such as Global Navigation Satellite Systems (GNSS), including the well-known Global Positioning System (GPS), have been successfully applied to virtually all areas of geophysical sciences. In this work, a method is described where Geographical Information Systems (GIS) are employed to build hourly ionospheric Total Electron Content (TEC) maps for 2011 over the southern Iberian Peninsula. The maps used GPS-derived geometryfree linear combinations attained from station data from the Algarve, Alentejo (Portugal), Andalusia, Murcia and Valencia (Spain) regions. Following the construction of the ionospheric maps, it was possible to relate these results to natural phenomena. The observed phenomena included diurnal and seasonal variations: daytime TEC maxima, nighttime TEC peaks, summer TEC value decreases, and spring and fall TEC maxima. After validation of these periodic phenomena, detection of non-periodic changes, such as solar flares and tectonic interactions with the ionosphere were attempted. The results showed a TEC increase following a selected solar flare event and a potential TEC build-up prior to the 2011 Lorca earthquake. Further studies could open up the possibility of building early warning systems. The presented methods, based on available software packages, are also of value in monitoring the effect of the ionosphere on radio signals, satellite and mobile communication, power grids, and for accurate GNSS navigation.


Introduction
The ionosphere is the subject of extensive scienti c research mainly due to the fact that it a ects the propagation of electromagnetic (EM) signals. Thus its presence must be taken into account whenever a satellite is sending radio signals through the atmosphere, and it is important to consider for radio communication and navigation systems in general. Today it is not only possible to correct for ionospheric refraction, but also to build ionospheric maps by computing the delay (or advance in the case of carrier phase) it causes in the GNSS signal, allowing GNSSderived data to be used for model building.
The ionosphere can be mapped regionally or globally, depending on the study area and the available ground GNSS stations. Prior literature such as in Camargo et al. (2000), Georgiadou (1994), Leick (1995), Orús et al. (2003), Ping et al. (2002 and, and Taylor et al. (2006) provides a comprehensive theoretical background for GPSderived ionospheric TEC modeling. Jin et al. (2011) argue that GPS satellites are capable of providing even more accurate detection of ionospheric parameters than traditional ionospheric detection methods such as ionosondes, scatter radars, topside sounders, onboard satellites, and in situ rockets, and at a cheaper user cost. Additionally, they developed a piece of software called Regional Ionospheric Mapping and Tomography (RIMT), which can monitor 2D TEC and map 3D ionospheric electron density distribution using GPS measurements. Having up to 1Hz GNSS data availability, high temporal resolution mapping can be carried out. Orús et al. (2003) and Ping et al. (2003), used data from regional and global GPS networks which were utilized to build Global Ionospheric Maps (GIMs) and Regional Ionospheric Maps (RIMs) over Europe and Japan.
Studies concerning speci c ionospheric behavior above the study area of southern Iberia have not yet been conducted. However, the technological background exists for such research because Spain and Portugal support networks of permanent GNSS stations.
One potential application of ionospheric mapping is the detection of earthquake signatures and other natural phenomena. The most seismologically active areas of the Iberian Peninsula are located in southern Iberia (Gibbons et al. 2003), more speci cally along the Baetic System (Andalusia and Murcia, Spain), and between the Goringe seamount and Cape Saint Vincent (Algarve, Portugal). These regions are also the most hazardous, especially along their coastline (due to possible tsunami event). During the year 2011, the time span of the data applied in the present work, there were two noteworthy seismic events inside the study area. The more energetic of these two occurred in southern Spain, later named the 2011 Lorca earthquake. This was a moderate magnitude 5.1 M W shallow-focus 1 km event on May 11, 2011 at 06:47:25 local time. Its coordinates were 37.699°N 1.673°W, 50 km southwest of Murcia, situated near a major fault, the Alhama de Murcia fault. Ouzounov et al. (2011) retrospectively analyzed spatial and temporal variations of physical parameters, such as GPS derived TEC, of the M W = 9 Tohoku Japanese earthquake of March 11, 2011. They could characterize the state of the ionosphere several days before the onset of the seismic event. The GPS-derived TEC values indicated an increase in electron density which had reached its maximum 3 days before the earthquake. They had found a positive correlation between the ionospheric anomalies and the Tohoku earthquake. Perrone et al. (2009) conducted research regarding Italian earthquakes, including the April 6, 2009 earthquake in L'Aquila. They investigated whether evidence for an ionospheric precursor can be found. Although they only considered 6.0> M W > 5.5 events, which are signi cantly less energetic than the Tohoku earthquake, empirical dependencies for the seismo-ionospheric disturbances relating to the earthquake magnitude and the epicenter distance were obtained and they were shown to be similar to the ones obtained for the Tohoku earthquake. They concluded that the similarity of these processes in di erent parts of the world suggests a uniformity of the processes during the earthquake preparation period both for powerful and moderate earthquakes. Lognonné et al. (2009) argue that tsunami waves propagating across long distances in the open ocean can induce atmospheric gravity waves by dynamic coupling at the surface, and thus can be detected from TEC measurements. In that work, ionospheric TEC monitoring of Eu-rope, California, and Japan was used for several investigations, in addition to the study of post-seismic signals. It is important to keep in mind that all of these studies were performed using post-processing data. Earthquake forecast, as a scienti c method, is still far from being accepted by the scienti c community and there are signicant doubts that it can ever be done reliably. Nevertheless, TEC monitoring, as it was demonstrated by the aforementioned articles, shows promising and important results.
As one nal example of an application of ionospheric mapping, Jin et al. (2011) describe a case study of GPS ionospheric mapping and tomography in a geomagnetic storm caused by solar activity. They conclude that a strong increase can be shown in the ionospheric observables during the storm, although they suggest further investigations. Thus great potential also exists for detecting solar e ects from GNSS-derived ionospheric TEC models. The present work describes a series of techniques that can be employed to build ionospheric maps, TEC time-series, and also describes a number of attempts to link these maps with other environmental variables.
There are two primary goals in this work. Firstly, this study aims to develop a GIS application which is able to build ionospheric maps from GNSS-derived TEC values, and to demonstrate how these vary in time above the selected study area. Daily, monthly, and seasonal changes of 2011 were studied. It is crucial to understand what natural phenomena can a ect GNSS accuracy, how they affect it, and what can be done to compensate for these errors. These can be partially achieved by building accurate ionospheric models. It is worth noting that solar ares that disturb the ionosphere can cause position errors in GPSbased navigation systems, such as the Wide Area Augmentation System (WAAS), as large as 50 meters if no corrections are applied. This can present an extreme hazard for aerial and other navigation systems (Hofmann-Wellenhof, 1994). Secondly, the results of this work can assist atmospheric and geophysical scientists to better understand the regional behavior of the ionosphere. A number of global and regional ionospheric mapping projects already exist, as described above, but this is the rst study that concentrates speci cally on the southern Iberian Peninsula.

. Study Area and GNSS resources
The study area is located in the southern portion of the Iberian Peninsula (Fig. 1) Figure 2 shows the distribution of all stations involved in this work. The study region overlaid with the Voronoi decomposition of the study area using the selected stations as sites. Note that on this map only the Priority 1 (P1) stations are considered. P1 stations had good enough quality data to be included in the TEC map building, while P2 stations had to be rejected.

. GNSS ground station selection
After selection of GNSS ground stations in the study area, RINEX data les were collected with 30 second sampling rates. The four regions of the study area contain a relatively high density of permanent GPS stations, but they belong to four di erent networks. Typically they are less than 100 km from each other; in some areas this is as short as 30-60 km. Due to the large amount of data storage, computing power, and computing time that TEC map computing requires, selection of a representative sub-network was required. This was done by selecting no more than 2 stations closer than approximately 60 km from each other, and with no geographical point further than 90 km from a station. The network should be as uniform as possible and dense enough that the interpolation between the intersection points of the lines receiver-satellites with the ionospheric layer, the called ionospheric pierce points (IPP), provide acceptable values, and also taking into account that the line-of-sight TEC data quality is proportional to 1/sin(ε), where ε is the elevation angle above horizon, thus the highest quality probing is accomplished near the stations zenith. There should also be stations near the boundary of the study area. The scheme that was used is shown in Fig. 1, with station selections in Table 1. This consists of a total of 25 stations. Although not completely ideal (up to a 118 km exists between two stations, while two other stations are 62 km from each other), this was the most optimal network that could be created from the available stations. The GIS analysis of the GNSS ground station network was performed in ArcGIS (version 9.3; ESRI Inc., Redlands, CA, http://www.esri.com). Precise station antenna coordinates are typically provided by the GNSS network data processing -these coordinates can be entered one by one into a point geometry shape le, resulting in the map shown in Fig. 1. An additional coastline contour (polygon geometry) shape le was also created. Based on this station shape le and the coastline contour shape le various spatial analyses can be conducted. Note the Priority column in Table 1. Priority 2 stations had to be eliminated from the nal subnetwork mainly due to occasional de cient and/or sporadic RINEX data provided by these stations.

. Determination of a representative GNSS ground network based on GIS
One of the basic assumptions behind our GNSS network selection was that an ideal GNSS monitoring network Table 1. The complete stations list. Including 4 character station names and other relevant information. "P" stands for Portuguese and "S" for Spanish stations. The data continuity of the Priority 2 stations was sporadic for the selected days; therefore they were eliminated from the nal processing.

St.# -ID Site Name Region Priority
Latitude Longitude should consist of uniformly distributed ground stations with stations at least 60 km from each other. This 60 km distance is double the maximal baseline (30 km) between a single station and a mobile receiver inside which the ionospheric errors are considered still acceptable (during normal ionospheric circumstances), so the relative TEC error will be in the range of 1 to 2% (Astafyeva et al. 2008). However, twice that distance was allowed, considering the probing locations are the IPPs for each receiver-satellite line-of-sight. For higher space and time resolution, it is necessary to process the data from more GNSS stations; however this is not always possible due to limitations in data processing and limitations in the availability of stations.
In order to mathematically describe the requirements for an optimal network, geostatistical methods can be applied to determine the best possible combination of stations. At each epoch, the maximum number of measurements equals the number of GNSS stations multiplied by the number of visible satellites from all stations. The ionospheric TEC values outside the measurement points can be estimated by interpolation. The GNSS stations situated in the study area are permanent stations operated by various institutes and companies independent of this work. Therefore the approach to mathematically describe, evaluate, and identify the optimal sampling scheme for this irregular structure was to rst decompose the study area through Voronoi tessellation ( Fig. 1) (Atsuyuki et al. 2000). In this analysis, there are n stations (Table 1) with appropriate geographic coordinates, that is s , s , . . . , sn, and each station/site s i has its corresponding Voronoi cell C i . C i contains all the points in the study area which are closer to s i than to any other station, s j , and where i ≠ j (Boots, 1986). Paláncz et al. (2006) have stated that the area of a Voronoi polygon generated by a site may be considered as the Region of Attraction (RA) of this site, simply because these points are closer to the measurements of this station than to any other. Therefore, it is expected that more accurate TEC mapping will be achieved for smaller polygon areas and less accurate mapping for larger ones.
It is now possible to rephrase the rst assumption of this section by stating that if stations are uniformly  (Dubois, 2000). For instance, isolated stations will have larger corresponding polygon surfaces. The distribution of the GNSS stations in the study area is non-uniform and denser than it was described in the beginning of this section. It is desired that the RA of the optimal sub-network have a homogeneous distribution and relatively low standard deviation. Dubois (2000) has recommended an additional measure that takes into account the distance of each station to its nearest neighbor and the surface of the Voronoi polygon. Let S V be the surface of the Voronoi polygon belonging to s i station and Sm the mean polygon surface de ned as: Where S total is the area of the whole study region and n is the number of Voronoi polygons in it. N i d is de ned as the distance between a site and its nearest neighbor. Using these variables, a quantity (for each station/polygon) called Coe cient of Representativity (CR) can be de ned: In an ideal case CR = 1. CR > 1 means that the station is isolated (i.e. it is relatively far from the other stations), while CR < 1 means that the station is close to a neighboring station (Paláncz et al. 2006). It can be concluded that the statistics of the polygon surface areas and the CR values provide tools to determine the regularity of a given sampling scheme. The CR should help identifying undersampled areas (Dubois, 2000). ArcGIS software was utilized to aid this geostatistical analysis. During the station distribution analysis the polygon areas had to be determined and the ArcToolbox/Analysis Tools set was employed, resulting in the RA. Additionally, the CR was calculated based on the nearest neighbor station. A MATLAB (version 2010a; MathWorks, http://www.mathworks.com) script was written to compute and plot the statistics of the values acquired through the geospatial analyses.
. Obtaining TEC values from Bernese GPS Software 5.0 Bernese GPS Software (version 5.0; Astronomical Institute, University of Bern, http://www.bernese.unibe.ch/) can build three types of ionospheric models based on spherical harmonic expansions. These are global, regional, and station speci c. For the purpose of this work the regional model was selected. For the processing strategy, precise point positioning (PPP) was selected, which is a special case of zero-di erence processing. The main purpose of performing PPP in Bernese using GPS observations from a number of stationary ground receivers is to obtain a set of station coordinates at a cm resolution level. From these PPP results the extraction of TEC information is possible. As opposed to other processing strategies, the satellite clock corrections are not estimated but are assumed to be known. They are introduced in the processing together with orbit information and Earth orientation parameters. These parameters can be obtained from the International GNSS Service (IGS) database. The parameters left to estimate are station clock corrections, coordinates, and troposphere parameters. The parameter estimation is done sequentially by executing a user-de ned list of Bernese scripts. These scripts employ data reading, preprocessing, conversion, synchronizing, and parameter estimation. Ultimately, the PPP solutions are computed, and optionally, ionospheric models can be generated (Dach et al. 2007).
The quality of the introduced information should be as good as possible for PPP, as all errors directly propagate to the position of the station. The regional TEC model may be written as: where β is the geographic latitude of the intersection point of the line receiver-satellite with the ionospheric layer (this point is called ionospheric pierce point (IPP)); s is the sun-xed longitude of the IPP; nmax is the maximum degree of the spherical harmonic expansion, Pnm are the normalized associated Legendre functions of degree n and order m based on a normalization function, and anm and bnm are the unknown TEC coe cients of the spherical harmonics, i.e. the regional ionosphere model parameters to be estimated (Dach et al. 2007). We used the phase observables for the geometry-free linear combination (L4) to estimate the ionosphere. A 10 degree cuto angle and 1/sin(ε) elevation-dependent weighting were selected. The type of the temporal modeling was static employing the geomagnetic reference frame. After PPP processing, the IONosphere map EXchange format (IONEX) les can be obtained (Schaer et al. 1998). The TEC map section contains coordinates and values in TEC units (TECU); the rootmean-square error (RMS) map section contains the RMS values for each location. Dach et al. (2007) recommends degree and order values around n=6 and m=6 for regional modeling. After running test processes with several degree and order values, n=5 and m=5 were selected. This was due to the facts that with the selected GNSS network, the higher values did not reveal additional ionospheric structures, but the processing time was increased signi cantly. The data les were processed from the selected stations for a 450 km single layer height. Latitudes were between 35°a nd 40°, with 0.5°steps, and longitudes were between 12°W and 2°E, with 0.5°steps.

. Interpolation and Visualization
Interpolation algorithms estimate the value at a given location as a weighted sum of data values at surrounding locations, with weights assigned according to functions that give a decreasing weight with increasing separation distance (Cressie 1990). In order to process the interpolation in ArcGIS and Surfer (version 9; Golden Software, https: //www.ssg-surfer.com/), the TEC values and their corresponding time parameters and coordinates were read from the IONEX les and transformed into an XYZ array using a custom MATLAB script, which was then imported into ArcGIS or Surfer. Once all the epochs were imported into ArcGIS or Surfer, the interpolations were made using inverse distance weighted (IDW), kriging, natural neighbor, spline, and trend methods (Davis, 1975). All methods were tested, and the kriging method was chosen (see Results). ArcGIS was used for calculating interpolation statistics, and Surfer was used for visualization. Composite images were rst created in MATLAB using a reference image consisting of the rst epoch or day of the series, and images from the remaining epochs or days. The reference image values (layer or data matrix) were subtracted from the referent image values using pixel-by-pixel subtraction.

. Spatial analysis of the GPS ground station distribution
The vertical TEC (VTEC) value is estimated in the IPP from the line-of-sight (integrated) TEC. Since these IPPs can be above sub-ionospheric points (SIPs) outside the study area, the ionosphere can be mapped beyond the area delineated in Fig. 1, as long as the elevation angle is not smaller than 15° (Dach et al. 2007). However, it is still desired that the GNSS sub-network be uniformly distributed.
In order to obtain ideal GNSS station coverage we require that none of the 30 km bu ers around each GNSS station intersect another 30 km bu er (an unnecessarily dense distribution), and that the areas not covered by 60 km bu ers be minimized. Additionally, di culties are presented by polygons extending over the ocean/sea, therefore the coastline is arti cially extended by a 30 km bu er zone. This new delimitation provides a reduced study area but an extended virtual land surface over which more reliable results are expected.
In the Methods section, CR was de ned to provide an additional statistical measure about a GNSS network distribution based on Voronoi polygon areas, but also taking nearest station distances into account. Performing the RA and the CR calculations with the help of GIS, it is now possible to visualize the computed CR results for each polygon (Fig. 2).

. TEC Map Interpolation
TEC values and their RMS were obtained from Bernese GPS 5.0 software as described in section 2.4. Typical RMS values are in the range of 0.1 -0.2 TECU, which corresponds to about 1% of the corresponding TEC values. Various interpolation methods were tested, which are illustrated in Fig. 3, with the TEC values shown on the z-axis. It can be visually observed that kriging and natural neighbor interpolation do not contain visible artifacts and that they are not over-smoothing the topography. Therefore they are optimal methods for our TEC data. Statistical measures were also used to mathematically assess the di erent interpolation methods. RMS deviation values can be seen in Table 2. RMS deviation is the standard deviation of the residuals between interpolated values and estimated values. Kriging interpolation has the lowest RMS deviation value; therefore it was selected as the interpolation method in the remainder of this study.

. Diurnal and seasonal changes in the ionosphere
Only certain ionospheric layers are signi cantly in uenced by solar radiation. Because this work presents a single layer model, these e ects and others are all integrated into a 2D surface above the study area (at an altitude of 450 km), thus representing a superimposition of the di erent layers. It is expected that the TEC maps are providing insights into speci c ionospheric temporal characteristics and quantitative relationships with in uencing factors (Afraimovich et al. 2009). Based on Kom-  Figure 3. Speci c values are based on a particular epoch; however they are representative of the entirety of the data. Kriging interpolation has the best results. jathy (1997), it is expected that there will be constant nonzero ionospheric TEC values from at least the D layer, even during nighttime, due to excitation by cosmic radiation.

IDW(p= ) IDW(p= ) Kriging Natural Neighbor
Since the E layer is strongly dependent on the zenith angle of the Sun, there is also an expected additional diurnal variation component. Furthermore, F1 is expected to add short term changes which are only present during day time. F2 itself does not have solar zenith angle dependence (diurnal variance), but it demonstrates a signi cant seasonal dependence. The global spatial distribution of the F2 layer also reveals a fundamental geomagnetic dependence (Komjathy 1997), which could indicate that this layer is affected by seismo-ionospheric coupling.
To investigate diurnal and seasonal ionospheric variations over southern Iberia, average TEC values over the study area for 1 hour increments were determined. In Fig. 4, hourly average TEC values and daily average TEC values are shown for every 15 th day of 2011 (±2 days when station data was of poor quality or missing). Diurnal changes are apparent, with a noon peak, a night peak, and a possible secondary daytime peak. Maximum TEC values peak around April to May (day 100 to 130) and in September to November (day 266 to 326). Maximum vales are lowest in the Winter, and are also reduced during the Summer months.

. Relation to solar activity
A solar are is an eruption on the Sun that occurs when energy stored in twisted magnetic elds is suddenly released. Flares produce a burst of radiation across the electromagnetic spectrum, from radio waves to X-rays and gamma rays. Solar ares are divided into classes, the most powerful of which are the X-class ares. There were two X-class ares during 2011, one of which is denoted X 6.9 which occurred on August 9, 2011. Fig. 5 displays the hourly averaged TEC values over the study area from August 10 th through August 12 th with the corresponding hourly averaged TEC values from August 9 th subtracted from each value. These days were selected because a coronal mass ejection is expected to reach Earth within one to four days. During the 72 hour period, large TEC maxima were observed centered around hour 20 and hour 45. . Possible relation to seismic activity TEC maps were analyzed in an attempt to detect a precursor or an ionospheric anomaly prior to the May 11, 2011 Lorca earthquake. Data from ve days before the earthquake to three days after the earthquake were processed. To better quantify any unusual activity, Fig. 6 displays the hourly averaged TEC values over the study area from May 7 th through May 14 th , with the corresponding hourly averaged TEC values from May 6 th subtracted from each value. To subtract out any longer-term trend in the data, a cubic polynomial equation was determined to reasonably t to the data, and the residuals are depicted. An unusually large residual is observed beginning at 80 hours, approximately 30 hours before the earthquake. Unusually high TEC values are observed over the study area one day before the earthquake. Also from this data, composite maps for the same epoch were constructed by subtracting out the TEC values from May 6 th in order to compare local variations (Fig. 7). Note that due to the very ne resolution, each map in Fig. 7 employs a di erent color scale. On the day of the earthquake and one day after, a large negative TEC anomaly is observed over the epicenter. Such negative anomalies have been observed in multiple other studies around the time of other earthquakes (Pulinets et al. 2004).

Discussion
In the rst part of this study, a methodology was presented for selecting an optimal GNSS sub-network in order to obtain GPS-derived data for TEC mapping. We can conclude, after analyzing Fig. 2, that an acceptable, regular GNSS ground network has been selected. All but four isolated stations are border stations of the sub-network. One possibility to improve station distribution uniformity would be to include GNSS stations from another network north of the study area. As a nal validation of the selected GNSS network, one can analyze the RMS maps provided in each IONEX le. These values are in the range of 0.1 TECU in most les and higher values only exist if there is a problem with a station data for a given day. This GNSS station selection methodology for ionospheric TEC mapping could be expanded or employed for di erent areas of Earth using the same techniques and principles. Problems may arise when it is applied on island-like geography or coastal areas, where uniformity cannot be guaranteed due to geographic constraints. At high latitudes and near the equatorial region there exist other ionospheric phenomena that require further studies, e.g. scintillations (see for instance Takashi et al. 2014). This work is limited to the mid-latitude ionosphere. Interpolated data were analyzed statistically (Table 2) and visually (Fig. 3). In this case, the only two viable interpolation methods based on statistical measures were kriging and natural neighbor. Most importantly, the RMS deviation values were the lowest when using kriging and natural neighbor, with kriging being slightly lower. Once the optimal interpolation was found, the next step was to construct ionospheric maps and time series.
In the second part of the study, TEC models were built using the selected GNSS sub-network in order to study seasonal and diurnal changes, and also to study how a selected coronal mass ejection (CME) and earthquake affected the ionosphere over the study area. The results suggest that the seasonal variability of TEC shows a semiannual cycle with higher values near the equinox and lower values near the solstice. This would be expected due to changes in the incidence angle of sunlight resulting in variations in excitation of the ionosphere from solar radiation. Similar observations were made by Huang et al. (2006). It is also notable that the yearly cycle in Fig. 4 does not return to its initial value (from 12 TECU to 24 TECU). This suggests that there is a correlation between the averaged TEC values and sunspot number progression (11-Year Solar Cycle), as solar activity increased until 2013 within Solar Cycle 24. Diurnal changes were also observed, with the daytime peak beginning when the solar terminator reaches the study area. During the nighttime the TEC value contours are nearly in horizontal lines (east-west or 90°azimuth), but as the terminator enters the area they can reach an azimuth of 45°or sometimes even down to 0°(north-south or vertical lines). TEC values rise rapidly during the day and only begin to drop during the late afternoon. This is an expected behavior since the E layer is almost entirely dependent on solar activity and the zenith angle of the Sun. Therefore these changes are likely re ecting variations in the E layer. Note that the E layer can practically vanish during the night; therefore the expected observation was a sudden decrease in TEC values after sunset and a constant low until sunrise. Contrarily, a nighttime peak was observed. This was consistent with previous studies such as in Horvath et al. (2000) and Huang et al. (2006). Mid-latitude nighttime TEC increases have been attributed to anomalous increases in the F2 region (Horvath et al. 2000). The physical processes underlying its formation are still not properly understood.
It has been shown that solar activity largely controls seasonal and diurnal cyclic TEC variability. But other ionospheric disturbances can also lead to considerable day-today variations in TEC. One of the sources of these disturbances is CMEs. In this work, the most energetic solar are of 2011 was investigated. It was found that the e ect of the CME on the ionosphere can be identi ed in the TEC data. The strongest in uence (a peak with a maximum of 7 TECU above the August 8 value) was on the day following the CME (August 9) at around 21:00 (UTC) hours, with even the third and the fourth day showing approximately +4 TECU maximum anomalies. In addition to the positive peaks, there were two negative peaks that were 3 to 4 TECU below the August 8 value for the same hour. It should be noted that these values may not be the real maximum and minimum anomalies due to the 1 hour time resolution, thus any short-term variation in TEC could remain undetected. The results show that the maximum values occur at di erent times every day after the solar are event and eventually diminish. Finally, a connection was attempted to be made between TEC trends and the 2011 Lorca earthquake. An attempted was made to nd a correlation between the seismic event and ionospheric anomalies appearing prior to or after its occurrence. The interaction between earthquakes and ionospheric disturbances has been studied for many years (Pulinets et al. 2004). Ionospheric disturbances following major earthquakes have been detected using GPS-derived TEC maps (Pulinets et al. 2004). At the same time, the existence of ionospheric precursors is still not widely accepted by the scienti c community. Therefore, the present study focused on days both before and after the Lorca earthquake. The reference day was May 6th and the last processed day was May 14th. Each of the nine processed days demonstrated typical diurnal variations, although a secondary daytime peak was less apparent on the day of the earthquake and the day after. The daily maximum value was the highest the day before the earthquake, and then dropped signi cantly the following days. An hourly average di erence analysis (Fig. 6) showed cyclic positive/negative uctuations. The highest residuals are the day before the earthquake, with a maximum of approximately 8 TECU. This positive anomaly lasted for nearly 20 hours, however there is also a negative anomaly on May 8. As a comparison, Ouzounov et al. (2011) observed that TEC values reached their maximum 3 days before the Tohoku earthquake. In order to verify that these anomalies are not triggered by solar activity, NOAA solar event reports were analyzed for May 7, May 8, and May 9, and no signi cant solar activity was reported, likely excluding solar ares as a cause of the positive anomaly. Fig. 7 reveals the 2D structure of the local ionospheric anomalies. Such images provide a good sense of the constantly changing spatial structure of the ionosphere above the study area. While solar are interactions with the ionosphere were unlikely, there are still other possible sources for such anomalies which cannot be excluded as yet. Prior results suggest that because phenomena other than solar activity, such as magnetic processes, predominantly a ect nighttime mid-latitude TEC (Horvath et al. 2000), potentially improved precursor detection could be achieved by concentrating only on the nighttime peak. Furthermore, a more detailed study analyzing the correlation of TEC with sunspot numbers should be carried out in order to eliminate this e ect. The time resolution also should be increased from the current 1 hour to a much smaller time interval, possibly even down to 30 seconds.