1.

Introduction

The Haiyuan fault zone, lying in the northeastern margin of the Tibetan plateau, is a major zone of tectonic discontinuity between active blocks and an important seismic belt in western China. It begins south of Xiaokou, Guyuan county, Ningxia Autonomous Region in southeast, and extends northwestward to Xingquanpu, Jingtai county, and Gansu province, belonging to the eastern section of the North Qilian fault system. This zone comprises three en echelon left-slip faults, overall trending in the N280 to 300°W area, with a total length $\sim 240\mathrm{km}$ (Fig. 1). Left slip has continued to characterize its motion since the late Quaternary (describes a period in geologic history) with a cumulative displacement up to 10 to 15 km.¹ The Haiyuan M8.5 earthquake of 1920 took place on the eastern segment of this fault zone, producing a left-lateral strike slip of up to 10.5 m, as well as a 237-km long surface rupture belt.^1–3

Fig. 1

(a) Tectonic setting of the Haiyuan fault zone. White lines are active faults. Thick white line denotes the Haiyuan fault composed of three segments. Red solid dots of different sizes are historical earthquakes of magnitude 6 to 6.9, 7 to 7.9 and greater than 8, respectively. Color background expresses elevation. Two gray boxes are the coverage of ENVISAT/ASAR images from track 2290 to 2018. Two gray dashed small boxes show the study areas of this work. (b) and (c) Average SAR intensity images for the enlarged gray dashed boxes in (a) that are central and eastern segments, respectively.

The analysis of long term and short term active behaviors of the Haiyuan fault zone is very important not only for estimating and predicting seismic hazards, but also for understanding the crustal motion in the northeastern margin of the Tibetan plateau. There have been a lot of studies on the slip rate of the Haiyuan fault zone. Using geological methods and investigations of the features of the terrain around the fault, and using geological dating, the Holocene slip rates of various segments of the Haiyuan fault zone are estimated to be 3 to $10\mathrm{mm}/\text{year}$.^4–7 Most of the recent GPS measurements suggest that this fault zone is in a left-slip state at rates of 3 to $8\mathrm{mm}/\text{year}$,^8–11 though individual work reports slip rates of only 1.1 to $1.7\mathrm{mm}/\text{year}$.¹² The discrepancy in these geodetic estimations may be associated with differences in measurement sites, station spacing, and analysis methods.

In recent years, InSAR technology has been applied to the measurement of crustal motion around faults and especially to the observation of co-seismic deformation fields. Although in theory, this method can detect displacement of mm orders of magnitude on the surface, interseismic crustal motion or fault slip is so slow and tiny that its signal can be blanketed by errors from decoherence factors such as satellite orbits atmospheric delay, and terrain conditions. To solve this problem, some improved InSAR techniques have been proposed, such as permanent or persistent scatterers InSAR (PSI), small baseline subsets (SBAS), and interferograms stacking. All these methods have been applied to detection of surface deformation.^13–42 Of them, the PSI uses highly coherent point targets, called “PS,” on the time series synthetic aperture radar (SAR) image with a high temporal stability of the backscattered signal, such as houses, bridges, big stones, and mountain crests, as objects for phase analysis, so that noncoherent image elements with low signal-to-noise ratio are avoided.^13,14 Once a PS data set is determined based on a criterion and algorithm, the analysis of the the phase time series for each PS allows one to detect near-field crustal motion of a fault.

In this work, we attempt to use the PSI method to conduct an experimental study of crustal deformation measurements on the Haiyuan fault zone. Our purpose is to examine whether this approach can reduce decorrelation of phase in repeat SAR images caused by some natural conditions in western China, and help better detect the near-field motion of active faults like the Haiyuan fault zone.

2.

InSAR Data and Processing

2.1.

Study Area and InSAR Data

