JOINT SIMULTANEOUS RECONSTRUCTION OF REGULARIZED BUILDING SUPERSTRUCTURES FROM LOW-DENSITY LIDAR DATA USING ICP

There are many applications for 3D city models, e.g., in visualizations, analysis, and simulations; each one requiring a certain level of detail to be effective. The overall trend goes towards including various kinds of anthropogenic and natural objects therein with ever increasing geometric and semantic details. A few years back, the featured 3D building models had only coarse roof geometry. But nowadays, they are expected to include detailed roof superstructures like dormers and chimneys. Several methods have been proposed for the automatic reconstruction of 3D building models from airborne based point clouds. However, they are usually unable to reliably recognize and reconstruct small roof superstructures as these objects are often represented by only few point measurements, especially in low-density point clouds. In this paper, we propose a recognition and reconstruction approach that overcomes this problem by identifying and simultaneously reconstructing regularized superstructures of similar shape. For this purpose, candidate areas for superstructures are detected by taking into account virtual sub-surface points that are assumed to lie on the main roof faces below the measured points. The areas with similar superstructures are detected, extracted, grouped together, and registered to one another with the Iterative Closest Point (ICP) algorithm. As an outcome, the joint point density of each detected group is increased, which helps to recognize the shape of the superstructure more reliably and in more detail. Finally, all instances of each group of superstructures are modeled at once and transformed back to their original position. Because superstructures are reconstructed in groups, symmetries, alignments, and regularities can be enforced in a straight-forward way. The validity of the approach is presented on a number of example buildings from the Vaihingen test data set.


INTRODUCTION
For several years, 3D city models assume a central role in urban and regional planning, surveying, navigation, and telecommunications and also allow in the environmental field precise analyses and simulations of pollutant, flood and noise propagation, and solar potential.The application areas are therefore quite diverse and lead at the present time to an increased demand for such models.Depending on the application, different requirements (e.g.accuracy, level of detail, etc.) for 3D city models have to be met.For this purpose, several methods for the reconstruction of buildings have been introduced; an overview is, e.g., given in (Haala and Kada, 2010).The automation rate has nowadays increased so that also large urban areas can be automatically reconstructed from airborne point cloud data.Depending on the data quality (e.g.density, accuracy, completeness, etc.), modern reconstruction approaches are generally suited to generate realistic looking 3D building models with detailed roof and façade structures.
However, problems with the existing automatic reconstruction approaches often occur if the given input data does not meet the desired or required quality criteria.In this case, only course and generalized models without details can be created.Measured points from roof superstructures (e.g.chimney, dormers, etc.) are thereby usually considered as noise and are therefore ignored during the reconstruction process.But for a detailed 3D city model and many analysis, complex superstructures on building roofs need to be reconstructed as well.Several approaches have been developed therefore that offer users the possibility to add them in a semi-automatic way or that utilize additional data sources like images (Rottensteiner and Briese, 2003;Habib et al., 2010).But most of them are not applicable in practice for a largescale reconstruction process due, for example, to missing additional data sources.
To solve this issue, we present in this paper a new approach for the fully-automatic reconstruction of regularized roof superstructures from low-density LIDAR data, which might feature partially occluded areas.The basic idea of the approach is to take advantage of the fact that often more than one instance of the same superstructure occurs on the same roof as e.g.mentioned in (Oude Elberink, 2008).Instead of reconstructing each superstructure independently from one another, our proposed approach detects first all instances of a roof superstructure and reconstructs them afterwards at once.For this purpose, a point cloud registration technique based on the Iterative Closest Point (ICP) algorithm (see section 3) is presented.It considers basic building knowledge (e.g.local directions, symmetries, repetitive structures, etc.) and increases locally the point density so that a more accurate model can be generated.For evaluation purposes, the proposed approach is tested on several buildings of residential areas in low-density point clouds (1.5-6 points/m²).It is shown that our method helps to find and reconstruct small rooftop details that would otherwise be missed.

