Source Parameters of the 2016 – 2017 Central Italy Earthquake Sequence from the Sentinel-1 , ALOS-2 and GPS Data

In this study, joint inversions of Synthetic Aperture Radar (SAR) and Global Position System (GPS) measurements are used to investigate the source parameters of four Mw > 5 events of the 2016–2017 Central Italy earthquake sequence. The results show that the four events are all associated with a normal fault striking northwest–southeast and dipping southwest. The observations, in all cases, are consistent with slip on a rupture plane, with strike in the range of 157◦ to 164◦ and dip in the range of 39◦ to 44◦ that penetrates the uppermost crust to a depth of 0 to 8 km. The primary characteristics of these four events are that the 24 August 2016 Mw 6.2 Amatrice earthquake had pronounced heterogeneity of the slip distribution marked by two main slip patches, the 26 October 2016 Mw 6.1 Visso earthquake had a concentrated slip at 3–6 km, and the predominant slip of the 30 October 2016 Mw 6.6 Norcia earthquake occurred on the fault with a peak magnitude of 2.5 m at a depth of 0–6 km, suggesting that the rupture may have reached the surface, and the 18 January 2017 Mw 5.7 Campotosto earthquake had a large area of sliding at depth 3–9 km. The positive static stress changes on the fault planes of the latter three events demonstrate that the 24 August 2016 Amatrice earthquake may have triggered a cascading failure of earthquakes along the complex normal fault system in Central Italy.


Introduction
A sequence of earthquakes from 2016 to 2017 occurred in Central Italy, including ten Mw ≥ 5.0 earthquakes and their aftershocks (Table 1 and Figure 1).The starting event was an Mw 6.2 earthquake at 1:36:36.2 on 24 August 2016 (Greenwich Mean Time, GMT) [1].Approximately 1 h later, an Mw 5.6 earthquake occurred 12 km northwest from the epicenter of the first main event [1].The second main event occurred two months later at 19:18:11 on 26 October 2016 (Mw 6.1) with an Mw 5.5 foreshock located approximately 8 km away [1].Only four days later, the largest earthquake (Mw 6.6) in this sequence occurred approximately between the location of the 24 August and 26 October 2016 events [1].On 18 January 2017, the fourth main event ruptured to the southeast of the 24 August 2016 event with an Mw 5.4 foreshock and two Mw > 5 aftershocks [1].According to the Global Centroid Moment Tensor (GCMT) [1], there are ten earthquakes greater than or equal to Mw 5 during the 2016-2017 Central Italy earthquake sequence.Previous studies have named the first three events after different cities near the epicenter, such as Amatrice or Norcia for the 24 August 2016 event [2][3][4][5][6][7], Visso or Ussita earthquake for the 26 October 2016 event [6,7], and Norcia earthquake for the 30 October 2016 event [6,7].The 2017 event is not named in open publications so far.For ease of analysis, we selected Amatrice (Event A), Visso (Event B), and Norcia (Event C) to name the first three events.For the 2017 event (Event D), we named it the 2017 Campotosto earthquake (Table 1).From the focal mechanism solutions released by different institutions, we can see that these four events all occurred as the result of shallow normal faulting in the Central Apennines (Tables 3−5).Near-field geodetic measurements with good coverage can provide independent constraints on the source parameters with high accuracy.Such high-precision source parameters can enhance our understanding of the tectonics in Central Italy and the internal relationships of these four events.
From seismological and GPS data, it was confirmed that the first event was a NNW-SSE trending, west dipping normal fault event, consistent with the aftershock sequence [2,3].The rupture history shows bilateral propagation and relatively high rupture velocity [2,3].Strong ground shaking was widely felt, and there were substantial landslides on the northwest and northeast sides of the epicenter [4].Lavecchia et al. [5] proposed a two-fault model to explain the coseismic geodetic measurements, although the improvement of data fitting for this two-fault geometry is not significant [4,5].The 26 October 2016 Mw 6.1 Visso earthquake was regarded as a double event, because its foreshock (Mw 5.5) was notably close in time and only 4 km away [6].The 30 October 2016 Norcia earthquake was investigated by different researchers, and there are consistent results in the magnitude of slip ranging from 2.6 m to 3.1 m [3, 6,7].The slip distribution for the 2017 Campotosto earthquake based on modeling of geodetic data is presented here for the first time.In this study, we combined Sentinel-1 Terrain Observation with Progressive Scans (TOPS), Advanced Land Observing Satellite-2 (ALOS-2) and GPS data to determine the source mechanisms for the four earthquakes that occurred in the Central Apennines in 2016-2017 (Table 1).The fault geometry and slip distribution were retrieved from a combination of nonlinear and linear inversions together with the Helmert variance component estimation (HVCE) weighting method [8].The interaction between the earthquakes was also investigated through Coulomb failure stress analysis.The spatiotemporal migration of the Central Italy earthquake sequence was also explored.Finally, we discuss the relationship among the rupture geometry, the fault slip patterns and the internal links of the earthquake sequence.

Tectonic Background
The Mediterranean region is seismically active due to the northward convergence at 4-10 mm/yr of the African plate with respect to the Eurasian plate along a complex plate boundary [9].One of the highest rates of seismicity in the Mediterranean region is found along the Calabrian subduction zone of Southern Italy and Central Italy [9].The 2016-2017 earthquake sequence occurred as the result of shallow normal faulting on a NW-SE oriented fault in the Central Apennines.The Apennines are a mountain range that runs from the Gulf of Taranto in the south to the southern edge of the Po basin in Northern Italy.This region is tectonically and geologically complex, involving both the subduction of the Adria microplate beneath Eurasia and subduction of the Apennines from east to west [9].A dense array of active normal fault systems accommodate approximately 3 mm/yr of extension across the chain [10].Among the main active tectonic structures in the area are the Vettore-Bove Mountains, the Laga Mountains, and the Norcia fault system [11].The central Apennine region has experienced several significant earthquakes in recorded history.The largest instrumentally recorded earthquake within 100 km of the 2016-2017 events was the 13 January 1915 M6.7 earthquake [9].Additionally, this region is the locus of other moderate-sized damaging earthquakes that struck Central Italy in modern times: the Mw 5.8, 1979 Norcia to the west, the Mw 6.0, 1997 Umbria-Marche earthquake sequence to the north, and the Mw 6.1, 2009 L'Aquila sequence to the south.The L'Aquila earthquake resulted in significant landslides in the local area, and was also followed by a vigorous aftershock sequence.The preliminary location of the 2016-2017 earthquake sequence appears to be in a gap between the aftershock sequences of the 1997 and 2009 events.
Since the occurrence of the 24 August 2016 event, approximately more than 20,000 earthquakes of M 2~5 were recorded by Institute of Geophysics and Volcanology (INGV) within ten months

