A Method for Unwrapping Large Gradient Phases in InSAR of Glacial Regions

Abstract

The present invention provides a method for unwrapping large gradient phases in InSAR in glacier areas, comprising the following steps: generating training samples through an interferogram simulation technique based on the differential value of DEM in glacier areas, training the training samples through a symmetric phase unwrapping network model constructed by mixing a convolutional layer and a Transformer, and then inputting the predicted phase discontinuity information as a probability mass map into the maximum flow / minimum cut algorithm to perform phase wrapping count estimation and complete phase unwrapping. The present invention utilizes the powerful learning and summarization ability and data mining ability of deep learning, and gets rid of the dependence of phase unwrapping on the assumption of phase continuity.

Description

Technical Field

[0001] The present invention belongs to the technical field of spaceborne synthetic aperture radar interferometry, and particularly relates to a method for unwrapping large gradient phases in InSAR in a glacier area. Background Technique

[0002] Under the background of global warming, glaciers are in a state of rapid mass loss and gradually becoming unstable. Due to the huge volume and mass of glacier materials, large-scale sudden detachment of glaciers (also known as glacier collapse) can trigger serious mountain disasters. Recent studies have found that glaciers often show obvious signs of surging before large-scale sudden detachment. In addition, glacier surging itself can also cause some mountain disasters, such as debris flows caused by meltwater generated by the rapid ablation of surging ice bodies, flash floods in dammed lakes caused by surging ice bodies blocking valleys, and ice lake outburst floods caused by surging bodies destroying glacier lakes. Therefore, monitoring glacier surging has important guiding significance for glacier disaster perception. Since glacier surging will cause significant changes in the surface elevation of glaciers, monitoring the change of ice surface elevation is one of the important means to study glacier surging. Spaceborne synthetic aperture radar interferometry (InSAR) has the ability of all-weather observation and can provide key data sources for continuous, timely and comprehensive monitoring of the change of glacier surface elevation. When generating the ice surface DEM, the single-pass dual-receive SAR image pair is not affected by temporal decorrelation on the glacier surface, glacier movement and atmospheric changes. Therefore, after the launch of the single-pass dual-receive image data of TanDEM-X in Germany in 2011, the InSAR technology has been widely applied in the field of monitoring the change of glacier surface elevation.

[0003] Phase unwrapping is the core step for InSAR technology to obtain the surface elevation. Since the original interferometric phase fringes in mountainous areas generated by general SAR images are very dense, the success rate of direct unwrapping is extremely low. It is necessary to subtract the terrain phase simulated by the external DEM from the original phase, and then unwrap the differential phase. The unwrapped differential phase is converted into an elevation difference, and a new DEM can be obtained by adding the external DEM. Traditional phase unwrapping methods are almost all based on the assumption of phase continuity, that is, the phase difference between adjacent pixels does not exceed π. However, when a glacier surges or collapses, the change in surface elevation may cause the phase difference between adjacent pixels on the InSAR differential phase map to exceed π. Using the traditional method for unwrapping will cause obvious jumps in the unwrapped phase, and finally lead to gross errors in the generated InSAR DEM in the key target area, making it impossible to be used to estimate the change of glacier surface elevation.

[0004] Almost all traditional InSAR phase unwrapping methods are based on the assumption of phase continuity, that is, the phase difference between adjacent pixels does not exceed. When the condition of the phase continuity assumption is damaged by the drastic phase change or noise caused by the significant elevation change in the glacier area, the traditional phase unwrapping method cannot obtain the accurate unwrapped phase. In addition, even if the condition of the phase continuity assumption can be satisfied, the path tracking method is prone to form a large number of unwrappable closed regions under low signal-to-noise ratio conditions, resulting in the "island" phenomenon; the optimization-based method is prone to spread the unwrapping error in the region with a large phase gradient to the entire unwrapping region; the statistics-based method has low efficiency in processing interferometric data and requires high computer hardware equipment.

