It is difficult to extract the boundary of complex planar points with nonuniform distribution of point density, concave envelopes, and holes. To solve this problem, an algorithm is proposed in this paper. Based on Delaunay triangulation, the maximum boundary angle threshold is introduced as the parameter in the extraction of the rough boundary. Then, the point looseness threshold is introduced, and the fine boundary extraction is conducted for the local areas such as concave envelopes and holes. Finally, the complete boundary result of the whole point set is obtained. The effectiveness of the proposed algorithm is verified by experiments on the simulated point set and practical measured point set. The experimental results indicate that it has wider applicability and more effectiveness in engineering applications than the state-of-the-art boundary construction algorithms based on Delaunay triangulation.

1. Introduction

The point set is basic spatial data for representing the shape of objects. GNSS survey, 3D laser scanning, aerial photography, and aerospace remote sensing technology are producing massive point set data. The accurate and efficient processing of these data is a challenge to building information modeling (BIM), geographic information (GIS), remote sensing (RS), computer information technology (IT), computer-aided design/manufacturing (CAD/CAM), and other related fields [1]. The boundary information of the point set is composed of discrete points representing the original contour features of the measured object. The rapid and efficient construction of boundary information from the discrete point set is fundamental and crucial for spatial data processing, which plays an important role in GIS-related fields such as geographic boundary determination [2, 3], building contour extraction [4, 5], map generalization [6, 7], statistics of the plot area, and the computation of road earthwork volume [8, 9]. The physical shape of the real object is diverse and can be simple or complex. This indicates that the point set boundary used to digitally represent the physical shape of the object is also diverse and complex. The complexity of this boundary is manifested by the following four problems: (1) the effective organization of the spatial relations (geometry, topology) of the points; (2) nonuniform distribution of the point density; (3) accurate expression of the information of the boundary shape (concavity and convexity); and (4) internal multihole problem. These are the main problems to be faced with the boundary construction of a complex planar point set.

At present, two kinds of algorithms are mainly exploited to construct the boundary of planar point set: the algorithm of point set based on convex hull and the algorithm of point set based on Delaunay triangulation. The former [1012] starts from the convex hull polygon that constitutes the point set and uses different contraction algorithms to obtain the boundary characteristics of the point set, which can handle the point set with nonuniform distribution and the point set with concave envelope region under certain conditions. The biggest drawback of this kind of algorithm is that it cannot deal with the internal multihole issues. The latter must first construct the Delaunay triangulation. Then, it extracts the corresponding boundary triangles, boundary line segments, or boundary points based on the geometric features of the triangulation to represent the contour characteristics of real objects. The Sculpture algorithm [13] and Alpha Shapes algorithm [14] are the main representatives of this kind of algorithm. The key point of the Sculpture algorithm and the improved algorithms such as algorithm [15] and algorithm [16] is to gradually delete triangles that are not needed from the outside to the inside of the triangulation, until a boundary point set that meets the given conditions is finally generated. In essence, the Sculpture algorithm is a triangulation stripping algorithm. The major drawback of this algorithm is that it is not stable enough to deal with a complex concave boundary and the boundary with internal holes, especially long and narrow holes. The Alpha Shapes algorithm was first proposed by Edelsbrunner. Its main feature is to give a strict mathematical definition of “point set shape.” On this basis, the algorithm constructs a two-dimensional point set boundary based on the Delaunay triangulation. The original Alpha Shapes algorithm is suitable for dealing with the boundary of the convex hull point set with uniform distribution, but it does not work well for a complex planar point set with concave envelopes and holes. Later, some improved versions of the Alpha Shapes algorithm, such as [17] and [18] are proposed. Based on local metric or global metric, these algorithms adopt different rules to solve the problems of constructing complex planar shape boundaries such as nonuniform density of point set and holes. However, some problems need to be addressed, such as too many required parameters and unsatisfactory extraction of the boundary of long and narrow holes.

