Method for screening of interferometric image pairs and use thereof

Abstract

The application relates to the technical field of image processing, in particular to an interferometric image pair screening method and application thereof. The method comprises the following steps: taking the collection time of a SAR image in an interferometric image pair as a node, and taking the connection line between the interferometric image pairs as an edge to construct an interferometric baseline network undirected graph; taking the sequence of the collection time as a sequence, and sequentially eliminating all target edges smaller than a preset threshold value on each graph node in the interferometric baseline network undirected graph, wherein the preset threshold value is the average value of the coherence coefficients of all edges connected to the graph node; and taking the remaining interferometric image pairs in the interferometric baseline network undirected graph after the target edges are eliminated as screened target interferometric image pairs. The application aims to solve the problem of how to obtain high-precision InSAR.

Description

Technical Field

[0001] This application relates to the field of image processing technology, and in particular to a method for screening interference image pairs and its application. Background Technology

[0002] InSAR (Interferometric Synthetic Aperture Radar) technology is an active microwave remote sensing system. Its basic working principle is to use a radar sensor to emit electromagnetic waves into the ground, and then image the backscattered echo signals from ground objects. Based on the phase information between the SAR sensor and the ground target recorded by the SAR (Synthetic Aperture Radar) images, SAR images of the same area with repeated orbits are interferometrically superimposed, and the phase differences are compared to infer the surface elevation or surface deformation information. InSAR can stack SAR images from multiple times to form multiple interferograms, obtain time-series deformation information, and capture the temporal evolution of the surface. It has wide applications in landslide identification, glacier displacement monitoring, mining area subsidence monitoring, and urban deformation monitoring.

[0003] Changes in environmental factors such as surface vegetation growth, human activities, water vapor variations, and landform changes can introduce decoherent noise between SAR images, affecting the quality of interferometric pairs. The quality of interferometric pairs is crucial for high-precision inversion of surface deformation. Selecting a poor-quality or fewer interferometric pairs will reduce the accuracy of InSAR deformation inversion, while selecting a larger number of interferometric pairs will increase computational costs and may not necessarily yield high-precision surface deformation results.

[0004] Therefore, this application aims to propose a selection method for filtering high-quality interferometric pairs from a large amount of interferometric image pair data, thereby obtaining high-precision InSAR results. Summary of the Invention

[0005] The main objective of this application is to provide a method for screening interferometric image pairs, aiming to solve the problem of how to obtain high-precision InSAR.

[0006] To achieve the above objectives, this application provides a method for screening interference image pairs, the method comprising:

[0007] S10, using the acquisition time of SAR images in the interferometric image pair as nodes and the lines connecting the interferometric image pairs as edges, constructs an undirected graph of the interferometric baseline network;

[0008] S20, in the order of the acquisition time, remove all target edges less than a preset threshold on each graph node in the undirected graph of the interferometric baseline network, wherein the preset threshold is the average coherence coefficient of all edges connected to the graph node.

[0009] S30, the remaining interferometric pairs in the undirected graph of the interferometric baseline network after removing the target edges are taken as the filtered target interferometric pairs.

[0010] Optionally, before step S10, the following steps are further included:

[0011] S40, Calculate the atmospheric delay phase value in the initial interferometric image pair;

[0012] S50, the atmospheric delay phase value is subtracted from the initial interferometric image to obtain the interferometric image pair as the undirected graph for constructing the interferometric baseline network.

[0013] Optionally, in S40, the expression for calculating the atmospheric delay phase value is:

[0014]

[0015] In the formula, ZTDk is the atmospheric delay phase value, representing the total zenith delay of the combination of the vertical stratification component and the horizontal turbulence component; k is the coordinate position; T represents the turbulence signal; x k It is the station coordinate vector in the local geocentric coordinate system; L0 represents the stratification component delay of sea level; This represents the remaining unmodeled residuals, including unmodeled stratification and turbulence signals; For height scale;

[0016] Among them, height scale The calculation expression is:

[0017]

[0018] In the formula, Indicates height scale, Indicates the lowest altitude. Indicates the highest altitude.

[0019] Optionally, the coherence coefficient The calculation expression is:

[0020]

[0021] In the formula, This represents the amplitude and phase of the main image. Indicates the amplitude and phase of the secondary image. express The complex conjugate of is used to calculate coherence.

[0022] Furthermore, to achieve the above objectives, this application also provides a method for inverting surface deformation, which includes the following steps:

[0023] Obtain the target interferometric image pair after screening based on the aforementioned interferometric image pair screening method;