The Haiyuan fault zone lies in a loess area of western China laced with numerous gullies, with an average elevation of 1700 to 2000 m and a maximum elevation of about 2500 m [Fig. 1(a)]. A dry climate and rare vegetation characterize this region, most of which is intact mountains. The exposed ridges and slopes are natural point targets for PS selection, so this area is suitable for measurements of fault activity using PSI. At present, due to limitations of the technology, the PSI can be used to detect crustal motion of small-sized areas, so we choose two small experimental areas for this work. One crosses the central section (104.75°E, 36.78°N) of the Haiyuan fault zone, as shown in Fig. 1(b). On the south and north sides of this area are mountainous lands with steep relief and little vegetation. Its central portion, where the Haiyuan fault zone runs, is relatively flat with towns, villages, and croplands. The Yellow River and flood land are present in the southwest corner. The northern side of the eastern section (105.60°E, 36.50°N) with relatively flat terrain, is where Haiyuan county is located as shown in Fig. 1(c). Here, on the south side of the Haiyuan fault zone are steep mountains, while on the north side is flat land with a town. Our tests show that the experimental area covering both sides of the fault zone with large terrain differences can yield bad results, so we have to select an area on the flat northern side.

The SAR data used in this work are ENVISAT/ASAR from the European Space Agency, with image mode, IS2 incidence mode, incidence 19 to 26 deg with a 23 deg average, C band and VV polarization, in a descending orbit from track 2290 to 2018. We have 17 scenes from 20040119 to 20100308 in track 2290, and 21 scenes from 20030813 to 20091209 in track 2018. Considering the spatial perpendicular baseline, temporal baseline, Doppler center frequency difference, and atmospheric errors, we selected one master image and form interferometric pairs (including one automatic pair). The ASAR data and interferometric parameters are shown in Tables 1 and 2 and Fig. 2, respectively. For most pairs, the perpendicular baseline is $<400\mathrm{m}$, which is helpful for maintaining good coherence in the mountainous area. Here, it should be mentioned that although there should be no spatial decorrelation occurring on the point target due to its dominated echo in the located pixel, in fact its phase is still affected by geometric decorrelation to some extent because of the limited number real point targets.

Table 1

Data and interferometric parameters of image pairs in the central segment of the Haiyuan fault.

No.	Master	Slave	Perpendicular baseline (m)	Time interval (d)	Registration error in range (pixel)	Registration error in azimuth (pixel)
1	20050523	20040119	$-331.8467$	$-490$	0.1733	0.1577
2	20050523	20040607	471.8301	$-350$	0.1475	0.1292
3	20050523	20041025	542.0414	$-210$	0.1824	0.1727
4	20050523	20050523	0.0000	0	0.1843	0.1590
5	20050523	20070806	$-59.9810$	805	0.0775	0.0998
6	20050523	20080303	$-346.9153$	1015	0.0278	0.0838
7	20050523	20080512	$-18.7991$	1085	0.1458	0.1237
8	20050523	20080825	140.3449	1190	0.0874	0.0901
9	20050523	20080929	$-220.2607$	1225	0.0	0.0
10	20050523	20081103	103.4028	1260	0.1488	0.1467
11	20050523	20090601	149.6371	1470	0.0674	0.0975
12	20050523	20090706	$-4.0218$	1505	0.0792	0.2458
13	20050523	20090914	420.5076	1575	0.0697	0.1457
14	20050523	20091123	307.9203	1645	0.0598	0.1697
15	20050523	20091228	$-244.5750$	1680	0.1236	0.1569
16	20050523	20100201	375.2442	1715	0.0698	0.1597
17	20050523	20100308	$-104.3174$	1750	0.1567	0.1388

Table 2

Data and interferometric parameters of image pairs in the eastern segment of the Haiyuan fault.

