# Point Cloud Registration Methods for Long-Span Bridge Spatial Deformation Monitoring Using Terrestrial Laser Scanning

Authors: Hongwei Zhang; Yanjie Zhu; Wen Xiong; C. S. Cai
Journal: Structural Control and Health Monitoring
 (2023)

---

## Abstract

In recent years, efforts have been devoted to utilizing terrestrial laser scanning for bridge spatial performance inspection, but they are still restricted to small or medium-span bridges, like some historical heritages. Due to the large-scale dimensional features of long-span bridges, applications of 3D point cloud techniques still remain challenging, such as the extra-long scan range and extreme-small incidence angle when scanning a bridge with a span over 1000 m. Moreover, rare attempts can be found for the performance evaluation of point cloud registration methods for long-span bridges as well, which is a critical basis for further spatial deformation recognition on the point cloud data. Hence, in this study, a cross-evaluation of three iterative closest point (ICP) registration methods is conducted for long-span suspension bridges, namely, traditional ICP, kd-tree-based ICP, and feature point-based ICP algorithms. We conducted field laser scanning on the Ma’anshan Yangtze Bridge, an 1880 m long suspension bridge located in China. The results show that the feature point-based ICP algorithm outperforms the other two in terms of convergence rate and execution time for a single iteration due to the smaller number of registration points compared to the other two algorithms. Moreover, it also gives more precise values in terms of bridge tower spatial deformation identification. Meanwhile, due to efficient data organization, the kd-tree-based ICP algorithm takes less time for a single iteration than the traditional ICP algorithm. Finally, two suggestions for algorithm improvement in terms of efficiency optimization and accuracy improvement of long-span bridge deformation analysis are proposed based on the assessment results.


---

## Body

## 1. Introduction

Bridge deformation is one of the most crucial indexes that can reflect its service condition; hence, research efforts have been devoted to identifying the spatial deformation of a bridge for further performance assessment. In terms of spatial deformation, targets of interest mainly include the vertical deflection and torsion of the main girder or the arch rib, the vertical and horizontal deviation of the pier or bridge tower, etc. To identify such deformations, strain-based solutions are the most popular and traditional ones. By investigating the relationship between deformation and strain, the former can be easily deduced from the latter. Strain gauges are mostly utilized because of their high accuracy and low cost [1, 2]. However, several disadvantages of strain gauges cannot be ignored, such as the complex installation process, short service life, and poor stability and durability [3].With the rapid development of surveying and mapping technology, optical instruments such as levels and total stations have been widely used for bridge deformation monitoring [4]. The popularization of these instruments realizes the noncontact remote measurement of bridge deformation. Such noncontact remote measuring technology was popularized for single-point testing with high precision but was limited by its poor continuity, which therefore makes it inapplicable for long-span bridge monitoring. Besides, the global positioning system (GPS) is another remote deflection measuring technique for bridge structures [5]. In contrast to the optical instruments, GPS can provide dynamic measurements in an automatic way, which greatly facilitates the practical applications. However, it has a more distinct measurement error than the aforementioned two methods; thus, GPS is more suitable for the deformation measurement of large-scale structures. In addition, the above noncontact methods have the same weakness that the finite single-point coordinates cannot well reflect the deformation of the whole bridge.To overcome the shortcomings of traditional noncontact inspection, terrestrial laser scanning (TLS), a rising 3D imaging technology, has emerged in recent years [6, 7]. The laser scanner emits and receives high-density laser beams at the target and calculates the angle and distance of each beam based on the principle of laser ranging to obtain the 3D coordinates of the target point. The object within the scanning range is established as a 3D coordinate model, that is, a 3D point cloud. Therefore, the laser scanning methods can not only overcome the challenges of traditional contact measuring solutions but also possess the merits of the most existing noncontact detection approaches, such as long-range measurement capabilities and a high degree of automation. At the same time, it can also model the massive 3D point coordinates of the measured target surface rather than the form of the little single points [8]. Due to the convenience of the data acquisition and the robustness of the point cloud data interpretation, this technique has been widely employed in landslide monitoring [6, 9], as-built modeling [10], surface defect detection [11–13], construction morphological identification [14, 15], medium and small span bridge health monitoring [16, 17], and bridge component recognition [18]. However, rare practices can be found for long-span bridges due to their large-scale dimensions, which results in the challenge of an extra-long scan range and an extreme-small incidence angle for laser scanners. For instance, according to our practices for using TLS on several long-span bridges, a bridge with a span over 1000 m requires a laser scanner with a scanning range of at least 500 m to cover the whole bridge girder by setting two stations on both sides. However, in practice, due to the space limitations underneath the bridge deck and the straight-line shape of the girder, the laser incidence angle is normally less than 10°, resulting in the low reflectivity of lasers and the reduction of scanning range accordingly. Thus, the point cloud conditions that are missing are inevitable, which will seriously affect the precision and robustness of further deflection identification. To overcome this challenge, substantial efforts should be given to point cloud interpretation methods; however, rare attempts can be found for long-span bridges.In terms of 3D point cloud interpretation methods, point cloud registration is the most fundamental and critical step to align point sets properly that were collected in a different coordinate system. Approaches to point cloud alignment have been available in a variety of fields [19], such as 3D model reconstruction [20], solar energy surveys [21], and cultural heritage management [22]. Again, practice on long-span bridges is still insufficient. It maximizes the overlap of several point sets, named source or target point clouds, by a certain standard. The establishment of correspondence and maximization of the coincidence between the source point cloud and the target point cloud are the key tasks of point cloud registration. Currently, point cloud registration can be divided into two steps: coarse registration and fine registration. The objective of coarse registration is to reduce the excessive translation and rotation difference between source point clouds and target point clouds. For instance, the handcrafted feature extraction-based method calculates the transformation matrix between point clouds through the geometric features extracted, such as intersection lines, geometric primitives, etc. [23–25]. The four-point congruent set (4pc)-based approach [26, 27] utilizes multiple groups of four coplanar points to compute the spatial transformation between point clouds, which improves the robustness of the algorithm. Based on the result of coarse registration, fine registration refines the approximate direction and translates between point clouds [19]. Instead of feature extraction, the fine registration directly obtains the correspondence from the 3D coordinate information of each point in the point cloud and finally determines the transformation of the point clouds. The normal distribution transform algorithm and the iterative closest point (ICP) algorithm are probably the two most commonly applied fine registration algorithms currently [9, 19]. By representing the point cloud as a set of Gaussian distributions with different probability density functions, the normal distributions transform algorithm achieves alignment based on the similarity of the probability density functions of the reference and transformed point clouds [28, 29]. This algorithm runs speedy but does not have sufficient accuracy, leading some scholars to apply it even in coarse registration [19]. Due to its conceptual simplicity, the ICP algorithm is the most commonly used approach for fine registration. This algorithm establishes the correspondence between point clouds based on the principle of Euclidean distance between the nearest points and then gradually reduces the sum of squares of Euclidean distances between the corresponding points through the transformation matrix at each iteration until the result converges [30].The ICP algorithm was first proposed by Besl and McKay [30] and Chen and Medioni [31]. They established point-to-point and point-to-face correspondence between point clouds, respectively, and finally both obtained transformation matrices that satisfied the error requirements. The iterative process of determining the corresponding point set and computing the best rigid body transformation makes the ICP algorithm good at aligning point clouds. However, due to its registration process, the ICP algorithm is sufficiently accurate but has a large convergence cost in most cases, and it also has strict initial attitude requirements. Subsequent research has optimized and improved the algorithm in various aspects. Invariant features were combined into the distance calculation of point pairs, and it made the method converge to the global optimum more easily than traditional methods [32]. The branch-and-bound scheme was integrated into ICP, and thus the restriction on the initial position was lifted [33]. Other studies updated the algorithm in terms of efficiency improvement, noise filtering, and error optimization [34, 35], which dwindled the operation time and iteration numbers.The virtues of ICP algorithms have been noticed by construction and bridge engineers in recent years. Cabaleiro et al. achieved the registration of the old and new bridge point clouds with the ICP algorithm [36]. Points exceeding the threshold value determined by the distance between the nearest points were considered to be deformation points in the bridge. This method can obtain the overall deformation of the bridge, unlike traditional monitoring, which can only acquire the change value of a limited number of single points. Vezocnik et al. [14] proposed an ICP-based registration, which utilized all points or points that have been randomly sampled point clouds to calculate the geometry of components for engineering quality control. The overall rotation and translation of several single-pier highway bridges were identified by the ICP registration algorithm [37], and the compression and torsion states of the main girders were obtained. On this basis, Kalasapudi et al. [38] proposed an improved ICP algorithm that used certain rules to eliminate the points in the bridge point cloud that undergo position changes, and obtained a more reasonable transformation matrix with deformation recognition results. The transformation relations were obtained by these ICP algorithms and their variants by registering old and new point sets and finally acquiring the deformation of the structures [9]. However, while the above studies mainly focus on point cloud registration and deformation analysis of bridges with a small or medium span, the effectiveness of ICP registration algorithms for long-span bridge point clouds still lacks in-field verification. In fact, due to the large dimensional features of the long-span bridges, their deformation patterns are more complex than those of small and medium-span bridges, which makes the accuracy of the registration algorithm more rigorous. Meanwhile, high precision inevitably results in high computational costs. Thus, it is a challenge and worth investigating an ICP registration algorithm to achieve a balance between registration accuracy and computation efficiency for long-span bridges.To identify the performance of various ICP-related algorithms used for long-span bridges, we conducted the evaluation of the traditional ICP algorithm, the kd-tree-based ICP algorithm, and the feature point-based ICP algorithm for bridge deformation analysis in this study. The kd-tree is ak-dimensional binary index tree developed from the basic tree structure, where each node represents a point in the k-dimensional space [39]. The idea behind the feature point-based approach is that existing registration algorithms cannot distinguish between rigid body motion and the deformation of structural elements [40]. The kd-tree-based ICP algorithm and the feature point-based ICP algorithm are respectively optimizations of the traditional ICP algorithms for point cloud data organization and registered point sets. Then, the Ma’anshan Yangtze River Bridge, a long-span cable-supported bridge, is used for the evaluation of the three algorithms. The assessment of the algorithms involves the convergence rate, execution time for a single iteration, root mean square (RMS) error between point clouds, transformation matrix, and precision of measured bridge deformation.The rest of this paper is organized as follows: Section2 introduces the theoretical background of three ICP-related algorithms. The long-span bridge profiles utilized for method assessment are presented in Section 3. Then, Section 4 illustrates the preprocessing before registration and the registration results of three ICP algorithms. The evaluation of the bridge deformation analysis results is then discussed in Section 5. This section also gives the optimal registration algorithm among the three algorithms and future improvements to the ICP algorithm for bridge tower deformation analysis followed by a conclusion in Section 6.