In summary, Delaunay triangulation can easily record the spatial information of an irregular planar point set in a specific data structure. Also, it has strong extensibility, and its triangulation generation algorithm is mature and stable. Therefore, compared with other types of algorithms, the Delaunay triangulation is relatively robust for constructing the boundary of a planar point set. Based on Delaunay triangulation, a new boundary construction algorithm is proposed in this paper to solve the existing problems in the boundary construction of the planar point set. The proposed algorithm first processes the convex hull boundary of the triangulation as a whole to form a rough boundary and then refines the processed results through the corresponding mathematical model to achieve a complete construction of the complex planar point set boundary. The results of the simulated comparative experiments and practical measurement experiments indicate that the proposed algorithm is effective and feasible.

2. Relevant Definitions

For the convenience of algorithm description, some basic concept definitions and symbols that will be used in this paper are presented here.

Let be the planar point set and be the Delaunay triangulation generated based on . Let be the triangle vertices, triangle edges, and triangles in, respectively:(1)Boundary edge: in the triangle set of , the edge belonging to only one triangle is denoted as ; the result of the algorithm is a set of boundary edges, i.e., .(2)Boundary point: in , the vertex connecting the boundary edge is denoted as . The boundary points form the boundary point set, i.e., .(3)Interior point: in , the nonboundary point is the interior point and is denoted as . The interior points form the interior point set, i.e., .(4)Boundary triangle: in the triangle set of , a triangle containing at least one boundary edge is denoted as . The boundary triangles form the boundary triangle set, i.e., .(5)Boundary angle: in a boundary triangle, the interior angle corresponding to the boundary edge is denoted as .(6)Concave envelope [19]: given a polygon covering a point set , Conv represents the convex hull of . represents the open area between and , and each open area is represented by a closed curve . The area formed by these closed curves is referred to as the concave envelope. In , the edges in the concave envelope are referred to as concave edges, and each concave envelope consists of a visible boundary edge and several concave edges.(7)Holes: a concave envelope with only concave edges and no boundary edges is a special case.

3. Algorithm Overview

The proposed algorithm consists of two steps: extracting the rough boundary and extracting the fine boundary. For preprocessing, a stable algorithm is exploited to generate the Delaunay triangulation . Based on this, the convex hull of the point set is found, and the basic data is prepared for the smooth progress of the algorithm. Firstly, according to the parameter selection rule of the boundary angle , a rough boundary is extracted by filtering inwardly from the initial edge of the convex hull. Then, taking the maximum density and edge length of the point as parameters, the local areas such as the concave envelopes and holes are further refined to form the final complete boundary.

3.1. Extracting Rough Boundaries
3.1.1. Algorithm Flow

The Delaunay triangulation of the point set needs to be generated as the input data. The algorithm of generating Delaunay triangulation is very mature, and an incremental insertion algorithm is exploited to generate the Delaunay triangulation of a point set by the proposed algorithm. Triangle point, triangle edge, and triangle are three kinds of geometric entities that constitute Delaunay triangulation. As listed in Table 1, a data structure with topological relations for Delaunay triangulation is designed to save the storage space of the device and improve the operation efficiency of the proposed algorithm.

It can be seen from Table 1 that the topological relationship of the geometric entities in Delaunay triangulation is defined in terms of triangle edge. After the start point and end point of the triangle edge are defined, the triangle edge becomes a vector and can conveniently record its adjacent triangles.

According to the Gestalt proximity principle [20], the boundary of a planar point set is formed by connecting points with similar distances from their peripheral points. As for a boundary triangle, if the boundary edge is very long and exceeds the given threshold, it violates the principle of proximity and needs to be deleted. Therefore, the maximum boundary edge length can be used as a filtering condition to extract the boundary. However, the length threshold is a global metric. Once set, it is fixed for the entire process of processing the point set, which is not effective in dealing with the nonuniform distribution point set. In this case, using the maximum boundary angle as the boundary filtering condition is considered. Since the angle in the triangle is determined by the relative length of the three sides constituting the triangle, it is a local metric and is suitable for processing the nonuniform distribution of points. Thus, a maximum boundary angle threshold is chosen as the boundary filtering condition to extract the rough boundary of the point set. In the extraction process, to ensure the uniqueness of the results, the boundary edges need to be arranged in descending order according to the size of the boundary angle. The extraction process is illustrated in Figure 1.

