Rock mass trace automatic extraction method and system based on point cloud

Abstract

The invention provides a rock mass trace automatic extraction method and system based on point cloud, and belongs to the technical field of geotechnical engineering / data processing. The method comprises the following steps: acquiring regular body point cloud data and rock mass point cloud data, taking the regular body point cloud data as a source domain based on an induced transfer learning thought, training a universal trace feature point extraction model by using an artificial neural network (ANN) technology, and then transferring the trained trace feature point extraction model to a target domain, namely actual rock mass point cloud. And carrying out coarse extraction on the rock mass point cloud trace feature points, and then sequentially carrying out fine extraction, feature point connection and trace combination processing on the rock mass trace feature points to obtain a final rock mass trace. According to the method, for the problem that threshold selection of different rock mass data densities is difficult, differences of different rock mass point clouds are considered, a trace feature point extraction method combining a multi-scale surface change index and machine learning is provided, threshold setting in the point clouds is not needed, and trace extraction precision and efficiency are improved.

Description

Technical Field

[0001] The present application relates to the field of geotechnical engineering / data processing technology, and in particular to a method and system for automatically extracting rock mass traces based on point clouds. Background Art

[0002] In recent years, the country has made steady progress in engineering construction, continuously increased investment in water conservancy and hydropower, coal mining and transportation, and the number of various rock engineering projects has continued to increase. Engineering safety is the core issue of engineering projects, and the safety factor of rock engineering is largely affected by the mechanical properties of rock. Rock trace is one of the important factors affecting the mechanical properties of rock. The International Society of Rock Mechanics recommends using trace length as one of the parameters for describing the rock structure surface (ISRM, 1978). Therefore, rapid and accurate identification of trace information is of great significance for rock stability analysis. Traditional methods usually use contact-type manual field measurement methods to collect rock data, such as measurement methods using tape measures and geological compasses, scanning line methods or window sampling methods. These methods have multiple problems. First, they are inefficient and require a lot of time and human resources to complete the measurement task. Second, the measurement results are prone to large deviations, and the accuracy of the data cannot be guaranteed. In addition, during field measurements, surveyors need to work in a complex and dangerous environment, which greatly increases the risk of operation. Finally, the measurement range of these methods is limited, and they cannot cover a large area of ​​rock, and cannot meet the needs of modern rock engineering for large-scale data collection. Therefore, the application of traditional methods in modern rock engineering is greatly restricted, and there is an urgent need to find more efficient, safe and accurate alternatives.

[0003] With the development of LiDAR (Light Detection and Ranging) technology and photogrammetry technology, high-precision three-dimensional point clouds of rock masses have become the main data source for extracting rock mass geometric parameter information. At present, in the technology of extracting traces from rock mass point clouds, thresholds are usually set in the color space and / or geometric space of the point cloud to preliminarily extract feature points (coarse extraction), and then the feature points are finely extracted and connected to obtain traces. However, the point cloud data density of different rock masses is usually quite different. If a unified threshold is used to extract feature points, the accuracy is often insufficient and it is difficult to meet the actual application requirements. In this regard, some technologies use the method of setting thresholds for different rock mass point clouds, which to a certain extent solves the problem of insufficient accuracy caused by improper threshold selection. However, on the one hand, since there is no unified method for selecting thresholds, even for the same point cloud data, the threshold selection results obtained by different methods still have great uncertainty, and it is difficult to guarantee the effect of feature point extraction. On the other hand, finding suitable thresholds for different point cloud data will undoubtedly greatly increase the workload.

[0004] Therefore, it is necessary to provide an improved technical solution to address the above-mentioned deficiencies in the prior art. Summary of the invention

[0005] The purpose of this application is to provide a method and system for automatic extraction of rock mass traces based on point clouds, so as to solve or alleviate the problems existing in the above-mentioned prior art.

[0006] In order to achieve the above objectives, this application provides the following technical solutions:

[0007] In a first aspect, the present application provides a method for automatically extracting rock mass traces based on point clouds, comprising:

[0008] Obtaining regular body point cloud data and rock mass point cloud data, and extracting a sample point set from the regular body point cloud data, wherein each sample point in the sample point set contains label information;

[0009] Calculate the multi-scale surface change index of each sample point in the sample point set to obtain a first feature set of the sample point set; calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data; wherein the total number of features of the first feature set and the total number of features of the second feature set are both equal to the scale n, the i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale;

[0010] Taking the first feature set as input features, combined with the label information contained in each sample point in the sample point set, a training data set is constructed; and using the training data set, based on the idea of ​​induced transfer learning, using artificial neural network technology ANN to train a general trace feature point extraction model, and then inputting the second feature set into the trained trace feature point extraction model to perform rough extraction of rock mass point cloud trace feature points to obtain a rough extracted feature point set;

[0011] Based on the coarsely extracted feature point set, the rock mass trace feature point fine extraction, feature point connection and trace merging processing are performed in sequence to obtain the final rock mass trace.

[0012] In combination with the first aspect, in some possible implementations, the regular body point cloud data or the rock body point cloud data is used as the input point cloud P, and any sample point in the sample point set or any point in the rock body point cloud data is used as the current point p; then the calculation steps of the multi-scale surface change index of the current point p are as follows:

[0013] Among all the points in the input point cloud P, search for the k points closest to the current point p to form a neighbor point set N k (p);

[0014] Calculate the set of nearest neighbor points N k (p) Any point p i The Euclidean distance to the current point p, and determine the set of neighboring points N k The maximum distance d from all points in (p) to the current point p k (p);

[0015] With d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, wherein the number of radii contained in the multi-scale search radius set is equal to the scale n;

[0016] Taking the current point p as the center and taking each radius in the multi-scale search radius set as the search radius, search to obtain the neighboring point set of the current point p at each scale;

[0017] Based on the neighbor point sets of the current point p at each scale, the surface change index values ​​of the current point p at each scale are calculated.

[0018] In combination with the first aspect, in some possible implementations, based on the neighboring point sets of the current point p at each scale, the surface change index values ​​of the current point p at each scale are calculated, specifically:

[0019] For the neighbor point set of the current point p at each scale, construct a neighbor point matrix of the current point p at each scale, the number of rows of the neighbor point matrix of the current point p at each scale is equal to the number of points included in the neighbor point set of the current point p at the scale, and the columns of the neighbor point matrix include three columns, namely, x coordinate, y coordinate and z coordinate;

[0020] Based on the neighbor point matrix of the current point p at each scale, the covariance matrix C of the current point p at each scale is calculated;

[0021] Performing eigenvalue decomposition on the covariance matrix C of the current point p at each scale by principal component analysis, and arranging the generated eigenvalues ​​in descending order to determine the maximum eigenvalue and the top three eigenvalues ​​of the current point p at each scale;

[0022] The ratio of the maximum eigenvalue of the current point p at each scale to the sum of the eigenvalues ​​arranged in the first three positions is calculated to obtain the surface change index of the current point p at each scale.

