Point cloud registration method for topographic change area driven by spatial similarity and multi-threshold optimization

Abstract

The invention belongs to the technical field of point cloud registration, and discloses a topographic change area point cloud registration method driven by spatial similarity and multi-threshold optimization. The method comprises the following steps: firstly, carrying out initial registration and voxelization segmentation on two stages of ground point cloud data; secondly, calculating spatial geometric features of each voxel, constructing a GMM-based point cloud spatial similarity measurement model, and preliminarily classifying all voxels into deformable voxels and stable voxels; carrying out recursive subdivision on the deformable voxels to screen stable voxels in the deformable voxels, and defining all marked stable voxels as initial stable voxels; thirdly, for points contained in the initial stable voxel, combining distance distribution characteristics from the points to a fitting plane, and further screening a stable region by adopting an adaptive dynamic threshold distance constraint strategy; and finally, implementing fine registration based on the extracted stable region point cloud. According to the method, the point cloud registration difficulty caused by the unknown position of the terrain stable area in surface deformation monitoring is effectively solved, and reliable technical support is provided for geological disaster monitoring.

Description

Technical Field

[0001] This invention belongs to the field of point cloud registration technology, specifically relating to a point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization. Background Technology

[0002] Three-dimensional laser point clouds, with their advantages of comprehensive, high-precision, and rich detail, have become a core data source for monitoring surface deformation. However, due to limitations such as differences in equipment performance, data acquisition errors, and local independent coordinate systems, multi-temporal point clouds often suffer from problems such as inconsistent references and systematic biases. Therefore, strict registration of multi-temporal point clouds is necessary before topographic change monitoring.

[0003] In recent years, scholars both domestically and internationally have conducted extensive research on point cloud registration, which can be categorized into coarse registration and fine registration based on the registration process. Coarse registration (also known as global registration) starts from an arbitrary initial pose and calculates the initial transformation parameters between two point clouds, thus providing reasonable initial values ​​for subsequent fine registration. Existing coarse registration methods mainly include registration based on geometric primitives such as points, lines, and surfaces, as well as 4PCS. Fine registration (also known as local registration), on the other hand, further optimizes the transformation parameters based on coarse registration to obtain high-precision registration results. Representative methods include ICP, NDT, and their various variants.

[0004] However, in various geological disaster scenarios, the location of stable terrain zones is unknown beforehand, making point cloud registration algorithms that ignore surface deformation prone to getting trapped in local optima. Therefore, accurate identification of stable terrain zones is crucial. The currently common strategy is to first divide the source and target point clouds into several voxels, and then identify stable zones based on indicators such as the difference in information entropy, centroid distance, and point-to-surface distance between corresponding voxels. However, these methods all ignore the local surface deformation characteristics within voxels when dealing with non-uniform surface deformation, easily misclassifying some deformed areas as stable zones.

[0005] In summary, given the challenge that existing research struggles to capture stable surface regions, leading to suboptimal registration accuracy, there is an urgent need to propose a new point cloud registration method for topographic change zones. Summary of the Invention

[0006] The purpose of this invention is to provide a point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization, which can accurately capture stable areas of the land surface to achieve precise registration.

[0007] To achieve the above objectives, the present invention adopts the following technical solution: A point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization includes the following steps: I. Perform initial registration and voxel segmentation on the source point cloud and target point cloud ground point cloud data; II. Calculate the spatial geometric features of each voxel, and calculate the feature similarity between voxels in the source point cloud and voxels in the corresponding target point cloud, and estimate the similarity threshold by combining the Gaussian Mixture Model (GMM). Then, based on the similarity threshold, all voxels were initially classified into deformable voxels and stable voxels; The deformable voxels are then recursively subdivided to screen for stable voxels. Finally, all labeled stable voxels are defined as initial stable voxels; III. Based on the initial stable voxels extracted in step II, first perform plane fitting on each initial stable voxel in the target point cloud to obtain the voxel fitting plane for each initial stable voxel in the target point cloud. Then, the distance from each point in each initial stable voxel in the source point cloud to the voxel fitting plane in the target point cloud determined by the corresponding initial stable voxel is calculated, and an adaptive dynamic threshold distance constraint strategy is used to further filter the stable region. IV. Based on the stable region point cloud extracted in step III, perform fine registration from the source point cloud to the target point cloud.