[0024] Select the target interferogram generated by superimposing SAR images corresponding to two different acquisition times from the target interferogram pair, and calculate the target interferogram pair corresponding to the target interferogram. The target interferogram pair is composed of the cumulative deformation of the radar line of sight, the residual terrain phase in the differential interferogram, the atmospheric delay phase, and the sum of the decoherent noise.

[0025] Calculate the velocity vector based on the target interferometric image pair and the acquisition time;

[0026] The minimum norm solution of the velocity vector is calculated using the singular value decomposition method. The minimum norm solution of the velocity vector is then integrated to obtain the deformation of the target interferogram.

[0027] Based on the numerical range of the deformation, the overlapping area, shadow area, perspective contraction area and undeformed area in the SAR image are determined. The overlapping area and the shadow area are masked to complete the surface deformation inversion.

[0028] Optionally, the expression for calculating the velocity vector is:

[0029]

[0030] In the formula, This represents the i-th target interferogram. This indicates the acquisition time corresponding to the i-th target interferogram. This represents the (i-1)th target interferogram. The acquisition time corresponding to the (i-1)th target interferogram.

[0031] Optionally, the step of calculating the minimum norm solution of the velocity vector using singular value decomposition and integrating the minimum norm solution of the velocity vector to obtain the deformation corresponding to the target interferogram includes:

[0032] Let the interference image pair value of the i-th interferogram be... The representation is:

[0033]

[0034] In the formula, T A Indicates the acquisition times A and T. B Indicates the acquisition time B, v pIt represents the rate of surface deformation.

[0035] The integrals of each time period over the time intervals of the master and slave images are represented by an M×N matrix:

[0036]

[0037] The generalized inverse matrix of matrix B is obtained by using singular value decomposition, and the minimum norm solution of the velocity vector is calculated based on the generalized inverse matrix.

[0038] Integrating the minimum norm solution of the velocity vector yields the deformation of the target interferogram.

[0039] In addition, to achieve the above objectives, this application also provides an application of the interferometric image pair screening method described above in surface deformation inversion.

[0040] In addition, to achieve the above objectives, this application also provides a computer system, the computer system comprising: a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the computer program is executed by the processor, it implements the steps of the interferometric image pair screening method or the surface deformation inversion method as described in any of the preceding claims.

[0041] In addition, to achieve the above objectives, this application also provides a computer-readable storage medium storing a computer program, which, when executed by the processor, implements the steps of the interferometric image pair screening method or the surface deformation inversion method as described in any of the preceding claims.

[0042] This application has at least the following beneficial effects:

[0043] 1. An undirected graph of the interferometric baseline network is constructed based on graph theory, and all target edges less than a preset threshold on each graph node in the undirected graph of the interferometric baseline network are removed. The preset threshold is the average coherence coefficient of all edges connected to the graph node. This method not only removes low-quality interferometric pairs but also ensures the integrity of the interferometric baseline network.

[0044] 2. Correct the atmospheric delay phase in the interferometric image pair data to improve the accuracy of the InSAR data structure;

[0045] 3. For the selected interferometric image pairs, surface deformation information is inverted based on SBAS-InSAR technology. During the inversion process, geometric distortion in the SAR image is taken into account. Areas with geometric distortion are visualized, marked, and then masked to ensure the accuracy of the InSAR results. Attached Figure Description

[0046] Figure 1 This is a flowchart illustrating the first embodiment of the method for screening interference image pairs in this application;

[0047] Figure 2 This is a schematic diagram of an undirected graph involved in an embodiment of this application;

[0048] Figure 3 This is a schematic diagram of a directed simple graph involved in an embodiment of this application;

[0049] Figure 4 This is a schematic diagram of the undirected graph of the interference baseline network before the target edge is removed, as described in the embodiments of this application.

[0050] Figure 5 This is a schematic diagram of the undirected graph of the interference baseline network after removing the target edges, as described in the embodiments of this application.

[0051] Figure 6 This is a schematic diagram of the phase standard deviation before and after GACOS correction in the embodiments of this application;

[0052] Figure 7 This is a radar visibility geometric distortion identification diagram related to the embodiments of this application;

[0053] Figure 8 This is an interferometric image pair connection diagram based on the fully connected method involved in the embodiments of this application;

[0054] Figure 9 This is an interferometric image pair connection diagram based on the average coherence coefficient threshold method involved in the embodiments of this application;

[0055] Figure 10 This is an interferometric image pair connection diagram based on the small baseline set method involved in the embodiments of this application;

[0056] Figure 11 This is an interferometric image pair connection diagram based on the interferometric image pair screening method proposed in this application, which is involved in the embodiments of this application.

[0057] Figure 12 This is a schematic diagram of the hardware operating environment of the computer system involved in the embodiments of this application.

