Solid Earth DInSAR Coseismic Deformation of the May 2011 M w 5 . 1 Lorca Earthquake ( southeastern Spain )

The coseismic superficial deformation at the region of Lorca (Murcia, southeastern Spain) due to the M w 5.1 earthquake on 11 May 2011 was characterized by a multidisciplinary team, integrating information from DInSAR, GPS and numerical modelling techniques. Despite the moderate magnitude of the event, quantitative information was obtained from the interferometric study of a pair of TerraSAR-X images. The DinSAR results defined the trace of the fault plane and evidenced uplift of the hanging wall block in agreement with the estimated deformation obtained through an elastic rupture dislocation numerical model. Meanwhile for the footwall block, interferometric results showed that tectonic deformation is masked by an important subsidence related to groundwater extraction previously identified at the area of study. Horizontal crustal deformation rates and velocity vectors, obtained from GPS stations existent at the area, were also coherent with the tectonic setting of the southern margin of the Iberian Peninsula and with the focal mechanism calculated for the Lorca event. The analysis of a continuous GPS site in Lorca showed good agreement with the horizontal N– S direction component relative to the numerical model and tectonics of the region. This is the first time at this seismic active area that a multitechnique analysis has been performed immediately after the occurrence of a seismic event, comparing the existing deformation data with a theoretical numerical model based on estimated seismic rupture dislocation.


Introduction
On 11 May 2011, two shallow moderate magnitude earthquakes occurred at less than 5 km northeast of the city of Lorca (Murcia, southeastern Spain).The first event (Mw 4.5) took place at 15:05 (UTC), and had a maximum intensity of VI in the European Macroseismic Scale (EMS).The second and main event (Mw 5.1) occurred at 16:47 (UTC), with an epicentre of coordinates 37.69 • N, 1.67 • W and a depth of 2 km (IGN, 2011), Fig. 1.The main event, which was assigned a maximum intensity of VII (IGN, 2011), did not cause surface rupture, but caused nine human deaths and extensive damages to dwelling buildings, schools and monuments (Irizarry et al., 2011).
Around the area of study there have occurred damaging earthquakes both in the historical and the instrumental periods.The city of Lorca had suffered two damaging earthquakes in 1579 and 1674, with intensities estimated to be VII (EMS) and VIII (EMS), respectively (Martínez-Solares and Mezcua, 2002).In 1829, an earthquake with a maximum intensity of X (EMS) and an estimated magnitude Ms 6.9 (Muñoz and Udías, 1991;Buforn et al., 2006) occurred near the town of Torrevieja (Alicante Province).
On 6 June 1977, an earthquake of magnitude m bLg 4.2 (Mezcua et al., 1984) was registered at 10 km SW from the city of Lorca and it had a maximum intensity of VI (EMS).In the region of Murcia, another damaging earthquake occurred near the town of Mula, 45 km to the NE of Lorca, on 2 February 1999.Its magnitude was M w 5.1 (Buforn et al.,Fig. 1. Study area in SE of Spain, showing the main fault traces in red lines; the GPS stations with their associated horizontal displacement vectors; the epicentral location of the main seismic event (red star) and aftershocks (yellow spots) (IGN, 2011); focal mechanism of the main earthquake (Delouis, 2011); the plate convergence between Eurasia-Nubia plates based on the MORVEL model (DeMets et al., 2010) is shown as a red vector on the lower right corner; the principal strain rates are calculated for the four stations connected by the green dashed lines (see Table 1 for numerical values).The black contour delimits the area for Figs. 3 and 4. AGV: Alto Guadalentín Valley; FAM: Alhama de Murcia Fault; TS: Sierra de la Tercia.
2005), and according to IGN (1999), the maximum intensity reached a value of VI (EMS).On 6 August 2002, an earthquake occurred 40 km west of Mula near the town of Bullas.It had a magnitude of M w 4.6 and a maximum intensity of V (EMS) (Buforn et al., 2005).In the same town on 29 January 2005, there was another earthquake near Bullas (Buforn et al., 2006) with a magnitude of M w 4.7; the IGN (2005) assigned a maximum intensity of VII (EMS).This last earthquake was more damaging than the one in 2002, probably due to the weakening effect produced by the 2002 previous shock on some structures.Though those events had mainly moderate magnitudes (M w ≤ 5.5), they produced damage to structures and noticeable alarm in the population (Gaspar-Escribano et al., 2008).Nevertheless, it should be pointed out that the Lorca event on 11 May 2011 is the first one during the instrumental period that has caused human deaths in the region.
The earthquakes on 11 May 2011 took place in the eastern part of the Betic Cordillera, along the Alhama de Murcia fault (FAM) (Bousquet, 1979).It is a high seismogenic potential strike-slip reverse fault, with a strike between N45 • E and N65 • E, and a movement of about 4-5 mm yr −1 since Neogene times (Masana et al., 2004;Vissiers and Meijninger, 2011).It is located close to the convergent plate limit between the Eurasian and African plates.The convergence direction of this fault has remained constant since late Miocene to present day (Martínez-Díaz, 2002).
According to Delouis (2011), the Lorca main event focal mechanism (Mw 5.1) showed a reverse sinistral motion, compatible with geological and GPS observations.One of the calculated fault planes coincides with the same orientation of the FAM (Fig. 1).
In accordance with IGN (2011), the number of aftershock events decreased quickly through time to less than five events per day in five days.Curiously, these events were not located along the inferred FAM fault plane, dipping towards NNW, but rather they were spread in seemingly non-linear fashion towards the SE from the main shock into the Alto Guadalentín Valley (AGV), Fig. 1.This supposed mislocation might have been due to the generation of aftershocks in a zone with high concentration of the static Mohr-Coulomb stresses out of the FAM plane.