[0008] Furthermore, based on the aforementioned point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization, this invention also proposes a corresponding point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization. Both are based on the same inventive concept and adopt the following technical solutions: A point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization includes the following modules: The preprocessing module is used to perform initial registration and voxel segmentation on the source point cloud and target point cloud ground point cloud data. The initial stable region screening module is used to calculate the spatial geometric features of each voxel, calculate the feature similarity between voxels in the source point cloud and voxels in the corresponding target point cloud, and estimate the similarity threshold by combining the Gaussian mixture model (GMM). Then, based on the similarity threshold, all voxels were initially classified into deformable voxels and stable voxels; The deformable voxels are then recursively subdivided to screen for stable voxels. Finally, all labeled stable voxels are defined as initial stable voxels; The stable region selection module is used to select the initial stable voxels extracted by the initial stable region selection module. First, it performs plane fitting on each initial stable voxel in the target point cloud to obtain the voxel fitting plane for each initial stable voxel in the target point cloud. Then, the distance from each point in each initial stable voxel in the source point cloud to the voxel fitting plane in the target point cloud determined by the corresponding initial stable voxel is calculated, and an adaptive dynamic threshold distance constraint strategy is used to further filter the stable region. And a fine registration module, used to perform fine registration from source point cloud to target point cloud on the stable region point cloud extracted by the stable region selection module.

[0009] Compared with the prior art, the present invention has the following advantages: As described above, this invention discloses a point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization. This method constructs a point cloud spatial similarity measurement model, which quantifies the similarity of the spatial structure features of points within voxels in the source and target point clouds, and dynamically estimates the classification threshold using a Gaussian mixture model to achieve preliminary voxel category discrimination. The method proposes a deformable voxel recursive subdivision algorithm, which accurately captures the boundaries between terrain deformation areas and stable areas through adaptive recursive subdivision of deformable voxels, improving the accuracy of local stable area boundary recognition. The method designs an adaptive dynamic threshold distance constraint strategy, which accurately extracts stable terrain areas through dynamically decaying multi-threshold optimization, further refining the range of stable areas and ensuring their reliability. Attached Figure Description

[0010] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the accompanying drawings used in the embodiments will be briefly described below.

[0011] Figure 1 This is a flowchart illustrating the point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to the present invention. Figure 2 This is a schematic diagram of the point cloud spatial similarity measurement model in an embodiment of the present invention; Figure 3 This is a schematic diagram of the corresponding voxel distance constraint in an embodiment of the present invention; Figure 4 These are point cloud experimental datasets for four different geological disaster areas in this embodiment of the invention; where (a), (b), (c) and (d) are Data1, Data2, Data3 and Data4, respectively. Figure 5 This is a schematic diagram comparing the accuracy of the method of the present invention with other registration methods on the experimental dataset in an embodiment of the present invention; where (a), (b), (c) and (d) are the accuracy comparisons of Data1, Data2, Data3 and Data4, respectively; Figure 6This is a global schematic diagram of the registration results of each registration method on Data1 in the embodiments of the present invention; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparison methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively; Figure 7 for Figure 6 Enlarged views of the red rectangular areas in the registration results of each method; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparative methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively; Figure 8 The diagram shows the registration results of each registration method on Data2 in the embodiments of the present invention; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparison methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively. Figure 9 for Figure 8 Enlarged views of the red rectangular areas in the registration results of each method; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparative methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively; Figure 10 The diagram shows the registration results of each registration method in the embodiments of the present invention on Data3; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparison methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively. Figure 11 for Figure 10 Enlarged views of the red rectangular areas in the registration results of each method; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparative methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively; Figure 12 The diagram shows the registration results of each registration method in the embodiments of the present invention on Data4; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparison methods ICP, tr-ICP, P2P-ICP, GICP, ICPux, VM and SV-M, respectively. Figure 13 for Figure 12 Enlarged views of the red rectangular areas in the registration results of each method; where (a) is the reference registration result; (b) is the registration result of the SMTO method of the present invention; (c), (d), (e), (f), (g), (h) and (i) are the registration results of the comparative methods ICP, tr-ICP, P2P-ICP, GICP, ICProx, VM and SV-M, respectively. Detailed Implementation

[0012] 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 a part of the embodiments of the present invention, and not all of them. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.

[0013] Example 1 This embodiment 1 describes a point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization, such as... Figure 1 As shown, it includes the following steps: I. Data Preprocessing.

[0014] First, initial registration is performed on the source and target point cloud data. To avoid getting stuck in local optima during fine-tuning, this embodiment uses the 4PCS algorithm for initial registration of the two point clouds. Specifically: First, in the source point cloud, a set of four coplanar points is selected, and the intersection point of the lines between the four points is calculated and denoted as e. The lengths of the two line segments are d1 and d2, and the ratios r1 and r2 of the two line segments divided by the intersection point are calculated for each line. In the target point cloud, all point pairs with distances of approximately d1 and d2 are searched, and the ratios are denoted as R1 and R2, respectively. Based on r1, the potential intersection points of the point pairs in R1 are calculated, and a range tree structure is constructed using the intersection points. Based on r2, the potential intersection points of the midpoint pairs in R2 are calculated, and for each intersection point, the nearest intersection point is found in the range tree structure, thus determining a coplanar four-point set. A transformation matrix is ​​calculated for each coplanar four-point set and applied to the source point cloud. Then, the number of points whose distance from the midpoint of the source point cloud to the midpoint of the target point cloud is less than a threshold is calculated; the optimal transformation is the one that brings the most points closer. Furthermore, to improve computational efficiency, the initially registered source and target point clouds are voxelized, with the voxel size determined based on experimental data.

