Background Point Filtering of Low-Channel Infrastructure-Based LiDAR Data Using a Slice-Based Projection Filtering Algorithm

A light detection and ranging (LiDAR) sensor can obtain richer and more detailed traffic flow information than traditional traffic detectors, which could be valuable data input for various novel intelligent transportation applications. However, the point cloud generated by LiDAR scanning not only includes road user points but also other surrounding object points. It is necessary to remove the worthless points from the point cloud by using a suitable background filtering algorithm to accelerate the micro-level traffic data extraction. This paper presents a background point filtering algorithm using a slice-based projection filtering (SPF) method. First, a 3-D point cloud is projected to 2-D polar coordinates to reduce the point data dimensions and improve the processing efficiency. Then, the point cloud is classified into four categories in a slice unit: Valuable object points (VOPs), worthless object points (WOPs), abnormal ground points (AGPs), and normal ground points (NGPs). Based on the point cloud classification results, the traffic objects (pedestrians and vehicles) and their surrounding information can be easily identified from an individual frame of the point cloud. We proposed an artificial neuron network (ANN)-based model to improve the adaptability of the algorithm in dealing with the road gradient and LiDAR-employing inclination. The experimental results showed that the algorithm of this paper successfully extracted the valuable points, such as road users and curbstones. Compared to the random sample consensus (RANSAC) algorithm and 3-D density-statistic-filtering (3-D-DSF) algorithm, the proposed algorithm in this paper demonstrated better performance in terms of the run-time and background filtering accuracy.


Introduction
As a type of active vision sensor, light detection and ranging (LiDAR) is a 3-D point cloud imaging vision sensor with the advantages of insensitivity to external light changes, strong adaptability to complex environments, wide coverage, and is informative. Currently, the most typical application of LiDAR in intelligent transportation systems is to detect road and traffic information for automatic driving vehicles [1][2][3][4]. Despite the rapid development of technologies related to self-driving cars in recent years, we cannot realize automatic driving by merely using on-board detectors. Fatal accidents may occur if the detectors of an autonomous vehicle make errors. For example, Tesla and Uber have experienced fatal accidents during their self-driving tests [5,6]. Therefore, it is very important for (2) The point cloud is classified into four categories in a slice unit: Valuable object points (VOPs), worthless object points (WOPs), abnormal ground points (AGPs), and normal ground points (NGPs). Based on the point cloud classification results, the traffic objects and their surrounding information would be easily identified from an individual frame of point cloud. (3) A machine learning model based on a backpropagation artificial neural network (BP-ANN) is used to dynamically fit the road gradient or employing the inclination of the LiDAR sensor, as to improve the adaptability of the algorithm in a complex environment.
The rest of the paper is structured as follows: In Section 2, we introduce the LiDAR sensor and the type classification of data points. In Section 3, we give a detailed description of the assumptions and definitions proposed in this paper. In Section 4, we illustrate the detailed methods for background filtering. The experiments and validation analysis are presented in Section 5. In Section 6, we summarize the contributions of our work, as well as the perspectives on future work.

LiDAR Sensor and Points Classification
The working principle of LiDAR is to obtain the distance information of objects by firing and receiving laser beams. There are four important parameters: The number of laser channels, rotating frequency, vertical angle resolution, and horizontal angle resolution. These parameters determine the amount of points and the speed of the LiDAR sensor output point cloud. In this paper, we develop the point cloud filtering algorithm based on a 360 • Velodyne VLP-16 LiDAR sensor that contains 16 laser channels. This sensor can detect ±15 • in the vertical range, and the vertical angle resolution is 2 • . We set the rotating frequency of LiDAR to 600 RPM and the corresponding horizontal angle resolution approximately 0.2 • . Figure 1 demonstrates a selected frame of a raw 3-D LiDAR point cloud from the VeloView (Velodyne LiDAR, San Jose, CA, USA) software and field photo corresponding to LiDAR firing.
Sensors 2020, 20, x FOR PEER REVIEW 4 of 25 (2) The point cloud is classified into four categories in a slice unit: Valuable object points (VOPs), worthless object points (WOPs), abnormal ground points (AGPs), and normal ground points (NGPs). Based on the point cloud classification results, the traffic objects and their surrounding information would be easily identified from an individual frame of point cloud. (3) A machine learning model based on a backpropagation artificial neural network (BP-ANN) is used to dynamically fit the road gradient or employing the inclination of the LiDAR sensor, as to improve the adaptability of the algorithm in a complex environment.
The rest of the paper is structured as follows: In Section 2, we introduce the LiDAR sensor and the type classification of data points. In Section 3, we give a detailed description of the assumptions and definitions proposed in this paper. In Section 4, we illustrate the detailed methods for background filtering. The experiments and validation analysis are presented in Section 5. In Section 6, we summarize the contributions of our work, as well as the perspectives on future work.

LiDAR Sensor and Points Classification
The working principle of LiDAR is to obtain the distance information of objects by firing and receiving laser beams. There are four important parameters: The number of laser channels, rotating frequency, vertical angle resolution, and horizontal angle resolution. These parameters determine the amount of points and the speed of the LiDAR sensor output point cloud. In this paper, we develop the point cloud filtering algorithm based on a 360° Velodyne VLP-16 LiDAR sensor that contains 16 laser channels. This sensor can detect ±15° in the vertical range, and the vertical angle resolution is 2°. We set the rotating frequency of LiDAR to 600 RPM and the corresponding horizontal angle resolution approximately 0.2°. Figure 1 demonstrates a selected frame of a raw 3-D LiDAR point cloud from the VeloView (Velodyne LiDAR, San Jose, CA USA) software and field photo corresponding to LiDAR firing. The LiDAR sensor outputs 3-D point cloud in spherical coordinates, which can be expressed as ( , , ), where is the distance between a data point and the LiDAR, is the vertical angle fixed and given by the Laser ID (identity document number), and is the horizontal rotation angle provided by the data packet of LiDAR. In this paper, we deal with the point cloud after projection in XYZ coordinates. The conversion relationship of the data point from spherical coordinates to XYZ coordinates is shown a s Figure 2 [29]. The LiDAR sensor outputs 3-D point cloud in spherical coordinates, which can be expressed as (d, ω, α), where d is the distance between a data point and the LiDAR, ω is the vertical angle fixed and given by the Laser ID (identity document number), and α is the horizontal rotation angle provided by the data packet of LiDAR. In this paper, we deal with the point cloud after projection in XYZ coordinates. The conversion relationship of the data point from spherical coordinates to XYZ coordinates is shown a s Figure 2 [29]. The spherical to XYZ conversion formulas are [30]:

LiDAR
A good filtering algorithm should retain as much valuable information as possible and remove as much worthless information as possible. That is, a good filtering algorithm should not simply filter all ground points. Instead, it should also retain all points with useful information for ground points, such as steps and curbs, which can also be used for road users' path prediction. There are also worthless information points in object points, such as leaves, which have little impact on traffic. Therefore, we divided the points in a point cloud into four categories in this paper: Valuable object points (VOPs), worthless object points (WOPs), abnormal ground points (AGPs), and normal ground points (NGPs).
VOPs contain points that can affect traffic, such as tree trunks, road users, and lighting poles. As for WOPs, they mainly include points that cannot be used to judge objects and have no impact on traffic, such as noise and leaves. For NGPs, road surface points that can be normally used by vehicles are included. The algorithm proposed in this paper aims to remove WOPs and NGPs from the point cloud as much as possible. According to AGPs and VOPs, we can analyze pedestrians, vehicles, and their surrounding information (such as nearness to curbstones, the relative location of the road surface, adjacent position to a light pole), which are valuable inputs for microcosmic traffic behavior analysis, pedestrian zebra-crossing behavior analysis, driver or pedestrian intention prediction, trajectory tracking in traffic accident prevention, and other various intelligent transportation applications. Figure 3 shows the schematic diagram of LiDAR scanning objects in a data frame of the point cloud. In the point cloud, the point and are the points firing on the ground by LiDAR. The point , , and are the points firing on the front object . The point , , and are within the same vertical block and form a slice within object . Point and are the points firing on the rear object . Point and are within the same vertical block and form a slice within object . The point , , , and are the points of the point , , , and mapping to the horizontal plane, respectively. As shown in Figure 3a,b, one object is composed of several slices in the same vertical plane. A vertical plane contains all the objects' slices, which are in the same direction, as shown in Figure 3c. We define that a slice is a point set containing the points firing on the same object and within the same vertical block. The spherical to XYZ conversion formulas are [30]:

Definition of Slice
A good filtering algorithm should retain as much valuable information as possible and remove as much worthless information as possible. That is, a good filtering algorithm should not simply filter all ground points. Instead, it should also retain all points with useful information for ground points, such as steps and curbs, which can also be used for road users' path prediction. There are also worthless information points in object points, such as leaves, which have little impact on traffic. Therefore, we divided the points in a point cloud into four categories in this paper: Valuable object points (VOPs), worthless object points (WOPs), abnormal ground points (AGPs), and normal ground points (NGPs).
VOPs contain points that can affect traffic, such as tree trunks, road users, and lighting poles. As for WOPs, they mainly include points that cannot be used to judge objects and have no impact on traffic, such as noise and leaves. For NGPs, road surface points that can be normally used by vehicles are included. The algorithm proposed in this paper aims to remove WOPs and NGPs from the point cloud as much as possible. According to AGPs and VOPs, we can analyze pedestrians, vehicles, and their surrounding information (such as nearness to curbstones, the relative location of the road surface, adjacent position to a light pole), which are valuable inputs for microcosmic traffic behavior analysis, pedestrian zebra-crossing behavior analysis, driver or pedestrian intention prediction, trajectory tracking in traffic accident prevention, and other various intelligent transportation applications. objects' slices, which are in the same direction, as shown in Figure 3c. We define that a slice is a point set containing the points firing on the same object and within the same vertical block.

Definition of Slice
Sensors 2020, 20, x FOR PEER REVIEW 6 of 25 N is defined as the amount of points contained in a slice, and if N > 1. If N = 1, this indicates that it is a noise point or a ground point of LiDAR. The slice is an invalid slice that does not have enough information for further processing. As shown in Figure 4a, the point and are the LiDAR data point in the same slice. The point and are the points of the point and mapping to the XY horizontal plane, respectively. Δ is the spatial distance between point and . Δ is the distance between point and in the XY horizontal plane. is the inclination angle of the object to the XY horizontal plane. and are the distance of point and to the ordinate origin, respectively. Their relationship can be formed as: Generally speaking, objects, such as people, vehicles, and buildings, are vertically aligned with the ground. It means that the object's inclination angle goes to 90°. Then: As shown in Figure 4b, ∆ is the minimum distance between object and object . There always have a social distance, safety distance, or viewing distance between pedestrians, vehicles, trees, and buildings. It can be deduced that: Therefore, a threshold value ∆ can be found to distinguish whether the point and belong to the same slice: Furthermore, the threshold value setting of ∆ will affect the accuracy of slice extraction within the object, particularly for the object, which is not perpendicular to the ground. For example, a too small ∆ may give rise to classification of the point and into different slices, when the point and actively belong to the same slice. On the contrary, a too large ∆ will result in the point and being classified into a slice, when the point and should be in different slices. From the experimental experiences, the threshold value of ∆ showed high accuracy and no obvious effect on the experiment results within the range of [0.1, 0.5]. In this paper, we set the threshold value as the median value of the range [0.1, 0.5], ∆ = 0.3 m. N is defined as the amount of points contained in a slice, and if N > 1. If N = 1, this indicates that it is a noise point or a ground point of LiDAR. The slice is an invalid slice that does not have enough information for further processing.
As shown in Figure 4a, the point α and β are the LiDAR data point in the same slice. The point α and β are the points of the point α and β mapping to the XY horizontal plane, respectively. ∆D β α is the spatial distance between point α and β. ∆R β α is the distance between point α and β in the XY horizontal plane. θ is the inclination angle of the object to the XY horizontal plane. R α and R β are the distance of point α and β to the ordinate origin, respectively. Their relationship can be formed as: Sensors 2020, 20, x FOR PEER REVIEW 7 of 25 (a) (b)

Definition of Abnormal and Normal Ground Points
We assumed that there is a base ground plane (BGP), which is a reference plane and absolutely smooth. The curbstones and road obstacles can be regarded as the static objects on the BGP, and a pothole is also a type of object. As shown in Figure 5, if there is a curbstone on the BGP, the LiDAR  Generally speaking, objects, such as people, vehicles, and buildings, are vertically aligned with the ground. It means that the object's inclination angle θ goes to 90 • . Then: As shown in Figure 4b, ∆D O2 O1 is the minimum distance between object O 1 and object O 2 . There always have a social distance, safety distance, or viewing distance between pedestrians, vehicles, trees, and buildings. It can be deduced that: Therefore, a threshold value ∆R th can be found to distinguish whether the point α and β belong to the same slice: ∆R Furthermore, the threshold value setting of ∆R th will affect the accuracy of slice extraction within the object, particularly for the object, which is not perpendicular to the ground. For example, a too small ∆R th may give rise to classification of the point α and β into different slices, when the point α and β actively belong to the same slice. On the contrary, a too large ∆R th will result in the point α and β being classified into a slice, when the point α and β should be in different slices. From the experimental experiences, the threshold value of ∆R th showed high accuracy and no obvious effect on the experiment results within the range of [0.1, 0.5]. In this paper, we set the threshold value as the median value of the range [0.1, 0.5], ∆R th = 0.3 m.

Definition of Abnormal and Normal Ground Points
We assumed that there is a base ground plane (BGP), which is a reference plane and absolutely smooth. The curbstones and road obstacles can be regarded as the static objects on the BGP, and a pothole is also a type of object. As shown in Figure 5, if there is a curbstone on the BGP, the LiDAR sensor obtains the laser point as A; if there is no curbstone on the BGP, the LiDAR sensor obtains the laser point as B. According the classification of the point cloud in Section 2, point A belongs to AGPs and point B belongs to NGPs.