Tectonic Background
The Mediterranean region is seismically active due to the northward convergence at 4-10 mm/year of the African plate with respect to the Eurasian plate along a complex plate boundary [9].One of the highest rates of seismicity in the Mediterranean region is found along the Calabrian subduction zone of Southern Italy and Central Italy [9].The 2016-2017 earthquake sequence occurred as the result of shallow normal faulting on a NW-SE oriented fault in the Central Apennines.The Apennines are a mountain range that runs from the Gulf of Taranto in the south to the southern edge of the Po basin in Northern Italy.This region is tectonically and geologically complex, involving both the subduction of the Adria microplate beneath Eurasia and subduction of the Apennines from east to west [9].A dense array of active normal fault systems accommodate approximately 3 mm/year of extension across the chain [10].Among the main active tectonic structures in the area are the Vettore-Bove Mountains, the Laga Mountains, and the Norcia fault system [11].The central Apennine region has experienced several significant earthquakes in recorded history.The largest instrumentally recorded earthquake within 100 km of the 2016-2017 events was the 13 January 1915 M6.7 earthquake [9].Additionally, this region is the locus of other moderate-sized damaging earthquakes that struck Central Italy in modern times: the Mw 5.8, 1979 Norcia to the west, the Mw 6.0, 1997 Umbria-Marche earthquake sequence to the north, and the Mw 6.1, 2009 L'Aquila sequence to the south.The L'Aquila earthquake resulted in significant landslides in the local area, and was also followed by a vigorous aftershock sequence.The preliminary location of the 2016-2017 earthquake sequence appears to be in a gap between the aftershock sequences of the 1997 and 2009 events.
Since the occurrence of the 24 August 2016 event, approximately more than 20,000 earthquakes of M 2~5 were recorded by Institute of Geophysics and Volcanology (INGV) within ten months (Figure 1) [12].This data set provides an opportunity to track the spatiotemporal evolution of seismicity and to investigate the interactivity of the earthquake sequence.In this study, we only show the general distribution and propagation pattern of the aftershocks.Therefore, we divided the aftershock sequence into four periods, and the rupture times of the four main events were regarded as the time division point.In the preliminary stage of this period (Figure S1a), the aftershocks exhibited a clear bilateral migration, and they were restricted near the epicenters of the 2016 Visso and Norcia earthquakes (green star and blue star).For the second period (Figure S1b), the migration pattern was systematic and formed an ellipse.Although some aftershocks were clustered near the epicenter of the 2016 Norcia earthquake, it did not rupture until four days later.In the third period (Figure S1c), the initial aftershocks covered the largest area and broke the migration limit that occurred in the first period.With the rupture of the 2017 event, the earthquake sequence marched to its ending.The aftershocks were immediately clustered around the epicenter and finally extended for a length of approximately 80 km in the NNW-SSE direction (Figure S1d).
There were some background earthquakes (ML < 2.5) before 24 August 2016, indicating that this area is not completely silent (Figure S2a).The magnitudes of the aftershocks were mainly concentrated near ML 2-3.The seismicity rate before the 24 August 2016 event was very low, and four peaks appeared, corresponding to the four main events (Figure S2b).The aftershocks of the four events all decayed following the Omori-like mode [13] (the aftershocks of the 2016 Visso and Norcia earthquakes mixed together, and it is difficult to separate them, so we did not fit the decay history of the 2016 Visso earthquake period).The decay rate for the 2016 Norcia earthquake period is obviously faster than the other periods (Figure S3).Considering the different decay history of the four periods, we therefore regard these four events as the independent shock type, rather than the main shock-aftershock type.

Sentinel-1 TOPS, ALOS-2 Interferograms, and GPS Data
We collected SAR images covering the 2016-2017 Central Italy earthquake sequence from the Sentinel-1 and ALOS-2 satellite (Table 2), operating on the C band and L band.As part of the Copernicus mission, Sentinel-1 is composed of a constellation of two satellites (A and B) carrying a C-band SAR instrument to assure data continuity provided by the ERS-1/2 and Envisat satellites.Sentinel-1A, launched on 3rd April 2014, operates in four image modes with different resolution and coverage: the Strip Map (SM), Interferometric Wide Swath (IWS), Extra Wide Swath (EWS), and Wave (WV) modes.During these modes, both the IWS and EWS modes employ the novel Terrain Observation with Progressive Scans (TOPS) acquisition mode [14,15].With the TOPSAR technique, in addition to steering the beam in range as in ScanSAR, the beam is electronically steered from backward to forward in the azimuth direction for each burst, avoiding scalloping and resulting in homogeneous image quality throughout the swath [14,16].The original TOPS data are composed of three sub-swaths, with resolutions of 5 m and 29 m along the range and azimuth directions, respectively.ALOS-2, the successor of the ALOS satellite, was launched on 24 May 2014 carrying an L-band PALSAR-2 instrument.The ALOS-2 satellite supports several image models, including the Spotlight mode with a width of 25 km, the Fine mode (Strip Map) with a swath width of 70 km, and ScanSAR (Wide Swath), with a width of 350 km.Compared with the C-band SAR, the PALSAR-2 preserves high coherence in interferometry over regions of rugged topography with high rainfall and dense vegetation, because of its longer wavelength [14,17].The data were processed using the GAMMA InSAR processing software [18].All of the interferograms were generated from the Single Look Complex (SLC) products.The multi-look ratio between the range and azimuth direction was 10:2 for the Sentinel-1 and 8:8 for ALOS-2 data.ALOS-2 Fine mode interferometry is the same as the traditional interferometry.For the TOPS interferometry, an accuracy of several thousandths of a pixel co-registration in the azimuth direction is required [14,19].Otherwise, there would be phase jumps between subsequent bursts.To ensure very high co-registration accuracy, a method for considering the effect of the scene topography and a spectral diversity method [20] considering the interferometric phase of the burst overlap region were used.After the high-quality co-registration between the TOPS SLC data, the topography effects were removed from the interferograms using precise orbits together with the 90 m resolution Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) [21].Slight height-related signals still remained in the interferograms.The interferograms were then filtered using a power spectrum filter [22] to reduce the effect of phase noise, and unwrapped using the branch cut method [23].The interferograms were finally geocoded to the WGS84 geographic coordinates with 90 m resolution.A linear function among the location (x, y), elevation (h), and error phase was estimated with observations away from the deformed areas to remove both the residual orbital errors and topography-correlated atmospheric delays.The decorrelated regions (white areas in the line-of-sight (LOS) displacement maps) extend approximately 1-2 km off the rupture, presumably indicating a zone of high damage and ground disruption around the fault [14].
The coseismic GPS measurements for the first three main shocks were provided by the INGV [24], and the measurements of the 2017 event were provided by the Nevada Geodetic Laboratory.There are 106 stations with three-component coseismic offset estimates for the 24 August 2016 Amatrice earthquake [24] and additional 22 stations released by Cheloni et al. [7], 127 stations for the 26 October 2016 earthquake [24], 114 stations for the 30 October 2016 earthquake [24] and 4 stations for the 18 January 2017 earthquake [25].The coseismic measurements were estimated from 3 days of continuous GPS station positions before and after the main shock, and the stations are around the epicentral area (up to 100 km).Because the GPS offsets are estimated from 3 days of continuous GPS station positions, we should note that the migration of aftershocks, and postseismic (aseismic) slip may affect the GPS offset.

Data Reduction and Weighting
In comparison to GPS and any other conventional geodetic means, the InSAR technique can provide measurements covering epicentral areas at a spatial resolution of several tens of meters.Millions of observations can usually be available for earthquake source modeling, significantly increasing computational cost [14].To reduce the number of points and improve computational efficiency, we adopted the resolution-based quadtree sampling approach [26] that can avoid oversampling in the far fields, due to the effect of observation error.The sampled data density depends on both the sensitivity of data to the fault model and the displacement gradients.
Because several different datasets were used for modeling, an appropriate determination of the weighting for each dataset was required during the inversions.Although the interferogram empirical errors derived from the 1-D covariance function can be used to weight the different datasets directly, Xu et al. [8] suggested that the HVCE method has various advantages, especially when the number of redundant observations is relatively large in a joint inversion.In this study, the HVCE method was employed to weight the different observations, and derive the coseismic slip distribution of the four main Central Italy earthquakes.

