REGISTRATION OF MMS LiDAR POINTS AND PANORAMIC IMAGE SEQUENCE USING RELATIVE ORIENTATION MODEL

We propose using the relative orientation model (ROM) of panoramic to register the MMS LiDAR points and panoramic image sequence, which has the wide applicability. The feature points, extracted and matched from panoramic image pairs, are used to solve the relative position and attitude parameters in the ROM, then, combining the absolute position and attitude parameters of the initial panoramic image, the MMS LiDAR points and panoramic image sequence are registered. First, we propose the position/attitude ROM (PA-ROM) and attitude ROM (A-ROM) of panoramic images respectively, which are apply to the position/attitude parameters both unknown and only the attitude parameters unknown. Second, we automatically extract and match feature points from panoramic image pairs using the SURF algorithm, as these mismatching points will affect the registration accuracy, the RANSAC algorithm and ROM were used to choose the best matching points automatically. Finally, we select the feature points manually from MMS LiDAR points and panoramic image sequence as the checkpoints, and then compare the registration accuracy of continuous/discontinuous panoramic image pairs. The results show that MMS LiDAR points and panoramic image sequence are registered accurately based on ROM (7.36 and 3.75 pixels is suitable for more road scenes.


INTRODUCTION
Mobile mapping systems (MMS) have been extensively used for spatial data acquisition, which can obtain optical image sequence (2D) and LiDAR points (3D) simultaneously. In order to capture optical images with large field of view, the panoramic camera is commonly used in MMS. The LiDAR points and optical image are different types of data, the registration of MMS LiDAR points and panoramic image sequence is important for many applications, e.g. texturing 3D models (Shan et al., 2019) (Yang and Dong, 2019), object extraction  (Yang et al., 2017), and point cloud classification  (Barnea and Filin, 2013).
There are many registration methods for 2D and 3D data, which can be divided into 2D-2D (Miled et al., 2016) (Plötz and Roth, 2017), 2D-3D (Liu and Stamos, 2012) (Shao et al., 2017) and 3D-3D (Corsini et al., 2013) (Zheng et al., 2013). The panoramic image is a special optical image (Zhu et al., 2018) (Zhu, 2019), the existing methods for registration of MMS LiDAR points and panoramic image sequence mostly rely on GPS/IMU (Original registration method) (Yao et al., 2017). The Original registration method directly obtains the position and orientation parameters from GPS/IMU, but the accuracy is affected by dropouts, etc., which lead to unreliable registration result (Miled et al., 2016). Another common method is Point-based registration method, it has the accurate precision, but automatically extract and match feature points from LiDAR points and panoramic image sequence is difficult, and the manual selection in image sequence is unrealistic. * Corresponding author. Bisheng Yang, E-mail addresses: bshyang@whu.edu.cn In addition, Wang et al. (Wang et al., 2012) proposed an automatic registration method for MMS LiDAR points and panoramas based on mutual information; however, only part of the panoramic image was used. Cui et al. (Cui et al., 2017) extracted line features from LiDAR points automatically and corresponding mono-images (1616×1232 pixels) semi-automatic, the registration accuracy was 4.244 pixels in panoramic camera; however, this method requires manual intervention. Li et al. (Li et al., 2018) extracted vehicles from panoramic images (2048×4096 pixels) via Faster-RCNN and took it as primitive pair, then particle swarm optimization was utilized to refine the translation parameters, the registration error less than 3 pixels; however, this method relies on the parked vehicles, the application is limited. Zhu et al. (Zhu et al., 2018) proposed using skyline to automatically register LiDAR points and panoramic image sequence, the brute force optimization was used to match skyline pixels and points. The mean registration error of 5 panoramic images (4000 × 8000 pixels) is approximately 10 pixels; however, the skyline pixels could represent very distant objects, while skyline points are missing in corresponding area, it may cause the mismatch of skyline pixels and points. Zhu et al.  proposed utilizing the feature points of road lamp and lane to register MMS LiDAR points and panoramic image sequence, which are easily extracted from LiDAR points and images, the mean registration accuracy is 5.84 pixels calculated by 31 panoramic images; however, this method will failure in the scene without road lamp and lane.
Above all, we propose utilizing the relative orientation model (ROM) of panoramic to register MMS LiDAR points and panoramic image sequence. The panoramic image sequences have high overlap, the feature points are automatically extracted and matched from panoramic image pairs, which used to solve the relative position and attitude parameters. Our registration method just tackle the image pairs (uninvolved LiDAR points), so it has larger applicability. First, the Position/Attitude ROM (PA-ROM) and Attitude ROM (A-ROM) are proposed respectively. Second, automatic extraction and matching of feature points from panoramic image pairs using SURF algorithm. Finally, we use the MMS LiDAR points and panoramic images in different road scenes for experiments. Figure 1 shows the flow chart of our registration method. This paper is organized as follows. Section 2 introduces the materials and methods, which includes two different road scenes dataset and ROM. In section 3, the comparative experiments are conducted to verify the effectiveness of our registration method. Finally, the conclusions are presented in section 4.

MATERIALS AND METHODS
We utilize the ROM of panoramic to register MMS LiDAR points and panoramic image sequence. First, the PA-ROM and A-ROM are proposed respectively. Second, automatic extraction and matching of feature points from panoramic image pairs by SURF algorithm. Finally, the elimination of mismatching feature points and registration accuracy evaluation were proposed. Fig.  2 and 3 show two different road scenes dataset, the panoramic cameras are Ladybug5 (Zhu et al., 2018). Figure 2 shows 10 panoramic image sequence (4000×8000 pixels), the interval between adjacent images is approximately 7 m. The length and height difference of the road are approximately 250 m and 4.6 m, which includes 3.8 million points. The MMS LiDAR points correspond to the panoramic image sequence at the same place, the initial values of the camera position (XS YS ZS) and attitude (rx ry rz) are known. Figure 3 shows 7 panoramic image sequence (4096×8192 pixels), the interval between adjacent images is about 20 m. The length and height difference of the road are approximately 350 m and 2.2 m, which includes 8.2 million points. These MMS LiDAR points and panoramic images are in the same place, the initial values of (XS YS ZS) is known, but (rx ry rz) is unknown. , plane map of the MMS LiDAR points (1600×1600, 320m×320m) and the panoramic image sequence, the red dots and number express 10 panoramic images (No.1-10).

The ROM of panoramic image
The ROM is often used to solve the exterior orientation elements in aerial photogrammetry, while the panoramic image adopts the spherical imaging model, with the field angle range of 360 °× 180 °. The ROM of panoramic image is based on the coplanarity condition, but the solution and mismatching feature points elimination are different. Firstly, we propose PA-ROM and A-ROM, there are suitable for the Position/Attitude initial parameters unknown and only position initial parameters unknown. Secondly, the feature points are extracted and matched from panoramic image pairs by SURF algorithm. Finally, we research the elimination method about mismatching points and the registration accuracy evaluation.

PA-ROM:
The panoramic camera is composed of multiple lenses, the images from different lenses are stitched together based on a fixed projection (Zhu, 2019). Eq. (1) shows the panoramic camera model using the spherical projection, the Y axis is perpendicular to the image plane (XOZ) (Zhu et al., 2018). (1) where, ( ) are the row and column coordinate of panoramic image, , is the rotation angle matrix, (x y z) are the coordinate of an object, D is used to convert the radians and image sizes (unit: pixels/rad).
The same name ray is expressed as: ( 2) where, , ( ) are the vertical and horizontal angles, , , row and col are the rows and columns of the panoramic image.
The baseline between two panoramic images is expressed as: ( 3) where, ( ) and ( ) are the imaging center of two panoramic images, , , .
Then, the coplanar equation constructed by the same name ray and baseline is shown in Eq.4.
( 4) where, the left and right epipolar line: ( ), ( ), and are the rotation angle matrix of two panoramic images.
Finally, remove the constant terms ( and ) that do not affect determinant calculation, Eq.4 is reduced to Eq.5. Eq.5 is the PA-ROM of panoramic image. In the model solving, separate the known quantities ( ), then Eq.5 is converted to: Then, the solution model composed of m points is expressed as: (7) where, , , .
The unknown matrix is calculated by Eq.8.
The relative position parameters and attitude parameter matrix ( ) can be calculated from by the following steps.
1) satisfy the Eq.9, then are calculated by .
2) is calculated by the characteristics of orthogonal matrix ( ), as shown in Eq.10.
Then, or (11) 3) Take and into Eq.6, so Multiply 2 sides by and , respectively, and then add them, thus, According to Eq.13, can be expressed as the function of , where, , Finally, is calculated by , as shown in Eq.15, then are calculated by Eq.12 and Eq.14.
The relative position parameters and attitude parameters matrix ( ) in the PA-ROM are obtained by above solution. It is assumed that the distance (s) between two panoramic images is known, which can be obtained directly by GPS, then the absolute position parameters and attitude parameter matrix ( ) of the next panoramic image are calculated by Eq.16.

