Automatic ground extraction method for laser point cloud

Abstract

The present invention provides a method for automatic ground extraction from laser point clouds, comprising: inputting a single-frame point cloud in a laser point cloud map package, performing box filtering on the single-frame point cloud, determining an area of ​​interest and extracting a single-frame area of ​​interest point cloud, and then judging whether the cumulative number of frames of the single-frame area of ​​interest point cloud is not less than a first threshold; superimposing all extracted single-frame area of ​​interest point clouds to obtain a superimposed area of ​​interest point cloud, and then removing noise from the superimposed area of ​​interest point cloud using a regional growing method to obtain a target area of ​​interest point cloud; calculating the ground undulation of the target area of ​​interest point cloud; and for the target point cloud map, using cloth filtering to complete ground extraction according to the ground undulation. The present invention improves the automation level of CSF ground extraction; and eliminates most irrelevant points, such as vegetation, buildings, etc., except for the noise points of some vehicles on the road.

Description

Technical Field

[0001] The present invention relates to the technical field of remote sensing and photogrammetry, and in particular, to an automated ground extraction method for laser point cloud, and more particularly, to an automated ground extraction method and system for laser radar point cloud by integrating region growing method and cloth filtering. Background Art

[0002] LiDAR point cloud map is an indispensable part of the autonomous driving technology stack, and extracting relevant elements from LiDAR point cloud map is the key to map application. Among them, extracting ground from point cloud is used in many scenarios, such as loop detection, sensor calibration, etc. Commonly used LiDAR point cloud ground extraction methods mainly include Random sample consensus (RANSAC), Cloth Simulation Filter (CSF), cloth filter, region growing method and ground extraction method based on deep learning. Among them, RANSAC and region growing method are sensitive to noise, and point cloud needs to be preprocessed to achieve better results, but the robustness is poor; CSF has better robustness in extracting ground, but it is necessary to set appropriate parameters according to the degree of ground undulation to obtain a better ground extraction effect, and manual participation is required to judge the degree of ground undulation; deep learning method also has better effect in extracting ground, but it requires a large amount of data to train the model, and the model is not universal. Generally, new models need to be trained for different data to achieve better results. The present invention proposes a universal automatic LiDAR point cloud ground extraction method that integrates CSF and region growing method.

[0003] Patent document CN113052881A (application number: CN202110316363.6) discloses an automatic registration method for indoor three-dimensional point clouds by extracting poles, which is based on avoiding the construction of mesh and normal direction estimation of point cloud data, and aims to solve the automatic registration of multi-station ground laser scanning point cloud data in building indoor scenes. It mainly includes pole extraction and point cloud registration. The method of the present invention is analyzed and verified using real building indoor point cloud data, and the registration parameters of multiple groups of point cloud data are solved. The registration results are compared with the results calculated by manually selecting targets, and the accuracy and reliability of the method of the present invention are further verified. However, the invention does not determine the degree of undulation of the ground by using the regional growing method to preliminarily extract the road surface in the lidar point cloud. Summary of the invention

[0004] In view of the defects in the prior art, the purpose of the present invention is to provide a method for automatic ground extraction of laser point clouds.

[0005] A laser point cloud automatic ground extraction method provided by the present invention comprises:

[0006] Step S1: Input a single-frame point cloud in the laser point cloud map package, perform box filtering on the single-frame point cloud, determine the region of interest and extract the single-frame region of interest point cloud, and then determine whether the cumulative number of frames of the single-frame region of interest point cloud is greater than or equal to the first threshold N a If yes, go to step S2, if no, repeat step S1;

[0007] Step S2: superimposing all the extracted single-frame ROI point clouds to obtain a superimposed ROI point cloud, and removing noise from the superimposed ROI point cloud using a region growing method to obtain a target ROI point cloud;

[0008] Step S3: Calculating the ground relief degree of the point cloud of the target area of ​​interest;

[0009] Step S4: for the target point cloud map, according to the ground undulation, cloth filtering is used to complete ground extraction.

[0010] Preferably, the box filtering of the single-frame point cloud in step S1 includes:

[0011] Pre-set the size parameters;

[0012] Along the road direction, a cuboid with preset dimensions is set as a box, and the area of ​​interest is inside the box;

[0013] The laser points within the area of ​​interest are retained, and the laser points outside the area of ​​interest are filtered out.