No.	Master	Slave	Perpendicular baseline (m)	Time interval (d)	Registration error in range (pixel)	Registration error in azimuth (pixel)
1	20041215	20030813	$-698.0954$	$-490$	0.1633	0.1577
2	20041215	20030917	511.2700	$-455$	0.1275	0.1192
3	20041215	20031022	869.0337	$-420$	0.1924	0.1927
4	20041215	20031231	876.1930	$-350$	0.1943	0.1390
5	20041215	20040519	185.0897	$-210$	0.0615	0.0928
6	20041215	20040623	$-53.5531$	$-175$	0.0198	0.0738
7	20041215	20040901	476.9821	$-105$	0.1098	0.0897
8	20041215	20041006	279.4629	$-70$	0.0794	0.0991
9	20041215	20041215	0.0000	0	0.0	0.0
10	20041215	20050330	$-660.0841$	105	0.1406	0.1443
11	20041215	20050921	203.8936	280	0.0674	0.0875
12	20041215	20080213	$-261.6804$	1155	0.0792	0.2058
13	20041215	20080319	110.1835	1190	0.0683	0.1446
14	20041215	20080528	68.9444	1260	0.0541	0.1625
15	20041215	20080806	302.2617	1330	0.1062	0.1481
16	20041215	20081015	151.2705	1400	0.0600	0.1557
17	20041215	20090617	440.9376	1645	0.1502	0.1378
18	20041215	20090722	13.1335	1680	0.0391	0.0204
19	20041215	20090930	$-51.5754$	1750	0.0601	0.0509
20	20041215	20091104	165.5058	1785	0.0709	0.1278
21	20041215	20091209	$-257.8017$	1820	0.0870	0.1579

Fig. 2

Perpendicular baselines and time intervals of interferometric pairs. (a) and (b) are for the central segment and eastern segment of the Haiyuan fault zone, respectively.

2.2.

Data Processing

We use the interferometric point target analysis (IPTA) module of the software Gamma to process the ASAR data. The whole flowchart is shown in Fig. 3. PS candidates (PSCs) are determined based on the co-registered single look complex (SLC) images. For each PSC, the SLC values are extracted and written to a point data stack. Next, the differential interferograms of the point data stack are calculated. From the differential phases of good PSCs, deformation, height correction, and atmospheric delay components are estimated using an iterative procedure. That is to say, the interferograms are only generated and interpreted for the selected points in IPTA. For efficiency and data storage reasons, vector format data structures are used instead of the raster data format used in conventional interferometry.⁴³

Fig. 3

Flowchart of data processing including geocoding procedure.

Before interferometric processing, we make a radiometric calibration of each SAR image used in this study acquired at different times, so that their intensity values are comparable, and good to the followed PSCs selection [convenient for selection of permanent scatterer candidates (PSCs) from many SAR images]. We also make perpendicular baseline corrections to the interferograms based on the accurate orbit data from Delft University. The initial topographic phase component is removed using the 3-arc-sec-SRTM DEM. Considering the mountainous terrain of the study area, we adopted a new co-registration algorithm based on a lookup table, which can greatly improve co-registration accuracy and increase coherence of the image pairs in the mountain area. As shown in Tables 1 and 2, most image pairs have a co-registration accuracy $<0.2\text{pixel}$, which is the precondition for success of the chosen processing methods. The detailed processing procedures are described below.

2.2.1.

Generation of PSCs

In processing the time-series InSAR analysis, one usually adopts amplitude stability and coherence stability or either of them (i.e., correlation) as pixel selection criteria. The choice of the selection criterion depends on the individual application and the point target processing method to be used. In this study, we use Gamma software programs to select the PSCs automatically.⁴³ Both amplitude stability (low temporal variability) and spectral phase diversity are utilized as selection criteria of PSCs. We first identify point target candidates (PSC1) based on the lower temporal variability of the backscattering (ratio of mean to standard deviation) in the average amplitude image of all co-registered SLCs, which is better for large SLC data stacks. Then we detect point target candidates (PSC2) based on low spectral phase diversity in a single co-registered SLC, which is better for small SLC data stacks. Finally, the two point target candidates are merged into a single one in which only those with common positions (column/pixel) are kept (Fig. 3). In the followed regression analysis, the quality of these PSCs is evaluated and some noisy PSCs are dropped.

In the experimental area of the central fault segment, we select 69143 PSCs from point targets in bare mountainous terrain. These PSCs are distributed nonuniformly, being relatively dense in the south and north mountainous sections, but sparse in the middle and southwest flat parts. Field investigations show sandstone, mudstone, and shale exposures in the northern and southern sections of the study area with steep terrains and linear ridges result in strong and stable backscattering, providing a lot of point targets. In the middle of the experimental area, some villages and towns are present in a flat terrain covered with farm fields and vegetation, thus few PS points can be identified. In the other experimental area, located at the eastern fault segment, we select 1204 PSCs that are concentrated around Haiyuan county. This implies that the distribution and density of PS points depend on ground surface conditions. As the quality of these PSCs would be evaluated and those of poor quality would be eliminated in the subsequent regressive analysis, the threshold for PSC selection set in this stage is not very high.