A-ROM:
The influence of position parameters on registration is less than attitude parameters, the position parameters obtained by GPS can be used directly, so only the ( ) in rotation angle matrix (PA-ROM) are all unknown, considering that the interval between panoramic image sequence is small, the attitude changes little, we adopt the small angle matrix. There are only three attitude parameters ( ) are unknown in the A-ROM.
Finally, the solution model composed of t points is expressed as Eq.19, the unknown matrix can be calculated by Eq.8.
The relative attitude parameters matrix ( ) in the A-ROM is obtained by , and the absolute attitude matrix ( ) of the next panoramic image can be obtained by Eq.16.

The elimination of mismatching feature points and registration accuracy evaluation
The feature points, extracted and matched automatically from panoramic image pairs, have mismatching condition, it will affect the registration accuracy. Therefore, we use the RANSAC algorithm and ROM to eliminate mismatching points and optimize the registration result. Assuming that the number of matching feature points is m, and n are randomly selected to solve the model, there are combinations. We take the average error ( ) as the optimization criterion to select the best one group, the calculation of in PA-ROM and A-ROM are analyzed below.
The distance ( ) from ( ) in current panoramic image to the core line, where the point is located, is expressed as: The distance ( ) from ( ) in next panoramic image to the core line, where the point is located, is expressed as: (22) Thus, the kernel error ( ) of the feature point is shown in Eq.23.
Finally, the average error ( ), taken as the evaluation indexes, is calculated by each .
where, , , The distance ( ) from ( ) in current panoramic image to the core line, where the point is located, is expressed as: (26) The distance ( ) from the feature point ( ) in next panoramic image to the core line, where the point is located, is expressed as: Finally, the kernel error ( ) of the feature point and the average error ( ) are the same with Eq.23 and Eq.24.
In above section, we give the optimal index ( ) of PA-ROM and A-ROM, traverse n of m feature points, and select one group with the least , so as to ensure the solution accuracy. The panoramic image has large scale (4000×8000 pixels or 4096×8192 pixels in this paper) and rich texture, it is necessary to keep only the high reliability matching points for calculation efficiency.