[0014] Preferably, the first threshold N in step S1 a It is determined according to the number of single-frame point cloud frames required to form the target point cloud map.

[0015] Preferably, in step S2, the noise removal is performed on the superimposed point cloud of the region of interest using a region growing method, comprising:

[0016] Step T1: Estimate the normal vector of each laser point, and use a method based on local surface fitting to search for other K points p adjacent to each laser point. i (i=1,2,…,K), these K neighboring points fit the plane P, and the normal vector of plane P is the normal vector of the laser point. Considering that plane P passes through the centroid p0 of the K neighboring points, the covariance matrix M of p0 is:

[0017]

[0018] Perform eigenvalue decomposition, where the eigenvector corresponding to the minimum eigenvalue is the normal vector of plane P, that is, the normal vector of the laser point;

[0019] Step T2: Select the seed point set from each laser point, and perform eigenvalue decomposition on the covariance matrix of each laser point to obtain three eigenvalues ​​λ0, λ1 and λ2. If their magnitude relationship is λ0≤λ1≤λ2, the curvature of the laser point is:

[0020]

[0021] The smaller the curvature, the flatter the neighborhood. Select the curvature less than the threshold δ th The laser points are taken as seed points to obtain a seed point set;

[0022] Step T3: Use the seed point set to remove noise points. For each seed point, use the KD-Tree method to search for all neighboring points of the seed point, and calculate the normal vector difference between each neighbor point and the seed point. If the difference is greater than or equal to the second threshold, it is determined that the neighbor point does not belong to the ground point. If the difference is not greater than the second threshold, then calculate whether the curvature of the neighbor point is less than the third threshold. If not, it is determined that the neighbor point does not belong to the ground point. If so, it is determined that the neighbor point is a ground point. After removing all noise points, the remaining ground points constitute the point cloud of the target area of ​​interest.

[0023] Preferably, the step S3 of calculating the ground relief degree of the target area of ​​interest point cloud includes:

[0024] Step Y1: Divide the point cloud of the target area of ​​interest into N blocks, where N is a positive integer and its initial value is arbitrarily selected;

[0025] Step Y2: For each point cloud, calculate its height mean H and height variance V:

[0026]

[0027]

[0028] Among them, R is the number of points in each point cloud, h i (i=1,2,…,R) is the height of the i-th point;

[0029] Step Y3: Determine whether the height variance of each point cloud is greater than the fourth threshold. If the height variance of a certain point cloud is greater than the fourth threshold, increase the N value, re-divide the point cloud of the target area of ​​interest, and repeat step Y2; if the height variance of all N point clouds is not greater than the fourth threshold, calculate the ground undulation S of the point cloud of the target area of ​​interest:

[0030]

[0031] Among them, H j(j=1,2,…,N) is the average height of the jth point cloud, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, θ is and Angle.

[0032] Preferably, in step S4, for the target point cloud map, cloth filtering is used to complete ground extraction according to the ground undulation, including:

[0033] Step Z1: According to the ground undulation degree S, determine whether the ground is flat, slightly undulating or steep, and then determine the Rigidity parameter of the cloth filter. The Rigidity parameter is the RI value, which is used to determine the movement times and movement distance of the grid particles. The flatter the ground, the larger the RI value, and the steeper the ground, the smaller the RI value;

[0034] Step Z2: Invert all laser points in the target point cloud map, initialize the cloth grid and set the grid size;

[0035] Step Z3: Project all laser points and grid particles onto the same horizontal plane, find the nearest neighbor of each grid particle, and record its elevation before projection;

[0036] Step Z4: for each movable grid particle, update its height according to the displacement caused by the gravity, and then if the height of the grid particle is less than or equal to the elevation IHV of its nearest neighbor before projection, set the height of the grid particle to IHV and set the grid particle to an immovable point; if the height of the grid particle is greater than IHV, further update its height according to the displacement caused by the internal driving factor of the grid particle; repeat step Z4 until the number of iterations reaches the preset value or the maximum height change of all grid particles is less than the preset value;

[0037] Step Z5: Calculate the height difference between each laser point and the corresponding grid particle. If the difference is less than or equal to the fifth threshold, the laser point is judged to be a ground point. If the difference is greater than the fifth threshold, the laser point is judged to be a non-ground point.

