3D Indoor Building Environment Reconstruction using Polynomial Kernel, Least Square Adjustment, Interval Analysis and Homotopy Continuation

: Nowadays, municipalities intend to have 3D city models for facility management, disaster management and architectural planning. Indoor models can be reconstructed from construction plans but sometimes, they are not available or very often, they differ from ‘as-built’ plans. In this case, the buildings and their rooms must be surveyed. One of the most utilized methods of indoor surveying is laser scanning. The laser scanning method allows taking accurate and detailed measurements. However, Terrestrial Laser Scanner is costly and time consuming. In this paper, several techniques for indoor 3D building data acquisition have been investigated. For reducing the time and cost of indoor building data acquisition process, the Trimble LaserAce 1000 range finder is used. The proposed approache use relatively cheap equipment: a light Laser Rangefinder which appear to be feasible, but it needs to be tested to see if the observation accuracy is sufficient for the 3D building modelling. The accuracy of the rangefinder is evaluated and a simple spatial model is reconstructed from real data. This technique is rapid (it requires a shorter time as compared to others), but the results show inconsistencies in horizontal angles for short distances in indoor environments. The range finder horizontal angle sensor was calibrated using a least square adjustment algorithm, a polynomial kernel, interval analysis and homotopy continuation.


INTRODUCTION
In 3D GIS, 3D spatial modelling is one of the most important aspects (Chen at al., 2008).3D spatial modelling involves the definition of spatial objects and spatial data models for visualization, interoperability and standards (Chen at al., 2008).Due to the complexity of the real world, 3D spatial modelling leads towards different approaches in different GIS applications.According to Chen et al. (2008), there is not a universal 3D spatial model that can be used in and shared between different applications.Different disciplines according to their input and output use different spatial data models.
The automatic reconstruction of urban 3D models has been a research area of photogrammetry for the past two decades (Haala and Kada, 2010).According to Habib et al. (2010), digital 3D modelling of complex buildings has been a challenge until now with photogrammetry technology.This leads towards semi-automated construction of complex 3D building models.Difficulties of interpretation of photogrammetric images for 3D city modelling, especially for complex buildings, motivated increasing demands for 3D point cloud technologies such as LiDAR (light detection and ranging), which can facilitate automated 3D building models.
According to Surmann et al. (2003), rapid characterization and quantification of complex environments with increasing demands has created a challenge for 3D data analysis.This crucial demand comes from different fields such as industrial automation, architecture, agriculture, construction and mine and tunnel maintenance.Precise 3D data is needed for facility management, factory design and regional and urban planning.Considering all the issues affecting fully automated construction of complex 3D building models, 3D indoor modelling is another aspect in the field of 3D city modelling which can make the current situation more complex.According to Deak et al. (2012), indoor location determination has become a crucial aspect in many different applications but unfortunately, a lack of standard is one of the challenges and there are more challenges encountered in this field.According to Donath and Thurow (2007), considering many fields of applications for building surveying and resulting different demands, representation of the building geometry is the most crucial aspect of a building survey.Due to the complexity of indoor environment, this field needs to be more researched.
Recently, there has been more interest for 3D building modeling based on LiDAR data, but extracting buildings from huge LiDAR datasets is difficult and time consuming and requires experienced technicians.Laser scanning technology started in the 1990s (Amato et al., 2003) and it can measure a 3D object surface with a high speed pulse.This technology is considered as a tool for remote and rapid data collection and it can be used in many different applications from urban and regional planning to architecture.A scanner can directly measure distance and reflection intensity of 3D object surfaces and automatically store collected data in a spatial database.Recent TLS technology can collect more than 500,000 points in a second with an accuracy of ±6 mm (Dongzhen et al., 2009).
Currently, a growing interest for indoor building surveying has been observed in GIS and BIM research communities (Isikdag and Underwood, 2010).The first community is interested in modelling existing building structures for emergency response and disaster management systems.Indoor building surveying is vital especially when other data sources such as research plans and architecture models are not available.The second community is interested in models with 'as-built' conditionsconstruction plans are often different from the final constructed building and it is rare that appropriate plans are available to the builders.In this case, the buildings and their rooms must be surveyed to obtain the locations of walls, edges, corners and their relationship with adjacent spaces (i.e.rooms, corridors).Unfortunately, many methods used for land surveying cannot be easily applied due to lack of GPS signal from satellites in indoor building environment, limited working area inside buildings especially in office space and very detailed environment with furniture and installations.There are four approaches that seem to be suitable for indoor surveying including: 1. Laser scanning which is expensive, time consuming and requires considerable modelling effort in order to fit sections of the surveyed point clouds to basic features such as walls.This results in extensive post-data collection manual work and there is no easy way to integrate individual scan results with the model of a complete complex building.
2. Traditional surveying with Total Station or equivalent is also possible, but conversion of captured data points into a building model is very complex.
3. Third approach, using a light rangefinder which integrates azimuth (from a digital compass) and inclination appears to be the most feasible for surveying indoor building environment, although it has a lower level of accuracy than Total Station and Laser Scanner based surveying approaches (Jamali et al.,2015).

