Point Ordering with Natural Distance Based on Brownian Motion

Abstract

We consider a curve reconstruction problem from unorganized point clouds with noise. In general, the result of curve reconstruction depends on how to select and order the representative points to resemble the shape of the clouds. We exploit a natural distance based on a property of one-dimensional Brownian motion to order sample points, which simultaneously reflect smoothness and nearness of points, so that our algorithm is able to reconstruct not only simple curves but also nonsimple curves. Numerous examples show that this algorithm is effective. The natural distance proposed in this paper is able to play an important role in a variety of fields of measuring the distance of points with considering direction.

1. Introduction

In many engineering fields such as the reverse engineering of geometric models and image processing of medical images, great attention has been shown to the problems of reconstructing a shape from sample points. The curve reconstruction problem we want to discuss in this article plays an especially important role in shape reconstruction problems (see, e.g., [1–4] and references therein) and has also been instrumental in the recent realizations of many algorithms.

Curve reconstruction problems can be regarded as the problem of fitting a curve or a family ℱ of curves to an irregularly spaced set 𝒫 of unorganized data points. Most existing techniques for solving the problem can be divided into two categories depending on whether the sample points of 𝒫 are distributed exactly on the family ℱ of curves without noise or are distributed around ℱ with some noise (possibly a point cloud around ℱ). The techniques in the first category can also be classified depending on situations where the curves are either closed or open, either smooth or nonsmooth, and the sample points are distributed either uniformly or non-uniformly and so on. In any case, the techniques in the first category are well developed by many researchers (e.g., see [5–10] and the references therein) and most of these techniques give mathematical guarantees. There are a few algorithms working for uniformly sampled curves such as alpha shape [9], 𝑟-regular shape [6], and EMST [11]. There are other algorithms that work for nonuniformly sampled curves such as Crust [5], Nearest Neighbor [8], Conservative Crust [12], Traveling Salesman Path [13], and Gathan [14, 15]. The algorithms provide a sampling condition with respect to geometric properties, which may result in a sample set not natural to human vision.

Recently, there are two algorithms which use two important properties of human vision: proximity and smoothness [16, 17]. Although the sample points are non-uniformly distributed, they can well reconstruct open simple and nonsmooth curves because they use several local features such as the mean of the distances and the mean of the angles of already-obtained edges. The creative vision function of VICUR [16] is a mathematical representation with experimental foundation. They considered only proximity and smoothness of the sample points on the Delaunay triangulation. That is, the candidates for the next point are limited to only points adjacent to the last point although the other points may produce smoother reconstructions. Thus, this is not able to reconstruct nonsimple curves and solve the curve reconstruction problem with noise.

Even though the second category, handling data points with noise, is closer to practical situations, relatively little research has been carried out on the second category (e.g., see [10, 18–20] and the references therein). Fang and Gossard [10] used a method based on spring energy minimization to approximate an unorganized point set with a curve, which requires a good initial guess to find the solution. Taubin and Ronfard [20] reconstructed a planar curve from unorganized points using an implicit simplicial curve, which is defined by a planar triangular mesh and the values at the vertices of the mesh. Lee [19] proposed an algorithm based on a modified least-squares method and Euclidean minimum spanning tree. This method is not appropriate for curves with self-intersections. Most curve reconstruction algorithms mentioned above are based on the Euclidean distance to compute the proximity and the adjacency without consideration for smoothness, so they do not solve nonsimple curve reconstruction.

Assume that the set 𝒫 of data points we want to handle fits into the second category. More precisely, let 𝒫={𝑞_𝑘=(𝑥_𝑘,𝑦_𝑘)∈ℝ²∶𝑘=1,2,…,𝑛}. As indicated by Dedieu and Favardin in [7], one of the main steps for solving the curve reconstruction problem is to give an order on a subset of 𝒫 to fit a curve interpolating the points in the subset. If we denote the ordered points as 𝑝₁,𝑝₂,𝑝₃, and so on, and we assume that 𝑝₁,…,𝑝_𝑖 are already ordered, then the central problem is to find a way of reasonably selecting a next successive point 𝑝_𝑖+1 from a subset 𝒞_𝑖⊂𝒫⧵{𝑝₁,…,𝑝_𝑖}, where we assume that the known facts are only the coordinates of points. For the case of the set 𝒫 in the first category, it is sufficient to estimate the angle between the lines for three consequence points, or estimate the Euclidean distance between nested points [7]. However, a separate usage of these two factors (angle and Euclidean distance) is not sufficient in the presence of noise.

