InSAR Phase Unwrapping Method Based on Terrain Gradient Sensing Network
By constructing a terrain gradient sensing network and combining gradient adaptive decoding and neural embedded post-processing, the problem of unifying noise suppression and detail preservation in InSAR phase unwrapping is solved, improving unwrapping accuracy and robustness, and adapting to complex terrain scenarios.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- SHIJIAZHUANG TIEDAO UNIV
- Filing Date
- 2026-02-02
- Publication Date
- 2026-07-31
AI Technical Summary
Existing InSAR phase unwrapping techniques are insufficient in terms of unwrapping accuracy and robustness in complex terrain gradient changes, noise interference, and dense phase jump scenes. They are difficult to balance noise suppression and detail preservation, and have low automation and poor scene adaptability.
An InSAR phase unwrapping method based on a terrain gradient sensing network is adopted. By constructing a terrain gradient sensing coding module, a gradient-guided dynamic weight multi-scale hollow pooling module, a gradient adaptive decoding module, and a neural embedded post-processing module, noise suppression and detail preservation are unified, thereby improving unwrapping accuracy and robustness.
This method addresses the issues of insufficient adaptation to complex terrain and lack of detail recovery by using end-to-end terrain gradient guidance and physical constraint optimization, thereby improving high accuracy, robustness, and scene adaptability.
Smart Images

