Distributed scatterer phase rapid estimation method based on improved PM

Abstract

The invention discloses a distributed scatterer phase rapid estimation method based on improved PM, and belongs to the technical field of observation data processing, and the method comprises the steps: recognizing distributed scatterers DS according to the amplitude information of a registered SAR image, and marking the homogeneous pixel SHP of each DS; calculating an interference phase mean value according to the SHP of the DS; constructing a sample covariance matrix C; iteratively estimating a main characteristic value of C and a corresponding characteristic vector by adopting matrix vector multiplication of PM; and constructing an adaptive termination strategy for the PM, and adaptively terminating iteration according to the DS quality and obtaining the optimal estimated value of the DS phase. According to the improved PM-based distributed scatterer phase rapid estimation method, the number of iterations is reduced, the convergence speed is improved, the algorithm is ensured to rapidly converge to an optimal solution for high-quality DS and timely terminate iteration for low-quality DS by setting a self-adaptive termination strategy, the calculation efficiency is improved, and the method is suitable for large-scale popularization and application. The problem that calculation is too time-consuming is effectively solved, calculation is simple, and the method is particularly suitable for parallel calculation.

Description

Technical Field

[0001] The present invention relates to the technical field of observation data processing, and in particular to a distributed scatterer phase fast estimation method based on improved PM. Background Art

[0002] Time-series InSAR technology has become one of the commonly used geodetic technologies due to its advantages such as high precision, wide coverage, and all-weather. The objects of time-series InSAR analysis are a series of special coherent targets, mainly including permanent scatterers (PS) and distributed scatterers (DS). Typical PS are mainly artificial buildings, so they are suitable for urban areas with dense artificial buildings. In non-urban areas where there are no artificial buildings, PS targets are very rare. DS is a natural scatterer widely distributed in nature, such as bare soil, deserts, and sparse low vegetation areas. Therefore, the combination of PS and DS can greatly increase the density of observation points in non-urban areas. However, the backscattering ability of each basic scatterer in the DS pixel is relatively balanced, and there is no dominant strong scattering unit. DS is greatly affected by spatiotemporal incoherence, and the signal coherence is poor and the signal-to-noise ratio is low. To this end, it is necessary to re-estimate the phase of DS according to some methods. A common method is to construct the covariance matrix of DS, and then use the eigendecomposition method to obtain the eigenvector corresponding to the maximum eigenvalue and use it as the DS estimated phase.

[0003] With the launch of new satellites, SAR earth observation has entered an era of high temporal and spatial resolution, multi-polarization, and wide bandwidth, and has also generated unprecedented amounts of data. Taking the Sentinel-1 SAR satellite as an example, a standard image has nearly 1 billion pixels, and the conventional feature decomposition method is too time-consuming to meet the actual production application. PM has the characteristics of simple calculation and small memory usage, but the calculation efficiency and convergence are restricted by the iteration initial value and iteration termination strategy. Summary of the invention

[0004] The purpose of the present invention is to provide a distributed scatterer phase fast estimation method based on improved PM, taking the interference phase mean on homogeneous pixels SHP as the iterative initial value of PM, avoiding PM divergence, reducing the number of iterations and improving the convergence speed, and by constructing an adaptive termination strategy, ensuring that the PM algorithm can quickly obtain solutions for different distributed scatterers DS, greatly improving the computational efficiency and solving the problem of time-consuming computation.

[0005] To achieve the above object, the present invention provides a distributed scatterer phase fast estimation method based on improved PM, comprising the following steps: S1. Using the two-sample hypothesis test method, according to the amplitude information of the registered SAR image, the distributed scatterer DS is identified and the homogeneous pixel SHP of each distributed scatterer DS is marked; S2, using homogeneous pixels SHP to calculate the interference phase mean of the distributed scatterer DS; S3, construct the sample covariance matrix C using homogeneous pixels SHP; S4, using PM's matrix-vector multiplication to iteratively estimate the principal eigenvalues ​​and corresponding eigenvectors of the sample covariance matrix C; S5, construct an adaptive termination strategy for PM; S6, obtain the DS phase of the distributed scatterer, through Iterations, the result obtained at the end of the iteration That is the optimal estimation of the DS phase of the distributed scatterer ,Right now: .

[0006] Preferably, in S1, the two-sample hypothesis testing method is one of Kolmogorov-Smirnov and Anderson-Darling.

