Coupling global graph-based difference enhancement SAR image change detection method

Abstract

The invention discloses a difference-enhanced SAR image change detection method based on a coupling global graph, and the method comprises the steps: firstly, dividing SAR images at two moments into pixel blocks, and selecting K pixel blocks with the highest similarity as the global adjacent points of the corresponding pixel blocks through the block similarity calculation, the method comprises the following steps: firstly, dividing an image block into a plurality of image blocks, constructing two corresponding global adjacent point sets # imgabs0 # and # imgabs1 #, then setting the divided image blocks as vertexes of a global graph, connecting each vertex with the respective global adjacent point set to construct edges of a global graph model, and then endowing each edge with a corresponding weight so as to complete the establishment of global graph models G1 and G2 at two moments; taking a union set according to the adjacent point sets of the two images, constructing a coupling adjacent point set # imgabs2, performing optimization processing on the coupling adjacent point set, and eliminating interference points in the point set to obtain an optimized point set psi 2K (vi); and then constructing an adjacent matrix according to the optimized coupling adjacent point set, and carrying out second-order neighborhood information convergence on the difference graph DI to obtain a difference graph DInew with higher robustness.

Description

Technical Field

[0001] The present invention belongs to the field of remote sensing image processing, which mainly involves change detection of synthetic aperture radar (SAR) images. The research can be applied to the fields of environmental monitoring, disaster assessment, urban planning, and flood monitoring. Background Art

[0002] SAR image change detection is an important technology that processes two SAR images of the same area at different times to obtain information about changes in the area. SAR sensors detect by actively emitting high-frequency microwave signals and are not affected by factors such as light intensity and weather. They have the advantage of all-day and all-weather observation. This makes the use of SAR sensors more popular among scholars than optical sensors in the field of change detection. Since SAR image change detection has been widely used in environmental monitoring, disaster assessment, urban planning, flood monitoring and other fields, the development of a high-precision and high-efficiency SAR image change detection method is of great significance for promoting the application of image change detection technology in actual production and life.

[0003] SAR image change detection can be divided into two categories: intensity change detection and coherence change detection. Intensity change detection usually uses the amplitude or logarithmic image of the SAR sensor, and uses some operators or methods to measure the change of the SAR image at two moments, aiming to identify the change of the average backscattered power in the scene, which is more beneficial for detecting areas with large change scales; while coherence change detection usually uses the real part, imaginary part and phase image of the SAR sensor, models the statistical characteristics of the image, and then judges the coherence information of the pixel points at two moments according to the probability density function derived from the model to measure the change, which is more advantageous for areas with subtle change scales; as a method of intensity change detection, the graph model is popular because of its high flexibility and strong robustness, and can take into account the detail information while maintaining effective detection of large-scale changes, so the present invention is designed based on the graph model method.

[0004] Currently, the most commonly used method for SAR image intensity change detection is the ratio method, which mainly performs ratio operations using the intensity values of pixel points or the intensity means of pixel points and their local neighborhoods to obtain a difference map. The ratio method has been widely used due to its high efficiency and simplicity. However, the ratio method only utilizes the information of the local neighborhood and the utilization method is relatively single. To address this issue, Wang et al. proposed a SAR image change detection method based on a heterogeneous graph with multi-attributes and multi-relationships (reference: J. Wang and A. Zhang, "SAR Image Change Detection Based on Heterogeneous Graph With Multi-attributes and Multi-relationships," in IEEE Access, vol. 10, pp. 44347-44361, 2022). This method can comprehensively utilize the intensity information and structural information of SAR images, and constructs a heterogeneous graph model of intensity-intensity, structure-structure, and intensity-structure through a graph model to describe the change information of pixel points, overcoming the problem that traditional methods only use intensity information. In fact, the above change detection methods can all effectively measure the change information of SAR images and have a good suppression effect on speckle noise. However, since the change detection results of these methods are directly generated after a single operation, misdetected and undetected pixel points in the generated results cannot be corrected. In view of the limitations of the existing technology in change detection in terms of measurement performance, the present invention is committed to researching and implementing a method for optimizing the difference map to achieve a more accurate and efficient change detection effect. Summary of the Invention

