CANOPY LIDAR POINT CLOUD DATA K-MEANS CLUSTERING WATERSHED

Airborne laser LiDAR has widely applied in the accurate extraction of single tree canopy for inventory of precision forestry. Due to the over-segmentation phenomenon occurring in the traditional watershed single-wood segmentation, this paper presents a method, called K − means clustering watershed for single tree segmentation. This method consists of four aspects: The first step is to filter the point cloud to eliminate the interference factors such as ground elevation and other factors that interfere with the LiDAR point cloud segmentation; The second step is to optimize the generation of CHM , generate a CMM based on CHM variable window detection, and obtain the treetop position to provide the pixel center position for subsequent K − means cluster segmentation; The third step is to use the K − means clustering algorithm to perform initial cluster segmentation to extract the target pixels of interest. At this time, the local maximum value detected by the variable window in the second step is used as the center pixel of the cluster; In the fourth step, an improved watershed algorithm based on the similarity of 4 neighborhoods is proposed. The improved watershed algorithm is applied to the K − means initial clustering image to segment the target area, and the oversegmentation results are subsequently processed, and the over-segmentation blocks are combined according to certain criteria. Identify the contour of single canopy from the CHM images of the experimental forest data. The experimental results show that the proposed algorithm can effectively solve the over-segmentation problem happening the traditional watershed algorithm. The accuracy of F, R and P parameters can be improved by 7.1%, 11% and 9.8%.


INTRODUCTION:
Due to China's regular forest resource inventory system, there is an urgent need for accurate and effective control and management of forest information. Traditional passive remote sensing technology methods, such as multi-spectral imaging technology, are affected by weather and surface coverage (Gougeon et al., 2001;Sajdak et al., 2014;Unger et al., 2014 ) . Airborne laser LiDAR technology has unique advantages over passive remote sensing technology, providing theoretical basis and technical support for the accurate extraction of single tree canopy and the development of precision forestry (Zhou et al., 2004;Zhou et al., 2014;Zhou et al., 2018;Mu et al., 2020). At present, there are roughly three types of existing single-wood segmentation algorithms (Jakubowski et al., 2013;Lu et al., 2014 ): ①Pre-process the point cloud data, such as denoising, filtering, classification to obtain vegetation point cloud and ground point cloud, and then use Triangulated Irregular Network interpolation (TIN) to generate digital surface model (DSM), interpolating the ground points to generate a digital elevation model (DEM), calculating the difference between the two to obtain a canopy height model ( ), and then according to the method of computer graphics to the The model is divided into regions. ② Use the spatial structure relationship between LiDAR point cloud data and the true characteristics of trees to set up classification rules and directly segment the three-dimensional point cloud . Lin Yi et al. (2017) Based on the theory of circle  Corresponding author: Guoqing Zhou; E-mail: gzhou@glut.edu.cn detection in 2017, obtained local extreme points as the center point, and judged the distance from other point clouds to the center point for clustering to extract information such as the position, height and DBH(Diameter at Breast Height) of a single tree, Use the aerial imagery of the same position to extract the single wood position. In 2014, the iterative watershed method of Coops et al. (2007) was used to better identify the lower canopy in forests dominated by coniferous forests, but it is not good for small and dense canopies of the same height. ③Combining the characteristics of the gray image and the spatial structure characteristics of the real point cloud of the tree, segment the point cloud and single tree . Li Xiang (2015) used LiDAR data, orthophotos, and second-class survey data in 2015 to conduct single tree canopy extraction studies using coniferous forests and broad-leaved forests with high canopy density using the marker-controlled area growth method. It was found that this method can extract the canopy of the forest area very well. By controlling the threshold, the algorithm is more flexible and the processing speed is faster. Guo Yushan et al. (2016) used the mark control watershed segmentation method to extract the canopy of the sparse forest area and dense forest area respectively. Khosravipour et al. Proposed a new algorithm in 2014, as Part of the laser pulse passes through the canopy and is directly reflected from the ground, forming a "close pits" in the canopy，which uses a subset of LiDAR points to close pits to generate a pitless CHM grid, and then uses the variable window local maximum method for single wood recognition .
Among these methods, the watershed transformation algorithm combined with morphological algorithm and regional growth idea has the advantages of good closed edge extraction effect and high target segmentation accuracy. Watershed transformation has been widely used in image segmentation and processing. It has the advantages of fast, simple and direct. It also has better segmentation effect for images with lower contrast, which can simplify subsequent processing procedures, but this method is more sensitive to noise. There is a problem of image oversegmentation. To this end, this paper proposes a − clustering watershed algorithm based on grid ( Canopy maximum model) to segment the single wood based on the problems of over-segmentation and noise sensitivity of the watershed segmentation algorithm.

