Differential tomographic sar staging sparse four-dimensional imaging method with anti-collapse rollback mechanism

Abstract

The application belongs to the technical field of radar signal processing, and relates to a differential layering SAR hierarchical sparse four-dimensional imaging method with a collapse-preventing rollback mechanism. OMP Fast rough positioning of scatterers is realized; the number of scatterers is determined based on CFAR+GLRT of complex degrees of freedom correction; dynamic dictionary clipping and FR L1 super-resolution fine adjustment is performed, a collapse-preventing rollback safety net taking amplitude ratio as a criterion is constructed, and when the algorithm degenerates, the safety net automatically rolls back to a robust rough solution; finally, a high-precision three-dimensional point cloud and a deformation rate field are generated. The application solves the problems of easy collapse of a traditional layering SAR under a pathological matrix, point cloud void and many false peaks, and takes into account imaging efficiency, super-resolution accuracy and system robustness, and is suitable for four-dimensional imaging scenes of urban buildings, terrain subsidence, disaster monitoring and the like of satellite-borne / machine-borne SAR.

Description

Technical Field

[0001] This invention belongs to the field of high-dimensional synthetic aperture radar imaging processing technology, and relates to a differential tomographic SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism. It is a method for differential tomographic synthetic aperture radar (D-TomoSAR) four-dimensional imaging of complex urban scenes, and is particularly suitable for accurate three-dimensional reconstruction and deformation rate inversion of multi-layered overlapping targets in complex urban scenes (such as large man-made buildings). Background Technology

[0002] Differential analytical synthetic aperture radar (D-TomoSAR) can resolve multiple scatterers distributed along the elevation direction within a single resolution cell by jointly processing SAR images with multiple baselines and multiple temporal phases, and simultaneously estimate the deformation rate of each scatterer, thus achieving four-dimensional imaging (three-dimensional space + deformation rate).

[0003] Currently, the mainstream processing framework is a hierarchical approach represented by CS-GLRT (Compressed Sensing-Generalized Likelihood Ratio Test). The typical process is: first, a coarse estimate is performed using compressed sensing algorithms (such as L1 norm minimization); then, the generalized likelihood ratio test (GLRT) is used to determine the number of scatterers (model order selection); and finally, the least squares (LS) method is used for fine parameter estimation. However, this framework and its derivatives have the following four significant drawbacks:

[0004] (i) Existing technologies are prone to numerical degradation and lack recovery mechanisms under locally ill-conditioned matrices.

[0005] In the fine estimation stage of LS (Laser-Solved Point Cloud), when the local observation matrix becomes ill-conditioned due to uneven baseline distribution or phase errors, it can lead to numerical solution collapse or severe distortion of amplitude estimation. Current technologies lack any detection or recovery mechanism for this, directly resulting in the complete loss of information for the corresponding pixels and the formation of irreversible data holes in the final 3D point cloud. This is the most fatal functional deficiency in existing frameworks.

[0006] (ii) Statistical mismatch: The selection of the GLRT model has theoretical flaws.

[0007] Existing GLRT methods neglect the fundamental characteristic that SAR observation data are complex signals when calculating test statistics. Their F-distribution degrees of freedom parameters are empirically set using the real domain, which is mismatched with the Hermitian properties of the complex signal covariance matrix. This leads to inaccurate decision thresholds, uncontrollable false alarm probabilities at low signal-to-noise ratios, and a tendency to generate spurious scatterers.

[0008] (iii) Homogeneous relay: Fine estimation lacks super-resolution capability

[0009] The fine estimation part of the existing framework uses the LS algorithm, which itself does not have the ability to exceed the Rayleigh resolution limit. When multiple scatterers are close in elevation, LS cannot effectively resolve them. At the same time, the coarse estimation of L1 minimization and the fine estimation of LS are methodologically homogeneous and cannot take advantage of the complementary advantages of heterogeneous algorithms in the coarse and fine stages.

[0010] (iv) Waste of resources: The per-pixel full-data processing strategy is inefficient