Considering the rules for deleting boundary edges, two special cases need to be focused on, which are shown in Figure 2: (1) chain edge and suspended edge, which are not part of a triangle; (2) the intersection point of more than two boundary edges, namely. When deleting a boundary edge, it is necessary to perform a topology check on the new boundary. If (1) and (2) appear, the deletion operation should be undone.

3.1.2. Selecting Maximum Boundary Angle Threshold

In the process of extracting the rough boundary, the threshold value of the maximum boundary angle cannot be determined by a definite calculation formula. Instead, it needs to be determined by experiments. In the proposed algorithm, the norm [15] calculated by formula (1) is used to evaluate the accuracy of the extraction results:

In formula (1), represents the polygon of the original graph, and represents the polygon of the graph surrounded by the boundary that is extracted by the algorithm. The norm represents the boundary extraction effect through the proportional relationship between the polygon area of the original graph and the polygon area of the graph surrounded by the extracted boundary. The smaller the , the higher the coincidence degree between the extracted boundary and the original graphic boundary, and the better the effect.

In this paper, the experiment was conducted through three groups of symbol shapes, including isoline A, F, and K. The number of experimental points was, respectively, 100, 150, 200, and 250, representing different point densities. Meanwhile, the points were uniformly distributed in the shape, and it is ensured that the point set within each shape could cover the shape boundary. The range of the maximum boundary angle was [0, 180] with an interval of 10 degrees. Figure 3 illustrates the boundary extraction effect of the shape F, and Figure 4 shows the corresponding change curve. The experimental results indicate that the threshold value of the maximum boundary angle falls within the range of [60, 100], and the algorithm can achieve satisfactory boundary extraction results.

3.2. Extraction of Fine Boundary
3.2.1. Algorithm Flow

The boundary obtained by the rough boundary extraction can represent the contour shape characteristics of the point set to a certain extent. However, the expression of the complex concave envelopes is not precise enough, and the boundary of holes cannot be processed completely. It can be seen from the red triangle in Figure 5(a) that there may be very small boundary angles when irregular triangles are in the concave envelope region. As a result, the boundary edges may be ignored in the process of deletion, thus affecting the quality of boundary results. To delete such boundary edges, the threshold of the maximum boundary angle must be lowered. However, this will cause the boundary to shrink excessively and become broken. In this case, the ideal boundary cannot be obtained either. Considering the boundary of holes shown in Figure 6, since its boundary edge is not the edge of the boundary triangle in , the rule of extracting the boundary edge through the maximum boundary angle is invalid. Therefore, it is necessary to design special rules for concave envelopes and holes to extract fine boundaries.

Since a hole is a special case of a concave envelope, the concave envelope described below contains the characteristics of the hole unless otherwise specified. It can be seen from Figures 5 and 6 that in the point set containing concave envelopes, the point density of the concave envelope region is significantly lower than that of the normal region. From the normal region to the concave envelope region or vice versa, the point density changes greatly, and the change can be expressed by a parameter called point looseness. Let be any point in the Delaunay triangulation generated by the given point set , and is the looseness of the point, which can be calculated by the following formula:where and are the average length and the standard deviation of the line segment between point and the connecting point, respectively.

The smaller the value of , the more subtle the difference between the edge length of point and the connecting point, and the less likely is to be a boundary point; otherwise, the more likely is to be a boundary point. As shown in Figures 5(b) and 6(b), the looseness of the normal boundary points and concave envelope boundary points is large, while that of the interior points is small. Therefore, a looseness threshold can be set as a parameter to detect the concave envelope region, and the looseness set of the concave envelope region is presented in the following formula:

The boundary point set of concave envelopes can be obtained by . Specifically, a large looseness of the concave envelope boundary points indicates that the length of the line segments between the boundary points and the connecting points is greatly different. Filtering these line segments to remove the unnecessary long edges is essential to the following processing. An edge length constraint parameter is obtained by the following formula:where is the average length of the line segment between point and the connecting point and is the average value of all points in the given point set .

For any point in the point set processed by the rough boundary extraction, if the edge length of its connecting points is greater than , it should be deleted according to the boundary edge deletion rules. Otherwise, the retained connecting points and edges are processed following the rough boundary extraction method to generate the required boundary results. The fine boundary extraction process is shown in Figure 7.