4.
A fourth approach based on photogrammetry technique, uses non-calibrated, non-metric cameras to extract 3D information from photographs.For indoor surveying, it is as simple as taking pictures.Additionally, images can be used for texture extractiontextures can be attached to walls, floors, and ceilings in the model which would increase the realism of visualization.
In this research, we provide a comparative analysis of 3D reconstruction and indoor survey of a building done using the Leica scanstation C10 and the Trimble LaserAce 1000 rangefinder (see Figure 1).The Trimble LaserAce 1000 has been used for outdoor mapping and measurements, such as forestry measurement and GIS mapping (Jamali et al., 2013).A rangefinder can be considered as a basic mobile Total Station with limited functionality and low accuracy.The Trimble LaserAce 1000 is a three-dimensional laser rangefinder with point and shoot workflow.This rangefinder includes a pulsed laser distance meter and a compass, which can measure distance, horizontal angle and vertical angle up to 150 meter without a target and up to 600 meter with a reflective foil target.In this research, we propose this device for indoor mapping and try to validate this technique in an indoor environment.Trimble LaserAce 1000 will decrease time and cost of surveying process ( Jamali et al. 2014).Following this introduction, in Section 2, the requirements of indoor building surveying are discussed.In Section 3, the range finder is calibrated using a least square adjustment algorithm.In Section 4, the range finder is calibrated using a polynomial kernel algorithm.In Section 5, the range finder is calibrated using interval analysis and homotopy continuation in order to control the uncertainty of the calibration and of the reconstruction of the building.Section 6 presents conclusions and future research.

REQUIREMENTS OF AN INDOOR BUILDING SURVEYING METHOD
In this research, the surveyor according to his experience and knowledge defined several requirements for indoor building surveying as follows: i. Number of control points and their positions To get better results with less shape deformation (e.g.intersection and gap between two rooms due to low accuracy of Trimble LaserAce 1000), for each door, there should be a control point inside the corridor.For example, if a room has two doors, there should be two control points in the corridor to access that certain room by its two doors (see Figure 2).To avoid narrow and wide angle propagation which will lead to shape deformation, Trimble LaserAce 1000 needs to be stationed in the center of each room.In reality, there will be some furniture in rooms which avoids putting Trimble LaserAce 1000 in the center of a room.Thus, for a fully furnished room, there should be at least two control points.

