MT-InSAR deformation monitoring method in combination with reservoir bank landslide deformation spatial-temporal characteristics

Abstract

The invention discloses an MT-InSAR deformation monitoring method in combination with reservoir bank landslide deformation spatial-temporal characteristics. The method comprises the following steps: S1, selecting SAR data; step S2, constructing small baseline set InSAR data, including SAR image registration, interferogram combination, differential interference processing, phase filtering, coherence estimation and interferogram optimization combination; s3, based on the dual-polarized SAR intensity image, performing earth surface classification in combination with random forest classification, gradient and direction constraint and a deep neural network; then adaptive homogeneous filtering is carried out on the interferometric phase based on a secondary classification result; s4, combining a plurality of Sigmoid functions and a linear creep and periodic function to obtain an adaptive deformation model of the reservoir bank landslide; s5, performing deformation parameter nonlinear search by using a genetic algorithm to determine an optimal parameter set; on the basis of low-frequency deformation, atmospheric delay and decoherence noise are separated through low-pass filtering of a space domain and high-pass filtering of a time domain, and finally a deformation time sequence and an average deformation rate are recovered.

Description

Technical Field

[0001] The invention relates to the technical field of reservoir bank landslide monitoring, and in particular to an MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation. Background Art

[0002] The strain softening of the reservoir bank rock mass caused by the reservoir impoundment can easily lead to the rebirth and resurgence of reservoir bank landslides. In addition, the dynamic water pressure generated by the rapid rise or fall of the reservoir water level will significantly reduce the stability of the landslide, which poses a potential threat to the safety of the reservoir.

[0003] Traditional monitoring methods can only obtain sparse measurement results at a few locations of reservoir bank landslides, such as GNSS measurement, leveling measurement, and crack meter monitoring. This poses a challenge to understanding the spatiotemporal evolution characteristics of landslide movement. Synthetic aperture radar interferometry (InSAR) is a ground observation technology developed in recent decades. With its high imaging resolution, non-contact and ability to trace historical deformation, this method has received special attention from researchers in reservoir bank landslide deformation monitoring. However, in addition to the deformation component, the interference phase also contains components such as terrain residuals, atmospheric delays, and de-correlated noise. This will interfere with the deformation phase and make it impossible to extract accurate deformation signals. To this end, many researchers combine SAR images with time-series InSAR (MT-InSAR) methods to suppress the de-coherence noise and atmospheric delay of the interferogram, so as to obtain the deformation time series of the study area. This is of great significance for understanding the evolution process and deformation mechanism of reservoir landslides. The commonly used MT-InSAR methods include permanent scatterer InSAR (PS-InSAR), small baseline set InSAR (SBAS-InSAR) technology and temporary coherence point InSAR (TCP-InSAR) method.

[0004] However, in the process of using InSAR technology to monitor reservoir bank landslides, there are still some bottleneck problems that have not been fully solved. On the one hand, reservoir bank landslides are generally distributed in areas with fewer buildings, resulting in a limited number of permanent scatterers in these areas. If PS-InSAR technology is used for monitoring, it will be difficult to obtain sufficiently dense monitoring points, thus losing the advantages of InSAR technology. In order to overcome this problem, some scholars have proposed adding distributed scatterers (DS) as monitoring targets on the basis of permanent scatterers to increase the density of landslide monitoring points. This method first uses the statistical characteristics of SAR intensity information in time series to determine the homogeneous points of the monitoring target, and then filters the phase of the homogeneous points to improve the phase quality of the DS target. Early studies mainly used the neighborhood pixels of the reference pixel for homogeneous filtering, such as the region growing method. Recently, some scholars have proposed various statistical hypothesis testing methods to select homogeneous pixels, including non-parametric tests and parametric tests. The non-parametric test method identifies homogeneous pixels by comparing the differences between the cumulative distribution functions of time samples. For example, Ferretti et al. used the Kolmogorov-Smirnov test to select statistically homogenous pixels. On this basis, other researchers have proposed methods such as the Baumgartner-Wei-Schindler test, the Kullback-Leibler test, and the Anderson-Darling test. However, these methods require the calculation and comparison of the cumulative distribution functions of all pixels, resulting in low computational efficiency. The parameter test directly compares the statistical characteristics of two samples. However, in the presence of deformation and atmospheric phase, these methods are difficult to make reliable assumptions about the overall distribution. In addition, they ignore the temporal properties of the samples, thereby reducing the reliability of DS selection. SAR image classification can automatically determine pixels as a series of meaningful homogeneous types based on backscattering characteristics, providing a promising idea for the selection of homogeneous pixels in urban areas. However, unlike cities, the scattering intensity of pixels in reservoir bank landslides tends to be consistent, which easily makes the entire landslide area divided into the same homogeneous area. Obviously, the homogeneous points selected by this method are unreasonable, which can easily lead to poor phase optimization effects.

