Permanent scatterer point extraction method, apparatus, device and medium

Abstract

The application discloses a permanent scatterer point extraction method, device, equipment and medium, the method is through obtaining multiple SAR images of different time phases, and the SAR images of different time phases are multi-level registration; wherein, the precision registration algorithm in the multi-level registration includes a coherent coefficient fitting method, and the coherent coefficient fitting direction is the tangent vector direction of the linear region in the SAR image; the linear region where the registered main image is located is subjected to image segmentation, and a super image is obtained; finally, the selection of the PS candidate point in each super pixel is carried out until the selected PS candidate point meets the relevant index parameter threshold condition and the density of the selected PS candidate point meets the preset PS point density threshold value, and finally, the PS point is obtained, the method can effectively overcome the defects that the distribution of the PS points extracted in the linear region is uneven and the density is low in the traditional PS point extraction algorithm, and the ground surface deformation detection precision is significantly improved. Accordingly, the application also provides a permanent scatterer point extraction device, a terminal equipment and a medium.

Description

Technical Field

[0001] This invention relates to the field of surface deformation detection technology, and in particular to a method, apparatus, equipment and medium for extracting permanent scattering points. Background Technology

[0002] Persistent scatterer interferometry (PSI) is an emerging remote sensing technology with the capability of monitoring surface deformation over large areas, long periods, and with high precision. The number and density of permanent scatterer points (PS points, defined here as PS points) extracted through time-series analysis are crucial to the accuracy of surface deformation detection in Persistent Scatterer Interferometric Synthetic Aperture Radar (PS-InSAR) technology. Existing methods for PS point extraction and densification mainly fall into the following categories:

[0003] (1) Phase simulation method: The method of extracting PS points by simulating the interferometric phase of the time series to determine the phase stability of ground targets.

[0004] (2) Phase analysis method: Starting from differential interferometric phase information, after a series of filtering and iterative operations, and combined with coherence index parameters, such as coherence coefficient, temporal coherence coefficient and comprehensive coherence coefficient, the stability of ground object coherence is evaluated to achieve the purpose of PS point selection.

[0005] The above-mentioned PS point extraction algorithm can effectively extract PS points in SAR images. However, when there are linear regions in the SAR image (the so-called linear region refers to the region in the image where the region of interest is a continuous linear object such as a subway or road), the above method has shortcomings in extracting PS points in linear regions. The main problem is that the obtained PS points have a low probability of falling on the linear region, resulting in low PS point density and thus low PSInSAR deformation detection accuracy in linear regions. Summary of the Invention

[0006] This invention provides a method, apparatus, device, and medium for permanent scatterer point extraction in various aspects, which can effectively overcome the shortcomings of traditional PS point extraction algorithms in terms of uneven distribution and low density of PS points in linear regions, and significantly improve the accuracy of subsequent PSInSAR for surface deformation detection in linear regions.

[0007] The first aspect of this invention provides a method for extracting permanent scatterer points, comprising:

[0008] Acquire multiple SAR images of the study area taken at different times, and identify the main image and auxiliary image from the multiple SAR images;

[0009] The auxiliary image and the main image are subjected to multi-level SAR image registration to register the registration result onto the main image; wherein, the multi-level registration includes coarse registration and fine registration, and the registration algorithm used for fine registration includes the coherence coefficient fitting method, and the coherence coefficient fitting direction is the tangent vector direction of the linear region in the SAR image;

[0010] The linear region where the registered main image is located is segmented to obtain homogeneous regions of SAR image;

[0011] The PS candidate points in each superpixel of the homogeneous region of the SAR image are selected until the selected PS candidate points meet the preset relevant index parameter threshold conditions and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points.

[0012] A second aspect of the present invention provides a permanent scattering point device, comprising:

[0013] The acquisition module is used to acquire multiple SAR images of different time phases taken in the study area, and to determine the main image and auxiliary image among the multiple SAR images;

[0014] The registration module is used to perform multi-level registration of the auxiliary image and the main image using SAR images, so as to register the registration result onto the main image; wherein, the multi-level registration includes coarse registration and fine registration, and the registration algorithm used for fine registration includes the coherence coefficient fitting method, and the coherence coefficient fitting direction is the tangent vector direction of the linear region in the SAR image;

[0015] The segmentation module is used to segment the linear region where the registered master image is located to obtain homogeneous regions of the SAR image.

[0016] The selection module is used to select PS candidate points in each superpixel in the homogeneous region of the SAR image until the selected PS candidate points meet the preset relevant index parameter threshold conditions and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points.

[0017] A third aspect of the present invention provides a terminal device including a processor, a memory, and a computer program stored in the memory and configured to be executed by the processor, wherein the processor executes the computer program to implement the permanent scatterer point extraction method as provided in the first aspect above.

[0018] A fourth aspect of the present invention provides a storage medium including a stored computer program, wherein, when the computer program is executed, it controls the device where the storage medium is located to perform the permanent scatterer point extraction method as provided in the first aspect above.

