A NOVEL APPROACH TO REGISTER MULTI-PLATFORM POINT CLOUDS FOR ROCKFALL MONITORING

Point cloud produced from technologies such as terrestrial laser scanning (TLS) and photogrammetry (terrestrial and aerial) are widely used in rockfall monitoring applications due to the wealth of data they provide. In such applications, the acquisition and registration of multi-epoch point clouds is necessary. In addition, point clouds can be derived from different sensors (e.g., lasers versus digital cameras) and different platforms (terrestrial versus aerial). Therefore, registration methods should be able to support multi-platform datasets. Currently, registration of multi-platform datasets is done with manual intervention, and automatic registration is difficult. While registration of TLS point clouds can be achieved by targets that are not on the rock surface, this is not the case for photogrammetric methods, as ground control points (GCPs) should be located on the rock surface. Such GCPs can be lost or destroyed with time, and re-establishing them is difficult. Automated registration often relies on feature-based algorithms with refinement using the iterative closest point (ICP) algorithm. This paper presents a novel registration approach of multi-epoch and multi-platform point clouds to support rockfall monitoring applications. The registration method is based on edges that are detected in the different datasets using α-molecules. The paper shows application examples of the novel approach at different rock slopes in Colorado. Results demonstrate that the developed method in many cases performs better than the well-known ICP method and can be used to register point clouds and support rockfall monitoring.


INTRODUCTION
Point cloud collection methods such as terrestrial laser scanning (TLS), terrestrial photogrammetry (TP), and aerial photogrammetry from small unmanned aerial systems (sUASs) are widely being used in rockfall monitoring applications, as they offer datasets with high accuracy and spatial resolution. For monitoring applications and multi-epoch comparisons, their registration into a common coordinate system is important. In most cases, a single remote sensing method is being used, but when practical, two or more methods can be combined.
Registration of point clouds often relies on targets and/or ground control points (GCPs) or registration algorithms. For photogrammetrically derived point clouds, placement of targets on the rock surface is often necessary. This is dangerous and difficult, and targets may get damaged with time, requiring their re-establishment. For TLS, methods registration of multiple scans can be achieved using targets. These targets do not need to be on the rock surface, which is beneficial. Using point cloud registration algorithms is not trivial, as multi-epoch datasets might present differences in data densities, data gaps and overlap, and quality. The Iterative Closest Point (ICP) (Besl and McKay 1992) is often used, which offers several advantages. However, the input point clouds need to have a good initial registration, as the algorithm can become trapped in local minima (Attia and Slama 2017; Li et al. 2020). To address this issue, several ICP variants have been developed over the years in order to enhance registration performance (e.g., Bae and Lichti 2008;Wujanz et al. 2018;Kromer et al. 2019;Li et al. 2020). It is questionable whether such methods facilitate multiplatform registration, and often they cannot be applied in an automatic way. In addition, there have been algorithms developed specifically for TLS data, taking advantage of accurate instrument levelling with electronic compensators; therefore, rotations need only be solved in the horizontal plane (e.g., Cai et al. 2019;Zang et al. 2019). Methods that were designed for TLS setups are not suitable for multi-platform points clouds such as those derived from TLS, TP, and sUAS methods, as in the general case the input datasets will not be levelled. Thus, the availability of point clouds from several sensors necessitates the use of registration algorithms suitable to automatically register multi-platform point clouds. Featurebased registration methods that rely on distinct geometric features (points, lines, planes) are more suitable for the built environment, as finding such features is easy (e.g., building corners and outlines, road segments and sidewalks), but their application in complex rock surfaces is difficult (Abellan et al. 2014). For registration of rock-surfaces, feature-based methods build feature vectors using the 3D coordinates of points and describing the local surface properties in the vicinity of a point. Based on those feature vectors point correspondences are identified, which are then used to estimate the transformation between two input point clouds. Some examples of featurebased algorithms are the scale-invariant feature transform (Lowe 1999), intrinsic shape signatures (Zhong 2009), and point feature histograms (Rusu et al. 2009); however, these methods are often used for initial registration, followed by more precise fine registration with the ICP algorithm.
Discontinuity traces are distinct rock features, which can be considered as edges i.e., rapid changes of elevation or brightness if an image is used. They play an important role in slope stability analysis, as they define weak planes in a rock mass where rocks can detach and fail (Battulwar et al. 2021). Several methods have been developed to extract such discontinuity traces from rock surface, and a review of these different methods and approaches can be found in Battulwar et al. (2021). In addition, edges are often used in image analysis and image matching (Gruen 2012). The wavelet transform allows multi-scale analysis; however, wavelets are not directionally sensitive and limit their edge detection ability. For this reason, newer methods of contourlets, curvelets, and shearlets were developed. A recently developed method named α-molecules (Grohs et al. 2016), offers a unified framework for multi-scale transform, including wavelets, curvelets, and shearlets (Reisenhofer and King 2019).
A new automatic registration method, based on detected rock discontinuity traces (edges), was developed in Bolkas et al. (2021) to support multi-epoch and multi-platform registration for rockfall monitoring applications. This paper presents the main algorithm steps, and we demonstrate its application in a new and more complex rock slope dataset to demonstrate the robustness of the algorithm. The following sections present and discuss the main algorithm steps, the datasets used in this paper, and the results related to rockfall detection, finishing with conclusions.

