Source complexity of the May 20 , 2012 , MW 5 . 9 , Ferrara ( Italy ) event

A Mw 3.9 foreshock on May 19, 2012, at 23:13 UTC, was followed at 02:03 on May 20, 2012, by a Mw 5.9 earthquake that hit a densely populated area in the Po Plain, west of the city of Ferrara, Italy (Figure 1). Over the subsequent 13 days, six Mw >5 events occurred; of these, the most energetic was a Mw 5.8 earthquake on May 29, 2012, 12 km WSW of the main shock. The tragic balance of this sequence was 17 casualties, hundreds of injured, and severe damage to the historical and cultural heritage of the area. From a seismological point of view, the 2012 earthquake was not an outstanding event in its regional context. The same area was hit in 1996 by a Mw 5.4 earthquake [Selvaggi et al. 2001], and previously in 1986 and in 1967 (DBMI11) [Locati et al. 2011]. The most destructive historical event was the 1570, Imax 8 event, which struck the town of Ferrara [Guidoboni et al. 2007, Rovida et al. 2011]. The 2012 seismic sequence lasted for several weeks and probably developed on a well-known buried thrust fault [Basili et al. 2008, Toscani et al. 2009, DISS Working Group 2010], at depths between 2 km and 10-12 km. […]


Introduction
A M W 3.9 foreshock on May 19, 2012, at 23:13 UTC, was followed at 02:03 on May 20, 2012, by a M W 5.9 earthquake that hit a densely populated area in the Po Plain, west of the city of Ferrara, Italy (Figure 1).Over the subsequent 13 days, six M W >5 events occurred; of these, the most energetic was a M W 5.8 earthquake on May 29, 2012, 12 km WSW of the main shock.The tragic balance of this sequence was 17 casualties, hundreds of injured, and severe damage to the historical and cultural heritage of the area.
From a seismological point of view, the 2012 earthquake was not an outstanding event in its regional context.The same area was hit in 1996 by a M W 5.4 earthquake [Selvaggi et al. 2001], and previously in 1986and in 1967 (DBMI11) [Locati et al. 2011].The most destructive historical event was the 1570, I max 8 event, which struck the town of Ferrara [Guidoboni et al. 2007, Rovida et al. 2011].
The 2012 seismic sequence lasted for several weeks and probably developed on a well-known buried thrust fault [Basili et al. 2008, Toscani et al. 2009, DISS Working Group 2010], at depths between 2 km and 10-12 km.Although the earthquake locations available from the catalog of the Italian Istituto Nazionale di Geofisica e Vulcanologia (INGV; National Institute for Geophysics and Volcanology) (ISIDe: http://iside.rm.ingv.it)do not image a well-defined fault plane, the most reliable interpretation is that of a south-dipping fault, coinciding with the known buried thrust.
Quasi real-time moment tensor estimates [Scognamiglio et al. 2009] of the most energetic events indicate an analogous deformation style, with the compression axis oriented towards the NNE.The moment tensor solution of the foreshock shows a similar reverse faulting style.Figure 1 shows the solutions for the May 20, 2012, M W 5.9 mainshock and the foreshock that occurred a few hours before.The mainshock accelerograms were recorded by the network managed by the Italian Civil Protection Department, and these indicate a peak ground acceleration of 0.303g on the vertical components at the MRN station, located 13 km west of the epicenter.
This study focuses on the analysis of the source mechanism of the May 20, 2012, mainshock.We use data from two distant arrays and four broadband seismic stations surrounding the epicentral area, to infer the characteristics of the finite rupture through multichannel and empirical Green's function (EGF) analysis, respectively.In particular, by exploiting the multichannel coherency of the array data, we obtain constraints on the duration and the gross shape of the source time function.We then focus our attention on the relative source time functions (RSTFs), the time histories and azimuthal variations of which demonstrate a complex pattern of slip and clear directivity effects.

Array analysis
For this part of the study, we use data from two distant clusters of stations.The first array consists of a subnet of the INGV Rete Sismica Nazionale (National Seismic Network), located in the Alta Val Tiberina, central Italy (Figure 1).The second array is composed by seven temporary stations that were deployed in the Larderello geothermal field (Tuscany, Italy), in the framework of the GAPSS (Geothermal Area Passive Seismic Source) experiment [Moretti et al. 2012, this volume].Hereinafter, these two arrays are referred to as the AVT and LRD arrays, respectively.The AVT array has an aperture of the order of 90 km, and average inter-station spacing of about 20 km.The LRD array has an aperture of about 50 km, and an average station spacing <10 km (Figure 1).Both of these arrays are homogeneous in terms of sensor types: Nanometrics Trillium 40s for the AVT array, and Nanometrics Trillium 120s Compact for the LRD deployment.Despite the large aperture of these arrays, the first few

Article history
Received July 17, 2012;accepted August 20, 2012.Subject classification: Empirical Green's functions, Array analysis, Earthquake source dynamics.1).(c) ZLCC as a function of propagation azimuth, for a constant horizontal slowness of 0.12 s/km as measured from the first P-wave arrival (see Table 1).Dashed lines in (b) and (c), theoretical horizontal slowness and propagation azimuth expected for the one-dimensional velocity structure of Li et al. [2007].Figure 3 (right).As for Figure 2, from the LRD array.See Table 1 for the theoretical and measured propagation parameters.

EMILIA EARTHQUAKES
seconds of the P-wave are coherent across the array elements, as depicted in panels (a) of Figures 2 and 3.
For each array, our analysis begins by aligning to the P-wave onset the 0.1 Hz to 1 Hz band-pass filtered velocity seismograms.For this procedure, we adjusted the manually picked arrival times using the inter-station time differences estimated via cross-correlation [Shearer 1997].The propagation parameters obtained from these adjusted arrival times (horizontal slowness and propagation azimuth) are listed in Table 1, together with the theoretical ones for the one-dimensional Earth model of Li et al. [2007].Then, we performed a systematic grid search over the azimuth-horizontal slowness plane.We used a grid extending by ±45˚ and ±0.1 s/km around the theoretical azimuth and horizontal slowness, respectively.For each trial azimuth-slowness pair, we aligned the seismograms to the array barycenter according to a plane-wave model [e.g., Equation (5) in Rost and Thomas 2002], and calculated the array-averaged zero-lag correlation coefficient (ZLCC) [Saccorotti and Del Pezzo 2000] over a 1.5-s-long time window sliding along the time-shifted seismograms with increments of 0.1 s.After this procedure, the ZLCC is expressed as a function of the two propagation parameters.The ZLCC maxima will thus correspond to a coherent, plane-wave crossing the array.
Figures 2 and 3 show the results from this analysis.At both arrays, the most striking feature that emerges from the ZLCC maps is the presence of at least two separate, well-correlated pulses that occur within the first 5 s after the P-wave onset.Both of these pulses are associated with similar propagation parameters, which are also very close to the theoretical ones (see also Table 1).Therefore, we infer that these two pulses are representative of the radiation from the source.On this basis, the total duration of the source time function would be of the order of 5 s and 4 s for arrays AVT and LRD, respectively.In addition, we observe that at AVT, the time separation between the two pulses is larger than that observed at LRD.This evidence appears therefore to suggest an overall westward source directivity, since the apparent duration of the source radiation observed from the SE (array AVT) is longer than that observed from the SW (array LRD).

Empirical Green's functions
EGF deconvolution is widely used for the analysis of ground-motion spectra [e.g., Mori 2003, and references therein], to retrieve source parameters (seismic moment, source radius, rupture duration, rise time and stress drop).In this method, the path, site and instrument effects are efficiently removed by deconvolving the waveform of a lower magnitude event from the main event waveform, with both recorded at the same site, and when the two events are collocated.The main condition to be satisfied by the mainshock and the associated EGF is the similarity of the respective waveforms, which implies closely located hypocenters and a similar focal mechanism.In addition, the EGF should represent a point-like source once compared to the main event.To satisfy this condition, the magnitude difference between the two earthquakes has to be at least one unit [Velasco et al. 1994].
When the source duration of the EGF is sufficiently small as compared to that of the main earthquake, then the former can be approximated to a delta function, and consequently, its waveform approximates the response of the medium between the focus and the measuring site.In principle, the closer the EGF source is to the mainshock, the larger is the frequency band which can be resolved by the method [Patton 1980].
From the theoretical point of view, the integral of the source time function represents the total seismic moment released by the earthquake.Due to rupture directivity, a station located along the fault strike in front of the rupture propagation will record an apparent source time function with a short duration and a large amplitude.Conversely, since the scalar moment has to be constant everywhere, the apparent source time function retrieved behind the rupture propagation will have a long duration and a low amplitude.
Therefore, for mainly unilateral ruptures, the assessment of the fracture propagation direction is quite straightforward, and there is no need for a large number of stations [e.g., Pino and Mazza 2000], providing these latter stations offer sufficient azimuthal coverage of the focal sphere.On the other hand, well-constrained results can also be obtained in the case of more complex seismic sources, even when just a few RSTFs are available [e.g., Pino and Di Luccio 2009].
In our EGF deconvolution, we used data from four stations that are located at different source-to-station directions (Figure 1), to ensure suitable azimuthal coverage.Their distances from the epicenter are between 150 km and 250 km.
To efficiently retrieve the best EGF, we performed a matched-filter analysis by sliding the mainshock recording along the continuous data streams spanning the May 19-22, 2012, time interval.Both time series were first band-pass filtered over the 0.02 Hz to 1 Hz frequency band.For each temporal position of the template wavelet (i.e., the mainshock FERRARA EARTHQUAKE SOURCE COMPLEXITY waveform), we calculated the cross-correlation function, and retained its maximum value.
The procedure is iterated for all of the four stations reported in Figure 1.The M W 3.9 event that occurred on May 19, 2012, a few hours before the mainshock, scores the maximum cross-correlation index (0.85) and we assume it as the EGF for the mainshock.
EGF deconvolution is performed in the frequency domain, and a simple water-level condition was adopted to stabilize the spectral divisions.The results of this process are shown in Figure 4.The RSTFs for the four analyzed stations are reported and located at the corresponding azimuth around the M W 5.9 mainshock focal mechanism.

Discussion and conclusions
Our preliminary results obtained through two independent datasets and methods show that the May 20, 2012, M W 5.9 event shows a relatively complex pattern of rupture.
Array analysis clearly shows the presence of at least two separate pulses that are separated by time intervals of the order of 2 s and 1.5 s at the AVT and LRD arrays, respectively.Although less clear, a third peak is also visible at both arrays some 5 s after the first pulse.Two closely spaced pulses marking the early stages of the rupture process are also shown by the RSTFs (Figure 4, arrows).Although these 'wiggles' are simultaneously observed at all of the RSTFs, thus ruling out the hypothesis that they might result from an artifact of the deconvolution process.The RSTF from station MAGA, located NW of the epicenter, also shows a third peak, the timing of which is consistent with the third pulse observed on the ZLCC maps.Moving towards the NE (station SABO), the envelope of the first two pulses becomes broader, which indicates an anti-directive direction.The same consideration holds for station ASQU, which is located in the SE quadrant.At both of these last two sites, the observation of the third pulse is probably inherited by the widening of the RSTFs.
These preliminary analyses lead us to infer that a great part of the energy was radiated by a source propagating towards the WSW.If we assume that the mainshock was associated with the south-dipping plane, this would correspond to an oblique, down-dip rupture propagation.However, the azimuthal distribution of the RSTF amplitudes suggests a more complex pattern of rupture propagation.Indeed, the overall amplitude of the RSTFs obtained at SABO and ASQU is much larger than what would be expected in the case of unilateral rupture propagation towards the WSW.Therefore, we hypothesize the development of a secondary rupture front, which propagated towards the east, roughly parallel to the fault strike.
DInSAR data [Salvi et al. 2012, this volume] show a main pattern of fringes that correspond to surface deformation compatible with a shallow slip patch located east of the epicenter.However, Salvi et al. [2012] also show a secondary area of ground uplift that is elongated in the WSW direction from the epicenter of the mainshock.Although mapping the superimposed effects of several M W >5 events in the sequence, such ground deformation patterns thus appear to confirm our hypothesis of multi-lateral rupture propagation.Hopefully, further analyses that include additional earthquakes and stations will provide tighter constraints on the dynamics of the rupture processes involved in the seismic sequence.

Figure 1 .
Figure 1.Map of the study area.Red triangles, locations of the four INGV stations used to calculate the relative source time functions; black triangle, location of the MRN accelerometer; black stars, locations of the two arrays, the geometries for which are illustrated in the two insets; black dots, aftershocks located during the first week following the May 20, 2012, mainshock; beach balls, focal solutions for the mainshock and the foreshock, as reported by the Time-Domain Moment Tensor catalog [Scognamiglio et al. 2009; see also http://cnt.rm.ingv.it/tdmt.html].The foreshock signal is used as the EGF (see text for details).

Figure 2
Figure 2 (left).(a) Vertical-component velocity seismograms from array AVT.Traces are aligned to the P-wave arrival times, and adjusted using cross-correlation.(b) ZLCC as a function of horizontal slowness, for a constant propagation azimuth equal to that measured for the first pulse (151.7˚;see Table1).(c) ZLCC as a function of propagation azimuth, for a constant horizontal slowness of 0.12 s/km as measured from the first P-wave arrival (see Table1).Dashed lines in (b) and (c), theoretical horizontal slowness and propagation azimuth expected for the one-dimensional velocity structure ofLi et al. [2007].Figure3(right).As for Figure2, from the LRD array.See Table1for the theoretical and measured propagation parameters.

Figure 4 .
Figure 4. Relative source time functions calculated through EGF deconvolution for the May 20, 2012, M W 5.9 earthquake.Arrows, the two primary pulses that appear to correspond to the correlation peaks shown by the ZLCC maps of Figures 2 and 3.
Figure 4. Relative source time functions calculated through EGF deconvolution for the May 20, 2012, M W 5.9 earthquake.Arrows, the two primary pulses that appear to correspond to the correlation peaks shown by the ZLCC maps of Figures 2 and 3.

Table 1 .
Li et al. [2007]een the propagation parameters of the first P-wave pulse obtained from array analyses and those predicted from theLi et al. [2007]velocity structure.