Identification and Analysis of Deformation Areas in the Construction Stage of Pumped Storage Power Station Using GB-InSAR Technology

To meet the needs of the rapid development of new energy sources, China is currently accelerating the construction of pumped storage power stations (PSPS). However, the complex geological environment and large-scale construction can easily cause deformation areas, such as collapse and sliding, on the surface of PSPS. This article adopts ground-based interferometric synthetic aperture radar (GB-InSAR) technology to monitor surface deformation using 206 ground-based radar (GBR) data collected from 2022/07/31 10:24 to 2022/08/06 14:57 with a PSPS in the construction stage in western China. In the selection of permanent scatterers (PSs) for GBR data processing, statistical homogeneous pixels (SHPs) information is considered in this article. Based on the fast SHP selection (FaSHPS) algorithm, we introduce intensity deviations (${D}_I$) to exclude pixels with unstable intensity performance in the temporal dimension and propose the ${D}_I + FaSHPS$ algorithm, which is experimentally proven to have higher SHPs selection accuracy and computational efficiency. By combining with the traditional PSs selection algorithm, the distribution density of monitoring points is effectively improved. In addition, based on the deformation results and on-site survey, this article analyzes and discusses the identified deformation areas in detail from the aspects of natural geography and man-made activities. The results show that the surface is relatively stable overall, except for the effects of man-made activities. This study has demonstrated the remarkable reliability and superiority of the GB-InSAR technology in identifying and analyzing deformation areas in the construction stage of PSPS.


I. INTRODUCTION
A PUMPED storage power station (PSPS) uses electrical energy from the low power load to pump water to the upper reservoir and then releases water to the lower reservoir to generate electricity during peak power load periods [1], [2]. Due to its role in grid peak regulation, voltage regulation, energy storage, and power stability control, PSPS have become an important type of hydropower projects around the world [3]. With some of the world's richest water resources, PSPS construction in China has been growing rapidly since 1968. Entering a new era, China has solemnly announced to the world that it aims to achieve "peak carbon emission" by 2030 and "carbon neutrality" by 2060. As a renewable and clean energy, the construction scale of PSPS will undoubtedly continue to increase [4]. The monitoring object is the first PSPS project in western China, which started construction in 2016. After the project is put into operation, it will be able to achieve joint and coordinated operation with thermal power, wind power, and solar power generation in Northwest China, playing an important role in achieving the "double carbon" goal.
Traditional means of deformation monitoring (e.g., level, total station, navigation, and positioning equipment) are limited to a small number of ground points and require manual operation or on-site deployment of instruments, with strict requirements for the terrain environment [5], [6]. 3-D laser scanning technology captures dense, irregular data, requiring larger data storage space and longer data processing times [7]. The main advantage of space-borne interferometric synthetic aperture radar (InSAR) systems is its wide spatial coverage, which is suitable for surface deformation monitoring in large study areas, such as cities, large slopes, and mining areas [8], [9], [10], [11]. However, the limitations of fixed orbit, long revisit time, low data resolution, and insensitivity to deformation in the north-south direction have hindered the wide application of space-borne InSAR systems in site-specific deformation monitoring. Compared with space-borne InSAR systems, groundbased (GB) InSAR systems can be deployed near the monitoring target and can also adjust the imaging period and observation direction, with higher data resolution, and their flexibility is a significant advantage in deformation monitoring [12], [13]. GB-InSAR systems have been developed over two decades, and they mainly consist of the following two types: synthetic This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ aperture radar (e.g., IBIS [14], FastGBSAR [15]) and real aperture radar (e.g., gamma portable radar interferometer, GPRI and GPRI-II [16]). In addition, GB-InSAR technology can measure the position changes of various natural objects [17], [18], [19], [20] such as slopes, mining areas, mountains and glaciers, as well as large buildings [21], [22], [23], such as bridges, dams, tall buildings, and hydropower stations with submillimeter accuracy.
Currently, the most commonly used data processing method for GB-InSAR technology is multi-temporal InSAR (MT-InSAR) [24], which mainly includes persistent scatterer InSAR (PS-InSAR) [25] and small baseline subset [26]. MT-InSAR can more accurately obtain the surface deformation rate and time series by modeling and analyzing the selected reflection phase stabilization points [27], [28], [29]. Coherence point selection plays a very important role in MT-InSAR, and its quality and density directly affect the extraction accuracy of deformation parameters [30], [31]. Generally, amplitude deviation index and coherence threshold are used as the selection criteria to select high-quality permanent scatterers (PSs) [32]. PSs return significantly more energy than other scatterers within a resolution pixel, thus reducing the decorrelation. Subsequently, the SqueeSAR method, which aims to improve the distribution density of monitoring points, has been proposed to accompany the increased demands of surface deformation monitoring [33]. This second generation PS-InSAR method uses the temporal amplitude values of two pixels within a window to construct a probabilistic distribution to identify statistically homogeneous pixels (SHPs). With the continuous emergence of algorithms for SHPs selection, the accuracy of its selection has been greatly improved [34], [35], [36]. In addition, as an important way to select PSs, SHPs effectively enriches the content of GB-InSAR technology [37].
In reality, the natural geological and geomorphological conditions of PSPS projects are often bad, with high and steep banks and complex topography and dense vegetation [38]. It is difficult to select a sufficient number of monitoring points for surface deformation analysis in such areas. In addition, geological hazards, such as slope failures and surface landslides, occasionally occur during the construction process, requiring timely, flexible safety monitoring of hazard areas. Furthermore, the identification and analysis of deformation areas during the construction stage of a PSPS require a large scanning angle (>180°) to obtain the wide area surface deformation information. Unfortunately, there is not much research related to GB-InSAR technology for such monitoring situations. In the face of the above predicament, this article has carried on the exploration and research full of practical significance. The rest of this article is organized as follows. First, this article presents an overview of the study area as well as the experimental design and data in Section II. Then, the article elaborates the MT-InSAR method in detail and focuses on the PSs selection method based on SHPs information in Section III. Subsequently, according to the PSs selection and deformation monitoring results, this article presents a detailed analysis and discussion in Section IV. Finally, Section V concludes this article..

A. Study Area
The PSPS studied in this article is located in the western region of China, approximately 128 km from downtown Xi'an in Shaanxi Province [39]. The study area is located at the eastern foot of the Taibai Mountains, the main peak of the Eastern Qinling Mountains, and is a medium-high mountainous area, adjacent to the Weihe Plain in the north, with a north-south tilted fault block structure. The study area is characterized by high mountains and deep ravines, dense vegetation, and complex topography. When fully completed, the PSPS will mainly serve the Shaanxi Power Grid, undertaking tasks such as peak regulation, valley filling, frequency modulation, and emergency backup. The entire PSPS project is mainly composed of upper reservoir, lower reservoir, water transmission system, underground plant, and switching station, while this study only includes the upper reservoir area of the PSPS [see Fig. 1(a)]. Among them, the normal storage level of the upper reservoir is 1392 m, with a corresponding storage capacity of 8.96 million m 3 , the dead water level is 1367 m, with a dead storage capacity of 460 000 m 3 , and the effective storage capacity is 8.5 million m 3 . The upper reservoir dam is a concrete panel rockfill dam with a maximum height of 125.90 m, and the elevation and length of the dam crest are 1396 m and 363 m, respectively. At dawn on 31 July 2022, a localized collapse occurred on the south side of the upper reservoir bank. Fortunately, there was no damage to people or equipment, but there is an urgent need for deformation monitoring of the collapse and its surrounding area to assess its surface stability.

1) Site Selection for Instrument Installation:
In order to monitor the deformation of the upper reservoir and especially the collapse area, it is fundamental to ensure that the data acquisition is through-viewed and the instrument installation is stable. The location of the instrument needs to meet the following conditions: 1) no obstacles on the ground affecting the installation and operation of the instrument, with a clear view; 2) avoid areas subject to ground vibrations from site construction; 3) avoid areas with unstable geological structures; 4) avoid areas susceptible to flooding or large changes in the water table; 5) keep the station as close to the surrounding atmosphere as possible to reduce the impact of atmospheric errors. As this is a short-term rapid monitoring, there is no need to build a fixed monitoring pier, i.e., the tripod is screwed to the ground for instrument installation. We organized a professional team to conduct a site survey and completed the site selection [see Fig. 1(a)] and instrument installation [see Fig. 1(b)] based on the actual site environment and installation requirements.
2) Data Acquisition Settings: The GPRI-II instrument with 360°rotational monitoring capability was chosen to meet the  needs of the monitoring range. Before formal monitoring, 1-3 images are first acquired quickly for inspection and their radar intensity should be clear on target and less than 60% of points less than 25 dB. Then, the measurement parameters are adjusted in time according to the quality of the data. The specific measurement parameters for the acquisition of radar data are given in Table I.

C. Experimental Data
Considering the harsh site environment of the upper reservoir, which is under construction, it is difficult to guarantee the absolute safety of instruments and personnel when conducting continuous monitoring for long periods of time. Therefore, from 31 July 2022 onwards, radar data collection will only be carried out during daytime hours. At the same time of radar data collection, professional engineers carried out field investigation, always paying attention to the surface morphology of the reservoir bank, and no obvious collapse and slip area was found again. In order not to interfere with the normal construction schedule of the PSPS, the radar data collection work was completed on August 6, 2022. A total of 227 images of GBR data were acquired. The quality of each image was assessed visually and then images with severe noise effects such as atmospheric and ambient noise were removed, resulting in a final selection of 206 images for subsequent data processing (see Table II). The image acquired near the center of the temporal baseline at 2022/08/03 12:23 was used as the master image and the remaining slave images were registered.
In the polar coordinate system, the intensity map generated by the master image is shown in Fig. 2(a), and the average intensity map generated by 206 temporal radar images after registration is shown in Fig. 2 Fig. 2(a) shows a clear image with less noise influence, indicating a high quality of data collection; Fig. 2(b) contains some noise signals, indicating that the data acquisition process is inevitably interfered by the complex surrounding environment but the overall image quality is satisfactory. In summary, the collected temporal radar data can meet the requirements of subsequent data processing and deformation inversion. In addition, we also used UAV survey technology to conduct field verification and auxiliary analysis of the study area, and generated digital surface model (DSM) and digital orthophoto maps (DOM) according to the aerial photographs. DSM is used to geocode radar data; DOM is mainly used for visual analysis of the study area.

III. METHODS
The GBR data in this article is Ku band, which has a short wavelength and is susceptible to spatial and temporal decoherence. Therefore, we chose the most suitable MT-InSAR method [25]. The MT-InSAR first generates a multimaster interferogram stack from adjacent radar data to reduce the effects of temporal  decoherence and atmospheric noise, then selects PSs based on the inverse of the amplitude deviation and the temporal coherence threshold, filters the interferometric phase corresponding to the PSs to remove atmospheric and noisy phases, and finally performs phase unwrapping to reconstruct the line-of-sight (LOS) directional deformation rate and temporal deformation based on minimum norm theory (see Section III-A). Uniquely in this article, the spatial homogeneity characteristics of the pixels are also taken into account when selecting the PSs, that is, PSs are selected based on the SHPs selected from the temporal amplitude values of the GBR data (see Section III-B).

A. MT-InSAR
The data processing workflow of MT-InSAR is shown in Fig. 3. Considering the short sampling interval of the 206 GBR data and the small variation of the atmosphere in a short period of time, this article generates multimaster interferometric stacks with the interval between two adjacent images as the temporal baseline threshold and also ensures redundant observations. In the end, a total of 409 interferometric stacks are generated. In addition, when the GB-InSAR system conducts deformation monitoring, the radar center position is fixed, so the spatial baseline is almost zero, and the baseline error can be ignored. Assuming, N radar data generates M wrapping interferograms, the unwrapping phase of the ith wrapping interferogram is expressed as [13] where ϕ def , ϕ atm , ϕ noi are the deformation phase, the atmospheric phase, and the noise phase, respectively; k is the ambiguity of the whole cycle. While generating the interferogram stack, PSs need to be selected where the radar scattering characteristics remain stable over the time series. In the classical PS-InSAR method, Ferretti et al. [25] use information on the amplitude deviation of the pixels in the temporal radar data to select the PSs. The amplitude deviation D A means by calculating the amplitude mean m A and standard deviation σ A for each pixel, which is expressed as The smaller the amplitude deviation, the more stable the amplitude information. The amplitude deviation only uses the amplitude information of the pixel, which has the advantages of small calculation and easy extraction. Another way of selecting PSs is based on the average coherence of the interferogram stack [40], [41]. The coherence of all pixels (x, y) in each interferogram needs to be calculated as follows: where I, J are the sliding neighborhood window sizes, respectively; m, s are the master and slave images, respectively; and * represents the complex conjugate.
The meanγ is then calculated from the temporal coherence value γ, and a suitable threshold is set to select the PSs. The higher the mean value of the temporal coherence, the more stable the radar scattering characteristics and the less the ground features are disturbed by noise. To maintain uniformity with the professional GAMMA software, PSs were selected in this article based on the inverse of the amplitude deviation. The threshold value was set to 1.5, and a total of 76 269 PSs were selected, as shown in Fig. 4(a). The threshold value of the temporal coherence meanγ was set to 0.8 and a total of 499 478 PSs were selected, as shown in Fig. 4(b). Generally, the PSs selected by the two methods are combined together for subsequent processing.
Commonly, the sum of the segmental linear deformations is used to represent the deformation values over the entire monitoring time. The deformation phase corresponding to PS (x, y) in the ith wrapping interferogram is Considering that the study area has a wide range of observations and significant height variations, the atmospheric phase is corrected by a linear regression model related to slant range and height. If (x , y ) is set as the reference point, the atmospheric phase can be expressed as where a 0 , a 1 , a 2 are coefficients; r is the slant range between the radar and the point target; and Δh x,y,x ,y is the height difference between the point target and the reference point. The noise phase can be removed by filtering in both the spatial and temporal dimensions. After removing the atmospheric and noise phases, the differential phase between each PS pixel and the reference pixel is expressed as According to (6) and (7), the deformation phase is a linear deformation, which is expressed as The matrix form is where A is an M × N matrix. Singular value decomposition A can get the least squares solution in the sense of the minimum norm, and then the deformation phase of the model can be obtained [42]. ω i x,y,x ,y in (6) is the residual phase, which includes atmospheric residual, nonlinear deformation, and noise. Atmospheric residual and noise can be removed by the filtering method. Then, the nonlinear deformation in the residual phase is summed with the linear deformation to obtain the total deformation of all PSs. Finally, the cumulative deformation is obtained by converting the interferogram stack into single master image form by least squares.