METHODS
The automatic registration algorithm is based on edge detection and consists of five main processing steps (1) data preparation and coarse registration (2) edge detection (3) edge matching (4) point correspondences identification, (5) and fine registration. More information about the registration algorithm, including a thorough assessment of the individual steps, can be found in Bolkas et al. (2021). Note that the algorithm was developed for outcrops that generally have a 2.5D shape and not for fully 3D point clouds.
The input datasets are the source and target point clouds. The target point cloud is used as a reference for registration, which can be geo-referenced or in a local reference frame. The source point cloud is the dataset that needs to be transformed to match the coordinate system of the target point cloud. In many cases, this process requires finding three rotations, three translations, and one scale. The scale parameter can be neglected in the case where laser scanner point clouds are used, but it should be included when photogrammetric point clouds are used, either terrestrial or aerial. Numerical simulations with various scenarios of rotation, translation, and scale have shown that the algorithm can register point clouds within few mm to their original root mean square error (RMSE), even when the source point cloud has a scale of 0.5 with respect to the target point cloud (Bolkas et al. 2021). For instance, the original RMSE difference between two georeferenced sUAS and TP datasets from a rockfall site in Idaho Springs, Colorado was 2.1 cm. We applied a 500-m offset, 45° rotation, and 0.5 scale to the sUAS dataset, and after application of the proposed algorithm, the RMSE is 2.3 cm (for more information and more scenarios see Bolkas et al. 2021).
The main algorithm steps are summarized in Figure 1. The point clouds are rotated through principal component analysis (PCA) to make the outcrop span in the x-and y-axis directions. This is possible because point clouds of rock outcrops are typically 2.5D. This achieves a coarse registration between the two datasets. However, depending on the spatial extent of the two input datasets, there can be some significant translations in the x-and y-axis and a rotation about the z-axis. Coarse registration is enhanced using a phase correlation transformation (one rotation along the z-axis and two translations) (Dimitrievski et al. 2016). Approximate cropping of the two input datasets to the same spatial extent can facilitate the coarse registration performed in this step.
With a user input grid size value, the point clouds are gridded, and the z-axis direction is converted to a grayscale. Next, edge detection takes place using symmetric α-molecules (Reisenhofer and King 2019). The tuning parameters for the symmetric α-molecules are found empirically with an approach similar to Bolkas et al. (2018). The goal in this step is to detect similar / identical edges in both datasets. One of the most important input parameters is the maximum feature length, as low values (e.g., 2 m) can result in highly fragmented edges, which can increase the number of falsely matched edges. For both datasets in this paper, we have found that a maximum feature length of 4 m is suitable.
Edge matching is then performed using the discrete Fréchet distance (DFD) (Eiter and Mannila 1994). The DFD is a measure of similarity between two polygonal curves, and measures the shortest "leash" or coupling distance that is needed to connect and traverse the two curves. The computation of the DFD can consider curves with a different number of points, which is beneficial for the problem of edge matching in rock surfaces. Because the DFD is measured in the metric space, an initial registration at the level of 1-2 meters is needed. To speed up processing, the DFD is computed for source edges that are within a radius of the reference edges (e.g., 10 m). For the computation of the radius, we use the edge midpoints. Fine tuning the DFD threshold value allows the user to keep strong edge matches. Figure 2 shows an example of edge matching and the DFD distance computation. Edges are identified in both the target and source datasets. It is expected that, in general, similar edges will be detected. Figure 2 shows two examples of the DFD computations, the first between a target edge and the corresponding source edge, and the second between the same target edge and a different source edge. It is seen that strong matches will have a low DFD value, while poor matches will have a higher DFD value; thus, with simple thresholding we can keep the strong matches.
Point correspondences in the matched edges are extracted using the sum of squared differences. This step is implemented using the Matlab function "matchFeatures" (MathWorks, Inc. 2020). The function offers the use of threshold to keep strong matches between corresponding points. The registration is achieved through a rigid body transformation using a Procrustes analysis (Kendall 1989). A threshold of two times the standard deviation of the residuals is used to remove outliers, and the Procrustes analysis is repeated (Bolkas et al. 2021). Although this is a simple approach for outlier removal, results obtained thus far suggest that it is sufficient. In the future, a more robust approach for outlier detection will be added to enhance robustness against outliers. After this step, it is considered that the two input point clouds are registered to a common coordinate system. The datasets can be transformed back to the target coordinate system using the PCA rotation values of the target dataset.
Estimations of multi-epoch difference between point cloud datasets are determined using the model-to-model cloud comparison (M3C2) algorithm (Lague et al. 2013), as implemented in Cloud Compare (CloudCompare 2015) using a point spacing of 10 cm. From the M3C2 distances between datasets we compute RMSE values, which are used for assessment of registration. For comparison and validation of the developed registration method we provide registration results using the ICP algorithm as implemented in Cloud Compare. In the computation of the RMSE values, we have excluded points in vegetated areas using a manually made vegetation mask, which was created using the photogrammetric datasets.

DATASETS
Datasets for this paper originate from two rockfall sites along interstate 70 in Colorado. The first site is located in Idaho Springs, Colorado ) and the second site is located in Glenwood Canyon, Colorado. The Idaho Springs site is part of the Front Range of the Colorado Rocky Mountains and is characterized by jointed biotite gneiss. The Glenwood Canyon site is located in the White River Uplift and is characterized by jointed coarse-grained granite. The sites are prone to rockfalls due to their steep geometry and weathering during freeze-thaw cycles. The size of rockfalls is typically at the 1 m 3 level or smaller, although larger events (between 1 m 3 and 10 m 3 ) occur every few years Weidner et al. 2019).
The Idaho Springs site has points clouds from TLS, TP, sUAS photogrammetry. The site has a spatial extent of approximately 55 m width by 35 m height. Datasets are available for 2018 and 2019, with the exception of the sUAS datasets, where only a 2019 data acquisition was conducted. There are 13 ground GCPs on the rock surface that facilitated geo-referencing of the photogrammetric datasets and continuous monitoring . Their coordinates were established using total station observations. For the purposes of this paper the TP and sUAS datasets are used (Figure 3). The Idaho Springs site utilizes a system of five fixed station cameras established by Kromer et al. (2019). All five stations use a Canon 5DSR digital camera, which offer a high resolution sensor of 50 megapixels. The stations are about 10-15 m apart and positioned parallel to the slope. The stations are located across the slope at a distance of approximately 70 m. The system automatically captures images two times a day, at noon and at 12:30 pm for redundancy. Using Agisoft Metashape and python scripting, a point cloud is automatically created and registered to previously generated point clouds without needing GCPs on the rock surface ). The average point spacing of the resulting dense point cloud is 1.5 cm. sUAS imagery was collected using an Aibot X6, which carried a Sony Alpha 6000 digital camera with 24.3 megapixels. The sUAS also had accurate Global Navigation Satellite System The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2021 XXIV ISPRS Congress (2021 edition) (GNSS) -Real Time Kinematic (RTK), which allowed cmlevel positioning, which also reduces the requirements for the number of GCPs to about 4-6 (Bolkas 2019). A 30 mm focal length flight was conducted on June 4, 2019 with 35 images collected from a distance of 70-80 m from the rock surface. The images were processed in Agisoft Metashape and the average point spacing of the resulting dense point cloud is 2.2 cm.
The Glenwood Canyon site has TLS and TP datasets from October 13 of 2018 and January 19 of 2019. The site has a spatial extent of 150 m width and 50 m height (Figure 4). TLS and photogrammetry data were captured concurrently for both dates. A Faro Focus x330 phase shift scanner was used to collect data for this site. To ensure adequate coverage of the complex slope surface, TLS datasets were collected from three different stations at a distance of approximately 40 m from the slope and 50 m between stations. The three component datasets were then aligned and merged for each date using a semiautomated feature based and ICP registration algorithm (Schovanec, 2020). The final merged point clouds have a point spacing of approximately 1-2 cm and registration error between component datasets of less than 1 cm. Photographs at the Glenwood site were taken with a handheld 18-megapixel Canon EOS Rebel T2i, and at approximately equal spacing between photos along the same 150 m width as the lidar stations. The 2018 dataset had 15 photos, and the 2019 dataset had 14 photos. Similar to the Idaho Springs site, the point cloud was created using Agisoft Metashape. To scale the model and reduce artifacts on the large slope, eight distinctive points were selected from the TLS scan and used as "pseudo" GCPs in Metashape. The TP dataset has an average point spacing of about 5 cm. Of note is that in the second data acquisition in 2019, there is also snow covering parts of the rock surface. These are depicted with white color in Figure 4c for the photogrammetric dataset and with blue color in Figure  4d for the laser scanning dataset.