[0011] Existing methods all perform coarse-to-fine processing on all pixels. However, in urban scenes, the number of complex pixels that overlap is usually no more than 15%, and this "one-size-fits-all" strategy results in a significant waste of computational resources. Summary of the Invention

[0012] The purpose of this invention is to overcome the shortcomings of existing technologies and provide a differential SAR hierarchical sparse four-dimensional imaging method with an anti-collapse rollback mechanism. This invention aims to solve the following technical problems:

[0013] (1) Solve the problem of point cloud holes: Establish an active fault tolerance mechanism that can automatically detect and recover when the fine-tuning algorithm crashes, thus avoiding data loss;

[0014] (2) Solve the problem of statistical inaccuracy: Correct the degree of freedom parameters of GLRT to make them strictly match the statistical characteristics of complex SAR data, so as to achieve accurate false alarm control;

[0015] (3) Solve the resolution bottleneck: Introduce heterogeneous algorithms with super-resolution capabilities in the fine-tuning stage to improve the ability to separate adjacent overlapping targets;

[0016] (4) Solve the problem of computational efficiency: Establish a selective fine-tuning triggering strategy to concentrate high-precision calculations on only a few complex pixels.

[0017] To achieve the objectives of this invention, the following technical solutions are adopted.

[0018] A hierarchical sparse four-dimensional imaging method for differential SAR with anti-collapse rollback mechanism is proposed. This method constructs the entire imaging process into three functional layers plus an independent safety mechanism: a global coarse scan layer, a model selection layer, a local fine-tuning layer, and an anti-collapse safety net layer. The anti-collapse safety net layer serves as an independent feedback branch connected in parallel with the local fine-tuning layer. The specific steps are as follows:

[0019] Step 1: Global coarse scan—preliminary localization of fast scatterers based on vectorized CS_OMP

[0020] For the preprocessed multi-baseline SAR high coherence point dataset, a joint elevation-deformation rate phase dictionary matrix is ​​constructed for each pixel. The vectorized compressed sensing orthogonal matching pursuit (CS_OMP) algorithm is used to quickly reconstruct all pixels one by one, extracting the coarse elevation position, initial deformation rate, amplitude vector, and iterative residual signal of the candidate scatterers in each pixel.

[0021] Step 2: Model Selection—Determination of the Number of Scatterers in CFAR and GLRT Based on Complex Degrees of Freedom Correction

[0022] This step differs from the conventional GLRT approach of constructing the likelihood ratio using the original observation vectors. Instead, the residual signal generated in step one is input into the model selection module. In this module, the Generalized Likelihood Ratio Test (GLRT) under constant false alarm rate (CFAR) control is used to determine the number K of actual overlapping scatterers within each pixel.

