CN114399603A - Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface - Google Patents

Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface Download PDF

Info

Publication number
CN114399603A
CN114399603A CN202210089365.0A CN202210089365A CN114399603A CN 114399603 A CN114399603 A CN 114399603A CN 202210089365 A CN202210089365 A CN 202210089365A CN 114399603 A CN114399603 A CN 114399603A
Authority
CN
China
Prior art keywords
point cloud
point
points
coal mine
curved surface
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202210089365.0A
Other languages
Chinese (zh)
Inventor
欧阳奇
林燕花
张兴兰
范月欣
杨伟静
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Chongqing University of Technology
Original Assignee
Chongqing University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Chongqing University of Technology filed Critical Chongqing University of Technology
Priority to CN202210089365.0A priority Critical patent/CN114399603A/en
Publication of CN114399603A publication Critical patent/CN114399603A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/08—Projecting images onto non-planar surfaces, e.g. geodetic screens
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
    • G06T3/4038—Image mosaicing, e.g. composing plane images from plane sub-images
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/60—Rotation of whole images or parts thereof
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/40—Image enhancement or restoration using histogram techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10028—Range image; Depth image; 3D point clouds
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing
    • G06T2207/10044—Radar image

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Computer Graphics (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • Processing Or Creating Images (AREA)

Abstract

The invention discloses a three-dimensional curved surface reconstruction method of a coal mine tunnel arch surface, which comprises the following steps: turning the laser radar by 90 degrees, inversely installing the laser radar on a guide rail of a mechanical arm base, and driving a laser radar collector through a hydraulic motor moving device to collect point cloud data of the arch surface of the coal mine tunnel according to a plurality of preset sampling points; removing outliers of point cloud data and simplifying the data by using a plurality of filters in a PCL library; estimating a point cloud normal by using a PCA algorithm; splicing and fusing the multiple point cloud images by using an SAC-IA coarse splicing and ICP fine splicing algorithm; the method comprises the following steps of adopting a Greedy PT algorithm to realize the reconstruction of the three-dimensional curved surface of the arch surface of the coal mine tunnel, and carrying out hole identification and online repair on a reconstructed triangular grid, thereby completing the reconstruction of the arch surface within a specified range, wherein the algorithm realization of the whole process is automatically completed by a computer; the Greedy PT algorithm is simple in principle, strong in real-time performance and high in calculation speed, and can meet the practical application of an engineering scene. The method has the advantages of high point cloud splicing precision and good reconstruction effect of the tunnel arch surface of the coal mine.

Description

Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface
Technical Field
The invention relates to the technical field of three-dimensional curved surface reconstruction, in particular to a three-dimensional curved surface reconstruction method for an arch surface of a coal mine tunnel.
Background
In recent years, the deep development of two-dimensional fusion enables the three-dimensional reconstruction technology based on the engineering environment to be widely applied to many engineering fields such as human-computer interaction, unmanned driving, cultural relic protection and the like, and the three-dimensional reconstruction technology based on the PCL (PCL) (Point Cloud library) library relates to multidisciplinary fusion such as computer graphics, geometric computation, sensors and the like. Because the traditional coal mine tunnel guniting work mostly adopts a manual spraying mode, the automatic anchor spraying robot is researched due to the problems of low working efficiency, frequent safety accidents and the like; in order to meet the actual guniting requirement of the tunnel arch surface of the coal mine, the tunnel arch surface needs to be three-dimensionally reconstructed; the method comprises the following steps that an explosion-proof three-dimensional laser radar is adopted to carry out data acquisition on a tunnel arch surface, however, due to the tunnel environment interference and the limitation of factors of the laser radar, point cloud data of the tunnel arch surface of a coal mine cannot be obtained all at one time, and part of point cloud information can be diluted by the tunnel arch surface, so that the subsequent three-dimensional reconstruction of the tunnel arch surface of the coal mine cannot be completed efficiently, and the complex tunnel environment has higher requirements on the installation position of the laser radar and the processing mode of the point cloud data; the existing three-dimensional curved surface reconstruction algorithm is easy to be interfered by complex environment factors, so that the three-dimensional curved surface reconstruction effect of the arch surface of the coal mine tunnel is poor, and the algorithm robustness is low.
According to the defects of the prior art, a method for accurately, quickly and robustly reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel is needed, wherein the method can effectively acquire the data of the arch surface of the coal mine tunnel in real time.
Disclosure of Invention
In order to solve the problems, the invention provides a three-dimensional curved surface reconstruction method for a tunnel arch surface of a coal mine, which can realize multi-station real-time acquisition of point cloud data of the tunnel arch surface, point cloud data preprocessing and point cloud splicing, and can also quickly and accurately finish the three-dimensional curved surface reconstruction of the tunnel arch surface of the coal mine.
In order to achieve the technical purpose, the technical scheme adopted by the invention is as follows:
a three-dimensional curved surface reconstruction method for an arch surface of a coal mine tunnel is characterized in that a laser radar is adopted to scan the arch surface of the coal mine tunnel in real time to obtain point cloud data, and the three-dimensional curved surface reconstruction is carried out on the collected point cloud data according to the following steps:
pretreatment: turning the laser radar by 90 degrees and inversely installing the laser radar on a guide rail of a mechanical arm base; acquiring tunnel arch surface data, and arranging a plurality of data acquisition stations; driving a laser radar to move on a guide rail at a constant speed at certain time intervals through a hydraulic motor moving device to obtain complete tunnel arch surface point cloud data, and respectively setting an initial point and a moving speed of the laser radar on the guide rail; setting the time interval of point cloud data collected by each observation station in the laser radar moving process;
s1: initializing, moving a laser radar on a guide rail of a mechanical arm base to an initial point, and establishing communication between a computer and the laser radar;
s2: the mechanical arm base guide rail drives the laser radar collector to sequentially move to the data sampling points on the guide rail at a certain speed, the laser collector collects real-time data of the arch surface of the coal mine tunnel according to a preset time interval, and transmits the collected point cloud data to the computer;
s3: the computer analyzes the point cloud data, judges whether noise, outliers and invalid points exist in the point cloud data, and performs point cloud data preprocessing on the point cloud data by adopting a filter;
s4: the computer adopts a principal component analysis method to carry out plane fitting on N points in K neighborhoods of each point according to the result of the data preprocessing of the last step so as to estimate normal vectors and calculate curvatures;
s5: completing point cloud image splicing on a plurality of point cloud images acquired by different sampling points on a guide rail according to key point extraction, SAC-IA rough splicing based on FPFH (field programmable gate flash) feature descriptors and ICP (inductively coupled plasma) fine splicing sequence, so that intersected areas are perfectly overlapped;
s6: the point cloud shape acquired by the arch surface of the coal mine tunnel is non-closed point cloud data, and a greedy triangulation algorithm is adopted to realize three-dimensional curved surface reconstruction according to the characteristics of the point cloud data;
s7: after the reconstruction of the three-dimensional curved surface of the tunnel arch surface of the coal mine is completed, performing online repair on a polygonal cavity generated by diluting part of point cloud data from the tunnel arch surface by adopting a cavity repair algorithm based on a triangular mesh;
s8: and the guide rail of the base of the mechanical arm drives the laser radar to move, so that the three-dimensional curved surface reconstruction of the arch surface of the coal mine tunnel is completed, the three-dimensional curved surface is returned to the initial point, and the next arch surface reconstruction is carried out after the guniting is completed.
By adopting the design, the computer collects the point cloud data of the tunnel arch surface of the coal mine through the laser radar collector, and performs point cloud data preprocessing according to the point cloud data obtained by each survey station to remove noise, discrete points and invalid points generated by interference of factors such as environment, the laser radar and the like, thereby reducing the data volume of the point cloud; after a plurality of simplified point cloud images are obtained, splicing the plurality of point clouds by adopting a point cloud data splicing technology to obtain complete point cloud data in a certain range of the arch surface of the coal mine tunnel; the computer carries out three-dimensional reconstruction of the tunnel arch surface according to the complete point cloud data and carries out online repair on the polygonal cavity generated after reconstruction; the method can realize the automatic process when used for three-dimensional reconstruction of the complex coal mine tunnel arch surface, can repair the triangular mesh which does not meet the requirements after reconstruction, prevents the conditions of uneven thickness, missing spraying and the like caused by direct guniting after reconstruction, reduces the re-spraying rate, and further improves the speed and the efficiency of the guniting work of the anchor spraying robot.
As a preferred scheme of the present invention, in steps S1 and S2, the initial point of movement of the mechanical arm base guide rail is (x0, y0), the movement speed is S, the time interval of the laser radar collector collecting point cloud data is T, and the collector establishes communication between the computer and the laser radar through the IP address;
the guide rail of the mechanical arm base drives a laser radar collector to collect point cloud data at a certain speed and time interval through a hydraulic motor moving device, so that the distance d between different data collection stations and the tunnel arch scanning width w are obtained, and a plurality of point cloud images can be generated in the whole guide rail range; the mechanical arm is arranged on the trolley of the mobile platform.
As a preferred embodiment of the present invention, in step S3, the specific steps of the point cloud data preprocessing are as follows:
s3-1: adopting a statistical filter to perform statistical analysis on K neighborhood points of each point in the point cloud, and calculating the average distance from the K neighborhood points to all the neighborhood points; then calculating the average value mu and the standard deviation sigma of the average distance of each point to determine a distance threshold thresh _ d; removing the point cloud with the average neighborhood distance of a certain point lower than or higher than the threshold range of the certain point according to the point cloud density distribution;
assume an original point cloud dataset P ═ { P ═ Pi(xi，yi，zi) I is more than or equal to 1 and less than or equal to n, the distance threshold is:
thresh_d＝μ+m·σ
the point cloud data set after the outlier filtering is as follows:
P′＝(u-m*σ，u+m*σ)
wherein m is a multiple of standard deviation (generally 1-3);
s3-2: importing P' into a straight-through filter, setting a plurality of dimension directions and different point cloud thresholds, and filtering points outside a parameter range to obtain point cloud data in a specified range;
s3-3: and importing the point cloud data obtained in the step S3-2 into a voxel filter or a conditional filter, and setting a plurality of voxel grid thresholds with different specifications and other conditions to simplify the point cloud data to obtain the point cloud preprocessing result, thereby greatly reducing the calculation amount of a computer.
By adopting the scheme, the computer filters the noise, outliers and invalid points in the point cloud data set, simplifies the original point cloud data, calculates the point cloud data in a specific range, effectively reduces the number of the point clouds and improves the calculation efficiency.
As a preferred scheme of the present invention, in step S4, the specific steps of the PCA algorithm for estimating the point cloud normal are as follows:
s4-1: for any point p in the point cloudi(1 is more than or equal to i and less than or equal to n) inquiring k field points p of the n field pointsij(j is more than or equal to 1 and less than or equal to k), calculating piCentroid of k neighborhood points
Figure BDA0003488538040000031
Figure BDA0003488538040000032
S4-2: constructing a covariance matrix of local features of the point cloud:
Figure BDA0003488538040000033
s4-3: calculating an eigenvalue λ of the covariance matrix C0，λ1，λ2，λ0≤λ1≤λ2And a feature vector v0，v1，v2λ and v are in one-to-one correspondence;
wherein v is0，v1，v2Orthogonal, eigenvectors v1，v2Determining a point piAn optimum tangent plane of0Orthogonal to the tangent plane, so point piNormal n toiCan be approximated by a feature vector v0Represents; lambda [ alpha ]0，λ1，λ2As its degree of variation in the direction of the respective eigenvector;
s4-4: calculating the point piOf (d) curvature τpi：
Figure BDA0003488538040000041
By adopting the scheme, the computer carries out normal and curvature estimation in a k neighborhood on the preprocessed point cloud to obtain the geometric characteristics of the point cloud, so that the feature point detection, the point cloud splicing and the curved surface reconstruction can be better carried out at a later stage.
As an improved scheme of the invention, in step S5, the computer completes point cloud image stitching on a plurality of point cloud images acquired by different sampling points on the guide rail according to the key point extraction, SAC-IA rough stitching based on FPFH feature descriptors and ICP fine stitching sequence, and the specific steps are as follows:
s5-1: extracting key points of the point cloud image by adopting an internal morphology descriptor algorithm:
set point cloud P ═ Pi(xi，yi，zi) I is more than or equal to 1 and less than or equal to n, and for each point piEstablishing a local coordinate system and setting a search radius prDetermining piIs the center of a sphere, prAll points within the radius sphere;
calculating the weight ωij：
Figure BDA0003488538040000042
|pi-pj|＜pr
Calculate each point piCovariance matrix of (2):
Figure BDA0003488538040000043
feature values of covariance matrix
Figure BDA0003488538040000044
Arranging in descending order; setting a threshold η1And η2Satisfy the following requirements
Figure BDA0003488538040000045
And
Figure BDA0003488538040000046
the points of (1) are key points, and iteration is carried out until all the key points are found;
s5-2: calculating the fast point feature histogram features of the point cloud to be spliced and the target point cloud to obtain each calculation point MpRelative relation with all the neighborhood points is establishedA simple point feature histogram; calculating the FPFH signature from the SPFH signature, denoted as F (M)P)：
Figure BDA0003488538040000047
Wherein d isiIs the Euclidean distance of the corresponding point pair;
s5-3: splicing the model point cloud and the target point cloud by using a sampling consistency initial registration algorithm; the algorithm comprises the following steps: (1) selecting m feature points to be registered from the model point cloud; (2) searching points similar to the FPFH characteristics of the model point cloud in the target point cloud, and selecting the point closest to the target point cloud as a corresponding relation point; (3) calculating a rotation and translation matrix of the corresponding point pair; the Huber function is used to represent the distance error sum function after the corresponding point pair is translated by rotation, and is recorded as
Figure BDA0003488538040000048
Figure BDA0003488538040000051
Wherein r is a distance threshold, | miI represents the Euclidean distance between the corresponding points of the ith group after transformation, minH (m)i) The corresponding is the optimal transformation matrix after the rough splicing;
s5-4: performing ICP fine splicing by taking the source point cloud P' and the target point cloud G after the rough splicing as input point clouds, and accelerating corresponding point pair search by adopting an optimal node to preferentially optimize a KD tree; all points P ' of point cloud P ' to be registered 'iSearching the corresponding point G with the nearest distance in the target point cloud GiForming corresponding point pairs; the rotation matrix R and the translational vector T are calculated such that the root mean square error M between corresponding pairs of pointskMinimum:
Figure BDA0003488538040000052
finally, a threshold value alpha (i.e. M) is setk-Mk+1< alpha) and maximum number of iterations Nmax。
By adopting the scheme, the computer carries out point cloud splicing on a plurality of point cloud images acquired by different sampling points on the guide rail of the mechanical arm base to obtain complete point cloud data within the range of 1.6m of the arch surface of the coal mine tunnel; therefore, the quick and efficient point cloud splicing technology is realized.
As an improved scheme of the present invention, the specific steps of implementing three-dimensional curved surface reconstruction by using greedy projection triangulation algorithm in step S6 are as follows:
as an improved scheme of the invention, a greedy projection triangulation algorithm is adopted in the step S6 to realize the three-dimensional reconstruction of the coal mine tunnel arch surface, and the method specifically comprises the following steps:
s6-1: adopting the recursive algorithm idea of dynamic programming problem to make three-dimensional points pass through the normal vector obtained in step 4
Figure BDA0003488538040000057
Determining a local tangent plane equation of the point O according to the coordinates of the target point O to obtain a two-dimensional tangent plane formed by adjacent points; if O is ═ x0，y0，z0)，
Figure BDA0003488538040000053
The tangent plane equation for the crossing point O is as follows:
A(x-x0)+B(y-y0)+C(c-c0)＝0
s6-2: projecting three-dimensional points to a two-dimensional tangent plane, storing the projected points in a projection matrix for rotation transformation calculation, and projecting the projection transformation matrix
Figure BDA0003488538040000054
The calculation process of (2) is as follows:
Figure BDA0003488538040000055
wherein, TMatrixFor the translation transformation matrix:
Figure BDA0003488538040000056
in the formula x0,y0，z0The translation amount of each coordinate axis;
Rxexpressed as a rotation matrix of a degrees around the x-axis:
Figure BDA0003488538040000061
Ryexpressed as a rotation matrix of θ degrees around the y-axis:
Figure BDA0003488538040000062
combining the tangent plane equation of step S6-1 and the above equation can obtain any point P (x)i，yi，zi) At tangent plane OПProjection of (2):
Figure BDA0003488538040000063
s6-3: performing in-plane triangulation on the point cloud obtained by projection by using a Delaunay-based spatial region growing algorithm, wherein the triangulation process meets the empty circumcircle characteristic and the maximum and minimum angle criterion of the Delaunay algorithm, the algorithm obtains the connection relation of each point and forms a complete triangular mesh curved surface by selecting a sample triangular plate as an initial curved surface and continuously expanding the boundary of the curved surface, and finally determines the topological connection among the original three-dimensional points according to the connection relation of the projection point cloud, and the obtained triangular mesh is a reconstructed curved surface model;
by adopting the scheme, the computer triangulates the point cloud data by adopting a greedy projection triangulation algorithm and splices triangular meshes established on a projection plane, so that the three-dimensional curved surface reconstruction of the tunnel arch surface of the coal mine is realized.
As a preferred scheme of the present invention, in step S7, a triangular mesh-based hole patching algorithm is used to perform online patching on a polygonal hole generated by diluting part of the point cloud from the arch surface, and the method specifically includes the following steps:
s7-1: detecting boundary points on the tangent plane, if a certain point P in the point cloud data is a boundary characteristic point, the k neighborhood point of the point P is located at one side of the boundary characteristic point; if P is an interior point, the k neighborhood points of P should fall around the P point; judging the boundary feature points by measuring the distribution uniformity condition of k neighborhoods in the point cloud, and adopting the maximum angle difference as the measurement of the distribution uniformity:
definition of OiThe projection point of the P neighborhood point on the tangent plane is taken as (i ═ 1, 2.. k), and the P nearest neighbor O is taken asiForm a line segment with P
Figure BDA0003488538040000064
To be provided with
Figure BDA0003488538040000065
As a reference, respectively calculate
Figure BDA0003488538040000066
Rotate clockwise to
Figure BDA0003488538040000067
Angle of (2)
Figure BDA0003488538040000068
Form an angle sequence
Figure BDA0003488538040000069
Ordering angles to obtain a new sequence of angles
Figure BDA00034885380400000610
Defining the angular sequence difference as
Figure BDA00034885380400000611
From L ═ L1，L2，...，Lk) To find the maximum angle sequence difference LmaxAs a basis for judging the boundary feature points; setting a threshold value when LmaxIf the value is larger than the threshold value, P is a boundary characteristic point, otherwise P is an internal point;
s7-2: after the boundary characteristic points are solved, the disordered characteristic points are ordered, the nearest neighbor points of the characteristic points are selected as connection points, and iteration is repeated until a closed boundary line is formed;
s7-3: after the boundary line is obtained, online repairing is carried out on the polygonal cavity according to the size of an included angle alpha between two adjacent edges in the cavity boundary;
(1) when alpha is less than or equal to 90 degrees, only constructing a triangular patch;
(2) when alpha is more than 90 degrees and less than or equal to 135 degrees, two triangular patches are constructed, filling points Q meet PiEqual fraction of Q & lt Pi-1PiPi+1；
(3) When alpha is more than 135 degrees and less than or equal to 200 degrees, three triangular patches are constructed at the moment, filling points Q, R are filled, and P is satisfiediQ and PiTrisection of R & lt Pi-1PiPi+1；
(4) When α > 200 °, three triangular patches are constructed at this time, filling point Q, R, and Δ QP is satisfiedi+1PiAnd Δ PPiRi-1Is a regular triangle;
s7-4: for the above hole repairing method, the newly added triangular patch needs to be subjected to validity check; can be calculated by
Figure BDA0003488538040000071
To
Figure BDA0003488538040000072
Direction of (delta)1And
Figure BDA0003488538040000073
to
Figure BDA0003488538040000074
Direction of (delta)2To judge if delta2＞δ1If so, the newly constructed triangular patch is unreasonable; the judgment can be carried out by checking whether the newly constructed triangular patch intersects with the original cavity polygon, and if the newly constructed triangular patch intersects with the original cavity polygon, the generated triangular patch is unreasonable.
By adopting the scheme, the computer identifies the holes of the reconstructed triangular mesh model, repairs the existing polygonal holes on line and checks the legality of the newly added triangular patch; finally, the reconstruction precision of the three-dimensional curved surface of the arch surface of the coal mine tunnel is improved, and the reconstruction effect is better; and finally, repeating the steps of the three-dimensional curved surface reconstruction algorithm to solve the engineering problem of reconstructing the three-dimensional curved surfaces within a plurality of ranges of 1.6m divided by the arch surface of the whole tunnel.
The method is characterized in that a three-dimensional reconstruction is performed on the arch surface of the coal mine tunnel once, the installation position of a laser radar is located above a guide rail, the laser radar is driven to move on the guide rail at a constant speed through a hydraulic motor moving device to acquire point cloud data of the arch surface of the coal mine tunnel, and the point cloud data of the arch surface of the coal mine tunnel is completely acquired and then subjected to three-dimensional curved surface reconstruction.
The invention has the technical effects that: the computer can acquire point cloud data of the arch surface of the coal mine tunnel in real time, preprocess the data, extract normal vectors and splice point clouds; the method has the advantages that the Greeny PT algorithm can be used for rapidly realizing the three-dimensional curved surface reconstruction of the arch surface of the coal mine tunnel and carrying out online repair on the generated polygonal cavity; the Greedy PT algorithm enables the reconstruction process to be more stable and the precision to be higher, is simple in principle, strong in real-time performance and high in calculation speed, and can meet the practical application of engineering scenes.
Drawings
FIG. 1 is a flow chart of a method for reconstructing a three-dimensional curved surface of an arch surface of a coal mine tunnel;
FIG. 2 is a flow chart of the pre-processing of point cloud data;
FIG. 3 is a feature area diagram;
FIG. 4 is a non-characteristic region diagram;
FIG. 5 is a flow chart for completing the stitching of multiple point cloud data;
FIG. 6 is a target neighborhood point map;
FIG. 7 is a diagram of a Bowyer-Watson point cloud data triangulation process;
FIG. 8 is a schematic diagram of the recognition of P point as a boundary feature point;
FIG. 9 is a schematic diagram of the identification of point P as an interior point;
FIG. 10 is a schematic view of an angle sequence L;
FIG. 11 is a schematic view of hole repair;
FIG. 12 is a feature view of a broken void polygon;
FIG. 13 is a diagram of the intersection of a newly constructed triangular patch with a void polygon;
FIG. 14 is a schematic structural diagram of a simulated coal mine tunnel soffit reconfiguration assembly.
Detailed Description
The technical solution of the present invention is further explained with reference to the drawings and the embodiments.
A three-dimensional curved surface reconstruction method for an arch surface of a coal mine tunnel is disclosed, the flow of the method is shown in figure 1, the method adopts a 16-line laser radar collector to scan the arch surface of the coal mine tunnel in real time to obtain point cloud data, and the three-dimensional curved surface reconstruction is carried out on the collected point cloud data according to the following steps:
pretreatment: turning the laser radar by 90 degrees and inversely installing the laser radar on a guide rail with the length of 1.6m of a base of the mechanical arm; in order to obtain complete tunnel arch surface data within the range of 1.6m at one time, a plurality of data acquisition stations are arranged; driving a laser radar to move on a guide rail at a constant speed at certain time intervals through a hydraulic motor moving device to obtain complete tunnel arch surface point cloud data, and respectively setting an initial point and a moving speed of the laser radar on the guide rail; setting the time interval of point cloud data collected by each observation station in the laser radar moving process;
s1: and initializing, moving the laser radar on the guide rail of the mechanical arm base to an initial point, and establishing communication between the computer and the laser radar.
S2: the mechanical arm base guide rail drives the laser radar collector to sequentially move to the pre-designed data sampling points on the guide rail at a certain speed, the laser collector collects real-time data of the arch surface of the coal mine tunnel according to a preset time interval, and the collected point cloud data are transmitted to the computer.
S3: the computer analyzes the point cloud data, judges whether noise, outlier and invalid point exist, and performs point cloud data preprocessing by using filters such as stateicalOutlierRemoval, PassThrough, VoxelGrid and the like.
S4: and (3) the computer estimates a normal vector and simultaneously calculates the curvature by performing plane fitting on N points in the K neighborhood of each point by adopting a Principal Component Analysis (PCA) according to the result of the data preprocessing in the previous step.
S5: and completing point cloud image splicing on a plurality of point cloud images acquired by different sampling points on the guide rail according to the key point extraction, SAC-IA rough splicing based on FPFH (field programmable gate flash) feature descriptors and ICP (inductively coupled plasma) fine splicing sequence, so that the intersected areas are perfectly overlapped.
S6: the point cloud shape obtained by the tunnel arch surface of the coal mine is non-closed point cloud data, and a greedy projection triangulation algorithm is adopted to realize three-dimensional curved surface reconstruction according to the characteristics of the point cloud data.
S7: after the three-dimensional curved surface of the tunnel arch surface of the coal mine is reconstructed, a hole repairing algorithm based on a triangular mesh is adopted to perform online repairing on a polygonal hole generated by diluting partial point cloud data from the tunnel arch surface.
S8: and the mechanical arm base guide rail drives the laser radar to move to a position of 1.6m, namely the three-dimensional curved surface reconstruction of the arch surface of the coal mine tunnel within a range of 1.6m is completed, the three-dimensional curved surface reconstruction returns to the initial point, and the arch surface reconstruction within the next range of 1.6m is carried out after the guniting within the range is completed.
In the embodiment, the Linux-ubuntu18.04 operating system, the ROS as a development platform, PandarXT-16 as a laser radar, and the PCL as an open source library are used in the steps.
As can be seen from fig. 1, the initial moving point of the robot base guide rail in step S1 and step S2 is (x0, y0), the moving speed is S, the time interval of the point cloud data collected by the lidar collector is T, and the collector enables the computer and the lidar to establish communication through the IP address.
The mechanical arm base is a guide rail 1.6m long, a laser radar collector is driven by a hydraulic motor moving device to collect point cloud data at a certain speed and time intervals, the distance d between different data collection stations and the tunnel arch scanning width w are obtained, and a plurality of point cloud images can be generated within the whole range of 1.6 m.
As can be seen from fig. 1 and fig. 2, in step S3, filters such as PassThrough, VoxelGrid, stateticaloutlierremove, and the like are used to perform point cloud data preprocessing, where the point cloud preprocessing includes noise and outlier removal, data interception, and data reduction, and specifically includes:
s3-1: firstly, point cloud data of a pcd file collected by a collector is loaded as input, and a REMOVENFOROmPoint cloud function is used for removing NAN (invalid) points; because the running speed of the computer is reduced due to too much data, a statistical outlierRemoval filter is adopted to remove outliers, the neighborhood point is set to be 50, and the distance threshold is set to be 1.0;
s3-2: secondly, importing the point cloud into a PassThrough filter, mainly operating in the z direction, and setting a threshold value to be (-8, 8), thereby filtering out points outside a parameter range to obtain point cloud data in a specified range;
s3-3: the data reduction mainly uses a VoxelGrid grid method, a voxel grid threshold value is set to be 1.0 x 1.0, and the centroid value of all data points in each grid is calculated to replace all points in the grid, so that the calculation amount of a computer is greatly reduced.
As can be seen from fig. 1, in step S4, the computer creates a normal estimation object through a normals estimation function, establishes a kd-tree data structure, and sets the number of K search nearest neighbor points, and the PCA algorithm estimates the normal vector of the point cloud surface specifically includes:
s4-1: for any point p in the point cloudi(1 is more than or equal to i and less than or equal to n) inquiring k field points p of the n field pointsij(j is more than or equal to 1 and less than or equal to k), calculating piCentroid of k neighborhood points
Figure BDA0003488538040000091
Figure BDA0003488538040000092
S4-2: constructing a covariance matrix of local features of the point cloud:
Figure BDA0003488538040000101
s4-3: computer protocolEigenvalues λ of the variance matrix C0，λ1，λ2(λ0≤λ1≤λ2) And a feature vector v0，v1，v2λ and v are in one-to-one correspondence;
wherein v is0，v1，v2Orthogonal, eigenvectors v1，v2Determining a point piAn optimum tangent plane of0Orthogonal to the tangent plane, so point piNormal n toiCan be approximated by a feature vector v0Represents; lambda [ alpha ]0，λ1，λ2As is the degree of variation in the direction of the respective eigenvector.
S4-4: computer computing point piOf (d) curvature τpi：
Figure BDA0003488538040000102
S4-5: and (3) describing the normal vector features:
as shown in fig. 3, the included angle of the normal vector of the point cloud in the local area is larger, which indicates that the geometric feature of the area is more prominent; and if the angle transformation of the normal vector of fig. 4 is not large, it indicates that the region is relatively smooth. Thus, a point p of the point cloud is definediArithmetic mean of angles with normal vectors of its neighboring points:
Figure BDA0003488538040000103
in the formula, thetaijAs a point cloud piAngle to the adjacent normal vector.
According to point piExtracting characteristic points from the included angle between the normal vector of the adjacent point, selecting proper threshold value delta, and when f isiAt > δ, point piHas a large curvature, so piIs a characteristic point; when f isiAt < delta, point piHas a small curvature, so piIs a non-characteristic point.
In order to ensure the consistency of the normal vector directions of the point cloud surface and more accurately extract the feature points, the direction of the feature points needs to be adjusted to satisfy the formula:
ni*nj＜0，(i≠j)
in the formula, niIs a cloud normal vector of origin, njAnd obtaining a target point cloud normal vector.
As can be seen from fig. 1 and fig. 5, in step S5, the computer first selects the point cloud data of the tunnel arch surface of the two stations adjacent to each other to be spliced, and then sequentially splices the point cloud data with the point cloud of the next survey station, and the specific steps are as follows:
s5-1: extracting the key points of the point cloud image by adopting an internal morphology descriptor algorithm ISSKeypoint 3D: set point cloud P ═ Pi(xi，yi，zi) I is more than or equal to 1 and less than or equal to n, and for each point piEstablishing a local coordinate system and setting a search radius prDetermining piIs the center of a sphere, prAll points within the radius sphere, as shown in fig. 6.
Calculating the weight ωij：
Figure BDA0003488538040000111
|pi-pj|＜pr
Calculate each point piCovariance matrix of (2):
Figure BDA0003488538040000112
feature values of covariance matrix
Figure BDA0003488538040000113
Arranging in descending order; setting a threshold η1And η2Satisfy the following requirements
Figure BDA0003488538040000114
And
Figure BDA0003488538040000115
the points of (1) are key points, and iteration is carried out until all the key points are found;
s5-2: creating an FPFH feature descriptor, inputting point cloud data, calculating to obtain the descriptors by establishing a K-d tree data structure and setting the number of K search nearest points, and storing the descriptors in an FPFHSignature33 function for preliminary matching in the registration process.
S5-3: in order to solve the problem that ICP is easy to fall into local optimum, a SAC-IA algorithm is used for roughly splicing a source point cloud P and a target point cloud G. The algorithm comprises the following steps: (1) selecting m feature points to be registered from the P; (2) searching points similar to the source point cloud FPFH characteristics in the G, and selecting the point closest to the source point cloud FPFH characteristics as a corresponding relation point; (3) calculating a rotation and translation matrix of the corresponding point pair; the Huber function is used to represent the distance error sum function after the corresponding point pair is translated by rotation, and is recorded as
Figure BDA0003488538040000116
Figure BDA0003488538040000117
Wherein r is a distance threshold, | miI represents the Euclidean distance between the transformed i-th group of corresponding points, min H (m)i) And correspondingly, the optimal transformation matrix after the rough splicing is performed, so that the initial splicing is completed.
S5-4: performing ICP fine splicing by taking the source point cloud P' and the target point cloud G after rough splicing as input point clouds, and establishing an optimal node First (Best Bin First, BBF) optimized KD tree data structure by adopting an IterativeClosestPoint function to accelerate corresponding point search; all points P ' of point cloud P ' to be spliced 'iSearching the corresponding point G with the nearest distance in the target point cloud GiForming corresponding point pairs; the rotation matrix R and the translational vector T are calculated such that the root mean square error M between corresponding pairs of pointskMinimum:
Figure BDA0003488538040000118
the rotation matrices R and T found for the first iteration are:
Figure BDA0003488538040000121
T＝[-4.44089e-16 6.66134e-16 -2.98023e-08]
finally, a threshold value alpha (i.e. M) is setk-Mk+1< alpha) and maximum number of iterations NmaxTo complete the stitching of multiple point cloud images so that the overlapping regions between them are perfectly registered.
As can be seen from fig. 1, in step S6, since the point cloud shape obtained from the coal mine tunnel arch surface is non-closed point cloud data, a greedy projection triangulation algorithm is adopted according to the characteristics of the point cloud data to reconstruct a three-dimensional curved surface of the coal mine tunnel arch surface, and the specific steps are as follows:
s6-1: adopting the recursive algorithm idea of dynamic programming problem to make three-dimensional points pass through the normal vector obtained in step 4
Figure BDA0003488538040000122
The coordinates of the positive target point O determine the local tangent plane equation of O to solve the two-dimensional tangent plane formed by the adjacent points. If O is ═ x0，y0，z0)，
Figure BDA0003488538040000123
The tangent plane equation for the crossing point O is as follows:
A(x-x0)+B(y-y0)+C(c-c0)＝0
s6-2: projecting three-dimensional points to a two-dimensional tangent plane, storing the projected points in a projection matrix for rotation transformation calculation, and projecting the projection transformation matrix
Figure BDA0003488538040000124
The calculation process of (2) is as follows:
Figure BDA0003488538040000125
wherein, TMatrixFor the translation transformation matrix:
Figure BDA0003488538040000126
in the formula x0，y0，z0The amount of translation of each coordinate axis.
RxExpressed as a rotation matrix of a degrees around the x-axis:
Figure BDA0003488538040000127
Ryexpressed as a rotation matrix of θ degrees around the y-axis:
Figure BDA0003488538040000128
combining the tangent plane equation of step S6-1 and the above equation can obtain any point P (x)i，yi，zi) At tangent plane OПProjection of (2):
Figure BDA0003488538040000129
s6-3: as shown in fig. 7, a point cloud obtained by projection is triangulated in a plane by adopting a Delaunay-based spatial region growing algorithm, the triangulation process of the point cloud should meet the empty circumcircle characteristic and the maximum and minimum angle criterion of the Delaunay algorithm, the algorithm obtains the connection relation of each point and forms a complete triangular mesh curved surface by selecting a sample triangular plate as an initial curved surface and continuously expanding the boundary of the curved surface, and finally determines the topological connection between each original three-dimensional point according to the connection relation of the projection point cloud, and the obtained triangular mesh is a reconstructed curved surface model;
as can be seen from fig. 1, after the reconstruction of the three-dimensional curved surface in step S6 is completed, the computer performs online repair on the polygonal hole generated after the reconstruction of the arch surface of the coal mine tunnel by using a hole repair algorithm based on a triangular mesh in step S7, and the specific steps are as follows:
s7-1: boundary point detection: the boundary feature points are determined by measuring the distribution uniformity of k neighborhoods in the point cloud, and the maximum angle difference is used as the measurement of the distribution uniformity, as shown in fig. 8 and 9.
Definition of OiThe projection point of the P neighborhood point on the tangent plane is taken as (i ═ 1, 2.. k), and the P nearest neighborhood point O is taken asiForm a line segment with P
Figure BDA0003488538040000131
To be provided with
Figure BDA0003488538040000132
As a reference, respectively calculate
Figure BDA0003488538040000133
Rotate clockwise to
Figure BDA0003488538040000134
Angle of (2)
Figure BDA0003488538040000135
Form an angle sequence
Figure BDA0003488538040000136
Obtaining by sorting the angles
Figure BDA0003488538040000137
Defining the angular sequence difference as
Figure BDA0003488538040000138
As shown in fig. 10.
From L ═ L1，L2，...，Lk) To find the maximum angle sequence difference LmaxSetting a threshold value as a basis for judging the boundary characteristic point, and taking the threshold value as LmaxIf the value is larger than the threshold value, P is a boundary characteristic point, otherwise P is an internal point;
s7-2: after the boundary characteristic points are solved, the disordered characteristic points are ordered, the nearest neighbor points of the characteristic points are selected as connection points, and iteration is repeated until a closed boundary line is formed;
s7-3: after the boundary line is obtained, the cavity is repaired on line according to the size of the included angle α between two adjacent edges in the cavity boundary, as shown in fig. 11:
(1) when alpha is less than or equal to 90 degrees, only one triangular patch is constructed, as shown in FIG. 11 (a);
(2) when alpha is more than 90 degrees and less than or equal to 135 degrees, two triangular patches are constructed, filling points Q meet PiEqual fraction of Q & lt Pi-1PiPi+1As shown in FIG. 11 (b);
(3) when alpha is more than 135 degrees and less than or equal to 200 degrees, three triangular patches are constructed at the moment, filling points Q, R are filled, and P is satisfiediQ and PiTrisection of R & lt Pi-1PiPi+1As shown in fig. 11 (c).
(4) When α > 200 °, three triangular patches are constructed at this time, filling point Q, R, and Δ QP is satisfiedi+1PiAnd Δ PPiRi-1Is a regular triangle, as shown in fig. 11 (d).
S7-4: and (3) checking the validity of the newly added triangular patch: as fig. 12 destroys the features of the hole polygon, the newly constructed triangular patch of fig. 13 intersects the hole polygon.
For the test of FIG. 12, by calculation
Figure BDA0003488538040000139
To
Figure BDA00034885380400001310
Angle delta of1And
Figure BDA00034885380400001311
to
Figure BDA00034885380400001312
Angle delta of2To judge if delta2＞δ1If so, the newly constructed triangular patch is unreasonable; in fig. 13, it is determined by checking whether the newly constructed triangular patch intersects the original cavity polygon, and if so, the generated triangular patch is not reasonable.
As can be seen from fig. 14, the area for performing one-time three-dimensional reconstruction on the arch surface of the coal mine tunnel needs to be 1.6m, so that the length of the guide rail of the mechanical arm base is designed to be 1.6m, the installation position of the laser radar is located above the guide rail, and the laser radar is driven to move on the guide rail at a constant speed by the hydraulic motor moving device to acquire point cloud data of the arch surface of the coal mine tunnel; the experiment uses 16-line explosion-proof laser radar, and due to the limit of polar line constraint of the laser radar, the vertical field of view can only acquire point cloud data within the range of (-15 degrees to +15 degrees) at a time, and cannot acquire point cloud within the range of 1.6m at a time, so a plurality of sampling points are designed on the guide rail, the point cloud data of each sampling point are subjected to point cloud splicing to acquire complete point cloud within the range of 1.6m of the arch surface of the coal mine tunnel, and the three-dimensional curved surface reconstruction of the arch surface of the coal mine tunnel is realized on the basis.
Finally, the above embodiments are only for illustrating the technical solutions of the present invention and not for limiting, although the present invention has been described in detail with reference to the preferred embodiments, it should be understood by those skilled in the art that modifications or equivalent substitutions may be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention, and all of them should be covered in the claims of the present invention.

