A low-overlap-rate point cloud registration method based on topology-metric decoupling

Abstract

A low overlap rate point cloud registration method based on topology-metric decoupling belongs to the field of computer vision and three-dimensional reconstruction technology, the three-dimensional laser scanning acquisition terminal scans the article and collects two point clouds at different angles as source point cloud and target point cloud, and realizes high robustness registration in the difficult scene with extremely low overlap rate; the present application solves the zero solution problem under low overlap rate by one-way matching and set through the asymmetric association strategy, which significantly improves the recall rate; at the same time, through the DGTG geometric pruning and probability manifold optimization mechanism, the anisotropic covariance and Lie algebra are used to eliminate the voxel quantization error, which greatly reduces the translation error and rotation error, and effectively breaks through the precision bottleneck caused by voxelization. The present application not only can guarantee the existence of solution under the condition of extremely low overlap rate, but also can eliminate the voxel quantization error through the manifold optimization based on covariance, and realize the high-precision and high-efficiency registration of low overlap rate point cloud.

Description

Technical Field

[0001] This invention belongs to the field of computer vision and 3D reconstruction technology, and in particular relates to a low overlap rate point cloud registration method based on topology-metric decoupling. Background Technology

[0002] Currently, well-known 3D point cloud registration techniques mainly employ Random Sample Consensus (RANSAC) and its variants, or feature matching methods based on deep learning. Traditional geometric methods typically use a "reciprocal consistency check" strategy to filter correspondences, requiring points in the source and target point clouds to be nearest neighbors to be considered valid. However, in scenarios with an overlap rate below 30%, due to the physical truncation of the sensor's field of view, over 95% of the true corresponding points only satisfy a one-way nearest neighbor relationship. The reciprocal check incorrectly discards these one-way valid geometric cues, causing the algorithm to fall into a "zero solution" dilemma in low-overlap areas because it cannot find enough corresponding points. Furthermore, to balance computational efficiency, existing techniques generally perform voxel downsampling on point clouds. This coordinate discretization inevitably introduces centimeter-level quantization errors, causing the final registration accuracy to be locked at the voxel resolution level, exhibiting a "flat-top accuracy effect," which cannot meet the high-precision alignment requirements at the sub-voxel level. While deep learning-based methods have improved feature extraction, their inference process heavily relies on high-performance GPUs and large amounts of training data, resulting in high power consumption and poor cross-scenario generalization ability, making them difficult to apply on edge computing devices with only low-power CPUs.

[0003] Therefore, there is an urgent need for a new technical solution to address this problem. Summary of the Invention

[0004] A low-overlap point cloud registration method based on topology-metric decoupling includes the following steps, which are performed sequentially:

[0005] Step 1: Use a 3D laser scanning acquisition terminal to scan and obtain point cloud data of the scanned object in two selected different directions, and record them as the source point cloud and the target point cloud respectively;

[0006] Step 2: Perform coordinate mapping and matrix transformation on the source point cloud and the target point cloud to obtain the source matrix and the target matrix respectively;

[0007] Mesh denoising is performed on the source matrix and the target matrix respectively to obtain two sets of sparse skeleton node matrices: the source skeleton coordinate matrix and the target skeleton coordinate matrix, and a KD-Tree index is generated.

[0008] Geometric encoding is performed on each skeleton node, resulting in two sets of region convolutional fusion FPFH feature matrices: the source feature matrix and the target feature matrix.

[0009] Step 3: Traverse each row of the source feature matrix in one direction to find the best matching row of the target feature matrix, form a matching pair, and then generate an initial rough correspondence set representing the index correspondence relationship, and establish a preliminary logical connection between the source point cloud and the target point cloud.

[0010] Step 4: Take two random matching pairs from the initial rough correspondence set, and use the maximum clique consistency principle of rigid constraints and DGTG dynamic geometric topology gating to remove outliers and obtain a correspondence set with consistent geometric structure, i.e., the maximum clique structure. Then, filter again and output the fine correspondence set.

[0011] Step 5: Perform a global attraction domain anchoring operation on the fine correspondence set to obtain the anchoring transformation matrix as the initial value for the next step of probability plane sliding;

[0012] Step 6: For the largest clique target skeleton node in the fine correspondence set, construct a virtual sliding plane. By minimizing the probability sliding loss function, drive the largest clique source node to slide along the tangent plane to the position that coincides with the target node, thereby achieving registration between the source point cloud and the target point cloud.

[0013] In step two, the grid noise reduction uses voxel filtering technology, and the specific method is as follows:

[0014] 1) Set the voxel side length For the source point cloud floating-point coordinate matrix and the target point floating-point coordinate matrix Spatial downsampling is performed separately to unify data density, where, Dimensions , Dimensions , The total number of points in the source point cloud. The total number of points in the target point cloud;

[0015] Set any point in the coordinate matrix Its spatial components in three-dimensional coordinates are Calculate the three-dimensional index vector of its voxels. :

[0016] ;

[0017] 2) Clustering: Retrievals falling within the same index of For each point, calculate its geometric centroid as the unique skeleton node. :

[0018] ;

[0019] In the formula, This indicates the first voxel in the source or target point cloud that falls within the current voxel grid. One point, ;

[0020] 3) Obtain the coordinate matrix of the source skeleton after dimensionality reduction. Dimension and the target skeleton coordinate matrix Dimension ,in, These represent the number of skeleton nodes in the source point cloud and the target point cloud after dimensionality reduction, respectively.

[0021] The method for generating the FPFH feature matrix in step two is as follows:

[0022] Step 1: Obtain the normal vector;

[0023] 1) Use KD-Tree to search for the first skeleton coordinate matrix. Nodes of Calculate the local centroid from each neighboring node. And construct the covariance matrix :

[0024] ;

[0025] In the formula, Represents a node The 1 neighboring node, ;

[0026] This means The geometric centroid of each neighboring node;

[0027] 2) Perform eigenvalue decomposition:

[0028] For covariance matrix Perform orthogonal decomposition and solve the equations Three eigenvalues ​​were obtained. and its corresponding eigenvectors ,and ;

[0029] in, This represents the principal direction with the largest variance in the data distribution, while This represents the direction with the smallest variance in the data distribution; geometrically, the direction perpendicular to the local plane has the smallest data variation, therefore, Define as a node normal vector ;

[0030] Step 2: Perform relative angle encoding (SPFH);

[0031] To eliminate the effects of coordinate system rotation, at the nodes Construct a local orthogonal coordinate system for the Dabu frame. ,in, , indicating that the normal vector is axis; Determined by the cross product of the line and the normal. axis; Determined by the right-hand rule axis;

[0032] The purpose of the Dabu frame is to: construct a dynamic local reference, ensuring rotational invariance of angular features; and calculate the angular deviation of neighboring normal vectors within this coordinate system. Generate a 33-dimensional histogram relative angle code. ;

[0033] Calculate the angular deviation of the neighbor normal vector in this coordinate system The specific calculation formula is as follows:

[0034] ;

[0035] ;

[0036] ;

[0037] In the formula, Represents the normal vector of the neighbor node and Projection angle deviation between axes; Represents the normal vector of the center node and , The angle of inclination between the vectors connecting two points; Represents the normal vector of the neighbor node exist Orientation rotation angle deviation on the plane;

[0038] The values ​​of the three angular deviations mentioned above are divided into 11 equally wide intervals for projection counting;

[0039] Step 3: Perform region convolution fusion FPFH, and concatenate the three 11-dimensional sub-histograms to construct the final 33-dimensional SPFH feature vector.

[0040] Perform distance-weighted convolution to transform the nodes The relative angle encoding SPFH and its neighboring nodes The relative angle encoding SPFH fusion is performed using the following formula:

[0041] ;

[0042] In the formula, Represents the merged nodes The final 33-dimensional feature vector; For nodes Relative angle encoding; For nodes with neighboring nodes Point-to-point distance, For nodes Relative angle encoding;

[0043] Step 4: Construct the feature matrix and its logical expression:

[0044] For the source skeleton coordinate matrix All nodes Each node is traversed through steps one through three to obtain the corresponding 33-dimensional FPFH feature vector;

[0045] This The vectors are stacked in row-matrix order in memory according to their point indices to construct the source feature matrix. ;

[0046] Similarly, for the target skeleton coordinate matrix By traversing and stacking, the target feature matrix is ​​constructed. ;

[0047] in, ;

[0048] Each row in the above matrix structure represents the geometric fingerprint of the corresponding node:

[0049] Dimensions 1 to 11: These mainly characterize whether a local surface is close to a planar structure. The higher the value in this range, the smaller the deviation of the normal direction in the neighborhood.

[0050] Dimensions 12 to 22: These mainly characterize whether a local region exhibits a cylindrical or rod-like structure. The higher the value in this range, the more significant the curvature change is in the direction of the line connecting the points.

[0051] Dimensions 23 to 33: mainly characterize whether a local region contains right angles or obvious edge structures. The higher the value in this range, the closer the angle between the neighborhood normals is to orthogonal.

[0052] At this point, the two sets of source point clouds and target point clouds, which only have coordinates, have been transformed into two sets of enhanced data structures that have both skeleton coordinates and geometric features.

[0053] The specific method for establishing the preliminary logical connection between the source point cloud and the target point cloud in step three is as follows:

[0054] 1) The source feature matrix Dimension Target feature matrix Dimension ; Using the source skeleton coordinate matrix With the target skeleton coordinate matrix As an auxiliary index, it is used to subsequently restore the feature index to the physical coordinate index;

[0055] 2) High-dimensional index construction:

[0056] In order to solve the problem To address the inefficiency of linear search within a 33-dimensional vector, a high-dimensional index is constructed, as follows:

[0057] target feature matrix Considered a For a high-dimensional dataset of data points, an index structure for a 33-dimensional feature space is constructed using the FLANN fast nearest neighbor search algorithm or a high-dimensional KD-Tree. This is used to divide the 33-dimensional feature space into several hyperplane regions, so that similar geometric fingerprints are adjacent to each other in memory address, thereby accelerating the retrieval.

[0058] 3) Asymmetric nearest neighbor search

[0059] Traversing the source feature matrix For each row of feature vectors in the dataset, use it as a query vector and execute a one-way matching strategy. The search algorithm finds the nearest target feature vector and initially determines that the two feature vectors are identical, belonging to the same node or source point. The distance metric used is... Norm, to obtain the best matching target index found in the search. Establish a preliminary logical connection between the source point cloud and the target point cloud;

[0060] ;

[0061] In the formula, The source characteristic matrix is ​​the first The eigenvectors of the row, ; For the target feature matrix and The nearest row feature vector, ;

[0062] 4) Result mapping and set generation:

[0063] For each source index traversed Record the best target index found in the search. This forms an index-correspondence tuple, or simply an index pair. ;

[0064] Then put all The pairs converge to form the initial rough correspondence set. : ;

[0065] At this time It is An integer matrix, where each row represents a pair of predicted origin points, and the total number of rows in the set is . , .

[0066] The specific steps for step four are as follows:

[0067] 1) Geometric pruning:

[0068] Through the initial rough correspondence set Any two points in and As a matching pair, the coordinates are obtained from the source skeleton matrix via indexing. Read source node and The coordinates, from the target skeleton coordinate matrix Read target node and The coordinates of the source node are calculated. and distance Calculate and obtain the target node and distance ;

[0069] The absolute value of the difference between the lengths of the two segments is calculated using the geometric residual formula. :

[0070] ;

[0071] Set the tolerance threshold for geometric consistency determination , Then determine and These two matching pairs are geometrically compatible, preserving logical connections; If it is not a match, it is considered an incorrect match and is cut off. and The logical connection between these two matching pairs completes the geometric pruning;

[0072] 2) Construct a compatibility matrix:

[0073] Initial rough correspondence set All possible matching pair combinations, total Group, construct a dimension as Correlation matrix Iterate through all matching pairs and fill the matrix elements according to the results of the geometric pruning described above. :

[0074] ;

[0075] 3) Further filtering to generate a refined correspondence set. :

[0076] For matrix Perform eigenvalue decomposition or exponentiation to obtain the principal eigenvector corresponding to the largest eigenvalue. , The dimension is ;

[0077] vector elements in For the first The scalar value of the row also represents the first row. The geometric consistency scores of each matching pair are sorted from highest to lowest, and the top scores are retained. Each matching pair generates a fine correspondence set. Fine correspondence set This constitutes the maximal clique in geometric graph theory, meaning that all preserved point pairs are rigidly consistent in spatial structure.

[0078] The specific steps for step five are as follows:

[0079] 1) Traverse the fine correspondence set Each index pair in, are the original physical coordinates of the source node and the target node in the c-th matching pair, respectively, and both are... Column vectors From the source skeleton coordinate matrix Obtained from the middle, From the target skeleton coordinate matrix Obtained from the middle;

[0080] 2) Geometric centroid calculation: The mean of all extracted node coordinate vectors is calculated to obtain the physical geometric centroids of the largest source skeleton and the target skeleton.

[0081] ;

[0082] In the formula, The physical geometric centroid of the largest source skeleton; The physical geometric centroid of the target skeleton of the largest cluster. The total number of matching pairs in the largest clique; For the iterative index of the summation operation, ;

[0083] 3) Decentralized mapping: Perform explicit vector subtraction to translate all points to a local coordinate system with the centroid as the origin, generating relative coordinate vectors:

[0084] ;

[0085] In the formula, Let be the relative coordinate vector of the source node in the c-th matching pair after the decentralization of the maximum clinker source skeleton. The relative coordinate vector of the target node in the c-th matching pair after the decentralization of the target skeleton of the maximum clique target;

[0086] 4) Relative coordinate set encapsulation: Stack all the obtained relative coordinate vectors column-wise to construct a relative coordinate set matrix. and :

[0087] ;

[0088] In the formula, The matrix of relative coordinates of the largest clique source skeleton, dimension ; The relative coordinate set matrix of the skeleton of the maximum clique target, dimension ;

[0089] 5) For the relative coordinate set matrix and Perform global attraction domain anchoring rotation analysis:

[0090] Covariance matrix construction: Calculate the point cloud reflecting the geometric correlation between the source skeleton point cloud and the target skeleton point cloud of the largest clique using relative coordinate sets. Attraction domain cross-covariance matrix :

[0091] ;

[0092] in, This represents the matrix transpose operation;

[0093] right Perform singular value decomposition (SVD) to decompose it into three orthogonal components:

[0094] ;

[0095] In the formula, Let represent the left singular vector matrix, the singular value diagonal matrix, and the transpose of the right singular vector matrix obtained from SVD decomposition, respectively.

[0096] Closed-form anchoring of rotation matrix: Based on the principle of orthogonal Pluke analysis, the optimal rotation matrix can be directly derived analytically. :

[0097] ;

[0098] To introduce a determinant term, used to force a correction to the reflection matrix generated by SVD, ensuring Satisfy rigid body rotation constraints in a right-handed coordinate system;

[0099] 6) For rotation matrices The physical geometric centroids of the maximum source skeleton and the target skeleton Perform translation calculations and matrix assembly:

[0100] Back-substitution for translation: using the transformed... and Given overlapping constraints, solve for the global anchoring translation vector. :

[0101] ;

[0102] Final matrix assembly: The rotation matrix and translation vector are combined into a standard homogeneous transformation matrix as the anchoring transformation matrix. :

[0103] ;

[0104] Will As the initial rigid transformation command.

[0105] The specific steps for step six are as follows:

[0106] 1) Construct a virtual sliding plane;

[0107] Fine correspondence set The Middle The original physical coordinates of the source node and the target node in the matching pair are respectively used as follows: This indicates that the corresponding coordinates are also used to represent the target node. Feature vectors are extracted from the neighborhood point set to construct a local orthogonal basis. ,in, For the target node Surface normal vector at that location; Target nodes Two orthogonal tangent vectors on the tangent plane;

[0108] Constructing anisotropic covariance matrices using the spectral decomposition theorem A planar structure is constructed, and extremely large tangential uncertainty and extremely small normal uncertainty are artificially set:

[0109] ;

[0110] In the formula, Let the variances of the uncertainties along the normal and tangential directions be defined respectively. ;

[0111] 2) Solve the probabilistic sliding optimization using Lie algebra:

[0112] Constructing the sliding loss function , The sum of Mahalanobis distances for matching pairs in all fine correspondence sets is used with the goal of minimizing the loss.

[0113] ;

[0114] in, The total number of matching pairs in the largest clique; Let be the position residual vector. ;

[0115] because variance in the tangential direction The maximum value is its inverse matrix. Since the eigenvalue weights in the tangential direction approach 0, the tangential residual is eliminated in the position residual vector formula, and only the normal residual term is retained; thus, each source node gains the degree of freedom to slide freely on the target surface until the normal error is reduced to zero.

[0116] 3) Dynamic evolution is performed through Gauss-Newton iteration, and the optimal perturbation vector that minimizes the loss is solved in the Lie algebra space. ;

[0117] For each matching pair, calculate the partial derivative of the position residual with respect to the perturbation vector, i.e., the Jacobian matrix. :

[0118] ;

[0119] in, Indicates the first A number of matching pairs Jacobian matrix; Let represent a 6-dimensional perturbation vector in the Lie algebra space, which is used as an independent variable in the differentiation here, and is also the pose update increment to be solved in this iteration; express identity matrix; Represents the coordinates of the skeleton nodes in the source point cloud; This indicates an antisymmetric operator, used here to transform a 3D vector Mapped to Antisymmetric matrix; This indicates that the current rigid body transformation matrix is ​​applied to the source point cloud skeleton node. The resulting transformed coordinate vector ; This represents the rotation matrix to be solved, used to describe the pose change of the source point cloud; This represents the translation vector to be solved, used to describe the spatial displacement of the source point cloud;

[0120] Construct and solve the linear increment equation to obtain the optimal Lie algebra update. :

[0121] ;

[0122] in, This represents the inverse matrix of the anisotropic covariance matrix constructed for the corresponding target skeleton node, used to weight the errors in different directions during optimization;

[0123] 4) Pose manifold update:

[0124] The update value obtained by using exponential mapping Restore to a rigid body transformation matrix and apply to the current pose. Update to obtain a new pose. :

[0125] ;

[0126] in, These represent the updated target pose matrix and the current pose matrix before the update, respectively. This represents the matrix exponential function used to map the Lie algebra increments in the tangent space back to the Lie group manifold.

[0127] Through the above design scheme, the present invention can bring the following beneficial effects:

[0128] This invention constructs a data flow-driven 3D point cloud processing system. At the physical level, the system consists of a 3D laser scanning acquisition terminal and a host computer connected to it via a high-speed data transmission interface. At the logical operation level, the host computer's memory is configured with a processing link consisting of a cascaded discretization feature mapping module, an asymmetric topology construction module, a geometrically robust screening module, and a probability metric refinement module.

[0129] This invention can achieve highly robust registration in challenging scenarios with extremely low overlap (such as the 3DLoMatch dataset). It solves the zero-solution problem under low overlap by using an asymmetric association strategy and a one-way matching union. Compared with the scheme that enables reciprocity constraints, it significantly improves the recall rate and solves the survival problem in extreme environments.

[0130] Meanwhile, as shown in Table 1, the comparative experimental results demonstrate that by employing DGTG geometric pruning and probabilistic manifold optimization mechanisms, and utilizing anisotropic covariance and Lie algebra iteration to eliminate voxel quantization errors, the translation and rotation errors are significantly reduced compared to the traditional point-to-point ICP method. In the standard scene (3DMatch dataset), a recall rate of 84.84% is achieved, and the average rotation error is reduced to 1.29 degrees, and the translation error to 5.88 cm. This accuracy metric approaches the physical limit of 5 cm voxel resolution for the input data, proving that the manifold optimization mechanism can effectively overcome the accuracy bottleneck caused by voxelization.

[0131] Furthermore, this method has been optimized so that it can run efficiently in parallel through multiple processes on a regular quad-core CPU. A single frame can be processed in about 19.00 seconds, without relying on GPU hardware, and has extremely high engineering practical value.

[0132] This method not only guarantees the existence of solutions even with extremely low overlap rates, but also eliminates voxel quantization errors through covariance-based manifold optimization, achieving high-precision and high-efficiency registration without the need for a GPU. Attached Figure Description

[0133] The present invention will be further described below with reference to the accompanying drawings and specific embodiments:

[0134] Figure 1 A diagram comparing the loss of correspondence caused by traditional reciprocal consistency constraints with the principle of the asymmetric association strategy of this invention;

[0135] Figure 2 This is a schematic diagram of the maximum clique topology screening and pruning process based on geometric side length consistency constraints in an embodiment of the present invention.

[0136] Figure 3 This is a comparison diagram of the error of traditional point-to-point rigid connection and the probabilistic plane sliding optimization principle based on covariance of the present invention. Detailed Implementation

[0137] A low overlap rate point cloud registration method based on topology-metric decoupling is proposed. The method utilizes a low overlap rate point cloud registration system based on topology-metric decoupling. The system includes a 3D laser scanning acquisition terminal and a host computer connected to it. The host computer includes a manually set port, a discretization feature mapping module, an asymmetric topology construction module, a geometric robust screening module, and a probability metric refinement module.

[0138] The three-dimensional laser scanning acquisition terminal is used to scan and obtain point cloud data of the scanned object in two selected different directions, and is recorded as the source point cloud and the target point cloud respectively.

[0139] The discretization feature mapping module performs coordinate mapping and matrixing on the source point cloud and the target point cloud respectively, obtaining the source matrix and the target matrix respectively; then, the source matrix and the target matrix are subjected to grid denoising to obtain two sets of sparse skeleton node matrices, the source skeleton coordinate matrix and the target skeleton coordinate matrix, and a spatial directory is constructed to form a KD-Tree index; then, the FPFH operator is used to perform geometric encoding on each skeleton node, generating two sets of region convolution fusion FPFH feature matrices that do not depend on absolute coordinates, namely the source feature matrix and the target feature matrix, thereby realizing the mapping from geometric space to feature space;

[0140] The asymmetric topology construction module is equipped with a one-way nearest neighbor indexer. Through the asymmetric topology construction module, each row of the source feature matrix is ​​traversed in one direction to find the best matching row of the target feature matrix, forming a matching pair, and then generating an initial rough correspondence set representing the index correspondence relationship, establishing a preliminary logical connection between the source point cloud and the target point cloud.

[0141] The geometric robust screening module is equipped with a RANSAC controller for iterative sampling. The geometric robust screening module takes two random matching pairs from the initial coarse correspondence set. Using rigid constraints and the maximum clique consistency principle of DGTG dynamic geometric topology gating, outliers are removed, i.e., erroneous connections are removed from the noisy index table. The geometrically consistent correspondence set, i.e. the maximum clique structure, is obtained. The fine correspondence set is then screened again and output, realizing topology pruning and retaining the skeleton nodes with correct indexes.

[0142] The probability metric refinement module performs a global attraction domain anchoring operation on the fine correspondence set to obtain the anchoring transformation matrix as the initial value for the next step of probability plane sliding. Then, for the maximum clique target skeleton node in the fine correspondence set, an anisotropic probability field model based on the covariance matrix is ​​constructed. The correctly indexed skeleton node is modeled as a probability manifold extending along the tangent plane, and a virtual sliding plane is constructed to restore the metric accuracy lost due to discretization. Then, through the generalized iterative nearest point GICP optimizer, a probability sliding loss function is constructed to minimize the Mahalanobis distance, driving the maximum clique source node to slide along the tangent plane to the coincident position, realizing sub-voxel level fusion of the source point cloud and the target point cloud, and outputting the final high-precision transformation matrix to achieve point cloud registration. Example

[0143] This embodiment is built based on the Python programming environment and the Open3D open-source library, and runs on a host computer equipped with an Intel Core i5-8265U (quad-core, 1.60GHz). The implementation steps of this invention are described in detail below with reference to the accompanying drawings:

[0144] Step 1: Use a 3D laser scanning acquisition terminal to scan and obtain point cloud data of the scanned object in the selected 0° direction and the corresponding 90° direction, and record them as the source point cloud and the target point cloud respectively.

[0145] Step 2: Topology-anchored dual-channel dimensionality reduction and region convolution fusion for FPFH feature embedding;

[0146] The core task of this step is to clean and characterize the data. The discretization feature mapping module, for the source point cloud and target point cloud to be registered, sequentially performs coordinate mapping, mesh denoising, and spatial topology construction, transforming the massive and noisy physical measurement points (i.e., the source and target point clouds) into two sets of sparse skeleton node matrices. Simultaneously, the region convolution fusion FPFH operator is used to geometrically encode each node, generating two sets of high-dimensional feature matrices that do not depend on absolute coordinates, establishing a mathematical foundation for subsequent asymmetric association.

[0147] in:

[0148] 1. Target of processing:

[0149] Input data stream: Contains two independent sets of raw point cloud data packets:

[0150] Source Cloud Mobile data to be registered;

[0151] Target point cloud : Reference data used as a benchmark.

[0152] 2. Perform coordinate analysis and matrix transformation:

[0153] Objective: To analyze the source and target data separately and transform unstructured signals into mathematical matrices.

[0154] Operation: Read two raw data streams and extract spatial components. .

[0155] Matrix construction: Allocate space in memory separately and construct the matrix. Source point cloud floating point type coordinate matrix and Target point floating-point coordinate matrix At this point, each row in the matrix represents the precise spatial location of a physical point.

[0156] 3. Perform mesh simplification and downsampling:

[0157] Objective: To utilize voxel filtering technology to... and Spatial downsampling is performed separately to unify data density.

[0158] Operation: Set uniform voxel side length For any point in a floating-point coordinate matrix Calculate the three-dimensional index vector of its voxels. :

[0159] ;

[0160] 2) Clustering: Retrievals falling within the same index of For each point, calculate its geometric centroid as the unique skeleton node. :

[0161] ;

[0162] In the formula, This indicates the first voxel in the source or target point cloud that falls within the current voxel grid. One point, ;

[0163] Results: The amount of data was drastically reduced, resulting in two sets of sparse skeleton node matrices: the source skeleton coordinate matrix. Dimension and the target skeleton coordinate matrix Dimension ,in, These represent the number of skeleton nodes in the source point cloud and the target point cloud after dimensionality reduction, respectively.

[0164] 4. Construct a KD-Tree index:

[0165] Objective: To construct spatial catalogs for the two sets of skeleton data respectively.

[0166] Operation: For each and Construct using recursive binary search method and .

[0167] Function: To ensure that subsequent operations, whether calculating neighbor lookups within the source point cloud for FPFH geometric features or finding corresponding points in the target point cloud during point cloud registration, can be performed smoothly. Level of rapid response, This indicates that the retrieval time increases logarithmically with the number of nodes, where This represents the number of skeleton nodes in the point cloud.

[0168] 5. Compute region convolution to fuse FPFH geometric features:

[0169] Objective: To separately and For each point in the graph, a "geometric ID card" is calculated, generating a feature descriptor that does not depend on absolute coordinates.

[0170] Step 1: Normal Vector Estimation

[0171] Operation: Use a KD-Tree search to find the first element in the skeleton coordinate matrix. Nodes of Calculate the local centroid from each neighboring node. And construct the covariance matrix :

[0172] ;

[0173] In the formula, Represents a node The 1 neighboring node, ;

[0174] This means The geometric centroid of each neighboring node;

[0175] Detailed explanation of eigenvalue decomposition:

[0176] Mathematical process: For real symmetric covariance matrices Perform orthogonal decomposition and solve the equation:

[0177] ;

[0178] Three eigenvalues ​​were obtained and its corresponding eigenvectors ,and .

[0179] Physical meaning: This step is essentially principal component analysis (PCA). This represents the principal direction with the largest variance in the data distribution, while This represents the direction with the smallest variance in the data distribution.

[0180] Conclusion: Geometrically, the direction perpendicular to the local plane exhibits the least change in data; therefore, [the following is a possible interpretation:] Define as a node normal vector ;

[0181] Step 2: Relative Angle Encoding SPFH:

[0182] Operation: To eliminate the effects of coordinate system rotation, at the nodes... Construct a local orthogonal coordinate system for the Darboux frame. ,in, , indicating that the normal vector is axis; Determined by the cross product of the line and the normal. axis; Determined by the right-hand rule axis;

[0183] The purpose of the Dabu frame is to: construct a "follow-up" local reference, ensuring rotational invariance of angular features; and calculate the angular deviation of neighboring normal vectors within this coordinate system. Generate a 33-dimensional histogram relative angle code. ;

[0184] Calculate the angular deviation of the neighbor normal vector in this coordinate system The specific calculation formula is as follows:

[0185] ;

[0186] ;

[0187] ;

[0188] In the formula, Represents the normal vector of the neighbor node and Projection angle deviation between axes; Represents the normal vector of the center node and , The angle of inclination between the vectors connecting two points; Represents the normal vector of the neighbor node exist Orientation rotation angle deviation on the plane.

[0189] The values ​​of the three angular deviations mentioned above are divided into 11 equally wide intervals for projection counting;

[0190] Step 3: Region Convolutional Fusion (FPFH):

[0191] Operation: Perform distance-weighted convolution to transform the nodes The relative angle encoding SPFH and its neighboring nodes The relative angle encoding SPFH fusion is performed using the following formula:

[0192] ;

[0193] In the formula, Represents the merged nodes The final 33-dimensional feature vector; For nodes Relative angle encoding; For nodes with neighboring nodes Point-to-point distance, For nodes Relative angle encoding;

[0194] The purpose of convolution is to expand the receptive field of features. By using distance weights, it allows previously isolated point features to incorporate contextual information from the macroscopic environment.

[0195] Step 4: Matrix Assembly and Final Output (MatrixAssembly) Assembly Logic:

[0196] For the source skeleton coordinate matrix All nodes Each node is traversed step one through step three, and the corresponding 33-dimensional FPFH feature vector is calculated for each node.

[0197] This The vectors are stacked in memory in row-wise stacking order according to their index, thereby constructing the source feature matrix. ;

[0198] Similarly, for the target skeleton coordinate matrix By traversing and stacking, the target feature matrix is ​​constructed. ;

[0199] in, ;

[0200] Each row in the above matrix structure represents the geometric fingerprint of the corresponding node:

[0201] Dimensions 1 to 11: These mainly characterize whether a local surface is close to a planar structure. The higher the value in this range, the smaller the deviation of the normal direction in the neighborhood.

[0202] Dimensions 12 to 22: These mainly characterize whether a local region exhibits a cylindrical or rod-like structure. The higher the value in this range, the more significant the curvature change is in the direction of the line connecting the points.

[0203] Dimensions 23 to 33: mainly characterize whether a local region contains right angles or obvious edge structures. The higher the value in this range, the closer the angle between the neighborhood normals is to orthogonal.

[0204] In summary, the two point clouds containing only coordinates have been transformed into two enhanced data structures with both skeletal coordinates and geometric features, thus fully preparing for feature matching in step two.

[0205] Step 3: Initial construction of the topological skeleton based on the ACE asymmetric association expansion strategy;

[0206] The purpose of this step is to establish a preliminary logical connection between the source and target point clouds using an asymmetric topology construction module. This utilizes the FPFH feature matrix generated in step one, no longer based on three-dimensional Euclidean space. Instead of performing a retrieval, it performs a nearest neighbor search in the 33-dimensional feature space.

[0207] This embodiment uses, as follows: Figure 1 The asymmetric strategy shown, namely unidirectional traversal of source features to find the best match of target features, aims to retain as many potential homologous point pairs as possible in low overlap scenarios, providing sufficient candidate samples for subsequent topology screening.

[0208] 1. Processing the input feature matrix of the object:

[0209] Source feature matrix Dimension ;

[0210] Target feature matrix Dimension ;

[0211] Auxiliary index: The source skeleton coordinate matrix generated in step one With the target skeleton coordinate matrix This is used to subsequently restore the feature index to the physical coordinate index.

[0212] 2. Construct a high-dimensional index

[0213] Purpose: In order to solve the problem in The problem of low efficiency in linear search within a 33-dimensional vector necessitates the construction of a high-dimensional index.

[0214] Operation: Transfer the target feature matrix Considered a A high-dimensional dataset with [number] data points. An index structure is constructed using the FLANN (FastLibrary for Approximate Nearest Neighbors) algorithm or a high-dimensional KD-Tree. . It is an index tree for a 33-dimensional feature space, used to divide the 33-dimensional feature space into several hyperplane regions, so that similar geometric fingerprints are adjacent to each other in memory address, thereby accelerating retrieval.

[0215] 3. Perform asymmetric nearest neighbor search;

[0216] Objective: To implement a one-way matching strategy of "finding the target from the source" and establish an initial connection.

[0217] Traversal logic: Traverse the source feature matrix Each row of feature vectors in the dataset is used as the query vector.

[0218] Distance metric: in The system searches for the target feature vector that is closest to the feature vector in the current row. It initially determines that the two feature vectors are identical, belonging to the same node, also known as a common source point. The distance metric used is... Norm, to obtain the best matching target index found in the search. :

[0219] ;

[0220] In the formula, The source characteristic matrix is ​​the first The eigenvectors of the row, ; For the target feature matrix and The nearest row feature vector, .

[0221] Combination Figure 1 Explanation of technical principles and effects:

[0222] Figure 1 (a) shows a traditional point cloud registration algorithm (such as ICP or RANSAC), which uses a "reciprocity check" mechanism, requiring that A to B and B to A be the nearest neighbors. However, in low-overlap scenarios, the target point cloud may lack corresponding overlapping regions, causing the reverse search to "drift," resulting in a large number of correct one-way matches, such as... Figure 1 In (a), during forward matching, the red dot points to the blue square, which is the correct edge. However, during reverse matching, the blue square points to the red cross. The red cross and the red dot do not coincide, which causes the bidirectional check to fail. Consequently, the correct edge is incorrectly deleted (FalseNegative), resulting in a broken skeleton.

[0223] Figure 1 (b) is a schematic diagram illustrating the use of an asymmetric association strategy instead of reverse consistency constraints in this embodiment. Only unidirectional matching is performed, requiring only the features of the source point cloud. The target point is considered in the feature space. If it is the best match, it will be retained.

[0224] Conclusion: Although this relaxation strategy introduces some noise, it preserves the topological skeleton of low-overlap regions to the maximum extent and avoids premature removal due to differences in perspective.

[0225] 4. Result Mapping and Set Generation:

[0226] Objective: To convert mathematical matching in feature space into physical indexes of point clouds.

[0227] Operation: For each source index traversed Record the best target index found in the search. This forms an index-correspondence tuple, or simply an index pair. .

[0228] Collection construction: All The pairs converge to form the initial rough correspondence set. : .

[0229] Data structure state: at this time It is An integer matrix, where each row represents a pair of predicted origin points, and the total number of rows in the set is . Under asymmetric full matching, .

[0230] Output status:

[0231] A set of hypothetical connections has been established: Although we do not yet know how the two sets of point clouds are rotated and translated, the algorithm has already determined that the two points are from the same source. For example, "point 105 in the source point cloud (table corner)" may correspond to "point 2048 in the target point cloud (table corner)".

[0232] Conclusion: This step achieved a metric association that "only recognizes shape, not location," laying the foundation for the next step of eliminating erroneous connections through topological relationships.

[0233] Step 4: Maximum clique consistency screening based on DGTG dynamic geometric topology gating;