Definition of Abnormal and Normal Ground Points
We assumed that there is a base ground plane (BGP), which is a reference plane and absolutely smooth. The curbstones and road obstacles can be regarded as the static objects on the BGP, and a pothole is also a type of object. As shown in Figure 5, if there is a curbstone on the BGP, the LiDAR sensor obtains the laser point as A; if there is no curbstone on the BGP, the LiDAR sensor obtains the laser point as B. According the classification of the point cloud in Section 2, point A belongs to AGPs and point B belongs to NGPs. In Figure 5, h1 is the height of the curbstone, r is the horizontal distance between the curbstone and LiDAR sensor, and h2 is the installation height of the LiDAR. ∆ is the horizontal distance between point A and B in the BGP. According to Equation (2), ∆ = Δ = Δ × cos( ).

Definition of Valuable and Worthless Object Points
In general, no matter the moving road users or the fixed tree trunks, they are all connected to the ground. We call them valuable objects. Road users are the main body of traffic, and other valuable objects will affect the trajectory of road users. Slices connected with or close to the ground are called valuable slices, and others are worthless slices. The corresponding points in slices are VOPs and WOPs. We define an area from the BGP to a certain height as a key region, as shown in Figure 6. The VOPs in the key region are called key valuable points. In Figure 5, h1 is the height of the curbstone, r is the horizontal distance between the curbstone and LiDAR sensor, and h2 is the installation height of the LiDAR. ∆r is the horizontal distance between point A and B in the BGP. According to Equation (2), ∆r = ∆R B A = ∆D B A × cos(ϕ).

Definition of Valuable and Worthless Object Points
In general, no matter the moving road users or the fixed tree trunks, they are all connected to the ground. We call them valuable objects. Road users are the main body of traffic, and other valuable objects will affect the trajectory of road users. Slices connected with or close to the ground are called valuable slices, and others are worthless slices. The corresponding points in slices are VOPs and WOPs. We define an area from the BGP to a certain height as a key region, as shown in Figure 6. The VOPs in the key region are called key valuable points. In Figure 6a, the orange area is the key region, the red points are the key valuable points, and the blue points are other types of points. In Figure 6b, the red slices are the valuable slices and the gray slices are worthless slices. Most of the points in the key region are key valuable points, which are denser than other points in the point cloud within the key region. Based on this assumption, we can extract all the key valuable points. Then, all the object points can be classified into VOPs and WOPs.

Slice-Based Projection Filtering Algorithm
In this section, we propose the detailed steps of the slice-based projection filtering (SPF) algorithm, which makes up four major parts: Projection, extraction of valid slices, AGPs and NGPs separation, and VOPs and WOPs separation. Figure 7 shows the flowchart of the proposed method.

Projection
Compared with the 2-D plane, background filtering using 3-D space coordinates will take advantage of the rich characteristic information of the object. However, it will increase the In Figure 6a, the orange area is the key region, the red points are the key valuable points, and the blue points are other types of points. In Figure 6b, the red slices are the valuable slices and the gray slices are worthless slices. Most of the points in the key region are key valuable points, which are denser than other points in the point cloud within the key region. Based on this assumption, we can extract all the key valuable points. Then, all the object points can be classified into VOPs and WOPs.

Slice-Based Projection Filtering Algorithm
In this section, we propose the detailed steps of the slice-based projection filtering (SPF) algorithm, which makes up four major parts: Projection, extraction of valid slices, AGPs and NGPs separation, and VOPs and WOPs separation. Figure 7 shows the flowchart of the proposed method. In Figure 6a, the orange area is the key region, the red points are the key valuable points, and the blue points are other types of points. In Figure 6b, the red slices are the valuable slices and the gray slices are worthless slices. Most of the points in the key region are key valuable points, which are denser than other points in the point cloud within the key region. Based on this assumption, we can extract all the key valuable points. Then, all the object points can be classified into VOPs and WOPs.

Slice-Based Projection Filtering Algorithm
In this section, we propose the detailed steps of the slice-based projection filtering (SPF) algorithm, which makes up four major parts: Projection, extraction of valid slices, AGPs and NGPs separation, and VOPs and WOPs separation. Figure 7 shows the flowchart of the proposed method.

Projection
Compared with the 2-D plane, background filtering using 3-D space coordinates will take advantage of the rich characteristic information of the object. However, it will increase the

Projection
Compared with the 2-D plane, background filtering using 3-D space coordinates will take advantage of the rich characteristic information of the object. However, it will increase the calculation dimensions and computation complexity in dealing with laser points. For example, in 3-D-DSF background filtering, Wu [25] divided the 3-D space into continuous small cubes with the cube size at 0.1 m × 0.1 m × 0.1 m. For a 40 m × 40 m × 8 m scanning range, it needed to be divided into 400 × 400 × 80 cubes. There is only 400 × 400 cells to be analyzed if it taken into the 2-D plane. The loop count of computation with the 2-D plane is 1/80 of the 3-D space method. Therefore, we attempted to project the 3-D point cloud into the 2-D plane. The 2-D point cloud data can be represented by rectangular coordinates (x, y) and polar coordinates (r, α), as shown in Figure 2. These two coordinate systems take the LiDAR sensor as the origin of coordinates. The data point can be expressed as: where x and y are the coordinate position within the X-axis and Y-axis, respectively, in rectangular coordinates; r is the distance between the point and LiDAR in the XY horizontal plane; and α is the horizontal rotation angle. The process of projection in rectangular coordinates includes the multiplication of four trigonometric functions, while the projection in the polar coordinate system only requires the multiplication of one trigonometric function. Therefore, projection in polar coordinates can effectively reduce the calculation time and computation resources. Furthermore, the scanning trajectory will be a closed circle or ellipse after being projected when a laser beam in LiDAR is scanned on the ground plane.
There will be two y values on a certain x value when rectangular coordinates are used to represent this trajectory curve. If polar coordinates are used, r and α are a one-by-one correspondence. Hence, polar coordinates can be used to represent the trajectory curve with a simpler equation, which is more in line with the characteristics of the LiDAR data. In this paper, the 3-D point cloud is projected in a polar coordinate plane for analysis, as shown in Figure 8.  Figure 2. These two coordinate systems take the LiDAR sensor as the origin of coordinates. The data point can be expressed as: where and are the coordinate position within the X-axis and Y-axis, respectively, in rectangular coordinates; is the distance between the point and LiDAR in the XY horizontal plane; and is the horizontal rotation angle. The process of projection in rectangular coordinates includes the multiplication of four trigonometric functions, while the projection in the polar coordinate system only requires the multiplication of one trigonometric function. Therefore, projection in polar coordinates can effectively reduce the calculation time and computation resources. Furthermore, the scanning trajectory will be a closed circle or ellipse after being projected when a laser beam in LiDAR is scanned on the ground plane.
There will be two values on a certain value when rectangular coordinates are used to represent this trajectory curve. If polar coordinates are used, and are a one-by-one correspondence. Hence, polar coordinates can be used to represent the trajectory curve with a simpler equation, which is more in line with the characteristics of the LiDAR data. In this paper, the 3-D point cloud is projected in a polar coordinate plane for analysis, as shown in Figure 8.   (Figure 8b). We find that the height information of the point cloud is hidden after the projection. However, this will not affect the processing results, as the purpose of the projection in this paper is to extract the valid slices in the point cloud.   (Figure 8b). We find that the height information of the point cloud is hidden after the projection. However, this will not affect the processing results, as the purpose of the projection in this paper is to extract the valid slices in the point cloud.