[0023] In some possible implementations, the maximum acquisition range is determined by the following formula:

[0024]

[0025] In the formula, r max Indicates the maximum acquisition range, d bb Represents the diagonal length of the bounding box of the input point cloud P, and |N| represents the number of points in the input point cloud P.

[0026] In some possible implementations, d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, including:

[0027] With d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. A fixed step size is used to set n equidistant breakpoints within the upper and lower limits of the search radius to form the multi-scale search radius set;

[0028] or,

[0029] With d k (p) is the lower limit of the search radius, and n breakpoints within the upper and lower limits of the search radius are obtained by adaptive iteration to form the multi-scale search radius set;

[0030] The adaptive iteration includes:

[0031] Assume d k (p) is the starting radius r0 of the iteration, the current point p to the neighboring point set N k The mean of the Euclidean distances of all points in (p) is used as the step size Δr0 of the first iteration, and the radius value r1 generated by the first iteration is:

[0032] r1=r0+Δr0,

[0033] Calculate the adaptive step size of the k1th iteration and the breakpoint value of each iteration:

[0034] r k1+1 =r k1 +Δr k1 ,

[0035] Δr k1 =Δr0·tanh(I(R*r k1 ))),

[0036]

[0037] I(R(rk1 ))=-log P(R(r k1 )),

[0038] Where k1 is the current iteration number, r k1+1 、r k1 are the breakpoint values ​​at the k1+1th and k1th iterations, Δr k1 is the step size of the k1th iteration, R(r k1 ) is the normalization function, I(R(r k1 )) represents the normalized r k1 The amount of information, P(R(r k1 )) is R(r k1 )’s probability distribution.

[0039] In some possible implementations, the rock mass trace feature point precise extraction includes:

[0040] The feature points in the coarsely extracted feature point set are stripped from the rock mass point cloud data to obtain a trace feature point cloud set.

[0041] The trace feature point cloud set is randomly sampled using a random sampling technique to obtain a sampling point set;

[0042] Based on the trace feature point cloud set and the sampling point set, an optimization target equation is constructed using the L1 median skeleton extraction algorithm. The expression of the optimization target equation is as follows:

[0043]

[0044] In the formula, x i is the i-th point in the sampling point set, p j is the jth point in the trace feature point cloud set, I and J are the number of points in the sampling point set and the trace feature point cloud set, θ is the weight function, γ i is the balance coefficient, x i′ is the i′th point in the sampling point set, and i′≠i, σ i It is point x i The linearity of the location;

[0045] The optimization objective equation is iteratively solved to obtain precisely extracted rock mass trace feature points.

[0046] In some possible implementations, connecting the characteristic points of the rock mass trace includes the following steps:

[0047] Eliminating external points of the precisely extracted rock mass trace feature points;

[0048] A seed point p0 is randomly selected from the precisely extracted rock mass trace feature points as a center point, and a neighbor point set of the point p0 is calculated;

[0049] For each neighbor point in the neighbor point set of the point p0, calculate its weight;

[0050] According to the calculated weight result, select the point with the largest weight as a new seed point, and connect the new seed point and the center point;

[0051] Take the new seed point as the new center point, find the neighbor point set of the center point and calculate the weight of each neighbor point, repeat the above steps until the trace is connected and get the connected trace.

[0052] In some possible implementations, the rock mass trace merging process is specifically performed by using a distance-based and collinearity-driven trace merging method, including the following steps:

[0053] The connected traces are smoothed using a B-spline interpolation method to obtain a smoothed trace;

[0054] Collecting all endpoints of the smoothed trace to construct an endpoint set;

[0055] Take any point in the endpoint set as the current center point p i , perform the following steps:

[0056] Search for the neighboring points of the center point and obtain the neighboring point set of the center point

[0057] Define the center point p i The vector to the jth neighbor point in its neighbor point set is Center point p i The vector pointing to the adjacent point on its trajectory is Defining vectors and The angle between them is θ j , according to the center point p i The Euclidean distance and angle between the jth neighbor point and the jth neighbor point are θ j Calculate the adjacent trace endpoints to the current center point p i The Euclidean distance d of the trace tt ;

[0058] If the angle is θ j is smaller than the preset angle threshold, and the distance d tt If the distance is less than the preset threshold, the jth neighbor point is assigned to the candidate point set and the weight of the jth neighbor point is calculated;

[0059] The point with the largest weight in the candidate point set is used as the connection point and connected to the center point p i connect;

[0060] The above steps are repeated until all endpoints in the endpoint set are calculated.

[0061] In some possible implementations, the weight of the jth neighbor point is expressed as follows:

[0062]

[0063] In the formula, ω j is the weight of the jth neighbor point, λ0, λ1, λ2 are the weights of the center point p i The characteristic value of pp is the center point p i The Euclidean distance to its jth nearest neighbor, d tt It is the adjacent trace endpoint to the current center point p i Euclidean distance of the trace.

[0064] In a second aspect, this embodiment provides a point cloud-based rock mass trace automatic extraction system, comprising:

[0065] A data acquisition module is configured to acquire regular body point cloud data and rock body point cloud data, and extract a sample point set from the regular body point cloud data, wherein each sample point in the sample point set contains label information;

[0066] A calculation module is configured to calculate the multi-scale surface change index of each sample point in the sample point set to obtain a first feature set of the sample point set; calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data; wherein the total number of features of the first feature set and the total number of features of the second feature set are both equal to the scale n, the i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale;

[0067] The rough extraction module is configured to use the first feature set as input features, and to construct a training data set in combination with label information contained in each sample point in the sample point set; and to use the training data set, in combination with the idea of ​​induced transfer learning, to train a general trace feature point extraction model using artificial neural network technology ANN, and then to input the second feature set into the trained trace feature point extraction model, to perform rough extraction on the rock mass point cloud trace feature points, and to obtain a rough extraction feature point set;

[0068] The post-processing module is configured to perform, based on the coarsely extracted feature point set, precise extraction of rock mass trace feature points, feature point connection, and trace merging processing in sequence to obtain a final rock mass trace.

[0069] The technical solution of the embodiment of the present application has the following beneficial effects:

[0070] In this embodiment, regular body point cloud data and rock body point cloud data are obtained, and a sample point set is extracted from the regular body point cloud data; multi-scale surface change indicators of each sample point in the sample point set and all points in the rock body point cloud data are calculated respectively to obtain a first feature set and a second feature set of the sample point set; then the first feature set is used as an input feature, and a training data set is constructed in combination with the label information contained in each sample point in the sample point set; and the training data set is used to train a general trace feature point extraction model based on the induced transfer learning idea using artificial neural network technology ANN, and then the second feature set is input into the trained trace feature point extraction model to perform rough extraction of rock body point cloud trace feature points to obtain a rough extracted feature point set; based on the rough extracted feature point set, rock body trace feature point fine extraction, feature point connection, and trace merging processing are sequentially performed to obtain a final rock body trace. This method aims to solve the problem of difficult threshold selection caused by different data densities of different rock masses. It fully considers the density differences of point clouds of different rock masses and proposes a trace feature point extraction (rough extraction) method combining multi-scale surface change index (MsSVI) with machine learning. There is no need to set a threshold in the point cloud to extract feature points, which improves the accuracy and efficiency of trace extraction. BRIEF DESCRIPTION OF THE DRAWINGS