[0015] II. Initial stable region screening.

[0016] Surface deformation typically includes rigid and non-rigid deformation. However, existing research often focuses on rigid deformation identification while neglecting the differences in local spatial features within voxels, thus limiting the accuracy of stable region identification. To address this, this embodiment constructs a point cloud spatial similarity metric model to initially screen stable regions, such as... Figure 2 As shown, firstly, the spatial features of each voxel are described, and the feature similarity between corresponding voxels of the source point cloud and the target point cloud is calculated; then, the feature similarity threshold is adaptively estimated by combining a Gaussian Mixture Model (GMM). The voxels are then classified; finally, deformable voxels are recursively subdivided to achieve preliminary screening of stable regions. The specific steps are as follows: II.1. Determination of stable voxels.

[0017] Due to factors such as equipment performance, observation angle, and environmental changes, multi-temporal point clouds often exhibit significant differences in data density and local distribution characteristics. Therefore, this embodiment uses the distance from the in-voxel point to the centroid as the weight to establish a weighted covariance matrix, and calculates its eigenvalues ​​λ1, λ2, and λ3, as well as the corresponding eigenvectors n1, n2, and n3. Based on this, a spatial feature descriptor F is designed, which is based on eigenvalues ​​and eigenvectors. The formula for F is: In the formula, Linearity measures whether the point cloud within a voxel exhibits a distinct linear distribution. For flatness, characterize the local flatness of voxels; Discreteness reflects the degree of dispersion of the spatial distribution of point clouds within a voxel; The variance is the total variance, representing the degree of dispersion of the point cloud within a voxel in each direction; It is anisotropic, reflecting the non-uniformity of point clouds within voxels in different directions; Feature entropy describes the uncertainty of local structure within a voxel; For traces, the strength of the geometric structure within a voxel is quantified; Verticality represents the degree of deviation of the local normal vector of a voxel from the vertical direction; and These represent the maximum and minimum curvatures, respectively, characterizing the maximum and minimum degree of bending of the voxel's inner surface shape.

[0018] Then, cosine similarity is used. The feature similarity between voxels in the source point cloud and their corresponding target point cloud is measured by the following formula: In the formula, and These represent the spatial features of the corresponding voxels in the source and target point clouds, respectively.

[0019] II.2. Adaptive estimation of feature similarity threshold using GMM .

[0020] Gaussian Mixture Models (GMMs) are linear combinations of multiple Gaussian distribution functions used to address situations where data from a unified dataset has multiple different distributions. They are used for adaptive estimation of feature similarity thresholds. The specific steps are as follows: The core theorem of GMM is that the probability density function is... , Expressed as: In the formula, x is the cosine similarity. K represents the number of components in the mixture; in this embodiment, K is set to 2. The Gaussian mixture weights are the weights of the k-th Gaussian distribution; It is the k-th Gaussian distribution with a mean of The covariance matrix is By learning the above parameters, GMM can achieve accurate threshold estimation of spatial feature similarity.

[0021] Gaussian Mixture Models (GMMs) typically employ the expectation maximization (EM) algorithm, an iterative estimation algorithm designed to maximize the log-likelihood function of the observed data. The EM algorithm consists of an E-step and an M-step; the E-step involves calculating the log-likelihood function for each data point given a parameter θ. Posterior probability of belonging to each Gaussian class k The posterior probability is also called the responsibility weight. The formula is: In the formula, Represents data point x i The probability density function value under the k-th Gaussian component Represents data point x i The probability density function value under the j-th Gaussian component; The M-step, also known as the maximization step, is based on the posterior probability obtained from the E-step. Update GMM parameters.

[0022] Repeat the E-step and M-step until convergence is achieved, obtaining the optimal parameters and the Gaussian mixture model. , The formula is: In the formula, Let represent the k-th Gaussian distribution.

[0023] When the GMM is fitted with two components, one representing "stability" and the other representing "deformation", the similarity threshold is... Take the intersection of the probability density functions of two Gaussian distributions: The intersection point x value obtained by solving equation (6) is the optimal similarity threshold. .

[0024] according to Perform voxel similarity determination when When the voxels in the source point cloud and the corresponding target point cloud have similar spatial feature distributions, it means that the corresponding surface has not been deformed; otherwise, it has been deformed. Therefore, according to Perform similarity determination between voxels in the source point cloud and corresponding voxels in the target point cloud. For similarity... Voxels that are stable are identified and retained at the current level; otherwise, they are marked as deformable voxels.