Claims (8)

1. A three-dimensional curved surface reconstruction method for an arch surface of a coal mine tunnel is characterized in that a laser radar is adopted to scan the arch surface of the coal mine tunnel in real time to obtain point cloud data, and the acquired point cloud data is subjected to three-dimensional curved surface reconstruction according to the following steps:
pretreatment: turning the laser radar by 90 degrees and inversely installing the laser radar on a guide rail of a mechanical arm base; acquiring tunnel arch surface data, and arranging a plurality of data acquisition stations; driving a laser radar to move on a guide rail at a constant speed at certain time intervals through a hydraulic motor moving device to obtain complete tunnel arch surface point cloud data, and respectively setting an initial point and a moving speed of the laser radar on the guide rail; setting the time interval of point cloud data collected by each observation station in the laser radar moving process;
s1: initializing, moving a laser radar on a guide rail of a mechanical arm base to an initial point, and establishing communication between a computer and the laser radar;
s2: the mechanical arm base guide rail drives the laser radar collector to sequentially move to the data sampling points on the guide rail at a certain speed, the laser collector collects real-time data of the arch surface of the coal mine tunnel according to a preset time interval, and transmits the collected point cloud data to the computer;
s3: the computer analyzes the point cloud data, judges whether noise, outliers and invalid points exist in the point cloud data, and performs point cloud data preprocessing on the point cloud data by adopting a filter;
s4: the computer adopts a principal component analysis method to carry out plane fitting on N points in K neighborhoods of each point according to the result of the data preprocessing of the last step so as to estimate normal vectors and calculate curvatures;
s5: completing point cloud image splicing on a plurality of point cloud images acquired by different sampling points on a guide rail according to key point extraction, SAC-IA rough splicing based on FPFH (field programmable gate flash) feature descriptors and ICP (inductively coupled plasma) fine splicing sequence, so that intersected areas are perfectly overlapped;
s6: the point cloud shape acquired by the arch surface of the coal mine tunnel is non-closed point cloud data, and a greedy triangulation algorithm is adopted to realize three-dimensional curved surface reconstruction according to the characteristics of the point cloud data;
s7: after the reconstruction of the three-dimensional curved surface of the tunnel arch surface of the coal mine is completed, performing online repair on a polygonal cavity generated by diluting part of point cloud data from the tunnel arch surface by adopting a cavity repair algorithm based on a triangular mesh;
s8: and the guide rail of the base of the mechanical arm drives the laser radar to move, so that the three-dimensional curved surface reconstruction of the arch surface of the coal mine tunnel is completed, the three-dimensional curved surface is returned to the initial point, and the next arch surface reconstruction is carried out after the guniting is completed.
2. The method for reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel according to claim 1, wherein in the steps S1 and S2, the initial moving point of the guide rail of the mechanical arm base is (x0, y0), the moving speed is S, the time interval of the point cloud data collected by the laser radar collector is T, and the collector enables the computer and the laser radar to establish communication through the IP address;
the guide rail of the mechanical arm base drives a laser radar collector to collect point cloud data at a certain speed and time interval through a hydraulic motor moving device, so that the distance d between different data collection stations and the tunnel arch scanning width w are obtained, and a plurality of point cloud images can be generated in the whole guide rail range; the mechanical arm is arranged on the trolley of the mobile platform.
3. The method for reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel according to claim 1, wherein in the step S3, the point cloud data is preprocessed by the specific steps of:
s3-1: adopting a statistical filter to perform statistical analysis on K neighborhood points of each point in the point cloud, and calculating the average distance from the K neighborhood points to all the neighborhood points; then calculating the average value mu and the standard deviation sigma of the average distance of each point to determine a distance threshold thresh _ d; removing the point cloud with the average neighborhood distance of a certain point lower than or higher than the threshold range of the certain point according to the point cloud density distribution;
assume an original point cloud dataset P ═ { P ═ Pi(xi，yi，zi) I is more than or equal to 1 and less than or equal to n, the distance threshold is:
thresh_d＝μ+m·σ
the point cloud data set after the outlier filtering is as follows:
P′＝(u-m*σ，u+m*σ)
wherein m is a multiple of standard deviation;
s3-2: importing P' into a straight-through filter, setting a plurality of dimension directions and different point cloud thresholds, and filtering points outside a parameter range to obtain point cloud data in a specified range;
s3-3: and importing the point cloud data obtained in the step S3-2 into a voxel filter or a conditional filter, and setting a plurality of voxel grid thresholds with different specifications and other conditions to simplify the point cloud data to obtain the point cloud preprocessing result, thereby greatly reducing the calculation amount of a computer.
4. The method for reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel according to claim 1, wherein: in step S4, the PCA algorithm estimates the point cloud normal by the specific steps of:
s4-1: for any point p in the point cloudi(1 is more than or equal to i and less than or equal to n) inquiring k field points p of the n field pointsij(j is more than or equal to 1 and less than or equal to k), calculating piCentroid of k neighborhood points
Figure FDA0003488538030000021
Figure FDA0003488538030000022
S4-2: constructing a covariance matrix of local features of the point cloud:
Figure FDA0003488538030000023
s4-3: calculating an eigenvalue λ of the covariance matrix C0，λ1，λ2，λ0≤λ1≤λ2And a feature vector v0，v1，v2λ and v are in one-to-one correspondence;
wherein v is0，v1，v2Orthogonal, eigenvectors v1，v2Determining a point piAn optimum tangent plane of0Orthogonal to the tangent plane, so point piNormal n toiCan be approximated by a feature vector v0Represents; lambda [ alpha ]0，λ1，λ2As its degree of variation in the direction of the respective eigenvector;
s4-4: calculating the point piOf (d) curvature τpi：
Figure FDA0003488538030000031
5. The method for reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel according to claim 1, wherein in the step S5, the computer completes the point cloud image mosaic of a plurality of point cloud images acquired by different sampling points on the guide rail according to the sequence of key point extraction, SAC-IA rough mosaic based on FPFH feature descriptors and ICP fine mosaic, and the method comprises the following specific steps:
s5-1: extracting key points of the point cloud image by adopting an internal morphology descriptor algorithm:
set point cloud P ═ Pi(xi，yi，zi) I is more than or equal to 1 and less than or equal to n, and for each point piEstablishing a local coordinate system and setting a search radius prDetermining piIs the center of a sphere, prAll points within the radius sphere;
calculating the weight ωij；
Figure FDA0003488538030000032
Calculate each point piCovariance matrix of (2):
Figure FDA0003488538030000033
feature values of covariance matrix
Figure FDA0003488538030000034
Arranging in descending order; setting a threshold η1And η2Satisfy the following requirements
Figure FDA0003488538030000035
And
Figure FDA0003488538030000036
the points of (1) are key points, and iteration is carried out until all the key points are found;
s5-2: calculating the fast point feature histogram features of the point cloud to be spliced and the target point cloud to obtain each calculation point MpEstablishing a simple point feature histogram according to the relative relation between the simple point feature histogram and all the neighborhood points; calculating the FPFH signature from the SPFH signature, denoted as F (M)P)：
Figure FDA0003488538030000037
Wherein d isiIs the Euclidean distance of the corresponding point pair;
s5-3: splicing the model point cloud and the target point cloud by using a sampling consistency initial registration algorithm; the algorithm comprises the following steps: (1) selecting m feature points to be registered from the model point cloud; (2) searching points similar to the FPFH characteristics of the model point cloud in the target point cloud, and selecting the point closest to the target point cloud as a corresponding relation point; (3) calculating a rotation and translation matrix of the corresponding point pair; the Huber function is used to represent the distance error sum function after the corresponding point pair is translated by rotation, and is recorded as
Figure FDA0003488538030000038
Figure FDA0003488538030000039
Wherein r is a distance threshold, | miI represents the Euclidean distance between the corresponding points of the ith group after transformation, min H (m)i) The corresponding is the optimal transformation matrix after the rough splicing;
s5-4: performing ICP fine splicing by taking the source point cloud P' and the target point cloud G after the rough splicing as input point clouds, and accelerating corresponding point pair search by adopting an optimal node to preferentially optimize a KD tree; all points P ' of point cloud P ' to be registered 'iSearching the corresponding point G with the nearest distance in the target point cloud GiForming corresponding point pairs; the rotation matrix R and the translational vector T are calculated such that the root mean square error M between corresponding pairs of pointskMinimum:
Figure FDA0003488538030000041
finally setting a threshold value alpha and the maximum iteration number Nmax，Mk-Mk+1＜α。
6. The method for reconstructing the three-dimensional curved surface of the coal mine tunnel arch surface according to claim 1, wherein a greedy projection triangulation algorithm is adopted in step S6 to realize the three-dimensional reconstruction of the coal mine tunnel arch surface, and the method comprises the following specific steps:
s6-1: adopting the recursive algorithm idea of dynamic programming problem to make three-dimensional points pass through the normal vector obtained in step 4
Figure FDA0003488538030000046
Determining a local tangent plane equation of the point O according to the coordinates of the target point O to obtain a two-dimensional tangent plane formed by adjacent points; if there is
Figure FDA0003488538030000042
The tangent plane equation for the crossing point O is as follows:
A(x-x0)+B(y-y0)+C(c-c0)＝0
s6-2: projecting three-dimensional points to a two-dimensional tangent plane, storing the projected points in a projection matrix for rotation transformation calculation, and projecting the projection transformation matrix
Figure FDA0003488538030000048
The calculation process of (2) is as follows:
Figure FDA0003488538030000047
wherein, TMatrixFor the translation transformation matrix:
Figure FDA0003488538030000043
in the formula x0，y0，z0The amount of translation of each coordinate axis.
RxExpressed as a rotation matrix of a degrees around the x-axis:
Figure FDA0003488538030000044
Ryexpressed as a rotation matrix of θ degrees around the y-axis:
Figure FDA0003488538030000045
combining the tangent plane equation of step S6-1 and the above equation can obtain any point P (x)i，yi，zi) At tangent plane O∏Projection of (2):
Figure FDA0003488538030000049
s6-3: the point cloud obtained by projection is triangulated in a plane by using a Delaunay-based space region growing algorithm, the trilateralization process of the point cloud meets the criteria of the characteristics of an empty circumcircle and the maximum minimum angle of the Delaunay algorithm, the algorithm obtains the connection relation of all points and forms a complete triangular mesh curved surface by selecting a sample triangular plate as an initial curved surface and continuously expanding the boundary of the curved surface, and finally, the topological connection between all original three-dimensional points is determined according to the connection relation of the projection point cloud, and the obtained triangular mesh is the reconstructed curved surface model.
7. The method for reconstructing the three-dimensional curved surface of the arch surface of the coal mine tunnel according to claim 1, wherein in the step S7, a triangular mesh-based hole repairing algorithm is adopted to perform online repairing on a polygonal hole generated by diluting part of point cloud from the arch surface, and the specific steps are as follows:
s7-1: detecting boundary points on the tangent plane, if a certain point P in the point cloud data is a boundary characteristic point, the k neighborhood point of the point P is located at one side of the boundary characteristic point; if P is an interior point, the k neighborhood points of P should fall around the P point; judging the boundary feature points by measuring the distribution uniformity condition of k neighborhoods in the point cloud, and adopting the maximum angle difference as the measurement of the distribution uniformity:
definition of OiThe projection point of the P neighborhood point on the tangent plane is taken as (i ═ 1, 2.. k), and the P nearest neighbor O is taken asiForm a line segment with P
Figure FDA0003488538030000051
To be provided with
Figure FDA0003488538030000052
As a reference, respectively calculate
Figure FDA0003488538030000053
Rotate clockwise to
Figure FDA0003488538030000054
Angle of (2)
Figure FDA0003488538030000055
Form an angle sequence
Figure FDA0003488538030000056
Ordering angles to obtain a new sequence of angles
Figure FDA0003488538030000057
Defining the angular sequence difference as
Figure FDA0003488538030000058
From L ═ L1，L2，...，Lk) To find the maximum angle sequence difference LmaxAs a basis for judging the boundary feature points; setting a threshold value when LmaxIf the value is larger than the threshold value, P is a boundary characteristic point, otherwise P is an internal point;
s7-2: after the boundary characteristic points are solved, the disordered characteristic points are ordered, the nearest neighbor points of the characteristic points are selected as connection points, and iteration is repeated until a closed boundary line is formed;
s7-3: after the boundary line is obtained, online repairing is carried out on the polygonal cavity according to the size of an included angle alpha between two adjacent edges in the cavity boundary;
(1) when alpha is less than or equal to 90 degrees, only constructing a triangular patch;
(2) when alpha is more than 90 degrees and less than or equal to 135 degrees, two triangular patches are constructed, filling points Q meet PiEqual fraction of Q & lt Pi-1PiPi+1；
(3) When alpha is more than 135 degrees and less than or equal to 200 degrees, three triangular patches are constructed at the moment, filling points Q, R are filled, and P is satisfiediQ and PiTrisection of R & lt Pi-1PiPi+1；
(4) When α > 200 °, three triangular patches are constructed at this time, filling point Q, R, and Δ QP is satisfiedi+1PiAnd Δ PPiRi-1Is a regular triangle;
s7-4: for the above hole repairing method, the newly added triangular patch needs to be subjected to validity check; can be calculated by
Figure FDA0003488538030000059
To
Figure FDA00034885380300000510
Run δ 1 and
Figure FDA00034885380300000511
to
Figure FDA00034885380300000512
Direction of (delta)2To judge if delta2＞δ1If so, the newly constructed triangular patch is unreasonable; the judgment can be carried out by checking whether the newly constructed triangular patch intersects with the original cavity polygon, and if the newly constructed triangular patch intersects with the original cavity polygon, the generated triangular patch is unreasonable.
8. The method for reconstructing the three-dimensional curved surface of the coal mine tunnel arch surface according to claim 1, wherein in an area where the coal mine tunnel arch surface is subjected to one-time three-dimensional reconstruction, the installation position of the laser radar is located above the guide rail, the laser radar is driven to move on the guide rail at a constant speed through a hydraulic motor moving device to acquire point cloud data of the coal mine tunnel arch surface, and the three-dimensional curved surface reconstruction is performed after the point cloud data of the complete coal mine tunnel arch surface is obtained.
CN202210089365.0A 2022-01-25 2022-01-25 Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface Pending CN114399603A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210089365.0A CN114399603A (en) 2022-01-25 2022-01-25 Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210089365.0A CN114399603A (en) 2022-01-25 2022-01-25 Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface

Publications (1)

Publication Number Publication Date
CN114399603A true CN114399603A (en) 2022-04-26

Family

ID=81232081

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210089365.0A Pending CN114399603A (en) 2022-01-25 2022-01-25 Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface

Country Status (1)

Country Link
CN (1) CN114399603A (en)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115690184A (en) * 2022-10-24 2023-02-03 西南交通大学 Tunnel face displacement measurement method based on three-dimensional laser scanning
CN116793245A (en) * 2023-08-24 2023-09-22 济南瑞源智能城市开发有限公司 Tunnel detection method, equipment and medium based on track robot
CN117961197A (en) * 2024-04-01 2024-05-03 贵州大学 Self-adaptive deviation rectifying method of unmanned turbine blade micropore electric machining unit

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115690184A (en) * 2022-10-24 2023-02-03 西南交通大学 Tunnel face displacement measurement method based on three-dimensional laser scanning
CN115690184B (en) * 2022-10-24 2024-02-06 西南交通大学 Tunnel face displacement measurement method based on three-dimensional laser scanning
CN116793245A (en) * 2023-08-24 2023-09-22 济南瑞源智能城市开发有限公司 Tunnel detection method, equipment and medium based on track robot
CN116793245B (en) * 2023-08-24 2023-12-01 济南瑞源智能城市开发有限公司 Tunnel detection method, equipment and medium based on track robot
CN117961197A (en) * 2024-04-01 2024-05-03 贵州大学 Self-adaptive deviation rectifying method of unmanned turbine blade micropore electric machining unit

