Coseismic displacement analysis of the 12 November 2017 Mw 7.3 Sarpol-e Zahab (Iran) earthquake from SAR Interferometry, burst overlap interferometry and offset tracking

Interferometric wide-swath mode of Sentinel-1, which is implemented by Terrain Observation by Progressive Scan (TOPS) technique, is the main mode of SAR data acquisition in this mission. It aims at global monitoring of large areas with enhanced revisit frequency of 6 days at the expense of reduced azimuth resolution, compared to classical ScanSAR mode. TOPS technique is equipped by steering the beam from backward to forward along the heading direction for each burst, in addition to the steering along the range direction, which is the only sweeping direction in standard ScanSAR mode. This leads to difficulty in measuring along-track displacement by applying the conventional method of multi-aperture interferometry (MAI), which exploits a double difference interferometry to estimate azimuth offset. There is a possibility to solve this issue by a technique called “Burst Overlap Interferometry” which focuses on the region of burst overlap. Taking advantage of large squint angle diversity of ~1° in burst overlapped area leads to improve the accuracy of ground motion measurement especially in along-track direction. We investigate the advantage of SAR Interferometry (InSAR), burst overlap interferometry and offset tracking to investigate coseismic deformation and coseismic-induced landslide related to 12 November 2017 Mw 7.3 Sarpol-e Zahab earthquake in Iran.

InSAR method fails to derive 3D component of displacement field as it only resolves the across-track component of the ground deformation, which is aligned in Line-of-Sight (LOS) direction from ground to the satellite. One solution to derive 3D deformation is to incorporate multiple across-track acquisitions acquired in right and left looking, i.e. descending and ascending tracks (Fialko, 2004;Froger et al., 2004;Funning et al., 2005). However, by this technique only two across-track components in right and left-looking sides can be inferred. The precise measurement of the along-track displacement would be helpful in resolving 3D displacements from InSAR method. This can be done either with amplitude images, through cross correlation technique (Strozzi et al., 2002), or with phase data, through double-difference interferometric phase (Bechor and Zebker, 2006). The first method is also called offset tracking, and the second is named multiple-aperture interferometry (MAI).
MAI technique is based on the split-beam InSAR processing. In this method the sub-aperture processing is applied on the raw SAR data to produce forward and backward looking SLC pairs with different squint angles, i.e. an offset angle with respect to the zero Doppler direction. Multi-aperture interferograms are then constructed by conjugate multiplication of two different-* Corresponding author looking pairs (forward-looking with backward-looking interferogram).
Any source of decorrelation either due to different squint angle of backward and forward looking interferograms or temporal and/or geometric decorrealation, which affects the coherence of the interferogram, also affects the performance of the MAI method. MAI has the sensitivity to resolve along-track displacement with the accuracy of around one tenth of azimuth pixel size. It means that it is able to map displacements greater than a few decimeter, which restricts its applicability to derive large deformation fields (Grandin et al., 2016).
TOPS imaging mode is a novel wide-swath imaging onboard the Sentinel-1 satellite. To acquire the wide-swath coverage, TOPS sensors use burst-by-burst scanning system, meaning that the TOPS antenna is rotated from backward to forward with a defined rotation rate to scan all the successive bursts. After the double scanning of the first burst, in both forward and backward look, the antenna rotates to illuminate the next burst in a different look angle, with an angular separation of ~1 degree compared to the previous burst. The subsequent bursts have an overlap region corresponding to ~10 percent of the burst length, which ensures that the final image is devoid of any gap.
The fast scanning of the antenna in TOPS imaging causes decreasing of dwell time, the time that an antenna steers an object, leading to the reduction of azimuth resolution (Fattahi et al., 2017). There is a trade-off between the swath width and the azimuth resolution. The width of images in TOPS mode increases by a factor of 3 or 5 compared to the ScanSAR mode, consequently yielding to a coarser azimuth resolution. Therefore, accuracy of the along-track measurement using offset tracking method, as influenced by ~3% of the azimuth resolution (Jung et al., 2014), decreases for Sentinel-1 TOPS wide swath mode. Spectral Diversity is a method proposed to modify the coregistration accuracy by using spectral properties of the complex SAR signal. Instead of cross correlation of the amplitude images, which is an essential part of the common coregistration processing, this method needs phase information of the different spectral look separated by different squint angles (Scheiber and Moreira, 2000). Unlike the MAI method, which only senses along-track shift, SD technique considers phase variation in both azimuth, or along-track, and range, across-track, with respect to Doppler rate behavior over range and azimuth directions. Therefore, it is able to extract the relative misregistration parameters in both range and azimuth directions, as cross correlation analysis does. The same methodology is also employed in Sentinel-1 processing in order to apply precise coregistration. TOPS scan mode takes the advantage of the burst overlap areas, as the pixels located within the overlap area are observed in two different squint angles. The variation in squint angles or equivalently the variation in Doppler centroid frequency is a key parameter to refine SD analysis, as its accuracy depends on a squint angular separation. This technique is referred to as Burst Overlap Interferometry (Grandin et al., 2016). It is also called "Enhanced Spectral Diversity" (ESD) in the studies focusing on refining the coregistration processing (Prats-Iraola et al., 2012) In this study we exploit Sentinel-1 data in IW TOPS mode acquired in both ascending and descending tracks to assess coseismic deformation associated with 12 November 2017, Mw 7.3, SarPol-e Zahab, Iran, earthquake. Interferometry results from ascending and descending orbits are used to extract displacement in both along-track and across-track directions. We use a combination of burst overlap interferometry and offset tracking to overcome the deficiency of offset tracking algorithm in dealing with wide TOPS images. We complement this SAR analysis by applying Bayesian inversion for a dislocation model in an elastic half-space to infer source parameters and slip model of the earthquake (Okada, 1985). We also use offset tracking displacement results in both along-track and across-track direction to evaluate displacement corresponding to a big landslide in the region that was triggered after this event.