Method
We jointly invert SAR and GPS data with a two-step approach: a nonlinear inversion to constrain the fault geometry with uniform slip followed by a linear inversion to retrieve the slip distribution on the fault plane.The observed data are modeled with a finite dislocation fault in an elastic and homogeneous half-space [27].The nonlinear inversion algorithm based on multi-peak particle swam optimization (M-PSO) [28] was used to invert nine fault parameters, including location, strike, dip, rake, length, depth, width, and slip, by minimizing the misfit between the observations and the model predictions, assuming a Poisson ratio of 0.25 and a shear modulus of 3.3 × 10 10 N/m 2 .The interferograms are subsampled to produce a set of points suitable for the inversion.The nonlinear joint inversion of the SAR and GPS data is therefore performed together with the standard deviation for each parameter, as derived from the uncertainty analysis by Monte Carlo simulations.
Under the determined fault geometry, the slip on the rupture plane shows a linear relationship with surface displacement based on classic linear-elastic dislocation theory.Fixing the fault geometry for the optimal fault plane determined in the uniform slip modeling enables the fault length and width to be extended along the strike and down-dip direction, respectively, and then discretized into small patches with an appropriate size.Next, a constrained least-squares method is used to solve the system: where d InSAR and d GPS are the InSAR and GPS observations, G is Green's matrix relating unit slip at the patches to the predicted displacement, S is the slip vector on each patch, ∇ 2 is the second-order finite difference approximation of the Laplacian operator used to avoid unphysical oscillating slip distribution, and k is the smoothing factor.First, we determined the optimal values for smoothness of the model and relative weighting between GPS and InSAR LOS data using the HVCE method as described by Xu et al. [8].In the HVCE method, the relative weight ratio is obtained according to the posterior misfit information from each dataset.The initial inversions were performed using the fault geometry determined in the uniform slip modeling.During inversion, the starting value for the relative weight ratio among the datasets was given as the same weight.After several iterations, the variance of unit weight (standard deviation) became almost uniform for the different datasets.Next, model solutions from 100 simulations perturbed with noise from the statistical properties, based on previous 1-D covariance functions of interferograms and variance of GPS data, were used to estimate the standard deviation from their distribution.

The 2016 Amatrice Earthquake
The Mw 6.2 2016 Amatrice earthquake occurred at 01:36 (UTC) between the towns of Norcia and Amatrice.A significant aftershock (Mw 5.6) followed the main shock, located 12 km to the NW.
The focal mechanism solutions from the United States Geological Survey (USGS) [9], GCMT [1] and the INGV [12] indicated the 2016 Amatrice earthquake was a primarily normal faulting event.
The data used in the inversion include vector displacements measured at 128 GPS stations and LOS displacements (Figure 2) derived from SAR data from four tracks of the Sentinel-1 and ALOS-2 satellite (Figure 1).The ground deformation retrieved from the four unwrapped interferograms is characterized by two NNW-SSE striking deformation lobes, with a maximum negative LOS displacement value of 18.5 cm for Sentinel-1 ascending track 117, 22.2 cm for Sentinel-1 descending track 022, 18.6 cm for ALOS-2 ascending track 197, and 30.4 cm for ALOS-2 descending track 092.

The 2016 Amatrice Earthquake
The Mw 6.2 2016 Amatrice earthquake occurred at 01:36 (UTC) between the towns of Norcia and Amatrice.A significant aftershock (Mw 5.6) followed the main shock, located 12 km to the NW.The focal mechanism solutions from the United States Geological Survey (USGS) [9], GCMT [1] and the INGV [12] indicated the 2016 Amatrice earthquake was a primarily normal faulting event.
The data used in the inversion include vector displacements measured at 128 GPS stations and LOS displacements (Figure 2) derived from SAR data from four tracks of the Sentinel-1 and ALOS-2 satellite (Figure 1).The ground deformation retrieved from the four unwrapped interferograms is characterized by two NNW-SSE striking deformation lobes, with a maximum negative LOS displacement value of 18.5 cm for Sentinel-1 ascending track 117, 22.2 cm for Sentinel-1 descending track 022, 18.6 cm for ALOS-2 ascending track 197, and 30.4 cm for ALOS-2 descending track 092.The inverted causative fault is located at a depth of 2.7 km (refer to the upper center of the rectangle fault plane).This fault is a normal fault dipping approximately 43.6° to the southwest, and the errors of the source parameters are relatively small.The slip is dominantly dip slip, with part leftlateral slip.The inverted geodetic moment magnitude is Mw 6.2, consistent with the GCMT [1] and USGS [9], but larger than the value from INGV [12] (Table 3).Fixing the fault geometry for the optimal fault plane determined in the uniform slip modeling enables the fault length and width to be extended to 36 km along the strike, and 18 km along the down-dip direction, respectively, and then discretized into 648 patches with a size of 1 km × 1 km.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weight ratio converged to 7.2: 2.7:1:1.5:0.8:0.003(Figure S4).The inverted causative fault is located at a depth of 2.7 km (refer to the upper center of the rectangle fault plane).This fault is a normal fault dipping approximately 43.6 • to the southwest, and the errors of the source parameters are relatively small.The slip is dominantly dip slip, with part left-lateral slip.The inverted geodetic moment magnitude is Mw 6.2, consistent with the GCMT [1] and USGS [9], but larger than the value from INGV [12] (Table 3).Fixing the fault geometry for the optimal fault plane determined in the uniform slip modeling enables the fault length and width to be extended to 36 km along the strike, and 18 km along the down-dip direction, respectively, and then discretized into 648 patches with a size of 1 km × 1 km.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weight ratio converged to 7.2: 2.7:1:1.5:0.8:0.003(Figure S4).The total released geodetic moment is approximately 2.16 × 10 18 Nm (Mw 6.2), and the majority of slip occurred on the fault with a peak magnitude of 1.26 m at a depth of 5-6 km (Figure 3).The best-fitting slip model consists of two distinct slip patches located beneath the two observed lobes of deformation, in agreement with the previous studies [2][3][4][5][6][7].The average slip errors are approximately 2.6 cm, except for the area down-dip from the left asperity, where the error increases to 14 cm, which can be attributed to reduced measurement constraints.The data are well reproduced but for some areas near the modeled fault (Figures S5 and S6).This pattern may be caused by the effect of the remaining topographic errors, atmospheric signals in the SAR measurements, and some postseismic deformation in the SAR measurement and GPS data.Nm (Mw 6.2), and the majority of slip occurred on the fault with a peak magnitude of 1.26 m at a depth of 5-6 km (Figure 3).The best-fitting slip model consists of two distinct slip patches located beneath the two observed lobes of deformation, in agreement with the previous studies [2][3][4][5][6][7].The average slip errors are approximately 2.6 cm, except for the area down-dip from the left asperity, where the error increases to 14 cm, which can be attributed to reduced measurement constraints.The data are well reproduced but for some areas near the modeled fault (Figures S5 and S6).This pattern may be caused by the effect of the remaining topographic errors, atmospheric signals in the SAR measurements, and some postseismic deformation in the SAR measurement and GPS data.

