COMPARISON OF POINT CLOUD REGISTRATION ALGORITHMS FOR BETTER RESULT ASSESSMENT – TOWARDS AN OPEN-SOURCE SOLUTION

: Terrestrial and airborne laser scanning, photogrammetry and more generally 3D recording techniques are used in a wide range of applications. After recording several individual 3D datasets known in local systems, one of the first crucial processing steps is the registration of these data into a common reference frame. To perform such a 3D transformation, commercial and open source software as well as programs from the academic community are available. Due to some lacks in terms of computation transparency and quality assessment in these solutions, it has been decided to develop an open source algorithm which is presented in this paper. It is dedicated to the simultaneous registration of multiple point clouds as well as their georeferencing. The idea is to use this algorithm as a start point for further implementations, involving the possibility of combining 3D data from different sources. Parallel to the presentation of the global registration methodology which has been employed, the aim of this paper is to confront the results achieved this way with the above-mentioned existing solutions. For this purpose, first results obtained with the proposed algorithm to perform the global registration of ten laser scanning point clouds are presented. An analysis of the quality criteria delivered by two selected software used in this study and a reflexion about these criteria is also performed to complete the comparison of the obtained results. The final aim of this paper is to validate the current efficiency of the proposed method through these comparisons.


INTRODUCTION
Terrestrial laser scanning and photogrammetry are two techniques that have proved their efficiency not only for documentation purposes, but also for the monitoring of natural, industrial or man-made sites.Whatever the final use of the produced 3D datasets, the raw products arising from these techniques are most of the time point clouds.While considering raw data that have not been directly georeferenced in a global reference frame, each individual point cloud acquired is expressed in its own local system.To build a consistent global point cloud, all the datasets first need to be grouped in a unique reference frame during the step known as registration.
Due to the large amount of sensors dedicated to the acquisition of 3D data, it becomes common practice in more and more projects to combine datasets from various sensors within a same model.This makes it possible to complete missing parts of a project later on, to work at different scales and levels of detail, or to reach particular areas that may be hardly recorded if using a unique acquisition method.As regards the existing software solutions dedicated to the registration of 3D datasets, two main lacks have been observed: first, even if various data formats may be imported within the program, the heterogeneous properties of datasets acquired using different sensors are not clearly managed.Differences in terms of point density or a priori precision for example are probably not considered during registration process, and this traceability gets lost in the finally delivered model.The second point is related to the way quality criteria are computed and precision is assessed.Even when working with a unique dataset, this may vary from one tool to another, which makes the result comparison between two methods tricky.Moreover, these criteria and their computations are most of the time not clearly explained.This 'black-box' effect can lead to misunderstandings or bad analysis of the result quality, which may be problematic for surveyors or people who have to ensure reliable results.
In order to find a solution to manage potential data heterogeneity, but also to overcome every lack of transparency encountered either during processing or quality assessment, it has been decided to develop an open-source global registration algorithm using the Matlab (MathWorks) language.The final aim of these developments is to simultaneously register multiple point clouds that may come from different sources.A particular attention will also be given to the result assessment at every step of the implemented solution.Since there is more than a unique right way to express quality for registration and georeferencing, the clear definition of useful criteria that makes sense for result interpretation is challenging.The two key elements in this project thus lie in the development of an appropriate method to integrate data heterogeneity, but also in the definition of meaningful criteria for the a posteriori quality assessment.The developed algorithm called MORPHεD (Multisensor glObal Registration based on Points for HEterogeneous 3D Data,pronounced 'morphed' or 'morph3D') is introduced in this paper, together with its extension Geo-MORPHεD for georeferencing.At its current development stage, the approach taking into account data heterogeneity is still under development and will be considered in a future contribution.The aim of this paper is thus to validate the first step of the algorithm, which deals with the global registration of several point clouds acquired with a unique instrument.
After a quick introduction to registration basic principle and a review of related work as well as existing solutions in Section 2, the proposed methodology is developed in Section 3. The problem formulation and the required data inputs are discussed in this section, and the computed residuals and quality criteria are explained.Section 4 is finally dedicated to the analysis of first results achieved for the simultaneous registration of ten laser scanning point clouds.Based on the same input dataset and criteria, a thorough comparison of these results with the ones obtained using two existing software solutions is finally carried out for the validation of the current version of the algorithm.

General formulation of the 3D transformation
Registration consists in a three-dimensional similarity transformation (also referred to as 3D conformal coordinate transformation) which is applied to the 3D coordinates of one or multiple point clouds to transform them into a reference coordinate system.A 6 or 7 parameter similarity transformation model is used, depending on the presence or not of a scale factor which is the same in all directions.This scale factor may be useful while working with non-scaled image-based point clouds for example.The remaining six parameters are three rotation angles and a 3-components translation vector.
Considering a point i observed from two systems A and B, the general formulation for such a transformation is: (1) where are 3D coordinates in A or B is a 3 × 3 orthogonal rotation matrix is a translation vector.

Related work
Registration is one of the first steps involved in point cloud processing and 3D model generation.As explained by (Bae and Lichti, 2008) or (Rabbani et al., 2007), this stage is crucial for the quality of the final products since registration errors propagate and accumulate with the individual point uncertainties.Registration errors also accumulate differently depending on the registration strategy applied.According to the scene configuration as well as to the digitization project, the registration may involve only two point clouds (pairwise registration) or multiple point clouds (global registration), that have been collected using one unique tool or several different devices.The number of point clouds composing the project and their source differences increase the complexity and the duration of the computation.
When the original reference frames are located far from each other, the registration is performed in two steps known as coarse and fine alignment of the point clouds.Whereas pre-aligned point clouds can be finely registered using the widespread ICP algorithm firstly proposed by (Besl and McKay, 1992) or one of its multiple variants as exposed in (Rusinkiewicz and Levoy, 2001), the first step consisting in coarse alignment of the point clouds present various solutions.If the most implicit methods are based on the use of artificial targets evenly distributed in the recorded scene as described in (Franaszek et al., 2009), the literature presents a large amount of other methods that tend to improve the level of automation or the pre-processing time.The use of natural features visible in the scene such as planar features for (Theiler and Schindler, 2012) or cylinders for (Hullo et al., 2012) has been considered, as well as the use of 2D information coming from range or intensity images in the work of (Weinmann et al., 2011).If the registration of point clouds using images is based on the detection of 2D keypoints in these images, more recent researches focus on the use of 3D keypoints directly detected in the 3D datasets, as proposed by (Theiler et al., 2014).At this stage, the proposed method is based on the use of artificial or natural common points as described in Section 3.2.If promising results are achieved, then the introduction of a 3D keypoint detection algorithm will be considered in the future, in order to increase the level of automation.

Software and tools
Mainly based on the principles mentioned before, software solutions are available to deal with the registration of laser scanning point clouds.When developed by the device manufacturers, such as FARO Scene or Trimble RealWorks among other companies, these software enable the automatic registration of the acquired datasets and provide visually good results.Unfortunately their algorithms and especially the way the accuracy is assessed are not clearly outlined in these commercial solutions.Moreover, they are initially not designed to automatically allow the combination with exterior datasets, which can be point clouds captured using image-based techniques or low-cost sensors.To consider the registration of point clouds acquired using different techniques, open-source solutions such as CloudCompare (Girardeau-Montaut, 2018) and also software from the academic domain such as 3DVEM (Fabado et al., 2013) can be found.In the latter case, even if the registration methodology implemented is not clearly detailed, the quality criteria applied are easier to interpret by the user than those provided by commercial solutions.
To confront the results delivered by the algorithm proposed in this paper with reference results, two software solutions have been selected to reprocess the same test dataset presented in Section 4.1.Among the commercial solutions, FARO Scene has been chosen since the considered dataset has been acquired using a FARO terrestrial laser scanner.As regards tools taken from the academic domain, 3DVEM has been preferred to CloudCompare which is limited to pairwise registration, because the aim of this paper is to assess a global registration methodology.The two selected tools are quickly described in the next subsections.

