Tailing pond dry beach length calculation method based on point cloud data
By acquiring 3D point cloud data of tailings dams using UAV lidar, and processing and generating facet pixel maps using the octree algorithm, the problems of time-consuming, labor-intensive, and inaccurate monitoring of dry beach length in existing tailings dams have been solved, achieving efficient and reliable monitoring results.
Patent Information
- Application Number
- CN202411734481.9
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2024-11-29
- Publication Date
- 2025-11-18
- Estimated Expiration
- 2044-11-29
AI Technical Summary
Existing methods for monitoring the dry beach length of tailings ponds rely on manual monitoring and image processing, which are time-consuming, labor-intensive, and lack accuracy and reliability.
A multi-rotor aerial surveying drone equipped with a lidar was used to acquire three-dimensional point cloud data. A three-dimensional point cloud map was generated using image processing software. The point cloud data was then processed using an octree algorithm to generate a surface voxel map, which was used for water surface segmentation and dam crest straight line determination to calculate the length of the dry beach.
It achieves highly reliable, accurate, and efficient calculation of tailings dam dry beach length, reducing errors caused by human intervention and improving the accuracy and efficiency of monitoring.
Smart Images

Figure CN119559350B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of mining technology, specifically relating to a method for calculating the length of dry beach in tailings dams based on point cloud data. Background Technology
[0002] Tailings dams are sites used to store tailings or other industrial waste discharged after ore beneficiation in metal or non-metal mines. They are essential infrastructure and environmental protection projects for mining enterprises. The high potential energy of tailings dams poses a potential danger of man-made debris flows at mining sites. A tailings dam failure could potentially trigger a major safety accident; therefore, monitoring tailings dams has always been a top priority for mining enterprises.
[0003] Tailings dams primarily intercept solid tailings and water from the tailings slurry through tailings dams (tailings dams consist of an initial dam and a stockpiling dam; the initial dam is the dam used during infrastructure construction to support the later tailings stockpile; the stockpiling dam is the dam constructed by filling and piling tailings above the initial dam during production). As the tailings slurry is gradually discharged, the water in the slurry gradually sinks and accumulates, while the solid tailings in the slurry form sediments on the interception surface of the tailings dam due to hydraulic impact. In the field of tailings dam safety technology, the surface portion of this sedimentary body exposed above the water surface is usually called the sedimentary beach or dry beach; the highest point of the dry beach is located at the intersection of the dry beach and the outer slope of the tailings dam, which is usually called the beach crest; the horizontal distance between the dry beach crest and the waterline inside the tailings dam is usually called the dry beach length. The dry beach length is an important indicator for measuring the safe operation of a tailings dam: if the dry beach length exceeds the standard, it may lead to safety accidents such as water level overflow and tailings dam failure. Therefore, the dry beach length of tailings dams must be monitored.
[0004] Currently, the monitoring methods for the dry beach length of tailings dams still rely on manual monitoring or simple image monitoring. Manual monitoring involves personnel arriving on-site to manually measure the dry beach length; however, this method is not only time-consuming and labor-intensive but also has poor accuracy. Image monitoring uses drones to acquire images of the dry beach portion of the tailings dam and then uses simple image processing to determine the dry beach length; however, this method relies solely on images, which are greatly affected by weather and site conditions, resulting in low accuracy and reliability. Summary of the Invention
[0005] The purpose of this invention is to provide a highly reliable, accurate, and efficient method for calculating the dry beach length of tailings ponds based on point cloud data.
[0006] The method for calculating the dry beach length of a tailings dam based on point cloud data provided by this invention includes the following steps:
[0007] S1. Obtain the three-dimensional point cloud data information of the target tailings dam;
[0008] S2. Process the 3D point cloud data obtained in step S1 to generate a 3D point cloud map and mark the straight line on the dam crest;
[0009] S3. Generate a voxel map based on the 3D point cloud map obtained in step S2;
[0010] S4. On the face pixel map obtained in step S3, perform water surface segmentation to obtain water surface outline and water surface height data;
[0011] S5. Using the straight line data of the dam crest determined in step S2, determine the dam crest plane and water surface outline, and calculate the dry beach length of the target tailings dam.
[0012] Step S1, which involves obtaining the three-dimensional point cloud data information of the target tailings dam, specifically includes the following steps:
[0013] A multi-rotor aerial surveying drone, equipped with lidar, was used to collect data on the target tailings dam and obtain three-dimensional point cloud data information of the target tailings dam.
[0014] Step S2 involves processing the 3D point cloud data obtained in step S1 to generate a 3D point cloud map and marking the straight lines on the dam crest. Specifically, this includes the following steps:
[0015] Using image processing software, the 3D point cloud data acquired at each moment in step S1 is registered according to the position and attitude of the UAV to generate a 3D point cloud map.
[0016] The dam crest plane is determined based on the straight line at the dam crest, the direction of gravity, and the direction perpendicular to the straight line at the dam crest; the straight line at the dam crest is determined by point O on the straight line and the direction vector. express.
[0017] Step S3, which involves generating a voxel map based on the 3D point cloud map obtained in step S2, includes the following steps:
[0018] A. Traverse the point cloud list;
[0019] B. Handling surface features and non-surface features;
[0020] C. Initialize all octrees;
[0021] D. Merging the planes yields the final facet matrix map.
[0022] Step A, which involves traversing the point cloud list, specifically includes the following steps:
[0023] The 3D point cloud map obtained in step S2 is a list of distorted point clouds in the world coordinate system;
[0024] Iterate through each point PV in the list of distorted point clouds:
[0025] Calculate the voxel position: Divide each coordinate of point pv by the maximum voxel size and round the result to the nearest integer to obtain the position of point pv in the voxel gradient, VoxeLoc.
[0026] Step B, which involves processing surface features and non-surface features, specifically includes the following steps:
[0027] The position of point pv obtained in step A in the voxel gradient is determined as follows:
[0028] If point pv is a face feature, then perform the following steps:
[0029] Check if the key-value pair set includes the position VoxeLoc: if it does, cache the point pv in the root node of the corresponding octree; if it does not include the position VoxeLoc, create a new octree structure and load the point pv into the root node of the new octree.
[0030] If point pv is a non-surface feature, then perform the following steps:
[0031] If an anomaly is detected in the input list of distorted point clouds, the anomaly information is recorded, the anomaly points are marked as invalid, and:
[0032] Record the coordinates, timestamps, and corresponding attribute information of the anomalies;
[0033] If the algorithm can continue running, skip the abnormal point and continue processing the next point in the point cloud list;
[0034] If the number of outliers exceeds the set threshold, or the algorithm cannot continue to run, perform the following operations: send a prompt signal, output an error message, and terminate the current algorithm;
[0035] Step C, which involves initializing all octrees, specifically includes the following steps:
[0036] After processing all points with both surface and non-surface features, octree initialization is performed on all octree structures in the key-value pair set to ensure that the internal nodes and leaf nodes of the octree are processed and that the octree correctly represents the spatial structure in the voxel map.
[0037] After the octree is initialized, the initialized face pixel map is obtained.
[0038] The octree initialization specifically includes the following steps:
[0039] Check the number of points in the input point cloud: Check if the number of points in the input point cloud exceeds a first preset threshold: If the number of points exceeds the first preset threshold, perform plane fitting to find the optimal plane containing the point cloud; if the number of points does not exceed the first preset threshold, proceed with the following steps:
[0040] Directly assign to leaf nodes: The point cloud is directly assigned to the corresponding octagonal leaf nodes, and the corresponding leaf nodes are marked as occupied to indicate that the corresponding leaf nodes contain point cloud data;
[0041] Record sparse regions: Record the currently processed region as a sparse region in the log file, and record the region information, including the coordinate range and number of points of the region;
[0042] Merging Neighboring Nodes: Check if the neighboring nodes of the current leaf node are also sparse nodes. If the neighboring nodes are equally sparse, merge them into a larger node to reduce the total number of nodes in the branch tree and optimize the spatial structure. If the neighboring nodes are not sparse, keep the current leaf node unchanged and record the state information of the corresponding neighboring nodes, including the node's data density, access frequency, and space occupation. Subsequently, check each of the other neighboring nodes of the current leaf node to evaluate their sparsity and merging potential. If the other neighboring nodes are equally sparse, merge the current leaf node with the other sparse neighboring nodes. If the other neighboring nodes are not sparse, mark them as "not to be merged" and monitor the corresponding data changes. At the same time, adjust the sparsity judgment criteria globally.
[0043] Apply default parameters: Apply default voxel parameters, including voxel size and density, to the sparse point cloud region; ensure that the voxel map of the sparse point cloud region is consistent with the voxel maps of other regions. Figure 1 To;
[0044] Anomaly Handling: If the number of points is less than the set threshold, it indicates that the point cloud cannot form a voxel structure that meets the set requirements. The corresponding point cloud is then marked as an anomaly point; detailed information about the anomaly point is recorded.
[0045] Continue the initialization process: After processing the current node, continue to initialize other octree nodes; ensure that all nodes have been processed to guarantee the integrity and accuracy of the octree structure;
[0046] Perform plane fitting: If plane fitting is successful and the number of points does not exceed the second set threshold, then stop further segmentation; if plane fitting is successful and the number of points exceeds the second set threshold, then the corresponding voxel will not be updated; if plane fitting fails, then perform octree segmentation until a plane is successfully generated or the set number of layers is reached.
[0047] Finally, the root node of the octree is output, which includes the data information of all point clouds.
[0048] The octree partitioning described herein specifically includes the following steps:
[0049] The input data for octree partitioning is: the root node of the octree to be partitioned;
[0050] The output data of octree partitioning is: an array of eight child nodes generated after partitioning, each child node representing a voxel;
[0051] The segmentation process is as follows:
[0052] Obtain point cloud data: Extract point cloud data from the input root node;
[0053] Check if the bottom layer has been reached: Determine if the root node has reached the maximum depth of the octree: If the maximum depth has been reached, no more nodes will be created; if the maximum depth has not been reached, proceed with the next steps.
[0054] Initialize child node array: Initialize a child node array of size 8, where each element in the child node array corresponds to a child node of the octree;
[0055] Assign point cloud data to child nodes: Traverse the point cloud data items in the extracted point cloud data, determine the position of the point cloud in the octree based on the voxel center and data item, and assign the point cloud data to the corresponding child node array;
[0056] Plane fitting for child nodes: For each child node in the child node array, determine whether the number of points in each child node is greater than the third set threshold. If it is greater than the third set threshold, then plane fitting is performed. If plane fitting is successful, then the current child node will no longer be octree split and the octree split of the next child node will begin. If plane fitting fails, then the current child node will continue to be octree split.
[0057] Returns an array of child nodes: After processing all child nodes, the final array of child nodes is obtained.
[0058] The aforementioned plane fitting specifically includes the following steps:
[0059] Input for plane fitting: point cloud data; the point cloud data includes a set of three-dimensional coordinate points, representing the point cloud data for which plane fitting is required;
[0060] The output of plane fitting is the parameters of the fitted plane; these parameters include the normal vector and centroid of the fitted plane.
[0061] The fitting process is as follows:
[0062] Calculate the centroid: Calculate the centroid of the input point cloud data: the coordinates of the centroid in each dimension, which is the average of the coordinates of all point cloud data in the corresponding dimension;
[0063] Calculate the covariance matrix: Calculate the covariance matrix of the input point cloud data;
[0064] Eigenvalue decomposition: Perform eigenvalue decomposition on the obtained covariance matrix to obtain the corresponding eigenvalues and eigenvectors corresponding to the eigenvalues;
[0065] Selecting the normal vector: Among the eigenvalues, select the eigenvector corresponding to the smallest eigenvalue as the normal vector of the fitting plane;
[0066] Successful Fitting Determination: The minimum eigenvalue is evaluated: if the minimum eigenvalue is less than the fourth set threshold, the plane fitting is considered successful, and the corresponding normal vector and centroid are returned; if the minimum eigenvalue is greater than or equal to the fourth set threshold, the plane fitting is considered unsuccessful, and the following steps are performed:
[0067] Record failure information: Record detailed information about the current point cloud data, including the number of points, centroid coordinates, covariance matrix, and eigenvalues; display a message indicating fitting failure;
[0068] Dimensionality reduction: The point cloud data is dimensionality reduced to reduce the impact of noise; the covariance matrix and eigenvalues of the dimensionality-reduced point cloud data are recalculated, and plane fitting is performed again;
[0069] Point cloud segmentation: Divide the current point cloud data into several subsets, and perform planar fitting on each subset; if a subset is successfully fitted, the fitting result of that subset is retained; if all subsets fail to fit, continue with subsequent processing.
[0070] Apply alternative fitting algorithms: Use the set alternative plane fitting algorithms to increase the probability of successful fitting; record the fitting results of the alternative plane fitting algorithms and compare them with the current results;
[0071] Mark as an abnormal region: If the number of fitting failures exceeds the set number, the point cloud data of the corresponding region will be marked as an abnormal region; the corresponding abnormal region will be marked in the voxel map;
[0072] Continue processing other regions: Skip the point cloud data that failed to fit at the moment and continue to perform plane fitting on the point cloud data of other regions; ensure that the algorithm process is not affected by the failure of fitting.
[0073] Generate an error report: The error report includes the point cloud data that failed to fit, the attempted fitting methods, and an analysis of the reasons for the failure; send the error report;
[0074] Return failure status: Returns a status indicating that the fit has failed at the current stage of the algorithm.
[0075] The merging of planes described in step D specifically includes the following steps:
[0076] The disjoint-set data structure algorithm is used to merge the planes:
[0077] initialization:
[0078] Create a disjoint-set data structure. The number of elements in the disjoint-set data structure is equal to the number of points on the plane.
[0079] Initialize each point as a separate set;
[0080] Perform edge scanning:
[0081] Traverse all edges in the plane;
[0082] For each edge, find the representative node of the set containing the two connected points, and determine: if the representative nodes are different, merge the two sets; if the representative nodes are the same, proceed with the following steps:
[0083] Ignore duplicate edges: Record the current edge as processed to avoid duplicate calculations; skip the current edge and continue processing the next edge;
[0084] Update edge information: The edge information includes the edge length and the number of connected vertices; ensure that the information is not updated repeatedly due to repeated processing of the same edge;
[0085] Maintain the internal structure of the set: If the internal structure of the set needs to be maintained, check and maintain the internal structure of the set; ensure that the internal structure of the set after the merge operation meets the set requirements;
[0086] Mark processing status: Record that the current edge is an internal edge and ensure that no more merging operations will be performed; update the merging status;
[0087] Optimize data structure: If the number of edges representing the same node exceeds a set value, optimize the data structure to reduce the number of search and merge operations;
[0088] To determine whether two points are in the same set, follow these steps:
[0089] Obtain the representative nodes for each of the two points, and then determine whether the representative nodes are the same: if the representative nodes are the same, it means that the two points are in the same set; if the representative nodes are different, it means that the two points are not in the same set.
[0090] Step S4 involves segmenting the water surface on the face pixel map obtained in step S3 to obtain the water surface outline and water surface height data. Specifically, this includes the following steps:
[0091] On the voxel map, all voxels are coarsely classified:
[0092] For each voxel, based on the normal of the internal fitting plane, determine the adjacent voxels reached by the plane extension from all adjacent voxels, and determine the occupancy status of the adjacent voxels: if all adjacent voxels are occupied, the voxel is located inside the point cloud map corresponding to the surface voxel map; otherwise, the voxel is located at the boundary of the point cloud map corresponding to the surface voxel map.
[0093] For each point in the voxel, search for K nearest neighbors, project the K nearest neighbors along the normal of the fitted plane, and make a judgment: if the normal projections are all biased to one side, the point is determined to be a boundary point; if the normal projections are evenly distributed, the point is determined to be an interior point.
[0094] Finally, the set of water surface boundary points is obtained, the water surface is segmented, and the corresponding water surface contour and water surface height data are obtained.
[0095] Step S5, which uses the straight line data of the dam crest determined in step S2 to determine the dam crest plane and water surface outline, and calculates the dry beach length of the target tailings dam, specifically includes the following steps:
[0096] Based on the dam crest straight line data determined in step S2, determine the dam crest plane;
[0097] Project the set of water surface boundary points obtained in step S4 onto the dam crest plane and calculate the projection distance;
[0098] The minimum projected distance is used as the minimum dry beach length of the target tailings dam to complete the calculation of the dry beach length of the tailings dam based on point cloud data.
[0099] The tailings dam dry beach length calculation method provided by this invention obtains the three-dimensional point cloud data of the target tailings dam and performs facet pixel map generation, water surface segmentation, projection and calculation on the obtained three-dimensional point cloud data. This method not only realizes the calculation of the tailings dam dry beach length based on point cloud data, but also has higher reliability, better accuracy and higher efficiency. Attached Figure Description
[0100] Figure 1 This is a schematic diagram of the method flow of the present invention. Detailed Implementation
[0101] like Figure 1 The diagram shown is a flowchart of the method of the present invention: The method for calculating the dry beach length of a tailings dam based on point cloud data disclosed in this invention includes the following steps:
[0102] S1. Obtain the 3D point cloud data of the target tailings dam; specifically including the following steps:
[0103] A multi-rotor drone equipped with lidar can be used to collect data on the target tailings dam and obtain three-dimensional point cloud data information of the target tailings dam.
[0104] S2. Process the 3D point cloud data obtained in step S1 to generate a 3D point cloud map and mark the straight lines on the dam crest; specifically including the following steps:
[0105] Using image processing software, the 3D point cloud data acquired at each moment in step S1 is registered according to the position and attitude of the UAV to generate a 3D point cloud map.
[0106] The dam crest plane is determined based on the straight line at the dam crest, the direction of gravity, and the direction perpendicular to the straight line at the dam crest; the straight line at the dam crest is determined by point O on the straight line and the direction vector. express;
[0107] S3. Generate a voxel map based on the 3D point cloud map obtained in step S2; including the following steps:
[0108] A. Traverse the point cloud list; specifically including the following steps:
[0109] The 3D point cloud map obtained in step S2 is a list of distorted point clouds in the world coordinate system;
[0110] Iterate through each point PV in the list of distorted point clouds:
[0111] Calculate the voxel position: Divide each coordinate of point pv by the maximum voxel size and round the result to the nearest integer to obtain the position of point pv in the voxel gradient, VoxeLoc.
[0112] B. Processing surface features and non-surface features; specifically including the following steps:
[0113] Surface features refer to certain points in a point cloud that represent a portion of a plane or curved surface;
[0114] The position of point pv obtained in step A in the voxel gradient is determined as follows:
[0115] If point pv is a face feature, then perform the following steps:
[0116] Check if the key-value pair set includes the position VoxeLoc: If the position VoxeLoc is already included, it means that the current voxel already has a corresponding octree structure, so the point pv is cached in the root node of the corresponding octree. This is for subsequent octree construction and processing; if the position VoxeLoc is not included, it means that a new voxel has been encountered, so a new octree structure is created and the point pv is loaded into the root node of the new octree.
[0117] If point pv is a non-surface feature, then perform the following steps:
[0118] If the input list of distorted point clouds is found to be abnormal, it indicates that the input point cloud contains abnormal features that do not conform to the expected data type or conditions. The abnormal information is recorded, and the abnormal points are marked as invalid. Furthermore:
[0119] Record the coordinates, timestamps, and other relevant attribute information of the anomalies to the log file;
[0120] If the algorithm can continue running, skip the abnormal point and continue processing the next point in the point cloud list;
[0121] If the number of outliers exceeds a preset threshold, or the algorithm is unable to process so much outlier data, the following actions will be taken: send a warning signal; output an error message to standard error output or the user interface; or terminate the current algorithm execution flow.
[0122] C. Initialize all octrees; this includes the following steps:
[0123] An octree is a data structure commonly used for partitioning three-dimensional spatial data. When processing three-dimensional point cloud data, an effective data structure is often needed to manage and organize these data points in order to perform efficient searching, retrieval, and spatial analysis. An octree is such a data structure that divides three-dimensional space into eight subspaces, each of which can be further subdivided until a predetermined condition is met.
[0124] After processing all points with both surface and non-surface features, octree initialization is performed on all octree structures in the key-value pair set to ensure that the internal nodes and leaf nodes of the octree are processed and that the octree correctly represents the spatial structure in the voxel map.
[0125] After the octree is initialized, the initialized face pixel map is obtained;
[0126] The octree initialization specifically includes the following steps:
[0127] Check the number of points in the input point cloud: Check if the number of points in the input point cloud exceeds a first preset threshold: If the number of points exceeds the first preset threshold, perform plane fitting to find the optimal plane containing the point cloud; if the number of points does not exceed the first preset threshold, proceed with the following steps:
[0128] Directly assign to leaf nodes: The point cloud is directly assigned to the corresponding octagonal leaf nodes, and the corresponding leaf nodes are marked as occupied to indicate that the corresponding leaf nodes contain point cloud data;
[0129] Record sparse regions: Record the currently processed region as a sparse region in the log file, and record the region information, including the region's coordinate range and number of points; if the system supports it, a warning message for sparse regions can be displayed in the user interface;
[0130] Merging Neighboring Nodes: The process checks if the neighboring nodes of the current leaf node are also sparse. If they are, they are merged into a larger node to reduce the total number of nodes in the tree and optimize the spatial structure. If they are not sparse, the current leaf node remains unchanged, and its state information, including data density, access frequency, and space usage, is recorded in detail. Subsequently, other neighboring nodes of the leaf node are checked one by one to assess their sparsity and merging potential. If other neighboring nodes are found to be equally sparse, they are prioritized for merging. Non-sparse neighboring nodes are marked as "not to be merged" and their data changes are continuously monitored. Simultaneously, the sparsity judgment criteria are dynamically adjusted globally to more flexibly address the merging needs of different nodes. Throughout the process, the results of each check and evaluation are meticulously recorded for reference and analysis in subsequent global optimization phases. This approach ensures flexibility in processing the current leaf node while providing ample data support for future structural optimization.
[0131] Apply default parameters: Apply default voxel parameters to the sparse point cloud region, including voxel size and density, etc.; ensure that the voxel map of the sparse point cloud region is consistent with the voxel maps of other regions. Figure 1 To;
[0132] Anomaly Handling: If the number of points is less than the set threshold, it indicates that the point cloud cannot form a voxel structure that meets the set requirements. The corresponding point cloud is then marked as an anomaly. Detailed information about the anomaly is recorded (e.g., recorded in the anomaly handling log). Consider whether further processing or filtering of these anomalies is required.
[0133] Continue the initialization process: After processing the current node, continue to initialize other octree nodes; ensure that all nodes have been processed to guarantee the integrity and accuracy of the octree structure;
[0134] Plane fitting is performed: If the plane fitting is successful and the number of points does not exceed the second set threshold, it means that these point clouds are roughly located on a plane. At this point, further segmentation is stopped, which means that the so-called tree base has been reached, and it also means that the child nodes of this node will no longer be created. If the plane fitting is successful and the number of points exceeds the second set threshold, the corresponding voxels will no longer be updated in terms of plane. At the same time, in order to save memory and computing resources, the cache corresponding to the temporary point cloud used for plane fitting can be released. If the plane fitting fails, it means that the point cloud distribution is too complex and it is not suitable to be represented by a plane. Therefore, octree segmentation is performed until a plane is successfully generated or the set number of layers is reached.
[0135] Finally, the root node of the octree is output, which includes the data information of all point clouds;
[0136] The complexity of octree initialization depends on several factors, including the number and distribution of the input point cloud, and a predetermined threshold. The following is a further analysis of this complexity:
[0137] Time complexity: The time complexity mainly consists of two parts: plane fitting and octree partitioning. Plane fitting usually involves calculating the geometric properties of the point cloud, such as centroid and normals. The time complexity of these operations is usually O(n). Octree partitioning involves recursive partitioning and the creation of child nodes. The time complexity can reach O(n log n), where n is the number of points in the point cloud.
[0138] Space complexity: The space complexity mainly depends on the depth of the octree and the number of child nodes of each node; in the worst case, if every node needs to be split, the space complexity may be O(n), where n is the number of points in the point cloud.
[0139] Octree partitioning specifically includes the following steps:
[0140] The input data for octree partitioning is: the root node of the octree to be partitioned;
[0141] The output data of octree partitioning is: an array of eight child nodes generated after partitioning, each child node representing a voxel;
[0142] The segmentation process is as follows:
[0143] Obtain point cloud data: Extract point cloud data from the input root node;
[0144] Check if the bottom layer has been reached: Determine if the root node has reached the maximum depth of the octree: If the maximum depth has been reached, no more nodes will be created; if the maximum depth has not been reached, proceed with the next steps.
[0145] Initialize child node array: Initialize a child node array of size 8, where each element in the child node array corresponds to a child node of the octree;
[0146] Assign point cloud data to child nodes: Traverse the point cloud data items in the extracted point cloud data, determine the position of the point cloud in the octree based on the voxel center and data item, and assign the point cloud data to the corresponding child node array;
[0147] Plane fitting for child nodes: For each child node in the child node array, determine whether the number of points in each child node is greater than the third set threshold. If it is greater than the third set threshold, then plane fitting is performed. If plane fitting is successful, then the current child node will no longer be octree split and the octree split of the next child node will begin. If plane fitting fails, then the current child node will continue to be octree split.
[0148] Returns an array of child nodes: After processing all child nodes, the final array of child nodes is obtained;
[0149] The complexity of the octree partitioning process mainly depends on the quantity and distribution of the point cloud data, as well as the maximum depth of the octree. In terms of space complexity, the algorithm needs to store the information of all child nodes, which may reach O(n) in the worst case, where n is the number of points in the point cloud. In terms of time complexity, the algorithm involves multiple traversals and recursive partitioning of the point cloud data, and its complexity may reach O(n log n) or higher.
[0150] Plane fitting specifically includes the following steps:
[0151] When processing 3D point cloud data, plane fitting is used to find a plane. It aims to estimate a planar model from a set of scattered 3D point cloud data. This planar model is typically defined by a normal vector and a point (usually the centroid).
[0152] Input for plane fitting: point cloud data; the point cloud data includes a set of three-dimensional coordinate points, representing the point cloud data for which plane fitting is required;
[0153] The output of plane fitting is the parameters of the fitted plane; these parameters include the normal vector and centroid of the fitted plane.
[0154] The fitting process is as follows:
[0155] Calculate the centroid: Calculate the centroid of the input point cloud data: the coordinates of the centroid in each dimension, which is the average of the coordinates of all point cloud data in the corresponding dimension;
[0156] Calculate the covariance matrix: Calculate the covariance matrix of the input point cloud data; the covariance matrix is a symmetric matrix that describes the distribution characteristics of point cloud data, and it provides information on the variance and covariance of the point cloud in various directions;
[0157] Eigenvalue decomposition: Perform eigenvalue decomposition on the obtained covariance matrix to obtain the corresponding eigenvalues and eigenvectors; eigenvalue decomposition usually produces three eigenvalues.
[0158] Selecting the normal vector: Among the eigenvalues, select the eigenvector corresponding to the smallest eigenvalue as the normal vector of the fitting plane; the smallest eigenvalue usually represents the direction in which the point cloud data is most dispersed, that is, the variance is the smallest in that direction;
[0159] Successful Fitting Determination: The minimum eigenvalue is evaluated. If the minimum eigenvalue is less than the fourth set threshold, the plane fitting is considered successful, and the corresponding normal vector and centroid are returned. If the minimum eigenvalue is greater than or equal to the fourth set threshold, the plane fitting is considered unsuccessful. This may be because the point cloud data is too scattered to be approximated by a single plane, or the distribution pattern of the point cloud data differs too much from that of the plane. Therefore, the following steps are performed:
[0160] Record failure information: Record (e.g., to a log file) detailed information about the current point cloud data, including the number of points, centroid coordinates, covariance matrix, and eigenvalues; display a fitting failure message so that users or developers can understand the current status;
[0161] Dimensionality reduction: Perform dimensionality reduction on the point cloud data, such as using PCA (Principal Component Analysis) or other dimensionality reduction algorithms, to reduce the impact of noise; recalculate the covariance matrix and eigenvalues of the dimensionality-reduced point cloud data, and try plane fitting again;
[0162] Point cloud segmentation: Divide the current point cloud data into several subsets, and perform planar fitting on each subset; if a subset is successfully fitted, the fitting result of that subset is retained; if all subsets fail to fit, continue with subsequent processing.
[0163] Apply alternative fitting algorithms: Use the set alternative plane fitting algorithms, such as RANSAC (Random Sample Consensus) algorithm, to increase the probability of successful fitting; record the fitting results of the alternative plane fitting algorithms and compare them with the current results;
[0164] Mark as an abnormal region: If the number of fitting failures exceeds the set number, the point cloud data of the corresponding region will be marked as an abnormal region; the corresponding abnormal region will be marked in the voxel map for subsequent processing or manual intervention.
[0165] Continue processing other regions: Skip the point cloud data that failed to fit at the moment and continue to perform plane fitting on the point cloud data of other regions; ensure that the algorithm process is not affected by the failure of fitting.
[0166] Generate an error report: The error report includes the point cloud data that failed to fit, the fitting methods attempted, and an analysis of the reasons for the failure; send the error report to the system administrator, etc., for further analysis and problem solving;
[0167] Returning a failure status: Returning a status indicating that the fit has failed at the current stage of the algorithm, so that the calling function or module can take appropriate action based on the failure status; if possible, it can consider trying to fit again in a later stage or adopting other alternative strategies;
[0168] The fourth threshold is used to determine whether the feature value is small enough, that is, whether the point cloud data is sufficiently dispersed, so that these points can be considered to be roughly located on a plane.
[0169] The complexity of the plane fitting process mainly depends on the number of points N in the point cloud; the complexity of calculating the centroid is O(N), and the complexity of calculating the covariance matrix is O(N). 2 The complexity of feature decomposition depends on the specific algorithm implementation, and is typically O(N). 3 Therefore, the time complexity of the entire plane fitting algorithm is likely O(N). 3 );
[0170] D. Merging the planes to obtain the final facet pixel map; specifically including the following steps:
[0171] The disjoint-set data structure algorithm is used to merge the planes:
[0172] The Union-Find algorithm is a highly efficient data structure that supports two main operations: Find and Union. In dealing with planar geometry problems, it is used to determine whether two regions are connected and to merge them. The core of the Union-Find is an undirected graph consisting of nodes and edges connecting them. Each node represents a set, and if there is an edge between two nodes, then the two sets can be merged.
[0173] Among them, Find: finds the representative node of the set to which a given node belongs. This operation can be used to determine whether two nodes belong to the same set.
[0174] Union: Merges the sets containing two nodes into a single set; to improve efficiency, union-set data structures typically employ two optimization strategies:
[0175] Path compression: In a search operation, the direct parent node of a node is set as the root node, thereby reducing the depth of subsequent search operations;
[0176] Merge by rank: In the merge operation, the tree with the smaller rank (height of the tree) is connected to the tree with the larger rank to prevent the tree from becoming too tall;
[0177] The steps for plane merging include:
[0178] initialization:
[0179] Create a disjoint-set data structure. The number of elements in the disjoint-set data structure is equal to the number of points on the plane.
[0180] Initialize each point as a separate set;
[0181] Perform edge scanning:
[0182] Traverse all edges in the plane;
[0183] For each edge, find the representative node of the set containing the two connected points, and determine: if the representative nodes are different, merge the two sets; if the representative nodes are the same, proceed with the following steps:
[0184] Ignore duplicate edges: Record the current edge as processed to avoid duplicate calculations; skip the current edge and continue processing the next edge;
[0185] Update edge information: The edge information includes the edge length and the number of connected vertices; ensure that the information is not updated repeatedly due to repeated processing of the same edge;
[0186] Maintain the internal structure of the set: If the internal structure of the set needs to be maintained (for example, each set may have parent-child relationships or other structures), then check and maintain the internal structure of the set; ensure that the internal structure of the set after the merge operation meets the set requirements;
[0187] Mark processing status: Record the current edge as an internal edge and ensure that no more merging operations will be performed; if visualization or user interaction is possible, update the visualization interface to reflect the current merging status;
[0188] Optimize data structure: If the number of edges representing the same node exceeds a set value, optimize the data structure to reduce the number of search and merge operations;
[0189] Maintaining algorithm efficiency: Optimize the disjoint-set data structure, for example, by using path compression and rank merging to improve efficiency;
[0190] To determine whether two points are in the same set, follow these steps:
[0191] Obtain the representative nodes for each of the two points, and then determine whether the representative nodes are the same: if the representative nodes are the same, it means that the two points are in the same set; if the representative nodes are different, it means that the two points are not in the same set.
[0192] Step S4 involves segmenting the water surface on the face pixel map obtained in step S3 to obtain the water surface outline and water surface height data. Specifically, this includes the following steps:
[0193] First, due to the low reflectivity of water, lidar cannot acquire point cloud data of the water surface, thus the water surface area appears as a void in the point cloud map. Observing the point cloud map reveals that the water surface area is typically a completely or semi-closed void. Extracting the boundaries of this void allows for the estimation of the water surface's extent and height. Therefore, the following steps are used for water surface segmentation:
[0194] On the voxel map, all voxels are coarsely classified:
[0195] For each voxel, based on the normal of the internal fitting plane, determine the adjacent voxels reached by the plane extension from all adjacent voxels, and determine the occupancy status of the adjacent voxels: if all adjacent voxels are occupied, the voxel is located inside the point cloud map corresponding to the surface voxel map; otherwise, the voxel is located at the boundary of the point cloud map corresponding to the surface voxel map.
[0196] For each point in the voxel, search for K nearest neighbors, project the K nearest neighbors along the normal of the fitted plane, and make a judgment: if the normal projections are all biased to one side, the point is determined to be a boundary point; if the normal projections are evenly distributed, the point is determined to be an interior point.
[0197] Finally, the water surface boundary point set is obtained, the water surface is segmented, and the corresponding water surface contour and water surface height data are obtained;
[0198] Step S5, which uses the straight line data of the dam crest determined in step S2 to determine the dam crest plane and water surface outline, and calculates the dry beach length of the target tailings dam, specifically includes the following steps:
[0199] Based on the dam crest straight line data determined in step S2, determine the dam crest plane;
[0200] Project the set of water surface boundary points obtained in step S4 onto the dam crest plane and calculate the projection distance;
[0201] The minimum projected distance is used as the minimum dry beach length of the target tailings dam to complete the calculation of the dry beach length of the tailings dam based on point cloud data.
Claims
1. A method for calculating the dry beach length of a tailings dam based on point cloud data, comprising the following steps: S1. Obtain the 3D point cloud data of the target tailings dam; specifically including the following steps: A multi-rotor aerial surveying drone, equipped with lidar, was used to collect data on the target tailings dam and obtain three-dimensional point cloud data information of the target tailings dam. S2. Process the 3D point cloud data obtained in step S1 to generate a 3D point cloud map and mark the straight lines on the dam crest; specifically including the following steps: Using image processing software, the 3D point cloud data acquired at each moment in step S1 is registered according to the position and attitude of the UAV to generate a 3D point cloud map. The dam crest plane is determined based on the straight line at the dam crest, the direction of gravity, and the direction perpendicular to the straight line at the dam crest; the straight line at the dam crest is formed by points on the straight line at the dam crest. and direction vector express; S3. Generate a voxel map based on the 3D point cloud map obtained in step S2; S4. On the face pixel map obtained in step S3, perform water surface segmentation to obtain water surface contour and water surface height data; specifically including the following steps: On the voxel map, all voxels are coarsely classified: For each voxel, based on the normal of the internal fitting plane, determine the adjacent voxels reached by the plane extension from all adjacent voxels, and determine the occupancy status of the adjacent voxels: if all adjacent voxels are occupied, the voxel is located inside the point cloud map corresponding to the surface voxel map; otherwise, the voxel is located at the boundary of the point cloud map corresponding to the surface voxel map. For each point in the voxel, search for K nearest neighbors, project the K nearest neighbors along the normal of the fitted plane, and make a judgment: if the normal projections are all biased to one side, the point is determined to be a boundary point; if the normal projections are evenly distributed, the point is determined to be an interior point. Finally, the water surface boundary point set is obtained, the water surface is segmented, and the corresponding water surface contour and water surface height data are obtained; S5. Using the dam crest straight line data determined in step S2, determine the dam crest plane and water surface outline, and calculate the dry beach length of the target tailings dam; specifically including the following steps: Based on the dam crest straight line data determined in step S2, determine the dam crest plane; Project the set of water surface boundary points obtained in step S4 onto the dam crest plane and calculate the projection distance; The minimum projected distance is used as the minimum dry beach length of the target tailings dam to complete the calculation of the dry beach length of the tailings dam based on point cloud data.
2. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 1, characterized in that... Step S3, which involves generating a voxel map based on the 3D point cloud map obtained in step S2, includes the following steps: A. Traverse the point cloud list; B. Handling surface features and non-surface features; C. Initialize all octrees; D. Merging the planes yields the final facet matrix map.
3. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 2, characterized in that... Step A, which involves traversing the point cloud list, specifically includes the following steps: The 3D point cloud map obtained in step S2 is a list of distorted point clouds in the world coordinate system; Iterate through each point in the list of distorted point clouds : Calculate voxel positions: Place the points Divide each coordinate by the largest voxel size and round the result to the nearest integer to obtain the point. Position in the voxel gradient ; Step B, which involves processing surface features and non-surface features, specifically includes the following steps: For the points obtained in step A Determine the position of the facet in the voxel gradient: If point If it is a face feature, then the following steps are performed: Check if the key-value pair set includes positions. If location is already included Then the point The cache is located in the root node of the corresponding octree; if the position is not included... Then create a new octree structure and put the points Load it into the root node of the new octree; If point If it is a non-surface feature, then the following steps are performed: If an anomaly is detected in the input list of distorted point clouds, the anomaly information is recorded, the anomaly points are marked as invalid, and: Record the coordinates, timestamps, and corresponding attribute information of the anomalies; If the algorithm can continue running, skip the abnormal point and continue processing the next point in the point cloud list; If the number of outliers exceeds the set threshold, or the algorithm cannot continue to run, the following actions will be taken: send a prompt signal, output an error message, and terminate the current algorithm.
4. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 3, characterized in that... Step C, which involves initializing all octrees, specifically includes the following steps: After processing all points with both surface and non-surface features, octree initialization is performed on all octree structures in the key-value pair set to ensure that the internal nodes and leaf nodes of the octree are processed and that the octree correctly represents the spatial structure in the voxel map. After the octree is initialized, the initialized face pixel map is obtained.
5. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 4, characterized in that... The octree initialization specifically includes the following steps: Check the number of points in the input point cloud: Check if the number of points in the input point cloud exceeds a first preset threshold: If the number of points exceeds the first preset threshold, perform plane fitting to find the optimal plane containing the point cloud; if the number of points does not exceed the first preset threshold, proceed with the following steps: Directly assign to leaf nodes: The point cloud is directly assigned to the corresponding octagonal leaf nodes, and the corresponding leaf nodes are marked as occupied to indicate that the corresponding leaf nodes contain point cloud data; Record sparse regions: Record the currently processed region as a sparse region in the log file, and record the region information, including the coordinate range and number of points of the region; Merging neighboring nodes: Check if the neighboring nodes of the current leaf node are also sparse nodes: If the neighboring nodes are also sparse, merge the neighboring nodes into a larger node to reduce the total number of nodes in the branch tree and optimize the spatial structure; if the neighboring nodes are not sparse nodes, keep the current leaf node unchanged and record the state information of the corresponding neighboring nodes, including the node's data density, access frequency, and space occupation; then, check the other neighboring nodes of the current leaf node one by one to evaluate the sparsity and merging potential of the other neighboring nodes: if the other neighboring nodes are also sparse, merge the current leaf node with the other sparse neighboring nodes; if the other neighboring nodes are not sparse, mark the other neighboring nodes as "not to be merged for now" and monitor the corresponding data changes; at the same time, adjust the sparsity judgment criteria globally. Apply default parameters: Apply default voxel parameters to the sparse point cloud region, including voxel size and density; ensure that the voxel map of the sparse point cloud region is consistent with the voxel map of other regions. Anomaly Handling: If the number of points is less than the set threshold, it indicates that the point cloud cannot form a voxel structure that meets the set requirements. The corresponding point cloud is then marked as an anomaly point; detailed information about the anomaly point is recorded. Continue the initialization process: After processing the current node, continue to initialize other octree nodes; ensure that all nodes have been processed to guarantee the integrity and accuracy of the octree structure; Perform plane fitting: If the plane fitting is successful and the number of points does not exceed the second set threshold, then stop further segmentation; if the plane fitting is successful and the number of points exceeds the second set threshold, then the corresponding voxel will no longer be updated; if the plane fitting fails, then perform octree segmentation until a plane is successfully generated or the set number of layers is reached. Finally, the root node of the octree is output, which includes the data information of all point clouds.
6. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 5, characterized in that... The octree partitioning described herein specifically includes the following steps: The input data for octree partitioning is: the root node of the octree to be partitioned; The output data of octree partitioning is: an array of eight child nodes generated after partitioning, each child node representing a voxel; The segmentation process is as follows: Obtain point cloud data: Extract point cloud data from the input root node; Check if the bottom layer has been reached: Determine if the root node has reached the maximum depth of the octree: If the maximum depth has been reached, no more nodes will be created; If the maximum depth has not been reached, proceed to the next step; Initialize child node array: Initialize a child node array of size 8, where each element in the child node array corresponds to a child node of the octree; Assign point cloud data to child nodes: Traverse the point cloud data items in the extracted point cloud data, determine the position of the point cloud in the octree based on the voxel center and data item, and assign the point cloud data to the corresponding child node array; Plane fitting for child nodes: For each child node in the child node array, determine whether the number of points in each child node is greater than the third set threshold. If it is greater than the third set threshold, then plane fitting is performed. If plane fitting is successful, then the current child node is no longer split into octrees, and the octree splitting of the next child node begins. If plane fitting fails, continue to perform octree partitioning on the current child node; Returns an array of child nodes: After processing all child nodes, the final array of child nodes is obtained.
7. The method for calculating the length of the dry beach of a tailings dam based on point cloud data according to claim 6, characterized in that... The aforementioned plane fitting specifically includes the following steps: Input for plane fitting: point cloud data; the point cloud data includes a set of three-dimensional coordinate points, representing the point cloud data for which plane fitting is required; The output of plane fitting is the parameters of the fitted plane; these parameters include the normal vector and centroid of the fitted plane. The fitting process is as follows: Calculate the centroid: Calculate the centroid of the input point cloud data: the coordinates of the centroid in each dimension, which is the average of the coordinates of all point cloud data in the corresponding dimension; Calculate the covariance matrix: Calculate the covariance matrix of the input point cloud data; Eigenvalue decomposition: Perform eigenvalue decomposition on the obtained covariance matrix to obtain the corresponding eigenvalues and eigenvectors corresponding to the eigenvalues; Selecting the normal vector: Among the eigenvalues, select the eigenvector corresponding to the smallest eigenvalue as the normal vector of the fitting plane; Successful Fitting Determination: The minimum eigenvalue is evaluated: if the minimum eigenvalue is less than the fourth set threshold, the plane fitting is considered successful, and the corresponding normal vector and centroid are returned; if the minimum eigenvalue is greater than or equal to the fourth set threshold, the plane fitting is considered unsuccessful, and the following steps are performed: Record failure information: Record detailed information about the current point cloud data, including the number of points, centroid coordinates, covariance matrix, and eigenvalues; display a message indicating fitting failure; Dimensionality reduction: The point cloud data is dimensionality reduced to reduce the impact of noise; the covariance matrix and eigenvalues of the dimensionality-reduced point cloud data are recalculated, and plane fitting is performed again; Point cloud segmentation: Divide the current point cloud data into several subsets, and perform planar fitting on each subset; if a subset is successfully fitted, the fitting result of that subset is retained; if all subsets fail to fit, continue with subsequent processing. Apply alternative fitting algorithms: Use the set alternative plane fitting algorithms to increase the probability of successful fitting; Record the fitting results of the alternative plane fitting algorithms and compare them with the current result; Mark as an abnormal region: If the number of fitting failures exceeds the set number, the point cloud data of the corresponding region will be marked as an abnormal region; the corresponding abnormal region will be marked in the voxel map; Continue processing other regions: Skip the point cloud data that failed to fit at the moment and continue to perform plane fitting on the point cloud data of other regions; ensure that the algorithm process is not affected by the failure of fitting. Generate an error report: The error report includes the point cloud data that failed to fit, the attempted fitting methods, and an analysis of the reasons for the failure; send the error report; Return failure status: Returns a status indicating that the fit has failed at the current stage of the algorithm.
8. The method for calculating the length of dry tailings dam beach based on point cloud data according to claim 7, characterized in that... The merging of planes described in step D specifically includes the following steps: The disjoint-set data structure algorithm is used to merge the planes: initialization: Create a disjoint-set data structure. The number of elements in the disjoint-set data structure is equal to the number of points on the plane. Initialize each point as a separate set; Perform edge scanning: Traverse all edges in the plane; For each edge, find the representative node of the set containing the two connected points, and determine: if the representative nodes are different, merge the two sets; if the representative nodes are the same, proceed with the following steps: Ignore duplicate edges: Record the current edge as processed to avoid duplicate calculations; skip the current edge and continue processing the next edge; Update edge information: The edge information includes the edge length and the number of connected vertices; ensure that the information is not updated repeatedly due to repeated processing of the same edge; Maintain the internal structure of the set: If the internal structure of the set needs to be maintained, check and maintain the internal structure of the set; ensure that the internal structure of the set after the merge operation meets the set requirements; Mark processing status: Record that the current edge is an internal edge and ensure that no more merging operations will be performed; Update merge status; Optimize data structure: If the number of edges representing the same node exceeds a set value, optimize the data structure to reduce the number of search and merge operations; To determine whether two points are in the same set, follow these steps: Obtain the representative nodes for each of the two points, and then determine whether the representative nodes are the same: if the representative nodes are the same, it means that the two points are in the same set; if the representative nodes are different, it means that the two points are not in the same set.
Citation Information
Patent Citations
Tailing pond dry beach important parameter extraction method and system based on three-dimensional laser point cloud
CN113218310A
Dry beach detection method and system and storage medium
CN114120119A