A method for monitoring three-dimensional deformation of a building

Abstract

The application provides a kind of housing three-dimensional deformation monitoring method.The monitoring method specifically includes the following steps:S1.multitrack SAR data acquisition, and the pretreatment is carried out to the multitrack SAR data collected;S2.the multitrack SAR image data after pretreatment in S1 step is respectively carried out PS-InSAR time sequence processing, and the LOS deformation and terrain error of each track PS point are inverted;S3.InSAR point cloud standardization extraction;S4.multisource InSAR point cloud fusion;S5.roof point cloud segmentation;S6.three-dimensional deformation decomposition, finally obtain the east direction, north direction, vertical direction three-dimensional deformation rate and time sequence of the whole roof of house, realize the accurate monitoring of three-dimensional deformation of house.The application solves the problems of single track InSAR data integrity, low multisource point cloud fusion precision, and inability to obtain three-dimensional deformation of house, and realizes large-scale, high-precision three-dimensional deformation monitoring of dense housing group in city and town.

Description

[0001] Technical Field

[0002] This invention relates to the field of surveying and remote sensing technology, specifically to a method for monitoring three-dimensional deformation of buildings based on InSAR point clouds.

[0003] Background Technology

[0004] Housing is a core urban infrastructure for people's livelihood, and its structural deformation is directly related to the safety of people's lives and property. With the acceleration of urbanization in my country, the number of old houses continues to increase, and housing deformation accidents caused by factors such as geological disasters, underground engineering construction, and extreme weather occur frequently. Therefore, monitoring the structural safety of houses has become a core task of urban safety operation and maintenance.

[0005] Traditional methods for monitoring building deformation mainly include precise leveling, total station polar coordinate measurement, static GNSS measurement, and ground-based lidar scanning. Precise leveling and total station measurements offer high accuracy, but their spatial coverage is discrete, limiting monitoring to a limited number of points. Furthermore, they are labor-intensive, inefficient, and difficult to implement for large-scale continuous monitoring. Static GNSS measurements enable all-weather automated monitoring, but their high equipment cost, low spatial resolution, and inability to acquire detailed building deformation characteristics are significant drawbacks. Ground-based lidar scanning can obtain three-dimensional point clouds of buildings, but data processing is complex, revisit cycles are long, and dynamic monitoring is difficult. None of these methods can meet the dynamic, automated, large-scale, and detailed monitoring needs of densely populated urban housing clusters in modern cities.

[0006] Synthetic Aperture Radar Interferometry (InSAR) technology, with its advantages of all-weather operation, wide coverage, high spatial resolution, non-contact operation, and low cost, has gradually become a core technology for monitoring surface deformation and building structures. Permanent Scatterer InSAR (PS-InSAR) technology can accurately extract deformation information from rigid permanent scatterers such as roofs, walls, and structures, generating high-density InSAR deformation point clouds, providing a new technical approach for building deformation monitoring. However, existing InSAR-based building deformation monitoring technologies have the following core technical shortcomings:

[0007] 1) Insufficient data integrity of single orbit: Due to the influence of observation geometry, terrain occlusion, shadow overlapping, and temporal decoherence, single orbit SAR satellites have incomplete PS point cloud extraction of buildings, missing data of key facades and corner areas, and cannot reflect the overall deformation shape of buildings.

[0008] 2) Low accuracy of multi-source point cloud fusion: No adaptive registration and fusion algorithm for multi-track InSAR point clouds has been established. The coordinate reference and time reference of point clouds from different tracks are not consistent. Noise and redundant data are not effectively removed, and the accuracy and reliability of deformation observation are greatly reduced.

[0009] 3) Lack of three-dimensional deformation decomposition capability: Traditional technology can only obtain one-dimensional deformation in the radar LOS direction, and cannot calculate the real three-dimensional deformation components of the building in the east, north and vertical directions, making it difficult to support the assessment of safety indicators such as building structure tilt and uneven settlement. Summary of the Invention

[0010] To address the shortcomings of existing technologies, this invention proposes a method for monitoring three-dimensional deformation of buildings based on InSAR point clouds. This method can solve problems such as poor integrity of single-track InSAR data, low accuracy of multi-source point cloud fusion, and inability to obtain three-dimensional deformation of buildings, thereby achieving large-scale, high-precision three-dimensional deformation monitoring of dense urban housing clusters.

[0011] To address the aforementioned technical problems, this invention provides a method for monitoring three-dimensional deformation of buildings, characterized by the following specific steps:

[0012] S1. Multi-track SAR data acquisition; Acquire time-series SAR image data with multiple orbits and multiple incident angles for the monitored buildings to ensure the integrity of PS point cloud extraction; and preprocess the acquired multi-source data to provide high-quality and stable multi-track SAR image data for subsequent PS-InSAR processing;

[0013] S2. PS-InSAR temporal data processing: For the preprocessed multi-orbit SAR image data in step S1, PS-InSAR temporal processing is carried out to invert the LOS deformation and terrain error of each orbit PS point;

[0014] S3. Standardized extraction of InSAR point clouds: The multi-track PS point clouds are processed in a unified manner to generate standardized InSAR deformed point clouds, which are used to prepare for subsequent fusion.

[0015] S4. Multi-source InSAR point cloud fusion: Construct an adaptive fusion algorithm for multi-source InSAR point clouds based on weighted iterative nearest point to achieve high-precision registration and fusion of point clouds from different orbits;

[0016] S5. Rooftop point cloud segmentation; Based on the fused complete InSAR point cloud, intelligent segmentation of the rooftop point cloud is achieved, and non-rooftop point clouds such as ground, vegetation, and structures are removed.

[0017] S6. Three-dimensional deformation decomposition: Based on the deformation observations of the multi-track LOS direction, a weighted least squares three-dimensional deformation decomposition model is constructed to solve the three-dimensional deformation components of the building in the east, north and vertical directions. Finally, the three-dimensional deformation rate and time series of the roof of the entire building in the east, north and vertical directions are obtained to achieve accurate monitoring of the three-dimensional deformation of the building.

[0018] The preferred technical solution of the present invention is as follows: In step S1, high-resolution SAR satellite data in C-band and X-band are selected, including COSMO, TerraSAR-X, and Gaofen-3. Multi-orbit data with a combination of ascending and descending orbits are preferred. The incident angle covers 20°~45°, the time span is not less than 6 months, and the number of time-series images is not less than 15 scenes.

[0019] The preferred technical solution of this invention is as follows: In step S1, the multi-track data covers the same monitoring area, with a spatial resolution better than 3m, a revisit period of less than 30 days, and no large-scale cloud or atmospheric interference, which meets the requirements of PS-InSAR time series processing; the preprocessing of multi-source data is to perform radiometric calibration, multi-view processing, and noise reduction on the original SAR image to eliminate image radiometric errors and noise.

[0020] A further technical solution of the present invention: The specific process of step S2 is as follows:

[0021] S201. Image registration and interferometric pair generation: Based on the main image, complete the accurate registration of multi-track temporal SAR images with a registration accuracy better than 1 / 2 pixel; select the image combination with the best spatiotemporal baseline to generate interferometric pairs and eliminate the phase error introduced by the spatiotemporal baseline;

[0022] S202. Selection of high-stability PS points: Based on the amplitude deviation index and phase stability threshold, rigid permanent scatterers of the building are screened, high-stability PS points of the building are selected, and a high-density PS point set is constructed. The high-stability PS points include the roof concrete surface, walls, and metal components.

[0023] S203. Interferometric phase decomposition and deformation inversion; Separate topographic phase, deformation phase, atmospheric phase and noise phase, and construct a PS point arc segment phase model; With the goal of maximizing the temporal coherence coefficient, complete atmospheric phase correction and phase unwrapping by weighted least squares inversion of single-point DEM error and deformation rate;

[0024] S204. Precision assessment and coarse screening: Calculate the deformation error of each PS point, remove low-precision points with errors greater than 1mm, and retain high-precision PS points.

[0025] Further technical solutions of the present invention:

[0026] The formula for the phase model of the PS point arc segment in step S203 is as follows:

[0027] (1)

[0028] In the formula, The total interference phase; For terrain phase; For deformed phase; Atmospheric phase; Noise phase;

[0029] In step S203, with the goal of maximizing the temporal coherence coefficient, the atmospheric phase correction and phase unwrapping are completed by weighted least squares inversion of single-point DEM errors and deformation rates, as follows: For the PS point network arc segment , No. Observation phase model of amplitude difference interferogram:

[0030] (2)

[0031] In the formula, For arc segment In the The observed phase of the amplitude difference plot; The radar wavelength; The radar incident angle; Slope distance; For the first Vertical baseline of the map; For the first Map time baseline; For arc segment Relative DEM error; For arc segment Relative linear deformation rate; Atmospheric phase; Noise phase;

[0032] The optimal solution for the arc parameters of the phase model is to maximize the time-domain coherence coefficient. Take the maximum value Objective function:

[0033] (3)

[0034] In the formula, The total number of differential interferograms; Imaginary unit .

[0035] Based on this, an overdetermined equation is constructed using the relative parameters of the arc segments, and the single-point DEM error is inverted using least squares. With deformation rate :

[0036] (4)

[0037] in, This is the arc-point correlation matrix (column sum is 0); Represents a single-point parameter vector; This is the optimal parameter vector for the arc segment.

[0038] The weighted least squares solution of the above formula can be expressed as:

[0039] (5)

[0040] Wherein, is a weight matrix, determined by , the higher the coherence coefficient, the greater the weight.

[0041] The preferred technical solution of the present invention: the specific steps of standardized extraction of InSAR point cloud in said S3 are as follows:

[0042] S301. Temporal gross error elimination; the criterion is used to eliminate gross error points in the deformation time series, and the judgment formula is:

[0043] (6)

[0044] wherein is the deformation value of the -th PS point of the -th orbit at time; is the mean value of the deformation sequence of the PS point; is the standard deviation of the deformation sequence of the PS point;

[0045] if the deformation value at a certain moment satisfies the above formula, it is determined as a gross error and eliminated;

[0046] S302. Coordinate datum unification; converting all orbital PS point clouds to the same geodetic coordinate system and projected coordinate system:

[0047] (7)

[0048] In the formula, ( ) are spatial rectangular coordinates; are latitude and longitude; is the radius of curvature of the prime vertical; is the eccentricity of the earth ellipsoid;

[0049] S303. Time datum unification; based on UTC time, interpolation processing is performed on deformation time sequences of different orbits, so that the time stamps of all point clouds correspond one to one; a cubic spline interpolation method is adopted, and the interpolated deformation sequence is expressed as:

[0050] (8)

[0051] Wherein, is the length of the original time series; is a cubic spline basis function;

[0052] S304. Point Cloud Recovery and Attribute Integration: The terrain phase error of each PS point is added to the external reference DEM, the DEM is height corrected, and the PS points are geocoded again to obtain InSAR point cloud results. Each PS point is assigned three-dimensional coordinates, LOS deformation, deformation rate, deformation mean square error, and orbital parameter attributes to generate a standardized InSAR deformed point cloud.

[0053] A further technical solution of the present invention: The specific process of multi-source InSAR point cloud fusion in step S4 is as follows:

[0054] S401. Initial Registration: Manually select more than 3 feature points from the multi-source InSAR deformed point cloud, calculate global geometric features such as weighted centroid and elevation statistical mean, without relying on point-by-point matching; complete the coarse registration of the multi-track point cloud by aligning the global geometric features, and eliminate the overall coordinate offset;

[0055] S402. Fine Registration: Construct a WICP registration objective function, using the deformation error of PS points as weights, and iteratively optimize the rotation matrix R and translation vector T to achieve fine alignment of the point cloud with a registration error of less than 0.2m; the function formula is:

[0056] (9)

[0057] In the formula, It is a rotation matrix; It is a translation vector; Source point cloud coordinates; The target point cloud coordinates;

[0058] The weighting coefficient is determined by the mean square error of the PS point deformation; the smaller the mean square error, the greater the weight.

[0059] S403. High-precision registration of point pairs is achieved by minimizing the objective function, realizing fine alignment and coordinate fusion of multi-source InSAR point clouds.

[0060] The preferred technical solution of this invention: The rooftop point cloud segmentation in step S5 specifically includes the following steps:

[0061] S501. Building outline extraction: Based on high-resolution optical images, extract the vector outline of buildings in the monitoring area, and use the vector outline to perform individual segmentation of the large-scale InSAR point cloud to obtain the InSAR point cloud of each building.

[0062] S502. Point cloud normal vector calculation; By calculating the point cloud normal vector and filtering by direction, non-roof points can be accurately eliminated. Among them, the roof point cloud normal vector is approximately vertically upward, the wall normal vector is horizontal, and the vegetation / noise point cloud normal vector is disordered and scattered.

[0063] S503. Roof point cloud segmentation: Remove wall point clouds based on normal vector direction, remove ground point clouds using height threshold (h≥0.5m), remove scattered points with normal vector variance greater than the threshold, retain roof connected regions with concentrated normal vector directions, delete small isolated points through connected component analysis, and finally obtain a clean roof InSAR deformed point cloud.