[0005] On the other hand, although the MT-InSAR method can effectively suppress decoherence noise, the interferogram contains terrain residuals, atmospheric delay and other components in addition to the deformation phase. Spatiotemporal filtering can suppress the influence of atmospheric delay on the interference phase. However, the terrain residual is difficult to be accurately separated from the deformation signal. Its influence on the interference phase is proportional to the spatial baseline. Usually, in the process of MT-InSAR data processing, we can assume that the deformation time series satisfies a priori deformation model. Based on this, the functional relationship between the interference phase and the deformation model parameters and the terrain residual will be established. In the deformation modeling strategy, the linear model is the most commonly used priori deformation model, which can be used to describe the deformation characteristics that move slowly in time (such as tectonic movement, ground subsidence, landslide creep, etc.). For specific deformations with significant periodic characteristics such as frozen soil, a simple linear model is prone to lose seasonal deformation signals. We can add periodic functions on the basis of the linear model to invert the deformation time series. In addition, the polynomial model can be used to fit the deformation of active volcanoes and sudden deformation in mining areas. It should be noted that the accuracy of the terrain residual and the deformation time series depends on the reliability of the priori deformation model. Generally, the deformation of reservoir landslides is mainly driven by rapid changes in hydraulic conditions (e.g., water level drop, water level rise, heavy rainfall), and shows significant step deformation characteristics in the time series. This step deformation undergoes a cyclic evolution process of stability, accelerated deformation, and stability. Obviously, the commonly used mathematical models are difficult to accurately describe the deformation law of reservoir landslides. If these models are used to invert the deformation of reservoir landslides, the deformation signal will be misunderstood as a terrain residual signal due to the low degree of fit between the true deformation and the empirical model, thereby reducing the accuracy of the deformation field. Summary of the invention

[0006] In view of the above technical problems, the present invention provides a MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation.

[0007] In order to solve the above technical problems, the technical solution proposed by the present invention is:

[0008] A MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation specifically includes the following steps:

[0009] Step S1, according to the demand for monitoring the deformation of reservoir bank landslide, multi-track high-resolution SAR data covering the monitoring area can be selected;

[0010] Step S2, constructing a small baseline set of InSAR data based on multi-period SAR images, including SAR image registration, interferogram combination, differential interferometry processing, phase filtering, coherence estimation and interferogram optimization combination;

[0011] Step S3, based on the dual-polarization SAR intensity image, combining random forest classification, slope and aspect constraints and deep neural network to achieve fine surface classification of the target area; then adaptive homogeneity filtering of the interference phase is performed based on the secondary classification result;

[0012] Step S4, combining multiple Sigmoid functions and linear creep and periodic functions to obtain an adaptive deformation model of the reservoir bank landslide;

[0013] Step S5, first use the genetic algorithm to perform nonlinear search for deformation parameters to determine the optimal parameter set; based on the low-frequency deformation, separate the atmospheric delay and incoherence noise through low-pass filtering in the spatial domain and high-pass filtering in the temporal domain, and finally restore the deformation time series and the average deformation rate.

[0014] As a further improvement of the above technical solution:

[0015] Preferably, the step S2 includes the following steps:

[0016] S2-1, SAR image registration: First, it is necessary to select the main image from the SAR time-series images, perform the registration process of the single main image, and crop the registered SAR image to obtain the SAR data of the monitoring area;

[0017] S2-2, Interference pattern combination: Estimate the spatiotemporal baseline between any two images based on the orbital data of the SAR image, select the spatiotemporal baseline threshold to combine the interferogram; construct a short baseline subset in the form of multiple main images, and combine the interferogram with a single subset;

[0018] S2-3, differential interferometry processing: interferometric processing is performed on the SAR image with a short temporal and spatial baseline to obtain an interferometric phase map;

[0019] S2-4, phase filtering: adaptive filtering to suppress the influence of noise;

[0020] S2-5, coherence estimation and interference pattern optimization combination: The interference quality of the interference pattern is evaluated by coherence analysis; the coherence coefficient of the interference phase describes the similarity between the same-name regions of the primary and secondary images.

[0021] Preferably, the step S3 includes the following steps:

[0022] S3-1, secondary image classification: firstly, random forest classification is performed using dual-polarization intensity texture maps, then the random forest classification results are optimized based on the slope aspect map constraints of the reservoir bank landslide, and the optimized classification results are used as training labels to run secondary classification using a deep neural network;

[0023] S3-2, connected domain construction: using SAE convolutional network to learn the features of ground object categories and their spatial distribution;

[0024] S3-3, homogeneous pixel selection: add a classifier after feature extraction and use a set of labeled samples for supervised training; obtain a SAE model with classification function to classify all SAR image pixels; based on the classification results, identify homogeneous pixels of each pixel in a predefined neighborhood;

[0025] S3-4, Interferometric Phase Optimization: In a set of homogeneous pixels, the optimal interferometric phase value for each pixel is reconstructed from the covariance matrix by maximizing the joint probability density function.

[0026] Preferably, in the construction of the connected domain, SAE constructs a deep structure by stacking multiple autoencoders, and trains each network using a greedy hierarchical learning method; the feature representation of the previous network will be used as the input of the latter; suppose a normalized n-dimensional input data set x, hidden layer y i (i∈1,…m) is expressed as:

[0027] y i =f i (wx+b)=f i (θ,x) (1)

[0028] Among them, f i is the autoencoder function, w represents the weight matrix, b represents the bias, and θ is the parameter set of w and b. The output layer x′ can be obtained from the m-dimensional hidden layer and the decoder:

[0029] x′≈g(f(θ,x))(2)

[0030] Where g is the decoding function. In the formula, the learning parameter set θ obtains the feature representation;

[0031] The optimal training parameters are determined by minimizing the loss function L:

[0032] L(x,x′)=‖xx′‖ 2 (3).

[0033] Preferably, in the interference phase optimization, the optimal interference phase value is:

[0034]

[0035] In the formula, Λ Ω is the complex covariance matrix of the target pixel, represents the optimized interference phase, is the optimized interference phase.

[0036] Preferably, the step S4 includes the following steps:

[0037] S4-1, Sigmoid function parameter adjustment: adjust the shape parameters of each Sigmoid function in the function group;