[0005] The existing deep learning phase unwrapping methods still have the following problems: (1) Most of the existing methods directly conduct experiments on simulated data, which makes the model insufficiently consider the spatial distribution characteristics of the real InSAR interferometric phase, and is prone to the problem of unbalanced phase classification under complex terrains; (2) The networks designed by the existing methods are relatively shallow, making it difficult to fit the phase data with more residuals, and almost all the deep learning models adopted use convolutional layers as the basic units, and the model performance is limited by problems such as the convolutional layer being insensitive to the global position of features and difficult to track the long-distance dependence relationships in the image; (3) Since the phase winding count estimation is based on the criterion of minimum norm, it is prone to spread the mispredicted phase gradient information to the entire region.

[0006] Therefore, there is a need in the art for a method for InSAR large-gradient phase unwrapping in glacier areas. Summary of the Invention

[0007] The purpose of the present invention is to provide a method for InSAR large-gradient phase unwrapping in glacier areas to solve the problems that almost all traditional InSAR phase unwrapping methods in the background art are based on the assumption of phase continuity and all have different disadvantages.

[0008] The technical solution of the present invention is a method for InSAR large-gradient phase unwrapping in glacier areas, including the following steps:

[0009] Step 1, generating training samples based on the interferogram simulation technology of the DEM difference value in the glacier area,

[0010] Step 2, training the training samples through a symmetric phase unwrapping network model constructed by mixing a convolutional layer and a Transformer,

[0011] Step 3, inputting the predicted phase discontinuity information as a probability mass map into the maximum flow / minimum cut algorithm for phase winding count estimation to complete phase unwrapping.

[0012] In a specific embodiment, in step 1, the difference value between the SRTM (Shuttle Radar Topography Mission) DEM data and the COP-DEM data of the mountainous area covered by glaciers is selected for training sample simulation.

[0013] In a specific embodiment, the specific process of step 1 includes: unifying the elevation benchmarks of the two-phase DEMs, registering the two-phase DEMs, differentiating the two-phase DEMs, simulating the terrain differential phase diagram based on the DEM differential diagram, then adding simulated atmospheric turbulence noise, Gaussian noise, and large-amplitude noise with random shapes. The large-amplitude noise is used to simulate the severe decorrelation caused by water bodies and mountain shadows, estimating the coherence of the simulated interference phase after adding noise, and finally wrapping the simulated absolute phase to obtain the simulated wrapped phase.

[0014] In a specific embodiment, in step 2, based on the simulated wrapped phase diagram, a residual diagram is calculated according to the "phase continuity hypothesis", the glacier label diagram is cropped according to the position of the simulated phase diagram, and finally the above two factors and the simulated wrapped phase diagram are used as the features input to the three channels of the deep learning model; the phase discontinuity points are calculated according to the simulated absolute phase diagram, and the phase discontinuity point diagram is used as the output of the deep model.

[0015] In a specific embodiment, in step 2, a symmetric phase unwrapping architecture GTPU-Net that combines CNN and Transformer is constructed, including an encoder, a decoder, and skip connections with attention mechanisms; the basic unit of the model encoder is the Vision Transformer.

[0016] In a specific embodiment, step 3 specifically includes: defining an energy function in the Markov random field, using the phase discontinuity information as the prior knowledge of the phase change used in the wrapped count estimation, constructing a Markov directed graph, calculating the weights of each edge of the directed graph according to the energy function, and solving the wrapped count k that minimizes the energy function based on the maximum flow / minimum cut algorithm to complete phase unwrapping.

[0017] The beneficial effects of the present invention include:

[0018] 1. The present invention utilizes the powerful learning and summarization ability and data mining ability of deep learning, getting rid of the dependence of phase unwrapping on the phase continuity hypothesis.

[0019] 2. The present invention uses the DEM difference value in the glacier area to simulate the bistatic InSAR differential phase, and adds various noises according to the phase characteristics in the glacier area, making the obtained InSAR phase training samples more in line with the characteristics of the glacier surface elevation change monitoring scenario.

[0020] 3. The present invention constructs a symmetric phase unwrapping network model that combines a convolutional layer and a Transformer. By means of skip connections, it fuses the local information extracted by the convolutional module and the global information extracted by the Transformer to achieve high-precision segmentation of interference images, solving the problems that the convolutional layer is insensitive to the global position of image features and difficult to track long-range dependence relationships in images.

