Rock RQD detection method and system based on multi-view point cloud fusion
Patent Information
- Application Number
- CN202611275596.5
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-08-21
- Publication Date
- 2026-09-25
AI Technical Summary
然而,这类方法难以从透视变换后的二维投影中精确恢复各断裂面片的空间位置关系,因而无法准确定位完整岩芯段的起止边界并测量其轴向真实长度,导致RQD计算结果不准确
1、本发明通过改进的RANSAC平面拟合与多因子置信度评估,能够随机采样对复杂不规则断口进行有效建模和准确筛选,通过对每个局部点云簇采用随机采样、候选内点筛选及优化拟合的迭代优化策略,并结合内点数量与均方根误差的双指标优化,能够准确拟合阶梯状断口中的多个断裂面片,另外,本发明引入空间分布因子、方向一致性因子、内点占比因子及误差因子计算综合置信度,能够有效剔除拟合质量差或几何不合理的伪平面,进而提高了不规则断口建模的鲁棒性。
Smart Images

Figure CN122820705A_ABST
Abstract
Description
Technical Field
[0001] This application relates to the field of rock mass quality evaluation technology, specifically to a rock RQD detection method and system based on multi-view point cloud fusion. Background Technology
[0002] Rock Quality Designation (RQD) is a core parameter for rock mass quality grading and engineering stability evaluation. It is defined as the percentage of the sum of the lengths of intact core segments with a length of not less than 10 cm in a single drilling run to the total footage reached in that run. In traditional measurements, operators remove the core from the core box and visually identify discontinuities such as fractures, breaks, and fracture zones to divide the core into several segments. Then, the length of each segment is measured along its lateral midline, and the lengths of intact segments ≥10 cm are accumulated. Finally, the sum is divided by the total footage reached in the drilling run to obtain the RQD value.
[0003] However, during core drilling and coring, the fracture surfaces at both ends of rock cores often exhibit irregular geometric shapes, commonly including oblique end faces and stepped fractures. In recent years, researchers have proposed various automatic RQD measurement algorithms, such as the core image segmentation method based on CascadeMask R-CNN. However, these methods struggle to accurately recover the spatial relationships of each fracture surface from the perspective-transformed two-dimensional projection, thus failing to accurately locate the start and end boundaries of a complete core segment and measure its true axial length, leading to inaccurate RQD calculation results. Summary of the Invention
[0004] This application provides a rock RQD detection method and system based on multi-view point cloud fusion, which can accurately locate the start and end boundaries of a complete rock core segment and measure its true axial length, thereby improving the reliability of subsequent RQD calculation results.
[0005] A first aspect of this application provides a rock RQD detection method based on multi-view point cloud fusion, the method comprising: Multiple laser scanning probes at different angles are used to collect multi-view point clouds of the rock core segment, and the multi-view point clouds are registered to obtain the fused point cloud of the rock core segment. Based on the fused point cloud of the core segment, extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction; Based on the first end point cloud set and the second end point cloud set, the fracture plane is fitted to obtain the first end plane set and the second end plane set, respectively. A confidence assessment is performed on the first end plane set and the second end plane set. Based on the confidence assessment results, the end planes in the end plane set are filtered to obtain the filtered first end plane set and the second end plane set. Based on the filtered first end plane set and second end plane set, calculate the equivalent length of the filtered first end plane set and second end plane set; Repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling run, and calculate the rock RQD of the current drilling run based on the equivalent length of each complete core segment.
[0006] In one possible implementation, the step of extracting a first end point cloud set and a second end point cloud set of the core segment along the main axis direction based on the merged point cloud of the core segment includes: Principal component analysis was performed on the fused point cloud of the core segment to obtain the principal axis direction of the core segment; Project each point in the core segment fusion point cloud onto the principal axis direction to obtain one-dimensional projected coordinates, and determine the minimum and maximum projected coordinates; The end cut-off width of the core segment is determined based on the minimum and maximum projection coordinates. Based on the end cut-off width, the core segment fused point cloud is screened once to obtain a first candidate point set and a second candidate point set; Calculate the normal vector of each point in the first candidate point set and the second candidate point set respectively; Based on the angle between the normal vector and the principal axis direction, the first candidate point set and the second candidate point set are further filtered to obtain the first end point cloud set and the second end point cloud set.
[0007] In one possible implementation, the step of fitting the fracture plane based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set, respectively, includes: The first end point cloud set is segmented using a region growing algorithm to obtain multiple local point cloud clusters; Spatial plane fitting is performed on each of the local point cloud clusters to obtain the candidate plane, the set of interior points of the corresponding plane, and the fitting error; Based on the set of interior points and the fitting error, the candidate planes are filtered to obtain the first set of end planes; Perform the same steps as described above on the second end point cloud set to obtain the second end plane set.
[0008] In one possible implementation, the step of performing spatial plane fitting on each of the local point cloud clusters to obtain candidate planes, the set of interior points of the corresponding planes, and the fitting error includes: For each local point cloud cluster, three points are randomly selected multiple times to determine an initial plane; Points with a distance less than a preset distance threshold are selected as candidate interior points based on the initial plane; The plane is refitted based on the candidate interior point set using the least squares method to obtain the optimized initial plane. Calculate the distance from all points in the current local point cloud cluster to the optimized initial plane, mark points whose distance is less than a preset distance threshold as interior points, and count the number of interior points; Calculate the root mean square error of all interior points relative to the optimized initial plane, and use it as the fitting error of the optimized plane; Repeat the above steps until the preset number of iterations is reached; Among the optimized planes obtained from all iterations, the plane with the most interior points is selected as the candidate plane. The set of interior points corresponding to the candidate plane and the fitting error are output as the fitting result of the local point cloud cluster.
[0009] In one possible implementation, the step of performing a confidence assessment on the first end plane set and the second end plane set, and filtering the end planes in the end plane set based on the confidence assessment result to obtain the filtered first end plane set and the second end plane set, includes: For candidate planes in the first end plane set and the second end plane set, calculate the inward point proportion factor based on the proportion of the number of inward points corresponding to the current plane to the total number of points in the local point cloud cluster. The error factor is calculated based on the ratio of the current plane fitting error to the preset maximum allowable error. Principal component analysis is performed on the set of interior points in the current plane to obtain the extension length of the interior points in the two principal directions in the plane. The product of the extension lengths is calculated as the coverage area of the interior points. The ratio of the coverage area to the total area of the end region of the local point cloud cluster is used as the spatial distribution factor. Obtain the angle between the normal vector of the current plane and the principal axis direction of the core segment, and calculate the orientation consistency factor; The confidence level of the current candidate plane is calculated based on the in-point proportion factor, error factor, spatial distribution factor, and orientation consistency factor. Calculate the confidence scores of all candidate planes in the first and second end plane sets, and filter the end planes in the end plane sets based on the confidence score calculation results to obtain the filtered first and second end plane sets.
[0010] In one possible implementation, the angle between the normal vector of the current plane and the principal axis direction of the core segment is obtained, and the orientation consistency factor is calculated, including: Calculate the angle between the normal vector of the current plane and the principal axis direction of the core segment; When the included angle value is within a preset first angle range, the direction consistency factor is determined to be a first preset value; When the included angle value is within a preset second angle range, the direction consistency factor is determined to be a second preset value, which is less than the first preset value; When the included angle value is within a preset third angle range, the direction consistency factor is determined to be a third preset value, which is less than the second preset value.
[0011] In one possible implementation, calculating the equivalent lengths of the filtered first and second end-plane sets based on the filtered first and second end-plane sets includes: Obtain the principal axis direction of the core segment; Calculate the center point of each plane in the first end plane set after filtering, and project each center point onto the main axis direction to obtain the first projection coordinate set; Calculate the center point of each plane in the filtered second end plane set, and project each center point onto the principal axis direction to obtain the second projection coordinate set; The equivalent length is calculated based on the first set of projected coordinates and the second set of projected coordinates.
[0012] In one possible implementation, calculating the equivalent length based on the first set of projected coordinates and the second set of projected coordinates includes: Select the minimum value from the first set of projected coordinates and the maximum value from the second set of projected coordinates; The difference between the maximum value and the minimum value is taken as the equivalent length.
[0013] In one possible implementation, registering the multi-view point cloud to obtain the core segment fused point cloud includes: Based on the fixed angles and relative positions between multiple laser scanning probes, the multi-view point clouds collected by each probe are transformed into the same coordinate system to complete coarse registration; The iterative nearest point algorithm is used to perform fine registration on the coarsely registered multi-view point clouds so that the overlapping areas of the point clouds from different views reach the preset registration accuracy. After removing overlapping and redundant points from the point cloud after fine registration, the fused point cloud of the core segment is obtained.
[0014] This example provides a method and system for RQD detection of rocks based on multi-view point cloud fusion. First, multiple laser scanning probes at different angles acquire multi-view point clouds of the core segment. Coarse registration and ICP fine registration are then performed to obtain a fused point cloud, eliminating blind spots from single views. Next, principal component analysis is performed on the fused point cloud to obtain the principal axis direction. One-dimensional projection is used to extract candidate regions at both ends, and secondary screening is performed using the angle between the normal vector and the principal axis to eliminate interference from side points, resulting in a clean set of end point clouds. Then, region growing is used to segment the end point clouds to obtain multiple local point cloud clusters. An improved RANSAC algorithm is used to fit multiple candidate planes, and confidence scores are calculated based on four factors to select high-quality end planes. Finally, the center points of the selected end planes are projected onto the principal axis to obtain the equivalent length. The effective lengths of all complete core segments are summed and divided by the total advance footage to obtain the accurate RQD value. This invention has the following beneficial effects: 1. This invention, through improved RANSAC plane fitting and multi-factor confidence assessment, can effectively model and accurately screen complex irregular fracture surfaces by random sampling. By adopting an iterative optimization strategy of random sampling, candidate interior point screening, and optimized fitting for each local point cloud cluster, and combining dual index optimization of interior point number and root mean square error, it can accurately fit multiple fracture surfaces in stepped fracture surfaces. In addition, this invention introduces spatial distribution factor, orientation consistency factor, interior point proportion factor, and error factor to calculate comprehensive confidence, which can effectively eliminate pseudo-planes with poor fitting quality or unreasonable geometry, thereby improving the robustness of irregular fracture surface modeling.
[0015] 2. This invention uses the end-plane projection range method to calculate the equivalent length, which can accurately measure the true axial length of the core segment in complex scenarios such as inclined end faces and stepped fractures. By projecting the center points of each plane in the selected first and second end-plane sets onto the core main axis, the minimum value of the first end projection coordinate and the maximum value of the second end projection coordinate are taken respectively, and the difference between the two is taken as the equivalent length of the complete core segment. This method can cover the outermost span of multi-plane fractures and is not affected by a single inclined fracture or local concavity and convexity, thereby improving the accuracy of subsequent RQD value calculation.
[0016] 3. This invention improves the data integrity of the end point cloud by multi-view point cloud fusion and normal secondary screening. This invention uses multiple laser scanning probes at mutually angled angles to collect data synchronously, and uses geometric coarse registration and ICP fine registration to eliminate occlusion blind spots, which can completely obtain the three-dimensional morphology of the inclined end face and stepped fracture. Furthermore, after truncating the end candidate point set, the angle between the normal vector of each point and the principal axis direction is calculated to eliminate side interference points, which improves the input quality of the point cloud data for subsequent plane fitting.
[0017] A second aspect of this application provides a rock RQD detection system based on multi-view point cloud fusion, the system comprising: The first module is used to acquire multi-view point clouds of the core segment using multiple laser scanning probes at different angles, and to register the multi-view point clouds to obtain a fused point cloud of the core segment. The second module is used to extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction based on the fused point cloud of the core segment. The third module is used to perform fracture plane fitting based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set, respectively. The fourth module is used to perform confidence assessment on the first end plane set and the second end plane set, and to filter the end planes in the end plane set according to the confidence assessment result to obtain the filtered first end plane set and the second end plane set. The fifth module is used to calculate the equivalent length of the filtered first end plane set and the filtered second end plane set based on the filtered first end plane set and the second end plane set. The sixth module is used to repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling cycle, and to calculate the rock RQD of the current drilling cycle based on the equivalent length of each complete core segment.
[0018] A third aspect of this application provides a terminal including a processor, an input device, an output device, and a memory, wherein the processor, input device, output device, and memory are interconnected, wherein the memory is used to store a computer program, the computer program including program instructions, and the processor is configured to invoke the program instructions to execute the step instructions as described in the multi-view point cloud fusion-based rock RQD detection method in the first aspect of this application.
[0019] A fourth aspect of this application provides a computer-readable storage medium storing a computer program for electronic data interchange, wherein the computer program causes a computer to perform some or all of the steps described in the multi-view point cloud fusion-based rock RQD detection method of the first aspect of this application.
[0020] A fifth aspect of this application provides a computer program product, comprising a non-transitory computer-readable storage medium storing a computer program operable to cause a computer to perform some or all of the steps described in the multi-view point cloud fusion-based rock RQD detection method of the first aspect of this application. The computer program product may be a software installation package. Attached Figure Description
[0021] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.
[0022] Figure 1 This application provides a schematic diagram of the overall process for a rock RQD detection method based on multi-view point cloud fusion. Figure 2 This application provides a schematic diagram of the structure of a multi-view point cloud acquisition device for a rock RQD detection method based on multi-view point cloud fusion. Figure 3 This application provides a schematic diagram of the overall structure of a rock RQD detection system based on multi-view point cloud fusion. Figure 4 This application provides a schematic diagram of the structure of a terminal. Figure label: Multi-view point cloud acquisition device-1, laser scanning probe-2, first module-3, second module-4, third module-5, fourth module-6, fifth module-7, sixth module-8. Detailed Implementation
[0023] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, and not all embodiments. Based on the embodiments of this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.
[0024] The terms "first," "second," etc., in the specification, claims, and accompanying drawings of this application are used to distinguish different objects, not to describe a specific order. Furthermore, the terms "comprising" and "having," and any variations thereof, are intended to cover non-exclusive inclusion. For example, a process, method, system, product, or apparatus that includes a series of steps or units is not limited to the listed steps or units, but may optionally include steps or units not listed, or may optionally include other steps or units inherent to these processes, methods, products, or apparatuses.
[0025] In this application, the reference to "embodiment" means that a specific feature, structure, or characteristic described in connection with an embodiment may be included in at least one embodiment of this application. The appearance of this phrase in various places throughout the specification does not necessarily refer to the same embodiment, nor is it a mutually exclusive, independent, or alternative embodiment. It will be explicitly and implicitly understood by those skilled in the art that the embodiments described in this application can be combined with other embodiments.
[0026] The rock RQD detection method based on multi-view point cloud fusion is applied to the rock RQD detection system based on multi-view point cloud fusion. Figure 1 A schematic diagram of the overall process of a rock RQD detection method based on multi-view point cloud fusion is shown. Figure 2 A schematic diagram of a multi-view point cloud acquisition device for rock RQD detection based on multi-view point cloud fusion is shown. Figure 1 and Figure 2 As shown, it includes: S1. Multiple laser scanning probes at different angles are used to collect multi-view point clouds of the rock core segment, and the multi-view point clouds are registered to obtain the fused point cloud of the rock core segment.
[0027] Step S1 includes the following sub-steps: S101. Based on the fixed angle and relative position relationship between multiple laser scanning probes, the multi-view point clouds collected by each probe are transformed into the same coordinate system to complete the coarse registration.
[0028] Among them, such as Figure 2 As shown, the multi-view point cloud acquisition device 1 includes multiple laser scanning probes 2. There are three laser scanning probes 2, with the middle probe installed perpendicular to the core axis, and the left and right probes respectively forming fixed angles of ±30° with the middle probe. The distances between the three probes and the core segment are kept equal. The rotation matrix and translation vector between the three probes are obtained in advance through calibration. Using these fixed geometric relationships, the point clouds acquired by the left and right probes are transformed into the coordinate system of the middle probe, respectively, to achieve the initial alignment of the multi-view point clouds.
[0029] S102. The iterative nearest point algorithm is used to perform fine registration on the multi-view point cloud after coarse registration, so that the overlapping area of the point cloud from different views reaches the preset registration accuracy.
[0030] In this process, the point cloud of the middle probe is used as the reference point cloud, and the point clouds of the left and right probes after coarse registration are used as the point clouds to be registered. The ICP algorithm is then used for fine registration. Specifically, the ICP algorithm repeatedly finds the closest point pairs between the two sets of point clouds, calculates the optimal rotation matrix and translation vector to minimize the average sum of squared distances between the closest point pairs, and stops iterating when the error change between two adjacent iterations is less than the registration accuracy threshold or the maximum number of iterations is reached, thus obtaining the finely registered point cloud.
[0031] S103. Remove overlapping and redundant points in the point cloud after fine registration to obtain the fused point cloud of the core segment.
[0032] In the finely registered multi-view point cloud, there were a large number of duplicated points in the overlapping areas. A voxel filtering method was used to divide the entire point cloud space into a three-dimensional grid. Within each grid, only the point closest to the grid center was retained, and all other duplicate points were removed. After redundancy removal, a single, uniform, and non-overlapping core segment fused point cloud was obtained.
[0033] S2. Based on the fused point cloud of the core segment, extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction.
[0034] Step S2 includes the following sub-steps: S201. Perform principal component analysis on the fused point cloud of the core segment to obtain the principal axis direction of the core segment.
[0035] In this example, the core segment fusion point cloud includes Each point, the core segment fusion point cloud is recorded as... ,in, , To determine the 3D coordinates of the i-th point in the merged point cloud of the core segment, first calculate the centroid of the point cloud: (1) In the formula, The coordinates of the centroid of the fused point cloud of the core segment are given. This represents the total number of points in the merged point cloud of the core segment. Let be the three-dimensional coordinates of the i-th point in the core segment fusion point cloud.
[0036] Furthermore, construct the covariance matrix: (2) In the formula, The covariance matrix of the fused point cloud of the core segment; This represents the total number of points in the merged point cloud of the core segment. Let i be the three-dimensional coordinates of the i-th point in the core segment fusion point cloud; The coordinates of the centroid of the fusion point cloud of the core segment; This is the matrix transpose symbol.
[0037] Furthermore, regarding the covariance matrix Eigenvalue decomposition yields three eigenvalues. ( The largest eigenvalue corresponds to the principal axis direction; , (for minor eigenvalues) and their corresponding eigenvectors ( To correspond to the eigenvalues respectively , , eigenvectors), maximum eigenvalue corresponding feature vector The direction is the principal axis direction of the core segment (pointing from one end to the other), denoted as the unit vector of the principal axis direction. ,in, This is the unit vector along the principal axis of the rock core. For feature vectors The length of the module.
[0038] S202. Project each point in the core segment fusion point cloud onto the main axis direction to obtain one-dimensional projection coordinates, and determine the minimum and maximum projection coordinates.
[0039] For each point Calculate each point Relative to the center of mass In the main axis direction Projected coordinates on: (3) In the formula, The point represents the one-dimensional projected coordinates of the point along the principal axis of the core. Let i be the three-dimensional coordinates of the i-th point in the core segment fusion point cloud; The coordinates of the centroid of the fusion point cloud of the core segment; This is the unit vector along the principal axis of the rock core.
[0040] The projected coordinates of all points form a set. Then find the minimum projected coordinates. and maximum projected coordinates ,in, It is the minimum value in the set of projected coordinates. This is the operator that takes the minimum value of all projected coordinates. For the first The one-dimensional projection coordinates of a point on the principal axis The maximum value in the set of projected coordinates. The operator is used to retrieve the maximum value of all projected coordinates, which correspond to the two ends of the core segment along the principal axis.
[0041] S203. Determine the end cut width of the core segment based on the minimum and maximum projection coordinates.
[0042] First, the total length of the core segment along the principal axis is calculated. According to the preset proportional coefficient In this example (The value range can be adjusted from 5% to 10% according to the actual situation), fixed length (The value can range from 0.5 cm to 2 cm). The end cut-off width is specified. Take the larger of the two values: (4) In the formula, The width of the cut-off section at the end of the rock core. This is a preset proportional coefficient. The total length of the core segment along the main axis. The preset fixed end cutting length.
[0043] The width is used to extract a sufficient number of point clouds from both ends to ensure that a stable end point cloud can be obtained even for short core segments.
[0044] S204. Based on the end cut-off width, the core segment fusion point cloud is screened once to obtain the first candidate point set and the second candidate point set.
[0045] Wherein, the projected coordinates satisfy The points are assigned to the first candidate point set. , will satisfy The points are assigned to the second candidate point set. . and The original point clouds correspond to the two ends of the core segment, including points on the actual fracture surface and possibly a small number of points on the side of the core.
[0046] S205. Calculate the normal vector of each point in the first candidate point set and the second candidate point set respectively.
[0047] Among them, for the first candidate point set (Second candidate point set) Similarly, each point in The normal vector is estimated by fitting a local plane using its neighborhood point set. Specifically, the normal vector is estimated by... Find points using the nearest neighbor search method neighbor set Then, the covariance matrix of the neighborhood point set is calculated: (5) In the formula, For point The covariance matrix of the neighborhood point set, For the nearest neighbor set The center of mass, For point The set of neighboring points, For the nearest neighbor set The number of points, For the nearest neighbor set The three-dimensional coordinates of a single neighboring point.
[0048] Furthermore, regarding Eigenvalue decomposition is performed, and the eigenvector corresponding to the smallest eigenvalue is the point. normal vector And by unifying the sign of the dot product of the vector pointing from the centroid to that point, its direction is adjusted to ensure that it points to the outside of the core.
[0049] S206. Based on the angle between the normal vector and the principal axis direction, perform a second screening on the first candidate point set and the second candidate point set to obtain the first end point cloud set and the second end point cloud set.
[0050] Among them, for the first candidate point set Each point in Calculate each point normal vector With the main axis direction The included angle: (6) In the formula, Let be the angle between the normal vector of the point and the principal axis of the core. For point The normal vector, It is the direction vector of the main axis.
[0051] Furthermore, the normal vectors of points on the fracture surface of the core fracture should be approximately perpendicular to the principal axis, while the normal vectors of points on the side of the core should be approximately parallel to the principal axis. Therefore, an angle threshold is set. , make the included angle Points that are considered side points are discarded, and points that are retained are discarded. point.
[0052] Furthermore, after a second round of screening, the first candidate point set... The remaining points constitute the first end point cloud set. Similarly, the second end point cloud set is obtained. .
[0053] S3. Perform fracture plane fitting based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set respectively.
[0054] Step S3 includes the following sub-steps: S301. The first end point cloud set is segmented using a region growing algorithm to obtain multiple local point cloud clusters.
[0055] The input is the first end point cloud set obtained in step S2. Because core fracture surfaces often exhibit stepped or irregular shapes, meaning the fracture surface may consist of multiple small, spatially separated planar fragments with different orientations, it is necessary to... The point cloud in the image is divided into several independent local point cloud clusters, and each cluster corresponds to a potential fracture surface.
[0056] Specifically, the implementation of the region growing algorithm includes: for Each point in Using the same method as step S205 The nearest neighbor method calculates the normal vector. The curvature is calculated by the ratio of the minimum eigenvalue of the neighborhood point set covariance matrix to the sum of all eigenvalues. The smaller the curvature value, the flatter the local area where the point is located, and the more suitable it is as a seed point.
[0057] Furthermore, All the point curvature Sort the points from smallest to largest, select the point with the smallest curvature as the first seed point, and mark it as "ungrown".
[0058] Furthermore, initialize an empty point cloud cluster. Add the current seed point And mark it as "grown".
[0059] Furthermore, search the neighborhood points of the current seed point, for each neighborhood point... If the following two conditions are met, then add it. And mark it as grown: For the angle between the normal vectors, the following condition is satisfied: (7) In the formula, To preset the threshold for the angle between normal vectors, ensure that the normal vectors of all points within the grown point cloud cluster have the same direction. For neighborhood points The normal vector, This is the normal vector of the seed point.
[0060] For the distance condition, the following is satisfied: (8) In the formula, The three-dimensional coordinates of the neighboring points The three-dimensional coordinates of the seed point. This is a preset distance threshold.
[0061] Furthermore, the newly added points are used as new seed points, and the above search and growth process is repeated until no new points meet the conditions, resulting in a complete point cloud cluster. .
[0062] Furthermore, from Removed from the list For the point, select the unclimbed point with the smallest curvature from the remaining points as the new seed point, and repeat the region growth step until... There are no remaining points.
[0063] Furthermore, set the minimum number of points in the point cloud cluster. The number of deleted points is less than Clusters are considered noise.
[0064] Finally, the first end point cloud set It is divided into multiple local point cloud clusters. Each cluster represents an independent fracture surface, and these clusters are used as the basic unit for subsequent planar fitting.
[0065] In this example, a region growing algorithm is used to segment the end point cloud, which can effectively handle multiple fracture surfaces with different orientations in a stepped fracture. By using the dual constraints of the normal vector angle and spatial distance, it is ensured that the points inside each point cloud cluster belong to the same geometric plane, providing reasonable input units for subsequent multi-plane fitting. At the same time, selecting seed points based on curvature can preferentially start growing from flat areas, improving the stability and accuracy of segmentation.
[0066] S302. Perform spatial plane fitting on each of the local point cloud clusters to obtain the candidate plane, the set of interior points of the corresponding plane, and the fitting error.
[0067] Step S302 includes the following sub-steps: S3021. For each local point cloud cluster, randomly select three points multiple times to determine an initial plane.
[0068] Among them, the maximum number of iterations is set. In each iteration In the middle, from the current local point cloud cluster Three non-collinear points are randomly selected from the data. The initial plane parameters determined by these three points include: Two direction vectors: (9) (10) In the formula, , , The three-dimensional coordinates of three non-collinear points randomly selected from a local point cloud cluster are given. , These are two direction vectors within the initial plane.
[0069] Plane normal vector (unnormalized): (11) In the formula, This is the unnormalized normal vector of the initial plane.
[0070] Plane unit normal vector: (12) In the formula, Let be the unit normal vector of the initial plane. Unnormalized normal vector The length of the module.
[0071] The plane equation is expressed as: (13) In the formula, Let be the coordinates of any point in space.
[0072] In this example, randomly selecting three points is a classic sampling strategy of the RANSAC algorithm. Three points are sufficient to uniquely determine a spatial plane and are not sensitive to outliers. Through multiple random samplings, the hypothesis consisting of all valid interior points can be covered.
[0073] S3022. Select points whose distance is less than a preset distance threshold as candidate interior point set according to the initial plane.
[0074] Specifically, for the current initial plane, the local point cloud cluster is calculated. All points Distance to the plane: (14) in, The unit normal vector calculated in step S3021, This is one of the three points used to define the plane. Because... Since it is already a unit vector, with a denominator of 1, we get: (15) In the formula, For point Distance to the initial plane, Let be the unit normal vector of the initial plane. Let be the three-dimensional coordinates of the i-th point in the local point cloud cluster; This is the reference point on the initial plane.
[0075] Furthermore, a distance threshold is set. (In this embodiment, the accuracy can be adjusted according to the point cloud acquisition, and can be 1 to 2 times the average spacing of the point cloud). All those that meet the requirements... The points are included in the candidate interior point set. And record the number of candidate interior points. .
[0076] S3023. The plane is refitted based on the candidate interior point set using the least squares method to obtain the optimized initial plane.
[0077] Among them, for the candidate interior point set The points in the matrix are subjected to least-squares plane fitting to obtain more accurate plane parameters. The specific steps are as follows: Let the candidate interior point set have One point, denoted as .
[0078] Calculate the centroid : (16) In the formula, For candidate interior point set The coordinates of the centroid, For candidate interior point set The total number of interior points, For candidate interior point set The three-dimensional coordinates of the j-th interior point in the matrix.
[0079] Construct the covariance matrix: (17) In the formula, For candidate interior point set The covariance matrix, For candidate interior point set The total number of interior points, For candidate interior point set The three-dimensional coordinates of the j-th interior point.
[0080] right Eigenvalue decomposition yields three eigenvalues. and the corresponding feature vectors Minimum eigenvalue corresponding feature vector That is, the optimized plane normal vector .
[0081] The optimized plane equation is obtained: (18) In the formula, To optimize the unit normal vector of the plane, Let the three-dimensional coordinates of any point in space be... For candidate interior point set The coordinates of the centroid.
[0082] S3024. Calculate the distance from all points in the current local point cloud cluster to the optimized initial plane, mark points whose distance is less than a preset distance threshold as interior points, and count the number of interior points.
[0083] Among them, the optimized plane parameters (normal vector) are used. Center of mass Recalculate the local point cloud clusters. All points Distance to the plane: (19) In the formula, For point Distance to the optimized plane To optimize the unit normal vector of the plane, Let i be the three-dimensional coordinates of the i-th point in the local point cloud cluster. The reference point on the optimized plane (i.e., the centroid of the candidate interior point set). for The length of the module.
[0084] Among them, due to It may be a non-unit vector, so it needs to be divided by its magnitude. A distance threshold should be used. (Same as step S3022), will satisfy The points are marked as interior points, thus obtaining the final set of interior points. Count the number of interior points .
[0085] S3025. Calculate the root mean square error of all interior points relative to the optimized initial plane, and use it as the fitting error of the optimized plane.
[0086] Among them, the final set of interior points is calculated. The root mean square error (RMSE) of the distance from all points to the optimized plane: (20) In the formula, To fit the root mean square error, For the final set of interior points The total number of interior points, for The point in the middle, For point The directed distance to the optimized plane, the root mean square error is the fitting error of the optimized plane, and it quantitatively describes the degree of dispersion of the interior points to the plane.
[0087] S3026. Repeat the above steps to reach the preset number of iterations.
[0088] After completing one iteration (S3021-S3025), the optimized plane (normal vector) obtained in this iteration is recorded. Center of mass ) and the number of its interior points The fitting error RMSE is then calculated. The process then returns to S3021 to begin the next iteration, continuing until the preset maximum number of iterations is reached.
[0089] S3027. Among all the optimized planes obtained through iterations, select the plane with the most interior points as the candidate plane.
[0090] Among them, complete all After the iteration, select the number of interior points from all the optimized planes recorded. The largest plane; if multiple planes have the same (or similar) number of interior points, select the plane with the smallest fitting error RMSE, and denote the selected plane as . Its interior set is The number of interior points is The fitting error is .
[0091] S3028. Output the set of interior points corresponding to the candidate plane and the fitting error as the fitting result of the local point cloud cluster.
[0092] Among them, the candidate planes selected in step S3027 and its corresponding set of interior points Fitting error As a current local point cloud cluster The final fitting result. Output triplet. This is for use in subsequent steps.
[0093] S303. Based on the set of interior points and the fitting error, the candidate planes are filtered to obtain the first set of end planes.
[0094] In this process, the candidate planes for each local point cloud cluster obtained in step S302 are initially screened, and planes with obviously substandard quality are removed. Specific screening rules include: S3031, Number of Interior Points Planes with fewer than the preset minimum number of interior points are discarded.
[0095] S3032. Planes with a fitting error RMSE greater than the preset maximum permissible error threshold are discarded.
[0096] S304. Perform the same steps as above on the second end point cloud set to obtain the second end plane set.
[0097] Among them, the second end point cloud set obtained in step S2 Repeat all operations from S301 to S303 to obtain the second end plane set. .
[0098] S4. Calculate the confidence level of the first end plane set and the second end plane set, and filter the end planes in the end plane set according to the confidence level assessment results to obtain the filtered first end plane set and the second end plane set.
[0099] Step S4 includes the following sub-steps: S401. For candidate planes in the first end plane set and the second end plane set, calculate the inward point proportion factor based on the proportion of the number of inward points corresponding to the current plane to the total number of points in the local point cloud cluster.
[0100] Among them, the first end plane set obtained after filtering in step S303 Each candidate plane in The corresponding local point cloud cluster is known. (Segmentation results from step S301) Total number of points is its interior set The points (from step S3028) are Interior point percentage factor Defined as the ratio of the number of interior points to the total number of points in the local point cloud cluster: (twenty one) In the formula, The percentage factor for interior points; This represents the number of interior points corresponding to the plane. This represents the total number of points in the local point cloud cluster corresponding to this plane. The larger the value, the more points in the local point cloud cluster corresponding to the plane are identified as inliers; that is, the higher the proportion of the point cloud that the plane can explain, and the better the fit between the plane and the actual fracture surface. Conversely, if... If the value is very small, it means that the plane only covers a few points in the local point cloud, resulting in poor fitting quality.
[0101] S402. Calculate the error factor based on the ratio of the fitting error of the current plane to the preset maximum allowable error.
[0102] For the same candidate plane The set of interior points has been calculated in step S3025. Root mean square error relative to the plane Preset maximum allowable fitting error Error factor Defined as: (twenty two) In the formula, For error factor; This represents the root mean square error of the fit to the plane. The maximum allowable fitting error is preset; This is a function that takes the minimum value.
[0103] First, calculate and The ratio is taken, truncated to no more than 1, and then 1 is subtracted from the ratio. When hour, ;when hour, . The range of values is also The larger the value, the smaller the fitting error and the higher the plane quality.
[0104] S403. Perform principal component analysis on the set of interior points in the current plane to obtain the extension length of the interior points in the two principal directions in the plane. Calculate the product of the extension lengths as the coverage area of the interior points. Use the ratio of the coverage area to the total area of the end region of the local point cloud cluster as the spatial distribution factor.
[0105] Among them, for candidate planes its interior set These interior points are distributed in the vicinity of the plane in three-dimensional space. To assess whether these interior points uniformly and comprehensively cover the actual area of the fracture surface (rather than just clustering in a small area), a spatial distribution factor is introduced. The specific calculation process is as follows: S4031, due to The points in the diagram are not strictly located on the plane. Above (where there is a small distance), first project them onto a plane. Above, the set of projection points is obtained. For any interior point Its plane (by normal vector) and center of mass Projection point on (definition) The calculation formula is: (twenty three) In the formula, interior point In plane The three-dimensional coordinates of the projection point on the surface; Candidate plane The corresponding set of interior points The three-dimensional coordinates of any interior point in the array; For plane The centroid (i.e., the candidate interior point set) (center of mass) The plane obtained in step S3023 The unit normal vector.
[0106] S4032, Project the point set Consider the points as a set on a two-dimensional plane, calculate the covariance matrix of these points, and perform eigenvalue decomposition to obtain two eigenvalues. , which correspond to the variances of the point set in the principal and secondary directions of its extension in the plane, respectively.
[0107] S4033, Define the main direction extension length as... The extension length in the secondary direction is This definition uses twice the square root of the eigenvalues as the scattering range of the point set in the corresponding direction.
[0108] S4034, Coverage area of interior points It is approximately the product of two extension lengths: (twenty four) In the formula, Let be the area of the rectangular outer envelope occupied by the interior point in the plane; The length extending in the main direction; This refers to the length extended in the secondary direction; and These are two eigenvalues obtained from principal component analysis of the projection set of interior points. .
[0109] S4035. Calculate the total area of the end region of this local point cloud cluster. Defined as the projection of all original points (including interior and non-interior points) in the cluster onto the plane. The area of the convex hull on the surface.
[0110] Among them, this local point cloud cluster All points projected onto the plane The projection point set is obtained. .
[0111] Furthermore, the convex hull algorithm is used to calculate... The 2D convex hull is used to obtain the convex hull polygon.
[0112] Furthermore, calculate the area of the convex hull polygon, denoted as . .
[0113] S4036. Calculate the spatial distribution factor, defined as the ratio of the area covered by the inner points to the total area of the terminal regions, truncated to no more than 1: (25) In the formula, Spatial distribution factor; The area covered by the interior point; The total area of the end region, The range of values is The closer this factor is to 1, the more the interior points cover most of the fracture surface, meaning the fitted plane has good spatial representativeness. If the factor is very small, it indicates that the interior points are only clustered in a small corner of the patch, and the fitted plane may only reflect local distortion rather than the true fracture surface orientation.
[0114] S404. Obtain the angle between the normal vector of the current plane and the principal axis direction of the core segment, and calculate the direction consistency factor.
[0115] Step S404 includes the following sub-steps: S4041. Calculate the angle between the normal vector of the current plane and the principal axis direction of the core segment.
[0116] For each candidate plane obtained after the aforementioned screening Its unit normal vector The unit vector of the principal axis direction of the core segment has been determined in step S3023. The included angle has already been obtained in step S201. Calculated using the dot product formula of the normal vector and the principal axis direction vector: (26) In the formula, The angle between the normal vector of the candidate plane and the principal axis direction of the core segment; Let be the unit normal vector of the candidate plane; This is the unit vector along the principal axis of the rock core. and These are the magnitudes of the two vectors, respectively.
[0117] included angle The smaller the value, the closer the fracture surface is to the core axis; the larger the value, the more inclined or even parallel the fracture surface is to the axis.
[0118] S4042. When the included angle value is within a preset first angle range, the direction consistency factor is determined to be a first preset value.
[0119] Based on prior physical knowledge of core fracture, the first angle range is pre-set as follows: This range corresponds to the ideal situation where the fracture surface is basically perpendicular to the axis and the end face is relatively flat. When the calculated included angle... When falling into this interval, the direction consistency factor will be adjusted. The first preset value is determined to be the value that indicates that the geometric relationship between the direction of the current plane and the core axis is completely in line with expectations and contributes the most to the subsequent confidence calculation.
[0120] S4043. When the included angle value is within a preset second angle range, the direction consistency factor is determined to be a second preset value, where the second preset value is less than the first preset value.
[0121] The preset second angle range is: This range corresponds to a situation where the end face has a certain degree of inclination but still falls within the normal fracture morphology. When the included angle... When the plane falls within this range, the directional consistency factor is set to the second preset value. The second preset value is less than the first preset value, which means that the direction of the current plane deviates from the ideal state, but is still within the acceptable range and contributes moderately to the confidence level.
[0122] S4044. When the included angle value is within a preset third angle range, the direction consistency factor is determined to be a third preset value, which is less than the second preset value.
[0123] The preset third angle range is: This range corresponds to unreasonable situations where the end face is severely tilted or the fracture surface is almost parallel to the axis, usually due to lateral noise points or pseudo-planes caused by fitting non-principal fracture surfaces. When the included angle... When the plane falls within this interval, the orientation consistency factor is set to the third preset value. If the third preset value is less than the second preset value, it indicates that the orientation of the current plane is highly inconsistent with physical expectations and contributes little to the confidence level. Through the above piecewise mapping, the orientation consistency factor becomes a value that ranges from... The dimensionless coefficients between the candidate plane and the core axis directly reflect the degree of geometric fit between the candidate plane's normal and the core axis, thus effectively eliminating geometrically unreasonable pseudo-planes.
[0124] S405. Calculate the confidence level of the current candidate plane based on the in-point proportion factor, error factor, spatial distribution factor, and orientation consistency factor.
[0125] For each candidate plane The interior point proportion factor has been calculated in steps S401 to S404 respectively. Error factor Spatial distribution factor and direction consistency factor The range of values for all four factors is [missing information]. The fitting quality of the plane was characterized from different dimensions. The overall confidence level of the candidate plane was calculated using a weighted summation method. The calculation formula is as follows: (27) In the formula, The overall confidence level of the candidate plane; The percentage factor for interior points; For error factor; Spatial distribution factor; For directional consistency factor; The preset weight coefficients for each item satisfy the following conditions: The weights for each factor are not restricted in this example, but can be adjusted based on the specific application scenario and point cloud quality. The calculated values... The range of values is The larger the value, the higher the overall quality of the candidate plane.
[0126] S406. Calculate the confidence of all candidate planes in the first end plane set and the second end plane set, and filter the end planes in the end plane set according to the confidence calculation results to obtain the filtered first end plane set and the second end plane set.
[0127] Among them, for the first end plane set For each candidate plane, repeat steps S401 to S405 to calculate the overall confidence level of each plane. Similarly, for the second end plane set... Each candidate plane in the algorithm undergoes the same confidence calculation process to obtain its own overall confidence score.
[0128] After obtaining the confidence scores of all candidate planes, the confidence scores are then determined according to a preset confidence threshold. Filter the planes in the end-plane set: The overall confidence level is... The flat surface is retained and considered a high-quality end surface; The planes that are not found are discarded and considered unreliable planes. After screening, the first set of end planes is complete. The planes retained in the middle constitute the first end plane set after filtering. Second end plane set The planes retained in the middle constitute the second end plane set after filtering. .
[0129] If, after screening, a certain end plane set becomes an empty set (i.e., the confidence of all candidate planes is below the threshold), it indicates that the point cloud quality of that end is poor or there is no effective fracture surface. In this case, as a degradation processing strategy, we revert to using only the outermost edge of the convex hull of the end point cloud for length estimation, or directly mark the core segment as an invalid segment (not included in the RQD statistics).
[0130] In this embodiment, if or If the value is empty, the length calculation for that core segment will be abandoned, and the operator will be prompted to manually verify it.
[0131] S5. Based on the filtered first end plane set and second end plane set, calculate the equivalent length of the filtered first end plane set and second end plane set.
[0132] Step S5 includes the following sub-steps: S501. Obtain the main axis direction of the core segment.
[0133] Among them, the main axis direction of the core section The direction, already obtained through principal component analysis in step S201, describes the axial direction of the core segment from the first end to the second end. In this step, the principal axis direction vector is directly reused, eliminating the need for recalculation.
[0134] S502. Calculate the center point of each plane in the first end plane set after filtering, and project each center point onto the main axis direction to obtain the first projection coordinate set.
[0135] Among them, the first end plane set after filtering Includes There are candidate planes. Each plane in ( Its plane equation is derived from the normal vector. and center of mass Confirmed. The center point of this plane is the centroid. .
[0136] Each center point Projected onto the principal axis direction Above, obtain the projected coordinates. The calculation formula is: (28) In the formula, Let be the projected coordinates of the center point of the first end plane along the principal axis. For the first end The center point (centroid) of each candidate plane. The reference origin (specifically, the centroid of the fusion point cloud of the core segment) ), This is the unit vector along the principal axis of the rock core.
[0137] in, As a reference origin, this embodiment uses the cloud centroid of the core segment fusion point calculated in step S201. Using the origin as a reference, the projected coordinates of all the first end planes constitute the first projected coordinate set. .
[0138] S503. Calculate the center point of each plane in the filtered second end plane set, and project each center point onto the main axis direction to obtain the second projection coordinate set.
[0139] Among them, the filtered second end plane set Includes There are candidate planes. Each plane in ( Its center point is also taken as the centroid of the plane. Using the same reference origin as in step S502. Project each center point onto the principal axis direction Above, obtain the projected coordinates. : (29) In the formula, For the second end plane set The Middle The projected coordinates of the center point of each plane along the principal axis; For the second end plane set The Middle a plane The center point; Using the origin as a reference, The unit vector along the principal axis of the core is the projected coordinates of all the second end planes, which constitute the second projected coordinate set. .
[0140] S504. Calculate the equivalent length based on the first projection coordinate set and the second projection coordinate set.
[0141] Step S504 includes the following sub-steps: S5041. Select the minimum value from the first set of projected coordinates and the maximum value from the second set of projected coordinates.
[0142] Among them, due to the main axis direction of the core section It points from the first end to the second end, therefore from the first projected coordinate set. Select the minimum value This value represents the outermost fracture surface at the first end (i.e., the fracture surface closest to the first end). From the second projected coordinate set... Select the maximum value from the middle This value represents the location of the outermost fracture surface at the second end. When there are multiple fracture surfaces at the end (e.g., a stepped fracture), the innermost and outermost planes together determine the effective length of the core segment.
[0143] S5042. The difference between the maximum value and the minimum value is taken as the equivalent length.
[0144] Among them, equivalent length The calculation formula is: (30) In the formula, This is the equivalent length of the current complete core segment. From the second projected coordinate set The maximum value selected (representing the outermost fracture surface position at the second end). From the first set of projected coordinates The minimum value selected (representing the position of the fracture surface closest to the inside of the core at the first end) is positive, indicating the axial distance along the main axis from the innermost fracture surface at the first end to the outermost fracture surface at the second end. If the difference is negative due to the coordinate system orientation, the absolute value is taken. This is the effective length of the current complete core segment, used for length selection in subsequent RQD calculations.
[0145] S6. Repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling run, and calculate the rock RQD of the current drilling run based on the equivalent length of each complete core segment.
[0146] Specifically, for each complete core segment pre-divided in the current drilling cycle, steps S1 to S5 are repeated to obtain the equivalent length of each segment. Then filter out Calculate the total length of the segments. Obtain the total drilling footage per drilling cycle. (This information can be read from the drilling records), and the rock quality indicators are calculated using the formula: (31) In the formula, As a rock quality indicator, This refers to the total length of all complete core segments with an equivalent length of not less than 10cm in the current drilling cycle. This represents the total footage length of the current drilling cycle.
[0147] This example provides a method and system for RQD detection of rocks based on multi-view point cloud fusion. First, multiple laser scanning probes at different angles acquire multi-view point clouds of the core segment. Coarse registration and ICP fine registration are then performed to obtain a fused point cloud, eliminating blind spots from single views. Next, principal component analysis is performed on the fused point cloud to obtain the principal axis direction. One-dimensional projection is used to extract candidate regions at both ends, and secondary screening is performed using the angle between the normal vector and the principal axis to eliminate interference from side points, resulting in a clean set of end point clouds. Then, region growing is used to segment the end point clouds to obtain multiple local point cloud clusters. An improved RANSAC algorithm is used to fit multiple candidate planes, and confidence scores are calculated based on four factors to select high-quality end planes. Finally, the center points of the selected end planes are projected onto the principal axis to obtain the equivalent length. The effective lengths of all complete core segments are summed and divided by the total advance footage to obtain the accurate RQD value. This invention has the following beneficial effects: 1. This invention, through improved RANSAC plane fitting and multi-factor confidence assessment, can effectively model and accurately screen complex irregular fracture surfaces by random sampling. By adopting an iterative optimization strategy of random sampling, candidate interior point screening, and optimized fitting for each local point cloud cluster, and combining dual index optimization of interior point number and root mean square error, it can accurately fit multiple fracture surfaces in stepped fracture surfaces. In addition, this invention introduces spatial distribution factor, orientation consistency factor, interior point proportion factor, and error factor to calculate comprehensive confidence, which can effectively eliminate pseudo-planes with poor fitting quality or unreasonable geometry, thereby improving the robustness of irregular fracture surface modeling.
[0148] 2. This invention uses the end-plane projection range method to calculate the equivalent length, which can accurately measure the true axial length of the core segment in complex scenarios such as inclined end faces and stepped fractures. By projecting the center points of each plane in the selected first and second end-plane sets onto the core main axis, the minimum value of the first end projection coordinate and the maximum value of the second end projection coordinate are taken respectively, and the difference between the two is taken as the equivalent length of the complete core segment. This method can cover the outermost span of multi-plane fractures and is not affected by a single inclined fracture or local concavity and convexity, thereby improving the accuracy of subsequent RQD value calculation.
[0149] 3. This invention improves the data integrity of the end point cloud by multi-view point cloud fusion and normal secondary screening. This invention uses multiple laser scanning probes at mutually angled angles to collect data synchronously, and uses geometric coarse registration and ICP fine registration to eliminate occlusion blind spots, which can completely obtain the three-dimensional morphology of the inclined end face and stepped fracture. Furthermore, after truncating the end candidate point set, the angle between the normal vector of each point and the principal axis direction is calculated to eliminate side interference points, which improves the input quality of the point cloud data for subsequent plane fitting.
[0150] For those consistent with the above, please refer to Figure 3 , Figure 3 This application provides a schematic diagram of the structure of a rock RQD detection system based on multi-view point cloud fusion, as an embodiment of the present application. Figure 3 As shown, the system includes: The first module 3 is used to acquire multi-view point clouds of the core segment using multiple laser scanning probes at different angles, and to register the multi-view point clouds to obtain a fused point cloud of the core segment. The second module 4 is used to extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction based on the fused point cloud of the core segment. The third module 5 is used to perform fracture plane fitting based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set respectively. The fourth module 6 is used to perform confidence assessment on the first end plane set and the second end plane set, and to filter the end planes in the end plane set according to the confidence assessment result to obtain the filtered first end plane set and the second end plane set. Module 7 is used to calculate the equivalent length of the filtered first end plane set and the filtered second end plane set based on the filtered first end plane set and the filtered second end plane set. Module 8 is used to repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling cycle, and to calculate the rock RQD of the current drilling cycle based on the equivalent length of each complete core segment.
[0151] For examples consistent with the above embodiments, please refer to... Figure 4 , Figure 4 This is a schematic diagram of the structure of a terminal provided in an embodiment of this application, such as... Figure 4 As shown, it includes a processor, an input device, an output device, and a memory, which are interconnected. The memory is used to store a computer program, which includes program instructions. The processor is configured to call the program instructions. The program includes instructions for performing the following steps. Multiple laser scanning probes at different angles are used to collect multi-view point clouds of the rock core segment, and the multi-view point clouds are registered to obtain the fused point cloud of the rock core segment. Based on the fused point cloud of the core segment, extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction; Based on the first end point cloud set and the second end point cloud set, the fracture plane is fitted to obtain the first end plane set and the second end plane set, respectively. A confidence assessment is performed on the first end plane set and the second end plane set. Based on the confidence assessment results, the end planes in the end plane set are filtered to obtain the filtered first end plane set and the second end plane set. Based on the filtered first end plane set and second end plane set, calculate the equivalent length of the filtered first end plane set and second end plane set; Repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling run, and calculate the rock RQD of the current drilling run based on the equivalent length of each complete core segment.
[0152] The above mainly describes the solutions of the embodiments of this application from the perspective of the method execution process. It is understood that, in order to achieve the above functions, the terminal includes the corresponding hardware structure and / or software modules for executing each function. Those skilled in the art should readily recognize that, in conjunction with the units and algorithm steps of the various examples described in the embodiments provided herein, this application can be implemented in hardware or a combination of hardware and computer software. Whether a function is executed in hardware or by computer software driving hardware depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application.
[0153] This application embodiment can divide the terminal into functional units according to the above method example. For example, each function can be divided into a separate functional unit, or two or more functions can be integrated into one processing unit. The integrated unit can be implemented in hardware or as a software functional unit. It should be noted that the unit division in this application embodiment is illustrative and only represents one logical functional division. In actual implementation, there may be other division methods.
[0154] This application also provides a computer storage medium storing a computer program for electronic data interchange, which causes a computer to perform some or all of the steps of any of the rock RQD detection methods based on multi-view point cloud fusion as described in the above method embodiments.
[0155] This application also provides a computer program product, which includes a non-transitory computer-readable storage medium storing a computer program that causes a computer to perform some or all of the steps of any of the rock RQD detection methods based on multi-view point cloud fusion as described in the above method embodiments.
[0156] It should be noted that, for the sake of simplicity, the foregoing method embodiments are all described as a series of actions. However, those skilled in the art should understand that this application is not limited to the described order of actions, as some steps may be performed in other orders or simultaneously according to this application. Furthermore, those skilled in the art should also understand that the embodiments described in the specification are preferred embodiments, and the actions and modules involved are not necessarily essential to this application.
[0157] In the above embodiments, the descriptions of each embodiment have different focuses. For parts not described in detail in a certain embodiment, please refer to the relevant descriptions in other embodiments.
[0158] In the several embodiments provided in this application, it should be understood that the disclosed apparatus can be implemented in other ways. For example, the apparatus embodiments described above are merely illustrative; for instance, the division of units is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be through some interfaces; the indirect coupling or communication connection between devices or units may be electrical or other forms.
[0159] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.
[0160] Furthermore, the functional units in the various embodiments of this application can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software program module.
[0161] If the integrated unit is implemented as a software program module and sold or used as an independent product, it can be stored in a computer-readable storage device (CMD). Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, or all or part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a memory and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods described in the various embodiments of this application. The aforementioned memory includes various media capable of storing program code, such as USB flash drives, read-only memory (ROM), random access memory (RAM), portable hard drives, magnetic disks, or optical disks.
[0162] Those skilled in the art will understand that all or part of the steps in the various methods of the above embodiments can be implemented by a program instructing related hardware. The program can be stored in a computer-readable storage device, which may include: a flash drive, a read-only memory, a random access memory, a magnetic disk, or an optical disk, etc.
[0163] The embodiments of this application have been described in detail above. Specific examples have been used to illustrate the principles and implementation methods of this application. The description of the above embodiments is only for the purpose of helping to understand the method and core ideas of this application. At the same time, for those skilled in the art, there will be changes in the specific implementation methods and application scope based on the ideas of this application. Therefore, the content of this specification should not be construed as a limitation of this application.
Claims
1. A rock RQD detection method based on multi-view point cloud fusion, characterized in that, include: Multiple laser scanning probes at different angles are used to collect multi-view point clouds of the rock core segment, and the multi-view point clouds are registered to obtain the fused point cloud of the rock core segment. Based on the fused point cloud of the core segment, extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction; Based on the first end point cloud set and the second end point cloud set, the fracture plane is fitted to obtain the first end plane set and the second end plane set, respectively. A confidence assessment is performed on the first end plane set and the second end plane set. Based on the confidence assessment results, the end planes in the end plane set are filtered to obtain the filtered first end plane set and the second end plane set. Based on the filtered first end plane set and second end plane set, calculate the equivalent length of the filtered first end plane set and second end plane set; Repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling run, and calculate the rock RQD of the current drilling run based on the equivalent length of each complete core segment.
2. The rock RQD detection method based on multi-view point cloud fusion according to claim 1, characterized in that, The step of extracting the first end point cloud set and the second end point cloud set of the core segment along the main axis direction based on the merged point cloud of the core segment includes: Principal component analysis was performed on the fused point cloud of the core segment to obtain the principal axis direction of the core segment; Project each point in the core segment fusion point cloud onto the principal axis direction to obtain one-dimensional projected coordinates, and determine the minimum and maximum projected coordinates; The end cut-off width of the core segment is determined based on the minimum and maximum projection coordinates. Based on the end cut-off width, the core segment fused point cloud is screened once to obtain a first candidate point set and a second candidate point set; Calculate the normal vector of each point in the first candidate point set and the second candidate point set respectively; Based on the angle between the normal vector and the principal axis direction, the first candidate point set and the second candidate point set are further filtered to obtain the first end point cloud set and the second end point cloud set.
3. The rock RQD detection method based on multi-view point cloud fusion according to claim 1, characterized in that, The step of fitting the fracture plane based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set respectively includes: The first end point cloud set is segmented using a region growing algorithm to obtain multiple local point cloud clusters; Spatial plane fitting is performed on each of the local point cloud clusters to obtain the candidate plane, the set of interior points of the corresponding plane, and the fitting error; Based on the set of interior points and the fitting error, the candidate planes are filtered to obtain the first set of end planes; Perform the same steps as described above on the second end point cloud set to obtain the second end plane set.
4. The rock RQD detection method based on multi-view point cloud fusion according to claim 3, characterized in that, The step of performing spatial plane fitting on each of the local point cloud clusters to obtain candidate planes, the set of interior points of the corresponding planes, and fitting error includes: For each local point cloud cluster, three points are randomly selected multiple times to determine an initial plane; Points with a distance less than a preset distance threshold are selected as candidate interior points based on the initial plane; The plane is refitted based on the candidate interior point set using the least squares method to obtain the optimized initial plane. Calculate the distance from all points in the current local point cloud cluster to the optimized initial plane, mark points whose distance is less than a preset distance threshold as interior points, and count the number of interior points; Calculate the root mean square error of all interior points relative to the optimized initial plane, and use it as the fitting error of the optimized plane; Repeat the above steps until the preset number of iterations is reached; Among the optimized planes obtained from all iterations, the plane with the most interior points is selected as the candidate plane. The set of interior points corresponding to the candidate plane and the fitting error are output as the fitting result of the local point cloud cluster.
5. The rock RQD detection method based on multi-view point cloud fusion according to claim 1, characterized in that, The step of performing a confidence assessment on the first and second end plane sets, and filtering the end planes in the end plane sets based on the confidence assessment results to obtain the filtered first and second end plane sets includes: For candidate planes in the first end plane set and the second end plane set, calculate the inward point proportion factor based on the proportion of the number of inward points corresponding to the current plane to the total number of points in the local point cloud cluster. The error factor is calculated based on the ratio of the current plane fitting error to the preset maximum allowable error. Principal component analysis is performed on the set of interior points in the current plane to obtain the extension length of the interior points in the two principal directions in the plane. The product of the extension lengths is calculated as the coverage area of the interior points. The ratio of the coverage area to the total area of the end region of the local point cloud cluster is used as the spatial distribution factor. Obtain the angle between the normal vector of the current plane and the principal axis direction of the core segment, and calculate the orientation consistency factor; The confidence level of the current candidate plane is calculated based on the in-point proportion factor, error factor, spatial distribution factor, and orientation consistency factor. Calculate the confidence scores of all candidate planes in the first and second end plane sets, and filter the end planes in the end plane sets based on the confidence score calculation results to obtain the filtered first and second end plane sets.
6. The rock RQD detection method based on multi-view point cloud fusion according to claim 5, characterized in that, Obtain the angle between the normal vector of the current plane and the principal axis direction of the core segment, and calculate the orientation consistency factor, including: Calculate the angle between the normal vector of the current plane and the principal axis direction of the core segment; When the included angle value is within a preset first angle range, the direction consistency factor is determined to be a first preset value; When the included angle value is within a preset second angle range, the direction consistency factor is determined to be a second preset value, which is less than the first preset value; When the included angle value is within a preset third angle range, the direction consistency factor is determined to be a third preset value, which is less than the second preset value.
7. The rock RQD detection method based on multi-view point cloud fusion according to claim 1, characterized in that, The step of calculating the equivalent lengths of the filtered first and second end plane sets based on the filtered first and second end plane sets includes: Obtain the principal axis direction of the core segment; Calculate the center point of each plane in the first end plane set after filtering, and project each center point onto the main axis direction to obtain the first projection coordinate set; Calculate the center point of each plane in the filtered second end plane set, and project each center point onto the principal axis direction to obtain the second projection coordinate set; The equivalent length is calculated based on the first set of projected coordinates and the second set of projected coordinates.
8. The rock RQD detection method based on multi-view point cloud fusion according to claim 7, characterized in that, The step of calculating the equivalent length based on the first projection coordinate set and the second projection coordinate set includes: Select the minimum value from the first set of projected coordinates and the maximum value from the second set of projected coordinates; The difference between the maximum value and the minimum value is taken as the equivalent length.
9. The rock RQD detection method based on multi-view point cloud fusion according to claim 1, characterized in that, The process of registering the multi-view point cloud to obtain the core segment fused point cloud includes: Based on the fixed angles and relative positions between multiple laser scanning probes, the multi-view point clouds collected by each probe are transformed into the same coordinate system to complete coarse registration; The iterative nearest point algorithm is used to perform fine registration on the coarsely registered multi-view point clouds so that the overlapping areas of the point clouds from different views reach the preset registration accuracy. After removing overlapping and redundant points from the point cloud after fine registration, the fused point cloud of the core segment is obtained.
10. A rock RQD detection system based on multi-view point cloud fusion, characterized in that, include: The first module is used to acquire multi-view point clouds of the core segment using multiple laser scanning probes at different angles, and to register the multi-view point clouds to obtain a fused point cloud of the core segment. The second module is used to extract the first end point cloud set and the second end point cloud set of the core segment along the main axis direction based on the fused point cloud of the core segment. The third module is used to perform fracture plane fitting based on the first end point cloud set and the second end point cloud set to obtain the first end plane set and the second end plane set, respectively. The fourth module is used to perform confidence assessment on the first end plane set and the second end plane set, and to filter the end planes in the end plane set according to the confidence assessment result to obtain the filtered first end plane set and the second end plane set. The fifth module is used to calculate the equivalent length of the filtered first end plane set and the filtered second end plane set based on the filtered first end plane set and the second end plane set. The sixth module is used to repeat the above steps to obtain the equivalent length of all complete core segments in the same drilling cycle, and to calculate the rock RQD of the current drilling cycle based on the equivalent length of each complete core segment.