[0038] S4-2, observation equation of adaptive deformation model: according to the components of interference phase and the acquisition time of the main and auxiliary images of the interference pattern, the observation equation is established; the equation includes deformation phase, terrain residual phase, atmospheric delay phase and incoherence noise;

[0039] S4-3, hypothesis test: test the significance of the adaptive deformation model based on the phase observation value of the small baseline set, and then test the significance of each parameter. The parameters retained after the hypothesis test are the parameters to be determined of the adaptive deformation model.

[0040] Preferably, the Sigmoid function parameter adjustment method is:

[0041] Assume that the monitoring area is at time t0, t1, t2, …, t N-1 A total of N SAR images were acquired; the images were combined into M interferometric pairs according to the temporal and spatial baseline thresholds. According to the D-InSAR principle, the phase in the interferogram generated by any two SAR images The composition is expressed as:

[0042]

[0043] in, is the deformation phase of the radar line of sight (LOS) at any time; is the atmospheric delay phase; is random noise; represents the phase associated with the terrain residual, which is:

[0044]

[0045] Among them, B ⊥ is the vertical baseline of the interference pair, Δ6 is the terrain residual of the pixel, θ and λ are the radar incident angle and wavelength respectively, and R0 is the slant distance from the SAR sensor to the ground; in addition, the MT-InSAR method represents the deformation phase in the interference diagram as:

[0046]

[0047] Among them, Δt is the time increment between two SAR images, v represents the deformation rate in the LOS direction, and a and f are the parameters of the periodic sine and cosine functions respectively;

[0048] The Sigmoid function is expressed as:

[0049]

[0050] Among them, the shape parameters b and w represent the center position and slope of the curve respectively; the inflection point moment and time span of the step characteristic accelerated deformation are described by adjusting the size of the parameters.

[0051] Preferably, the observation equation of the adaptive deformation model is:

[0052] The deformation model is expressed as:

[0053]

[0054] Where d is the cumulative deformation at time t, and c is the parameter for adjusting the amplitude of the Sigmoid function. Each Sigmoid function in the formula represents a deformation acceleration signal, which has a year as the period according to the deformation law of reservoir landslide.

[0055] For t A and t B Any high coherence pixel in the differential interferogram generated by the SAR image at any moment can construct the observation equation:

[0056]

[0057] Among them, B ⊥ is the interferometric vertical baseline, Δh is the pixel terrain residual, R0 is the slant distance from the SAR sensor to the ground, t B The deformation, t A The deformation, is the residual phase.

[0058] Preferably, the step S5 specifically includes the following steps:

[0059] S5-1, GA algorithm parameter estimation: GA algorithm is used to encode the problem to be solved into a string through a genetic coding method such as binary coding or tree coding, and an initial population is randomly generated in the solution space. The genes in the population include linear rate, terrain residual, periodic function parameters and step function coefficients;

[0060] S5-2, high-frequency deformation extraction: Calculate the low-frequency signal of each interference pair using the solved deformation parameters and terrain residuals; Subtract the modeled deformation phase and terrain residual phase from the interference phase observation equation to obtain the residual phase; The residual phase contains atmospheric delay, decoherence noise and high-frequency deformation phase; Perform spatial domain low-pass filtering and temporal domain high-pass filtering on the residual phase to restore the high-frequency deformation phase;

[0061] S5-3, deformation time series estimation: add the low-frequency deformation phase and the high-frequency deformation phase to obtain the complete deformation phase in the interference phase; solve the deformation time series according to the conversion relationship between deformation and phase;

[0062] S5-4, average deformation rate: The deformation time series is used to reveal the temporal evolution characteristics of the reservoir bank landslide; the spatial pattern of deformation is described by the average deformation rate.

[0063] Preferably, in step S5, the low-frequency deformation phase and the high-frequency deformation phase are added to obtain a complete deformation phase in the interference phase; M unwrapped interference phases The phase time series of N SAR images The relationship between them is:

[0064]

[0065] Where B represents the design matrix between the interferogram and the SAR image; the least squares method is used to solve the time series of the deformation phase:

[0066]

[0067] Through the conversion relationship between phase and deformation, the temporal deformation of reservoir bank landslide is restored.

[0068] The MT-InSAR deformation monitoring method combined with the spatiotemporal characteristics of reservoir bank landslide deformation provided by the present invention has the following advantages compared with the prior art:

[0069] (1) The MT-InSAR deformation monitoring method of the present invention, which combines the spatiotemporal characteristics of reservoir bank landslide deformation, can achieve fine surface classification of the target area based on dual-polarization SAR intensity images, combined with random forest classification, slope and aspect constraints, and deep neural networks. Subsequently, the interference phase is adaptively homogeneously filtered, which can improve the phase quality while ensuring the accuracy of the InSAR deformation monitoring results.

[0070] (2) The MT-InSAR deformation monitoring method of the present invention combines the spatiotemporal characteristics of reservoir bank landslide deformation, and proposes an interference phase optimization method combining random forest and deep neural network to address the problems of poor interference phase quality in the reservoir bank landslide study area and difficulty in balancing the accuracy and density of monitoring point targets. And by adjusting the amplitude, boundary and center position of multiple Sigmoid functions, this method has been widely used in the field of reservoir landslide displacement prediction.

[0071] (3) The MT-InSAR deformation monitoring method of the present invention, which combines the temporal and spatial characteristics of reservoir bank landslide deformation, can perform refined InSAR deformation monitoring for complex reservoir bank landslide scenarios. BRIEF DESCRIPTION OF THE DRAWINGS

[0072] Figure 1 It is the overall flow chart of the MT-InSAR monitoring method for reservoir bank landslide of the present invention.

