HyP3 cloud platform-based earth surface deformation monitoring method, medium and electronic equipment

Abstract

The invention provides a surface deformation monitoring method based on a HyP3 cloud platform, a medium and electronic equipment. The method comprises the following steps: acquiring an unwrapping graph set; an optimized unwrapping graph set is obtained through automatic selection of unwrapping graphs; a corrected unwrapping graph set is obtained through unwrapping graph error correction; and through network inversion, atmospheric errors, orbital errors and DEM errors are removed, and a final data processing result, namely a deformation time sequence, is obtained. According to the method, a user is supported to carry out SAR data processing by using data distributed by the HyP3 cloud platform and own demands, and two steps are added on the basis of a traditional SBAS method, so that the precision of a deformation result based on HyP3 cloud platform data processing is equivalent to the precision of directly using original data, and the efficiency is greatly improved.

Description

Technical Field

[0001] The present invention relates to the field of geographic information technology, and in particular to a surface deformation monitoring method, medium and electronic equipment based on a HyP3 cloud platform. Background Art

[0003] Since the 1990s, the differential synthetic aperture radar interferometry (D-InSAR) technology has been widely used in geological disaster deformation monitoring and has become one of the important technologies for geological disaster deformation monitoring. In the 21st century, with the successive launch of satellites equipped with SAR sensors, the amount of SAR data has increased rapidly, and multi-temporal InSAR (MT-InSAR) technology has emerged. Compared with traditional D-InSAR technology, MT-InSAR technology can obtain surface displacement data with millimeter-level accuracy and realize long-term deformation monitoring, which has important application value. In addition, the application scope of MT-InSAR technology is also constantly expanding. It has been applied to the monitoring of surface displacement caused by natural disasters such as ground subsidence, earthquakes, landslides and volcanoes, as well as groundwater mining, coal mining and building deformation monitoring.

[0004] Although MT-InSAR technology has significant advantages in surface deformation monitoring and is widely used, its data processing requires solid radar-related theoretical knowledge and high-performance hardware resources. Therefore, some organizations provide SAR earth observation systems and online InSAR data processing services based on cloud computing. Among them, LiCSAR is easy to use, has an adapted time series analysis package, and can quickly obtain large-scale deformation data with relatively reliable accuracy. However, the spatial resolution of LiCSAR products is only 100m, which is more suitable for large-scale structural deformation monitoring and difficult to meet the needs of small and medium-scale surface deformation monitoring, such as ground subsidence deformation monitoring. In addition, LiCSAR products have no available data in some areas and are less customizable.

[0005] In contrast, the HyP3 cloud platform is a cloud computing-based online service platform for Sentinal-1SAR data processing, which supports users to process SAR data according to their own needs. However, there are relatively few studies on the application of this cloud platform, and the HyP3 cloud platform product does not have the phase closure property and other problems, resulting in low accuracy of the results obtained based on the product, which limits its application and popularization. Summary of the invention

[0006] The purpose of the present invention is to overcome the defects of the prior art and provide a surface deformation monitoring method, medium and electronic equipment based on the HyP3 cloud platform, which supports users to process SAR data according to their own needs and has high deformation monitoring accuracy.

[0007] To achieve the above object, the present invention adopts the following technical solutions:

[0008] In a first aspect of the present application, a surface deformation monitoring method based on the HyP3 cloud platform is provided, the method comprising:

[0009] Get the disentangled graph collection;

[0010] Obtain optimized disentanglement graph set through automatic selection of disentanglement graphs;

[0011] Obtain a corrected unwrapping graph set by correcting the unwrapping graph error;

[0012] After network inversion, atmospheric error removal, orbit error removal and DEM error removal, the final data processing result is obtained, and the data processing result includes a deformation time series.

[0013] Optionally, the automatic selection of the disentanglement graph includes:

[0014] Use the minimum spanning tree algorithm to obtain a set of disentangled graphs, called the minimum disentangled graph set;

[0015] The coherence map is used to calculate the average spatial coherence of each disentangled image, and a set of high-coherence disentangled images is obtained by setting a spatial coherence threshold.