[0071] Figure 1 A schematic flow chart of a method for automatic extraction of rock traces based on point clouds according to some embodiments of the present application.

[0072] Figure 2 A schematic diagram of the structure of a point cloud-based rock mass trace automatic extraction system provided according to some embodiments of the present application. DETAILED DESCRIPTION

[0073] The terms "first", "second", "third" and "fourth" etc. in the specification and claims of the present application and the drawings are used to distinguish different objects, rather than to describe a specific order. In addition, the terms "including" and "having" and any variations thereof are intended to cover non-exclusive inclusions. For example, a process, method, system, product or device that includes a series of steps or units is not limited to the listed steps or units, but optionally includes steps or units that are not listed, or optionally includes other steps or units inherent to these processes, methods, products or devices.

[0074] Reference to "embodiments" herein means that a particular feature, structure, or characteristic described in conjunction with the embodiments may be included in at least one embodiment of the present application. The appearance of the phrase in various locations in the specification does not necessarily refer to the same embodiment, nor is it an independent or alternative embodiment that is mutually exclusive with other embodiments. It is explicitly and implicitly understood by those skilled in the art that the embodiments described herein may be combined with other embodiments.

[0075] "Multiple" means two or more. "And / or" describes the relationship between related objects, indicating that there can be three relationships. For example, A and / or B can mean: A exists alone, A and B exist at the same time, and B exists alone. The character " / " generally indicates that the related objects are in an "or" relationship.

[0076] The embodiments of the present application are described below in conjunction with the accompanying drawings.

[0077] Method Embodiment

[0078] LiDAR technology obtains high-precision three-dimensional terrain data by emitting laser beams and measuring the reflection time. It not only provides extremely high spatial resolution, but also can penetrate vegetation cover to obtain rock mass information under the surface. Photogrammetry technology generates high-precision digital models of rock mass through three-dimensional reconstruction of multi-view photos. These technologies are far superior to traditional contact measurement methods in terms of data acquisition accuracy. In addition, LiDAR and photogrammetry technologies also have significant advantages in safety. Since they are non-contact measurement technologies, operators do not need to conduct on-site measurements in dangerous rock mass environments, which greatly reduces personal risks. These technologies are also not limited by complex terrain and harsh climatic conditions, and can be measured in a wider range of environments, providing continuous and comprehensive rock mass data. In summary, the application of LiDAR and photogrammetry technology not only improves the efficiency and accuracy of rock mass data acquisition, but also significantly improves safety and applicability, becoming an indispensable tool in modern rock mass engineering. With the continuous advancement and popularization of technology, these advanced measurement methods will play an increasingly important role in the field of rock mass engineering.

[0079] Therefore, in view of the limitations of the traditional trace measurement method, this embodiment proposes a method for automatically extracting rock traces based on point clouds based on LiDAR technology. The execution subject of this method is a personal computer, server, tablet and other electronic devices, such as Figure 1 As shown, the method includes:

[0080] Step S101: Obtain regular body point cloud data and rock mass point cloud data, and extract a sample point set from the regular body point cloud data, wherein each sample point in the sample point set contains label information.

[0081] Point cloud data is a collective data structure used to represent objects or environments in three-dimensional space. Point clouds are composed of many discrete points, each with three-dimensional coordinates (x, y, z). Point clouds are composed of points on the surface or volume of these geometric bodies. In some specific application scenarios, point cloud data can also include other attributes, such as color (RGB value), intensity, normal vector, etc.

[0082] A regular body point cloud refers to a collection of points obtained by sampling a regular geometric body in three-dimensional space. A regular body usually refers to a geometric body with a regular shape, such as a cube, a sphere, a cylinder, etc. A rock mass point cloud is point cloud data of an actual rock mass obtained through laser scanning, photogrammetry, and other methods.

[0083] In this embodiment, extracting traces from point clouds may include the following steps: rough extraction of feature points, fine extraction, characteristic connection, and trace merging. Among them, the rough extraction of feature points is implemented based on the idea of ​​Inductive Transfer Learning. Specifically, trace extraction includes two aspects: extracting traces based on regular body point clouds and extracting traces based on actual rock point cloud data. The former is called the source task and the latter is called the target task. The purpose of this embodiment is to use the knowledge learned from regular body point cloud data (source task) to improve the performance of extracting traces from actual rock point cloud data (target task), and to improve the generalization ability of the model by sharing knowledge between related tasks.

[0084] In order to learn the knowledge in the regular body point cloud, this embodiment first extracts a sample point set from the regular body point cloud data, and each sample point in the sample point set contains label information, where the label refers to the category label or identification to which the sample point belongs, such as the sample point is a plane point, or the sample point is a sharp point.

[0085] In one example, manual extraction can be used to collect plane points and edge points (sharp points) in the regular volume point cloud, including 7000 plane points and 3148 sharp points, and a total of 10148 sample points (feature points) are collected.

[0086] Step S102: Calculate the multi-scale surface change index of each sample point in the sample point set to obtain a first feature set of the sample point set; calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data.

[0087] The goal of step S102 is to calculate the feature information of the rock mass point cloud. Among them, the surface variation index refers to the parameter used to evaluate the morphological changes of the rock mass surface, which is used to describe the geometric characteristics, deformation, etc. of the rock mass surface. The multi-scale surface variation index (Multi-scale Surface Variation Index, MsSVI for short) refers to calculating the geometric characteristics of the rock mass surface at each point at multiple scales to form a set of surface variation indicators of the point at different scales, that is, a feature set.

[0088] Preferably, in one embodiment, the surface change index can be calculated by a principal component analysis method, specifically, the principal component analysis method is used to quantify the main direction of the surface change at each point of the rock mass, and the quotient of the first eigenvalue and the sum of all eigenvalues ​​is calculated to quantify the contribution of the deformation in this main direction to the overall deformation. The specific calculation method will be described in detail later.

[0089] It should be noted that, in addition to the calculation method of the principal component analysis method mentioned above, the surface change index can also adopt traditional surface change indicators, such as the following parameters: amplitude parameters (such as the mean square value M and root mean square R of the rock joint boss height), slope parameters (such as the first-order derivative of the rock joint boss height (slope Z2), the second-order derivative of the boss height (curvature Z3)), mixed parameters (such as structure function S, autocorrelation function AC and spectral density function), etc., or other parameters that describe the surface morphological characteristics of the rock mass, which are not limited to this in this embodiment. In other words, several achievable forms of surface change indicators are given above, but in actual applications they are not limited to the above-mentioned methods. In short, parameters obtained based on rock point cloud computing and used to evaluate changes in rock surface morphology also fall within the scope of protection of this application.

