A wide-area deformation adaptive InSAR network adjustment method

Abstract

The application relates to the technical field of InSAR and discloses a wide-area deformation adaptive InSAR network adjustment method, which comprises the following steps: extracting an effective pixel coordinate set; based on the effective pixel coordinate set, a multi-target optimization function is constructed, the optimal model structure and the regularization parameter of each overlapping area are determined; a mapping model is constructed, the observation weight of each pixel is dynamically adjusted to suppress the influence of abnormal values, and the mapping model parameters of each image are obtained; the parameter feasible region of the curved surface fitting model is adaptively determined according to the geometric shape of the overlapping area, and the spatial curved surface parameters are obtained; the optimal absolute correction parameters of the global image are solved, and the consistent splicing result with the minimum global error and taking the central image as the reference is obtained. Through the introduction of the adaptive window mechanism driven by the time series coherence stability index and the terrain roughness, the adaptability of the robust identification of seasonal incoherent pixels and the extraction of effective observation points under complex terrain is enhanced.

Description

Technical Field

[0001] This invention relates to the field of InSAR technology, and more specifically, to a wide-area deformation adaptive InSAR network adjustment method. Background Technology

[0002] Temporal synthetic aperture radar interferometry (TS-InSAR) technology can acquire high-precision and time-stable surface deformation information within a certain range, and has been widely used in deformation-prone areas such as tectonic movements, geological hazards, and engineering construction. With the in-depth development of Earth system science, the need for monitoring surface deformation has expanded from point-like or local deformation objects (such as landslides and mining subsidence) to regional and even intercontinental scale geodynamic processes, directly driving the transformation of InSAR technology towards wide-area deformation inversion research.

[0003] However, due to limitations in satellite orbit design, imaging swath width, and observation geometry, a single satellite platform can hardly cover a large, continuous area in a single transit. In practical applications, it is usually necessary to stitch and fuse InSAR deformation fields from multiple orbits and frames. This stitching process faces the following core challenges:

[0004] (1) Track and line-of-sight (LOS) geometric deviation: Different tracks have independent imaging geometry and processing references. The deviation between track and line-of-sight (LOS) geometry will induce cross-track long-wavelength errors, which will weaken the spatial comparability of multi-track data. The spatial correlation characteristics of propagation medium delay will further superimpose regional system deviations, resulting in a gradual drift in the large-scale deformation field. If the relative references between tracks and frames are not finely corrected, it will cause a jump from the overlapping area to the non-overlapping area, making the cross-regional deformation comparison lose a unified scale.

[0005] (2) Insufficient pixel quality control: Most existing methods do not explicitly introduce coherence constraints in the pixel selection of overlapping areas, resulting in low-quality observations being directly included in the adjustment system. Even with a coherence mask using a fixed threshold, the rigid selection strategy is difficult to adapt to the significant spatial heterogeneity within the overlapping area. In narrow-band overlapping areas, coherence often exhibits a directional structure that is lower along the track direction and higher perpendicular to the strip direction, making the observation points in some areas too sparse and weakening the spatial constraint capability of the adjustment model.

[0006] (3) Insufficient adaptability of model structure: The difference fitting in the overlapping area generally relies on empirical models with fixed structures, which are difficult to cope with the complex geometric constraints in wide-area stitching. In the same track narrow band region, local noise can easily lead to overfitting of high-order models and produce cross-frame pseudo tilt; while in the cross-track wide band region, spatial heterogeneity can easily lead to underfitting of low-order models and cannot effectively absorb trend errors.

[0007] (4) Insufficient robustness of adjustment: Due to the influence of low coherence noise and local phase unwrapping jumps, outliers are often present in the same observations. If the adjustment system lacks an effective robust estimation mechanism, these local gross errors will spread and accumulate along the adjustment network to the whole domain, eventually destroying the consistency of the deformation field.

[0008] (5) Error propagation: Relying solely on pairwise registration between adjacent images, errors accumulate as the stitching path extends, and spatial closure loop contradictions are easily triggered in areas with multiple overlaps.

[0009] Existing methods addressing the aforementioned issues have developed into three main technical approaches: benchmark constraints, geometric consistency adjustment of overlapping regions, and error field modeling. While methods based on Global Navigation Satellite System (GNSS) benchmark constraints provide a high-precision absolute deformation benchmark, the sparse and spatially uneven distribution of GNSS stations still results in significant inconsistencies between orbital boundaries and low-coherence regions. Uncontrolled adjustment methods based on geometric consistency of overlapping regions can alleviate dependence on GNSS stations to some extent. Their core lies in utilizing the relative constraints provided by the rate differences in overlapping regions to unify the offsets between orbits. Building upon this, researchers have successively proposed adaptive grid adjustment models, K-fold cross-validation models, and methods incorporating physical parameter corrections to accurately model the error field. Summary of the Invention

[0010] The technical problem to be solved by this invention is to solve the problem of deformation rate field fragmentation in wide-area InSAR multi-track and multi-frame stitching in the prior art, and to achieve high-precision wide-area deformation field consistency stitching without external control points.

[0011] To address the aforementioned technical problems, this invention proposes a wide-area deformation adaptive InSAR network adjustment method, comprising:

[0012] S1. Obtain the coherence values ​​of each pixel in the overlapping area of ​​adjacent image pairs, use a multi-level adaptive coherence mask to perform layered screening of each pixel in the overlapping area, and determine the effective pixels based on the consistency condition between the coherence values ​​and the statistical features of the local neighborhood, and extract the set of coordinates of the effective pixels.

[0013] S2. Based on the effective pixel coordinate set, construct a multi-objective optimization function that includes weighted difference residuals, in-orbit fitting residuals and regularization terms, and jointly optimize the preset candidate model structure and regularization intensity to determine the optimal model structure and regularization parameters for each overlapping area.

[0014] S3. Under the constraints of the optimal model structure and regularization parameters, a mapping model is constructed based on the set of effective pixel coordinates. The mapping model parameters are iteratively updated using an alternating least squares framework, and the observation weights of each pixel are dynamically adjusted to suppress the influence of outliers, thus obtaining the mapping model parameters for each image.

[0015] S4. Calculate the residual value of the overlapping area using the mapping model parameters, perform spatial surface fitting on the residual value of the overlapping area using random consistency sampling, and adaptively determine the feasible region of the parameters of the surface fitting model according to the geometric shape of the overlapping area to obtain the spatial surface parameters.

[0016] S5. Integrate the mapping model parameters and spatial surface parameters into a relative observation vector, use the observation quality evaluation index of the overlapping area as the observation edge weight, and incorporate all observation edges into the sparse matrix adjustment framework to solve for the optimal absolute correction parameters of the global image, thereby obtaining a consistent stitching result with the central image as the reference and the global error minimized.

[0017] The beneficial effects of this invention are as follows:

[0018] 1. This invention, through MLACM, introduces local statistical consistency and neighborhood coverage constraints, combined with a hierarchical screening strategy, effectively avoiding the problems of sparsity of effective points or noise amplification caused by fixed threshold screening. MLACM can suppress accidental high-coherence points while retaining true stable pixels in weakly coherent regions, forming an effective coherent dataset with more robust statistical characteristics and more balanced spatial distribution. Furthermore, by introducing a temporal coherence stability index and an adaptive window mechanism driven by terrain roughness, it enhances the robustness of identifying seasonally decoherent pixels and the adaptability of effective observation point extraction under complex terrain.

[0019] 2. This invention constructs a multi-objective optimization function, which can adaptively determine the optimal model structure and regularization parameters, balancing model complexity and inter-track consistency, and effectively preventing overfitting in cross-frame regions within the same track and underfitting in cross-track regions. Furthermore, through DBSCAN density clustering and consensus distance optimization strategy, the robustness of selecting the optimal parameter combination in multi-objective optimization is improved, avoiding the problem that a single Euclidean distance criterion is prone to getting trapped in local optima when parameter evaluation indicators are close.

[0020] 3. This invention uses DRWF-based ALS to dynamically adjust observation weights according to the residual distribution, automatically dividing observations into retention zones, deweighting zones, and elimination zones, effectively suppressing the interference of abnormal observation points on parameter estimation and ensuring the robustness of iterative convergence.

[0021] 4. This invention uses geometrically constrained adaptive RANSAC spatial residual surface correction to automatically determine the feasible region of the fitting model parameters based on the aspect ratio of the overlapping area. It avoids false slopes in the weak observation direction and fully characterizes the two-dimensional spatial trend in the square region, effectively eliminating the directional spatial bias caused by noisy structures and geometric shapes.

[0022] 5. This invention uses global consistency adjustment to unify all pairwise difference constraints into the sparse matrix adjustment framework, eliminating the contradictions of error propagation and spatial closure loops in long strip or large-scale cross-track splicing, and obtaining velocity field splicing results with good closure and consistent trends. Attached Figure Description

[0023] The invention will now be further described with reference to the accompanying drawings.

[0024] Figure 1 This is a flowchart illustrating a wide-area deformation adaptive InSAR network adjustment method as an example of the present invention;

[0025] Figure 2 This is a comparison chart of the multi-frame deformation rate results before and after correction in three time periods (March 2017 - December 2018, January 2019 - December 2020, January 2021 - December 2023) of the test area of ​​this invention example;

[0026] Figure 3 This is a comparison chart of the accuracy statistics (RMSE, STD) of the original results and the adaptive InSAR network adjustment results in the example of this invention.

[0027] Figure 4 for Figure 2 One of the detailed consistency analysis diagrams of the mid-section line;

[0028] Figure 5 for Figure 2 Part Two: Detail Consistency Analysis of the Mid-Section Line;

[0029] Figure 6 for Figure 2 The third in a series of detailed consistency analysis diagrams of the mid-section lines. Detailed Implementation

[0030] The present invention will now be described in detail with reference to the accompanying drawings, which will make the technical approach and operation steps of the present invention clearer.

[0031] According to an embodiment of the present invention, a wide-area deformation adaptive InSAR network adjustment method is provided.

[0032] The present invention will now be further described in conjunction with the accompanying drawings and specific embodiments, such as... Figure 1 As shown, according to an embodiment of the present invention, a wide-area deformation adaptive InSAR network adjustment method is provided, comprising:

[0033] S1. Obtain the coherence values ​​of each pixel in the overlapping area of ​​adjacent image pairs, use a multi-level adaptive coherence mask to perform layered screening of each pixel in the overlapping area, and determine the effective pixels based on the consistency condition between the coherence values ​​and the statistical features of the local neighborhood, and extract the set of coordinates of the effective pixels.

[0034] It should be explained that the effective homonym extraction based on the Multi-Level Adaptive Coherence Mask (MLACM) involves introducing local statistical consistency, neighborhood coverage constraints, and a hierarchical screening strategy. The coverage and local median of high-coherence pixels within the sliding convolution window are calculated, and three screening modes—strict, moderate, and lenient—are set. When the coherence value of a single pixel and the coverage and median of its neighborhood simultaneously meet the consistency condition, it is considered a valid observation. The effective pixels in the overlapping area are extracted as the basis for the adjustment network.

[0035] It should be explained that the consistency condition refers to the comprehensive constraint between the temporal coherence value of a single pixel and the statistical characteristics of its local neighborhood. Specifically, it requires the simultaneous satisfaction of the following three sub-conditions: 1) The temporal coherence value of the pixel itself. ≥ , , That is, not lower than the coherence threshold set for the current level; 2) in a local sliding window centered on the pixel. Within, the proportion of pixels that meet the high coherence threshold ≥ 3) Median of coherence within the window ≥ This is to rule out the possibility that the pixel is an isolated noise peak. Only when all three conditions mentioned above are met is the pixel considered a valid observation.

[0036] Specifically, based on the initial threshold, the calculation is performed for each pixel. Coverage of high-coherence pixels within the center-based sliding convolution window With local median Set three-level screening thresholds. Initially, a strict mode is used. When the strict mode fails to meet the minimum effective point threshold in narrow band regions along the same track, it switches to a medium mode. If the effective observation set still does not meet the size requirement, it is further relaxed to a relaxed mode to ensure that basic spatial constraints are maintained in low-coherence regions. Based on this, the effective pixel coordinate set for subsequent network adjustment is extracted. The specific formula is as follows:

[0037] ;

[0038] ;

[0039] ;

[0040] in, The optimal level to satisfy the point requirements (0 is the most stringent, 2 is the most lenient); These are the cell coordinates; To effectively estimate the minimum number of points required; For the first Binary masks generated by the hierarchical strategy; For the first Image in the center of the window Coherence value at the location; This represents the median of coherence within the window. This represents the minimum proportion of highly coherent points within the window. , , These are the strict mode threshold, the medium mode threshold, and the lenient mode threshold, respectively. The percentage of pixels within the window that meet the high threshold requirement; This represents the total number of pixels within the window. This is an indicator function (1 if the condition is met, 0 otherwise); Therefore A local sliding window centered on the target area; For the first Image in neighboring pixels The coherence value.

[0041] In the initial stage of applying MLACM, high coherence constraints are preferred, i.e. , , The size is 3×3; in the narrow band region of the same track, strict mode may not meet the minimum number of points threshold requirement, so it switches to relaxed mode, i.e. If the effective observation set still does not meet the size requirement, the threshold will be further relaxed, i.e. This is to ensure that basic spatial constraints can still be maintained in regions of low coherence.