EXPERIMENTS AND ANALYSIS
In order to verify the applicability of the ROM registration method (hereinafter referred to as N), datasets I and II are used in the experiments, the feature point registration method (M) is used to compare. First, we register 10 panoramic images in dataset I and 7 panoramic images in dataset II with M. Second, we use the SURF algorithm to extract and match feature points of panoramic image sequence. Finally, the accuracy of M, continuous registration of N (N1) and discontinuous registration of N (N2) are compared.
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B1-2020, 2020 XXIV ISPRS Congress (2020 edition)

The registration with control points
(XS YS ZS) and of the starting panoramic image need to know before using ROM to register panoramic image sequence, it has a great influence on the subsequent images. In addition, method M is taken as a reference for comparison of method N1 and N2.
We select the feature points manually from MMS LiDAR points and panoramic image sequence, which are also the checkpoints for accuracy evaluate. The feature points satisfy two conditions: 1) appear on most panoramic images; 2) the spatial distribution is uniform. As shown in Figure 4(a, b), 38 and 17 feature points are selected from points may be occluded in some panoramic images, so it is not all points involved in the registration. The panoramic image sequence are registered well by these feature points. As shown in Table 1. and 2, n express the number of control points in the registration, express the root mean square error (RMSE) of checkpoints in each image. According to the registration error distribution, the top/bottom of panoramic image, and the objects with close distance are prone to large errors, we remove the feature points at these areas. There are at least 34 points (No.6) used in dataset I, and 12 points (No.7) used in equals 4.69 and 2.83 pixels respectively, which meets the accuracy requirements of the starting image.