RELATED WORK
Up to now, not much research has so far explicitly dealt with the automatic reconstruction of small roof objects.One reason might be the fast development of sensors that provide nowadays denser and more accurate point clouds than before.Consequently, the processing of low-density point clouds becomes in many cases obsolete.However, if for example an urban area from a certain time has to be reconstructed, these new sensors cannot be utilized.Therefore, there will be always a demand for automatic reconstruction approaches that are able to handle and to reconstruct as many details as possible from low-density point clouds.
An early reconstruction approach that tackles this problem for flat superstructures on top of flat building roofs is described in (Stilla and Jurkiewicz, 1999).It creates a histogram and makes use of the peaks to segment the height data and to recognize areas of possible superstructures.Minor peaks with a certain gap to the main and to other minor peaks are then examined regarding their extent and compactness.Afterwards, points of an accepted minor peak are approximated by a prismatic object.This approach is in practice only applicable to a small number of building roofs in an urban area and strongly depends on predefined parameters.
Therefore, most research effort has been put to the following strategy: first, an initial model that consists of the main building components is reconstructed and then the best fitted parametric shapes from a predefined library are added during the (modeldriven) reconstruction of small roof objects.An example is given in (Brédif et al., 2007), in which the problems of reconstructing small objects on rooftops are discussed and which can be extended to point clouds.It presents a parametric roof superstructure reconstruction that makes use of a Minimum Description Length (MDL) energy minimization technique.Further examples that define different dormer shapes as superstructures, which are then placed on top of a base model, are presented in (Milde and Brenner, 2009) and (Vosselman and Dijkman, 2001).These model-driven approaches are generally better suited than data-driven methods due to the low number of input points.But on the one hand, they fail if the point density is too low and on the other hand the reconstruction capability is always limited to a predefined library.To weaken the latter disadvantage, the number of shapes in the predefined library can be extended causing an increased computational time.
To reduce the impact of the exhaustive search on the computational time, some approaches integrate a coarse superstructure detection step to detect the shape type before a shape refinement step is carried out.For example, in (Satari, 2012) a support vector machine (SVM) is used for the recognition of three predefined dormer types.It utilizes the gradient and azimuth values of the normal vectors of the superstructure and the underlying roof plane as a discriminating factor.Afterwards, the initial superstructure model is refined and added.Tested on a data set with an estimated density of 9 points/m², the overall qualitative and quantitative assessments confirmed to a certain extent the potential of this strategy.Another example for a coarse superstructure detection step is given in (Dornaika and Brédif, 2008).
For point clouds with a very high density, a generative statistical top-down approach for the reconstruction of superstructures is presented in (Huang et al., 2011).It searches for superstructures in the area above a base roof by using two simple parametric primitives (flat and gable roofs).The selection and the parameter estimation are driven by the Reversible Jump Markov Chain Monte Carlo technique.Tested on a data set with an average point density of about 25 points/m², roof superstructures have been reconstructed with an average z-error of 0.05 m.The success of the method is, however, highly dependent on the point density.
Summarizing the previous research, it deals mainly with modeldriven approaches because they are usually more stable for sparse point clouds than data-driven approaches.To overcome the computational time caused by an exhaustive search technique, a coarse detection of the superstructure type can be performed.However, all approaches reconstruct superstructures independently, wherefore a more or less high point density is required to get an accurate model.

REGISTRATION OF TWO POINT SETS BASED ON ITERATIVE CLOSEST POINT
During the last decades, several registration approaches have been introduced to find a transformation that best represents the relative displacement between two surfaces with an overlapping area; examples are given in (Silva et al., 2005).In general, they can be roughly divided into coarse and fine registration techniques.In coarse registration methods no prior knowledge of the relative spatial position is required.For this purpose, they extract certain features from one point set and assign them to the extracted features of the other point set.Depending on the number of features and the assignment certainty, they usually provide only rough alignments.For a more precise alignment, prior knowledge regarding an initial transformation is used in fine registration approaches.One of the most common method is summarized in the following paragraphs.
Iterative Closest Point is an algorithm employed to register two (partially) overlapping point sets in a common coordinate system (Besl and McKay, 1992).For this purpose, it iteratively minimizes the mean squared error (MSE) between the point sets by applying rigid transformations.Compared to other fine registration techniques, an exact point-to-point correspondence between the point sets is not required.Nowadays, there exist many variants of ICP.An overview is given in (Rusinkiewicz and Levoy, 2001), (Grün and Akca, 2005), and (Pomerleau et al., 2015).
Essentially, for two given point sets  and , which are already roughly registered, each iteration in ICP is composed of the following steps: 1. Find for each point in  the closest point in .
2. Estimate a transformation that minimizes the MSE of the point pairs from the previous step.3. Transform all points in  by applying the obtained transformation.4. Start the next iteration with the first step again if a predefined stopping criteria is not fulfilled.
In the first step, a correspondence search between points of  to  is performed by a nearest-neighbor search.To improve the speed of ICP, several approaches have been developed to increase the performance of computing corresponding points, because this is the most time consuming part of ICP.Usually k-d trees (Bentley, 1975) are used to accelerate the search.Several variants of the k-d tree have been developed to improve the performance like a cached (Nüchter et al., 2007) or a GPGPU (Wu et al., 2015) version.
For the second step of each iteration, a rotation and translation is computed that minimize the MSE between both point sets.The objective function that has to be minimized for two given point sets  = {  } and  = {  } with  = {1, … , } and  = {1, … , } is therefore: where   is the closest point in  to the point   ∈ ,  the rotation matrix and  the translation vector.For the minimization of the objective function and the calculation of the rigid transformation, there exist already several methods; for example, singular value decomposition (SVD) is used in (Arun et al., 1987), orthonormal matrices in (Horn et al., 1988) and dual quaternions in (Walker et al., 1991).
Afterwards, the resulting rotation and translation can be directly applied in the third step to all points in .Then, the algorithm starts again from the beginning if a lower MSE value than in the iteration before is received.Otherwise the algorithm terminates and a sequence of transformations is given that can be merged to a single rigid transformation, which minimizes the MSE with the final point correspondences at once.
Due to the greedy nature of the ICP algorithm, it converges always to a local minimum but not necessarily to the global minimum.Therefore, to get useful results, it requires a good coarse a priori alignment of the given point clouds.Furthermore, the method as presented above is limited to rigid transformations (rotations and translations), so that it is for example not capable of handling scaling.

RECONSTRUCTION PROCESS
In this section the workflow of our fully automatic reconstruction approach for regularized roof superstructures is outlined.To perform it simultaneously with the reconstruction process of the building, it does not rely on an initial building model but it relies on a common segmentation result.Our reconstruction approach consists roughly of the following three steps: (1) detection of appropriate roof superstructure candidate points; (2) grouping of similar roof superstructures; (3) roof superstructure modeling and construction.
To determine appropriate candidate points from the input data, our approach makes use of virtual points which are located below measured surface points.In conjunction with those surface points that belong to already detected segments, candidate areas of superstructures and their surface points are identified (see section 5).In the next step, a registration technique based on the ICP algorithm is carried out to transform candidate points, which hint to the same superstructure type, into a common coordinate system so that their areas overlap each other.This results in a higher point density for all superstructure types that possess more than one instance in the real building (see section 6).Afterwards, only one superstructure model is created for each superstructure type using feature-driven modeling techniques.Instances of it are then placed at the candidate areas.Finally, the whole building model is checked for overlapping instances and, if necessary, corrected (see section 7).The basic idea of this workflow is visualized in Figure 1.To proof the applicability, the proposed approach is integrated in the fully automatic building reconstruction approach presented in (Wichmann and Kada, 2014).It consists of the following four main steps: (1) segmentation of planar regions; (2) rule-based recognition of building features; (3) adjustment of half-spaces to emphasize regularities that are often present in real buildings; (4) construction of the final building model utilizing Boolean operations.An overview of the whole building reconstruction process is shown in Figure 2.
Figure 2. Overview of the proposed reconstruction approach of regularized superstructures in low-density point clouds that has been integrated in the building reconstruction approach presented in (Wichmann and Kada, 2014).

DETECTION OF APPROPRIATE SUPERSTRUCTURE CANDIDATE POINTS
A crucial point for the automatic reconstruction of superstructures is the detection of those points that belong to a superstructure.This task is challenging because it has to deal with a small number of points so that their surface shape is often not clearly recognizable during the reconstruction process of the base roof.Therefore, many problems might occur if the given point cloud is sparse and noisy.To improve the visualizing capabilities, the presented examples in this section are shown in a higher point density than in reality.
For the detection of appropriate candidate points, our proposed approach utilizes virtual points as they are defined in (Kada and Wichmann, 2012).They have been originally introduced in their sub-surface segmentation method to enlarge small segments.As a result, those segment patches are merged together that share the same or similar plane properties but are disconnected due to the presence of superstructures.For our purpose only virtual points below surface points are taken into account.Since these virtual points are located always below measured surface points, they can be used as an indicator for areas where superstructures are present in the x-y plane.Some examples are given in Figure 3. Figure 3. Segmented surface points enriched with virtual points (red).Unsegmented surface points are colorized in black.However, virtual points as defined above cannot be used alone as a reliable criterion for the presence of superstructures as shown in Figure 4.In fact their related surface points can be considered as a superset of points that are located on a superstructure.Therefore, a selection of these points is mandatory.For this purpose, a connected component analysis (CCA) is carried out to separate the virtual points into different sets so that neighboring virtual points belong to the same point set.The number of point sets can be considered as the maximum number of superstructures on a base roof.Afterwards, a minimal bounding rectangle is calculated in the x-y plane for each point set.It represents the approximate area and the location of the potential superstructure, wherefore we use instead in the following the term candidate area.Point sets with a small candidate area are treated as outliers and discarded at this point.
In order to verify whether a candidate area of a point set comprises a superstructure, the candidate area is slightly enlarged, inversely to the point density of the given input data.Then, all those surface points are selected that are located within the enlarged candidate area.Based on the selected surface points, another minimal bounding rectangle is calculated in the x-y plane and slightly enlarged.A candidate area is only considered as the location of a superstructure if the enlarged bounding rectangle of the surface points encloses all the virtual points of its point set.All other candidate areas and their point sets are discarded.Examples of the selection workflow are given in Figure 5 and Figure 6.In the first example, all virtual points are enclosed by the enlarged minimal bounding rectangle of those surface points that are located within the enlarged minimal bounding rectangle of the virtual points.Therefore, the candidate area of this point set is considered as a location of a superstructure.In contrast, the point set in the second example is discarded because the enlarged bounding rectangle of the selected surface points does not enclose all virtual points in the candidate area.