[0016] The optimized disentangled graph obtained by merging the minimum disentangled graph set and the high coherence disentangled graph set is called the optimized disentangled graph set.

[0017] Optionally, before correcting the unwrapping image error, a phase unwrapping error correction algorithm based on compressed sensing is used to correct the erroneous pixel phase.

[0018] Optionally, the phase unwrapping error correction algorithm includes:

[0019] Perform reference phase recovery, the calculation formula is as follows:

[0020] φ=φ raw +φ ref (1)

[0021] where φ raw is the original unwrapped phase, φ ref is the reference phase of the unwrapping graph, φ is the phase after the reference phase is restored, and the reference phase of the HyP3 unwrapping graph is obtained through the distributed parameter file;

[0022] Phase unwrapping error correction based on compressed sensing, phase unwrapping error correction is to establish constraints on the unwrapping graph using the phase closure property. The time dimension constraint of the unwrapping phase is expressed as:

[0023]

[0024] where φ α，β，γ is the phase closure error at the coherent point, k α , k β , k γ are the phase unwrapping error correction numbers of α, β, and γ, respectively. α, β, and γ are three unwrapping images that form a closed loop in the time domain. Filtering and multi-view operations will destroy the phase closure property of the unwrapping image, and it is necessary to perform rounding operations on the closed phase. round[·] represents the integer operator;

[0025] Assuming that the phase unwrapping error in the unwrapping graph set is sparse, that is, for a pixel set specified in the time domain, the number of pixels with phase unwrapping error is much smaller than the number of elements in the set. Under this assumption, compressed sensing is used to calculate the phase unwrapping error and recover the phase, L 0 The optimization criteria are:

[0026] argmin||k||0 (3)

[0027] Where k is the vector of the phase unwrapping error correction number to be estimated. Solving formula (3) is an NP-hard problem. In order to quickly solve the phase unwrapping error, we can replace formula (3) with L 0 The problem is transformed into L 1 The norm problem is solved, and the constraints are as shown in formula (2), which is L 1 The objective function of the norm is as follows:

[0028] argmin||k||1 (4)

[0029] In addition, the filtering and multi-look steps in D-InSAR data processing will break the phase closure property of the unwrapped image in the time domain, namely:

[0030] φ α +φ β +φ γ ≠0 (5)

[0031] φ α ,φ β ,φ γ are the unwrapped phases of α, β, and γ respectively;

[0032] The rounding operator in formula (2) will lead to incorrect constraints and thus make it impossible to obtain a feasible solution. Therefore, a slack variable is added to the integer linear programming problem, and formula (2) becomes

[0033]

[0034] The optimization function becomes:

[0035]

[0036] where k i and i are the target variables of the optimization function, respectively representing the phase unwrapping error correction number and relaxation variable at the pixel point, wk i and wy i are their weights respectively; K is the vector of phase unwrapping error at the specified pixel point, expressed as K = [k1, ..., k s ] T , Y is expressed as Y=[y1,...,y M ] T , L is the Mx1 vector of closure errors calculated for all constraints, A is the M×S design matrix, B is the identity matrix, S and M are the number of unwrapped graphs and the number of constraints, respectively. Note that the slack variable y needs to be given here i The value is much larger than k i The weight of; Formula (7) is a nonlinear integer programming problem, which is transformed into an integer linear programming problem by adding slack variables and solved using the interior point method;

[0037] To ensure that vector K is sparse, we use L 1 Norm instead of L 0 The premise of solving the compressed sensing problem is to use the norm, that is, when the number of phase unwrapping errors is significantly smaller than the number of interference patterns, then L 0 The problem is effectively reduced to 1 To solve the minimization problem, we need to evaluate the phase unwrapping error correction result for each pixel. The sum of the absolute values ​​of the sparse vectors is represented as cost. When the cost is less than a given threshold, a good result is obtained through compressed sensing. The formula is as follows:

[0038]

[0039] cost<th (9)

[0040] Where th is a given threshold; in order to meet the prerequisite of sparsity, th is usually taken as 10% of the number of disentangled graphs.