[0234] Purpose of this step: The core task of this step is to use the geometrically robust screening module to shift the registration focus from the "high-dimensional feature space" back to the "three-dimensional physical space". Figure 2 This is a schematic diagram of the maximum clique topology screening and pruning process based on geometric edge length consistency constraints in an embodiment of the present invention. The initial rough correspondence set that has already been obtained is, for example... Figure 2As shown in (a), this involves a large number of assumptions with similar features but incorrect physical locations. Therefore, utilizing the physical axiom of "Euclidean distance between two points remains unchanged" in rigid body transformation, a topologically consistent graph model is constructed. By calculating the geometric residual matrix and performing eigenvalue decomposition, all outliers violating rigid constraints are eliminated, ultimately locking a set of fine correspondences with strictly locked geometric structures, as shown in (a). Figure 2 The maximal clique structure shown in (c) is a large clique structure.

[0235] 1. Target of processing:

[0236] Input data stream:

[0237] Initial rough correspondence set :Include For candidate matching relationships, corresponding Figure 2 The messy connections in (a).

[0238] Source skeleton coordinate matrix and target skeleton coordinate matrix Used to query the index corresponding to each index Physical coordinates.

[0239] Physical Euclidean distance between two points in the source cloud The calculation formula is as follows: The calculation results correspond to Figure 2 (b) ,in, and Represent the source skeleton coordinate matrix respectively The index is and The three-dimensional coordinate vector corresponding to the row; This represents the straight-line distance between these two skeleton nodes in the physical space of the source point cloud.

[0240] The corresponding target point cloud Physical Euclidean distance between two points The calculation formula is as follows: The calculation results correspond to Figure 2 (b) ,in, and These represent the target skeleton coordinate matrices respectively. The index is and The three-dimensional coordinate vector corresponding to the row; This represents the straight-line distance between these two skeleton nodes in the physical space of the target point cloud.

[0241] 2. Perform geometric residual calculations and pruning:

[0242] Object to be processed: Initial coarse correspondence set Any two matching pairs in and .

[0243] Operation: Execute the "virtual ranging" operation, such as... Figure 2 As shown in (b). First, from the source skeleton coordinate matrix... Read source node and Calculate the distance from the coordinates. Then from the target skeleton coordinate matrix Read target node and Calculate the distance from the coordinates. and distance .

[0244] The absolute value of the difference between the lengths of the two segments is calculated using the geometric residual formula. :

[0245] ;

[0246] Perform logical judgments, corresponding to Figure 2 (b):

[0247] Set the tolerance threshold for geometric consistency determination ;

[0248] like ,For example Figure 2 (b) In the diagram, Then determine and This matching pair is geometrically compatible and preserves logical connections.

[0249] like ,For example Figure 2 (b) Showing Then determine and This matching pair must contain an error, and its connection needs to be severed to complete geometric pruning;

[0250] 3. Construct a compatibility matrix:

[0251] Processing object: All possible matching pair combinations, total Group.

[0252] Operation: Construct a dimension of Correlation matrix Iterate through all combinations and fill the matrix elements according to the results of the geometric pruning described above. :

[0253] ;

[0254] Product generated: A "geometric compatibility map" was generated, such as Figure 2 As shown in (c).

[0255] Physical meaning: This is essentially an adjacency matrix of an undirected graph. For example... Figure 2 In (c), nodes are "matching pairs," and edges represent "mutual support." Correct matches form a tight "clique," while incorrect matches... Figure 2 The isolated node in (c).

[0256] 4. Perform eigenvalue decomposition and maximum clique extraction:

[0257] Object to be processed: The matrix that was just generated .

[0258] operate:

[0259] Eigenvalue decomposition: the real symmetric incidence matrix Perform eigenvalue decomposition or exponentiation to find the principal eigenvector corresponding to its largest eigenvalue. Dimension .

[0260] Weighted sorting: vector Each element in For the first The scalar value of the row also represents the first row. The geometric consistency score of a pair of matches, also known as the authority score (Eigenvector Centrality). A higher score indicates that the pair is more geometrically consistent. Figure 2 (c) shows the center of the global suction basin.

[0261] Truncate the output: Sort the scores from highest to lowest, keeping only the top scores. A number of matching pairs.

[0262] Product generated: Fine correspondence set .

[0263] Physical meaning: such as Figure 2 As shown in (c), this set constitutes the maximal clique in geometric graph theory, meaning that all preserved point pairs are rigidly consistent in spatial structure.

[0264] 5. Summary of physical meaning and data flow:

[0265] Input state: a large number of assumptions, such as Figure 2 As shown in (a), the data contains a mixture of real structure and a large number of random noise lines, such as... Figure 2 The line between the two blue dots on the right in (a) is the noise line, also known as the initial chaotic connection, where the signal-to-noise ratio is extremely low.

[0266] Core action: Geometric invariant check, such as Figure 2 As shown in (b). Utilizing the physical constraint of rigid body transformation: the invariance of side lengths, each pair of assumptions is... Quick complexity check:

[0267] Preservation logic: such as Figure 2 (b) As shown in the upper part, when the distance between two points in the source point cloud is visually consistent with the distance between the corresponding two points in the target point cloud (e.g., the difference is within the tolerance), the connection is determined to meet the rigid constraint. The connection is marked with a bold black equal sign "=" in the figure and labeled "Reserved".

[0268] Pruning logic: such as Figure 2 (b) As shown in the lower part, when the lengths of the source connection (such as a short line) and the target connection (such as a long line) are significantly mismatched (e.g., 0.5m vs. 1.5m), the connection is determined to violate physical constraints and is forcibly cut off. In the figure, the connection is marked with a bold black inequality sign “≠” and labeled with the word “pruning”.

[0269] This process cuts off all erroneous connections with mismatched lengths, among which... This indicates constant-time complexity, meaning that the computation time is fixed and extremely low, and does not increase with the size of the data.

[0270] Output state: maximally clique structure, such as Figure 2 As shown in (c), after the above pruning, only the stable connection structure centered on the black node remains in the figure. Thus, we have obtained a set of point pairs that are completely identical in spatial geometry, namely the maximal clique in the figure. This set of point pairs no longer depends on fuzzy feature similarity, but rather on strict physical distance constraints.

[0271] Function: It provides an absolutely reliable mathematical anchor for solving the unique transformation matrix in the fifth step.

[0272] Step 5: Closed Pose Calculation Based on Global Basein Anchoring

[0273] This step is the initial pose establishment stage of the registration framework of this invention. Through the probabilistic metric refinement module, a global attraction domain anchoring operation is performed on the fine correspondence set, i.e., the maximal clique structure. This step completely white-boxes the pose solution process by explicitly constructing a decentralized relative coordinate system and utilizing the singular value decomposition (SVD) of the attraction domain cross-covariance matrix to non-iteratively anchor the globally optimal rotation matrix in one step. With translation vector This establishes the central location of the global attraction basin, providing a high-confidence initial value for subsequent subvoxel refinement.

[0274] 1. Macro-level processing objects:

[0275] Overall input stream:

[0276] Fine correspondence set This is a set of "maximum cluster skeleton indexes" that have been preserved after geometric pruning and strictly satisfy rigid constraints in spatial structure.

[0277] Source skeleton coordinate matrix and target skeleton coordinate matrix : As a physical coordinate lookup table, it is used to extract the specific three-dimensional values ​​of each skeleton node in the source point cloud and the target point cloud based on the index.

[0278] Overall output stream:

[0279] Anchoring transformation matrix The unique rigid transformation solution output after global anchoring calculation describes the pose of the center of the globally optimal attraction domain.

[0280] 2. Perform white-box construction of the relative coordinate set to achieve centroid alignment and decoupling:

[0281] The purpose of this step is to explicitly eliminate the interference of translation components on rotation estimation and construct a linearly independent pure rotation feature field. Traditional iterative algorithms (such as ICP) are like groping in a fog at low overlap, easily getting trapped in local optima. This step uses a pure skeleton filtered by DGTG to directly calculate the analytical solution of the energy function. The pose is then "pinned" to the lowest point of the Global Attraction Basin.

[0282] Detailed operating steps:

[0283] Physical coordinate extraction: Traversing the fine correspondence set Each index pair in the matrix is ​​derived from the source skeleton coordinate matrix. and target skeleton coordinate matrix Extract the original physical coordinates of the source node and the target node from the c-th matching pair respectively. and .

[0284] Geometric centroid calculation: The mean of all extracted node coordinate vectors is calculated to obtain the physical geometric centroids of the maximum source skeleton and the target skeleton.

[0285] ;

[0286] In the formula, The physical geometric centroid of the largest source skeleton; The physical geometric centroid of the target skeleton of the largest cluster. This represents the total number of matching pairs in the largest clique. For the iterative index of the summation operation, ;

[0287] Decentralized mapping: Performs explicit vector subtraction to translate all points into a local coordinate system with the centroid as the origin, generating relative coordinate vectors.

[0288] ;

[0289] In the formula, Let be the relative coordinate vector of the source node in the c-th matching pair after the decentralization of the maximum clinker source skeleton. The relative coordinate vector of the target node in the c-th matching pair after the decentralization of the target skeleton of the maximum clique target;

[0290] Relative coordinate set encapsulation (SetAssembly): Stacks all calculated relative coordinate vectors column-wise to construct a relative coordinate set matrix. and :

[0291] ;

[0292] In the formula, The matrix of relative coordinates of the largest clique source skeleton, dimension ; The relative coordinate set matrix of the skeleton of the maximum clique target, dimension .