GROUPING OF SIMILAR SUPERSTRUCTURES
Once all appropriate point sets are detected, all those sets are assigned to the same group that belong to the same type of superstructure.However, due to the low point density and the small of the candidate areas, an automatic extraction of certain features or a shape fitting algorithm is usually not feasible to complete this task.Therefore, our proposed approach utilizes the ICP principles as explained in section 3 to increase the point density for each superstructure type.For this purpose, all remaining virtual points from the previous section are replaced by their related surface points and those candidate areas that feature a similar size are initially assigned to the same group.To start the ICP iterations for a group, an initial rough transformation in form of a rigid transformation is carried out for each of its point sets.The translations ensure that the centroids of the candidate areas of a group meet at one point.Furthermore, the rotations accomplish that the x-y normal vectors of the base roof segments, on which the point sets of a group are located, point afterwards in the same direction.Thereby certain angles are preferred for two point sets of the same group if they are located on the same base roof; e.g.0 degree if both point sets belong to the same segment and 180 degree if the point sets belong to opposite segments.Furthermore, if two point sets are not located on the same base roof, additionally the direction of the base roof is taken into account.An example of the initial rough transformation process is shown in Figure 7.The fine registration process starts thereafter in each group to estimate the point correspondences of the point set with the highest and the second highest number of surface points.Then, a transformation that minimizes the MSE of the point correspondences is calculated and applied.But instead of a rigid transformation, our method restricts ICP to use only translations.Rotations are not considered during the fine registration process.The objective function that has to be minimized in each iteration for two given point sets  = {  } and  = {  } with  = {1, … , } and  = {1, … , } is defined as where   is the closest point in  to the point   ∈  and  the translation vector.Both point sets are merged if the final MSE value is lower than a predefined threshold and the fine registration process starts again.Otherwise the group is split and each of the point set is assigned to one group before the process starts again.In the latter case, the remaining point sets of the original group have to be tested for both groups but only the registration with the lower MSE value is taken.Once all point sets of a group are merged, the final groups whose point sets originated from the same group are registered and merged if the MSE value is lower than a predefined value.This threshold is reduced if the group merging supports the symmetric distribution of the original point sets on the base roof.An example, which separates the given point sets into the two groups G1,1 and G1,2, is given in Figure 8. First, all point sets in the given example are assigned to the same group because their candidate areas feature a similar size.Then, the surface points of the yellow and the purple point sets are tested if they belong to the same superstructure type.Due to their final MSE value, they are separated into different groups.Therefore, the remaining pink point set has to be tested for both groups and is finally assigned to the group of the yellow point set.This is at the same time the final result of the grouping process because the registration of both groups results in a too high MSE value.