[0019] Compared with existing technologies, the permanent scatterer point extraction method provided in this invention acquires multiple SAR images of the study area taken at different times, and determines the main image and auxiliary image from these images. Multi-level registration of the auxiliary and main images is then performed to register the registration result onto the main image. This multi-level registration includes coarse and fine registration. The fine registration uses a coherence coefficient fitting algorithm, with the fitting direction being the tangent vector direction of the linear region in the SAR image. Next, the linear region containing the registered main image is segmented to obtain homogeneous SAR image regions. Finally, candidate PS points are selected from each superpixel in the homogeneous SAR image regions until the selected candidate PS points meet preset threshold conditions for relevant index parameters and the density of the selected candidate PS points meets a preset PS point density threshold, thus obtaining the final PS points. This invention effectively overcomes the shortcomings of traditional PS point extraction algorithms in terms of uneven PS point distribution and low density in linear regions, significantly improving the accuracy of subsequent PSInSAR surface deformation detection in linear regions. Accordingly, embodiments of the present invention also provide a permanent scatterer point extraction device, a terminal device, and a medium. Attached Figure Description

[0020] Figure 1 This is a schematic flowchart of the permanent scatterer point extraction method provided in an embodiment of the present invention;

[0021] Figure 2 This is a schematic diagram of a multi-level registration process provided in an embodiment of the present invention;

[0022] Figure 3 This is a schematic diagram illustrating the sub-pixel level registration offset obtained by a coherence coefficient fitting method according to an embodiment of the present invention.

[0023] Figure 4 This is a schematic diagram of the fitting direction during the correlation coefficient fitting process provided in the embodiments of the present invention;

[0024] Figure 5 This is a comparison diagram of image search performed by the super-image segmentation algorithm provided in this embodiment of the invention and the standard K-means algorithm;

[0025] Figure 6 This is a schematic diagram of the PS point filtering and extraction process provided in an embodiment of the present invention;

[0026] Figure 7This is a flowchart of a method for removing flat-ground phases from interference phases based on ephemeris parameters and DEM data, provided by an embodiment of the present invention. Detailed Implementation

[0027] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0028] The following explanations are provided for the relevant terms appearing in this application:

[0029] Image registration: Given a point in the master image, find the corresponding point in the auxiliary image.

[0030] Same point: The image point formed by the same point on the ground in different images.

[0031] Direction and distance: Direction refers to the row direction, and distance refers to the column direction.

[0032] Coherence refers to the degree of coherence between two images. Coherence is an indicator of the processing effect of InSAR images. The level of coherence determines the quality of the interference fringes, and thus the result of the interferometric processing. SAR image interferometric processing is based on the coherence between images. Therefore, in this application, the coherence index is used to select interferometric image pairs and to determine the feasibility of the selected candidate PS points.

[0033] Sub-pixel level: This refers to the subdivision between two adjacent pixels. The input value is typically 1 / 2, 1 / 3, or 1 / 4. This means that each pixel will be divided into smaller units, and interpolation algorithms will be applied to these smaller units. For example, if 1 / 4 is chosen, it is equivalent to treating each pixel as four pixels horizontally and vertically. Therefore, if a 5x5 pixel image is processed with 1 / 4 sub-pixel precision, it is equivalent to creating a 20x20 discrete pixel array, which is then interpolated.

[0034] See Figure 1 , Figure 1 This is a schematic flowchart of the permanent scatterer point extraction method provided in this embodiment of the invention. The permanent scatterer point extraction method provided in this embodiment of the invention includes steps S11 to S14:

[0035] S11. Acquire multiple SAR images of the study area taken at different times, and determine the main image and auxiliary image among the multiple SAR images.

[0036] S12. Perform multi-level SAR image registration on the auxiliary image and the main image to register the registration result onto the main image; wherein, the multi-level registration includes coarse registration and fine registration, and the registration algorithm used for fine registration includes the coherence coefficient fitting method, and the coherence coefficient fitting direction is the tangent vector direction of the linear region in the SAR image.

[0037] S13. Perform image segmentation on the linear region where the registered main image is located to obtain homogeneous regions of SAR image;

[0038] S14. Select PS candidate points in each superpixel in the homogeneous region of the SAR image until the selected PS candidate points meet the preset relevant index parameter threshold conditions and the density of the selected PS candidate points meets the preset PS point density threshold, and obtain the final PS points.

[0039] Based on the technical solution provided in the above embodiments, it acquires multiple SAR images of different time phases taken in the study area, and determines the main image and auxiliary image from the multiple SAR images; then, it performs multi-level registration of the auxiliary image and the main image to register the registration result onto the main image; wherein, the multi-level registration includes coarse registration and fine registration, and the registration algorithm used for fine registration includes the coherence coefficient fitting method, and the coherence coefficient fitting direction is the tangent vector direction of the linear region in the SAR image; next, the linear region where the registered main image is located is segmented to obtain homogeneous regions of the SAR image; finally, PS candidate points are selected in each superpixel of the homogeneous region of the SAR image until the selected PS candidate points meet the preset relevant index parameter threshold conditions and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points. The embodiments of the present invention can effectively overcome the defects of uneven distribution and low density of PS points extracted in linear regions by traditional PS point extraction algorithms, and significantly improve the accuracy of subsequent PSI nSAR for surface deformation detection in linear regions.

[0040] In one optional implementation, the fine registration specifically includes:

[0041] Based on coarse registration of SAR images, oversampling is performed on the SAR images;

[0042] The oversampled SAR images are fitted with coherence coefficients to achieve fine registration.

[0043] In one optional implementation, the coherence coefficient fitting method specifically includes:

[0044] Based on the coarse registration of SAR images, the correlation coefficients of each pixel in its neighborhood are fitted with a surface, with the registered pixel location as the center. The fitting direction of the surface is the tangent vector direction of the linear region.

[0045] Find the maximum coherence coefficient along the current tangent vector direction of the linear region, and obtain the position corresponding to the maximum coherence coefficient as the sub-pixel level matching point to complete the fine registration.

[0046] Specifically, the fitting direction refers to the normal of the fitted surface. The fitting direction determines the quality of the fitted surface. In this application, the characteristics of the linear region are considered, and the fitting direction of the surface is set to the tangent vector direction (and tangent direction) of the linear region to obtain the optimal surface. Then, the sub-pixel coordinates can be obtained from the peak of the fitted surface to obtain the sub-pixel level offset, so as to complete the fine registration.

[0047] In this embodiment, multi-level registration includes a configuration strategy from coarse registration to sub-pixel-level fine registration. The multi-level registration process can be found in [link to documentation]. Figure 2 The coarse registration stage primarily aims to determine the offset between the main and auxiliary images. This can be achieved by having the user select a corresponding point on each of the two images, obtaining the coordinates of that point on each image, and then subtracting the coordinates to obtain the offset value between the two images. After coarse registration, pixel-level registration can be performed, followed by sub-pixel-level fine registration.

[0048] The following provides a detailed explanation of each registration stage.

[0049] The coarse registration stage is mainly to obtain the offset between the main image and the auxiliary image, which can be calculated using the following formula in one implementation:

[0050] p sla (m,n)=p mas (m,n)+offset(m,n) (1)

[0051] In the formula, p mas (m,n) represents a point on the main image, p sla (m,n) represents the corresponding points of the auxiliary image, and offset(m,n) represents the offset between the main image and the auxiliary image.

[0052] The offset accuracy requirement is within 30 pixels. You can select a corresponding point on the two images to be registered, and then obtain the coordinates of the point on the main and auxiliary images respectively. Subtracting the coordinates will give you the offset value between the two images.

[0053] In another implementation, the coarse registration stage can also automatically calculate the offset between two images using the imaging trajectory information. The implementation process of this method is as follows:

[0054] Step A: Using the RD imaging geometric model (Equation (2)), first calculate the ground coordinates corresponding to the center point of the main image, then find the coordinates of the corresponding points of the ground points on the auxiliary image, and solve the offset between the main image and the auxiliary image from Equation (1). The relevant parameters in the equation are obtained through the image header file.

[0055]

[0056] Step B, obtaining the registration offset of the main image center point, is as follows:

[0057] B-1. The geographic coordinates (X, Y) of the corresponding ground point are directly calculated based on the image coordinates (row and column numbers known) of the center point of the main image and the SAR image header information. s ,Y s Z s ).

[0058] B-2, Iterative Solution (X) s ,Y s Z s ), corresponding to the azimuth time T of the auxiliary image imaging point a The solution for the azimuth time is mainly obtained by using an iterative approximation method. The basis of the iteration is the Doppler equation E2 in equation (2). The ratio of the Doppler equation and its differential equation is used as the azimuth time T. a The correction value is obtained by iterative accumulation to obtain the accurate azimuth time. The specific steps include: 1) T a 1) Initialize to the mean value of the SAR imaging time; 2) Calculate the relative position vector S between the satellite and the ground point. sp 3) Calculate the satellite velocity vector V S ;4) Calculate the satellite acceleration vector A a 5) Calculate the Doppler equation F t =S sp V S and for F t Find the partial differential f with respect to t. t =V S V S +A a S sp , get T a The correction is ΔT = -F t / f t 6) If ΔT is greater than a certain threshold (usually set to 10) e-15 If the maximum number of iterations has not been reached (which can be set to 10), the iteration continues, and ΔT is corrected to T. a 6) Once the iteration termination condition is met, the accurate T can be obtained. a .

[0059] B-2, (X) s ,Y s Z s ) Calculation of the distance-time axis T corresponding to the imaging points of the auxiliary image r The distance equation is mainly based on the precise satellite position vector and geodetic position vector, as well as the distance equation E1 in equation (2). The actual satellite-to-ground slant distance R is obtained using the distance equation, and then the distance-to-time vector T is obtained from the relationship between the electromagnetic wave propagation speed C and time. r That is, T r =R / C.

[0060] B-3, Obtain the azimuth time T a and distance to time T r Then, based on the imaging parameters of the SAR image header file, the imaging coordinates of the ground point corresponding to the center point of the main image in the auxiliary image can be obtained, and finally the offset between the coordinates of the main image and the auxiliary image can be obtained.

[0061] The pixel-level registration stage mainly involves obtaining the pixel-level registration offset, which can be achieved through oversampling. In this embodiment, pixel-level registration based on moving window search is used. This method primarily involves selecting a search window on the auxiliary image based on the coarse matching results, then calculating the similarity between the two windows according to registration evaluation metrics. Gross errors are eliminated through moving the search window and a consistency check to obtain the pixel-level registration offset. Specifically, this includes:

[0062] Step a: Selection of registration points

[0063] Based on the required number of points to be matched and the size of the overlapping area obtained from the coarse matching, windows are evenly divided on the main image, with the center point of the window as the point to be registered.

[0064] Step b: Calculation of registration evaluation indicators

[0065] For each point to be registered, registration evaluation parameters are calculated within the search window based on the search range (typically 30*30 pixels) and the size of the matching window. To speed up the calculation, the correlation coefficient is usually chosen as the registration evaluation parameter; the point with the highest correlation coefficient is the matching point. The correlation coefficient is calculated as follows:

[0066]

[0067] In the formula, corr is the correlation coefficient; g i,j ,g' i+r.j+c , respectively, represent the amplitude intensity at corresponding points in the main image and the auxiliary image; m, n are the sizes of the matching window, and r, c are the offsets for coarse matching.

[0068] Step c: Matching point consistency check

[0069] Using windows for correlation estimation often results in biased values, leading to discrepancies in registration. Furthermore, the presence of speckle noise in SAR images means that the point with the highest correlation may not be the optimal match. Therefore, a consistency check is necessary to ensure that the offsets of most points remain consistent, thus eliminating significant errors.

[0070] A typical consistency check can be performed as follows: A coherence coefficient threshold is set. Matching points with a coherence coefficient greater than this threshold are clustered using the minimum distance method to obtain cluster center values. Simultaneously, the "clustering degree" is calculated by counting all points smaller than a certain radius, with the cluster center as the center of a sphere. If the "clustering degree" is less than a given threshold (e.g., half the total number of valid matching points), registration has failed; otherwise, points larger than the radius are discarded as gross errors.

[0071] The sub-pixel-level fine registration stage aims to obtain offset values ​​with lower precision. In this embodiment of the invention, sub-pixel-level registration employs a correlation coefficient fitting method. Considering the characteristics of linear regions, the direction of the tangent vector in the linear region is introduced as the direction for coherence coefficient fitting, which improves the coherence of the linear region. The coherence coefficient fitting method does not oversample the pixel itself, but rather correlates the coherence coefficient values ​​of surrounding pixels with an adjustment function to obtain the maximum value of the function as the sought corresponding point. Specifically, it includes the following steps:

[0072] Step (1): Based on the oversampled SAR image pixels, take the registered pixel location as the center and obtain the correlation coefficient at the location of the center of each pixel in its neighborhood.

[0073] Step (2) involves surface fitting of the correlation coefficients within the neighborhood range, such as... Figure 3 As shown, there are 5 coherence coefficients within the neighborhood. These 5 coherence coefficients are fitted using a quadratic parabolic equation. Let the general equation of the parabola be:

[0074] f(S) = A + B·S + C·S 2 (4)

[0075] In the formula, S is the independent variable, and A, B, and C are parameters that can be obtained using the indirect adjustment method. The address of the vertex k of the parabola should then be:

[0076] k=iB / (2·C) (5)

[0077] When the coherence coefficients of three adjacent pixels are used for parabolic fitting, the following system of equations can be obtained:

[0078]

[0079] Where ρi-1 , ρ i , ρ i+1 , where is the coherence coefficient. With the coordinate system shifted to point i, from the above equation, we get:

[0080] A = ρ i

[0081] B=(ρ i+1 -ρ i-1 ) / 2

[0082] C=(ρ i+1 -2ρ i +ρ i-1 ) / 2 (7)

[0083] but

[0084]

[0085] Step (3) The fine registration method based on coherence coefficient fitting is essentially a two-dimensional extension of the above theory. Two-dimensional curve fitting is a mathematical solution process, which is often achieved through interpolation in specific implementation. Therefore, the basic idea of ​​the fine registration method based on coherence coefficient fitting is to interpolate the coherence coefficient values ​​in the two-dimensional neighborhood range to achieve two-dimensional curve fitting of the coherence coefficient, thereby obtaining the best matching position.

[0086] The specific formula for calculating the interpolation coherence coefficient within the neighborhood range is as follows:

[0087]

[0088] In the formula, z ij Here, z(x,y) represents the coherence coefficient value, and z(x,y) represents the interpolated correlation coefficient value. The coherence coefficients are interpolated at a certain pixel spacing (e.g., 0.001 pixels). After the coherence coefficient interpolation is complete, the maximum coherence coefficient (i.e., the peak value of the fitted surface) is found along the current tangent vector direction (i.e., the tangent direction) of the linear region. The position corresponding to the maximum coherence coefficient is the more accurate registration position. This method can effectively improve the registration accuracy in the vector direction, reduce the influence of spatial decorrelation, and thus potentially obtain more PS candidate points.

[0089] For a linear region, the direction of its tangent vector Θ is defined by the centerline of the linear region, such as... Figure 4 ,but:

[0090]

[0091] In the formula, (x i y i ), (x i+1 y i+1 ) are the coordinates of the current point P and the next point Q on the center line of the linear region, respectively.

[0092] In one optional implementation, the step of segmenting the registered master image to obtain homogeneous SAR image regions specifically includes:

[0093] The superpixel segmentation algorithm (also known as the SLIC algorithm) is used to segment the registered master image to obtain homogeneous regions of the SAR image.

[0094] In this embodiment, the specific steps for performing super-image segmentation on the registered main image using a superpixel segmentation algorithm include:

[0095] S131. Initialize seed points (cluster centers) and distribute them evenly within the image according to the set number of superpixels. Assuming the image has a total of N pixels, pre-segmented into K superpixels of the same size, then the size of each superpixel is N / K, and the distance (step size) between adjacent seed points is approximately S = sqrt(N / K).

[0096] S132. Reselect a seed point within the n*n neighborhood of the seed point (generally n=3). The specific method is: calculate the gradient values ​​of all pixels within this neighborhood, and move the seed point to the location with the smallest gradient within that neighborhood. This is done to avoid the seed point falling on the contour boundary with a large gradient, which could affect subsequent clustering results.

[0097] S133. Assign a class label (i.e., which cluster center each pixel belongs to) to each pixel within its neighborhood. Unlike the standard K-means algorithm, which searches the entire image, SLIC limits its search range to 2S*2S, which can accelerate the algorithm's convergence. See below. Figure 5 At this point, the desired superpixel size is S*S, but the search range is 2S*2S.

[0098] S134. Distance metric. This includes color distance and spatial distance. For each searched pixel, calculate its distance to the seed point. The distance calculation method is as follows:

[0099]

[0100] Where, d c Representing color distance, for SAR images, the color channel component (a i b i ) and (a j b j ) represent the amplitude intensity and phase value of the current pixel and the seed point, respectively. d s N represents spatial distance. s It is the maximum spatial distance within a class, defined as Applies to each cluster. Maximum color distance N cThe distance metric varies depending on both the image and the clustering method. A fixed constant m (ranging from [1,40], typically 10) is used instead. The final distance metric D' is as follows:

[0101]

[0102] S135, Iterative optimization. The above steps are iterated until the error converges (which can be understood as the cluster center of each pixel no longer changing).

[0103] S136. Enhance Connectivity. The above iterative optimization may result in the following defects: multiple connectivity issues, excessively small superpixel sizes, and a single superpixel being fragmented into multiple discontinuous superpixels. These issues can be resolved by enhancing connectivity. The main idea is: create a new label table with all elements set to -1. Following a "Z" shaped pattern (from left to right, from top to bottom), discontinuous superpixels and excessively small superpixels are reassigned to neighboring superpixels. Pixels that have been traversed are assigned their corresponding labels until all pixels have been traversed. After traversal, the SAR image segmentation is complete.

[0104] In one optional implementation, the relevant index parameters include coherence coefficient and amplitude deviation. The step of selecting PS candidate points in each superpixel of the homogeneous region of the SAR image until the selected PS candidate points meet the preset threshold conditions of the relevant index parameters and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points, specifically includes:

[0105] Set thresholds for coherence coefficient and amplitude deviation;

[0106] Pixels with an amplitude deviation greater than the amplitude deviation threshold and a coherence coefficient greater than the coherence threshold within the superpixel are selected as PS candidate points;

[0107] When the density of PS candidate points selected in the superpixel is greater than or equal to the preset PS point density threshold, all PS candidate points currently selected for each superpixel are taken as the final PS points.

[0108] When the density of PS candidate points selected in the superpixel is less than the PS point density threshold, the amplitude deviation threshold and the coherence coefficient threshold are modified, and PS candidate points are selected according to the modified amplitude deviation threshold and the modified coherence coefficient threshold until the density of selected PS candidate points meets the PS point density threshold.

[0109] In this embodiment, the implementation process of the above steps can be found in [reference needed]. Figure 6 ,include:

[0110] S41. Based on the registration results of the main image and the auxiliary image, generate a differential interferogram, including:

[0111] Multi-temporal differential interferometry (DInSAR) processing is performed on the main image and the auxiliary image, including:

[0112] DEM resampling: Acquire external DEM data (available data sources include GOTOPO30 (1km resolution) and SRTM (3 radian resolution)) and resample the DEM to the size of the interferometric image pixels;

[0113] Removal of the flat phase effect: Based on Doppler, slant range, and ellipsoid equations, DEM data is transformed into an interferometric coordinate system, and the flat phase is simulated and removed from the interferometric phase. The specific flowchart is as follows. Figure 7 ;

[0114] S42. Set the amplitude deviation threshold; based on the generated differential interferogram, calculate the amplitude deviation of the pixels in each superpixel of the interferogram, and select the pixels with amplitude deviation greater than the amplitude deviation threshold as the initial PS selection points.

[0115] S43. The interference phase model is established based on the initial PS selection points as follows:

[0116] Δφ i =φ topo_ε_i +φ def_i +φ atm_i +φ noise_i (14)

[0117] The interference phase model includes the three main parts of the interference phase, namely the deformation phase φ. def_i DEM error phase φ topo_ε_i and atmospheric delay phase φ atm (x,y), the remaining part is the noise phase φ noise_i By setting a certain threshold, the selected PS candidate points are connected into a spatially interconnected sparse network. Assuming that the adjacent points are x and y, the phase gradient model between the two points can be expressed as:

[0118] Δφ diff (x,y)=Δφ topo_ε (x,y)+Δφ linear (x,y)+Δφ non-linear (x,y)+Δφ atm (x,y)+Δφ noise (x,y) (15)

[0119] In the formula, Δφ topo_ε_i ,Δφ atm (x,y),Δφ noise_iLet Δφ represent the DEM error phase gradient, atmospheric delay phase gradient, and noise phase gradient, respectively. linear (x,y),Δφ non-linear (x, y) represent the deformation phase gradient Δφ, respectively. def_i The linear deformation phase component and the nonlinear deformation phase component.

[0120] S44. Solving for the deformation phase, DEM error phase, and atmospheric delay phase of the initial PS points, i.e., solving equation (15). Solving equation (15) is a very complex linear process, and it is difficult to directly obtain the deformation phase between adjacent points. The ultimate goal of using interferometry to obtain deformation is to obtain the deformation phase, which can be called the effective phase, then:

[0121] Δφ atm (x,y)+Δφ noise (x,y)=Δφ diff (x,y)-Δφ topo_ε (x,y)-Δφ linear (x,y)-Δφ non-linear (x,y)(16)

[0122] For the initial PS points adjacent to each other, Δφ non-linear (x,y) is very small and has little impact on the overall result; it can be considered noise. Therefore, the above equation can be transformed into:

[0123] Δω(x,y)=Δφ diff (x,y)-(Δφ topo_ε (x,y)+Δφ linear (x,y)) (17)

[0124] When the residual phase between adjacent points is sufficiently small, it is generally considered that most phases of adjacent points in a time-series differential interferogram set (i.e., differential interferograms of multiple different time phases) can satisfy (-π, π], and phase unwrapping can then be performed. The optimal estimate of the DEM error phase and linear deformation phase is obtained by iteratively searching across all time-series differential interferograms, achieving the highest overall coherence across all interferograms. Using a weighted least squares method, phase unwrapping can be performed by selecting an initial PS point as a reference point and traversing the sparse network. From this, the coherence coefficient γ can be defined. s As a weighting indicator:

[0125]

[0126] In the formula, γ s The value of γ ranges from [0, 1] and can be used to evaluate phase dispersion. For a PS candidate point, when its phase dispersion across all interferograms is low, γ...s Approaching 1. Therefore, by maximizing γ s To obtain the linear deformation rate and DEM error.

[0127] S45. Filtering and Separating the Atmospheric Delay Phase of the Initial PS Point: The residual phase after subtracting the deformation phase and elevation error phase from the differential interferometric phase includes the atmospheric delay phase and the nonlinear deformation and noise phase. By filtering the temporal and spatial signal attributes of the atmospheric effect phase and the nonlinear deformation, the two can be separated. The atmospheric effect phase is random in the time dimension and belongs to the high-frequency signal, while it is continuous in the spatial dimension and belongs to the low-frequency signal; the nonlinear deformation phase is a low-frequency signal in the time dimension. High-frequency and low-frequency filtering in both time and space dimensions separates the two. For a point target, the mean value of the residual phase in N images... Atmospheric phase estimation as the main image. The low-pass portion in the time dimension is used as an estimate of the nonlinear deformation.

[0128] S46. Iteratively selecting initial PS points to obtain PS candidate points: The atmospheric delay phase of the PS candidate points is restored to the resolution scale of the original radar image using an interpolation algorithm, becoming the atmospheric phase screen. After subtracting the atmospheric phase screen from the differential interferometric phase, all pixels are analyzed using a time-series differential interferometry model. The above steps are repeated until the pixel's γ... s When a value meets the set coherence threshold, it is determined as a candidate point for PS, and the PS point set is obtained.

[0129] S47. Based on the number and density of the PS point set obtained in step S46, continuously adjust the amplitude deviation threshold in step S42 and the coherence threshold γ in step S44. s This causes the PS point set density to exceed a set threshold in each superpixel.

[0130] Accordingly, embodiments of the present invention also provide a permanent scatterer point extraction device, comprising:

[0131] The acquisition module is used to acquire multiple SAR images of different time phases taken in the study area, and to determine the main image and auxiliary image among the multiple SAR images;

[0132] The registration module is used to perform multi-level registration of the auxiliary image and the main image using SAR images, so as to register the registration result onto the main image; wherein, the multi-level registration includes coarse registration and fine registration, and the registration algorithm used for fine registration includes the coherence coefficient fitting method, and the coherence coefficient fitting direction is the tangent vector direction of the linear region in the SAR image;

[0133] The segmentation module is used to segment the linear region where the registered master image is located to obtain homogeneous regions of the SAR image.

[0134] The selection module is used to select PS candidate points in each superpixel in the homogeneous region of the SAR image until the selected PS candidate points meet the preset relevant index parameter threshold conditions and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points.

[0135] In one optional implementation, the coherence coefficient fitting method specifically includes:

[0136] Based on the coarse registration of SAR images, the correlation coefficients of each pixel in its neighborhood are fitted with a surface, with the registered pixel location as the center. The fitting direction of the surface is the tangent vector direction of the linear region.

[0137] Find the maximum coherence coefficient along the current tangent vector direction of the linear region, and obtain the position corresponding to the maximum coherence coefficient as the sub-pixel level matching point to complete the fine registration.