Valid Slices Extraction
The purpose of the valid slice extraction operation is to find the object points in the point cloud based on the properties of the slice. We used the matrix C 0 = [D, W, A] to represent the raw data of the LiDAR rotation output, where D, W, and A represent the column vectors that contain d, ω, and α, respectively, and each row corresponds to a point collected by the LiDAR sensor. We used matrix C 1 = [R, W, A] to represent the matrix after projecting matrix C 0 . We used C 1 [ f lag] to obtain all row vectors in the C 1 satisfying the condition of the flag.
Suppose P x = (r x , ω x , α x ) is a certain row of C 1 and we need to find other points in the same slice with P x , it is the simplest way to traverse all points in C 1 [α = α x ] and calculate the distance from P x . However, there are many points in the point cloud, and it will take a long time to traverse all points, so the implementation efficiency of this method is very low. For this reason, this paper proposed a solution based on matrix operations to extract effective slices efficiently, as shown in Algorithm 1.
Get the difference value between two elements in R sequences, obtain delta 4: According to the definition and attributes of valid slice, the valid slices can be extracted based on two conditions: Condition 1: If points belong to the same slice, these points have the same horizontal rotation angle. Condition 2: The horizontal distance of two points ∆R β α in the same slice is less than the threshold value of ∆R th . Therefore, we suppose P1, P2, and P3 are three points that are all contained in C 1 . If P1 and P3 belong to the same slice, P1 and P2 also belong to the same slice, then it is clear that P2 and P3 also belong to the same slice. Thus, we only need to determine the distance between two adjacent laser points within the setting horizontal rotation angle.
We can obtain the matrix where two adjacent rows are also adjacent in real space after sorting the rows of matrix C 1 by A and then by W. Then, it is easy to get the distance difference of all adjacent points to LiDAR. Finally, we can reasonably set the threshold value of ∆R th to extract all points in the point cloud that belong to the valid slices. In algorithm 1, matrix F1 is used to mark all points, where true indicates that the point belongs to a certain slice and false indicates that the point is a single point.
There are one or more single points between the different slices. Based on the differences, the slice can be differentiated from and labeled. In the proposed algorithm, we used the matrix F2 to mark different slices. Within the matrix F2, 0 indicates that the point is a single point; if the number is greater than 0, it indicates the index of the slice that the point belongs to.
By the experiments, the runtime of the proposed algorithm is about 18 ms compared to the traversing method, which needed 1.5 s to extract the valid slice. Therefore, the proposed algorithm of this paper is more efficient.

Abnormal and Normal Ground Points Separation
The surface, slope, curbstones, road obstacles, and other information of the ground are also very useful for trajectory predictions of road users with the consideration of the edge of different function areas of the road. The purpose of the AGPs and NGPs separation operation is to abstract all the ground points that could impact the road users.

Ground Point Extraction
After extracting the valid slices from the point cloud, the remaining points are a single point set. The single point set is composed of the ground points, noise points, and edge points of objects. In the LiDAR scanning range, the ground points are closer to the BGP. Therefore, we can obtain the ground points by setting the condition as follows: where h X is the height of point X; Z X is the coordinates of the point X on the Z-axis, where the horizontal plane is the BGP; and H gp is the threshold of ground points. In this paper, the threshold value H gp was based on the height between the point cloud and BGP. We divided 0-1 m into 100 parts on average and counted the amount of points in each part. The height distribution of the point cloud is shown in Figure 9.
As shown in Figure 9, there is a critical point (0.26, 40) in Figure 9. The point cloud with 0-1 m is divided into two parts at the critical point (0.26, 40). The height of most ground points is less than 0.26 m. The height of most objects' points is higher than 0.26 m. Thus, we set H gp = 0.26 in this paper.
where ℎ is the height of point ; is the coordinates of the point on the -axis, where the horizontal plane is the BGP; and is the threshold of ground points. In this paper, the threshold value was based on the height between the point cloud and BGP. We divided 0-1 m into 100 parts on average and counted the amount of points in each part. The height distribution of the point cloud is shown in Figure 9.

Separation Method
In the absence of AGPs, the projected NGPs should be a closed curve track, noted as L1. When the projection plane of the LiDAR point cloud is parallel to the ground, the trajectory is a circle. The trajectory is a curve similar to an ellipse when LiDAR has an angle from the ground. When there is an AGP at a horizontal angle α, the point should deviate from the curve L1. Therefore, we can fit the trajectory curve L1, and use it to decide whether a point P in the ground points belongs to the NGPs or AGPs according to the distance from the curve L1. The judgement condition is: where R P is the actual distance between the projected point P and LiDAR; R Predict |α = α P is the distance predicted by the curve L1 at the horizontal roatation angle α P ; Abs( * ) is the absolute function; and P th is the threshold value to differentiate AGPs and NGPs; the range of P th is: where K is a coefficient, and K = h1/(h2 − h1); h1 is the height of the point; and h2 is the installation height of the LiDAR. We selected the data of laser channel 0 with the smallest vertical angle to train the trajectory curve L1, as the points where the laser of this channel fell on BGP were the closest to the LiDAR, and contained the most ground information. In addition, we selected a frame of the point cloud including more ground information, so as to fit the curve more accurately. After fitting the curve L1, the curves can be applied to the point cloud data of all other frames.
The α and r of all points obtained by a single channel laser beam in a scanning period are one-by-one corresponding. As shown in Figure 10, with α taken as the horizontal axis and r as the vertical axis, the processes of the fitting curve can be described as the following steps: Step 1: We used a quadratic equation to fit the ground points detected by the LiDAR laser channel 0 by least squares, firstly: f it f unction1 = λ ·r 2 + b· r + c, where r is the distance between the point and LiDAR in the XY horizontal plane as mentioned above; and λ, b, and c are regression coefficients, which can be obtained by least squares regression. As shown in Figure 10a, the green scatter is the ground point detected by the laser channel 0, noted as GP0. The orange curve is the fitting curve of the ground points, noted as L0. It can be seen that there are some points far below the fitting curve L0. These are noise points or AGPs in GP0.
Step 2: Base on the fitting curve L0 fitted by the laser channel 0, we used L0 to filter the points detected by LiDAR laser channel 1 to channel 15. With the same horizontal angle α, if the points detected by the LiDAR laser channel 1 to channel 15 are far below the corresponding point in L0, these are removed before stepping into the next step.
Step 3: After filtering out all the noise points or AGPs in step 1 and step 2, the remaining points are the NGPs. We used the remaining points to fit the final curve L1 by the Sinusoid function as: where α is the horizontal rotation angle of the point as mentioned above; and A, δ, ϑ, and γ are the regression coefficients, which can be obtained by least squares regression.
where is the actual distance between the projected point and LiDAR; | = is the distance predicted by the curve L1 at the horizontal roatation angle ; ( * ) is the absolute function; and Ρ is the threshold value to differentiate AGPs and NGPs; the range of Ρ is: where is a coefficient, and = ℎ1/(ℎ2 − ℎ1); h1 is the height of the point; and h2 is the installation height of the LiDAR. We selected the data of laser channel 0 with the smallest vertical angle to train the trajectory curve L1, as the points where the laser of this channel fell on BGP were the closest to the LiDAR, and contained the most ground information. In addition, we selected a frame of the point cloud including more ground information, so as to fit the curve more accurately. After fitting the curve L1, the curves can be applied to the point cloud data of all other frames.
The and of all points obtained by a single channel laser beam in a scanning period are one-by-one corresponding. As shown in Figure 10, with taken as the horizontal axis and as the vertical axis, the processes of the fitting curve can be described as the following steps: Step 1: We used a quadratic equation to fit the ground points detected by the LiDAR laser channel 0 by least squares, firstly: where is the distance between the point and LiDAR in the XY horizontal plane as mentioned above; and , , and are regression coefficients, which can be obtained by least squares regression.
As shown in Figure 10a, the green scatter is the ground point detected by the laser channel 0, noted as GP0. The orange curve is the fitting curve of the ground points, noted as L0. It can be seen that there are some points far below the fitting curve L0. These are noise points or AGPs in GP0. As shown in Figure 10b, the green scatters are the NGPs. The orange curve is the final fitting curve L1. It can be seen that some blue scatters are far below the fitting curve L1. These are noise points or AGPs within the raw data. The remaining points nearly around the fitting curve L1 map well with the NGPs. Most of the blue points nearly around the fitting curve L1 overlap with the corresponding NGPs. It means that the fitting curve L1 can effectively distinguish the AGPs and NGPs.