2.2.2.

Geocoding of SAR data

For InSAR data processing, orthorectification means geocoding of SAR data. This is the coordinate transformation between the coordinates of range-Doppler coordinates of the SAR and orthonormal map coordinates. This type of geocoding takes the image in the radar geometry (range azimuth) and re-samples the data into a map projection, such as UTM with an associated datum. Depending on the ground surface considered for geocoding, it can be distinguished between ellipsoid corrected and terrain corrected geocoding. The former does not require a digital elevation model (DEM) and solely considers a constant height, while the latter needs the DEM and takes into account the surface topography which definitely has an effect on the look vector.⁴³

Due to the slope terrain of this study area, we use the terrain corrected geocoding method. The geocoding procedure includes three steps: First, it builds an initial geocoding lookup table between the SAR and map geometry based on the orbital data and the parameters of the master image to be geocoded as well as that of the DEM. Second, it refines the lookup table through simulating the SAR intensity image based on the DEM and the SAR imaging geometry, and calculates the registration offsets between the simulated SAR image and the true SAR image. Third, it implements the transform from radar geometry into map geometry based on the refined geocoding lookup table, with the advantage of just a single re-sampling and interpolation (Fig. 3).

2.2.3.

Calculation of initial differential phases

We extracted SLC data of PSCs and generated initial interferograms by interferometric processing focused on only the point data stack, as well as extracted elevation values of PSCs from 3-arc-sec SRTM DEM and calculated the simulated topographic phases (initial phase model) of the PSCs. Then by subtracting the initial phase model from the point interferograms, the initial differential phases are determined (Fig. 4). Figure 4 shows 17 differential interferograms of PSCs with different perpendicular baselines in the central segment of the Haiyuan fault zone. It is clear that for the interferometric pairs with a short perpendicular baseline $<300\mathrm{m}$, the initial differential interferograms are relatively smooth, implying good coherence. For those pairs with relatively long perpendicular baselines $>300\mathrm{m}$, however, they are rough in variation and disordered in color, indicative of high noise level. It suggests that the perpendicular baseline is still a very important influence factor for point interferometric processing, so we adopt an appropriate strategy of regression analysis according to perpendicular baselines of all point interferograms.

Fig. 4

Initial differential interferograms of PSCs of 17 image pairs in the central segment of the Haiyuan fault zone. The interferograms with perpendicular baselines longer than 300 m are very noisy. Blue asterisk: reference point.

2.2.4.

Regression analysis of point differential phases

The initial differential phases of PSCs may contain various phase error components, such as orbital error, topographic residual, atmospheric delay, and other noise. These error phases can be estimated and removed through regression analyses iteratively based on a two-dimensional linear phase model, which indicates the topographic phase is linearly dependent on the perpendicular baseline component of the interferometric pairs, and the crustal motion rate is linearly dependent on the time interval between the image pairs. Thus the related slopes of regressions correspond to relative terrain height corrections and relative linear deformation rates. During regression analyses, a Delaunay triangle grid is constructed and phase unwrapping of differential phases of PSCs is completed by the minimum cost flow (MCF) algorithm. The phase standard deviation and residual phases derived from the regression are used as a measure for the quality of phase unwrapping and regressive fitting. Low phase standard deviation and smooth residual phases are indicative of correct unwrapping, whereas large phase standard deviation and uneven residual phases exhibiting block-looking variation mean wrong unwrapping and iterative calculation needs to continue.⁴³