SUPERSTRUCTURE MODELING AND CONSTRUCTION
The resulting groups of the previous section are considered during the subsequent modeling process as different superstructure types and the number of merged point sets in a group as the number of instances of it.Instead of reconstructing the superstructure of each original point set independently, all instances of a group are constructed at once.Due to the previous merging process, a higher point density than in the input data is already achieved if more than one instance of a superstructure type has been detected.Moreover, in this case also single gaps in the surface points of an instance, which might occur due to partially occluded areas, are automatically closed by the surface points of other instances.
Therefore, several features like ridge lines can now be directly extracted in each group.Taking into account the extracted features, points of planar regions are segmented and planes estimated.An adjustment of these planes is carried out to emphasize on the one side rectangularity, if the direction of the normal vector projected onto the x-y plane is close to a certain angle, and on the other side symmetries, which is for example often present in the rooftop slope of a saddleback superstructure.Afterwards, half-space modeling, which is a form of solid modeling (Mäntylä, 1988;Foley et al., 1990), is carried out for the construction of the superstructure model.For this, the adjusted planes are formulated as hyperplanes whose normal vectors have positive z-values and each of these hyperplanes are used to create half-spaces.Further half-spaces are added by taking the existing half-spaces into account so that their limiting hyperplanes correspond to the vertical surfaces of the superstructure.Equal eave heights in the final model are favored.Then, instances of the so far defined model are created and placed in each candidate area of the original point set of the group.The location of an instance is defined by the inverse rigid transformation of the previous section, which consists of the initial rough transformation and the translations during the iterations.In conjunction with the surface of the base roof, the final superstructure model is created.Special care has to be taken that the superstructure models do not overlap each other.In such a case, the facade hyperplanes of the overlapping models are translated accordingly.

RESULTS
For evaluation purposes, the proposed approach is tested both on generated data (1.5-4 points/m²) and on several buildings located in residential parts of the Vaihingen test data set (4-6 points/m²), which is provided by the German Society for Photogrammetry, Remote Sensing and Geoinformation (DGPF) (Cramer, 2010).At the moment, our implementation is limited to those superstructures that are located on a single roof plane.For this dominant kind of superstructures, the approach presented in this paper is in general suitable for their automatic reconstruction in low-density point clouds.Some results are presented in Figure 9 and Figure 10.As shown, it is generally able to reconstruct superstructures that would otherwise be missed.In particular, our approach provides accurate models if many instances of the same superstructure type are present.In this case, also partially occluded superstructures can be reconstructed, as shown in Figure 11.These models are furthermore always topological correct and reflect regularities that are usually present in man-made objects.However, special care has to be taken, so that a superstructure model is not placed in an area of point outliers.The size of a candidate area is often a quite sufficient indication but it cannot deal with an accumulation of outliers that might occur in reality.

CONCLUSION
In this paper, we present a fully automatic reconstruction approach for building superstructures from low density point clouds.In contrast to already existing approaches, which usually reconstruct superstructures individually and independently in a model-driven way, our proposed approach helps to find and reconstruct several rooftop details that would otherwise be missed.For this purpose, surface points that belong to the same superstructure type are transformed using an ICP based registration technique so that they overlap each other.The resulting point density is, depending on the number of instances of a superstructure type, higher than in the input data so that features are more likely to be recognized.Therefore, the presented approach is especially suited to large buildings with many superstructures of the same type.In this case, superstructures can be reconstructed in point clouds with a point density of ~1.5 points/m².Furthermore, symmetries, alignments, and regularities can be enforced in a straight-forward way because superstructures are reconstructed in groups.
However, further research is needed to extend the reconstruction capabilities in low density point clouds.For example, the proposed approach is limited to superstructures that belong to a single roof segment.Superstructure point sets that do not comply with this condition cannot be registered at the moment.Furthermore, in some cases a low rotation during the iterations of the grouping process could help to get other point correspondences so that only a reduced number of iterations is needed.But the conditions for a rotation and the allowed rotation circumference have still to be investigated.

Figure 1 .
Figure 1.From top to bottom: segmentation result (unsegmented points in black); segmentation result after adding virtual points (red), and estimation of candidate areas; transformation result of the candidate points; roof superstructure reconstruction result; reconstruction result of the whole building.

Figure 4 .
Figure 4. Result of the subsurface segmentation.Virtual points are colorized in red.

Figure 5 .
Figure 5. From top to bottom: building overlaid with surface points and the CCA result; bounding rectangles (virtual points); bounding rectangles (surface points).

Figure 6 .
Figure 6.From top to bottom: building overlaid with surface points and the CCA result; bounding rectangles (virtual points); bounding rectangles (surface points).

Figure 7 .
Figure 7. Example of the initial grouping and the initial rough transformation process.

Figure 8 .
Figure 8. Example workflow of the fine registration process.

Figure 9 .
Figure 9. Left: Top view of the segmentation results.Unsegmented points are colorized in black.Right: The resulting building model reconstructed by our proposed approach.

Figure 10 .
Figure 10.Some reconstructed buildings of the Vaihingen test data set.

Figure 11 .
Figure 11.Left: Top view of the segmentation result.Unsegmented points are colorized in black.Right: The resulting building model reconstructed by our proposed approach.