A 2-point ransac algorithm with additional multiple constraints

Abstract

The application discloses a 2-point RANSAC (random sample consensus) algorithm with additional multiple constraints. In the point cloud registration based on local shape descriptors, the pairing of two point clouds is established through the descriptors, but there are a large number of false pairs in the pairing, and the 2-point RANSAC algorithm with additional multiple constraints can exclude the influence of the false pairs and calculate accurate conversion parameters. The algorithm randomly selects two pairs in each iteration, and if the two pairs meet the multiple constraint conditions, the conversion parameters are calculated by using the two pairs, and the number of inliers is calculated, and if the multiple constraint conditions are not met, the next iteration is jumped to. Finally, the conversion parameters corresponding to the maximum number of inliers are the best conversion parameters.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of three-dimensional point cloud data processing, in particular to a 2-point RANSAC algorithm with additional multiple constraints. BACKGROUND

[0002] In many practical applications, such as scene reconstruction, 3D object recognition, position recognition, point cloud registration plays a crucial role, and therefore receives extensive attention. The purpose of point cloud registration is to transform point clouds in different coordinate systems into the same coordinate system to obtain complete point clouds. In the field of computer vision and computer graphics, there have been a large number of studies on point cloud registration based on local shape descriptors. Because the matching pairs established by local shape descriptors contain a large number of incorrect matching pairs, it is a challenging task to calculate accurate transformation parameters in this case. In view of this problem, a large number of algorithms have been proposed.

[0003] Existing algorithms can be divided into two categories: RANSAC (random sample consensus) type methods and other algorithms. Other algorithms include geometric constraint cluster (GCC), game theoretic matching (GTM), local and global voting (LGV), etc. GCC uses geometric distance constraints to find matching pairs that satisfy the constraint conditions as correct matching pairs. GTM also uses geometric distance constraints to find matching pairs that are compatible with each other, and excludes matching pairs that are not compatible. LGV adopts two stages to find correct matching pairs. In the local voting stage, a compatibility score is calculated for each matching pair, and matching pairs with high compatibility scores are extracted and input into the global voting stage. In the global voting stage, a transformation parameter is calculated for each extracted matching pair, and then global voting is calculated. This type of method still has poor accuracy, and many incorrect matching pairs are still selected and retained, which leads to the inability to calculate accurate transformation parameters.

[0004] Many scholars improve the RANSAC algorithm and propose many improved methods, which are collectively referred to as RANSAC-like algorithms. The RANSAC algorithm randomly selects three pairs each time, calculates the transformation parameters, and then converts the original point cloud to the target point cloud according to the transformation parameters, and calculates the number of points in the overlapping area of the two point clouds as the confidence score. The 1-point RANSAC algorithm randomly selects one pair, calculates the transformation parameters, and because the rotation matrix is calculated using the local reference frame (LRF) of the key point, only one pair is enough to calculate the transformation parameters. The 2SAC-GC (2-point based sample consensus with global constraint) algorithm randomly selects two pairs to calculate the transformation parameters, and uses a geometric distance constraint and an angle constraint to reject incorrect pairs and speed up the calculation. The algorithm uses the normal vector of the key point to calculate the rotation matrix, so two pairs are enough to calculate the transformation parameters. This type of algorithm is greatly affected by the proportion of incorrect pairs. When the proportion of incorrect pairs is large, more iterations are required, resulting in a large amount of calculation time, so the RANSAC-like algorithm has poor calculation efficiency. SUMMARY

[0005] The purpose of the present application is to overcome the shortcomings of the prior art and provide a 2-point RANSAC algorithm with additional multiple constraints for point cloud registration based on local shape descriptors. By adding multiple constraint conditions to reject incorrect pairs during the random sampling process, when the randomly selected two pairs do not meet the constraint conditions, the current iteration is skipped, the accurate transformation parameters can be calculated, and the calculation efficiency is good.

[0006] To achieve the above purpose, the technical scheme is as follows.

[0007] A 2-point RANSAC algorithm with additional multiple constraints for point cloud registration, after obtaining an initial pair set, the original point cloud is converted to the target point cloud by the transformation parameters, and automatic point cloud registration is realized, comprising the following steps:

[0008] Step S1, randomly select two pairs from the initial pair set;

[0009] Step S2, if the selected two pairs meet the multiple constraint conditions, calculate the transformation parameters;

[0010] Step S3, according to the transformation parameters calculated in step S2, convert the original point cloud to the target point cloud, and calculate the number of points in the overlapping area of the two point clouds as the number of inliers;