[0038] Preferably, the target point cloud map is composed of N a Single-frame point clouds are superimposed to avoid insufficient features. a is the first threshold.

[0039] Compared with the prior art, the present invention has the following beneficial effects:

[0040] 1. The present invention uses the region growing method to preliminarily extract the road surface in the laser radar point cloud to determine the undulation of the ground, thereby improving the automation level of CSF ground extraction;

[0041] 2. The present invention utilizes the fact that vehicles equipped with laser radars will travel on the road. By limiting a long and narrow area along the road direction for each frame of radar point cloud data, the potential area of ​​the road is obtained, thereby eliminating most of the irrelevant points except for the noise points of some vehicles on the road, such as vegetation, buildings, etc. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] Other features, objects and advantages of the present invention will become more apparent from the detailed description of non-limiting embodiments made with reference to the following drawings:

[0043] Figure 1 It is a flow chart of the present invention;

[0044] Figure 2 Flowchart for box filtering;

[0045] Figure 3 Flowchart for removing noise using the region growing method;

[0046] Figure 4 Flow chart for ground relief calculation;

[0047] Figure 5 This is the fabric filtering flow chart;

[0048] Figure 6 To overlay the point cloud of the area of ​​interest;

[0049] Figure 7 Point cloud of the area of ​​interest;

[0050] Figure 8 Extract ground points of single frame point cloud for CSF;

[0051] Fig. 9 Extract non-ground points of single-frame point cloud for CSF;

[0052] Fig.10 The point cloud of the target area of ​​interest after removing noise using the region growing method;

[0053] Fig.11 Extract ground points of local point cloud map for CSF under mild undulation;

[0054] Fig.12 Extract non-ground points of the local point cloud map for CSF under mild undulation conditions;

[0055] Fig.13 Extract ground points of local point cloud map for CSF in flat state;

[0056] Fig.14 Extract ground points of local point cloud map for CSF in flat state. DETAILED DESCRIPTION

[0057] The present invention is described in detail below in conjunction with specific embodiments. The following embodiments will help those skilled in the art to further understand the present invention, but are not intended to limit the present invention in any form. It should be noted that, for those of ordinary skill in the art, several changes and improvements can also be made without departing from the concept of the present invention. These all belong to the protection scope of the present invention.

[0058] Embodiment 1:

[0059] According to the present invention, a laser point cloud automatic ground extraction method is provided. Figure 1-Figure 14 As shown, including:

[0060] Step S1: Input a single-frame point cloud in the laser point cloud map package, perform box filtering on the single-frame point cloud, determine the region of interest and extract the single-frame region of interest point cloud, and then determine whether the cumulative number of frames of the single-frame region of interest point cloud is greater than or equal to the first threshold N a If yes, go to step S2, if no, repeat step S1;

[0061] Specifically, the box filtering of the single-frame point cloud in step S1 includes:

[0062] Along the road direction, a cuboid with preset dimensions is set as a box, and the area of ​​interest is inside the box;

[0063] The laser points within the area of ​​interest are retained, and the laser points outside the area of ​​interest are filtered out.

[0064] Specifically, the first threshold N in step S1 a It is determined according to the number of single-frame point cloud frames required to form the target point cloud map.

[0065] Step S2: superimposing all the extracted single-frame ROI point clouds to obtain a superimposed ROI point cloud, and removing noise from the superimposed ROI point cloud using a region growing method to obtain a target ROI point cloud;

[0066] Specifically, in step S2, the noise removal is performed on the superimposed point cloud of the region of interest by using the region growing method, including:

[0067] Step T1: Estimate the normal vector of each laser point, and use a method based on local surface fitting to search for other K points p adjacent to each laser point. i(i=1,2,…,K), these K neighboring points fit the plane P, and the normal vector of plane P is the normal vector of the laser point. Considering that plane P passes through the centroid p0 of the K neighboring points, the covariance matrix M of p0 is:

[0068]

[0069] Perform eigenvalue decomposition, where the eigenvector corresponding to the minimum eigenvalue is the normal vector of plane P, that is, the normal vector of the laser point;

[0070] Step T2: Select the seed point set from each laser point, and perform eigenvalue decomposition on the covariance matrix of each laser point to obtain three eigenvalues ​​λ0, λ1 and λ2. If their magnitude relationship is λ0≤λ1≤λ2, the curvature of the laser point is:

[0071]

[0072] The smaller the curvature, the flatter the neighborhood. Select the curvature less than the threshold δ th The laser points are taken as seed points to obtain a seed point set;

[0073] Step T3: Use the seed point set to remove noise points. For each seed point, use the KD-Tree method to search for all neighboring points of the seed point, and calculate the normal vector difference between each neighbor point and the seed point. If the difference is greater than the second threshold, it is determined that the neighbor point does not belong to the ground point. If the difference is not greater than the second threshold, then calculate whether the curvature of the neighbor point is less than the third threshold. If not, it is determined that the neighbor point does not belong to the ground point. If so, it is determined that the neighbor point is a ground point. After removing all noise points, the remaining ground points constitute the point cloud of the target area of ​​interest.

[0074] Step S3: Calculating the ground relief degree of the point cloud of the target area of ​​interest;

[0075] Specifically, the step S3 of calculating the ground relief degree of the target area of ​​interest point cloud includes:

[0076] Step Y1: Divide the point cloud of the target area of ​​interest into N blocks, where N is a positive integer and its initial value is arbitrarily selected;

[0077] Step Y2: For each point cloud, calculate its height mean H and height variance V:

[0078]

[0079]

[0080] Among them, R is the number of points in each point cloud, h i (i=1,2,…,R) is the height of the i-th point;

[0081] Step Y3: Determine whether the height variance of each point cloud is greater than the fourth threshold. If the height variance of a certain point cloud is greater than the fourth threshold, increase the N value, re-divide the point cloud of the target area of ​​interest, and repeat step Y2; if the height variance of all N point clouds is not greater than the fourth threshold, calculate the ground undulation S of the point cloud of the target area of ​​interest:

[0082]

[0083] Among them, H j (j=1,2,…,N) is the average height of the jth point cloud, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, θ is and Angle.

[0084] Step S4: for the target point cloud map, according to the ground undulation, cloth filtering is used to complete ground extraction.

[0085] Specifically, in step S4, for the target point cloud map, cloth filtering is used to complete ground extraction according to the ground undulation, including:

[0086] Step Z1: According to the ground undulation degree S, determine whether the ground is flat, slightly undulating or steep, and then determine the Rigidity parameter of the cloth filter. The Rigidity parameter is the RI value, which is used to determine the movement times and movement distance of the grid particles. The flatter the ground, the larger the RI value, and the steeper the ground, the smaller the RI value;

[0087] Step Z2: Invert all laser points in the target point cloud map, initialize the cloth grid and set the grid size;

[0088] Step Z3: Project all laser points and grid particles onto the same horizontal plane, find the nearest neighbor of each grid particle, and record its elevation before projection;

[0089] Step Z4: for each movable grid particle, update its height according to the displacement caused by the gravity, and then if the height of the grid particle is less than or equal to the elevation IHV of its nearest neighbor before projection, set the height of the grid particle to IHV and set the grid particle to an immovable point; if the height of the grid particle is greater than IHV, further update its height according to the displacement caused by the internal driving factor of the grid particle; repeat step Z4 until the number of iterations reaches the preset value or the maximum height change of all grid particles is less than the preset value;

[0090] Step Z5: Calculate the height difference between each laser point and the corresponding grid particle. If the difference is less than or equal to the fifth threshold, the laser point is judged to be a ground point. If the difference is greater than the fifth threshold, the laser point is judged to be a non-ground point.

[0091] Specifically, the target point cloud map is composed of N a Single-frame point clouds are superimposed to avoid insufficient features. a is the first threshold.

[0092] Embodiment 2:

[0093] Embodiment 2 is a preferred example of Embodiment 1, and is used to illustrate the present invention in more detail.

[0094] The core point of the present invention is to determine the undulation of the ground by preliminarily extracting the road surface in the laser radar point cloud using the region growing method, thereby improving the automation level of CSF ground extraction.

[0095] The core point of the present invention is to utilize the characteristic that vehicles equipped with laser radars will travel on the road, and by limiting a rectangular area along the road direction for each frame of radar point cloud data, obtain an open road area of ​​interest, thereby eliminating most of the irrelevant points except for the noise points of some vehicles on the road, such as vegetation, buildings, etc., and this road area of ​​interest will reflect the undulation of the ground.