[0058] The realization of the purpose, functional features and advantages of this application will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0059] To better understand the above technical solutions, exemplary embodiments of this disclosure will be described in more detail below with reference to the accompanying drawings. While exemplary embodiments of this disclosure are shown in the drawings, it should be understood that this disclosure can be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of this disclosure to those skilled in the art.

[0060] First Embodiment

[0061] To reduce the impact of decoherence noise on InSAR deformation inversion, extensive research has been conducted on the optimal selection of InSAR interferometric pairs. Currently, there are three main types of methods for optimizing InSAR interferometric pair selection: The first type relies on the past experience of InSAR data processors, manually removing low-quality interferometric pairs through visual interpretation, retaining only the more ideal pairs. This method is time-consuming, labor-intensive, and inefficient, making it difficult to meet the needs of processing large volumes of long-term InSAR data. Furthermore, this method depends on the prior knowledge of the data processors, introducing a degree of subjectivity. The second type uses the average coherence coefficient of all interferometric pairs as a threshold for selection, i.e., the average coherence coefficient threshold method. Existing research indicates that this method achieves relatively ideal results in urban areas where coherence variation is minimal. However, in mountainous areas with complex terrain, coherence often varies seasonally; summer rainfall and high vegetation cover result in lower coherence, while winter is drier with lower vegetation cover, resulting in higher coherence. Due to the large difference in coherence, averaging with this type of method may lead to breaks in the interferometric baseline network, making subsequent surface deformation information inversion impossible. The third type of method involves threshold constraints on temporal and spatial baselines, i.e., the short-temporal baseline method. The effectiveness of this method depends on the setting of the temporal and spatial baseline thresholds. Setting the threshold too high increases temporal and spatial decorrelation, producing redundant low-coherence image pairs; setting the threshold too low reduces the number of connected interferometric image pairs, lowering the accuracy of surface deformation information inversion. Therefore, this type of method has significant instability; different regions may have different optimal thresholds, and the setting of temporal and spatial baseline thresholds relies heavily on past experience.

[0062] In this embodiment, to address the shortcomings of the aforementioned InSAR interferometric image pair optimization selection method, a graph theory-based interferometric image pair selection method is proposed:

[0063] Reference Figure 1 In this embodiment, the method for screening interference image pairs includes the following steps:

[0064] S10, using the acquisition time of SAR images in the interferometric image pair as nodes and the lines connecting the interferometric image pairs as edges, constructs an undirected graph of the interferometric baseline network;

[0065] In this embodiment, the interferometric image pair consists of at least two SAR images of the same area at different times and / or different perspectives. Each SAR image is recorded with a timestamp when it is acquired. The timestamp is used as the acquisition time of the SAR image and is used as a graph node. The lines connecting the interferometric image pairs are used as edges to construct an undirected graph of the interferometric baseline network.

[0066] It should be noted that in an undirected graph of an interferometric baseline network, a target node that is connected to another node is considered an interferometric pair. The connection between a target node and other nodes can be one or more. For example, if nodes B and C are connected to node A, then nodes A and B are an interferometric pair, and nodes A and C are also an interferometric pair.

[0067] It should be noted that graph theory is used in this embodiment to screen interferometric pairs. Graph theory is a branch of mathematics that uses graphs to study the objective world. A graph in graph theory consists of some points and lines connecting these points. Points represent things or objects, and the lines connecting points are called edges, which represent the relationships between things or objects (in this embodiment, they are used to represent interferometric pairs). Each edge has a corresponding value called the edge weight, which represents the importance between two points (i.e., the importance of the interferometric pair).

[0068] In graph theory, a graph G=(V,E) refers to two sets (V,E), where V is the set of vertices and E is the set of edges. A graph is a data structure composed of the relationships between the vertex set V and the edge set E. (See also...) Figure 2 The diagram illustrates an undirected graph. Assume there is a line connecting x and y in graph G, denoted as x→y, where point x is called the starting point of the edge and point y is called the ending point. If both x→y and y→x exist simultaneously, these two edges are merged into one, and graph G is then called an undirected graph.

[0069] Additionally, refer to Figure 3 The diagram shown is a schematic of a directed simple graph. Suppose that in graph G, there is only one edge from any starting point x to y, and there are no edges from itself to itself. Then graph G is also called a simple graph or a directed simple graph.

[0070] S20, in chronological order of the acquisition time, remove all target edges less than a preset threshold on each graph node in the undirected graph of the interferometric baseline network, wherein the preset threshold is the average coherence coefficient of all edges connected to the graph node.

[0071] In this embodiment, in the constructed undirected graph of the interferometric baseline network, each graph node corresponds to an interferometric image pair, and each interferometric image pair has its corresponding acquisition time. The graph nodes are sorted according to the order of acquisition time, and the graph node associated with the earliest interferometric image pair is marked as the first node. Starting from the first node, the average weight of all edges connected to the node is used as a threshold to remove all edges on the node that are below the threshold for optimization and filtering.