[0073] Figure 2 It is a schematic diagram of adjusting the shape parameters of the Sigmoid function to describe the step deformation according to the present invention.

[0074] Figure 3 It is the InSAR average deformation rate map of a typical reservoir bank landslide in the existing technology.

[0075] Figure 4 This is a GNSS accuracy verification diagram of the InSAR time-series deformation of a typical reservoir bank landslide in the existing technology. DETAILED DESCRIPTION

[0076] The specific embodiments of the present invention are described in detail below. It should be understood that the specific embodiments described herein are only used to illustrate and explain the present invention, and are not used to limit the present invention.

[0077] like Figure 1 As shown, the present invention proposes an MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation, which specifically includes the following steps:

[0078] Step S1, SAR data acquisition

[0079] According to the needs of reservoir bank landslide deformation monitoring, multi-orbit high-resolution SAR data covering the monitoring area can be selected.

[0080] Step S2: Small baseline set InSAR interferometry processing

[0081] The data product provided by the SAR satellite platform is a single-look complex (SLC) image, which contains amplitude and phase information. The phase is used to measure the distance between the observed target and the sensor. Deformation is the change in distance between the time intervals between the acquisition of two images. Interferometric processing of two SAR images of repeated orbits to restore deformation phase information is the basis of all time-series InSAR processing techniques. Usually, in addition to the deformation phase, a single interference pattern also includes other components such as terrain phase, atmospheric delay and incoherence noise. The present invention constructs a small baseline set InSAR data based on multi-time period SAR images, including SAR image registration, interference pattern combination, differential interference processing, phase filtering, coherence estimation and interference pattern optimization combination.

[0082] The interferometric phase processing flow before MT-InSAR parameter estimation and temporal deformation inversion in the present invention comprises the following steps:

[0083] S2-1, SAR image registration

[0084] There are certain spatial and temporal differences in the SAR images of repeated orbits, and these geometric deviations will affect the interferometric phase. First, it is necessary to select the main image from the SAR time-series images and perform the registration process of the single main image. The image registration accuracy reaches the sub-pixel level, generally reaching 0.0001 pixels. The registered SAR image is cropped to obtain the SAR data of the monitoring area.

[0085] S2-2, Interference pattern combination

[0086] The temporal and spatial baselines between any two images are estimated based on the orbital data of the SAR images. In order to effectively suppress the influence of decoherence and terrain residuals on the interferometric phase, the temporal and spatial baseline thresholds are selected to combine the interferograms. In the reservoir bank landslide scenario, short baseline subsets are constructed in the form of multiple main images. And in order to avoid rank deficiency during the inversion of temporal deformation, the interferograms are combined in a single subset.

[0087] S2-3, Differential Interference Processing

[0088] Interferometric processing is performed on SAR images with short temporal and spatial baselines to obtain interferometric phase maps. 30m resolution SRTM DEM data can be used for geocoding of SAR images and terrain phase removal. Polynomial fitting can be used for orbit error removal.

[0089] S2-4, phase filtering

[0090] The phase of the incoherent noise in the interference pattern is random, and filtering can effectively suppress the influence of noise. Adaptive (ADF) filtering takes the minimum mean square error as the criterion and can automatically adjust the autocorrelation function between noise and phase to optimize the filter effect. In addition, the non-local filtering method can maintain good deformation edge characteristics.

[0091] S2-5, Coherence estimation and interference pattern optimization combination

[0092] The interference quality of the interferogram is evaluated by coherence analysis. The coherence coefficient of the interferometric phase is used to describe the similarity between the regions with the same name in the primary and secondary images, and is between [0,1]. 0 indicates complete loss of coherence, and 1 indicates that the scattering characteristics of the two images are completely consistent. Interferograms with low coherence in the deformation area should be eliminated to avoid gross errors in subsequent parameter estimation and temporal deformation.

[0093] Step S3, neural network interference phase optimization

[0094] In the reservoir bank landslide area, the presence of decorrelation factors such as noisy speckle noise, large surface motion and spatial phase gradient, and strong atmospheric effects reduces the quality of the interferometric phase, making it difficult to balance the accuracy and density of the extracted monitoring point targets.

[0095] In view of the poor quality of the interferometric phase in the reservoir bank landslide area, and the difficulty in taking into account the target accuracy and density of the monitoring points, in order to increase the number of effective monitoring points, this invention will adopt a phase optimization method based on a neural network for the InSAR interferometric phase. First, based on the dual-polarization SAR intensity image, the fine surface classification of the target area is achieved by combining random forest classification, slope and aspect constraints, and deep neural networks; then, based on the secondary classification results, the interferometric phase is adaptively homogeneously filtered to improve the phase quality while ensuring the accuracy of the InSAR deformation monitoring results.

[0096] The specific implementation plan is:

[0097] S3-1, secondary image classification

[0098] Firstly, random forest classification is performed using the dual-polarization intensity texture map. Then, the random forest classification results are optimized based on the slope aspect map constraints of the reservoir bank landslide. The optimized classification results are used as training labels, and a deep neural network is used to run secondary classification.

[0099] Random forest classification is a typical supervised learning algorithm, that is, an ensemble algorithm based on decision trees. Its "randomness" is reflected in two aspects: first, back-sampling from the original sample set; second, reverse sampling from the original features, and then selecting the best features from them. The subsets of samples and features form a decision tree, and multiple decision trees form a "forest". In this "forest", all trees are integrated together to vote on the classification results. Finally, the category with the highest number of votes is output as the classification result.

[0100] Deep neural networks have many layers and wide widths, and can be mapped to any function, with better classification results. In actual application, 30% of the labeled samples are randomly selected for model training, and the rest are used as validation sets. During training, the monitoring area is divided into multiple categories according to the requirements of actual homogeneous area selection and the backscattering intensity characteristics of typical land cover.

