A New Method for Registration of 3D Point Sets with Low Overlapping Ratios

Available online at www.sciencedirect.com

ScienceDirect Procedia CIRP 27 (2015) 202 – 206

13th CIRP conference on Computer Aided Tolerancing

A new method for registration of 3D point sets with low overlapping ratios Y. F. Wua, W. Wang*b, K. Q. Lub, Y. D. Weia, Z. C. Chena a

Department of Mechanical Engineering, Zhejiang University, Hangzhou 310027, China;

b

School of Mechanical Engineering, Hangzhou Dianzi University, Hangzhou 310018, China.

* Corresponding author. Tel.: +86-571-8691-9155; fax: +86-571-8691-9155. E-mail address: [email protected]

Abstract The registration of multi-view point sets is often used in surface reconstruction for complex parts, and many researchers are striving to improve its accuracy and efficiency. A popular algorithm for registration of 3D point sets is the Iterative Closest Point (ICP). With the increase of overlapping ratios, the convergence accuracy and speed of ICP algorithm improves. In most cases of 3D scanning systems, the overlapping regions of multi-view point sets often locate at the boundaries. Based on this characteristic, a new method for registration of multi-view point sets with low overlapping ratios is proposed. Firstly, two kinds of sampling strategies are introduced, and a great number of point pairs are constructed. Through surface classification, curvature and neighborhood curvature matching of point pairs, the overlapping regions are selected. Then a Levenberg-Marquardt ICP (LM-ICP) method is used to search optimal registration parameters with the points in the overlapping regions. Finally, experiments for registration of multi-view point sets are performed, and results show that the method can get good registration accuracy for point sets with low overlapping ratios. © 2015 2015 The The Authors. Authors. Published Published by by Elsevier Elsevier B.V B.V.This is an open access article under the CC BY-NC-ND license © Peer-review under responsibility of the organizing committee of 13th CIRP conference on Computer Aided Tolerancing. (http://creativecommons.org/licenses/by-nc-nd/4.0/). Peer-review under responsibility of the organizing committee of 13th CIRP conference on Computer Aided Tolerancing Keywords: registration; 3D point sets; low overlapping ratio; curvature matching; LM-ICP

1. Introduction The registration of 3D data sets is an important work for geometrical reconstruction in reverse engineering. It is also widely applied in depth data fusion [1], medical image registration [2], computer vision [3], mold manufacture and other fields. There are many different methods for surface digitization in reverse engineering, which can be mainly divided into two groups: contact and non-contact [4]. For the object with large size and complex shape, due to the limitation of scan range and the presence of occlusions [5], it is essential to measure multi-view point sets to cover the entire object and register them to the same coordinate system. The registration of point sets is typically performed pairwise, and a well-known algorithm to solve the registration problem is the Iterative Closest Point (ICP) algorithm proposed by Besl and McKay [6]. ICP is a non-linear local search algorithm, and its robustness is influenced by initial values and search strategies. Many researches involving sampling initial points, searching closest point, rejecting

wrong point pairs and searching optimal transformation parameters are reviewed in the literatures [7-10]. The ICP algorithm is suitable for handling point sets with high overlapping ratios. When overlapping regions occupy a small part of point sets, results become less accurate. In practice, a 50% overlapping ratio of the region is critical [7]. The increase of overlapping ratio can help to improve the registration accuracy of two point sets [11]. However, the total number of multi-view point sets increases, which may bring a large amount of computation and global registration accuracy problems. In industrial applications, graphic markers such as concentric circles, cross lines, ect., are widely used to ensure the registration efficiency and accuracy of point sets with low overlapping ratios. In practice, for most 3D laser scanning systems, the overlapping regions of multi-view point sets often locate at the boundaries. Based on this characteristic, a new method for registration of multi-view point sets with low overlapping ratios is proposed, and it doesn’t need any auxiliary graphic marker. Section 2 begins by describing two kinds of sampling strategies to gain point pairs which are used to search

2212-8271 © 2015 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). Peer-review under responsibility of the organizing committee of 13th CIRP conference on Computer Aided Tolerancing doi:10.1016/j.procir.2015.04.067