ii. Number of repetitions of measurements
Based on the surveyor's skills, for a certain distance and angle measurement, there should be at least three observations.Trimble LaserAce 1000 shows inconsistency in horizontal angle in an indoor building environment.By increasing number of repetitions for a measurement, the probability of a mistake or blunder occurring will be decreased.In this paper, we use least squares adjustment and polynomial kernel as statistical methods and interval analysis and homotopy continuation as mathematical methods to calibrate our low accuracy surveying equipment in an indoor building environment.Thus, a higher number of measurements will led to a more accurate average using least square adjustment and more accurate measurement intervals using interval analysis and homotopy continuation.
iii.Time Time is one of the most important factors in this research.The cost of a surveying project will be decreased by reducing the time of data collection.There should be a balance between the number of repetitions and the time to collect data of a room.The time for collecting data of a room with three repetitions (horizontal angle, vertical distance, horizontal distance and vertical distance) and one control point using Trimble LaserAce 1000 is between 5 to 10 minutes based on our experiments.As can be seen in Figure 3, Trimble LaserAce 1000 is handy and easy to use, which can decrease the time of data collection.Based on the surveyor experience, Rangefinder is 2 times quicker compared to a traditional Total Station (e.g.Leica 307 TCR) and 10 times quicker compared to a terrestrial laser scanner (e.g.Leica Scan Station C10).

RANGEFINDER CALIBRATION
Coordinates measured by rangefinder are not as precise as laser scanner or total station measurements.As seen in Figures 5 and 6, results of Trimble LaserAce 1000 shows deformation of building geometry.Figure 7 shows a 3D point cloud collected by Leica scanstation C10.
According to device specifications, the accuracies of the Leica scanstation C10, Trimble LaserAce 1000 are as shown in Table 1.called the problem data.The vector y is called the measurement vector, and the matrix A the design or input matrix.The vector  ∶=  −  is called as the residual error vector.The input matrix A is written as = [  1 ⋯   ], where  ∈   is the  ℎ column of A, j=1, …, n.
The problem is to compute the x that minimizes the sum of the squares of the residuals We are trying to find the best approximation of y in terms of a linear combination of the columns of A. Thus, the least square problem amounts to project (find the minimum Euclidean distance) the vector y on the span of the vectors  's (that is to say the range of A).The 3D building measured by the Trimble LaserAce 1000 can be calibrated and reconstructed from the Leica scanstation C10 based on the least square adjustment algorithm, in the form of absolute orientation.Least square adjustment is a well-known algorithm in surveying engineering which is used widely by engineers to get the best solution in the sense of the minimization of the sum of the squares of the residuals, which is obtained as in the following normal equations, which express that the total differential of the sum of squares of residuals is zero.Least square adjustment for a linear system is shown in Equation ( 1).
Where L = observations X = unknowns A = coefficient of unknowns W=observation's weight N = (A T W A) Considering two points, Pa= (XA, YA, ZA) from the Leica C10 and Pc= (XC, YC, ZC) from the Trimble LaserAce 1000, the absolute orientation problem can be defined as the transformation between two coordinates systems (Leica C10 and Trimble LaserAce 1000).The relationship between measuring devices can be solved by using absolute orientation.Absolute orientation can be found by a set of conjugate pairs: {(Pc,l, Pa,l), (Pc,2 Pa,2), ... , (Pc,n, Pa,n)}.For a pair of common points in both (camera coordinates and absolute coordinates) systems; rotation, scale and translation components can be calculated by Equations 2 to 4, where the matrix R with coefficients RXX, RXY, RXZ, RYX, RYY, RYZ, RZX, RZY and RZZ, is the matrix of a linear transformation combining a 3D rotation (that can be decomposed into the combinations of 3 rotations along the x, y and z axes) and a scaling, and its determinant is the scaling parameter (since the determinant of a rotation matrix must equal 1).
XA=RXX XC + RXY YC + RXZ ZC + TX (2) YA=RYX XC + RYY YC + RYZ ZC + TY (3) ZA=RZX XC + RZY YC + RZZ ZC + TZ (4) Twelve unknown parameters, including nine linear transformation (combined rotation and scaling) parameters and three translations components need to be solved.Each conjugate pair yields three equations.The minimum number of required points to solve for the absolute orientation is thus four common points.Practically, to get better results with higher accuracy, a higher number of points need to be used.The coordinates of the points measured by the rangefinder can be adjusted, or their maximum error can be minimized, by adjusting the coefficients of the rotation matrix R and the translation vector (see adjustment results in Table 2).Room number nine has been selected by the researcher to calculate its absolute orientation parameters.Absolute orientation can be found by computing the rotation matrix R and the translation vector for any given point.Any points measured by the rangefinder can be transferred or absolutely oriented by using the corresponding matrix A arrays.Results from calibrating the Trimble LaserAce 1000 based on the least square adjustment (Absolute orientation) using the Leica C10 data were calculated (see Table 3).8).
Since the highest measurement uncertainties are those of the horizontal angles measured by the magnetometer of the rangefinder, we focus on the calibration of these magnetometer measurements.
In order to check the applicability of least squares to range finder calibration, we have applied the least squares linear model resulting from the calibration of horizontal angles of room 1 to rooms 10 and 9.The results are shown in Tables 4, 5 and 6.We can observe that the results for room 9 are not satisfactory at all.The least squares methods used in this section assume a linear statistical model of propagation of the errors and a normal probability distribution function of the measurements.However, in any real measurement experiment, one can observe that no probability distribution function actually fits the data set to any desired degree of accuracy.In Section 5, we will see how we can relax these assumptions.