Similar Documents

Publication Publication Date Title
CN114399603A (en) Three-dimensional curved surface reconstruction method for coal mine tunnel arch surface
CN111080627B (en) 2D +3D large airplane appearance defect detection and analysis method based on deep learning
CN107123164B (en) Three-dimensional reconstruction method and system for keeping sharp features
Wang et al. Fully automated generation of parametric BIM for MEP scenes based on terrestrial laser scanning data
CN111986115A (en) Accurate elimination method for laser point cloud noise and redundant data
CN112927360A (en) Three-dimensional modeling method and system based on fusion of tilt model and laser point cloud data
CN111508074A (en) Three-dimensional building model simplification method based on roof contour line
CN113628263A (en) Point cloud registration method based on local curvature and neighbor characteristics thereof
Das et al. 3D scan registration using the normal distributions transform with ground segmentation and point cloud clustering
Qin et al. Automated reconstruction of parametric bim for bridge based on terrestrial laser scanning data
CN113327276B (en) Mobile measurement-oriented general mass point cloud data registration method
CN115290097B (en) BIM-based real-time accurate map construction method, terminal and storage medium
Bassier et al. Comparison of Wall Reconstruction algorithms from Point Cloud Data for as-built BIM
CN113192200A (en) Method for constructing urban real scene three-dimensional model based on space-three parallel computing algorithm
CN111489416A (en) Tunnel axis fitting method and application in calculation of over-under excavation square measure
CN115222883A (en) Electric power tower reconstruction method based on foundation LiDAR point cloud
CN113593038A (en) Tunnel point cloud central line automatic extraction and triangulation network construction method
CN115222884A (en) Space object analysis and modeling optimization method based on artificial intelligence
Yuan et al. 3D point cloud recognition of substation equipment based on plane detection
CN117162098B (en) Autonomous planning system and method for robot gesture in narrow space
Zhao et al. Completing point clouds using structural constraints for large-scale points absence in 3D building reconstruction
Zhu et al. Feature line based building detection and reconstruction from oblique airborne imagery
Zhao et al. Design of 3D reconstruction system on quadrotor Fusing LiDAR and camera
Gollub et al. A partitioned approach for efficient graph–based place recognition
Lu et al. Bolt 3D point cloud segmentation and measurement based on DBSCAN clustering

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination