IDENTIFYING CORRESPONDING SEGMENTS FROM REPEATED SCAN DATA

It has been demonstrated that surface changes in the order of millimeters are detectable using terrestrial laser scanning. In practice however, it is still virtually impossible to detect such small changes from for example repeated scans of a complex industrial scene. The three main obstructions are, ﬁrst, a priori uncertainty on what objects are actually changing, second, errors introduced by registration, and third, difﬁculties in the identiﬁcation of identical object parts. In this paper we introduce a method enabling efﬁcient identiﬁcation that can also be applied to evaluate the quality of a registration. The method starts with a pair of co-registered point clouds, at least partially representing the same scene. First, both point clouds are segmented, according to a suited homogeneity criterion. Based on basically orientation and location, corresponding segment parts are identiﬁed, while lack of correspondence leads to the identiﬁcation of either occlusions or large changes. For the corresponding segments, subtle changes at the millimeter level could be analyzed in a next step. An initial version of the new method is demonstrated on repeated scan data of a metro station experiencing heavy construction works.


INTRODUCTION
Detecting changes in scenes sampled by laser scanning data is a topic with a rapidly increasing number of applications, (Vosselman and Maas, 2010), chapter 7. Administrations, ranging from individual communities to complete countries are initiating airborne laser ranging archives.In The Netherlands for example, the second version of a nation wide archive is almost completed, while plans have been initiated for a third version, (AHN, 2000).Meanwhile, car based laser mobile mapping systems, (Toth, 2009), are now available in many countries and several researches are initiated to investigate their usage in both built-up, (Haala et al., 2008) and natural environments, (Barber and Mills, 2007).As a consequence, methodology aiming at the detection of changes from such data sources is also under active research, (Champion et al., 2009).
The most easy way to obtain laser point clouds from multiple epochs is by repeatively scanning from a fixed position using a terrestrial laser scanner, e.g.(Little, 2006, Lindenbergh andPfeifer, 2005).It has been demonstrated, (Gordon and Lichti, 2007), that in a well controlled, close range laser experiment it is possible to measure vertical displacements with an accuracy of ± 0.29 mm.To obtain such accuracies from an analysis of repeated scan data of an arbitrary complex scene is still very challenging however.One major issue is that in general it is on forehand unknown what and to what extent objects in the scene were changed.
A computationally efficient method to identify changes is by means of octrees, (Barber et al., 2008, Girardeau-Montaut et al., 2005), that divide the 3D space in small cubical cells, each containing part of the point cloud.Differences in the octree division for different point clouds correspond to changes in the point clouds.Another efficient way is to reduce the 3D problem to 2D by comparing the relative visibility of different point clouds from one, possibly virtual, viewpoint, (Zeibak and Filin, 2007).By com-bining the use of octrees with Triangulated Irregular Networks, it becomes possible to analyze changes for individual scan points, (Kang and Lu, 2010).
A disadvantage of these methods is that they are not object oriented.An alternative is to first identify objects and in a next step consider if these objects are changed.One way to identify objects in a point cloud is by means of a segmentation method, (Matikainen et al., 2009).In this paper we present a new approach that systematically evaluates segmentations of point clouds of different epochs.The method aims at identifying corresponding parts of segments, sampling the same object.Segments that have no correspondence in other epochs can be linked to occlusions or large changes well above the signal to noise ratio.What remains are segments corresponding to the same object surface.In a final step it can be decided using statistical testing methods, whether locally this object surface has changed in a significant way, (Levin andFilin, 2010, Van Gosliga et al., 2006).

IDENTIFICATION ALGORITHM
Below we discuss all ingredients of the proposed method in detail.We focus on the identification of parts of planar segments that represent the same object.

Algorithm overview
For convenience we consider here the case of just two point clouds, PC1 and PC2, of approximately the same scene.The method is easily extended in case more point clouds are given.For the moment we also only consider a segmentation into planar segments.Algorithm 1 summarizes the major steps.

Registration and Segmentation
First both point clouds should be registered in a common coordinate system, using any available method.Then both point clouds end if 25: end for 26: return INT, list of segment intersections are segmented into planar patches, for example using the methods described by (Rabbani et al., 2006).It should be noted that often the segmentation of two different point clouds representing the same scene will lead to different results, even when the same parameter settings are used: Distinct scanners that operated from distinct scan positions will result in point clouds with locally different point densities and noise levels.Both point density and noise level will affect the segmentation.

Segment Plane location
In the following, a first match of segments from both point clouds will be performed, based on two segment properties.These two properties together parameterize the position of the plane containing the local planar segment relative to the chosen coordinate system.The first property is the segment orientation, given by the unit normal vector n = (a, b, c) of the segment.The second property is the distance to the origin of the plane containing the segment.Let c S = (cx, cy, cz) denotes the center of gravity of the points in the segment.The best fitting plane in the least squares sense of the points in the segment passes through c S and the distance dS of the plane containing S, to the origin is simply given as where T denotes transposition.To increase numerical stability it is convenient to take the center of gravity of, say, the first point cloud as the origin of the local coordinate system.

Segment Plane matching
The combination of normal vector and distance to the origin are used to determine if two segments, each of one point cloud, are in a common plane, compare Algorithm 1, line 12.A more formal formulation is that the following two conditions should hold for two segments S 1 i and S 2 j from the first and second point cloud respectively: where d(S 1 i ) and d(S 2 j ) denote the distance to the origin from the first and second segment plane, while n(S 1 i ) and n(S 2 j ) denote normal vectors of the first and second segments respectively.Here, ∆d is a preset distance threshold and ∆αn a preset angular threshold, that should be specified by the user, depending on e.g. the quality of the registration and the local noise level of the scan points.The angle (u, v) between two vectors u and v is computed as where in addition the two possible antipodal orientations of a normal should be taken into account.

Reduction to 2D
If for a pair of segments (S 1 i , S 2 j ) both Equation 2 and 3 hold, both segments are approximately in the same plane but are not necessarily intersecting.They could for example both represent a different part of one big wall.Therefore the area of intersection Sij should be determined.
As by now we know that S 1 i and S 2 j are approximately in the same plane, we can reduce the problem at first approximation to a 2D problem as follows.Let c i be the center of gravity of the points in segment S 1 i .We consider the projections of the segments S 1 i and S 2 j on the plane through c i perpendicular to the normal vector of the segment S1.In practice such projection can be obtained as follows, compare (Lay, 2003).From the centralized segment with C 1 i the matrix with each row equal to the center of gravity of the points in segment S 1 i , the eigen vector matrix E 1 i of the variance-covariance matrix of S 1 i is determined: The first eigen vector points in the direction of maximal variation in segment S 1 i , second eigen vector in the remaining orthogonal direction of maximal variation, while the last points in the remaining orthogonal direction that has least variation, see also Figure 1.As a consequence, the first two eigen vectors span (in most cases) the least squares plane through the (centralized) points of segment S 1 i .Therefore the desired projection is obtained by expressing the (centralized) segment points of both segments with respect to this eigenvector basis for the first, reference segemnt: In this new coordinate system, the first two coordinates express the location of a point in the segment plane of S 1 i , while the third coordinate gives the distance of the point to the plane.
In Figure 2 the results of such projection are illustrated.

2D Segment intersection
What remains, see also Figure 2, is to find the intersection of two planar segments s 1 i and s 2 j .To increase efficiency, this is done in two steps.First, for both segments s 1 i and s 2 j , the bounding Cartesian coordinates are determined, where e.g.x m i denotes the minimal x-coordinate of the points in segment s 1 i and y M j the maximal y-coordinate of the points in segment s 2 j .If the bounding rectangular boxes for both segments do not intersect, the segments themselves certainly do not intersect.If the bounding boxes do intersect, the second step only has to consider the points of both segments in the intersection.
In the second step the points from each segment (in the intersection) are organized in a raster, with a preset width, of, say, 5 cm.If a raster cell contains a segment point, it is given the value 1, otherwise it gets the value 0. The intersection of both segments is now determined from overlaying the two rasters.Only if two corresponding raster cells both have the value 1, it belongs to the intersection of the two segments.This method is flexible in the sense that it also deals with segments with holes, corresponding e.g. to a window in a door.This method can be further refined by using quad-trees instead of simple rasters: the quad-tree structure can be used to more precisely define the boundary of the segments.Each raster cell containing segment points, but 4-adjacent to an empty raster cell, is subdivided into four smaller cells, until a minimal cell size is reached.The intersection of two such quad-trees, each corresponding to a segment, basically works the same as intersecting two rasters.

Implementation
The described method is relatively easy to implement.Using least squares adjustment, (Teunissen, 2000), the normal of each segment can be derived.The same least squares adjustment can be used to project candidate corresponding segments to the least squares plane of the first segment.What remains is constructing a raster or quad-tree for a relative small number of possibly parts

DISCUSSION
In this section, the potential of the method described in Section 2 is illustrated for a case It should be noted however that the method used to obtain the results below was an initial method that differences at details from the method described in Section 2.

Data description
For the case study we consider two point clouds representing a part of the Rotterdam Central metro station.The selected area an area of 8 × 3 m of the wall on the southern side of the tunnel, see Figure 3.During the first epoch, acquired at April 24, 2006, approximately 350 thousand points of this specific part of the tunnel were sampled from one standpoint using the Z+F Imager 5003 phase based scanner.During the third epoch, d.d. November 20, 2007, 3.5 million points were obtained, from two different standpoints, compare also Figure 3, left, using the Faro LS880 phase based scanner.From one of these standpoints, two scans were acquired, resulting in three available scans in this epoch.
Registration of individual scans was performed using the software Cyclone from Leica.Both control points and object matching were used.In Epoch 1, the reported accuracy was 3 mm, in Epoch 3, it was 6 mm.Control points were additionally measured with a total station, and both point clouds were georeferenced to the same absolute coordinate system RDNAP.The two resulting point clouds are shown in Figure 4.The effect of grafitti on the inmeasurements is clearly visible in both epochs.In Epoch 3, on the right, the footbridge is visible, indicated by a red ellipse, which was not yet present in Epoch 1.

Segmentation
Point clouds from both scans were segmented into planar patches using the region growing method described in (Rabbani et al., 2006).The method requires three parameters, a number, k, of neighbors for each point, an angular threshold, αp, to compare local surface normals, and a distance threshold, dp, that considers the distance between the point at hand and a growing segment.Here, the parameter values k = 30, αp = 30 • and dp = .03mwere used.These setting resulted in 20 segments for Epoch 1, and 494 segments for Epoch 3. The increase in the number of segments from Epoch 1 to Epoch 3 already indicates that more different objects are in the scene in Epoch 3.

Corresponding segments
For deciding if two segments from two different epochs intersect, some additional parameter values have to be set, as described in detail in Section 2.3.First, it is decided if two segments are in the same plane.For this purpose, the maximal difference, αS, in orientation between two segments is set at 10 • .Moreover, the distance, dS, between segments, as defined in Section 2.3, should not be more than 0.10 m.To determine if two segments in the same plane intersect, no additional parameters were needed.
The resulting corresponding segments for the case study are indicated by color correspondences in Figure 5.The dark blue segment pair, representing the main wall, contains 57 % of the points in Epoch 1, versus 45 % of the points in Epoch 3. The cyan colored pair of segments, representing the lower part of the wall contains 18 % and 13 % of the points in Epoch 1 and Epoch 3, respectively.The only remaining significant segment pair is the gray one, representing a small band of wall between the blue and the cyan segments, containing 5.4 % and 2.4 % of the data in 1 and Epoch 3, resp.
2: Top: Two intersecting segments, with blue and red points resp., in a common segment plane; Middle: Same, but zoomed in on the red points; Right: Two non-intersecting segments in the same segment plane.Axes are labeled in centimeter.

Possible large changes
Points that have no correspondence are shown in Figure 6.Points that are marked by a capital 'A' are located at the edges of the data sets and have no correspondence because the two data sets were cropped in a slightly different way; Points marked by 'B' are situated behind the location where later the footbridge would be placed in Epoch 1, and on the footbridge in Epoch 3. Similarly, points marked 'C' correspond to points behind and on water pipes that were added to the scene between epochs.

Possible changes in intersecting segments
Figure 7 illustrates the distances between points from corresponding segment parts.The method proposed in this paper allows to efficiently identify such correspondences.Still very interesting differences may occur between points from corresponding segments.A first type of differences may be caused by registration errors.Errors caused by registration are expected to be locally stable, that is, for a point cloud covering several tenths of meters, the registration error between two point clouds is expected to be more or less the same over a stretch of one meter.
To mitigate the effect of registration errors, the mean difference between points belonging to corresponding segments can be removed, which is how Figure 7 is obtained.What remains are differences that can be caused by either real local differences in the surface elevation or by local scanner artifacts.In Figure a circular pattern is visible corresponding to some systematic differences in the order of a few millimeter.
There are two possible explanations for this particular pattern.The first is hardware related.In Figure 3 it can be seen that in both cases the scanner was located directly in front of the region of interest.This may have resulted in some saturation effect that could explain the circular pattern.Another possible explanation is physical.The signal may be the result of ground water outside the tunnel.The footbridge is placed at this location to protect a water inlet.Below the footbridge, a water tap is constructed to inject water in the ground around the tunnel.The water tap is connected to the water pipes that run through the tunnel.The ground around the tunnel presses on the exterior of the tunnel.When the water is injected in the ground, the pressure on the tunnel increases.The increase is the largest around the injection point, below the footbridge.Further away the pressure decreases, which may explain the circular pattern.

CONCLUSIONS AND FUTURE WORK
In this paper a new method has been presented that in an efficient way identifies corresponding parts of different point clouds representing the same planar feature.This method allows to decompose repeated point clouds in components corresponding to large changes or occlusions on one hand and components sampling the same object surface on the other hand.The method is computationally efficient, given a segmentation into planar segments.The method can be extended to segmentations including other geometric shapes, like cylinders or spheres.As a consequence it would be suitable for analyzing repeated point clouds of most human made objects.In combination with the methods of e.g.(Zeibak and Filin, 2007) a complete work flow can be obtained considering both possible changes close to the signal to noise ratio, occlusions and large changes.

Figure 1 :
Figure 1: Basis consiting of eigenvectors of the variance covariance matrix of the points in a segment.

Figure 3 :Figure 4 :
Figure 3: Left: Schematic representation of standpoints in both epochs and area of interest.Right: Photograph of the footbridge taken in January 2010.