- Email: [email protected]
Assessing planar asymmetries in shipbuilding from point clouds
Accepted Manuscript Assessing planar asymmetries in shipbuilding from point clouds Javier Roca-Pardiñas, Celestino Ordó ñez, Carlos Cabo, Agusín MenéndezDíaz PII: DOI: Reference:
S0263-2241(16)30746-1 http://dx.doi.org/10.1016/j.measurement.2016.12.048 MEASUR 4510
To appear in:
Measurement
Received Date: Accepted Date:
27 October 2015 26 December 2016
Please cite this article as: J. Roca-Pardiñas, C. Ordó ñez, C. Cabo, A. Menéndez-Díaz, Assessing planar asymmetries in shipbuilding from point clouds, Measurement (2016), doi: http://dx.doi.org/10.1016/j.measurement.2016.12.048
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Assessing planar asymmetries in shipbuilding from point clouds Javier Roca-Pardi˜ nasa , Celestino Ord´ on ˜ezb , Carlos Cabob , Agus´ın Men´endez-D´ıazc a Department
of Statistics and Operations Research, University of Vigo, 36208 Vigo, Spain of Mining Exploitation and Prospecting, University of Oviedo, 33600 Oviedo, Spain c Department of Construction and Manufacturing Engineering, University of Oviedo, 33203 Gij´ on, Spain Author for the correspondence: Celestino Ord´ on ˜ez, E-mail address: [email protected], Telephone number: +34985458027, Fax number: +34985458188 b Department
Abstract This paper introduces a statistical methodology to identify planar asymmetries in surfaces. An algorithm based on bivariate kernel smoothers is used to estimate non-parametrically the surface from a point cloud. Also bootstrap-based procedures are proposed for testing a null hypothesis establishing that the surface has planar symmetry. After checking the validity of the proposed methodology with simulated data, it was applied to study symmetries in a boat hull. Keywords: symmetry, bootstrap, kernel smoothing, laser scanner, point cloud
1. Introduction Symmetry is a property present in nature but also in artificial objects, such as architectural models and industrial products [1]. It plays an important role in human visual perception and aesthetics [2] but also in the economy and efficiency of structures and animate or inanimate organisms [3]. From a mathematical point of view, a symmetry is a set of transformations applied to an object that preserve the properties of the structure [4]. Translations, rotations and reflections are typical geometric transformations corresponding to symmetries. Although it is widely associated to geometric transformations [5], symmetry concerns other disciplines, such as music [6]. The importance of symmetry has originated a number of studies concerning this topic throughout time. The authors in [7] describe a variety of applications of symmetry in digital geometry processing and shape analysis. They also analyse the strengths and limitations of some algorithms for symmetry detection. Some examples of symmetry detection in images can be seen in [8–10]. In the architectural and industrial fields, symmetry is also a prime characteristic of
Preprint submitted to Measurement
December 31, 2016
some objects that can be used in reverse engineering. In [11] a number of aspects of the geometry of mechanical parts, such as geometric restrictions, slots, holes or symmetries, which should be considered to improve boundary representation models, are described. The authors in [12] analyse the cognitive process involved in human symmetry detection to develop a computational method to detect symmetries. In [13] the authors developed a methodology to detect symmetries on a vessel from a point cloud measured with a terrestrial laser scanner that has some similarites with that presented in this paper. However, in this work we extend the analysis to the whole surface of the vessel, without restricting it to cross-sections. In addition, it has other significant differences related to the bootstrapping method used, the application of binning techniques to reduce computational time, and the analysis of the effect of the standard deviation in the error term. Symmetries are also fundamental in the design of objects. Thus, symmetries were used to detect the design intent in approximate CAD models [14]. In [15] symmetries in CAD models are detected using Boolean combinations of the object features. Symmetry is also a basic concept when parts are assembled in one piece [16, 17]. This paper addresses the problem of contrasting the symmetry of an object surface function. There are several methods in the literature that tackle the symmetry of objects or functions. Some of these methods are centered on determining a symmetry axis [18, 19] or [20], among others. In most of these papers it is admitted that the object surface function is observed without stochastic errors. Other authors analyse if an image or data are symmetric with respect to a known or unknown axis. Recently, tests for reflection symmetry and rotational invariance of two-dimensional functions or images based on Zernike polynomials were proposed [21]. A symmetry test for a regression function with a bivariate predictor based on the L2 distance between the original function and its reflection was proposed in [22]. In this work we propose a slightly different vision of the regression itself, where errors are considered in three dimensions, not just in the Z component. Our research aims to provide a methodology to detect planar symmetries in objects from point clouds (obtained from laser scanners, close-range photogrammetry or other measurement techniques). The method is based on fitting a surface to the point cloud using a non-parametric orthogonal regression method, sectioning the surface into two halves and comparing these using confidence intervals for the difference between each half obtained by bootstrapping. Thus, it can be determined whether each section can be considered symmetric in relation to a hypothetical plane of symmetry. The hypothesis test must be able to accept the null hypothesis when it is true while it has enough power to reject that hypothesis when it is false. This work also includes linear binning as an efficient technique to manage point clouds of several million points, since these massive point clouds are quite frequent nowadays. The rest of the paper is structured as follows: Section 2 introduces the problem of detecting planar symmetry in shipbuilding. The methodology used to adjust a surface to the point cloud is explained in Section 3. Section 4 shows 2
the methododolgy developed to detect symmetries in the surface. In Section 5 a simulation study to validate the method is carried out. In Section 6 the proposed methodology is used to detect symmetries in a yacht hull from a point cloud obtained using a terrestrial laser scanner. Finally the most relevant conclusions of the research are collected. 2. Testing symmetries in a boat: problem statement The methodology proposed in this paper was used to determine planar symmetry in a yacht measuring 12 m long, with a beam of 4.5 m and a draft of 2.4 m (see Figure 1). Symmetry along a vertical plane containing the fore-aft axis is important for correct sailing behavior, as asymmetry can cause a shift in the volume of displaced water and in a different thrust on the port and starboard sides of the hull, causing the vessel to list. Another possible drawback of asymmetries is a lateral thrust that causes the vessel to veer from a straight-line course. The relevance of these effects depends on the planned use of the vessel, its characteristics and dimensional requirements. The use of point clouds obtained from laser scanners or close-range photogrammetry for reverse engineering problems, such as the detecion of symmetries, is currently common in industry applications apart from shipbuilding. [23].
Figure 1: Images of the yacht hull.
The yacht hull was scanned using a Leica Cyrax 2500 time-of-flight terrestrial laser scanner (TLS) able to record up to 1000 points per second with a nominal accuracy of ±6mm(1σ) in position and ±4mm(1σ) in distance. The measurent range is 1.5 to 50 m. The scan density (minimum point-to-point measurement spacing) can be modified, being its maximum value less than 1 mm. Terrestrial laser scanners determine the position of the points calculating the distances and the horizontal and vertical angles between the center of the scanner and the points. As it name indicates, time-of-flight scanners measured distances to the points determining the time needed for a laser-beam to travel from the laser emissor to the object and back. Most TLS devices use a rotating mirror located inside the instrument to reflect the emited laser-beam in different 3
directions, allowing to emit a great number of laser-beams in a very short time. Consequently, thousand of points can be measured in a few seconds. This type of laser scanners are able to measure distances of up to several hundred meters under ideal circumstances. For its part, phase-based scanners operate in a different way: they estimate the phase difference between transmitted and reflected waves. These scanners range is shorter than that of the time-of-flight scanners, but, on the contrary, have greater precision and are also much much faster [24]. The scanner was placed in 12 different locations to measure the entire boat. All the scans were registered in the same coordinate system using magnetic reflective targets that were automatically recognized. Using a scan density about one point per square centimeter, aproximatelly 1,600,000 points were necessary to cover the entire yatch hull. Each of these points {P1• , . . . , Pn• } with Pi• = (Xi• , Yi• , Zi• ), being n the number of observations, can be expressed as Pi• = Pi + i
(1)
where ! X YPi Z= (xi , yi , zi ) represents the point of the object surface and i = εi , εi , εi the three dimensional measurement error on the ith point [25]. The only restriction for the distribution of the error term is that it has to have zero mean. Being (x, y, z) the spatial coordinates of each point on the object surface, coordinate z is obtained adjusting a smooth surface z = m (x, y) to the data. Given the point cloud, the first step in our method is to reconstruct the ship hull obtained through the estimate m ˆ of that surface. Once the boat hull surface has been reconstructed, the next step is to detect asymmetries in this surface. Specifically, we are interested in determining if the ”left” and ”right” parts of the surface are symmetric in relation to some unknown plane given by x = x0 . Testing planar symmetry can considered as a hypothesis test where the null hypothesis is H0 : ∃ x0 so that m(x0 + x, y) = m(x0 − x, y) ∀(x, y)
(2)
This equation is valid for symmetries with respect to a vertical plane, and it should be modified to tackle non-vertical symmetry planes. Planar symmetry is present in many man-made objects; particularly, planar symmetry is important in ship construction where, except in rare cases, a ship’s hull is symmetrical about a fore and aft plane. Symmetry affects speed navigation, righting position, stability and maneuverability. 3. Non-parametric estimation algorithm. This section describes the algorithm used to fit the unknown function m to the point cloud {P1• , . . . , Pn• }. For any given point (x0 , y0 ) , a smoothed version of the orthogonal regression (or Deming regression [26]) is proposed. In this type of regression the estimated surface constitutes the best fit to the set of the
4
observed sample points {P1• , . . . , Pn• } by minimizing the perpendicular distances from the points to the surface. The proposed estimators are based on the fact that surface (x, y, m(x, y)) can be approximated by the plane ax + by + cz + d = 0 in values (x, y) near (x0 , y0 ). Thus, the estimated coefficients (ˆ a, ˆb, cˆ) correspond to the third local principal component [27, 28]. The estimate of the true surface m(x0 , y0 ) at each point (x0 , y0 ) is obtained as follows: First, compute the sample ˆ obtained from the point cloud (P • , . . . , Pn• ) weighted 3 × 3 covariance matrix Σ 1 by (W1 , . . . , Wn ) with
Wi =
Wh (Xi•
−
x0 , Yi•
(Xi• − x0 )2 + (Yi• − y0 )2 − y0 ) = exp − h
(3)
ˆ as and obtain (ˆ a, ˆb, cˆ) as the third principal component of the PCA using Σ the input covariance matrix. Then, calculate a plane normal to this vector and passing through (X¯• , Y¯• , Z¯• ) a ˆ(x − X¯• ) + ˆb(y − Y¯• ) + cˆ(z − Z¯• ) = 0 Finally, evaluating that plane in (x0 , y0 ) and solving in z, an estimation of m(x0 , y0 ) is obtained m(x ˆ 0 , y0 ) = Z¯• −
a ˆ(x0 − X¯• ) + ˆb(y0 − Y¯• ) cˆ
(4)
Moreover the vector
a ˆ, ˆb, cˆ
u ˆ(x0 , y0 ) = p a ˆ2 + ˆb2 + cˆ2
(5)
represents the estimated orthogonal vector to the surface z = m(x, y) at the point (x0 , y0 ). Note that the kernel weight Wi in (10) depends on the euclidian distance between (x0 , y0 ) and (Xi• , Yi• ) and, in addition, contains a smoothing parameter h. See [29–31] for information about kernel smoothing. Instead of kernel smoothers, other methods such as penalized and plate regression splines could have been used[32, 33], respectively. 3.1. Bandwidth selection According to the definition of Wh , the influence on the estimate m(x ˆ 0 , y0 ) depends on the distance of (Xi• , Yi• ) to (x0 , y0 ). The smoothing parameter h in (10), generally called the bandwidth, is a critical factor that controls the extension of that influence. When h is small the estimate m(x ˆ 0 , y0 ) tends to adjust too much to the data while it produces a very high variance. On the contrary, when h is large the estimated surface tends to approximate a plane, so the the adjustment to the object is quite poor, being the variance small. 5
Determining an adequate bandwidth is a classical problem in kernel smoothing. Sometimes it is given manually according to the knowledge of the problem. Other times it is obtained automatically by means of cross-validation. In the latter case the bandwidth h is selected by minimizing
CV (h) =
n X i=1
min x,y
q
(Xi• − x)2 + (Yi• − y)2 + (Zi• − m ˆ −i (x, y))2
(6)
where m ˆ −i (x, y) indicates the fit of m(x, y) leaving out the i-sample data (Xi• , Yi• , Zi• ). 4. Testing symmetry In the previous section the mechanism to allow the nonparametric reconstruction of the smooth surface z = m(x, y) has been proposed. In this section we discuss the problem of testing for symmetry of the function m. For testing the null hypothesis in (2) we propose to use the test statistic given by Z 2 T = min (m(x ˆ 0 − x, y) − m(x ˆ 0 + x, y)) dxdy) (7) x0
Note that if H0 is verified, the value of T value should be close to zero. For a given significance α, the decision rule based on T consists of rejecting the null hypothesis if T > T α , where T α is the (1 − α)-percentile of T under the null hypothesis. Nevertheless, within a non-parametric regression context, the theory for determining such percentiles is not closed, and resampling methods such as bootstrap introduced by Efron [39] can be applied instead. Here, we propose an ”orthogonal” version of the traditional wild bootstrap [40, 41]. The steps of the proposed procedure are as follows: ˆ Step 1: Given the original sample P1• , . . . , Pn• obtain the estimated surface m as explained above. Then obtain the symmetric estimated surface m ˆ S under the null hypothesis of symmetry. In particular, for a given point (x, y) the estimated symmetric surface m ˆ S (x, y) can be obtained as the mean of the ”left” and ”right” estimated curves. Explicitly, the expression of m ˆ S (x, y) is given by m ˆ S (x, y) = 0.5m(x, ˆ y) + 0.5m(2ˆ ˆ x0 − x, y) where x ˆ0 is the value of x0 that minimizes the expression in (7). Step 2: Obtain the orthogonal distances to the symmetric surface q 2 2 2 Di = (Xi• − x0i ) + (Yi• − yi0 ) + (Zi• − m ˆ S (x0i , yi0 )) ! where x0i , yi0 is obtained as !
q 2 ˆ S (x, y)) x0i , yi0 = arg min (Xi• − x)2 + (Yi• − y)2 + (Zi• − m x,y
6
(8)
! ! Note that Pi0 = x0i , yi0 , m ˆ S x0i , yi0 is the point on the estimated symmetric surface m ˆ S whose distance to the measured point (Xi• , Yi• , Zi• ) is minimal. Moreover, it is easy to show that the orthogonal vectors to the surface at points (x, y) and (2x0 − x, y) are related, under the null hypothesis, by means of a 90 degree rotation about y axis us (x, y) = T · u(2x0 − x, y) being
0 T= 0 1
0 1 0
-1 0 0
Bearing this in mind, for each (x0i , yi0 ) the estimate of the true u(x0i , yi0 ) under H0 can be obtained as follows u ˆS (x0i , yi0 ) = u ˆ(x0i , yi0 ) + T · u ˆ(2ˆ x0 − x0i , yi0 ) Step 3: For b = 1, . . . , B (e.g. P1•,b , . . . , Pn•,b in accordance with
B=1000) generate the bootstrap samples
Pi•,b = Pi0 +
u ˆS (x0i , yi0 ) · Ri ||ˆ uS (x0i , yi0 )||
where || · || represents the euclidean norm, and Ri being ( √ Di · (1−2√5) with probability p = Ri = Di · (1+2 5) with probability p =
√ 5+ 5 10 √ 5− 5 10
that is, Ri follows a distribution that resembles the Mammen distribution [42]. Then, compute T •b with the generated bootstrap in the same way that T had been obtained with the original sample. Step 4: Finally, the test rule consists of rejecting the null hypothesis if T > Tˆα , where Tˆα is the empirical (1 − α)-percentile of values T •1 , . . . , T •B obtained before. 5. Computational aspects Non-parametric regression modelling can be a very time consuming process when, as in our case, the number of data n is very large. Accordingly, some thecniques aimed to reduce this computing time, such a binning, has been proposed. Information regarding bining thecniques for density estimation can be found in [34–37]. Specifically, in this work we use the binning algorithm for multivariate kernel estimation proposed by [38]. ˜ • < ... < X ˜ • , Y˜ • < ... < Y˜ • and Z˜ • < ... < Z˜ • being M equidistant Let X 1 1 1 M M M points along the X, Y and Z directions respectively. The linear binning assigns 7
˜ • , Y˜ • , Z˜ • a weight that depends on the relative distances to each grid point X j k l between the observations and the eight closest grid points according to the expression
˜ •i W jkl
• • ˜ • • ˜• ˜ • − Y − Z Xi − X Y Z j i i k l 1 − 1 − = 1 − δX δY δZ
+
+
(9) +
with x+ = max {0, x} and δX , δY and δZ the distances between two neighbouring knots on X, Y , Z, respectively. The binning version of the orginal estimate m(x ˆ 0 , y0 ), which expression is given in (4), is obtained as previously explained, ˆ from the M 3 grid points P˜ • = but computing now the covariance matrix Σ jkl • ˜ • ˜• ˜ X , Y , Z with associated weights j
k
l
• ˜ jkl (x0 , y0 ) W
˜ • − x0 )2 + (X ˜ • − y0 )2 (X j k = exp − h (
)
• ˜ jkl W
(10)
In the same way, the bandwidth selection procedure the CV criterion in (6) can be approximated by
CV (h) ≈
M X M X M X
• ˜ jkl W
j=1 k=1 l=1
· min x,y
r
2 ˜ • − x)2 + (Y˜ • − y)2 + Z˜ • − m (X ˆ −(j,k,l) (x, y) j k k
where m−(j,k,l) (x, y) indicates now the binning approximation fit at (x, y) leaving out the (j, k, l)-th binning data point. The number of grid points selected reflects a balance between the degree of adjustment and the computer time. The distribution of the covariates and the regularity of the surface m, also have influence on the suitable number of grid points. Particularly, we propose to choose a number of nodes M0 which meets that not significant differences exists between the estimation obtained with this value and that obtained with a greater number of nodes (M > M0 ). 6. Simulation Study In order to assess the behavior of our proposed algorithms it was applied to simulated data. The procedure is performed in four steps, which are described as follows. 1. The first step consists in the definition of a theoretical mathematical surface Z = m(X, Y ). In our case we have defined a polynomial surface of degree 2 in (X, Y ) Z = m(X, Y ) = a · X 2 + Y 2 where a is a real constant, being X and Y independent variables with uniform distribution on [−1, 1], 8
2. In the second step each point measured with the scanning is simulated according to ! Y Z (Xi• , Yi• , Zi• ) = (Xi , Yi , Zi ) + εX i , εi , εi ; i = 1, ..., n ! Y Z where the error term εX is supposed to be uniformly distributed i , εi , εi on a sphere of a random radius distributed in accordance with a N (0, σ). Several values of σ (from σ = 0.01 to σ = 0.05) were considered in order to obtain different error variances. 3. The proposed estimation algorithm is applied to the point cloud, obtaining an approximation m(x, ˆ y) to the theoretical surface in each point (x, y). To evaluate the behaviour of the estimates we consider a grid of 400 equidistant points in the interval [−1, 1], given by 2(i − 1) 2(j − 1) (ui , vj ) = −1 + , −1 + f or 1 ≤ i, j ≤ 20 19 19 At each grid point (ui , vj ) the nearest point in the true surface is computed q ! 0 0 2 ˆ y)) ui , vj = arg min (ui − x)2 + (vj − y)2 + (m(ui , vj ) − m(x, x,y
Finally, we consider the mean euclidean error (MEE) given by M EE = 400−1
20 q 20 X X ! 2 ! ! 2 2 (ui − u0i ) + vj − vj0 + m(ui , vj ) − m ˆ u0i , vj0 i=1 j=1
Table 1 summarizes numerical average results of the considered error (MEE) over 1000 replicated samples for two different sample sizes (n = 1 × 106 and n = 10 × 106 ) as well as the computational time in brackets (using a PC with a processor 2.53 GHz and RAM 4 GB.), for an increasing number of grid points M . In this table an error variance σ 2 = 0.0062 has been considered. Moreover, the MEE and the limits for the 95% simulation interval are shown in 2 for different σ, considering n = 1 × 106 and M = 100. From Table 1 we can see that the choice of M is a trade-off between the computational time and the error when estimating the true surface. As M increases, the M EEs decrease, but the computing time may increase substantially. On the other hand, the error reduction is very small from M = 100 on. Moreover, it can be appreciated that errors for n = 1 × 106 are similar to those obtained for n = 10 × 106 . That is, there is not much gain in increasing the scan density after a certain level. Finally, as can be seen in Table 2, M EE increases with σ, as was expected. Figure 2 depicts a section of the estimated surface (center) and zoom of different regions (right and left). The continuous black line represents the theoretical surface, the red line the estimation and the gray line each of the bootstrap repetitions for different values of M (M =40, 100 y 200). Finally, we explore the validity of the test considering the null hypothesis H0 of symmetry given in (2). Different values were considered for a, ranging 9
n = 1 × 106 n = 10 × 106
20 42.6 (2) 42.6 (13)
40 10.3 (8) 10.2 (20)
60 4.5 (18) 4.4 (32)
Number M of binning grid points 80 100 120 140 160 3.1 2.9 2.2 1.8 1.6 (38) (67) (117) (172) (285) 2.7 2.5 2.2 1.7 1.5 (55) (91) (139) (204) (376)
180 1.5 (395) 1.4 (444)
Table 1: Simulation-based averages (based on 1000 replications) for M EE (and computational times in seconds in brackets), depending on the sample size n and the number M of grid points.
σ 0.003 0.006 0.009 0.012 0.015
MEE (mm) 2.6 (1.7,3.1) 2.9 (2.2,3.2) 3.9 (2.9,4.5) 4.2 (3.5,6.0) 4.3 (3.5,6.5)
Table 2: Simulation-based averages (based on 1000 replications) for M EE (and the limits for the 95% simulation interval in brackets), depending on σ for n = 1 × 106 and M = 100
from 0 to 10 mm. It should be noted that the value a = 0 corresponds to the null hypothesis (m is symmetric with respect to x0 = 0), while if a > 0 the null hypothesis H0 is not true. In the test 1000 bootstrap samples were generated for calculating type I error and 200 bootstrap samples for calculating the power under the alternative hypothesis. Remember that error type I provides the decision of rejecting the null hypothesis when it is true, and that it is fixed in advance by the significance level. Then, the contrast evidences a good test perfomance when error type I matches the significance level considered each time. On the other hand, the power of the test reflects its capacity to reject the null hypothesis when it is false. The type I error (expressed in %) for different significance levels and sample sizes is shown in Table 3 for different error variances σ. As can be observed from this table, the statistical test shows a good perfomance given that estimated type I errors are close to their nominal values. Similarly, Table 4 shows the rejection percentages for a not null values, that is the estimated power of the test. The power of the test allows studying the trade-off between point density, given by n, and the quality of the point cloud, given by σ . For instance, the power obtained for n = 1000 and σ=0.01 is similar to that for n = 10000 and σ=0.03. This result suggests that although we have a point cloud of low quality it is possible to get a good power if we have a sufficient number of points. Power curves are represented in Figure 3. As can be appreciated, the test produces satisfactory power curves, given that the probability of rejection increases whith the values of a. Also, it can be seen that for a standard deviation σ=0.01 (a usual value for TLS equipment), n = 10, 000 and a significance level of 5%, that probability of rejecting the null hypothesis is 99.5% when a = 2mm, that is, when the paraboloid separates a maximum of 2 mm from the symmetric. Note 10
200 1.3 (518) 1.2 (557
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Y (m)
−0.73
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.32 0.26
0.0
0.06
0.1
0.12
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=40
−0.5
0.0
0.5
X (m)
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Y (m)
−0.73
0.32 0.26
0.0
0.06
0.1
0.12
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=100
−0.5
0.0
0.5
X (m)
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Z (m)
−0.73
0.34 0.28
0.0
0.05
0.1
0.08
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=200
−0.5
0.0
0.5
X (m)
Figure 2: Cross-sections of the estimated curve (center) and zooms of different zones (right and left) for M =40, 100 and 200. The solid black line represents the theoretical, the red lines the estimation and the gray line each repetition boostrap.
also, from Figure 3, the increase of the power test when the number of points in the cloud point increases. That is, the more points in the simulation, the better the fit to the theoretical surface. 7. Results and discussion Applying the method for rebuilding the yacht hull from the point cloud, as described in Sections 3,4 and 5, we obtained the surface shown in Figure 4 .
11
n
1000
3000
5000
10000
σ 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05
1% 0.4 1.1 1.1 1.5 1.6 0.7 0.7 1.3 1.3 1.0 0.7 1.0 1.4 1.6 1.1 1.7 0.9 1.2 2.1 1.1
Significance level 5% 10% 15% 20% 3.5 7.1 11.9 17.3 4.8 8.6 13.7 18.6 5.2 10.4 15.3 20.2 5.1 10.7 15.7 21.3 4.4 10.1 15.6 20.7 5.0 11.1 15.8 20.7 5.5 12.4 16.7 21.6 5.8 12.2 17.6 22.1 5.1 12.2 18.5 24.3 5.8 11.3 17.3 23.0 4.8 9.1 14.3 19.8 5.1 9.3 15.4 21.7 5.1 9.4 16.1 21.8 5.6 10.9 16.1 22.9 4.3 9.6 15.5 21.6 4.6 9.9 15.3 21.0 4.3 11.0 15.9 19.7 5.2 10.9 16.6 21.5 5.4 11.4 16.8 22.0 3.8 9.9 14.7 19.5
Table 3: Estimated type I error (in percent) for the test.
a n
1000
10000
σ 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05
2 20.3 8.0 4.8 4.0 3.3 99.5 65.5 33.5 21.0 9.5
4 76.0 27.5 16.0 9.5 5.3 100 100 89.8 70.3 38.0
6 98.8 55.0 29.3 18.0 12.5 100 100 99.8 97.3 75.0
8 100 82.5 48.3 28.8 21.5 100 100 100 99.8 95.0
10 100 94.5 68.3 45.0 29.0 100 100 100 100 99.5
Table 4: Powers obtained for a significance level 5%.
From this surface we obtained 20 cross-sections perpendicular to the central fore-aft line, one for each 0.25 m of the yacht hull, as shown in Figure 5. In order to study the symmetry of the hull the statistical test proposed in Section 3 was applied. The p-value obtained is lower than 0.01, reaching statistical significance for the lack of symmetry. Figure 6 depicts several cross-
12
section from the bow of the yacht (solid line) and its symmetrical half (dotted line) that allows a visual comparison of the symmetry on both sides. As mentioned previously, the contrast established that there is no planar symmetry, but in Figure 6 it can be appreciated that differences between estimated and symmetric curves are quite small (just a few millimeters). As most boats are not 100% symmetrical, it is reasonable to consider certain tolerance depending on the dimensional requirements of the hull. If a tolerance L is allowed, then it is interesting to contrast the null hypothesis: H0 : ∃ x0 so that m(x0 + x, y) − m(x0 − x, y)
= 0∀(x, y)
(11)
For testing H0 a statistical test similar to the one exposed in Section 3 has been used but considering now the following statistical test. Z 2 ˆ 0 − x, y) − m(x T = min (min {|m(x ˆ 0 + x, y)| − L, 0}) dxdy (12) x0
For this kind of hull it is known that a symmetry not bigger than 5mm does not significantly affect the operation of the vessel. Considering L=4mm a p-value = 0.47 was obtained, resulting in a non significative contrast, so that the null hypothesis is accepted. Then the yacht can be considered symmetric inside the limits of tolerance. 8. Conclusions In this paper a statistical methodology is proposed to detect planar asymmetries from point clouds. First, it allows to recuperate an object surfaced from a point cloud using a non-parametric bivariate orthogonal regression based on kernel smoothers. Then, a global statistical test is described to verify that a real scanned object has planar symmetry. An important aspect of the proposed methodology is its capacity to control error type I and, at the same time, to reject the null hypothesis when it is false. The results of a simulation with a paraboloid carried out with a maximum of 1,000 points showed an increase of the test power with the scan density. Then, we conclude that symmetries can be better distinguished with a greater scan density. However, the simulation showed that detection of errors type I does not improve significantly after a specific density level. The proposed methodology can be easily extended to the study of planar symmetry establishing a tolerance, given that a perfect symmetry is not expected in many situations. This approach was applied to study planar symmetry on a vessel. The results obtained show that the vessel can be considered symmetric when a tolerance of 5 mm is considered. The test was applied to the whole vessel but it could be applied to parts of the scanned object in the same way.
13
Acknowledgements The authors gratefully acknowledge the financial support from projects FC15-GRUPIN14-033 of the Fundaci´on para el Fomento en Asturias de la Investigaci´on Cient´ıfica Aplicada y la Tecnolog´ıa (FICYT) (Spain) (FEDER support included), and MTM2011-23204 of the Spanish Ministry of Science and Innovation (FEDER support included). References [1] N.J. Mitra, M. Pauly, Symmetry for architectural design,. Adv. Archit Geom. (2008) 13-16 . [2] C. Ducottet, J. Daniere, M. Moine, J.P. Schon, M. Courbon, Localization of objects with circular symmetry in a noisy image using wavelet transforms and adapted correlation, Pattern Recogn. 27(3) (1994) 351-364. [3] H. Weyl, Symmetry, Princeton University Press, 1952. [4] I. Stewart, Why Beauty is Truth: A History of Symmetry, 2007. [5] E.H. Lockwood, R.H. MacMillan, Geometric Symmetry. Cambridge University Press, 1978. [6] W. Hodges, The geometry of music. Music and Mathematics: Pythagoras to Fractals, Oxford University Press, 2003.
From
[7] N.J. Mitra, M. Pauli, M. Wand, D. Ceylan, Symmetry in 3D geometry: Extraction and Applications. EUROGRAPHICS 2012. [8] G. Marola, IEEE T. Pattern Anal. 11(1) (1989) 104-108. [9] C. Sun, J. Sherrah, 3d symmetry detection using the extended gaussian image. IEEE T. Pattern Anal. 19 (2) (1997), 164-168. [10] Xuan Yang, Jihong Pei, Weixin Xie, Rotation Registration of Medical Images Based on Image Symmetry, Lect. Notes in Comput. Sc. 3644 (2005) 68-76. [11] B. Mills, F. Langbein, A. Marshall, R. Martin, Estimate of frequencies of geometric regularities for use in reverse engineering of simple mechanic components, Tech. Rep. GVG 2001-1. [12] S.J. Tate, G.E.M. Jared, Recognising symmetry in solid models, Comput. Aided D. 35(7) (2003) 673-692. ´ [13] J. Roca-Pardi˜ nas, F. L´opez -Alvarez, C. Ord´on ˜ez, A. Men´endez, Terrestrial laser scanning used to detect asymmetries in boat hulls, Opt. Eng. 51(1) (2012) 92-97.
14
[14] M. Lia, F.C. Langbein, R.R. Martin, Detecting design intent in approximate CAD models using symmetry, Comput. Aided D. 42(3) (2010) 183201. [15] J. Jiang, Z. Chen, K. He, A feature-based method of rapidly detecting global exact symmetries in CAD models, Comput. Aided D. 45(8-9)(2013) 1081-1094. [16] L.M. Rosario, W.A. Knight, Design for Assembly Analysis: Extraction of Geometric Features from a CAD System Data Base, CIRP Annals Manufacturing Technology 38(1) (1989), pp:13-16. [17] S.J. Tate , G.E.M. Jared, K.G. Swift, Detection of symmetry and primary axes in support of proactive design for assembly, SMA ’99 Proceedings of the fifth ACM symposium on Solid modeling and applications, 1999 pp. 151-158. [18] M.J. Atallah, On symmetry detection. IEEE T. Comput, C-34 (1985) 663666. [19] N. Kiryathi, Y. Gofman, Detecting symmtery in grey level images: the global optimiza- tion approach, Int. J. Comput Vision 29 (1998) 29-45. [20] Z. Xiao, Z. Hou, C. Miao, J. Wang, Using phase information for symmetry detection, Pattern Recogn. Lett. 26 (2005) 1985-1994. [21] B. Nicolai, H. Holzmann, M. Pawlak, Testing for Image Symmetries with Application to Confocal Microscopy. IEEE T. Inform. Theory, 55(4) (2009) 1841-1855. [22] M. Birke, H. Dette, K Stahljans, Testing symmetry of a nonparametric bivariate regression function. J. Nonparametr Stat. 23(2) (2011), 547-565. [23] M. Chang, S.C. Park, Reverse engineering of a symmetric object.Comput & Ind. Eng., 55(2) (2008) 301-320. [24] T. Schulz, H. Ingensand, , Terrestrial Laser Scanning - Investigation and Applications for High Precision Scanning.FIG Working Week, Athens, Greece, 2004, pp. 1-15. [25] F. de As´ıs L´opez, C. Ord´on ˜ez, J. Roca - Pardi˜ nas, S. Garc´ıa - Cort´es, Point cloud comparison under uncertainty. Application to beam bridge measurement with terrrestrial laser scanning, Measurement 51 (2014) 259-264. [26] W.E. Deming, Statistical adjustment of data. Wiley, NY, 1985 [27] T. Hastie,W. Stuetzle, Principal Curves, J. Am. Stat. Assoc. 84 (1989) 502-516.
15
[28] F. de As´ıs - L´opez, C. Ord´on ˜ez, J. Roca - Pardi˜ nas, S. Garc´ıa - Cort´es, A statistical method for geometry inspection from point clouds, App. Math. Comput. 242 (2014) 562-568. [29] R.L. Eubank, Spline smoothing and nonparametric regression, New York, 1988. [30] W. Hardle, Applied nonparametric regression, Cambridge University Press, Cambridge, 1990. [31] P. Sarda, P. Vieu, Smoothing kernel regression, In: Smoothing and Regression: Approaches, Computation, and Application, New York: Wiley, 1999. [32] D. Ruppert D, M.P. Wand, R.J. Carroll, Semiparametric regression, Cambridge University Press, Cambridge, 2003. [33] S.N. Wood, Thin plate regression splines, J. R. Stat. Soc. B 65 (2003) 95-114. [34] B.W. Silverman, Density Estimation For Statistics And Data Analysis, Chapman and Hall, London, 1986. [35] D.W. Scott, S.J. Sheather, Kernel density estimation with binned data, Commun Stat - Theor M 14 (1985) 1353-1359. [36] P. Hall , M.P. Wand, On the accuracy of binned kernel density estimators, J. Multivariate Anal. 56(2)(1996) 165-184 [37] J. Fan, J.S. Marron, Fast implementation of nonparametric curve estimators, J. Comput. Graph. Stat. 3 (1994) 35-56. [38] M.P Wand, Fast Computation of Multivariate Kernel Estimators, J. Comp. Graph. Stat. 3 (1994) 433-45. [39] B. Efron Bootstrap methods: another look at the jacknife, Ann. Stat. 7 (1979) 1-26. [40] C.F.J. Wu, Jackknife, Bootstrap and Other Resampling Methods in Regression Analysis, Ann. Stat 14(4 )(1986) 1261-1295. [41] R.Y. Liu, Bootstrap Procedures under some Non-I.I.D. Models. The Annals of Statistics 16(4) (1988) 1696-1708. [42] E. Mammen, Bootstrap and Wild Bootstrap for High Dimensional Linear Models, Ann. Stat. 21(1) (1993) 255-285.
16
n= 1000
●
80
σ=0.01 σ=0.02 σ=0.03 σ=0.04 σ=0.05
60
60
80
100
100
n= 1000
● ●
40
40
●
●
●
●
20
20
● ● ●
● ●
0
0
●
4
6
8
10
0
2
4
6
a (mm)
a (mm)
n= 3000
n= 3000
8
100
2
100
0
10
●
●
80
80
● ●
60
60
●
●
40
40
●
●
20
20
● ●
●
0
0
●
0
2
4
6
8
10
0
2
4
a (mm)
6
8
10
a (mm)
100
n= 5000
100
n= 5000 ●
● ●
80
80
●
●
60
60
●
●
40
40
●
20
20
●
●
●
0
0
●
0
2
4
6
8
10
0
2
4
a (mm)
6
8
10
●
●
8
10
a (mm)
●
100
n= 10000
100
n= 10000 ●
●
80
80
●
60
60
●
40
40
●
●
20
20
●
●
0
0
●
0
2
4
6
8
10
0
a (mm)
2
4
6
a (mm)
Figure 3: Percentage rejection for the test for increasing a nominal levels of 5% (left) and 10%(right),and sample sizes of n = 1000 (top row), n = 3000 (second row), n = 5000 (third row) and n = 10000 (bottom row).
17
Figure 4: Yatch hull surface estimated from the point cloud.
0.5
1.0
1.5
Z (m)
2.0
2.5
3.0
Cross−sections
−2
−1
0
1
2
Y (m)
Figure 5: Cross-sections considered in the analysis of asymmetries.
18
1.50
1.60
2.0
1.70
1.80
X= −2.605
1.15
1.20
1.25
0.70
1.5 0.5
0.55
0.60
0.65
1.0
Z (m)
1.10
−1.5
−1.0
−0.5
0.0
0.5
1.0
1.5
0.00
0.05
0.10
0.15
0.20
1.00
1.05
1.10
1.15
1.20
Y (m)
0.86
0.7
0.90
0.8
0.94
0.9
Z (m)
1.0
1.1
1.06
1.10
1.2
1.14
X= 0.553
−1.5
−1.0
−0.5
0.0
0.5
1.0
1.5
0.50
0.55
0.60
0.65
Y (m)
Figure 6: Cross-sections of the estimated surface (left) and zooms (rigth), representing the estimation (solid line) and the symmetric (dotted line).
19
S0263-2241(16)30746-1 http://dx.doi.org/10.1016/j.measurement.2016.12.048 MEASUR 4510
To appear in:
Measurement
Received Date: Accepted Date:
27 October 2015 26 December 2016
Please cite this article as: J. Roca-Pardiñas, C. Ordó ñez, C. Cabo, A. Menéndez-Díaz, Assessing planar asymmetries in shipbuilding from point clouds, Measurement (2016), doi: http://dx.doi.org/10.1016/j.measurement.2016.12.048
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Assessing planar asymmetries in shipbuilding from point clouds Javier Roca-Pardi˜ nasa , Celestino Ord´ on ˜ezb , Carlos Cabob , Agus´ın Men´endez-D´ıazc a Department
of Statistics and Operations Research, University of Vigo, 36208 Vigo, Spain of Mining Exploitation and Prospecting, University of Oviedo, 33600 Oviedo, Spain c Department of Construction and Manufacturing Engineering, University of Oviedo, 33203 Gij´ on, Spain Author for the correspondence: Celestino Ord´ on ˜ez, E-mail address: [email protected], Telephone number: +34985458027, Fax number: +34985458188 b Department
Abstract This paper introduces a statistical methodology to identify planar asymmetries in surfaces. An algorithm based on bivariate kernel smoothers is used to estimate non-parametrically the surface from a point cloud. Also bootstrap-based procedures are proposed for testing a null hypothesis establishing that the surface has planar symmetry. After checking the validity of the proposed methodology with simulated data, it was applied to study symmetries in a boat hull. Keywords: symmetry, bootstrap, kernel smoothing, laser scanner, point cloud
1. Introduction Symmetry is a property present in nature but also in artificial objects, such as architectural models and industrial products [1]. It plays an important role in human visual perception and aesthetics [2] but also in the economy and efficiency of structures and animate or inanimate organisms [3]. From a mathematical point of view, a symmetry is a set of transformations applied to an object that preserve the properties of the structure [4]. Translations, rotations and reflections are typical geometric transformations corresponding to symmetries. Although it is widely associated to geometric transformations [5], symmetry concerns other disciplines, such as music [6]. The importance of symmetry has originated a number of studies concerning this topic throughout time. The authors in [7] describe a variety of applications of symmetry in digital geometry processing and shape analysis. They also analyse the strengths and limitations of some algorithms for symmetry detection. Some examples of symmetry detection in images can be seen in [8–10]. In the architectural and industrial fields, symmetry is also a prime characteristic of
Preprint submitted to Measurement
December 31, 2016
some objects that can be used in reverse engineering. In [11] a number of aspects of the geometry of mechanical parts, such as geometric restrictions, slots, holes or symmetries, which should be considered to improve boundary representation models, are described. The authors in [12] analyse the cognitive process involved in human symmetry detection to develop a computational method to detect symmetries. In [13] the authors developed a methodology to detect symmetries on a vessel from a point cloud measured with a terrestrial laser scanner that has some similarites with that presented in this paper. However, in this work we extend the analysis to the whole surface of the vessel, without restricting it to cross-sections. In addition, it has other significant differences related to the bootstrapping method used, the application of binning techniques to reduce computational time, and the analysis of the effect of the standard deviation in the error term. Symmetries are also fundamental in the design of objects. Thus, symmetries were used to detect the design intent in approximate CAD models [14]. In [15] symmetries in CAD models are detected using Boolean combinations of the object features. Symmetry is also a basic concept when parts are assembled in one piece [16, 17]. This paper addresses the problem of contrasting the symmetry of an object surface function. There are several methods in the literature that tackle the symmetry of objects or functions. Some of these methods are centered on determining a symmetry axis [18, 19] or [20], among others. In most of these papers it is admitted that the object surface function is observed without stochastic errors. Other authors analyse if an image or data are symmetric with respect to a known or unknown axis. Recently, tests for reflection symmetry and rotational invariance of two-dimensional functions or images based on Zernike polynomials were proposed [21]. A symmetry test for a regression function with a bivariate predictor based on the L2 distance between the original function and its reflection was proposed in [22]. In this work we propose a slightly different vision of the regression itself, where errors are considered in three dimensions, not just in the Z component. Our research aims to provide a methodology to detect planar symmetries in objects from point clouds (obtained from laser scanners, close-range photogrammetry or other measurement techniques). The method is based on fitting a surface to the point cloud using a non-parametric orthogonal regression method, sectioning the surface into two halves and comparing these using confidence intervals for the difference between each half obtained by bootstrapping. Thus, it can be determined whether each section can be considered symmetric in relation to a hypothetical plane of symmetry. The hypothesis test must be able to accept the null hypothesis when it is true while it has enough power to reject that hypothesis when it is false. This work also includes linear binning as an efficient technique to manage point clouds of several million points, since these massive point clouds are quite frequent nowadays. The rest of the paper is structured as follows: Section 2 introduces the problem of detecting planar symmetry in shipbuilding. The methodology used to adjust a surface to the point cloud is explained in Section 3. Section 4 shows 2
the methododolgy developed to detect symmetries in the surface. In Section 5 a simulation study to validate the method is carried out. In Section 6 the proposed methodology is used to detect symmetries in a yacht hull from a point cloud obtained using a terrestrial laser scanner. Finally the most relevant conclusions of the research are collected. 2. Testing symmetries in a boat: problem statement The methodology proposed in this paper was used to determine planar symmetry in a yacht measuring 12 m long, with a beam of 4.5 m and a draft of 2.4 m (see Figure 1). Symmetry along a vertical plane containing the fore-aft axis is important for correct sailing behavior, as asymmetry can cause a shift in the volume of displaced water and in a different thrust on the port and starboard sides of the hull, causing the vessel to list. Another possible drawback of asymmetries is a lateral thrust that causes the vessel to veer from a straight-line course. The relevance of these effects depends on the planned use of the vessel, its characteristics and dimensional requirements. The use of point clouds obtained from laser scanners or close-range photogrammetry for reverse engineering problems, such as the detecion of symmetries, is currently common in industry applications apart from shipbuilding. [23].
Figure 1: Images of the yacht hull.
The yacht hull was scanned using a Leica Cyrax 2500 time-of-flight terrestrial laser scanner (TLS) able to record up to 1000 points per second with a nominal accuracy of ±6mm(1σ) in position and ±4mm(1σ) in distance. The measurent range is 1.5 to 50 m. The scan density (minimum point-to-point measurement spacing) can be modified, being its maximum value less than 1 mm. Terrestrial laser scanners determine the position of the points calculating the distances and the horizontal and vertical angles between the center of the scanner and the points. As it name indicates, time-of-flight scanners measured distances to the points determining the time needed for a laser-beam to travel from the laser emissor to the object and back. Most TLS devices use a rotating mirror located inside the instrument to reflect the emited laser-beam in different 3
directions, allowing to emit a great number of laser-beams in a very short time. Consequently, thousand of points can be measured in a few seconds. This type of laser scanners are able to measure distances of up to several hundred meters under ideal circumstances. For its part, phase-based scanners operate in a different way: they estimate the phase difference between transmitted and reflected waves. These scanners range is shorter than that of the time-of-flight scanners, but, on the contrary, have greater precision and are also much much faster [24]. The scanner was placed in 12 different locations to measure the entire boat. All the scans were registered in the same coordinate system using magnetic reflective targets that were automatically recognized. Using a scan density about one point per square centimeter, aproximatelly 1,600,000 points were necessary to cover the entire yatch hull. Each of these points {P1• , . . . , Pn• } with Pi• = (Xi• , Yi• , Zi• ), being n the number of observations, can be expressed as Pi• = Pi + i
(1)
where ! X YPi Z= (xi , yi , zi ) represents the point of the object surface and i = εi , εi , εi the three dimensional measurement error on the ith point [25]. The only restriction for the distribution of the error term is that it has to have zero mean. Being (x, y, z) the spatial coordinates of each point on the object surface, coordinate z is obtained adjusting a smooth surface z = m (x, y) to the data. Given the point cloud, the first step in our method is to reconstruct the ship hull obtained through the estimate m ˆ of that surface. Once the boat hull surface has been reconstructed, the next step is to detect asymmetries in this surface. Specifically, we are interested in determining if the ”left” and ”right” parts of the surface are symmetric in relation to some unknown plane given by x = x0 . Testing planar symmetry can considered as a hypothesis test where the null hypothesis is H0 : ∃ x0 so that m(x0 + x, y) = m(x0 − x, y) ∀(x, y)
(2)
This equation is valid for symmetries with respect to a vertical plane, and it should be modified to tackle non-vertical symmetry planes. Planar symmetry is present in many man-made objects; particularly, planar symmetry is important in ship construction where, except in rare cases, a ship’s hull is symmetrical about a fore and aft plane. Symmetry affects speed navigation, righting position, stability and maneuverability. 3. Non-parametric estimation algorithm. This section describes the algorithm used to fit the unknown function m to the point cloud {P1• , . . . , Pn• }. For any given point (x0 , y0 ) , a smoothed version of the orthogonal regression (or Deming regression [26]) is proposed. In this type of regression the estimated surface constitutes the best fit to the set of the
4
observed sample points {P1• , . . . , Pn• } by minimizing the perpendicular distances from the points to the surface. The proposed estimators are based on the fact that surface (x, y, m(x, y)) can be approximated by the plane ax + by + cz + d = 0 in values (x, y) near (x0 , y0 ). Thus, the estimated coefficients (ˆ a, ˆb, cˆ) correspond to the third local principal component [27, 28]. The estimate of the true surface m(x0 , y0 ) at each point (x0 , y0 ) is obtained as follows: First, compute the sample ˆ obtained from the point cloud (P • , . . . , Pn• ) weighted 3 × 3 covariance matrix Σ 1 by (W1 , . . . , Wn ) with
Wi =
Wh (Xi•
−
x0 , Yi•
(Xi• − x0 )2 + (Yi• − y0 )2 − y0 ) = exp − h
(3)
ˆ as and obtain (ˆ a, ˆb, cˆ) as the third principal component of the PCA using Σ the input covariance matrix. Then, calculate a plane normal to this vector and passing through (X¯• , Y¯• , Z¯• ) a ˆ(x − X¯• ) + ˆb(y − Y¯• ) + cˆ(z − Z¯• ) = 0 Finally, evaluating that plane in (x0 , y0 ) and solving in z, an estimation of m(x0 , y0 ) is obtained m(x ˆ 0 , y0 ) = Z¯• −
a ˆ(x0 − X¯• ) + ˆb(y0 − Y¯• ) cˆ
(4)
Moreover the vector
a ˆ, ˆb, cˆ
u ˆ(x0 , y0 ) = p a ˆ2 + ˆb2 + cˆ2
(5)
represents the estimated orthogonal vector to the surface z = m(x, y) at the point (x0 , y0 ). Note that the kernel weight Wi in (10) depends on the euclidian distance between (x0 , y0 ) and (Xi• , Yi• ) and, in addition, contains a smoothing parameter h. See [29–31] for information about kernel smoothing. Instead of kernel smoothers, other methods such as penalized and plate regression splines could have been used[32, 33], respectively. 3.1. Bandwidth selection According to the definition of Wh , the influence on the estimate m(x ˆ 0 , y0 ) depends on the distance of (Xi• , Yi• ) to (x0 , y0 ). The smoothing parameter h in (10), generally called the bandwidth, is a critical factor that controls the extension of that influence. When h is small the estimate m(x ˆ 0 , y0 ) tends to adjust too much to the data while it produces a very high variance. On the contrary, when h is large the estimated surface tends to approximate a plane, so the the adjustment to the object is quite poor, being the variance small. 5
Determining an adequate bandwidth is a classical problem in kernel smoothing. Sometimes it is given manually according to the knowledge of the problem. Other times it is obtained automatically by means of cross-validation. In the latter case the bandwidth h is selected by minimizing
CV (h) =
n X i=1
min x,y
q
(Xi• − x)2 + (Yi• − y)2 + (Zi• − m ˆ −i (x, y))2
(6)
where m ˆ −i (x, y) indicates the fit of m(x, y) leaving out the i-sample data (Xi• , Yi• , Zi• ). 4. Testing symmetry In the previous section the mechanism to allow the nonparametric reconstruction of the smooth surface z = m(x, y) has been proposed. In this section we discuss the problem of testing for symmetry of the function m. For testing the null hypothesis in (2) we propose to use the test statistic given by Z 2 T = min (m(x ˆ 0 − x, y) − m(x ˆ 0 + x, y)) dxdy) (7) x0
Note that if H0 is verified, the value of T value should be close to zero. For a given significance α, the decision rule based on T consists of rejecting the null hypothesis if T > T α , where T α is the (1 − α)-percentile of T under the null hypothesis. Nevertheless, within a non-parametric regression context, the theory for determining such percentiles is not closed, and resampling methods such as bootstrap introduced by Efron [39] can be applied instead. Here, we propose an ”orthogonal” version of the traditional wild bootstrap [40, 41]. The steps of the proposed procedure are as follows: ˆ Step 1: Given the original sample P1• , . . . , Pn• obtain the estimated surface m as explained above. Then obtain the symmetric estimated surface m ˆ S under the null hypothesis of symmetry. In particular, for a given point (x, y) the estimated symmetric surface m ˆ S (x, y) can be obtained as the mean of the ”left” and ”right” estimated curves. Explicitly, the expression of m ˆ S (x, y) is given by m ˆ S (x, y) = 0.5m(x, ˆ y) + 0.5m(2ˆ ˆ x0 − x, y) where x ˆ0 is the value of x0 that minimizes the expression in (7). Step 2: Obtain the orthogonal distances to the symmetric surface q 2 2 2 Di = (Xi• − x0i ) + (Yi• − yi0 ) + (Zi• − m ˆ S (x0i , yi0 )) ! where x0i , yi0 is obtained as !
q 2 ˆ S (x, y)) x0i , yi0 = arg min (Xi• − x)2 + (Yi• − y)2 + (Zi• − m x,y
6
(8)
! ! Note that Pi0 = x0i , yi0 , m ˆ S x0i , yi0 is the point on the estimated symmetric surface m ˆ S whose distance to the measured point (Xi• , Yi• , Zi• ) is minimal. Moreover, it is easy to show that the orthogonal vectors to the surface at points (x, y) and (2x0 − x, y) are related, under the null hypothesis, by means of a 90 degree rotation about y axis us (x, y) = T · u(2x0 − x, y) being
0 T= 0 1
0 1 0
-1 0 0
Bearing this in mind, for each (x0i , yi0 ) the estimate of the true u(x0i , yi0 ) under H0 can be obtained as follows u ˆS (x0i , yi0 ) = u ˆ(x0i , yi0 ) + T · u ˆ(2ˆ x0 − x0i , yi0 ) Step 3: For b = 1, . . . , B (e.g. P1•,b , . . . , Pn•,b in accordance with
B=1000) generate the bootstrap samples
Pi•,b = Pi0 +
u ˆS (x0i , yi0 ) · Ri ||ˆ uS (x0i , yi0 )||
where || · || represents the euclidean norm, and Ri being ( √ Di · (1−2√5) with probability p = Ri = Di · (1+2 5) with probability p =
√ 5+ 5 10 √ 5− 5 10
that is, Ri follows a distribution that resembles the Mammen distribution [42]. Then, compute T •b with the generated bootstrap in the same way that T had been obtained with the original sample. Step 4: Finally, the test rule consists of rejecting the null hypothesis if T > Tˆα , where Tˆα is the empirical (1 − α)-percentile of values T •1 , . . . , T •B obtained before. 5. Computational aspects Non-parametric regression modelling can be a very time consuming process when, as in our case, the number of data n is very large. Accordingly, some thecniques aimed to reduce this computing time, such a binning, has been proposed. Information regarding bining thecniques for density estimation can be found in [34–37]. Specifically, in this work we use the binning algorithm for multivariate kernel estimation proposed by [38]. ˜ • < ... < X ˜ • , Y˜ • < ... < Y˜ • and Z˜ • < ... < Z˜ • being M equidistant Let X 1 1 1 M M M points along the X, Y and Z directions respectively. The linear binning assigns 7
˜ • , Y˜ • , Z˜ • a weight that depends on the relative distances to each grid point X j k l between the observations and the eight closest grid points according to the expression
˜ •i W jkl
• • ˜ • • ˜• ˜ • − Y − Z Xi − X Y Z j i i k l 1 − 1 − = 1 − δX δY δZ
+
+
(9) +
with x+ = max {0, x} and δX , δY and δZ the distances between two neighbouring knots on X, Y , Z, respectively. The binning version of the orginal estimate m(x ˆ 0 , y0 ), which expression is given in (4), is obtained as previously explained, ˆ from the M 3 grid points P˜ • = but computing now the covariance matrix Σ jkl • ˜ • ˜• ˜ X , Y , Z with associated weights j
k
l
• ˜ jkl (x0 , y0 ) W
˜ • − x0 )2 + (X ˜ • − y0 )2 (X j k = exp − h (
)
• ˜ jkl W
(10)
In the same way, the bandwidth selection procedure the CV criterion in (6) can be approximated by
CV (h) ≈
M X M X M X
• ˜ jkl W
j=1 k=1 l=1
· min x,y
r
2 ˜ • − x)2 + (Y˜ • − y)2 + Z˜ • − m (X ˆ −(j,k,l) (x, y) j k k
where m−(j,k,l) (x, y) indicates now the binning approximation fit at (x, y) leaving out the (j, k, l)-th binning data point. The number of grid points selected reflects a balance between the degree of adjustment and the computer time. The distribution of the covariates and the regularity of the surface m, also have influence on the suitable number of grid points. Particularly, we propose to choose a number of nodes M0 which meets that not significant differences exists between the estimation obtained with this value and that obtained with a greater number of nodes (M > M0 ). 6. Simulation Study In order to assess the behavior of our proposed algorithms it was applied to simulated data. The procedure is performed in four steps, which are described as follows. 1. The first step consists in the definition of a theoretical mathematical surface Z = m(X, Y ). In our case we have defined a polynomial surface of degree 2 in (X, Y ) Z = m(X, Y ) = a · X 2 + Y 2 where a is a real constant, being X and Y independent variables with uniform distribution on [−1, 1], 8
2. In the second step each point measured with the scanning is simulated according to ! Y Z (Xi• , Yi• , Zi• ) = (Xi , Yi , Zi ) + εX i , εi , εi ; i = 1, ..., n ! Y Z where the error term εX is supposed to be uniformly distributed i , εi , εi on a sphere of a random radius distributed in accordance with a N (0, σ). Several values of σ (from σ = 0.01 to σ = 0.05) were considered in order to obtain different error variances. 3. The proposed estimation algorithm is applied to the point cloud, obtaining an approximation m(x, ˆ y) to the theoretical surface in each point (x, y). To evaluate the behaviour of the estimates we consider a grid of 400 equidistant points in the interval [−1, 1], given by 2(i − 1) 2(j − 1) (ui , vj ) = −1 + , −1 + f or 1 ≤ i, j ≤ 20 19 19 At each grid point (ui , vj ) the nearest point in the true surface is computed q ! 0 0 2 ˆ y)) ui , vj = arg min (ui − x)2 + (vj − y)2 + (m(ui , vj ) − m(x, x,y
Finally, we consider the mean euclidean error (MEE) given by M EE = 400−1
20 q 20 X X ! 2 ! ! 2 2 (ui − u0i ) + vj − vj0 + m(ui , vj ) − m ˆ u0i , vj0 i=1 j=1
Table 1 summarizes numerical average results of the considered error (MEE) over 1000 replicated samples for two different sample sizes (n = 1 × 106 and n = 10 × 106 ) as well as the computational time in brackets (using a PC with a processor 2.53 GHz and RAM 4 GB.), for an increasing number of grid points M . In this table an error variance σ 2 = 0.0062 has been considered. Moreover, the MEE and the limits for the 95% simulation interval are shown in 2 for different σ, considering n = 1 × 106 and M = 100. From Table 1 we can see that the choice of M is a trade-off between the computational time and the error when estimating the true surface. As M increases, the M EEs decrease, but the computing time may increase substantially. On the other hand, the error reduction is very small from M = 100 on. Moreover, it can be appreciated that errors for n = 1 × 106 are similar to those obtained for n = 10 × 106 . That is, there is not much gain in increasing the scan density after a certain level. Finally, as can be seen in Table 2, M EE increases with σ, as was expected. Figure 2 depicts a section of the estimated surface (center) and zoom of different regions (right and left). The continuous black line represents the theoretical surface, the red line the estimation and the gray line each of the bootstrap repetitions for different values of M (M =40, 100 y 200). Finally, we explore the validity of the test considering the null hypothesis H0 of symmetry given in (2). Different values were considered for a, ranging 9
n = 1 × 106 n = 10 × 106
20 42.6 (2) 42.6 (13)
40 10.3 (8) 10.2 (20)
60 4.5 (18) 4.4 (32)
Number M of binning grid points 80 100 120 140 160 3.1 2.9 2.2 1.8 1.6 (38) (67) (117) (172) (285) 2.7 2.5 2.2 1.7 1.5 (55) (91) (139) (204) (376)
180 1.5 (395) 1.4 (444)
Table 1: Simulation-based averages (based on 1000 replications) for M EE (and computational times in seconds in brackets), depending on the sample size n and the number M of grid points.
σ 0.003 0.006 0.009 0.012 0.015
MEE (mm) 2.6 (1.7,3.1) 2.9 (2.2,3.2) 3.9 (2.9,4.5) 4.2 (3.5,6.0) 4.3 (3.5,6.5)
Table 2: Simulation-based averages (based on 1000 replications) for M EE (and the limits for the 95% simulation interval in brackets), depending on σ for n = 1 × 106 and M = 100
from 0 to 10 mm. It should be noted that the value a = 0 corresponds to the null hypothesis (m is symmetric with respect to x0 = 0), while if a > 0 the null hypothesis H0 is not true. In the test 1000 bootstrap samples were generated for calculating type I error and 200 bootstrap samples for calculating the power under the alternative hypothesis. Remember that error type I provides the decision of rejecting the null hypothesis when it is true, and that it is fixed in advance by the significance level. Then, the contrast evidences a good test perfomance when error type I matches the significance level considered each time. On the other hand, the power of the test reflects its capacity to reject the null hypothesis when it is false. The type I error (expressed in %) for different significance levels and sample sizes is shown in Table 3 for different error variances σ. As can be observed from this table, the statistical test shows a good perfomance given that estimated type I errors are close to their nominal values. Similarly, Table 4 shows the rejection percentages for a not null values, that is the estimated power of the test. The power of the test allows studying the trade-off between point density, given by n, and the quality of the point cloud, given by σ . For instance, the power obtained for n = 1000 and σ=0.01 is similar to that for n = 10000 and σ=0.03. This result suggests that although we have a point cloud of low quality it is possible to get a good power if we have a sufficient number of points. Power curves are represented in Figure 3. As can be appreciated, the test produces satisfactory power curves, given that the probability of rejection increases whith the values of a. Also, it can be seen that for a standard deviation σ=0.01 (a usual value for TLS equipment), n = 10, 000 and a significance level of 5%, that probability of rejecting the null hypothesis is 99.5% when a = 2mm, that is, when the paraboloid separates a maximum of 2 mm from the symmetric. Note 10
200 1.3 (518) 1.2 (557
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Y (m)
−0.73
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.70
0.71
0.72
0.73
0.74
0.75
0.50
0.52
0.54
0.56
0.58
0.60
0.32 0.26
0.0
0.06
0.1
0.12
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=40
−0.5
0.0
0.5
X (m)
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Y (m)
−0.73
0.32 0.26
0.0
0.06
0.1
0.12
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=100
−0.5
0.0
0.5
X (m)
−0.30
−0.26
−0.22
0.56 0.50
0.4
−0.71
Z (m)
−0.73
0.34 0.28
0.0
0.05
0.1
0.08
0.2
−0.75
0.3
0.50
0.5
0.56
0.6
M=200
−0.5
0.0
0.5
X (m)
Figure 2: Cross-sections of the estimated curve (center) and zooms of different zones (right and left) for M =40, 100 and 200. The solid black line represents the theoretical, the red lines the estimation and the gray line each repetition boostrap.
also, from Figure 3, the increase of the power test when the number of points in the cloud point increases. That is, the more points in the simulation, the better the fit to the theoretical surface. 7. Results and discussion Applying the method for rebuilding the yacht hull from the point cloud, as described in Sections 3,4 and 5, we obtained the surface shown in Figure 4 .
11
n
1000
3000
5000
10000
σ 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05
1% 0.4 1.1 1.1 1.5 1.6 0.7 0.7 1.3 1.3 1.0 0.7 1.0 1.4 1.6 1.1 1.7 0.9 1.2 2.1 1.1
Significance level 5% 10% 15% 20% 3.5 7.1 11.9 17.3 4.8 8.6 13.7 18.6 5.2 10.4 15.3 20.2 5.1 10.7 15.7 21.3 4.4 10.1 15.6 20.7 5.0 11.1 15.8 20.7 5.5 12.4 16.7 21.6 5.8 12.2 17.6 22.1 5.1 12.2 18.5 24.3 5.8 11.3 17.3 23.0 4.8 9.1 14.3 19.8 5.1 9.3 15.4 21.7 5.1 9.4 16.1 21.8 5.6 10.9 16.1 22.9 4.3 9.6 15.5 21.6 4.6 9.9 15.3 21.0 4.3 11.0 15.9 19.7 5.2 10.9 16.6 21.5 5.4 11.4 16.8 22.0 3.8 9.9 14.7 19.5
Table 3: Estimated type I error (in percent) for the test.
a n
1000
10000
σ 0.01 0.02 0.03 0.04 0.05 0.01 0.02 0.03 0.04 0.05
2 20.3 8.0 4.8 4.0 3.3 99.5 65.5 33.5 21.0 9.5
4 76.0 27.5 16.0 9.5 5.3 100 100 89.8 70.3 38.0
6 98.8 55.0 29.3 18.0 12.5 100 100 99.8 97.3 75.0
8 100 82.5 48.3 28.8 21.5 100 100 100 99.8 95.0
10 100 94.5 68.3 45.0 29.0 100 100 100 100 99.5
Table 4: Powers obtained for a significance level 5%.
From this surface we obtained 20 cross-sections perpendicular to the central fore-aft line, one for each 0.25 m of the yacht hull, as shown in Figure 5. In order to study the symmetry of the hull the statistical test proposed in Section 3 was applied. The p-value obtained is lower than 0.01, reaching statistical significance for the lack of symmetry. Figure 6 depicts several cross-
12
section from the bow of the yacht (solid line) and its symmetrical half (dotted line) that allows a visual comparison of the symmetry on both sides. As mentioned previously, the contrast established that there is no planar symmetry, but in Figure 6 it can be appreciated that differences between estimated and symmetric curves are quite small (just a few millimeters). As most boats are not 100% symmetrical, it is reasonable to consider certain tolerance depending on the dimensional requirements of the hull. If a tolerance L is allowed, then it is interesting to contrast the null hypothesis: H0 : ∃ x0 so that m(x0 + x, y) − m(x0 − x, y)
= 0∀(x, y)
(11)
For testing H0 a statistical test similar to the one exposed in Section 3 has been used but considering now the following statistical test. Z 2 ˆ 0 − x, y) − m(x T = min (min {|m(x ˆ 0 + x, y)| − L, 0}) dxdy (12) x0
For this kind of hull it is known that a symmetry not bigger than 5mm does not significantly affect the operation of the vessel. Considering L=4mm a p-value = 0.47 was obtained, resulting in a non significative contrast, so that the null hypothesis is accepted. Then the yacht can be considered symmetric inside the limits of tolerance. 8. Conclusions In this paper a statistical methodology is proposed to detect planar asymmetries from point clouds. First, it allows to recuperate an object surfaced from a point cloud using a non-parametric bivariate orthogonal regression based on kernel smoothers. Then, a global statistical test is described to verify that a real scanned object has planar symmetry. An important aspect of the proposed methodology is its capacity to control error type I and, at the same time, to reject the null hypothesis when it is false. The results of a simulation with a paraboloid carried out with a maximum of 1,000 points showed an increase of the test power with the scan density. Then, we conclude that symmetries can be better distinguished with a greater scan density. However, the simulation showed that detection of errors type I does not improve significantly after a specific density level. The proposed methodology can be easily extended to the study of planar symmetry establishing a tolerance, given that a perfect symmetry is not expected in many situations. This approach was applied to study planar symmetry on a vessel. The results obtained show that the vessel can be considered symmetric when a tolerance of 5 mm is considered. The test was applied to the whole vessel but it could be applied to parts of the scanned object in the same way.
13
Acknowledgements The authors gratefully acknowledge the financial support from projects FC15-GRUPIN14-033 of the Fundaci´on para el Fomento en Asturias de la Investigaci´on Cient´ıfica Aplicada y la Tecnolog´ıa (FICYT) (Spain) (FEDER support included), and MTM2011-23204 of the Spanish Ministry of Science and Innovation (FEDER support included). References [1] N.J. Mitra, M. Pauly, Symmetry for architectural design,. Adv. Archit Geom. (2008) 13-16 . [2] C. Ducottet, J. Daniere, M. Moine, J.P. Schon, M. Courbon, Localization of objects with circular symmetry in a noisy image using wavelet transforms and adapted correlation, Pattern Recogn. 27(3) (1994) 351-364. [3] H. Weyl, Symmetry, Princeton University Press, 1952. [4] I. Stewart, Why Beauty is Truth: A History of Symmetry, 2007. [5] E.H. Lockwood, R.H. MacMillan, Geometric Symmetry. Cambridge University Press, 1978. [6] W. Hodges, The geometry of music. Music and Mathematics: Pythagoras to Fractals, Oxford University Press, 2003.
From
[7] N.J. Mitra, M. Pauli, M. Wand, D. Ceylan, Symmetry in 3D geometry: Extraction and Applications. EUROGRAPHICS 2012. [8] G. Marola, IEEE T. Pattern Anal. 11(1) (1989) 104-108. [9] C. Sun, J. Sherrah, 3d symmetry detection using the extended gaussian image. IEEE T. Pattern Anal. 19 (2) (1997), 164-168. [10] Xuan Yang, Jihong Pei, Weixin Xie, Rotation Registration of Medical Images Based on Image Symmetry, Lect. Notes in Comput. Sc. 3644 (2005) 68-76. [11] B. Mills, F. Langbein, A. Marshall, R. Martin, Estimate of frequencies of geometric regularities for use in reverse engineering of simple mechanic components, Tech. Rep. GVG 2001-1. [12] S.J. Tate, G.E.M. Jared, Recognising symmetry in solid models, Comput. Aided D. 35(7) (2003) 673-692. ´ [13] J. Roca-Pardi˜ nas, F. L´opez -Alvarez, C. Ord´on ˜ez, A. Men´endez, Terrestrial laser scanning used to detect asymmetries in boat hulls, Opt. Eng. 51(1) (2012) 92-97.
14
[14] M. Lia, F.C. Langbein, R.R. Martin, Detecting design intent in approximate CAD models using symmetry, Comput. Aided D. 42(3) (2010) 183201. [15] J. Jiang, Z. Chen, K. He, A feature-based method of rapidly detecting global exact symmetries in CAD models, Comput. Aided D. 45(8-9)(2013) 1081-1094. [16] L.M. Rosario, W.A. Knight, Design for Assembly Analysis: Extraction of Geometric Features from a CAD System Data Base, CIRP Annals Manufacturing Technology 38(1) (1989), pp:13-16. [17] S.J. Tate , G.E.M. Jared, K.G. Swift, Detection of symmetry and primary axes in support of proactive design for assembly, SMA ’99 Proceedings of the fifth ACM symposium on Solid modeling and applications, 1999 pp. 151-158. [18] M.J. Atallah, On symmetry detection. IEEE T. Comput, C-34 (1985) 663666. [19] N. Kiryathi, Y. Gofman, Detecting symmtery in grey level images: the global optimiza- tion approach, Int. J. Comput Vision 29 (1998) 29-45. [20] Z. Xiao, Z. Hou, C. Miao, J. Wang, Using phase information for symmetry detection, Pattern Recogn. Lett. 26 (2005) 1985-1994. [21] B. Nicolai, H. Holzmann, M. Pawlak, Testing for Image Symmetries with Application to Confocal Microscopy. IEEE T. Inform. Theory, 55(4) (2009) 1841-1855. [22] M. Birke, H. Dette, K Stahljans, Testing symmetry of a nonparametric bivariate regression function. J. Nonparametr Stat. 23(2) (2011), 547-565. [23] M. Chang, S.C. Park, Reverse engineering of a symmetric object.Comput & Ind. Eng., 55(2) (2008) 301-320. [24] T. Schulz, H. Ingensand, , Terrestrial Laser Scanning - Investigation and Applications for High Precision Scanning.FIG Working Week, Athens, Greece, 2004, pp. 1-15. [25] F. de As´ıs L´opez, C. Ord´on ˜ez, J. Roca - Pardi˜ nas, S. Garc´ıa - Cort´es, Point cloud comparison under uncertainty. Application to beam bridge measurement with terrrestrial laser scanning, Measurement 51 (2014) 259-264. [26] W.E. Deming, Statistical adjustment of data. Wiley, NY, 1985 [27] T. Hastie,W. Stuetzle, Principal Curves, J. Am. Stat. Assoc. 84 (1989) 502-516.
15
[28] F. de As´ıs - L´opez, C. Ord´on ˜ez, J. Roca - Pardi˜ nas, S. Garc´ıa - Cort´es, A statistical method for geometry inspection from point clouds, App. Math. Comput. 242 (2014) 562-568. [29] R.L. Eubank, Spline smoothing and nonparametric regression, New York, 1988. [30] W. Hardle, Applied nonparametric regression, Cambridge University Press, Cambridge, 1990. [31] P. Sarda, P. Vieu, Smoothing kernel regression, In: Smoothing and Regression: Approaches, Computation, and Application, New York: Wiley, 1999. [32] D. Ruppert D, M.P. Wand, R.J. Carroll, Semiparametric regression, Cambridge University Press, Cambridge, 2003. [33] S.N. Wood, Thin plate regression splines, J. R. Stat. Soc. B 65 (2003) 95-114. [34] B.W. Silverman, Density Estimation For Statistics And Data Analysis, Chapman and Hall, London, 1986. [35] D.W. Scott, S.J. Sheather, Kernel density estimation with binned data, Commun Stat - Theor M 14 (1985) 1353-1359. [36] P. Hall , M.P. Wand, On the accuracy of binned kernel density estimators, J. Multivariate Anal. 56(2)(1996) 165-184 [37] J. Fan, J.S. Marron, Fast implementation of nonparametric curve estimators, J. Comput. Graph. Stat. 3 (1994) 35-56. [38] M.P Wand, Fast Computation of Multivariate Kernel Estimators, J. Comp. Graph. Stat. 3 (1994) 433-45. [39] B. Efron Bootstrap methods: another look at the jacknife, Ann. Stat. 7 (1979) 1-26. [40] C.F.J. Wu, Jackknife, Bootstrap and Other Resampling Methods in Regression Analysis, Ann. Stat 14(4 )(1986) 1261-1295. [41] R.Y. Liu, Bootstrap Procedures under some Non-I.I.D. Models. The Annals of Statistics 16(4) (1988) 1696-1708. [42] E. Mammen, Bootstrap and Wild Bootstrap for High Dimensional Linear Models, Ann. Stat. 21(1) (1993) 255-285.
16
n= 1000
●
80
σ=0.01 σ=0.02 σ=0.03 σ=0.04 σ=0.05
60
60
80
100
100
n= 1000
● ●
40
40
●
●
●
●
20
20
● ● ●
● ●
0
0
●
4
6
8
10
0
2
4
6
a (mm)
a (mm)
n= 3000
n= 3000
8
100
2
100
0
10
●
●
80
80
● ●
60
60
●
●
40
40
●
●
20
20
● ●
●
0
0
●
0
2
4
6
8
10
0
2
4
a (mm)
6
8
10
a (mm)
100
n= 5000
100
n= 5000 ●
● ●
80
80
●
●
60
60
●
●
40
40
●
20
20
●
●
●
0
0
●
0
2
4
6
8
10
0
2
4
a (mm)
6
8
10
●
●
8
10
a (mm)
●
100
n= 10000
100
n= 10000 ●
●
80
80
●
60
60
●
40
40
●
●
20
20
●
●
0
0
●
0
2
4
6
8
10
0
a (mm)
2
4
6
a (mm)
Figure 3: Percentage rejection for the test for increasing a nominal levels of 5% (left) and 10%(right),and sample sizes of n = 1000 (top row), n = 3000 (second row), n = 5000 (third row) and n = 10000 (bottom row).
17
Figure 4: Yatch hull surface estimated from the point cloud.
0.5
1.0
1.5
Z (m)
2.0
2.5
3.0
Cross−sections
−2
−1
0
1
2
Y (m)
Figure 5: Cross-sections considered in the analysis of asymmetries.
18
1.50
1.60
2.0
1.70
1.80
X= −2.605
1.15
1.20
1.25
0.70
1.5 0.5
0.55
0.60
0.65
1.0
Z (m)
1.10
−1.5
−1.0
−0.5
0.0
0.5
1.0
1.5
0.00
0.05
0.10
0.15
0.20
1.00
1.05
1.10
1.15
1.20
Y (m)
0.86
0.7
0.90
0.8
0.94
0.9
Z (m)
1.0
1.1
1.06
1.10
1.2
1.14
X= 0.553
−1.5
−1.0
−0.5
0.0
0.5
1.0
1.5
0.50
0.55
0.60
0.65
Y (m)
Figure 6: Cross-sections of the estimated surface (left) and zooms (rigth), representing the estimation (solid line) and the symmetric (dotted line).
19