K-MEANS CLUSTERING IMPROVED WATERSHED ALGORITHM
Watershed is a region-based widespread technique for image segmentation. This method originated from mathematical morphology (Serra, 1982), and a break-through in applicability was achieved by Vincent and Soille (1991), since it can effectively extract the area of interest (AOI) of the imagery. Moreover, the watershed algorithm has the characteristics of less time-consuming, highly accurate locating at imagery edges, simple formulas and operation, and can carry out the parallel processing using computer relative to morphological segmentation algorithms. In order to simulate the flood process, it is necessary to first calculate the complement set of CHM and CHM (Figure 1), which is similar to the two catchment basins (Chen et al., 2006). The basic idea of the watershed algorithm is that an imaged LiDAR data set is regarded as a topographic landscape with ridges and valleys (Preim and Botha 2013). The elevation values of the LiDAR point cloud data are quantified by the gray values of the respective pixels. Based on such a gray (3D) representation, the watershed transform decomposes the imagery into catchment basins. For each of local minimum values, a catchment basin comprises all points whose path of steepest descent terminates at this minimum value. The watershed transform algorithm separates the basins from each other, and decomposes an imagery completely and thus assigns each pixel either to a region or a watershed. Due to noisy image data, the watershed algorithm could cause a large number of small regions segmentations. This phenomenon is known as the "oversegmentation" problem.
In order to overcome these shortcomings above, many scholars have proposed many algorithms to improve the watershed algorithm. The improvement methods can be roughly divided into: image pre-processing before watershed transformation and region merging after watershed transformation. The preprocessing before segmentation generally focuses on the noise reduction of the image. Its purpose is to filter out the irrelevant small targets in the image, ignore the small details inside the target, and then perform a watershed transformation to obtain a better segmentation. Most of the current improved watershed algorithms are the method of the pre-processing before segmentation. This article presents the K-means clustering algorithm to improve the over-segmentation caused by the watershed algorithm, the median filter is used to remove noise from the image, and the improved K-means clustering algorithm is used to perform multi-threshold segmentation on the gradient image to obtain the continuous boundary contour of the target area, thereby improving the over-segmentation problem of the watershed algorithm.

Coarse Segmentation of − Clustering
The − clustering algorithm is an unsupervised clustering algorithm. It performs initial clustering on areas with similar gray levels in the image to remove noise (Liu et al., 2012;Ren et al., 2013 ). The highest point of the local pixel maximum detection combined with the crown maximum model is taken as the highest point of the tree, and then each local highest point is used as the seed point of − clustering for clustering (Morsdorf et al., 2004) . The method is: Suppose ( , ) is the pixel coordinate of the digital image, and ( , ) is its gray scale function. Q ( ) represents the -th region after the -th clustering, and m ( +1) represents the average value of the -th category (cluster center) after the + 1th clustering.
Step 1. Suppose the initial clustering center is: 1 , 2 ,…, , The first clustering centers 1 , 2 ,…, represent the pixels at the treetop located by the variable window; Step 2. For each point ( , ), find the closest clustering center according to a given criterion, and assign it to the area of this type, namely: In the formula, m ( ) ( = 1,2, … , ) represents the average value of the -th class after the i-th cluster, then ( , ) ∈ Q ( ) . After the completion, update the clustering center after i + 1 clusters to get a new clustering center m ( +1) ; Step 3. The new clustering center of Q ( ) obtained in step 2 is (2) Let the value of 2 be the smallest, = 1, 2, … , k. Then: is the number of samples in Q ( ) ; Step 4. For all = 1，2， … ，k , if m ( +1) = m ( ) , the iteration is terminated, otherwise the above process loop is continued.
After − divides the image cluster into K categories, the tree canopy contour is generally an irregular and complete area under the influence of . Therefore, the maximum connected component ( = 1,2, . . , k) is taken as the singlewood area. Let be a subset of pixels in the image. If there is a path between all pixels of , the pixels are considered to be connected. In order to extract the target of interest to the maximum, this paper adopts the four-neighbored method. In the binary image, the contour line is extracted by judging the gray levels of other pixels in the four neighborhoods of each pixel of the CHM. There are:

Watershed Segmentation Based on K-means Clustering
In order to effectively solve the problem of over-segmentation in the traditional watershed algorithm for image segmentation. Continue to use the watershed algorithm on the basis of the above K-means algorithm processing, the specific steps are as follows: Step 1. Let K-means segment the image into blocks, Let each block be , ( , ) ∈ , ( , ), the gray value is ( , ), and the size of each block is , 0 < < . Calculate the gray variance of each block 2 .
Step 2. Let and be the two adjacent blocks, then define the gray-scale variance of the two adjacent blocks 2 as: Step 3. Let and be adjacent segment blocks, ( , ) ∈ , ( , ) ∈ , ( , ) are adjacent points on the edge of the segment block, ( , ) and ( , ) satisfy the fourneighborhood relationship.
is the number of pixels on the edge of the divided block that satisfy the four-neighborhood relationship, and the mean value of the standard deviation of the gray value of the four-neighbored pixels on the edge of the divided block is: Step 4. Define the similarity between two adjacent partitions as: Suppose that a threshold is initially selected to maximize the variance between 1 and 2 , 1 is the pixel with gray level {0,1, … , } , 2 is a group of pixels with gray level [ + 1, … , − 1] . The ratio of the between-class variance to the overall gray-scale variance of the image can divide the image gray into two categories. Generally, the value that can divide the image into two categories on the gray scale is called the optimal threshold, as follows: The value range is [0,1], and the generated threshold T is used to threshold the neighboring divided blocks.

Summary of The Study Area
The Saihanba forest area is located at the northernmost tip of Hebei Province, and is divided into two geomorphic units: the upper plateau of the dam and the mountainous area of northern Hebei. Saihanba belongs to a temperate semi-humid and semiarid monsoon climate zone. It is 58.6km long from north to south and 65.6km wide from east to west. The geographical coordinates are longitude l16 ° 51 '~ 117 ° 39' east, latitude 42 ° 02 '~ 42 ° 36' north, the forest area of Saihanba is 68708.05m2, and the forest coverage is 72.3 %. The experimental area selected for this study is 760 m2, the average flight speed when collecting data is 61.7m / s, the flight altitude is 150-200m, the zonal overlap is 67%, the data includes four echoes, and the point cloud density is 50 ( pts / m2). The tree species in the study area selected in this experiment are mainly pine trees, a total of 175 trees, with an average height of 20-30 m.

Point Cloud Filtering
When acquiring point cloud data, while receiving accurate echoes, due to the influence of equipment accuracy and environmental factors, there will be some errors. Such error-free points are generally significantly higher or lower than the feature points. The presence of noise points will affect the accuracy of subsequent work such as segmentation and parameter extraction. In this paper, the radius filter is used to denoise the point cloud. A circle is drawn around a point to calculate the number of points that fall outside the circle. Jiao Yitao chose 20m as the average tree height of forest stand based on the airborne LiDAR point cloud The radius, combined with the actual situation in this paper to determine 15m as a radius to draw a circle, when the number is less than the given value, then delete the point, the number is less than the given value to retain the point (Jiao et al., 2015 ). Since this method relies on empirical values, the error is relatively large.

Optimized Generation of
is a surface model that expresses the height of the forest vegetation canopy surface from the ground.In the process of acquiring and generating data, it is easy to produce pits, that is, pixels in the local range of the gray image are smaller than the pixels in the neighborhood. The formation of pits is related to many factors, such as:

1.
In the process of LiDAR data acquisition, inaccurate stitching of the flight belts is likely to cause elevation errors between points of the same name in the overlapping area. Using this point cloud data interpolation to generate is prone to pits; 2.
Due to the existence of the canopy gap, the laser pulse is penetrative to it, so the echo signals of the shrubs and ground under the forest canopy can be obtained, and the grid size when processing such laser point cloud interpolation to generate Too small can easily cause pits; 3.
In addition, because the values of DSM and DEM are different from the interpolation method, improper interpolation method is likely to cause pits in DSM, and then there are pits in the generated .

Variable Window Detection Treetop
Popescu and Wynne (2003) proposed a variable window method for detection of a tree height and crown size from LiDAR point cloud data. This Section presents the tree height is manually measured; and the crown size is measured through their average crown diameter along two vertical profiles. The window size can be calculated by: = 29.77 − 1.46 × + 0.03 × 2 (10) where TC is the crown diameter of each individual tree, TH is the average tree height of individual single tree.
The window size for tree crown detection is predicted by Equation (11) and is used to accurately detects the position of the treetop. In order to obtain canopy surface model, "on" operation to the image generated by the above operation is carried out. That is, the image is first eroded and then expanded, which is used to eliminate small objects, smooth the shape boundary, and do not change its area. It can remove small particle noise and break the adhesion between objects, in which each grid unit value is replaced by the height of the point cloud data in the neighborhood. Get the largest model of the canopy, that is, form the CMM model, and is expressed by The experimental result for model generation using Eq. (11) is depicted in Figure 4b. A comparison analysis with model is conducted and the result is depicted in Figure 4a. The experimental results demonstrate that the false treetops can be significantly eliminated, and the phenomenon of the missing treetops is largely reduced.

