A LiDAR data processing method for complex forest environments

Abstract

The present invention provides a laser radar data processing method for a complex forest environment, which obtains original point cloud data of a target forest area through an airborne laser radar system or a backpack laser radar system, obtains preprocessed point cloud data through preprocessing, obtains non-ground point cloud data by filtering out ground points in the preprocessed point cloud data, and classifies the ground objects in the non-ground point cloud data to obtain ground object classification point cloud data; based on the ground object classification point cloud data, a three-dimensional canopy point cloud model is constructed, and based on a single tree parameter extraction method, the single tree spatial position and key forestry parameters of individual trees in the target forest area are obtained, and the single tree parameter extraction method includes crown segmentation and allometric growth equations. The present invention obtains ground object vertical structure information, removes ground point data and classifies ground objects through a radar system, so as to achieve high-precision acquisition of single tree and forest stand level parameters in a complex forest environment.

Description

Technical Field

[0001] The invention belongs to the technical field of complex forest data processing, and in particular relates to a laser radar data processing method for a complex forest environment. Background Art

[0002] In existing forest management, it is often necessary to analyze and calculate various parameters of trees in the forest, including tree height, diameter at breast height, canopy density, tree species, etc. These data are of great significance to forest resource management, ecological protection and sustainable development. Through LiDAR technology, various forest data can be obtained. LiDAR obtains three-dimensional information of trees by emitting laser pulses, including tree height, canopy shape and distribution; by processing the point cloud data obtained by LiDAR scanning, a detailed three-dimensional model can be generated and various parameters of trees can be extracted. It is more efficient and accurate than traditional ground survey methods, and is particularly suitable for large-scale forest surveys and dynamic monitoring. However, the current processing of complex forest data extracted by LiDAR has the following problems:

[0003] (1) The canopy segmentation and parameter estimation are mainly based on the canopy height model (CHM). The canopy surface height relative to the ground is obtained by extracting the difference between the normalized digital surface model (DSM) and the digital elevation model (DEM). However, in complex forest stands with high canopy density, CHM is often difficult to accurately reconstruct the canopy surface morphology due to factors such as canopy occlusion and edge loopholes, and cannot effectively express the vertical structure information of the canopy, resulting in insufficient accuracy in subsequent tree positioning and parameter estimation.

[0004] (2) When processing lidar data, most of them rely on simple threshold segmentation or classification algorithms based on single features. The classification accuracy is usually low, resulting in confusion between ground points and object points, which in turn affects the accuracy of subsequent tree parameter extraction;

[0005] (3) Due to problems such as tree density, obstruction, and measurement accuracy limitations, there are large errors in the estimation of the precise location of individual trees and forestry parameters. The positioning accuracy of individual trees is insufficient, resulting in inaccurate estimates of key indicators such as forest standing volume and forest density.

[0006] Therefore, we need to develop a lidar data processing method for complex forest environments that can obtain vertical structure information of forests, high-precision point cloud processing, and accurate single tree parameter estimation to provide more reliable data support for forest management. Summary of the invention

[0007] The purpose of the present invention is to provide a lidar data processing method for complex forest environments to solve the problems of insufficient accuracy of canopy structure analysis, low accuracy of terrain and ground object classification, and low accuracy of single tree positioning and parameter estimation in the existing forest data processing schemes mentioned in the above background technology.

[0008] To achieve the above object, the present invention provides a laser radar data processing method for a complex forest environment, the method is specifically as follows:

[0009] Acquire original point cloud data of the target forest area by using an airborne laser radar system or a backpack laser radar system, and preprocess the original point cloud data to obtain preprocessed point cloud data;

[0010] Filtering out ground points in the preprocessed point cloud data by an adaptive terrain filtering method to obtain non-ground point cloud data, and classifying ground objects in the non-ground point cloud data by a multi-feature fusion classification method to obtain ground object classification point cloud data;

[0011] The adaptive terrain filtering method includes dividing the pre-processed point cloud data into ground points and non-ground points, and filtering out the ground points in the pre-processed point cloud data to obtain non-ground point cloud data; based on the non-ground point cloud data, further identifying the ground object category to which the non-ground points belong;

[0012] The adaptive terrain filtering method comprises: gridding the pre-processed point cloud data, constructing a Delaunay triangulation network, and calculating the normal vector and roughness coefficient of each triangular face in the Delaunay triangulation network;

[0013] The multi-feature fusion classification method includes predicting the feature vectors of the laser points of the non-ground point cloud data by a trained classifier to obtain a category label of each of the laser points;