The 2016 Visso Earthquake
There is only one Sentinel-1 ascending track 044 interferograms covering the epicentral area of this event (Figure 4a).The maximum observed LOS change in the entire epicentral area reached −21.9 cm.One ALOS-2 ascending track 197 interferograms recorded the cumulative displacement pattern relative to the 2016 Amatrice and Visso earthquake (Figure 4b).Although the ALOS-2 interferogram is spanning both the Amatrice and Visso earthquake, the deformation signal was not mixed together.We then tried to use the deformation part relative to the Visso earthquake.The maximum observed LOS change in the part of 2016 Visso earthquake reached −25.8 cm.
For the 2016 Visso earthquake, the occurrence of clear surface fault scarp can be taken as a model constraint, as the trace of the Vettore fault.We therefore fixed the dip as 40°, and ran the inversion with all other parameters free (Table 4).Using the fault geometry determined with the uniform slip modeling, we extend the fault plane into 1 km by 1 km elements for a total of 540 subfaults.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weight ratio converged to 69.4:1.4:1:0.006(Figure S7).The comparison of observed and modeled SAR and GPS data shows good consistency but for the area near the modeled fault trace (Figures S8 and S9).This may be due to an over-simplification of the model for the complexity of the fault system.
The total geodetic moment in the 2016 Visso earthquake is approximately Nm (Mw 6.14), which is larger than those of USGS [9], GCMT [1], and INGV [12].The slip is almost entirely distributed on a patch of 15 × 6 km 2 , with a peak magnitude of 0.9 m at a depth of 4-5 km (Figure 5).The rake angle averaged over the whole slipping area is −89.7°, and it is primarily normal faulting.

The 2016 Visso Earthquake
There is only one Sentinel-1 ascending track 044 interferograms covering the epicentral area of this event (Figure 4a).The maximum observed LOS change in the entire epicentral area reached −21.9 cm.One ALOS-2 ascending track 197 interferograms recorded the cumulative displacement pattern relative to the 2016 Amatrice and Visso earthquake (Figure 4b).Although the ALOS-2 interferogram is spanning both the Amatrice and Visso earthquake, the deformation signal was not mixed together.We then tried to use the deformation part relative to the Visso earthquake.The maximum observed LOS change in the part of 2016 Visso earthquake reached −25.8 cm.
For the 2016 Visso earthquake, the occurrence of clear surface fault scarp can be taken as a model constraint, as the trace of the Vettore fault.We therefore fixed the dip as 40 • , and ran the inversion with all other parameters free (Table 4).Using the fault geometry determined with the uniform slip modeling, we extend the fault plane into 1 km by 1 km elements for a total of 540 subfaults.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weight ratio converged to 69.4:1.4:1:0.006(Figure S7).The comparison of observed and modeled SAR and GPS data shows good consistency but for the area near the modeled fault trace (Figures S8 and S9).This may be due to an over-simplification of the model for the complexity of the fault system.
The total geodetic moment in the 2016 Visso earthquake is approximately 1.81 × 10 18 Nm (Mw 6.14), which is larger than those of USGS [9], GCMT [1], and INGV [12].The slip is almost entirely distributed on a patch of 15 × 6 km 2 , with a peak magnitude of 0.9 m at a depth of 4-5 km (Figure 5).The rake angle averaged over the whole slipping area is −89.7 • , and it is primarily normal faulting.

The 2016 Norcia Earthquake
The Sentinel-1 ascending track 044 recorded this largest event.The ascending data show major motion towards the satellite, with the LOS value ranging from −50.1 cm to 26.7 cm.One ALOS-2 interferogram also spans the Norcia earthquake, and shows LOS values ranging from −66 cm to −22 cm (Figure 6).The GPS site VETT has the maximum horizontal displacement of 36 cm and 13 cm for the east and north components, respectively.
For the 2016 Norcia earthquake, the occurrence of clear surface fault scarp can be taken as a model constraint, as the trace of the Vettore fault.We therefore fixed the strike as 160°, and ran the inversion with all other parameters free (Table 5).Similar to the procedure described in Section 4.1, we use the fault geometry of the parameters from the best-fitting uniform slip model (Table 5) to estimate the slip distribution for the 2016 Norcia earthquake.The relative weight ratios were determined by the HVCE method mentioned in Section 3.2 and finally converged to 13.4:0.7:1:0.001(Figure S10).The comparison of observed and model SAR and GPS data show good consistency (Figures S11 and S12).

The 2016 Norcia Earthquake
The Sentinel-1 ascending track 044 recorded this largest event.The ascending data show major motion towards the satellite, with the LOS value ranging from −50.1 cm to 26.7 cm.One ALOS-2 interferogram also spans the Norcia earthquake, and shows LOS values ranging from −66 cm to −22 cm (Figure 6).The GPS site VETT has the maximum horizontal displacement of 36 cm and 13 cm for the east and north components, respectively.
For the 2016 Norcia earthquake, the occurrence of clear surface fault scarp can be taken as a model constraint, as the trace of the Vettore fault.We therefore fixed the strike as 160°, and ran the inversion with all other parameters free (Table 5).Similar to the procedure described in Section 4.1, we use the fault geometry of the parameters from the best-fitting uniform slip model (Table 5) to estimate the slip distribution for the 2016 Norcia earthquake.The relative weight ratios were determined by the HVCE method mentioned in Section 3.2 and finally converged to 13.4:0.7:1:0.001(Figure S10).The comparison of observed and model SAR and GPS data show good consistency (Figures S11 and S12).

The 2016 Norcia Earthquake
The Sentinel-1 ascending track 044 recorded this largest event.The ascending data show major motion towards the satellite, with the LOS value ranging from −50.1 cm to 26.7 cm.One ALOS-2 interferogram also spans the Norcia earthquake, and shows LOS values ranging from −66 cm to −22 cm (Figure 6).The GPS site VETT has the maximum horizontal displacement of 36 cm and 13 cm for the east and north components, respectively.
For the 2016 Norcia earthquake, the occurrence of clear surface fault scarp can be taken as a model constraint, as the trace of the Vettore fault.We therefore fixed the strike as 160 • , and ran the inversion with all other parameters free (Table 5).Similar to the procedure described in Section 4.1, we use the fault geometry of the parameters from the best-fitting uniform slip model (Table 5) to estimate the slip distribution for the 2016 Norcia earthquake.The relative weight ratios were determined by the HVCE method mentioned in Section 3.2 and finally converged to 13.4:0.7:1:0.001(Figure S10).The comparison of observed and model SAR and GPS data show good consistency (Figures S11 and S12).The rupture fault is dominated by normal motion with a slight strike-slip component.The total geodetic moment in the 2016 Norcia earthquake was 8.97 × 10 18 Nm (Mw 6.6), consistent with the USGS [9] and GCMT [1] solutions, but slightly larger than INGV [12] (Table 5).The predominant slip occurred on the fault with a peak magnitude of 2.5 m at a depth of 0-6 km, which suggests that the rupture reached the surface (Figure 7).The result is consistent with other papers [3,6,7].
Remote Sens. 2017, 9, 1182 14 of 21 The rupture fault is dominated by normal motion with a slight strike-slip component.The total geodetic moment in the 2016 Norcia earthquake was 18 10 97 .8 × Nm (Mw 6.6), consistent with the USGS [9] and GCMT [1] solutions, but slightly larger than INGV [12] (Table 5).The predominant slip occurred on the fault with a peak magnitude of 2.5 m at a depth of 0-6 km, which suggests that the rupture reached the surface (Figure 7).The result is consistent with other papers [3,6,7].

The 2017 Campotosto Earthquake
Despite the copious snow cover and adverse weather conditions, the Sentinel-1 and ALOS-2 satellite enables us to retrieve the Earth's surface deformation induced by the 2017 event.By processing the images acquired across the main event by the Sentinel-1 and ALOS-2 satellite, we generated the coseismic interferogram, although most parts of the image have a low signal-to-noise ratio (Figure 8).The maximum observed LOS change in the entire epicentral area reaches up to −9.6 cm for Sentinel-1 and −12.2 cm for ALOS-2 (Figure 8).5).The predominant slip occurred on the fault with a peak magnitude of 2.5 m at a depth of 0-6 km, which suggests that the rupture reached the surface (Figure 7).The result is consistent with other papers [3,6,7].