[0096] The core point of the present invention is to measure the ground undulation by estimating the normal vector and height index S of each laser point. The calculation formula is:

[0097]

[0098] Among them, H j (j=1,2,…,N) is the average height of the jth point cloud, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, θ is and Angle.

[0099] The present invention inputs a single-frame point cloud in a point cloud map package, and then performs box filtering on the point cloud to extract the area of ​​interest. Box filtering is a filter that filters point clouds in a limited area. When used, the size of the box filter needs to be preset. The area of ​​interest is within the box, and the laser points within the area of ​​interest are retained, while the laser points outside the area of ​​interest are filtered out. The size of the box filter used in the present invention requires a rectangular parallelepiped that meets the size requirements along the road direction. The purpose is to remove a large amount of information that is irrelevant to the ground, leaving only a part of the ground and a small number of noise points. All the extracted single-frame area of ​​interest point clouds are superimposed. When the number of superimposed point cloud frames is sufficient to generate a superimposed area of ​​interest point cloud, the regional growing method is executed to remove noise points from the superimposed area of ​​interest point cloud to obtain the target point cloud of interest. The reason for superimposing point clouds is that the ground features of single-frame point clouds are not obvious enough, and the effect of directly using filtering methods to extract the ground on single-frame point clouds is poor. After removing noise points by region growing method, the main body in the area of ​​interest is the ground. Then, according to the ground undulation S defined in the present invention, the ground undulation degree of the area of ​​interest is measured and divided into three states: flat, slightly undulating, and steep. Then, it is transferred to CSF ​​to extract the ground from the local point cloud map. Each part will be described in detail below.

[0100] The present invention first uses box filtering as a method to remove most of the noise points in the point cloud. Because the present invention only needs to obtain the undulation of the ground instead of extracting the entire ground before using CSF, and the car of the laser radar point cloud data is collected on the road, the undulation of the ground can reflect the undulation of the ground of the point cloud data. Therefore, the key to this step is to set the box filter size. A rectangular block of the same size should be preset in the laser radar coordinate system along the road direction. The purpose is to retain the laser points in the area of ​​interest after removing irrelevant points outside the area of ​​interest through box filtering. At this time, the main body in the area of ​​interest is the ground.

[0101] The present invention uses the region growing method to remove the noise points outside the road surface from the region of interest extracted by box filtering. The first step is to estimate the normal vector of each laser point, and then use the method based on local surface fitting to search for other K points p adjacent to each laser point. i (i=1,2,…,K), and then use the least squares method to calculate the K nearest points to fit the plane P. The normal vector of the plane P is the normal vector of the laser point. This plane can be expressed as:

[0102]

[0103] In the formula, is the normal vector of plane P, d is the distance from P to the origin of the coordinate system, and Pi is the i-th point among K points;

[0104] Considering that the plane P passes through the centroid p0 of K neighboring points, the covariance matrix M of p0 is:

[0105]

[0106] Perform eigenvalue decomposition, where the eigenvector corresponding to the minimum eigenvalue is the normal vector of plane P, that is, the normal vector of the laser point;

[0107] Secondly, the seed point set is selected from each laser point. For each laser point, the covariance matrix is ​​decomposed to obtain three eigenvalues ​​λ0, λ1 and λ2. Assuming that the eigenvalue relationship is λ0≤λ1≤λ2, the curvature of the laser point is:

[0108]

[0109] Among them, δ is the curvature, the smaller it is, the flatter the neighborhood is, and the seed point is the point whose curvature is less than the threshold δ th The laser points are obtained to obtain the seed point set, λ0 is the first eigenvalue of the matrix M, λ1 is the second eigenvalue of the matrix M, and λ2 is the third eigenvalue of the matrix M.

[0110] Then, the seed point set is used to remove noise points. For each seed point, KD-Tree is used to search for all neighboring points of the seed point, and the normal difference between the neighboring point and the seed point is calculated. If the difference is greater than the second threshold, the neighbor point is judged not to be a ground point. If the difference is not greater than the second threshold, then it is calculated whether the curvature of the neighbor point is less than the third threshold. If not, the neighbor point is judged not to be a ground point. If so, the neighbor point is judged to be a ground point. After removing all noise points, the remaining ground points constitute the point cloud of the target area of ​​interest.