[0014] Based on the classified point cloud data of the ground objects, a three-dimensional canopy point cloud model is constructed, and based on the single tree parameter extraction method, the single tree spatial position and key forestry parameters of individual trees in the target forest area are obtained. The single tree parameter extraction method includes crown segmentation and allometric growth equation.

[0015] Based on the above-mentioned scheme, the airborne lidar system is mounted on a flying platform to scan the target forest area within a preset relative flight altitude to obtain a large range of the original point cloud data; the backpack lidar system is integrated into a portable backpack device, and the close-range high-density original point cloud data is collected by operators walking in the forest area.

[0016] Based on the above scheme, the pre-processed point cloud data is gridded, specifically:

[0017] The pre-processed point cloud data is divided into grids in the XY direction. The size is d×d, and the initial grid set is constructed , the initial grid set All of the grids are included ;

[0018] The construction of the Delaunay triangulation network is specifically as follows:

[0019] In the grid Select the one with the lowest elevation value. points as initial ground seed points , construct a seed point set , the seed point set All the initial ground seed points are included ;

[0020] Based on the seed point set Construct a Delaunay triangulation and calculate the normal vector of each triangle in the Delaunay triangulation and roughness coefficient ;

[0021] Based on the normal vector and the roughness coefficient , adaptively calculate the elevation threshold of each triangle patch and angle threshold ;

[0022] Traverse each of the triangular facets, and find the non-initial ground seed points in the neighborhood of the triangular facets. The seed point is recorded as a non-seed point, and the height difference between the non-seed point and the triangle is calculated. Angle between the normal vector and , if the non-seed point satisfies and , then the non-seed point is marked as a candidate ground point and added to the seed point set .

[0023] Based on the above scheme, the multi-feature fusion classification method specifically includes:

[0024] Based on the local neighborhood information of the point cloud, the geometric features of each laser point in each of the non-ground point cloud data are calculated, and the geometric features include: normal vector , curvature and characteristic entropy ;

[0025] Performing distance correction and filtering denoising preprocessing on the reflection intensity information of the laser points in the non-ground point cloud data;

[0026] Based on the scale characteristics of different objects, the geometric features and radiation features of the laser point at multiple different neighborhood scales are extracted, and the geometric features and the radiation features are combined into the final feature vector of the laser point. ;

[0027] The feature vector is analyzed by using Z-score or Max-Min method. The dimensions of the features in are scaled and normalized to eliminate the differences in the dimensions of different features and obtain the normalized feature vector ;

[0028] In the non-ground point cloud data, typical ground object points are manually marked, and the typical ground object points include ground, tree trunks, tree crowns, buildings, etc., and feature vectors of the typical ground object points are extracted as training samples ,in is the category label of the typical feature point, representing the feature category;

[0029] Select a classifier and pass the training samples Inputting into the classifier for model learning, and establishing a mapping relationship between the feature vector and the ground object category in the classifier;

[0030] The normalized feature vector of each laser point is obtained by the trained classifier. Prediction is performed to obtain the category label of each laser point and obtain the ground object classification point cloud data.

[0031] Based on the above scheme, the construction of a three-dimensional canopy point cloud model includes: the input of the three-dimensional canopy point cloud model is the ground object classification point cloud data, the point cloud data in the ground object classification point cloud data is divided into a two-dimensional grid in the XY plane, and then the voxels in the Z direction are divided into three-dimensional voxels, the point cloud density inside each of the three-dimensional voxels is quantified and normalized, the canopy surface height of the two-dimensional grid is extracted, and the three-dimensional canopy point cloud model is output.

[0032] Furthermore, the construction of the three-dimensional canopy point cloud model is specifically as follows:

[0033] The point cloud data is: , ,in Represents the i-th point cloud data point; Representing point cloud data The category label of the XY plane is , the voxel height in the Z direction is ; Project the point cloud data onto the XY plane, according to Divide the point cloud data range into grids to obtain a two-dimensional grid , , ,in and are the number of two-dimensional grids in the X and Y directions respectively, and satisfy: , ;

[0034] In each of the two-dimensional grids Based on the Z-axis direction The voxel division in the Z direction is performed for the step size, and the two-dimensional grid is expanded into three-dimensional voxels ,in, is the number of voxels in the Z direction and satisfies: ;

[0035] Traverse each of the three-dimensional voxels , and count the point cloud data inside it Category label , if exists (i.e., tree category), the three-dimensional voxel Mark as , otherwise marked as ;

[0036] By property value Quantify each of the Internal point cloud density: ,in, is an indicative function, when the point cloud data Category label When it is 1, it is 1, otherwise it is 0; for the three-dimensional voxel that does not contain the point , attribute value Recorded as 0;

[0037] Regarding the The attribute value of Normalize in the range 0-1: ,in and Respectively all the of The minimum and maximum values ​​of

[0038] In each of the two-dimensional grids Within, extract the maximum of , the The Z value of the two-dimensional grid The canopy surface height : ,make , based on all the two-dimensional grids The canopy surface height Construct the canopy height model, denoted as ;

[0039] The three-dimensional canopy point cloud model is output as: , and the canopy height model: .

[0040] Based on the above scheme, the tree crown segmentation includes: Search for the local maximum point in the neighborhood, extract all the pixels with the largest local height and mark them as , the search neighborhood size is larger than the single tree crown size, based on the Eliminate points whose relative height difference in the search neighborhood is less than the first preset threshold, and retain points whose height is greater than the second preset threshold as treetop points , as described For marking, Perform marker-controlled watershed segmentation.

[0041] Compared with the prior art, the present invention has the following advantages and effects:

[0042] (1) By using airborne laser radar to penetrate the forest canopy to obtain vertical structural information, a three-dimensional canopy point cloud model (3D-CPM) was constructed, which not only accurately depicts the canopy surface morphology, but also further expresses the vertical distribution structure inside the canopy, greatly improving the accuracy and detail of forest tree canopy structure analysis;

[0043] (2) The point cloud processing strategy that combines adaptive terrain filtering with multi-feature fusion classification can filter out ground points in complex forest environment data with high precision. At the same time, multi-feature fusion is used to perform semantic classification of ground objects, which significantly improves the accuracy and robustness of point cloud classification and ensures the data accuracy of subsequent tree parameter extraction to a high degree.

[0044] (3) By precisely locating the spatial position of each tree and combining it with airborne lidar data, it is possible to accurately calculate the key forestry parameters of each tree and obtain tree stand level indicators per hectare, greatly improving the accuracy of forest resource assessment and providing more reliable data support for forest management. BRIEF DESCRIPTION OF THE DRAWINGS

[0045] In order to more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the following is a brief introduction to the drawings required for the specific embodiments or the description of the prior art. In all the drawings, similar elements or parts are generally identified by similar reference numerals. In the drawings, the elements or parts are not necessarily drawn according to the actual scale.

[0046] Figure 1 It is a flow chart of a laser radar data processing method for a complex forest environment provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0047] In order to more clearly explain the purpose, technical solutions and advantages of the present invention, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all of the embodiments. The example implementation methods can be implemented in various forms and should not be understood as being limited to the examples described herein. On the contrary, these implementation methods are provided to make the present invention more comprehensive and complete, and to fully convey the concepts of the example implementation methods to those skilled in the art.

[0048] In addition, the described features, structures or characteristics may be combined in one or more embodiments in any suitable manner. In the following description, many specific details are provided to provide a full understanding of the embodiments of the present invention. However, those skilled in the art will appreciate that the technical solution of the present invention can be practiced without one or more of the specific details, or other methods, components, devices, steps, etc. may be adopted. In other cases, known methods, devices, implementations or operations are not shown or described in detail to avoid blurring various aspects of the present invention.

[0049] The block diagrams shown in the accompanying drawings are merely functional entities and do not necessarily correspond to physically independent entities. That is, these functional entities may be implemented in software form, or in one or more hardware modules or integrated circuits, or in different networks and / or processor devices and / or microcontroller devices.

[0050] The flowcharts shown in the accompanying drawings are only exemplary and do not necessarily include all the contents and operations / steps, nor must they be executed in the order described. For example, some operations / steps can be decomposed, and some operations / steps can be combined or partially combined, so the actual execution order may change according to actual conditions.

[0051] The present invention is described in detail below in conjunction with specific embodiments:

[0052] As attached Figure 1As shown, an embodiment of the present invention provides a laser radar data processing method for a complex forest environment, and the specific steps of the method are as follows:

[0053] Step S1, obtaining original point cloud data of a target forest area through an airborne laser radar system or a backpack laser radar system, and preprocessing the original point cloud data to obtain preprocessed point cloud data;

[0054] Step S2, filtering out ground points in the preprocessed point cloud data by an adaptive terrain filtering method to obtain non-ground point cloud data, and classifying the ground objects in the non-ground point cloud data by a multi-feature fusion classification method to obtain ground object classification point cloud data;

[0055] Step S3: construct a three-dimensional canopy point cloud model based on the classified point cloud data, and obtain the single tree spatial position and key forestry parameters of individual trees in the target forest area based on a single tree parameter extraction method.