203

Y.F. Wu et al. / Procedia CIRP 27 (2015) 202 – 206

overlapping regions. Then section 3 presents the surface classification, curvature matching, neighborhood curvature matching of point pairs and the construction of overlapping regions. Section 4 briefly illustrates the LM-ICP method used in registration of point sets. As an example, a Stanford bunny is used to examine the key parameters in section 5. Finally, a conclusion is drawn in section 6.

used in Section 3.

2. Sampling strategies

3. Searching overlapping regions

We define A, B the two point sets to be registered. L denotes a loop list containing sample points from A. U denotes the sample points from B. The purpose of sampling points is to construct point pairs {L, U} which are used to search overlapping regions of A. P and Q are the search results of overlapping regions (PA, QB). V denotes the real overlapping regions.

3.1. Surface classification

2.1. Sampling points in A Bounding boxes of point sets A, B are respectively constructed, and the center of each box is estimated. The Principal Component Analysis (PCA) [12] method is used to obtain the orthogonal components of point sets A. As shown in Fig. 1, the local coordinate system is established by using the center O (O is the nearest point from the center of the bounding box) and the principal components. The normal vector n at the center O is estimated by using covariance analysis [13], and the axis OZ can be easily separated from the principal components. A plane (dotted line in Fig. 1) is constructed through axis OZ, then two edge points (hollow points) on the plane are obtained, and the sample points (solid points) are determined as follows: 1) The distance d between O and edge point is calculated, and the sample point is determined with 90% of d from O on the plane; 2) The plane is rotated with a fixed angle, and the remaining sample points can be obtained; 3) All sample points are stored in L.

2.2. Sampling points in B The overlapping region of point sets B is unknown, so it is impossible to sample points in specified regions. Uniform sampling is used to construct point sets U. Plenty of point pairs {L, U} are constructed, which are

According to Table 1, the surface type can be classified by Gaussian curvature and mean curvature. The moving least square (MLS) method [14] is used to estimate the Gaussian curvature and the mean curvature of the local surface. If the surface types of each point pair in {L, U} are different, the point pair will be removed. The plane type influences the results of searching overlapping regions, so point pairs with plane type are also removed. When | K | < 0.01 and | H | < 0.01, K and H are regarded to be zero. Through experimental tests, 3/4 of {L, U} are filtered out by surface classification. If K, H of some points in L is 0, these points should be resampled, and 90% of d is reduced to 85% of d. Table 1. K, H classification

H<0 H=0 H>0

K<0

K=0

K>0

Hyperbolic concave Hyperbolic symmetric Hyperbolic convex

Cylindrical concave

Elliptical concave

Planar

Impossible

Cylindrical convex

Elliptical convex

3.2. Curvature matching Let (p, q) be one of the remaining point pairs in {L, U}. If (p, q) meets Equation (1), (p, q) will be kept in {L, U}. ­ k1 ( p)  k1 (q) / ( k1 ( p)  k1 (q) )  H1 ° (1) ® ° ¯k2 ( p )  k2 ( q ) / ( k2 ( p )  k2 ( q ) )  H 2 Where, k1, k2 are the principle curvature estimated by MLS method [14], and H1, H2 denote the curvature threshold. 3.3. Neighborhood curvature matching Through surface classification and curvature matching, point pairs in {L, U} are drastically reduced. Zero mean Normalized Cross Correlation (ZNCC) is employed to denote the curvature neighborhood matching degree. ZNCC is widely used due to its robustness in template matching, motion analysis, stereo vision and industrial inspection [15]. The expressions of ZNCC are as follows: zncc1 ( p, q)

¦ (k ( p )  k ( p))(k (q )  k (q)) ¦ (k ( p )  k ( p)) ¦ (k (q )  k (q)) ¦ (k ( p )  k ( p))(k (q )  k (q)) ¦ (k ( p )  k ( p)) ¦ (k (q )  k (q))

zncc2 ( p, q)

i

1

1

i

i

1

1

i

1

2

i

i

2

2

1

i

2

i

i