Gaussian Filtering
Gaussian filtering to smooth the generated , thereby obtaining . Gaussian filtering is a common method for suppressing insignificant local maxima in treetop detection, and it can flexibly apply standard deviation.The two-dimensional isotropic Gaussian function can be calculated using formula (12), where x is the distance from the horizontal axis to the origin, y is the distance from the vertical axis to the origin, and σ is the deviation of the standard Gaussian distribution.
The calculated by equation (13) is used as the basic data for the subsequent extraction of single tree structure parameters. Using to extract single tree crown vertices can effectively reduce false detection of tree crown vertices.

Accuracy Assessment and Comparative Analysis
In order to quantitatively evaluate the accuracy of the proposed K-means watershed algorithm for single tree extraction from LiDAR point cloud data, the tree recognition rate (noted as, R), recognition accuracy rate (noted as, P), and F value (F-score) were used as indexes, and are expressed by In the above formula, TP is the number of correctly divided trees; FN is the number of trees that have not been divided, but have been mistakenly divided into other trees, also known as missing points; FP is the number of trees that are more divided.
Where the extents of the F, R and P values are between 0 and 1. According to the F, R, P scores of the segmentation results, Precision quantitative analysis of the traditional watershed algorithm and K-means watershed algorithm. It can be seen from Table 1 that the F value based on − combined with the improved watershed segmentation method can reach 0.9 in the experimental area, which is better than the traditional watershed segmentation accuracy. Compared with the traditional watershed segmentation method, the improved watershed segmentation method combined with − clustering has higher tree recognition rate and segmentation accuracy rate. In addition, we can find that the recognition rate of traditional watershed segmentation is close to the improved watershed segmentation method, but its accuracy is low. Therefore, based on the data in the experimental area, from the three aspects of the F-score, Recall, and Precision segmentation accuracy judgment parameters, the following conclusions can be drawn: − combined with improved watershed segmentation can more effectively achieve crown segmentation, and effectively solve the emergence of traditional segmentation The problem of over-segmentation.
The − watershed algorithm proposed in This article is compared with the traditional watershed algorithm. − clustering method is used to segment initial coarse gray image, the result is shown in Figure 5(d); while the traditional watershed algorithm does and the result is depicted in Figure 5(c); the improvement based on − clusteringThe watershed algorithm divides the , as shown in Figure 5(f).
Finally, a single-wood segmentation result display based on Kmeans clustering watershed segmentation algorithm is formed (see Figure 6). The number at the top of the tree represents the tree ID. Partially enlarged 3D display

CONCLUSION AND ANALYSIS
In order to overcome the problem of over-segmentation and noise encountered in traditional watershed algorithm, This article presents − watershed algorithm for extraction of single treetop from airborne LiDAR point cloud data. In this method, − initializes clustering, then threshold segmentation is further conducted, finally the watershed transformation is applied to segment the single treetop regions. This method was confirmed in the airborne LiDAR data of theSaihanba Forest Proving Ground, which is 58.6 kilometers long from north to south and 65.6 kilometers wide from east to west.
The experimental results show that the P and R values calculated by the traditional watershed algorithm are lower than the K-means cluster watershed algorithm, which is the main factor affecting the comprehensive segmentation accuracy index F. The traditional watershed algorithm has a lot of oversegmentation and omission-segmentation and K-means clustering watershed algorithm can solve this problem to a large extent. Through quantitative precision analysis, it can be found that compared with the traditional watershed algorithm, the improved watershed algorithm has F, R, and P parameters increased by 7.1%, 11%, and 9.8%, respectively. It can be analyzed by the formula that the number of unsplit trees and wrongly split trees are reduced. Therefore, the segmentation result of the improved watershed algorithm is more accurate. Due to the high degree of forest canopy, the internal adhesion of CHM images is caused, and the clarity of the CHM image itself is limited, resulting in some trees not being well recognized. This will also be the next direction that this research needs to explore. It will consider improving the model from the internal enhancement and external contour recognition of the CHM image to improve the research of the internal segmentation line trend.