[0056] Preferably, the airborne laser radar system in step S1 is mounted on a flying platform, and scans the target forest area within a preset relative flying altitude to obtain the raw point cloud data in a large range; the backpack laser radar system is integrated into a portable backpack device, and the workers collect the close-range high-density raw point cloud data by walking in the forest area; the data collection output of the target forest area includes: the raw point cloud data, a large number of discrete three-dimensional points and their attribute information and trajectory data, and the carrier platform position and posture information when the raw point cloud data is obtained.

[0057] Specifically, the flight platform includes a fixed-wing aircraft or a helicopter, and the preset relative flight altitude includes 700-1000m.

[0058] Preferably, the preprocessing of the original point cloud data in step S1 includes noise point removal, point cloud registration and point cloud segmentation, specifically:

[0059] The noise point removal includes statistical filtering, distance filtering, etc., combined with the distribution characteristics of noise points in complex forest scenes, to adaptively adjust the filtering threshold;

[0060] The point cloud registration achieves accurate registration by combining corner point features and normal vector features, thereby achieving a multi-view point cloud registration effect;

[0061] The point cloud segmentation adopts the octree space partitioning method to store the global point cloud in blocks, and divides the point cloud into point cloud data blocks with consistent scale and hierarchical organization. Each point cloud data block has a unique octree index number for fast retrieval and calling of the point cloud data block.

[0062] Preferably, in step S2, the adaptive terrain filtering method specifically includes:

[0063] Step S201: Grid division is performed on the pre-processed point cloud data in the XY direction. The size is d×d, and the initial grid set is constructed , the initial grid set All of the grids are included ; The d×d is preferably 5m×5m;

[0064] Step S202: In the grid Select the one with the lowest elevation value. points as initial ground seed points , construct a seed point set , the seed point set All the initial ground seed points are included ; The n is preferably 10;

[0065] Step S203: Based on the seed point set Construct a Delaunay triangulation and calculate the normal vector of each triangle in the Delaunay triangulation and roughness coefficient ;

[0066] Specifically, the normal vector for: ,in , are the two edge vectors of the triangular patch; the roughness coefficient for: ,in , are the elevation values ​​and the mean elevation values ​​of the three vertices of the triangle;

[0067] Step S204: based on the normal vector and the roughness coefficient , adaptively calculate the elevation threshold of each triangle patch and angle threshold ;

[0068] Specifically, the elevation threshold for , the angle threshold for ,in is an empirical parameter, which is determined according to the different types of objects and their physical characteristics in the triangular facet; in this embodiment, it is preferred ;

[0069] More specifically, since the physical characteristics of different types of ground objects, such as height, density, leaf shape, etc., have different effects on the pre-processed point cloud data, in this embodiment, the empirical parameters can also be selected as follows:

[0070] Large trees (such as pine and fir): a=1.5, b=0.3, c=π / 6, d=π / 4;

[0071] Shrubs (such as bushes): a=1.2, b=0.2, c=π / 8, d=π / 6;

[0072] Herbs (such as grassland, wheat): a=0.8, b=0.1, c=π / 12, d=π / 8.

[0073] Step S205: traverse each of the triangular facets, and the non-initial ground seed points in the neighborhood of the triangular facets The seed point is recorded as a non-seed point, and the height difference between the non-seed point and the triangle is calculated. Angle between the normal vector and , if the non-seed point satisfies and , then the non-seed point is marked as a candidate ground point and added to the seed point set ; The size of the neighborhood is preferably 0.5m;

[0074] Step S206, iteratively execute steps S203 to S205 until the seed point set When the number of iterations is no longer increased or reaches the preset number, the seed point set The output is the filtered ground point; the preset number of iterations is preferably 10 times;

[0075] It should be noted that, through the adaptive terrain filtering method, the pre-processed point cloud data is divided into ground points and non-ground points, and the ground points in the pre-processed point cloud data, i.e., the seed point set, are filtered out. , and obtain non-ground point cloud data; on the basis of obtaining the ground points, further identify the ground object categories to which the non-ground points belong, such as trees, buildings, etc.

[0076] Preferably, the multi-feature fusion classification method in step S2 specifically includes:

[0077] Step S207: Calculate the geometric features of each laser point in each non-ground point cloud data based on the local neighborhood information of the point cloud, wherein the geometric features include: normal vector , curvature and characteristic entropy ; The neighborhood information is preferably a KNN neighborhood or a voxel neighborhood;

[0078] Step S208, performing distance correction and filtering denoising preprocessing on the reflection intensity information of the laser points in the non-ground point cloud data;