[0005] The object of the present invention is to address the problems existing in the existing change detection methods, and propose a difference-enhanced SAR image change detection method based on a coupled global graph, which strengthens the utilization of the neighborhood information of pixels and optimizes the generated initial difference map to obtain a difference map with higher separability, thereby improving the anti-noise ability and detection accuracy of the change detection method.

[0006] To achieve the above object, the implementation steps adopted by the present invention are as follows:

[0007] (1) Obtain two-time SAR images after preprocessing:

[0008] Input two SAR images X1 and X2 to be detected that have been preprocessed (including rectification and registration).

[0009] (2) Image block generation:

[0010] According to the image mean, boundary filling is performed on SAR images at two moments with N pixel points respectively, and the filled images are divided into pixel blocks to generate an overlapping image block tensor with a dimension size of 3×3×N.

[0011] (3) Construct the global neighbor set:

[0012] For a single SAR image, K most similar image blocks are assigned to each image block according to the intensity distance between blocks to construct the global neighbor set The method for obtaining this neighbor set is as follows:

[0013]

[0014] Among them, dist(v i , v j ) is the intensity distance function between blocks, and its calculation method is as follows:

[0015]

[0016] In the formula, η i represents the vector reconstructed by the image block v i by columns, and η j represents the reconstructed vector of the image block v j . (4) Global graph model construction:

[0017] The divided image blocks are set as the vertices of the global graph, and the global neighbor set of each vertex is obtained according to formula (1); and each vertex is connected to its global neighbor set to construct the edges of the graph model. Then, corresponding weights are assigned to each edge to complete the establishment of the global graph model. Define the global graph as G = {V, E, ω}, and it can be expressed in the following form:

[0018] V = {v1,..., v n}} (3)

[0019]

[0020] ω(v i , v j ) = exp(-dist 2 (v i , v j )) (5)

[0021] In the formula, dist 2 (v i , v j ) is the square value of the intensity distance between image blocks, and ω is the global edge weight.

[0022] (5) Change metric based on the global graph:

[0023] Obtain the global neighbor sets on the SAR images at two moments respectively and And build the global graph models G1 and G2 based on these, then the problem of image change measurement can be transformed into the attribute matching task of the global graph models. This attribute matching can be obtained from the vertices of the global graph and the subgraph structures connected to them, and the corresponding attribute calculation methods are as follows:

[0024] DI(v i ) = DI v (v i ) + DI sg (v i ) (6)

[0025] In the formula, DI v (v i ) represents the direct change information that occurs at the corresponding vertex, and DI sg (v i ) represents the indirect change information brought by the subgraph structure connected to the vertex. The specific calculation is as follows:

[0026]

[0027] In the formula is the block mean of the image block represented by vertex in the global graph G1. Correspondingly, is the block mean of the image block represented by vertex in the global graph G2.

[0028] (6) Initial Otsu threshold calculation:

[0029] Use the Otsu threshold method to calculate the threshold for the change measurement result DI in (5) to obtain the threshold for the binary classification of the current difference map DI.

[0030] (7) Coupled adjacent point set construction:

[0031] For the global adjacent point sets and constructed in step (5), take the union to obtain the coupled adjacent point set Its calculation method is as follows:

[0032]

[0033] (8) Coupled adjacent point set optimization and coupled global graph construction:

[0034] Since there are some interfering vertices in the coupled adjacent point set, in order to enhance the difference map, optimization needs to be carried out within the coupled neighborhood. For this purpose, a neighborhood optimization rule based on the adjacent difference mean is proposed, and the optimization can be obtained by the following formula:

[0035]

[0036] where is the number of the coupled adjacent point sets of vertex v i According to the average value of the adjacency differences, the coupled adjacent point sets of vertex v i are optimized, and the adjacent points with vertex differences less than the average value of the adjacency differences are retained, so as to construct the optimized adjacent point sets.

[0037] For the optimized point set ψ 2K (v i ), the weight W(v i , v j ) is set and calculated as follows:

[0038]