Separation Results
After obtaining the trajectory equations of all laser channels, we can distinguish AGPs and NGPs in all ground points by setting the threshold value P th reasonably. In practical applications, the values of h1 and h2 are fixed. According Equation (9), the max(P th ) will change as the actual distance R P is changing. As to distinguish the AGPs and NGPs accurately, the threshold value P th should change as the as the actual distance R P changes. In the experiment, we set h1 equal to 0.26 m based on the height of the curbstones and h2 equal to 1.8 m based on the installation height of the LiDAR. According K = h1/(h2 − h1) in Equation (9), K ≈ 0.1688. Thus, 0 < P th < R P × 0.1688. The experiments results showed that a too small threshold value P th will give rise to recognition of the NGPs as AGPs when the points are near the LiDAR sensor. On the contrary, a too large threshold value P th will result in the AGPs being identified as NGPs when the points are far away from the LiDAR sensor. In this paper, we used the median value of the range (0, R P × 0.1688), setting P th = R P × 0.0844. The experiment result is shown in Figure 11.
As shown in Figure 11a, the blue scatters are NGPs, and the orange scatters are AGPs. It can be effectively distinguished by the proposed method. If there are a large number of AGPs and a small number of NGPs in a certain area, it can be deduced that this area is a non-driving area due to too many obstacles within this area. The gray areas are non-driving areas as shown in Figure 11a. Figure 11b shows the experimental location in the Baidu Map and the field scene scanned by the LiDAR sensor in the experiment. the LiDAR. According K = ℎ1/(ℎ2 − ℎ1) in Equation (9), K ≈ 0.1688 . Thus, 0 < Ρ < × 0.1688. The experiments results showed that a too small threshold value Ρ will give rise to recognition of the NGPs as AGPs when the points are near the LiDAR sensor. On the contrary, a too large threshold value Ρ will result in the AGPs being identified as NGPs when the points are far away from the LiDAR sensor. In this paper, we used the median value of the range (0, × 0.1688), setting Ρ = × 0.0844. The experiment result is shown in Figure 11. As shown in Figure 11a, the blue scatters are NGPs, and the orange scatters are AGPs. It can be effectively distinguished by the proposed method. If there are a large number of AGPs and a small

Valuable and Worthless Object Points Separation
We can obtain the object points in the point cloud after extracting the valid slice operation. However, not all object points are useful for traffic behavior analysis. Therefore, the separating VOPs and WOPs operation is required to find the points from all object points, which are useful for traffic behavior analysis.

Key Region Optimization
If the projection plane of the LiDAR point cloud is parallel to the BGP, the Z information of the point cloud can be directly used to obtain the key region. In fact, there may be a slope in the BGP. In this situation, it is hard to ensure the projection plane remains parallel to the BGP. Therefore, if there is an inclination angle between the BGP and the projection plane, it will be difficult to obtain the accurate key region by using only the Z information of the point cloud. The error increases with the distance from the LiDAR sensor, as shown in Figure 12. number of NGPs in a certain area, it can be deduced that this area is a non-driving area due to too many obstacles within this area. The gray areas are non-driving areas as shown in Figure 11a. Figure 11b shows the experimental location in the Baidu Map and the field scene scanned by the LiDAR sensor in the experiment.

Valuable and Worthless Object Points Separation
We can obtain the object points in the point cloud after extracting the valid slice operation. However, not all object points are useful for traffic behavior analysis. Therefore, the separating VOPs and WOPs operation is required to find the points from all object points, which are useful for traffic behavior analysis.