As the original differential interferograms have distinct coherences, we adopt the following procedures for regression analyses. First, regression processing and phase unwrapping are done on those interferometric pairs (layers) with a perpendicular baseline $<300\mathrm{m}$ in order to search the layers of high quality. Using these good layers, the DEM error model and linear deformation rate model for some PS points are reliably established. Second, more layers and points are unwrapped based on these models, while those points of poor quality are removed. Third, perpendicular baseline correction is performed by the least square algorithm based on stable points with motion rates $<2\mathrm{mm}/\text{year}$. The regression processing is conducted again using the optimized baseline, and the DEM error model and linear deformation rate model are renewed. Finally, the atmospheric delay phase is estimated and eliminated by spatial filtering to residual phases and a new linear deformation model after atmospheric correction as well as the uncertainty of parameter estimates and phase standard deviation is generated again. These resultant images are shown in Figs. 5Fig. 6Fig. 7Fig. 8Fig. 9Fig. 10Fig. 11–12. It should be mentioned that the reliability of estimation of the linear motion rate is related to the quality of PS points, while the quality of PS points is associated with their density. A trade-off is needed between the quality and the density of PS points.

Fig. 5

Motion rate maps of the study area on the central segment of the Haiyuan fault zone using different PS point amounts. The relative motion rate across the fault is about $5\mathrm{mm}/\text{year}$. Red-dashed line denotes the location of the Haiyuan fault zone. Blue asterisk: reference point. (a) $\mathrm{PS}=67665$, (b) $\mathrm{PS}=59584$, (c) $\mathrm{PS}=39993$ and (d) $\mathrm{PS}=12333$.

Fig. 6

Phase standard deviations corresponding to different PS amounts of the central segment of the Haiyuan fault zone. With decreasing PS points, phase standard deviations decline. White-dashed line is the trace of the Haiyuan fault zone. Blue asterisk: reference point. (a) $\mathrm{PS}=67665$, (b) $\mathrm{PS}=59584$, (c) $\mathrm{PS}=39993$ and (d) $\mathrm{PS}=12333$.

Fig. 7

PS point deformation rates versus topographic elevations in the case of PS point amount 59584. (a) PS point deformation rate map. (b) PS point elevations of the study area. (c) Estimate uncertainties of PS point deformation rates. (d) Estimate uncertainties of PS point elevation correction. Red-dashed line donates the Haiyuan fault. Black-dashed line is the elevation 1700 m. Black-dotted double lines in (a) denote profiles of deformation rates. In the mountainous area with elevation above 1700 m, the reliability of rate estimate drops. Blue asterisk: reference point.

Fig. 8

Deformation rate profiles across the fault. Locations of profiles are shown in Fig. 7(a). Profiles extend from north to south for about 22 km. Gray thick dashed lines denote the Haiyuan fault. Red-dashed lines show the fit curves. (a) Profile $\mathrm{A}-{\mathrm{A}}^{\prime }$ and (b) profile $\mathrm{B}-{\mathrm{B}}^{\prime }$.

Fig. 9

Cumulative displacement-time curves of four PS points north of the Haiyuan fault. Abscissas are time in days. Ordinates are displacement in meter.

Fig. 10

Time series of cumulative displacement images of PS points for all interferograms respect to the reference scene 20050523 in the central segment study area of the Haiyuan fault. Blue asterisk: reference point.

Fig. 11

(a) Deformation rate map of PS points in the eastern segment study area of the Haiyuan fault. (b) Corresponding phase standard deviations. This study area is in the north side of the Haiyuan fault. Asterisk: reference point.

Fig. 12

Time series of cumulative displacement images of PS points for all interferograms with respect to the reference scene 20041215 in the eastern segment study area of the Haiyuan fault. Asterisk: reference point.

3.

Results and Analysis

Figure 5 shows the resultant linear deformation rates in four conditions with different number of PS points in the central segment of the Haiyuan fault zone, which are 67665, 59584, 39993, and 12333, respectively [Figs. 5(a)–5(d)]. On the whole, the motion rates corresponding to the four PS point amounts are similar, each rate map exhibits even variation with smooth colors, devoid of local patches, but with an apparent motion rate gradient across the Haiyuan fault zone, which just agrees with the feature of fault motion. On the other hand, some colored spots with low quality can be found in the condition of more PS points, such as Figs. 5(a), especially in the northern mountain land of the study area. If the number of PS points is reduced by setting a lower threshold of phase standard deviation in the regression analysis, these low-quality points could decrease considerably or even disappear [Figs. 5(b)–5(d)]. It implies that both the quality of PS points and reliability of motion rate estimates are enhanced with a decrease of PS points. However, the observation from sparse and few PS points would lose the advantage of revealing the spatial variation trend of the motion rate around the fault zone and motion differences on either sides. Thus, a trade-off between the number of PS points and the quality of deformation rate estimates is required, as in Fig. 5(b), which should be the optimal result in a relative sense.