POLYNOMIAL KERNEL
In this section, polynomial kernel as a non-parametric non-linear method is used to calibrate Trimble LaserAce 1000 horizontal angle.In kernel methods, data is embedded into a vector space and then linear and non-linear relations between data will be investigated.Basically, kernels are functions that return inner products between of the images of data points in vector space.Each kernel k(x,z) has an associated feature mapping  which takes input x∈X (input space) and maps it to F (feature space).Kernel k(x,z), takes two inputs and gives their similarity in F space.
F needs to be a vector space with a dot product to define it, also called a Hilbert space (Browder and Petryshyn, 1967;Akhiezer and Glazman, 2013).For k to be a kernel function, a Hilbert space F must exist and k must define a dot product and be positive definite.

∫ 𝑑𝑥∫ 𝑑𝑧𝑓(𝑥)𝑘(𝑥, 𝑧)𝑓(𝑧) > 0
There are several kernel function including linear kernel and Radial Basis Function (RBF) kernel In this section, we use polynomial kernel to find similarities between the rangefinder horizontal angles and the total station horizontal angles and then, to use calculated parameters to calibrate the rangefinder horizontal angles.Table 7 represents results of calibrated horizontal angle of Trimble LaserAce 1000 using different degrees of polynomial kernel, while Table 8 represents the analog results for room 10 using least square parameters calculated for room 1.As can be seen in Table 8, these results are not satisfactory.The results for room 9 are even less satisfactory (thus, omitted).