What we try to do in this paper is to expound a new persuasive measure (called a natural distance) giving a smoothness and closeness between points, and applicable to order the points in the second category by combining naturally two facts: angle and Euclidean distance. The natural distance to pick up the next consecutive point 𝑝_𝑖+1 from 𝒞_𝑖 is defined by𝑓_𝑖(𝑞)=𝑡_𝑖𝑠(𝑞)+𝐾^𝑖(𝑞)𝑡_𝑖(𝑞),𝑞∈𝒞_𝑖, (1.1) where 𝑡_𝑖(𝑞)=|𝑝_𝑖𝑞|cos𝜃,𝑠_𝑖(𝑞)=|𝑝_𝑖𝑞|sin𝜃,𝜃=𝜃(𝑞)=cos⁻¹(𝑝_𝑖−1𝑝_𝑖⋅𝑝_𝑖𝑞|𝑝_𝑖−1𝑝_𝑖||𝑝_𝑖𝑞|). (1.2) Hereafter, we call 𝑡_𝑖(𝑞) a time distance, 𝑠_𝑖(𝑞)/𝑡_𝑖(𝑞) a standardized probability distance, and 𝐾 a subjective weight. By taking into account simple properties of Brownian motion, we will interpret a naturalness about the distance 𝑓_𝑖(𝑞) ordering the points in a subset of 𝒫. In particular, we will show that, the time distance 𝑡_𝑖(𝑞) is a diffusing time measuring the closeness between points, while the standardized probabilistic distance 𝑠_𝑖(𝑞)/𝑡_𝑖(𝑞) is a function of measuring the smoothness between points. Finally, we will show that, under an appropriate lower bound of the subjective weight 𝐾, the natural distance 𝑓_𝑖(𝑞) provides a plausible answer in ordering the points of 𝒫. Also we will provide remarkable results through several figures of numerical experiments.

2. Natural Distance

In this section, we will discuss a naturalness for the distance 𝑓_𝑖(𝑞) defined in (1.1). In order to do this, we must first define a naturalness for a measure between points applicable to order the points. For chosen consecutive points 𝑝_𝑖−1 and 𝑝_𝑖 which are already ordered, let 𝒞_𝑖 be a subset of 𝒫 such that for any 𝑞∈𝒞_𝑖 the angle between 𝑝_𝑖−1𝑝_𝑖 and 𝑝_𝑖𝑞 is acute. If a measure 𝜇 prefers to pick up the next point 𝑝_𝑖+1∈𝒞_𝑖 satisfying the following conditions, then we say that the measure has a naturalness and call it a natural distance. (C.1) Suppose that 𝑝_𝑖 and 𝑞₁,𝑞₂,…,𝑞_𝑚∈𝒞_𝑖 are on the same straight line. Then the measure prefers the closest one to the previous point 𝑝_𝑖 among 𝑞_𝑗, 𝑗=1,2,…,𝑚.(C.2) Suppose that 𝑞₁,𝑞₂,…,𝑞_𝑚∈𝒞_𝑖 are on the same circle with center 𝑝_𝑖 and the angle between 𝑝_𝑖−1𝑝_𝑖 and 𝑝_𝑖𝑞_𝑘, 1≤𝑘≤𝑚, is the smallest angle among those between 𝑝_𝑖−1𝑝_𝑖 and 𝑝_𝑖𝑞_𝑗 for all 𝑗=1,2,…,𝑚,. Then the measure prefers 𝑞_𝑘 among 𝑞_𝑗, 𝑗=1,…,𝑚 as the next point of 𝑝_𝑖.

Two conditions (C.1) and (C.2) are similar to the first and second criterions of [7], respectively. To show a naturalness of 𝑓_𝑖(𝑞), we will use a simple property of one-dimensional Brownian motion introduced in the following subsection.

An easier way to understand Brownian motion is to think of it as a continuous version of the random walk process. As documented in an enormous amount of literature on stochastic processes such as [21, 22], Brownian motion has several properties which explain the irregular motion observed in nature; for example, in two-dimensional space the motion of pollen grains suspended in liquid which does not flow. Among those, we use only the following property throughout this paper. (P.B) For given time 𝑇, 𝐵(𝑇,⋅) starting at 0 is a random variable normally distributed with mean 0 and variance 𝑇, where 𝐵(𝑇,⋅),𝑇≥0 is the Brownian motion process.

Roughly speaking, (P.B) means that the set of original points on the real line at time 0, which move randomly, becomes the set which is normally distributed with mean 0 and variance 𝑇 after 𝑇 time.

2.1. A Motivation of the Distance 𝑓_𝑖(𝑞)