GPS data
The GPS data considered in this study consist of the analysis of one continuous GPS (CGPS) station located in Lorca (LORC) and six campaign data belonging to the CuaTe-Neo (Cuantificación de la Tectónica actual y Neotectónica) geodetic network, which was established in 1996 to quantify the current rates of crustal deformation in the eastern part of the Betic Cordillera (Colomina et al., 1999).It consists of 15 geodetic markers covering an area approximately 6000 km 2 in Murcia and Almeria.The presented data are based on the analysis of four campaigns performed in 1997, 2002, 2006(Khazaradze et al., 2008) ) and 2009.During each campaign every site was observed for three consecutive days (in some cases up to five days) in sessions of at least 8 hours.Some of the sites were observed throughout the campaign in a continuous mode, since it was possible to leave the instruments unattended.
Horizontal velocities of the six CuaTeNeo sites falling within the study area are shown in Fig. 1.These velocities represent an inter-seismic phase of deformation averaged over the 12 years of observations.The stations on the SE side of the FAM (PURI and GANU) show oblique compression with left-lateral direction of motion, 1.9 ± 0.5 mm yr −1 , relative to the stations on the NW (MELL and TERC) in accordance to geological observations and calculated focal mechanism of the Lorca earthquake (Fig. 1).The decomposition of the station PURI vector into fault parallel and perpendicular components, give velocity estimates of 0.29 mm yr −1 and 1.87 mm yr −1 , respectively.This decomposition was performed by projecting the GPS calculated vector of PURI to the trace of the FAM fault, deduced from the focal mechanism and used in the numerical model (N65 • E).In consequence, the relative weight of compression vs. extension is significantly higher than expected from the focal mechanism.On the contrary, when we used an alternative N45 • E orientation for the FAM fault, more in accordance to the geology (see Fig. 1), the fault parallel and perpendicular velocities of PURI became 0.92 mm yr −1 and 1.66 mm yr −1 , respectively.These values are in better agreement with approximately 30 % strike slip and 70 % compression component of the Lorca earthquake focal mechanism (Delouis, 2011).
We also calculated strain rates for the four stations surrounding the earthquake epicentre (TERC, MELL, PURI and GANU) using the SSPX software (Cardozo and Allmendinger, 2009).The obtained results (Fig. 1 and Table 1) are also in good agreement with the calculated focal mechanism of the main event.Specifically, the orientations of ε1 (∼N77 • E) and ε2 (∼N167 • E), deduced from the GPS velocities, agree well with the orientations of the T-axis (N111 • E) and P-axis (N178 • W) (Delouis, 2011).Although this agreement is not always guaranteed or necessary, since the two measurements represent two different physical parameters, strain and stress, which are not always coincident (see Wang (2000) for further details).The NNW orientation of the GPS velocities of the three stations located closer to the coast (MONT, PURI and GANU) indicates that the main driving force behind their motion is related to the convergence between the Nubia and Eurasia plates (see the convergence vector in Fig. 1).This, convergence according to the MORVEL model (DeMets et al., 2010), equals to 5.6 ± 0.3 mm yr −1 and is oriented N322.1 • E ± 1.8 • .These values show that roughly one-third of the overall convergence between the two tectonic plates is taking place on the southern margin of the Iberian Peninsula, at the Betic Cordillera.
These values also agree with other published GPS results in a wider region of SE Spain (e.g.Koulali et al., 2011;Pérez-Peña et al., 2010), which provide velocities ranging between 0 to 2 mm yr −1 oriented in the NNW direction.However, a detailed comparison with our results presented in this paper is not possible, since these studies do not include observations that fall within the geographic area considered here.
Two weeks after the occurrence of the Lorca earthquake, the UB group organized a special post-event campaign of the CuaTeNeo sites located near the epicentre of the earthquake.The preliminary results indicated no signs of detectable coseismic deformation.Nevertheless, our analysis of a continuous GPS site LORC, belonging to the Meristemum public network (Garrido et al., 2011) and located in the city of Lorca (Figure 1), showed a horizontal co-seismic jump of 5 to 6 mm (5.5 ± 0.2 mm) towards the north and statistically marginal subsidence in the vertical direction of −2.0 ± 0.9 mm (Fig. 2, Table 2).It must be pointed out that before the occurrence of the Lorca earthquake, our analysis of the LORC data from 2008 had already indicated a highly anomalous behaviour at this site, which was most likely due to the instability of the terrain surrounding the building where it is located (Echeverria et al., 2011).Specifically, we had estimated that the LORC CGPS station was subsiding with a rate of