[0090] In order to calculate the indicators at different scales, it is necessary to first clarify the scale range of the analysis (scale n), that is, at how many scales to calculate the surface change indicators at different scales. In this embodiment, the scale n can be determined according to the research purpose and the availability of data, or it can be manually specified based on the domain knowledge of experts. This embodiment does not limit the specific method for determining the scale n.

[0091] Exemplarily, the scale n can be 5, 6 or 7.

[0092] For the convenience of description, the surface change index of the sample points of the regular body point cloud at different scales is called the first feature set, and the surface change index of the rock mass point cloud data at different scales is called the second feature set. Among them, the total number of features of the first feature set and the total number of features of the second feature set are equal to the scale n. The i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale. Thus, a matrix with different data points as rows and surface change indexes at different scales as columns (features) is formed.

[0093] Step S103: Using the first feature set as input features and combining the label information contained in each sample point in the sample point set, a training data set is constructed; and using the training data set, based on the idea of ​​induced transfer learning, an artificial neural network technology ANN is used to train a general trace feature point extraction model, and then the second feature set is input into the trained trace feature point extraction model to perform rough extraction of the rock point cloud trace feature points to obtain a rough extracted feature point set.

[0094] It should be noted that in the training process of the general trace feature point extraction model, in order to speed up the convergence of the model, M feature points can be collected in the rock point cloud data to make a test data set. For example, 3300 plane feature points and 2668 sharp feature points can be collected in the rock point cloud, for a total of 5968 feature points, and these points can be made into a test data set.

[0095] Step S104: Based on the coarsely extracted feature point set, the rock mass trace feature point fine extraction, feature point connection, and trace merging processing are performed in sequence to obtain the final rock mass trace.

[0096] Fine extraction of feature points is also called local optimization of feature points. The accuracy of the trace can be improved through local optimization. Specifically, after obtaining a set of roughly extracted feature points, local optimization of rock traces can be achieved through clustering or regional growing algorithms, or local optimization can be achieved based on machine learning technology according to the curvature change of the three-dimensional point cloud. The specific implementation of the above-mentioned fine extraction can be carried out with reference to the existing technology. For the sake of simplicity, this embodiment will not be described in detail.

[0097] The trace feature points after precise extraction are connected and merged to achieve trace integration and obtain the final rock trace.

[0098] In summary, the technical solution provided in this embodiment uses the model trained by transfer learning ideas and ANN technology to avoid manual threshold selection during trace feature extraction and improve the feature extraction accuracy.

[0099] In some embodiments, regular body point cloud data or rock body point cloud data is used as input point cloud P, and any sample point in the sample point set or any point in the rock body point cloud data is used as the current point p; after reading the input point cloud P, the characteristic information calculation steps of the current point p are as follows:

[0100] Step S121: Search for the k points closest to the current point p among all the points in the input point cloud P to form a neighbor point set N k (p).

[0101] That is to say, for any point p in the input point cloud P, a smaller number of neighboring points k is defined, and the set of k nearest neighbor points of point p is N k (p). For ease of description, we can define the neighbor point set N k (p) is called the first nearest neighbor point set.

[0102] Furthermore, in order to improve the neighbor point set N k The quality of (p) also includes: removing abnormal points according to the density index. Specifically, for p i ∈N k (p), define p i The Euclidean distance between p and i ,p), the density index D(p) of point p is calculated as:

[0103]

[0104] Where h is the weight function, which is used for the distribution of spatial weights; |N k (p)| is the number of neighboring points.

[0105] The density threshold is set in advance. If the density index D(p)>density threshold, point p is judged as an abnormal point, otherwise it is a non-abnormal point. By calculating the density index of each point in the point cloud and removing points with excessive density (usually noise or outliers), the quality of point cloud data can be effectively improved. In addition, this local density-based filtering method can adaptively process data density changes in different regions, thereby more accurately identifying and removing noise.

[0106] Step S122: Calculate the neighbor point set N k (p) Any point p i The Euclidean distance to the current point p, and determine the set of neighboring points N k The maximum distance d from all points in (p) to the current point p k (p).

[0107] It should be noted that the neighboring point set N of the outlier point after removal is k The number of points in (p) may be less than or equal to the set of neighboring points N before removal.k (p) is the number of points in the matrix. k For all non-abnormal points in (p), calculate the maximum distance d from them to the current point p k (p) (also called the kth distance):

[0108] d k (p) = max{d(p i ,p)|p i ∈N k (p)} (2)

[0109] Step S123: k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, where the number of radii contained in the multi-scale search radius set is equal to the scale n.

[0110] In some embodiments, the maximum acquisition range is determined by the following formula:

[0111]

[0112] In the formula, r max Indicates the maximum acquisition range, d bb Represents the diagonal length of the bounding box of the input point cloud P, and |N| represents the number of points in the input point cloud P.

[0113] Among them, the point cloud bounding box is an algorithm for solving the optimal bounding space of a discrete point set, which approximately replaces complex geometric objects with geometric bodies with slightly larger volumes and simpler characteristics (such as cubes or directed bounding boxes, etc.).

[0114] In formula (3), the diagonal length of the bounding box is a measure of the overall range of the point cloud, which can reflect the overall size of the point cloud in three-dimensional space. Using it as a basis to determine the maximum acquisition range can ensure that the acquisition range is reasonable in the global range, avoiding the deviation that may be caused by considering only local information; at the same time, the cube root of the number of points reflects the density of the point cloud to a certain extent. By multiplying by 3 and dividing by the cube root of the number of points in the point cloud, the change in point cloud density is taken into account while maintaining globality. By adjusting the maximum acquisition range in this way, a more appropriate acquisition range can be obtained in point cloud areas of different densities, achieving a comprehensive consideration of adaptability and globality.

[0115] Step S124: Taking the current point p as the center and taking each radius value in the multi-scale search radius set as the search radius, a search is performed to obtain a neighbor point set of the current point p at each scale.

[0116] Step S125: Calculate the surface variation index value of the current point p at each scale based on the neighbor point set of the current point p at each scale.

[0117] Among them, the neighbor point set of the current point p at each scale can be recorded as In the formula, Represents the neighboring point set of the current point p at the ni-th scale. For ease of description, we can is called the second nearest neighbor point set. The second nearest neighbor point set at all scales constitutes the NN representation. The surface change index values ​​of the current point p at various scales, such as Z2, Z3, AC, etc., can be calculated by a point cloud algorithm. The specific calculation method can refer to the prior art. For the sake of brevity, this embodiment will not be described in detail here.

[0118] From the above steps, it can be seen that by gradually expanding the search radius to obtain neighboring point sets at different scales, the trace feature point extraction model can analyze the local geometric characteristics of each point at multiple scales, revealing features at different levels through analysis at different scales, which is conducive to capturing the details of complex shapes and structures and improving the accuracy of trace feature point extraction.

[0119] On the basis of the above embodiments, in some optional embodiments, in addition to traditional indicators such as Z2, Z3, AC, etc., the values ​​of surface change indicators at various scales can also be taken by principal component analysis method to quantify the main direction of surface change at each point of the rock mass, and calculate the quotient of the first eigenvalue and the sum of all eigenvalues ​​to quantify the contribution of deformation in this main direction to the overall deformation.

[0120] The implementation method is described in detail below, and can be implemented specifically according to the following steps:

[0121] Step S1251: For the neighbor point set (second neighbor point set) of the current point p at each scale, construct a neighbor point matrix of the current point p at each scale. The number of rows of the neighbor point matrix of the current point p at each scale is equal to the number of points included in the neighbor point set of the current point p at the scale. The columns of the neighbor point matrix include three columns, namely the x-coordinate, the y-coordinate and the z-coordinate.

[0122] That is to say, in the neighbor point matrix of each scale, the number of neighbor points is the number of rows of the matrix (assuming k points). The columns are the three-dimensional coordinates (x, y, z) of each point, each coordinate is a column, a total of 3 columns. And different scales correspond to different neighbor point matrices.

[0123] Step S1252: Based on the neighbor point matrix of the current point p at each scale, calculate the covariance matrix C of the current point p at each scale.

[0124] Among them, the covariance matrix C is a 3×3 square matrix, denoted as: C 3×3 , and its calculation formula is as follows:

[0125]

[0126] In the formula, is the center point of the neighbor point set (the second neighbor point set) (i.e. the mean of the three-dimensional coordinates), p i ={p1,p2,…p k}, k is the number of points in the neighboring point set.

[0127] Step S1253: Perform eigenvalue decomposition on the covariance matrix C of the current point p at each scale through principal component analysis (PCA), and arrange the generated eigenvalues ​​in order from large to small to obtain the eigenvalues ​​of the matrix C λ0≤λ1≤λ2, thereby determining the maximum eigenvalue of the current point p at each scale and the eigenvalues ​​ranked in the top 3.

[0128] Step S1254: Calculate the ratio of the maximum eigenvalue of the current point p at each scale to the sum of the eigenvalues ​​ranked in the top three positions, and obtain the surface change index of the current point p at each scale.

[0129] The expression of surface variation index (SVI) is as follows:

[0130]

[0131] Where σ(p) is the surface variation index, and λ0, λ1, and λ2 are the eigenvalues ​​of the matrix C from large to small.

[0132] Since PCA can identify the main direction of change in the data, that is, the direction with the largest data variation, in the rock mass point cloud, the direction with the largest data variation usually corresponds to the main direction of rock mass surface deformation. By calculating the quotient of the first eigenvalue λ0 and the sum of all eigenvalues ​​(λ0, λ1, λ2), the contribution of the deformation in the main direction of the rock mass to the overall deformation can be quantified. The larger σ(p) is, the more the surface changes of the rock mass are mainly concentrated in this direction, which helps to identify the main deformation mode of the rock mass. In other words, using the quotient of the first eigenvalue λ0 and the sum of all eigenvalues ​​(λ0, λ1, λ2) as the surface change index can effectively simplify the large number of points and complex deformation information of the rock mass point cloud data, further simplify the analysis process, and only need to focus on the most important deformation direction, while ignoring those secondary directions with smaller contributions, so as to improve the accuracy and reliability of the analysis.

[0133] Repeat steps S1251 to S1254 until all points in the input point cloud P are calculated.

[0134] In some embodiments, d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, which can be achieved in the following two ways:

[0135] Method 1: Using d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. A fixed step size is used to set n equidistant breakpoints within the upper and lower limits of the search radius to form a multi-scale search radius set.

[0136] The search radius is gradually increased by a fixed step size, so that each increase in the search radius is controllable. Sensitivity analysis can be performed on this basis, that is, observing how small changes in the search radius affect the search results, which is conducive to accurately controlling the search range. At the same time, the fixed step size method is relatively simple and intuitive, easy to implement in algorithms or software, reducing the difficulty and cost of technical implementation, making this method more feasible in practical applications. However, due to the difference in density of different rock point clouds, the search method with a fixed step size increase may result in the inability to search for neighboring points under certain step sizes, thereby resulting in data loss. For this reason, in some exemplary embodiments using a fixed step size method, after step S124, it also includes: for the neighbor point set of the current point p at each scale, it is determined whether the number of points in the neighbor point set at each scale is greater than a preset number threshold. If so, step S125 is executed. If not, step S125 is not executed, but the surface change index value at the corresponding scale is assigned to a preset fixed value. In this way, the problems caused by missing data can be compensated. At the same time, the fixed step size method is used to gradually increase the search radius and calculate the MsSVI information. The rough extraction of extreme feature points combined with machine learning can, to a certain extent, solve the problem of difficult threshold selection caused by different data densities of different rock masses.

[0137] Method 2: Using d k (p) is the lower limit of the search radius. An adaptive iteration method is used to obtain n breakpoints within the upper and lower limits of the search radius and form a multi-scale search radius set.

[0138] The adaptive iteration method is a method proposed to solve the problem that some search radius data may be missing under a fixed step size. The specific steps of the adaptive iteration include:

[0139] Assume d k (p) is the starting radius r0 of the iteration, the distance from the current point p to the neighboring point set N kThe mean value d of the Euclidean distances of all points in (p) (the first nearest neighbor point set) mean (p) is used as the step size Δr0 of the first iteration, then the radius value r1 generated by the first iteration is:

[0140] r1=r0+Δr0 (7)

[0141] Among them, the average distance (mean of the Euclidean distance) of point p is d mean The calculation formula of (p) is as follows:

[0142]

[0143] The meanings of the parameters in the formula are the same as those in formula (1) and will not be repeated here.

[0144] Then, the adaptive step size of the k1th iteration and the breakpoint value of each iteration (i.e., the search radius value) are calculated:

[0145] r k1+1 =r k1 +Δr k1 ,

[0146] Δr k1 =Δr0·tanh(I(R(r k1 ))),

[0147]