[0138] In one optional implementation, the relevant index parameters include coherence coefficient and amplitude deviation. The step of selecting PS candidate points in each superpixel of the homogeneous region of the SAR image until the selected PS candidate points meet the preset threshold conditions of the relevant index parameters and the density of the selected PS candidate points meets the preset PS point density threshold, thus obtaining the final PS points, specifically includes:

[0139] Set thresholds for coherence coefficient and amplitude deviation;

[0140] Pixels with amplitude deviation greater than the amplitude deviation threshold and coherence coefficient greater than the coherence coefficient threshold within the superpixel are selected as PS candidate points. When the density of PS candidate points selected in the superpixel is greater than or equal to the preset PS point density threshold, all PS candidate points currently selected for each superpixel are taken as the final PS point.

[0141] When the density of PS candidate points selected in the superpixel is less than the PS point density threshold, the amplitude deviation threshold and the coherence coefficient threshold are modified, and PS candidate points are selected according to the modified amplitude deviation threshold and the modified coherence coefficient threshold until the density of selected PS candidate points meets the PS point density threshold.

[0142] This invention also provides a terminal device, which includes a processor, a memory, and a computer program stored in the memory and configured to be executed by the processor. When the processor executes the computer program, it implements steps S11 to S14 of the permanent scatterer point extraction method described above. Alternatively, when the processor executes the computer program, it implements the functions of each module in the above-described device embodiments, such as an acquisition module, a registration module, a segmentation module, and a selection module.

[0143] For example, the computer program may be divided into one or more modules, which are stored in the memory and executed by the processor to complete the present invention. The one or more modules may be a series of computer program instruction segments capable of performing specific functions, which describe the execution process of the computer program in the permanent scatterer point extraction device / terminal equipment.

[0144] The permanent scattering point extraction device / terminal equipment can be a computing device such as a desktop computer, laptop, handheld computer, or cloud server. The permanent scattering point extraction device / terminal equipment may include, but is not limited to, a processor and memory. Those skilled in the art will understand that the schematic diagram is merely an example of the permanent scattering point extraction device / terminal equipment and does not constitute a limitation on the device. It may include more or fewer components than illustrated, or combine certain components, or use different components. For example, the permanent scattering point extraction device / terminal equipment may also include input / output devices, network access devices, buses, etc.

[0145] The processor can be a Central Processing Unit (CPU), or other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. The general-purpose processor can be a microprocessor or any conventional processor. The processor is the control center of the permanent scattering point extraction device / terminal equipment, connecting all parts of the device via various interfaces and lines.

[0146] The memory can be used to store the computer program and / or modules. The processor implements various functions of the permanent scatterer point extraction device / terminal device by running or executing the computer program and / or modules stored in the memory and calling the data stored in the memory. The memory may mainly include a program storage area and a data storage area. The program storage area may store the operating system, at least one application program required for a function (such as sound playback function, image playback function, etc.), etc.; the data storage area may store data created according to the use of the mobile phone (such as audio data, phonebook, etc.). In addition, the memory may include high-speed random access memory, and may also include non-volatile memory, such as hard disk, memory, plug-in hard disk, smart media card (SMC), secure digital card (SD) card, flash card, at least one disk storage device, flash memory device, or other volatile solid-state storage device.

[0147] Wherein, if the modules / units integrated in the permanent scatterer point extraction device / terminal equipment are implemented as software functional units and sold or used as independent products, they can be stored in a computer-readable storage medium. Based on this understanding, all or part of the processes in the methods of the above embodiments of the present invention can also be implemented by a computer program instructing related hardware. The computer program can be stored in a computer-readable storage medium, and when executed by a processor, it can implement the steps of the various method embodiments described above. The computer program includes computer program code, which can be in the form of source code, object code, executable files, or certain intermediate forms, etc.

[0148] This invention also provides a storage medium comprising a stored computer program, wherein the computer program, when running, controls the device containing the storage medium to execute steps S11 to S14 of the above-described permanent scatterer point extraction method.

[0149] The storage medium may include: any entity or device capable of carrying the computer program code, recording media, USB flash drive, portable hard drive, magnetic disk, optical disk, computer memory, read-only memory (ROM), random access memory (RAM), electrical carrier signal, telecommunication signal, and software distribution medium, etc.

[0150] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications are also considered to be within the scope of protection of the present invention.

Claims