FARO
Scene is largely used to handle the 3D point clouds produced with the terrestrial laser scanners of the FARO Focus series, from data import to export through the registration and georeferencing of point clouds.Several methods are offered for the registration, based either on natural or artificial targets placed by the user in the environment or based on top views into the 3D data, which is especially suited for indoor environments.A cloud to cloud refinement method is also available.The detection and matching of artificial targets that are spheres or flat targets is automatic, unless the user wants to manually add some specific points or assign the correspondences.The version of FARO Scene used in this study is the currently most recent, namely Scene 7.1.1.81.

3DVEM -
Register GEO is developed by the research group GIFLE from the University of Valencia (Spain), and it is also dedicated to the registration and the georeferencing of 3D data.The registration is based on common points whose coordinates are given by the users, with or without pre-allocated identifiers for correspondences since the matching can be made automatically.Several point clouds can be handled in a same computation, thus the functionalities are close to those offered by the proposed algorithm.A report containing a summary about the performed computations, the final transformation parameters and the computed residuals can finally be generated.The way the residuals and quality criteria are handled in this report is discussed in Section 4.3 and they will be confronted to the method adopted in our approach.Note that the version used to perform the tests in this paper is 3DVEM 1.0.2.8.

GLOBAL REGISTRATION APPROACH
To perform a global registration of several point clouds in a unique computation, the basic 3D transformation principle described in Section 2.1 has to be extended.As a matter of fact, it is not sufficient to consider all the connections between adjacent point clouds in the dataset.A reference frame for the transformation also needs to be defined.The method evaluated in this paper is described in this section.

Basic principle and formulation
While considering a global registration approach, all the relationships between the point clouds being registered should be considered and expressed in a mathematical form.The sets of parameters that will be simultaneously estimated during adjustment computation are the transformation parameters that link each moving point cloud to the reference system.For a dataset composed of n point clouds, n-1 point clouds are moving point clouds and one cloud is taken as reference.Thus there will be as many individual sets of six or seven parameters (depending on the presence of a scale factor) to determine as there are moving clouds.Even if the unknown parameters to determine connect moving point clouds to the reference cloud, the links between adjacent moving data have to be taken into account in order to achieve a global approach.For this reason two kinds of relationships will be expressed using two different sets of equations, depending on the considered tie point.These relationships are illustrated in Figure 1 and are detailed in the following paragraphs.This idea is inspired by the aerotriangulation principle for independent models by simultaneous transformations, which is a compromise between the sequential construction of strip models and simultaneous bundle adjustment.This principle taken from the analytical photogrammetry domain is thoroughly described in (Wolf and Dewitt, 2010), and the equations given hereafter have been adapted from this book to fit the current problematic.

Figure 1. Direct and indirect connections between moving point clouds and the reference system
The first set of equations is written for a common point i which belongs to the reference point cloud and to a moving cloud A, as depicted in Figure 1.It has the same form as the general Equation (1) given in Section 2.1, using adapted subscripts: (2) In this equation the terms are as defined for Equation (1), except that the coordinates of the point i known in the system of A (Xi A ) are transformed in the reference system (Xi Ref ).Thus the transformation parameters contained in the rotation matrix and translation vector RA Ref and tA Ref are related to the reference system.To complete this first set of equations, a tie point j that links two adjacent point clouds A and B which are not the reference cloud is considered, as shown in Figure 1.To enforce that such a tie point must coincide between adjacent clouds, equations of the following form are written: (3)

Data inputs
Since registration is based on common points, files containing point coordinates have to be associated to each imported point cloud.If the correspondences between targets in each scan are known, pre-allocated point identifiers can be given by the user.However the automatic matching of common points is much comfortable.It has been developed through an iterated process based on the RANSAC paradigm initially introduced by (Fischler and Bolles, 1981).In this case, files containing only three-dimensional coordinates are necessary.This automatic approach based on an adapted version of RANSAC principle has the advantage of being robust to outliers, thus tie points poorly chosen or artificial targets displaced during the scanning process will be excluded from the correspondences.