[0148] I(R(r k1 ))=-log P(R(r k1 ) ((8)

[0149] In the formula, k1 is the current iteration number, k1 is a positive integer greater than 1, r k1+1 、r k1 are the breakpoint values ​​at the k1+1th and k1th iterations, i.e., the search radius values ​​at the k1+1th and k1th under multi-scale, Δr k1 is the step size of the k1th iteration, R(r k1 ) is the normalization function, I(R(r k1 )) represents the normalized r k1 The amount of information, P(R(r k1 )) is R(r k1 ), that is, the point cloud probability under the current search radius.

[0150] In formula (7), the first search radius is calculated based on the starting radius and the step size of the first iteration. For the iterations after the first one, the current search radius r is calculated first. k1The amount of information I under the condition is mapped to [-1,1] by the hyperbolic tangent function tanh, and multiplied by the step size of the initial iteration to obtain the adaptive step size Δr k1 , and then get the next search radius r k1+1 . This formula maps the amount of information to between [-1,1] through the hyperbolic tangent function tanh, and multiplies it by the step size of the initial iteration, so that the subsequent step size can grow according to the curve law, that is, the step size is adaptively adjusted according to the change in the amount of information. In this way, when the amount of information is large, the step size increases accordingly, thereby speeding up the search speed; when the amount of information is small, the step size decreases accordingly to improve the search accuracy. This adaptive adjustment helps to find a balance between search speed and accuracy. At the same time, this adaptive step size strategy helps to avoid falling into the local optimal solution during the multi-scale search process. When approaching the local optimum, the amount of information may decrease, resulting in a decrease in the step size, allowing the algorithm to explore the area more carefully and possibly find a better solution.

[0151] In summary, the adaptive step size is used to gradually increase the search radius, calculate the adaptive radius MsSVI information, and combine it with machine learning to perform rough extraction of extreme feature points, which solves the problem of difficult threshold selection caused by different data densities of different rock masses. At the same time, the adaptive step size method can automatically adapt to features of different scales and realize effective feature representation and recognition, making it more flexible and effective when processing point cloud data of different rock masses.

[0152] In some embodiments, the precise extraction of rock mass trace feature points includes:

[0153] Step S141: stripping the feature points in the coarsely extracted feature point set from the rock mass point cloud data to obtain a trace feature point cloud set Among them, J represents the trajectory feature point cloud set P trace The total number of midpoints, R is the set of real numbers, p j Represents the trajectory feature point cloud set P trace The jth feature point in .

[0154] Step S142: Using random sampling technology to collect the trace feature point cloud set P trace Perform random sampling to obtain a set of sampling points Where I represents the total number of points in the sampling point set X, x i Represents the i-th sampling point in the sampling point set.

[0155] Step S143: Based on the trace feature point cloud set and the sampling point set, an optimization target equation is constructed using an L1-skeleton extraction algorithm. The expression of the optimization target equation is as follows:

[0156]

[0157] In the formula, x i is the i-th point in the sampling point set X, p j is the jth point in the trace feature point cloud set, I and J are the number of points in the sampling point set and the trace feature point cloud set, θ is the weight function, γ i is the balance coefficient, x i′ is the i′th point in the sampling point set, and i′≠i, σ i It is point x i The degree of linearity at the location, ||.|| represents the L1 norm of the matrix, that is, the absolute value and total value of all elements of the matrix are calculated, i′∈I\{i} means that the range of i′ is the set of all values ​​of I minus i.

[0158] Among them, the expression of the weight function θ is as follows:

[0159]

[0160] Where r is the input variable and h1 is the parameter representing the range.

[0161] σ i The expression is as follows:

[0162]

[0163] In the formula, σ i It is point x i The linearity of the location, It's point p i The eigenvalues ​​are arranged from large to small.

[0164] In the optimization objective equation expressed by formula (9), the purpose of the first term is to move the points in the sampling point set X to the trace feature point cloud set P according to the distance weight (i.e., weight function θ) trace The second term is the penalty term, whose purpose is to control the points in the sampling point set X to be as close as possible to the trace feature point cloud set P trace The linearity is high, and the distance between the points in the sampling point set X is kept constant. In this way, the point positions in the sampling point set X can be continuously sampled so that they are distributed in the trace feature point cloud set P. trace on the skeleton position.

[0165] Step S144: iteratively solving the optimization target equation to obtain precisely extracted rock mass trace feature points.

[0166] Specifically, the solution can be obtained by the alternating direction multiplier method (ADMM) or other optimization algorithms. During the solution process, the point positions in the sampling point set X are iteratively updated to meet the constraints and minimize the objective function, thereby obtaining the precisely extracted rock mass trace feature points. The update formula for the point positions in the sampling point set X is as follows:

[0167]

[0168] In the formula, k represents the number of iterations, μ, are constant parameters, The remaining parameters are the same as those in formula (9) and will not be described in detail here.

[0169] According to formula (12), iterative solution is performed to obtain the precisely extracted rock mass trace feature points, which are denoted as P precise .

[0170] In some embodiments, a weighted line growing algorithm is used to automatically connect rock mass traces, including the following steps:

[0171] Step S151: extract the rock mass trace feature points P precise Perform outlier removal.

[0172] The steps for removing outliers are the same as those for removing abnormal points, and will not be described in detail here.

[0173] Step S152: randomly select a seed point p0 from the precisely extracted rock mass trace feature points as the center point, and calculate the neighbor point set of point p0. For the sake of distinction, it is called the third neighbor point set N p0 ,

[0174]

[0175] Step S153: Calculate the weight of each neighbor point in the neighbor point set of the center point. The calculation formula is as follows:

[0176]

[0177] Where pli represents the li-th seed point, σ(p li ) indicates directionality, d(p l ,p li ) represents the distance, s(p li ) represents the angle, λ s , t and λ d They are the control angles s(p li ), directionality σ(pli ) and the distance d(p l ,p li )'s balance coefficient.

[0178] In formula (13), by combining cosine similarity, directionality and distance to calculate the weight, the relationship between the precisely extracted rock mass trace feature points and their neighboring points can be considered more comprehensively, thereby improving the accuracy of trace connection.

[0179] Step S154: According to the calculated weight result, select the point with the largest weight as the next seed point, ie, the new seed point, and connect the new seed point and the center point.

[0180] Step S155: Take the new seed point as the new center point, find the neighbor point set of the center point and calculate the weight of each neighbor point, repeat the above steps S153 to S154 until the trace is connected and a connected trace is obtained.

[0181] In some embodiments, the rock mass trace merging process is specifically performed by using a trace merging method based on distance and collinearity drive, including the following steps:

[0182] Step S161: The connected traces are smoothed using the B-spline interpolation method to obtain smoothed traces.

[0183] Smoothing through B-spline interpolation can effectively remove or weaken the influence of noise points, making the trace smoother, thereby improving the accuracy of merging. At the same time, B-spline interpolation can maintain the high-order continuity of the data, making the smoothed point cloud data more continuous in shape and trend, and closer to the real trace.

[0184] Step S162: Collect all endpoints of the smoothed trace and construct an endpoint set P v ={p0,p1,p2,...}.

[0185] Step S163: At the endpoint set P v Take any point as the current center point p i , perform the following steps:

[0186] Search for the neighboring points of the center point and obtain the neighboring point set of the center point For the sake of distinction, it is referred to as the fourth nearest neighbor point set.

[0187] Define the center point p i The vector to the jth neighbor point in its neighbor point set is Center point p i The vector pointing to the adjacent point on its trajectory is Defining vectors and The angle between them is θ j , according to the center point p i The Euclidean distance and angle between the jth neighbor point and the jth neighbor point are θ j Calculate the adjacent trace endpoints to the current center point p i The Euclidean distance d of the trace tt , Euclidean distance d tt The calculation method is as follows:

[0188]

[0189] In the formula, is the center point p i The Euclidean distance to its jth nearest neighbor, also denoted as d pp .

[0190] If the angle is θ j Less than the preset angle threshold θ deg , and the distance d tt Less than the preset distance threshold θ dis , determine the jth neighbor point and the center point p i Is it connected? If it is connected, ignore the point; if it is not connected, assign the jth neighbor point to the candidate point set and calculate the weight of the jth neighbor point; the expression of the weight is:

[0191]

[0192] In the formula, ω j is the weight of the jth neighbor point, λ0, λ1, λ2 are the eigenvalues, d pp is the center point p i The Euclidean distance to its jth nearest neighbor, d tt It is the adjacent trace endpoint to the current center point p i Euclidean distance of the trace.

[0193] In the weight formula (15), the center point p is introduced i The distance between its neighboring points and the adjacent trace endpoints to the current center point p i The distance of the trace makes the weight calculation very sensitive to the change of spatial position, which helps to better capture the local and global spatial relationship when processing point cloud data. jThe introduction of makes the weight calculation directional, which can better reflect the method of neighboring points relative to the center point. By using the eigenvalues ​​λ0, λ1, and λ2, the relationship between the angle and the distance can be effectively adjusted, which can more comprehensively describe the relationship between the center point and the neighboring points, thereby improving the accuracy and robustness of the weight calculation. In summary, this embodiment adopts a trace merging method based on distance and collinearity drive to complete the merging of trace breakpoints. Combined with distance measurement and collinearity features, it can more accurately connect and reconstruct broken traces, thereby optimizing the trace extraction results and improving the accuracy of subsequent analysis.

[0194] Step S164: Take the point with the largest weight in the candidate point set as the connection point and connect it to the center point p i connect.

[0195] Step S165: Repeat the above steps until all endpoints in the endpoint set are calculated.

[0196] In summary, the main advantages of the method provided in this embodiment are: (a) the model trained by using the idea of ​​transfer learning and ANN technology can avoid artificial threshold selection during trace feature extraction, thereby improving the accuracy of feature extraction; (b) the L1 skeleton extraction algorithm based on random sampling is used, and no point cloud triangulation operation is required, so it can be directly used for point cloud data; (c) in view of the phenomenon that the traces connected by the line growth algorithm are mostly discontinuous and have jagged edges, B-spline interpolation is used for smoothing, and a trace merging algorithm driven by the angle and distance of the trace endpoints is proposed to reduce the phenomenon of many trace discontinuities. Through trace integration, the continuity of the trace results can be improved, which is closer to the real trace.

[0197] System Example

[0198] This embodiment provides a rock mass trace automatic extraction system based on point cloud. Figure 2 As shown, the system includes: a data acquisition module 201, a calculation module 202, a rough extraction module 203 and a post-processing module 204. Among them:

[0199] The data acquisition module 201 is configured to acquire regular body point cloud data and rock body point cloud data, and extract a sample point set from the regular body point cloud data, each sample point in the sample point set contains label information;

[0200] The calculation module 202 is configured to calculate the multi-scale surface change index of each sample point in the sample point set to obtain a first feature set of the sample point set; calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data; wherein the total number of features of the first feature set and the total number of features of the second feature set are both equal to the scale n, and the i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale;

[0201] The rough extraction module 203 is configured to use the first feature set as input features, and to construct a training data set in combination with label information contained in each sample point in the sample point set; and to use the training data set, in combination with the idea of ​​induced transfer learning, to train a general trace feature point extraction model using artificial neural network technology ANN, and then to input the second feature set into the trained trace feature point extraction model, to perform rough extraction on the rock mass point cloud trace feature points, and obtain a rough extraction feature point set;

[0202] The post-processing module 204 is configured to perform, based on the coarsely extracted feature point set, precise extraction of rock mass trace feature points, feature point connection, and trace merging processing in sequence to obtain the final rock mass trace.

[0203] The point cloud-based rock mass trace automatic extraction system provided in this embodiment can implement the steps and processes of the point cloud-based rock mass trace automatic extraction method provided in any of the above embodiments, and achieve the same technical effects, which will not be repeated here one by one.

[0204] The above description is only a preferred embodiment of the present application and is not intended to limit the present application. For those skilled in the art, the present application may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present application shall be included in the protection scope of the present application.

Claims

1. A method for automatic extraction of rock mass traces based on point cloud, characterized in that: include: Obtaining regular body point cloud data and rock mass point cloud data, and extracting a sample point set from the regular body point cloud data, wherein each sample point in the sample point set contains label information; Calculate the multi-scale surface change index of each sample point in the sample point set to obtain a first feature set of the sample point set; calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data; wherein the total number of features of the first feature set and the total number of features of the second feature set are both equal to the scale n, the i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale; Taking the first feature set as input features, combined with the label information contained in each sample point in the sample point set, a training data set is constructed; and using the training data set, based on the idea of ​​induced transfer learning, using artificial neural network technology ANN to train a general trace feature point extraction model, and then inputting the second feature set into the trained trace feature point extraction model to perform rough extraction of rock mass point cloud trace feature points to obtain a rough extracted feature point set; Based on the coarsely extracted feature point set, the rock mass trace feature point fine extraction, feature point connection and trace merging processing are performed in sequence to obtain the final rock mass trace.

2. The method according to claim 1, characterized in that The regular body point cloud data or the rock mass point cloud data is used as the input point cloud P, and any sample point in the sample point set or any point in the rock mass point cloud data is used as the current point p; the calculation steps of the multi-scale surface change index of the current point p are as follows: Among all the points in the input point cloud P, search for the k points closest to the current point p to form a neighbor point set N k (p); Calculate the neighbor point set N k (p) Any point p i The Euclidean distance to the current point p, and determine the set of neighboring points N k The maximum distance d from all points in (p) to the current point p k (p); With d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, wherein the number of radii contained in the multi-scale search radius set is equal to the scale n; Taking the current point p as the center and taking each radius in the multi-scale search radius set as the search radius, search to obtain the neighboring point set of the current point p at each scale; Based on the neighbor point sets of the current point p at each scale, the surface change index values ​​of the current point p at each scale are calculated.

3. The method according to claim 2, characterized in that Based on the neighboring point sets of the current point p at each scale, the surface change index values ​​of the current point p at each scale are calculated, specifically: For the neighbor point set of the current point p at each scale, construct a neighbor point matrix of the current point p at each scale, the number of rows of the neighbor point matrix of the current point p at each scale is equal to the number of points included in the neighbor point set of the current point p at the scale, and the columns of the neighbor point matrix include three columns, namely, x coordinate, y coordinate and z coordinate; Based on the neighbor point matrix of the current point p at each scale, the covariance matrix C of the current point p at each scale is calculated; Performing eigenvalue decomposition on the covariance matrix C of the current point p at each scale by principal component analysis, and arranging the generated eigenvalues ​​in descending order to determine the maximum eigenvalue and the top three eigenvalues ​​of the current point p at each scale; The ratio of the maximum eigenvalue of the current point p at each scale to the sum of the eigenvalues ​​arranged in the first three positions is calculated to obtain the surface change index of the current point p at each scale.

4. The method according to claim 2 or 3, characterized in that: The maximum acquisition range is determined by the following formula: In the formula, r max Indicates the maximum acquisition range, d bb Represents the diagonal length of the bounding box of the input point cloud P, and |N| represents the number of points in the input point cloud P.

5. The method according to claim 4, characterized in that With d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. Multiple breakpoints are determined within the upper and lower limits of the search radius to form a multi-scale search radius set, including: With d k (p) is the lower limit of the search radius, and the maximum acquisition range of the input point cloud P is the upper limit of the search radius. A fixed step size is used to set n equidistant breakpoints within the upper and lower limits of the search radius to form the multi-scale search radius set; or, With d k (p) is the lower limit of the search radius, and n breakpoints within the upper and lower limits of the search radius are obtained by adaptive iteration to form the multi-scale search radius set; The adaptive iteration includes: Assume d k (p) is the starting radius r0 of the iteration, the current point p to the neighboring point set N k The mean of the Euclidean distances of all points in (p) is used as the step size Δr0 of the first iteration, and the radius value r1 generated by the first iteration is: r1=r0+Δr0, Calculate the adaptive step size of the k1th iteration and the breakpoint value of each iteration: r k1+1 =r k1 +Δr k1 , Δr k1 =Δr0·tanh(i(R(r k1 ))), i(R(r k1 ))=-log P(R(r k1 )), Where k1 is the current iteration number, r k1+1 、r k1 are the breakpoint values ​​at the k1+1th and k1th iterations, Δr k1 is the step size of the k1th iteration, R(r k1 ) is the normalization function, i(R(r k1 )) represents the normalized r k1 The amount of information, P(R(r k1 )) is R(r k1 )’s probability distribution.

6. The method according to claim 1, characterized in that Precise extraction of rock mass trace feature points, including: The feature points in the coarsely extracted feature point set are stripped from the rock mass point cloud data to obtain a trace feature point cloud set. The trace feature point cloud set is randomly sampled using a random sampling technique to obtain a sampling point set; Based on the trace feature point cloud set and the sampling point set, an optimization target equation is constructed using the L1 median skeleton extraction algorithm. The expression of the optimization target equation is as follows: In the formula, x i is the i-th point in the sampling point set, p j is the jth point in the trace feature point cloud set, i and J are the number of points in the sampling point set and the trace feature point cloud set, θ is the weight function, γ i is the balance coefficient, x i′ is the i′th point in the sampling point set, and i′≠i, σ i It is point x i The linearity of the location; The optimization objective equation is iteratively solved to obtain precisely extracted rock mass trace feature points.

7. The method according to claim 6, characterized in that The rock mass trace feature point connection includes the following steps: Eliminating external points of the precisely extracted rock mass trace feature points; Randomly select a seed point p0 from the precisely extracted rock mass trace feature points as the center point, and search for a set of neighboring points of the center point; For each neighbor point in the neighbor point set of the center point, calculate its weight; According to the calculated weight result, select the point with the largest weight as a new seed point, and connect the new seed point and the center point; Take the new seed point as the new center point, find the neighbor point set of the center point and calculate the weight of each neighbor point, repeat the above steps until the trace is connected and get the connected trace.

8. The method according to claim 7, characterized in that The rock mass trace merging process is specifically to use a trace merging method based on distance and collinearity drive to merge the rock mass traces, including the following steps: The connected traces are smoothed using a B-spline interpolation method to obtain a smoothed trace; Collecting all endpoints of the smoothed trace to construct an endpoint set; Take any point in the endpoint set as the current center point p i , perform the following steps: Search for the neighboring points of the center point and obtain the neighboring point set of the center point Define the center point p i The vector to the jth neighbor point in its neighbor point set is Center point p i The vector pointing to the adjacent point on its trajectory is Defining vectors and The angle between them is θ j , according to the center point p i The Euclidean distance and angle between the jth neighbor point and the jth neighbor point are θ j Calculate the adjacent trace endpoints to the current center point p i The Euclidean distance d of the trace tt ; If the angle is θ j is smaller than the preset angle threshold, and the distance d tt If the distance is less than the preset threshold, the jth neighbor point is assigned to the candidate point set and the weight of the jth neighbor point is calculated; The point with the largest weight in the candidate point set is used as the connection point and connected to the center point p i connect; The above steps are repeated until all endpoints in the endpoint set are calculated.

9. The method according to claim 8, characterized in that The weight of the jth neighbor point is expressed as follows: In the formula, ω j is the weight of the jth neighbor point, λ0, λ1, λ2 are the weights of the center point p i The characteristic value of pp is the center point p i The Euclidean distance to its jth nearest neighbor, d tt It is the adjacent trace endpoint to the current center point p i Euclidean distance of the trace.

10. A point cloud-based rock mass trace automatic extraction system, characterized in that: include: A data acquisition module is configured to acquire regular body point cloud data and rock body point cloud data, and extract a sample point set from the regular body point cloud data, wherein each sample point in the sample point set contains label information; A calculation module, configured to calculate a multi-scale surface variation index of each sample point in the sample point set to obtain a first feature set of the sample point set; Calculate the multi-scale surface change index of all points in the rock mass point cloud data to obtain a second feature set of the rock mass point cloud data; wherein the total number of features of the first feature set and the total number of features of the second feature set are both equal to the scale n, and the i-th feature in the first feature set represents the surface change index value of each sample point in the sample point set at the i-th scale; the i-th feature in the second feature set represents the surface change index value of all points in the rock mass point cloud data at the i-th scale; The rough extraction module is configured to use the first feature set as input features, and to construct a training data set in combination with label information contained in each sample point in the sample point set; and to use the training data set, in combination with the idea of ​​induced transfer learning, to train a general trace feature point extraction model using artificial neural network technology ANN, and then to input the second feature set into the trained trace feature point extraction model, to perform rough extraction on the rock mass point cloud trace feature points, and to obtain a rough extraction feature point set; The post-processing module is configured to perform, based on the coarsely extracted feature point set, precise extraction of rock mass trace feature points, feature point connection, and trace merging processing in sequence to obtain a final rock mass trace.

Images