We propose a new distance measure between 2dimensional points to provide a total order for an entire point set and to reflect the correct geometric meaning of the naturalness of the point ordering. In general, there is no total order for 2dimensional point sets, so curve reconstruction algorithms do not solve the selfintersection problem because the distance used in the previous methods is the Euclidean distance. A natural distance based on Brownian motion was previously proposed to solve the selfintersection problem. However, the distance reflects the wrong geometric meaning of the naturalness. In this paper, we correct the disadvantage of the natural distance by introducing a polarnatural distance, and we also propose a new curve reconstruction algorithm that is based on the polarnatural distance. Our experiments show that the new distance adequately reflects the correct geometric meaning, so nonsimple curve reconstruction can be solved.
1. INTRODUCTION
The problems of reconstructing a shape from sample points appear in many scientific and engineering applications. Because of its practical importance, many algorithms have been proposed over the last two decades in the fields such as the reverse engineering of geometric models and image processing of medical images. Especially, curve reconstruction plays an important role in the shape reconstruction problems because most of boundaries of shapes are represented by the set of curves.
Curve reconstruction is the problem of computing a piecewise linear approximation to a curve from a set of sample points. So curve reconstruction consists of two steps: the first step is how to order the given points and the second step is how to interpolate the ordered points
[1]

[4]
. The first step is called the
point ordering problem
. Many approaches have been suggested for the point ordering problem. Edelsbrunner et al.
[5]
defined the αshapes of point sets as the underlying space of a simplicial complex. Attali proposed the
r
regular shape method, where the
r
regular shapes are characterized by requiring that any circle passing the points on the boundary has radius greater than
r
[6]
. These methods deal with only uniform sample points. For the nonuniform sample points, Amenta et al.
[7]
proposed the first algorithm to reconstruct a curve from nonuniform sample points with guarantee. This algorithm uses the Voronoi diagram and the Delaunay triangulation of the sample points. Dey et al.
[8]
proposed the nearest neighbor approach based on the properties of Voronoi diagrams. These algorithms work only under the assumption that the sample points are dense and do not work for nonsimple curve reconstruction because they reconstruct curves without the consideration of curve's orientation
[9]