[0021] 4. The present invention estimates the winding count based on the Markov random energy field and the maximum flow / minimum cut algorithm, solving the problem of the inconsistency between the unwrapped phase obtained by regression prediction based on deep learning and the original interference phase, and at the same time suppressing the propagation of phase gradient information with incorrect predictions in low-coherence regions, making the InSAR large-gradient phase unwrapping result in the glacier area more reliable.

[0022] In addition to the objectives, features, and advantages described above, the present invention has other objectives, features, and advantages. The following will refer to the drawings to further elaborate on the present invention in detail. BRIEF DESCRIPTION OF THE DRAWINGS

[0023] The drawings forming a part of this application are used to provide a further understanding of the present invention. The schematic embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute an improper limitation to the present invention. In the drawings:

[0024] Figure 1 It is a simulation comparison diagram of sample data; where: (a) is the unwrapped phase diagram; (b) is the interference diagram; (c) is the glacier boundary diagram; (d) is the discontinuity points calculated based on the phase continuity assumption; (e) is the coherence diagram; (f) is the discontinuity points calculated based on the unwrapped phase.

[0025] Figure 2 It is a comparison diagram of the interference diagram of a certain glacier in Xinluhai, the unwrapped phase generated by the minimum-cost flow method based on GAMMA software, and the unwrapped phase generated by the present invention; where the marked area 1 is the ice tongue area.

[0026] Figure 3 It is a comparison diagram of the interference diagram of another glacier in Xinluhai, the unwrapped phase generated by the minimum-cost flow method based on GAMMA software, and the unwrapped phase generated by the present invention; where the marked area 2 is the ice tongue area.

[0027] Figure 4 It shows the phase continuity of the InSAR interference diagram of the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0028] The following will elaborate on the embodiments of the present invention in detail with reference to the drawings. The specific embodiments described herein are merely used to explain the present invention and are not used to limit the present invention.

[0029] Embodiment 1

[0030] In order to make the phase sample information more conform to the scenario of differential interferogram in the glacier area and ensure the effectiveness of model training, based on the differential maps of multi-temporal real active glacier area DEMs, the present invention conducts InSAR interferometric phase simulation considering the significant elevation changes on the glacier surface and multi-source noises. The principle and process are as follows:

[0031] Considering that the effective information utilized by the model is phase gradient information, and the spatial resolution of the interferogram after multi-look processing is about 30m, the present invention selects the difference value between SRTM DEM data and COP-DEM data in the mountainous area covered by glaciers for training sample simulation. Using the orbital parameters of the TanDEM-X satellite, the difference value between SRTM DEM data and COP-DEM data is converted into the topographic phase in the radar coordinate system:

[0032]

[0033] where ψ(m, n) is the topographic phase of point (m, n), B ⊥ is the vertical baseline length, λ is the wavelength, R is the satellite orbital height, θ is the satellite incidence angle, and h(m, n) is the DEM difference map after inverse geocoding. Berlin noise is added to the topographic phase with a probability of 30% to simulate the influence of atmospheric turbulence on the interferometric phase, and a large gradient deformation phase simulated by a two-dimensional Gaussian function is added with a probability of 20% to obtain the initial simulated phase.

[0034] Random complex Gaussian noise is added to the initial simulated phase map, and the coherence of the interferogram is estimated based on the added noise, which is used as the phase quality index for traditional algorithm unwrapping. At the same time, considering that the terrain in the glacier area is complex, factors such as water bodies and mountain shadows may cause serious decorrelation, large-amplitude noise with a random shape is added to the simulated phase map with a probability of 10%, and finally the simulated wrapped phase is obtained through phase wrapping.

[0035] The present invention regards the phase unwrapping problem as a semantic segmentation problem. However, different from the existing deep learning phase unwrapping methods that classify images into multiple classes based on the wrap count, the present invention defines the phase discontinuity points as the points where the adjacent phase difference in the unwrapped phase map exceeds, and only performs binary classification on the phase map according to whether the pixel is a phase discontinuity point, reducing the stringent requirements for model performance. Using the interferogram, residual map, and glacier boundary map as input data, and the phase discontinuity point prediction map as output data. At the same time, in order to solve the problems that the convolutional layer is not sensitive to the global position of features and it is difficult to track the long-distance dependence relationship in the image and ensure the performance and stability of the model, the present invention adopts a symmetric phase unwrapping architecture GTPU-Net that combines a CNN and a Transformer, which is composed of an encoder, a decoder, and skip connections with attention mechanisms added.

[0036] To solve the problem of binary classification sample imbalance, the present invention selects a strategy of fusing two loss functions, namely binary cross-entropy loss function (BCELoss) and focal loss function (Focal Loss). At the same time, two model trainings are carried out for the horizontal and vertical discontinuities of pixel points, and the predicted phase discontinuity information is input into the unwrapping algorithm as a probability mass map to estimate the phase wrapping count.

[0037] The present invention estimates the phase wrapping count based on the maximum flow / minimum cut theory and phase discontinuity information to obtain the unwrapped phase. Figure 4 Shows the continuity of the pixel point (i, j) and its first-order neighborhood in the interferogram. Among them, (i, j) ∈ G0, h ij and v ij represent horizontal and vertical discontinuities respectively.

[0038] Define the energy function in the Markov random field:

[0039]

[0040] Among them, k is the wrapping count, ψ is the interference phase, V(.) is the clique potential, and represent the pixel horizontal and vertical differences respectively. At the same time, there are:

[0041]

[0042]

[0043]

[0044]

[0045] The goal of this method is to solve the wrapping count k that minimizes the energy function. Replace h ij and v ij with the phase continuity probabilities in the horizontal and vertical directions predicted by the GTPU-Net model as the prior knowledge of phase change used in the wrapping count estimation. For the convex potential V, the minimization of E(k|ψ) can be achieved through a series of binary optimizations. Each binary problem is mapped to a binary image, and the binary minimization of the objective function can be obtained by calculating the maximum flow / minimum cut on the image.

[0046] The interferogram simulation method in step 1 of the present invention is not a common numerical simulation and DEM inversion, but an interferometric phase simulation method that takes into account the characteristics of glacier surface elevation changes, based on the differential map of multi-period DEMs of real surging glacier areas and fuses multi-source noises.

[0047] The present invention designs a symmetric phase unwrapping architecture that combines a convolutional layer and a Transformer, and adds an attention mechanism to the skip connection, which can solve the problems that the convolutional layer is insensitive to the global position of image features and difficult to track long-distance dependencies in the image. At the same time, in order to enable the model to better learn and summarize the phase discontinuity information in the mountain glacier surge scenario, semantic information closely related to the distribution of phase discontinuity points is selected as the model input, and the prediction of phase discontinuity information is completed based on the encoding and decoding of phase information, thereby ensuring the robustness of the model.

[0048] The present invention does not use the L 1 norm minimum as the criterion for estimating the wrapping count. Instead, a Markov directed graph is constructed based on the interferogram, and the phase discontinuity information predicted by the deep learning model is used as prior information of phase change and input into the maximum flow / minimum cut model. This method can effectively suppress the propagation of the phase gradient information mispredicted by the deep learning model in the estimation of the phase wrapping count. Compared with the wrapping count estimation method based on the L 1 norm, this method is more robust.

[0049] The present invention concerns an InSAR phase training sample simulation technique that takes into account the characteristics of glacier surface elevation changes. Most of the existing interferogram simulation methods use simple numerical simulations or DEM back-calculations, lacking sufficient consideration of the spatial distribution characteristics of real InSAR differential interferometric phases. The present invention proposes an interferogram simulation method based on the DEM difference values in the glacier area, simulating the large-gradient phases caused by significant elevation changes in the glacier area with multi-temporal DEM difference values, and then fusing multi-source noises to fully restore the spatial distribution characteristics of InSAR differential interferometric phases in the real glacier surge scenario, improving the quality of training samples and ensuring the effectiveness of model training.

[0050] The present invention concerns a network model construction technique that takes into account the global features of InSAR interferograms. The performance of the network models used in existing deep learning phase unwrapping methods is restricted by problems such as the insensitivity of convolutional layers to the global position of features and the difficulty in tracking long-distance dependencies in images. The present invention constructs a symmetric phase unwrapping network model that combines a convolutional layer and a Transformer, and fuses the local information extracted by the convolutional module and the global information extracted by the Transformer through skip connections to achieve high-precision segmentation of interferometric images.