B. PSs Selection of SHPs
According to the above mentioned, we generally select PSs based on indicators such as amplitude deviation and temporal coherence. However, this traditional selection method only takes into account the temporal dimension of the radar data and does not make sufficient use of the spatial relationships between the pixels. Moreover, as the monitoring period grows, the temporal coherence decreases and so does the PSs density. Thus, this article incorporates content of the SHPs from the SqueenSAR method in the selection of the PSs and enhances it slightly [33]. Currently, the algorithms selected for SHPs are mainly classified into nonparametric statistical algorithms and parametric statistical algorithms. In this article, the algorithms in each class that are widely used and have high accuracy are separately experimented in the study area. In this article, the algorithms that are widely used and with high accuracy in each category are tested in the study area respectively.
The BWS test is one of the typical representatives of the nonparametric statistical algorithm, the basic idea of which is to weight the squared values of the empirical distribution differences [43]. It is assumed that statistical calculations are performed on reference pixel B ref and neighborhood pixel B neig in the temporal radar data, namely, as follows: where R i and H i are the rank of the reference pixel sample and the neighborhood pixel sample, respectively. In addition, the parametric statistical algorithm represented by fast SHP selection (FaSHPS) is also widely used [34]. The core principle of the FaSHPS is to transform the hypothesis testing problem into confidence interval estimation. The SHPs of the reference pixel are determined based on their similarity to neighborhood pixels in the time series. When a confidence interval is constructed at a given confidence level 1 − α, all values within the interval are considered as plausible values for the unknown parameter while values outside the interval are discarded. According to the central limit theorem, with the increase of radar data N , the average amplitudeĀ(p) of pixel p is closer to the normal distribution, namely: where P {·} is the probability, z 1−α/2 is the quantile at confidence level 1 − α, μ(p) is the expectation of pixel p, and Var(A(p)) is the true variance of the temporal amplitude of pixel p. According to the statistical theory of radar images, the temporal amplitudes in homogeneous regions obey the Rayleigh distribution. By the property that the coefficient of variation (CV) of the Rayleigh distribution is constant, for radar data with multilook L, the CV A of the temporal amplitude can be expressed as where E(·) is the expectation. Combined with the interval estimation theory, the FaSHPS takes the average amplitude value of the reference pixel as the true value and the average amplitude value of the neighborhood pixel as the pending estimation. Then (11) is substituted as whereμ ref (p) andμ neig are the average of the temporal amplitudes of the reference pixel p and its neighborhood pixels, respectively. Ifμ neig lies within the confidence interval of (13), it is considered as SHP.
Since GBR data have a shorter temporal baseline than spaceborne radar data, some noise with short time correlation will also be considered as backscatter-stabilized scatterers. In addition, the wide scanning area and complex spatial environment of the study area in this article will exacerbate the selection of noisy pixels. Considering the above actual situation, this article introduces intensity deviation, and first excludes the unstable pixels according to the intensity deviation, and then proceeds to SHPs selection. The intensity deviation D I is where σ I and m I are the standard deviation and average of the temporal intensity, respectively. The smaller the intensity deviation, the more stable the temporal backscattering of the pixel. We used a probability density function to count the probability distribution of intensity deviations in the radar data, as shown in Fig. 5. Fig. 5 shows the calculated intensity deviation values in the range 0-15. For intensity deviation values less than about 1.5, the probability is generally high, indicating that the pixels are stable in the temporal dimension, while for intensity deviation values greater than about 2, the probability is low and the curve is smooth, indicating that a small number of pixels are not stable in the temporal dimension. Unstable pixels may correspond to surface features with weak backscatter or may be noisy, which can seriously affect the accuracy of subsequent data processing. Therefore, in this article, we set the probability 0.1 as the threshold value, and pixels with probability of intensity deviation value lower than 0.1 (i.e., intensity deviation value greater than 1.79) are considered as noise and are excluded. The aim is to improve the quality of the selected SHPs.

A. PSs Selection Results and Analysis
In this article, the BWS test and FaSHPS are used for experiments, respectively, and the results of SHPs selection are shown in Fig. 6(a) and (c). After excluding noise pixels according to the intensity deviation threshold, D I + BWS and D I + FaSHPS were used to select SHPs. The results are shown in Fig. 6(e) and (g), respectively. Fig. 6 combined with the intensity map (see Fig. 2) shows that the pixels excluded by the intensity deviation threshold are basically located in regions without backscattering or noise; the SHPs obtained by FaSHPS does not distinguish well between noisy and strongly scattering regions and are worse than the BWS test; the SHPs obtained by D I + FaSHPS compensate well for the shortcomings of FaSHPS and do not differ significantly from D I + BWS. Based on the number of SHPs selected, the pixels can be divided into PSs and distributed scatterers [44], [45]. Previous studies and analyses have concluded that the maximum number of SHPs considered as PSs is set to 10 [37], [40]. Therefore, in this article, pixels with the number of SHPs less than ten were selected as PSs. The results of PSs obtained from BWS test, FaSHPS, D I + BWS, and D I + FaSHPS are shown in Fig. 6(b), (d), (f), and (h), respectively. Furthermore, in order to compare the advantages and deficiencies of each algorithm more objectively, we have evaluated them quantitatively, as given in Table III. Table III lists that after noise pixels are excluded by the intensity deviation threshold, the number of PSs is significantly reduced, with little difference between D I + BWS and   of 771 577 PSs were filtered, and then a mask was set up to filter out some noisy PSs again based on indicators such as phase standard deviation and atmospheric noise, resulting in the retention of 618 173 monitoring points [see Fig. 7(d)]. The number and density of monitoring points determine whether the deformation characteristics of the study area can be reflected. Therefore, the selection of PSs considering SHPs proposed in this article not only effectively excludes the noisy points with unstable temporal scattering but also significantly increases the number and distribution density of monitoring points, which can reflect the deformation characteristics of the study area in more detail and has excellent applicability.