[0023] The key technological innovation in this step lies in introducing a complex variable degree-of-freedom correction logic for complex signals: when calculating the GLRT test statistic (i.e., the F-test statistic corresponding to the residual power reduction; the F-test is a statistical test method based on the F-distribution, commonly used to compare variances, test the overall significance of regression models, or compare whether there are significant differences between multiple population means in analysis of variance (ANOVA), the traditional real-domain degree-of-freedom parameters are corrected to a form that considers the Hermitian properties of the complex covariance matrix— , Where N is the number of SAR image baselines, The number of scatterers before the current iteration assumption. The number of scatterers is added as a new hypothesis. Based on this adjusted degree of freedom, the F-test statistic is constructed as (mean square error of gain / mean square error of residuals):

[0024] ,

[0025] in: , This indicates the decrease in residual power. This represents the residual signal of the scatterer. Combined with a preset constant false alarm rate threshold... ( =1e -3 The threshold is determined by calculating the inverse function of the F-distribution. (The inverse function of the F-distribution is implemented using a built-in MATLAB function package) Function to calculate a given false alarm rate The decision threshold is set below. Precisely eliminate those caused by The false noise peaks detected are used to obtain the true number of scatterers K for each pixel.

[0026] Step 3: Local Fine-tuning—Local adjustments based on dynamic dictionary pruning High-resolution parameter estimation and crash-resistant security network

[0027] This step is the core of the invention, employing a strategy of "first judgment, then diversion, and main fine-tuning of the parallel safety network".

[0028] For pixels identified as having overlapping properties in step 2 (i.e., K≥2), the local fine-tuning module is triggered; for pixels with K<2, the settings from step 1 are directly applied. The result is output as the final solution.

[0029] The local fine-tuning module consists of two parallel branches:

[0030] (1) Main operation branch — High-resolution solution

[0031] Step 1 Centered on the output coarse elevation position, adjacent grids are dynamically extracted from the original joint phase dictionary matrix to form a dimension-reduced sub-dictionary. Then, the FISTA-based Reweighted L1 Minimization algorithm is invoked. The pixel is then subjected to super-resolution fine-tuning to obtain high-precision estimates of the scatterer elevation, deformation rate, and amplitude.

[0032] The reweighted L1 minimization algorithm based on the FISTA framework is an improved sparse reconstruction method built upon the standard L1 sparse minimization. Standard L1 sparse minimization introduces L1 norm constraints into the objective function to make the amplitude coefficient vector of the scatterer as sparse as possible, thereby selecting a small number of effective scatterers from multiple candidate units in the elevation-deformation rate joint parameter space. However, standard L1 sparse minimization applies the same level of penalty to all candidate scattering units, which can easily lead to the compression of the amplitude of the real strong scatterers, while also having limited ability to suppress weak noise components, sidelobes, and spurious scattering peaks. To address these issues, this invention introduces a reweighting mechanism on top of standard L1 sparse minimization, extending the unified L1 sparse constraint to a weighted L1 sparse constraint. The weights are adaptively updated based on the amplitude of the sparse solution obtained in each iteration, so that candidate scattering units with larger amplitudes receive smaller weights and weaker penalties, while candidate scattering units with smaller amplitudes receive larger weights and stronger penalties. Ultimately, the real scatterers can be effectively preserved, and noise components and spurious peaks can be further suppressed.

[0033] Meanwhile, this invention employs the FISTA framework for fast solution of the weighted L1 sparse optimization problem. The FISTA framework first updates the current estimation result using the gradient of the data fidelity term to reduce the error between the sparse reconstruction result and the observed signal; then, it applies sparse constraints to the gradient-updated result through soft thresholding shrinkage, suppressing candidate scattering results with smaller amplitudes; finally, it combines the results from both iterations through momentum acceleration to improve the iterative convergence speed. Therefore, The algorithm retains the characteristic of standard L1 sparsity minimization being applicable to sparse scatterer inversion, while improving the stability of scatterer peak extraction and the estimation accuracy of super-resolution parameters through reweighting mechanism and FISTA fast iterative framework.

[0034] The algorithm flow is as follows:

[0035] 1. Input the observed signal y, dictionary matrix A, iteration number iter and regularization parameter μ. First, organize the observed signal into column vectors, initialize the sparse solution x=0, and set the initial weight W=1. Since the initial weight of all candidate scattering units is 1, the first iteration is equivalent to the standard L1 sparse solution.

[0036] 2. Estimate the Lipschitz constant Lf of the gradient term based on the dictionary matrix A, and determine the gradient descent step size: This step size is used to control the stability of subsequent gradient updates;

[0037] 3. Under the current weight W, use FISTA to solve the L1 norm weighted sparse optimization problem, and obtain the current sparse solution x through gradient descent, soft threshold shrinkage and momentum acceleration;

[0038] 4. Update weights and output the final sparse solution: After completing one round of FISTA solving, update the weights according to the magnitude of the current sparse solution: δ represents a stabilization factor to avoid a denominator of zero. This update mechanism reduces the weight and penalty for strong scatterers, while increasing the weight and penalty for weak noise components. After multiple weighted iterations, the final sparse reconstruction result x is output, which is used for subsequent scatterer peak extraction, elevation estimation, and deformation rate estimation.

[0039] With CS OMP Methodologically, they belong to a heterogeneous relationship: CS OMP As a greedy algorithm, it provides a robust global coarse-scan base by tracking residuals through successive projections, while As a convex optimization algorithm, it achieves local super-resolution fine-tuning by leveraging FISTA acceleration and reweighted iteration. The two constitute a heterogeneous cascade system that balances efficiency, accuracy, and anti-collapse robustness.

[0040] (2) Bypass monitoring and rollback branch - anti-collapse safety net

[0041] ① Criterion Selection Principle

[0042] This mechanism uses the scatterer amplitude ratio as the criterion, instead of phase, complex correlation coefficient, or solution process parameters (such as the number of iterations or the magnitude of the residual), for the following reasons:

[0043] Physical truth invariance: The amplitude of a ground object scatterer possesses physical truth properties and does not change drastically with the iteration of solution algorithms. The amplitude of a real target should remain highly consistent from coarse to fine estimation.

[0044] Robust Benchmarking: Coarse Scan Algorithm (CS) OMP Although its accuracy is limited, the solution process is highly robust and always returns a physically robust solution. This invention uses this as a benchmark reference value to characterize the "physical truth" and to verify the output of the fine-tuning algorithm.

[0045] Typical characteristics of numerical degradation: The most significant symptom of numerical degradation in high-precision algorithms under ill-conditioned matrix conditions is that the amplitude of the solved scatterer deviates significantly from the true physical value, manifesting as amplitude collapse (too small) or abnormal amplification. The amplitude ratio can directly and quantitatively capture this anomaly, which is more sensitive and essential than indirectly monitoring the solution process parameters.

[0046] Avoiding phase and position sensitivity: When fine-tuning degrades, its phase information is usually randomly disordered, and elevation estimation may "jump". Using scalar amplitude for comparison on the same elevation grid naturally avoids the risk of misjudgment caused by phase misalignment and position drift of complex signals.

[0047] ② Monitoring and rollback execution logic

[0048] Reference extraction: CS in step 1 OMP After the solution is completed, the amplitude information of each scatterer within the pixel on different elevation grids is retained. In step 3 After the solution is completed, the maximum amplitude value corresponding to the scatterer solution output is extracted, and at the same time CS is extracted OMP The maximum amplitude value corresponding to the coarse-scan scatterer is used as the reference amplitude A OMP .

[0049] Degradation determination: Calculate the amplitude ratio , where: A FR is the maximum amplitude value of the FR L1 fine adjustment. The ratio R is compared with the preset safety threshold T

[0050] Threshold setting logic: The value range of the safety threshold T is set to 0.5 to 0.9, and its logical boundaries are as follows

[0051] Lower limit 0.5 (physical negation line): If the amplitude of the fine adjustment is less than half of the coarse-scan solution, it means that the target signal power loss exceeds 6 dB, which physically constitutes a self-negation of the "true existence" of the target. This boundary excludes any situation that can be interpreted as a normal error and is an absolute condition for triggering rollback

[0052] Upper limit 0.9 (normal tolerance line): Ideal fine adjustment will bring more concentrated energy of the scattering center, making the main peak amplitude slightly increase relative to the coarse estimate (the ratio is usually ≥1.0). 0.9 constitutes the lower boundary of the normal random fluctuation range. Below this value, it indicates that the signal energy of the high-precision algorithm is less than that of the low-precision coarse estimate, which is physically unreasonable and constitutes the first warning for determining degradation

[0053] Decision gray area: The range of 0.5 - 0.9 constitutes a physically interpretable and engineeringly operable transition zone for reliably distinguishing the normal fluctuations and substantial degradation of the algorithm

[0054] Automatic rollback: When R < T, it is determined that FR L1 has numerical degradation at this pixel, and the system automatically discards the FR​​​​​​​​​​

[0056] Step 4: Point Cloud Generation and Post-processing

[0057] Sidelobe suppression: When the elevation distance between any two scatterers within the same pixel is less than (0.5 * Rayleigh resolution of the tomographic direction), only the one with the larger amplitude is retained; the final solution of all pixels in the dataset (including the CS rolled back by the safety net) is collected. OMP Solution: CS not triggered for fine-tuning OMP Solution and normal output FR L1 (Fine adjustment) is performed to calculate the spatial three-dimensional coordinates of elevation values ​​and deformation rates. Then, intensity filtering is applied: the entire point cloud is sorted in descending order of amplitude, and a hard truncation is performed, retaining only the top 95% of strong scattering points. The final output is a clean 3D point cloud map of the building and a deformation rate distribution map.

[0058] Compared with the prior art, the present invention has the following advantages:

[0059] (1) This invention establishes a crash-prevention security network mechanism, which solves the problem of data gaps in a hierarchical framework.

[0060] An active fault-tolerance mechanism is introduced into the D-TomoSAR hierarchical processing framework. This is achieved by real-time monitoring of FR. L1 Fine-tuning and CS OMP The amplitude ratio of the coarse scan solution can instantly detect the numerical degradation of the high-precision algorithm and automatically roll back to a robust coarse scan solution. This mechanism provides a "soft landing" capability for the entire imaging process—even if the local ill-conditioned matrix causes the optimal algorithm to fail, the system can still output physically reasonable estimation results, ensuring the integrity and continuity of the 3D reconstruction of the building.

[0061] (2) This invention significantly improves the statistical rigor of model order selection and the accuracy of false alarm control.

[0062] By introducing a degree-of-freedom correction logic based on the Hermitian property of the complex covariance matrix, the F-distribution threshold in the GLRT test is made to accurately match the true statistical distribution of SAR complex observation data. Compared with the existing approach of directly using real-domain degrees of freedom, this method can strictly control the false alarm probability in low signal-to-noise ratio regions to a preset CFAR level (0.1%), fundamentally reducing the number of pseudo-scatterers in 3D point clouds.

[0063] (3) The heterogeneous algorithm of this invention achieves the optimal balance between computational efficiency and elevation resolution.

[0064] This invention employs a tiered triggering strategy: computationally intensive rendering is enabled only on a small number of complex pixels where model selection determines there is overlap (K≥2). L1 High-resolution algorithm; for the vast majority of simple pixels containing only a single scatterer, the computationally efficient CS algorithm is directly adopted.OMP Solution. Meanwhile, local fine-tuning employs dynamic dictionary pruning technology, further compressing the matrix size for high-precision solutions. This design ensures super-resolution stripping capability in complex overlapping areas such as building corners while avoiding the enormous time overhead of global high-precision calculations.

[0065] (4) This invention enhances the robustness and engineering practicality of the overall imaging process.

[0066] Through CS OMP Coarse sweep, complex degree of freedom correction GLRT, dynamic dictionary FR L1 Through the four-fold synergy of fine-tuning and anti-collapse safety net, this invention constructs a robust imaging framework that is highly immune to uneven baseline distribution, low signal-to-noise ratio, and phase error, significantly improving the engineering practicality and reliability of D-TomoSAR technology in complex urban scenarios. Attached Figure Description

[0067] Figure 1 This is a flowchart of a specific embodiment of the method described in this invention;

[0068] Figure 2 This is a cropped main image intensity thumbnail used in an embodiment of the present invention;

[0069] Figure 3 This is a distribution amplitude diagram of the high coherence points extracted in an embodiment of the present invention;

[0070] Figure 4 This is the final 3D point cloud map of the building generated in this embodiment of the invention;

[0071] Figure 5 This is the final building deformation rate diagram generated in this embodiment of the invention;

[0072] Figure 6 This is a spatial distribution diagram of scatterers superimposed on the main image according to an embodiment of the present invention. Detailed Implementation

[0073] The present invention will be further described in conjunction with the accompanying drawings and embodiments.

[0074] In this embodiment, as Figure 1 As shown, 18 images of TerraSAR-X / TanDEM-X rising orbit staring spotlight mode were used, and after cropping, they covered buildings in a typical urban center area, spanning approximately 26 months.

[0075] Step 1: Data Preprocessing and Dictionary Construction

[0076] Registration, cropping, multi-look processing, interferometry, coherent scattering point selection, differential interferometry, and phase error compensation were performed on 18 multi-baseline SAR images. The tenth image was selected as the master image, and high-coherence points were extracted as the dataset for subsequent processing, determining their azimuth-range boundary range. All preprocessing was performed using the professional SAR data processing software GAMMA.

[0077] The tomographic search range is set to -10 meters to 80 meters, with a search step size of 1 meter (91 grids in total). The known basic radar system parameters are: wavelength lambda = 0.0311 meters, main image proximity R = 613672.8765 meters, baseline span ΔBv = 454.9573 meters, and the Rayleigh resolution limit in the tomographic direction is Rayleigh. resolution =(lambda*R) / (2*ΔBv), the final calculated result is 20.97m, and the set step size is much smaller than the Rayleigh resolution; the deformation rate is set to search range from -2 mm / year to 2 mm / year, and the search step size is 0.2 mm / year (a total of 21 grids). Based on this, a joint phase dictionary matrix Φ is constructed for each high coherence pixel, with dimensions of 18 rows (number of baselines) × (91 × 21) columns. The main image intensity thumbnail used in this embodiment is as follows. Figure 2 As shown, the extracted high coherence points are as follows: Figure 3 As shown, high coherence points are displayed as white highlights, while the remaining pixels are black backgrounds.

[0078] Step 2: CS OMP Global coarse scan and CFAR+GLRT model order selection

[0079] For each pixel within the boundary, call vectorized CS. OMP Algorithm. During the algorithm iteration process, each time, the dictionary atom with the largest absolute value of the inner product with the current residual signal is selected, and its corresponding elevation and deformation rate index are recorded. The amplitude coefficient of the selected atom is updated using the least squares method until the residual power decreases no longer significantly or the preset maximum number of iterations is reached. Output: coarse position of the candidate scatterer, amplitude vector, and final residual signal for each pixel.

[0080] The residual signal is input into the model order selection module. For each iteration, it is assumed that the number of scatterers is K. old Increase to K new In this case, calculate the relative decrease in residual power and construct the F-test statistic:

[0081] ,

[0082] ,

[0083] ,

[0084] ,

[0085] The degrees of freedom parameters of the numerator and denominator are respectively adopted using complex-corrected parameters. and N=18 represents the number of baselines.

[0086] Constant False Alarm Rate (CFAR) detection follows the Neyman-Pearson criterion: maximizing the detection probability given a false alarm probability Pfa. In radar signal processing, the typical range for Pfa detection in engineering applications of CFAR is 1e. -3 up to 1e -6 The range. In urban building scenes, the number of overlapping targets is usually 2-4, the signal is relatively stable, and the signal-to-noise ratio is moderate. When the number of baselines is limited (such as N=18 in this embodiment), the constant false alarm rate Pfa=1e is set. -3 This is sufficient to distinguish between genuine overlays and noise spurious peaks, eliminating the need for a more stringent PFA. Setting the PFA too stringently (e.g., 1e) would be problematic. -4 1e -5 ), through the inverse function of the F distribution Calculating the decision threshold will significantly increase the threshold value, causing real weak scatterers to be misclassified as noise spurious peaks and rejected because they cannot cross the excessively high threshold. If the actual calculated F-statistic... If the value exceeds this threshold, the newly added candidate scatterer is considered statistically significant and retained; otherwise, it is judged as a noise spurious peak and the iteration is terminated. The final result is the number K of the actual overlapping scatterers within each pixel.

[0087] Step 3: Local FR_L1 fine-tuning and crash prevention safety net execution

[0088] Traverse all pixels in the image: For pixels with K=0 or K=1, directly set CS OMP The result serves as the final solution; for pixels with K ≥ 2, the local fine-tuning module is triggered.

[0089] The local fine-tuning module first uses CS OMP The given elevation positions of each scatterer are centered on CS. OMP The algorithm's super-resolution capability is limited. To avoid coarse scanning errors, a local sub-dictionary, centered on the elevation grid and containing 5 to 10 times the tomographic search step size (the search range is close to Rayleigh resolution), is extracted from the original joint phase dictionary matrix Φ. The tomographic search step size is 1 meter; to ensure computational efficiency, this example uses 5 times the step size, i.e., M=5. Furthermore, to prevent the search range from going out of bounds, the minimum search range of the dictionary is limited to the first grid, and the maximum search range is limited to the maximum number of grids.

[0090] This pruning operation significantly reduces the dictionary dimension to be solved from 18×(91×21) to 18×(11×21), thus substantially reducing the FR. L1 The computational burden.

[0091] FR L1 The algorithm is based on the FISTA acceleration framework, and its core iterative formula is:

[0092]

[0093] The algorithm mainly achieves signal separation and denoising through "inner and outer loops":

[0094] The inner loop (FISTA accelerated solution) uses gradient descent and soft thresholding to approximate the real signal and adds a Nesterov momentum term (i.e., each iteration uses the difference from the previous iteration as inertia compensation), which greatly improves the convergence speed.

[0095] Outer loop (dynamic reweighting): After each inner loop, the penalty weight formula is updated based on the obtained signal amplitude. The core idea is that the larger the amplitude of a point (real building), the smaller the penalty, and the smaller the amplitude of a point (background noise / sidelobe), the larger the penalty, eventually approaching zero.

[0096] Where y represents the input signal: the multi-baseline SAR observation vector after pixel differential processing; A represents the dictionary matrix: a joint phase grid (elevation-rate) containing spatial and temporal baselines; x represents the signal to be determined: the solved high-resolution backscatter profile (non-zero values ​​represent real buildings); μ represents the regularization parameter: controlling the denoising intensity. The larger the value, the more background clutter is filtered out, and the cleaner the inversion result; δ represents the stability factor: a very small constant (e.g., 0.05) to prevent the denominator from being zero when calculating the weight W, thus avoiding numerical collapse.

[0097] The specific execution logic of the crash prevention security network:

[0098] ① Instantiation of decision parameters

[0099] For each overlay pixel that triggers fine-tuning (K≥2), extract FR. L1 The maximum amplitude value A in the output scatterer profile FR and CS OMP The maximum amplitude value A of the coarse sweep solution OMP Calculate the ratio R.

[0100] ② Experimental selection for safety threshold

[0101] The logic for setting the safety threshold T has been explained in the aforementioned technical solution section. In this embodiment, to prioritize system robustness and sensitivity to the detection of local ill-conditioned matrix degradation, T=0.9 is preferably used as the safety threshold. A rollback is triggered when R < T, with CS... OMP The solution is the final output; when R≥T, FR is used. L1 High-resolution solution.

[0102] It should be noted that 0.5 ≤ T ≤ 0.9 are all within the acceptable safety threshold range of this invention. Lower thresholds focus more on reducing the probability of erroneous rollbacks, while higher thresholds focus more on improving degradation detection sensitivity and system stability. This embodiment uses T=0.9 as the preferred engineering parameter, which does not constitute a limitation on the scope of protection of this invention.

[0103] Step 4: Point Cloud Post-processing and Visualization

[0104] The scattering volume solutions from the full image output are processed in 3D coordinates. Sidelobe suppression is performed: if the elevation distance between any two scattering volumes within the same pixel is less than 10 grid steps (0.5 * Rayleigh resolution in the tomographic direction, i.e., 10 meters), only the one with the larger amplitude is retained. Subsequently, all effective scattering volumes are sorted in descending order of their amplitude values, and the points with the highest cumulative energy percentage (i.e., retaining strong scattering volumes and filtering out weak noise points) are selected to generate the final 3D building point cloud.

[0105] By mapping elevation values ​​and deformation rate values ​​using color bands, a 3D point cloud map and a deformation rate map of the building are obtained. The experimental results are as follows: Figure 4 , Figure 5 As shown. Figure 4 This represents the final 3D point cloud map of the building generated by this invention. X represents the azimuth coordinate, Y represents the range coordinate, and Z represents the elevation coordinate, with different elevation values ​​indicated by color. Figure 4 It can be seen that the building's maximum height is around 42 meters, which is not much different from the actual height. The overall outline structure is complete, and there are obvious elevation undulations at the building's edges and facade area. No obvious large-area voids appear in the point cloud. Figure 5 This indicates the final building deformation rate diagram generated by the present invention, with three-dimensional spatial coordinates and... Figure 4 Consistent with each other, the point cloud color represents the magnitude of the deformation rate; the deformation rate distribution at the bottom of the building is around 0, with no significant settlement, while significant deformation is observed in localized areas of the building facade. Single scatterers (K=1), two scatterers (K=2), and complex superimposed scatterers (K≥3) are overlaid on SAR images (red, green, blue) to obtain a scatterer distribution map. Experimental results are as follows... Figure 6 As shown. Figure 6This is a spatial distribution map of scatterers superimposed on the main image, used to show the spatial distribution characteristics of scatterer pixels of different complexities. It can be seen that about 80% of the pixels in this mid-rise building are single scatterer points, and only one effective target is encountered when the radar line of sight passes through them; double scatterers and complex superimposed scatterer points are mainly concentrated at the intersection of the building's edge facade and in areas with strong scattering.

Claims

1. A differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism, characterized in that, Includes the following steps: Step 1, Global Coarse Scan: Based on the preprocessed multi-baseline SAR high coherence point dataset, a joint phase dictionary matrix of elevation and deformation rate is constructed. A vectorized compressed sensing orthogonal matching pursuit algorithm is used for pixel-by-pixel reconstruction to extract the coarse location of candidate scatterers, initial deformation rate values, amplitude vectors, and iterative residual signals. Step 2, Model Selection: Using the iterative residual signal as input, the numerator and denominator degrees of freedom of the F-test statistic are obtained through complex degree-of-freedom correction. A decision threshold is calculated using the inverse function of the F-distribution, combined with a preset constant false alarm rate, to eliminate noise spurious peaks in the vectorized compressed sensing orthogonal matching pursuit algorithm and determine the number of true scatterers per pixel. Step 3, Local Fine-tuning: For K≥2... For each pixel, the dictionary is dynamically cropped centered on the coarse scan position, and super-resolution fine-tuning is performed using a reweighted L1 minimization algorithm based on the FISTA framework. Simultaneously, a crash prevention rollback safety net is enabled, and the algorithm is judged to have degraded based on the amplitude ratio. If degraded, it automatically rolls back to the coarse solution of the vectorized compressed sensing orthogonal matching pursuit algorithm. For pixels with K < 2, the result of the vectorized compressed sensing orthogonal matching pursuit algorithm is directly used as the final solution output. Step 4: Point cloud generation and post-processing: Sidelobe suppression is performed on the scatterer elevation solution extracted from the same pixel, and the three-dimensional coordinate calculation is completed. The final solutions of all pixels are collected, and intensity filtering is performed to output a three-dimensional point cloud and deformation rate distribution map.

2. The differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The preprocessing includes registration, cropping, multi-view, interferometry, coherent scattering point selection, differential interferometry processing, and phase error compensation.

3. The differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The complex degrees of freedom correction adopts the form of adapting the Hermitian properties of the complex covariance matrix, molecular degrees of freedom Denominator degrees of freedom Where: N is the number of SAR image baselines, Assuming the number of forward scatterers for iteration, To add the number of hypothetical scatterers.

4. The differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The preset constant false alarm rate Pfa is set to 1×10. -3 .

5. A differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The criterion for the anti-crash rollback safety net is the amplitude ratio: , where: A FR To achieve the maximum adjustment result, A OMP The maximum amplitude of the coarse scan result; when R < T, it is judged as degenerate, and T is the judgment threshold, T∈[0.5, 0.9].

6. The differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The dynamic cropping dictionary is a sub-dictionary of local grids centered on the elevation position obtained by global coarse scanning, which reduces the solution dimension and computational load.

7. A differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, The final solution includes the CS rolled back via the crash prevention rollback safety net. OMP Solution: CS not triggered for fine-tuning OMP Solution and normal output FR L1 Fine-tuning.

8. A differential SAR hierarchical sparse four-dimensional imaging method with anti-collapse rollback mechanism according to claim 1, characterized in that, Sidelobe suppression: When the elevation distance between any two scatterers in the same pixel is less than 0.5 * Rayleigh resolution of the tomographic direction, only the one with the larger amplitude is retained; Intensity filtering: The point cloud of the whole image is arranged in descending order of amplitude, and hard truncation is performed to retain only the strong scatterer points with the first 95% intensity.

Images