## 2. Theoretical Background of Registration Algorithms

As aforementioned, point cloud registration consists of the coarse registration and the fine registration. According to previous research, algorithms for coarse registration roughly align the two point clouds and have already provided good input for fine registration [38]. Hence, this study will focus on the algorithms for fine registration, covering the traditional ICP algorithm and two improved ICP algorithms, i.e., the kd-tree-based ICP algorithm and the feature point-based ICP algorithm.


### 2.1. Traditional ICP Algorithm

The ICP algorithm is a precise and efficient iterative algorithm for the registration of 3D shapes, including free-form curves and surfaces [30], which has been widely applied for 3D point cloud fine registration [19]. Given a sufficient set of initial rotations and translations, the mean square distance metric can be globally minimized over all six degrees of freedom of the point cloud to achieve geometric model alignment under unknown correspondence. The proposed method is applicable to the registration of 3D shapes, including point sets, implicit curves, and parametric surfaces and obtains good results [30].The Euclidean distanced (p1, p2) between two points p1=x1,y1,z1 and p2=x2,y2,z2 is denoted as follows: (1)dp1,p2=x2−x12+y2−y12+z2−z12.It is easy to understand that the minimum distance from a pointp to a 3D entity E can be represented by the minimum value of the distance from point p to a point pE in the 3D entity E. The minimum distance between p to E is written as follows:(2)dminp,E=mindp,pE.Based on this principle, this algorithm finds the closest points to the source point cloud in the target point cloud and forms a one-to-one correspondence between two point sets. For instance,pA in the source point cloud and pB in the target point cloud are a pair of corresponding points, satisfying the condition given by(3)dpB,pA=mindpB,A.The transformation matrix is a parameter that characterizes the differences in translation and rotation between the source point cloud and the target point cloud. In other words, in the ideal state where the two point clouds are morphologically identical, the translation and rotation discrepancies between the target point cloud and the source point cloud can be eliminated after applying the spatial transformation to the target point cloud by using the transformation matrix. The unit quaternion method is employed to calculate the spatial transformation matrix due to its excellent performance in the calculation of 3D spatial transformations. The unit quaternion is represented as a four-dimensional vector:qR=q0q1q2q3 for q0≥0 and q02+q12+q22+q32=1. Due to the arithmetic rule of the unit quaternion, any point pi in 3D space can be represented by the quaternion method, as given in the following equation where q0=0,q1=pi,x,q2=pi,y,q3=pi,z:(4)pi=0+pi,xi+pi,yi+pi,zi.According to the properties of quaternions, the rotation matrix can be given by(5)R=q02+q12−q22−q322q1q2−q0q32q1q2−q0q32q1q2−q0q3q02+q22−q12−q322q2q3−q0q12q1q3−q0q22q2q3+q0q1q02+q32−q12−q32.At the same time, three parameters are added to represent the translation vectorT as(6)T=q4q5q6.Thus, the overall transformation matrix is integrated as(7)TM=R,Tt.Let the points in the source point cloudA = {ai}and the points in the target point cloud B = {bi}already correspond with the same index. In practice, due to the possible difference in the spatial morphology of the two point clouds collected at two various times, there is still a deviation between A and B after applying TM to B. This deviation can be denoted as the error function in the following equation where Na is the number of points of A. The algorithm keeps the least squares iterating the process of computing and applying the transformation matrix until error function f is less than a given threshold that determines the accuracy of registration.(8)f=1Na∑i=1Naai−Rbi−T2.The least-squares iteration is also utilized to minimize the error functionf in the following equation [31]. The distance between the point and the surface is used as the error measure, which leads to the error function f being transformed into the form of(9)f=∑i=1Nads2TMbi,S,where ds is the distance from a point to a plane and S is the normal of the point cloud A at ai.The above two algorithms with highly similar frameworks are called traditional ICP algorithms.The flow of traditional ICP algorithm is given in Figure1.Figure 1
Flowchart of traditional ICP algorithm.


### 2.2. Kd-Tree-Based ICP Algorithm

The kd-tree-based ICP algorithm improves the efficiency of registration from the perspective of optimizing matching point pairs by exploiting the ability of the kd-tree to organize large amounts of data [39]. Specifically, the nearest point matching in the traditional ICP algorithm was replaced by using the kd-tree method, which can reduce the time cost of the algorithm.The kd-tree method, where each node can represent a point in thek-dimensional space, is a binary index tree-based data classification approach [41]. The discriminator of the kd-tree in each layer is utilized for a branch decision. The kd-tree takes the form of an empty tree or a binary tree with the following properties:(1)
Its left and right subtrees are also kd-tree(2)
If its left subtree is not empty, the value of thedth dimension of all nodes in the left subtree is less than the value of the dth dimension of its parent node, where d is the value of the discriminator in the parent node(3)
If its right subtree is not empty, the value of thedth dimension of all nodes in the right subtree is greater than or equal to the value of the dth dimension of its parent node, where d is the value of the discriminator in the parent nodeThe kd-tree not only inherits the advantage that the average lookup length of a binary tree is only1+4logn where n is the number of points [41] but also is more suitable for indexing three-dimensional spatial points or for organizing and managing data with multiple attributes due to its discriminator.It is necessary to consider the information of the neighboring points of each data point to simplify the processing of the point clouds due to the clutter of massive point cloud data [41]. The steps of the kd-tree-based point cloud algorithm that can execute the range search of the point cloud are listed as follows:(1)
DefineN (i) as the point cloud neighborhood of the initialized node, r as the search range of the point cloud neighborhood, and s (i, r) as the search function.(2)
The algorithm determines whether pointi is a leaf node of some search subtree, and if so, returns the node. If the lookup range includes the neighborhood range of the search, then all leaf nodes of that subtree are returned.(3)
IfL (left (i)) where left (i) is the left subtree of i and r have duplicate parts, then execute the search function s (left (i), r). If L (right (i)) where right (i) is the left subtree of i and r have duplicate parts, then execute the search function s (right (i), r) and finally return to the point cloud neighborhood.The traditional ICP registration algorithm has to calculate a large number of iterations of point data, causing low efficiency of the algorithm and point matching bias, eventually affecting the reliability of the registration results.Using kd-tree instead of the nearest point finding part of the traditional ICP algorithm can reduce the algorithm’s time cost. Then a more reliable correspondence between the source point cloud and the target point cloud is established, followed by the calculation of the transformation matrix and achievement of the point clouds registration. The flow of the kd-tree-based ICP algorithm is illustrated in Figure2:Figure 2
Flowchart of the kd-tree-based ICP algorithm.


### 2.3. Feature Point-Based ICP Algorithm

Considering that the previous ICP algorithm used all points for registration, the feature point-based ICP algorithm considered that there were interference points in the point cloud that needed to be removed before registration [38]. It is worth mentioning that these interference points are not noise points that were removed in the earlier point cloud preprocessing, but rather points that have changed position in the point cloud. If these interference points are involved in the process of calculating the transformation matrix, the registration results obtained may have unpredictable errors. Therefore, determining and excluding these interference points before registration by developing some rules is the kernel of this approach. The rules for identifying interference points in the feature point-based algorithm and the method of registering two point clouds using feature points will be introduced next.Principal component analysis (PCA) is a statistical method that converts a set of potentially correlated variables into a set of linearly uncorrelated variables through an orthogonal transformation [42]. The transformed set of variables is called the principal component. The 2D “major plane” that is near most of the points in the point cloud can be obtained from PCA. Specifically, PCA can be used as a plane fitting method to find a plane close to the most of points in the 3D point cloud and ignore a small portion of “anomalous” points far from the majority of the point cloud [43].After obtaining the major planes Ps and Pr in the source point cloud and the target point cloud, the interference points are determined by several key parameters in Table1.Table 1
Major parameters for interference points determination.
Distance between Ps and PrD (Ps, Pr)Distance from the point in the source point cloud to PsDistance from the point in the source point cloud to PsD (a (i), ps)D (b (i), ps)Change in the distance between corresponding pointsD (i) = |D (a (i), ps)−D (b (i), ps)|D (i) > D (ps, pr)D (ps, pr) <= D (i)Interference pointsFeature pointsAs presented in Table1, it is determined that the point is an interference point that needs to be excluded when the difference in the distance from the corresponding point pair to the major plane Ps is greater than the distance between the two major planes. The calculation of the distance from each point to the plane is based on the fitted major plane of the source point cloud to facilitate the calculation of the corresponding point pair changes D (i).The adoption of the feature points filtered by the above PCA-based rules as the set of alignment points is an enhancement of this algorithm. Due to the rigid body transformation principle of the traditional ICP algorithm, the point cloud rigid body consisting of these feature points allows the algorithm to calculate a transformation matrix that may be closer to the theoretical rigid body transformation matrix. The registration is completed according to the remaining steps of the ICP algorithm after the feature points are obtained, that is, after the spatial transformation matrix that makes the error function less than the threshold is acquired iteratively. The flow of this algorithm is shown in Figure3.Figure 3
Flowchart of feature point-based ICP algorithm.