[0101] S3-2, Connected Domain Construction

[0102] The present invention uses SAE convolutional network to learn the features of ground object categories and their spatial distribution. The model uses unsupervised mode to learn features from unlabeled raw data, which can better describe the features of raw images and improve the accuracy of image classification results.

[0103] Specifically, SAE builds a deep structure by stacking multiple autoencoders and trains each network using a greedy hierarchical learning method. The feature representation of the previous network will be used as the input of the latter. Given a normalized n-dimensional input dataset x, the hidden layer y i (i∈1,…m) can be expressed as:

[0104] yi =f i (wx+b)=f i (θ,x) (1)

[0105] Among them, f i is the autoencoder function, w represents the weight matrix, b represents the bias, and θ is the parameter set of w and b. The output layer x′ can be obtained from the m-dimensional hidden layer and the decoder:

[0106] x′≈g(f(θ,x))(2)

[0107] Where g is the decoding function. In the formula, the learning parameter set θ obtains the feature representation.

[0108] The optimal training parameters are determined by minimizing the loss function L:

[0109] L(x,x′)=‖xx′‖ 2 (3).

[0110] S3-3, homogeneous pixel selection

[0111] The SAE model after unsupervised training is just a feature extractor without classification capability. It is necessary to add a classifier after feature extraction and use a set of labeled samples for supervised training. The classifier uses SVM, random forest or softmax model.

[0112] This training process will use error back propagation to fine-tune the network parameters, and finally obtain a SAE model with classification function to classify all SAR image pixels. Based on the classification results of the SAE network, homogeneous pixels of each pixel can be identified within a predefined neighborhood. To ensure efficiency, all pixels can use the same rectangular window. Within this window, the neighbors within the eight-connected domain of the target pixel are marked as homogeneous pixels.

[0113] S3-4, Interferometric Phase Optimization

[0114] In a set of homogeneous pixels, the optimal interferometric phase value for each pixel can be reconstructed from the covariance matrix by maximizing the joint probability density function:

[0115]

[0116] In the formula, Λ Ω is the complex covariance matrix of the target pixel, represents the optimized interference phase, is the optimized interference phase.

[0117] The optimized interference phase effectively suppresses the incoherence noise of low-coherence monitoring point targets, thereby improving the phase reliability and point density.

[0118] Step S4: InSAR adaptive deformation model of reservoir bank landslide

[0119] In the MT-InSAR method, the rationality of the spatiotemporal evolution characteristics of deformation depends on the reliability of the prior deformation model. The MT-InSAR method suppresses the influence of spatiotemporal decorrelation and atmospheric delay on the quality of interference phase by combining interference patterns. In reservoir landslide scenes with few artificial buildings and severe decoherence, this method is usually based on multiple main images. The present invention adds a Sigmoid function that can better describe the significant step characteristics of reservoir bank landslides affected by reservoir water level and heavy rainfall on the basis of linear and periodic functions to construct an adaptive function group that describes complex deformation characteristics.

[0120] The invention combines a plurality of Sigmoid functions and traditional linear creep and periodic functions to obtain an adaptive deformation model of a reservoir bank landslide.

[0121] In the present invention, constructing the InSAR adaptive deformation model includes the following steps:

[0122] S4-1, Sigmoid function parameter adjustment

[0123] The Sigmoid function is a monotonically increasing, bounded and symmetrical S-shaped curve function between 0 and 1. It contains two nonlinear shape parameters, such as the inflection point of accelerated deformation and the time span in formula (6). In order to make the function group reasonably describe the stability, accelerated deformation and stable periodic motion characteristics of the reservoir bank landslide, the shape parameters of each Sigmoid function in the function group are adjusted in this embodiment. The specific parameter adjustment effect is as follows: Figure 2 .

[0124] Specifically, assume that the monitoring area is at time t0, t1, t2, …, t N-1 A total of N SAR images were acquired. According to the temporal and spatial baseline threshold, these images were combined into M interferometric pairs. According to the D-InSAR principle, the phase in the interferogram generated by any two SAR images The composition is expressed as:

[0125]

[0126] in, is the deformation phase of the radar line of sight (LOS) at any time; is the atmospheric delay phase; is random noise. represents the phase associated with the terrain residual, which is:

[0127]

[0128] Among them, B ⊥ is the vertical baseline of the interference pair, Δh is the terrain residual of the pixel, θ and λ are the radar incident angle and wavelength respectively, and R0 is the slant distance from the SAR sensor to the ground. In addition, in the interference pattern, the traditional MT-InSAR method usually uses linear deformation or periodic model to represent the deformation phase:

[0129]

[0130] Among them, Δt is the time increment of the two SAR images, v represents the deformation rate in the LOS direction, and a and f are the parameters of the periodic sine and cosine functions, respectively. Considering that N SAR images can be combined into M interference pairs, the observation equation of pixel by pixel or differential arc segment can be established based on high coherence points. The least squares method can be used to solve the terrain residual and deformation parameters. However, the deformation model is the key to parameter estimation and time series deformation inversion. If the deformation model is unreasonable, it will affect the estimation accuracy of deformation and DEM error. The Sigmoid function is a monotonically increasing, bounded and symmetric S-shaped curve function between 0-1, which can be expressed as:

[0131]

[0132] Among them, the shape parameters b and w represent the center position and slope of the curve respectively; by adjusting the size of these two parameters, it can be used to describe the turning point and time span of the step characteristic accelerated deformation. Therefore, this function can reasonably reflect the movement characteristics of reservoir landslides in time series.