Correspondences, reference system and initial transformations
To define the reference point cloud, the correspondences that have been automatically defined or pre-allocated using identifiers are considered.A direct relationship is established between two point clouds if they have at least three tie points in common.Nevertheless all the correspondences are listed for each point cloud even if they are indirect, which means based on less than three common points.The point cloud with the highest number of direct connections is chosen as the reference point cloud which will define the reference system for the registration.If several point clouds meet this condition, then the number of overall correspondences is observed to select as reference the point cloud with the most correspondences.If there are still several candidates after this second condition, the point cloud located in a central position regarding the list of files is arbitrarily defined as reference data, assuming that the data have been acquired successively in the spatial distribution.
Based on these direct and indirect relationships, a graph is built to represent the links between all stations, as shown in Figure 2. Nevertheless it may occur that some point clouds are not directly connected to the reference, and thus no initial transformation to the reference can be computed for these clouds.In this case a path that links such a cloud to the reference node is established using the generated graph, in order to compute an initial transformation based on a composition of several transformations.For direct connections in the graph, approximate values of the transformation parameters are computed regarding the solution proposed by (Dewitt, 1996).

Adjustment
To estimate in a unique computation all the individual sets of transformation parameters between moving point clouds and the reference system, the Generalized Least Squares (GLS) method (also called General Least Squares in (Ghilani, 2010)) is applied to perform the adjustment.The functional model which is composed of the two kinds of equations developed in Section 3.1 is written in a general form: (4) This model is composed by a set of n 0 equations expressing the relationships between the u 0 unknowns x (transformation parameters) and the n observations l (tie point coordinates).
Observation equations that link adjacent point clouds apart from the reference using Equation (3) already follow this general form.However the first set of equations following the notation in Equation (2) needs to be slightly modified to be used in the functional model, so that: (5) Since Equations (3) as well as (5) applied in the functional model are nonlinear, they have to be linearized using Taylor's series approximations and therefore initial approximations x0 of the unknown parameters are required.According to (Ghilani, 2010) but using a slightly different notation, the linearized observation equations for the GLS approach can be written in matrix form: (6) where dx denotes the corrections on the unknown parameters, v is the residuals vector related to observations, and K is called misclosure vector.In Equation ( 6), two Jacobian matrices (sometimes called design matrices) are involved: Jacobian matrix A contains the partial derivatives of the equations with respect to the unknown parameters, and B contains the partial derivatives with respect to the observations.While still employing matrix algebra, Equation ( 6) has the solution: ( 7) with ( 8) The square matrix M in Equation ( 8) has dimensions (n 0 × n 0 ) and is called equivalent weight matrix in (Ghilani, 2010).The equivalent weight matrix involves the (n × n) weight matrix W defined by the user, which contains on its diagonal the individual weights (w 1 ,…,w n ) associated to each observation.
This last aspect may justify the use of GLS instead of parametric approach to solve the adjustment in the proposed context.In the parametric approach of least squares, the weight matrix is a square matrix which dimensions correspond to the number of equations.However in the developed method, the number of equations n 0 is not equal to the number of observations n, since a single equation involves two sets of coordinate observations (one in each related point cloud).Moreover in Equations ( 3) and ( 5) applied in the functional model, observations and unknown parameters are hard to dissociate because of the multiplicative factors.Applying GLS method, though, a specific weight can be assigned to each of the n observations through a diagonal (n × n) weight matrix.This will be interesting in future work for weighting the observations depending on their source.

Potential indirect georeferencing
Once all point clouds have been registered into a same reference system, a second transformation of these grouped data may be necessary depending on the project.Georeferencing consists in the transformation of a dataset known in a local system into a global system, such as a national framework for example.One speaks about indirect georeferencing when the point clouds have first been registered, before the registered dataset is georeferenced.Indirect georeferencing can be carried out using the developed algorithm, if necessary.This step requires to have determined beforehand the three-dimensional coordinates of representative points in the global system, thanks to tacheometric observations for example.These points may be artificial targets if some have been placed in the scene, or natural points that are selected due to their representativeness.Points selected to perform the georeferencing must also be known in the local reference system anyway.Most of the time when point-based registration is first performed, these points are chosen among the tie points used for the registration.In this case, mean coordinates of the registered tie points observed from several stations are determined.
Direct georeferencing, which consists in the direct transformation of individual point clouds into the final global system, is not foreseen for the moment, but this option could be easily added to the algorithm.