## 2.1. Traditional ICP Algorithm

The ICP algorithm is a precise and efficient iterative algorithm for the registration of 3D shapes, including free-form curves and surfaces [30], which has been widely applied for 3D point cloud fine registration [19]. Given a sufficient set of initial rotations and translations, the mean square distance metric can be globally minimized over all six degrees of freedom of the point cloud to achieve geometric model alignment under unknown correspondence. The proposed method is applicable to the registration of 3D shapes, including point sets, implicit curves, and parametric surfaces and obtains good results [30].The Euclidean distanced (p1, p2) between two points p1=x1,y1,z1 and p2=x2,y2,z2 is denoted as follows: (1)dp1,p2=x2−x12+y2−y12+z2−z12.It is easy to understand that the minimum distance from a pointp to a 3D entity E can be represented by the minimum value of the distance from point p to a point pE in the 3D entity E. The minimum distance between p to E is written as follows:(2)dminp,E=mindp,pE.Based on this principle, this algorithm finds the closest points to the source point cloud in the target point cloud and forms a one-to-one correspondence between two point sets. For instance,pA in the source point cloud and pB in the target point cloud are a pair of corresponding points, satisfying the condition given by(3)dpB,pA=mindpB,A.The transformation matrix is a parameter that characterizes the differences in translation and rotation between the source point cloud and the target point cloud. In other words, in the ideal state where the two point clouds are morphologically identical, the translation and rotation discrepancies between the target point cloud and the source point cloud can be eliminated after applying the spatial transformation to the target point cloud by using the transformation matrix. The unit quaternion method is employed to calculate the spatial transformation matrix due to its excellent performance in the calculation of 3D spatial transformations. The unit quaternion is represented as a four-dimensional vector:qR=q0q1q2q3 for q0≥0 and q02+q12+q22+q32=1. Due to the arithmetic rule of the unit quaternion, any point pi in 3D space can be represented by the quaternion method, as given in the following equation where q0=0,q1=pi,x,q2=pi,y,q3=pi,z:(4)pi=0+pi,xi+pi,yi+pi,zi.According to the properties of quaternions, the rotation matrix can be given by(5)R=q02+q12−q22−q322q1q2−q0q32q1q2−q0q32q1q2−q0q3q02+q22−q12−q322q2q3−q0q12q1q3−q0q22q2q3+q0q1q02+q32−q12−q32.At the same time, three parameters are added to represent the translation vectorT as(6)T=q4q5q6.Thus, the overall transformation matrix is integrated as(7)TM=R,Tt.Let the points in the source point cloudA = {ai}and the points in the target point cloud B = {bi}already correspond with the same index. In practice, due to the possible difference in the spatial morphology of the two point clouds collected at two various times, there is still a deviation between A and B after applying TM to B. This deviation can be denoted as the error function in the following equation where Na is the number of points of A. The algorithm keeps the least squares iterating the process of computing and applying the transformation matrix until error function f is less than a given threshold that determines the accuracy of registration.(8)f=1Na∑i=1Naai−Rbi−T2.The least-squares iteration is also utilized to minimize the error functionf in the following equation [31]. The distance between the point and the surface is used as the error measure, which leads to the error function f being transformed into the form of(9)f=∑i=1Nads2TMbi,S,where ds is the distance from a point to a plane and S is the normal of the point cloud A at ai.The above two algorithms with highly similar frameworks are called traditional ICP algorithms.The flow of traditional ICP algorithm is given in Figure1.Figure 1
Flowchart of traditional ICP algorithm.


## 2.2. Kd-Tree-Based ICP Algorithm

The kd-tree-based ICP algorithm improves the efficiency of registration from the perspective of optimizing matching point pairs by exploiting the ability of the kd-tree to organize large amounts of data [39]. Specifically, the nearest point matching in the traditional ICP algorithm was replaced by using the kd-tree method, which can reduce the time cost of the algorithm.The kd-tree method, where each node can represent a point in thek-dimensional space, is a binary index tree-based data classification approach [41]. The discriminator of the kd-tree in each layer is utilized for a branch decision. The kd-tree takes the form of an empty tree or a binary tree with the following properties:(1)
Its left and right subtrees are also kd-tree(2)
If its left subtree is not empty, the value of thedth dimension of all nodes in the left subtree is less than the value of the dth dimension of its parent node, where d is the value of the discriminator in the parent node(3)
If its right subtree is not empty, the value of thedth dimension of all nodes in the right subtree is greater than or equal to the value of the dth dimension of its parent node, where d is the value of the discriminator in the parent nodeThe kd-tree not only inherits the advantage that the average lookup length of a binary tree is only1+4logn where n is the number of points [41] but also is more suitable for indexing three-dimensional spatial points or for organizing and managing data with multiple attributes due to its discriminator.It is necessary to consider the information of the neighboring points of each data point to simplify the processing of the point clouds due to the clutter of massive point cloud data [41]. The steps of the kd-tree-based point cloud algorithm that can execute the range search of the point cloud are listed as follows:(1)
DefineN (i) as the point cloud neighborhood of the initialized node, r as the search range of the point cloud neighborhood, and s (i, r) as the search function.(2)
The algorithm determines whether pointi is a leaf node of some search subtree, and if so, returns the node. If the lookup range includes the neighborhood range of the search, then all leaf nodes of that subtree are returned.(3)
IfL (left (i)) where left (i) is the left subtree of i and r have duplicate parts, then execute the search function s (left (i), r). If L (right (i)) where right (i) is the left subtree of i and r have duplicate parts, then execute the search function s (right (i), r) and finally return to the point cloud neighborhood.The traditional ICP registration algorithm has to calculate a large number of iterations of point data, causing low efficiency of the algorithm and point matching bias, eventually affecting the reliability of the registration results.Using kd-tree instead of the nearest point finding part of the traditional ICP algorithm can reduce the algorithm’s time cost. Then a more reliable correspondence between the source point cloud and the target point cloud is established, followed by the calculation of the transformation matrix and achievement of the point clouds registration. The flow of the kd-tree-based ICP algorithm is illustrated in Figure2:Figure 2
Flowchart of the kd-tree-based ICP algorithm.


## 2.3. Feature Point-Based ICP Algorithm

Considering that the previous ICP algorithm used all points for registration, the feature point-based ICP algorithm considered that there were interference points in the point cloud that needed to be removed before registration [38]. It is worth mentioning that these interference points are not noise points that were removed in the earlier point cloud preprocessing, but rather points that have changed position in the point cloud. If these interference points are involved in the process of calculating the transformation matrix, the registration results obtained may have unpredictable errors. Therefore, determining and excluding these interference points before registration by developing some rules is the kernel of this approach. The rules for identifying interference points in the feature point-based algorithm and the method of registering two point clouds using feature points will be introduced next.Principal component analysis (PCA) is a statistical method that converts a set of potentially correlated variables into a set of linearly uncorrelated variables through an orthogonal transformation [42]. The transformed set of variables is called the principal component. The 2D “major plane” that is near most of the points in the point cloud can be obtained from PCA. Specifically, PCA can be used as a plane fitting method to find a plane close to the most of points in the 3D point cloud and ignore a small portion of “anomalous” points far from the majority of the point cloud [43].After obtaining the major planes Ps and Pr in the source point cloud and the target point cloud, the interference points are determined by several key parameters in Table1.Table 1
Major parameters for interference points determination.
Distance between Ps and PrD (Ps, Pr)Distance from the point in the source point cloud to PsDistance from the point in the source point cloud to PsD (a (i), ps)D (b (i), ps)Change in the distance between corresponding pointsD (i) = |D (a (i), ps)−D (b (i), ps)|D (i) > D (ps, pr)D (ps, pr) <= D (i)Interference pointsFeature pointsAs presented in Table1, it is determined that the point is an interference point that needs to be excluded when the difference in the distance from the corresponding point pair to the major plane Ps is greater than the distance between the two major planes. The calculation of the distance from each point to the plane is based on the fitted major plane of the source point cloud to facilitate the calculation of the corresponding point pair changes D (i).The adoption of the feature points filtered by the above PCA-based rules as the set of alignment points is an enhancement of this algorithm. Due to the rigid body transformation principle of the traditional ICP algorithm, the point cloud rigid body consisting of these feature points allows the algorithm to calculate a transformation matrix that may be closer to the theoretical rigid body transformation matrix. The registration is completed according to the remaining steps of the ICP algorithm after the feature points are obtained, that is, after the spatial transformation matrix that makes the error function less than the threshold is acquired iteratively. The flow of this algorithm is shown in Figure3.Figure 3
Flowchart of feature point-based ICP algorithm.


## 3. Ma’anshan Yangtze River Bridge

### 3.1. Profile of the Bridge

In this section, the parameters and details of the bridge being modeled in 3D are provided. The Ma’anshan Yangtze River Bridge is a super-kilometer cable-supported bridge built in 2011 located in Anhui Province, China. The main bridge of the Ma’anshan Yangtze River Bridge is a three-tower, two-span suspension bridge with a symmetrical arrangement, and the span of the main bridge girder is 2 × 1080 m. The realistic view of the Ma’anshan Yangtze River Bridge is presented in Figure4.Figure 4
Ma’anshan Yangtze river bridge.In this study, the south bridge tower point cloud is utilized for the evaluation of the algorithm. The south bridge tower is a concrete portal frame structure consisting of four parts: the tower column, the lower crossbeam, the upper crossbeam, and the saddle hood at the top of the tower column. The tower column is a reinforced concrete structure with a height of 165.3 m and its widths in the vertical and horizontal directions are 6.0 m and 10.0 m, respectively.


### 3.2. 3D Point Cloud of the Bridge Tower

The towers of suspension bridges bear the important task of transferring the loads from the main cables to the foundation, which results in the spatial deformation of the towers being a significant indicator in structural health monitoring. Specifically, the longitudinal and lateral offset of the bridge towers are the key metrics for measurement. The total station-based inspection method derives the offset values by identifying the difference between the measurement points at the top and bottom of the bridge tower. The measurement points at the bottom are easy to determine due to their proximity and stability. However, the measurement points at the top of the bridge towers are not only far away but also have a tendency to shake due to the vibrations of the towers, which causes difficulties and errors in the measurements.Based on measurement experience and the working principle of the 3D laser scanner, four stations were placed on the ground in the vicinity of the bridge tower to acquire its 3D point cloud data. The Lecia P50 3D laser scanner was employed for point cloud data acquisition in this study. This device is one of the most accurate 3D scanning devices available, with a single point accuracy of 1.2 mm + 10 ppm.The establishment of the final integral bridge tower point cloud requires a process of stitching together in scanner-adapted postprocessing software, cyclone, which can directly edit the data from the scanner’s original format. Figure5(a) shows the point cloud of the bridge tower and the approximate location of the measuring stations, while Figure 5(b) presents the 3D laser scanner at work. It is noted that according to the scanner accuracy parameters and the real distance from the instrument to the top of the tower (up to approximately 200 m), the point cloud error at the top of the tower is approximated to be 1.4 mm at the maximum.Figure 5
Construction of the bridge tower point cloud: (a) bridge tower point cloud and measuring stations. (b) 3D laser scanner at work.
(a)(b)


## 3.1. Profile of the Bridge

In this section, the parameters and details of the bridge being modeled in 3D are provided. The Ma’anshan Yangtze River Bridge is a super-kilometer cable-supported bridge built in 2011 located in Anhui Province, China. The main bridge of the Ma’anshan Yangtze River Bridge is a three-tower, two-span suspension bridge with a symmetrical arrangement, and the span of the main bridge girder is 2 × 1080 m. The realistic view of the Ma’anshan Yangtze River Bridge is presented in Figure4.Figure 4
Ma’anshan Yangtze river bridge.In this study, the south bridge tower point cloud is utilized for the evaluation of the algorithm. The south bridge tower is a concrete portal frame structure consisting of four parts: the tower column, the lower crossbeam, the upper crossbeam, and the saddle hood at the top of the tower column. The tower column is a reinforced concrete structure with a height of 165.3 m and its widths in the vertical and horizontal directions are 6.0 m and 10.0 m, respectively.


## 3.2. 3D Point Cloud of the Bridge Tower

The towers of suspension bridges bear the important task of transferring the loads from the main cables to the foundation, which results in the spatial deformation of the towers being a significant indicator in structural health monitoring. Specifically, the longitudinal and lateral offset of the bridge towers are the key metrics for measurement. The total station-based inspection method derives the offset values by identifying the difference between the measurement points at the top and bottom of the bridge tower. The measurement points at the bottom are easy to determine due to their proximity and stability. However, the measurement points at the top of the bridge towers are not only far away but also have a tendency to shake due to the vibrations of the towers, which causes difficulties and errors in the measurements.Based on measurement experience and the working principle of the 3D laser scanner, four stations were placed on the ground in the vicinity of the bridge tower to acquire its 3D point cloud data. The Lecia P50 3D laser scanner was employed for point cloud data acquisition in this study. This device is one of the most accurate 3D scanning devices available, with a single point accuracy of 1.2 mm + 10 ppm.The establishment of the final integral bridge tower point cloud requires a process of stitching together in scanner-adapted postprocessing software, cyclone, which can directly edit the data from the scanner’s original format. Figure5(a) shows the point cloud of the bridge tower and the approximate location of the measuring stations, while Figure 5(b) presents the 3D laser scanner at work. It is noted that according to the scanner accuracy parameters and the real distance from the instrument to the top of the tower (up to approximately 200 m), the point cloud error at the top of the tower is approximated to be 1.4 mm at the maximum.Figure 5
Construction of the bridge tower point cloud: (a) bridge tower point cloud and measuring stations. (b) 3D laser scanner at work.
(a)(b)


## 4. Point Cloud Registration Results

In this section, three ICP algorithms, i,e., the traditional ICP algorithm, the kd-tree-based ICP algorithm, and the feature point-based ICP algorithm are evaluated by using the point cloud data of the Ma’anshan Yangtze River Bridge, which was collected in 2020 and 2022 as previously described in Section3.Table2 gives the parameters of the three algorithms in the iteration process. Firstly, in terms of the point matching method, traditional ICP and feature point-based ICP both use the nearest neighbor method to find the nearest point in the target point cloud to the registered point, while kd-tree-based ICP utilizes the kd-tree search principle to achieve point correspondence. Then, from the perspective of point sets participating in registration, traditional ICP and kd-tree-based ICP use all points for registration, while feature point-based ICP only selects the filtered feature points during the algorithm iterations. These two parameters constitute the principle of the algorithms and influence the results of the registration.Table 2
Major parameters of the ICP algorithms.
AlgorithmPoint matching methodRegistration point setMatrix calculationMinimization metricAccelerationTraditional ICPNearest neighbor searchingAll pointsUnit quaternionPoint-to-pointExtrapolateKd-tree-based ICPKd-tree-based searchingAll pointsUnit quaternionPoint-to-pointExtrapolateFeature point-based ICPNearest neighbor searchingFeature pointsUnit quaternionPoint-to-pointExtrapolateTo accurately represent the impact of the unique operating principles of the three algorithms on the registration results, the parameters in the iteration process are set to be consistent except for the two parameters mentioned above. As listed in Table2, all three algorithms employ the unit quaternion method for calculating transformation matrices, and the point-to-point method is chosen as the minimization metric between point clouds. Finally, all three algorithms utilize extrapolation for acceleration.In this section, the preprocessing of the point clouds is first performed to serve a good and consistent input for registration. Then, the registration results of the three algorithms are evaluated based on the point cloud data of the bridge towers acquired at different times in the following two aspects: (1) iteration process and RMS error, and (2) transformation matrix and measured deformation results.


### 4.1. Preprocessing for Registration

The three ICP algorithms are considered fine registration algorithms due to their pairwise point matching properties. Generally, the fine registration highly depends on the accuracy of the initial input. In other words, the degree of coarse registration between the source point cloud and the target point cloud may greatly affect the precision of the results of fine alignment. Therefore, many studies have been conducted on the automated coarse registration between two point cloud sets, and reliable registration results have been achieved.In this study, the 4PCS–based algorithm [27] is selected to roughly align the point clouds, which sample the four points with an approximately coplanar layout as a base set. This algorithm extracts the cross-sections at the same location at the base of the bridge tower in both point clouds, and the four corner points in the cross-section are treated as a registered point set. After applying the resulting transformation matrix to the bridge tower target point cloud, the initial pose for the fine registration is established. Figure 6(a) shows the original position of the bridge tower point clouds, and Figure 6(b) illustrates the two point clouds that have been roughly aligned. The x, y, and z axes represent the lateral, longitudinal, and vertical directions of the bridge, respectively.Figure 6
Bridge tower point clouds: (a) before the coarse registration and (b) after the coarse registration.
(a)(b)Besides, point cloud density is also a factor that affects registration results in preprocessing. Point cloud density is defined as the average point spacing of a point cloud, which determines the number of points that can be accommodated in a defined 3D space. This factor determines the operation time of the algorithm due to the pairwise calculation in each iteration. Meanwhile, the degree of similarity between the densities of the two point clouds has an impact on the matching of the corresponding points. In other words, the closer the density of the point clouds, the higher the accuracy of the point pair matching. Therefore, the authors adjust the density of the point clouds to be approximately the same in the generic postprocessing software by modifying the point spacing.The operations above of coarse registration and point cloud density uniformity are conducted in CloudCompare. CloudCompare is a 3D point cloud processing software based on the GPL open source license, which can realize point cloud import and export, view and edit, and other functions related to the processing of the point cloud.


### 4.2. Registration Results from Three ICP Algorithms

To assess the performance of algorithms during the iteration, the root mean square (RMS) error concerning the number of iterations is employed, which is defined in the following equation:(10)ek=1Na∑i=1Naaik−bik2,where ek is the RMS error after k iterations, aik and bik are the corresponding point pairs in the source point cloud A and target point cloud B of the kth iteration, respectively, and Na is the number of the participated points during the iteration.As the transformation matrix between two point clouds reflects the direct iteration result after algorithm convergences, it is also chosen for algorithm comparison. Besides, to verify the reliability of the measured bridge deformation, the displacement values of the top of the bridge tower in three algorithms are utilized. The program implementation of the three ICP algorithms and the acquisition of the displacement values of the bridge tower top is based on the MATLAB programming language.


#### 4.2.1. Iteration Process and RMS Error

Based on the convergence theorem [30], the ICP algorithm converges monotonically to a local minimum with the mean square distance objective function. This means that after each iteration, the least squares registration generically reduces the average distance between the corresponded point pairs, and the average distance is usually expressed by the RMS error in equation (10). The smaller RMS error for the same number of iterations reflects the greater ability of the algorithm to overlap point clouds. Meanwhile, the RMS error reduction rates and final values vary among different ICP algorithms due to their different iterative approaches. Figure 7 illustrates the registration point set and point matching in the iteration process of the three algorithms. The traditional ICP algorithm selects all points as the registration point set and matches the points with the nearest point searching (Figure 7(a)). The kd-tree-based ICP algorithm also utilizes all points in point cloud data for registration but uses the kd-tree-based approach as the point matching method (Figure 7(b)). The feature point-based ICP algorithm only employs the filtered feature points as the registration point set with the same point matching method as the traditional ICP algorithm (Figure 7(c)).Figure 7
The registration point set and the point matching of (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)To observe the convergence trend of the RMS error during the iteration process and the final value after convergence, this section visualizes how the RMS error varies with the iteration process, as shown in Figures8 and 9. Without loss of generality, sufficient iteration numbers are utilized for three algorithms. As shown in Figure 8, the RMS errors from three ICP algorithms that vary with the operation time of the algorithm under the same iteration number, in which the preprocessing time has already been eliminated. Each node in the curve represents the beginning of an iteration. For a given number of iterations (NoI), the total operation time (TOT) is 950 s, 29 s, and 0.88 s for the traditional ICP algorithm (Figure 8(a)), the kd-tree-based ICP algorithm (Figure 8(b)), and the feature point-based ICP algorithm (Figure 8(c)), respectively. The average execution time for a single iteration of the three algorithms is 31.67 s, 0.97 s, and 0.03 s, in that order.Figure 8
The operation time and convergence speed for (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)Figure 9
RMS errors with three ICP algorithms.More apparently, Figure8(c) indicates that the convergence operation time (COT) of the feature point-based ICP algorithm is only 0.15 s, while the kd-tree-based ICP algorithm takes nearly 23 s (Figure 8(b)), and the traditional ICP algorithm uses 782 s (Figure 8(a)). Therefore, it can be concluded that the feature point-based ICP algorithm has the shortest operation time and the fastest convergence speed. Although the other two algorithms converge at the same iteration step, the kd-tree-based ICP algorithm converges faster than the traditional ICP algorithm due to its shorter execution time for a single iteration.The converged RMS error represents the maximum overlap of two point clouds that the algorithms can achieve. Figure9 plots the RMS error variation curves of the three algorithms under the same coordinate system, with the number of iterations as the horizontal coordinate. To achieve comparability, like the other two algorithms, the feature point-based ICP algorithm employs all points to calculate the RMS error, even though it uses only the filtered feature points in the registration.It can be observed in Figure9 that the two curves of the traditional ICP algorithm and kd-tree-based ICP algorithm 2 completely coincide, which shows that the kd-tree method does not help to reduce the RMS error. Figure 9 also indicates that the feature point-based algorithm only achieves an RMS error of 0.134, compared to the other two algorithms that can reduce the RMS error to 0.075 after convergence. The feature point-based ICP algorithm is not as good as the other two algorithms at overlapping point clouds.


#### 4.2.2. Transformation Matrix and Measured Bridge Deformation

When the RMS error converges to a stable value, the transformation matrix can be obtained to represent the final registration result. Therefore, the translation and rotation values (Tx,Ty,Tz,Rx,Ry,Rz) reflected by the elements in the final transformation matrix are considered to be the ideal deformation measurement results. Based on the correlation formula written in equations (11) and (12), the values of translation and rotation that change with the number of iterations are shown in Figures 10 and 11, where TM is the final transformation matrix, Tx, Ty, and Tz are the distances that the point cloud travels in the positive direction along the three axes, Rx, Ry, and Rz are the angles that the point cloud rotates clockwise around the three axes. As anticipated, the spatial deformation results in the six dimensions of all three algorithms converge after a certain iteration step.Figure 10
Translation parameters for the number of iterations with three ICP algorithms: (a)Tx, (b)Ty, and (c)Tz.
(a)(b)(c)Figure 11
Rotation parameters for the number of iterations with three ICP algorithms: (a)Rx, (b)Ry, and (c)Rz.
(a)(b)(c)Unlike the monotonic form of change in RMS error, some fluctuations in the trend of deformation values can be observed before convergence, such asTy in Figure 10(b) and Tz in Figure 10(c) for feature point-based ICP algorithms, Ry in Figure 11(b) for the traditional ICP algorithm and kd-tree-based ICP algorithm. At the same time, there is also a difference in the number of iterations required for convergence of the deformation values of the same ICP algorithm, but none is greater than the number of convergence iterations of the RMS error in Section 4.2.1. The number of iterations to reach convergence for the parameters of the three algorithms is summarized in Table 3.(11)TM=a11,a12,a21,a22,a13,a14,a23,a24,a31,a32,0,0,a33,a34,0,1,(12)Tx=a14,Ty=a24,Tz=a34Rx=atan2a32,a33Ry=atan2−a31,a322+a332,Rz=atan2a21,a11.Table 3
The convergence numbers of major parameters of three algorithms.
Algorithm nameThe convergence numbersRMS errorTxTyTzRxRyRzTraditional ICP249102411910Kd-tree-based ICP249102411910Feature point-based ICP5555554Based on the conversion relations in equation (12) and the spatial transformation matrices at the convergence in the following equations, change values for the six dimensions reflected by the transformation matrix are calculated, as given in Table 4. These parameters represent the difference in spatial position between the two point clouds calculated by the three different algorithms.Table 4
Deformation for the six dimensions reflected by the transformation matrices.
AlgorithmTx (m)Ty (m)Tz (m)Rx (°)Ry (°)Rz (°)Traditional ICP0.61330.05380.0688−0.5900.002−0.001Kd-tree-based ICP0.61330.05380.0688−0.5900.002−0.001Feature point-based ICP0.23060.00640.0437−0.0690.034−0.034For the traditional ICP algorithm and kd-tree-based ICP algorithm, the entire point cloud of the bridge tower is registered as a rigid body. Therefore, the deformation values using the traditional ICP algorithm and kd-tree-based ICP algorithm in Table4 are the rigid body deformation results of the bridge tower calculated by these two algorithms. On the other hand, the feature point-based ICP algorithm, compared to the feature points in the source point cloud, the feature points in the target point cloud are only rigidly deformed by the overall deformation of the bridge but not internally deformed by external loads [38]. The overall bridge deformation is interpreted as a stress-free change in the relative position between the bridge and the environment. In contrast, internal bridge deformation is a stress deformation caused by external forces. In this case, the deformation results in Table 4 only denote the difference in spatial location of the rigid part in the bridge tower caused by the overall bridge deformation. Therefore, bridge tower deformation identification of the feature point-based ICP algorithm firstly requires the transformation matrix in equation (15) to eliminate the effect of the overall bridge deformation, and then calculate the distance between nonfeature points.(13)TMtraditionalICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(14)TMkd−tree−basedICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(15)TMfeaturepoint−basedICP=1.000000.00380−0.003801.000000.000500.230600.000700.00640−0.00050−0.000700.000000.000001.000000.043700.000001.00000.Since the point cloud data acquired by terrestrial laser scanning is perpendicular to a common horizontal plane in the vertical direction, the difference in tilt angle between any two point clouds in the vertical direction is the same in all cases, which gives comparability of deformation results between different algorithms regarding Rx and Ry. At the same time, the bridge tower top displacement in the lateral and longitudinal deflections is exactly the key index for structural health monitoring of cable-supported bridges; thus, this study selects the top displacement of the bridge tower identified by three algorithms as the final deformation result.The lateral and longitudinal deformations of the bridge tower identified by using three ICP-related algorithms are summarized in Table5 and are compared with the observation data measured by using total station. Considering the total station measurements as the baseline, slight difference can be observed of the lateral displacement of the bridge tower, i.e., 0.063 m∼0.101 m and 0.086 m identified by ICP-related methods and total station, respectively. However, a significant difference can be observed in terms of the longitudinal bridge tower top displacement. The displacement recognized by using the traditional ICP algorithm and kd-tree-based ICP algorithm is 1.752 m, which is approximately eight times greater than the feature point-based ICP algorithm, i.e., 0.203 m. Compared with the baseline measurement, i.e., 0.167 m recorded by using a total station, the outperformance of the featured point-based ICP can be demonstrated; for instance, only 0.036 m differs from the baseline measurement, while 1.585 m offset can be observed when using traditional ICP or kd-tree-based ICP.Table 5
Bridge tower spatial deformation.
MethodsBridge tower displacementLongitudinal direction (m)Lateral direction (m)Traditional ICP1.7520.063Kd-tree-based ICP1.7520.063Feature point-based ICP0.2030.101Total station measurement0.1670.086