B. Deformation Monitoring Results and Analysis
By combining the D I + FaSHPS algorithm with traditional PSs selection, this article uses MT-InSAR to invert the surface deformation of 618 173 monitoring points. In order to combine the analysis and discussion with the field landform, we carried out geocoding during data processing, and finally converted the deformation monitoring results from the polar coordinate system to the WGS84 geographic coordinate system for display. The LOS deformation rate calculated from 206 radar data is shown in Fig. 8. Fig. 8 shows that the overall surface deformation characteristics in the study area are stable, with deformation rates mostly in the range of −1.0-1.0 mm/d. There are several areas of significant surface deformation in the study area, which are distributed in the reservoir basin, reservoir banks, and surrounding vegetation areas of the PSPS, with a minimum deformation rate of −11.7 mm/d and a maximum rate of 11.8 mm/d.
In the reservoir basin of PSPS, the rugged ground surface and construction have resulted in some areas being decoherent and surface deformation cannot be monitored. On the reservoir banks of PSPS, where the view is open and the surface is mostly made of concrete, the coherence is higher, and thus, the monitoring points are relatively dense and the surface deformation can be well monitored. In the vegetated areas around PSPS, the short radar wavelengths and weak penetration lead to sparse monitoring points and large fluctuations in deformation rates, which to some extent affect the accuracy of the monitoring results.
The results in Fig. 8 illustrate that we can identify the deformation areas but are not yet able to capture the surface deformation characteristics for each day of the period from 2022/07/31 10:24 to 2022/08/06 14:57. Therefore, we depict the temporal cumulative deformation of the study area, as shown in Fig. 9. The temporal figure shows that several obvious deformation areas increase in extent and magnitude with time over the monitored time period. In the study area, the maximum and minimum values of cumulative deformation are 91.3 mm and −92.9 mm, respectively.