[12]
. Kim et al.
[3]
proposed a point ordering approach with a natural distance to reflect the orientation. The natural distance is based on properties of Brownian motion so that the approach enables us to give 1dimentional order to points on 2dimentional space. However, the approach missed important geometric property: the right direction of Brownian motion.
In this paper, we introduce a new natural distance based on polar coordinates so that the distance enables our algorithm to reflect the direction properly and resolve the nonsimple curve reconstruction.
2. POLARNATURAL DISTANCE
In this section, we review the properties of Brownian motion and the natural distance
[3]
. First of all, we will discuss naturalness for a directionbased distance. In order to do this, we define the naturalness for a distance between points which is suitable to order the points on a 2D plane. For chosen consecutive points
p
_{i1}
and
p_{i}
which are already ordered, we define the ordering direction vector as vector
. For the ordering direction, what will be the next point
p
_{i+1}
?
Let
P
be the sample point sets. We define the candidate set
C
_{i+1}
(
r_{i}
) for the next point
p
_{i+1}
as
as shown in
Fig. 1
, where
d
(
p
,
q
) is the Euclidean distance between
p
and
q
, and
r_{i}
is the radius of the ordering neighborhood. The radius
r_{i}
is computed by the average length of the consecutive points:
The candidate set C_{i+1} for the next point p_{i+1}.
If a measure μ prefers to pick up the next point
q
∈
C
_{i+1}
satisfying the following conditions, then we say that the measure has a naturalness and call it as a
natural distance
:
(C1) Suppose
p_{i}
and
q
_{1}
,
q
_{2}
,⋯ ,
q_{m}
∈
C
_{i+1}
are on the same line. Then the measure prefers to choose the closest one
q_{k}
to the previous point
p_{i}
among
q_{j}
’s :
(C2) Suppose
q
_{1}
,
q
_{2}
,⋯
q_{m}
∈
C
_{i+1}
are on the same circle with center
p_{i}
, and the angle between
is strictly smaller than the angle between
. Then
q_{k}
is more preferable to
p_{j}
:
The condition (C2) implies that, among the points with the same Euclidian distance from
p_{i}
, we want to choose the point which makes the arc, made by connecting the chosen points, most smooth as shown in
Fig. 2
. To show the naturalness of distance, we adopt a simple property of one dimensional Brownian motion.
Geometric conditions for naturalness of ordering: (a) (C1) condition and (b) (C2) condition
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 [17], [18], Brownian motion has several properties which explain the irregular motion observed in nature; for example, in twodimensional 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 T, Brownian motion starting at 0 is a random variable normally distributed with mean 0 and variance T.
Roughly speaking, (P.B) means that the set of original points on the real line at time 0, which move randomly, become the set which is normally distributed with mean 0 and variance T after T time. Now, we will introduce the polarnatural distance function by understanding the right geometric meaning of the Brownian property.
We assume that
p
_{i1}
and
p_{i}
are already chosen and
q
is an arbitrary point in
C
_{i+1}
. In order to determine the next point
p
_{i+1}
, Kim et al.
[3]
defined the natural distance in the following manner. They applied a 2dimension transformation so that
p
_{i1}
is on the negative
x
axis and
p_{i}
is the origin of
R
^{2}
. In order to help the geometric understanding, let
R
(
p
_{i1}
,
p_{i}
) be the ray with the initial point
p_{i}
and the directional vector
and
L
(
q
,
p
_{i1}
,
p_{i}
) be the perpendicular line to
R
(
p
_{i1}
,
p_{i}
) passing through the point q satisfying (C1) as shown in
Fig. 3
. Therefore, the arbitrary point
q
can be represented by the new coordinates q = (
q_{t}
,
q_{s}
) as follows:
where
New coordinates for point q [3]
The coordinate
q_{t}
is the Euclidian distance between
p_{i}
and the intersection of
R
(
p
_{i1}
,
p_{i}
) and
L
(q,
p
_{i1}
,
p_{i}
), and the coordinate
q_{s}
is the signed Euclidian distance between
q
and the intersection of
R
(
p
_{i1}
,
p_{i}
) and
L
(
q
,
p
_{i1}
,
p_{i}
). Then we can think the pair (
q_{t}
,
q_{s}
) as the new coordinates of
q
. With the property (P.B) of Brownian motion, we further think of the
q_{t}
axis as the time axis with the original point
p_{i}
, and
q_{s}
axis as the axis representing the value of Brownian motion.
Therefore, if
q_{t}
> 0,
q_{s}
is considered as a sample point of the random variable which has the normal distribution with mean 0 and variance q
_{t}
as shown in
Fig. 4
, then we consider the point
q
as one of the normally distributed points after q
_{t}
time, which is randomly diffused from the original point.
The distribution of Brownian motion [3]
The value of distribution cannot directly be applied to the natural distance because the distribution changes according to the position of
q
. So we transform the distribution into the standard normal distribution so that the random variable
q_{s}
is transformed to the standard random variable
. It is well known that the larger the value of the standard random variable is, the smaller the probability is. Therefore, the problem of determining the next point
p
_{i+1}
is transformed into that of computing the minimum of the standard random variables. Kim et al.
[3]
defined a natural distance by the following manner:
where
K
is a constant and plays a role of a subjective weight.
It is well known that the natural distance of Kim et al.
[3]
marked a turning point in the research on point ordering problem. Most of the distances used in the previous point ordering research cannot give a total order to the set of points on
R
^{2}
, whereas the natural distance can do that. So we could solve the problem of reconstructing curves having selfintersections. However, there are obvious shortcomings in the natural distance. The first is that the random variable
used in the natural distance cannot reflect the right geometric meaning of the conditions (C1) and (C2. The second is that the result of point ordering is very sensitive to the value of subjective weight. The condition (C1) shows that all points on the same line with
p_{i}
have the same value
q_{s}
, the condition (C2) shows that all points on the same circle with center
p_{i}
have the same value
q_{t}
. So the random variable
cannot resemble the right geometric meaning of conditions. In this paper, we introduce a new natural distance which may give a total order to twodimensional points and reflect the right geometric meaning of conditions (C1) and (C2). First of all, we define the polarnatural coordinates as (
r_{i}
(
q
),
θ_{i}
(
q
)), where
. By using this coordinates and understanding the geometric meaning as shown in
Fig. 5
, we define the
polarnatural distance
as
The new distribution of Brownian motion based on the polarnatural coordinates
3. POLARNATURAL DISTANCE BASED POINT ORDERING ALGORITHM
In this section, we will describe a practical algorithm for ordering points in an unorganized point cloud. Assume that there are arbitrarily scattered n points, called sample points, in
R
^{2}
. Denote the set of all sample points by P. The point ordering problem aims to choose a subset
SP
of
P
, where
SP
is a wellordered set satisfying the natural conditions (C1) and (C2) mentioned in Section 2. In the above theoretical section, we have assumed that the first point
p
_{1}
and the second point
p
_{2}
are known. This easily enabled us to initialize the algorithm at the first point with an initial direction
. In general, the result of the point ordering problem is dependent on the initial direction. But, in practice,
p
_{1}
and
p
_{2}
are not known. Moreover, the wellordered subset SP may have several connected components, and may be open or closed according to the distribution of the data points. So we need to analyze the distribution of the sample points in order to find out the property of the subset
SP
. In general, the sample points of the unorganized point cloud for curve reconstruction may cluster in several groups and the points of each group are densely located near meaningful trajectories. So, the shape of each cluster looks like the union of bands. The objective of our algorithm is to remove outliers in each cluster and find out a wellordered subset that plays a role of skeletons of clusters. Our algorithm consists of three major steps: data clustering, determination of the initial direction, and local point ordering.
 3.1 Data Clustering
Most partitioning methods cluster objects based on the Euclidean distance between objects. Such methods can find only sphericalshaped clusters and encounter difficulty at discovering clusters of arbitrary shapes. We cannot apply such clustering methods directly to our problem because the shape of the input point set P is unknown. In this paper, we adopt the DBSCAN clustering algorithm based on the density, which is useful to filter out outliers and discover clusters of an arbitrary shape. The general idea is to continue growing the given cluster as long as the density in the neighborhood exceeds some threshold; that is, for each data point within a given cluster, the neighborhood of a given radius
cr
, called
clustering radius
has to contain at least a minimum number of points. By exploiting DBSCAN algorithm, we may obtain
nc
clusters
CL
^{1}
CL
^{2}
,…,
CL^{nc}
, in P of a bandlike shape. For each cluster
CL^{j}
, we choose a point which may well represent the character of the shape of the cluster as the initial point, denoted by
cp^{j}
.
 3.2 Determination of Initial Direction
Now, we will try to determine the initial direction of the natural distance at the initial point
cp^{j}
for each cluster
CL^{j}
. For the sake of convenience, we set the initial point as
p
_{1}
. In order to analyze the distribution of points near
p
_{1}
, we compute the neighborhood
N
(
p
_{1}
,
r
) of
p
_{1}
with a given radius r, called
ordering radius
;
N
(
p
_{1}
,
r
) = {
p
∈
CL^{j}

pd
(
p
_{1}
,
p
) ≤
r
}, where
pd
(
p
,
q
) is the PolarNatural distance between two points pand q, which is defined in Section 2. And then we apply the PCA (Principal Components Analysis) to the neighborhood because the principal component analysis of a set of points gives us the mean, an orthogonal frame, and the standard deviation of the neighbors. Let
q
_{1}
,
q
_{2}
,⋯,
q_{m}
be the points in the neighborhood
N
(
p
_{1}
,
r
). We compute the mean
m
of the neighbors such that
and then we construct the covariance matrix We compute the eigenvector of the covariance matrix which corresponds to the maximal eigenvalue so that it plays a role of the major principal axis in the neighborhood. Let
be the unit vector on the principal axis. Then the neighborhood
N
(
p
_{1}
,
r
) may be subdivided into two half discs:
Where
If
HD
_{–}
= ø or
HD
_{+}
= ø, then the point
p
_{1}
is one of the two end points of SP. Otherwise the point is a middle point of
SP
so that we should apply the following oneway ordering algorithm to both
HD
_{+}
and
HD
_{–}
with the direction vectors
, respectively. Let
SP
_{+}
and
SP
_{–}
be the solutions of the oneway ordering problem with the initial directions
, respectively. Then the solution
SP
of the original point ordering problem may be obtained by merging the subsets
SP
_{+}
and
SP
_{–}
, i.e.,
SP
=
SP
_{+}
∪
SP
_{–}
.
 3.3 Local Point Ordering
Now, we are ready to explain the oneway point ordering algorithm with the first point
p
_{1}
and the initial direction
. First of all, we have to select the candidates of the next point
p
_{2}
in the neighborhood
N
(
p
_{1}
,
r
). It is well known that the candidate
C
_{2}
set is the same as the half disc
HD
_{+}
(
p
_{1}
,
r
). So, without loss of generality, we may put
C
_{i+1}
=
HD
(
p
_{1}
,
r
). Next, we select the second point
p
_{2}
in
N
(
p
_{1}
,
r
) satisfying the following condition:
f
_{1}
(
p
_{1}
) = min
_{p∈HD+(p1)}
f
_{1}
(
p
), where the function
f
_{1}
(
p
) is derived from the initial direction vector
. We have known that the above condition satisfies the properties of conditions (C. 1) and (C. 2). Then the point
p
_{2}
is added to the wellordered set so that
SP
_{+}
= {
p
_{1}
,
p
_{2}
}.
In order to find out the third point
p
_{3}
we compute the candidate set
C
_{3}
that is the half disc
HD
_{+}
(
p
_{2}
,
r
) with the direction vector
. The third point
p
_{3}
should satisfy the following condition:
f
_{2}
(
p
_{2}
) = min
_{p∈HD+(p2)}
f
_{2}
(
p
), where the function
f
_{2}
(
p
) is derived from the direction vector
. Then
SP
_{+}
= {
p
_{1}
,
p
_{2}
,
p
_{3}
}. We apply the above process until the candidate set
C
_{i+1}
is empty or the intersection of
C
_{i+1}
and
HD
_{–}
(
p
_{1}
) is not empty. When the algorithm meets the former condition, it generates one of the two end points for an open curve, so it has to be applied to the other region of the cluster with the direction
. If the latter is satisfied, then the new obtained point
p
_{n1}
is near the first point
p
_{1}
.We can consider the wellordered set is a closed piecewise curve so that there is no need to continue this oneway point ordering algorithm in the cluster. The solution of the oneway point ordering problem is
SP
_{+}
= {
p
_{1}
,
p
_{2}
, ⋯ ,
p
_{n1}
}. By similar process, we may get
SP
_{–}
= {
p
_{n2}
,
p
_{n2+1}
, ⋯ ,
p
_{1}
} as the solution of the backward oneway ordering algorithm. Now, we may obtain the first connected wellordered subset SP =
SP
_{+}
∪
SP
_{–}
. Therefore, the global point ordering problem can be solved by applying the above process to all of clusters in P. The outline of our natural point ordering algorithm is as follows.
 3.4 Time Complexity
The three major steps of this algorithm are to partition the point set P into several clusters, to find out the neighbor points in
N
(
p
,
r
), and to select one point who has the minimum natural distance among them. The computational complexity of DBSCAN algorithm is
O
(
n log n
), where n is the number of the input cloud points. It is wellknown that the algorithm is sensitive to the userdefined parameters. The second process of this algorithm runs in
O
(
n
) time. The process is to find out a neighborhood
N
(
p
,
r
) of a given point p and to apply the principal component analysis to the neighborhood in order to derive the initial direction. The third process is repeatedly to find out a candidate set
C
_{i+1}
, and to select one point who has the minimum natural distance among the candidates. It takes
O
(
kM
), where
k
is the number of the output points in
SP
and
. Our algorithm is outputsensitive. However, in the worst case, the total time complexity of this algorithm may be O(
n
^{2}
), where the most of candidate sets have only one point as neighbors.
4. EXPERIMENTS
Our experimental results demonstrate the superiority of our point ordering algorithm based on the polarnatural distance. First of all, we compare the results of curve reconstruction for a sampled point set which are obtained by the natural distance and by the polarnatural distance.
The point set is obtained by sampling from points on two circles with different radius.
Figure 6
presents two piecewise linear curves obtained by the natural point ordering method (a) and the polarnatural ordering method (b). These figures show quite different curves. The curve generated by our distance is more natural than that by the natural distance in
[3]
. Moreover, the polarnatural poing ordering method can be used to reconstruct a curve for the sample points set containing the noise points not on the curve as well as those on it.
The results of curve reconstruction for a sampled point set: the natural distance and the polarnatural distance
Fig. 7
shows the results of curve reconstruction for a point cloud set. The cloud poing set was obtained by mouse dragginganddropping so that the global figure of the set has several selfintersections.
Fig. 7
(a) does not overcome the selfintersection problem, whereas
Fig. 7
(b) well reflects the global configuration of the data set.
The results of curve reconstruction for a point cloud set: (a) the natural distance and (b) the polarnatural distance
Fig. 8
shows the process of our polarnatrual point ordering algorithm.
Fig. 8
(a) is the initial point set. The set was obtained by mouse dragginganddropping so that the global figure of the set is the union of two ‘C’type point sets.
Fig. 8
(b) shows the middle state which has two curves. Each curve has two components
SP
_{+}
and
SP
_{–}
. The final result (c) has two curves. Since our algorithm has a special key that is polarnatural distance it can solve the selfintersection problem. Therefore our algorithm can be applied to nonsimple curve reconstruction.
The process of our point ordering algorithm: (a) point set (initial state), (b) each curve has two components (middle state), and (c) two curves (final state)
5. CONCLUSION
In this paper, we proposed an improved curve reconstruction algorithm by exploiting a new distance. We already introduced a natural distance in
[3]
which is based on a useful property of Brownian motion. The distance enables us to provide a totally order for the set of arbitrarily scattered points. However, the work of
[3]
missed the geometric meaning of the naturalness conditions (C1) and (C2). In this paper we proposed the polarnatural distance that is based on the polar coordinates. The totally ordered set generated by the polarnatural distance is more natural than that by the natural distance. By connecting these ordered points, we can obtain desired or meaningful curves. Also the polarnatural point ordering method can be used in the reconstruction of smooth curves, even though the data set contains noise points. The several experimental results show the superiority of our point ordering algorithm based on the polarnatural distance.
BIO
HyoungSeok Kim
He received the M.S. and Ph.D. degrees in department of Mathematics from KAIST (Korea Advanced Institute of Science and Technology) in 1992 and 1998, respectively. During 19981999, he worked at ETRI as a postDoctorial researcher. He is a professor in dept. of Multimedia Engineering at DongEui University from 1999. His interesting research fields are computer graphics, computer game engine, and applied mathematics.
HoSook Kim
She received the M.S. and Ph.D. degrees in department of Computer Engineering from Ewha Womans University in 1999 and 2005, respectively. She was an associate professor in division of Advance Computer Information at DongEui Institute of Technology from 2001 to 2009. She is a teacher in dept. of mathematics and computer science at Korea Science Academy of KAIST from 2010. Her interesting research fields are database, data mining, GIS, and mobile technology.
Hoppe H.
,
DeRose T.
,
Duchamp T.
,
McDonald J.
,
Stuetzle W.
“Surface Reconstruction from unorganized points,“
Proc. SIGGRAPH’ 92
1992
71 
78
Boissonnat J. D.
,
Geiger B.
“Three cimentional reconstruction of complex shapes based on the Delaunay triangulation,“
Proc. of Biomedical Image Process, Biomed. Visualizaiton
1993
964 
975
Kim Philsu
,
Kim HyoungSeok
2010
“Point ordering with the natural distance based on Brownian motion,”
Mathematical problems in engineering
article ID 450460
2010
17 
Boissonnat D.
,
Cazals F.
2002
“Smooth surface reconstruction via natural neighbor interpolation of distance functions,”
Comp. Geom. Theo. Appl.
22
185 
203
DOI : 10.1016/S09257721(01)000487
Edelsbrunner H.
,
Kirkpatrick D. G.
,
Seidel R.
1983
“On the shape of a set of points in the plane,”
IEEE Trans. Information Theory
551 
559
DOI : 10.1109/TIT.1983.1056714
Attali D.
“rregular shape reconstruction from unorganized points,”
Proc. of 13th Ann. Sympos. Comput. Geom.
1997
248 
253
Amenta N.
,
Bern M.
,
Eppstein D.
1998
“The crust and the betaskeleton: combinatorial curve reconstruction,”
Graphical Models
60
125 
135
DOI : 10.1006/gmip.1998.0465
Dey T. K.
,
Kumar T. K.
“A simple provable algorithm for curve reconstruction,”
Proc. 10th. ACMSIAM Symp. Discr. Algorithms
1999
893 
894
Fang L.
,
Gossard D. C.
“Fitting 3D curves to unorganized data points using deformable curves,”
Proc. of CG International'92
1992
535 
543
Dedieu J. P.
,
Favardin C. H.
1994
“Algorithms for ordering unorganized points along parametrized curves,”
Numerical Algorithms
6
160 
200
DOI : 10.1007/BF02149768
Taubin G.
,
Ronfard R.
1996
“Implicit simplical models for adaptive curve reconstruction,”
IEEE Trans. on Pattern Analysis and Machine Intelligence
18
321 
325
DOI : 10.1109/34.485559
Nguyen Thanh An
,
Zeng Yong
2008
“VICUR: A humanvisionbased algorithm for curve reconstruction,”
Robotics and Computer Integrated manufacturing
24
824 
834
DOI : 10.1016/j.rcim.2008.03.007