[0011] Step S4: Repeat steps S1 to S3 until the maximum number of iterations is reached. Output the transformation parameters corresponding to the maximum number of interior points as the optimal transformation.

[0012] In step S1, the initial pairing set is obtained through local shape description. The specific process is as follows:

[0013] First, key points are extracted from the two point clouds to be registered. This can be done using ISS key point detectors, 3D Harris key point detectors, etc. Then, local shape descriptors of the key points are calculated, such as TOLDI descriptors, LoVS descriptors, TriLCI descriptors, etc. Finally, the similarity of the descriptors is compared using the nearest neighbor similarity ratio (NNSR) to establish a pair between the two point clouds and obtain an initial pairing set.

[0014] The multiple constraints mentioned in step S2 include:

[0015] a. Geometric distance constraints

[0016] Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, then p i and p j The distance between them is equal to q i and q j The distance between them, therefore the geometric distance constraint is defined as

[0017] |d(p i ,p j )-d(q i ,q j )|<δ1 (1)

[0018] In the above formula, |·| represents taking the absolute value, d represents the distance between two key points, δ1 represents the threshold, which is generally set to 2pr, and pr represents the point cloud resolution;

[0019] b. Local reference axis angle constraint

[0020] Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and qj are two key points of the original point cloud; if the two pairs satisfy multiple constraint conditions, then the normal included angle between p i and p j is equal to the normal included angle between q i and q j , so the local reference frame angle constraint is defined as

[0021]

[0022] In the above formula, arccos(·) represents the inverse cosine function, LRA represents the local reference axis of the key point, and δ2 represents a threshold value, which is generally set to 10 degrees;

[0023] c. Local reference frame angle constraint

[0024] Two pairs (p i , q i ) and (p j , q j ) are given, p i and p j are two key points of the original point cloud, and q i and q j are two key points of the target point cloud; if the two pairs satisfy multiple constraint conditions, then the local reference frame included angle between p i and p j is equal to the local reference frame included angle between q i and q j , and the local reference frame angle constraint is defined as

[0025]

[0026] In the above formula, trace(·) represents the trace of a matrix, V represents the local reference frame of the key point, and δ3 is a threshold value, which is generally set to 10 degrees;

[0027] d. Normal projection distance constraint

[0028] Two pairs (p i , q i ) and (p j , q j ) are given, p i and p j are two key points of the original point cloud, and q i and q j are two key points of the target point cloud; if the two pairs satisfy multiple constraint conditions, then the projection distance from one key point to the normal vector of the other key point is fixed, so for the two pairs (p i , q i ) and (p j,q j ), calculate its four projected distances

[0029]

[0030] In the above formula, n represents the normal vector of the key point;

[0031] Therefore, the normal projection distance constraint is defined as

[0032] |d1-d3|<δ4,|d2-d4|<δ4 (5)

[0033] In the formula, δ4 is a threshold, which is usually set to 2pr.

[0034] The calculation process for the transformation parameters (including the rotation matrix and translation vector) mentioned in step S2 is as follows:

[0035] Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, four points are calculated based on the key points and their normal vectors as follows:

[0036]

[0037] From point Construct the following matrix

[0038]

[0039] Calculate the covariance matrix Σ

[0040]

[0041] In the above formula, It is the average of A. It is the average of Y, and ones(4,1) represents a 4x1 vector with all elements equal to 1; by performing singular value decomposition on the covariance matrix, we obtain

[0042] USV T =svd(Σ) (9)

[0043] The rotation matrix is ​​then calculated as follows:

[0044] R = VDU T (10)