The 2017 Campotosto Earthquake
Despite the copious snow cover and adverse weather conditions, the Sentinel-1 and ALOS-2 satellite enables us to retrieve the Earth's surface deformation induced by the 2017 event.By processing the images acquired across the main event by the Sentinel-1 and ALOS-2 satellite, we generated the coseismic interferogram, although most parts of the image have a low signal-to-noise ratio (Figure 8).The maximum observed LOS change in the entire epicentral area reaches up to −9.6 cm for Sentinel-1 and −12.2 cm for ALOS-2 (Figure 8).

The 2017 Campotosto Earthquake
Despite the copious snow cover and adverse weather conditions, the Sentinel-1 and ALOS-2 satellite enables us to retrieve the Earth's surface deformation induced by the 2017 event.By processing the images acquired across the main event by the Sentinel-1 and ALOS-2 satellite, we generated the coseismic interferogram, although most parts of the image have a low signal-to-noise ratio (Figure 8).The maximum observed LOS change in the entire epicentral area reaches up to −9.6 cm for Sentinel-1 and −12.2 cm for ALOS-2 (Figure 8).
Considering the low signal-to-noise ratio of the SAR data and only four GPS site observations, we decided to use the geological structure and trace of the Gorzano fault as the model constraint.We then fixed the strike and dip as 160 • and 40 • at the same time, and ran the inversion with all other parameters free (Table 6).We also use the same strategy mentioned above to invert the fault geometry and then to detail the slip pattern on the fault.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weighting between the GPS, Sentinel ascending LOS data, ALOS-2 ascending LOS data and the smoothness factor are fixed at 613.5:47.2:1:0.08(Figure S13).The comparison of observed and model Sentinel-1 and GPS data shows good consistency, but the ALOS-2 displacement is not well reproduced, especially in the northeast of the model fault trace (Figures S14 and S15).This may be due to the low signal-to-noise ratio of the original SAR data.
The total released geodetic moment is approximately 1.04 × 10 18 Nm (Mw 5.97).The inferred optimal slip model suggests that the coseismic slip occurred within an area of approximately 15 × 9 km 2 , and the maximum was approximately 0.36 m (Figure 9).Because of reduced measurement constraints, the average slip errors are large around the main slip area.Considering the low signal-to-noise ratio of SAR data and only four GPS site observations, we decided to use the geological structure and trace of the Gorzano fault as the model constraint.We then fixed the strike and dip as 160° and 40° at the same time, and ran the inversion with all other parameters free (Table 6).We also use the same strategy mentioned above to invert the fault geometry and then to detail the slip pattern on the fault.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weighting between the GPS, Sentinel ascending LOS data, ALOS-2 ascending LOS data and the smoothness factor are fixed at 613.5:47.2:1:0.08(Figures S13).The comparison of observed and model Sentinel-1 and GPS data shows good consistency, but the ALOS-2 displacement is not well reproduced, especially in the northeast of the model fault trace (Figures S14 and S15).This may be due to the low signal-to-noise ratio of the original SAR data.
The total released geodetic moment is approximately The inferred optimal slip model suggests that the coseismic slip occurred within an area of approximately 15 × 9 km 2 , and the maximum was approximately 0.36 m (Figure 9).Because of reduced measurement constraints, the average slip errors are large around the main slip area.Considering the low signal-to-noise ratio of the SAR data and only four GPS site observations, we decided to use the geological structure and trace of the Gorzano fault as the model constraint.We then fixed the strike and dip as 160° and 40° at the same time, and ran the inversion with all other parameters free (Table 6).We also use the same strategy mentioned above to invert the fault geometry and then to detail the slip pattern on the fault.We use the HVCE method mentioned in Section 3.2 to determine the relative weight between the datasets and smooth factor.The relative weighting between the GPS, Sentinel ascending LOS data, ALOS-2 ascending LOS data and the smoothness factor are fixed at 613.5:47.2:1:0.08(Figures S13).The comparison of observed and model Sentinel-1 and GPS data shows good consistency, but the ALOS-2 displacement is not well reproduced, especially in the northeast of the model fault trace (Figures S14 and S15).This may be due to the low signal-to-noise ratio of the original SAR data.
The total released geodetic moment is approximately The inferred optimal slip model suggests that the coseismic slip occurred within an area of approximately 15 × 9 km 2 , and the maximum was approximately 0.36 m (Figure 9).Because of reduced measurement constraints, the average slip errors are large around the main slip area.

Static Coulomb Stress Change
A coseismic fault slip induces a stress change to the surrounding area that may affect the state of neighboring faults and seismicity rates.To test whether more recent earthquakes in the Central Italy earthquake sequence may have been triggered by previous ones, we calculated the progression of the Coulomb stress changes using the Coulomb 3.3 software [29] based on the distributed slip source.From the Coulomb failure law, the static Coulomb failure stress ∆CFS on a given receiver fault plane following an event is given by where ∆τ is the fault-parallel shear stress change on the receiver fault plane (positive when sheared in the direction of fault slip), ∆σ is the fault-normal stress change (positive when unclamping), and µ is the effective coefficient of friction.Increasing Coulomb stress, due to the stress readjustment following the event, tends to promote the failure of that plane.Decreasing the Coulomb stress, representing stress release, inhibits the failure.Although the Coulomb stress following an event is relatively smaller than the tectonic stress required for an earthquake, earthquake activity may be promoted when stress is increased on a fault by as little as 0.1 bar [30][31][32].The effective coefficients of friction and the shear modulus were set to be 0.4 and 3.3 × 10 10 N/m 2 , respectively.We examine the internal triggering effects on each earthquake of the 2016-2017 Central Italy earthquake sequence from preceding earthquakes.First, we calculated the static Coulomb stress change induced by the 2016 Amatrice earthquake on the fault plane of the 2016 Visso earthquake (Figure 10a, Table S1).Then, the joint effect of the first two large events (2016 Amatrice and 2016 Visso earthquake) on the 2016 Norcia earthquake was calculated (Figure 10b).Next, the joint effect of the first three events on the 2017 event was calculated (Figure 10c).We find that the latter three earthquakes lay in the stress triggering zones where the Coulomb stress changes are positive, including the 2016 Visso and Norcia, and the 2017 Campotosto earthquakes.Thus, these three mainshocks were promoted by their preceding mainshocks through coseismic Coulomb stress transfer (Figure 10).The 2016 Amatrice earthquake exerts an approximately 0.15 bar Coulomb stress increase (refer to the mean Coulomb stress increase on the main slip area of the Visso earthquake) on the fault plane of the 2016 Visso earthquake (Figure 10a).These two earthquakes exhibited similar strike and fault dip (Tables 3 and 4), suggesting that a connected normal fault system in this area was activated.The 2016 Norcia earthquake received the largest coseismic Coulomb stress change, approximately 6.02 bar, from the 2016 Visso earthquake.The overall coseismic Coulomb stress change that was imparted by the previous two earthquakes (the 2016 Amatrice and Visso earthquakes) was merely 6.08 bar, quantitatively indicating that the 2016 Visso earthquake played a much more important role in promoting the 2016 Norcia earthquake.The 2017 Campotosto earthquake received the largest coseismic Coulomb stress change, approximately 0.24 bar, from the 2016 Amatrice earthquake, and 0.05 bar from the 2016 Norcia earthquake.The overall Coulomb stress increase on the fault of the 2017 Campotosto earthquake exerted by the previous three events was appropriately 0.30 bar.This finding suggests that the 2017 Campotosto earthquake should have been brought closer to failure.