[0007] Preferably, in S2, the mean interference phase on the homogeneous pixel SHP is set to the initial estimate of the phase of the distributed scatterer DS , and used as the initial value of PM iteration to avoid PM divergence and reduce the number of iterations, the initial value of distributed scatterer DS phase iteration The calculation method is: ; in, Indicates The interference phase on the interference pair, M Indicates the number of homogeneous pixels SHP.

[0008] Preferably, in S3, constructing the covariance matrix C specifically includes: for each distributed scatterer DS, using its homogeneous pixel SHP to construct the covariance matrix C, that is: (2); Among them, the superscript represents the conjugate transpose of a matrix, is the complex reflectance of each homogeneous pixel SHP in the time series, Represents a set of homogeneous pixels SHP.

[0009] Preferably, S4 specifically includes the following steps: S41, the interference phase mean calculated by formula (1) As an initial estimate of the main eigenvector; S42, according to the matrix-vector multiplication of PM, iteratively update the main eigenvalue and eigenvector, after After iterations: (3); (4); in, Indicates The feature vector obtained by the iteration is Indicates The eigenvalues ​​obtained by iteration; S43, Normalized feature vector : (5); because , when the number of iterations As it approaches infinity, we have: (6); (7); in, is the theoretical principal eigenvalue of the sample covariance matrix C The corresponding feature vector.

[0010] Preferably, in S5, the number of iteration steps is It is impossible to take infinity, and the termination condition must be set for the iteration. The constructed adaptive iteration termination strategy is: (8); in, and Respectively represent Second and The main eigenvalue obtained by the iteration; Indicates the rate of change of the main eigenvalue. When the rate of change of the main eigenvalue is less than 1e-6, it means that the main eigenvalue and eigenvector have converged. Indicates that the maximum number of iterations is 50.

[0011] Therefore, the present invention adopts the above-mentioned distributed scatterer phase fast estimation method based on improved PM, and takes the interference phase mean on the homogeneous pixel SHP as the iterative initial value of PM, so as to avoid PM divergence, reduce the number of iterations and improve the convergence speed. By setting an adaptive termination strategy, it is ensured that the algorithm can quickly obtain the optimal solution for different distributed scatterers DS. Compared with the existing methods, the computational efficiency of the method of the present invention is greatly improved, which effectively solves the problem of excessively time-consuming calculations. Moreover, because the calculation is simple, it is particularly suitable for parallel calculations, which can further reduce the calculation time.

[0012] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0013] Figure 1 It is a flow chart of an embodiment of a distributed scatterer phase fast estimation method based on improved PM of the present invention; Figure 2 is a graph showing the relationship between the number of iterations and the DS quality of an embodiment of a distributed scatterer phase fast estimation method based on improved PM of the present invention; Figure 3 It is a comparison of the DS phase accuracy estimated by a distributed scatterer fast phase estimation method based on improved PM of the present invention and a conventional eigendecomposition method; Figure 4 It is a time-consuming comparison between a distributed scatterer phase fast estimation method based on improved PM of the present invention and a conventional eigendecomposition method. DETAILED DESCRIPTION

[0014] The technical solution of the present invention is further described below through the accompanying drawings and embodiments.

[0015] Unless otherwise defined, the technical terms or scientific terms used in the present invention should be understood by people with ordinary skills in the field to which the present invention belongs. The words "first", "second" and similar words used in the present invention do not indicate any order, quantity or importance, but are only used to distinguish different components. "Include" or "comprise" and similar words mean that the elements or objects appearing before the word include the elements or objects listed after the word and their equivalents, without excluding other elements or objects. "Connect" or "connected" and similar words are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect. "Up", "down", "left", "right" and the like are only used to indicate relative positional relationships. When the absolute position of the described object changes, the relative positional relationship may also change accordingly.

[0016] Embodiment 1 The present invention provides a distributed scatterer phase fast estimation method based on improved PM, the process of which is as follows: Figure 1 As shown, the following steps are included: S1. Using two-sample hypothesis test methods such as Kolmogorov-Smirnov and Anderson-Darling, the distributed scatterers DS are identified according to the amplitude information of the registered SAR image and the homogeneous pixels SHP of each distributed scatterer DS are marked.