Key Region Optimization
If the projection plane of the LiDAR point cloud is parallel to the BGP, the Z information of the point cloud can be directly used to obtain the key region. In fact, there may be a slope in the BGP. In this situation, it is hard to ensure the projection plane remains parallel to the BGP. Therefore, if there is an inclination angle between the BGP and the projection plane, it will be difficult to obtain the accurate key region by using only the Z information of the point cloud. The error increases with the distance from the LiDAR sensor, as shown in Figure 12. As shown in Figure 12a, point and point are formed by LiDAR laser beams falling on the BGP. The projected points are and , respectively. is the inclination angle between BGP and the projection plane. When another laser beam falls on point on the BGP, the Z-coordinate of point should be: where , , and are the horizontal distance from the point , , and to the LiDAR As shown in Figure 12a, point A and point B are formed by LiDAR laser beams falling on the BGP. The projected points are A and B , respectively. ∂ is the inclination angle between BGP and the projection plane. When another laser beam falls on point X on the BGP, the Z-coordinate of point X should be: where R A , R B , and R X are the horizontal distance from the point A, B, and X to the LiDAR sensor, respectively; and Z A , Z B , and Z X are the distance deviation in calculation the Z-coordinate of point A, B, and X, respectively. Clearly, if there is an inclination angle between the BGP and the projection plane, the ∂ is changing for different points within the point cloud, particularly for the object points. It requires at least two channel laser points to calculate the inclination angle between the BGP and the projection plane.
According to the LiDAR user manual and programming guide [29], the laser ID 0 and 2 are the two laser beam closer to the ground. In the experiment, these two laser beams generate more ground points than the other laser beam in LiDAR. In this paper, the NGPs generated by laser ID 0 and 2 are used to calculate the inclination angle ∂ between BGP and the projection plane. The NGPs are obtained based on the process results of Section 4.3, as shown in Figure 12b.
In the LiDAR data packet, we can obtain the XYZ coordinates of the point and its horizontal rotation angle. Here, A(x A , y A , z A ) is the coordinate information of point A obtained by the laser ID 0, B(x B , y B , z B ) is the coordinate information of point B obtained by the laser ID 2, and point A and point B have the same horizontal rotation angle. According to Equation (12) and Figure 12a, the inclination angle ∂ between BGP and the projection plane can be calculated as: All the point pairs in Figure 12a were used to calculate tan(∂), and the calculation results are shown in Figure 13a, where the green scatters is the value of tan(∂) at different horizontal rotation angles (0-2π).
As shown in Figure 13a, there is noise in the calculation results, which is caused by the errors of the sensor itself and the uneven ground. The continuous curve of the tan(∂) value cannot be obtained directly as the points used to calculate it are discrete, or there is a shied on the road.
To improve the flexibility of the algorithm, it is necessary to fully consider the inclination angle information between the BGP and the projection plane. In addition, the road surface condition also should be considered. To adapt to different road surface conditions and employ inclination, we attempted to use a machine learning algorithm to obtain accurately the inclination angle between the BGP and the projection plane and identify the key region. A machine learning model based on the backpropagation artificial neuron network (BP-ANN) was developed to fit the scatter of thte calculated tan(∂) to obtain an accurate and continuous actual tan(∂).
The BP-ANN is a multilayer feedforward neural network composed of an input layer, a hidden layer, and an output layer, as shown in Figure 13b. In practical applications, the road pavement conditions within the scanning range of LiDAR are different. Thus, in the different LiDAR scanning horizontal rotational angles α, the inclination angle ∂ may be different. In the BP-ANN, we input the horizontal rotational angle α to get the corresponding road pavement inclination ∂ with the output of tan(∂) by the BP-ANN-based machine learning model. It means that tan(∂) = Func(α) and Func( * ) is obtained by the BP-ANN model. The training data set of α and tan(∂) were the points on the BGP with the same horizontal rotational angle α, using Equation (12) or Equation (13) to calculate the corresponding tan(∂). Therefore, the input layer and output layer contained only one neuron each as only one variable was needed for each layer. The variables α and tan(∂) were used for the input layer and the output layer, respectively. The number of neurons of the hidden layer was set as 32. We used the activation function for neurons in the hidden layer for the function that could converge faster than the sigmoid function, while the neuron in the output layer did not use any activation function to avoid gradient explosion. During the model training period, the mean square error function was used for the loss function, the optimizer was Adam, and the batch size was set as 32. The fitting curve of tan(∂) is shown as the red line in Figure 13a.
In the LiDAR data packet, we can obtain the XYZ coordinates of the point and its horizontal rotation angle. Here, A( , , ) is the coordinate information of point A obtained by the laser ID 0, B( , , ) is the coordinate information of point B obtained by the laser ID 2, and point A and point B have the same horizontal rotation angle. According to Equation (12) and Figure 12a, the inclination angle between BGP and the projection plane can be calculated as: All the point pairs in Figure 12a were used to calculate tan( ), and the calculation results are shown in Figure 13a, where the green scatters is the value of tan( ) at different horizontal rotation angles (0-2π).

X Input Layer
Hidden Layer Output Layer As shown in Figure 13a, there is noise in the calculation results, which is caused by the errors of the sensor itself and the uneven ground. The continuous curve of the tan( ) value cannot be obtained directly as the points used to calculate it are discrete, or there is a shied on the road.
To improve the flexibility of the algorithm, it is necessary to fully consider the inclination angle information between the BGP and the projection plane. In addition, the road surface condition also should be considered. To adapt to different road surface conditions and employ inclination, we In addition, we also used three other algorithms to fit the tan(∂) curve. A machine learning model based on long short-term memory (LSTM) was used to fit the tan(∂) curve [31]. The complete network model is shown in Figure 13c, in which each LSTM cell contains 15 neurons. The time-step of the model was set to 2, due to the tan(∂) being predicted by keeping the memory of several previous position points when both the LiDAR and the BGP are fixed. Then, the memory of the changing tan(∂) was transferred between two adjacent batches until the internal state of the model was reset. During the model training period, the mean square error function was used for the loss function, the optimizer was Adam, and the batch size was set as 8. After the completion of the model, we set the time-step to 1 for the prediction, and the obtained continuous tan(∂) curve is shown as the blue line in Figure 13a.
The second method was the Ridge regression [32], where the coefficient of the penalty term was set as 0.5. The fitting curve of tan(∂) is shown as the yellow line in Figure 13a. The last method was support vector regression (SVR) [33], where the radial basis was used for the kernel function. The fitting result of the tan(∂) curve is shown as the orange line in Figure 13a.
From the fitting results, we can see that the fitting results based on LSTM are the best in the range of [1.5π, 2π]. The fitting curve of tan(∂) based on BP-ANN is similar to LSTM. The loss curve generated by the mean square error function with the LSTM model dropped faster than the BP-ANN model. However, the loss values were almost the same after 100 epochs, as shown in Figure 13d. In practical applications, the LSTM-based machine learning model is much more complex than the BP-ANN machine learning model. Therefore, we chose the simpler BP-ANN machine learning model to obtain the accurate inclination angle between the BGP and the projection plane.
After identification, there is an inclination angle between the BGP and the projection plane, and after obtaining the accurate tan(∂) in the corresponding horizontal rotation angle, it is easy to identify a point within the point cloud that belongs to the key region by: where z is the Z-coordinate output by the LiDAR sensor; H LiDAR is the installation height of the LiDAR sensor; H KR is the height of the key region; and Z X is the distance deviation, which can be calculated according to Equation (12). In China, the average height of male adults is about 1.71 m, and it is about 1.6 m for female adults [34]. The height of small cars is generally within 1.6 m and 1.8 m. Road users are the main research object in the field of transportation. Therefore, in the experiments, we set the value of H KR as 1.5 m, slightly lower than the general height of road users. Within this height range, most road users can be included, and most roadside leaves can be excluded.

Valuable Object Points Extraction
The key valuable points are more point intensive than other types of points after projecting all points in the key region. Therefore, the DBSCAN (density-based spatial clustering of applications with noise) [35] algorithm can be used to extract the key value points from the key region. With DBSCAN, one object can be joined in the clustering area if the point density is higher than or equal to a pre-defined minimum number of points (MinPts), and the distance to the center is less than a pre-defined epsilon (ε). In the experiments, we set Minpts = 3 and ε = 0.008 to cluster the projected points in the key region. The key value point cluster results are shown as Figure 14a.
Sensors 2020, 20, x FOR PEER REVIEW 17 of 25 as 1.5 m, slightly lower than the general height of road users. Within this height range, most road users can be included, and most roadside leaves can be excluded.