Result assessment
The aim of developing an internal algorithm for registration is not only to be able to complete it with further developments, but also to bring more transparency into the result assessment.This goes from the way the residuals are computed, to a reflexion about the choice of meaningful quality criteria.Since registration and georeferencing does not have the same purpose, residuals are computed differently.

Registration assessment:
To this end, residuals are computed based on tie points that have been transformed into the reference system.Since all the tie points used to perform the global registration are not necessarily defined in the reference point cloud, it is not always possible to have reference values for all common points.As a matter of fact, most probable values of these tie point coordinates are defined by computing the mean values of all the transformed coordinates in the reference system, for each common point.These mean values include the coordinates of tie points being present in the reference cloud, even if these coordinates are untransformed.
For each point cloud, the computed residuals are 3D deviations between the transformed coordinates of common points associated to this cloud, and the mean values computed before.
Even if the reference point cloud is static in its reference system, residuals are computed for this reference data since the raw coordinates of its tie points differ from the mean coordinates.For each registered point cloud, this makes it possible to determine a standard deviation based on the residuals of its tie points.These statistical values give an idea of the internal precision of the process and are good indicators to detect potential problems with the tie points.Indeed, a high value of this standard deviation indicates potential high residuals for the tie points of the considered point cloud.Thus tie points with high residuals that should probably be excluded from the global registration can be easily detected.

3.6.2
Georeferencing assessment is slightly different since external and independent reference coordinates are known and used to determine the transformation.For each reference point, a 3D deviation is computed between the reference value and the mean value of georeferenced coordinates, taken from transformed stations where the considered point appears.A Root Mean Square Error (RMSE) is then calculated based on these residuals.Note that it would be more correct to compute the deviations between reference values and the local coordinates used to achieve the georeferencing.However, since it is not known if FARO Scene and 3DVEM make use of mean local coordinates of tie points to process their georeferencing, the first approach using the mean georeferenced coordinates to compute the residuals is applied in this paper.
If all the reference points have been used during the georeferencing process, they can no more be considered as really independent to the computation.In this case, a RMS error can still be calculated but it will rather define the internal precision reached by the computation.To be able to give an absolute precision (or accuracy), residuals have to be computed based on reference points that have not been included in the georeferencing process.These can thus be used as external checkpoints.Both approaches are possible within the proposed algorithm, depending on the availability or not of an external reference dataset.Regardless of the approach, in the case of indirect georeferencing, the residual values are computed for the reference points.Thus reference points presenting high deviations can be excluded from the georeferencing.In contrast, announcing precision values for each point cloud makes more sense when a direct georeferencing is performed.

EXPERIMENTAL RESULTS AND COMPARISONS
First results obtained using the methodology described in this paper are discussed in this section, after a quick presentation of the test dataset.

Dataset and results
The dataset used in this context has been acquired in the interior part of the Saint-Pierre-le-Jeune Catholic Church in Strasbourg (France).Ten laser scanning stations have been set up using a FARO Focus X330 terrestrial laser scanner.Spheres have been evenly placed in the digitized environment for the automatic registration of the data using FARO Scene software.Among the 12 spheres that have been placed in the recorded environment, the centers of 11 spheres have been accurately measured with a total station for the later georeferencing of the project.The distribution of scanning stations and of the 12 spheres is depicted in Figure 3. Since the implemented algorithm as well as 3DVEM software require common point coordinates as input data, spheres have been automatically detected using FARO Scene and their local coordinates have been listed before performing the registration.These coordinates have been reused as input data with MORPHεD and 3DVEM, thus the three methods can be compared based on the same inputs.As regards the results delivered by MORPHεD, Figure 4a shows a top view of the ten unorganized point clouds before the registration, whereas Figure 4b gives an overview of the registered data.The considered dataset mainly covers the transept of the church between the nave and the chancel, even if the scanned points cover almost the whole building due to the scanning range of the device, as it can be seen in Figure 5.
Based on a simple visual inspection, it appears in Figure 4b that all the point clouds have been correctly placed in the reference system.This qualitative analysis will be completed with a quantitative analysis in the next section.