[0017] S2. Use homogeneous pixels SHP to calculate the mean interference phase of the distributed scatterer DS. The essence of the distributed scatterer DS phase estimator is a filter that aims to improve the coherence quality of the distributed scatterer DS phase. Spatial mean filtering is a simple and effective method to improve the quality of interference phase. Its biggest drawback is that the filtering error is large when the scattering types of each pixel in the filter window are different. This defect can be suppressed if mean filtering is performed on homogeneous pixels SHP. Therefore, the mean interference phase on homogeneous pixels SHP is set as the initial estimate of the distributed scatterer DS phase. , and used as the initial value of PM iteration to avoid PM divergence and reduce the number of iterations, the initial value of distributed scatterer DS phase iteration The calculation method is: = (1); in, Indicates The interference phase on the interference pair, M Indicates the number of homogeneous pixels SHP.

[0018] S3. Use homogeneous pixel SHP to construct the sample covariance matrix C. For each distributed scatterer DS, use its homogeneous pixel SHP to construct the covariance matrix C, that is: (2); Among them, the superscript represents the conjugate transpose of a matrix, is the complex reflectance of each pixel on the SHP of each homogeneous pixel in the time series, Represents a set of homogeneous pixels SHP.

[0019] S4, using PM's matrix-vector multiplication, iteratively estimate the main eigenvalues ​​and corresponding eigenvectors of the sample covariance matrix C. The specific process is: S41, the interference phase mean calculated by formula (1) As an initial estimate of the main eigenvector; S42, according to the matrix-vector multiplication of PM, iteratively update the main eigenvalue and eigenvector, after After iterations: (3); (4); in, Indicates The feature vector obtained by the iteration is Indicates The eigenvalues ​​obtained by iteration; S43, Normalized feature vector : (5); because , when the number of iterations As it approaches infinity, we have: (6); (7); in, is the theoretical principal eigenvalue of the sample covariance matrix C The corresponding feature vector.

[0020] S5. Set the adaptive termination strategy. Iteration steps It is impossible to take infinity, so the termination condition must be set for the iteration. The constructed adaptive iteration termination strategy is: (8); in, and Respectively represent Second and The main eigenvalue obtained by the iteration; Indicates the rate of change of the main eigenvalue. When the rate of change of the main eigenvalue is less than 1e-6, it means that the main eigenvalue and eigenvector have converged. It shows that the maximum number of iterations is 50. For high-quality distributed scatterers, DS can quickly converge to the exact solution, while for low-quality distributed scatterers, DS can terminate the iteration in time.

[0021] set up The basis is as follows: In practical applications, after estimating the phase of the distributed scatterer DS, the quality of the estimation needs to be evaluated: (9); in, and denote the observed phase and estimated phase of the distributed scatterer DS, respectively. N represents the number of SAR images, and Re(·) represents the real part operation. indicates the quality of the estimate, the larger the value, the better the estimate; When it is less than a certain value (such as ), indicating that the estimation quality is poor, and the distributed scatterer DS should be discarded. Actual SAR data show that, if Figure 2 Shown: For those For the distributed scatterer DS, only a few iterations (about 35 times) are needed to converge; The distributed scatterer DS is a low-quality target that should be discarded. Even if it does not converge, it will not affect the subsequent results.

[0022] S6, obtain the DS phase of the distributed scatterer, obtained when the iteration ends That is the optimal estimation of the DS phase of the distributed scatterer ,Right now: .

[0023] In order to verify the effectiveness and computational efficiency of the method described in this embodiment, the method described in this embodiment and the conventional eigendecomposition method are used to estimate the DS phase of the distributed scatterer, and the phase estimation accuracy and computational efficiency of the two methods are compared.

[0024] 1. Accuracy assessment The simulation data contains 100 scenes of SAR data. According to the imaging parameters of Sentinel-1A, the time interval is set to 12 days, the wavelength is set to 55.6 mm, the deformation rate is set to 1 mm / a, and the SHP is set to 300. The method described in this embodiment and the conventional feature decomposition method are used to estimate the DS phase of the distributed scatterer, and each is repeated 100,000 times, and then the root mean square error (RMSE) is calculated to evaluate the accuracy. The experimental results are shown in Figure 2. Figure 3 As shown, it can be found that the accuracy of the method described in this embodiment is completely consistent with the conventional method, which shows that the DS phase of the distributed scatterer estimated by the method described in the present invention is reliable.