## 4.1. Preprocessing for Registration

The three ICP algorithms are considered fine registration algorithms due to their pairwise point matching properties. Generally, the fine registration highly depends on the accuracy of the initial input. In other words, the degree of coarse registration between the source point cloud and the target point cloud may greatly affect the precision of the results of fine alignment. Therefore, many studies have been conducted on the automated coarse registration between two point cloud sets, and reliable registration results have been achieved.In this study, the 4PCS–based algorithm [27] is selected to roughly align the point clouds, which sample the four points with an approximately coplanar layout as a base set. This algorithm extracts the cross-sections at the same location at the base of the bridge tower in both point clouds, and the four corner points in the cross-section are treated as a registered point set. After applying the resulting transformation matrix to the bridge tower target point cloud, the initial pose for the fine registration is established. Figure 6(a) shows the original position of the bridge tower point clouds, and Figure 6(b) illustrates the two point clouds that have been roughly aligned. The x, y, and z axes represent the lateral, longitudinal, and vertical directions of the bridge, respectively.Figure 6
Bridge tower point clouds: (a) before the coarse registration and (b) after the coarse registration.
(a)(b)Besides, point cloud density is also a factor that affects registration results in preprocessing. Point cloud density is defined as the average point spacing of a point cloud, which determines the number of points that can be accommodated in a defined 3D space. This factor determines the operation time of the algorithm due to the pairwise calculation in each iteration. Meanwhile, the degree of similarity between the densities of the two point clouds has an impact on the matching of the corresponding points. In other words, the closer the density of the point clouds, the higher the accuracy of the point pair matching. Therefore, the authors adjust the density of the point clouds to be approximately the same in the generic postprocessing software by modifying the point spacing.The operations above of coarse registration and point cloud density uniformity are conducted in CloudCompare. CloudCompare is a 3D point cloud processing software based on the GPL open source license, which can realize point cloud import and export, view and edit, and other functions related to the processing of the point cloud.


## 4.2. Registration Results from Three ICP Algorithms

To assess the performance of algorithms during the iteration, the root mean square (RMS) error concerning the number of iterations is employed, which is defined in the following equation:(10)ek=1Na∑i=1Naaik−bik2,where ek is the RMS error after k iterations, aik and bik are the corresponding point pairs in the source point cloud A and target point cloud B of the kth iteration, respectively, and Na is the number of the participated points during the iteration.As the transformation matrix between two point clouds reflects the direct iteration result after algorithm convergences, it is also chosen for algorithm comparison. Besides, to verify the reliability of the measured bridge deformation, the displacement values of the top of the bridge tower in three algorithms are utilized. The program implementation of the three ICP algorithms and the acquisition of the displacement values of the bridge tower top is based on the MATLAB programming language.


### 4.2.1. Iteration Process and RMS Error

Based on the convergence theorem [30], the ICP algorithm converges monotonically to a local minimum with the mean square distance objective function. This means that after each iteration, the least squares registration generically reduces the average distance between the corresponded point pairs, and the average distance is usually expressed by the RMS error in equation (10). The smaller RMS error for the same number of iterations reflects the greater ability of the algorithm to overlap point clouds. Meanwhile, the RMS error reduction rates and final values vary among different ICP algorithms due to their different iterative approaches. Figure 7 illustrates the registration point set and point matching in the iteration process of the three algorithms. The traditional ICP algorithm selects all points as the registration point set and matches the points with the nearest point searching (Figure 7(a)). The kd-tree-based ICP algorithm also utilizes all points in point cloud data for registration but uses the kd-tree-based approach as the point matching method (Figure 7(b)). The feature point-based ICP algorithm only employs the filtered feature points as the registration point set with the same point matching method as the traditional ICP algorithm (Figure 7(c)).Figure 7
The registration point set and the point matching of (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)To observe the convergence trend of the RMS error during the iteration process and the final value after convergence, this section visualizes how the RMS error varies with the iteration process, as shown in Figures8 and 9. Without loss of generality, sufficient iteration numbers are utilized for three algorithms. As shown in Figure 8, the RMS errors from three ICP algorithms that vary with the operation time of the algorithm under the same iteration number, in which the preprocessing time has already been eliminated. Each node in the curve represents the beginning of an iteration. For a given number of iterations (NoI), the total operation time (TOT) is 950 s, 29 s, and 0.88 s for the traditional ICP algorithm (Figure 8(a)), the kd-tree-based ICP algorithm (Figure 8(b)), and the feature point-based ICP algorithm (Figure 8(c)), respectively. The average execution time for a single iteration of the three algorithms is 31.67 s, 0.97 s, and 0.03 s, in that order.Figure 8
The operation time and convergence speed for (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)Figure 9
RMS errors with three ICP algorithms.More apparently, Figure8(c) indicates that the convergence operation time (COT) of the feature point-based ICP algorithm is only 0.15 s, while the kd-tree-based ICP algorithm takes nearly 23 s (Figure 8(b)), and the traditional ICP algorithm uses 782 s (Figure 8(a)). Therefore, it can be concluded that the feature point-based ICP algorithm has the shortest operation time and the fastest convergence speed. Although the other two algorithms converge at the same iteration step, the kd-tree-based ICP algorithm converges faster than the traditional ICP algorithm due to its shorter execution time for a single iteration.The converged RMS error represents the maximum overlap of two point clouds that the algorithms can achieve. Figure9 plots the RMS error variation curves of the three algorithms under the same coordinate system, with the number of iterations as the horizontal coordinate. To achieve comparability, like the other two algorithms, the feature point-based ICP algorithm employs all points to calculate the RMS error, even though it uses only the filtered feature points in the registration.It can be observed in Figure9 that the two curves of the traditional ICP algorithm and kd-tree-based ICP algorithm 2 completely coincide, which shows that the kd-tree method does not help to reduce the RMS error. Figure 9 also indicates that the feature point-based algorithm only achieves an RMS error of 0.134, compared to the other two algorithms that can reduce the RMS error to 0.075 after convergence. The feature point-based ICP algorithm is not as good as the other two algorithms at overlapping point clouds.


### 4.2.2. Transformation Matrix and Measured Bridge Deformation