[0039] where represents the square of the intensity distance between the coupled adjacent point v j and point v i on the global graph G1. Similarly represents the square of the intensity distance between the coupled adjacent point v j and point v i on the global graph G2.

[0040] (9) Difference enhancement based on the coupled global graph

[0041] For the optimized coupled point set ψ 2K (v i ) constructed in step (8), using this point set and the newly assigned weight values, a coupled adjacency matrix P is constructed, and the difference enhancement of the change metric result DI is performed based on the coupled adjacency matrix.

[0042] DI new = unvec(P 2 .vec(DI) (15)

[0043] where vec(DI) represents flattening the difference graph matrix DI into a column vector, and unvec(.) represents the opposite operation, that is, restoring the column vector to the original specification matrix. The final difference graph DI is generated through the above formula new . In the formula, the adjacency matrix P is obtained by the following formula:

[0044]

[0045] For the strongly changing regions with significant changes in the original difference graph DI, the points in these regions are retained through neighborhood optimization, while the adjacent points with smaller changes are removed. Through the information aggregation of the adjacency matrix, the regions with a high degree of change are enhanced.

[0046] In contrast, for the weakly varying regions with insignificant changes, neighborhood optimization can preserve the points with slight changes and eliminate adjacent strongly varying points. Through the information aggregation of the adjacency matrix, the regions with low change degree are suppressed. By this method, the contrast between the strongly varying regions and the weakly varying regions is enhanced, and the separability of the difference map is improved.

[0047] (10) Output the change detection result.

[0048] Segment the final DI with the calculated initial threshold new . Obtain the non-changing class and the changing class, and obtain the change detection result map.

[0049] The present invention has the following advantages compared with the prior art:

[0050] 1) The present invention uses the average value of local image patches and the globally coupled adjacent point set to construct the edges of the graph model, optimizes the usage mode of neighborhood information, and significantly improves the robustness and accuracy in the change detection process.

[0051] 2) The present invention reduces the influence of interference points in the vertex coupled adjacent point set of the graph model by performing one update and optimization on the initial difference map, thereby achieving a smoother and more accurate change measurement result.

[0052] Experimental results show that the present invention has more accurate detection accuracy and better anti-noise ability compared with the existing change detection methods (ATDI, INLPG, NR, and WG). Description of the Drawings

[0053] Figure 1 is the implementation flowchart of the method of the present invention;

[0054] Figure 2 is the schematic diagram of constructing the coupled global neighborhood in the method of the present invention;

[0055] Figure 3 is the change detection result map of the present invention applied to two-time SAR images in the Ottawa area;

[0056] Figure 4 is the change detection result map of the present invention applied to two-time SAR images in the Yellow River Basin area. Detailed Embodiment

[0057] Refer to Figure 1 As shown, the specific implementation of the SAR image change detection method based on the coupled global graph includes the following steps:

[0058] Step 1. Obtain the two-time SAR images after preprocessing: Input the two SAR images X1 and X2 to be detected after preprocessing (including rectification and registration).

[0059] Step 2. Image block generation:

[0060] According to the image mean value, padding the boundaries of the SAR images at two moments with N pixel points respectively, and dividing the padded images into pixel blocks to generate an overlapping image block tensor with a dimension size of N*3*3.

[0061] Step 3. Constructing the global adjacent point set: For a single SAR image, assign K most similar image blocks to each image block according to the intensity distance between blocks to construct the global adjacent point set

[0062] Step 4. Global graph model construction: Set the divided image blocks as the vertices of the global graph, and obtain the global neighbor set of each vertex; then connect each vertex with its global neighbor set to construct the edges of the global graph model, and then assign corresponding weights to each edge to complete the establishment of the global graph model.

[0063] Step 5. Change metric based on the global graph: Perform change metric on the global graph models G1 and G2 established for the SAR images at two moments to generate the initial difference map DI.

[0064] Step 6. Initial Otsu threshold calculation: Use the Otsu threshold method to calculate the threshold for the change metric result DI to obtain the binary classification threshold of the current difference map DI.

[0065] Step 7 Coupled adjacent point set construction: For the global adjacent point set constructed in Step 5 and perform the union operation to obtain the coupled adjacent point set