It worth noting that the processing sequence of the fine boundary extraction cannot be reversed. The concave envelopes (without holes) should be processed before the holes are processed. As for the algorithm, it is necessary to filter the point looseness before the boundary edge. This can filter out the most irregular triangles. Besides, the point looseness traversal and edge length traversal of connecting points are arranged in a descending order to ensure the uniqueness of the boundary extraction result.

3.2.2. Determination of Looseness Threshold

Determination of the point looseness threshold is a crucial step in the process of fine boundary extraction. According to the definition of point looseness, the looseness of the boundary point is significantly greater than that of the interior point. The determination of the point looseness threshold can be regarded as a conic segmentation problem. The objective function shown in formula (5) is exploited to detect the point set and find the point in the transition region between the initial boundary and the final boundary to automatically determine . To this end, it is necessary to construct a set which contains the looseness of all points arranged in an ascending order:where is the th element of the set ; is the mean value of the set ; is the standard deviation of the set ; and is the adjustment parameter. The larger is, the larger the corresponding looseness threshold will be, and more boundary details will be retained.

After judged by the objective function, the looseness set of the point set is segmented into two parts. To achieve the best effect of conic segmentation (the minimum intraclass variance and the maximum interclass variance), the PBM index [21] is introduced to automatically solve the parameter . This index is a relative evaluation index, which can meet the clustering criteria of compactness and separation. The formulas for calculating this index are as follows:where represents the number of classification; represents the number of points in the classification ; represents the centroid of classification ; represents the internal distance of classification (the sum of distances between all points in the classification and the centroid of classification); represents the internal distance of the clustering of point set divided into classes; and represents the separation degree of the classification of the point set. The higher the PBM index, the higher the reliability of classification results.

4. Simulated Comparative Experiments

To verify the effectiveness of the proposed algorithm, this paper designed simulation experiments for three scenarios, including a nonuniform distribution point set, a uniform distribution point set with concave envelopes, and a uniform distribution point set with holes. Meanwhile, to verify the applicability of the proposed algorithm, the Delaunay triangulation-based algorithms and were taken as the reference. The experimental results and analysis were as follows.

4.1. Nonuniform Distribution Point Set

150 points were drawn manually in a circular area with a radius of 35 length units by the CAD software. As shown in Figure 8, the density of the points decreased from left to right. The two-dimensional coordinate data file of these points was exported from the CAD software and then imported into the algorithm verification program as the point set data. The results are shown in Figure 8 and Table 2.

The experimental results indicate that the three algorithms can successfully extract the boundary of the point set, and the roughness of the extracted boundary varies with the algorithm parameters. Since the algorithm used the threshold of the global measurement of the edge length as a parameter, the effect of boundary extraction was not stable. For a small , this algorithm might produce more sawtooth in the areas with smaller point densities. Both the algorithm and the proposed algorithm obtained better boundary results.

4.2. Uniform Distribution Point Set with Concave Envelopes

200 points were drawn manually with uniform distribution in a closed area with a width of about 20 length units. The closed area was drawn by the CAD software, and it was surrounded by the peripheral outline of the letter C in Arial Black font, as shown in Figure 9. The two-dimensional coordinate data file of these points was exported from the CAD software and then imported into the algorithm verification program as the point set data. The results are shown in Figure 9 and Table 3.

The experimental results indicate that the three algorithms can successfully extract the boundary of the point set with good effect. Compared with the algorithm and algorithm, the proposed algorithm provided a reasonable mathematical definition and detection method based on the fact that the conic segmentation index had a strong theoretical basis.

4.3. Uniform Distribution Point Set with Holes

280 points were drawn manually with uniform distribution in a circular area with a radius of 35 length units by the CAD software. Then, 40 points were deleted along the horizontal direction of the diameter to form a narrow hole, as shown in Figure 10. The two-dimensional coordinate data file of these points was exported from the CAD software and imported into the algorithm verification program as the point set data. The results are shown in Figure 10 and Table 4.

The experimental results indicate that the algorithm failed to extract the hole boundary, and the algorithm extracted the hole boundary only when is appropriate, while the proposed algorithm extracted the hole boundary with better effect than the algorithm when is appropriate. Although the proposed algorithm can extract the boundary of holes, it may also fail sometimes and the extraction result will be quite different from that of the actual situation if the parameter value is not reasonable. In this case, artificial discrimination is needed in practical applications.