[0072] In this embodiment, each edge of the undirected graph of the interferometric baseline network corresponds to a preset weight value, which is the coherence coefficient of the graph node.

[0073] The coherence coefficient is a standardized covariance function with a value range of [0,1]. 0 indicates complete decoherence, and 1 indicates complete coherence. The larger the coherence coefficient, the higher the coherence and the more similar the echo signals between two SAR images. The interferometric pairs in the resulting interferogram can accurately reflect the distance difference between the two echoes and can obtain more accurate information on surface deformation.

[0074] Alternatively, coherence reflects the degree of similarity between two radar echo signals. The coherence r of two zero-mean circular Gaussian complex random signals y1 and y2 can be defined as:

[0075]

[0076] Ideally, the expected value E{} in the above formula can be obtained by calculating the pooled average for each pixel of a large number of interferometric image pairs acquired simultaneously under the same conditions. However, in reality, this is usually difficult to achieve. Therefore, assuming that the random process is stationary and ergodic within a window of N pixels, the spatial average of these N pixels can replace the pooled average to obtain the coherence coefficient. :

[0077]

[0078] The threshold is obtained by averaging the correlation coefficients corresponding to each edge. .

[0079] S30, the remaining interferometric pairs in the undirected graph of the interferometric baseline network after removing the target edges are taken as the filtered target interferometric pairs.

[0080] As an example, refer to respectively Figure 4 and Figure 5 The diagrams shown depict the undirected graph of the interferometric baseline network before and after removing the target edge. Starting from node A, there are a total of 6 edges connected to node A, with each edge having a weight of C. AB =0.45, CAC =0.22, C AD =0.27, C AE =0.39, C AF =0.42, C AG =0.28, calculated C avg =0.29, therefore, edges AC, AD, and AG are removed, that is, the three interference pairs AC, AD, and AG are removed. The remaining AB, AE, and AF are all selected as target interference pairs, and then proceed to the next node B, until the last node G ends.

[0081] It is worth noting that the interferometric image pair screening process involved in this embodiment does not involve the removal of nodes (i.e., the SAR image itself is not removed). What is removed are the low-quality interferometric image pairs among the various interferometric image pairs that interfere with the node, so as to ensure the integrity of the interferometric baseline network.

[0082] In the technical solution provided in this embodiment, the lines connecting interferometric image pairs are used as edges of the graph, and the coherence coefficients of the interferometric image pairs are used as edge weights. Starting from the first node, the average weight of all edges connected to that node is used as a threshold to remove all edges on that node that are below the threshold, until the last node. In the method of this application, each node is connected by edges, which not only removes low-quality interferometric image pairs but also ensures the integrity of the interferometric baseline network.

[0083] Second Embodiment

[0084] Based on the first embodiment, in this embodiment, InSAR technology is similar to geodetic techniques such as Global Navigation Satellite System (GNSS) and Very Long Baseline Interferometry (VLBI). They all use radar signals to acquire information. When radar signals pass through the atmosphere, a delay occurs, causing the electromagnetic wave carrier phase to shift. This is called the atmospheric delay effect.

[0085] According to atmospheric stratification theory, InSAR atmospheric delay can be divided into two categories: ionospheric atmospheric delay and tropospheric atmospheric delay. Existing studies have shown that the ionospheric effect on the phase delay of electromagnetic waves is inversely proportional to wavelength. For example, L-band ALOS / PALSAR data (wavelength 23.6 cm) experiences a greater phase delay in the ionosphere than C-band Sentinel-1 data (wavelength 5.6 cm). Furthermore, using multi-temporal InSAR data processing can also mitigate the ionospheric delay effect to some extent. Therefore, this embodiment considers the atmospheric delay effect of the troposphere on InSAR but does not consider the ionospheric atmospheric delay effect.

[0086] Furthermore, the troposphere mainly consists of two components: dry air and water vapor. Dry air undergoes contraction under changes in atmospheric pressure and temperature, resulting in a delay effect on radar signals. However, the delay effect caused by dry air is relatively stable and its variation is relatively small, generally negligible. Therefore, water vapor, as the main factor causing atmospheric delay in InSAR, is highly variable and difficult to predict. Water vapor varies with time and space, and radar signals undergo refraction when passing through water vapor, causing changes in the signal propagation path and direction, thus producing phase delay errors and affecting the accuracy of InSAR results.

[0087] Therefore, this embodiment introduces the GACOS atmospheric model to correct the atmospheric delay phase in the data. Specifically:

[0088] S40, Calculate the atmospheric delay phase value in the initial interferometric image pair;

[0089] S50, the atmospheric delay phase value is subtracted from the initial interferometric image to obtain the interferometric image pair as the undirected graph for constructing the interferometric baseline network.

[0090] Furthermore, in S40, the expression for calculating the atmospheric delay phase value is:

[0091]

[0092] In the formula, ZTDk is the atmospheric delay phase value, representing the total zenith delay of the combination of the vertical stratification component and the horizontal turbulence component; k is the coordinate position; T represents the turbulence signal; x k It is the station coordinate vector in the local geocentric coordinate system; L0 represents the stratification component delay of sea level; This represents the remaining unmodeled residuals, including unmodeled stratification and turbulence signals; For height scale;

[0093] Among them, height scale The calculation expression is:

[0094]

[0095] In the formula, Indicates height scale, Indicates the lowest altitude. Indicates the highest altitude.

[0096] It is worth noting that not all GACOS products can effectively correct atmospheric delay errors for every interferometric image pair. On the contrary, GACOS correction may introduce additional noise. This is because GACOS relies on meteorological data from weather forecast models. If there is intense atmospheric activity during a certain period, the forecast model may struggle to capture all atmospheric changes, leading to the introduction of new noise in InSAR atmospheric delay correction. Therefore, the standard deviation (STD) of the interferometric image pair should also be used as an indicator to evaluate the effectiveness of atmospheric correction.

[0097] Exemplarily, in one specific embodiment, reference is made to Figure 6 The diagram showing the phase standard deviation before and after GACOS correction illustrates that, after GACOS correction on 435 interferometric pairs, the STD decreased in 330 pairs (75.8% of all pairs), while the STD increased in 105 pairs (24.2% of all pairs), with a maximum increase of 0.08 rad. Overall, GACOS effectively mitigates errors caused by atmospheric delay.

[0098] Third Embodiment

[0099] This embodiment provides a method for further inversion of surface deformation information on optimized interferometric image pairs, the steps of which are as follows:

[0100] Step S100: Obtain the target interference image pair after screening based on the aforementioned interference image pair screening method;

[0101] Step S200: Select the target interferogram generated by superimposing SAR images corresponding to two different acquisition times from the target interferogram pair, and calculate the target interferogram pair corresponding to the target interferogram. The target interferogram pair is composed of the cumulative deformation of the radar line of sight, the residual terrain phase in the differential interferogram, and the sum of decoherent noise.

[0102] Specifically, for the i-th (i=1,2,…,N+1) target interferogram generated by interferometric overlay of two SAR images at times TA and TB, the interferometric pair of pixels with azimuth coordinate x and range coordinate y can be represented as:

[0103]

[0104] In the formula, This represents the cumulative deformation along the radar line-of-sight direction. The residual topographic phase in the differential interferogram. For atmospheric delayed phase, To remove coherent noise.

[0105] It should be noted that the atmospheric delay phase removed in the second embodiment is a global atmospheric effect, but local atmospheric residuals still exist in the interferometric image pair. Therefore, the formula for this interferometric image pair still includes the atmospheric delay phase. .

[0106] Step S300: Calculate the velocity vector based on the target interferometric image pair and the acquisition time;

[0107] Optionally, in order to obtain physically meaningful time-series deformation information, the phase in the above equation is represented as a velocity vector between two acquisition times. The product of time and time, i.e.:

[0108]

[0109] The velocity vector is obtained by sorting. The calculation expression is:

[0110]

[0111] In the formula, This represents the i-th target interferogram. This indicates the acquisition time corresponding to the i-th target interferogram. This represents the (i-1)th target interferogram. The acquisition time corresponding to the (i-1)th target interferogram.

[0112] Step S400: The minimum norm solution of the velocity vector is calculated using the singular value decomposition method. The minimum norm solution of the velocity vector is integrated to obtain the deformation corresponding to the target interferogram.

[0113] Furthermore, this step specifically includes:

[0114] Step S401, let the interference image pair value of the i-th interferogram be... The representation is:

[0115]

[0116] In the formula, T A Indicates the acquisition times A and T. B Indicates the acquisition time B, v p It represents the rate of surface deformation.

[0117] Step S402: Represent the integral of each time period over the time interval between the master and slave images using an M×N matrix:

[0118]

[0119] Step S403: Use singular value decomposition to calculate the generalized inverse matrix of matrix B, and calculate the minimum norm solution of the velocity vector based on the generalized inverse matrix;

[0120] Step S404: Integrate the minimum norm solution of the velocity vector to obtain the deformation corresponding to the target interferogram.