[0041] Optionally, the set of unwrapping graphs uses the same reference point, and the phases of the unwrapping graphs are all phase gradients relative to the reference point, expressed as:

[0042] φ adj =φ corr -φ refnew (10)

[0043] where φ adj is the unwrapped phase after the unified reference, φ corr is the unwrapped phase after the phase unwrapping error is corrected, φ refnewThe unwrapped phase of the newly selected reference.

[0044] Optionally, the network inversion, atmospheric error removal, orbit error removal and DEM error removal are processed using open source MintPy software.

[0045] In a second aspect of the present application, a computer-readable storage medium is provided, in which a computer program is stored. When the computer program is loaded and executed by a processor, the method described in any one of the first aspects is adopted.

[0046] In a third aspect of the present application, an electronic device is provided, comprising a memory, a processor, and a computer program stored in the memory and capable of running on the processor, wherein when the processor loads and executes the computer program, the method described in any one of the first aspects is adopted.

[0047] In summary, due to the adoption of the above technical solution, the beneficial effects of the present invention are:

[0048] The present invention designs a time-series surface deformation InSAR monitoring technology route based on the HyP3 cloud platform product. This technology route aims at the problem of large phase unwrapping error. On the basis of the traditional small baseline set InSAR (SBAS-InSAR) technology route, it adds the steps of automatic selection of unwrapping images and correction of phase unwrapping errors to suppress the influence of phase unwrapping errors on deformation monitoring results, and has high reliability for surface deformation monitoring. BRIEF DESCRIPTION OF THE DRAWINGS

[0049] Figure 1 This is the MT-InSAR workflow diagram based on the Hyp3 product of the present invention;

[0050] Figure 2 This is an optical image of Yan'an New District in Yan'an City;

[0051] Figure 3 It is a spatiotemporal baseline network;

[0052] Figure 4 is the optimized spatiotemporal baseline network;

[0053] Figure 5 It is a temporal coherence comparison diagram;

[0054] Figure 6 is the temporal coherence statistical histogram;

[0055] Figure 7 is the phase closure error statistical histogram;

[0056] Figure 8 is the pixel point with conflicting conditions and the cost of each point y value.

[0057] Fig. 9 is the statistical histogram of deformation rate differences;

[0058] Fig.10 It is a time series comparison chart of point P using different methods. DETAILED DESCRIPTION

[0059] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, rather than all the embodiments. Based on the embodiments in the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0060] See also Figure 1 The present embodiment provides a surface deformation monitoring method based on the HyP3 cloud platform, including: obtaining a set of unwrapped maps; automatically selecting and obtaining an optimized set of unwrapped maps through unwrapped maps; obtaining a corrected set of unwrapped maps through unwrapped map error correction; and obtaining a final data processing result through network inversion, atmospheric error removal, orbit error removal, and DEM error removal, wherein the data processing result includes a deformation time series.

[0061] The automatic selection of disentangled graphs includes: in order to ensure the connectivity of the spatiotemporal baseline network, a set of disentangled graphs is obtained using the minimum spanning tree algorithm, which is called the minimum disentangled graph set; the average spatial coherence of each disentangled graph is calculated using the coherence graph, and a high coherence disentangled graph set is obtained by setting the spatial coherence threshold; the minimum disentangled graph set and the high coherence disentangled graph set are merged to obtain the optimized disentangled graph, which is called the optimized disentangled graph set.

[0062] The unwrapped graphs containing a large number of phase unwrapping errors are excluded by the unwrapped graph automatic selection algorithm, but there are still pixels with phase unwrapping errors in the unwrapped graphs. The phase unwrapping error correction algorithm based on compressed sensing can better identify and correct the wrong phase unwrapped pixels and show good results. Therefore, the phase unwrapping error correction algorithm based on compressed sensing is used to correct the wrong pixel phase before the unwrapped graph error correction, and three aspects are introduced: reference phase recovery, phase unwrapping error correction based on compressed sensing, and quality function. The unwrapped graph set obtained after the phase unwrapping error correction is performed on the optimized unwrapped graph set is called the corrected unwrapped graph set.