[0025] II.3. Recursive Subdivision of Deformable Voxels. Voxel subdivision based on geometric features tends to remove stable regions mixed in with deformable voxels, resulting in the loss of a large number of point clouds that can be used for fine registration. Therefore, this embodiment further recursively subdivides the deformable voxels obtained in step II.2 to filter stable regions, such as... Figure 2 As shown. Specifically, firstly, the deformed voxel octree obtained in step II.2 is subdivided; then, for each sub-voxel, feature similarity estimation and stable voxel identification are performed again according to steps II.1 to II.2. Steps II.1 and II.2 are repeated until the number of point clouds within the deformed voxel is less than a preset threshold or a predetermined number of layers is reached; finally, all labeled stable voxels are defined as initial stable voxels.

[0026] III. Refinement of the stable region. Due to factors such as voxel size and calculation deviations in voxel geometry, the initial stable voxels inevitably contain non-terrain stable points. To further eliminate rigid deformable voxels contained within the initial stable voxels, this embodiment designs an adaptive dynamic threshold distance constraint strategy to further screen stable regions. This strategy sets the number of iterations and simultaneously adjusts the distance threshold t. s Perform iterative updates and obtain the optimal distance threshold t through iterative updates. s Each iteration involves three key steps: defining the objective function, selecting the stable region, and updating the threshold. Specifically: III.1. Define the objective function.

[0027] To minimize the registration error between stable points in the source and target point clouds, this embodiment uses the root mean square error (RMSE) of the distance between corresponding points as the optimization objective function: In the formula, n represents the number of corresponding points that can be successfully matched. and Let R and t represent the coordinates of the k-th corresponding point between the source and target point clouds, respectively. R is the rotation matrix, and t is the translation vector. The position update of the source point cloud can be achieved by using the optimal transformation parameters R and t.

[0028] III.2. Screening for stable regions.

[0029] like Figure 3 As shown, this embodiment employs a method based on planar geometric consistency to further screen stable regions. For example... Figure 1As shown in the refinement section of the intermediate stability region, based on the initial stable voxels extracted in step II, firstly, plane fitting is performed on each initial stable voxel in the target point cloud to obtain the voxel fitting plane for each initial stable voxel in the target point cloud; then, the distance d from each point in each initial stable voxel in the source point cloud to the voxel fitting plane determined by the corresponding initial stable voxel in the target point cloud is calculated. i Then, the optimal distance threshold t is obtained according to the adaptive dynamic threshold distance constraint strategy. s If the distance d i Less than or equal to the optimal distance threshold t s If the value is positive, then mark that point as a stable point. Finally, calculate the proportion of stable points within each initial stable voxel. ,in, This indicates that the initial stable voxel in the source point cloud satisfies d i Less than or equal to the optimal distance threshold t s The total number of points, m represents the total number of points within the initial stable voxel in the source point cloud. If this ratio is not less than a set ratio η, the initial stable voxel is determined to be a stable voxel. In this embodiment, η is set to 50%. Compared to existing methods that rely solely on the geometric distance between centroids or boundary points, this strategy introduces geometric consistency constraints from points to the local fitting plane, which can effectively improve the reliability of stable region identification.

[0030] III.3. Threshold Update.

[0031] During each iteration, the distance threshold t from each point within each initially stable voxel in the source point cloud to the corresponding initially stable voxel fitting plane in the target point cloud is... s The threshold needs to be updated gradually as the computational accuracy improves. This is because the accuracy of stable region extraction largely depends on the threshold setting. While a smaller threshold helps eliminate outliers, it may also miss some stable regions; conversely, a larger threshold easily introduces unstable points, affecting subsequent registration accuracy. To address this, this embodiment introduces an adaptive dynamic threshold update strategy based on AdaDelta. This method uses changes in the overlap ratio of stable regions to drive the distance threshold t. s The adaptive adjustment gradually optimizes the transformation parameters to minimize the objective function (Equation (3)).

[0032] Specifically, firstly, the initial distance threshold t s0 The distance is set to the average distance from each point within each initially stable voxel in the source point cloud to the fitted plane of the corresponding initially stable voxel in the target point cloud: In the formula, This represents the coordinates of all points in the source point cloud voxels. This represents the parameters after fitting the target point cloud to the voxel plane, where N is the number of stable voxels.

[0033] Considering that the changing trend of the stable proportion in adjacent iterations can more reasonably adjust the distance threshold, this embodiment uses the difference in the stable proportion between two adjacent iterations to construct the descent gradient; at the same time, to avoid the distance threshold t s Dramatic changes occur during the iteration process, so a decay coefficient is introduced, and the distance threshold adjustment is calculated based on this coefficient. Subsequently, gradient descent is used to update the distance threshold t. s The formula is as follows: In the formula, λ is the proportionality coefficient, which is used to enhance the influence of the stable proportion on the threshold update; γ represents the stable voxel ratio in each iteration; γ is the decay coefficient; and ε is a minimal constant term to prevent the denominator from being zero. During this process, the fixed parameters are set as follows: , , .