Before proving the naturalness for 𝑓_𝑖(𝑞), we first explain a motivation of the derivation of 𝑓_𝑖(𝑞). To do this, we begin with the following interpretation for the set 𝒞_𝑖. For two ordered points 𝑝_𝑖−1 and 𝑝_𝑖, (𝑖≥2), which were previously chosen, let 𝑅(𝑝_𝑖−1,𝑝_𝑖) be the ray with initial point 𝑝_𝑖 and directional vector 𝑝_𝑖−1𝑝_𝑖 and let 𝐿(𝑞,𝑝_𝑖−1,𝑝_𝑖) be the perpendicular line to 𝑅(𝑝_𝑖−1,𝑝_𝑖) passing through the point 𝑞∈𝒞_𝑖+1 satisfying (C.1). Also define the signed angle 𝜃  (−𝜋/2<𝜃<𝜋/2) between 𝑝_𝑖−1𝑝_𝑖 and 𝑝_𝑖𝑞 by ⟨𝑝cos𝜃=^𝑖−1𝑝_𝑖,𝑝_𝑖𝑞⟩|𝑝_𝑖−1𝑝_𝑖||𝑝_𝑖,𝑞| (2.1) where the sign of 𝜃 is positive if 𝑞 is located on the left-hand side of 𝑅(𝑝_𝑖−1,𝑝_𝑖) and is negative on the right-hand side. Here, ⟨⋅,⋅⟩ denotes the inner product in 𝐑². Define two functions 𝑡_𝑖∶𝒞_𝑖+1→𝐑⁺ and 𝑠_𝑖∶𝒞_𝑖+1→𝐑 by 𝑡_𝑖(𝑞)=|𝑝_𝑖𝑞|cos𝜃,𝑠_𝑖(𝑞)=|𝑝_𝑖𝑞|sin𝜃. (2.2)

The function 𝑡_𝑖(𝑞) is the Euclidean distance between 𝑝_𝑖 and the intersection of 𝑅(𝑝_𝑖−1,𝑝_𝑖) and 𝐿(𝑞,𝑝_𝑖−1,𝑝_𝑖), and the function |𝑠_𝑖(𝑞)| is the Euclidean distance between 𝑞 and the intersection of 𝑅(𝑝_𝑖−1,𝑝_𝑖) and 𝐿(𝑞,𝑝_𝑖−1,𝑝_𝑖). Then we can consider the pair (𝑡_𝑖(𝑞),𝑠_𝑖(𝑞)) as the new coordinates representation of 𝑞 (see Figure 1).

Figure 1 

New coordinates for a scattered point $q$ with $\theta >0$ when $i=2.$

With the property (P.B) of Brownian motion, we further think of the 𝑡_𝑖-axis as the time axis with the original point 𝑝_𝑖, and 𝑠_𝑖-axis as the axis representing the value of Brownian motion. Therefore, if 𝑡_𝑖(𝑞)=𝑡>0, 𝑠_𝑖(𝑞) is considered as a sample point of the random variable which has the normal distribution with mean 0 and variance 𝑡 as shown in Figure 2, then we consider the point 𝑞 as one of the normally distributed points after 𝑡 time, which is randomly diffused from original point.

Figure 2 

The distribution of Brownian motion after time $t=t_2(q).$

With this setting, we want to find a function 𝐹 providing the distances between 𝑝_𝑖 and candidate points of 𝑝_𝑖+1, which prefer to choose a point 𝑞∈𝒞_𝑖+1 satisfying the conditions (C.1) and (C.2) where 𝐹 has the minimum value. Then what is suitable for variables of 𝐹? We suggest two factors as variables of 𝐹.

We adopt the function 𝑡_𝑖, called 𝑖th time distance, as its first variable. The condition (C.1) in the previous section requires that 𝐹 should be strictly increasing for 𝑡_𝑖. This implies that we prefer to choose the point taking the smaller diffusion time to 𝑝_𝑖+1.

Similarly thinking, one may take the function |𝑠_𝑖| as the second variable of 𝐹. As shown in Figure 3, the probability density 𝑝 satisfies that, for arbitrary 𝜖>0 and two points 𝑞,𝑞^′∈𝒞_𝑖+1 such that 𝑡_𝑖(𝑞)=𝑡_𝑖(𝑞^′)=𝑡 and |𝑠_𝑖(𝑞)|>|𝑠_𝑖(𝑞^′)|, _{𝑁(𝑠𝑖(𝑞),𝜖)}𝑝(𝑡,0,𝑠)𝑑𝑠<_{𝑁(𝑠𝑖(𝑞′),𝜖)}𝑝(𝑡,0,𝑠)𝑑𝑠. (2.3)

Figure 3 