[0111] The present invention proposes a ground undulation S to measure the degree of ground undulation, and the calculation formula is:

[0112]

[0113] Among them, H j (j=1,2,…,N) is the average height of the jth point cloud, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, θ is and Angle.

[0114] The calculation flow chart is as follows Figure 4 As shown, first, the point cloud of the target area of ​​interest after denoising by the growing region method needs to be divided into N blocks, where N is a positive integer and its initial value is arbitrarily selected. For each point cloud, its height average H and height variance V are calculated. The calculation method is as follows:

[0115]

[0116]

[0117] Among them, R is the number of points in each point cloud, h i (i=1,2,…,R) is the height of the i-th point;

[0118] Determine whether the height variance of each point cloud is greater than the fourth threshold. If the height variance of a certain point cloud is greater than the fourth threshold, increase the N value, re-divide the target area of ​​interest point cloud, and then repeat for each point cloud, calculate its height average H and height variance V; if the height variance of all N point clouds is not greater than the fourth threshold, calculate the ground undulation S of the target area of ​​interest point cloud.

[0119] Finally, the present invention uses cloth filtering to extract the ground of the local point cloud map. The basic formula of cloth filtering is

[0120]

[0121] Among them, m represents the mass of the particle, which is generally 1, X represents the position of the particle in the cloth at time t, and F ext (X, t) represents the external driving factors (gravity, collision, etc.), F inte (X,t) represents the internal driving factors.

[0122] First, assume that there is only one external driving factor F ext (X, t), and let m be 1, then the above formula can be expressed as

[0123] x(t+Δt)=2X(t)-X(t-Δt)+GΔt2

[0124] Among them, Δt is the time step and G is a constant.

[0125] Then consider the internal driving factor F inte (X, t), randomly select two adjacent particles. If both particles are movable, move them the same distance in opposite directions. If one particle is immovable, move the other particle. If both particles have the same height, do not move them. The displacement can be calculated by the following formula:

[0126]

[0127] in, is the displacement of the particle; when the particle is movable, b is equal to 1, and when it is immovable, b is equal to 0; pi is the adjacent particle of p0, It is to normalize the point to the unit vector (0,0,1) in the vertical direction. T .

[0128] The ground undulation is used to determine whether the ground is flat, slightly undulating or steep, and this parameter is passed to the cloth filter to adjust the parameter rigidity, RI, in the cloth filter. When the ground is flat, the RI value is large, and when the ground is steeper, the RI value is smaller.

[0129] The parameter rigidity, RI, is used to describe the number of times a particle moves. For example, when RI is 1, it moves once, and the distance moved is 1 / 2 of the height difference between the two particles; when RI is 2, it moves twice, and moves 3 / 4. Similarly, when RI is 3, it moves three times, and the distance moved is 7 / 8 of the height difference between the two particles. The larger the RI value, the harder the "cloth". Generally speaking, a larger RI value is required for flat land, and a smaller RI value is required for steep slopes and mountainous areas.

[0130] After that, the laser points in the input local point cloud map are inverted, the cloth grid is initialized and the grid size is set. All laser points and grid particles are projected onto the same horizontal plane, the nearest neighbor of each grid particle is found, and the elevation before projection (intersection height value, IHV) is recorded. For each movable grid particle, its height is updated according to the displacement caused by its gravity, and then if the height of the grid particle is less than or equal to the elevation IHV of its nearest neighbor before projection, the height of the grid particle is set to IHV and the grid particle is set to an immovable point; if the height of the grid particle is greater than IHV, the height of the grid particle is further updated according to the displacement caused by the internal driving factors; repeat the above steps until the number of iterations reaches the preset value or the maximum height change of all grid particles is less than the preset value. Calculate the height difference between each laser point and the corresponding grid particle. If the difference is less than or equal to the fifth threshold, the laser point is judged to be a ground point. If the difference is greater than the fifth threshold, the laser point is judged to be a non-ground point.

[0131] The present invention firstly a The single-frame point clouds are superimposed to form a superimposed point cloud map of the area of ​​interest, and each frame of the point cloud in the superimposed target point cloud map is box filtered in the laser radar coordinate system to obtain the point cloud of the target area of ​​interest. The superimposed point cloud of the area of ​​interest is as follows Figure 6 , the point cloud of the area of ​​interest is as follows Figure 7 .