[0079] It should be noted that, in this embodiment, considering that different ground objects have different reflection characteristics to laser energy, the reflection intensity information of the laser point is selected to distinguish the types of ground objects.

[0080] Step S209: based on the scale characteristics of different objects, extract the geometric features and radiation features of the laser point at multiple different neighborhood scales, and combine the geometric features and the radiation features into the final feature vector of the laser point. ; The field scales include 1m, 2m, and 5m;

[0081] Step S210: Use Z-score or Max-Min method to The dimensions of the features in are scaled and normalized to eliminate the differences in the dimensions of different features and obtain the normalized feature vector ;

[0082] Step S211: in the non-ground point cloud data, typical ground object points are manually marked, the typical ground object points include ground, tree trunks, tree crowns, buildings, etc., and feature vectors of the typical ground object points are extracted as training samples. ,in is the category label of the typical feature point, representing the feature category.

[0083] Step S212: Select a classifier and convert the training samples The data are input into the classifier for model learning, and a mapping relationship between the feature vector and the ground object category in the classifier is established; the classifier includes random forest, SVM, and neural network.

[0084] Step S213: The normalized feature vector of each laser point is obtained by the trained classifier. Predicting is performed to obtain the category label of each laser point, obtaining the classified point cloud data of the ground object, and post-processing and optimizing the classified point cloud data of the ground object by using methods such as graph cuts and conditional random fields;

[0085] It should be noted that through the above steps, the laser point cloud in complex forest scenes can be accurately divided into different semantic categories such as ground, trees, buildings, etc., with an overall classification accuracy of more than 90%; key parameters in the solution such as terrain filter threshold and classifier parameters need to be repeatedly debugged according to specific data sets and application requirements in order to achieve optimal performance.

[0086] Preferably, in step S3, a three-dimensional canopy point cloud model is constructed, including: the input of the three-dimensional canopy point cloud model is the ground object classification point cloud data, the point cloud data in the ground object classification point cloud data is divided into a two-dimensional grid in the XY plane, and then the voxels in the Z direction are divided into three-dimensional voxels, the point cloud density inside each three-dimensional voxel is quantified and normalized, the canopy surface height of the two-dimensional grid is extracted, and the three-dimensional canopy point cloud model is output, which is specifically as follows:

[0087] Step S301: the point cloud data is: , ,in Represents the i-th point cloud data point; Representing point cloud data The category label of the XY plane is , the voxel height in the Z direction is ; Project the point cloud data onto the XY plane, according to Divide the point cloud data range into grids to obtain a two-dimensional grid , , ,in and are the number of two-dimensional grids in the X and Y directions respectively, and satisfy: , ;

[0088] Step S302: In each of the two-dimensional grids Based on the Z-axis direction The voxel division in the Z direction is performed for the step size, and the two-dimensional grid is expanded into three-dimensional voxels ,in, is the number of voxels in the Z direction and satisfies: ;in, It means rounding up to ensure that the divided grids and voxels can cover the entire point cloud range;

[0089] Step S303: traverse each of the three-dimensional voxels , and count the point cloud data inside it Category label , if exists (i.e., tree category), the three-dimensional voxel Mark as , otherwise marked as ;

[0090] Step S304: By attribute value Quantify each of the Internal point cloud density: ,in, is an indicative function, when the point cloud data Category label When it is 1, it is 1, otherwise it is 0; for the three-dimensional voxel that does not contain the point , attribute value Recorded as 0;

[0091] Step S305: The attribute value of Normalize in the range 0-1: ,in and Respectively all the of The minimum and maximum values ​​of , which facilitates comparative analysis between different point cloud data sets;

[0092] Step S306: In each of the two-dimensional grids Within, extract the maximum of , the The Z value of the two-dimensional grid The canopy surface height : ,make , based on all the two-dimensional grids The canopy surface height Construct the canopy height model, denoted as ;

[0093] Step S307: Output the three-dimensional canopy point cloud model as follows: , and the canopy height model: ;

[0094] It should be noted that the The quality of the construction depends largely on whether the attribute information of the point cloud is complete. Ideally, the input point cloud should contain ground / non-ground classification information. , Tree / non-tree classification information And echo type information etc.; among them It can be obtained by the adaptive terrain filtering method mentioned above in the embodiment, and Reasonable annotation is required during the point cloud processing stage.

[0095] Preferably, the single tree parameter extraction method described in step S3 is based on crown segmentation and allometric growth equations to obtain the single tree spatial position and key forestry parameters of individual trees in the target forest area, specifically as follows:

[0096] Step S308: for the three-dimensional canopy point cloud model: , the canopy height model: And the digital elevation model: Enter data and set parameters as: Crown segmentation parameters , which are the minimum and maximum thresholds of the crown area, respectively, used to eliminate unreasonable segmentation results; smoothing threshold: SmoothThreshold, used to control the smoothness of the crown boundary; allometric growth model: HeightDBHModel, which is the tree height-diameter at breast height growth model.

[0097] Step S309: Search for the local maximum point in the neighborhood, extract all the pixels with the largest local height and mark them as , the search neighborhood size is larger than the single tree crown size, based on the Eliminate points whose relative height difference in the search neighborhood is less than the first preset threshold, and retain points whose height is greater than the second preset threshold as treetop points , as described For marking, Perform marker-controlled watershed segmentation;

[0098] Specifically, the values ​​of the first preset threshold and the second preset threshold depend on the forest stand density and canopy density. For example, if the height difference of trees in the forest is large, the scale of the single tree crown is about 5 meters, and the forest terrain is relatively flat, then:

[0099] The first preset threshold is set to 0.5 meters to 1 meter, which is used to remove some small-scale undulations (such as ground undulations or small shrubs) without removing the real treetops.

[0100] The second preset threshold is set to 2 meters; only those local maximum points that are significantly higher than the surrounding points are considered as treetop points, thus avoiding misidentification of low or obscured trees;

[0101] It should be noted that the first preset threshold and the second preset threshold are set for different levels of terrain features and tree height differences. The first preset threshold is used to remove noise points with local height differences less than the threshold, and the second preset threshold is used to screen out real treetop points. Through reasonable settings, the accuracy of treetop point extraction can be effectively improved, ensuring that the subsequent watershed segmentation effect is more accurate.

[0102] Step S310: Calculate the The gradient of the gradient map is Gradient, and the gradient reflects the local height difference: ,in , ; Treat the gradient map as a terrain surface, start from each mark, iterate along the direction of the steepest gradient descent until the marks converge to form a watershed boundary, and segment the gradient map; perform morphological filtering on the segmented gradient map to eliminate area or The area is smoothed and the boundary is obtained to obtain the final crown segmentation result, which is recorded as ,in represents the segmentation of the crown, ;

[0103] Step S311: In each of the segmented tree crowns The corresponding three-dimensional canopy point cloud model Search for the highest point in the Z direction in the three-dimensional voxel and extract the elevation value As the crown height; in the split crown In range Extract the average ground elevation , tree height Estimated to be: ; Based on the tree height-DBH measured sample data, the least squares regression method was used to establish the allometric growth equation suitable for the target tree species: ,in is the breast diameter (cm), is the tree height (m), and is the parameter to be fitted; Substitute the DBH into the calculation to estimate the DBH of each segmented crown. : ,in To estimate tree height;

[0104] Step S312: for each segmented tree crown Fit the minimum bounding rectangle, where the long side and short side Corresponding to the long and short diameters of the crown width, the average value is used as the estimated value of the crown width : ;

[0105] Step S313: based on the estimated tree height , estimate DBH , using the tree species volume model to estimate the volume of a single tree ; Exemplarily: ,in, is the volume of a single tree (m³), is the breast diameter (cm), is the tree height (m), is the parameter to be fitted, and the tree volume model is calibrated by the measured sample data. , The volume of each tree can be calculated: ;

[0106] Step S314: Calculate the number of trees per hectare , the sum of the cross-sectional areas of chest height , average breast diameter , average tree height Volume per hectare The specific calculation formula of the stand level index is as follows: , , , , where S is the plot area in hectares (ha), and T is the total number of split trees;

[0107] Step S315: Outputting a single tree parameter table and a forest stand parameter table, wherein the single tree parameter table , the stand parameter table ;

[0108] It should be noted that this embodiment realizes the parameter extraction of single tree and forest stand levels in a complex forest environment through the above-mentioned steps based on the three-dimensional canopy point cloud model 3D-CPM, and obtains important forest stand indicators such as the spatial position, tree height, breast diameter, crown width, and volume per hectare of each tree.