Figure 6 displays phase standard deviations corresponding to PS point amounts 67665, 59584, 39993, and 12333, respectively, which measure the quality of regression analyses and reliability of motion rate estimates. It shows that with a decrease of PS point number, the standard deviations decrease. For example, in Figs. 6(b)–6(d), these deviations are as small as $\le 1.0\mathrm{rad}$ in general, while in Fig. 6(a) with relatively more PS points, the deviation increases in the north and south parts and $\ge 1.0\mathrm{rad}$, especially in the north mountains, which suggests that there are some points of low quality there. The motion rates of PS points are also more heterogeneous and dispersed, indicative of poor reliability [Fig. 5(a)]. The reason for this may be radar shadows and overlapping in the mountain areas which lead to a decline of point quality. Therefore, in the authors’ opinion, it would be better to set different and proper thresholds of criteria for PS point selection according to specific conditions of topography and vegetation coverage, and perform regression analyses separately in distinct regions.

Figure 7 shows relationship between PS point deformation rates, estimate uncertainty and terrain for the case of PS point amount 59584 with standard deviations $<1.0\mathrm{rad}$, representing a reliable result. The deformation rates [Fig. 7(a)] are relatively smooth and exhibit a clear tendency of gradient differences across the fault in despite of the nonuniform distribution of PS points. Specifically, south of the fault is dominated by red-yellow points, whereas it transits into yellow-light blue north of the fault. From such a contrast of colors, the relative deformation rate in the line of sight (LOS) is estimated to be about $5\mathrm{mm}/\text{year}$ of left slip.

The estimate uncertainties of deformation rates of most PS points are below $0.6\mathrm{mm}/\text{year}$ [Fig. 7(c)], except for a few points in the north and other places where the uncertainties are $>0.6\mathrm{mm}/\text{year}$. The elevation correction uncertainties are of the same feature, mostly $<2\mathrm{m}$ except for individual points in the north [Fig. 7(d)]. The elevation of PS points in nearly the whole study area is below 1700 m, as shown as the light-blue color region [Fig. 7(b)], while altitudes are $>1700\mathrm{m}$ in the north mountainous portion shown by yellow and red colors. Comparison of Figs. 7(b), 7(c), and 7(d) shows that the uncertainties of both motion rate and elevation correction estimates are related to real elevations. The large uncertainties appear in mountainous portions with elevations over 1700 m, while those values are relatively small in lower terrains. It means that the accuracy of motion rate estimates using PSI is relatively poor in steep mountains.

From the profile $\mathrm{A}–{\mathrm{A}}^{\prime }$ [Figs. 7(a) and 8(a)], a considerable gradient of deformation rates is present on either side of the fault zone. In the north, the motion rates are $-6$ to $-7\mathrm{mm}/\text{year}$, while in the south they become $-1$ to $-2\mathrm{mm}/\text{year}$, implying a relative change of about $5\mathrm{mm}/\text{year}$ across the fault zone. On the profile $\mathrm{B}–{\mathrm{B}}^{\prime }$ [Figs. 7(a) and 8(b)], motion rate differences are also present on either side of the fault. A relative rate of $5\mathrm{mm}/\text{year}$ appears in the north of the Haiyuan fault zone, and it is around 0 in the south, also yielding a relative change of $5\mathrm{mm}/\text{year}$ across the fault zone (same as on the profile $\mathrm{A}–{\mathrm{A}}^{\prime }$).

In order to describe the motion history of the study area, we select four PS points on the north side of the central Haiyuan fault zone and plot their cumulative displacement–time curves (Fig. 9). They exhibit linear changes with time, and a largely consistent cumulative displacement of about 4 cm on average with a negative sense, i.e., increase of the distance between the satellite and the Earth. Assuming the ground moves along with the fault zone only horizontally, these changes mean westward motion of PS points on the north side in the mode of the right looking and descending orbit, which is consistent with the left slip of the Haiyuan fault zone.