[0132] The reason for using single-frame point cloud overlay to form an overlay point cloud map of the area of ​​interest is that the ground points of single-frame point cloud are sparse, resulting in unclear ground features, and the effect of using CSF extraction is poor. For details, see Figure 8 and Fig. 9 , while superimposing the point cloud map of the area of ​​interest is to superimpose dozens of frames of point clouds in the same coordinate system, which can avoid this problem.

[0133] The extracted overlay point cloud map of the area of ​​interest is used to remove some noise points on the ground using the regional growing method, such as Fig.10 ,and Figure 7 In comparison, most of the non-ground points are removed, reducing the impact of noise on the calculation of ground relief.

[0134] Calculate the normal vector of each laser point in the point cloud of the target area of ​​interest after removing all noise points, and bring it into the ground relief formula to calculate the ground relief S. The calculation formula is:

[0135]

[0136] Among them, H j (j=1,2,…,N) is the average height of the jth point cloud, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, θ is and Angle.

[0137] The requirement for dividing the point cloud of the target area of ​​interest into N blocks is that the variance of the Z coordinate of each point cloud should be kept within a small range. If the variance is large, it means that the number of blocks N is small, resulting in a large degree of undulation inside the segmented point cloud, and the average height value H is not representative. When the value of S is close to 0, the ground undulation state is flat. The larger the value of S, the greater the degree of ground undulation.

[0138] By calculation Fig.10 The ground undulation S of the point cloud of the target area of ​​interest in the middle shows that the ground undulation of the point cloud of the target area of ​​interest is large, which is a slight undulation state. It is passed to CSF ​​and extracted using CSF Figure 6 The ground in the local point cloud map can be obtained Fig.11 and Fig.12 result.

[0139] If a flat state is input to CSF, then extract Figure 6 The ground in the local point cloud map can be obtained Fig.13 and Fig.14 result.

[0140] It can be seen from the comparison of the results in the figure that the present invention can significantly improve the degree of automation of CSF extraction of point cloud ground. For ground surfaces with different undulations, better point cloud ground extraction effects can be obtained without artificially adjusting CSF parameters.

[0141] Those skilled in the art know that, in addition to implementing the system, device and its various modules provided by the present invention in a purely computer-readable program code, it is entirely possible to implement the same program in the form of logic gates, switches, application-specific integrated circuits, programmable logic controllers and embedded microcontrollers by logically programming the method steps. Therefore, the system, device and its various modules provided by the present invention can be considered as a hardware component, and the modules included therein for implementing various programs can also be considered as structures within the hardware component; the modules for implementing various functions can also be considered as both software programs for implementing the method and structures within the hardware component.

[0142] The above describes the specific embodiments of the present invention. It should be understood that the present invention is not limited to the above specific embodiments, and those skilled in the art can make various changes or modifications within the scope of the claims, which does not affect the essence of the present invention. In the absence of conflict, the embodiments of the present application and the features in the embodiments can be combined with each other arbitrarily.

Claims

1. A method for automatic ground extraction of laser point clouds, characterized in that: include: Step S1: Input a single-frame point cloud in the laser point cloud map package, perform box filtering on the single-frame point cloud, determine the region of interest and extract the single-frame region of interest point cloud, and then determine whether the cumulative number of frames of the single-frame region of interest point cloud is greater than or equal to the first threshold N a If yes, go to step S2, if no, repeat step S1; Step S2: superimposing all the extracted single-frame ROI point clouds to obtain a superimposed ROI point cloud, and removing noise from the superimposed ROI point cloud using a region growing method to obtain a target ROI point cloud; Step S3: Calculating the ground relief degree of the point cloud of the target area of ​​interest; Step S4: for the target point cloud map, according to the ground undulation, use cloth filtering to complete ground extraction; The step S1 of performing box filtering on the single-frame point cloud comprises: Pre-set size parameters; A cuboid with preset dimensions is set along the road direction as a box, and the area of ​​interest is inside the box; Keep the laser points within the area of ​​interest and filter out the laser points outside the area of ​​interest; The step S3 of calculating the ground relief degree of the point cloud of the target area of ​​interest includes: Step Y1: Divide the point cloud of the target area of ​​interest into N blocks, where N is a positive integer and its initial value is arbitrarily selected; Step Y2: For each point cloud, calculate its height mean H and height variance V: Among them, R is the number of points in each point cloud, h i is the height of the i-th point, i=1,2,…,R; Step Y3: Determine whether the height variance of each point cloud is greater than the fourth threshold. If the height variance of a certain point cloud is greater than the fourth threshold, increase the N value, re-divide the point cloud of the target area of ​​interest, and repeat step Y2; if the height variance of all N point clouds is not greater than the fourth threshold, calculate the ground undulation S of the point cloud of the target area of ​​interest: Among them, H j is the average height of the jth point cloud, j = 1, 2, ..., N, W is the number of points in the point cloud of the target area of ​​interest, is the vertical unit normal vector, is the unit normal vector of the kth point, k=1,2,…,W,θ is and Angle.