[0293] This stage produces: two sets of relative coordinate matrices that satisfy the zero-mean property. and .

[0294] 3. Global Basin Anchoring Rotation Resolution:

[0295] This step deals with the relative coordinate set matrix. and .

[0296] Detailed operating steps:

[0297] Covariance matrix construction: Calculating the covariance matrix reflecting the geometric correlation between two sets of point clouds using relative coordinate sets. Attraction domain cross-covariance matrix :

[0298] ;

[0299] in, This represents the matrix transpose operation;

[0300] Singular Value Decomposition (SVD): for Perform singular value decomposition (SVD) to decompose it into three orthogonal components:

[0301] ;

[0302] In the formula, Let represent the left singular vector matrix, the singular value diagonal matrix, and the transpose of the right singular vector matrix obtained from SVD decomposition, respectively.

[0303] Closed-form anchoring of rotation matrix: Based on the principle of orthogonal Procrustes analysis, the optimal rotation matrix is ​​directly obtained analytically. :

[0304] ;

[0305] Introducing determinant terms into the formula This is to forcibly correct the reflection matrix that SVD may produce, to ensure Strictly satisfy rigid body rotation constraints in the right-hand coordinate system.

[0306] This step produces: a mathematically globally optimal anchoring rotation matrix. .

[0307] 4. Translation calculation and matrix assembly (this step):

[0308] Object to be processed: Rotation matrix With geometric centroid .

[0309] 6) For rotation matrices The physical geometric centroids of the maximum source skeleton and the target skeleton Perform translation calculations and matrix assembly:

[0310] Back-substitution for translation: using detailed operation steps:

[0311] Back-substitution for translation: Using the constraint that "the transformed source centroid should coincide with the target centroid", solve for the translation vector. :

[0312] ;

[0313] Final matrix assembly: The rotation matrix and translation vector are combined into a standard homogeneous transformation matrix, which serves as the anchoring transformation matrix. :

[0314] ;

[0315] This stage outputs: Initial rigid transformation command. , as the initial rigid transformation command.

[0316] It is a mathematically interpretable and reproducible global optimal pose solution, which will be used as the initial value for the probabilistic plane sliding in step five.

[0317] Step Six: Probabilistic Plane Sliding Refinement Based on Uncertain Manifold

[0318] This step is the final accuracy rescue (MetricRescue) stage of the registration framework of this invention. In the probabilistic metric refinement module, a special "virtual sliding plane" is constructed to address the voxelization quantization error introduced by computational efficiency in the previous steps. Unlike traditional algorithms that force "hard connections" between points, this step uses the anisotropic covariance matrix to mathematically construct a "probabilistic lubrication layer." By minimizing the energy function based on Mahalanobis distance, the algorithm only penalizes normal deviations and ignores tangential sliding, thereby driving the source point cloud to automatically "slide" into the sub-voxel accuracy true position in the Lie group space.

[0319] 1. Macro-level processing objects:

[0320] Overall input stream:

[0321] Anchoring transformation matrix At this point, although the point cloud is roughly aligned, the feature points are locked in the center of the discrete voxel grid (e.g., 5cm resolution), resulting in significant quantization errors.

[0322] Fine correspondence set: including Matching index.

[0323] Overall output stream:

[0324] Refine the transformation matrix Breaking through the limitations of voxel resolution, achieving a final pose accuracy of less than 1 cm at the sub-voxel level.

[0325] 2. Perform anisotropic covariance modeling to construct a "virtual sliding plane";

[0326] The objects processed in this step are:

[0327] Input: The relative coordinate vectors of the target nodes after the decentralization of the maximum clique target skeleton and the geometric features of their spatial neighborhoods.

[0328] Output: The anisotropic covariance matrix corresponding to this node.

[0329] Detailed mathematical operations:

[0330] Basis vector extraction: fine correspondence set The Middle The original physical coordinates of the source node and the target node in the matching pair are respectively used as follows: This indicates that the corresponding coordinates are also used to represent the target node. Principal component analysis (PCA) is performed on the neighborhood point set to extract eigenvectors and construct a local orthogonal basis. ,in, For the target node Surface normal vector at that location; Target nodes Two orthogonal tangent vectors on the tangent plane.

[0331] Planar construction formula: Construct the anisotropic covariance matrix using the spectral decomposition theorem. In order to achieve mathematically Figure 3 (b) shows a “planar” effect, in which we artificially set extremely large tangential uncertainty and extremely small normal uncertainty:

[0332] ;

[0333] In the formula, Let the variances of the uncertainties along the normal and tangential directions be defined respectively. ;

[0334] Mathematical and physical implications: This matrix describes a "probability ellipsoid" that extends infinitely in the tangent plane direction and is extremely flat in the normal direction. This completes the mathematical transformation from "discrete points" to a "continuous plane".

[0335] 3. Perform "probability sliding" optimization, i.e., solve for Lie algebras:

[0336] The objects processed in this step are:

[0337] Input: Fine correspondence set The corresponding source skeleton nodes act as sliders; the constructed covariance field It serves as a sliding track.

[0338] Output: Optimal Lie algebra perturbation and the updated pose .

[0339] Detailed mathematical operations:

[0340] Constructing the sliding loss function: Defining the objective function The sum of Mahalanobis distances for matching pairs in all fine correspondence sets is used with the goal of minimizing the loss.

[0341] ;

[0342] in, The total number of matching pairs in the largest clique; Let be the position residual vector. ;

[0343] in, This represents the rotation matrix to be solved, used to describe the pose changes of the source point cloud; This represents the translation vector to be solved, used to describe the spatial displacement of the source point cloud;

[0344] because variance in the tangential direction The maximum value is its inverse matrix. Since the eigenvalue weights in the tangential direction approach 0, the tangential residual is eliminated in the position residual vector formula, and only the normal residual term is retained; thus, each source node gains the degree of freedom to slide freely on the target surface until the normal error is reduced to zero.

[0345] The dynamic evolution is performed using Gauss-Newton iteration, solving for the optimal perturbation that minimizes the loss in the Lie algebra space. ;

[0346] Jacobian matrix: Calculates the derivative of the residual with respect to the Lie algebra :

[0347] ;

[0348] in, Indicates the first A number of matching pairs Jacobian matrix; Let represent a 6-dimensional perturbation vector in the Lie algebra space, which is used as an independent variable in the differentiation here, and is also the pose update increment to be solved in this iteration; express identity matrix; Represents the coordinates of the skeleton nodes in the source point cloud; This indicates an antisymmetric operator, used here to transform a 3D vector Mapped to Antisymmetric matrix; This indicates that the current rigid body transformation matrix is ​​applied to the source point cloud skeleton node. The resulting transformed coordinate vector ;

[0349] Incremental Equations: Construct and solve linear incremental equations to obtain the optimal Lie algebra update amount. :

[0350] ;

[0351] in, This represents the inverse matrix of the anisotropic covariance matrix constructed from the corresponding target skeleton nodes, i.e., the information matrix, which is used to weight the errors in different directions during optimization.

[0352] The update value obtained by using exponential mapping Restore to a rigid body transformation matrix and apply to the current pose. Update to obtain a new pose. :

[0353] ;

[0354] in, These represent the updated target pose matrix and the current pose matrix before the update, respectively; This represents the matrix exponential function used to map the Lie algebraic increments in the tangent space back to the Lie group manifold. Final pose. No longer limited by a 5cm mesh, it has achieved a leap from discrete approximation to sub-voxel-level physical bonding.

[0355] 4. Summary of physical meaning and data flow:

[0356] Traditional state: voxel center lock-up effect, such as Figure 3 As shown in (a), the blue nodes represent the transformed source point cloud, the red nodes represent the discrete voxel centers of the target point cloud, and the black vertical lines on the right represent the real physical surfaces (such as walls). In traditional methods, the algorithm forcibly pulls the blue nodes towards the red nodes through rigid connections (the black diagonal lines in the figure). Since the red nodes are only the center of the mesh and not the actual wall locations, this "point-to-point" hard constraint means that the blue nodes can never reach the black vertical lines, resulting in an unavoidable quantization error.

[0357] Core action: Probabilistic plane sliding optimization, such as... Figure 3 As shown in (b), the target point (red node) is modeled as a probability field represented by a green ellipse using the anisotropic covariance matrix. Tangential relaxation: The black double-headed arrows drawn vertically inside the green ellipse vividly represent the significant uncertainty imposed in the tangential plane direction. This means the algorithm "tells" the source point that the positional difference along the wall direction (i.e., the arrow direction) is negligible. Normal Constraint: The green ellipse is compressed very flat in the horizontal direction, indicating a very strong constraint in the normal direction (perpendicular to the wall).

[0358] Output state: Sub-voxel level physical fit. Guided by the probability field, the blue node is no longer forcibly dragged to the red node, but instead, like a slider, "slides" along the tangential track of the green ellipse to the position that best conforms to the physical constraints. As a result, the blue node precisely fits the black vertical line (the real surface) in the horizontal direction, successfully breaking through the voxel resolution limit represented by the red node.

