Technische Universität Wien
Orientation and Processing of Airborne Laser Scanning data
Department of Geodesy and Geoinformation - Research Groups Photogrammetry and Remote Sensing
Performs orientation (co-registration) of multiple point cloud datasets using the Iterative Closest Point algorithm.
This module requires the installation of the Matlab runtime.
The Iterative Closest Point (ICP) algorithm (Chen & Medioni (1991), Besl & McKay (1992)) improves the alignment of two (or more) point clouds by minimizing iteratively the discrepancies within the overlap area of these point clouds. Nowadays the term ICP does not necessarily refer to the algorithm presented in the original publications, but rather to a group of surface matching algorithms which have in common the following aspects:
I: correspondences are established iteratively (i.e. not only once) C: as correspondence the closest point, or more generally, the corresponding point, is used P: correspondences are established on a point basis (instead of using e.g. interpolated surfaces or rasters) The following flow diagram illustrates the functionality of this module. The main iteration loop consists of 7 main steps, which are described below.
In the simplest case, the ICP algorithm is used to improve the orientation of a single loose point cloud with respect to a single fixed point cloud. The main steps of the algorithm are explained on the basis of this case. The initial orientation of the two point clouds is visualized in Figure 2. As can be seen, the point clouds are already roughly aligned, which is a main requirement of the ICP algorithm.
Step 1: Find overlap: The overlap area of the point clouds is determined by using so-called voxel hulls. Voxel hulls are a low resolution representation of the volume occupied by a point cloud. For the computation of the voxel hull the object space is subdivided into a voxel structure. The voxel hull of a point cloud consists of all voxels which contain at least one point of the point cloud. The parameter voxelSize defines the edge length of a single voxel.
Step 2: Selection: A subset of points is selected within the overlap area in one point cloud. For this, points are first selected uniformly in object space, which leads to a homogeneous distribution of the points within the overlap area (uniform sampling). The mean sampling distance in each coordinate direction is defined by the parameter samplingDist. The selected subset of points can be further refined with the strategies normal space sampling or maximum leverage sampling, which are described in section Advanced selection strategies.
Step 3: Matching: Find the closest points (= nearest neighbors) of the selected subset in the other point cloud. This step leads to a set of correspondences which possibly also includes some outliers.
Step 4: Rejection: Rejection of false correspondences (outliers) based on the compatibility of points. Correspondences are rejected if:
Since it is not guaranteed that with these strategies all outliers are rejected, a robust adjustment method is used in step 6 (Minimization) for the detection and deactivation of the remaining ones.
Step 5: Weighting: In this step a weight ( \(0 \le w \le 1\)) is assigned to each correspondence. Weights are based on:
Step 6: Minimization: Estimation of the transformation parameters (for the loose point cloud) by a least squares adjustment. The adjustment minimizes the sum of squared point-to-plane distances. The point-to-plane distance of two corresponding points is defined as the orthogonal distance of one point to the fitted plane of the other point.
Step 7: Transformation: Transformation of the loose point cloud with the estimated parameters. The transformation model can be chosen with the parameter trafoType.
After step 7, a convergence criterion is tested. If it is not met, the process restarts with a new iteration from step 1 (Figure 1). However, if subsets are defined by the parameter subsetRadius, the iteration loop restarts from step 3 and the first two steps are only performed once. The total number of iterations heavily depends on the initial orientation of the points clouds, but typically, between 3 and 10 iterations should be sufficient to reach the global minimum.
Note that opalsICP is not restricted to two point clouds, but can handle an arbitrary number of point clouds simultaneously.
The selection of points (step 2) can heavily affect the final alignment accuracy. Therefore, opalsICP offers two advanced selection strategies, which are especially useful for point clouds where one normal direction is predominating, but the data still includes some valuable features for the alignment. First, points are always selected with the uniform sampling strategy described above, where the mean sampling distance is defined by the parameter samplingDist. To select a subset of this selection, one of the two following strategies can be applied:
A comparison of the selection strategies offered by opalsICP is given in Figure 8. For most ICP variants, this ALS scene is rather difficult because only one feature - the ditch - can constrain the transformation at the finest level. Since normal space sampling and maximum leverage sampling consider the usefulness of points for the alignment process, the number of points can be dramatically reduced, without affecting the solubility of the adjustment. Further information can be found in Glira et al. 2015.
Note that the default parameter values were chosen to be well suited for ALS point clouds. For point clouds at other scales, the Examples section may serve as a reference.
Possible values: zShift .... compute z-shift only shifts .... compute 3D-shifts only rigid ..... compute rigid 3D transformation (3 rotations + 3 shifts = 6 parameters) helmert ... compute 3D Helmert transformation (3 rotations + 3 shifts + 1 scale = 7 parameters) full ...... compute full 3D affine transformation (12 parameters)Allows to specify the transformation model applied to the point clouds within the ICP algorithm.
The data used in the following examples can be found in the $OPALS_ROOT/demo/ directory.
In this example a small strip adjustment (3 strips) is calculated with opalsICP. Besides the mandatory parameters inFile and outDirectory, we also set the samplingDist to 5m and the reduction point (redPoint) since we are directly using LAS data:
Please note, that all values displayed in this visualization are also included in the command line output:
B | S : GLOBALICP > RUNICP > ICP ITERATION 5 OF 5 > RESULTS B | T : ------------------------------------------------------------------------------- B | T : GENERAL STATISTICS: B | T : iteration corresp. std(dp) mean(dp) norm(dx) B | T : 1 2422 0.04827 0.00189 6.42463 B | T : 2 2369 0.03621 0.00070 5.78050 B | T : 3 2363 0.03632 0.00073 0.94349 B | T : 4 2365 0.03605 0.00074 0.29409 B | T : new: 5 2362 0.03611 0.00075 0.11907 B | T : where dp = vector of all (signed) point-to-plane distances B | T : dx = vector of parameter increments B | T : ------------------------------------------------------------------------------- B | T : POINT CLOUD DEPENDANT STATISTICS: B | T : point_cloud corresp. std(dp) mean(dp) file B | T : [1] 1121 0.03689 0.00034 G111 B | T : [2] 2059 0.03478 0.00087 G112 B | T : [3] 1544 0.03726 0.00090 G113 B | T : where dp = vector of all (signed) point-to-plane distances associated to a single point cloud B | T : ------------------------------------------------------------------------------- B | E : GLOBALICP > RUNICP > ICP ITERATION 5 OF 5 > RESULTS
The following files are created in the folder specified by the parameter outDirectory:
icpALS\G111.las icpALS\G112.las icpALS\G113.las
opalsLog.xml): icpALS\ICPLog.txt
icpALS\ICP.mat
The following files are created in the temporary folder (not explicitly specified in the above command):
tempICP\G111.mat tempICP\G112.mat tempICP\G113.matThese point cloud files correspond to the LAS file format and can be opened in MATLAB (as objects of class pointCloud).
Please note that the transformed output files (icpALS\*.las) are exported in the same data format (LAS version, point data record format) as the input LAS files. Please use Module Import to convert the LAS files to OPALS Datamanager format for further use within OPALS.
In this example the orientation of six overlapping close range TLS scans is improved by opalsICP. The initial orientation of the scans is shown on the left-hand side of Figure 10 (coordinate units = mm). Please note: since the default parameter values of this module were chosen to be well suited for ALS point clouds, compared to the previous example, a few more parameters have to be specified:
The new OPALS Datamanager files are placed in the folder icpCloseRange.
The orientation of three TLS scans of a room is improved in this example. As a first step, the data are imported to ODM format:
The following call of opalsICP differs from the previous examples in two points:
icpTLS. Instead, the transformation parameters are written to the file icp_output.xml. This file can afterwards be used to directly transform the input data into any format supported by OPALS. 
Now, we can directly create the transformed LAS files with Module Export using Parameter Mapping to set the transformation parameters trafo:
Besl P., McKay N., 1992. A Method for Registration of 3-D Shapes. In: IEEE Trans. PAMI, Vol. 14, No. 2: 239-256.
Chen Y., Medioni, G., 1991. Object Modeling by Registration of Multiple Range Images. In: Proc. IEEE Conf. on Robotics and Automation, Sacramento, California: 2724-2729.
Glira P., Pfeifer N., Briese C., Ressl C., 2015: A correspondence framework for ALS strip adjustments based on variants of the ICP algorithm.
1.8.17