[0109] It should be noted that the existing methods are mainly based on the canopy height model (CHM) for crown segmentation and parameter estimation; CHM obtains the canopy surface height relative to the ground by extracting the difference between the normalized digital surface model (nDSM) and the digital elevation model (DEM), but in complex forest stands with high canopy density, CHM is often difficult to accurately reconstruct the canopy surface morphology due to factors such as canopy occlusion and edge loopholes, and cannot effectively express the vertical structure information of the canopy, resulting in insufficient accuracy in subsequent tree positioning and parameter estimation. In this embodiment, the advantages of airborne LiDAR data in penetrating the canopy and obtaining vertical structure information are fully utilized. By constructing a three-dimensional canopy point cloud model (3D-CPM), the canopy surface morphology is accurately depicted while further expressing the vertical distribution structure inside the canopy; on this basis, the spatial position of a single tree is obtained by combining the crown segmentation and the single tree parameter extraction method of the allometric equation, and at the same time, the key forestry parameters such as tree height, breast diameter, crown width, and timber volume of each tree are reliably estimated.

[0110] In this embodiment, the original point cloud data of the target forest area is obtained by an airborne laser radar system or a backpack laser radar system, and the preprocessed point cloud data is obtained by preprocessing. The ground points in the preprocessed point cloud data are filtered out to obtain non-ground point cloud data, and the ground objects in the non-ground point cloud data are classified to obtain ground object classification point cloud data; based on the ground object classification point cloud data, a three-dimensional canopy point cloud model is constructed, and based on the single tree parameter extraction method, the single tree spatial position and key forestry parameters of individual trees in the target forest area are obtained, and the single tree parameter extraction method includes crown segmentation and allometric growth equation; the advantages of airborne laser radar data are fully utilized, and innovative point cloud processing, modeling and parameter extraction methods are adopted to achieve high-precision acquisition of single tree and forest stand level parameters in complex forest environments, which can better solve the shortcomings of the existing CHM method, greatly improve the efficiency and accuracy of forest resource monitoring, and have important application value in forest resource management, forestry ecology research and other fields.

[0111] Those skilled in the art will readily appreciate other embodiments of the present invention after considering the specification and practicing the invention disclosed herein. This application is intended to cover any modification, use or adaptation of the present invention, which follows the general principles of the present invention and includes common knowledge or customary techniques in the art that are not disclosed by the present invention. The specification and examples are intended to be exemplary only, and the true scope and spirit of the present invention are indicated by the claims. It should be understood that the present invention is not limited to the precise structure described above and shown in the accompanying drawings, and various modifications and changes can be made without departing from the scope thereof. The scope of the present invention is limited only by the appended claims.

Claims

1. A laser radar data processing method for a complex forest environment, characterized in that: include: Acquire original point cloud data of the target forest area by using an airborne laser radar system or a backpack laser radar system, and preprocess the original point cloud data to obtain preprocessed point cloud data; Filtering out ground points in the preprocessed point cloud data by an adaptive terrain filtering method to obtain non-ground point cloud data, and classifying ground objects in the non-ground point cloud data by a multi-feature fusion classification method to obtain ground object classification point cloud data; The adaptive terrain filtering method includes dividing the pre-processed point cloud data into ground points and non-ground points, and filtering out the ground points in the pre-processed point cloud data to obtain non-ground point cloud data; Based on the non-ground point cloud data, further identifying the ground object category to which the non-ground point belongs; The adaptive terrain filtering method comprises: gridding the pre-processed point cloud data, constructing a Delaunay triangulation network, and calculating the normal vector and roughness coefficient of each triangular face in the Delaunay triangulation network; The multi-feature fusion classification method includes predicting the feature vectors of the laser points of the non-ground point cloud data by a trained classifier to obtain a category label of each of the laser points; Based on the classified point cloud data of the ground objects, a three-dimensional canopy point cloud model is constructed, and based on a single tree parameter extraction method, the single tree spatial position and key forestry parameters of individual trees in the target forest area are obtained, wherein the single tree parameter extraction method includes canopy segmentation and allometric growth equation; The construction of the three-dimensional canopy point cloud model comprises: the input of the three-dimensional canopy point cloud model is the ground object classification point cloud data, the point cloud data in the ground object classification point cloud data is divided into a two-dimensional grid in the XY plane, and then the voxels in the Z direction are divided into three-dimensional voxels, the point cloud density inside each three-dimensional voxel is quantified and normalized, the canopy surface height of the two-dimensional grid is extracted, and the three-dimensional canopy point cloud model is output; The construction of the three-dimensional canopy point cloud model is specifically as follows: The point cloud data is: , ,in Represents the i-th point cloud data point; Representing point cloud data The category label of the XY plane is , the voxel height in the Z direction is ; Project the point cloud data onto the XY plane, according to Divide the point cloud data range into grids to obtain a two-dimensional grid , , ,in and are the number of two-dimensional grids in the X and Y directions respectively, and satisfy: , ; In each of the two-dimensional grids Based on the Z-axis direction The voxel division in the Z direction is performed for the step size, and the two-dimensional grid is expanded into three-dimensional voxels ,in, is the number of voxels in the Z direction and satisfies: ; Traverse each of the three-dimensional voxels , and count the point cloud data inside it Category label , if exists , that is, the point whose category is tree, then the three-dimensional voxel Mark as , otherwise marked as .