[0066] Step 8 Coupled adjacent point set optimization and coupled global graph construction: Since there are some interference points in the coupled adjacent point set, in order to enhance the difference map, optimization needs to be carried out in the coupled neighborhood. For this purpose, a neighborhood optimization rule based on the mean of adjacent differences is proposed, and a coupled global graph model is constructed for the optimized coupled global neighborhood.

[0067] Step 9 Difference enhancement based on the coupled global graph: Perform difference enhancement on the change metric result DI through the random walk matrix P. At the same time, the second-order nearest neighbor relationship is introduced to obtain a more robust graph model.

[0068] Step 10 Output the change detection result: Segment the final DI with the calculated initial threshold new . Obtain the non-changing class and the changing class to obtain the change detection result map.

[0069] The effectiveness of the present invention is further illustrated by the following two groups of real SAR image experiments.

[0070] 1. Experimental setup:

[0071] Experimental environment: MATLAB R2022a, Intel(R) i5-9300 CPU 2.4GHz, Window10

[0072] Experimental data: Figure 3 SAR images before and after the flood in the Ottawa area obtained by the Radarsat-1 SAR satellite system, with a resolution of 10m and an image size of 290×350 pixels; Figure 4 SAR images of the Yellow River Basin area obtained by the ERS-2 SAR satellite system, with a resolution of 30m and an image size of 301×301 pixels.

[0073] Experimental parameters: Pixel block size ω = 3, global adjacent point adjacency rate K = 0.005.

[0074] Experimental comparison methods: The ATDI method, INLPG method, NR method, and WG method were respectively used.

[0075] Experimental result evaluation indicators: The number of missed detections, the number of false detections, the correct rate, and the Kappa coefficient were respectively selected to objectively evaluate the experimental results.

[0076] 2. Experimental results and analysis:

[0077] In Experiment 1, the ATDI method, INLPG method, NR method, WG method, and the method of the present invention were respectively applied to perform change detection on the two-time SAR images of the Ottawa area. The experimental results are as Figure 3 shown, where Figure 3 (a) is the SAR image of the Ottawa area at the previous moment when the flood disaster occurred, Figure 3 (b) is the SAR image of the Ottawa area at the later moment when the flood receded, Figure 3 (c) is the reference map of change detection by manual interpretation, Figure 3 (d) - (h) respectively correspond to the change detection results of the ATDI method, INLPG method, NR method, WG method, and the method of the present invention.

[0078] In Experiment 2, four change detection methods were respectively used to study the two-time SAR images of the Yellow River Basin area. The experimental results are as Figure 4 shown, where Figure 4 (a) is the SAR image before the event in the Yellow River Basin area, Figure 4 (b) is the SAR image after the event in the Yellow River Basin area, Figure 4 (c) is the reference map of change detection by manual interpretation, Figure 4(d)-(h) respectively correspond to the change detection results of the ATDI method, INLPG method, NR method, WG method and the method of the present invention.

[0079] From Figure 3 , Figure 4 and Table 1, it can be seen that the present invention is superior to other change detection methods in terms of detection accuracy and the ability to resist speckle noise. This is because the present invention uses the establishment of an optimized coupled adjacent point set method and the introduction of a difference map enhancement method, so that the method of the present invention can further improve the accuracy of the change detection method while effectively suppressing speckle noise, making the accuracy of the change detection significantly higher than other methods.

[0080] Table 1 Change detection results of SAR images in Ottawa and the Yellow River Basin

[0081]

Claims