[0064] The preferred technical solution of this invention: In step S6, for each pair of PS corresponding points, establish the mapping relationship between LOS deformation and deformation in the East (E), North (N), and Vertical directions; combine three or more ascending and descending orbit observations to construct a weighted least squares three-dimensional deformation decomposition model, and calculate the deformation rate and time series of the entire roof area of ​​each building in the East, North, and Vertical directions; LOS direction mapping relationship with three-dimensional components:

[0065] (10)

[0066] In the formula, For LOS deformation; The radar incident angle; This is the radar azimuth angle; These are the deformation components in the east, north, and vertical directions;

[0067] By combining the LOS deformation components of multiple orbits, a three-dimensional deformation decomposition model is calculated based on weighted least squares:

[0068] (11)

[0069] In the formula, It is a three-dimensional deformation component matrix; The coefficient matrix of the observation equation; This is the weight matrix; Let be the LOS-directed deformation observation vector matrix.

[0070] The present invention has the following beneficial effects:

[0071] (1) This invention improves the integrity of point cloud data: By jointly collecting multi-track SAR data, it effectively overcomes the problem of missing point clouds caused by terrain occlusion, shadows, and overlapping of single-track data, greatly improving the integrity of building PS point cloud extraction and fully reflecting the overall deformation characteristics of the building.

[0072] (2) This invention improves the accuracy of point cloud fusion: a multi-orbit InSAR point cloud adaptive fusion algorithm based on WICP is constructed, which takes into account the difference in observation accuracy of point clouds on different orbits, and realizes high-precision registration and fusion of point clouds on different orbits.

[0073] (3) The present invention realizes accurate calculation of three-dimensional deformation of houses: Based on the deformation observation values ​​of multi-track LOS direction, a weighted least squares three-dimensional deformation decomposition model is constructed, which can accurately obtain the three-dimensional deformation characteristics of houses in the east, north and vertical directions, and provide a reliable basis for the structural safety assessment of houses.

[0074] (4) The present invention has constructed a complete chain of technology system: forming a complete technical process of multi-track data acquisition - temporal InSAR processing - multi-source point cloud extraction and fusion - roof point cloud intelligent segmentation - three-dimensional deformation decomposition, realizing the fine monitoring of three-dimensional deformation of buildings. Attached Figure Description

[0075] Figure 1 This is a technical flowchart of the present invention;

[0076] Figure 2 This is an image showing the effect of InSAR point cloud fusion.

[0077] Figure 3 This is a framework diagram for 3D deformation monitoring of buildings based on InSAR point clouds. Detailed Implementation

[0078] 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.

[0079] The embodiment provides a method for monitoring 3D deformation of buildings based on InSAR point clouds, targeting buildings in a specific city. The process is as follows: Figure 1 Specifically, it includes the following steps:

[0080] S1. Multi-track SAR data acquisition; Acquire time-series SAR image data with multiple orbits and multiple incident angles for the monitored buildings to ensure the integrity of PS point cloud extraction; and preprocess the acquired multi-source data to provide high-quality and stable multi-track SAR image data for subsequent PS-InSAR processing;

[0081] A specific city was selected as the monitoring target, encompassing brick-concrete and frame-structured buildings, constructed 20-30 years ago, with high building density and complex facades. The city is also undergoing subway and other construction projects, which could easily cause uneven settlement and tilting deformations in the buildings, presenting typical monitoring needs. The preferred satellite data source was the X-band TerraSAR-X satellite, supplemented by high-resolution SAR data such as COSMO. An ascending and descending orbit combination was used for observation, with an incident angle coverage of 30°–38°. The time series parameters included a time span greater than 8 months, a cumulative number of valid images greater than 20, a revisit period of 11 days, and a spatial resolution less than 3m. Auxiliary data included high-resolution optical images of the same region, high-precision DEM data, and building vector contour maps. The raw SAR images underwent preprocessing: radiometric calibration, multi-look processing, speckle noise removal, and topographic radiometric correction were performed to eliminate radiometric errors and system noise, resulting in a high-quality time-series SAR image dataset.

[0082] Step 2. PS-InSAR Time Series Data Processing; For the preprocessed multi-orbit SAR image data from Step S1, PS-InSAR time series processing is performed to invert the LOS deformation and terrain error of each orbit PS point; the specific steps are as follows:

[0083] S201. Image registration and interferometric pair generation: Based on the master image, sub-pixel level registration of temporal SAR images is completed, with a registration accuracy better than 1 / 2 pixel; Based on spatiotemporal baseline threshold screening, the optimal single master image interferometric pair is generated to minimize the phase error introduced by the spatiotemporal baseline.

