Coseismic Slip Distribution and Coulomb Stress Change of the 2023 M W 7.8 Pazarcik and M W 7.5 Elbistan Earthquakes in Turkey

: On 6 February 2023, the M W 7.8 Pazarcik and the M W 7.5 Elbistan earthquakes occurred in southeastern Turkey, close to the Syrian border, causing many deaths and a great deal of property destruction. The Pazarcik earthquake mainly damaged the East Anatolian Fault Zone (EAFZ). The Elbistan earthquake mainly damaged the Cardak fault (CF) and the Do˘gan¸sehir fault (DF). In this study, Sentinel-1A ascending (ASC) and descending (DES) orbit image data and pixel offset tracking (POT) were used to derive surface deformation fields in the range and azimuth directions induced by the Pazarcik and Elbistan earthquakes (hereinafter referred to as the Turkey double earthquakes). Utilizing GPS coordinate sequence data, we computed the three-dimensional surface deformation resulting from the Turkey double earthquakes. The surface deformation InSAR and GPS results were combined to invert the coseismic slip distribution of the EAFZ, CF, and DF using a layered earth model. The results show that the coseismic ruptures of the Turkey double earthquakes were dominated by left-lateral strike-slips. The maximum slip was 7.76 m on the EAFZ and about 8.2 m on the CF. Both the earthquakes ruptured the surface. The Coulomb failure stress (CFS) was computed based on the fault slip distribution and the geometric parameters of all the active faults within 300 km of the M W 7.8 Pazarcik earthquake’s epicenter. The CFS change resulting from the Pazarcik earthquake suggests that the subsequent Elbistan earthquake was triggered by the Pazarcik earthquake. The Antakya fault experienced an increase in CFS of 8.4 bars during this double-earthquake event. Therefore, the M W 6.3 Uzunba˘g earthquake on 20 February 2023 was jointly influenced by the Turkey double earthquakes. Through stress analysis of all the active faults within 300 km of the M W 7.8 Pazarcik earthquake’s epicenter, the Ecemis segment, Camliyayla fault, Aadag fault, Ayvali fault, and Pula segment were all found to be under stress loading. Particularly, the Ayvali fault and Pula segment exhibited conspicuous stress loading, signaling a higher risk of future seismic activity