1. A method for detecting changes in difference-enhanced SAR images based on a coupled global graph, comprising the following steps: (1) Obtain preprocessed SAR images at two moments: Two pre-processed (including correction and registration) SAR images X1 and X2 to be detected are input. (2) Image block generation: According to the image mean, the boundaries of the two-time SAR images with N pixels are filled respectively, and the filled images are divided into pixel blocks to generate overlapping image block tensors with a dimension size of 3× 3× N. (3) Constructing a global neighbor point set: For a single SAR image, K most similar image blocks are assigned to each image block according to the intensity distance between blocks to construct a global neighbor point set The adjacent point set is obtained as follows: Among them, dist(v i , v j ) is the inter-block intensity distance function, which is calculated as follows: Where η i Represents the image block v i The vector reconstructed by column, η j Represents the image block v j The reconstructed vector. (4) Construction of global graph model: The divided image blocks are set as vertices of the global graph, and the global neighbor set of each vertex is obtained according to formula (1); and each vertex is connected to its own global neighbor set to construct the edge of the graph model. Then, each edge is assigned a corresponding weight to complete the establishment of the global graph model. The global graph is defined as G = {V, E, ω}, which can be expressed as follows: V={v1,...,v n } (3) ω(in i ,in j )=exp(-dist 2 (in i ,in j )) (5) Where dist 2 (v i , v j ) is the square of the intensity distance between image patches, and ω is the global edge weight. (5) Change measurement based on the global graph: Obtain the global neighbor set on the SAR images at two times and By building global graph models G1 and G2, the image change measurement problem can be transformed into an attribute matching task of the global graph model. The attribute matching can be obtained from the vertices of the global graph and its connected subgraph structure, and the corresponding attribute calculation method is as follows: DI(v i )=DI v (v i )+DI sg (v i ) (6) Where DI v (v i ) indicates the direct change information of the corresponding vertex, DI sg (v i ) represents the indirect change information brought by the subgraph structure connected to the vertex. The specific calculation is as follows: In the formula is a vertex in the global graph G1 The block mean of the image block represented by, correspondingly, is a vertex in the global graph G2 The block mean of the image block represented. (6) Initial Otsu threshold calculation: The Otsu threshold method is used to calculate the threshold of the change measurement result DI in (5) to obtain the threshold of the binary classification of the current difference map DI. (7) Construction of coupled neighbor point set: For the global neighbor point set constructed in step (5) and Perform a union to obtain the coupled adjacent point set It is calculated as follows: (8) Optimization of coupled neighbor point sets and construction of coupled global graph: Since there are some interfering vertices in the coupled adjacent point set, in order to enhance the difference graph, it is necessary to optimize within the coupled neighborhood. For this purpose, a neighborhood optimization rule based on the mean of the adjacent difference is proposed. The optimization can be obtained by the following formula: In the formula For vertex v i The number of coupled adjacent point sets of . According to the mean of adjacent difference, for vertex v i The coupled adjacent point set is optimized, and the adjacent points whose vertex difference is less than the mean of the adjacent difference are retained to construct the optimized adjacent point set. For the optimized point set ψ 2K (v i )’s weight W(v i , v j )Set the calculation method as follows: In the formula Represents the coupled adjacent point v j With point v i The squared value of the intensity distance on the global graph G1, similarly Represents the coupled adjacent point v j With point v i The squared intensity distance on the global graph G2. (9) Difference enhancement based on coupled global graph For the optimized coupling point set ψ constructed in step (8) 2K (v i ), use the point set and the newly assigned weight value to construct the coupling adjacency matrix P, and perform difference enhancement on the change measurement result DI based on the coupling adjacency matrix. OF new =unvec(P 2 .vec(DI)) (15) Where vec(DI) means flattening the difference map matrix DI into a column vector, and unvec(.) means the opposite operation, that is, restoring the column vector to the original matrix. The final difference map DI is generated by the above formula. new . The adjacency matrix P is obtained by the following formula: For the strongly changed areas with significant changes in the original difference map DI, the points in these areas are retained through neighborhood optimization, while the adjacent points with less changes are removed. After the information of the adjacency matrix is ​​aggregated, the areas with high degree of change are enhanced. On the contrary, for weakly changed areas where the changes are not obvious, neighborhood optimization can retain points with small changes and remove adjacent points with strong changes. After the information of the adjacency matrix is ​​aggregated, the areas with low changes are suppressed. In this way, the contrast between the areas with strong changes and the areas with weak changes is enhanced, and the separability of the difference map is improved. (10) Output the change detection results. Segment the final DI using the calculated initial threshold new . Obtain the non-changing class and the changing class, and obtain the change detection result graph.

Images