[0133] S4-2, observation equation of adaptive deformation model

[0134] According to the components of the interference phase and the acquisition time of the main and secondary images of the interference pattern, the observation equation is established. The equation includes deformation phase, terrain residual phase, atmospheric delay phase and decoherence noise. In the deformation phase, considering that the value range of the Sigmoid function is between 0-1, an amplitude parameter is set before each Sigmoid function to adjust the size of the step deformation. Due to the different MT-InSAR frameworks, the observation equations will be different. The pixel-by-pixel phase equation is established in the SBAS-InSAR framework based on the unwrapped phase. However, based on PS-InSAR and TCP-InSAR, the phase difference equation between arc segments is established, and the equation contains the difference in phase components between the pixels that constitute the arc segment.

[0135] Specifically, based on the traditional deformation model of MT-InSAR, the Sigmoid function is embedded as a function group to perform parameter estimation and time series deformation inversion. The new deformation model is expressed as:

[0136]

[0137] Where d is the cumulative deformation at time t, and c is the parameter for adjusting the amplitude of the Sigmoid function. Each Sigmoid function in the formula represents a deformation acceleration signal, which has a year as the period according to the deformation law of reservoir landslide.

[0138] Combining the basic principles of D-InSAR, for t A and t B Any high coherence pixel in the differential interferogram generated by the SAR image at any moment can construct the observation equation:

[0139]

[0140] Among them, B ⊥ is the interferometric vertical baseline, Δh is the pixel terrain residual, R0 is the slant distance from the SAR sensor to the ground, t B The deformation, t A The deformation, is the residual phase.

[0141] S4-3, Hypothesis Testing

[0142] Since the real deformation is complex and changeable, in order to accurately describe the deformation characteristics on the long time series, the present invention first tests the significance of the adaptive deformation model with the F distribution at a significance level of 0.01 based on the phase observation value of the small baseline set, and then tests the significance of each parameter in the segment with the Tu distribution at a significance level of 0.01 to eliminate the insignificant parameters to be determined. The parameters retained after the hypothesis test are the parameters to be determined of the adaptive deformation model.

[0143] Step S5: InSAR deformation estimation of reservoir bank landslide

[0144] InSAR deformation includes average deformation rate and time series, which are used to reveal the spatial and temporal evolution characteristics of reservoir bank landslides. Among them, deformation parameter estimation is a key step in MT-InSAR, which includes deformation parameters and terrain residuals. Considering that the Sigmoid function is a nonlinear function, the shape parameters in the observation equation cannot be directly solved by the least squares method.

[0145] The present invention aims at the problem of nonlinear parameters in the deformation model. First, a genetic algorithm is used to perform nonlinear search of deformation parameters to determine the optimal parameter set. On the basis of low-frequency deformation, atmospheric delay and incoherent noise are separated by low-pass filtering in the spatial domain and high-pass filtering in the time domain, and finally the deformation time series and average deformation rate are restored. The following steps are included:

[0146] S5-1, GA algorithm parameter estimation

[0147] The present invention utilizes GA algorithm to encode the problem to be solved into the form of a string through a genetic coding method of binary coding or tree coding, and randomly generates an initial population in the solution space. The genes in the population include linear rate, terrain residual, periodic function parameters and step function coefficients.

[0148] In the present invention, the minimum two norm of the residual phase is set as the fitness function to evaluate the convergence condition of the population. First, the size of each initial individual gene is set, that is, the solution space search range of each parameter to be determined; secondly, the fitness function value corresponding to each gene is calculated and it is determined whether the iteration termination convergence condition is met. If not, the population is optimized by selection, crossover and mutation operations to generate a new population, and the fitness function value is repeatedly calculated. Until the iteration termination condition of the fitness function is met, the optimization operation is terminated, and the population gene generated by the process is used as the global optimal solution for the deformation model estimation parameters.

[0149] Specifically, the derivative of the Sigmoid function with respect to a variable can be expressed by itself:

[0150] sigmoid′(t)=sigmoid(t)·(1-sigmoid(t)) (11)

[0151] The traditional method uses Taylor series expansion to linearize the observation equation and obtain the error equation of the interference phase:

[0152]

[0153] Among them, X (3n+2)×1 =[v,c1,b1,w1,…c n ,b n ,w n ,Δh] T are the parameters to be determined for the deformation model and DEM error; is the deformation model phase corresponding to the approximate value of the nonlinear parameter; A is the coefficient matrix of the deformation parameter in the observation equation, expressed as:

[0154]

[0155] Combining the basic principles of the Gauss-Markov model, the traditional method uses the least squares method to obtain the deformation parameters

[0156]

[0157]

[0158] The deformation parameters and terrain residuals in the adaptive deformation model are obtained.

[0159] Although the series expansion method can obtain the parameters to be solved by the least squares method, the real deformation in the reservoir bank landslide is complex and changeable, and it is difficult to obtain the approximate value of the nonlinear parameters of each pixel, and it is easy to introduce model errors in the parameter inversion process. The GA algorithm first encodes the problem to be solved into a string form through a genetic coding method such as binary coding or tree coding, and randomly generates an initial population in the solution space. The genes in the population include linear rate v, terrain residual Δh, periodic function parameters [a, f], and step function coefficients [c, b, w]. In the present invention, the least squares norm of the residual phase is set as the fitness function to evaluate the convergence condition of the population, as shown in the following formula:

[0160]