The probabilities of the $ϵ$-neighborhoods of two points $q,\hspace{0.17em}\hspace{0.17em}q\text{'}$ for $i=2.$

However, the function 𝑠_𝑖 is not enough as the second variable of 𝐹 because it is not independent of the first variable 𝑡_𝑖. Assume that 𝑞,𝑞^′∈𝒞_𝑖+1 are the points having different reaching times, that is, 𝑡_𝑖(𝑞)≠𝑡_𝑖(𝑞^′). Since 𝑠_𝑖(𝑞)(𝑠_𝑖(𝑞^′)) is interpreted as a sample point of the random variable with a normal distribution with variance 𝑡_𝑖(𝑞)(𝑡_𝑖(𝑞^′),resp.), we cannot compare the probability in the same way as above. For this reason, it is necessary to standardize random variables with normal distributions with the variance 𝑡_𝑖 to be the standardized random variable with normal distribution with the same mean and the variance 1. So we adopt the function |𝑠_𝑖|/𝑡_𝑖 as the second variable of 𝐹. We call this function 𝑖th standardized probability distance. Easily speaking, if we ignore the effect of the 𝑖th time distance, we think that the original point 𝑝_𝑖 is more naturally attainable to 𝑞^′ than 𝑞 whenever |𝑠_𝑖(𝑞^′)|/𝑡_𝑖(𝑞^′)>|𝑠_𝑖(𝑞)|/𝑡_𝑖(𝑞) (see Figure 4). The condition (C.2) in the previous section requires that 𝐹 should be strictly increasing for this second variable, though this is not a sufficient condition for (C.2).

Figure 4 

Standardizations of $s_i(q)$ and $s_i(q\text{'})$ when $i=2.$

Consequently, our problem is rewritten as the following, if 𝐹 is assumed to be 𝐶¹. Choose 𝑞∈𝒞_𝑖+1 satisfying the conditions (C.1) and (C.2) where the suitable (still undetermined) 𝐶¹-function 𝐹(𝑡_𝑖(𝑞),|𝑠_𝑖(𝑞)|𝑡_𝑖)(𝑞) (2.4) has the minimum value. Here the function 𝐹(𝑥,𝑦)∶𝐑⁺×𝐑⁺∪{0}→𝐑⁺ strictly increases for both variables, that is, 𝐹_𝑥>0 and 𝐹_𝑦>0.

There are plenty of such increasing function 𝐹 with two variables. Among all, the simplest one will be the linear function of the form 𝐹(𝑥,𝑦)=𝑥+𝐾𝑦, (2.5) for a positive constant 𝐾. This is mostly caused by the convenience of calculation, but the function has several good features as the solution for the point ordering problem. The reader can check these, exploring some examples in the next section.

In this spirit, we define the 𝑖th natural distance 𝑓_𝑖  (𝑖≥2) by, for 𝑞∈𝒞_𝑖+1 and a positive constant 𝐾,𝑓_𝑖(𝑞)=𝑡_𝑖(𝑞)+𝐾|𝑠_𝑖(𝑞)|𝑡_𝑖.(𝑞) (2.6) We call 𝐾 the subjective weight, since it represents how much weight is given to the second factor.

2.2. Naturalness of 𝑓_𝑖(𝑞)

To prove the naturalness of 𝑓_𝑖(𝑞), we need to make some technical assumptions. Choosing points on the same level curve is problematic, where the points have the same natural distance from 𝑝_𝑖. If points lie on the same level curve with a different Euclidean distance from 𝑝_𝑖, we choose the closest point from 𝑝_𝑖 to the Euclidean distance. Also, if two points lie on the same level curve with the same Euclidean distance from 𝑝_𝑖, we choose the point lying on the right-hand side of the direction 𝑝_𝑖−1𝑝_𝑖.

The 𝑖th natural distance 𝑓_𝑖 has the following properties, which can be easily checked.

(P.1) If 𝑞₁ lies in the interior of the level curve 𝑡_𝑖+𝐾(|𝑠_𝑖|/𝑡_𝑖)=𝑓_𝑖(𝑞), then 𝑓_𝑖(𝑞)>𝑓_𝑖(𝑞₁),(P.2) If 𝑓_𝑖(𝑞)<𝑡_𝑖(𝑞₂), then 𝑓_𝑖(𝑞)<𝑓_𝑖(𝑞₂),(P.3) If 2√3𝑓_𝑖(𝑞)^3/2/9𝐾<|𝑠_𝑖(𝑞₃)|, then 𝑓_𝑖(𝑞)<𝑓_𝑖(𝑞₃),(P.4) If 𝑞₄ lies on the level curve 𝑡_𝑖+𝐾(|𝑠_𝑖|/𝑡_𝑖)=𝑓_𝑖(𝑞), then 𝑓_𝑖(𝑞)=𝑓_𝑖(𝑞₄).

Some of these can be used to simplify numerical algorithms ordering the points. Here 𝑞,𝑞_𝑗∈𝒞_𝑖+1 for 𝑗=1,2,3,4 (see Figure 5). (P.1) implies that the level curves do not intersect each other in (𝑡_𝑖,𝑠_𝑖)-plane, so the natural distance 𝑓_𝑖(𝑞) satisfies (C.1). (P.3) is derived by the fact that the maximum value of the level curve 𝑡_𝑖+𝐾(|𝑠_𝑖|/𝑡_𝑖)=𝑓_𝑖(𝑞) is 2√3𝑓_𝑖(𝑞)^3/2/9𝐾. In the following theorem, we will show that the natural distance 𝑓_𝑖(𝑞) locally satisfies (C.2).

Figure 5 

The level curve passing through a point $q$ for $i=2.$

Theorem 2.1. Assume that 𝐾>(2783+48√3)^1/4√𝑟₀ (2.7) for some fixed positive 𝑟₀. Then the natural distance 𝑓_𝑖(𝑞) satisfies the condition (C.2).

Proof. To show (C.2), it is enough to show that, for any fixed 𝑐>0 and 0<𝑟<𝑟₀, the system of equations 𝑥²+𝑦²=𝑟²,𝑦𝑥+𝐾√𝑥=𝑐 (2.8) has at most one root in the first quadrant 𝑥>0,  𝑦>0.
The second equation of (2.8) gives the equation 𝑦²=𝑥(𝑥²−2𝑐𝑥+𝑐²)/𝐾². Hence from the first equation of (2.8), we have 𝑔(𝑥)=𝐾²𝑟²,where𝑔(𝑥)=𝑥(𝑥²+(𝐾²−2𝑐)𝑥+𝑐²). (2.9) Thus, it suffices to prove that two curves 𝑦=𝑔(𝑥) and 𝑦=𝐾²𝑟² meet at most one point in the first quadrant 𝑥>0,𝑦>0. Note that, if the cubic function 𝑔(𝑥) has three zeros, all of them must be less than zero because 𝑥²+(𝐾²−2𝑐)𝑥+𝑐²=(𝑥−𝑐)²+𝐾²𝑥>0, and hence 𝑔(𝑥)>0 provided 𝑥>0. In this case, the intersection of two curves 𝑦=𝑔(𝑥) and 𝑦=𝐾²𝑟² on the first quadrant 𝑥>0,𝑦>0 is at most one.
First, assume that 𝑔(𝑥)=0 has only one zero at 𝑥=0. If the quadratic equation 𝑔^′(𝑥)=0 does not have solution or has negative solutions, then 𝑔(𝑥) is a strictly increasing function for 𝑥>0, and hence the system has at most one root. Therefore the only case that has to be considered is the case 𝐾²<(2−√3)𝑐,orequivalently𝑐>(2+√3)𝐾². (2.10) Note that, in this case, the quadratic equation 𝑔^′(𝑥)=0 has two distinct positive solutions 𝑥₁=2𝑐−𝐾²−𝑋3,𝑥₂=2𝑐−𝐾²+𝑋3, (2.11) where 𝑋=𝑐²−4𝑐𝐾²+𝐾⁴, and the minimal extreme value of 𝑔(𝑥) is obtained at 𝑥₂. Using the relation 𝑔^′(𝑥₂)=0, we can derive 𝑔(𝑥₂2)=−9𝑥₂(𝐾⁴−4𝑐𝐾²+𝑐²1)−9𝑐²(𝐾²2−2𝑐)=−9𝑥₂𝑋²−19𝑐²(𝐾²=1−2𝑐)27(−2𝑋³+2𝑐³+𝐾²(15𝑐²−12𝐾²𝑐+2𝐾⁴)). (2.12) From (2.10) it can be easily verified that 𝑋 has the maximum value 𝑐. So 𝑔(𝑥₂) is bounded such as 𝑔(𝑥₂𝐾)≥²27(15𝑐²−12𝐾²𝑐+2𝐾⁴√)≥83+483𝐾27⁶. (2.13) The last inequality is derived from 𝑐>(2+√3)𝐾². Hence if the inequality 𝑟²𝐾²≤83+48√3𝐾27⁶,orequivalently𝐾>(2783+48√3)^1/4√𝑟, (2.14) holds, then the proof is complete.

From Theorem 2.1, we know that points lying on the same level curve with the same Euclidean distance to the points should be located symmetrically with respect to the 𝑡_𝑖-axis. Thus if there exist such points, then there are only two in 𝑁(𝑝_𝑖,𝑟₀).

Now, we need to make a practical remark concerned with Theorem 2.1. As seen in its proof, the lower bound for 𝐾 is not obtained optimally, so we conjecture that this bound can be smaller in most cases. For practical usage of our method, we think that it is not inevitably necessary to give a consideration for the bound, since the purpose of our point ordering is to construct (or reconstruct) meaningful shapes (or the original shape) after connecting all ordered points. For this reason, in the following experiments, we are not concerned with the lower bound for 𝐾 on using the natural point ordering method.

Before closing this section, we remark on the effect of the subjective weights 𝐾 on the point orderings. It is obvious that, if 𝐾 is very big, the 𝑖th natural distance 𝑓_𝑖 is more sensitive to magnitudes of the changes of the 𝑖th standardized probability distance |𝑠_𝑖|/𝑡_𝑖 than those of the 𝑖th time distance 𝑡_𝑖. This fact can be checked in Figure 6, which presents the level curves of (2.6) for several values of the subjective weight 𝐾. As 𝐾 increases, the width between the line 𝑠_𝑖=0 and level curve becomes narrower.

Figure 6 

Level curves for various subjective weights $K.$

3. Practical Algorithm

In this section, we will describe a practical algorithm for ordering points in an unorganized point cloud as well as a sample point set.

Let 𝒫 be a set of arbitrarily scattered 𝑛 points in 𝐑². The point ordering problem aims to choose a family 𝒮={𝒫_𝑗} of subsets of 𝒫, where each subset 𝒫_𝑗 is a well-ordered set satisfying the natural conditions (C.1) and (C.2). Assume that 𝑝^𝑗₁,𝑝^𝑗₂,…,𝑝^𝑗_𝑚 are already ordered so that we may denote that 𝒫_𝑗=[𝑝^𝑗₁,𝑝^𝑗₂,…,𝑝^𝑗_𝑚] with the last point 𝑝^𝑗_𝑚 and the direction 𝑝^𝑗_𝑚−1𝑝^𝑗_𝑚. We define the connectivity 𝑡^𝑗_𝑚 of an ordered subset as the mean of the lengths of edges of the ordered subset: 𝑡^𝑗_𝑚=(1/(𝑚−1))∑^𝑚−1_𝑙=1𝑑^𝑗_𝑙, where 𝑑^𝑗_𝑙 is the distance between 𝑝^𝑗_𝑙 and 𝑝^𝑗_𝑙+1.

We explain a process for computing the set of candidates for the next point and choosing the next point among the candidates, called NextPointFinder. Let 𝑁(𝑝^𝑗_𝑚,𝑟^𝑗_𝑚) be a neighborhood with center 𝑝^𝑗_𝑚 and radius 𝑟^𝑗_𝑚=𝑎𝑡^𝑗_𝑚 satisfying the following condition: 𝑁(𝑝^𝑗_𝑚,𝑟^𝑗_𝑚)={𝑞∈ℛ∶𝑑(𝑝^𝑗_𝑚,𝑞)≤𝑟^𝑗_𝑚}, (3.1) where ℛ is the remaining point set such that ℛ=𝒫⧵⋃^𝑗_𝑙=1𝒫_𝑙 and 𝑎 is a parameter. We define the candidate set 𝒞_𝑚 by 𝒞_𝑚={𝑞∈𝑁(𝑝^𝑗_𝑚,𝑟^𝑗_𝑚)∶𝑝^𝑗_𝑚−1𝑝^𝑗_𝑚⋅𝑝^𝑗_𝑚𝑞>0} (3.2) as shown in Figure 7. We find out the point 𝑝^𝑗_𝑚+1 among the elements of 𝒞_𝑚 such that 𝑓_𝑚(𝑝_𝑚+1)=min_{𝑞∈𝒞𝑚}{𝑓_𝑚(𝑞)}.

Figure 7 

The candidate points of $p_{m+1}^j$ are in the right half-disc ${𝒞}_m$ with center $p_m^j$ and radius $r_m^j.$

There is a main difference between our method and VICUR [16] in the process. Our method permits all of points 𝑞∈𝒞_𝑚 as a candidate for 𝑝_𝑚+1 whether [𝑝_𝑚,𝑞] is a Delaunay edge or not, whereas VICUR considers only the points 𝑞, where [𝑝_𝑚,𝑞] is a Delaunay edge. We think that it plays a role of a key to solve a nonsimple curve reconstruction problem.

In order to deal with the first category problem, we assume that there is no noise in the point set 𝒫.

NADIAS (NAtural DIstance based Algorithm for Sample points)
Step 1. Compute Delaunay triangulation DT(𝒫).Step 2. Let 𝒫_𝑗 be the 𝑗th local solution generated by NextPointFinder. Let 𝑒^𝑗=[𝑝^𝑗₁,𝑝^𝑗₂] be the shortest Delaunay edge in ℛ. Apply NextPointFinder with the last point 𝑝^𝑗₂ and the direction 𝑝^𝑗₁𝑝^𝑗₂ until the candidate set is empty. Let 𝒫_𝑗+ be the solution of this process. Apply NextPointFinder with the last point 𝑝^𝑗₁ and the direction 𝑝^𝑗₂𝑝^𝑗₁ until the candidate set is empty. Let 𝒫_𝑗− be the solution of this process. Merge the two subsets 𝒫_𝑗+ and 𝒫_𝑗− into 𝒫_𝑗, that is, 𝒫_𝑗=𝒫_𝑗+∪𝒫_𝑗−. Update 𝒮=𝒮∪{𝒫_𝑗}. Update ℛ=𝒫⧵⋃^𝑗_𝑙=1𝒫_𝑙.Step 3. If |ℛ|>1, then repeat Step 2.Step 4. If |𝒮|>1, then apply the following bridging process to a pair of 𝒫_𝑖 and 𝒫_𝑗 of 𝒮: for two ordered subsets 𝒫_𝑖=[𝑝^𝑖₁,𝑝^𝑖₂,…,𝑝^𝑖_𝑙] and 𝒫_𝑗=[𝑝^𝑗₁,𝑝^𝑗₂,…,𝑝^𝑗_𝑚], if the distance between 𝑝^𝑖₁ (or 𝑝^𝑖_𝑙) and 𝑝^𝑗₁ (or 𝑝^𝑗_𝑚) is less than the average of the connectivity 𝑡^𝑖_𝑙 and 𝑡^𝑗_𝑚, then our algorithm connects these two ordered subsets.

Step 2 says that our algorithm follows a type of an incremental greedy approach: at each step it chooses the next point from among the candidate points that is best under the natural distance, including “close to the last point” and “in the direction of current travel”. The bridging process is similar to the curve-curve connectivity of [16] and makes our method reconstruct curves with sharp corners.

The time complexity of our algorithm is 𝑂(𝑛log𝑛), where 𝑛 is the number of the input points because the dominant process of our algorithm is to compute Delaunay triangulation of them and it runs in 𝑂(𝑛log𝑛). The main process of our algorithm is NextPointFinder which is repeatedly to find out a candidate set 𝒞_𝑖, and to select one point that has the minimum natural distance among the candidates. It takes 𝑂(𝑘𝑀), where 𝑘 is the number of the output points in 𝑃 and 𝑀 is a constant.

For the second category problem, we outline our algorithm as follows.

NADIAC (NAtural DIstance based Algorithm for Cloud points)
Step 1. Compute Delaunay triangulation 𝐷𝑇(𝒫). Step 2. Choose a point 𝑝₁∈𝒫. Compute a neighborhood 𝑁(𝑝₁,𝑟_𝑐), where 𝑟_𝑐 is a constant radius. Apply Principal Component Analysis to 𝑁(𝑝₁,𝑟_𝑐) in order to obtain the initial direction ID which plays a role of the major principal axis in the neighborhood. Step 3. Apply NextPointFinder with the last point 𝑝₁ and the direction ID until the candidate set is empty. Let 𝑃₊ be a local solution generated by NextPointFinder. Step 4. Apply NextPointFinder with the last point 𝑝₁ and the direction −ID until the candidate set is empty. Let 𝑃₋ be a local solution generated by NextPointFinder. Step 5. 𝑃=𝑃₊∪𝑃₋.

4. Experimental Results

In this section, we tested several reconstruction problems where each problem is focused on the following facts: (1) nonsimple curve from sample points, (2) sharp corners from sample points, (3) the difference between Natural distance and Euclidean distance, (4) the effect of the subjective weight 𝐾, and (5) the treatment of cloud sets.

Figure 8 shows that our method can reconstruct a nonsimple curve from sample points. Figure 8(a) is the Delaunay triangulation of the sample points. In order to explain the difference between VICUR [16] and our method, we zoom in the Delaunay triangulation as shown in 8(b). The adjacent points to 𝑃 are 𝐴, 𝐵, and 𝐶, whereas 𝑄 is not adjacent to 𝑃. So VICUR does not permit 𝑄 as a candidate for the next point of 𝑃 although the edge 𝑃𝑄 is smoother to the already-ordered edge 𝑂𝑃 than the other Delaunay edges. On the other hand, our algorithm deals with all of four points as candidates for the next point of 𝑃 since they belong to the neighborhood of 𝑃 with radius 𝑎⋅𝑡, where 𝑡 is its connectivity and 𝑎=2.5 is set. Figure 8(d) shows a nonsimple curve reconstructed by our method. The red point in Figure 8(d) means the initial point of our algorithm.

(a)

(b)

(c)

(d)

Figure 8 

Nonsimple curve reconstruction: (a) Delaunay triangulation, (b) enlarged Delaunay triangulation at the dense area, (c) reconstructed curve generated by [16], and (d) reconstructed curve generated by our method.

Figure 9 shows another example for a nonsimple curve from sample points as shown in Figure 9(a). In this example, we choose the initial point as that having a minimum edge length between points. In Figure 9(b), the boundary point between red line and green line is the chosen initial point. Thus, our algorithm works in two opposite directions. Figures 8 and 9 are smooth nonsimple curves.

(a)

(b)

Figure 9 

Nonsimple curve reconstruction: (a) sample points and (b) result by NADIAS.

Figure 10 shows the process of our algorithm, called NADIAS (NAtural DIstance based Algorithm for Sample points). Figure 10(a) is an example for sample points on a piecewise smooth curve. For the sample data, NADIAS finds out the shortest Delaunay edge and then applies NextPointFinder process two times with opposite directions. Figure 10(b) shows the result which has two components of smooth curves. The boundary of the red line and the green line is the initial point of the first component, and the boundary of the black line and the blue line is one of the second component. Figure 10(c) is the result after our bridging process. This figure contains one component because the minimum distance between the end points of two components is smaller than the given threshold distance. Figure 11 shows that NADIAS is able to well extract the several components whether or not each component is open and whether or not its shape is sharp.

(a)

(b)

(c)

Figure 10 

The process of NADIAS: (a) sample points, (b) NextPointFinder, and (c) bridging.

(a)

(b)

(c)

Figure 11 

Examples with multicomponents.

There are two different results of point ordering for the same set of sample points by using Euclidean distance and Natural distance in Figure 12. The ordered set with Euclidean distance in (a) does not contain all sample points, contrary to that with natural distance in (b). It implies that there is an example in which parameters can be chosen for the natural point ordering that will reconstruct the curve, while parameters chosen for Euclidean distance may not.

(a) Euclidean distance [5]

(b) Natural distance

Figure 12 

Point orderings based on the different distance measures.

In Figure 13, the complicated examples for three different subjective weights are presented. For the set of scattered sample points in (a), the case 𝐾=0.5 of (c) produces a meaningful closed curve of `@' almost completely. One also finds the fact that, as 𝐾 increases, our method produces smoother shapes. Figures 13(b) and 13(c) are generated by the first algorithm which regards the point data as sample points, so there are two curves. On the other hand, Figure 13(d) shows the result of the second algorithm. 𝐾 is very large, so the smoothness is emphasized more than the proximity. Therefore the algorithm meets the stop condition before most of points are applied to the algorithm.

(a) sample data

(b) K=0.3

(c) K=0.5

(d) K=3.0

Figure 13 

A complicated example using the natural point ordering method.

Figure 14 shows a more meaningful example of the natural point ordering problem. The sample data are generated by drawing a spiral and then adding three small circles to the spiral. One may understand that the major configuration of the dataset is a spiral. Our natural point ordering method can cut out the minor configurations such as three small circles from the spiral, so the ordered set is only on the spiral.

(a)

(b)

Figure 14 

Extraction of the major configuration: sample points (a) and the result (b) with $K=0.8.$

Figure 15 shows the results of NADIAC for point clouds. The cloud data of a spiral type are generated by random process and the clouds of the right example are generated by drawing a nonsimple curve in several times. Although the set of sample points is unorganized, our natural point ordering algorithm can extract the ordered subset from the point cloud, which plays a role of representative in their local neighborhood.

(a)

(b)

Figure 15 

Natural point orderings for point clouds: a point cloud generated by random process (a), a point cloud generated by drawing with mouse (b).

5. Conclusion

We have presented a new method for constructing curves from unorganized point clouds with noise. In general, the result of curve reconstruction depends on how to select and order the representative points to resemble the shape of the clouds. In this paper, in contradiction to the previous curve reconstruction algorithms based on Euclidean distance, we exploit a natural distance to reflect orientation to the ordering of sample points, so that our algorithm is able to reconstruct not only simple curves but also nonsimple curves. Moreover, for unorganized point clouds, this method efficiently extracts the skeletons of the clouds by cutting out the outliers, even though the result by our method is sensitive to the initial point and the initial direction. The method may be extended in two ways. The first thing is to extend the method to the surface reconstruction problem in three-dimensional space. We think that one should use the properties of two-dimensional Brownian motion in order to get a natural distance similar to that in this paper. Second, it is possible to improve our method which can control the subjective weight 𝐾 to obtain the wanted shape. We conjecture that this can be fulfilled if we weaken the condition (C.1).

Acknowledgments

The authors would like to express their sincere gratitude to reviewers for many helpful comments and valuable suggestions on the first draft of this paper. This work was supported by Korea Research Foundation Grant KRF-2003-041-C00039.