Methods of comparison
The quantitative analysis of the proposed algorithm is based on a comparison of the obtained results with results achieved using FARO Scene and 3DVEM software.The comparison can be undertaken at different levels: for example, the finally adjusted transformation parameters as well as the available quality criteria can be compared to this end.
The risk while comparing the adjusted transformation parameters lies in the way the rotation matrix is expressed, because this may differ between the three approaches.This would result in different rotation angles that are not directly comparable.Moreover, a difference of some hundredth of radians for the computed rotation angles does not give a clear idea of the difference in terms of precision between the different registration solutions.For this reason, it has been decided to deal with the residuals computed based on tie points (or reference points) to provide a clear comparison of the precision reached after registration (or georeferencing) of the point clouds.That being said, the values to confront need to be computed the same way in order be comparable, which is visibly not the case for Scene and 3DVEM.The first problem to face is thus the understanding of the quality criteria delivered by the software, as explained in the next section.

Residuals and precision computation
The computation of residuals is performed very differently depending on the software used.Even the terminology applied varies from one software to another, and it is not always very explicit which values have been actually computed.A detailed analysis has been carried out for the two software solutions considered in the comparison, and the following are assumptions that have been verified through empiric observations after reprocessing of the values for the test dataset.

Residuals in FARO Scene:
In the current version of Scene used during this study, the calculated residuals are called target tensions and represent distance errors expressed in millimeters.These values are 3D deviations computed pairwise (between two point clouds) for each target.For example if a target is observed in four different scans, all in all six deviations are computed because six is the binomial coefficient to form pairs from a set of four points.Considering all possible combinations in a project, this list of residuals may be long.Besides, this approach can be questionable since these residuals relate to pairwise deviations, but not global deviations.To get an idea of the precision for each registered point cloud, a mean target distance error (in mm) is computed as being the mean value of all previously defined deviations where the considered point cloud is involved.After georeferencing is processed, residuals are computed between paired point clouds as described before using the transformed coordinates of the targets used.3D deviations between the reference coordinates and the transformed coordinates in the point clouds where these references appear are also added to the list of residuals.Thus absolute residuals with respect to references as well as relative values between coordinates from two georeferenced point clouds are grouped into a unique list.Finally, based on this list of 3D deviations a mean value is derived for each point cloud as after the registration, which is in our view not the most representative solution to investigate the result of indirect georeferencing.

Residuals in 3DVEM:
As regards the 3DVEM software, considering the report which is generated after registration, it seems that the coordinates of tie points available in the reference point cloud are taken as references to compute the residuals.As a matter of facts, residuals and 3D deviations related to the reference cloud are zero.For the remaining tie points that do not appear in the reference point cloud, coordinates are taken as references from one of the moving point clouds where they appear in a way is unfortunately not described.This assumption is made based on the fact that residuals are zero for each of the remaining tie points in one or several moving clouds.Thus the coordinates chosen as reference to compute residuals are not selected uniformly.Using these residuals, 3D deviations are computed for each tie point of each scan, but no unique value is synthetized per scan.While considering georeferencing, residuals are computed for each reference point used for the global transformation.These deviations are the differences between the reference coordinates, and one set of georeferenced coordinates for each point arbitrarily chosen among the list of georeferenced coordinates from all the point clouds.Once again, this approach is questionable and not clearly explained.

A thorough comparison
of the three methods only makes sense if comparable values are confronted.To this end, all the residuals and precision values have been reprocessed the same way for FARO Scene and 3DVEM, using the approach implemented in MORPHεD and detailed in Subsection 3.6.To reprocess the residuals, only the transformed coordinates which are available after registration as well as georeferencing in both tools are useful.The following analyses take only into account these reprocessed residuals and precisions.

Comparisons and assessment of the method
4.4.1 Registration: Assessing the registration is possible thanks to the residuals that have been computed for each tie point of each station.These residuals are not exactly the same for the three registration tools, varying from some millimeters but with the same order of magnitude.This is due to the fact that different transformation parameters and thus different transformed coordinates of tie points have been estimated in each case.delivered by the three methods.Having a look at this graph, one can assume that the developed computation methodology is closer to the one applied in 3DVEM than in Scene because of very close standard deviations.However the processing of further datasets is required to confirm such an assumption.Indeed, the small variations are quite insignificant and may be due to rounding errors during the computation process.That being said, the same trend is observed for each point cloud with the three methods.Based on this result, the registration offered by the proposed algorithm can be validated.This analysis strengthens the visual rendering presented in Figure 4b.