Discussion
Moderate-sized (approximately Mw 5.0-6.6)normal faulting earthquakes are common in the Central Apennines (Figure 1).Before the 2016-2017 Central Italy earthquake sequence, there were two large earthquakes in this region that occurred to the northwest (1997 Umbria-Marche seismic sequence) and the southeast (2009 L'Aquila earthquake).Our slip model shows that the overall slip area of the 2016-2017 earthquake sequence extends approximately 65 km along the strike of the geologically active faults (Figure 11b), and fills the seismic gap between the Umbria-Marche seismic sequence and the L'Aquila earthquake sequence (Figure 1).The inversion results show that the four events aligned to the Vettore and Gorzano fault, and show a variation that is consistent with the locus of the active faults.
For the 2016 Amatrice earthquake, the distributed slip model shows two main distinctive asperities due to the lateral heterogeneity of the friction coefficient on the fault [2][3][4].Our inverted result is consistent with previous studies [2][3][4][5][6][7].Because the GPS and InSAR data can provide independent constraints on the fault parameters in this study, the location of the fault rupture should be more accurate.Our results show a shallower slip patch compared with Chiaraluce et al.'s results, inverted with strong-motion data.For the 2016 Norcia earthquake, the distributed slip model shows that the rupture may reach to the surface.In addition to the main slip, the middle layer, 5-10 km, is found to have a near-uniform slip with a value of approximately 2 m, although there is a relative large uncertainty in the same area.The 2017 Campotosto earthquake ruptured the 10-15 km long segment that had remained unbroken after the previous events.The main shock occurred on the southern end of the 2016-2017 Central Italy earthquake sequence.This area may have been partially activated during the 2016 Amatrice earthquake.In this event, we found a main patch of slip (with peak slip of 0.36 m) located at a depth of 2-8 km.The slip does not reach the surface in the area.

Discussion
Moderate-sized (approximately Mw 5.0-6.6)normal faulting earthquakes are common in the Central Apennines (Figure 1).Before the 2016-2017 Central Italy earthquake sequence, there were two large earthquakes in this region that occurred to the northwest (1997 Umbria-Marche seismic sequence) and the southeast (2009 L'Aquila earthquake).Our slip model shows that the overall slip area of the 2016-2017 earthquake sequence extends approximately 65 km along the strike of the geologically active faults (Figure 11b), and fills the seismic gap between the Umbria-Marche seismic sequence and the L'Aquila earthquake sequence (Figure 1).The inversion results show that the four events aligned to the Vettore and Gorzano fault, and show a variation that is consistent with the locus of the active faults.
For the 2016 Amatrice earthquake, the distributed slip model shows two main distinctive asperities due to the lateral heterogeneity of the friction coefficient on the fault [2][3][4].Our inverted result is consistent with previous studies [2][3][4][5][6][7].Because the GPS and InSAR data can provide independent constraints on the fault parameters in this study, the location of the fault rupture should be more accurate.Our results show a shallower slip patch compared with Chiaraluce et al.'s results, inverted with strong-motion data.For the 2016 Norcia earthquake, the distributed slip model shows that the rupture may reach to the surface.In addition to the main slip, the middle layer, 5-10 km, is found to have a near-uniform slip with a value of approximately 2 m, although there is a relative large uncertainty in the same area.The 2017 Campotosto earthquake ruptured the 10-15 km long segment that had remained unbroken after the previous events.The main shock occurred on the southern end of the 2016-2017 Central Italy earthquake sequence.This area may have been partially activated during the 2016 Amatrice earthquake.In this event, we found a main patch of slip (with peak slip of 0.36 m) located at a depth of 2-8 km.The slip does not reach the surface in the area.The second and the fourth events all had an Mw > 5 foreshock.The Visso earthquake is approximately 30 km northwest of the 2016 Amatrice earthquake, and overlaps with the southern end of the 1997 sequence.A Mw 5.6 earthquake occurred less than two hours before the Visso earthquake, which was approximately just 8 km away.The Mw 5.6 foreshock occurred in close proximity to the area of the aftershocks distribution of the 2016 Amatrice earthquake.The migration of the aftershocks (the 2016 Amatic earthquake) may have triggered this Mw 5.6 earthquake and finally promoted the failure of the 2016 Visso earthquake.Interestingly, the aftershocks that followed the 2016 Visso earthquake also appeared at the location of the 2016 Norcia earthquake, but this area did not rupture until four days later.The 2017 Campotosto earthquake was located at the northern end of the 2009 L'Aquila earthquake sequence.The initial rupture, similar to the 2016 Visso earthquake, may also be related to the Mw 5.4 foreshock.To check if there is any triggering relationship between the four events, we used the slip distribution determined in Section 4.1 to calculate the static Coulomb stress perturbation on the fault plane of the 2016 Visso earthquake and found that there was low stress change at the rupture area of the 2016 Visso earthquake.We also used the slip distribution of the 2016 Amatrice and Visso earthquakes to calculate the static Coulomb stress perturbation on the fault plane of the 2016 Norcia earthquake, and found that the main rupture area of the 2016 Norcia earthquake was located in the area of increasing stress.The 2017 event was also tested, and showed stress increasing in the rupture area.Therefore, there may be two different triggering mechanisms between the four events.The migration of the aftershocks of the 2016 Amatrice earthquake resulted in the rupture of the Mw 5.6 earthquake (we consider this earthquake as the foreshock of the 2016 Visso earthquake).The foreshock promoted the failure of the 2016 Visso earthquake.In addition, the rupture of the 2016 Visso earthquakes led to the stress accumulation in the local region, and finally, resulted in the rupture of 2016 Norcia earthquake.The 2016 Amatrice and Norcia earthquakes together promoted the failure of the 2017 Campotosto earthquake.For the Coulomb stress change of the Central Italy earthquake sequence, Papadopoulos et al. [33] gave a different version.Papadopoulos et al. [33] calculated the regional stress change affected by the earthquakes, while we are to resolve stress imparted to specified faults.We think we calculated the stress change in different situations compared with Papadopoulos et al. [33].The second and the fourth events all had an Mw > 5 foreshock.The Visso earthquake is approximately 30 km northwest of the 2016 Amatrice earthquake, and overlaps with the southern end of the 1997 sequence.A Mw 5.6 earthquake occurred less than two hours before the Visso earthquake, which was approximately just 8 km away.The Mw 5.6 foreshock occurred in close proximity to the area of the aftershocks distribution of the 2016 Amatrice earthquake.The migration of the aftershocks (the 2016 Amatic earthquake) may have triggered this Mw 5.6 earthquake and finally promoted the failure of the 2016 Visso earthquake.Interestingly, the aftershocks that followed the 2016 Visso earthquake also appeared at the location of the 2016 Norcia earthquake, but this area did not rupture until four days later.The 2017 Campotosto earthquake was located at the northern end of the 2009 L'Aquila earthquake sequence.The initial rupture, similar to the 2016 Visso earthquake, may also be related to the Mw 5.4 foreshock.To check if there is any triggering relationship between the four events, we used the slip distribution determined in Section 4.1 to calculate the static Coulomb stress perturbation on the fault plane of the 2016 Visso earthquake and found that there was low stress change at the rupture area of the 2016 Visso earthquake.We also used the slip distribution of the 2016 Amatrice and Visso earthquakes to calculate the static Coulomb stress perturbation on the fault plane of the 2016 Norcia earthquake, and found that the main rupture area of the 2016 Norcia earthquake was located in the area of increasing stress.The 2017 event was also tested, and showed stress increasing in the rupture area.Therefore, there may be two different triggering mechanisms between the four events.The migration of the aftershocks of the 2016 Amatrice earthquake resulted in the rupture of the Mw 5.6 earthquake (we consider this earthquake as the foreshock of the 2016 Visso earthquake).The foreshock promoted the failure of the 2016 Visso earthquake.In addition, the rupture of the 2016 Visso earthquakes led to the stress accumulation in the local region, and finally, resulted in the rupture of 2016 Norcia earthquake.The 2016 Amatrice and Norcia earthquakes together promoted the failure of the 2017 Campotosto earthquake.For the Coulomb stress change of the Central Italy earthquake sequence, Papadopoulos et al. [33] gave a different version.Papadopoulos et al. [33] calculated the regional stress change affected by the earthquakes, while we are to resolve stress imparted to specified faults.We think we calculated the stress change in different situations compared with Papadopoulos et al. [33].

Conclusions and Outlook
The complete SAR Interferometry coverage and adequate GPS data provides a unique opportunity to determine the location, geometry, and focal mechanism of the seismogenic fault.The recorded Central Italy earthquake sequence since 24 August 2016 consists of four main events initiated by an Mw 6.2 earthquake along the Gorzano and Vettore faults.GPS and InSAR observations provide important constraints on the source parameters and slip distribution.The results indicate that the four events are all associated with a normal fault striking northwest-southeast and dipping southwest.The observations, in all cases, are consistent with slip on planar surfaces, with dips in the range of 39 to 44 • , and strikes in the range of 157 to 164 • , that penetrate the uppermost crust to a depth of 0 to 8 km.The main characteristics of these four events are that the 2016 Amatrice earthquake had pronounced heterogeneity of slip distribution marked by two slip patches, the 2016 Visso earthquake had a concentrated slip at 3-6 km, the predominant slip of the 2016 Norcia earthquake occurred on the fault with a peak magnitude of 2.5 m at a depth of 0-6 km, suggesting that the rupture have reached the surface, and the 2017 Campotosto earthquake had a concentrated slip patch at depth 3-9 km.A calculation of the static stress change on the fault planes demonstrates that the 2016 Visso earthquake may have been triggered by the 2016 Amatrice earthquake, the 2016 Norcia earthquake was mainly triggered by the 2016 Visso earthquake, and the 2017 Campotosto earthquake was mainly triggered by the 2016 Amatrice and Norcia earthquakes.
Our work demonstrates natural hazard response in characterizing surface deformation, as well as source slip distribution in a median size continental earthquake at shallow depth.With relatively short temporal and spatial baseline between SAR acquisitions, InSAR can reveal details of the near-to-far-field coseismic displacement, which is important in constraining source slip distribution, as well as assessing surface damage in a short period of time.It is very helpful for quick hazard responses and public access [4,[34][35][36].In recent years, a new generation of SAR satellite missions has been proposed or launched, with the promise to boost our observational capability through reduced revisit times; we think quick hazard responses will be more publicly accessible through the remote sensing.

Figure 1 .
Figure 1.Tectonic setting of the 2016-2017 Central Italy earthquake sequence.Yellow circles are the aftershocks with ML ≥ 2.0 from Institute of Geophysics and Volcanology (INGV) between 24 August 2016 and 30 June 2017.Black boxes outline the spatial coverage of the Sentinel-1 images and Advanced Land Observing Satellite-2 (ALOS-2) images for the 24 August 2016 Amatrice earthquake.Red lines are the active faults from INGV. Green triangle denotes GPS sites used in this paper.The focal mechanism plots represent the moment tensor solutions from Global Centroid Moment Tensor (GCMT).Red focal mechanisms denote the 2016-2017 Central Italy earthquake sequence (Mw ≥ 5); black focal mechanisms denote the history earthquake (1976-2015, Mw ≥ 5).

Figure 1 .
Figure 1.Tectonic setting of the 2016-2017 Central Italy earthquake sequence.Yellow circles are the aftershocks with ML ≥ 2.0 from Institute of Geophysics and Volcanology (INGV) between 24 August 2016 and 30 June 2017.Black boxes outline the spatial coverage of the Sentinel-1 images and Advanced Land Observing Satellite-2 (ALOS-2) images for the 24 August 2016 Amatrice earthquake.Red lines are the active faults from INGV. Green triangle denotes GPS sites used in this paper.The focal mechanism plots represent the moment tensor solutions from Global Centroid Moment Tensor (GCMT).Red focal mechanisms denote the 2016-2017 Central Italy earthquake sequence (Mw ≥ 5); black focal mechanisms denote the history earthquake (1976-2015, Mw ≥ 5).

Figure 3 .
Figure 3. (a) Surface projection of the coseismic slip model of 24 August 2016 Amatrice earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red lines represent the active fault (GF, Gorzano Fault; VF, Vettore Fault).Gray dots denote the aftershocks of ML > 3 from 24 August to 25 October 2016 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia and Amatrice.The red star corresponds to the epicenter of 24 August 2016 Italy earthquake, the green star denotes the epicenter of 26 October 2016 earthquake, the blue star denotes the epicenter of 30 October 2016 earthquake and the golden yellow star denotes the epicenter of 18 January 2017 earthquake [12].

Figure 3 .
Figure 3. (a) Surface projection of the coseismic slip model of 24 August 2016 Amatrice earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red lines represent the active fault (GF, Gorzano Fault; VF, Vettore Fault).Gray dots denote the aftershocks of ML > 3 from 24 August to 25 October 2016 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia and Amatrice.The red star corresponds to the epicenter of 24 August 2016 Italy earthquake, the green star denotes the epicenter of 26 October 2016 earthquake, the blue star denotes the epicenter of 30 October 2016 earthquake and the golden yellow star denotes the epicenter of 18 January 2017 earthquake [12].

Figure 4 .
Figure 4. Coseismic interferograms of 26 October 2016 Visso earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P197A (b).The interferogram are rewrapped with an interval of 6 cm.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake [12].The green star corresponds to the epicenter of 26 October 2016 Visso earthquake [12].The blue rectangle outlines the study area.

Figure 5 .
Figure 5. (a) Surface projection of the coseismic slip model of 26 October 2016 Visso earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red lines represent the active fault.Gray dots denote the aftershocks of ML > 3 from 26 October to 29 October 2016 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia.The green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 4 .
Figure 4. Coseismic interferograms of 26 October 2016 Visso earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P197A (b).The interferogram are rewrapped with an interval of 6 cm.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake [12].The green star corresponds to the epicenter of 26 October 2016 Visso earthquake [12].The blue rectangle outlines the study area.

Figure 4 .
Figure 4. Coseismic interferograms of 26 October 2016 Visso earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P197A (b).The interferogram are rewrapped with an interval of 6 cm.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake [12].The green star corresponds to the epicenter of 26 October 2016 Visso earthquake [12].The blue rectangle outlines the study area.

Figure 5 .
Figure 5. (a) Surface projection of the coseismic slip model of 26 October 2016 Visso earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red lines represent the active fault.Gray dots denote the aftershocks of ML > 3 from 26 October to 29 October 2016 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia.The green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 5 .
Figure 5. (a) Surface projection of the coseismic slip model of 26 October 2016 Visso earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red lines represent the active fault.Gray dots denote the aftershocks of ML > 3 from 26 October to 29 October 2016 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia.The green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 6 .
Figure 6.Coseismic interferograms of 30 October 2016 earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P196A (b).The blue star corresponds to the epicenter of 30 October 2016 Norcia earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 7 .
Figure 7. (a) Surface projection of the coseismic slip model of 30 October 2016 Norcia earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 30 October, 2016 to 17 January, 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia and Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake, the green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 6 .
Figure 6.Coseismic interferograms of 30 October 2016 earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P196A (b).The blue star corresponds to the epicenter of 30 October 2016 Norcia earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 6 .
Figure 6.Coseismic interferograms of 30 October 2016 earthquake obtained from the ascending Sentinel-1 track T044 (a), the ascending ALOS-2 track P196A (b).The blue star corresponds to the epicenter of 30 October 2016 Norcia earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 7 .
Figure 7. (a) Surface projection of the coseismic slip model of 30 October 2016 Norcia earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 30 October, 2016 to 17 January, 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia and Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake, the green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 7 .
Figure 7. (a) Surface projection of the coseismic slip model of 30 October 2016 Norcia earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 30 October 2016 to 17 January 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Norcia and Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake, the green star denotes the epicenter of 26 October 2016 Visso earthquake, the blue star denotes the epicenter of 30 October 2016 Norcia earthquake [12].

Figure 8 .
Figure 8. Coseismic interferograms of 18 January 2017 Campotosto earthquake obtained from the ascending Sentinel-1 track T117 (a), the ascending ALOS track P197 (b).The golden yellow star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 9 .
Figure 9. (a) Surface projection of the coseismic slip model of 18 January 2017 Campotosto earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 18 January 2017 to 30 June 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake and the golden star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].

Figure 8 .
Figure 8. Coseismic interferograms of 18 January 2017 Campotosto earthquake obtained from the ascending Sentinel-1 track T117 (a), the ascending ALOS track P197 (b).The golden yellow star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 8 .
Figure 8. Coseismic interferograms of 18 January 2017 Campotosto earthquake obtained from the ascending Sentinel-1 track T117 (a), the ascending ALOS track P197 (b).The golden yellow star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].The interferograms are rewrapped with an interval of 6 cm.The blue rectangle outlines the study area.

Figure 9 .
Figure 9. (a) Surface projection of the coseismic slip model of 18 January 2017 Campotosto earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 18 January 2017 to 30 June 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake and the golden star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].