[0084] S202. PS Point Optimization and Point Set Construction: Based on the dual thresholds of amplitude deviation index less than 0.25 and phase stability greater than 0.85, rigid permanent scatterers such as concrete slabs on building roofs and metal components are screened to construct a high-density PS point set; 68,000 PS points are extracted from track 1, 72,000 PS points are extracted from track 2, and 60,000 PS points are extracted from track 3, for a total of 200,000 initial PS points.

[0085] S203. Interferometric phase decomposition and deformation inversion; Separate topographic phase, deformation phase, atmospheric phase and noise phase, and construct a PS point arc segment phase model; With the goal of maximizing the temporal coherence coefficient, complete atmospheric phase correction and phase unwrapping by weighted least squares inversion of single-point DEM error and deformation rate;

[0086] The PS point interferometric phase consists of four parts: terrain phase, deformation phase, atmospheric phase, and noise phase. Its formula is:

[0087]

[0088] In the formula, The total interference phase; For terrain phase; For deformed phase; Atmospheric phase; Noise phase;

[0089] For PS point network arcs , No. Observation phase model of amplitude difference interferogram:

[0090]

[0091] In the formula, For arc segment In the The observed phase of the amplitude difference plot; The radar wavelength; The radar incident angle; Slope distance; For the first Vertical baseline of the map; For the first Map time baseline; For arc segment Relative DEM error; For arc segment Relative linear deformation rate; Atmospheric phase; Noise phase;

[0092] The optimal solution for the arc parameters of the phase model is to maximize the time-domain coherence coefficient. Take the maximum value Objective function:

[0093]

[0094] In the formula, The total number of differential interferograms; Imaginary unit .

[0095] Based on this, an overdetermined equation is constructed using the relative parameters of the arc segments, and the single-point DEM error is inverted using least squares. With deformation rate :

[0096]

[0097] in, This is the arc-point correlation matrix (column sum is 0); Represents a single-point parameter vector; This is the optimal parameter vector for the arc segment.

[0098] The weighted least squares solution of the above formula can be expressed as:

[0099]

[0100] Wherein, is a weight matrix, which is determined by , and the higher the coherence coefficient is, the larger the weight is.

[0101] S204. Accuracy evaluation and rough screening; calculate the mean square error of deformation of each PS point, eliminate low-precision points with error greater than 1 mm, and retain 180,000 high-precision PS points to ensure the quality of the initial point cloud.

[0102] S3. Standardized extraction of InSAR point cloud; perform unified processing on multi-track PS point clouds to generate standardized InSAR deformation point clouds, and prepare for subsequent fusion; the specific steps are as follows:

[0103] S301. Temporal gross error elimination; adopt the 3σ criterion to traverse the deformation time series of each PS point, eliminate 3200 abnormal gross error points with a gross error elimination rate of 1.7%, and purify the temporal deformation data and point cloud quality; the judgment formula is:

[0104]

[0105] Wherein is the deformation value of the -th PS point on the -th track at moment; is the mean value of the deformation sequence of this PS point; is the standard deviation of the deformation sequence of this PS point;

[0106] if the deformation value at a certain moment satisfies the above formula, it is determined as a gross error and eliminated;

[0107] S302. Coordinate reference unification; convert all orbital PS point clouds from the radar coordinate system to the WGS-84 geodetic coordinate system plus UTM projection coordinate system, complete the unification of three-dimensional coordinates, and eliminate coordinate deviation between tracks;

[0108]

[0109] In the formula, ( ) are space rectangular coordinates; are longitude and latitude; is the radius of curvature of the prime vertical; is the eccentricity of the earth ellipsoid;

[0110] S303. Time Standardization: Using UTC standard time as the standard, cubic spline interpolation is employed to resample the deformed sequence, ensuring complete alignment of the timestamps in the ascending and descending orbit point clouds, and achieving uniform time resolution. The deformed sequence after cubic spline interpolation is represented as follows:

[0111]

[0112] in, The length of the original time series; Cubic spline basis functions;

[0113] S304. Point cloud restoration and attribute integration; correct the PS point DEM error to the reference DEM and complete geocoding relocation; assign each PS point with attributes such as three-dimensional coordinates, LOS deformation, deformation rate, deformation error, orbit number, incident angle, and azimuth angle to generate standardized ascending and descending InSAR deformed point clouds.

[0114] S4. Multi-source InSAR point cloud fusion: Construct an adaptive fusion algorithm for multi-source InSAR point clouds based on weighted iterative nearest points to achieve high-precision registration and fusion of point clouds from different orbits, such as... Figure 2 The weighted iterative nearest point (WICP) algorithm is used to achieve high-precision registration and fusion of multi-track point clouds, as detailed below:

[0115] S401. Initial Registration: Manually select three or more feature points from the multi-source InSAR deformed point cloud, calculate global geometric features such as weighted centroid and elevation statistical mean, without relying on point-by-point matching. By aligning global geometric features, complete the coarse registration of the multi-track point cloud and eliminate overall coordinate offset;

[0116] S402. Fine Registration; Construct a WICP objective function, using the mean square error of PS point deformation as weight (the smaller the mean square error, the larger the weight), iteratively optimize the rotation matrix R and translation vector T to achieve fine alignment of the point cloud, with a registration error of less than 0.2m. The WICP objective function formula is:

[0117] (9)

[0118] In the formula, It is a rotation matrix; It is a translation vector; Source point cloud coordinates; The target point cloud coordinates;

[0119] The weighting coefficient is determined by the mean square error of the PS point deformation; the smaller the mean square error, the greater the weight.

[0120] S403. Point cloud fusion: Unify the coordinates and attributes of multi-source PS point clouds, and generate 176,800 complete InSAR deformed point clouds after fusion, with the point cloud density increased by 250% compared to a single orbit.

[0121] S5. Rooftop Point Cloud Segmentation; Based on the fused complete InSAR point cloud, intelligent segmentation of the rooftop point cloud is achieved, removing non-rooftop point clouds such as ground, vegetation, and structures; the specific steps are as follows:

[0122] S501. Monomerization and Segmentation

[0123] Based on the vector outline of the building, the large-scale point cloud is individually cropped to separate the point cloud of a single building and remove invalid points outside the area.

[0124] S502. Normal Vector Selection

[0125] Calculate the point cloud normal vector: the roof point cloud normal vector is approximately vertically upward, the wall point cloud normal vector is horizontal, and the vegetation / noise point cloud normal vector is scattered; remove 100,800 point clouds from the wall, ground, and vegetation according to the normal vector direction.

[0126] S504. Threshold and Connectivity Optimization

[0127] By setting a height threshold h greater than 0.5m, 22,000 ground points were removed. Small isolated points were deleted through connected component analysis, resulting in a rooftop PS point cloud of 54,000 points containing multiple orbits.

[0128] S6. Three-dimensional deformation decomposition; Based on multi-track LOS deformation observations, a weighted least squares three-dimensional deformation decomposition model is constructed to solve the three-dimensional deformation components in the east, north, and vertical directions of the building. Finally, the three-dimensional deformation rates and time series of the entire roof area in the east, north, and vertical directions are obtained, enabling accurate monitoring of the building's three-dimensional deformation. For example... Figure 3 Based on multi-track LOS deformation, the three-dimensional true deformation of a single building is calculated. The specific steps are as follows:

[0129] S601. Rooftop Matching with Same Name Points

[0130] By using the nearest neighbor method, PS corresponding points and their corresponding LOS direction observations of multiple orbits are selected. Pairs of corresponding points with large distances are removed, resulting in 2 sets of corresponding point pairs, with 3 points in each set.

[0131] S602. Construction of 3D Decomposition Model

[0132] For each pair of PS corresponding points, establish the mapping relationship between LOS deformation and deformation in the East (E), North (N), and Vertical (U) directions. Combine three or more ascending and descending orbit observations to construct a weighted least squares three-dimensional deformation decomposition model. LOS direction mapping relationship with three-dimensional components:

[0133]

[0134] In the formula, For LOS deformation; The radar incident angle; This is the radar azimuth angle; These are the deformation components in the east, north, and vertical directions;

[0135] S603. Three-dimensional deformation calculation;

[0136] By combining the LOS deformation components of multiple orbits, a three-dimensional deformation decomposition model is calculated based on weighted least squares:

[0137]

[0138] In the formula, It is a three-dimensional deformation component matrix; The coefficient matrix of the observation equation; This is the weight matrix; The LOS-directed deformation observation vector matrix;

[0139] The calculation yields the deformation rates and time series of the entire roof area of ​​each building in the east, north, and vertical directions, enabling precise monitoring of the three-dimensional deformation of the buildings.

[0140] The above description is merely one embodiment of the present invention, and while it is detailed and specific, it should not be construed as limiting the scope of the invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the appended claims.

Claims