[0034] Update the distance threshold t in each iteration s The corresponding stable region results were obtained. The root mean square error (RMSE) of the source and target point clouds was calculated based on the obtained stable region results for registration, as shown in Equation 7. In this embodiment, P2P-ICP is used to register the source and target point clouds of the obtained stable region results. The minimum RMS error corresponds to the distance threshold t. s That is, the optimal distance threshold t s And using this optimal distance threshold t s The refined results of the stable region were obtained.

[0035] IV. Point cloud precision registration.

[0036] Traditional ICP algorithms employ a point-to-point nearest neighbor search strategy, using the minimization of the distance between corresponding points in the source and target point clouds as the registration criterion. This approach is susceptible to interference from terrain undulations and changes in point cloud density, making accurate matching of point pairs difficult. In contrast, P2P-ICP uses the minimization of the distance between a point in the source point cloud and the plane containing the corresponding point in the target point cloud as the registration criterion. This effectively reduces the impact of multi-temporal point cloud density differences, better resists erroneous point pair matching, and improves iterative convergence efficiency. Therefore, this embodiment uses the P2P-ICP algorithm to perform fine registration of the source and target point cloud data from the stable region point cloud extracted in step III. The solved transformation parameters are then applied to the entire unregistered point cloud to achieve coordinate transformation from the source to the target point cloud for subsequent deformation analysis.

[0037] The effectiveness and robustness of the proposed Point Cloud Registration Method for TerrainChange Areas Based on Spatial Similarity And Multi-Threshold Optimization (SMTO) will be verified below in different terrain and disaster scenarios.

[0038] like Figure 4 As shown in Table 1, this embodiment selected four representative LiDAR point cloud datasets, including rockfall areas, glacier areas, landslide areas, and glacier front areas. Each dataset contains two phases of point cloud data, with collection time spans ranging from several months to several years, and each area exhibits significant topographic variations. The two phases of point cloud data for each area were rigorously registered manually, and the registration parameters were used as ground truth values ​​to verify the accuracy of each registration method.

[0039] Table 1. Statistical information of the experimental dataset in this invention.

[0040] To evaluate the performance of SMTO, it was compared with seven classic registration methods, including fine registration algorithms that do not consider the deformation region, such as ICP, tr-ICP, P2P-ICP, and GICP, and registration algorithms that take the deformation region into account, such as ICPux, VM, and SV-M. Rotation error (eR) and translation error (et) were introduced as quantitative indicators for accuracy evaluation, and the specific calculation formulas are as follows: In the formula, and These represent the calculation of the rotation matrix and translation vector, respectively. and These represent the reference transformation matrix and the translation vector, respectively. It is the trace of the matrix.

[0041] like Figure 5 As shown, the registration accuracy of each method varies across different datasets, but the overall trend is relatively consistent. Specifically, all methods achieve the highest registration accuracy on Data1 because this dataset has a complete point cloud structure, low deformation, and provides a higher proportion of stable regions, thus yielding higher accuracy registration results. In contrast, all methods have the largest errors on Data3, mainly due to significant deformation in this dataset. The results for Data2 and Data4 fall between the two, with Data2 showing a slightly lower eR but a larger et, while Data4 is generally more balanced. This indicates that the degree of deformation between multi-temporal point clouds is the main factor affecting registration accuracy. When the deformation is small, all methods can obtain relatively stable registration results; however, under severe deformation, the performance of each method is affected to varying degrees.

[0042] A comparison of the performance of various methods reveals that traditional ICP and its improved algorithms (tr-ICP, P2P-ICP, GICP) have lower accuracy. This is mainly because these methods do not consider surface deformation, making them prone to getting trapped in local optima under severe surface deformation conditions, as seen in Data2 and Data3. The ICProx method slightly outperforms GICP in Data2 because it introduces a stable region filtering mechanism; however, its reliance on point cloud distance variations results in insufficient stability. VM and SV-M show relatively stable performance across the datasets, with SV-M outperforming VM in translation accuracy, indicating that introducing stepwise constraints on stable regions can effectively reduce local mismatches. In contrast, SMTO consistently maintains the best performance across all four datasets, especially in Data3 and Data4 scenarios, where its translation error et is significantly lower than the comparative methods, validating its comprehensive advantages in both rotation and translation. Table 2 shows the average computational accuracy of SMTO compared with other registration methods. On average, SMTO has an eR of 0.035 rad, which is about 2.8% to 31.4% lower than other methods; and an et of 2.025 cm, which is 64.4% to 92.0% lower than other methods.

[0043] Table 2. Average computational accuracy of SMTO compared to other registration methods