C. Discussion of Reservoir Bank Results
Considering that the reservoir basin is the main site of construction, its surface deformation is mainly affected by construction, so it is of little significance to analyze it specifically. In addition, the vegetated areas around the PSPS are far from the center of construction operations and it is more practical to use aerial photography equipment, such as UAV, for regular inspections or manual field surveys. Therefore, this article focuses on the analysis of the deformation monitoring results of the reservoir bank. The reservoir bank has mostly been concreted and its top is relatively flat and solid, its slope is relatively steep, even nearly perpendicular to the ground in some places. Often debris and stones slip or even roll down the slope of the reservoir bank. The reservoir bank is also the area most prone to safety accidents, which can easily cause damage to workers and property. The collapse that occurred at dawn on 31 July 2022 was located in the reservoir bank [see Fig. 1(c)], which is one of the main reasons why this article focuses on the analysis of reservoir bank deformation.
In order to analyze the surface deformation monitoring results of the reservoir bank more clearly, we only kept the monitoring points covering the reservoir bank, as shown in Fig. 10. The deformation rate results from the reservoir bank show that there is no significant surface deformation in most areas, except for a few sporadic areas with small magnitudes of deformation. We selected five relatively obvious deformation areas for focus analysis, whose geographical location layout is shown in the black rectangle box in Fig. 10, and named them I-V. In addition, with the help of the DSM generated by the UAV images, the Areas I-V are displayed in 3-D visualization. The debris from the collapse in Area I is piled up as a slope body, with the upper part of the slope deformed more significantly and to a greater extent than the lower part. Area II is more stable overall, with deformation in some of the monitoring points at the bottom. Area III has a gentler slope and deformation is concentrated in the lower part of the slope. In Area IV, the slope is reinforced with lattice anchor retaining walls, but there is still surface deformation. Area V and Area III have similar topography, with a gentle slope and deformation sporadically distributed in the lower part of the slope. In summary, the 3-D visualization enables us to understand and study the deformation characteristics of the reservoir bank in a more intuitive and visual way.
We selected a suitable number of feature points at different geographical locations in Areas I-V. The geographical locations and numbers of all feature points are shown in Fig. 10, and the temporal deformation trends of all feature points are depicted in Fig. 11. In addition to the statistics of deformation rate and cumulative deformation of all feature points, we also calculated the Hurst index for all feature points using the rescaled range (R/S) analysis method [46], and the results are given in Table IV. The Hurst index of the temporal deformation values lies between 0 and 1. When the Hurst index is close to 0, it means that the temporal deformation trend is weaker and may have anticontinuity; when the Hurst index is 0.5, it means that the temporal deformation trend is random; when the Hurst index is close to 1, it means that the temporal deformation trend is more continuous [47]. In Area I: the deformation trend at P1 and P2 is moderately continuous, with a cumulative deformation value of approximately 15 mm; the deformation trend at P3, P4, and P5 is less continuous, and the surface at P4 and P5 is basically stable; overall, the deformation trend gradually decreases from the top to the bottom of the collapse body. In Area II: P6-P8 are at the bottom of the reservoir bank, where the safety risk is low, and their cumulative deformation values are around 10 mm and the deformation trend in the latter days tends to be stable and not strong in continuity. In Area III: P9 has a dramatic fluctuation in deformation trend until 2022/8/3, after which the deformation trend is more or less stable; P10 has very obvious deformation at 2022/8/3 and 2022/8/4, and is relatively flat for the rest of the time; overall, the surface almost tends to stabilize. In Area IV: P11 and P12 have obvious deformation trends, with cumulative deformation values of 49.23 mm and 33.07 mm respectively, and their Hurst indexes are greater than 0.5, at 0.60 and 0.57, respectively. It can be seen that although this area has been slope supported, there is still approximately 5 mm/d of surface movement, which is slightly continuous. Therefore, powerful measures for safety protection and focused inspection and monitoring of this area will be required in the future. In Area V: the deformation at P13 tends to converge; the deformation at P14 has an almost linear continuous trend, which requires on-site verification and protection. By analyzing the temporal deformation in Areas I-V, we can understand and grasp the deformation situation of the reservoir bank and also provide reliable data support for the relevant technical staff to identify and analyze geological hazard potential areas.

D. Deformation Factors Discussion
According to the results of the deformation rate and temporal deformation in the study area, deformation is found to occur to varying degrees in the reservoir basin, reservoir bank, and surrounding vegetation areas. The analysis of the results, combined with the field survey, revealed that the factors affecting surface deformation fall into two main categories: natural geography and man-made activities. Sudden large-scale deformation, such as collapse, is mainly influenced by natural geographical factors, while slow small-scale deformation, such as slip, is mainly influenced by man-made activities.