When the RMS error converges to a stable value, the transformation matrix can be obtained to represent the final registration result. Therefore, the translation and rotation values (Tx,Ty,Tz,Rx,Ry,Rz) reflected by the elements in the final transformation matrix are considered to be the ideal deformation measurement results. Based on the correlation formula written in equations (11) and (12), the values of translation and rotation that change with the number of iterations are shown in Figures 10 and 11, where TM is the final transformation matrix, Tx, Ty, and Tz are the distances that the point cloud travels in the positive direction along the three axes, Rx, Ry, and Rz are the angles that the point cloud rotates clockwise around the three axes. As anticipated, the spatial deformation results in the six dimensions of all three algorithms converge after a certain iteration step.Figure 10
Translation parameters for the number of iterations with three ICP algorithms: (a)Tx, (b)Ty, and (c)Tz.
(a)(b)(c)Figure 11
Rotation parameters for the number of iterations with three ICP algorithms: (a)Rx, (b)Ry, and (c)Rz.
(a)(b)(c)Unlike the monotonic form of change in RMS error, some fluctuations in the trend of deformation values can be observed before convergence, such asTy in Figure 10(b) and Tz in Figure 10(c) for feature point-based ICP algorithms, Ry in Figure 11(b) for the traditional ICP algorithm and kd-tree-based ICP algorithm. At the same time, there is also a difference in the number of iterations required for convergence of the deformation values of the same ICP algorithm, but none is greater than the number of convergence iterations of the RMS error in Section 4.2.1. The number of iterations to reach convergence for the parameters of the three algorithms is summarized in Table 3.(11)TM=a11,a12,a21,a22,a13,a14,a23,a24,a31,a32,0,0,a33,a34,0,1,(12)Tx=a14,Ty=a24,Tz=a34Rx=atan2a32,a33Ry=atan2−a31,a322+a332,Rz=atan2a21,a11.Table 3
The convergence numbers of major parameters of three algorithms.
Algorithm nameThe convergence numbersRMS errorTxTyTzRxRyRzTraditional ICP249102411910Kd-tree-based ICP249102411910Feature point-based ICP5555554Based on the conversion relations in equation (12) and the spatial transformation matrices at the convergence in the following equations, change values for the six dimensions reflected by the transformation matrix are calculated, as given in Table 4. These parameters represent the difference in spatial position between the two point clouds calculated by the three different algorithms.Table 4
Deformation for the six dimensions reflected by the transformation matrices.
AlgorithmTx (m)Ty (m)Tz (m)Rx (°)Ry (°)Rz (°)Traditional ICP0.61330.05380.0688−0.5900.002−0.001Kd-tree-based ICP0.61330.05380.0688−0.5900.002−0.001Feature point-based ICP0.23060.00640.0437−0.0690.034−0.034For the traditional ICP algorithm and kd-tree-based ICP algorithm, the entire point cloud of the bridge tower is registered as a rigid body. Therefore, the deformation values using the traditional ICP algorithm and kd-tree-based ICP algorithm in Table4 are the rigid body deformation results of the bridge tower calculated by these two algorithms. On the other hand, the feature point-based ICP algorithm, compared to the feature points in the source point cloud, the feature points in the target point cloud are only rigidly deformed by the overall deformation of the bridge but not internally deformed by external loads [38]. The overall bridge deformation is interpreted as a stress-free change in the relative position between the bridge and the environment. In contrast, internal bridge deformation is a stress deformation caused by external forces. In this case, the deformation results in Table 4 only denote the difference in spatial location of the rigid part in the bridge tower caused by the overall bridge deformation. Therefore, bridge tower deformation identification of the feature point-based ICP algorithm firstly requires the transformation matrix in equation (15) to eliminate the effect of the overall bridge deformation, and then calculate the distance between nonfeature points.(13)TMtraditionalICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(14)TMkd−tree−basedICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(15)TMfeaturepoint−basedICP=1.000000.00380−0.003801.000000.000500.230600.000700.00640−0.00050−0.000700.000000.000001.000000.043700.000001.00000.Since the point cloud data acquired by terrestrial laser scanning is perpendicular to a common horizontal plane in the vertical direction, the difference in tilt angle between any two point clouds in the vertical direction is the same in all cases, which gives comparability of deformation results between different algorithms regarding Rx and Ry. At the same time, the bridge tower top displacement in the lateral and longitudinal deflections is exactly the key index for structural health monitoring of cable-supported bridges; thus, this study selects the top displacement of the bridge tower identified by three algorithms as the final deformation result.The lateral and longitudinal deformations of the bridge tower identified by using three ICP-related algorithms are summarized in Table5 and are compared with the observation data measured by using total station. Considering the total station measurements as the baseline, slight difference can be observed of the lateral displacement of the bridge tower, i.e., 0.063 m∼0.101 m and 0.086 m identified by ICP-related methods and total station, respectively. However, a significant difference can be observed in terms of the longitudinal bridge tower top displacement. The displacement recognized by using the traditional ICP algorithm and kd-tree-based ICP algorithm is 1.752 m, which is approximately eight times greater than the feature point-based ICP algorithm, i.e., 0.203 m. Compared with the baseline measurement, i.e., 0.167 m recorded by using a total station, the outperformance of the featured point-based ICP can be demonstrated; for instance, only 0.036 m differs from the baseline measurement, while 1.585 m offset can be observed when using traditional ICP or kd-tree-based ICP.Table 5
Bridge tower spatial deformation.
MethodsBridge tower displacementLongitudinal direction (m)Lateral direction (m)Traditional ICP1.7520.063Kd-tree-based ICP1.7520.063Feature point-based ICP0.2030.101Total station measurement0.1670.086


## 4.2.1. Iteration Process and RMS Error

Based on the convergence theorem [30], the ICP algorithm converges monotonically to a local minimum with the mean square distance objective function. This means that after each iteration, the least squares registration generically reduces the average distance between the corresponded point pairs, and the average distance is usually expressed by the RMS error in equation (10). The smaller RMS error for the same number of iterations reflects the greater ability of the algorithm to overlap point clouds. Meanwhile, the RMS error reduction rates and final values vary among different ICP algorithms due to their different iterative approaches. Figure 7 illustrates the registration point set and point matching in the iteration process of the three algorithms. The traditional ICP algorithm selects all points as the registration point set and matches the points with the nearest point searching (Figure 7(a)). The kd-tree-based ICP algorithm also utilizes all points in point cloud data for registration but uses the kd-tree-based approach as the point matching method (Figure 7(b)). The feature point-based ICP algorithm only employs the filtered feature points as the registration point set with the same point matching method as the traditional ICP algorithm (Figure 7(c)).Figure 7
The registration point set and the point matching of (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)To observe the convergence trend of the RMS error during the iteration process and the final value after convergence, this section visualizes how the RMS error varies with the iteration process, as shown in Figures8 and 9. Without loss of generality, sufficient iteration numbers are utilized for three algorithms. As shown in Figure 8, the RMS errors from three ICP algorithms that vary with the operation time of the algorithm under the same iteration number, in which the preprocessing time has already been eliminated. Each node in the curve represents the beginning of an iteration. For a given number of iterations (NoI), the total operation time (TOT) is 950 s, 29 s, and 0.88 s for the traditional ICP algorithm (Figure 8(a)), the kd-tree-based ICP algorithm (Figure 8(b)), and the feature point-based ICP algorithm (Figure 8(c)), respectively. The average execution time for a single iteration of the three algorithms is 31.67 s, 0.97 s, and 0.03 s, in that order.Figure 8
The operation time and convergence speed for (a) the traditional ICP algorithm, (b) kd-tree-based ICP algorithm, and (c) feature point-based ICP algorithm.
(a)(b)(c)Figure 9
RMS errors with three ICP algorithms.More apparently, Figure8(c) indicates that the convergence operation time (COT) of the feature point-based ICP algorithm is only 0.15 s, while the kd-tree-based ICP algorithm takes nearly 23 s (Figure 8(b)), and the traditional ICP algorithm uses 782 s (Figure 8(a)). Therefore, it can be concluded that the feature point-based ICP algorithm has the shortest operation time and the fastest convergence speed. Although the other two algorithms converge at the same iteration step, the kd-tree-based ICP algorithm converges faster than the traditional ICP algorithm due to its shorter execution time for a single iteration.The converged RMS error represents the maximum overlap of two point clouds that the algorithms can achieve. Figure9 plots the RMS error variation curves of the three algorithms under the same coordinate system, with the number of iterations as the horizontal coordinate. To achieve comparability, like the other two algorithms, the feature point-based ICP algorithm employs all points to calculate the RMS error, even though it uses only the filtered feature points in the registration.It can be observed in Figure9 that the two curves of the traditional ICP algorithm and kd-tree-based ICP algorithm 2 completely coincide, which shows that the kd-tree method does not help to reduce the RMS error. Figure 9 also indicates that the feature point-based algorithm only achieves an RMS error of 0.134, compared to the other two algorithms that can reduce the RMS error to 0.075 after convergence. The feature point-based ICP algorithm is not as good as the other two algorithms at overlapping point clouds.


## 4.2.2. Transformation Matrix and Measured Bridge Deformation

