Speeding up coarse point cloud registration by threshold-independent baysac match selection

: This paper presents an algorithm for the automatic registration of terrestrial point clouds by match selection using an efficiently conditional sampling method -- threshold-independent BaySAC (BAYes SAmpling Consensus) and employs the error metric of average point-to-surface residual to reduce the random measurement error and then approach the real registration error. BaySAC and other basic sampling algorithms usually need to artificially determine a threshold by which inlier points are identified, which leads to a threshold-dependent verification process. Therefore, we applied the LMedS method to construct the cost function that is used to determine the optimum model to reduce the influence of human factors and improve the robustness of the model estimate. Point-to-point and point-to-surface error metrics are most commonly used. However, point-to-point error in general consists of at least two components, random measurement error and systematic error as a result of a remaining error in the found rigid body transformation. Thus we employ the measure of the average point-to-surface residual to evaluate the registration accuracy. The proposed approaches, together with a traditional RANSAC approach, are tested on four data sets acquired by three different scanners in terms of their computational efficiency and quality of the final registration. The registration results show the st.dev of the average point-to-surface residuals is reduced from 1.4 cm (plain RANSAC) to 0.5 cm (threshold-independent BaySAC). The results also show that, compared to the performance of RANSAC, our BaySAC strategies lead to less iterations and cheaper computational cost when the hypothesis set is contaminated with more outliers.


INTRODUCTION
Terrestrial laser scanning (TLS) is effectively used in urban mapping for applications like as-built documentation, change detection, 3D reconstruction of architectural details and building facades. However, one of the biggest problems encountered while processing the scans is the so-called registration step in which rigid transformation parameters (RTPs) are determined in order to align scans obtained from different scan positions into one common coordinate system. RTPs are computed from minimizing distances between common points, obtained by correspondence tracking, in scenes sampled from different viewpoints. In general two registration goals are distinguished, coarse and fine registration. A coarse registration of sufficient quality is typically needed as input for a successful fine registration (Vosselman & Maas, 2010). This paper is focused on coarse registration.
Coarse registration should be able to compute an initial estimate of the transformation between two point clouds that are quite far apart. A well-known approach is to compute motion invariant features (or descriptors) for both sets and consecutively match these features. The feature extraction and matching step is generally approached in three different ways: Manual pointbased matching, target-based matching and feature-based matching. However, existing fully automated approaches may only work in certain contexts.
Besides differing on the type of descriptors used, coarse registration algorithms also differ on the matching methods used to establish correspondences between scans (e.g. ordinary Least Squares adjustment, Genetic Algorithms (Brunnström & Stoddart, 1996), Principal Component Analysis ''PCA'' (Chung & Lee, 1998) or RANdom Sample And Consensus ''RANSAC'' (Fischler & Bolles, 1981)). Kang et al. (2009) proposed an iterative algorithm, using Lease Squares adjustment, for optimal transformation parameters estimation. Akca (2010) achieved the parameter estimation using a Generalized Gauss-Markov model. Grant et al. (2012) proposed a point-to-plane fine registration approach that also uses a General Least Squares adjustment model. In many cases RANSAC is used to automatically remove false matches. However, the original RANSAC algorithm assumes constant prior probabilities for all data points and chooses hypothesis sets at random, which leads to many iterations and high computational costs as many hypothesis sets will be contaminated with outliers.
Therefore, the scope of this paper is focused on improving the efficiency of coarse matching by improving on the RANSAC The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLI-B5, 2016 XXIII ISPRS Congress, 12-19 July 2016, Prague, Czech Republic paradigm, The SIFT-based reflectance image matching method (Kang et al., 2009) is employed to determine initial 3D correspondences. For robust estimation of transformation parameters from 3D correspondences contaminated by outliers, we implement a conditional sampling methodoptimized BaySAC (Kang et al., 2014) to always select the minimum number of data required with the highest inlier probabilities as hypothesis set. Moreover, LMedS method is used to construct the cost function that is used to determine the optimum model. We evaluate the results of the optimized BaySAC strategy compared to a traditional RANSAC approach both by considering the required computational efforts and by comparing the quality of the resulting registration. For the latter, we compare the RMSEs of the distances between matching points in the final inlier set and, in addition, evaluate the Euclidean distances between the transformed points in the analyzed scan to the surface fitted to the fixed scan.

DETERMINATION OF THE PRIOR INLIER PROBABILITIES OF DATA POINTS
Determining the proper prior inlier probability is important for the BaySAC algorithm. We implemented the statistical testing process proposed in (Kang et al., 2014), which is generic and can be applied to any BaySAC problem. The process is implemented using a histogram that illustrates the distribution of the discrete hypothesis model parameter sets that are computed during the different iterations and the degree of convergence of each candidate parameter set that describes how the other sets converge to it. The degree of convergence of a bin in the histogram is calculated as the number of parameter sets in that bin divided by the total number of parameter sets.
Our presented strategy is based on a histogram to dynamically evaluate the convergence of the hypothesis parameters sets during the hypothesis testing process ( Fig.1 shows a histogram to evaluate the convergence of the rigid transformation parameters. Different convergent cluster of the hypothesis sets will be presented in the histogram. We select the oldest parameter set as the reference point for each convergent cluster of parameter sets. The more the hypothesis sets converge to a cluster, the more possible is it that the reference parameters set of the cluster is correct. Therefore, we present the convergence degree of a cluster to evaluate the correct possibility. The convergence degree is a percentage which describes the number of sets converge to a parameters sets cluster. It is calculated by dividing the number of the hypothesis sets in the cluster by the total number of hypothesis sets.
In this paper, the rigid transformation model in point cloud registrationis used as the hypothesis models.
The hypothesis testing process starts with a RANSAC strategy which chooses initial data sets at random. After that, during each iteration, the distribution of the hypothesis parameters sets is updated by adding the parameters set corresponding to the newly evaluated hypothesis sets. If a newly added set doesn't converge to any existing cluster, it will be regarded as a new cluster in the distribution. Else, the number of solutions in an existing cluster is increase by one and the convergence degree of the different clusters is updated accordingly.
When the degree of convergence of a cluster in the distribution of parameter solutions reaches a predefined threshold, the first hypothesis set in that cluster is used to determine the prior inlier probabilities of the data points according Eq. (1): where P i denotes the prior probability of point i, D i is the distance between point i and the fitted primitive, and m represents the predefined threshold for outlier identification, which is set as five times the point precision.

PROBABILITY UPDATING
After determining the prior inlier probabilities of each corresponding pair, Eq.
(2) (Kang et al., 2014) is employed to update the inlier probabilities during consecutive iterations: where I is the set of all inliers, H t is the hypothesis set of n data points used in iteration t of the hypothesis testing process, P t-1 (i  I) and P t (i  I) denote the inlier probabilities for data point i during iterations t-1 and t, respectively, k is the number of points consistent with the model during a test, and D is the total number of data points.

CONSTRUCTION OF THE COST FUNCTION USING LMEDS
We applied the LMedS method (Rousseeuw & Leroy, 1987) to construct the cost function that is used for determining the optimum model to reduce the influence of human factors and improve the robustness of the model estimate.
We first randomly select a subset of corresponding points from the samples to calculate the model parameters and then calculate the deviation (i.e., distances between corresponding points) of all the other sample points with respect to the model. We use the equation as follows to calculate the deviation in the where i D is the distance from point i to the model, and N is the number of data points used to test the model parameter.
We used Med deviation to replace the number of inlier points complying with the hypothesis model as the criteria to determine the optimum model. After the iterative process of the hypothesis test is finished, we select the model with the smallest Med deviation as the optimum model.

RESULT EVALUATION
Accuracy is of great importance for point cloud registration. Therefore, various error metrics have been defined to measure the registration accuracy of point clouds (Simon, 1996;Maas, 2000;Rusinkiewicz & Levoy, 2001). These include the change in rotation angle and/or translation vector, the distances between corresponding points and the distances between points and their corresponding surfaces. Among them, point-to-point and point-to-surface error metrics are most commonly used, both of which are utilized in this paper. As the registration process is implemented using 3D corresponding points, we naturally use the st.dev of the distances between corresponding points in the first scan and their transformed matches from the second scan, to evaluate the quality of the registration results. However, the error evaluated as a point to point difference between two scans after registration in general consists of at least two components, a random measurement error component and a systematic error component as a result of a remaining error in the found rigid body transformation. Still, in most publications the quality of a registration is expressed in terms of standard deviations, which assumes an underlying normal distribution. Thus we employ the measure of the average pointto-surface residual, introduced by (Lindenbergh & Pfeifer, 2005) to evaluate the registration accuracy.

Distances between corresponding points
The first method to evaluate the quality of the resulting registration is by considering the distances at corresponding points as follows. Applying the methods on the initial set of corresponding pairs results in general different sets of correct correspondences. As indicated before, the correct correspondences are used to estimate a final rigid body transformation TR using least squares. The resulting transformation TR is applied to each inlying corresponding point in the second scan. Finally the st.dev of the distances between corresponding points in the first scan and their transformed matches from the second scan is determined. Advantage of this approach is that it is expected that matches are found throughout the overlap between the two scans.

Average point-to-surface residual
The error evaluated as a point to point difference between two scans after registration normally comprises at least two components, i.e. random measurement error and registration error. As errors in the computed transformation parameters lead to non-random differences between the fixed and transformed scans, we consider the registration error as systematic.
Random errors are also known as compensating errors, since they tend to partially cancel themselves in a series of measurements. In this paper, we implement the U test (Lehman & Romano, 2005) to estimate the normality of the distribution of the residuals. The formula to compute U value in terms of the skewness or kurtosis of the distribution is as follow.
If the random measurement error is reduced, the remaining difference is expected to approach the real registration error. First, the surface representations, e.g. plane and curve, are obtained from the first scan using the RANSAC algorithm. We then compute the distances between the transformed points (from the second scan) to the surface of the fixed reference, the first scan. For the purpose of reducing the effect of random measurement noise, the point-to-surface residuals are averaged to a regular grid spanned over the considered surface segment. The average point-to-surface residual per cell is expected to reduce the random error so that this residual approaches the real registration error.

EXPERIMENT AND ANALYSIS
The performances of the threshold-independent BaySAC algorithm was evaluated on real datasets in terms of both registration accuracy and computational efficiency.

Description of the test data
Four datasets are considered, acquired respectively by the 3D SwissRanger camera, by the LMS-Z620 and VZ-400 laser scanners from RIEGL, as well as by the FARO LS 880 laser scanner (Fig.2) Two robust strategies, i.e. RANSAC, BaySAC-LMedS, were respectively utilized to prune false correspondences.

Registration accuracy
As presented in Section 5, the registration accuracies of RANSAC and BaySAC-LMedS on four datasets were evaluated using two measures, the average distance between inlier correspondences after registration and the average point-tosurface residual. Fig.3 shows that a certain number of correspondences (e.g. 5 of Dataset 1) were identified for each dataset as check points. Tab.2 lists the evaluation of registration accuracies in terms of the average distance between the selected correspondences after registration. BaySAC-LMedS achieved overall higher accuracy (the average RMSs 8.1 mm) compared with the performance of the plain RANSAC (the average RMSs 22.9 mm).

Average transformed-point-to-analyzed-surface residual per cell
For the purpose of registration quality evaluation we selected ten planar segments from Datasets 1, 2, 3 and 4 and nine curved-surface segments from Dataset 2. These segments are visualized in Fig.5 To evaluate the registration error, we computed the Euclidean distances between the transformed points in the analyzed scan and the surface fitted to the fixed scan. Non-zero distances are mainly caused by three factors: measurement errors, surface fitting errors and registration errosr. The first factor is believed to be random, therefore the influence of the two other factors are analyzed in the following sections.    The average point-to-surface residuals are given in Tab 4. All signed average point-to-surface residuals are below 1 mm, which proves that the signed point-to-surface residuals can partially cancel themselves after averaging, compared with the RMS of the point-to-surface residuals.

2) Points to surface registration residuals
To evaluate the registration error, we computed the Euclidean distances between the transformed points in the analyzed scan to the surface fitted to the fixed scan. Fig.6 shows the histograms of these distances. Compared to Fig. 4, in Fig. 6 obvious shifts from the origin are visible, which we interpret as being caused by registration errors because the surface fitting error is proven to be random.
Afterwards, we superimposed a regular grid with an edge length of 10cm on the selected surface segment and calculated the average point-to-surface difference for each cell. Fig.7 ( Segment D1) shows that both of the magnitudes and signs of the average point-to-surface differences per cell are non-random and are expected to approach the real registration error. The spatial distribution of the residuals also appears systematic (Fig.7), for instance most larger differences occur in the bottom left part of Fig.7(a). Fig.7(b) shows a distribution of differences with smaller absolute values, compared to the distributions shown in Fig.7(a).
The average differences between the transformed points from the registered scan and the reference surface are shown in Tab

Computational efficiency
Hypothesis set evaluation is an iterative process, so the calculational efficiency of the proposed strategies was evaluated in terms of the number of iterations and computation time. As described in Section 2, the BaySAC-LMedS strategy consists of two parts, i.e. the random part (BaySAC-LMedS -Random) and the BaySAC part. The random part consists of the iterations implemented with a RANSAC strategy, during which the inlier probabilities of the different correspondences are determined as well.
As shown in Tab.6, the ratio between the number of inliers and the number of possible correspondences in Dataset 1 reaches 81%. The number of iterations for RANSAC varies from 6 to   18. In attempts No. 1,2,7,8,as presented in Section 2, the number of iterations for RANSAC is so small that the BaySAC-LMedS process failed in the determination of initial inlier probabilities of data points before the RANSAC process ended, therefore the BaySAC-LMedS degraded to a RANSAC process. The total number of iterations for BaySAC-LMedS varies from 2 to 11. The results show that, although the numbers of consumed iterations of BaySAC-LMedS is smaller than those of RANSAC as a whole, the differences are not so obvious due to the high percentage of inliers.
The percentage of inliers in Dataset 2 (e.g. point cloud pair I-II) reduces to 38% (Tab.6). The number of iterations for RANSAC varies between 283 and 437. The random part of BaySAC-LMedS takes between 4 and 26 iterations. As presented in Section 2, after the determination of the initial inlier probability for each correspondence, the strategy of hypothesis set evaluation changes from RANSAC to BaySAC. Therefore, the total number of iterations, compared to those for a plain RANSAC strategy, sharply reduces to the range from 14 to 76. The results show that, although the number of possible correspondences is only 61 (Tab.6), the decrease of percentage of inliers leads to large differences between the performance of RANSAC and BaySAC-LMedS. When the percentage of inliers for Dataset 3 and 4 reduces further to 23% and 25%, the computational efficiency of RANSAC becomes much lower compared to BaySAC-LMedS. Tab. 6 compares the different statistics for the evaluated strategies. The average consumed time of the three strategies is listed in Tab. 7. As shown in Tab. 7, the random parts of the time consumed by BaySAC-LMedS was to determine the prior inlier probabilities of data points, which are the core of BaySAC-CONV.
After the determination of the prior inlier probabilities, the use of BaySAC remarkably reduced the computational time needed to find a good model compared with the time consumed by plain RANSAC. For instance, in Tab.7 60 milliseconds were spent by BaySAC-LMedS to find a good transformation from Dataset 3, while the time costed by plain RANSAC is 982 milliseconds.

CONCLUSIONS
In this paper, we have proposed conditional sampling algorithm for point cloud registration using the threshold-independent BaySAC framework. The measure of the average point-tosurface residual was employed to evaluate the registration accuracy for the purpose of reducing the random measurement error and thus approaching the real registration error.
The proposed algorithms were implemented and evaluated on four TLS data sets in terms of their computational efficiency and the quality of the registration results. The results show that, compared to the performance of the original RANSAC framework, the computational efficiency of BaySAC-LMedS scales with the percentage of inliers. The smaller this percentage, the higher the computation efficiency gain of BaySAC-LMedS on RANSAC. When the percentage of inliers becomes very large (e.g. 81%), the BaySAC-LMedS strategy may fail and thus degrades to a RANSAC process. The registration results also indicate that, among the two different considered strategies, BaySAC-LMedS achieves the highest registration quality on all experimental datasets.
The proposed algorithm aims at robustly estimating a single transformation between two scans. Future work will consider how BaySAC-LMedS can be applied for the co-registration of multiple point clouds in the simultaneous estimation of multiple transformations.