[0121] In this embodiment, optionally, since InSAR technology uses side-looking active transmission and reception of microwave signals for imaging, the signals transmitted by the SAR sensor in side-looking imaging are inevitably affected by terrain undulations, resulting in a difference between the recorded situation in the image and the actual surface situation. This difference is due to imaging geometry and is called geometric distortion of SAR images, including overlay, shadows, and perspective contraction. During the SAR sensor imaging process, even slight changes in elevation can produce large-scale distortions in the image. Among these, the local incident angle and the radar line-of-sight incident angle have the most significant effects. Therefore, the deformation includes the local incident angle and the radar line-of-sight incident angle.

[0122] Step S500: Based on the numerical range of the deformation, determine the overlapping area, shadow area, perspective contraction area and undeformed area in the SAR image, and mask the overlapping area and the shadow area to complete the surface deformation inversion.

[0123] In this embodiment, the local incident angle and the radar line-of-sight incident angle are used as deformation variables, and the SAR region is determined based on the numerical range of the two.

[0124] Optionally, a local incident angle can be set. The radar line of sight is at the angle of incidence. The numerical range in which it lies can be as follows:

[0125]

[0126] For example, in one specific implementation, the interferometric baseline network, constructed from 30 fully connected Sentinel-1A images of the study area from January 10, 2022 to December 24, 2022, is considered as an undirected graph. Combined with external DEM data, and based on the LSM (Layover and Shadow Map) algorithm and the R-index, radar visibility geometric distortion identification is performed on the study area, resulting in the following... Figure 7 As shown.

[0127] from Figure 7It can be seen that the shadow and layered areas of the Sentinel-1A satellite in the study area are mainly concentrated in the eastern and northern regions, and are distributed in strips. The shadow accounts for 0.2% of the entire study area, and the layered areas account for 4.4%. In actual InSAR data processing, the deformation information obtained from the shadow and layered areas is unusable, while the foreshortening and good visibility areas are usable.

[0128] Fourth embodiment

[0129] Based on any of the foregoing embodiments, this embodiment verifies the effectiveness and accuracy of the aforementioned screening of interference image pairs.

[0130] In this embodiment, the method of this application is compared with the traditional fully connected method (without interference pair optimization), the average coherence coefficient threshold method, and the small baseline set method (temporal baseline set to 120 days, spatial baseline set to 50%, totaling 214 image pairs). The interference pair connectivity of the four methods is described in detail below. Figure 8-11 The diagrams shown are interferometric image pair connection diagrams based on the fully connected method, the average coherence coefficient threshold method, the small baseline set method, and the interferometric image pair connection diagram based on the interferometric image pair screening method proposed in this application.

[0131] from Figure 8-11 It can be seen that the average coherence coefficient thresholding method improves the overall coherence of the interferometric baseline network. However, due to the complex mountainous area of ​​the study area, the coherence varies significantly with the seasons. Interferometric image pairs in January, February, March, April, November, and December have higher coherence, while those in other months have lower coherence. After overall averaging, interferometric image pairs in May and June are almost completely eliminated, resulting in a broken interferometric baseline network that cannot be used for subsequent surface deformation inversion. The method proposed in this application treats the entire interferometric baseline network as a graph composed of nodes and edges from a graph theory perspective. The SAR image acquisition time is used as the node, the connection between interferometric image pairs is used as the edge, and the coherence coefficient of the interferometric image pair is used as the edge weight. Starting from the first node, the average weight of all edges connected to that node is used as the threshold to eliminate all edges below the threshold at that node, until the last node. In the method proposed in this application, each node is connected by edges, which both eliminates low-quality interferometric image pairs and ensures the integrity of the interferometric baseline network.

[0132] Since the average coherence threshold method cannot be used for subsequent surface deformation inversion, we continue to quantitatively compare the fully connected method, the small baseline set method, and the interferometric pair optimization method based on graph theory. We use surface deformation inversion coherence, effective interferogram ratio, and inversion error as indicators for comparison. InSAR technology captures and inverts surface deformation information by comparing the phase difference of radar images. Coherence quantifies the consistency of phase information between radar images; high surface deformation inversion coherence means good consistency of phase information, enabling the acquisition of more accurate surface deformation information. The effective interferogram ratio refers to the proportion of interferograms participating in the inversion calculation among all interferograms. A higher effective interferogram ratio means that more interferograms are involved in the inversion calculation within the monitoring area, which helps improve the reliability and accuracy of monitoring. RMSE refers to the error in the surface deformation inversion process; the lower the value, the higher the model fit and the more accurate the surface deformation information obtained. The average values ​​of all pixels for the three indicators of surface deformation inversion coherence, effective interferogram ratio, and inversion error RMSE for the three methods are shown in Table 1:

[0133] Table 1. Coherence of surface deformation inversion, effective interferogram ratio, and overall average RMSE for the three methods

[0134]

[0135] It is evident that the fully connected method exhibits the worst accuracy, with a surface deformation inversion coherence of 0.18, an effective interferogram ratio of 70.41%, and an RMSE of 3.52 rad. The fully connected method contains a large number of low-quality interferometric pairs, resulting in the lowest coherence and effective interferogram ratio, and the largest inversion error during surface deformation inversion. This demonstrates that more interferometric pairs do not necessarily guarantee reliable deformation inversion results and increase computational costs, highlighting the importance of optimized interferometric pair selection for InSAR data processing. Compared to the fully connected method, the small baseline set method improves surface deformation inversion coherence by 0.13, increases the effective interferogram ratio by 11.01%, and reduces the RMSE by 0.89 rad. This indicates that the small baseline set method effectively suppresses temporal and spatial decorrelation, thereby improving the accuracy of surface deformation inversion. The interferometric image pair optimization method based on graph theory proposed in this application achieves the highest accuracy, with a surface deformation inversion coherence of 0.40, an effective interferogram ratio of 88.13%, and an RMSE of only 2.31 rad. Compared to the method without interferometric image pair optimization, the surface deformation inversion coherence is improved by 0.22, the effective interferogram ratio is improved by 17.72%, and the RMSE is reduced by 1.21 rad. Compared to the small baseline set method, the surface deformation inversion coherence is improved by 0.09, the effective interferogram ratio is improved by 6.71%, and the RMSE is reduced by 0.32 rad. Furthermore, the number of interferometric image pairs involved in the inversion calculation is reduced by 40 pairs. In summary, the InSAR interferometric image pair optimization method based on graph theory proposed in this application can effectively suppress the influence of decoherence noise on InSAR deformation inversion, selecting as few high-quality interferometric image pairs as possible from a large amount of data to improve the reliability and accuracy of InSAR surface deformation inversion.

[0136] As one implementation scheme, Figure 12 This is a schematic diagram of the hardware operating environment of the computer system involved in the embodiments of this application.

[0137] like Figure 12As shown, the computer system may include: a processor 1001, such as a CPU; a memory 1005; a user interface 1003; a network interface 1004; and a communication bus 1002. The communication bus 1002 is used to enable communication between these components. The user interface 1003 may include a display screen or an input unit such as a keyboard; optionally, the user interface 1003 may also include a standard wired interface or a wireless interface. The network interface 1004 may optionally include a standard wired interface or a wireless interface (such as a Wi-Fi interface). The memory 1005 may be high-speed RAM or non-volatile memory, such as a disk drive. Optionally, the memory 1005 may also be a storage device independent of the aforementioned processor 1001.

[0138] Those skilled in the art will understand that Figure 12 The computer system architecture shown does not constitute a limitation on the computer system and may include more or fewer components than shown, or combine certain components, or have different component arrangements.

[0139] like Figure 12 As shown, the memory 1005, as a storage medium, may include an operating system, a network communication module, a user interface module, and computer programs. The operating system is a program that manages and controls the hardware and software resources of the computer system, as well as the operation of the computer programs and other software or programs.

[0140] exist Figure 12 In the computer system shown, the user interface 1003 is mainly used to connect to the terminal and communicate data with the terminal; the network interface 1004 is mainly used to communicate data with the backend server; and the processor 1001 can be used to call the computer program stored in the memory 1005.

[0141] In this embodiment, the computer system includes: a memory 1005, a processor 1001, and a computer program stored in the memory and executable on the processor, wherein:

[0142] When processor 1001 calls a computer program stored in memory 1005, it performs the following operations:

[0143] S10, using the acquisition time of SAR images in the interferometric image pair as nodes and the lines connecting the interferometric image pairs as edges, constructs an undirected graph of the interferometric baseline network;

[0144] S20, in chronological order of the acquisition time, remove all target edges less than a preset threshold on each graph node in the undirected graph of the interferometric baseline network, wherein the preset threshold is the average coherence coefficient of all edges connected to the graph node.

[0145] S30, the remaining interferometric pairs in the undirected graph of the interferometric baseline network after removing the target edges are taken as the filtered target interferometric pairs.

[0146] When processor 1001 calls a computer program stored in memory 1005, it performs the following operations:

[0147] S40, Calculate the atmospheric delay phase value in the initial interferometric image pair;

[0148] S50, the atmospheric delay phase value is subtracted from the initial interferometric image to obtain the interferometric image pair as the undirected graph for constructing the interferometric baseline network.

[0149] When processor 1001 calls a computer program stored in memory 1005, it performs the following operations:

[0150] Obtain the target interferometric image pair after screening based on the aforementioned interferometric image pair screening method;

[0151] Select the target interferogram generated by superimposing SAR images corresponding to two different acquisition times from the target interferogram pair, and calculate the target interferogram pair corresponding to the target interferogram. The target interferogram pair is composed of the cumulative deformation of the radar line of sight, the residual terrain phase in the differential interferogram, the atmospheric delay phase, and the sum of the decoherent noise.

[0152] Calculate the velocity vector based on the target interferometric image pair and the acquisition time;

[0153] The minimum norm solution of the velocity vector is calculated using the singular value decomposition method. The minimum norm solution of the velocity vector is then integrated to obtain the deformation of the target interferogram.

[0154] Based on the numerical range of the deformation, the overlapping areas, shadow areas, perspective contraction areas, and undeformed areas in the SAR image are determined. The overlapping areas and shadow areas are then masked to complete the surface deformation inversion. When processor 1001 calls the computer program stored in memory 1005, it performs the following operations:

[0155] Let the interference image pair value of the i-th interferogram be... The representation is:

[0156]

[0157] In the formula, T AIndicates the acquisition times A and T. B Indicates the acquisition time B, v p It indicates the rate of change in the shape of the Earth's surface.

[0158] The integrals of each time period over the time intervals of the master and slave images are represented by an M×N matrix:

[0159]

[0160] The generalized inverse matrix of matrix B is obtained by using singular value decomposition, and the minimum norm solution of the velocity vector is calculated based on the generalized inverse matrix.

[0161] Integrating the minimum norm solution of the velocity vector yields the deformation of the target interferogram.

[0162] Furthermore, those skilled in the art will understand that all or part of the processes in the methods of the above embodiments can be implemented by a computer program instructing related hardware. The computer program includes program instructions and can be stored in a storage medium, which is a computer-readable storage medium. The program instructions are executed by at least one processor in a computer system to implement the process steps of the embodiments of the above methods.

[0163] Therefore, this application also provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the steps of the interferometric image pair screening method or the surface deformation inversion method as described in the above embodiments.

[0164] The computer-readable storage medium can be any computer-readable storage medium capable of storing program code, such as a USB flash drive, portable hard drive, read-only memory (ROM), magnetic disk, or optical disk.

[0165] It should be noted that, since the storage medium provided in the embodiments of this application is the storage medium used to implement the methods of the embodiments of this application, those skilled in the art can understand the specific structure and variations of the storage medium based on the methods described in the embodiments of this application, and therefore will not be repeated here. All storage media used in the methods of the embodiments of this application fall within the scope of protection of this application.

[0166] Those skilled in the art will understand that embodiments of this application can be provided as methods, systems, or computer program products. Therefore, this application can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, this application can take the form of a computer program product embodied on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0167] This application is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of this application. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart... Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0168] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0169] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0170] It should be noted that any reference signs placed between parentheses in the claims should not be construed as limiting the claims. The word "comprising" does not exclude the presence of components or steps not listed in the claims. The word "a" or "an" preceding a component does not exclude the presence of a plurality of such components. This application can be implemented by means of hardware comprising several different components and by means of a suitably programmed computer. In a unit claim enumerating several means, several of these means may be embodied by the same item of hardware. The use of the words first, second, and third, etc., does not indicate any order. These words can be interpreted as names.

[0171] Although preferred embodiments of this application have been described, those skilled in the art, upon learning the basic inventive concept, can make other changes and modifications to these embodiments. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments as well as all changes and modifications falling within the scope of this application.

[0172] Obviously, those skilled in the art can make various modifications and variations to this application without departing from the spirit and scope of this application. Therefore, if such modifications and variations fall within the scope of the claims of this application and their equivalents, this application also intends to include such modifications and variations.

Claims

1. A method for screening interference image pairs, characterized in that, The method includes the following steps: S10, using the acquisition time of SAR images as nodes and the lines connecting interferometric image pairs as edges, constructs an undirected graph of the interferometric baseline network; S20, in chronological order of the acquisition time, remove all target edges less than a preset threshold on each graph node in the undirected graph of the interferometric baseline network, wherein the preset threshold is the average coherence coefficient of all edges connected to the graph node. S30, the remaining interferometric pairs in the undirected graph of the interferometric baseline network after removing the target edges are taken as the filtered target interferometric pairs.

2. The application of the interferometric image pair screening method as described in claim 1 in surface deformation inversion.

3. A computer system, characterized in that, The computer system includes: a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the computer program, when executed by the processor, implements the steps of the interferometric image pair screening method as described in claim 1.

4. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores a computer program that, when executed by a processor, implements the steps of the method for screening interferometric image pairs as described in claim 1.

Images