2

2

2

i

2

i

2

(2)

1

2

i

Fig. 1. Sampling points from point sets A

1

i

i

2

2

(3)

204

Y.F. Wu et al. / Procedia CIRP 27 (2015) 202 – 206

(4) zncc1 ( p, q)  zncc2 ( p, q) Where, zncc (p, q) is the comprehensive neighborhood curvature matching degree.Ck1,Ck2 are the mean value of the principle curvature of k-nearest neighborhood points. If some points in L have no corresponding points in U after curvature matching, the neighborhood curvature matching degree is set to zero. The neighborhood curvature matching degree of the remaining point pairs in {L, U} is calculated by using Equation (2) to (4), and the maximum value for each point in L is selected and stored orderly in a loop M. zncc( p, q)

3.4. Constructing overlapping regions Loop M is finally obtained. Values are extracted orderly from M with a fixed number. Chains of values are constructed, and the chain with the largest sum is used to search P. P (the dotted line in Fig. 2) is constructed by using the starting point, the ending point of and the coordinate system O-XYZ. The overlapping region of point sets B, namely, Q, is also obtained by using the method above.

break else k = k +1 end for The procedures of C1, C2 and C3 are implemented by the LM-ICP method proposed by Fitzgibbon [16]. This method replaces the Euclidean distance with the Chamfer distance and uses a Levenberg-Marquartd algorithm to compute Tk, which accelerates the convergence of data registration, while keeping high accuracy. Especially in the treatment of high overlapping ratios, the LM-ICP method is superior to the standard ICP method. 5. Experiments Shown in Fig. 3 are two pieces of point sets. The number of each point sets is 14,195 and 12,151 with overlapping ratios of 37.77% and 44.13% respectively. Searching overlapping regions and registering point sets were performed on MATLAB 2007b with the computer configuration of Core 2 2.0 GHz CPU and 2 GB memory.

Matching & Filtering Uniform sampling

O

Point sets A Overlapping region of point sets A, namely P

Point

sets

Fig. 2. Searching overlapping region of point sets A Fig. 3. Multi-view point sets (Stanford bunny)

4. LM-ICP method The ICP algorithm is used to estimate the transformation matrix between P and Q. D = {di | i =1, 2, … , Np } is the distance space, and di is calculated by Equation (5). di ( pi , Q) min q j  pi (5) x j X

H = {hi | hi Q, i = 1, 2, … , Np} is the nearest point sets to P, and C1 denotes the search procedure. Given (P, H), the transformation matrix T and distance space D can be computed, and C2, C3 denotes the computation procedure respectively. The coordinates of P are updated until the variation of D is less than threshold. The pseudocode for ICP algorithm is as follows: T initialize P0=P, T0=[1,0,0,0,0,0,0] D0 = C3 (P0, Q)

for k = 0 to N do

Hk = C1 (Pk, Q) Tk = C2 (Pk, Hk ) Pk+1 = Tk·Pk Dk+1 = C3 (Pk+1, Q) if |Dk -Dk+1| is less than τ then print Tk

5.1. Searching overlapping regions 36 points were sampled from point sets A and put into L in the clockwise direction. A continuous 20-point chain was used to search overlapping regions (In fact, 17 points in L locate on the overlapping region). Three uniform sample ratios for point sets B and six curvature thresholds were tested. Fig. 4 (a) shows that the recognition rate of overlapping regions (PˆV/ V) raises up to a peak then drops down with the increase of curvature threshold. Small curvature threshold makes the neighborhood curvature matching degree of {pi| pi L, pi V} become zero, and large curvature threshold increase the chance that {pi| pi L, pi V} gain high neighborhood matching degree, which may influence the search results. Fig. 4 (b) shows that the computing time is multiplied when the sampling rate increases. Based on the testing results, the curvature threshold of 0.2 and sampling ratio of 1/10 were used in the search algorithm, and the results are shown in Fig. 5 (a).

Y.F. Wu et al. / Procedia CIRP 27 (2015) 202 – 206

205

Fig. 5. (a) searching results of overlapping regions; (b) registration results of P and Q; (c) registration results of point sets A and B