INTERVAL ANALYSIS AND HOMOTOPY CONTINUATION
Interval analysis is a well-known method for computing bounds of a function, being given bounds on the variables of that function (E.Ramon Moore and Cloud, 2009).The basic mathematical object in interval analysis is the interval instead of the variable.The operators need to be redefined to operate on intervals instead of real variables.This leads to an interval arithmetic.In the same way, most usual mathematical functions are redefined by an interval equivalent.Interval analysis allows one to certify computations on intervals by providing bounds on the results.The uncertainty of each measure can be represented using an interval defined either by a lower bound and a higher bound or a midpoint value and a radius.
In this paper, we use interval analysis to model the uncertainty of each measurement of horizontal angle and horizontal distance done by the range finder.We represent the geometric loci corresponding to each surveyed point as functions of the bounds of each measurement.Thus, for distances observed from a position of the range finder, we represent the possible position of the surveyed point by two concentric circles centered on the position of the range finder and of radii the measured distance plus and minus the uncertainty on the distance respectively (see Figure 9).For horizontal angles observed from a position of the range finder, we represent the possible position of the surveyed point by two rays emanating from the position of the range finder and whose angles with respect to a given point or the North are the measured angle plus and minus the uncertainty on the horizontal angle respectively (see Figure 9).Therefore, the surveyed point must be within a region bounded by these 4 loci: in between 2 concentric circles and 2 rays.Proceeding in the same way for each room, we get the geometric loci for each room and for the union of the surveyed rooms (see Figure 9).
A homotopy is a continuous deformation of geometric figures or paths or more generally functions: a function (or a path, or a geometric figure) is continuously deformed into another one (Allgower and Georg, 1990).The use of homotopies can be tracked back to works of Poincaré (1881-1886), Klein (1882-1883), and Berstein (1910) (Allgower and Georg, 1990).The use of homotopies to solve non-linear systems of equations may be traced back at least to Lahaye (1934) (Allgower and Georg, 1990).A homotopy between two continuous functions f and f from a topological space X to a topological space Y is defined as a continuous map H: X × [0, 1] → Y from the Cartesian product of the topological space X with the unit interval [0, 1] to Y such that H(x, 0) = f0, and H(x, 1) = f1, where x ∈ X.The two functions f0 and f1 are called respectively the initial and terminal maps.The second parameter of H, also called the homotopy parameter, allows for a continuous deformation of f0 to f1 (Allgower and Georg, 1990).Two continuous functions f0 and f1 are said to be homotopic, denoted by f0 ≃ f1, if, and only if, there is a homotopy H taking f0 to f1.Being homotopic is an equivalence relation on the set C(X, Y) of all continuous functions from X to Y.
In this paper, we used homotopy to calibrate the range finder.The main idea is that assuming a linear model and a normal probability distribution function, we only assume that the calibration of the set of our range finder measurements with respect to the set of measurements of our total station can be done continuously, because there is no discontinuity in the n-dimensional space corresponding to the space of measurements performed using the range finder and the total station.Even though, not all real numbers are representable in a digital measurement device, we can assume that all the real numbers corresponding to measurements can be obtained physically, and it is just the fixed point notation used by the digital measurement device, that limits the set of representable real numbers to a discrete subset of the set of real numbers.Thus, we can compute the calibration of the range finder as a continuous function mapping our measurements obtained using our range finder to the measurements obtained using our total station.
The results of the linear homotopy continuation are presented in Figures 10 and Tables 10-13.One can calibrate the differences of horizontal angles observed with the rangefinder to the differences of horizontal angles observed with the theodolite.One can start from any point and point and assume that the measurement of the horizontal angle of that point by the rangefinder will not be changed by the calibration process.Without loose of generality, this point can be the first observed point.Now the idea for the calibration is that we are using each one of the intervals between measurements of horizontal angles made with the rangefinder, and we calibrate the new measurements of horizontal angles made by the rangefinder in each one of these intervals as a non-linear homotopy, where the homotopy parameter is the relative position of the measured horizontal angle in between the bounds of the enclosing interval of rangefinder horizontal angles.This homotopy calibration can be visualized as the continuous deformation of each sector (defined by the rangefinder horizontal angle intervals of room 1) of a plastic disk (corresponding to the old time theodolite graduated disk) to the corresponding sector of the total station's theodolite graduated disk.We used the intervals of total station horizontal angles of room 1 as reference intervals for all other horizontal angle measurements in rooms 1, 9 and 10.In the remainder of the paper, "the homotopy parameter of a horisontal angle measurement" is equivalent to "the relative position of the horizontal angle measurement in the corresponding reference interval".We fitted a polynopmial of degree 5 through the 4 points whose x-coordinates are the homotopy parameters of the horizontal angles measured by the rangefinder in room 10 and the corresponding homotopy parameters of the horizontal angles measured by the total station in room 10 and the points (0,0) and (1,1).We used this polynomial as the convex homotopy functlon that maps the uncalibrated homotopy paramter to the calibrated one.
The initial and terminal maps correspond respectively to the mappings between the uncalibrated and calibrated horizontal angles at the start point and the end point of the enclosing interval of horizontal angles measured by the range finder.We can observe that, contrary to the least squares calibration, the only limitation of this interval analysis and homotopy continuation based calibration is the precision of the fixed point arithmetic used by the computing device used for the calibration.The results shown in Tables 11 and 12 prove that homotopies give the best results both in terms of RMSE and the L ∞ metric.The proposed research is to demonstrate the feasibility of interior surveying for full 3D building modelling..The main objective of this research was to propose a methodology for data capturing in indoor building environment.A rangefinder was compared to a high accurate surveying device (Leica scanstation C10) using weighted least squares, polynomial kernel and a novel technique based on interval analysis and homotopy continuation.In an indoor environment, the Trimble LaserAce 1000 showed inconsistencies within the uncertainty ranges claimed by the manufacturer for short distances in the horizontal angle.Rangefinder data was calibrated by least square adjustment (absolute orientation) which shows a maximum error of 13 centimeters and a minimum error of 6 centimeters using the Leica scanstation C10 as a benchmark.By opposition, the combined interval analysis and homotopy continuation technique calibration obtained by continuous deformation of the function mapping the rangefinder measurements to the theodolite measurements allows a much better match, whose only limitation is the fixed point arithmetic of the computing device used to perform the computation.Results from polynomial kernel are not satisfactory (see Tables 8 and 11).