Valuable Object Points Extraction
The key valuable points are more point intensive than other types of points after projecting all points in the key region. Therefore, the DBSCAN (density-based spatial clustering of applications with noise) [35] algorithm can be used to extract the key value points from the key region. With DBSCAN, one object can be joined in the clustering area if the point density is higher than or equal to a pre-defined minimum number of points ( ), and the distance to the center is less than a pre-defined epsilon ( ). In the experiments, we set = 3 and = 0.008 to cluster the projected points in the key region. The key value point cluster results are shown as Figure 14a. As shown in Figure 14a, the green scatters are the key valuable points obtained from the point cloud in Figure 8. The red scatters are sparse points that belong to WOPs. According to Algorithm 1, we can obtain the matrix 2, which contains the index of all valid slices. Therefore, we can find all the index numbers contained in the key valuable points and then find all the VOPs by combining with the matrix 2. Figure 14b is generated by the NGPs that were obtained based on the process methods of Section 4.3, and has the same raw point cloud of Figure 14a. According the fitting curve of tan( ), we can infer that the height of the left side is higher than the right side in Figure 14. If only the Z coordinate is used for the identification of the key region, we will obtain less area for the left side and more in the right side. Table 1 shows the numbers of VOPs obtained in the left and right sides of Figure 14, with or without the fitting of the inclination angle between the BGP and the projection plane. Table 1. The number of VOP extraction with or without fitting of the inclination angle between the BGP and the projection plane. As shown in Figure 14a, the green scatters are the key valuable points obtained from the point cloud in Figure 8. The red scatters are sparse points that belong to WOPs. According to Algorithm 1, we can obtain the matrix F2, which contains the index of all valid slices. Therefore, we can find all the index numbers contained in the key valuable points and then find all the VOPs by combining with the matrix F2. Figure 14b is generated by the NGPs that were obtained based on the process methods of Section 4.3, and has the same raw point cloud of Figure 14a. According the fitting curve of tan(∂), we can infer that the height of the left side is higher than the right side in Figure 14. If only the Z coordinate is used for the identification of the key region, we will obtain less area for the left side and more in the right side. Table 1 shows the numbers of VOPs obtained in the left and right sides of Figure 14, with or without the fitting of the inclination angle between the BGP and the projection plane. Table 1. The number of VOP extraction with or without fitting of the inclination angle between the BGP and the projection plane.

Fitting the Inclination Angle
The As shown in Table 1, on the left side of the key region, the number of VOPs obtained after fitting the inclination angle was more than that without fitting, and it was an adverse result on the right side. This phenomenon demonstrates that greater accuracy in the key region can be found after fitting the angle. In addition, the number of VOPs with fitting was more than that without fitting, which shows that the greater the accuracy in the key region, the more VOPs that will be obtained. Figure 15 shows the AGPs and VOPs extraction results with or without the fitting of the inclination angle between the BGP and the projection plane.
Sensors 2020, 20, x FOR PEER REVIEW 18 of 25 As shown in Table 1, on the left side of the key region, the number of VOPs obtained after fitting the inclination angle was more than that without fitting, and it was an adverse result on the right side. This phenomenon demonstrates that greater accuracy in the key region can be found after fitting the angle. In addition, the number of VOPs with fitting was more than that without fitting, which shows that the greater the accuracy in the key region, the more VOPs that will be obtained. Figure 15 shows the AGPs and VOPs extraction results with or without the fitting of the inclination angle between the BGP and the projection plane. As shown in Figure 15a, the VOPs within the green ellipse region are missing in Figure 15b. Furthermore, the inclination angle between the BGP and the projection plane not only affects the extraction result of VOPs but also affects the extraction result of AGPs. As shown in Figure 15, the AGPs within the orange ellipse regions in Figure 15a are missing in Figure 15b, which does not have fitting of the inclination angle. This is due to the left side being higher than the right side as discussed above in Figure 14, and a large number of ground points are far from LiDAR on the left side being identified as individual single points, which belong to the objects in Figure 15b. In addition, as shown in Figure 15b, some single points within the red ellipse region, which belong to the objects on the right side, are identified as the AGPs. Therefore, it is necessary to fit the inclination angle to improve the accuracy in the identification of the category of points. As shown in Figure 15a, the VOPs within the green ellipse region are missing in Figure 15b. Furthermore, the inclination angle between the BGP and the projection plane not only affects the extraction result of VOPs but also affects the extraction result of AGPs. As shown in Figure 15, the AGPs within the orange ellipse regions in Figure 15a are missing in Figure 15b, which does not have fitting of the inclination angle. This is due to the left side being higher than the right side as discussed above in Figure 14, and a large number of ground points are far from LiDAR on the left side being identified as individual single points, which belong to the objects in Figure 15b. In addition, as shown in Figure 15b, some single points within the red ellipse region, which belong to the objects on the right side, are identified as the AGPs. Therefore, it is necessary to fit the inclination angle to improve the accuracy in the identification of the category of points.

Filtering Results
The first experiment was designed to illustrate the obtained filtering results. We evaluated the performance of the proposed method of this paper and compared it to a popular RANSAC algorithm by Behley [23] and the 3-D-DSF algorithm provided by Wu [25]. For the 3-D-DSF approach, the side length of the cube was 0.1 m.
The point cloud data were collected by the Velodyne VLP-16 LiDAR sensor at Central Avenue, Nanling campus, Jilin University. The location of the LiDAR sensor is latitude 43.85 n, longitude 125.33 w. Approximately 1200 frames of point cloud were collected and used as input data in the experiments. The background point filtering results of the three algorithms are shown in Figure 16. Figure 16a shows the raw point cloud before processing. Figure 16b shows that the RANSAC algorithm can only filter the ground points from the raw point cloud. More processing stages are required to remove the worthless points from the point cloud by other methods. As shown in Figure 16c, there are still some noise points left after the background point filtering. We needed to remove noise points using clustering methods, such as k-means and DBSCAN, which increases the execution time of the algorithm [36]. Furthermore, the 3-D-DSF algorithm ignores certain useful information generated by fixed objects, such as curbstones and tree trunks, which can affect the trajectory of a road user in summer in tropical areas. Figure 16d shows the filtering result by the proposed algorithm SPF. Different colors in the figure present different objects clustered by DBSCAN. From the perspective of filtering quality, the proposed algorithm was able to effectively extract useful points in the point cloud, including the AGPs that belong to the ground and the VOPs that belong to the objects. Through AGPs and VOPs, we can not only obtain the information of all road users in the point cloud but also the environment information that can affect the trajectory of the road users. Therefore, analyzing the data after filtering by the proposed algorithm can help to perceive the traffic information more comprehensively, accurately, and quickly.
To compare the filtering results in more detail, we used the amount of the point cloud before and after filtering to evaluate the algorithms. The experimental results are compared in Table 2 Figure 16a shows the raw point cloud before processing. Figure 16b shows that the RANSAC algorithm can only filter the ground points from the raw point cloud. More processing stages are required to remove the worthless points from the point cloud by other methods. As shown in Figure  16c, there are still some noise points left after the background point filtering. We needed to remove noise points using clustering methods, such as k-means and DBSCAN, which increases the execution time of the algorithm [36]. Furthermore, the 3-D-DSF algorithm ignores certain useful information generated by fixed objects, such as curbstones and tree trunks, which can affect the trajectory of a road user in summer in tropical areas. Figure 16d shows the filtering result by the proposed algorithm SPF. Different colors in the figure present different objects clustered by DBSCAN. From the perspective of filtering quality, the As shown in Table 2, the 3-D-DSF algorithm can filter out 96.7% of the points from the point cloud. The VOPs and AGPs excluded percentages were 84.5% and 96.8%, respectively, as shown in Figures 17c and 18c, respectively. However, the 3-D-DSF algorithm filtered out too much valuable information. Thus, the amount of valuable points (VOPs and AGPs) after background point filtering was only 956, which was the minimum. However, it still contained some noise after the background point filtering by 3-D-DSF. The percentage of noise point in the valuable points was about 4.29%.
filtering by the proposed algorithm can help to perceive the traffic information more comprehensively, accurately, and quickly.
To compare the filtering results in more detail, we used the amount of the point cloud before and after filtering to evaluate the algorithms. The experimental results are compared in Table 2, and Figures 17 and 18. The raw point cloud, raw VOPs, and raw AGPs in the experiments are shown in Figure 16a, 17a, and 18a, respectively.  The background point filtering percentage of RANSAC was only 15.9% in the raw point cloud filtering. The VOPs and AGPs excluded percentages were 11.7% and 7.0%, respectively, as shown in Figures 17b and 18b, respectively. This showed that the RANSAC algorithm retained too much worthless information. This is why the amount of valuable points after background point filtering was only 6321, which was the maximum. However, it needed to further exclude about 74% of the worthless points after the raw point cloud filtering. The percentage of valuable point was only 26% after background point filtering by RANSAC, which was the minimum. As shown in Table 2, the 3-D-DSF algorithm can filter out 96.7% of the points from the point cloud. The VOPs and AGPs excluded percentages were 84.5% and 96.8%, respectively, as shown in Figure 17c and 18c, respectively. However, the 3-D-DSF algorithm filtered out too much valuable information. Thus, the amount of valuable points (VOPs and AGPs) after background point filtering was only 956, which was the minimum. However, it still contained some noise after the background point filtering by 3-D-DSF. The percentage of noise point in the valuable points was about 4.29%.
The background point filtering percentage of RANSAC was only 15.9% in the raw point cloud filtering. The VOPs and AGPs excluded percentages were 11.7% and 7.0%, respectively, as shown in Figure 17b and 18b, respectively. This showed that the RANSAC algorithm retained too much worthless information. This is why the amount of valuable points after background point filtering was only 6321, which was the maximum. However, it needed to further exclude about 74% of the worthless points after the raw point cloud filtering. The percentage of valuable point was only 26% after background point filtering by RANSAC, which was the minimum. The SPF algorithm, the proposed method in this paper, filtered out 80.1% of the points from the raw point cloud. The VOPs and AGPs excluded percentages were 20.0% and 14.8%, respectively, as shown in Figures 17d and 18d, respectively. It could exclude all of the worthless points in the raw point cloud background filtering. The percentage of valuable point was almost 100% after background point filtering by SPF in the raw point cloud filtering, which was the maximum.
In terms of the filtering quality and accuracy for further processing in identifying the traffic object, the SPF algorithm was able to retain more valuable information and remove more worthless information for further processing in the stage of traffic object identification and trajectory tracking. In this respect, we concluded that the SPF algorithm performed better than RANSAC and 3-D-DSF.