DATA PROCESSING
We exploit the existing Sentine1-1 wide-swath data covering the coseismic period of the 12 November 2017, Mw 7.3, SarPol-e Zahab Iran earthquake. Pre-earthquake images include those acquired on 7 November and 11 November 2017 for descending and ascending tracks, respectively. Post-event acquisitions were acquired on 19 November and 17 November 2017 for descending and ascending tracks, respectively.
First conventional interferometric processing was done to derive coseismic across-track displacement map. The interferometric processing starts with interferogram generation that includes several steps including (1) multi-looking, (2) complex multiplication of coregistered images, (3) subtraction of reference phase due to the reference ellipsoid, (4) reduction of phase simulated by an external DEM for topographic correction, and (5) phase unwrapping and geocoding. We use 90-meter SRTM DEM as an external DEM for our processing. All the interferometric processing was done using GAMMA 1 SAR software.
The processing above is followed by burst overlapped interferometry and offset tracking analysis to achieve along-track and across-track deformation. First, the initial offsets are estimated based on ESD analysis in which offset vectors are only localized in burst overlap areas. The bilinear polynomial function is then used to determine offset values all over the image, which are in turn considered as an initial value in offset tracking processing in order to measure final displacement field. We refer to this method as combined offset tracking in the rest of the paper Finally, we invert all the results from interferometry and the combined offset tracking method using a nonlinear inversion procedure to infer the fault geometry parameters. We employ open-source software called Geodetic Bayesian Inversion Software (available at http://comet.nerc.ac.uk/gbis/) to apply non-linear inversion. The software uses Markov-chain Monte Carlo algorithm and Metropolis-Hastings algorithm (Hastings, 1970) to calculate the posterior probability distribution of each unknown parameter. Bayesian inversion allows the fault characteristics vary in geometry to minimize the residual between the modeled and observed deformation. To obtain the slip model we extend the fault width and length along the fault dip and fault azimuth direction, discretize the rectangular fault plane into a grid of 1 km by 1 km patches, and invert for variable dip-slip and strike-slip using the following cost function (Harris and Segall, 1987): Where G is the dislocation Green function, s is the slip, H is the finite difference approximation of the Laplacian operator and  is the smoothing factor that controls the trade-off between data misfit and model roughness. We estimate =0.15 which is determined by use of L-curve plot (Hansen, 1999). Figure 1 shows the across-track coseismic interferogram from InSAR processing for two different orbits corresponding to ascending (Fig. 1a) and descending tracks (Fig. 1b). The colour shows LOS changes in 1 meter, with color cycles blue to red indicating motion towards the satellite. The semielliptical pattern pictured in ascending interferogram (Fig.1 a) indicates the ground moved towards the satellite over an area of about 40 by 60 km located in the southwest of the earthquake epicenter with a maximum LOS motion of ~ 1 meter. In the descending interferogram (Fig. 1b) two different deformation patterns are resolved in the unwrapped interferogram. The across-track movement in both ascending and descending track are consistent with ENE oblique thrust faulting. Figure 2 shows the retrieved deformation field derived from the combined offset tracking method for ascending and descending data. The across-track deformation here (Fig. 2a, 2c) has generally the same pattern as that illustrated in Fig. 1. There are also some local deformation spots exceeding 1 m, which are not seen in interferometry results due to large displacement gradient, mainly related to localized slope instability and rock falls triggered by the earthquake. As seen in Figure 2, the coseismic displacement pattern is highlighted better in across-track observations than along-track ones. This is due to the ramp effect that dominates the along-track displacement, probably caused by ionosphere effects. Such effects cannot be separated from displacement signal through burst overlap interferometry and offset tracking analysis. Therefore, offset tracking results do not resolve much information in along-track direction for this earthquake. For dislocation modelling we only considered across-track measurements. It is worth noting that similar to across-track results, the alongtrack results also exhibit some localized pattern in areas of large displacement gradient, mainly related to landslides and rockfalls triggered by the earthquake. The model fits all the across-track interferograms with an overall RMS of 10 cm. The absence of significant signal in the residual map ( Fig. 3f and Fig. 4f) shows that the vast majority of the coseismic across-track displacement can be modelled through uniform slip solution.  (d-f) are the same as (a-c), but for unwrapped result. The outline of the data used for modelling is depicted in Fig.2 with a dash rectangle labelled A.

RESULTS
We next fix the orientation and location of the fault by those parameters that were obtained using the Bayesian and do the inversion to resolve a distributed slip model for varying strikeslip, dip-slip or rakes parameters. Figure 5 shows the slip distribution map, in which the slip varies smoothly over 70x30 km in the length and width of the fault. Slip distribution map shows a homogenous slip pattern with a maximum slip of approx. 5 meter occurring at a depth of 17 km on the fault plane.
The estimated geodetic moment of the model shown in Figure 5 is 1.06×〖10〗^20, considering a rigidity modulus of 30 Gpa. This is equivalent to a moment magnitude Mw=7.3, which is in agreement with seismic observations. Following the coseismic analysis, we investigate further the applicability of combined offset tracking method to detect large deformation induced by rock falls or landslides in the area. The temporal decorrelation due to the fast ground motion within the landslide area makes the InSAR measurement unfeasible to estimate the landside motion. The combined offset tracking method, however, has a lower sensitivity to temporal decorrelation. In this study we focus on a region depicted by a solid rectangle (rectangle B) in Fig. 2a. Figure 6 shows the displacement filed in both along-track and across-track direction in this region. The results show ground motion of up to 20 m occurred in this area. Field investigation performed following this observation confirmed that a big landslide happened in this region.

CONCLUSION
In this study we presented the applicability of Sentinel-1 wide swath TOPS mode to resolve coseismic displacement and the local ground motion triggered by 12 November 2017 Mw 7.3 Sarpol-e Zahab earthquake in Iran. A combination of repeat-pass InSAR, burst overlap interferometry and offset tracking is exploited to retrieve the along-track and across-track deformation fields. We then used Bayesian-based inversion to infer source parameters and slip model of the earthquake. Our results suggest that the earthquake was generated by a blind ENE oblique thrust faulting with maximum slip of approx. 5 m at a depth of 17 km. In addition we derived deformation field related to a large landslide motion triggered by this event. Due to the lack of the coherence within the landslide area caused by temporal decorrelation, InSAR analysis was not able to measure landslide motion. However, both along-track and across-track displacement fields derived from the combination of offset tracking and burst overlap interferometry could map the landslide pattern and measure a motion around ~20 meter in LOS direction.