[0063] The phase unwrapping error correction algorithm includes:

[0064] The unwrapping reference point of the InSAR unwrapping image is called the reference benchmark, and its phase is the reference phase. In the phase unwrapping step, the reference point has no deformation by default, and its phase is set to 0. Therefore, the step of setting the reference point phase to its reference phase is called reference phase recovery. In the phase unwrapping step, the HyP3 unwrapping image set has the problem of inconsistent reference benchmarks and lack of phase closure. Therefore, this paper needs to perform reference phase recovery before correcting the phase unwrapping error. The calculation formula is as follows:

[0065] φ=φ raw +φ ref (1)

[0066] where φ raw is the original unwrapped phase, φ ref is the reference phase of the unwrapping graph, φ is the phase after the reference phase is restored, and the reference phase of the HyP3 unwrapping graph is obtained through the distributed parameter file;

[0067] Phase unwrapping error correction based on compressed sensing, phase unwrapping error correction is to establish constraints on the unwrapping graph using the phase closure property. The time dimension constraint of the unwrapping phase is expressed as:

[0068]

[0069] where φ α，β，γ is the phase closure error at the coherent point, k α , k β , k γ are the phase unwrapping error correction numbers of α, β, and γ, respectively. α, β, and γ are three unwrapping images that form a closed loop in the time domain. Filtering and multi-view operations will destroy the phase closure property of the unwrapping image, and it is necessary to perform rounding operations on the closed phase. round[·] represents the integer operator;

[0070] Assuming that the phase unwrapping error in the unwrapping graph set is sparse, that is, for a pixel set specified in the time domain, the number of pixels with phase unwrapping error is much smaller than the number of elements in the set. Under this assumption, compressed sensing is used to calculate the phase unwrapping error and recover the phase, L 0 The optimization criteria are:

[0071] argmin||k||0 (3)

[0072] Where k is the vector of the phase unwrapping error correction number to be estimated. Solving formula (3) is an NP-hard problem. In order to quickly solve the phase unwrapping error, we can replace formula (3) with L 0 The problem is transformed into L 1 The norm problem is solved, and the constraints are as shown in formula (2), which is L 1 The objective function of the norm is as follows:

[0073] argmin||k||1 (4)

[0074] In addition, the filtering and multi-look steps in D-InSAR data processing will break the phase closure property of the unwrapped image in the time domain, namely:

[0075] φ α +φ β +φ γ ≠0 (5)

[0076] φ α ,φ β ,φ γ are the unwrapped phases of α, β, and γ respectively;

[0077] The rounding operator in formula (2) will lead to incorrect constraints and thus make it impossible to obtain a feasible solution. Therefore, a slack variable is added to the integer linear programming problem, and formula (2) becomes

[0078]

[0079] The optimization function becomes:

[0080]

[0081] where k i and i are the target variables of the optimization function, respectively representing the phase unwrapping error correction number and relaxation variable at the pixel point, wk i and wy i are their weights respectively; K is the vector of phase unwrapping error at the specified pixel point, expressed as K = [k1, ..., k s ] T , Y is expressed as Y=[y1,...,y M ] T , L is the M×1 vector of closure errors calculated for all constraints, A is the M×S design matrix, B is the identity matrix, S and M are the number of unwrapped graphs and the number of constraints, respectively. Note that the slack variable y needs to be given here i The value is much larger than k i The weight of; Formula (7) is a nonlinear integer programming problem, which is transformed into an integer linear programming problem by adding slack variables and solved using the interior point method;

[0082] To ensure that vector K is sparse, we use L 1 Norm instead of L 0 The premise of solving the compressed sensing problem is to use the norm, that is, when the number of phase unwrapping errors is significantly smaller than the number of interference patterns, then L 0 The problem is effectively reduced to1 To solve the minimization problem, we need to evaluate the phase unwrapping error correction result for each pixel. The sum of the absolute values ​​of the sparse vectors is represented as cost. When the cost is less than a given threshold, a good result is obtained through compressed sensing. The formula is as follows:

[0083]

[0084] cost<th (9)

[0085] Where th is a given threshold; in order to meet the prerequisite of sparsity, th is usually taken as 10% of the number of disentangled graphs.

[0086] The set of disentangled graphs uses the same reference point, and the phases of the disentangled graphs are all phase gradients relative to the reference point, expressed as:

[0087] φ adj =φ corr -φ refnew (10)

[0088] where φ adj is the unwrapped phase after the unified reference, φ corr is the unwrapped phase after the phase unwrapping error is corrected, φ refnew The unwrapped phase of the newly selected reference.

[0089] The subsequent MT-InSAR workflow includes network inversion, atmospheric error removal, orbit error removal, and DEM error removal. The final data processing results include deformation time series and annual average deformation rate. The above operations are processed using the open source MintPy software.

[0090] In addition, it should be noted that: an embodiment of the present application also provides a computer-readable storage medium, in which a computer program is stored. When the computer program is loaded and executed by a processor, the method described in the above embodiment is executed.

[0091] In addition, it should be noted that an embodiment of the present application also provides an electronic device, including a memory, a processor, and a computer program stored in the memory and capable of running on the processor. When the processor loads and executes the computer program, the method described in the above embodiment is executed.

[0092] This embodiment designs a time-series surface deformation InSAR monitoring technology route based on the HyP3 cloud platform product. This technology route aims at the problem of large phase unwrapping error. On the basis of the traditional small baseline set InSAR (SBAS-InSAR) technology route, it adds the steps of automatic selection of unwrapping images and correction of phase unwrapping errors to suppress the influence of phase unwrapping errors on deformation monitoring results, and has high reliability for surface deformation monitoring.

[0093] In order to evaluate the performance of the proposed technical route, this paper selects Yan'an New District, Yan'an City, Shaanxi Province as the experimental area. The optical images of the experimental area are as follows: Figure 2 As shown. The Yan'an New District project officially started construction on April 17, 2012, and the geotechnical engineering of the new district was fully completed in August 2018. Studies have observed that ground subsidence had occurred in the Yan'an New District since 2014, and it slowed down in 2018 when the construction was completed. In order to avoid phase decoherence caused by filling and excavation, this paper uses Sentinel-1A SAR data from 2020 to 2021 for experiments. 34 ascending SAR images with orbit number 84 were selected on the ASF website and submitted to HyP3 for D-InSAR processing. The parameters for data processing of the GAMMA software include a multi-view ratio of 10×2, a time baseline threshold of 60 days, and a spatial baseline threshold of 150m. HyP3 generated a total of 126 interference pairs, and the spatiotemporal baseline network of the data is shown as follows. Figure 3 shown.

[0094] The experiment cuts the disentanglement map, DEM and coherence map to the experimental area. Figure 2 The range shown in the optical image is consistent. In the experiment, we used the automatic selection of disentanglement graphs to optimize the disentanglement graph set. The results are shown in Figure 4 As shown in the figure. The dotted line in the figure represents the disentangled graph that needs to be discarded, and the color of the arc represents its corresponding average spatial coherence. It can be seen that the spatial coherence of the discarded disentangled graph is lower than 0.5, and there is no independent subset of the spatiotemporal baseline.

[0095] Phase unwrapping error is the main factor affecting temporal coherence. The more accurate the unwrapping phase is, the higher the corresponding temporal coherence is, and vice versa. This paper calculates the temporal coherence of the optimized unwrapping map set and quantitatively evaluates and analyzes the automatic selection algorithm of the unwrapping map. The temporal coherence of the pixel is defined as:

[0096]

[0097] where φ i is the unwrapped phase, and is the re-estimated winding phase.