1) Natural Geographical Factors:
The upper reservoir area is dominated by Middle Devonian strata and Indochinese intrusive granodiorite strata. The lithology of the upper section of the Middle Devonian Gudaoling Formation is mainly white and grayish marble and dark gray crystalline limestone [48]. The Indochinese intrusion lithology is gray-white granodiorite with semiautomorphic fine-grained structure and blocky tectonics. The lithology of the reservoir basin and reservoir bank is mainly crystalline limestone, and in local areas, there is a pattern of intercalation between marble and crystalline limestone, and basically no fracture structure is developed. In addition, the groundwater in the upper reservoir is mainly bedrock fracture water, with a small amount of carbonate karst fracture water locally. Fissure leakage is the main potential leakage pattern in the reservoir area. Especially in the thin parts of the mountain, in    the pass, and in the part that runs through the adjacent valley are the locations where seepage can occur. Fortunately, through field investigation, it is found that the key areas of the reservoir bank are reinforced with concrete, and support measures are taken to minimize the influence of natural geographical factors, and the surface is basically stable.
An example is the collapse of the reservoir bank, where the collapse area is adjacent to the mountain pass in a thin area. The optical images before and after the collapse are shown in Fig. 12. In Fig. 12(a), there are some striped black voids on the back edge of the reservoir bank before the collapse, which are more similar to ground fissures. This also warns the responsible personnel to promptly identify key locations where deformation may be induced according to the natural geomorphology and geological conditions, so as to ensure that the construction of the PSPS is safe and error-free.
2) Man-Made Activity Factors: There are two main types of man-made activity factors affecting surface deformation: on-site construction and unscheduled blasting. At the upper reservoir basin, where dozens of machines, including heavy trucks and excavators, work in tandem every day; on-site photography of one moment is shown in Fig. 13(a). In the area of construction operations, the ground is accompanied by slight vibrations. In addition, there will be unscheduled blasting during the construction of PSPS in order to facilitate excavation of soil and rock. The location of the blasting during data collection was located in the southwest corner of the upper reservoir, and a distant view just after blasting is shown in Fig. 13(b). On-site verification of the areas of significant deformation based on the results of the surface deformation identification revealed that deformation of the reservoir basin was mainly caused by construction. However, deformation of the reservoir bank can be  affected by both construction and blasting but the reservoir bank was relatively less affected by construction due to its distance from the construction site.
Two blasts occurred during the monitoring period on 3 August 2022 and 5 August 2022, respectively. In order to discuss the impact of blasting on surface stability, we take the focus Area I as an example and analyze the temporal deformation characteristics of its P1-P5 feature points in 3 August 2022 and 5 August 2022, as shown in Fig. 14. Fig. 14(a) shows that the blasting on 3 August 2022 resulted in the most obvious change in the deformation of P1 feature point, and the other feature points were also significantly affected. Fig. 14(b) shows that the effect of blasting on the deformation amplitude of each feature point on 5 August 2022 is generally less than that on 3 August 2022. This indicates that with the deposition of collapse debris over time, the surface stability is enhanced, and the influence of blasting is also reduced. Generally speaking, blasting causes short periods of ground shaking, resulting in ground deformation fluctuations, but causes little or no accelerated ground deformation over time. In addition to paying attention to the safety of personnel and equipment during blasting, the probability of secondary hazards caused by blasting is low.

V. CONCLUSION
This article uses the GB-InSAR technology, which has the advantage of large-angle rotational scanning to achieve surface deformation monitoring in the upper reservoir area of a PSPS in western China. We collected 206 GBR data during the monitoring time period from 2022/07/31 10:24 to 2022/08/06 14:57 and obtained temporal deformation information for the study area based on the MT-InSAR method. Some conclusions can be drawn as following.
1) In selecting PSs based on GBR data, this article takes into account information from SHPs. Based on the FaSHPS algorithm, we introduce intensity deviations to exclude pixels with unstable intensity performance in the temporal dimension and propose the D I + FaSHPS algorithm, which has higher accuracy in the selection of SHPs and computational efficiency. Therefore, in addition to selecting PSs using both the inverse of the amplitude deviation and the temporal coherence, this article incorporates the D I + FaSHPS algorithm screening of PSs. Experimental results show that this not only enriches the selection method of PSs but also effectively improves the distribution density of monitoring points. 2) This article identifies several obvious deformation areas, including the collapse area, from the deformation rate and temporal deformation results. We have analyzed and discussed in detail the deformation areas, especially those located on the reservoir bank, in relation to both physical geography and man-made activities in conjunction with the on-site survey, and found that the surface is relatively stable overall, except for the effects of man-made activities. In addition, we found that there are several deformation areas on the reservoir bank that have a tendency to continue to deform and require powerful measures for safety protection and prior inspection and monitoring. Throughout the study, the GB-InSAR technology is more portable and flexible than monitoring means, such as space-borne InSAR, and has its own unique application advantages in the identification and analysis of such small-scale surface deformation in PSPS.