Georeferencing accuracy:
To finally assess the global accuracy reached for this dataset after georeferencing using Geo-MORPHεD, external checkpoints unused during the transformation estimation but whose coordinates are known in the global system are necessary.To this end, georeferencing has been reprocessed for the considered dataset with the proposed algorithm only.Among the 11 targets that have been measured using a total station, six targets have been used to achieve this second georeferencing.The five remaining spheres can thus be defined as checkpoints, since they are not involved in the computation.3D deviations with respect to the reference values are then calculated for these checkpoints, and a RMS error of 2 mm can be derived from these residuals listed in Table 8.Table 8.Georeferencing residuals of five checkpoints The reached georeferencing accuracy of 2 mm is a satisfactory result regarding the a priori precision of the laser scanner used and the overall dimensions of the building.This is another indicator of the efficiency of the algorithm presented in this paper.Note that the identifiers of the spheres and stations are visible in Figure 3. Nevertheless, no detailed analysis of the individual results is made here, since the real accuracy assessment of the test dataset is not the topic of this paper.
To go further in the analyses and complete the quality criteria offered in the presented algorithm, one could imagine determining and plotting error ellipsoids for the stations or for the tie points (registration), but also for reference points (georeferencing).This would help for the visual inspection of reliable and non-reliable points, and will be considered in future work.

CONCLUSION AND FUTURE WORK
In this paper, an open source algorithm dedicated to the global registration of multiple point clouds has been presented.The mathematical model used to that end has been inspired by the aerotriangulation principle in photogrammetry domain.The global registration of the 3D data is based on the use of common points such as targets or manually selected representative points, which may be automatically and robustly matched if the initial correspondences are unknown.In this paper, the algorithm has been tested for the global registration of ten laser scanning point clouds representing the interior part of a church.
To validate the results obtained using MORPHεD introduced in this paper, the same dataset composed of raw point clouds and their associated common point files has been processed with one commercial and one academic software.To compare the results delivered by the three solutions, a thorough analysis of the quality criteria and a reflexion about the meaning of these values have been carried out.To achieve a rigorous comparison, the residuals and statistic criteria computed based on the common points have been reprocessed in the same way.The high similarity between the resulting values enabled to conclude about the efficiency of the proposed algorithm, which builds a strong basis for future developments.Besides, the advantage of using MORPHεD is the transparency of the quality assessment.
As regards the supported input 3D data, for the moment the algorithm is dedicated to the registration of point clouds from different data formats.A first way to improve the methodology will be the solution proposed to handle the heterogeneity of datasets acquired with various sensors.Even if this paper focused on the evaluation of the raw method, the idea of weighting the datasets depending on their source is under study.Besides, another point to be considered in future developments is the integration of other 3D data than point clouds, like for example surveying measurements from a total station or from GNSS.This would be a further step in the combination of heterogeneous data.

Figure 2 .
Figure 2. Graph presenting the connections between stations based on at least three tie pointsexample of a case study with 10 point clouds

Figure 3 .
Figure 3. Overview of the distribution of the scanning stations and of the registration targets (spheres)

Figure 4 .
Figure 4. Top view of the ten unorganized point clouds before registration (a) and after registration (b) using MORPHεD

Figure 6 .
Figure 6.Comparison of standard deviations (mm) computed for each point cloud after registration For georeferencing stage, only an internal precision can be assessed if all the external reference points have been used during the global transformation estimation.3D deviations are computed for each reference point while considering an indirect georeferencing, as described before.These values are shown on the graph in Figure7, based on the results delivered by each of the three compared solutions.As for the registration step before, the same trend is observed regarding each compared tool.Moreover, based on these residuals a same RMS error of 2.3 mm is calculated for each of the three tools.The results provided by the algorithm extension Geo-MORPHεD for georeferencing are thus validated.

Figure 7 .
Figure 7.Comparison of residuals (mm) computed for each reference point used during georeferencing