[0025] 2. Computational efficiency test The computation time is related to the number of SAR images. Five groups of SAR data sets with 100, 150, 200, 250 and 300 scenes were simulated. Then the DS phase was estimated using the method described in this embodiment and the conventional feature decomposition method. The calculation was repeated 10,000 times and the cumulative time was calculated. The experimental results are shown in Figure 2. Figure 4 As shown, it can be found that as the number of SARs increases, the computational efficiency of the method described in this embodiment is significantly higher than that of the conventional eigendecomposition method, and the computational time is only 1 / 10 of that of the conventional eigendecomposition method.

[0026] Therefore, the present invention adopts the above-mentioned distributed scatterer phase fast estimation method based on improved PM, and takes the interference phase mean on the homogeneous pixel SHP as the iterative initial value of PM, so as to avoid PM divergence, reduce the number of iterations and improve the convergence speed. By setting an adaptive termination strategy, it is ensured that the algorithm can quickly obtain the optimal solution for different distributed scatterers DS. Compared with the existing methods, the computational efficiency of the method described in the present invention is greatly improved, which effectively solves the problem of excessively time-consuming calculations, and is particularly suitable for parallel computing.

[0027] Finally, it should be noted that the above embodiments are only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that they can still modify or replace the technical solution of the present invention with equivalents, and these modifications or equivalent replacements cannot cause the modified technical solution to deviate from the spirit and scope of the technical solution of the present invention.

Claims

1. A distributed scatterer phase fast estimation method based on improved PM, characterized in that: The following steps are involved: S1. Using the two-sample hypothesis test method, according to the amplitude information of the registered SAR image, the distributed scatterer DS is identified and the homogeneous pixel SHP of each distributed scatterer DS is marked; S2, using homogeneous pixels SHP to calculate the interference phase mean of the distributed scatterer DS; S3, construct the sample covariance matrix C using homogeneous pixels SHP; S4, using PM's matrix-vector multiplication to iteratively estimate the principal eigenvalues ​​and corresponding eigenvectors of the sample covariance matrix C; S5, construct an adaptive termination strategy for PM; S6, obtain the DS phase of the distributed scatterer, through Iterations, the result obtained at the end of the iteration That is the optimal estimation of the DS phase of the distributed scatterer ,Right now: .

2. The method for fast phase estimation of distributed scatterers based on improved PM according to claim 1, characterized in that: In S1, the two-sample hypothesis testing method is one of Kolmogorov-Smirnov and Anderson-Darling.

3. The method for fast phase estimation of distributed scatterers based on improved PM according to claim 1, characterized in that: In S2, the mean interference phase on the homogeneous pixel SHP is set as the initial estimate of the distributed scatterer DS phase. , and used as the initial value of PM iteration to avoid PM divergence and reduce the number of iterations, the initial value of distributed scatterer DS phase iteration The calculation method is: ； in, Indicates The interference phase on the interference pair, M Indicates the number of homogeneous pixels SHP.

4. The method for fast phase estimation of distributed scatterers based on improved PM according to claim 3 is characterized in that: In S3, the covariance matrix C is constructed as follows: for each distributed scatterer DS, the covariance matrix C is constructed using its homogeneous pixel SHP, that is: ； Among them, the superscript represents the conjugate transpose of a matrix, is the complex reflectance of each homogeneous pixel SHP in the time series, Represents a set of homogeneous pixels SHP.

5. The method for fast phase estimation of distributed scatterers based on improved PM according to claim 4, characterized in that: S4 specifically includes the following steps: S41, the interference phase mean calculated by formula (1) As an initial estimate of the main eigenvector; S42, according to the matrix-vector multiplication of PM, iteratively update the main eigenvalue and eigenvector, after After iterations: (3)； (4)； in, Indicates The feature vector obtained by the iteration is Indicates The eigenvalues ​​obtained by iteration; S43, Normalized feature vector : (5)； because , when the number of iterations As it approaches infinity, we have: (6)； (7)； in, is the theoretical principal eigenvalue of the sample covariance matrix C The corresponding feature vector.

6. The method for fast phase estimation of distributed scatterers based on improved PM according to claim 5, characterized in that: In S5, the number of iterations k It is impossible to take infinity, and the termination condition must be set for the iteration. The constructed adaptive iteration termination strategy is: (8)； in, and Respectively represent Second and The main eigenvalue obtained by the iteration; Indicates the rate of change of the main eigenvalue. When the rate of change of the main eigenvalue is less than 1e-6, it means that the main eigenvalue and eigenvector have converged. Indicates that the maximum number of iterations is 50.

Images