[0098] Figure 5Figures a and b are the temporal coherence calculated for all disentangled graph sets and the temporal coherence calculated for the optimized disentangled graph set, respectively. The temporal coherence of the optimized disentangled graph set has been improved overall. Figure 5 The areas with low temporal coherence in the middle image (a) are concentrated in mountainous areas with vegetation. This is because mountainous areas with vegetation are prone to phase decoherence, resulting in phase unwrapping errors, which in turn lead to low temporal coherence. Figure 5 The temporal coherence of the image b in this area has been improved to a certain extent. This is because the experiment uses the automatic unwrapping selection algorithm to exclude the unwrapped images with low average spatial coherence, retain the unwrapped images with relatively high spatial coherence, and reduce the number of phase unwrapping errors in this area, thereby showing a certain degree of improvement in temporal coherence.

[0099] Figure 6 It is a statistical histogram of temporal coherence, indicating the proportion of the number of pixels with different temporal coherence. The lines in the figure respectively indicate the proportion of temporal coherence pixels of all disentangled image sets and optimized disentangled image sets. It can be seen that in the interval of temporal coherence from 0 to 0.93, the number of pixels of all disentangled image sets is more than that of the optimized disentangled image sets, while the number of pixels of the optimized disentangled image sets is more between 0.93 and 1, indicating that the temporal coherence of the optimized disentangled image sets has been improved to a certain extent, and the number of pixels with phase unwrapping errors has been reduced. This also proves that the disentangled image automatically selects to exclude the disentangled image with large phase unwrapping errors, reducing the phase unwrapping error of the overall disentangled image.

[0100] In order to make the disentangled graph set have the phase closure property, the experiment uses formula (2) to restore the reference phase of the disentangled graph, and the phase closure difference before and after the reference phase is restored is statistically calculated, as follows: Figure 7 shown. Figure 7 Middle a picture and Figure 7 Figure b in the middle shows the statistical histograms of the phase closure error before and after the reference phase recovery. It can be seen that the set of unwrapping images that do not perform the reference phase recovery step does not have the phase closure property. This is because the reference benchmark of each unwrapping image is different. In addition to the phase unwrapping error and phase closure error, there is also a systematic error caused by the different unwrapping image benchmarks in the closed phase.

[0101] Figure 8 The pixels in the figure a represent pixels where the constraints conflict, and are called conflicting pixels. Figure 8 Figure b shows the cost obtained using the slack variable NLP algorithm y , where

[0102]

[0103] It is found through statistics that the NLP algorithm with added slack variables can obtain the phase unwrapping error correction number of all pixels. However, the NILP algorithm without adding slack variables cannot obtain a feasible solution that meets the constraints on the conflicting pixels.

[0104] The experiment uses the quality function to evaluate the results of phase unwrapping error correction and masks pixels with cost greater than 10. The average temporal coherence after masking is improved by 0.03, because the quality function can identify pixels that do not meet the sparsity requirement.

[0105] Figure 5 Figure c in the figure shows the temporal coherence of the corrected unwrapped map set. It can be seen that the average temporal coherence of the corrected unwrapped map set is about 0.95, which is generally higher than the optimized unwrapped map set. This shows that the phase unwrapping error correction algorithm corrects part of the phase unwrapping error, thereby improving the temporal coherence. The improvement in areas with a lot of vegetation is not significant enough. This is because the phase unwrapping error of the pixels in this area is difficult to meet the sparsity in the time dimension and is masked in the quality function, so the improvement in average temporal coherence is not significant.

[0106] In order to further verify the reliability of the proposed technical route, the results of processing the corresponding SAR data using GAMMA software were compared. The parameter configuration adopted a multi-view ratio of 10×2, a spatial baseline threshold of ±150m, and a time baseline threshold of 60 days. All unwrapped maps used a unified reference benchmark. GAMMA software was used to remove atmospheric and orbital errors. High-quality unwrapped maps were selected through manual visual interpretation to calculate the deformation time series and deformation rate. This paper subtracts the deformation rate of the GAMMA results from the proposed technical route and statistically analyzes the results, such as Fig.10 As shown, it can be seen that the difference distribution is normally distributed, concentrated near 0, and the range is mainly concentrated between -10mm and 10mm. The proportion of differences falling in the interval (-10mm-10mm) is 85.7%.