[0051] The InSAR phase winding number estimation technology based on Markov random energy field and maximum flow / minimum cut algorithm. In the existing phase winding number estimation methods, the wrongly predicted winding number gradient information is easily propagated to the entire target area. The present invention designs an energy function considering the prior distribution information of phase discontinuity based on Markov random field, transforms the interferogram into a directed graph including a source point and a sink point, assigns weights to the edges of the directed graph by using the energy function, and optimizes the winding number estimation through the maximum flow / minimum cut algorithm, effectively suppressing the propagation of the wrongly predicted phase gradient information.

[0052] Figure 1 Shows the simulated data generated by the InSAR phase training sample simulation method considering the characteristics of glacier surface elevation change proposed by the present invention. Among them, the interferogram, glacier boundary and residual map are the input data of the deep learning model, and the phase discontinuity point prediction map is used as the output data. Subsequently, the discontinuity point prediction map is used as the InSAR quality map and input into the maximum flow / minimum cut algorithm for winding count estimation, and finally phase unwrapping is completed.

[0053] Figure 2 、 Figure 3 Shows a case of Xinluhai Glacier monitoring. The ice tongue area is the most active section of glacier action and also the ablation area of the glacier. In the unwrapping map generated by GAMMA software in the figure, the ice tongue areas 1 and 2 marked are not correctly unwrapped. However, using the unwrapping method of the present invention and combining with the actual morphology of Xinluhai Glacier, it can be seen that part of the glacier has been correctly unwrapped.

[0054] The above content is a further detailed description of the present invention in combination with specific preferred embodiments. It cannot be determined that the specific implementation of the present invention is only limited to these descriptions. For those of ordinary skill in the technical field to which the present invention belongs, without departing from the concept of the present invention, several simple deductions and substitutions can still be made, which should all be regarded as belonging to the protection scope of the present invention.

Claims

1. A method for unwrapping large gradient phases in InSAR in glacial regions, characterized in that, It includes the following steps: Step 1: Generate training samples by using the interferogram simulation technology based on the differential value of the glacier area DEM. In Step 1, the differential value of the SRTM (Shuttle Radar Topography Mission) DEM data and the COP-DEM data in the mountainous area covered by glaciers is selected for training sample simulation. The specific process of Step 1 includes: unifying the elevation benchmarks of the two-phase DEMs, registering the two-phase DEMs, differentiating the two-phase DEMs, simulating the terrain differential phase diagram based on the DEM differential map, then adding simulated atmospheric turbulence noise, Gaussian noise, and large-amplitude noise with random shapes. The large-amplitude noise is used to simulate the severe decorrelation caused by water bodies and mountain shadows. Estimate the coherence of the simulated interference phase after adding noise. Finally, wrap the simulated absolute phase to obtain the simulated wrapped phase diagram. Step 2: Train the training samples through a symmetric phase unwrapping network model constructed by mixing a convolutional layer and a Transformer. In Step 2, a symmetric phase unwrapping architecture GTPU-Net mixing CNN and Transformer is constructed, including an encoder, a decoder, and a skip connection with an attention mechanism added. The basic unit of the model encoder is the Vision Transformer. Calculate the residual map based on the simulated wrapped phase diagram according to the "phase continuity hypothesis", crop the glacier boundary according to the position of the simulated wrapped phase diagram, and finally use the residual map, glacier boundary, and interferogram as the features input to the three channels of the deep learning model. Calculate the phase discontinuity points based on the simulated absolute phase diagram, and use the phase discontinuity point map as the output of the deep learning model. Step 3: Input the predicted phase discontinuity information as a probability mass map into the maximum flow / minimum cut algorithm to estimate the phase wrapping count and complete phase unwrapping. The specific content of Step 3 includes: define the energy function in the Markov random field, use the phase discontinuity information as the prior knowledge of the phase change used in the wrapping count estimation, construct a Markov directed graph, calculate the weights of each edge of the directed graph according to the energy function, and solve the wrapping count k that minimizes the energy function based on the maximum flow / minimum cut algorithm to complete phase unwrapping.

Images