Figure CN122017835B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of synthetic aperture radar interferometry, and more particularly to an InSAR phase unwrapping method based on a terrain gradient sensing network. Background Technology
[0002] Synthetic Aperture Radar Interferometry (SAR) technology, with its advantages of all-weather, all-time capability and high spatial resolution, has become a core method for monitoring surface deformation and elevation inversion. Phase unwrapping, as a crucial step in data processing, directly determines the reliability of the final surface data in terms of its accuracy and scene adaptability. The core objective of phase unwrapping is to recover continuous phase information reflecting the true surface morphology and deformation from the wrapped phase within the range of [-π, π]. However, in practical applications, problems such as complex terrain gradient changes, noise interference, dense phase jumps, and uneven distribution of incoherent regions consistently constrain the improvement of unwrapping accuracy and robustness.
[0003] Current phase unwrapping techniques are mainly divided into two categories: traditional methods and deep learning methods. Traditional phase unwrapping methods have inherent drawbacks: path-following methods, minimum norm methods, and network flow methods all rely on manually designed physical models or mathematical rules, lacking adaptive adjustment capabilities. In low signal-to-noise ratio and dense residual point scenarios, they are easily affected by noise and measurement errors, leading to the accumulation of unwrapping errors. When facing steep terrain and highly decoherent regions, it is difficult to balance noise suppression and detail preservation, often resulting in unwrapped phase breaks, over-smoothing, or the global propagation of local errors. Furthermore, they suffer from low automation, and it is difficult to balance computational efficiency with scene adaptability.
[0004] Existing deep learning methods suffer from significant optimization bottlenecks. First, they lack adaptability to complex terrain. Most mainstream models employ fixed-structure multi-scale modules, lacking terrain gradient awareness mechanisms and failing to dynamically adjust the receptive field based on terrain steepness. This results in overly smooth details in flat areas and insufficient context capture in steep regions. Second, edge detail recovery lacks specificity. The decoder upsampling process does not effectively utilize phase gradient information, leading to inaccurate phase transition localization at abrupt terrain changes, resulting in blurred edges and discontinuous stripe breaks. Third, physical plausibility is disconnected from the processing flow. Some methods rely on traditional unwrapping frameworks or offline post-processing, lacking dedicated physical logic for phase unwrapping, and their loss functions are limited, leading to numerically accurate but physically distorted unwrapping results. Furthermore, some models rely on additional auxiliary data or large-scale labeled data, limiting their generalization ability. Summary of the Invention
[0005] The technical problem to be solved by this invention is to provide an InSAR phase unwrapping method that can achieve a balance between noise suppression and detail preservation, and improve the accuracy, robustness and scene adaptability of phase unwrapping.
[0006] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is: an InSAR phase unwrapping method based on a terrain gradient sensing network, the method comprising:
[0007] S1: Obtain the phase unwrapping dataset and divide the phase unwrapping dataset into a training dataset, a validation dataset, and a test dataset;
[0008] S2: Construct a terrain gradient sensing network, including a data loading and preprocessing module, a phase map terrain gradient change sensing encoding module, a gradient-guided dynamic weight multi-scale hollow pooling module, a gradient adaptive decoding module, and a neural embedded post-processing module.
[0009] S3: Input the training dataset into the terrain gradient perception network for iterative training to generate a trained terrain gradient perception network. The validation dataset is used to monitor the training effect in real time. The loss value is calculated through a multi-constraint loss function, and the model parameters are iteratively optimized through backpropagation. After training is completed, the optimal model weights are saved.
[0010] S4: Input the test dataset into the trained terrain gradient perception network and output the predicted unwrapped phase map.
[0011] A further technical solution is as follows: the phase map terrain gradient change sensing encoding module includes an encoder branch and a terrain gradient sensing module branch, which are connected in parallel with the input layer and processed synchronously. The encoder branch consists of a first residual unit, a first max pooling layer, a second residual unit, a second max pooling layer, a third residual unit, a third max pooling layer, a fourth residual unit, a fourth max pooling layer, a fifth residual unit, a fifth max pooling layer, and a sixth residual unit connected in series; each residual unit consists of two convolutional layers, batch normalization, and an activation function, and the input and convolutional output are added through residual connections; the gradient sensing branch consists of a horizontal gradient convolutional layer, a vertical gradient convolutional layer, a normalization and shaping unit, a gradient distribution map generation unit, and a dynamic weight generation unit; the gradient distribution map is used for spatial guidance in the decoding stage, and the dynamic weights are used for multi-scale response modulation.
[0012] A further technical solution is as follows: the gradient-guided dynamic weight multi-scale dilated pooling module includes a first dilated convolution branch with an inflation rate of 1, a second dilated convolution branch with an inflation rate of 3, a third dilated convolution branch with an inflation rate of 6, a fourth dilated convolution branch with an inflation rate of 12, and a global context branch. The outputs of the five branches are concatenated in the channel dimension and then fused with a 1×1 convolution layer, batch normalization, and activation function to output bottleneck enhancement features. Among them, the dynamic weights establish weight associations with the four dilated convolution branches, and adaptively modulate the responses of each branch to enhance the contextual expression of dense stripes and phase change regions, while maintaining details and stability in flat regions.
[0013] A further technical solution is that the gradient adaptive decoding module includes a first gradient adaptive upsampling unit, a second gradient adaptive upsampling unit, a third gradient adaptive upsampling unit, a fourth gradient adaptive upsampling unit, and a fifth gradient adaptive upsampling unit connected in series and respectively connected to the first to fifth connection layers; each gradient adaptive upsampling unit includes a transposed convolutional upsampling layer with a kernel size of 2×2 and a stride of 2; a gradient weight generation layer receives a gradient distribution map matching the current scale and generates spatial attention weights through a 1×1 convolution and a sigmoid function; an enhancement convolution includes two 3×3 convolutions and batch normalization; an activation function and a residual fusion layer add the enhancement branch to the main upsampling branch after weighting it by attention; the upsampling output of each level is concatenated with the corresponding connection layer output in the channel dimension, refined by the convolutional residual fusion layer, and then sent to the next upsampling unit, so that the resolution is restored from 8×8 to 256×256 sequentially, and key detail restoration is performed on the high gradient region throughout the entire process.
[0014] A further technical solution involves a neural embedded post-processing module comprising two parts: residual correlation correction and continuity optimization. The residual correlation correction part, based on the adjacent differences of the phase, constructs continuous residual indicator features through a sine map. It then generates a residual attention map using a shallow convolutional network, which is concatenated with the phase map. Another convolutional network predicts the correction amount, applying targeted correction only in high-confidence residual regions to avoid excessive intervention in normal regions. The continuity optimization part uses a fixed 3×3 Laplacian operator to calculate discontinuity responses, generating adaptive smoothing weights accordingly. Simultaneously, a small convolutional network predicts the smoothing correction amount. The combination of these two methods weighted refinement of the phase: enhancing global continuity in flat regions and preserving true transition details in abrupt regions. Finally, a 1×1 convolution outputs an unwrapped phase map.
[0015] The beneficial effects of adopting the above technical solution are as follows: The method is based on an encoder-decoder architecture, constructing a terrain gradient-aware encoding module and a gradient adaptive decoding module. It extracts terrain gradient information through an improved Sobel operator, generating dual-path guiding signals to drive dynamic weight ASPP multi-scale feature fusion and differential upsampling of the decoder, respectively. Simultaneously, it embeds a neural embedded differentiable post-processing module, transforming residual correction and phase continuity constraints into differentiable units for end-to-end collaborative training with the network. This method, through end-to-end terrain gradient guidance and physical constraint optimization, solves problems such as insufficient adaptation to complex terrain, lack of detail recovery, and physical distortion, achieving a balance between noise suppression and detail preservation, and improving the accuracy, robustness, and scene adaptability of phase unwrapping. Attached Figure Description
[0016] The present invention will now be described in further detail with reference to the accompanying drawings and specific embodiments.
[0017] Figure 1 This is a flowchart illustrating the method described in an embodiment of the present invention;
[0018] Figure 2 This refers to the partial training set wrapped phase data in the method described in this embodiment of the invention;
[0019] Figure 3 This refers to a portion of the training set unwrapped phase data in the method described in this embodiment of the invention.
[0020] Figure 4 This is a schematic diagram of the terrain gradient sensing network in the method described in the embodiments of the present invention;
[0021] Figure 5 This is a schematic diagram of the terrain gradient sensing module in the method described in this embodiment of the invention;
[0022] Figure 6 This is a schematic diagram of the neural embedded post-processing module in the method described in the embodiment of the present invention;
[0023] Figure 7 This refers to the winding phase test data in the method described in the embodiments of the present invention;
[0024] Figure 8 This refers to the unwrapping phase test data in the method described in the embodiments of the present invention. Detailed Implementation
[0025] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.
[0026] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and those skilled in the art can make similar extensions without departing from the spirit of the invention. Therefore, the invention is not limited to the specific embodiments disclosed below.
[0027] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be described in detail below with reference to the accompanying drawings and specific embodiments.
[0028] A schematic diagram of an InSAR phase unwrapping method based on a terrain gradient sensing network is shown in this embodiment. Figure 1 As shown, the following steps may be included:
[0029] S1: Obtain the phase unwrapping dataset and divide the phase unwrapping dataset into a training dataset, a validation dataset, and a test dataset.
[0030] Specifically, the phase unwrapping dataset used is the phase unwrapping dataset proposed by Zhou Lifan's team. This dataset is based on real terrain data of 30 meters / pixel across a certain area, covering various landform types, and includes a training set of 30,000 simulated phase pairs, a validation set of 1,000 simulated phase pairs, and a test set of 100 real InSAR phase pairs. Part of the training set unwrapped phase data in this embodiment of the invention is shown below. Figure 2 The training set unwrapping phase data shown in the figure and corresponding to each figure are as follows: Figure 3 As shown, they are used for core model training, standard simulated scenario untangling performance verification, and real data generalization ability evaluation, respectively.
[0031] S2: Construct a terrain gradient sensing network, which includes a data loading and preprocessing module, a phase map terrain gradient change sensing encoding module, a gradient-guided dynamic weight multi-scale hollow pooling module, a gradient adaptive decoding module, a neural embedded post-processing module, and an output.
[0032] Specifically, a schematic diagram of the terrain gradient sensing network structure in an embodiment of the present invention is shown below. Figure 4 As shown. The specific structure and data processing flow of each part are as follows:
[0033] To ensure the consistency and effectiveness of model input, the data loading and preprocessing module performs standardization and format conversion operations on the input wrapped phase map. This standardizes the numerical range of the input phase map to the [-π, π] interval, eliminating scale differences between different data sources; resizes the phase map to 256×256 pixels to ensure consistent data dimensions during batch training or prediction, adapting to model input requirements; converts the single-channel phase map to tensor format with dimensions [B, 1, 256, 256], where B is the batch size; and performs enhancement operations such as random flipping and rotation on the training dataset to improve the model's generalization ability and avoid overfitting.
[0034] The phase map terrain gradient change sensing encoding module receives preprocessed data, extracts multi-scale phase features, and generates terrain gradient guidance information. It consists of an encoder with six encoding layers and a terrain gradient sensing module forming two parallel branches: the encoder branch increases the number of feature channels layer by layer from 1 to 256, while the spatial resolution halves layer by layer from 256×256 to 8×8. Each encoding layer of the encoder uses two 3×3 convolutions as its core, with a stride of 1 to avoid losing boundary information. Batch normalization (BN) and ReLU activation are applied to accelerate training convergence and enhance nonlinear representation. Subsequently, the input x is added to the convolutionally transformed features through residual connections, forming y = x + F(x), where:
[0035] ;
[0036] and These are the first and second 3×3 convolutions, respectively. This is a batch normalization operation. The parallel terrain gradient sensing module calculates and normalizes the horizontal and vertical gradients of the input phase map, generating a gradient distribution map and dynamic weights. A schematic diagram of the terrain gradient sensing module in this embodiment of the invention is shown below. Figure 5 As shown. This unit is deployed in parallel with 6 coding layers, synchronously receiving the preprocessed wrapped phase map with dimensions [B, 1, 256, 256], where B is the batch size. The feature dimensions of the phase map itself are not changed throughout the process; only gradient features are extracted as the global guidance basis. It receives the original input phase map and generates gradient guidance information including a gradient distribution map and a dynamic dilation rate weight vector, providing prior information on the steepness of the terrain for subsequent modules. This includes four steps:
[0037] The first step is gradient calculation. Using the Sobel operator, adapted to the periodic characteristics of the phase map, the first-order partial derivatives of the input phase map are calculated along the x-axis and y-axis respectively, yielding the horizontal and vertical gradients. The calculation formulas are as follows:
[0038] ;
[0039] in For horizontal gradient, For vertical gradient, and A 3×3 adaptive Sobel convolution kernel is used; based on this, the gradient magnitude representing the steepness of the terrain and the gradient direction representing the terrain orientation are further calculated, and the calculation formula is as follows:
[0040] ;
[0041] in For horizontal gradient, The first step is to output a gradient magnitude map G and a gradient direction map θ, both with dimensions [B, 1, 256, 256], representing the vertical gradient. The second step is gradient normalization. To eliminate the scale differences in gradient values between different batches and different terrain phase maps, the maximum and minimum values of the gradient magnitude map are calculated for each batch, and a normalization operation is performed on each pixel. The calculation formula is as follows:
[0042] ;
[0043] in, This represents the maximum value of the gradient magnitude plot. To find the minimum value of the gradient magnitude map, output a normalized gradient magnitude map. And the gradient direction pattern θ that retains the original values.
[0044] The third step provides a basis for the gradient guidance requirements of the decoding layer. First, a 1×1 lightweight convolution operation is performed on the normalized gradient magnitude map to keep the number of channels at 1, generating a pixel-by-pixel attention weight map. Adaptive scaling is then applied to the feature responses at different spatial locations using the attention weight map, giving higher weights to high gradient and structurally abrupt regions to enhance edge and detail representation, and lower weights to low gradient and flat regions to suppress noise perturbations. Finally, a gradient distribution map with dimensions [B, 1, 256, 256] is generated, and its calculation formula is as follows:
[0045] ;
[0046] in Represents a 1×1 lightweight convolution operation. This represents a gradient-guided spatial weighting operation on the normalized gradient magnitude map. This is the gradient distribution map. The gradient distribution map is directly passed to each upsampling branch of the decoding layer, providing a basis for the gradient guidance and enhancement stage in the decoding phase.
[0047] The fourth step is the generation of dynamic expansion rate weights. This step aims to provide a weight basis for gradient adaptation for multi-scale feature fusion in the Hollow Spatial Pyramid Pooling (ASPP) module. Specifically, it first generates weights for the normalized gradient magnitude map. A global average pooling operation is performed to obtain the batch-level average gradient value. Then, based on this average gradient value, adaptive weights are assigned to the four different expansion rate branches of the Hollow Space Pyramid Pooling (ASPP) module. These weights are encapsulated as learnable parameters, ultimately outputting a dynamic expansion rate weight vector of dimension 4. The calculation formula is as follows:
[0048] ;
[0049] in This represents the global average pooling operation. This represents a branch weight assignment operation based on the average gradient value. This is a dynamic expansion rate weight vector, which will be directly passed to the Hollow Spatial Pyramid Pooling (ASPP) module for weighted fusion of features from different expansion rate branches, thereby achieving dynamic adaptation of multi-scale features to terrain gradients.
[0050] The gradient-guided dynamic weight multi-scale dilated pooling module receives features from the last layer of the encoder and dynamic dilation rate weights. It is responsible for dynamically fusing multi-scale contextual features, overcoming the adaptation limitations of the fixed multi-scale module. It includes four parallel basic branches with base dilation rates of 1, 3, 6, and 12, and one global average pooling branch. First, each of the parallel basic branches consists of a 3×3 convolution corresponding to the base dilation rate, a batch normalization (BN) layer, and a ReLU activation, extracting features at a specific scale. Based on the dynamic dilation rate weights, the output features of each basic branch are differentially weighted: high-gradient regions enhance the response of branches with large dilation rates, while flat regions enhance the response of branches with small dilation rates. The calculation formula is as follows:
[0051] ;
[0052] in, This is the original output of the i-th basic branch. For dynamically enhanced features, an enhancement factor of 0.5 is used to balance the contributions of each branch. Secondly, the global average pooling branch, composed of adaptive average pooling output size, 1×1 convolution, BN layer, and ReLU activation, is responsible for extracting global contextual information. The output features of the global average pooling branch are upsampled to an 8×8 size to maintain dimensionality consistent with the parallel base branch features. Finally, the four dynamically enhanced branches are concatenated with the global features along the channel dimension and fused into a 256-channel feature through a 1×1 convolution. The calculation formula is as follows:
[0053] ;
[0054] in, It is a fused feature map output by the void space pyramid pooling module, which is a feature map output by the coded residual block. After being processed by dilated convolutions with four different dilation rates, the dynamic weights generated by the terrain gradient sensing unit are... The weighted summation can be further optimized by 1×1 convolution and batch normalization as needed to unify the number of channels and improve feature stability.
[0055] The gradient adaptive decoding module is responsible for restoring spatial resolution and fusing encoder features. It consists of 5 gradient adaptive upsampling units and 5 transposed convolutional layers. The number of feature channels decreases layer by layer from 256 to 1, and the spatial resolution is restored layer by layer from 8×8 to 256×256. The adaptive upsampling unit introduces gradient-guided detail enhancement to avoid the edge blurring problem of simple transposed convolution. It first performs standard upsampling on the input features, and then doubles the spatial resolution through a 2×2 transposed convolution with a stride of 2. The calculation formula is as follows:
[0056] ;
[0057] in, Input features for the module, Features after standard upsampling The kernel size is [size]. The step size is set; then, the gradient distribution map output by the terrain gradient sensing unit is... Upsampling to Generate spatial attention weights for the same size For targeted enhancement of high gradient edge regions, the calculation formula is as follows:
[0058] ;
[0059] in, This is the gradient distribution map after upsampling. For channel splicing, A 3×3 convolutional layer is used to extract fused features; the enhanced features are then fed into an upsampled residual block for further refinement. In the upsampled residual block and skip-connection fusion stage, the output features are combined with the skip-connection features at the corresponding scale from the encoder. The data is concatenated along the channel dimension and fused using two 3×3 convolutional layers. Gradient propagation and feature reuse are stabilized using a residual approach. Simultaneously, the spatial attention weights obtained in the previous step are combined with the data. Before fusion, the splicing features are slightly adjusted to enhance the information transmission of high-gradient edges. The calculation formula is as follows:
[0060] ;
[0061] in, For the features of the corresponding layer of the encoder, This is for channel splicing operations. This is the final output of the decoding layer and serves as the input for the next decoding level.
[0062] The neural embedded post-processing module and output perform final refinement on the phase feature map of the last layer output of the decoder. It includes a differentiable branch cutting unit and a differentiable least squares unit, which are connected in series. A schematic diagram of the neural embedded post-processing module of this embodiment is shown below. Figure 6 As shown in the diagram, the differentiable branch-cutting unit simulates the residual correction logic of the traditional branch-cutting method, achieving accurate detection and targeted correction of residual points. First, the phase difference ΔI between adjacent pixels in the horizontal and vertical directions of the phase map is calculated. Residual points are identified using the numerical features of |sin(ΔI)|. Then, a shallow convolutional neural network is used to apply attention weights to the residual point region, generating a residual attention map. Attention graph of residuals After concatenating the input phase map I along the channel dimension, it is input into the CNN to generate a targeted correction that only applies to the residual region, avoiding interference with the normal phase region. The calculation formula is as follows:
[0063] ;
[0064] in This is the phase diagram after residual correction. A larger value indicates a higher probability that the corresponding pixel is a residual point. Differentiable least squares units are used to simulate the phase continuity constraint logic of traditional least squares methods, completing phase optimization while balancing global phase continuity and local detail preservation. A fixed-weight Laplacian operator is used, and its calculation formula is:
[0065] ;
[0066] Phase diagram after residual correction Convolution operations are performed to obtain the phase discontinuity response L. A larger value for L indicates a stronger phase discontinuity in the corresponding region. Adaptive smoothing weights are generated based on the discontinuity response L; for flat regions, a smaller value for L results in a smoother response. A value close to 1 strengthens the global continuity of the phase; for abrupt change regions, a large L value... Approaching 0, preserving local details, its calculation formula is:
[0067] ;
[0068] in, For adaptive smoothing weights, This is the phase discontinuity response. The residual-corrected phase map is input into a convolutional network to generate a smoothing correction. Combined with smoothing weights, phase continuity optimization is completed, and the final output is a globally continuous phase map with complete local details. Its calculation formula is as follows:
[0069] ;
[0070] in, The optimized phase map is the final output of the module. This is the phase diagram after residual correction. To achieve adaptive smoothing weights, the dual requirements of global phase continuity and local detail preservation can be met simultaneously, resulting in an output unwrapped phase tensor with a size of 1×256×256.
[0071] Step S3: Input the training dataset into the terrain gradient-aware deep learning phase unwrapping network for iterative training. Verify the dataset for real-time monitoring of training performance. Calculate the loss value using a multi-constraint loss function, iteratively optimize the model parameters via backpropagation, and save the optimal model weights after training. Specifically, the steps are as follows:
[0072] (1) Set the initial learning rate, maximum learning rate, number of training batches, and number of training rounds; the initial learning rate is 0.0001, the maximum learning rate is 0.01, the number of training batches is 4, and the number of training rounds is 100.
[0073] (2) During training, batch normalization and early stopping strategies are used to prevent network overfitting.
[0074] (3) The Adam gradient optimization algorithm is used to optimize the training of the terrain gradient perception network structure. During the optimization process, a composite function combining pixel-level precision loss and gradient consistency loss is used as the loss function.
[0075] The loss function, which integrates pixel-level accuracy and gradient structure consistency, guides the model to simultaneously optimize numerical accuracy and physical plausibility. Its calculation formula is as follows:
[0076] ;
[0077] Where 0.3 is the weight coefficient of the gradient loss. For pixel-level loss, For gradient consistency loss, the balance between the two constraints is verified through experiments. Pixel-level loss is calculated by taking the mean absolute error of the predicted phase versus the true phase pixel-by-pixel to ensure numerical accuracy. The calculation formula is as follows:
[0078]
[0079] Where B is the batch size, and H×W=256×256 is the image size. To predict the untangling phase, For the true untangling phase, This is a pixel-level loss. The gradient consistency loss constrains the gradient of the predicted phase to be consistent with the true phase gradient, avoiding edge blurring and physical distortion. Its calculation formula is:
[0080]
[0081] in, , These are the gradient operators in the x and y directions, respectively. This is the gradient consistency loss.
[0082] (4) Repeat steps (2) and (3) to perform iterative training until the terrain gradient sensing network structure converges, obtain the final terrain gradient sensing network model and weights for interferometric phase unwrapping, and save the optimal model weights.
[0083] Step S4: Input the test dataset into the trained terrain gradient perception network and output the predicted unwrapped phase map.
[0084] Specifically, the optimal model weights saved during the training phase are loaded into the terrain gradient perception network to complete initialization; then, the entangled phase data of the test set is input into the network for inference, wherein the entangled phase test data in this embodiment of the invention is as follows: Figure 7 As shown; after forward computation via the network, the corresponding interference phase unwrapping prediction result is directly obtained. The unwrapping phase test data of this embodiment of the invention is as follows: Figure 8 As shown.
[0085] It is worth noting that, for the phase map terrain gradient change sensing coding module, this invention designs an encoder branch and a terrain gradient sensing module branch, which are connected in parallel with the input layer and processed synchronously. The encoder branch consists of a first residual unit, a first max pooling layer, a second residual unit, a second max pooling layer, a third residual unit, a third max pooling layer, a fourth residual unit, a fourth max pooling layer, a fifth residual unit, a fifth max pooling layer, and a sixth residual unit connected in series. Each residual unit consists of two convolutional layers, batch normalization, and an activation function, and the input and convolutional output are added through residual connections. The terrain gradient sensing module branch consists of a horizontal gradient convolutional layer, a vertical gradient convolutional layer, a normalization and shaping unit, a gradient distribution map generation unit, and a dynamic weight generation unit. The gradient distribution map is used for spatial guidance in the decoding stage, and the dynamic weights are used for multi-scale response modulation.
[0086] For the gradient-guided dynamic weight multi-scale dilated pooling module, this invention designs a first dilated convolution branch with an inflation rate of 1, a second dilated convolution branch with an inflation rate of 3, a third dilated convolution branch with an inflation rate of 6, a fourth dilated convolution branch with an inflation rate of 12, and a global context branch. The outputs of the five branches are concatenated in the channel dimension and then fused with a 1×1 convolution layer, batch normalization, and activation function to output bottleneck enhancement features. Among them, the dynamic weights are associated with the four dilated convolution branches, and the responses of each branch are adaptively modulated to enhance the contextual representation of dense stripes and phase change regions, while maintaining details and stability in flat regions.
[0087] For the gradient adaptive decoding module, this invention designs a first gradient adaptive upsampling unit, a second gradient adaptive upsampling unit, a third gradient adaptive upsampling unit, a fourth gradient adaptive upsampling unit, and a fifth gradient adaptive upsampling unit connected in series and respectively connected to the first to fifth connection layers. Each gradient adaptive upsampling unit includes a transposed convolutional upsampling layer with a kernel size of 2×2 and a stride of 2. A gradient weight generation layer receives a gradient distribution map matching the current scale and generates spatial attention weights through a 1×1 convolution and a sigmoid function. An enhancement convolution includes two 3×3 convolutions and batch normalization. An activation function and a residual fusion layer add the enhancement branch to the main upsampling branch after weighting it by attention. The upsampling output of each level is concatenated with the corresponding connection layer output in the channel dimension, refined by the convolutional residual fusion layer, and then sent to the next upsampling unit, so that the resolution is restored from 8×8 to 256×256 sequentially. Key detail restoration is performed on high gradient regions throughout the entire process.
[0088] For the neural embedded post-processing module, this invention designs two parts: residual correlation correction and continuity optimization. The residual correlation correction part is based on the adjacent differences of the phase, constructing continuous residual indicator features through a sine map. A shallow convolutional network generates a residual attention map, which is then concatenated with the phase map. Another convolutional network predicts the correction amount, applying targeted correction only in high-confidence residual regions to avoid over-intervention in normal regions. The continuity optimization part uses a fixed 3×3 Laplacian operator to calculate discontinuity responses, generating adaptive smoothing weights accordingly. Simultaneously, a small convolutional network predicts the smoothing correction amount. The combination of these two parts performs weighted refinement of the phase: enhancing global continuity in flat regions and preserving true transition details in abrupt regions. Finally, a 1×1 convolution outputs an unwrapped phase map.
[0089] Specific embodiments of the present invention have been described above. It should be understood that the present invention is not limited to the specific embodiments described above, and those skilled in the art can make various changes or modifications within the scope of the claims, which do not affect the essence of the present invention. Unless otherwise specified, the embodiments and features described in this application can be arbitrarily combined with each other.
Claims
1. A method for InSAR phase unwrapping based on terrain gradient perception network, characterized in that, The method includes: S1: Obtain the phase unwrapping dataset and divide the phase unwrapping dataset into a training dataset, a validation dataset, and a test dataset; S2: Construct a terrain gradient sensing network, which includes a data loading and preprocessing module, a phase map terrain gradient change sensing encoding module, a gradient-guided dynamic weight multi-scale hollow pooling module, a gradient adaptive decoding module, and a neural embedded post-processing module. The phase map terrain gradient change sensing encoding module includes an encoder branch and a terrain gradient sensing module branch, which are connected in parallel with the input layer and processed synchronously. The terrain gradient sensing module branch includes a horizontal gradient convolutional layer, a vertical gradient convolutional layer, a normalization and shaping unit, a gradient distribution map generation unit, and a dynamic weight generation unit connected in series. The gradient distribution map is used for spatial guidance in the decoding stage, and the dynamic weight is used for multi-scale response modulation. The gradient-guided dynamic weight multi-scale dilated pooling module includes a first dilated convolution branch with an inflation rate of 1, a second dilated convolution branch with an inflation rate of 3, a third dilated convolution branch with an inflation rate of 6, a fourth dilated convolution branch with an inflation rate of 12, and a global context branch. The outputs of the five branches are concatenated along the channel dimension and then fused with a 1×1 convolution layer, batch normalization, and activation function to output bottleneck enhancement features. The dynamic weights are associated with the four dilated convolution branches, and the responses of each branch are adaptively modulated. The gradient adaptive decoding module is responsible for restoring spatial resolution and fusing encoder features, and contains a total of 5 gradient adaptive upsampling units and 5 transposed convolutional layers; The neural embedded post-processing module includes a residual correlation correction part and a continuity optimization part. The residual correlation correction part is based on the adjacent differences of the phase, constructs continuous residual indicator features through sine mapping, generates a residual attention map with the help of a shallow convolutional network, and then concatenates it with the phase map before another convolutional network predicts the correction amount. The continuity optimization part uses a fixed 3×3 Laplacian operator to calculate the discontinuity response and generates adaptive smoothing weights accordingly. At the same time, a small convolutional network predicts the smoothing correction amount, and the two are combined to perform weighted refinement of the phase. Finally, an unwrapped phase map is output through a 1×1 convolution. S3: Input the training dataset into the terrain gradient perception network for iterative training to generate a trained terrain gradient perception network. The validation dataset is used to monitor the training effect in real time. The loss value is calculated through a multi-constraint loss function, and the model parameters are iteratively optimized through backpropagation. After training is completed, the optimal model weights are saved. S4: Input the test dataset divided from the phase unwrapped dataset into the trained terrain gradient perception network, and output the predicted unwrapped phase map.
2. The terrain gradient-aware network-based InSAR phase unwrapping method of claim 1, wherein, In step S1: The phase unwrapping dataset includes a training set of 30,000 simulated phase pairs, a validation set of 1,000 simulated phase pairs, and a test set of 100 real InSAR phase pairs.
3. The InSAR phase unwrapping method based on terrain gradient sensing network as described in claim 1, characterized in that: The data loading and preprocessing module is used to standardize and convert the input wrapped phase map, standardizing the numerical range of the input phase map to the interval [-π, π]; resizing the phase map to 256×256 pixels to ensure consistent data dimensions during batch training or prediction, adapting to the model input requirements; converting the single-channel phase map to tensor format with dimensions [B, 1, 256, 256], where B is the batch size; and performing augmentation operations on the training dataset to improve the model's generalization ability.
4. The method of claim 1, wherein the terrain gradient-aware network is a convolutional neural network (CNN). The encoder branch sequentially includes a first residual unit, a first max pooling layer, a second residual unit, a second max pooling layer, a third residual unit, a third max pooling layer, a fourth residual unit, a fourth max pooling layer, a fifth residual unit, a fifth max pooling layer, and a sixth residual unit connected in series. Each residual unit includes two convolutional layers, batch normalization, and an activation function, and the input and convolutional output are added together through residual connections.
5. The InSAR phase unwrapping method based on terrain gradient sensing network as described in claim 1, characterized in that: The phase map terrain gradient change sensing encoding module receives preprocessed data, extracts multi-scale phase features, and generates terrain gradient guidance information. It consists of an encoder with six encoding layers and a terrain gradient sensing module forming two parallel branches: the encoder branch increases the number of feature channels layer by layer from 1 to 256, while the spatial resolution halves layer by layer from 256×256 to 8×8. Each encoding layer of the encoder uses two 3×3 convolutions as its core, with a stride of 1, and undergoes batch normalization (BN) and ReLU activation to accelerate training convergence and enhance nonlinear representation. Subsequently, the input x is added to the convolutionally transformed features through a residual connection to form... ,in: ; and These are the first and second 3×3 convolutions, respectively. For batch normalization operations; The parallel terrain gradient sensing module calculates and normalizes the horizontal and vertical gradients of the input phase map, generating a gradient distribution map and dynamic weights. The terrain gradient perception module is deployed in parallel with 6 coding layers. It synchronously receives the preprocessed wrapped phase map with dimensions [B,1,256,256], where B is the batch size. It receives the original input phase map and generates gradient guidance information including gradient distribution map and dynamic dilation rate weight vector, which provides prior information on the steepness of the terrain for subsequent modules. The terrain gradient sensing module's processing method includes four steps: 1) Gradient calculation: Using the Sobel operator adapted to the periodic characteristics of the phase map, the first-order partial derivatives of the input phase map are calculated along the x-axis and y-axis respectively to obtain the horizontal and vertical gradients. The calculation formula is as follows: ; in For horizontal gradient, For vertical gradient, and It uses a 3×3 adapted Sobel convolution kernel; The formula for calculating the gradient magnitude of each pixel, which represents the steepness of the terrain, and the gradient direction, which represents the terrain orientation, is as follows: ; in For horizontal gradient, For vertical gradients, the output is a gradient magnitude map with dimensions [B, 1, 256, 256]. and gradient direction pattern ; 2) Gradient normalization: Calculate the maximum and minimum values of the gradient magnitude map in batches, and perform normalization for each pixel. The calculation formula is as follows: ; in, This represents the maximum value of the gradient magnitude plot. To find the minimum value of the gradient magnitude map, output a normalized gradient magnitude map. and gradient direction plot that retains the original values ; 3) To provide a basis for gradient guidance in the decoding layer, a 1×1 lightweight convolution operation is first performed on the normalized gradient magnitude map to maintain the number of channels at 1, generating a pixel-wise attention weight map. Adaptive scaling is then applied to the feature responses at different spatial locations using this attention weight map, giving higher weights to high-gradient and structurally abrupt regions to enhance edge and detail representation, and lower weights to low-gradient and flat regions to suppress noise perturbations. Finally, a gradient distribution map with dimensions [B, 1, 256, 256] is generated, calculated using the following formula: , in Represents a 1×1 lightweight convolution operation. This represents a gradient-guided spatial weighting operation on the normalized gradient magnitude map. This is the gradient distribution map; the gradient distribution map is directly passed to each upsampling branch of the decoding layer, providing a basis for the gradient guidance enhancement stage in the decoding phase; 4) To generate dynamic inflation rate weights, first normalize the gradient magnitude map. A global average pooling operation is performed to obtain the batch-level average gradient value. Then, based on this average gradient value, adaptive weights are assigned to the four different inflation rate branches of the holed space pyramid pooling ASPP module. These weights are encapsulated as learnable parameters, ultimately outputting a dynamic inflation rate weight vector with dimension 4. The calculation formula is as follows: ; in This represents the global average pooling operation. This represents a branch weight assignment operation based on the average gradient value. This is a dynamic expansion rate weight vector, which is directly passed to the hollow space pyramid pooling ASPP module to weight and fuse features from different expansion rate branches, thereby achieving dynamic adaptation of multi-scale features to terrain gradients.
6. The InSAR phase unwrapping method based on terrain gradient sensing network as described in claim 1, characterized in that: The processing method of the gradient-guided dynamic weight multi-scale hollow pooling module includes: Each parallel basic branch includes a 3×3 convolution corresponding to the basic dilation rate, a BN batch normalization layer, and a ReLU activation function to extract features at a specific scale. The output features of each basic branch are differentiated and weighted based on dynamic inflation rate weights. In high-gradient regions, the response of branches with high inflation rates is enhanced, while in flat regions, the response of branches with low inflation rates is enhanced. The calculation formula is as follows: ; in, For the first The original output of each basic branch, The feature is dynamically enhanced, with 0.5 as the enhancement factor to balance the contributions of each branch; Secondly, the global average pooling branch includes adaptive average pooling output size, 1×1 convolution, BN layer and ReLU activation, which is responsible for extracting global context information; The output features of the global average pooling branch are upsampled to an 8×8 size to maintain the same dimension as the features of the parallel base branch. Finally, the four dynamic enhancement branches are concatenated with the global features along the channel dimension and fused into a 256-channel feature through a 1×1 convolution. The calculation formula is as follows: ; in, It is a fused feature map output by the void space pyramid pooling module, which is a feature map output by the coded residual block. After being processed by dilated convolutions with four different dilation rates, the dynamic weights generated by the terrain gradient sensing unit are... We obtain the weighted sum.
7. The InSAR phase unwrapping method based on terrain gradient sensing network as described in claim 1, characterized in that: The gradient adaptive decoding module includes a first gradient adaptive upsampling unit, a second gradient adaptive upsampling unit, a third gradient adaptive upsampling unit, a fourth gradient adaptive upsampling unit, and a fifth gradient adaptive upsampling unit connected in series, and respectively connected to the first to fifth connection layers. Each gradient adaptive upsampling unit includes a transposed convolutional upsampling layer with a kernel size of 2×2 and a stride of 2. A gradient weight generation layer receives a gradient distribution map matching the current scale and generates spatial attention weights through a 1×1 convolution and a sigmoid function. An enhancement convolution includes two 3×3 convolutions and batch normalization. An activation function and a residual fusion layer add the enhancement branch to the main upsampling branch after weighting it by attention. The upsampling output of each level is concatenated with the corresponding connection layer output in the channel dimension, refined by the convolutional residual fusion layer, and then sent to the next upsampling unit, so that the resolution is restored from 8×8 to 256×256 sequentially. In the entire process, key detail restoration is performed on the high gradient region.
8. The InSAR phase unwrapping method based on terrain gradient sensing network as described in claim 7, characterized in that, In the gradient adaptive decoding module: The gradient adaptive decoding module first performs standard upsampling on the input features, then doubles the spatial resolution through a 2×2 transposed convolution with a stride of 2. The calculation formula is as follows: ; in, Input features for the module, Features after standard upsampling The kernel size is [size]. The step size is set; then, the gradient distribution map output by the terrain gradient sensing unit is... Upsampling to Generate spatial attention weights for the same size For targeted enhancement of high gradient edge regions, the calculation formula is as follows: ; in, This is the gradient distribution map after upsampling. For channel splicing, A 3×3 convolutional layer is used to extract fused features; the enhanced features are then fed into an upsampled residual block for further refinement. In the upsampling residual block and skip connection fusion stage, the output features are compared with the skip connection features of the encoder at the corresponding scale. The data is concatenated along the channel dimension and then fused using two 3×3 convolutional layers to extract features. Gradient propagation and feature reuse are stabilized using a residual approach. Simultaneously, the spatial attention weights obtained in the previous step are combined with these features. Before fusion, the splicing features are slightly adjusted to enhance the information transmission of high-gradient edges. The calculation formula is as follows: ; in, For the features of the corresponding layer of the encoder, This is for channel splicing operations. This is the final output of the decoding layer and serves as the input for the next decoding level.