1. A method for monitoring three-dimensional deformation of a building, characterized in that, The specific steps are as follows: S1. Multi-track SAR data acquisition; Acquire time-series SAR image data with multiple orbits and multiple incident angles for the monitored buildings to ensure the integrity of PS point cloud extraction; and preprocess the acquired multi-source data to provide high-quality and stable multi-track SAR image data for subsequent PS-InSAR processing; S2. PS-InSAR temporal data processing: For the preprocessed multi-orbit SAR image data in step S1, PS-InSAR temporal processing is carried out to invert the LOS deformation and terrain error of each orbit PS point; S3. Standardized extraction of InSAR point cloud: The multi-track PS point cloud is processed in a unified manner to generate a standardized InSAR deformed point cloud, which is used to prepare for subsequent fusion. S4. Multi-source InSAR point cloud fusion: Construct an adaptive fusion algorithm for multi-source InSAR point clouds based on weighted iterative nearest point to achieve high-precision registration and fusion of point clouds from different orbits; S5. Rooftop point cloud segmentation; Based on the fused complete InSAR point cloud, intelligent segmentation of the rooftop point cloud is achieved, and non-rooftop point clouds such as ground, vegetation, and structures are removed. S6. Three-dimensional deformation decomposition: Based on the deformation observations of the multi-track LOS direction, a weighted least squares three-dimensional deformation decomposition model is constructed to solve the three-dimensional deformation components of the building in the east, north and vertical directions. Finally, the three-dimensional deformation rate and time series of the roof of the entire building in the east, north and vertical directions are obtained to achieve accurate monitoring of the three-dimensional deformation of the building.

2. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that: The data selection in step S1 uses high-resolution SAR satellite data in C-band and X-band, including COSMO, TerraSAR-X, and Gaofen-3. Priority is given to multi-orbit data with a combination of ascending and descending orbits, with an incident angle coverage of 20°~45°, a time span of no less than 6 months, and no fewer than 15 time-series images.

3. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that: In step S1, multi-track data covers the same monitoring area, with a spatial resolution better than 3m, a revisit period of less than 30 days, and no large-scale cloud cover or atmospheric interference, meeting the requirements for PS-InSAR time series processing. The preprocessing of multi-source data involves radiometric calibration, multi-view processing, and noise reduction of the original SAR images to eliminate image radiometric errors and noise.

4. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that: The specific process for step S2 is as follows: S201. Image registration and interferometric pair generation: Based on the main image, complete the accurate registration of multi-track temporal SAR images with a registration accuracy better than 1 / 2 pixel; select the image combination with the best spatiotemporal baseline to generate interferometric pairs and eliminate the phase error introduced by the spatiotemporal baseline; S202. Selection of high-stability PS points: Based on the amplitude deviation index and phase stability threshold, rigid permanent scatterers of the building are screened, high-stability PS points of the building are selected, and a high-density PS point set is constructed. The high-stability PS points include the roof concrete surface, walls, and metal components. S203. Interferometric phase decomposition and deformation inversion; Separate topographic phase, deformation phase, atmospheric phase and noise phase, and construct a PS point arc segment phase model; With the goal of maximizing the temporal coherence coefficient, complete atmospheric phase correction and phase unwrapping by weighted least squares inversion of single-point DEM error and deformation rate; S204. Precision assessment and coarse screening: Calculate the deformation error of each PS point, remove low-precision points with errors greater than 1mm, and retain high-precision PS points.

5. The method for monitoring three-dimensional deformation of a building according to claim 4, characterized in that, The formula for the phase model of the PS point arc segment in step S203 is as follows: （1） In the formula, The total interference phase; For terrain phase; For deformed phase; Atmospheric phase; Noise phase; In step S203, with the goal of maximizing the temporal coherence coefficient, the atmospheric phase correction and phase unwrapping are completed by weighted least squares inversion of single-point DEM errors and deformation rates, as follows: For the PS point network arc segment , No. Observation phase model of amplitude difference interferogram: （2） In the formula, For arc segment In the The observed phase of the amplitude difference plot; The radar wavelength; The radar incident angle; Slope distance; For the first Vertical baseline of the map; For the first Map time baseline; For arc segment Relative DEM error; For arc segment Relative linear deformation rate; Atmospheric phase; Noise phase; The optimal solution for the arc parameters of the phase model is to maximize the time-domain coherence coefficient. Take the maximum value Objective function: （3） In the formula, This represents the total number of differential interferograms. Imaginary unit ; Based on this, an overdetermined equation is constructed using the relative parameters of the arc segments, and the single-point DEM error is inverted using least squares. With deformation rate : （4） in, This is the arc-point correlation matrix (column sum is 0); Represents a single-point parameter vector; This represents the optimal parameter vector for the arc segment. The weighted least squares solution to the above equation can be expressed as: （5） in, The weight matrix is ​​formed by... It is confirmed that the higher the coherence coefficient, the greater the weight.

6. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that, The specific steps for InSAR point cloud normalization extraction in S3 are as follows: S301. Timing gross error removal; using... The criterion for removing gross errors from deformed time series is as follows: （6） in For the first The first orbit A PS point at The deformation value at time; This is the mean of the deformation sequence at point PS; The standard deviation of the deformation sequence at the PS point; If the deformation value at a certain moment satisfies the above formula, it is judged as a gross error and discarded; S302. Unify coordinate datum; transform all orbital PS point clouds to the same geodetic coordinate system and projected coordinate system: （7） In the formula, ( ) are spatial rectangular coordinates; are latitude and longitude; is the radius of curvature in the prime vertical; is the eccentricity of the earth ellipsoid; S303. Time standardization; Using UTC time as the reference, interpolation is performed on the deformation time series of different orbits to ensure a one-to-one correspondence of timestamps in all point clouds; cubic spline interpolation is used, and the interpolated deformation sequence is represented as follows: （8） in, The length of the original time series; Cubic spline basis functions; S304. Point Cloud Recovery and Attribute Integration: The terrain phase error of each PS point is added to the external reference DEM, the DEM is height corrected, and the PS points are geocoded again to obtain InSAR point cloud results. Each PS point is assigned three-dimensional coordinates, LOS deformation, deformation rate, deformation mean square error, and orbital parameter attributes to generate a standardized InSAR deformed point cloud.

7. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that, The specific process of multi-source InSAR point cloud fusion in step S4 is as follows: S401. Initial Registration: Manually select more than 3 feature points from the multi-source InSAR deformed point cloud, calculate global geometric features such as weighted centroid and elevation statistical mean, without relying on point-by-point matching; complete the coarse registration of the multi-track point cloud by aligning the global geometric features, and eliminate the overall coordinate offset; S402. Fine Registration: Construct a WICP registration objective function, using the deformation error of PS points as weights, and iteratively optimize the rotation matrix R and translation vector T to achieve fine alignment of the point cloud with a registration error of less than 0.2m; the function formula is: （9） In the formula, It is a rotation matrix; It is a translation vector; Source point cloud coordinates; The target point cloud coordinates; The weighting coefficient is determined by the mean square error of the PS point deformation; the smaller the mean square error, the greater the weight. S403. High-precision registration of point pairs is achieved by minimizing the objective function, realizing fine alignment and coordinate fusion of multi-source InSAR point clouds.

8. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that, Roof point cloud segmentation in step S5 specifically includes the following steps: S501. Building outline extraction: Based on high-resolution optical images, extract the vector outline of buildings in the monitoring area, and use the vector outline to perform individual segmentation of the large-scale InSAR point cloud to obtain the InSAR point cloud of each building. S502. Point cloud normal vector calculation; By calculating the point cloud normal vector and filtering by direction, non-roof points can be accurately eliminated. Among them, the roof point cloud normal vector is approximately vertically upward, the wall normal vector is horizontal, and the vegetation / noise point cloud normal vector is disordered and scattered. S503. Roof point cloud segmentation: Remove wall point clouds based on normal vector direction, remove ground point clouds using height threshold (h≥0.5m), remove scattered points with normal vector variance greater than the threshold, retain roof connected regions with concentrated normal vector directions, delete small isolated points through connected component analysis, and finally obtain a clean roof InSAR deformed point cloud.

9. The method for monitoring three-dimensional deformation of a building according to claim 1, characterized in that, In step S6, for each pair of PS corresponding points, a mapping relationship is established between LOS deformation and deformation in the East (E), North (N), and Vertical directions. Combining observations from three or more ascending and descending orbits, a weighted least squares three-dimensional deformation decomposition model is constructed to calculate the east, north, and vertical deformation rates and time series for the entire roof area of ​​each building. The mapping relationship between the LOS direction and the three-dimensional components is as follows: （10） In the formula, For LOS deformation; The radar incident angle; This is the radar azimuth angle; These are the deformation components in the east, north, and vertical directions; By combining the LOS deformation components of multiple orbits, a three-dimensional deformation decomposition model is calculated based on weighted least squares: （11） In the formula, It is a three-dimensional deformation component matrix; The coefficient matrix of the observation equation; This is the weight matrix; Let be the LOS-directed deformation observation vector matrix.

Images