Introduction
On 6 February 2023, the M W 7.8 Pazarcik earthquake occurred in central and southern Turkey.The earthquake's epicenter was located at (37.226 • N, 37.014 • E), with a focal depth of 17.5 km.The Pazarcik earthquake damaged the Amanos, Pazarcik, and Erkenek segments of the East Anatolian Fault Zone (EAFZ), with a fault length of 300 km.About 9 h later, the M W 7.5 Elbistan earthquake occurred, with the epicenter located at (38.011 • N, 37.196 • E), and a focal depth of 13.5 km.The Elbistan earthquake damaged the Cardak fault (CF) and Do ganşehir fault (DF) [1], with a rupture length of 180 km.Thousands of homes were damaged, and the death toll in Turkey and Syria exceeded 50,000 [2].A comprehensive understanding of the earthquake hazards, fault geometry, energy release, and deformation trends can offer guidance and serve as a reference for earthquake emergency response and post-disaster reconstruction [3][4][5].
The EAFZ is located at the boundary between the Anatolian Microplate and the Arabian Plates [6].The EAFZ is located in southeastern Turkey and is a 700 km-long strike-slip fault zone.Since the 20th century, seven major earthquakes (M W > 6.0) have occurred in the EAFZ (Figure 1 and Table A1).These earthquakes disrupted parts of the EAFZ structure and resulted in two seismic void zones (the Puturge segment and Pazarcik segments) [7].A 45 km rupture formed during the 2020 Mw 6.7 Sivrice earthquake in the Puturge segment [8].According to CFS analysis, the Pazarcik segment may generate destructive earthquakes with Mw ≥ 7.3 [9,10].According to some geological [7,[11][12][13] and geodetic data [14][15][16][17][18], the EAFZ runs northeast between the Pazarcik and Palu segments, and north-northeast on the Amanos segment.The strike-slip rate along the EAFZ exhibits a gradual decrease from the northeast to the southwest.Specifically, the Puturge-Palu segments have a rate ranging from 10 to 13 mm/year.The Pazarcik-Erkenek segments have a rate ranging from 4 to 7 mm/year.The Amanos segment has a rate ranging from 1 to 3 mm/year.The accumulated interseismic strain on the main fault is distributed to a number of large secondary faults [19,20].Although some fault segments in the southwest EAFZ, such as the Amanos, Pazarcik, and Erkenek segments, have low long-term slip rates, there is still the possibility of large earthquakes occurring.
Differential interferometry SAR (D-InSAR) has the advantages of a large range, a high accuracy, and a low cost, and it can monitor cm-level or mm-level deformation [21].However, the accuracy of this method depends on the coherence of the SAR images.When obtaining surface deformation data in areas with poor coherence, the unwrapping results exhibit interference.The pixel offset tracking (POT) method utilizes SAR image intensity information to derive azimuth and range surface deformation offsets with a sub-pixel-level accuracy through a multi-level registration process.Currently, it is widely applied to understand the ground displacements caused by earthquakes [22][23][24].The surface deformation accuracy obtained with the POT method can reach 1/10~1/100 of the resolution of SAR images [25,26].Although the accuracy of the POT method is lower than that of the D-InSAR method, the advantage of POT technology is that it can simultaneously obtain two-dimensional deformation in the azimuth and range directions and does not rely on the coherence of SAR images, which prevents decoherence introducing errors into the deformation results.The POT method does not require phase unwrapping and can avoid the errors caused by phase unwrapping.It has certain advantages when applied to detecting large gradient deformation fields, such as large landslides, earthquakes, and volcanic eruptions.The POT method is particularly suitable for application to the Turkey double earthquakes.
Three-dimensional coseismic surface deformation information about the Turkey double earthquakes has been obtained using the POT method and small-baseline subset persistent scatterer interferometry (SPC) methods [27].Based on the deformation field and coseismic slip distribution of the two earthquakes, some prior researchers have conducted an in-depth analysis of the mutual triggering relationships among the Pazarcik, Elbistan, and Uzunba g earthquakes, and the impact of these earthquakes on the geological faults near the epicenter of the Pazarcik earthquake [2,28].The existing research primarily focuses on the slip model, surface rupture, and seismic risk assessment of geological faults on the EAFZ.However, there has been limited exploration into the quantitative analysis of CFS change in the active faults within 300 km of the Pazarcik earthquake's epicenter.These CFS changes were triggered jointly by the Turkey double earthquakes.A comprehensive assessment of the associated seismic risk is also needed.
ments have a rate ranging from 10 to 13 mm/year.The Pazarcik-Erkenek segments have a rate ranging from 4 to 7 mm/year.The Amanos segment has a rate ranging from 1 to 3 mm/year.The accumulated interseismic strain on the main fault is distributed to a number of large secondary faults [19,20].Although some fault segments in the southwest EAFZ, such as the Amanos, Pazarcik, and Erkenek segments, have low long-term slip rates, there is still the possibility of large earthquakes occurring.The objective of this study is to accurately understand the deformation characteristics and seismogenic structures of the Turkey double earthquakes, utilizing the Turkish active fault database and the coseismic slip model to comprehensively assess the seismic hazards of the surrounding faults.In Section 2, the InSAR and GPS data employed in this study are introduced.In Section 3, the basic principles and algorithmic flow of the POT method, coseismic sliding model, and Coulomb rupture stress are introduced.In Section 4, we describe the computation of the coseismic slip model.Our process combines GPS and InSAR data to compute the model.Our analysis was conducted on the coseismic slip distribution and fault rupture of the fault.We used the focal mechanism solution and coseismic slip distribution results provided by the USGS to calculate the CFS change and study the stress-triggering relationship between the Turkey double earthquakes.We analyzed the seismic risk based on the CFS change in active faults within 300 km of the Pazarcik earthquake's epicenter, which provides a basis for earthquake surface rupture interpretation and the study of seismogenic structures.In Section 5, the work conducted in this study and the relevant results obtained are summarized.

InSAR Data
The SAR dataset includes Sentinel-1A satellite ascending and descending orbits data covering the 2023 Pazarcik and Elbistan earthquake areas before and after the Turkey double earthquake (Table 1) (website: https://scihub.copernicus,accessed on 12 April 2023).The time baselines of the ascending and descending orbit data are both 12 days.The image coverage is shown in Figure 2. The time baselines of the ascending and descending orbit data are both 12 days.The image coverage is shown in Figure 2.

GPS Data
There are 23 GPS stations in the earthquake area.The GPS data were obtained from the Nevada Geodetic Laboratory (website: h p://geodesy.unr.edu/,accessed on 6 March 2023).The data sampling rate is 5 min [29].This dataset includes a time series of deviations in the coordinates of each GPS station from its daily average coordinates.Based on the earthquake occurrence times of the Turkey double earthquakes announced by the United States Geological Survey (UGUS), the three-dimensional surface displacement caused by the Turkey double earthquakes were calculated.
The GPS data were utilized to assess the accuracy of InSAR surface deformation modeling.Since there are only 9 GPS stations within the Sentinel-1 image coverage, we only used data from 9 GPS stations for accuracy verification.The three-dimensional surface displacement of the GPS data was projected onto the range direction using the following equation: where is the distance deformation component; , , and are the three-dimensional surface deformation of GPS; is the satellite flight azimuth angle; is the incident angle; and is the noise.

GPS Data
There are 23 GPS stations in the earthquake area.The GPS data were obtained from the Nevada Geodetic Laboratory (website: http://geodesy.unr.edu/,accessed on 6 March 2023).The data sampling rate is 5 min [29].This dataset includes a time series of deviations in the coordinates of each GPS station from its daily average coordinates.Based on the earthquake occurrence times of the Turkey double earthquakes announced by the United States Geological Survey (UGUS), the three-dimensional surface displacement caused by the Turkey double earthquakes were calculated.
The GPS data were utilized to assess the accuracy of InSAR surface deformation modeling.Since there are only 9 GPS stations within the Sentinel-1 image coverage, we only used data from 9 GPS stations for accuracy verification.The three-dimensional surface displacement of the GPS data was projected onto the range direction using the following equation: where d ran is the distance deformation component; d n , d e , and d u are the three-dimensional surface deformation of GPS; φ is the satellite flight azimuth angle; θ is the incident angle; and δ d ran is the noise.

Pixel Offset Tracking
The surface deformation fields of the Turkey double earthquakes were obtained using the POT method.First, SAR images before and after the earthquakes were precisely registered using orbital and terrain data from the Shuttle Terrain Radar Mission.Then, cross-correlation calculations were performed comparing the SAR amplitude images to measure the image offset caused by the earthquakes.Since the Sentinel-1A SLC image (resolution: range × azimuth = 2.33 × 14 m) has different resolutions in the range and azimuth directions, a 64 × 16 window (range direction × azimuth direction, about 150 m × 224 m) was employed for image matching.The images were oversampled by a factor of 64 during the matching process to extract offsets with an accuracy of 1/64 pixel, but due to the influence of noise, the accuracy should not be above theoretical accuracy.In order to further reduce the noise in the offset, the offset results in the range and azimuth directions were subjected to median filtering, with a filter width of 500 m.The results were geocoded with 0.0002 • (~22 m) longitude and latitude sample spacing.

Coseismic Slip Distribution Inversion
The inversion of the coseismic slip distribution on the seismogenic fault was conducted by integrating GPS and InSAR data.Due to the 12-day time span of the InSAR deformation results, distinguishing the coseismic deformations caused by the Turkey double earthquakes was impractical.GPS data, characterized by a high-precision and continuous time series, can be employed to monitor surface deformations before and after earthquakes.This capability allows for the more effective differentiation of the surface deformations caused by the dual earthquakes in Turkey, thus providing accurate input data for coseismic slip inversion.However, typically, the number of GPS stations is limited, resulting in a lower spatial resolution.Additionally, the stations located in close proximity to the epicenter may be adversely affected by the earthquake, leading to data gaps or errors.Therefore, it is necessary to integrate the GPS and InSAR results to obtain a more comprehensive and precise assessment of coseismic slip.The accuracy of the azimuth deformation results from InSAR is limited; therefore, in the inversion process, only the range deformation results from the ascending and descending orbits were utilized.To mitigate spatial correlation in the InSAR deformation outcomes and reduce the influence of far-field noise on inversion results, while simultaneously enhancing computational efficiency, quadtree downsampling was implemented [30], with the variance threshold set at 0.008 2 m.The coseismic slip distribution was computed using the steepest descent method (SDM) [31].The model is as follows: where G is the Green's function based on the layered earth model; s is the slip of the underground fault plane; y is the surface deformation of InSAR and GPS; and ε is the residual between the model and the surface deformation results.In order to obtain highresolution results, the fault plane was discretized into 5 km × 5 km sub-fault blocks.Since the surface deformation results of InSAR and GPS are non-uniform and the Green's function matrix is usually rank-deficient or ill-conditioned, additional constraints are needed to ensure the stability of the inversion.The smoothing constraint we adopted is as follows: where y is the surface deformation of InSAR and GPS; s is the slip on each sub-fault block; G is the Green's function; H is the second-order Laplace point difference operator; α is the optimal smoothing factor; and W is the weight matrix of the observed values.
Figure 3 shows a compromise curve illustrating the trade-off between the roughness of the coseismic slip model and the degree of data fit under varying smoothing factors.The optimal smoothing factor was selected at the inflection point [32].It was conclusively determined that the optimal smoothing factor for the Turkey double earthquakes was 0.15.The layered earth model was used to study the impact of layering of the earth's medium on the coseismic deformation fields.Following the existing research [33][34][35], the parameters we established for the layered earth model are outlined in Table A2.These parameters are incorporated into the layered earth model to calculate the Green's function G.The weight matrix W was determined based on the ratio of the sampling points in the InSAR to GPS data.
Remote Sens. 2024, 16, x FOR PEER REVIEW 6 of 18 determined that the optimal smoothing factor for the Turkey double earthquakes was 0.15.The layered earth model was used to study the impact of layering of the earth's medium on the coseismic deformation fields.Following the existing research [33][34][35], the parameters we established for the layered earth model are outlined in Table A2.These parameters are incorporated into the layered earth model to calculate the Green's function G.The weight matrix W was determined based on the ratio of the sampling points in the InSAR to GPS data.

Coulomb Failure Stress
An earthquake rupture is a process of stress release and readjustment.Positive and negative stress changes correspond to the loading and unloading of fault loads, respectively, influencing the occurrence of earthquakes.Analyzing the CFS change in strong earthquakes is an effective method for predicting future earthquake risks in the vicinity of an epicenter.Utilizing the Coulomb 3.4 software package [36], and based on the published Coulomb rupture accuracy and slip distribution parameters for Turkey double earthquakes by UGUS, we calculated the coseismic static CFS change.The formula is as follows: where ∆ is the CFS change; ∆ is the shear stress change in the fault sliding direction, which can be obtained according to the stress change tensor; ∆ is the shear stress change; and is the effective friction coefficient, and its typical value is 0.4 [37].

Coseismic Surface Deformation
The coseismic horizontal GPS deformation caused by the Turkey double earthquakes is shown in Figure 4a,b.The most significant horizontal displacement during the Pazarcik earthquake was observed at the ANTP station, where a substantial 0.4 m movement occurred in the northeast direction.The stations situated to the north of the EAFZ primarily exhibited southwestward motions, while the stations to the south predominantly shifted in a northeastern direction.These observations suggest that the fault's primary mode of motion may be characterized by strike-slip faulting in the northeast direction.The most significant horizontal displacement during the Elbistan earthquake was observed at the EKZ1 station, where a substantial 4.4 m movement occurred in the west direction.The northern side of the earthquake's epicenter predominantly shifted westward, while the southern side primarily moved eastward; therefore, the fault is mainly a left-lateral strike-

Coulomb Failure Stress
An earthquake rupture is a process of stress release and readjustment.Positive and negative stress changes correspond to the loading and unloading of fault loads, respectively, influencing the occurrence of earthquakes.Analyzing the CFS change in strong earthquakes is an effective method for predicting future earthquake risks in the vicinity of an epicenter.Utilizing the Coulomb 3.4 software package [36], and based on the published Coulomb rupture accuracy and slip distribution parameters for Turkey double earthquakes by UGUS, we calculated the coseismic static CFS change.The formula is as follows: where ∆ σ f is the CFS change; ∆ τ is the shear stress change in the fault sliding direction, which can be obtained according to the stress change tensor; ∆ σ n is the shear stress change; and µ ′ is the effective friction coefficient, and its typical value is 0.4 [37].

Coseismic Surface Deformation
The coseismic horizontal GPS deformation caused by the Turkey double earthquakes is shown in Figure 4a,b.The most significant horizontal displacement during the Pazarcik earthquake was observed at the ANTP station, where a substantial 0.4 m movement occurred in the northeast direction.The stations situated to the north of the EAFZ primarily exhibited southwestward motions, while the stations to the south predominantly shifted in a northeastern direction.These observations suggest that the fault's primary mode of motion may be characterized by strike-slip faulting in the northeast direction.The most significant horizontal displacement during the Elbistan earthquake was observed at the EKZ1 station, where a substantial 4.4 m movement occurred in the west direction.The northern side of the earthquake's epicenter predominantly shifted westward, while the southern side primarily moved eastward; therefore, the fault is mainly a left-lateral strike-slip along the east-west direction.The cumulative horizontal deformation caused by the two earthquakes according to 23 GPS stations near the epicenter is shown in Figure 4c.Based on the deformation measurements from both the ascending (T014) and descending (T021) orbits in the range and azimuth directions (see Figure 5), the deformation directions on either side of the seismogenic fault for the two earthquakes are opposing.The opposite surface deformations on both sides of the faults indicate that the Turkey double earthquakes were left-lateral strike-slip earthquakes.The range displacement results are significantly be er than the azimuth displacement results.This is because the azimuth resolution of the Sentinel-1A data is 14 m and the range resolution is 2.33 m.The range resolution is much higher than the azimuth resolution.The azimuth deformation characteristics of the Elbistan earthquake obtained using the POT method are not prominent due to the seismogenic fault's E-W orientation, which is nearly perpendicular to the flight direction.
The comparison of range deformation between GPS and InSAR is shown in Table 2.The error statistics in Table 3 were derived based on the residuals between the InSAR and GPS results in Table 2.The difference between the InSAR and GPS deformation results is at the sub-decimeter level, which is consistent with the theoretical accuracy of the POT method [25,26].There is no obvious systematic bias in the InSAR results.Based on the deformation measurements from both the ascending (T014) and descending (T021) orbits in the range and azimuth directions (see Figure 5), the deformation directions on either side of the seismogenic fault for the two earthquakes are opposing.The opposite surface deformations on both sides of the faults indicate that the Turkey double earthquakes were left-lateral strike-slip earthquakes.The range displacement results are significantly better than the azimuth displacement results.This is because the azimuth resolution of the Sentinel-1A data is 14 m and the range resolution is 2.33 m.The range resolution is much higher than the azimuth resolution.The azimuth deformation characteristics of the Elbistan earthquake obtained using the POT method are not prominent due to the seismogenic fault's E-W orientation, which is nearly perpendicular to the flight direction.
The comparison of range deformation between GPS and InSAR is shown in Table 2.The error statistics in Table 3 were derived based on the residuals between the InSAR and GPS results in Table 2.The difference between the InSAR and GPS deformation results is at the sub-decimeter level, which is consistent with the theoretical accuracy of the POT method [25,26].There is no obvious systematic bias in the InSAR results.

Coseismic Slip Distribution
The Pazarcik earthquake mainly damaged the EAFZ (Amanos, Pazarcik and Erkenek segments), and the Elbistan earthquake mainly damaged the CF and DF.Both the earthquakes ruptured to the surface, and the fault activity was mainly controlled by left-lateral strike-slip motion.
The slip of the EAFZ is primarily distributed along the EAFZ (S1) (see Figure 6).The fault is relatively long and almost vertical, and the surface rupture length is about 360 km.The EAFZ (S1) has two obvious slip zones.The first slip zone is situated within 150~220 km of the fault strike, with a depth of 0~15 km and a maximum slip of 7.76 m.The second slip zone is situated within the range of 250~290 km along the fault strike, with a depth of 0~20 km and the maximum slip of 7.0 m.The fault slip is mainly a left-lateral strike-slip with a small amount of thrust slip.Aftershocks in the EAFZ are mainly distributed within the main slip zone, indicating that the stress caused by coseismic rupture may promote the occurrence of aftershocks.However, the aftershock activity at the base of the EAFZ is relatively weak, which may be related to the enhanced friction rate characteristics below the locking depth [38].It is worth noting that although the Pazarcik earthquake occurred in the EAFZ (S3), the fault rupture is weak, with an average slip length of less than 1 m. Figure 7 shows the observed deformation and simulated deformation fields and residual result diagrams from left to right.The simulated deformation field shows good consistency with the observed values, with the average residuals for the ascending and descending orbits being 7.22 cm and 6.35 cm, respectively.This validates the reliability of the coseismic slip model.The presence of residuals may be related to observational noise or the simplification of the fault geometry [28].Additionally, the inconsistent time span between the GPS and SAR images may also contribute to the residual errors.The surface rupture length of the CF and DF caused by the Elbistan earthquake is about 220 km, which is much shorter than that of EAFZ.The slip of the CF is primarily concentrated in the CF (S2), situated within the range of 90~140 km along the fault strike, with a depth of 0~15 km and a maximum slip of 8.2 m.The CF and DF are relatively short, the rupture is predominantly concentrated in the middle of the CF, and the slip is relatively concentrated.
Table 4 compares the fault rupture depth and maximum slip amount inverted by different institutions.The rupture depth obtained in this paper is generally consistent with the results inverted by various institutions, but the maximum slip value is lower.Possible reasons for this include the spatial distribution of the observational data, the constraints, and the smoothing conditions, any or all of which could affect the size of the rupture surface in the coseismic slip distribution model, thereby influencing the magnitude of the maximum slip.This reflects the inherent complexities and variabilities in modeling such geological phenomena, necessitating a nuanced interpretation of the data and results.Figure 7 shows the observed deformation and simulated deformation fields and residual result diagrams from left to right.The simulated deformation field shows good consistency with the observed values, with the average residuals for the ascending and descending orbits being 7.22 cm and 6.35 cm, respectively.This validates reliability of the coseismic slip model.The presence of residuals may be related to observational noise or the simplification of the fault geometry [28].Additionally, the inconsistent time span between the GPS and SAR images may also contribute to the residual errors.

Coseismic Coulomb Stress Disturbance and Regional Seismic Risk Assessment
The Coulomb stress changes were calculated using both the coseismic slip model presented in this paper and the one provided by USGS.Specifically, 81.1% of the aftershocks occurred in the regions of increased Coulomb stress, as determined using the USGSestimated coseismic slip model, while 79.3% occurred in the regions calculated using the coseismic slip model presented in this paper.Additionally, as indicated in Table 4, there is a high degree of agreement between the coseismic slip model presented in this paper and the USGS-estimated coseismic slip model.Consequently, the USGS-estimated coseismic slip model was utilized for the subsequent analysis of Coulomb stress changes, as shown in Figure 8.The CFS change near the epicenter of the Elbistan earthquake was 3.7 bar, which is sufficient to induce a significant fault rupture [41,42], indicating that the onset of the Elbistan earthquake was driven by the Pazarcik earthquake.The aftershocks after the Pazarcik earthquake primarily occurred in regions where the CFS change was positive; therefore, they were largely promoted by the Pazarcik earthquake.

Coseismic Coulomb Stress Disturbance and Regional Seismic Risk Assessment
The Coulomb stress changes were calculated using both the coseismic slip model presented in this paper and the one provided by USGS.Specifically, 81.1% of the aftershocks occurred in the regions of increased Coulomb stress, as determined using the USGS-estimated coseismic slip model, while 79.3% occurred in the regions calculated using the coseismic slip model presented in this paper.Additionally, as indicated in Table 4, there is a high degree of agreement between the coseismic slip model presented in this paper and the USGS-estimated coseismic slip model.Consequently, the USGS-estimated coseismic slip model was utilized for the subsequent analysis of Coulomb stress changes, as shown in Figure 8.The CFS change near the epicenter of the Elbistan earthquake was 3.7 bar, which is sufficient to induce a significant fault rupture [41,42], indicating that the onset of the Elbistan earthquake was driven by the Pazarcik earthquake.The aftershocks after the Pazarcik earthquake primarily occurred in regions where the CFS change was positive; therefore, they were largely promoted by the Pazarcik earthquake.The CFS changes in the faults induced by the Turkey double earthquakes are illustrated in Figure 9.The strike angles and dip angles of the receiving faults were obtained from the active faults database in Turkey [43].The slip angles of the left-lateral, right-lateral, and thrust faults were determined as 0°, 180° and 90°, respectively, simplifying the geometry of the receiving faults.The results show that the Turkey double earthquakes caused significant stress disturbances in surrounding active fault regions.The further The CFS changes in the faults induced by the Turkey double earthquakes are illustrated in Figure 9.The strike angles and dip angles of the receiving faults were obtained from the active faults database in Turkey [43].The slip angles of the left-lateral, right-lateral, and thrust faults were determined as 0 • , 180 • and 90 • , respectively, simplifying the geometry of the receiving faults.The results show that the Turkey double earthquakes caused significant stress disturbances in surrounding active fault regions.The further away from the epicenter the location, the smaller the CFS change was.The cumulative CFS changes in each fault induced by the Turkey double earthquakes are presented in Table A3.The faults that experienced stress accumulation during this event include the Palu and Puturge segments of the EAFZ, the Ecemis segment of the Central Anatolian Fault Zone, the NAFZ, the SATZ, and the Camliyayla, Aadag, Ayvali, Malatya, Heltepe, Pulumur, and Antakya faults.The CFS change caused by this double earthquake had a relatively small impact on faults far from the epicenters, such as the NAFZ, SATZ, Heltepe, and Pulumur, and the impact of the change on them can be ignored.

Conclusions
The surface deformation caused by the 2023 Pazarcik and Elbistan earthquakes was determined using Sentinel-1A images and GPS data in this study.The deformation characteristics, fault slip, and potential impact on surrounding faults of this earthquake event were revealed, and the main conclusions are as follows: (1) The geometry of the ruptured faults in the Turkey double earthquakes is very complex, with ground fault lengths of 360 km and 220 km, respectively.Fault slips occurred at depths of 0-15 km.Both the earthquakes were left-lateral strike-slip earthquakes.The peak sliding value was situated near the surface, at approximately 8.2 m.Along the main fault, three conspicuous main slip zones were observed, two of which extended to the surface.Based on the distribution of aftershocks and the Mw 6.0+ earthquakes in previous years, as reported by UGUS (Figure 1), the M W 6.7 Sivrice earthquake on the Puturge segment in 2020 alleviated the accumulated stress in this segment [44].Therefore, the earthquake risk in this segment is low.
Although the Palu segment released some stress during the 2010 M W 6.1 Karakoçan earthquake, it did not rupture.The number of aftershocks in this segment is small, indicating that the possibility of a large earthquake still exists.
In the Ecemis segment of CAFS, the Camliyayla, and Aadag faults, no medium or large earthquake has occurred in last 100 years.Moreover, stress accumulation may have occurred during the 1998 Adana earthquake [45].The CFS increase along the Ayvali fault during this double earthquake is as high as 4.99 bars.Therefore, the Turkey double earthquakes may have accelerated the rupture of the Ecemis segment of CAFS and the Camliyayla, Aadag, and Ayvali faults.
The Antakya fault experienced an increase in CFS of 8.4 bars during this doubleearthquake event.Therefore, the Turkey double earthquakes collectively facilitated the occurrence of the Mw 6.3 Uzunba g earthquake on 20 February 2023.According to the aftershock distribution in the study area (Figure 2), it is evident that there is a higher concentration of aftershocks near the Malatya and Antakya faults.The Malatya and Antakya faults experienced ruptures during the Elbistan and Uzunba g earthquakes, releasing accumulated stress [2,28]; therefore, the Malatya and Antakya faults have a low seismogenic potential.

Conclusions
The surface deformation caused by the 2023 Pazarcik and Elbistan earthquakes was determined using Sentinel-1A images and GPS data in this study.The deformation characteristics, fault slip, and potential impact on surrounding faults of this earthquake event were revealed, and the main conclusions are as follows: (1) The geometry of the ruptured faults in the Turkey double earthquakes is very complex, with ground fault lengths of 360 km and 220 km, respectively.Fault slips occurred at depths of 0-15 km.Both the earthquakes were left-lateral strike-slip earthquakes.The peak sliding value was situated near the surface, at approximately 8.2 m.Along the main fault, three conspicuous main slip zones were observed, two of which extended to the surface.(3) The Turkey double earthquakes subjected the Ecemis segment of CAFS and the Camliyayla, Aladag, and Ayvali faults to stress loading.The Ayvali fault exhibited a conspicuous CFS-loading condition, indicating a higher risk of future earthquakes, necessitating ongoing monitoring and risk assessment.The Pula fault released some stress during the 2010 Mw 6.1 Karakoçan earthquake.However, there was no significant fault rupture on this fault during this double-earthquake event, and the number of aftershocks in this segment was limited.Consequently, there might be a substantial accumulation of stress in this segment, suggesting the potential for significant earthquakes in the future.

Figure 1 .Figure 1 .
Figure 1.Structural map of the study area.The black beach balls are the focal mechanism solution for MW 6.0+ earthquakes before 2023; The different color lines represent the extent of each fault segment; The red beach balls are the focal mechanism solution of the Turkey double earthquakes and the MW 6.0+ aftershocks; The yellow dots in the figure represent earthquakes (MW > 4.5) recorded byFigure 1. Structural map of the study area.The black beach balls are the focal mechanism solution for M W 6.0+ earthquakes before 2023; The different color lines represent the extent of each fault segment; The red beach balls are the focal mechanism solution of the Turkey double earthquakes and the M W 6.0+ aftershocks; The yellow dots in the figure represent earthquakes (M W > 4.5) recorded by the United States Geological Survey (USGS) solution.(All earthquake information in the figure comes from the USGS.website: https://scihub.copernicus,accessed on 16 June 2023).

Figure 2 .
Figure 2. Ground coverage of SAR images used for the Turkey double earthquakes; Green dots represent aftershocks; Red stars indicate the epicenters of the Turkey double earthquakes.

Figure 2 .
Figure 2. Ground coverage of SAR images used for the Turkey double earthquakes; Green dots represent aftershocks; Red stars indicate the epicenters of the Turkey double earthquakes.

Figure 3 .
Figure 3. Compromise curve between slip model roughness and data fit for Pazarcik earthquake (a) and Elbistan earthquake (b).

Figure 3 .
Figure 3. Compromise curve between slip model roughness and data fit for Pazarcik earthquake (a) and Elbistan earthquake (b).
Remote Sens. 2024, 16, x FOR PEER REVIEW 7 of 18 slip along the east-west direction.The cumulative horizontal deformation caused by the two earthquakes according to 23 GPS stations near the epicenter is shown in Figure 4c.

Figure 4 .
Figure 4. GPS coseismic horizontal deformation caused by the Turkey double earthquakes.The blue and red arrows in (a,b) indicate the Turkey double earthquakes coseismic horizontal displacement; (c) shows the cumulative coseismic horizontal displacement caused by the Turkey double earthquakes and aftershocks; Red stars indicate the epicenters of the Turkey double earthquakes.

Figure 4 .
Figure 4. GPS coseismic horizontal deformation caused by the Turkey double earthquakes.The blue and red arrows in (a,b) indicate the Turkey double earthquakes coseismic horizontal displacement; (c) shows the cumulative coseismic horizontal displacement caused by the Turkey double earthquakes and aftershocks; Red stars indicate the epicenters of the Turkey double earthquakes.

Figure 5 .
Figure 5. Surface deformation of InSAR observations of Pazarcik earthquake and Elbistan earth quake obtained by POT method.(a,b) are the surface deformations in the range and azimuth direc tions of the ascending orbit (T014); (c,d) are the surface deformations in the range and azimuth di rections of the descending orbit (T021).The black triangles are the location of the GPS stations; Th red stars indicates the epicenters of the Turkey double earthquakes.

Figure 5 .
Figure 5. Surface deformation of InSAR observations of Pazarcik earthquake and Elbistan earthquake obtained by POT method.(a,b) are the surface deformations in the range and azimuth directions of the ascending orbit (T014); (c,d) are the surface deformations in the range and azimuth directions of the descending orbit (T021).The black triangles are the location of the GPS stations; The red stars indicates the epicenters of the Turkey double earthquakes.

Figure 6 .
Figure 6.Coseismic slip distribution of the Turkey double earthquakes.(a,b) show the slip distribution map of the CF; (c) shows the slip distribution map of the DF; (d-f) show the coseismic slip distribution of the EAFZ; (g) shows the three-dimensional distribution of coseismic slip distribution; the black dots represent aftershocks; The yellow stars indicates the epicenters of the Turkey double earthquakes..

Figure 6 .
Figure 6.Coseismic slip distribution of the Turkey double earthquakes.(a,b) show the slip distribution map of the CF; (c) shows the slip distribution map of the DF; (d-f) show the coseismic slip distribution of the EAFZ; (g) shows the three-dimensional distribution of coseismic slip distribution; the black dots represent aftershocks; The yellow stars indicates the epicenters of the Turkey double earthquakes.

Figure 6 .
Figure 6.Coseismic slip distribution of the Turkey double earthquakes.(a,b) show the slip distribution map of the CF; (c) shows the slip distribution map of the DF; (d-f) show the coseismic slip distribution of the EAFZ; (g) shows the three-dimensional distribution of coseismic slip distribution; the black dots represent aftershocks; The yellow stars indicates the epicenters of the Turkey double earthquakes..

Figure 7 .
Figure 7. Observed values, simulated values and residuals of coseismic deformation fields.(a,d) are the observation results of the POT method; (b,e) are the simulation results of the coseismic sliding model; (c,f) are the prediction results residuals from simulation results; solid black lines indicate fault segments; yellow stars indicate the epicenters of the Turkey double earthquakes.

Figure 7 .
Figure 7. Observed values, simulated values and residuals of coseismic deformation fields.(a,d) are the observation results of the POT method; (b,e) are the simulation results of the coseismic sliding model; (c,f) are the prediction results residuals from simulation results; solid black lines indicate fault segments; yellow stars indicate the epicenters of the Turkey double earthquakes.

Figure 8 .
Figure 8. CFS change of the Pazarcik earthquake.The receiving fault mechanism is the focal mechanism solution of the Elbistan earthquake.Circles represent aftershocks that preceded the Elbistan earthquake.Yellow stars indicate the epicenters of the Turkey double earthquakes.

Figure 8 .
Figure 8. CFS change of the Pazarcik earthquake.The receiving fault mechanism is the focal mechanism solution of the Elbistan earthquake.Circles represent aftershocks that preceded the Elbistan earthquake.Yellow stars indicate the epicenters of the Turkey double earthquakes.

( 2 )
According to the CFS change, the Pazarcik earthquake caused a CFS change of 3.7 bars near the center of the Elbistan earthquake, which propelled the Elbistan earthquake.The Pazarcik and Elbistan earthquakes increased the CFS change (8.4 bars) in the Antakya fault, which facilitated the occurrence of the MW 6.3 Uzunbağ earthquake on 20 February 2023.

( 2 )
According to the CFS change, the Pazarcik earthquake caused a CFS change of 3.7 bars near the center of the Elbistan earthquake, which propelled the Elbistan earthquake.The Pazarcik and Elbistan earthquakes increased the CFS change (8.4 bars) in the Antakya fault, which facilitated the occurrence of the M W 6.3 Uzunba g earthquake on 20 February 2023.

Table 2 .
Comparison between the range component of coseismic deformation obtained by the POT method and the GPS displacement projected to the range direction.

Table 2 .
Comparison between the range component of coseismic deformation obtained by the PO method and the GPS displacement projected to the range direction.

Table 3 .
The error statistics for InSAR and GPS surface deformations.("Residual Mean" refers to the average residual for the ascending and descending deformation results, respectively, while "Overall Residual Mean" represents the overall average residual for both ascending and descending results.)

Table 4 .
Comparison of fault rupture depth and maximum slip inverted by different institutions.