Cumulative displacement images of all 17 interferometric pairs are shown in order of time (Fig. 10). Each of them reveals displacements of PS points over the interval between two acquisition times of one image and the reference image. The colors of these serial displacement images change slowly and stably with time, implying a linear variation of PS point motion. They should reflect a long-term slow interseismic motion of the crust rather than local changes caused by the atmosphere or other factors. It is noted that four PS point interferograms before 2008 (20040119, 20040607, 20041025, and 20070806) have a relatively conspicuous color contrast with respect to those after 2008, perhaps indicating a change of crustal motion around 2008. The other 12 interferograms after 2008 (20080303 to 20100201) have a small color contrast, probably indicating a stable and slow motion process. It is likely that such a phenomenon is associated with the Wenchuan great earthquake of 2008, or implies an accelerated motion on the Haiyuan fault zone even before that event. It is reported that an anomalous deformation was observed at several sites on the Qilian Shan and Haiyuan faults $>600\mathrm{km}$ to the north of Wenchuan,^44,45 consistent with the above results. Unfortunately, there is little ENVISAT/ASAR data from 2006 and 2007 available for the study area, limiting further exploration of the crustal motion anomalies prior to the 2008 Wenchuan great earthquake.

The resultant motion rate maps of PS points on the eastern segment of the Haiyuan fault zone are shown in Figs. 11 and 12, respectively, where the PS points are mainly concentrated in Haiyuan county in the center of the experimental area. In this test area, the fault motion rates revealed by most PS points are about 5 to $6\mathrm{mm}/\text{year}$ [Fig. 11(a)], and the corresponding phase standard deviations are mostly $<0.7\mathrm{rad}$ [Fig. 11(b)]. A time series of cumulative displacement images of PS points (Fig. 12), similar to those in the central segment of the study area, display a gradually stable variation trend, and a total displacement up to 4 cm from August 2003 to December 2009. However, the cumulative displacement images before 2008 show more visible color differences than those after 2008, which is a very interesting phenomenon. We shall study it based on more data in the future.

4.

Discussion

The radar satellite measures range variation in the LOS between two flights over the same target. In theory, such measurements can be converted into strike slip and dip slip of a fault by a geometric relationship (Fig. 13). For a fault like the Haiyuan fault zone, we assume that Dstrike is a strike slip of the fault (southeastward is positive), Ddip is the dip slip (upward is positive), Dlos is displacement in LOS (toward satellite is positive), $\theta$ is the incidence angle of the radar, $\alpha$ is the orbit azimuth of the satellite (clockwise is positive), and $\beta$ is the strike of the fault (including the angle to the north). Then, we have the following relationship:

Fig. 13

Sketch map for conversion of LOS deformation to strike slip along with the Haiyuan fault with ENVISAT/ASAR data on descending orbit. (a) Deformation project from ground range to slant range. (b) Deformation project from fault strike to ground range. $\theta$: incidence angle of radar pulse, $\alpha$: azimuth angle of radar flight orbit, and $\beta$: fault strike.

In this work, we use the images ENVISAT/ASAR of a descending pass with an average incidence of 23 deg and a satellite orbit azimuth 192 deg. The Haiyuan fault zone tends in 280 to 300 deg, so we take 290 deg. If we assume that this fault is of a pure strike-slip motion, the calculation by the relationship above shows that a motion rate of $5\mathrm{mm}/\text{year}$ in LOS corresponds to $12.9\mathrm{mm}/\text{year}$ along the fault. If we assume that the fault also has a $1\mathrm{mm}/\text{year}$ of dip slip in addition to its strike-slip motion, a strike-slip rate of $10.5\mathrm{mm}/\text{year}$ can be derived from the $5\mathrm{mm}/\text{year}$ in LOS.

GPS measurements (3 to $8\mathrm{mm}/\text{year}$) are the average values on a large scale. In this work, the PSI analysis attempts to reveal motions in a small area, likely with a localized feature. Besides, the results from this approach might include the dip slip component of the fault, while GPS is only a horizontal motion. Thus, we think the relatively high value in this work ($10.51\mathrm{mm}/\text{year}$ with dip slip or $12.9\mathrm{mm}/\text{year}$ pure strike slip, see Table 3) is acceptable, though it does have some error. Of course, this issue should be studied further.