Figure 2 :
Figure 2: Position of control points by Trimble LaserAce 1000

Figure 3 :
Figure 3: Time of data collection and 3D building modeling

Figure 4 :
Figure 4: Equipment comparison based on the cost.

Figure 8 :
Figure 8: Model calibrated and reconstructed based on the least square adjustment; rangefinder (red dash lines), total station (blue lines) and non-calibrated rangefinder (black lines).

Figure 9 :
Figure 9: The geometric loci of each corner of a room as a function of all the measurements

Table 2 :
Coefficient of unknowns including rotation, scale and translation parameters (matrix A).

Table 3 :
LaserAce 1000 calibration based on the least square adjustment (Absolute orientation).Considering the Leica C10 data as absolute coordinates, differences between two coordinates can be referred as the Trimble LaserAce 1000 accuracy.The model calibrated and reconstructed using the Leica C10 is shown in Figure.8.Model in black lines represents model reconstructed from raw data of Trimble LaserAce 1000 and model in blue lines represents model reconstructed from Leica C10.Calibrated model of Trimble LaserAce 1000 based on the least square adjustment algorithm from Leica C10 data can be seen as red dash line model (see Figure

Table 4 :
Calibration of room 1 rangefinder horizontal angle measurements using total station horizontal angle measurements by least squares.

Table 5 :
Calibration of room 10 rangefinder horizontal angle measurements using total station horizontal angle measurements by least squares calculated from room 1.

Table 6 :
Calibration of room 9 rangefinder horizontal angle measurements using total station horizontal angle measurements by least squares calculated from room 1.

Table 7 :
Calibration of room 1 rangefinder horizontal angle measurements using total station horizontal angle measurements by polynomial kernels.

Table 8 :
Calibration of room 10 rangefinder horizontal angle measurements using total station horizontal angle measurements by polynomial kernel calculated from room 1.

Table 9 :
Calibration of room 1 rangefinder horizontal angle measurements by homotopies.

Table 10 :
Calibration of room 10 rangefinder horizontal angle measurements by homotopies

Table 11 :
Calibration of room 10 rangefinder horizontal angle measurements using total station horizontal angle measurements by homotopies, least squares and polynomial kernel.

Table 12 :
Calibration of room 9 rangefinder horizontal angle measurements using total station horizontal angle measurements by homotopies, least squares and polynomial kernel.