2. The method for automatic ground extraction of laser point cloud according to claim 1, characterized in that: The first threshold N in step S1 a It is determined according to the number of single-frame point cloud frames required to form the target point cloud map.

3. The method for automatic ground extraction of laser point cloud according to claim 1, characterized in that: In step S2, the noise removal is performed on the superimposed point cloud of the region of interest using a region growing method, including: Step T1: Estimate the normal vector of each laser point, and use a method based on local surface fitting to search for other K points p adjacent to each laser point. l , l=1,2,…,K, these K neighboring points fit the plane P, and the normal vector of plane P is the normal vector of the laser point. Considering that plane P passes through the centroid p0 of the K neighboring points, the covariance matrix M of p0 is: Perform eigenvalue decomposition, where the eigenvector corresponding to the minimum eigenvalue is the normal vector of plane P, that is, the normal vector of the laser point; Step T2: Select the seed point set from each laser point, and perform eigenvalue decomposition on the covariance matrix of each laser point to obtain three eigenvalues ​​λ0, λ1 and λ2. If their magnitude relationship is λ0≤λ1≤λ2, the curvature of the laser point is: The smaller the curvature, the flatter the neighborhood. Select the curvature less than the threshold δ th The laser points are taken as seed points to obtain a seed point set; Step T3: Use the seed point set to remove noise points. For each seed point, use the KD-Tree method to search for all neighboring points of the seed point, and calculate the normal vector difference between each neighbor point and the seed point. If the difference is greater than the second threshold, it is determined that the neighbor point does not belong to the ground point. If the difference is not greater than the second threshold, then calculate whether the curvature of the neighbor point is less than the third threshold. If not, it is determined that the neighbor point does not belong to the ground point. If so, it is determined that the neighbor point is a ground point. After removing all noise points, the remaining ground points constitute the point cloud of the target area of ​​interest.

4. The method for automatic ground extraction of laser point cloud according to claim 1, characterized in that: In step S4, for the target point cloud map, cloth filtering is used to complete ground extraction according to the ground undulation, including: Step Z1: According to the ground undulation degree S, determine whether the ground is flat, slightly undulating or steep, and then determine the Rigidity parameter of the cloth filter. The Rigidity parameter is the RI value, which is used to determine the movement times and movement distance of the grid particles. The flatter the ground, the larger the RI value, and the steeper the ground, the smaller the RI value; Step Z2: Invert all laser points in the target point cloud map, initialize the cloth grid and set the grid size; Step Z3: Project all laser points and grid particles onto the same horizontal plane, find the nearest neighbor of each grid particle, and record its elevation before projection; Step Z4: for each movable grid particle, update its height according to the displacement caused by the gravity, and then if the height of the grid particle is less than or equal to the elevation IHV of its nearest neighbor before projection, set the height of the grid particle to IHV and set the grid particle to an immovable point; if the height of the grid particle is greater than IHV, further update its height according to the displacement caused by the internal driving factor of the grid particle; repeat step Z4 until the number of iterations reaches the preset value or the maximum height change of all grid particles is less than the preset value; Step Z5: Calculate the height difference between each laser point and the corresponding grid particle. If the difference is less than or equal to the fifth threshold, the laser point is judged to be a ground point. If the difference is greater than the fifth threshold, the laser point is judged to be a non-ground point.

5. The method for automatic ground extraction of laser point cloud according to claim 2, characterized in that: The target point cloud map is composed of N a Single frame point clouds are superimposed to avoid insufficient features. a is the first threshold.

Images