[0042] 1) To address the issue that some pixels may experience sudden decoherence only during certain time periods due to seasonal vegetation changes, abrupt changes in land cover, or short-term strong disturbances, this invention further introduces a temporal coherence stability mechanism. Specifically, for each pixel, its coherence time series across all interferometric pairs is extracted, and a temporal stability index (TSI) based on MAD is calculated. TSI measures the robust dispersion of the coherence time series curve from the median level at each time period, effectively identifying pseudo-high coherence points that only show sudden abnormal drops within individual interferometric windows, while remaining unaffected by excessive interference from a few extreme jump values. It is important to emphasize that for TSI exceeding a preset threshold... For pixels, this invention does not employ a direct elimination strategy, but instead uses a downgrade retention process: the pixel is downgraded one level within the currently determined three-level screening hierarchy (strict / medium / loose). Pixels that originally met the strict mode are downgraded to the medium mode, and pixels in the medium mode are downgraded to the loose mode. The downgraded pixels still participate in the subsequent adjustment network as valid observations, thus maintaining the continuity of spatial constraints in the overlapping area; however, their weighting coefficients are correspondingly reduced in each step S2-S5, so that the final adjustment result is dominated by highly coherent and stable pixels, while the downgraded pixels play an auxiliary role in spatial constraints.

[0043] 2) Considering the significantly uneven spatial distribution of InSAR coherence under complex terrain conditions such as mountainous and urban areas, a fixed-size sliding convolution window cannot simultaneously ensure both statistical robustness and spatial resolution. This invention utilizes a DEM to extract the local terrain roughness at each pixel within the overlapping area. (Using the standard deviation of neighborhood elevation as a metric), the size of the sliding window is adaptively adjusted accordingly: in In flat areas, a window size of 9×9 is used to accommodate more statistical samples, enhancing the robustness of coverage and median estimation; For steep slopes, valleys, or densely built-up areas, the window size should be 3×3 to reduce the aliasing effect of heterogeneous object pixels and maintain the accuracy of local constraints. When the value falls between these two extremes, the window side length is determined by a linear mapping. This adaptive mechanism makes the effective pixel selection no longer dependent on a single fixed spatial constraint rule, and is more in line with the monitoring conditions of uneven coherence distribution under actual surface conditions.

[0044] Finally, a binary mask of effective pixels in the overlapping area is generated through a multi-level strategy, and the set of effective pixel coordinates for subsequent network adjustment is extracted from this mask. .

[0045] Specifically, S1 includes:

[0046] The local terrain roughness of each pixel coordinate is extracted using a digital elevation model. The sliding convolution window size of each pixel is adaptively adjusted according to the local terrain roughness. The largest window size is used in flat areas where the terrain roughness is less than or equal to a first preset threshold, and the smallest window size is used in steep slopes or densely built-up areas where the terrain roughness is greater than or equal to a second preset threshold. When the terrain roughness is between the first and second preset thresholds, the window size is determined by linear mapping.

[0047] Based on a preset initial high coherence threshold, the coverage and local median of high coherence pixels within the sliding convolution window centered on each pixel are calculated.

[0048] Three levels of screening are set: strict screening mode (i.e., strict mode above), medium screening mode (i.e., medium mode above), and lenient screening mode (i.e., lenient mode above). Each screening mode corresponds to a different coherence threshold, coverage threshold, and local median threshold. In the initial stage, the threshold combination corresponding to the strict screening mode is used to determine each pixel in the overlapping area pixel by pixel.

[0049] A cell is marked as a valid cell when its coherence value is greater than or equal to the coherence threshold of the current strict filtering mode, the coverage within the window containing the cell is greater than or equal to the coverage threshold of the current strict filtering mode, and the local median within the window containing the cell is not lower than the local median threshold of the current strict filtering mode.

[0050] The system counts the number of valid pixels marked in the current strict filtering mode. If the number of valid pixels reaches the preset threshold, the filtering is terminated; otherwise, it switches to the lenient filtering mode.

[0051] The overlapping area is re-evaluated pixel by pixel using the threshold combination corresponding to the lenient filtering mode until the number of valid pixels reaches the preset valid point threshold or all filtering modes have been traversed.

[0052] Extract the coherence time series of each pixel marked as valid by each screening mode on all interference pairs, and calculate the time series stability index measured by median absolute deviation.

[0053] All pixels that have been selected through a three-level stratified screening process and simultaneously meet the coherence threshold, coverage threshold, and local median threshold conditions, and have been adjusted for temporal stability constraints, will be considered as valid pixels for subsequent network adjustment, and the set of valid pixel coordinates will be output.

[0054] S2. Based on the effective pixel coordinate set, construct a multi-objective optimization function that includes weighted difference residuals, in-orbit fitting residuals and regularization terms, and jointly optimize the preset candidate model structure and regularization intensity to determine the optimal model structure and regularization parameters for each overlapping area.

[0055] It needs to be explained that the multi-objective joint optimization of the adjustment model and regularization parameters involves: for each pair of adjacent images, a candidate model set and multiple preset regularization intensities are defined. A multi-objective optimization function vector containing weighted difference residuals, in-orbit fitting residuals, and regularization terms is constructed. This vector is input into alternating least squares iteratively, in the order of constant model and linear model. Pareto optimization is used to eliminate inferior solutions across all evaluation dimensions to construct a set of non-dominated solutions. Preliminary selection is performed by calculating the normalized Euclidean distance of each solution to the ideal origin. Furthermore, the globally optimal model structure and regularization parameters are determined using DBSCAN density clustering (i.e., spatial clustering) and consensus distance. This invention defines a candidate model set. and in multiple preset regularization strengths Next, construct the multi-objective optimization function vector:

[0056] ;

[0057] In the formula, Represents the set of candidate models; Indicates transpose; This represents a regularization term, which serves as a penalty for the complexity of the linear model to prevent overfitting. This represents the weighted difference residual, used to measure the difference between two images. Consistency of alignment within overlapping areas; high-coherence pixels are given higher weights during calculation to ensure they dominate differential consistency evaluation. Indicates the first In-orbit fitting residuals of the image; Indicates the first The in-orbit fitting residuals of the image are used to measure the performance of the current model structure and (The search space is set based on orbital differences) The combination of latent variables and parameters... , The interpretability of the original velocity fields of two adjacent images; This represents a vector of multi-objective optimization functions.

[0058] Specifically, the formulas for calculating the weighted difference residuals and in-orbit fitting residuals of the multi-objective optimization function vector are as follows:

[0059] ;

[0060] ;

[0061] ;

[0062] in, and Images , In pixels The coordinates determine the fitting rate of the model (constant, linear). The general structure is defined as follows: , . The total weight of the two images in that pixel, i.e. ; , Images , In pixels The weight of the coordinates is the square of the coherence value; , For images , In pixels The deformation rate value of the coordinate. This is the regularization term. Model selection refers to: in the multi-objective optimization process, both constant and linear models are evaluated sequentially. By comparing the normalized Euclidean distances of the optimal solutions in the Pareto front, the final model to be used for the image pair is determined, along with the corresponding optimal regularization strength. .