5. Experimental Study on the Strip Terrain Point Set of Mountain Highway

To verify the effectiveness of the algorithm for the practical measured point set, 17 km long highway strip terrain data composed of more than 8000 GPS points were selected as the experimental data. The practical measured point set was obtained from a strip topographic map survey task undertaken by the author in May 2016, and the highway strip was located in a section of the rural road from Guquan town to Shuiyang town in Xuancheng city, Anhui province. The purpose of the task was to provide topographic maps for local road planning. The technical requirements of the task were as follows: (1) The topographic map should be in digital format and the scale should be 1 : 1000. (2) The width of the topographic map was 50 m along both sides of the road centerline. (3) The mapping content consisted of vegetations, river networks, roads, buildings, and structures along the centerline. The general situation of the task site is shown in Figure 11.

The survey task in May 2016 was not specifically performed to verify the effectiveness of the proposed algorithm. As for a mountainous highway, there are many belt-shaped terrain bends, dense collection points in high and low undulating areas, and sparse collection points in flat areas. This is characterized by a uniform distribution of the point set, many concave envelopes (as shown in the red circle in Figure 12(a), and some holes such as the pit ponds in some areas (as shown in Figures 12(c) and 12(d)). Thus, the GPS terrain point set obtained by the survey task is suitable to be used as the data source to verify the effectiveness of the proposed algorithm.

8156 GPS terrain points were used in the experiment, and the length of the road was about 17 km. The experiment was conducted on a laptop equipped with Intel® Core™ i7-10750H [email protected] GHz and 16 GB memory. The operating system was 64-bit Windows 10, and the experimental program of the proposed algorithm was implemented by C# program language based on Microsoft Visual Studio 2010 IDE. The experimental results are shown in Figure 12 and Table 5. Figure 12(b) to Figure 12(d) are the enlarged images of the boundaries extracted by the boundary angle thresholds at the highway curve delineated with red circle in Figure 12(a).

The experimental results indicate that the boundary of the concave envelopes in the northwest can be extracted well. Meanwhile, the extracted boundary results change from rough to fine as the maximum boundary angle threshold changes from large to small. This result is consistent with that of the simulation experiment. The number of holes extracted by the proposed algorithm changes from less to more as the maximum boundary angle threshold α changes from large to small. However, the correct number of holes in accordance with the actual position does not increase. For a larger , some holes with smaller area are ignored, and a smaller number of holes are extracted. For a smaller , the tiny holes are retained, and the number of holes increases. It should be noted that the extracted holes are not completely correct. Some of them are “pseudo holes” caused by the nonuniform distribution of practical measured points and need artificial discrimination.

6. Conclusion

Based on the Delaunay triangulation and two constraint parameters, i.e., the threshold of maximum boundary angle and the threshold of point looseness , an algorithm is proposed in this paper for constructing the boundary of the complex planar point set with nonuniform point density distribution, concave envelopes, and holes. The proposed algorithm consists of rough boundary extraction and fine boundary extraction. Compared with the existing algorithms, the proposed algorithm has two main advantages:(1)The proposed algorithm has wide applicability, and it can be applied to construct complex planar boundary with nonuniform point density distribution, concave envelopes, holes, etc.(2)A reasonable mathematical definition and a detection method based on conic segmentation index are given for processing concave envelopes, which provides a strong theoretical basis for the proposed algorithm.

This paper has not given a quantitative analysis on the execution efficiency of the proposed algorithm, namely, the time complexity and the space complexity. At present, the application of this algorithm is only limited to the two-dimensional planar point set. Further research should be conducted on the quantitative analysis of efficiency and three-dimensional surface reconstruction.

Data Availability

The data used to support the findings of this research are generated from experiments.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.


This research was funded by the Nature Science Research Project of Anhui Education Department (KJ2018JD04), the Science and Technology Project of Housing and Urban and Rural Construction in Anhui Province (2020-YF23), and the Natural Science Research Project of Anhui University (KJ2020JD01). The authors would like to thank MJ Editor (http://www.mjeditor.com) for its linguistic assistance during the preparation of this manuscript.