[0161] First, the size of each initial individual gene is set, that is, the search range of the solution space of each parameter to be determined; secondly, the fitness function value corresponding to each gene is calculated and it is determined whether the iterative termination condition of formula (15) is met. If not, the population is optimized by selection, crossover and mutation to generate a new population, and the fitness function value is calculated repeatedly. When the iterative termination condition of the fitness function is met, the optimization operation is terminated, and the population gene generated by this process is used as the global optimal solution of the deformation model estimation parameters.

[0162] S5-2, high frequency deformation extraction

[0163] In this embodiment, the low-frequency signal of each interference pair is calculated using the solved deformation parameters and terrain residuals. The modeled deformation phase and terrain residual phase are subtracted from the interference phase observation equation to obtain the residual phase. The residual phase contains atmospheric delay, incoherence noise and high-frequency deformation phase. Because the atmospheric delay and incoherence noise of the high-frequency deformation (non-modeled deformation) signal have different characteristics in terms of temporal and spatial characteristics.

[0164] In this embodiment, the residual phase is first decomposed into the SAR time series based on the least squares method. High-frequency deformation shows low-frequency characteristics in both space and time domains. Combined with the characteristics of atmospheric delay, low-frequency space and high-frequency time, and decoherence noise, both of which are high-frequency signals in space and time, the residual phase on the SAR image time series is low-pass filtered in the space domain and high-pass filtered in the time domain to restore the high-frequency deformation phase.

[0165] S5-3, Deformation time series estimation

[0166] By adding the low-frequency deformation phase and the high-frequency deformation phase, the complete deformation phase in the interference phase can be obtained. Based on the combination of interference patterns, the functional relationship between the time-series SAR image and the interference pattern is established. All interference patterns belong to the same baseline subset, and the least squares method can be used to directly solve the time series of the deformation phase. According to the conversion relationship between deformation and phase, the deformation time series can be solved.

[0167] The low-frequency deformation phase and the high-frequency deformation phase are added together to obtain the complete deformation phase in the interference phase. The phase time series of N SAR images The relationship between can be modeled as:

[0168]

[0169] Among them, B represents the design matrix between the interferogram and the SAR image. All interferograms belong to the same baseline subset, and the least squares method can directly solve the time series of the deformation phase:

[0170]

[0171] Through the conversion relationship between phase and deformation, the temporal deformation of reservoir bank landslide can be restored.

[0172] S5-4, average deformation rate

[0173] Deformation time series can be used to reveal the temporal evolution characteristics of reservoir bank landslides. However, the spatial pattern of deformation is usually described by the average deformation rate.

[0174] In this embodiment, the least squares solution is performed on the cumulative deformation of the time series SAR images element by element to obtain the average deformation rate in the time series to reveal the spatial evolution of the landslide.

[0175] Finally, the temporal and spatial evolution characteristics of reservoir bank landslides were obtained.

[0176] The above implementation cases are only preferred embodiments of the present invention and are not intended to limit the present invention in any form. Although the present invention has been disclosed as above with preferred embodiments, they are not intended to limit the present invention. Therefore, any simple modification, equivalent changes and modifications made to the above embodiments according to the technical essence of the present invention without departing from the content of the technical solution of the present invention shall fall within the scope of protection of the technical solution of the present invention.

Claims

1. A MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation, characterized in that: The specific steps include: Step S1, according to the demand for monitoring the deformation of reservoir bank landslide, multi-track high-resolution SAR data covering the monitoring area can be selected; Step S2, constructing a small baseline set of InSAR data based on multi-period SAR images, including SAR image registration, interferogram combination, differential interferometry processing, phase filtering, coherence estimation and interferogram optimization combination; Step S3, based on the dual-polarization SAR intensity image, combining random forest classification, slope and aspect constraints and deep neural network to achieve fine surface classification of the target area; then adaptive homogeneity filtering of the interference phase is performed based on the secondary classification result; Step S4, combining multiple Sigmoid functions and linear creep and periodic functions to obtain an adaptive deformation model of the reservoir bank landslide; Step S5, first use the genetic algorithm to perform nonlinear search for deformation parameters to determine the optimal parameter set; based on the low-frequency deformation, separate the atmospheric delay and incoherence noise through low-pass filtering in the spatial domain and high-pass filtering in the temporal domain, and finally restore the deformation time series and the average deformation rate.

2. The MT-InSAR deformation monitoring method according to claim 1, characterized in that: The step S2 includes the following steps: S2-1, SAR image registration: First, it is necessary to select the main image from the SAR time-series images, perform the registration process of the single main image, and crop the registered SAR image to obtain the SAR data of the monitoring area; S2-2, Interference pattern combination: Estimate the spatiotemporal baseline between any two images based on the orbital data of the SAR image, select the spatiotemporal baseline threshold to combine the interferogram; construct a short baseline subset in the form of multiple main images, and combine the interferogram with a single subset; S2-3, differential interferometry processing: interferometric processing is performed on the SAR image with a short temporal and spatial baseline to obtain an interferometric phase map; S2-4, phase filtering: adaptive filtering to suppress the influence of noise; S2-5, coherence estimation and interference pattern optimization combination: The interference quality of the interference pattern is evaluated by coherence analysis; the coherence coefficient of the interference phase describes the similarity between the same-name regions of the primary and secondary images.

3. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 2 is characterized in that: The step S3 includes the following steps: S3-1, secondary image classification: firstly, random forest classification is performed using dual-polarization intensity texture maps, then the random forest classification results are optimized based on the slope aspect map constraints of the reservoir bank landslide, and the optimized classification results are used as training labels to run secondary classification using a deep neural network; S3-2, connected domain construction: using SAE convolutional network to learn the features of ground object categories and their spatial distribution; S3-3, homogeneous pixel selection: add a classifier after feature extraction and use a set of labeled samples for supervised training; obtain a SAE model with classification function to classify all SAR image pixels; based on the classification results, identify homogeneous pixels of each pixel in a predefined neighborhood; S3-4, Interferometric Phase Optimization: In a set of homogeneous pixels, the optimal interferometric phase value for each pixel is reconstructed from the covariance matrix by maximizing the joint probability density function.

4. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 3 is characterized in that: In the construction of the connected domain, SAE constructs a deep structure by stacking multiple autoencoders and trains each network using a greedy hierarchical learning method; the feature representation of the previous network will be used as the input of the latter; suppose a normalized n-dimensional input data set x, hidden layer y i (i∈1,…m) is expressed as: y i =f i (wx+b)=f i (θ,x) (1) Among them, f i is the autoencoder function, w represents the weight matrix, b represents the bias, and θ is the parameter set of w and b. The output layer x′ can be obtained from the m-dimensional hidden layer and the decoder: x′≈g(f(θ,x)) (2) where g is the decoding function. In the formula, the learning parameter set θ obtains the feature representation; The optimal training parameters are determined by minimizing the loss function L: L(x,x′)=‖x-x′‖ 2 (3)。 5. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 4 is characterized in that: In the interference phase optimization, the optimal interference phase value is: In the formula, Λ Ω is the complex covariance matrix of the target pixel, represents the optimized interference phase, is the optimized interference phase.

6. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 3 is characterized in that: The step S4 includes the following steps: S4-1, Sigmoid function parameter adjustment: adjust the shape parameters of each Sigmoid function in the function group; S4-2, observation equation of adaptive deformation model: according to the components of interference phase and the acquisition time of the main and auxiliary images of the interference pattern, the observation equation is established; the equation includes deformation phase, terrain residual phase, atmospheric delay phase and incoherence noise; S4-3, hypothesis test: test the significance of the adaptive deformation model based on the phase observation value of the small baseline set, and then test the significance of each parameter. The parameters retained after the hypothesis test are the parameters to be determined of the adaptive deformation model.

7. The MT-InSAR deformation monitoring method according to claim 6, characterized in that: The Sigmoid function parameter adjustment method is: Assume that the study area is at time t0, t1, t2, …, t N-1 A total of N SAR images were acquired; the images were combined into M interferometric pairs according to the temporal and spatial baseline thresholds. According to the D-InSAR principle, the phase in the interferogram generated by any two SAR images The composition is expressed as: in, is the deformation phase of the radar line of sight (LOS) at any time; is the atmospheric delay phase; is random noise; represents the phase associated with the terrain residual, which is: Among them, B ⊥ is the vertical baseline of the interference pair, Δh is the terrain residual of the pixel, θ and λ are the radar incident angle and wavelength respectively, and R0 is the slant distance from the SAR sensor to the ground; in addition, the MT-InSAR method represents the deformation phase in the interference diagram as: Among them, Δt is the time increment between two SAR images, v represents the deformation rate in the LOS direction, and a and f are the parameters of the periodic sine and cosine functions respectively; The Sigmoid function is expressed as: Among them, the shape parameters b and w represent the center position and slope of the curve respectively; the inflection point moment and time span of the step characteristic accelerated deformation are described by adjusting the size of the parameters.

8. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 7 is characterized in that: The observation equation of the adaptive deformation model is: The deformation model is expressed as: Where d is the cumulative deformation at time t, and c is the parameter for adjusting the amplitude of the Sigmoid function. Each Sigmoid function in the formula represents a deformation acceleration signal, which has a year as the period according to the deformation law of reservoir landslide. For t A and t B Any high coherence pixel in the differential interferogram generated by the SAR image at any moment can construct the observation equation: Among them, B ⊥ is the interferometric vertical baseline, Δh is the pixel terrain residual, R0 is the slant distance from the SAR sensor to the ground, t B The deformation, t A The deformation, is the residual phase.

9. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 6 is characterized in that: The step S5 specifically includes the following steps: S5-1, GA algorithm parameter estimation: GA algorithm is used to encode the problem to be solved into a string through a genetic coding method such as binary coding or tree coding, and an initial population is randomly generated in the solution space. The genes in the population include linear rate, terrain residual, periodic function parameters and step function coefficients; S5-2, high-frequency deformation extraction: Calculate the low-frequency signal of each interference pair using the solved deformation parameters and terrain residuals; Subtract the modeled deformation phase and terrain residual phase from the interference phase observation equation to obtain the residual phase; The residual phase contains atmospheric delay, decoherence noise and high-frequency deformation phase; Perform spatial domain low-pass filtering and temporal domain high-pass filtering on the residual phase to restore the high-frequency deformation phase; S5-3, deformation time series estimation: add the low-frequency deformation phase and the high-frequency deformation phase to obtain the complete deformation phase in the interference phase; solve the deformation time series according to the conversion relationship between deformation and phase; S5-4, average deformation rate: The deformation time series is used to reveal the temporal evolution characteristics of the reservoir bank landslide; the spatial pattern of deformation is described by the average deformation rate.

10. The MT-InSAR deformation monitoring method combining the temporal and spatial characteristics of reservoir bank landslide deformation according to claim 3 is characterized in that: In step S5, the low-frequency deformation phase and the high-frequency deformation phase are added to obtain a complete deformation phase in the interference phase; M unwrapped interference phases The phase time series of N SAR images The relationship between them is: Where B represents the design matrix between the interferogram and the SAR image; the least squares method is used to solve the time series of the deformation phase: Through the conversion relationship between phase and deformation, the temporal deformation of the reservoir bank landslide is restored.

Images