[0063] Specifically, the range normalization of the three components of the multi-objective optimization function vector is used to eliminate dimensional differences, and the non-dominated solution set is calculated. Normalized Euclidean distance from each solution to the ideal origin (i.e., the point where all components are local minima). After constructing the Pareto non-dominated solution set and calculating the normalized Euclidean distance, this invention further performs density-based spatial clustering (DBSCAN) on the non-dominated solutions, grouping candidate solutions that are close in distance in the target space. Grouping solutions into the same cluster. After clustering, each cluster has two inherent properties: the cluster representative solution (the solution with the smallest average normalized distance to all members in the cluster) and the cluster size (the number of members in the cluster). Based on this, this invention uses consensus distance instead of the single normalized Euclidean distance as the selection criterion:

[0064] ;

[0065] In the formula, For the first The normalized Euclidean distance from the representative solution of each cluster to the ideal origin. For the first The size of a cluster, This is the consensus strength adjustment factor (default value is 0.5). The physical meaning of this criterion is: the larger the cluster size, the more non-dominated solutions are present in that cluster. A consensus reached near a given set indicates that the set is less sensitive to data disturbances and parameter fine-tuning, thus exhibiting stronger robustness; consensus distance is achieved through... For the original distance Regularization is applied to ensure that solutions that are slightly farther apart in large clusters are selected before solutions that are slightly closer apart in small clusters (especially isolated point clusters). The final consensus distance is then chosen. The smallest cluster represents the solution as the global optimum. Specifically, when all non-dominated solutions are grouped into a single cluster, consensus optimization automatically degenerates into classical Euclidean distance minimization, maintaining backward compatibility with the original method. Furthermore, an absolute accuracy tolerance is set: when the weighted difference residual of the constant term model under optimal regularization strength is better than this threshold, no linear term is introduced.

[0066] Considering the ideal situation and It should approach zero, therefore, through calculation Normalized Euclidean distance from each solution to the ideal origin (minimum point) Find the globally optimal combination with the shortest distance, which is the optimal solution for each overlapping region:

[0067] ;

[0068] This invention sets the absolute accuracy tolerance at 0.008m / a, that is, when the constant model is in When the data quality in the overlapping region is already better than the threshold, it indicates that the data quality is good and the error is gradual, and no further linear terms are introduced. To effectively suppress overfitting, The search space is set according to the differences in the orbits: in the same orbit mode, the system error is simple, so it is set to 0.2-1.2 (step size 0.1) to emphasize stability; in the cross-orbit mode, it is set to 0.05-1.0 (step size 0.05) to more finely adapt to complex noise. This represents the weighted difference residual, which measures the difference in deformation rate between two adjacent images under the selected model and regularization intensity in the overlapping area. , Representing the impact and The in-orbit fitting residuals; Indicates the first The in-orbit fitting residuals of the image after range normalization; Indicates the first The in-orbit fitting residuals of the image after range normalization; This represents the weighted difference residuals after range normalization; This represents the regularization term after range normalization.

[0069] Specifically, S2 includes:

[0070] For each pair of adjacent images, a candidate model set and multiple preset regularization intensities are defined, and a multi-objective optimization function vector containing weighted difference residuals, in-orbit fitting residuals, and regularization terms is constructed; the candidate model set includes constant models and linear models;

[0071] The constant model and the linear model are combined with each preset regularization intensity to form candidate parameter combinations. Each candidate parameter combination is then input into the alternating least squares framework for iterative solution to obtain the weighted difference residuals and in-track fitting residuals corresponding to each candidate parameter combination.

[0072] The weighted difference residuals and in-orbit fitting residuals corresponding to each candidate parameter combination, together with the corresponding regularization term values, are used as the values ​​of each component of the multi-objective optimization function vector. The values ​​of each component of the multi-objective optimization function vector are normalized by range to eliminate dimensional differences. Pareto optimization is used to eliminate inferior solutions in all evaluation dimensions to construct a set of non-dominated solutions.

[0073] Calculate the normalized Euclidean distance from each solution in the non-dominated solution set to the ideal origin, and perform spatial clustering on the non-dominated solutions. Combine candidate models with regularization strengths that are close in distance in the target space into the same scheme cluster. Use consensus distance instead of the single normalized Euclidean distance as the selection criterion, and select the cluster with the smallest consensus distance as the global optimal model structure and regularization parameter.

[0074] An absolute precision tolerance is set. When the weighted difference residual of the constant term model under the optimal regularization strength is better than the absolute precision tolerance, it is determined that the data quality of the current overlapping area is good and the error is gradual. No more linear terms are introduced, and the constant model is used as the final selected optimal model structure.

[0075] S3. Under the constraints of the optimal model structure and regularization parameters, a mapping model is constructed based on the set of effective pixel coordinates. The mapping model parameters are iteratively updated using an alternating least squares framework, and the observation weights of each pixel are dynamically adjusted to suppress the influence of outliers, thus obtaining the mapping model parameters for each image.

[0076] It should be explained that the systematic differences between different image pairs in the overlapping area caused by orbital geometry, line-of-sight deviation, and time baseline are explicitly modeled as shared latent variables to reflect the potential true rate of their combined influence. A mapping model is constructed, including parameters such as scale scaling factor, constant offset term, and spatial linear term. An alternating least squares (ALS) framework is used to alternately update the latent variables and mapping model parameters. Simultaneously, a dynamic robust weight function (DRWF) is constructed. After each iteration, the observation residuals are calculated to obtain standardized residual statistics, which automatically divide the observation weights into retention, deweighting, and elimination zones, and robustly solve for the scale and offset terms.

[0077] Based on the obtained This invention explicitly models the systematic differences between different image pairs in the overlapping region caused by orbital geometry, line-of-sight deviation, and time baseline as shared latent variables. This is used to characterize the potential true surface deformation rate of the region. In the adjustment model (i.e., the mapping model), the first... Image in pixels The original observation rate at the location and The mapping relationship between them is determined by the scale scaling factor, constant offset, and spatially dependent polynomial coefficients. Therefore, for any... Construct the following weighted least squares objective function:

[0078] ;

[0079] Among them, the mapping model The general structure is defined as follows:

[0080] ;

[0081] ;

[0082] In the formula, This represents the weighted objective function; Represents polynomial parameters, ; Represents a pixel Observation weights of coordinates; Indicates the first Image in pixels The original observed deformation rate of the coordinates; Indicates the pixel to be estimated Latent variables of coordinates; Indicates the scale factor; This represents the system constant offset term; Indicates the coefficient of the spatial linear term; Represents spatial coordinate components; Represents the set of valid cell coordinates; This represents the cell index.

[0083] and The scale consistency and system offset compensation between images were jointly determined. The strong coupling between these two factors makes it highly susceptible to ill-conditioned normal equations when using ordinary least squares (ALS) for solution. Therefore, ALS was adopted as the core framework for solving the adjustment model. ALS, given an initial... Under the conditions, solve for Minimized and in After the update, the solution is obtained through joint least squares inverse solution. This alternating update ensures that each step reaches conditional optimum under the current weight structure, thereby effectively mitigating [the problem]. Instability arises from the coupling between cross-track and cross-frame parameters. In ALS iterative solutions, due to factors such as unwrapping errors, complex cross-track geometric differences, and incompletely modeled atmospheric phases, the actual observation residuals often deviate from a Gaussian distribution and exhibit heavy-tailed characteristics. To address this, this invention constructs a DRWF based on the residual distribution, dynamically adjusting the observation weights in each iteration through the statistical distribution of the residuals.