[0359] Function: Improves registration accuracy from "centimeter level" (voxel grid limit) to "millimeter level" (sub-voxel physical limit).

[0360] To verify the effectiveness of the registration method based on topology-metric decoupling proposed in this invention, a quantitative comparative experiment was conducted on an internationally recognized standard indoor scene dataset. The experimental results are shown in Table 1.

[0361] Table 1. Statistical table of comparative experimental results

[0362]

[0363] (1) Experimental setup and dataset introduction: This experiment selects two representative benchmark datasets to test the algorithm performance in normal and extreme scenarios respectively:

[0364] 3DMatch Dataset (High Overlap Rate): This dataset contains scan pairs with an overlap rate of over 30%, representing typical indoor 3D reconstruction scenes, and is used to test the accuracy and stability of algorithms in standard environments.

[0365] 3DLoMatch dataset (low overlap rate): This dataset is the core test field of this invention, with a point cloud overlap rate of only 10% to 30%. Due to the extremely small overlap area, traditional algorithms are prone to failure. This dataset is specifically used to verify the ability of this invention to solve the "low overlap rate registration problem".

[0366] (2) Evaluation Indicators: Table 1 uses three core indicators to quantify registration performance:

[0367] RR (Registration Recall): Measured as a percentage. It represents the probability that the algorithm successfully registers data. A higher RR indicates stronger robustness and less likelihood of failure.

[0368] RE (RotationError): Unit is degrees. It represents the deviation of the predicted attitude from the true attitude in terms of rotation angle. The lower the value, the higher the accuracy.

[0369] TE (Translation Error): The unit is centimeters. It represents the deviation of spatial displacement; the lower the value, the higher the accuracy.

[0370] (3) Comparative analysis of experimental results: Table 1 compares the performance differences between traditional methods (such as RANSAC, TEASER++) and deep learning methods (such as DGR, PointDSC) and the present invention (marked as GeoRescue(Ours) in the table).

[0371] On the 3DLoMatch dataset with low overlap (core advantage): the recall (RR) of traditional mainstream algorithms RANSAC-1M and RANSAC-4M is only 0.67% and 0.45% respectively, which means that traditional methods are almost completely ineffective (close to zero solution) in low overlap scenarios. In contrast, the recall of this invention reaches 41.27%. This huge performance improvement (nearly 100 times) strongly proves that the "asymmetric association strategy" in step three and the "maximum clique geometric pruning" in step four of this invention can effectively preserve sparse correspondences and solve the survival problem in extreme environments.

[0372] On the 3DMatch dataset with high overlap (accuracy advantage): Although it is a traditional non-learning method, the accuracy metrics of this invention surpass many deep learning models. As shown in the table, the rotation error (RE) of this invention is reduced to 1.29°, and the translation error (TE) is reduced to 5.88cm. Considering that the voxel grid resolution of the input data is 5cm, the error level of this invention is close to the physical limit. This directly verifies the effectiveness of the "probabilistic plane sliding based on covariance" mechanism in step six, proving that the mechanism successfully overcomes the voxel center locking effect and achieves sub-voxel level high-precision alignment.

[0373] In summary, the data in Table 1 objectively demonstrates that the present invention achieves both high robustness (i.e., high recall) in low overlap scenarios and high precision (i.e., low error) in standard scenarios without the need for GPU hardware acceleration, representing a significant technological advancement.

Claims

1. A low-overlap point cloud registration method based on topology-metric decoupling, characterized in that: Includes the following steps, And the following steps are performed in sequence: Step 1: Use a 3D laser scanning acquisition terminal to scan and obtain point cloud data of the scanned object in two selected different directions, and record them as the source point cloud and the target point cloud respectively; Step 2: Perform coordinate mapping and matrix transformation on the source point cloud and the target point cloud to obtain the source matrix and the target matrix respectively; Mesh denoising is performed on the source matrix and the target matrix respectively to obtain two sets of sparse skeleton node matrices: the source skeleton coordinate matrix and the target skeleton coordinate matrix, and a KD-Tree index is generated. Geometric encoding is performed on each skeleton node, resulting in two sets of region convolutional fusion FPFH feature matrices: the source feature matrix and the target feature matrix. Step 3: Traverse each row of the source feature matrix in one direction to find the best matching row of the target feature matrix, form a matching pair, and then generate an initial rough correspondence set representing the index correspondence relationship, and establish a preliminary logical connection between the source point cloud and the target point cloud. Step 4: Take two random matching pairs from the initial rough correspondence set, and use the maximum clique consistency principle of rigid constraints and DGTG dynamic geometric topology gating to remove outliers and obtain a correspondence set with consistent geometric structure, i.e., the maximum clique structure. Then, filter again and output the fine correspondence set. Step 5: Perform a global attraction domain anchoring operation on the fine correspondence set to obtain the anchoring transformation matrix as the initial value for the next step of probability plane sliding; Step 6: For the largest clique target skeleton node in the fine correspondence set, construct a virtual sliding plane. By minimizing the probability sliding loss function, drive the largest clique source node to slide along the tangent plane to the position that coincides with the target node, thereby achieving registration between the source point cloud and the target point cloud.

2. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 1, characterized in that: In step two, the grid noise reduction uses voxel filtering technology, and the specific method is as follows: 1) Set the voxel side length For the source point cloud floating-point coordinate matrix and the target point floating-point coordinate matrix Spatial downsampling is performed separately to unify data density, where, Dimensions , Dimensions , The total number of points in the source point cloud. The total number of points in the target point cloud; Set any point in the coordinate matrix Its spatial components in three-dimensional coordinates are Calculate the three-dimensional index vector of its voxels. : ； 2) Clustering: Retrievals falling within the same index of For each point, calculate its geometric centroid as the unique skeleton node. : ； In the formula, This indicates the first voxel in the source or target point cloud that falls within the current voxel grid. One point, ; 3) Obtain the coordinate matrix of the source skeleton after dimensionality reduction. Dimension and the target skeleton coordinate matrix Dimension ,in, These represent the number of skeleton nodes in the source point cloud and the target point cloud after dimensionality reduction, respectively.

3. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 1, characterized in that: The method for generating the FPFH feature matrix in step two is as follows: Step 1: Obtain the normal vector; 1) Use KD-Tree to search for the first skeleton coordinate matrix. Nodes of Calculate the local centroid from each neighboring node. And construct the covariance matrix : ； In the formula, Represents a node The 1 neighboring node, ; This means The geometric centroid of each neighboring node; 2) Perform eigenvalue decomposition: For covariance matrix Perform orthogonal decomposition and solve the equations Three eigenvalues ​​were obtained. and its corresponding eigenvectors ,and ; in, This represents the principal direction with the largest variance in the data distribution, while This represents the direction with the smallest variance in the data distribution; geometrically, the direction perpendicular to the local plane has the smallest data variation, therefore, Define as a node normal vector ; Step 2: Perform relative angle encoding (SPFH); To eliminate the effects of coordinate system rotation, at the nodes Construct a local orthogonal coordinate system for the Dabu frame. ,in, , indicating that the normal vector is axis; Determined by the cross product of the line and the normal. axis; Determined by the right-hand rule axis; The purpose of the Dabu frame is to: construct a dynamic local reference, ensuring rotational invariance of angular features; and calculate the angular deviation of neighboring normal vectors within this coordinate system. Generate a 33-dimensional histogram relative angle code. ; Calculate the angular deviation of the neighbor normal vector in this coordinate system The specific calculation formula is as follows: ； ； ； In the formula, Represents the normal vector of the neighbor node and Projection angle deviation between axes; Represents the normal vector of the center node and , The angle of inclination between the vectors connecting two points; Represents the normal vector of the neighbor node exist Orientation rotation angle deviation on the plane; The values ​​of the three angular deviations mentioned above are divided into 11 equally wide intervals for projection counting; Step 3: Perform region convolution fusion FPFH, and concatenate the three 11-dimensional sub-histograms to construct the final 33-dimensional SPFH feature vector. Perform distance-weighted convolution to transform the nodes The relative angle encoding SPFH and its neighboring nodes The relative angle encoding SPFH fusion is performed using the following formula: ； In the formula, Represents the merged nodes The final 33-dimensional feature vector; For nodes Relative angle encoding; For nodes with neighboring nodes Point-to-point distance, For nodes Relative angle encoding; Step 4: Construct the feature matrix and its logical expression: For the source skeleton coordinate matrix All nodes Each node is traversed through steps one through three to obtain the corresponding 33-dimensional FPFH feature vector; This The vectors are stacked in row-matrix order in memory according to their index, thereby constructing the source feature matrix. ; Similarly, for the target skeleton coordinate matrix By traversing and stacking, the target feature matrix is ​​constructed. ; in, ; Each row in the above matrix structure represents the geometric fingerprint of the corresponding node: Dimensions 1 to 11: Characterize whether the local surface is close to a planar structure. The higher the value in this range, the smaller the deviation of the normal direction in the neighborhood. Dimensions 12 to 22: Characterize whether a local region exhibits a cylindrical or rod-like structure. The higher the value in this range, the more significant the curvature change is in the direction of the line connecting the points. Dimensions 23 to 33: Characterize whether a local region contains right angles or obvious edge structures. The higher the value in this range, the closer the angle between the neighborhood normals is to orthogonal. At this point, the two sets of source point clouds and target point clouds, which only have coordinates, have been transformed into two sets of enhanced data structures that have both skeleton coordinates and geometric features.

4. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 3, characterized in that: The specific method for establishing the preliminary logical connection between the source point cloud and the target point cloud in step three is as follows: 1) The source feature matrix Dimension Target feature matrix Dimension ; Using the source skeleton coordinate matrix With the target skeleton coordinate matrix As an auxiliary index, it is used to subsequently restore the feature index to the physical coordinate index; 2) High-dimensional index construction: In order to solve the problem To address the inefficiency of linear search within a 33-dimensional vector, a high-dimensional index is constructed, as follows: target feature matrix Considered a For a high-dimensional dataset of data points, an index structure for a 33-dimensional feature space is constructed using the FLANN fast nearest neighbor search algorithm or a high-dimensional KD-Tree. This is used to divide the 33-dimensional feature space into several hyperplane regions, so that similar geometric fingerprints are adjacent to each other in memory address, thereby accelerating the retrieval. 3) Asymmetric nearest neighbor search Traversing the source feature matrix For each row of feature vectors in the dataset, use it as a query vector and execute a one-way matching strategy. The search algorithm finds the nearest target feature vector and initially determines that the two feature vectors are identical, belonging to the same node or source point. The distance metric used is... Norm, to obtain the best matching target index found in the search. Establish a preliminary logical connection between the source point cloud and the target point cloud; ； In the formula, The source characteristic matrix is ​​the first The eigenvectors of the row, ; For the target feature matrix and The nearest row feature vector, ; 4) Result mapping and set generation: For each source index traversed Record the best target index found in the search. This forms an index-correspondence tuple, or simply an index pair. ; Then put all The pairs converge to form the initial rough correspondence set. : ; At this time It is An integer matrix, where each row represents a pair of predicted origin points, and the total number of rows in the set is . , .

5. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 4, characterized in that: The specific steps for step four are as follows: 1) Geometric pruning: Through the initial rough correspondence set Any two points in and As a matching pair, the coordinates are obtained from the source skeleton matrix via indexing. Read source node and The coordinates, from the target skeleton coordinate matrix Read target node and The coordinates of the source node are calculated. and distance Calculate and obtain the target node and distance ; The absolute value of the difference between the lengths of the two segments is calculated using the geometric residual formula. : ； Set the tolerance threshold for geometric consistency determination , Then determine and These two matching pairs are geometrically compatible, preserving logical connections; If it is not a match, it is considered an incorrect match and is cut off. and The logical connection between these two matching pairs completes the geometric pruning; 2) Construct a compatibility matrix: Initial rough correspondence set All possible matching pair combinations, total Group, construct a dimension as Correlation matrix Iterate through all matching pairs and fill the matrix elements according to the results of the geometric pruning described above. : ； 3) Further filtering to generate a refined correspondence set. : For matrix Perform eigenvalue decomposition or exponentiation to obtain the principal eigenvector corresponding to the largest eigenvalue. , The dimension is ; vector elements in For the first The scalar value of the row also represents the first row. The geometric consistency scores of each matching pair are sorted from highest to lowest, and the top scores are retained. Each matching pair generates a fine correspondence set. Fine correspondence set This constitutes the maximal clique in geometric graph theory, meaning that all preserved point pairs are rigidly consistent in spatial structure.

6. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 5, characterized in that: The specific steps for step five are as follows: 1) Traverse the fine correspondence set Each index pair in, are the original physical coordinates of the source node and the target node in the c-th matching pair, respectively, and both are... Column vectors From the source skeleton coordinate matrix Obtained from the middle, From the target skeleton coordinate matrix Obtained from the middle; 2) Geometric centroid calculation: The mean of all extracted node coordinate vectors is calculated to obtain the physical geometric centroids of the largest source skeleton and the target skeleton. ； In the formula, The physical geometric centroid of the largest source skeleton; The physical geometric centroid of the largest cluster target skeleton. The total number of matching pairs in the largest clique; For the iterative index of the summation operation, ; 3) Decentralized mapping: Perform explicit vector subtraction to translate all points to a local coordinate system with the centroid as the origin, generating relative coordinate vectors: ； In the formula, Let be the relative coordinate vector of the source node in the c-th matching pair after the decentralization of the maximum clinker source skeleton. The relative coordinate vector of the target node in the c-th matching pair after the decentralization of the target skeleton of the maximum clique target; 4) Relative coordinate set encapsulation: Stack all the obtained relative coordinate vectors column-wise to construct a relative coordinate set matrix. and : ； In the formula, The matrix of relative coordinates of the largest clique source skeleton, dimension ; The relative coordinate set matrix of the skeleton of the maximal clique target, dimension ; 5) For the relative coordinate set matrix and Perform global attraction domain anchoring rotation analysis: Covariance matrix construction: Calculate the point cloud reflecting the geometric correlation between the source skeleton point cloud and the target skeleton point cloud of the largest clique using relative coordinate sets. Attraction domain cross-covariance matrix : ； in, This represents the matrix transpose operation; right Perform singular value decomposition (SVD) to decompose it into three orthogonal components: ； In the formula, Let represent the left singular vector matrix, the singular value diagonal matrix, and the transpose of the right singular vector matrix obtained from SVD decomposition, respectively. Closed-form anchoring of rotation matrix: Based on the principle of orthogonal Pluke analysis, the optimal rotation matrix can be directly derived analytically. : ； To introduce a determinant term, used to force a correction to the reflection matrix generated by SVD, ensuring Satisfy rigid body rotation constraints in a right-handed coordinate system; 6) For rotation matrices The physical geometric centroids of the maximum source skeleton and the target skeleton Perform translation calculations and matrix assembly: Back-substitution for translation: using the transformed... and Given overlapping constraints, solve for the global anchoring translation vector. : ； Final matrix assembly: The rotation matrix and translation vector are combined into a standard homogeneous transformation matrix as the anchoring transformation matrix. : ； Will As the initial rigid transformation command.

7. The low overlap rate point cloud registration method based on topology-metric decoupling according to claim 6, characterized in that: The specific steps for step six are as follows: 1) Construct a virtual sliding plane; Fine correspondence set The Middle The original physical coordinates of the source node and the target node in the matching pair are respectively used as follows: This indicates that the corresponding coordinates are also used to represent the target node. Feature vectors are extracted from the neighborhood point set to construct a local orthogonal basis. ,in, For the target node Surface normal vector at that location; Target nodes Two orthogonal tangent vectors on the tangent plane; Constructing anisotropic covariance matrices using the spectral decomposition theorem A planar structure is constructed, and extremely large tangential uncertainty and extremely small normal uncertainty are artificially set: ； In the formula, Let the variances of the uncertainties along the normal and tangential directions be defined respectively. ; 2) Solve the probabilistic sliding optimization using Lie algebra: Constructing the sliding loss function , The sum of Mahalanobis distances for matching pairs in all fine correspondence sets is used with the goal of minimizing the loss. ； in, The total number of matching pairs in the largest clique; Let be the position residual vector. ; because variance in the tangential direction The maximum value is its inverse matrix. Since the eigenvalue weights in the tangential direction approach 0, the tangential residual is eliminated in the position residual vector formula, and only the normal residual term is retained; thus, each source node gains the degree of freedom to slide freely on the target surface until the normal error is reduced to zero. 3) Dynamic evolution is performed through Gauss-Newton iteration, and the optimal perturbation vector that minimizes the loss is solved in the Lie algebra space. ; For each matching pair, calculate the partial derivative of the position residual with respect to the perturbation vector, i.e., the Jacobian matrix. : ； in, Indicates the first A number of matching pairs Jacobian matrix; Let represent a 6-dimensional perturbation vector in the Lie algebra space, which is used as an independent variable in the differentiation here, and is also the pose update increment to be solved in this iteration; express identity matrix; Represents the coordinates of the skeleton nodes in the source point cloud; This indicates an antisymmetric operator, used here to transform a 3D vector Mapped to Antisymmetric matrix; This indicates that the current rigid body transformation matrix is ​​applied to the source point cloud skeleton node. The resulting transformed coordinate vector ; This represents the rotation matrix to be solved, used to describe the pose changes of the source point cloud; This represents the translation vector to be solved, used to describe the spatial displacement of the source point cloud; Construct and solve the linear increment equation to obtain the optimal Lie algebra update. : ； in, This represents the inverse matrix of the anisotropic covariance matrix constructed for the corresponding target skeleton node, which is used to weight the errors in different directions during optimization. 4) Pose manifold update: The update value obtained by using exponential mapping Restore to a rigid body transformation matrix and apply to the current pose. Update to obtain a new pose. : ； in, These represent the updated target pose matrix and the current pose matrix before the update, respectively; This represents the matrix exponential function used to map the Lie algebra increments in the tangent space back to the Lie group manifold.

Images