Figure 9 .
Figure 9. (a) Surface projection of the coseismic slip model of 18 January 2017 Campotosto earthquake.(b) Standard deviations in slip from Monte Carlo estimation with 100 perturbed datasets.Red line represents the active fault.Gray dots denote the aftershocks of ML > 3 from 18 January 2017 to 30 June 2017 [12].Color shows the slip magnitude in meters, and arrows correspond to the slip directions.White square denotes the city of Amatrice.The red star corresponds to the epicenter of 24 August 2016 Amatrice earthquake and the golden star corresponds to the epicenter of 18 January 2017 Campotosto earthquake [12].

Figure 10 .
Figure 10.(a) Coseismic Coulomb stress change on the fault plane of the 26 October 2016 event triggered by the 24 August 2016 event; green contour lines: the 26 October 2016 event coseismic slip contours every 0.2 m; magenta contour shows the 24 August 2016 event rupture area; gray dots denote the aftershocks of ML > 3 from 24 August to 25 October 2016 [12].(b) Stress change induced on the 30 October 2016 event triggered by 24 August and 26 October 2016 events; black contour lines show the 30 October 2016 event coseismic slip contours every 0.5 m; green and magenta contour show the 26 October and the 24 August 2016 events rupture area, respectively; gray dots denote the aftershocks of ML > 3 from 26 October 2016 to 17 January 2017 [12].(c) Stress change induced on the fault of the 2017 event triggered by the first three events; golden yellow contour lines show the 18 January 2017 event coseismic slip contours every 0.2 m; magenta, green, and black contours show the 24 August, 26 October, and 30 October 2016 events rupture area, respectively; gray dots denote the aftershocks of ML > 3 from 18 January to 30 June 2017 [12].