[0084] Specifically, step a: fix By using weighted least squares inverse solution, Minimization of the estimated value Step b: Fix Solve Minimized Update values; Step c: Calculate the observation residuals under the latest parameters, and update the observation weights using DRWF; Repeat steps a-c until the parameter change rate is lower than the preset convergence threshold or the maximum number of iterations is reached.

[0085] Specifically, in the After one round of iteration completes a parameter update, the calculation is performed. Observation residuals of each pixel ,Right now To avoid a few extreme gross errors contaminating the overall error scaling estimate, MAD is introduced to calculate robust standard deviation, and a standardized residual statistic is constructed:

[0086] ;

[0087] In the formula, Represents standardized residual statistics; Represents the robust standard deviation estimate of the residual sequence in the overlapping region; Represents the median operator; represents the normality consistency factor, which makes MAD an unbiased estimator of the standard deviation under the Gaussian distribution assumption; here it is set to 0.6745. This represents the observation residual.

[0088] Based on the magnitude of the standardized residuals, the observed data are automatically divided into a retention zone, a deweighted zone, and a rejection zone:

[0089] ;

[0090] In the formula, As a weighting factor; and These are the preservation threshold and the exclusion threshold, respectively. When the residual is less than When the residual is between two values, it is considered normal observation noise, and the weight remains unchanged; when the residual is between two values, it is considered a suspicious value, and the weight is reduced using a quadratic function; when the residual is greater than two values, it is considered a suspicious value. When this occurs, it is considered a gross error, and the weight is reset to zero. This invention sets the threshold to... and .based on The observation weights of the current pixels are dynamically updated to obtain robust equivalent weights. Furthermore, the maximum number of iterations is set to 100, when the parameter change rate ( The convergence threshold is less than the preset convergence threshold (set to). If the number of iterations remains stable after more than 10 iterations, it is considered to have converged.

[0091] Specifically, S3 includes:

[0092] Based on the effective set of pixel coordinates, a mapping model containing scale factor, offset term and spatial first-order term parameters is constructed for each pixel in the overlapping area of ​​each image. The latent variable is used to characterize the potential real surface deformation rate of the overlapping area, and the observation weight is used to measure the reliability of the observation value of each pixel. A weighted least squares objective function is constructed.

[0093] An iterative solution is performed using an alternating least squares framework. Latent variables are fixed, and the estimated values ​​of the mapping model parameters are obtained by weighted least squares. The updated values ​​of the latent variables are obtained by fixing the mapping model parameters. The observation residuals under the current mapping model parameters are calculated, and the observation weights of each pixel are updated according to the dynamic robust weight function.

[0094] Specifically, updating the observation weights of each pixel according to the dynamic robust weight function includes:

[0095] After the mapping model parameters are updated in the current iteration, the difference between the original observed deformation rate of each pixel in the overlapping area and the predicted value of the current mapping model is calculated to obtain the observation residual of each pixel.

[0096] Robust standard deviation estimates of the residual sequences in the overlapping region are calculated by combining the median absolute deviation, and standardized residual statistics are constructed.

[0097] Based on the magnitude of the standardized residuals, the observation data is automatically divided into a retention zone, a deweighting zone, and a removal zone. In the retention zone, the observation weights remain unchanged. In the deweighting zone, the observation weights are reduced. In the removal zone, the observation weights are reset to zero, resulting in the updated observation weights for each pixel.

[0098] Repeat the process of solving for the estimated values ​​of the mapping model parameters, solving for the updated values ​​of the latent variables, and calculating the observation residuals until the rate of change of the mapping model parameters is lower than the preset convergence threshold or the maximum number of iterations is reached, and output the mapping model parameters of each image.

[0099] S4. Calculate the residual value of the overlapping area using the mapping model parameters, perform spatial surface fitting on the residual value of the overlapping area using random consistency sampling, and adaptively determine the feasible region of the parameters of the surface fitting model according to the geometric shape of the overlapping area to obtain the spatial surface parameters.

[0100] It should be explained that after ALS convergence, Random Sample Consensus (RANSAC) is used to fit a low-frequency surface to the spatial residuals. The aspect ratio of the overlapping region is defined, and an aspect ratio threshold is set. The magnitude of the aspect ratio of the overlapping region and the threshold are compared to determine which directions of fitting parameters to retain. A lower limit proportion of inlier points participating in the fitting and a residual tolerance are set to ensure robustness of the fit. The low-frequency spatial surface parameters are obtained by performing conventional least-squares refit using the selected set of high-confidence inlier points.

[0101] RANSAC is used to find a spatially coherent residual structure from the residual point set, that is, to find a set of parameters to be estimated. To make it within the error threshold By maximizing the number of data points (inliers), robust trend parameters can be estimated.

[0102] ;

[0103] ;

[0104] in, For the optimal set of surface coefficients; The set of parameters to be estimated, i.e. ; For pixels The residuals of the coordinates after ALS convergence; Fit a surface to the spatial trend; These are the coordinates of the pixel in the local planar coordinate system. The residual threshold for RANSAC to determine interior points; 、( , ), ( , , ) are the constant term coefficients of the surface fitting model, respectively. , The coefficient of the first term in the direction, in , The coefficients of the quadratic terms in the direction and its intersection direction; This represents the final, residual minute error; The model order is used. In this invention, a pixel refers to its coordinates. .

[0105] Specifically, this invention designs an adaptive spatial trend modeling strategy driven by geometric constraints. The overlapping region is defined in... , Boundary extrema in direction and Then the geometric aspect ratio Further based on and The feasible region of adaptive model parameters was constructed:

[0106] ;

[0107] in, It is a three-dimensional or six-dimensional vector space. The dimension shown is... ; , All are aspect ratio thresholds, where the following settings are provided: =2.0, =0.5: At that time, the overlapping area is The axial extension is much greater than The algorithm forcibly removes all axes that are related to the axis. The relevant parameters; conversely, At that time, only retain In the case of directional parameters, all parameters are retained; otherwise, all parameters are retained. In The set of parameters to be estimated, i.e. , The geometric aspect ratio, i.e. , The feasible region for adaptive model parameters, that is, through Determine which coefficient to substitute into the parameter to be estimated, and then... Together they form the spatial trend fitting surface as parameters. For residual surfaces Perform fitting.

[0108] Based on the above process, a lower limit for the proportion of interior points is set at 60%, and the residual tolerance is set at 0.005m / a to ensure the robustness of the fit. RANSAC maximizes the cardinality of interior points under the guidance of directional geometric priors, effectively suppressing interference caused by geometric degradation in weakly observed directions. Conventional least-squares refitting is then performed using the selected high-confidence interior point set to obtain... This provides a reliable relative benchmark difference vector for the subsequent construction of a global absolute adjustment network.

[0109] Specifically, S4 includes:

[0110] Based on the difference between the original observed deformation rate and the predicted value of the mapping model in the overlapping area of ​​each image, the residual value of each pixel in the overlapping area is calculated.

[0111] The random consistency sampling method is adopted to randomly select the smallest sample subset from the residual values ​​for spatial surface fitting, and the number of interior points of each candidate surface within the preset residual tolerance is counted. The surface coefficient set with the largest number of interior points is taken as the optimal surface coefficient set.

[0112] Calculate the aspect ratio of the overlapping area, and adaptively determine the feasible region of the parameters of the spatial surface fitting model based on the comparison between the aspect ratio and the preset aspect ratio threshold.

[0113] By using the selected set of interior points and constrained by the parameter feasible region, the parameters of the spatial surface are obtained through least-squares refitting.

[0114] S5. Integrate the mapping model parameters and spatial surface parameters into a relative observation vector, use the observation quality evaluation index of the overlapping area as the observation edge weight, and incorporate all observation edges into the sparse matrix adjustment framework to solve for the optimal absolute correction parameters of the global image, thereby obtaining a consistent stitching result with the central image as the reference and the global error minimized.

[0115] It should be explained that all pairwise difference constraints are incorporated into the sparse matrix adjustment framework, integrating the ALS and RANSAC parameters of each image into a unified relative observation vector. The product of the RANSAC in-point ratio in the overlapping area and the average coherence of MLACM is used as the observation edge weight to construct a global relative observation matrix. The spatial topological center image is selected as the absolute reference benchmark, and constraints are added. The optimal absolute correction parameters for all images in the global domain are obtained by constructing and inverting the augmented method equation, resulting in a consistent stitching result with the center image as the benchmark and minimizing global error. This invention employs a two-step stitching strategy from local to global: first, robust estimation is used to obtain the relative spatial benchmark differences between adjacent image pairs; then, global consistency adjustment is introduced to solve for the optimal absolute correction parameters.

[0116] Before constructing the global network, the obtained local parameters are integrated into unified relative observations. ALS provides the systematic migration, scale variation, and one-dimensional linear terms for each image, while RANSAC obtains the low-frequency spatial surface parameters. The definition of the first... ALS parameter vector of the image For any image pair with effective overlap ( The relative observation vector is composed of the difference between the ALS parameters and the RANSAC surface parameters of the two, i.e. .

[0117] All Image is defined as adjustment network The set of nodes in Overlapping image pairs are defined as the set of valid edges. , ,Right now Corresponding image pairs ( For any edge Its relative observation equation and weights can be established as follows:

[0118] ;

[0119] in, Let be the adjustment residual vector of this local observation edge; For the first The weights of each observation edge. To make the adjustment network more responsive to high-quality observation constraints, It is defined as the product of the reliability of the spatial geometric fitting of the overlapping region and the quality of the observed phase. and These represent the number of interior points in the overlapping region during the RANSAC surface correction process and the total number of pixels involved in the fitting process, respectively. This represents the average coherence of the effective point set extracted by the MLACM algorithm within the overlapping region.

[0120] By combining all valid observation edges of the network, a global relative observation matrix is ​​constructed: .

[0121] In the formula, For the whole region The absolute correction parameter matrix of the image to be estimated; For the reason A vector of constant terms formed by concatenation; To include all The global relative residual vector; This is the network topology correlation matrix. This represents the transpose operation of a matrix or vector.

[0122] because Essentially, it is a free network whose normal equation coefficient matrix is ​​not of full rank, exhibiting rank deficiency in the corresponding parameter dimensions. In the absence of globally accurate external control points, In the (adjustment network), the spatial topological center image is selected as the absolute reference datum (i.e., constraints are added). This is equivalent to setting the correction parameters of the reference image to 0, and using it as the benchmark for the correction parameters of all other images. A weighted least squares objective function with additional constraints is constructed as follows:

[0123] ;

[0124] in, This is the solution that minimizes the objective function; The relative observation weight vector in the weighted least squares function; The function is a diagonal matrix indicator. For the first The weights of each observation edge; For the first The weights of each observation edge, and Let be the total number of observed edges; by constructing the augmented method equation and finding its inverse, the optimal absolute correction parameters for all images in the entire domain can be obtained. Ultimately, a consistent stitching result is obtained with the central image as the reference and the global error minimized.

[0125] Specifically, S5 includes:

[0126] For each pair of adjacent images with effective overlap, calculate the difference between the mapping model parameters of the two images and the difference between the spatial surface parameters of the two images, and use the combination of the difference between the mapping model parameters and the difference between the spatial surface parameters as the relative observation vector of the observation edge.

[0127] Calculate the observation edge weight for each overlapping region. The observation edge weight is the product of the proportion of random consistent sampling points in the overlapping region and the average coherence of effective pixels.

[0128] Define all images as adjustment network nodes, define image pairs with effective overlap as effective edges, and solve the relative observation equations of all effective observation edges to construct the global relative observation matrix.

[0129] In the absence of external control points, the spatial topological center image is selected as the absolute reference benchmark and constraints are added. A weighted least squares objective function with additional constraints is constructed. By constructing the augmented method equation and inverting it, the optimal absolute correction parameters for all images in the entire domain are obtained, resulting in a consistent stitching result with the center image as the benchmark and the global error minimized.

[0130] According to another embodiment of the present invention, an adaptive InSAR network adjustment system for monitoring surface deformation is also provided, the system comprising:

[0131] The effective homonym extraction module is used to obtain the coherence values ​​of each pixel in the overlapping area of ​​adjacent images. It uses a multi-level adaptive coherence mask to perform hierarchical screening of each pixel in the overlapping area, and determines the effective pixels based on the consistency condition between the coherence values ​​and the statistical features of the local neighborhood, and extracts the set of effective pixel coordinates.

[0132] The joint optimization module for adjustment model and regularization parameters is used to construct a multi-objective optimization function based on the effective set of pixel coordinates, which includes weighted difference residuals, in-orbit fitting residuals and regularization terms, and to jointly optimize the preset candidate model structure and regularization intensity to determine the optimal model structure and regularization parameters for each overlapping area.

[0133] The dynamic robust alternating solution module is used to construct a mapping model based on the effective pixel coordinate set under the constraints of the optimal model structure and regularization parameters. It uses an alternating least squares framework to iteratively update the mapping model parameters and dynamically adjusts the observation weight of each pixel to suppress the influence of outliers, thereby obtaining the mapping model parameters of each image.

[0134] The spatial residual surface correction module is used to calculate the residual value of the overlapping area using the mapping model parameters, perform spatial surface fitting on the residual value of the overlapping area using random consistency sampling, and adaptively determine the feasible region of the parameters of the surface fitting model according to the geometric shape of the overlapping area to obtain the spatial surface parameters.

[0135] The global consistency adjustment module integrates the mapping model parameters and spatial surface parameters into a relative observation vector. It uses the observation quality evaluation index of the overlapping area as the weight of the observation edge and incorporates all observation edges into the sparse matrix adjustment framework to solve for the optimal absolute correction parameters of the global image, thereby obtaining a consistent stitching result with the central image as the reference and the global error minimized.

[0136] It should be noted that this invention uses a province in central China as a test area to verify the adaptive InSAR network adjustment method proposed in this invention.

[0137] For example, a certain province in central China is located in the core zone of the Loess Plateau (34.34°N-40.44°N, 110.14°E-114.33°E) on the east bank of the middle reaches of the Yellow River and the west side of the North China Plain. The Taihang Mountains and Lüliang Mountains together form a regional uplifted mountain system on the east and west sides, while a goose-shaped fault basin has developed in the central part. As an important energy base, its coalfields are distributed in a belt-like pattern along the fault zones, longer in the north and south and narrower in the east and west. Long-term high-intensity mining has induced various types of complex deformation fields, including goaf subsidence, ground fissures, and landslides, exhibiting significant spatiotemporal heterogeneity.

[0138] The data used in this invention include:

[0139] (1) Sentinel-1A SAR dataset: 1817 C-band Sentinel-1A (Sentinel-1A) up-orbit SLC data provided by the European Space Agency, spanning from March 2017 to December 2023, operating in wide-swath interferometry (IW) mode, with VV polarization, covering 11 frames across three orbits: Path11, Path113, and Path40. The data is organized into an interferometric pair network using a strategy of connecting three interferograms of adjacent dates, with a multi-look ratio of range:azimuth = 8:2.

[0140] (2) Auxiliary data: 30m resolution SRTM DEM data, used to remove the terrain phase in time-series differential interferometry; European Space Agency precise orbit determination ephemeris data, used to correct satellite orbit errors; GACOS atmospheric correction data (spatial resolution 90m, temporal resolution 1 minute), used to correct atmospheric delay errors.

[0141] Since 2021, some orbital imagery has been missing, making it impossible to calculate the annual deformation rate. Therefore, the imagery is processed by combining images from three time periods, such as... Figure 2 As shown: from March 2017 to December 2018, from January 2019 to December 2020, and from January 2021 to December 2023, three sets of initial time-series InSAR deformation rate results were obtained, where (a1)-(c1) are the uncorrected results for the corresponding time periods, and (a2)-(c2) are the corrected results of the present invention for the corresponding time periods. Figure 3 By comparing the uncorrected results with the accuracy statistics of the method of the present invention in terms of root mean square error (RMSE) and standard deviation (STD), it can be seen that the method of the present invention has a significant improvement in accuracy. Figures 4-6 All are in Figure 2 Results of local area detail consistency analysis at the profile line in three time periods ( Figure 4 The corresponding period is from March 2017 to December 2018; Figure 5 The corresponding period is from January 2019 to December 2020; Figure 6 The corresponding time period is January 2021 to December 2023, where (a1)-(i1) represent the deformation rate of adjacent images at the profile line after correction; (a2)-(i2) represent the kernel density results of the deformation rate of adjacent images before and after correction; and (a3)-(i3) represent the histogram distribution of the pixel-level rate difference between adjacent images before and after correction. It can be seen that after correcting the original temporal InSAR deformation rate results for the three time periods, the deformation rate of adjacent images has good consistency and is significantly improved compared with the original uncorrected results.

[0142] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A wide-area deformation adaptive InSAR network adjustment method, characterized in that, include: S1. Obtain the coherence values ​​of each pixel in the overlapping area of ​​adjacent image pairs, use a multi-level adaptive coherence mask to perform layered screening of each pixel in the overlapping area, and determine the effective pixels based on the consistency condition between the coherence values ​​and the statistical features of the local neighborhood, and extract the set of coordinates of the effective pixels. S2. Based on the effective pixel coordinate set, construct a multi-objective optimization function that includes weighted difference residuals, in-orbit fitting residuals and regularization terms, and jointly optimize the preset candidate model structure and regularization intensity to determine the optimal model structure and regularization parameters for each overlapping area. S3. Under the constraints of the optimal model structure and regularization parameters, a mapping model is constructed based on the set of effective pixel coordinates. The mapping model parameters are iteratively updated using an alternating least squares framework, and the observation weights of each pixel are dynamically adjusted to suppress the influence of outliers, thus obtaining the mapping model parameters for each image. S4. Calculate the residual value of the overlapping area using the mapping model parameters, perform spatial surface fitting on the residual value of the overlapping area using random consistency sampling, and adaptively determine the feasible region of the parameters of the surface fitting model according to the geometric shape of the overlapping area to obtain the spatial surface parameters. S5. Integrate the mapping model parameters and spatial surface parameters into a relative observation vector, use the observation quality evaluation index of the overlapping area as the observation edge weight, and incorporate all observation edges into the sparse matrix adjustment framework to solve for the optimal absolute correction parameters of the global image, thereby obtaining a consistent stitching result with the central image as the reference and the global error minimized.

2. The wide-area deformation adaptive InSAR network adjustment method according to claim 1, characterized in that, S1 includes: The local terrain roughness of each pixel coordinate is extracted using a digital elevation model. The sliding convolution window size of each pixel is adaptively adjusted according to the local terrain roughness. The largest window size is used in flat areas where the terrain roughness is less than or equal to a first preset threshold, and the smallest window size is used in steep slopes or densely built-up areas where the terrain roughness is greater than or equal to a second preset threshold. When the terrain roughness is between the first and second preset thresholds, the window size is determined by linear mapping. Based on a preset initial high coherence threshold, the coverage and local median of high coherence pixels within the sliding convolution window centered on each pixel are calculated. Three levels of screening are set: strict screening mode, medium screening mode and lenient screening mode. Each screening mode corresponds to different coherence threshold, coverage threshold and local median threshold. In the initial stage, the threshold combination corresponding to the strict screening mode is used to determine each pixel in the overlapping area pixel by pixel. A cell is marked as a valid cell when its coherence value is greater than or equal to the coherence threshold of the current strict filtering mode, the coverage within the window containing the cell is greater than or equal to the coverage threshold of the current strict filtering mode, and the local median within the window containing the cell is not lower than the local median threshold of the current strict filtering mode. The system counts the number of valid pixels marked in the current strict filtering mode. If the number of valid pixels reaches the preset threshold, the filtering is terminated; otherwise, it switches to the lenient filtering mode. The overlapping area is re-evaluated pixel by pixel using the threshold combination corresponding to the lenient filtering mode until the number of valid pixels reaches the preset valid point threshold or all filtering modes have been traversed. Extract the coherence time series of each pixel marked as valid by each screening mode on all interference pairs, and calculate the time series stability index measured by median absolute deviation. All pixels that have been selected through a three-level stratified screening process and simultaneously meet the coherence threshold, coverage threshold, and local median threshold conditions, and have been adjusted for temporal stability constraints, will be considered as valid pixels for subsequent network adjustment, and the set of valid pixel coordinates will be output.

3. The wide-area deformation adaptive InSAR network adjustment method according to claim 1, characterized in that, S2 includes: For each pair of adjacent images, a candidate model set and multiple preset regularization intensities are defined, and a multi-objective optimization function vector containing weighted difference residuals, in-orbit fitting residuals, and regularization terms is constructed; the candidate model set includes constant models and linear models; The constant model and the linear model are combined with each preset regularization intensity to form candidate parameter combinations. Each candidate parameter combination is then input into the alternating least squares framework for iterative solution to obtain the weighted difference residuals and in-track fitting residuals corresponding to each candidate parameter combination. The weighted difference residuals and in-orbit fitting residuals corresponding to each candidate parameter combination, together with the corresponding regularization term values, are used as the values ​​of each component of the multi-objective optimization function vector. The values ​​of each component of the multi-objective optimization function vector are normalized by range to eliminate dimensional differences. Pareto optimization is used to eliminate inferior solutions in all evaluation dimensions to construct a set of non-dominated solutions. Calculate the normalized Euclidean distance from each solution in the non-dominated solution set to the ideal origin, and perform spatial clustering on the non-dominated solutions. Combine candidate models with regularization strengths that are close in distance in the target space into the same scheme cluster. Use consensus distance instead of the single normalized Euclidean distance as the selection criterion, and select the cluster with the smallest consensus distance as the global optimal model structure and regularization parameter. An absolute precision tolerance is set. When the weighted difference residual of the constant term model under the optimal regularization strength is better than the absolute precision tolerance, it is determined that the data quality of the current overlapping area is good and the error is gradual. No more linear terms are introduced, and the constant model is used as the final selected optimal model structure.

4. The wide-area deformation adaptive InSAR network adjustment method according to claim 3, characterized in that, The formula for constructing a multi-objective optimization function vector containing weighted difference residuals, in-orbit fitting residuals, and regularization terms is as follows: ； In the formula, Represents the set of candidate models; Indicates transpose; Represents the regularization term; This represents the weighted difference residual; Indicates the first In-orbit fitting residuals of the image; Indicates the first In-orbit fitting residuals of the image; This represents a vector of multi-objective optimization functions.

5. The wide-area deformation adaptive InSAR network adjustment method according to claim 1, characterized in that, S3 includes: Based on the effective set of pixel coordinates, a mapping model containing scale factor, offset term and spatial first-order term parameters is constructed for each pixel in the overlapping area of ​​each image. The latent variable is used to characterize the potential real surface deformation rate of the overlapping area, and the observation weight is used to measure the reliability of the observation value of each pixel. A weighted least squares objective function is constructed. An iterative solution is performed using an alternating least squares framework. Latent variables are fixed, and the estimated values ​​of the mapping model parameters are obtained by weighted least squares. The updated values ​​of the latent variables are obtained by fixing the mapping model parameters. The observation residuals under the current mapping model parameters are calculated, and the observation weights of each pixel are updated according to the dynamic robust weight function. Repeat the process of solving for the estimated values ​​of the mapping model parameters, solving for the updated values ​​of the latent variables, and calculating the observation residuals until the rate of change of the mapping model parameters is lower than the preset convergence threshold or the maximum number of iterations is reached, and output the mapping model parameters of each image.

6. The wide-area deformation adaptive InSAR network adjustment method according to claim 5, characterized in that, The formula for constructing the weighted least squares objective function is as follows: ； In the formula, This represents the weighted objective function; Represents polynomial parameters; Represents a pixel Observation weights of coordinates; Indicates the first Image in pixels The original observed deformation rate of the coordinates; Indicates the pixel to be estimated Latent variables of coordinates; Represents the mapping model; This represents the set of valid cell coordinates.

7. The wide-area deformation adaptive InSAR network adjustment method according to claim 5, characterized in that, The process of updating the observation weights of each pixel according to the dynamic robust weight function includes: After the mapping model parameters are updated in the current iteration, the difference between the original observed deformation rate of each pixel in the overlapping area and the predicted value of the current mapping model is calculated to obtain the observation residual of each pixel. Robust standard deviation estimates of the residual sequences in the overlapping region are calculated by combining the median absolute deviation, and standardized residual statistics are constructed. Based on the magnitude of the standardized residuals, the observation data is automatically divided into a retention zone, a deweighting zone, and a removal zone. In the retention zone, the observation weights remain unchanged. In the deweighting zone, the observation weights are reduced. In the removal zone, the observation weights are reset to zero, resulting in the updated observation weights for each pixel.

8. The wide-area deformation adaptive InSAR network adjustment method according to claim 7, characterized in that, The formula for constructing the standardized residual statistic is as follows: ； In the formula, Represents standardized residual statistics; Represents the robust standard deviation estimate of the residual sequence in the overlapping region; Represents the median operator; Represents the normality uniformity factor; This represents the observation residual.

9. The wide-area deformation adaptive InSAR network adjustment method according to claim 1, characterized in that, S4 includes: Based on the difference between the original observed deformation rate and the predicted value of the mapping model in the overlapping area of ​​each image, the residual value of each pixel in the overlapping area is calculated. The random consistency sampling method is adopted to randomly select the smallest sample subset from the residual values ​​for spatial surface fitting, and the number of interior points of each candidate surface within the preset residual tolerance is counted. The surface coefficient set with the largest number of interior points is taken as the optimal surface coefficient set. Calculate the aspect ratio of the overlapping area, and adaptively determine the feasible region of the parameters of the spatial surface fitting model based on the comparison between the aspect ratio and the preset aspect ratio threshold. By using the selected set of interior points and constrained by the parameter feasible region, the parameters of the spatial surface are obtained through least-squares refitting.

10. The wide-area deformation adaptive InSAR network adjustment method according to claim 1, characterized in that, S5 includes: For each pair of adjacent images with effective overlap, calculate the difference between the mapping model parameters of the two images and the difference between the spatial surface parameters of the two images, and use the combination of the difference between the mapping model parameters and the difference between the spatial surface parameters as the relative observation vector of the observation edge. Calculate the observation edge weight for each overlapping region, where the observation edge weight is the product of the proportion of random consistent sampling points within the overlapping region and the average coherence of effective pixels. Define all images as adjustment network nodes, define image pairs with effective overlap as effective edges, and solve the relative observation equations of all effective observation edges to construct the global relative observation matrix. In the absence of external control points, the spatial topological center image is selected as the absolute reference benchmark and constraints are added. A weighted least squares objective function with additional constraints is constructed. By constructing the augmented method equation and inverting it, the optimal absolute correction parameters for all images in the entire domain are obtained, resulting in a consistent stitching result with the center image as the benchmark and the global error minimized.

Images