LORC GPS
Numerical model S-N direction 5.5 ± 0.2 mm 8.6 mm W-E direction 0.0 ± 0.2 mm 1.0 mm Up direction −2.0 ± 0.9 mm −2.5 mm 98.5 ± 1.9 mm yr −1 (Table 3) and moving horizontally with a rate of 12.7 ± 1.2 mm yr −1 in the SE direction with respect to stable Eurasia.

DInSAR analysis
The use of Synthetic Aperture Radar Differential Interferometry (DInSAR) technique for quantifying coseismic deformations has been previously used at the Africa-Eurasian plate boundary at the western Mediterranean area, e.g. in Morocco (Belabbes et al., 2009;Akoglu et al., 2006); and in the Iberian Peninsula, also at the same seismogenic area in the Betic Cordillera (González et al., 2009).Nevertheless, this is the first time in this area that a processing has been performed immediately after the occurrence of a seismic event and it has been compared to theoretical numerical modelled vertical elastic deformation based on estimated seismic rupture dislocation.
A DInSAR processing (Hanssen, 2001;Mora et al., 2007) of one pre-and one post-event stripmap TerraSAR-X image (25 July 2008 and 14 May 2011) was performed for the Lorca event.In order to reduce temporal decorrelation and avoid non-seismic deformation phenomena, a shorter temporal baseline would be desirable.Unfortunately, that was the only pre-event image available in the TerraSAR-X archive for the study zone.Topography was cancelled employing an interpolated SRTM DTM.Atmospheric effects were considered non-significant, as the detected fringe spatial gradient does not correspond to the typical atmospheric pattern (Hanssen, 2001).
Figure 3a shows the filtered deformation fringes of the differential phase in radians.Each color cycle is equivalent to a deformation of 1.55 centimetres along the line of sight (LOS) of the radar (∼35 • of incidence angle).The well marked fringe pattern is aligned along the trace of the FAM, showing a defined deformation gradient perpendicular to its trace.Further to the S and SE of the epicentral area, the quality of the signal (measured by the coherence parameter) is lower outside the urban areas, mostly associated to agricultural fields in the AGV, showing a concentric low coherence fringe pattern.
A vertical displacement map was generated by unwrapping the phase of the differential interferogram (Costantini, 1998).A high coherence pixel with a zero deformation value according to the numerical model (see next section) was employed to fix the solution.A median filter was applied to the deformation map to reduce the impact of the low coherence pixels (Fig. 3b).
The northern (hanging wall) block of the fault has a maximum upward movement of about 3 cm that agrees with the reported focal mechanism, while the southern (footwall) block of the fault shows a maximum downward movement of 18 cm (Fig. 3b).There is a remarkable difference of order of magnitude between the displacements in each one of the blocks.The limit of these movement tendencies coincides clearly with the FAM trace (Fig. 3a), reflecting also the local change in the strike of the fault from N35 • E to N60 • E and the geological contrasts between the sediments of AGV and the Tertiary rocks of Sierra de la Tercia (TS).
The maximum downward movement of 18 cm in the southern block would represent a constant rate of movement of about 64.2 mm yr −1 .This motion is most likely caused by non-tectonic deformation (Table 3; González et al., 2011) and similar to the terrain instability reported earlier for the LORC CGPS station.It must be pointed out that the accuracy of the calculated deformation may be lower in the southern block than in the northern block.This is justified by the large deformation gradient (compared to the DInSAR signal's wavelength) in this area.As not all the pixels have enough quality to participate in the unwrapping, some fringes may be skipped and the overall deformation may be underestimated.Assuming this, we consider compatible both DIn-SAR's 64.2 mm yr −1 constant rate and the inter-seismic subsidence rate of 98.5 mm yr −1 detected at the LORC CGPS station (Table 3 and Fig. 2).

Numerical model of coseismic vertical deformation
The surface deformation numerical model produced by the M w 5.1 earthquake was generated using the method of Wang et al. (2003), which considers an elastic deformation field.In a first step, the Green's functions are computed for a number of source depths and distances given a layered half-space velocity crustal model.For the Lorca area, the chosen crustal model was taken from Dañobeitia et al. (1998) which consists of 7 layers between 0 and 35 km in depth.
A second step sets the source rectangular rupture surface by defining six fault parameters: slip, length, width, strike, dip and rake of the dislocation.For the present analysis, the first three parameters were set at 15 cm, 4 km and 2 km respectively, attending the mean values for a M w 5.1 given by Wells and Coppersmith (1994) and in agreement with the seismic moment (4.9 × 10 23 dyn × cm) calculated by Delouis (2011).The moment tensor inversion, calculated by the same author, was used to determine the orientation and slip  of the fault, i.e. 245 • for the strike, 65 • for the dip and 58 • for the rake.According to the hypocenter depth of 2 km, and a fault width of 2 km, we assumed that the top of the rupture stops at a depth of 1 km (Fig. 3b).
Considering the geometry described above, the model predicts a maximum vertical deformation of 1 cm subsidence 2 km SE from the epicentre and around 4 cm uplift near the epicentre (Fig. 3b).As shown in Table 3, even if there is a good agreement between DInSAR and numerical model results in the northern block, there is a noticeable discrepancy in the southern one.
It should be pointed out that slip, length and width of the rupture area have been chosen considering the average values proposed by Wells and Coppersmith (1994); the estimation of those parameters assumes a standard value for the shear modulus (µ= 30 GPa) that might not fit in this case study given the zone peculiarities, such as the existence of prefractured rocks and the shallowness of the focal depth.
A lower value for the shear modulus might be considered and a parametric analysis was performed to compute a new vertical surface displacement by considering a lower value of µ with the subsequent different compatible rupture and slip dimensions at diverse hypocentral depths.Finally, we compared the results with those obtained from numerical model and DInSAR results.
Considering the current velocity models used for observational seismology (Delouis et al., 2009) we took a S-wave velocity v s of 1.9 km s −1 and a density ρ of 2 g cm −3 to obtain a lower limit for the shear modulus value µ = ρ v 2 s = 7.22 GPa, which is approximately 6 times lower than the µ considered in the first calculations presented above.Therefore, to maintain the same value of seismic moment tensor, Mo, it is necessary to multiply by 6 the factor sliptimes rupture surface.The model had been constructed by considering two alternative assumptions on the amount of slip: the first one with an average slip of 15 cm and the second with an average slip of 30 cm.For each case, three different focal depths were chosen.For the case which considers an average slip of 15 cm, a hypocentral depth of about 4 km provides a maximum surface displacement of 4.5 cm.And for an average slip of 30 cm and a hypocentral depth of 6 km, maximum surface displacement of 4 cm is obtained.Thus, the results from the µ variations show similar vertical displacement configurations at land surface of the northern block.Therefore, even by considering a lower value of µ, which implies different slip values (between 15 and 30 cm) and hypocentral depths (between 2 and 6 km), we obtained in the northern block synthetic displacements in accordance with the first assumption.

Solid
Considering again the former results, once they were validated by reviewing other µ values, we present two vertical displacement profiles, perpendicular to the strike of the fault plane used for the numerical model and parallel to the schematic section A-A (Fig. 3c).The traces of both profiles are shown in Fig. 3b.Profile I-I trace passes through the maximum of DInSAR vertical movement, mainly concentrated to the northwest of the epicentre.Profile II-II trace passes through the maximum value of the numerical model, very close to the epicentre.In Fig. 3c, computed and DIn-SAR vertical displacements (U z -axe) were compared along the profiles.Positive values are located at the northern block of the fault while the ones in the southern block are negative.Maximum values in the northern block are of the same order of DInSAR and numerical model results, very close to 4.0 cm, though differences are reflected in their distribution: along profile I-I , DInSAR measured movements are constant decreasing rapidly close to the fault trace; while along profile II-II , DInSAR shows more irregularities, decreasing constantly towards the fault trace.
At the southern block differences are more accentuated.The DInSAR movement increases constantly almost at the same rates in both profiles until reaching more than 10.0 cm, while the model results tendency is for a maximum of 0.5 cm and going back zero as it gets away from the fault trace.Clearly in this block there is the superimposition of other vertical movement source added to the earthquake deformation.Unfortunately, separation in the DInSAR "signal" of coseismic deformation component from other terrain instabilities is not possible.Moreover, as the numerical model magnitude is of the order of the DInSAR uncertainty, it is not possible to differentiate coseismic deformation from other type of subsidence.
Due to its acquisition geometry, the sensitivity of the DIn-SAR signal to horizontal displacements is significantly lower than to vertical displacements.In fact it is almost blind to the N-S component.Therefore, even if it is possible to project LOS deformation in its horizontal component due to its high degree of uncertainty, we have not performed any comparison between DInSAR and the modelled horizontal movement components.Therefore, only GPS data has been compared to the horizontal model component.
The horizontal deformation predicted by the numerical model (Fig. 4 and Table 2) agrees well with the co-seismic jump observed at the LORC CGPS station from the analysis of the N-S, E-W and vertical GPS time-series.The examination of the modelled horizontal motion makes easier to appreciate why we were not able to detect any co-seismic motion even at the closest CuaTeNeo station TERC after performing a post-earthquake campaign.

Discussion and conclusions
Using the numerical model as reference for the coseismic displacement, we found a good agreement between the DIn-SAR measurements (3 cm) and the model estimated values (4 cm) on the northern hanging wall block of the fault.This match, as well as the distribution of the vertical movement gradient along the FAM trace, allows to state that the numerical model is a good approximation of the coseismic deformation (Fig. 3b).The largest difference is the areal extent of the deformation and concentration of displacements on the north-western sector of the study area.This difference might be due to changes of geology and definition of local tectonic blocks (Martínez-Díaz, 2002) in the area, and to a possible heterogeneous rupture process not considered in our uniform dislocation model.
On the other hand, there are differences for the southern block with both DInSAR and CGPS (LORC station) results (Table 3).González and Fernández (2011) report important subsidence rates (Table 3) at the AGV sedimentary basin, of about 100 mm yr −1 , due to intensive groundwater extraction, which could be responsible of the large differences obtained by numerical and field techniques (DInSAR and GPS).Considering possible DInSAR movement magnitude underestimation in this area due to low coherence pixels, both DInSAR and LORC measurements show reasonable correspondence.
This intensive groundwater extraction might generate changes of the stress field that could contribute to the generation of the aftershocks (IGN, 2011) out of the FAM trace, but further away within AGV (Fig. 1).While model predicted coseismic vertical deformation is less than few centimetres, as mentioned earlier, subsidence related deformation is of the order of tens of centimetres.Unfortunately, as the model magnitude is of the order of the DInSAR, uncertainty it is not possible to differentiate coseismic deformation from the groundwater extraction related subsidence.The coseismic jump of the LORC station is of the same order as that obtained by the numerical model, and only remarkable in the north direction.
In synthesis, by combining remote sensing measurements (DInSAR), in-situ field measurements (CGPS station LORC) and numerical models of fault rupture, we were able to characterize the coseismic deformation for the 11 May 2011, Mw 5.1 earthquake.This is the first time that results from interferometry technique are obtained and confirmed by a multi-technique and multi-disciplinary study for an earthquake in Spain.

Fig. 2 .
Fig. 2. Time-series of LORC continuous GPS station belonging to the Meristemum network (http://gps.medioambiente.carm.es/).The analysis was performed by the UB group using the GAMIT-GLOBK software from MIT. Reference frame is ITRF08.Vertical green lines depict dates of hardware changes and/or earthquake occurrence.The 1st jump in 2009 is due to an antenna change and the 2nd jump is due to a Lorca earthquake of 11 May 2011.The N-S component, where the co-seismic offset of 5 to 6 mm can be seen, includes a zoom of the time-series.
Fig. 3. (a) Filtered deformation fringes of the differential interferometric phase in radians.Each color cycle is equivalent to deformation along LOS of about 1.55 cm.The limit of the Lorca urban area is shown for reference; (b) Vertical displacement map.The colour scale shows the displacement measured by DInSAR.The black curves show the isovalues obtained by the numerical model in cm.The red rectangle and cross-section line A-A show the dimensions and location of the rupture plane considered in the numerical model.White lines, I-I and II-II , are the traces of vertical displacement profiles shown on Fig. 3c.The limit of the Lorca urban area is shown for reference; (c) Vertical displacement profiles comparing measured DInSAR vs. numerical model results.

Fig. 4 .
Fig. 4. Horizontal coseismic displacement field.A projection of the modelled fault is shown as a rectangle.Green star depicts an epicentre of the 11 May 2011 Mw 5.1 Lorca earthquake.LORC CGPS station includes measured (green) and modelled (red) displacement vectors.TERC CuaTeNeo station location is shown to illustrate a low level of expected co-seismic deformation at this site.

Table 1 .
Strain rate estimates (maximum and minimum principal axis, ε1 and ε2 respectively) based on the GPS velocities of the TERC, MELL, PURI and GANU stations.Negative sign means compression.

Table 2 .
Comparison of coseismic displacement components between the LORC CGPS station and numerical model.Negative vertical displacement means subsidence.

Table 3 .
Comparison of the numerical model vertical displacement results with reported values of maximum displacements and calculated rates from different GPS and DInSAR studies, at both FAM blocks.Negative sign means subsidence.
* Measurements made before the Lorca earthquake on 11 May 2011.