6. Conclusions This paper proposed a new method for registration of multi-view point sets with low overlapping ratios. Two kinds of sampling strategies are used to construct plenty of point pairs. Through surface classification, curvature and neighborhood curvature matching of point pairs, the overlapping regions are selected. Then a LevenbergMarquardt ICP (LM-ICP) method is used to search optimal registration parameters with the points in the overlapping regions. Through comparison of different parameters, the search result of overlapping regions is greatly influenced by the curvature threshold and sample ratio. With curvature threshold of 0.2 and sample ratio of 0.1, our method can get good registration accuracy of point sets with low overlapping ratios. An important application in reverse engineering is to measure geometrical deviations between CAD model and the real product. Our method is suitable to register freeform surfaces without any graphic markers. Only accurate registration can guarantee the accurate computation of deviations. Further research will mainly focus on improving computational efficiency with feature recognition, and more complex objects will be studied. Fig. 4. (a) search results under different curvature thresholds and different sample ratios; (b) comparison of computing time

5.2. Registration of point sets The LM-ICP method was used to register two pieces of point sets shown in Fig. 5 (a), and the result is shown in Fig. 5 (b) with the RMS error of 0.266 mm and the computation time of 0.89 s. Fig. 5 (c) shows the registration result of A and B.

Acknowledgements This research has been supported by the Natural Science Foundation of China (Grant No. 51105332, 51275465), the Scientific Research Foundation and the projects of teambuilding for academic backbone of Hangzhou Dianzi University, China (Grant No.2012JXR001). References [1] Werghi N, Fisher R, et al. Object reconstruction by incorporating geometric constraints in reverse engineering. Comput Aided Design 1999;6(31):363-399. [2] Lester H, Arridge SR. A survey of hierarchical non-linear medical image registration. Pattern Recogn 1999;1(32):129-149. [3] Wendt A. A concept for feature based data registration by simultaneous consideration of laser scanner data and photogrammetric images. ISPRS J Photogramm 2007;2(62):122. [4] Varady T, Martin, RR, Cox J. Reverse engineering of geometric models An introduction. Comput Aided Design 1997;4(29):255-268. [5] Gelfand N, Ikemoto L, et al. Geometrically stable sampling for the ICP algorithm. 3DIM 2003; 260-7. [6] Besl PJ, McKay ND. A method for registration of 3-D shapes. IEEE T Pattern Anal 1992;2(14):239-256. [7] Salvi J, Matabosch C, et al. A review of recent range image registration methods with accuracy evaluation. Image Vision Comput 2007;5(25):578-596. [8] Rusinkiewicz S, Levoy M. Efficient Variants of the ICP Algorithm. 3DIM 2001;145-152. [9] Rodrigues M, Fisher R, Liu Y. Special Issue on Registration and Fusion of Range Images. Comput Vis Image Und 2002;1-3(87):1-7. [10] Liu YH. Automatic registration of overlapping 3D point clouds using closest points. Image Vision Comput 2006;7(24):762-781.

206

Y.F. Wu et al. / Procedia CIRP 27 (2015) 202 – 206

[11] Chetverikov D, Svirko D, Stepanov D, et al. The trimmed iterative closest point algorithm. Proceedings - International Conference on Pattern Recognition 2002;3:545-8. [12] Kim C, Son H, Kim C. Fully automated registration of 3D data to a 3D CAD model for project progress monitoring. Automat Constrp 2013;35:587-594. [13] Pauly M, Gross M, Kobbelt LP. 2002. Efficient simplification of pointsampled surfaces. Proceedings of the IEEE Visualization Conference 2002:163-170.

[14] Yang PH, Qian XP. Direct computing of surface curvatures for point-set surfaces. SPBG 2007:29-36. [15] Di Stefano L, Mattoccia S, Tombari F. ZNCC-based template matching using bounded partial correlation. Pattern Recogn Lett 2005;14(26): 21292134. [16] Fitzgibbon A W. Robust registration of 2D and 3D point sets. Image Vision Comput 2003;21:1145-1153.