[0045] In the above formula, D = diag(1,1,det(UV)T diag(·) represents a diagonal matrix, and det(·) represents taking the determinant of a matrix;

[0046] The translation vector is calculated as

[0047]

[0048] In the above formula, is the average of q i and q j is the average of p i and p j .

[0049] Step S3: According to the conversion parameters calculated in step S2, the original point cloud is converted to the target point cloud, and the number of points in the overlapping area of the two point clouds is calculated as the number of inliers, and the specific calculation process is as follows:

[0050] Step S31: According to the calculated conversion parameters, the original point cloud is converted to the target point cloud, and for each converted point, find its nearest point in the target point cloud, if the distance between this point and its nearest point is less than 3pr, then this point is considered to be an inlier, and the total number of inliers is counted.

[0051] Compared with the prior art, the beneficial effects of the present application are:

[0052] 1. The method of the present application uses multiple constraint conditions to determine whether the currently selected pair is a correct pair, and if it is determined to be an incorrect pair, the current iteration is skipped, which can reduce the calculation time and the proposed algorithm has high calculation efficiency.

[0053] 2. Experimental analysis shows that: compared with the 1-point RANSAC and OSAC algorithms, the method of the present application has significantly better calculation efficiency and better calculation accuracy, and can obtain accurate conversion parameters. BRIEF DESCRIPTION OF DRAWINGS

[0054] Figure 1 is a flowchart of the 2-point RANSAC algorithm with additional multiple constraints of the present application;

[0055] Figure 2 is the registration result of the two pairs of point clouds. DETAILED DESCRIPTION

[0056] ​With reference to the drawings and the embodiments of the present application, the technical solutions in the embodiments of the present application will be described clearly and completely. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all the embodiments of the present application. Based on the embodiments of the present application, all other embodiments obtained by those skilled in the art without creative efforts belong to the scope of the present application.

[0057] Embodiment 1

[0058] As shown in the following, in this embodiment, a 2-point RANSAC algorithm with additional multiple constraints is provided for point cloud registration. After obtaining an initial matching set, the original point cloud is converted to the target point cloud by conversion parameters, and automatic point cloud registration is realized, including the following steps: Figure 1

[0059] Step S1, two matches are randomly selected from the initial matching set;

[0060] Step S2, if the selected two matches satisfy the multiple constraint conditions, the conversion parameters are calculated;

[0061] Step S3, according to the conversion parameters calculated in step S2, the original point cloud is converted to the target point cloud, and the number of points in the overlapping area of the two point clouds is calculated as the number of inliers;

[0062] Step S4, steps S1-S3 are repeated until the maximum number of iterations is reached, and the conversion parameters corresponding to the maximum number of inliers are output as the best conversion.

[0063] In step S1, the initial matching set is obtained by local shape description, and the specific process is as follows:

[0064] First, key points are extracted from two point clouds to be registered, which can use ISS key point detector, 3D Harris key point detector, etc., then the local shape descriptor of the key points is calculated, such as TOLDI descriptor, LoVS descriptor, TriLCI descriptor, etc., finally, the similarity of the descriptors is compared by using the nearest neighbor similarity ratio (NNSR), and the matching between the two point clouds is established to obtain the initial matching set.

[0065] The multiple constraint conditions in step S2 include:

[0066] a. Geometric distance constraint

[0067] Two matches (p i ,q i ) and (p j ,q j ) are given, p i and p j are two key points of the original point cloud, q i and q​j are two key points of the original point cloud; if the two pairs satisfy multiple constraint conditions, the normal angle between p i and p j is equal to the normal angle between q i and q j , so the geometric distance constraint is defined as

[0068] |d(p i ,p j )-d(q i ,q j )|<δ1 (1)

[0069] In the above formula, |·| represents taking the absolute value, d represents the distance between two key points, δ1 represents the threshold value, which is generally set to 2pr, pr represents the point cloud resolution;

[0070] b. Local reference axis angle constraint

[0071] Given two pairs (p i , q i ) and (p j , q j ), p i and p j are two key points of the original point cloud, q i and q j are two key points of the target point cloud; if the two pairs satisfy multiple constraint conditions, the normal angle between p i and p j is equal to the normal angle between q i and q j , so the local reference frame angle constraint is defined as

[0072]

[0073] In the above formula, arccos(·) represents the inverse cosine function, LRA represents the local reference axis of the key point, and δ2 represents the threshold value, which is generally set to 10 degrees;

[0074] c. Local reference frame angle constraint

[0075] Given two pairs (p i , q i ) and (p j , q j ), p i and p j are two key points of the original point cloud, q i and q j are two key points of the target point cloud; if the two pairs satisfy multiple constraint conditions, the normal angle between p i and p jThe local reference frame included angle between them is equal to q i and q j The local reference frame included angle between them, the local reference frame angle constraint is defined as

[0076]

[0077] In the above formula, trace(·) represents the trace of a matrix, V represents the local reference frame of the key point, δ3 is a threshold value, generally set to 10 degrees;

[0078] d. Normal projection distance constraint

[0079] Two pairs (p i , q i ) and (p j , q j ) are given, p i and p j are two key points of the original point cloud, and q i and q j are two key points of the target point cloud; if the two pairs satisfy the multiple constraint conditions, the projection distance on the normal vector from one key point to another key point is fixed, so for the two pairs (p i , q i ) and (p j , q j ) given, the four projection distances d1, d2, d3 and d4 are calculated

[0080]

[0081] In the above formula, n represents the normal vector of the key point;

[0082] Therefore, the normal projection distance constraint is defined as

[0083] |d1-d3|<δ4,|d2-d4|<δ4 (5)

[0084] In the formula, δ4 is a threshold value, generally set to 2pr.

[0085] The calculation process of the conversion parameters (including the rotation matrix and the translation vector) in step S2 is as follows:

[0086] Two pairs (p i , q i ) and (p j , q j ) are given, p i and p j are two key points of the original point cloud, and q i and q jare two key points of the target point cloud; if two pairs satisfy multiple constraint conditions, four points are calculated according to the key points and their normal vectors as follows

[0087]

[0088] The points are constructed as follows

[0089]

[0090] The covariance matrix Σ is calculated

[0091]

[0092] In the above formula, is the average of A, is the average of Y, and ones(4, 1) represents a 4-row 1-column vector with all elements being 1; singular value decomposition is performed on the covariance matrix to obtain

[0093] USV T = svd(Σ) (9)

[0094] Then the rotation matrix is calculated as

[0095] R = VDU T (10)

[0096] In the above formula, D = diag(1, 1, det(UV T )), diag(·) represents a diagonal matrix, and det(·) represents taking the determinant of a matrix;

[0097] The translation vector is calculated as

[0098]

[0099] In the above formula, is the average of q i and q j , is the average of p i and p j .

[0100] Step S3: According to the conversion parameters calculated in step S2, the original point cloud is converted to the target point cloud, and the number of points in the overlapping region of the two point clouds is calculated as the number of inliers, and the specific calculation process is as follows:

[0101] Step S31: According to the calculated conversion parameters, the original point cloud is converted to the target point cloud, and for each converted point, its nearest point in the target point cloud is found, and if the distance from this point to its nearest point is less than 3pr, it is considered that this point is an inlier, and the total number of inliers is counted.

[0102] Example 2

[0103] In the following, the algorithm provided by the present application is compared with the 1-point RANSAC and OSAC algorithms through a specific experiment to verify the superiority of the algorithm of the present application.

[0104] The experimental process of the present embodiment is as follows: Apartment and City scene point clouds are used to implement the experiment, a pair of point clouds of the Apartment scene is obtained by a FARO Focus 3D X330 HDR scanner, and a pair of point clouds of the City scene is obtained by a Leica C10 laser scanner. For each pair of point clouds, a 3D Harris detector is used to extract key points, a TriLCI descriptor is calculated, and then the similarity of the descriptors is compared by using the NNSR algorithm to establish the pairing between the two point clouds. Then, the conversion parameters are calculated by using the method of the present application, the original point cloud is converted to the target point cloud according to the conversion parameters, and the point cloud registration is realized.

[0105] The algorithm of the present application is compared with the 1-point RANSAC and OSAC algorithms, the number of iterations of the 1-point RANSAC is set to 500 times, the number of iterations of the OSAC algorithm is set to 2000 times, the effect of the algorithm is evaluated from the registration accuracy and the calculation time, the registration accuracy is measured by the rotation error and the translation error, and the calculation is as follows

[0106]

[0107] error t ||t true -t|| (13)

[0108] In the above formula, R represents the estimated rotation matrix, R true represents the true value of the rotation matrix; t represents the estimated translation vector, t true represents the true value of the translation parameter.

[0109] The true rotation matrix and translation vector of the two pairs of point clouds are calculated by artificial coarse registration and ICP fine registration. Each of the three algorithms is run 100 times, the average rotation error, the average translation error and the average calculation time are calculated, and the experimental results are shown in Table 1 below.

[0110] Table 1 Registration error and calculation time of three algorithms

[0111]

[0112]

[0113] As can be seen from Table 1, the average time of the 2-point RANSAC algorithm proposed in the present application is obviously much less than that of the 1-point RANSAC and OSAC algorithms, indicating that the algorithm proposed in the present application has high calculation efficiency. In the Apartment scene, the OSAC algorithm often fails to register, resulting in large average rotation error and translation error.

[0114] As shown in Figure 2 The registration results of the two pairs of point clouds are shown, indicating that the two pairs of point clouds have been well registered together, which is sufficient to provide a good initial pose to the fine registration algorithm, avoiding the fine registration from falling into a local optimal solution.

[0115] The above describes the preferred implementation method of the present application, but the present application is not limited to the embodiments described, and those skilled in the art can make various equivalent modifications or replacements without departing from the spirit of the present application. These equivalent modifications or replacements are all included in the scope defined by the claims of the present application.

Claims

1. A 2-point RANSAC method with additional multiple constraints for point cloud registration, which, after obtaining an initial pairing set, transforms the original point cloud to the target point cloud through transformation parameters to achieve automatic point cloud registration, characterized in that... Includes the following steps: Step S1: Obtain the initial pairing set from the initial pairing set through local shape description, and randomly select two pairs; Step S2: If the two selected pairs satisfy multiple constraints, then calculate the transformation parameters; Multiple constraints include: a. Geometric distance constraints Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, then p i and p j The distance between them is equal to q i and q j The distance between them, therefore the geometric distance constraint is defined as |d(p i ,p j )-d(q i ,q j )|<δ1 (1) In the above formula, |·| represents taking the absolute value, d represents the distance between two key points, δ1 represents the threshold, set to 2pr, and pr represents the point cloud resolution; b. Local reference axis angle constraint Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, then p i and p j The normal angle between them is equal to q i and q j The angle between the normals, therefore the local reference frame angle constraint is defined as... In the above formula, arccos(·) represents the inverse cosine function, LRA represents the local reference axis of the key point, and δ2 represents the threshold, which is set to 10 degrees. c. Local reference frame angle constraints Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, then p i and p j The angle between the local reference frames is equal to q. i and q j The angle between the local reference frames is defined as follows: In the above formula, trace(·) means to find the trace of the matrix, V represents the local reference frame of the key point, and δ3 is a threshold, set to 10 degrees. d. Normal projection distance constraint Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, the projection distance on the normal vector from one key point to the other key point is constant. Therefore, for the given two pairs (p... i ,q i ) and (p j ,q j ), calculate its four projected distances In the above formula, n represents the normal vector of the key point; Therefore, the normal projection distance constraint is defined as |d1-d3|<δ4,|d2-d4|<δ4 (5) In the formula, δ4 is a threshold value, set as 2pr; Step S3: Based on the conversion parameters calculated in step S2, convert the original point cloud to the target point cloud, and calculate the number of points in the overlapping area of ​​the two point clouds as the number of internal points; Step S4: Repeat steps S1 to S3 until the maximum number of iterations is reached. Output the transformation parameters corresponding to the maximum number of interior points as the optimal transformation.

2. The 2-point RANSAC method with additional multiple constraints according to claim 1, characterized in that, In step S1, the initial pairing set is obtained through local shape description. The specific process is as follows: First, key points are extracted from two point clouds to be registered using the ISS key point detector and the 3D Harris key point detector. Then, the local shape descriptors of the key points are calculated. Finally, the nearest neighbor similarity ratio (NNSR) is used to compare the similarity of the descriptors and establish a pairing between the two point clouds to obtain the initial pairing set.

3. The 2-point RANSAC method with additional multiple constraints according to claim 1, characterized in that, The transformation parameters mentioned in step S2 include the calculation process of the rotation matrix and translation vector, which is as follows: Given two pairs (p) i ,q i ) and (p j ,q j ), p i and p j These are the two key points of the origin cloud, q i and q j These are two key points in the target point cloud; if the two pairs satisfy multiple constraints, four points are calculated based on the key points and their normal vectors as follows: From point Construct the following matrix Calculate the covariance matrix Σ In the above formula, It is the average of A. It is the average of Y, and ones(4,1) represents a 4x1 vector with all elements equal to 1; by performing singular value decomposition on the covariance matrix, we obtain USV T =svd(Σ) (9) The rotation matrix is ​​then calculated as follows: R=VDU T (10) In the above formula, D = diag(1,1,det(UV) T )), diag( · ) represents a diagonal matrix, and det(·) represents taking the determinant of the matrix; The translation vector is calculated as follows In the above formula, It is q i and q j On average, It is p i and p j The average.

4. The 2-point RANSAC method with additional multiple constraints according to claim 1, characterized in that, Step S3 involves converting the original point cloud to the target point cloud based on the conversion parameters calculated in step S2, and calculating the number of points in the overlapping area of ​​the two point clouds as the number of interior points. The specific calculation process is as follows: Step S31: Based on the calculated transformation parameters, transform the original point cloud to the target point cloud. For each transformed point, find its nearest point in the target point cloud. If the distance from this point to its nearest point is less than 3pr, then this point is considered an interior point. Count all interior points.

Images