Peak horizontal vibrations from GPS response spectra in the epicentral areas of the 2016 earthquake in central Italy

ABSTRACT GPS measurements recorded at active sites in central Italy during the seismic events of greatest magnitude (24 August, 26 and 30 October 2016) were processed in kinematic mode according to the Precise Point Positioning (PPP) technique. The resulting data were the displacements and, by derivation with respect to time, the velocities and instantaneous horizontal accelerations. Elastic response spectra along the orthogonal walls of the site (if the GPS antenna was fixed to a building) or along the geographical directions (if the antenna was fixed to the ground) were obtained from the derived accelerations. The maximum amplitudes, i.e. “peak vibrations”, were then extracted from the response spectra. These peaks, unlike the co-seismic movements, represent the maximum instantaneous vibrations recorded following the “shock” produced by the seismic waves and thus are representative of both the discomfort perceived by the populations and the structural damage. This study shows that GPS is becoming an increasingly important tool to measure and monitor the dynamic responses of a structure. The results also provide a complete picture of the displacements induced by the seismic sequences in the earthquake-affected areas, leaving unresolved some questions concerning the localization of the phenomena and the causes of the structural deformations.


Introduction
In rigid body dynamics, the elastic period of vibration depends on the mass and stiffness (Chopra 1995), while the amplitude of vibration depends on the external force and the damping of the system. From this, it is possible to extend the concept of elastic vibration of a rigid body damped to a single degree of freedom to sites hosting GPS (Global Positioning System) receivers, especially when they are located near the epicentre of strong earthquakes. The earthquake's intensity can be used as an external force to calculate the maximum displacement amplitudes or periods of vibration. This type of study was carried out for the earthquake in Emilia (Italy) in 2012 (Gatti 2018).
In recent years, GPS has proved to be an excellent tool for precise instantaneous positioning, also in dynamic conditions. This has greatly increased its areas of application, especially in structural monitoring (Schaal and Larocca 2009;Yi et al. 2009;Yi et al. 2010aYi et al. , 2010bMoschas and Stiros 2011;Yi et al. 2013a;Breuer et al. 2015;Yigit 2016;Yigit and Gurlek 2017). Hence, it was decided to exploit this ability for GPS stations in central Italy during the 2016 earthquake.
The instantaneous displacements of the phase centres of the antennas (fixed to buildings or to the ground) were calculated from the GPS measurements (Shi et al. 2010;Kobori et al. 2015;Psimoulis et al. 2015): each displacement, or rather its second derivative with respect to time, was considered the external force of a rigid body from which to extract the elastic response spectra of the site and from them the values of the maximum amplitude along two orthogonal directions coinciding either with the walls of the building (antenna fixed to a building) or with the geographical directions (antenna fixed to the ground). These amplitudes represent the peak vibrations, i.e. the main source of the discomfort perceived by the populations as well as the cause of the structural damage and deformations.
The study involved a regional census of GPS recording sites operating in the area of the epicentres of the three above-mentioned seismic events. Hourly recordings centred on the time of each earthquake were extracted and processed in kinematic mode according to the precise point positioning (PPP) technique. The resulting data were the displacements, velocities and instantaneous horizontal accelerations, and subsequently the peak amplitudes of the vibrations (corresponding to the maximum value on the ordinate of the elastic response spectra). The peak vibrations were then compared with the peak ground displacements (PGD) extracted from the recordings of the national accelerometer network (RAN) stations.
Finally, the positions of the centre of gravity of the peak vibrations were determined for each event in order to assess the direction of rotation of the seismic force.
The sample of processed measurements and resulting peak values is particularly significant in view of the variability of the earthquakes' intensity, duration and acceleration (the maximum recorded horizontal components of peak ground acceleration reached 0.5 g in the North-South direction and 0.6 g in the East-West directionsource Istituto Nazionale di Geofisica e Vulcanologia; INGV 2016) as well as the geological characteristics and structural type of the sites (geometry, mass, height, structural constraints, etc.), as reported by (Stewart et al. 2017).
The work is divided as follows: (1) a description of the GPS recording sites; (2) a summary of the literature on the co-seismic movements recorded during the 2016 earthquake in central Italy; (3) processing of the GPS measurements in kinematic mode according to the PPP technique; (4) calculation of the elastic response spectra; (5) the results; (6) comparison with the data provided by the RAN; (7) conclusions.
Unlike the co-seismic movements, which provide the magnitude of superficial deformations of the earth's crust after the seismic event, the values reported in this study quantify the maximum instantaneous vibrations recorded at the sites following the shock produced by the seismic waves. Hence, they are descriptive of both the discomfort perceived by the populations and the structural damage.
Each site can monitor an area of ca. 20 £ 20 km, with a distance from the three epicentres ranging from 3 to 59 km. Figure 1 shows the geographical location of the sites and the location of the epicentres of the 24 August, 26 October and 30 October 2016 earthquakes.
The GPS measurements were recorded by means of a high-performance choke ring antenna (Figure 2(a)) firmly fixed to a building at its highest point (Figure 2 than two storeys high, with the exception of two buildings of four and five storeys each. In other cases, the antenna is fixed to the ground via a reinforced concrete pillar no taller than one metre above ground (Figure 2(c)). The antenna is connected to a dual-frequency receiver with a high sampling rate (Figure 2(d)). Power is supplied either by the mains or by photovoltaic devices. The recordings can be sent remotely or via LAN or wireless: there is always a supplementary power supply unit.  (Macerata) (3.1 cm North-East) and Castello di Campi (Perugia) (2.7 cm South-West), while the Savelli station (Perugia) recorded a 1.7 cm subsidence; for 30 October, the maximum horizontal displacements were observed at Monte Vettore and San Pellegrino (Perugia), with displacements of 38.3 cm North-East and 26 cm South-West, respectively, while the maximum vertical displacements were observed at Arquata del Tronto (Ascoli Piceno), Rifugio Perugia (Perugia) and San Pellegrino (Perugia), with subsidences of 44.6, 26.1 and 17.1 cm, respectively. The Monte Vettore station recorded an uplift of 5.5 cm.

GPS measurement processing
The hourly files recorded during the most intense seismic events having a sampling interval greater than 1 Hz were isolated from the data acquired at the GPS sites, according to the indications of Li et al. (2006), Yi et al. (2013b) and Stiros (2015a, 2015b).
Hourly recordings for the three days were not present at all the sites, either because the receiver had not yet been installed or because the recording was interrupted due to technical problems or the intensity of the event itself. In other cases, the hourly recording was incomplete. Therefore, of the 24 recording sites shown in Figure 1, only 17 had complete records during the 24 August event, 19 sites during the 26 October event and 20 sites during the 30 October event.
The GPS measurements consist of phase measurements and code measurements: they represent the input of the calculation model known in the literature as PPP model (Zumberge et al. 1997;Bertiger et al. 2010). There are various scientific codes that resolve the PPP model: BERNESE (http:// www.bernese.unibe.ch), GAMIT (http://www-gpsg.mit.edu/simon/gtgk), GIPSY (http://gipsy.jpl.nasa. gov/orms/goa) and Coulomb US Geological Survey. For this study, we used GIPSY OASIS II developed by the Jet Propulsion Laboratory (JPL) as single-receiver ambiguity resolution in kinematic positioning (function of time) or PPP (Webb and Zumberge 1996;Yi et al. 2011Yi et al. , 2012Yi et al. , 2013cCheng et al. 2017), RINEX 2.11 input file format, measurement types: dual frequency P code and phase; JPL's precise orbit and clock products in the ITRF08 reference system, cutoff 5 , average PDOP less than 2, IGS standards satellite antenna phase centre offset, antenna type in RINEX input file, tropospheric gradients (Bar-Sever et al. 1998), second-order ionospheric delay (Kedar et al. 2003).
The output of the least squares solution of this model consists of the instantaneous coordinates of the phase centre of the antenna in ECEF ITRF08. For convenience, we chose to express them in the North-East local system whose origin coincides with the antenna's phase centre and the axes coordinated, respectively, with the directions tangent to the meridian and to the parallel passing through the origin. In this way, the instantaneous coordinates subsequent to time t = 0 represent the instantaneous displacements of the antenna, i.e. the vibrations of the building to which it is rigidly fixed: based on the literature they can be considered affected by an error of 0.5 cm/sec (Psimoulis et al. 2008), although recent studies (Xu et al. 2013;Geng et al. 2017) have demonstrated that by combining the GPS measurements with those of the GNSS system (e.g. GPS+GLONASS) one can obtain higher precisions for high sampling rates. For a considerable number of the stations, it was necessary to rotate the North and East instantaneous coordinates along the directions parallel to the building's external walls: the alignment was achieved with a simple rotation equal to the directional angle of one wall with respect to true North. Finally, the coordinates-instantaneous displacements were reduced to zero mean: the latter version or its second derivative with respect to time represented the external force for calculation of the response spectra (discussed in the next section). By way of example, Figure 2 reports the second derivative with respect to time of the instantaneous displacement of the Norcia-NRCI (Perugia) site along the transverse and longitudinal walls following the seismic event of 24 August UTC 01:36:32. By convention, we indicate the walls of the GPS site of shorter length as transverse and those of greater length as longitudinal. The calculation performed for all the stations is shown in Figure 1.

Calculation of the response spectra and extraction of the amplitude
The single building hosting the GPS site was considered to be an elastic oscillator damped to a single degree of freedom (Chopra 1995 ) subjected to the external forces numerically defined in the previous section. Therefore, the numerical calculations were carried out in the frequency domain. First, we performed the Fourier transform of the external forces, reducing the instantaneous displacements or their second derivatives with respect to time from 3600 to 512-1024 s (the latter extracted around the peak instantaneous displacement), and then determined their ratio with the transfer function of the oscillator: the latter is represented by a complex number that depends on the sampling frequency of the GPS measurements, the period T and the damping.
Setting a damping value of 0.05 and a T value, we extracted the maximum value of the ratio at different frequencies. Extraction of the maximum value was carried out n times, each time increasing T by 0.01 s, up to a total of 3 s. The curve obtained by plotting the T values from zero to three on the abscissa and the corresponding maximums (amplitudes) on the ordinate is the response spectrum. The response spectra were calculated in terms of displacement, pseudovelocity and pseudoacceleration from which it was possible to estimate the period of vibration T and the highest amplitude values: the last value represents the maximum vibration at the site during the earthquake. For this calculation (Elnashai and Di Sarno 2008), we developed a dedicated code in the Matlab environment. By way of example, Figure 3 shows the plots of the response spectra along the transverse and longitudinal walls and their mean at the Norcia-NRCI site obtained from accelerations derived from the GPS displacements ( Figure 4). The maximum amplitude of the spectrum calculated along the transverse wall is 0.08 m, corresponding to a period of vibration in the same direction of 0.17 s; the maximum amplitude of the spectrum calculated along the longitudinal wall is 0.06 m, corresponding to a period of vibration of 0.16 s; the mean amplitude is ca. 0.07 m and corresponds to a period of vibration of 0.17 s. North-South and East-West directions (antenna fixed to the ground) for each of the three seismic events.

Results of the measurements
The vectors are centred on the phase centres of the antennas: their direction coincides with the longitudinal and transverse directions for antennas fixed to the buildings, with the North-South and East-West directions for antennas fixed to the ground; since these are absolute values, by convention the direction was chosen as outgoing.
Tables 1-3 report the mean peak vibration for each site calculated in correspondence of the maximum value on the ordinate of the mean spectrum of each site: the order depends on the distance of the site from the epicentre of the considered earthquake.
Since there are two experimental determinations of T (along the two orthogonal directions of the walls of the sites), we made the comparison homogeneous by extrapolating a single experimental value of T equal to the mean of the two spectra.
The same tables also show the mean elastic period of vibration, which is virtually constant for the same site. Table 4 reports, for each seismic event, the comparison between the mean of the instantaneous vibrations from the GPS measurements and the mean of the PGD from the RAN station recordings. For homogeneity of comparison, the values were extracted from the RAN stations located in the same area as that of the GPS stations used in this study.
As seen in Table 4, the maximum difference is less than one centimetre, the minimum difference is two millimetres and the mean of the differences in the East component is zero.   The geographical position of the centre of vibration was also calculated for each earthquake (Figure 8). There is a clear inversion of position between the first and the subsequent events. In the case of the 24 August earthquake, the centre is located to the left of and higher than the epicentre, demonstrating a clockwise rotation induced by the seismic force (red arrow in Figure 8): in this case, a compression of the immersed ground surface was produced to the right of the epicentre and a traction to the left of it. During the October events, the positions were reversed, passing to the right of and lower than the epicentre, demonstrating an anti-clockwise rotation (yellow and green arrows  in Figure 8) induced by the force: consequently the compressed zones were reversed with respect to the stretched ones. Finally, the mean value of the mean periods of vibration calculated at the sites with the antenna fixed to the ground was estimated. This value (0.21 sec) is interesting because it coincides with the value of the first mode of vibration of a two-storey building (T = 0.2 n seconds), where n is the number of storeys with an intrados of 3.00 m.

Conclusions
The present study confirms that the GPS technique is a valid tool to measure or monitor the dynamic responses of structures: new receivers are equipped with boards of up to 100 sampling rates, while the measuring accuracy combined with potent computational algorithms allows estimation of the vibrations even with minimal peaks (ca. 1 cm) such as those found on the buildings most distant from the epicentre. Therefore, GPS station networks can supplement accelerometer networks, especially in areas where the latter are most deficient, or be collocated with them in order to set up real-time warning systems, as advised by Bock et al. (2011) and Geng et al. (2013).
The resulting numerical values (in full agreement with those measured by the RAN stations and with magnitude of up to 9 cm) complete the picture of the vibrations induced by the seismic sequences in the earthquake-affected areas. Unlike the co-seismic movements (which provide only the magnitude of the deformations suffered superficially by the earth's crust after the seismic event), they fully explain the discomfort suffered by the populations and the great structural damage.
Finally, the T value of 0.21 sec is interesting because it coincides with the resonance value of a two-storey building.