Runtime
The execution time of the algorithm plays an important role in practical applications, especially in area traffic behavior and safety analysis. Therefore, the experiment was designed to evaluate the running performance of the three algorithms. The three algorithms were implemented in Python3.5 and were deployed on a desktop computer equipped with an Intel Core i7-7500U CPU (2.70 GHz), with 8 GB of RAM, and with a Windows 10 operating system. The raw point cloud was collected by the Velodyne VLP-16 LiDAR sensor located at Central Avenue, Nanling campus, Jilin University. The raw point cloud including about 1200 data frames was used as the data input of the three algorithms. The performance results are shown in Table 3. As shown in Table 3, the average background point filtering time was about 51 ms per frame by the SPF algorithm. The running times of RANSAC and 3-D-DSF were 198 and 113 ms per frame, respectively. In addition, the 3-D-DSF algorithm needed multi-frames to perform the algorithm. It would not be suitable for on-line intelligent transportation applications.

Conclusions and Discussion
A LiDAR sensor has advantages in obtaining the high-resolution fine-grained position, direction, speed, and trajectory data of vehicles and pedestrians. However, the point cloud generated by a LiDAR scan contains not only road user points but also buildings, trees, and ground points. To speed up the micro-level traffic data extraction from the LiDAR sensor, removing the background points using a suitable background filtering algorithm is necessary. In this paper, we developed a novel method of background point filtering for low-channel infrastructure-based LiDAR.
First, we projected a 3-D point cloud into the 2-D polar coordinate points for filtering to reduce the data processing dimensions and computing time. Secondly, we divided the points in the point cloud into four categories: AGPs, NGPs, VOPs, and WOPs. Then, we used the slice as a unit to analyze the point attributes, and to separate and extract the VOPs and AGPs from a single frame of the point cloud. To dynamically fit different road conditions and situations, such as different road gradients and different LiDAR scanning angles, a machine learning method based on the artificial neural network was used to improve the adaptability of the algorithm in complex environments.
Finally, compassion was done to verify the efficiency and accuracy of the method in this paper. Compared to the methods of RANSAC and 3-D-DSF, the proposed algorithm of this paper retained more valuable points, removed more worthless points, and required fewer computational sources. In addition, there were less parameters required for calibration in the proposed algorithm. We have the following major findings: (1) Projecting 3 = D point cloud into 2-D polar coordinates can reduce the computation complexity and memory resource in the point cloud loop processing during the background point filtering.
The processing time using 2-D polar coordinates is about one quarter of the processing time using the 3-D point cloud directly. (2) Using the vertical slice unit, which is different from the horizontal slice unit mentioned in the the state-of-the-art literature, can effectively classify the point cloud into VOPs, WOPs, AGPs, and NGPs. It can further remove the worthless points in the background point filtering.
The accurate and excluding efficiency is about 5 times that is not based on the vertical slice unit.
(3) In practical applications, road pavement conditions and LiDAR installation statuses need to be considered, such as the road gradient and LiDAR sensor inclination. Therefore, a BP-ANN machine learning model was proposed to fit the inclination angle between the BGP and the projection plane to improve the flexibility of the background point algorithm in a complex environment.
However, the present study and its methodological approach are not free from limitations. The first limitation is that the slice-based background point filtering algorithm did not consider the object far away from the LiDAR sensor, which can be obtained from only one single point. Based on the vertical slice unit rule, these points presenting the objects are considered as noise in the algorithm and filtered out. Therefore, it needs to identify the effective filtering range of the proposed algorithm. A second limitation is that the non-driving area far away from the LiDAR sensor cannot be detected accurately, due to a lack of reference points, and the threshold used to differentiate point types increases as the distance increases in the present study. However, the non-driving area is important for pedestrian trajectory tracking and intention prediction in intelligent transportation applications. The third limitation is that the present study did not consider the weather conditions. The proposed algorithm only evaluated during the daytime on a sunny day. However, the bright day and the black night have a different impact on the laser refection intensity of the LiDAR sensor. Snow and rain falling on the ground will influence the accuracy of ground point identification. A fourth limitation is the efficiency and accuracy of traffic object identification and traffic object trajectory tracking from the VOPs and AGPs. The effectiveness and necessity to retain AGPs for assistant identification of the traffic object and trajectory should be further evaluated. The final limitation is that the present study used low-density point cloud from the Velodyne VLP-16 low-channel LiDAR sensor to evaluate the efficiency and accuracy of the proposed algorithm. The high-density point cloud data from a higher channel LiDAR sensor, such as Velodyne VLP-32C and HDL-64, or other manufacturers' LiDAR product models with the same or higher channel, should be used to evaluate the proposed algorithm to test and validate the efficiency, accuracy, and portability of the proposed algorithm in different LiDAR products. For further studies, such problems and factors should be considered and added in the modeling, testing, and validation of the model.

Conflicts of Interest:
The authors declare no conflict of interest.