Figure 10 .
Figure 10.(a) Coseismic Coulomb stress change on the fault plane of the 26 October 2016 event triggered by the 24 August 2016 event; green contour lines: the 26 October 2016 event coseismic slip contours every 0.2 m; magenta contour shows the 24 August 2016 event rupture area; gray dots denote the aftershocks of ML > 3 from 24 August to 25 October 2016 [12].(b) Stress change induced on the 30 October 2016 event triggered by 24 August and 26 October 2016 events; black contour lines show the 30 October 2016 event coseismic slip contours every 0.5 m; green and magenta contour show the 26 October and the 24 August 2016 events rupture area, respectively; gray dots denote the aftershocks of ML > 3 from 26 October 2016 to 17 January 2017 [12].(c) Stress change induced on the fault of the 2017 event triggered by the first three events; golden yellow contour lines show the 18 January 2017 event coseismic slip contours every 0.2 m; magenta, green, and black contours show the 24 August, 26 October, and 30 October 2016 events rupture area, respectively; gray dots denote the aftershocks of ML > 3 from 18 January to 30 June 2017 [12].

Figure 11 .
Figure 11.(a) Three-dimensional views of the fault models of the 2016-2017 Central Italy earthquake sequence (red for the 2016 Amatrice earthquake, green for the 2016 Visso earthquake, blue for the 2016 Norcia earthquake and golden yellow for the 2017 Campotosto earthquake); (b) The coseismic rupture region of the four main events of the 2016-2017 Central Italy earthquake sequence.The purple contours represent the coseismic rupture regions (slip ≥ 0.2 m) of the 2016 Amatrice earthquake.The green contours display the coseismic rupture region (slip ≥ 0.2 m) of Visso earthquake.The blue contours show the coseismic rupture region (slip ≥ 0.5 m) of Visso earthquake.The golden yellow contour denotes the coseismic rupture region (slip ≥ 0.2 m) of 2017 Campotosto earthquake.The stars and gray dots show the main shock and large aftershocks (ML ≥ 3) of the 2016-2017 Central Italy earthquake sequence.

Figure 11 .
Figure 11.(a) Three-dimensional views of the fault models of the 2016-2017 Central Italy earthquake sequence (red for the 2016 Amatrice earthquake, green for the 2016 Visso earthquake, blue for the 2016 Norcia earthquake and golden yellow for the 2017 Campotosto earthquake); (b) The coseismic rupture region of the four main events of the 2016-2017 Central Italy earthquake sequence.The purple contours represent the coseismic rupture regions (slip ≥ 0.2 m) of the 2016 Amatrice earthquake.The green contours display the coseismic rupture region (slip ≥ 0.2 m) of Visso earthquake.The blue contours show the coseismic rupture region (slip ≥ 0.5 m) of Visso earthquake.The golden yellow contour denotes the coseismic rupture region (slip ≥ 0.2 m) of 2017 Campotosto earthquake.The stars and gray dots show the main shock and large aftershocks (ML ≥ 3) of the 2016-2017 Central Italy earthquake sequence.

Table 2 .
Details of satellite synthetic aperture radar (SAR) images used in this study.

Table 3 .
Source parameters of the 24 August 2016 earthquake.
[5]e: Lavecchia et al.[5]has two fault planes for the Amatrice earthquake.The Depth in the Uniform slip model and Distributed slip model denotes the top center of the fault plane.S1: Sentinel-1; CSK: Cosmo-SkyMed.

Table 4 .
Source parameters of the 26 October 2016 earthquake.

Table 5 .
Source parameters of the 30 October 2016 Norcia earthquake.

Table 6 .
Source parameters of the 18 January 2017 earthquake.