[0107] This paper extracts point P (position as Figure 2 The deformation time series is shown in Fig.10 As shown. The two lines in the figure represent the deformation time series of the two technical routes. It can be seen that the settlement trends of the two sets of monitoring results are almost the same, but the settlement curves do not completely overlap. There are three reasons for this: First, the set of disentanglement graphs obtained by manual visual interpretation during the GAMMA processing is not completely consistent with the set of disentanglement graphs automatically selected by the disentanglement graph; secondly, the parameters in the GAMMA software processing process cannot be completely consistent with the processing parameters of HyP3; finally, the different software used in the time series processing will also cause deviations.

[0108] in conclusion

[0109] (1) The phase unwrapping error of the HyP3 unwrapping graph is large and it does not have the phase closure property.

[0110] (2) Both the automatic selection of unwrapping images and the correction of phase unwrapping errors can suppress the influence of phase unwrapping errors on deformation monitoring results to a certain extent. Based on the HyP3 product and the technical route proposed in this paper, MT-InSAR monitoring results that are approximately consistent with those based on original SAR data and GAMMA software can be obtained, proving that this technical route can be used for surface deformation monitoring.

[0111] (3) The reference phase recovery step can restore the reference benchmark of the unwrapped graph set and solve the problem that the HyP3 product does not have the phase closure property.

[0112] (4) Experiments show that compared with the traditional nonlinear integer programming model, the NLP model with added slack variable constraints can also obtain feasible solutions on conflicting pixels, proving that the proposed method can solve the constraint conflict problem in the phase unwrapping error correction algorithm based on compressed sensing.

[0113] (5) The experiment found that compared with the case where no evaluation was performed on the monitoring results, the average temporal coherence of the results after the evaluation function masking was improved by 0.03. This proves that the use of this evaluation function can reduce the number of pixels whose phase unwrapping errors cannot be corrected to a certain extent, reduce the overall phase unwrapping error of the unwrapping image, and improve the accuracy of deformation monitoring.

[0114] The above contents are further detailed descriptions of the present invention in combination with specific preferred embodiments, and it cannot be determined that the specific implementation of the present invention is limited to these descriptions. For ordinary technicians in the technical field to which the present invention belongs, several simple deductions or substitutions can be made without departing from the concept of the present invention, which should be regarded as falling within the protection scope of the present invention.

Claims

1. A surface deformation monitoring method based on HyP3 cloud platform, characterized in that: The method comprises: Get the disentangled graph collection; Obtain optimized disentanglement graph set through automatic selection of disentanglement graphs; Obtain a corrected unwrapping graph set by correcting the unwrapping graph error; After network inversion, atmospheric error removal, orbit error removal and DEM error removal, the final data processing result is obtained, and the data processing result includes a deformation time series.

2. A surface deformation monitoring method based on HyP3 cloud platform as claimed in claim 1, characterized in that: The automatic selection of the disentanglement graph includes: Use the minimum spanning tree algorithm to obtain a set of disentangled graphs, called the minimum disentangled graph set; The coherence map is used to calculate the average spatial coherence of each disentangled image, and a set of high-coherence disentangled images is obtained by setting a spatial coherence threshold. The optimized disentangled graph obtained by merging the minimum disentangled graph set and the high coherence disentangled graph set is called the optimized disentangled graph set.

3. The surface deformation monitoring method based on the HyP3 cloud platform according to claim 1, characterized in that: Before the unwrapping image error correction, a phase unwrapping error correction algorithm based on compressed sensing is used to correct the erroneous pixel phase.

4. The surface deformation monitoring method based on the HyP3 cloud platform as claimed in claim 3, characterized in that: The phase unwrapping error correction algorithm includes: Perform reference phase recovery, the calculation formula is as follows: f=f raw +φ ref (1) where φ raw is the original unwrapped phase, φ ref is the reference phase of the unwrapping graph, φ is the phase after the reference phase is restored, and the reference phase of the HyP3 unwrapping graph is obtained through the distributed parameter file; Phase unwrapping error correction based on compressed sensing, phase unwrapping error correction is to establish constraints on the unwrapping graph using the phase closure property. The time dimension constraint of the unwrapping phase is expressed as: where φ α，β，γ is the phase closure error at the coherent point, k α , k β , k γ are the phase unwrapping error correction numbers of α, β, and γ, respectively. α, β, and γ are three unwrapping images that form a closed loop in the time domain. Filtering and multi-view operations will destroy the phase closure property of the unwrapping image, and it is necessary to perform rounding operations on the closed phase. round[·] represents the integer operator; Assuming that the phase unwrapping error in the unwrapping graph set is sparse, that is, for a pixel set specified in the time domain, the number of pixels with phase unwrapping error is much smaller than the number of elements in the set. Under this assumption, compressed sensing is used to calculate the phase unwrapping error and recover the phase, L 0 The optimization criteria are: argmin||k||0 (3) Where k is the vector of the phase unwrapping error correction number to be estimated. Solving formula (3) is an NP-hard problem. In order to quickly solve the phase unwrapping error, we can replace formula (3) with L 0 The problem is transformed into L 1 The norm problem is solved, and the constraints are as shown in formula (2), which is L 1 The objective function of the norm is as follows: argmin||k||1 (4) In addition, the filtering and multi-look steps in D-InSAR data processing will break the phase closure property of the unwrapped image in the time domain, namely: f α +φ β +φ γ ≠0 (5) φ α ,φ β ,φ γ are the unwrapped phases of α, β, and γ respectively; The rounding operator in formula (2) will lead to incorrect constraints and thus make it impossible to obtain a feasible solution. Therefore, a slack variable is added to the integer linear programming problem, and formula (2) becomes The optimization function becomes: where k i and i are the target variables of the optimization function, respectively representing the phase unwrapping error correction number and relaxation variable at the pixel point, wk i and wy i are their weights respectively; K is the vector of phase unwrapping error at the specified pixel point, expressed as K = [k1, ..., k s ] T , Y is represented by Y=[y1,...,y M ] T , L is the M×1 vector of closure errors calculated for all constraints, A is the M×S design matrix, B is the identity matrix, S and M are the number of unwrapped graphs and the number of constraints, respectively. Note that the slack variable y needs to be given here. i The value is much larger than k i The weight of; Formula (7) is a nonlinear integer programming problem, which is transformed into an integer linear programming problem by adding slack variables and solved using the interior point method; To ensure that the vector K is sparse, we use L 1 Norm instead of L 0 The premise of solving the compressed sensing problem is to use the norm, that is, when the number of phase unwrapping errors is significantly smaller than the number of interference patterns, then L 0 The problem is effectively reduced to 1 To solve the minimization problem, we need to evaluate the phase unwrapping error correction result for each pixel. The sum of the absolute values ​​of the sparse vectors is represented as cost. When the cost is less than a given threshold, a good result is obtained through compressed sensing. The formula is as follows: cost <th (9) Where th is a given threshold; in order to meet the prerequisite of sparsity, th is usually taken as 10% of the number of disentangled graphs.

5. The surface deformation monitoring method based on HyP3 cloud platform according to claim 1, characterized in that: The set of disentangled graphs uses the same reference point, and the phases of the disentangled graphs are all phase gradients relative to the reference point, expressed as: f adj =φ corr -f refnew (10) where φ adj is the unwrapped phase after the unified reference, φ corr is the unwrapped phase after the phase unwrapping error is corrected, φ refnew The unwrapped phase of the newly selected reference.

6. The surface deformation monitoring method based on the HyP3 cloud platform according to claim 1, characterized in that: The network inversion, atmospheric error removal, orbit error removal and DEM error removal are processed using the open source MintPy software.

7. A computer-readable storage medium, wherein a computer program is stored in the computer-readable storage medium, characterized in that: When the computer program is loaded and executed by a processor, the method according to any one of claims 1 to 6 is adopted.

8. An electronic device comprising a memory, a processor, and a computer program stored in the memory and capable of running on the processor, characterized in that: When the processor loads and executes the computer program, the method according to any one of claims 1 to 6 is adopted.

Images