When the RMS error converges to a stable value, the transformation matrix can be obtained to represent the final registration result. Therefore, the translation and rotation values (Tx,Ty,Tz,Rx,Ry,Rz) reflected by the elements in the final transformation matrix are considered to be the ideal deformation measurement results. Based on the correlation formula written in equations (11) and (12), the values of translation and rotation that change with the number of iterations are shown in Figures 10 and 11, where TM is the final transformation matrix, Tx, Ty, and Tz are the distances that the point cloud travels in the positive direction along the three axes, Rx, Ry, and Rz are the angles that the point cloud rotates clockwise around the three axes. As anticipated, the spatial deformation results in the six dimensions of all three algorithms converge after a certain iteration step.Figure 10
Translation parameters for the number of iterations with three ICP algorithms: (a)Tx, (b)Ty, and (c)Tz.
(a)(b)(c)Figure 11
Rotation parameters for the number of iterations with three ICP algorithms: (a)Rx, (b)Ry, and (c)Rz.
(a)(b)(c)Unlike the monotonic form of change in RMS error, some fluctuations in the trend of deformation values can be observed before convergence, such asTy in Figure 10(b) and Tz in Figure 10(c) for feature point-based ICP algorithms, Ry in Figure 11(b) for the traditional ICP algorithm and kd-tree-based ICP algorithm. At the same time, there is also a difference in the number of iterations required for convergence of the deformation values of the same ICP algorithm, but none is greater than the number of convergence iterations of the RMS error in Section 4.2.1. The number of iterations to reach convergence for the parameters of the three algorithms is summarized in Table 3.(11)TM=a11,a12,a21,a22,a13,a14,a23,a24,a31,a32,0,0,a33,a34,0,1,(12)Tx=a14,Ty=a24,Tz=a34Rx=atan2a32,a33Ry=atan2−a31,a322+a332,Rz=atan2a21,a11.Table 3
The convergence numbers of major parameters of three algorithms.
Algorithm nameThe convergence numbersRMS errorTxTyTzRxRyRzTraditional ICP249102411910Kd-tree-based ICP249102411910Feature point-based ICP5555554Based on the conversion relations in equation (12) and the spatial transformation matrices at the convergence in the following equations, change values for the six dimensions reflected by the transformation matrix are calculated, as given in Table 4. These parameters represent the difference in spatial position between the two point clouds calculated by the three different algorithms.Table 4
Deformation for the six dimensions reflected by the transformation matrices.
AlgorithmTx (m)Ty (m)Tz (m)Rx (°)Ry (°)Rz (°)Traditional ICP0.61330.05380.0688−0.5900.002−0.001Kd-tree-based ICP0.61330.05380.0688−0.5900.002−0.001Feature point-based ICP0.23060.00640.0437−0.0690.034−0.034For the traditional ICP algorithm and kd-tree-based ICP algorithm, the entire point cloud of the bridge tower is registered as a rigid body. Therefore, the deformation values using the traditional ICP algorithm and kd-tree-based ICP algorithm in Table4 are the rigid body deformation results of the bridge tower calculated by these two algorithms. On the other hand, the feature point-based ICP algorithm, compared to the feature points in the source point cloud, the feature points in the target point cloud are only rigidly deformed by the overall deformation of the bridge but not internally deformed by external loads [38]. The overall bridge deformation is interpreted as a stress-free change in the relative position between the bridge and the environment. In contrast, internal bridge deformation is a stress deformation caused by external forces. In this case, the deformation results in Table 4 only denote the difference in spatial location of the rigid part in the bridge tower caused by the overall bridge deformation. Therefore, bridge tower deformation identification of the feature point-based ICP algorithm firstly requires the transformation matrix in equation (15) to eliminate the effect of the overall bridge deformation, and then calculate the distance between nonfeature points.(13)TMtraditionalICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(14)TMkd−tree−basedICP=1.000000.01000−0.010301.000000.000040.613300.000020.05380−0.00004−0.000020.000000.000001.000000.068800.000001.00000,(15)TMfeaturepoint−basedICP=1.000000.00380−0.003801.000000.000500.230600.000700.00640−0.00050−0.000700.000000.000001.000000.043700.000001.00000.Since the point cloud data acquired by terrestrial laser scanning is perpendicular to a common horizontal plane in the vertical direction, the difference in tilt angle between any two point clouds in the vertical direction is the same in all cases, which gives comparability of deformation results between different algorithms regarding Rx and Ry. At the same time, the bridge tower top displacement in the lateral and longitudinal deflections is exactly the key index for structural health monitoring of cable-supported bridges; thus, this study selects the top displacement of the bridge tower identified by three algorithms as the final deformation result.The lateral and longitudinal deformations of the bridge tower identified by using three ICP-related algorithms are summarized in Table5 and are compared with the observation data measured by using total station. Considering the total station measurements as the baseline, slight difference can be observed of the lateral displacement of the bridge tower, i.e., 0.063 m∼0.101 m and 0.086 m identified by ICP-related methods and total station, respectively. However, a significant difference can be observed in terms of the longitudinal bridge tower top displacement. The displacement recognized by using the traditional ICP algorithm and kd-tree-based ICP algorithm is 1.752 m, which is approximately eight times greater than the feature point-based ICP algorithm, i.e., 0.203 m. Compared with the baseline measurement, i.e., 0.167 m recorded by using a total station, the outperformance of the featured point-based ICP can be demonstrated; for instance, only 0.036 m differs from the baseline measurement, while 1.585 m offset can be observed when using traditional ICP or kd-tree-based ICP.Table 5
Bridge tower spatial deformation.
MethodsBridge tower displacementLongitudinal direction (m)Lateral direction (m)Traditional ICP1.7520.063Kd-tree-based ICP1.7520.063Feature point-based ICP0.2030.101Total station measurement0.1670.086


## 5. Discussion

The results in Section4.2.1 indicate that the traditional ICP algorithm and the kd-tree-based ICP algorithm can achieve a greater point cloud overlap than the feature point-based ICP algorithm. The potential reason for this finding is that the traditional ICP algorithm and the kd-tree-based ICP algorithm take all points in the point cloud as alignment targets, while the feature point-based ICP algorithm only targets the filtered feature points that are only part of the overall point cloud, which ignores the alignment of nonfeature points. The smaller number of registration points is also the reason that the feature point-based ICP algorithm converges much faster than the other two algorithms. The relationship between the RMS error and the number of iterations for the traditional ICP algorithm and the kd-tree-based ICP algorithm is identical due to the same registration object. However, the operation time of the kd-tree-based ICP algorithm is much smaller than the kd-tree-based ICP algorithm because of the better organization of the point cloud data.Another observation is that the parameters of the spatial transformation matrix do not always converge in the form of a monotonic smooth curve, whereas the RMS error iteration curve remains monotonically decreasing until it converges. A convergence characteristic of the ICP algorithm is summarized as follows: monotonic convergence for the mean square distance function and not necessarily monotonic convergence for all the independent spatial transformation parameters.The difference between the feature point-based ICP algorithm and the other two algorithms regarding the deformation values on top of the bridge towers can be observed in Table5. According to the response on a suspension bridge with a main span of 1080 m under common live loads based on the finite element method (FEM) [44], the concrete tower top longitudinal offset should range from 2.9 mm to 258.6 mm for a height of approximately 170 m. Therefore, the deformation identified by using the feature point-based ICP algorithm is acceptable, i.e., 0.203 m, while the results identified using the other two algorithms are considered to be unreasonable results, i.e., 1.752 m. Meanwhile, the lateral and longitudinal displacement values of the feature point-based algorithm are closer to the data measured by using the total station when compared with the other two ICP-related algorithms. Thus, it can be inferred that the results of the feature point-based algorithm have a high confidence level.The reason that the feature point-based ICP outperforms the other two is dependent on the set of points selected for registration, which will be treated as rigid bodies. The characteristic of a rigid body is that it can only generate rigid deformation without internal forces. When dealing with the bridge tower point data, the featured point-based ICP extracted a certain set of points for registration, while the other two utilized the whole bridge tower points for registration. In fact, due to the constraint effect of the pile cap and the pulling effect of the unequal horizontal combined force of the left and right main cables, the longitudinal deformation of the tower is similar to that of a bent cantilever beam, as shown in Figure12, where FL and FR are the unequal horizontal forces provided by the left and right main cables, respectively. Unlike rigid deformation, this form of deformation generates internal forces in the structure, and it is more significant that only the segment within a limited length from the base of the bridge tower can be approximated as a rigid body. Therefore, only the portion of the point cloud that conforms to the motion characteristics of the rigid body is extracted by a point filtering rule in the feature point-based ICP algorithm. This particular registration is more in line with the deformation characteristics of the bridge tower, resulting in the obtained tower top displacement values that are closer to the actual value. On the contrary, registering the entire bridge tower as a rigid body in the traditional ICP algorithm and the kd-tree-based ICP algorithm may lead to unreliable tower top displacements.Figure 12
Longitudinal deformation of the tower under the horizontal force.It can be observed that the lateral displacement results of the three algorithms are less different from the total station measured results, which is because the lateral displacement effect of the bridge towers is generally not significant [44], resulting in the bridge tower point cloud being registered more like a rigid body in the lateral direction. It also reflects the importance of exploiting the rigid part in the registration of point clouds for long-span bridges.Previous studies of the ICP algorithm have been devoted to minimizing the RMS error of the two point clouds, which is one of the most critical metrics for the evaluation of the algorithm. However, in this study, the feature point-based ICP algorithm with a higher RMS error yields more reasonable deformation recognition results, while the other two algorithms with a lower RMS error give results that exceed the normal range. Hence, in addition to the concern for RMS error, the improvement of the ICP algorithm should focus on the deformation principle of the bridge structure. In other words, determining the points in the point cloud of bridge components that can be used as input to the ICP algorithm is an issue that must be taken into account when programming the algorithm.


## 6. Conclusions

This study evaluated the performance of different ICP algorithms for point cloud registration of long-span bridges. The traditional ICP algorithm and its two variants, the kd-tree-based ICP algorithm and the feature point-based ICP algorithm, were assessed in terms of convergence rate, execution time for a single iteration, RMS error between point clouds, transformation matrix, and precision of measured bridge deformation. Several conclusions can be highlighted, as follows:(1)
The convergence of the feature point-based ICP algorithm is the fastest among the three algorithms due to fewer registration points, while the other two algorithms converge at the same iteration step because of the same registration point set.(2)
Due to the efficient data organization capability, the execution time for a single-step iteration of the kd-tree-based algorithm is reduced compared to the traditional ICP algorithm. And the feature point-based ICP algorithm has the shortest execution time for a single-step iteration because of the dramatically reduced amount of point cloud data used for registration.(3)
In terms of the value of the RMS error between point clouds, the feature point-based algorithm is inferior to the other two algorithms in its ability to overlap two point clouds. The reason is that the feature point-based algorithm only registers the filtered feature points and ignores the overlap of nonfeature points.(4)
By comparing with the simulative results, the longitudinal displacement at the top of the tower should range from 2.9 mm to 258.6 mm. The longitudinal offset value of 203 mm from the tower top is within this range, thus, the feature point-based algorithm is considered to give an acceptable result. On the other hand, the longitudinal displacement results of the other two algorithms exceed the upper limit of the range by 577.5%, which are considered unreasonable results. Moreover, compared with the total station measurement, the outperformance of the feature point-based algorithm can also be demonstrated. The reason is that the feature point-based algorithm is more consistent with the deformation law of long-span bridges, which exploits the rigid part of the point cloud for registration, while other algorithms directly use all points in the point cloud to participate in the registration.From the assessment of the three algorithms for bridge tower deformation analysis, the two suggestions can be concluded as follows: (1) further efforts can be devoted to find the optimal point cloud management solutions to save the computational cost and (2) in terms of the reliability of the long-span bridge deformation analysis, the automatic and accurate selection of registration points that conform to the laws of the rigid body motion need to be taken into account in future research.



---

*Source: 2629418-2023-02-09.xml*