1. A permanent scatterer point extraction method characterized by, The method comprises the following steps: acquire a plurality of SAR images of different time phases taken by shooting a study area, and determine a main image and an auxiliary image in the plurality of SAR images; perform multi-level registration of the SAR images on the auxiliary image and the main image to register the registration result to the main image; wherein the multi-level registration comprises coarse registration and fine registration, the fine registration adopts a registration algorithm comprising a coherent coefficient fitting method, and the coherent coefficient fitting direction is a tangent vector direction of a linear region in the SAR image; perform image segmentation on a linear region where the registered main image is located to obtain a SAR image homogeneous region; perform selection of a PS candidate point in each superpixel in the SAR image homogeneous region until the selected PS candidate point meets a preset related index parameter threshold condition and the density of the selected PS candidate point meets a preset PS point density threshold to obtain a final PS point; wherein the related index parameters comprise a coherent coefficient and an amplitude deviation, and the performing selection of a PS candidate point in each superpixel in the SAR image homogeneous region until the selected PS candidate point meets a preset related index parameter threshold condition and the density of the selected PS candidate point meets a preset PS point density threshold to obtain a final PS point specifically comprises: setting a threshold of the coherent coefficient and the amplitude deviation; selecting a pixel point with an amplitude deviation greater than an amplitude deviation threshold and a coherent coefficient greater than a coherent coefficient threshold in the superpixel as a PS candidate point; when the density of the selected PS candidate point in the superpixel is greater than or equal to a preset PS point density threshold, taking all the PS candidate points currently selected in each superpixel as the final PS point; when the density of the selected PS candidate point in the superpixel is less than the PS point density threshold, modifying the amplitude deviation threshold and the coherent coefficient threshold, and performing selection of a PS candidate point according to the modified amplitude deviation threshold and the modified coherent coefficient threshold until the density of the selected PS candidate point meets the PS point density threshold.

2. The permanent scatterer point extraction method of claim 1, wherein, The coherent coefficient fitting method specifically comprises: on the basis of coarse registration of the SAR image, fitting a surface to the correlation coefficients at the positions of the centers of the pixel points in the neighborhood range of the pixel position obtained by registration, wherein the fitting direction of the surface is a tangent vector direction of the linear region; finding a maximum value of the coherent coefficient along the current tangent vector direction of the linear region, obtaining the position corresponding to the maximum value of the coherent coefficient as a matching point at the sub-pixel level to complete fine registration.

3. The permanent scatterer point extraction method of claim 1, wherein, The performing image segmentation on a linear region where the registered main image is located to obtain a SAR image homogeneous region specifically comprises: performing segmentation on the linear region where the registered main image is located by using a superpixel segmentation algorithm to obtain a SAR image homogeneous region.

4. The permanent scatterer point extraction method of claim 1, wherein, The fine registration specifically comprises: on the basis of coarse registration of the SAR image, performing oversampling on the SAR image; performing coherent coefficient fitting on the SAR image after oversampling to complete fine registration.

5. A permanent scatterer point extraction apparatus characterized by, The method comprises the following steps: an acquisition module configured to acquire a plurality of SAR images of different time phases taken by shooting a study area, and determine a main image and an auxiliary image in the plurality of SAR images; The registration module is configured to perform multi-stage registration of the auxiliary image and the main image to obtain a registration result, and to register the registration result to the main image; wherein the multi-stage registration comprises coarse registration and fine registration, and the fine registration employs a registration algorithm comprising a coherence coefficient fitting method, and the fitting direction of the coherence coefficient fitting method is a tangent vector direction of a linear region in the SAR image; The segmentation module is configured to perform image segmentation on the linear region in which the main image is located to obtain a SAR image homogeneous region. The selection module is configured to select a PS candidate point in each superpixel in the SAR image homogeneous region until the selected PS candidate point meets a preset threshold condition of a related index parameter and the density of the selected PS candidate points meets a preset PS point density threshold, and to obtain a final PS point. The related index parameter comprises a coherence coefficient and an amplitude deviation, and the selection of the PS candidate point in each superpixel in the SAR image homogeneous region until the selected PS candidate point meets the threshold condition of the related index parameter and the density of the selected PS candidate points meets the PS point density threshold to obtain the final PS point specifically comprises: setting a threshold of the coherence coefficient and the amplitude deviation; selecting a pixel point in the superpixel as the PS candidate point when the amplitude deviation of the pixel point is greater than an amplitude deviation threshold and the coherence coefficient of the pixel point is greater than a coherence coefficient threshold; when the density of the selected PS candidate points in the superpixel is greater than or equal to a preset PS point density threshold, regarding all the PS candidate points currently selected in each superpixel as the final PS point; 6. The permanent scatterer point extraction apparatus of claim 5, wherein, when the density of the selected PS candidate points in the superpixel is less than the PS point density threshold, modifying the amplitude deviation threshold and the coherence coefficient threshold, and selecting the PS candidate point according to the modified amplitude deviation threshold and the modified coherence coefficient threshold until the density of the selected PS candidate points meets the PS point density threshold. The coherence coefficient fitting method specifically comprises: on the basis of coarse registration of the SAR image, fitting a surface to the correlation coefficients at positions of centers of pixel points in a neighborhood range of a pixel location obtained by registration, wherein the fitting direction of the surface is a tangent vector direction of the linear region; 7. A terminal device, characterized by, finding a maximum value of the coherence coefficient along the tangent vector direction of the linear region to obtain a position corresponding to the maximum value of the coherence coefficient as a matching point at a sub-pixel level, and completing fine registration.

8. A storage medium, characterized by The storage medium comprises a stored computer program, wherein the computer program controls a device in which the storage medium is located to perform the permanent scatterer point extraction method according to any one of claims 1 to 4 when the computer program is running. The storage medium comprises a stored computer program, wherein the computer program controls a device in which the storage medium is located to perform the permanent scatterer point extraction method according to any one of claims 1 to 4 when the computer program is running.

Images