Table 3

Strike-slip rates of the Haiyuan fault in experimental areas converted from LOS motion rates.

The PSI analysis reveals a gradient contrast of motion rates across the central segment of the Haiyuan fault zone. This is not what we expected if the upper crust portion is presumably locked above the seismogenic layer. Nevertheless, such a phenomenon has been observed on many active faults in the world.^20,23 Using the interferogram stacking method that converts average phase change into velocity, Cavalié et al.²² found a zone of high velocity gradient across the western Haiyuan fault zone which is a few kilometers wide. They suggest a screw dislocation model in an elastic half-space with locking at depth or an alternative model with shallow creep and deep aseismic slip to explain the InSAR observations. The current PSI cannot resolve the discrimination between these two models by itself. We expect to combine PSI, GPS, and seismic observations to yield a better interpretation on the motion behavior of the Haiyuan fault in the future.

This work shows that the PS point cumulative displacements of all interferograms in the two experimental areas exhibit conspicuous changes before and after 2008, likely indicating the effect of the 2008 Wenchuan great earthquake on the Haiyuan fault zone. Because only a few SAR scenes before the shock are available, with data from 2006 and 2007 absent and only three scenes available for 2005 in two tracks, our effort to study the variation of motion rates on the Haiyuan fault zone before and after the Wenchuan event is limited. We hope to acquire more relevant data to fill this gap in our following work.

This work demonstrates that the reliability of the processing result of PSI is associated with the number and quality of PS points that are related to surface coverage of the study area. PSI can work well in mountainous regions with moderate vegetation where PS points are usually dense and stable. But in such settings, the distribution of PS points are sometimes not uniform due to the concentration of partial PS points on ridges and slopes that are facing the incidence of radar wave fronts. In the regions with heavy vegetation, PS points are sparse and unstable, thus the usage of PSI is limited. A possible approach to solve this problem is to combine PS and CR (corner reflector) techniques, i.e., to increase observational point targets by the deployment of corner reflectors in the places lacking PS points.

This work was conducted on a small area of the central and eastern Haiyuan fault zone, thus cannot describe the variation of near-field motion rates along the entire fault zone as well as a wider range across the fault zone. Meanwhile, the SAR data used a span only about 6 year. With the increasing spatial scale of a study area, the phase errors of PS points from baselines, topography, and atmospheric delay would become larger, probably leading to failure of phase unwrapping. This is an issue in PSI that remains to be solved. One possible approach is to perform processing on many small areas separately and then stitch them together so that these small areas are continuous and comparable in motion rates. It will be a focused topic of our future work. On the other hand, we will use other time-series InSAR methods such as SBAS, and stacking to research larger regions.

5.

Conclusions

In this work, we attempt to detect near-field crustal motion of the Haiyuan fault zone based on 38 ENVISAT/ASAR images of two neighbor descending tracks through a PSI analysis. The 3-arc-sec SRTM DEM and IPTA module of the software Gamma are used in the data processing. The results show that the motion rates of PS points in the two experimental areas of the fault zone are all smooth and consistent in spatial distribution, exhibiting a feature of steady and linear variance in time. Meanwhile, a conspicuous contrast of motion rates is present on opposite sides of the Haiyuan fault zone, and a relative motion rate in LOS about of about $5\mathrm{mm}/\mathrm{a}$ is derived, which can be converted into the left-lateral strike slip of $13\mathrm{mm}/\text{year}$ of the fault under the assumption of pure strike-slip faulting. This value is larger than 3 to 10 and 3 to $8\mathrm{mm}/\text{year}$ determined by geological and GPS observations, respectively. A possible explanation for this discrepancy is that the estimate from geology is a long-term average (e.g., since Holocene) of crustal motion or fault slip, and estimates from GPS are averaged on a large scale, while the result from PSI represents a near-field motion picture of a fault with a localized character. Besides, the motion rates from PSI analysis likely contain dip-slip displacement of the fault in addition to the horizontal slip. Of course, this issue requires further studies with more data and modeling.