2. The laser radar data processing method for a complex forest environment according to claim 1, characterized in that: The airborne laser radar system is mounted on a flying platform, and scans the target forest area within a preset relative flight altitude to obtain the original point cloud data in a large range; the backpack laser radar system is integrated into a portable backpack device, and the high-density original point cloud data of the close range is collected by operators walking in the forest area.

3. The laser radar data processing method for a complex forest environment according to claim 1, characterized in that: The grid division of the pre-processed point cloud data is specifically as follows: The pre-processed point cloud data is divided into grids in the XY direction. The size is d×d, constructing the initial grid set , the initial grid set All of the grids are included ; The construction of the Delaunay triangulation network is specifically as follows: In the grid Select the one with the lowest elevation value. points as initial ground seed points , construct a seed point set , the seed point set All the initial ground seed points are included ; Based on the seed point set Construct a Delaunay triangulation and calculate the normal vector of each triangle in the Delaunay triangulation and roughness coefficient ; Based on the normal vector and the roughness coefficient , adaptively calculate the elevation threshold of each triangle patch and angle threshold .

4. The laser radar data processing method for a complex forest environment according to claim 3, characterized in that: The adaptive terrain filtering method further comprises: Traverse each of the triangular facets, and find the non-initial ground seed points in the neighborhood of the triangular facets. The seed point is recorded as a non-seed point, and the height difference between the non-seed point and the triangle is calculated. Angle between the normal vector and , if the non-seed point satisfies and , then the non-seed point is marked as a candidate ground point and added to the seed point set .

5. The laser radar data processing method for a complex forest environment according to claim 1, characterized in that: The multi-feature fusion classification method specifically includes: Based on the local neighborhood information of the point cloud, the geometric features of each laser point in each of the non-ground point cloud data are calculated, and the geometric features include: normal vector , curvature and characteristic entropy ; Performing distance correction and filtering denoising preprocessing on the reflection intensity information of the laser points in the non-ground point cloud data; Based on the scale characteristics of different objects, the geometric features and radiation features of the laser point at multiple different neighborhood scales are extracted, and the geometric features and the radiation features are combined into the final feature vector of the laser point. ; The feature vector is analyzed by using Z-score or Max-Min method. The dimensions of the features in are scaled and normalized to eliminate the differences in the dimensions of different features and obtain the normalized feature vector ; In the non-ground point cloud data, typical ground object points are manually annotated, and the typical ground object points include ground, tree trunks, tree crowns, and buildings, and feature vectors of the typical ground object points are extracted as training samples. ,in is the category label of the typical feature point, representing the feature category.

6. The laser radar data processing method for a complex forest environment according to claim 5, characterized in that: The multi-feature fusion classification method also includes: Select a classifier and pass the training samples Inputting into the classifier for model learning, and establishing a mapping relationship between the feature vector and the ground object category in the classifier; The normalized feature vector of each laser point is obtained by the trained classifier. Prediction is performed to obtain the category label of each laser point and obtain the ground object classification point cloud data.

7. The laser radar data processing method for a complex forest environment according to claim 1, characterized in that: The constructing of the three-dimensional canopy point cloud model also includes: By property value Quantify each of the Internal point cloud density: ,in, is an indicative function, when the point cloud data Category label When it is 1, it is 1, otherwise it is 0; for the three-dimensional voxel that does not contain the point , attribute value Recorded as 0; Regarding the The attribute value of Normalize in the range 0-1: ,in and Respectively all the of The minimum and maximum values ​​of In each of the two-dimensional grids Within, extract the maximum of , the The Z value of the two-dimensional grid The canopy surface height : ,make , based on all the two-dimensional grids The canopy surface height Construct the canopy height model, denoted as .

8. The laser radar data processing method for complex forest environment according to claim 7, characterized in that: The tree crown segmentation comprises: In the Search for the local maximum point in the neighborhood, extract all the pixels with the largest local height and mark them as , the search neighborhood size is larger than the single tree crown size, based on the Eliminate points whose relative height difference in the search neighborhood is less than the first preset threshold, and retain points whose height is greater than the second preset threshold as treetop points , as described For marking, Perform marker-controlled watershed segmentation.

Images