Idaho Springs site
For the first site, we are using the TP dataset from 2018 as the reference. For the TP registration, 13 GCPs were used. To show how the proposed algorithm can be utilized in multi-platform datasets, we are using the sUAS dataset acquired in 2019. For georeferencing, the sUAS dataset uses 0 GCPs, relying solely on GNSS. This scenario resembles a case where all GCPs were lost in the second data acquisition and one must rely on a registration algorithm to connect the two epochs. We assume that an accurate pre-calibration for the sUAS is available (as it is in this case); therefore, most of the error in the computed RMSE values will reflect misregistration. Figure 5 shows the matched edges for the TP and sUAS datasets in the registration using the developed algorithm. In general, the figure shows good matches. Some mismatches are found in the upper part of the site caused by vegetation; however, poor point correspondences have limited impact on the registration algorithm, as they are removed through thresholding in the Procrustes step. For this dataset, about 1,700 points were used for the final registration with Procrustes analysis.

Figure 5.
Matched edges in the Idaho Springs dataset for the TP and sUAS datasets Table 1 shows the RMSE statistics for these scenarios. Without a registration algorithm applied, the RMSE is 35.2 cm, which reflects direct georeferencing relying on GNSS. When the ICP The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2021 XXIV ISPRS Congress (2021 edition) is applied, the RMSE reduces to 3.3 cm, and when the developed algorithm is applied the RMSE reduces to 2.5 cm. Note that the RMSE using all 13 GCPs for the sUAS dataset is 2.1 cm; therefore, the 2.5 cm achieved shows that the developed algorithm adequately registered the two point clouds. Figure 6 depicts the change histograms for each case, showing the higher accuracy achieved by the developed algorithm (higher peak and smaller tails).   Rockfalls are detected at the 95% confidence level with respect to the RMSE values. Thus, in Figure 7a, changes higher than 3 cm are considered rockfalls, and in Figure 7b, changes higher than 5 cm are considered rockfalls. In general, we find that similar rockfalls are detected, showing how the developed algorithm can be used to support accurate rockfall detection using photogrammetric datasets when an accurate camera precalibration is available.

Glenwood Canyon site
For the second site, the TLS registration between the 2018 and 2019 epochs was already at the 1 cm level and focus is placed on the TP dataset registration. In this example, the TP dataset are registered with respect to the TLS datasets, thus comprising a multi-platform registration scenario. Alignment of the TP dataset was achieved using the TLS-derived GCPs, which might be insufficient. The registration using the TLS-derived GCPs shows an agreement at the 12.4 cm level for both the TP 2018 and 2019 datasets when using the TLS from the 2018 data acquisition ( Table 2). The ICP manages to reduce these values to 9.5 cm and 10.0 cm, respectively. When the developed algorithm is being implemented instead of the ICP, the RMSE values drop to 8.2 cm and 8.3 cm for the 2018 and 2019 TP datasets, respectively.  Figure 8 shows the matched edges that were used for this registration. Note that the same α-molecules input parameters with the Idaho Springs site were used, showing that after some initial training and experience minimal tuning is needed. For the The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2021 XXIV ISPRS Congress (2021 edition) 2018 comparison (Figure 8a) about 3,300 points were then identified from the matched edges to be used in the Procrustes analysis, and for the 2019 comparison (Figure 8b) about 3,000 points were identified. Figure 9 shows the histograms for the TP (2019) comparisons, depicting the smaller tails and higher peak in the histogram when the developed algorithm is used. These results further demonstrate the contribution of the developed algorithm to register multi-platform datasets and support rockfall monitoring from multiple point cloud data sources.   Figure 10 shows the detected changes at the 95% confidence level. Figure 10b shows the detected changes for the TLS (2018) versus TLS (2019) comparison, using 3 cm as the 95% confidence level, and Figure 10c shows the detected changes for the TLS (2018) versus TP (2019) comparison using the developed algorithm for registration and 17 cm as the 95% confidence level. Most of the detected changes in Figure 10b are due to snow coverage; therefore, few or no rockfalls are detected in the remaining parts of the rock surface (see also Figure 10a for the snow areas). Even though the photogrammetric datasets have a lower accuracy, they also confirm that no changes take place in the rock surface, and most changes are due to the presence of snow in the second epoch.

CONCLUSIONS
Accurate registration of point clouds is important for change detection and rockfall monitoring. Existing registration methods are often based on the ICP and its variants, while feature-based methods are often used for initial alignment only. This paper presented the main components of a new automatic algorithm (Bolkas et al. 2021) and demonstrated its application in a new and more complex rock slope dataset. The application of the algorithm in a new site validates the robustness of the developed approach. Application examples from two test sites in Colorado showed how the developed algorithm can be used to register multi-epoch and multi-platform point clouds and support rockfall monitoring applications.