Extraction and matching of feature points from panoramic image sequence
We use the SURF algorithm to extract and match feature points of panoramic image pairs, which consider the accuracy and efficiency. In the SURF algorithm, the extraction parameter (Q) and matching parameter (T) determine the number of extracted and matched feature points. In panoramic images, the left and right sides are objects on both sides of road, the top and bottom are occupied by the sky and road, respectively. The rich texture are mainly located in the middle area of panoramic image, we divide the panoramic image into left and right parts, then extract and match the feature points respectively.
In the matching, we exchange the matched and to be matched image, only retain the feature points matched successfully in both directions. Then, the trajectory is approximately a straight line, the object in image sequence gradually moves to both sides, this direction filtering can eliminate some mismatching points. In addition, the number and distribution of feature points will affect the registration efficiency and accuracy, we use the coordinate difference of matching points in different images (external filtering) and the same image (internal filtering) to optimize the matching points respectively. As shown in Figure 5(a, b), white lines are the connection of feature points before filtering, black dots and black lines are the feature points and their connection after filtering, the blue, red and green lines are the feature points connection in different filtering processes. The trajectory locate on the right side of road, the objects on right side are more inclined to distribute at the bottom of panoramic image, so the left and right sides images start from line 1000 of panoramic image, and the left side image are 1250×4000 pixels, the right side image are 1450×4000 pixels.
According to above method, set Q=2000, as shown in the upper left of Figure 5( Figure 5(b) shows the feature points matching on the right side of panoramic image, the parameters are the same as the left image in the extraction, matching and filtering process. Finally, the matching points from the left and right images make up the matching points of panoramic image pair No.1-No.2. No.9-10 No.1-5 No.1-10 Figure 6. Feature points matching of 3 pairs in dataset I Then, we extract and match the feature points of 10 panoramic images in dataset I. In addition, the discontinuous panoramic image pairs also can be matched, as the similarity is weakened with the increase of image pair spacing, the parameters in discontinuous pairs are relaxed to ensure the number of matching points, here, set Q=1000 and T=0.5. If the matching points after filtering are insufficient, the feature points are directly used. Limited to space, Figure 6 only shows the matching results of 3 image pairs, No.9-10 is a continuous pair, the similarity of these 2 images is well, there are 14 and 5 points are successfully matched from the left and right sides. No.1-5 and No.1-10 are discontinuous pairs, there are 7 and 1 points are matched from the left and right sides in No.1-5, but 2 and 0 points in No.1-10, which doesn't meet the number requirement, so the feature points after coincidence matching are directly used in No.1-10.

Dataset
: The feature point extraction and matching of panoramic image pairs in dataset II is the same as dataset I. The size of these images is 4096×8192, the trajectory is close to the right side of road, so the images on the left and right sides both start from the line 1500, the left image is 800×4096 pixels, the right image is 1000×4096 pixels. In the experiment, set Q=2000, T=0.4, the external filtering parameter is 10 pixels, and the internal filtering parameter is 20 pixels. Figure 7 shows The feature points of 7 panoramic images in datas extracted and matched. Similar with dataset I, the discontinuous pairs also can be used, here, set Q=1000 and T=0.5. If the matching points after filtering are insufficient, the feature points are used directly. Figure 8 shows the matching results of three image pairs, No.6-7 is a continuous pair, which have high similarity and easy to match. No.1-7 has the largest spacing and worst similarity, the road lane and moving platform points are mismatching points in the left and right sides, even if the er is guarantee, the parameters in ROM cannot be solved accurately.

The registration of panoramic image sequence using ROM
The A-ROM is used to register the MMS LiDAR points and panoramic image sequences in datasets I and II. First, calculate the relative attitude parameters by the matching feature points; Second, the continuous (N1) and discontinuous (N2) registration methods based on the A-ROM are experimented, they are compared with the feature point registration method (M).

3.3.1
: We select 6 points from all matching points to calculate the relative attitude parameters of panoramic image pairs. By traversing all combinations, the minimum average parallax of both images is taken as the optimal result, and then the registration of next pair is carried out. Table 3 shows nine continuous pairs in dataset I, q, t and c are the number of matched feature points, combinations and iterative optimization respectively, v the registration error (N1).  Table 3. The registration of continuous pairs in dataset I As shown in Table 3, the combinations ( ) increases sharply with the number of feature points (q), so we only keep the necessary matching points. We optimize parameters by minimizing the average parallax (v), and then traverse all combinations, the final registration results are obtained after 34 times (c) optimization at most. The average parallax (v) of nine panoramic image pairs in dataset I is 1.19~6.21 pixels, and the registration error ( ) is 4.93~10.67 pixels. These nine pairs are optimized different times, the least is No.7-8, which has 9 matching points and 84 combinations, the final result is obtained after 2 times of optimization. The most of which are No.2-3 and No.9-10, these pairs both have 19 matching points and 27132 combinations, the results are obtained after 13 and 31 times of optimization respectively. Figure 9 shows the registration results of No.2-No.10 panoramic images using the A-ROM. In above experiments, we register the continuous image pairs (adjacent images), which is called the continuous registration method based on the ROM (N1). In addition, the discontinuous image pairs composed of the image to be registered and the starting image also can be registered, which is called the discontinuous registration method based on the ROM (N2). The extraction and matching of feature points in discontinuous image pairs (No.1-3, No.1--10) have been analyzed above, other process about optimization and accuracy evaluating are the same as N1, here, we only give the final registration results. As shown in Table 4, No.1 is the starting image both in method N1 and N2, which registered by method M and the accuracy is 4.69 pixels. We compare the accuracy of these three registration methods, the registration error of No.10 is 62.20 pixels with method N2, the reason is that the distance between No.1 and No.10 is large, there are no enough correct matching points. Therefore, excluding No.1 and No.10, the average registration errors of No.2-No.9 panoramic images are 5.98, 7.36 and 7.84 pixels respectively. This experiment shows method N1 is slightly better than method N2. In addition, the registration error of method N1 is accumulated, except for No.6, which increases gradually from No.1 to No.10. It shows that method N1 has limitations, the registration error will accumulate to unacceptable level with the increase of images. Method N2 may be failure when the image pairs spacing is large. For these situations, we can reselect and recalculate the starting image with method M, then, method N1 and N2 are used to register the subsequent images.
3.3.2 Dataset : Using the same processing as dataset I, we  Table 5. The comparison of M, N1 and N2 in dataset (pixels) As shown in Table 5, No.1 is the starting image, the registration accuracy is 2.83 pixels using method M. The number of matching points in No.1-6 pair less than 6, so this registration is failure with method N2. The registration error of No.7 is 55.69 pixel with method N2, which is much large than method M (3.44 pixel) and N1 (4.33 pixel), the reason is that No.1-7 pair spacing is too large and no enough correct matching points. Therefore, excluding No.1,No.6 and No.7, the average registration errors of No.2-No.5 panoramic images are 3.34, 3.75 and 4.65 pixels, respectively. This experiment shows method N1 is close to method M, both better than method N2. In addition, method N1 is better than method M in No.3,No.4 and No.6.

CONCLUSIONS
In this paper, we propose the registration method of MMS LiDAR points and panoramic image sequence based on ROM, which includes the PA-ROM and A-ROM of panoramic image, extract and match feature points, and the registration accuracy evaluation. In the experiments, panoramic image sequences of different datasets are used to verify the effectiveness of our registration method. The continuous and discontinuous registration methods based on ROM are compared, the results show that the continuous is superior to the discontinuous registration, but there are both lower than the control points registration method. Our registration method transforms the registration of panoramic image and LiDAR points into the matching of feature points between images, which avoids the tedious work of extracting feature points from images and LiDAR points, so it is suitable for different road scenes. In addition, the model used in continuous and discontinuous registration is A-ROM, which adopts the small angle attitude matrix, it is suitable for the small attitude change between panoramic image pairs, for the large attitude change, such as the road bending and the distance pairs, the PA-ROM can be used to solve.