[0044] like Figure 6 , Figure 8 , Figure 10 and Figure 12As shown, overall, the performance of different methods varies significantly across different datasets, but SMTO maintains high accuracy and stability across all datasets. Specifically, as... Figure 6 and Figure 7 As shown, in Data1, the overall registration results of all methods are relatively good, with SMTO achieving the highest accuracy and its results being basically consistent with the reference registration results. ICProx and VM are next, while the SV-M method, although considering the stability region constraint, still exhibits deviations at the boundaries. Furthermore, the four registration methods that do not consider deformation show more severe deviations. For example... Figure 8 and Figure 9 As shown in Data2, due to the instability of most areas within the study region caused by glacial melting, the performance differences among the registration methods are significant. ICP, tr-ICP, P2P-ICP, and GICP failed to identify stable regions, ultimately introducing severe registration errors and resulting in large deviations from the reference results. In comparison, while ICPox, VM, and SV-M achieved relatively better registration results, their overall errors remained large. In contrast, SMTO registration results showed high consistency with the reference results. Figure 10 and Figure 11 As shown, the Data3 scenario exhibits severe landslide deformation, but SMTO, compared to other methods, can accurately distinguish between the deformation zone and the stable zone, and completes registration only based on the stable zone. The final registration result is essentially consistent with the reference result. Figure 12 and Figure 13 As shown, in Data4, SMTO maintains good registration performance, with results closely matching the reference registration. Overall, SMTO demonstrates stable and excellent registration results across all four datasets, and its comprehensive performance significantly outperforms other classic methods.

[0045] To further verify the contribution of SMTO's innovations to point cloud registration performance, this embodiment conducted ablation experiments on four datasets, including: ① establishing a point cloud spatial similarity metric model versus direct distance constraints; ② recursive subdivision of deformable voxels versus no recursive subdivision; ③ dynamic thresholding with layer-by-layer optimization versus fixed thresholding. As shown in Table 3, the results show that SMTO's improvement in point cloud registration performance in terrain-changing areas is mainly reflected in the translation error et. Specifically, compared with the method of direct distance constraints, establishing a point cloud spatial similarity metric model reduces the translation error et by 61.3%, 7.7%, 50.9%, and 56.2% in the four datasets, respectively. This indicates that the model helps to identify and eliminate deformable voxels, reduce mismatches, and thus improve the reliability of point cloud registration in terrain-changing areas. Compared with the method without recursive subdivision, SMTO reduces the translation error et by 70.9%, 13.5%, 50.8%, and 45.6% in the four datasets, respectively. This indicates that the strategy can accurately capture the boundaries of stable areas through the subdivision of deformable voxels, thereby improving matching accuracy. Compared with the registration method with a fixed threshold, the dynamic threshold optimized layer by layer reduced the translation error et by 41.7%, 18.8%, 70.2% and 75.4%, respectively. This shows that the mechanism can effectively alleviate the interference of local rigid deformation on the registration process, thereby enhancing the convergence stability and robustness of point cloud registration in areas with different terrain changes.

[0046] Table 3 Comparison of the accuracy of the method of the present invention in ablation experiments

[0047] To address the problem of low point cloud registration accuracy caused by inaccurate identification of stable terrain zones in traditional methods, this invention proposes a point cloud registration method for terrain change zones driven by spatial similarity and multi-threshold optimization. This method first uses a point cloud spatial similarity metric model to initially identify stable and deformable voxels. Then, it recursively subdivides deformable voxels and couples a layer-by-layer optimized threshold dynamic update strategy to further subdivide and identify deformable voxels, thereby improving the accuracy of stable voxel identification. Validation using datasets from four different geological disaster scenarios shows that the proposed method performs best in terms of average rotation error eR and average translation error et, and the overall registration effect is closer to the actual registration effect. This indicates that the proposed method has high stable zone identification capability and good point cloud registration accuracy. In summary, the proposed method effectively solves the problem of high-precision point cloud registration in geological disaster scenarios where the location of deformable zones is unknown, providing reliable technical support for geological disaster monitoring.

[0048] Example 2 This embodiment 2 describes a point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization. This point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization is based on the same inventive concept as the point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization in embodiment 1 above.

[0049] The point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization in this embodiment includes the following modules: The preprocessing module is used to perform initial registration and voxel segmentation on the source point cloud and target point cloud ground point cloud data.

[0050] The initial stable region screening module is used to calculate the spatial geometric features of each voxel, calculate the feature similarity between voxels in the source point cloud and voxels in the corresponding target point cloud, and estimate the similarity threshold by combining the Gaussian mixture model (GMM). Then, based on the similarity threshold, all voxels are initially classified into deformable voxels and stable voxels. The deformable voxels are then recursively subdivided to screen the stable voxels. Finally, all the marked stable voxels are defined as initial stable voxels.

[0051] The stable region selection module is used to select the initial stable voxels extracted by the initial stable region selection module. First, it performs plane fitting on each initial stable voxel in the target point cloud to obtain the voxel fitting plane of each initial stable voxel in the target point cloud. Then, it calculates the distance from each point in each initial stable voxel in the source point cloud to the voxel fitting plane determined by the corresponding initial stable voxel in the target point cloud, and uses an adaptive dynamic threshold distance constraint strategy to further select stable regions.

[0052] And a fine registration module, used to perform fine registration from source point cloud to target point cloud on the stable region point cloud extracted by the stable region selection module.

[0053] It should be noted that any content not mentioned in the above-mentioned functional modules of the spatial similarity and multi-threshold optimization-driven point cloud registration system for terrain change areas described in this embodiment can be referred to the step description of the corresponding method in the above embodiment 1, and will not be elaborated here.

[0054] The present embodiment has now been described in detail with reference to the accompanying drawings. Based on the above description, those skilled in the art should have a clear understanding of the present invention. Of course, the above description is only a preferred embodiment of the present invention and does not limit the patent scope of the present invention. Any equivalent structural transformations made using the content of the present invention's specification and drawings under the inventive concept of the present invention, or direct / indirect applications in other related technical fields, are included within the patent protection scope of the present invention and should be protected by the present invention.

Claims

1. A point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization, characterized in that, Includes the following steps: I. Perform initial registration and voxel segmentation on the source point cloud and target point cloud ground point cloud data; II. Calculate the spatial geometric features of each voxel, and calculate the feature similarity between voxels in the source point cloud and voxels in the corresponding target point cloud, and estimate the similarity threshold by combining the Gaussian Mixture Model (GMM). Then, based on the similarity threshold, all voxels were initially classified into deformable voxels and stable voxels; The deformable voxels are then recursively subdivided to screen for stable voxels. Finally, all labeled stable voxels are defined as initial stable voxels; III. Based on the initial stable voxels extracted in step II, first perform plane fitting on each initial stable voxel in the target point cloud to obtain the voxel fitting plane for each initial stable voxel in the target point cloud. Then, the distance from each point in each initial stable voxel in the source point cloud to the voxel fitting plane in the target point cloud determined by the corresponding initial stable voxel is calculated, and an adaptive dynamic threshold distance constraint strategy is used to further filter the stable region. IV. Based on the stable region point cloud extracted in step III, perform fine registration from the source point cloud to the target point cloud.

2. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 1, characterized in that, Step II specifically involves: II.

1. Calculate the feature similarity between voxels in the source point cloud and the corresponding voxels in the target point cloud; First, a weighted covariance matrix is ​​established using the distance from the inlier point to the centroid as the weight, and its eigenvalues ​​λ1, λ2, and λ3 are calculated, along with the corresponding eigenvectors n1, n2, and n3. Based on this, a spatial feature descriptor F is designed, which is based on eigenvalues ​​and eigenvectors. The formula for F is: In the formula, Linearity measures whether the point cloud within a voxel exhibits a distinct linear distribution. For flatness, characterize the local flatness of voxels; Discreteness reflects the degree of dispersion of the spatial distribution of point clouds within a voxel; The variance is the total variance, representing the degree of dispersion of the point cloud within a voxel in each direction; It is anisotropic, reflecting the non-uniformity of point clouds within voxels in different directions; Feature entropy describes the uncertainty of local structure within a voxel; For traces, the strength of the geometric structure within a voxel is quantified; Verticality represents the degree of deviation of the local normal vector of a voxel from the vertical direction; and These are the maximum and minimum curvatures, representing the maximum and minimum degree of bending of the voxel's inner surface shape; Then, cosine similarity is used. The feature similarity between voxels in the source point cloud and their corresponding target point cloud is measured by the following formula: In the formula, and These represent the spatial features of the corresponding voxels in the source and target point clouds, respectively. II.

2. Adaptive estimation of feature similarity threshold using GMM ,when When the time is right, it means that the corresponding surface has not deformed; otherwise, it has deformed. Thus, all voxels are divided into two categories: deformable voxels and stable voxels. II.

3. Further recursively subdivide the deformable voxels obtained in step II.2; first, subdivide the deformable voxel octree into sub-voxels; then, perform feature similarity estimation and stable voxel identification on each sub-voxel again according to steps II.1 to II.2; Repeat steps II.1 and II.2 until the number of point clouds within the deformed voxels is less than a preset threshold or a predetermined number of layers is reached; finally, define all marked stable voxels as initial stable voxels.

3. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 2, characterized in that, In step II.2, the feature similarity threshold is adaptively estimated using GMM. The steps are as follows: First, use the probability density function Describe the data distribution. Expressed as: In the formula, x is the cosine similarity. K represents the number of components in the mixture. The Gaussian mixture weights are the weights of the k-th Gaussian distribution; It is the k-th Gaussian distribution with a mean of The covariance matrix is ; Then, the expectation-maximization algorithm is used for iterative estimation; the expectation-maximization algorithm includes an E-step and an M-step. Step E is as follows: Given the parameter θ, calculate x for each data point. i Posterior probability of belonging to each Gaussian class k , The formula is: In the formula, Represents data point x i The probability density function value under the k-th Gaussian component Represents data point x i The probability density function value under the j-th Gaussian component; The M-step then uses the posterior probability obtained from the E-step. Update GMM parameters; Repeat the E-step and M-step until convergence is achieved, obtaining the optimal parameters and the Gaussian mixture model. , The formula is: In the formula, This represents the k-th Gaussian distribution; When the GMM is fitted to two components, the similarity threshold is... Take the intersection of the probability density functions of two Gaussian distributions: The intersection point x value obtained by solving this equation is the optimal similarity threshold. .

4. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 1, characterized in that, Step III specifically involves: Based on the initial stable voxels extracted in step II, the points contained within each initial stable voxel in the target point cloud are first fitted to a plane. The distance d from each point within each initial stable voxel in the source point cloud to the fitted plane of the corresponding initial stable voxel in the target point cloud is then calculated. i Then, the optimal distance threshold t is obtained according to the adaptive dynamic threshold distance constraint strategy. s ; If the distance d i Less than or equal to the optimal distance threshold t s If the initial stable voxel is not less than a set ratio, then the initial stable voxel is marked as a stable voxel. Finally, the proportion of stable voxels within each initial stable voxel is calculated. If the proportion is not less than a set ratio, then the initial stable voxel is determined to be a stable voxel.

5. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 4, characterized in that, The adaptive dynamic threshold distance constraint strategy is as follows: First, define the objective function; Root mean square error is chosen as the objective function for optimization. In the formula, n represents the number of corresponding points that can be successfully matched. and Let R and t represent the coordinates of the k-th corresponding point between the source point cloud and the target point cloud, respectively. R is the rotation matrix and t is the translation vector. Then, the initial distance threshold t s0 The distance is set to the average distance from each point within each initially stable voxel in the source point cloud to the fitted plane of the corresponding initially stable voxel in the target point cloud: In the formula, This represents the coordinates of all points in the source point cloud voxels. This represents the parameters after fitting the plane to the voxel corresponding to the target point cloud, where N is the number of stable voxels. Finally, set the number of iterations for the distance threshold t. s Perform iterative updates; Each distance threshold t s Each iterative update yields a corresponding stable region result. The root mean square error (RMSE) of the source and target point clouds is then calculated based on the obtained stable region results, with the minimum RMS error corresponding to the distance threshold t. s As the optimal distance threshold t s And using this optimal distance threshold t s The refined results of the stable region were obtained.

6. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 5, characterized in that, The distance threshold t s The iterative update is specifically as follows: The descent gradient is constructed using the difference in the stable proportions between two adjacent iterations; a decay coefficient is introduced, and the distance threshold t is calculated based on this coefficient. s Adjust the amount; then update the distance threshold t using gradient descent. s : In the formula, λ is the scaling factor, used to enhance the stability of the scaling factor relative to the distance threshold t. s The impact of the update; γ is the proportion of stable voxels in each iteration; γ is the decay coefficient; ε is a constant term to prevent the denominator from being zero.

7. The point cloud registration method for terrain change areas driven by spatial similarity and multi-threshold optimization according to claim 1, characterized in that, Step IV specifically involves: Based on the stable area point cloud extracted in step III, the fine registration of the two phases of ground point cloud data is completed, and the solved transformation parameters are applied to the entire unregistered point cloud to realize the coordinate transformation from the source point cloud to the target point cloud.

8. A point cloud registration system for terrain change areas driven by spatial similarity and multi-threshold optimization, characterized in that, Includes the following modules: The preprocessing module is used to perform initial registration and voxel segmentation on the source point cloud and target point cloud ground point cloud data. The initial stable region screening module is used to calculate the spatial geometric features of each voxel, calculate the feature similarity between voxels in the source point cloud and voxels in the corresponding target point cloud, and estimate the similarity threshold by combining the Gaussian mixture model (GMM). Then, based on the similarity threshold, all voxels were initially classified into deformable voxels and stable voxels; The deformable voxels are then recursively subdivided to screen for stable voxels. Finally, all labeled stable voxels are defined as initial stable voxels; The stable region selection module is used to select the initial stable voxels extracted by the initial stable region selection module. First, it performs plane fitting on each initial stable voxel in the target point cloud to obtain the voxel fitting plane for each initial stable voxel in the target point cloud. Then, the distance from each point in each initial stable voxel in the source point cloud to the voxel fitting plane in the target point cloud determined by the corresponding initial stable voxel is calculated, and an adaptive dynamic threshold distance constraint strategy is used to further filter the stable region. And a fine registration module, used to perform fine registration from source point cloud to target point cloud on the stable region point cloud extracted by the stable region selection module.

Images