A method and system for locating underground hidden dangers using ground penetrating radar based on fusion positioning

Abstract

The present invention provides a method and system for locating underground hidden dangers of a ground-penetrating radar based on fusion positioning, and the method comprises the following steps: step S1, synchronously collecting and recording ground-penetrating radar data, GNSS positioning data and Lidar point cloud data; step S2, using a point cloud registration algorithm to register Lidar point cloud data frame by frame, calculating the pose of each frame of point cloud in the point cloud coordinate system, and obtaining point cloud pose data; step S3, fusing GNSS positioning data and point cloud pose data, and establishing a mutual conversion relationship; step S4, performing pose conversion on each frame of point cloud according to the fused pose, and constructing a point cloud map; step S5, searching for the corresponding point cloud according to the GNSS coordinates of the underground hidden danger, calculating the relative coordinates of the point cloud in the point cloud map, and converting the relative coordinates into global coordinates. The present invention can more accurately realize the detection and positioning of underground hidden dangers, and has strong anti-interference ability.

Description

Technical Field

[0001] The present invention relates to an underground hidden danger locating method, in particular to an underground hidden danger locating method of a ground penetrating radar based on fusion positioning, and to an underground hidden danger locating system of a ground penetrating radar using the underground hidden danger locating method of a ground penetrating radar based on fusion positioning. Background Art

[0002] Ground penetrating radar is a non-destructive electromagnetic wave detection technology. Ground penetrating radar transmits high-frequency electromagnetic waves underground through antennas and receives echoes from underground targets to detect the material characteristics and distribution patterns of underground media. Ground penetrating radar has the advantages of fast detection speed, high detection efficiency, high resolution and easy operation, and is widely used in the field of underground hidden danger detection.

[0003] Ground penetrating radar needs to be used in conjunction with the global satellite navigation system (GNSS). While scanning underground targets, ground penetrating radar collects and records GNSS data to obtain the location information of underground targets. Commonly used GNSS systems include GPS global positioning system and Beidou satellite navigation system, and are often combined with RTK real-time differential technology and INS inertial navigation technology to improve the accuracy of satellite positioning.

[0004] However, in urban road environments, due to the presence of many buildings, viaducts and trees, satellite signals are blocked or reflected. In many cases, even with RTK and INS, accurate positioning information cannot be obtained, and the positioning error can reach several meters or even tens of meters. Since ground penetrating radar only scans underground targets and lacks environmental reference objects above the surface, GNSS positioning deviation makes it difficult to accurately determine the specific location of hidden dangers during underground hidden danger re-survey, which greatly troubles the underground hidden danger re-survey work. Summary of the invention

[0005] The technical problem to be solved by the present invention is to provide a ground penetrating radar underground hidden danger positioning method based on fusion positioning that can more accurately locate underground hidden dangers, and further provide a ground penetrating radar underground hidden danger positioning system that adopts the ground penetrating radar underground hidden danger positioning method based on fusion positioning.

[0006] In this regard, the present invention provides a method for locating underground hidden dangers using a ground penetrating radar based on fusion positioning, comprising the following steps:

[0007] Step S1, synchronously collecting and recording ground penetrating radar data, GNSS positioning data of the global satellite navigation system, and Lidar point cloud data of the laser radar, and using a timestamp generated by a unified clock source when recording the GNSS positioning data and the Lidar point cloud data;

[0008] Step S2, using a point cloud registration algorithm to register the Lidar point cloud data frame by frame, calculating the pose of each frame of the point cloud in the point cloud coordinate system, and obtaining point cloud pose data;

[0009] Step S3, fusing the GNSS positioning data and the point cloud pose data to establish a mutual conversion relationship;

[0010] Step S4, performing posture conversion on each frame of point cloud according to the fused posture to construct a point cloud map;

[0011] Step S5, for the underground hidden dangers interpreted and output by the ground penetrating radar data, the corresponding point cloud is searched according to its GNSS coordinates, the relative coordinates of the point cloud corresponding to the underground hidden danger in the point cloud map are calculated through the point cloud registration algorithm, and the relative coordinates are converted into global coordinates to realize the positioning of the underground hidden dangers.

[0012] A further improvement of the present invention is that in step S1, the clock of the GNSS positioning data is used as a unified clock source to generate a timestamp mark, the GNSS positioning data and the Lidar point cloud data are recorded, and the timestamp of the Lidar point cloud data is obtained through the following sub-steps:

[0013] Step S101, obtaining clock data of the first GNSS positioning data, wherein the corresponding GNSS clock is recorded as T1, and the corresponding local clock is recorded as C1;

[0014] Step S102, obtaining clock data of the second GNSS positioning data, wherein the corresponding GNSS clock is recorded as T2, and the corresponding local clock is recorded as C2;

[0015] Step S103, obtaining the Lidar point cloud data of the current frame, and the corresponding local clock is recorded as C;

[0016] Step S104, by formula Record the timestamp T of the Lidar point cloud data of the current frame.

[0017] A further improvement of the present invention is that step S2 includes the following sub-steps:

[0018] Step S201, preprocessing each frame of point cloud, filtering points whose distance is less than a first threshold or greater than a second threshold, and performing voxel downsampling;

[0019] Step S202, taking the first frame point cloud as an initial reference point cloud, and taking the initial reference point cloud as the origin of the point cloud coordinate system;

[0020] Step S203, performing point cloud registration on the current frame point cloud and the reference point cloud using a normal distribution transformation algorithm or an iterative closest point algorithm to obtain the position and posture of the current frame point cloud relative to the reference point cloud;

[0021] Step S204, adding the point cloud that is successfully registered and whose movement posture is greater than the third threshold to the reference point cloud after posture conversion, and marking it as a key frame;

[0022] Step S205, discarding point clouds whose poses with the current point cloud are greater than a fourth threshold from the reference point cloud;

[0023] Step S206, repeating steps S203 to S205 until the last frame.

[0024] A further improvement of the present invention is that step S3 includes the following sub-steps:

[0025] Step S301, taking the UTM coordinates of the first GNSS positioning data as the origin of the GNSS coordinate system, and calculating the offset of the UTM coordinates of each GNSS positioning data relative to the UTM coordinates of the origin;

[0026] Step S302, sampling the point cloud coordinate trajectory so that the point cloud coordinate number is consistent with the GNSS coordinate number, and obtaining the point cloud coordinate trajectory and its coordinates corresponding to the same timestamp for each GNSS positioning data;

[0027] Step S303, aligning the point cloud coordinate trajectory and the GNSS coordinate trajectory using Proctor analysis to obtain a transformation matrix from the point cloud coordinate system to the GNSS coordinate system, and establishing a mutual conversion relationship between the GNSS positioning data and the point cloud pose data.

[0028] A further improvement of the present invention is that in step S302, the process of sampling the point cloud coordinate trajectory is as follows: for each coordinate point in the GNSS coordinate trajectory, two adjacent coordinate points are found from the point cloud coordinate trajectory, and the two coordinate points are respectively recorded as coordinates P1 (X1, Y1, Z1) and coordinates P2 (X2, Y2, Z2), and the time is t1 and t2 respectively. i By (X i ,Y i ,Z i ), i is a natural number representing the cyclic sequence, and the formula (X t ,Y t ,Z t )=(X1+S x ×Δt,Y1+S y ×Δt,Z1+S z ×Δt) calculate the point cloud coordinates P at time t t (X t ,Y t ,Z t ), where Δt = t-t1, which indicates the time from t1 to t; the speed at t

[0029] A further improvement of the present invention is that in step S303, the process of aligning the point cloud coordinate trajectory and the GNSS coordinate trajectory using Prospectus analysis is as follows: by G = {G1, G2, ... G n} represents the GNSS coordinate trajectory, P = {P1, P2, ... P n} represents the point cloud coordinate trajectory, and the transformation parameters are represented by R as the rotation matrix and T′ as the translation matrix. Perform Protts analysis, where D represents the Protts distance between the point cloud coordinate track and the GNSS coordinate track, and n represents the number of coordinate points.

[0030] A further improvement of the present invention is that step S4 includes the following sub-steps:

[0031] Step S401, for each key frame point cloud, using the transformation matrix between the point cloud coordinate system and the GNSS coordinate system, converting its point cloud coordinate system pose into the GNSS coordinate system pose;

[0032] Step S402, for each key frame point cloud, adding it to the point cloud map after being converted into a GNSS coordinate system pose;

[0033] Step S403, performing voxel down-sampling on the point cloud map to form a final GNSS coordinate system point cloud map.

[0034] A further improvement of the present invention is that step S5 includes the following sub-steps:

[0035] Step S501, searching and calculating the corresponding timestamp from the GNSS positioning data according to the GNSS coordinates of the underground hidden danger;

[0036] Step S502, searching the Lidar point cloud data for a point cloud with a timestamp closest to the timestamp according to the timestamp;

[0037] Step S503, performing point cloud registration on the point cloud corresponding to the underground hidden danger and the GNSS coordinate system point cloud map by using a normal distribution transformation algorithm or an iterative closest point algorithm, and obtaining the position and posture of the underground hidden danger in the GNSS coordinate system point cloud map;

[0038] Step S504, adding the relative coordinates of the underground hidden danger in the GNSS coordinate system point cloud map to the UTM coordinates of the origin of the GNSS coordinate system to obtain the global UTM coordinates of the hidden danger, and converting them into longitude and latitude coordinates.

[0039] A further improvement of the present invention is that in step S501, the process of searching and calculating the corresponding timestamp from the GNSS positioning data is as follows: searching for the two GNSS coordinates closest to the underground hidden danger, the timestamps of these two GNSS coordinates are T1′ and T2′, and the distances to the underground hidden danger are D1 and D2, respectively, using the formula Calculate the timestamp T′ corresponding to the underground hidden danger.

[0040] The present invention also provides a ground penetrating radar underground hidden danger positioning system based on fusion positioning, which adopts the ground penetrating radar underground hidden danger positioning method based on fusion positioning as described above, and includes a ground penetrating radar, a global satellite navigation system, a laser radar and a collection host, wherein the ground penetrating radar, the global satellite navigation system and the laser radar are respectively directly connected to the collection host or connected through a network; the ground penetrating radar is mounted on the tail of a ground penetrating vehicle, and the global satellite navigation system and the laser radar are installed on the top of the ground penetrating vehicle.

[0041] Compared with the prior art, the beneficial effects of the present invention are: for the special use environment of underground hidden dangers, it integrates ground penetrating radar data, GNSS positioning data of the global satellite navigation system and Lidar point cloud data of the laser radar to form a combination of underground hidden danger data and surface environment data, which can more accurately realize the detection and positioning of underground hidden dangers, greatly improve the positioning reliability and data validity of underground hidden danger data, and has strong anti-interference ability. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] Figure 1 It is a schematic diagram of the working process of an embodiment of the present invention;

[0043] Figure 2 is a schematic diagram of positioning trajectory fusion according to an embodiment of the present invention;

[0044] Figure 3 It is a structural schematic diagram of another embodiment of the present invention. DETAILED DESCRIPTION

[0045] The preferred embodiments of the present invention will be further described in detail below in conjunction with the accompanying drawings.

[0046] like Figure 1 As shown, this example provides a method for locating underground hidden dangers using a ground penetrating radar based on fusion positioning, comprising the following steps:

[0047] Step S1, synchronously collecting and recording ground penetrating radar data, GNSS positioning data of the global satellite navigation system, and Lidar point cloud data of the laser radar, and using a timestamp generated by a unified clock source when recording the GNSS positioning data and the Lidar point cloud data; wherein GNSS is the abbreviation of the global satellite navigation system 2, and Lidar is the abbreviation of the laser radar 3;

[0048] Step S2, using the point cloud registration algorithm to register the Lidar point cloud data frame by frame, calculate the pose of each frame of the point cloud in the point cloud coordinate system, and obtain the point cloud pose data; the pose in this example represents position + pose (direction), and the pose data of the point cloud is also required when synthesizing the point cloud map. Therefore, this example uses pose data, that is, position and pose data;

[0049] Step S3, fusing the GNSS positioning data and the point cloud pose data to establish a mutual conversion relationship;

[0050] Step S4, performing posture conversion on each frame of point cloud according to the fused posture to construct a point cloud map;

[0051] Step S5, for the underground hidden dangers interpreted and output by the ground penetrating radar data, the corresponding point cloud is searched according to its GNSS coordinates, the relative coordinates of the point cloud corresponding to the underground hidden danger in the point cloud map are calculated through the point cloud registration algorithm, and the relative coordinates are converted into global coordinates to realize the positioning of the underground hidden dangers.

[0052] In step S1 of this example, it is preferred to use the GNSS clock as a unified clock source to generate a timestamp mark to record the GNSS positioning data and the Lidar point cloud data. Since the frame rates of the GNSS clock data and the Lidar point cloud data are different, this example uses the following sub-steps to obtain the timestamp of the Lidar point cloud data shown:

[0053] Step S101, obtaining clock data of the first GNSS positioning data, wherein the corresponding GNSS clock is recorded as T1, and the corresponding local clock is recorded as C1;

[0054] Step S102, obtaining clock data of the second GNSS positioning data, wherein the corresponding GNSS clock is recorded as T2, and the corresponding local clock is recorded as C2;

[0055] Step S103, obtaining the Lidar point cloud data of the current frame, and the corresponding local clock is recorded as C;

[0056] Step S104, by formula Record the timestamp T of the Lidar point cloud data of the current frame.

[0057] In step S2 described in this example, the point cloud registration algorithm is used to register the point cloud data frame by frame, and the pose of each frame of point cloud in the point cloud coordinate system is calculated; in practical applications, it is preferably calculated through sub-steps S201 to S206.

[0058] In step S201 described in this example, each frame of point cloud is preprocessed, points that are too close and too far are filtered out, points with a distance less than the first threshold or greater than the second threshold are filtered, and voxel downsampling is performed; the first threshold is a preprocessing low threshold set by the user, and the first threshold in this example is preferably 1 meter; the second threshold is a preprocessing high threshold set by the user, and the second threshold in this example is preferably 200 meters; during the voxel downsampling process, the voxel size is preferably 1×1×1 meter; of course, in practical applications, the first threshold, the second threshold, and the voxel size can be customized and adjusted according to the actual situation and needs.

[0059] The distance D of the point is calculated using the formula When the distance D of the point satisfies the condition min_d < D < max_d, this point is retained, otherwise it is removed: min_d is the first threshold, max_d is the second threshold, and X, Y, and Z represent the point cloud coordinates; voxel downsampling means dividing the three-dimensional point cloud space into small grids (i.e., voxels) with fixed length, width, and height, each point in the point cloud is included in one of the cells (voxels), and then each voxel is sampled, for example, using the centroid of the voxel as the sampling point to replace all points in the entire voxel, that is, voxel downsampling is achieved.

[0060] The length, width, and height of the voxel are respectively represented as L, W, and H, and the voxel coordinates of the point are calculated using the following method: The centroid sampling point of the voxel is calculated using the following method:

[0061] In step S202 described in this example, the first frame of point cloud is used as the initial reference point cloud, and the initial reference point cloud is used as the origin (0,0,0)(0,0,0) of the point cloud coordinate system.

[0062] In step S203 described in this example, the current frame of point cloud and the reference point cloud are registered using the Normal Distribution Transform (NDT) algorithm or the Iterative Closest Point (ICP) algorithm. Preferably, the voxel size is 1×1×1 meter, the convergence value is selected as 0.01, the step size is selected as 0.1, and the maximum number of iterations is selected as 32 to obtain the pose of the current frame of point cloud relative to the reference point cloud.

[0063] Among them, the NDT algorithm first divides the scanned space into cells (voxels) and divides the sampled point cloud into different cells (voxels). Assuming that the point cloud surface conforms to the normal distribution, the normal distribution PDF (probability density function) parameters of each square are calculated, and the PDF parameters are used as an approximate expression of the local point cloud surface; the direction and smoothness of the local surface can be reflected by the eigenvalues ​​and eigenvectors of the covariance matrix.

[0064] For each cell (voxel), the sample set is {V i =[V i1 ,V i2 ,V im ] T |1≤i≤n} represents the point, n is the number of samples, and the PDF parameters of the cell (voxel) are calculated as follows, including the mean and the covariance matrix C, the mean Covariance matrix

[0065] When using NDT registration, the algorithm goal is to find the posture of the current scan frame to maximize the possibility that it is located on the reference point cloud surface. The parameters that need to be optimized are the transformation of the current scan point cloud, including rotation and translation. The NDT algorithm uses Newton iteration to optimize parameters and find the maximum likelihood posture transformation to perform point cloud positioning.

[0066] With P = {P1, P2, ... P n} represents the set of points to be registered in the current scan frame. The transformation parameters are represented by R as the rotation matrix and T′ as the translation matrix. For each point P in the current scan frame, i , through the formula X i =R×P i +T′ to transform the posture: point X after the posture transformation i , through the formula Calculate the probability distribution function f(X i ).

[0067] The NDT algorithm makes the objective function to calculate, θ represents the likelihood function; since the NDT algorithm does not use point features for calculation and matching during the registration process, it has a fast calculation speed, stable time consumption, and has a certain tolerance for small differences.

[0068] This example can also be implemented using the ICP algorithm, which is essentially an optimal registration method based on the least squares method. The ICP algorithm repeatedly selects the corresponding closest point pair and calculates the optimal rigid body transformation until the convergence condition of the registration is met. Use R′={R′1,R′2,…R′ n} represents the reference point cloud, P = {P1, P2, ... P n} represents the set of points to be registered, and the ICP algorithm makes the following objective function T′ represents the translation matrix.

[0069] The ICP algorithm first searches for the nearest corresponding point in the reference point cloud for each point in the point cloud to be registered, forming a pair of matching points. Then the sum of the Euclidean distances of all matching point pairs is used as the error objective function, and SVD (Singular Value Decomposition) is used to calculate the rotation matrix R and translation matrix T' to minimize the error. After the point cloud to be registered is transformed according to the rotation matrix R and translation matrix T', the matching point pairs are searched again and the new rotation matrix R and translation matrix T' are solved. This process is repeated until the convergence condition is reached or the maximum number of iterations is reached.

[0070] The ICP algorithm has the advantage of high registration accuracy, but since each iteration requires re-finding matching point pairs, its operating efficiency is lower than that of the NDT algorithm. Therefore, in this example, the NDT algorithm and the ICP algorithm are preferably combined, and the two algorithms are automatically switched according to actual application requirements. For example, the NDT algorithm is used by default. When encountering special circumstances or extreme usage environments, if the registration in step S204 is unsuccessful, the ICP algorithm is returned to be re-registered.

[0071] In step S204 of this example, the point cloud that is successfully registered (convergence flag is 1, and registration score is less than 50, indicating successful registration) and whose mobile pose is greater than the third threshold is discarded, so as to limit the size of the reference point cloud to reduce memory consumption and improve registration efficiency; the reference point cloud is added after pose conversion and marked as a key frame. The third threshold is a preset mobile pose threshold, which can be customized and modified according to actual needs. The third threshold is preferably 5 meters by default.

[0072] In step S205 of this example, point clouds whose pose difference with the current point cloud is greater than a fourth threshold are discarded from the reference point cloud. The fourth threshold is a pre-set point cloud pose threshold, which can be customized and modified according to actual needs. It is used to discard point clouds whose pose difference with the current point cloud is too large, and to limit the size of the reference point cloud to reduce memory consumption and improve alignment efficiency.

[0073] Step S206, repeating steps S203 to S205 until the last frame.

[0074] In this example, step S3 preferably includes sub-steps S301 to S303, wherein step S301 uses the UTM coordinates of the first GNSS positioning data as the origin of the GNSS coordinate system, and calculates the offset of the UTM coordinates of each GNSS positioning data relative to the origin UTM coordinates, which is its coordinate in the GNSS coordinate system. The UTM coordinates are Universal Transverse Mercator Grid System, which is a commonly used plane rectangular geographic coordinate system. The origin UTM coordinates are expressed as (E o ,N o ,H o ), the UTM coordinates of the current GNSS are expressed as (N i ,E i ,H i ), then the coordinates in the GNSS coordinate system are: (X i ,Y i ,Z i )=(E i -E o ,N i -N o ,H i ), it is worth noting that the H i What is needed is the absolute height from the ground, not the relative height from the origin, so there is no need to compare it with H. o Perform differential processing.

[0075] In step S302 of this example, the point cloud coordinate trajectory is sampled so that the number of point cloud coordinates is consistent with the number of GNSS coordinates, and for each GNSS positioning data, the point cloud coordinate trajectory and its coordinates corresponding to the same timestamp are obtained; in this example, for each coordinate point in the GNSS coordinate trajectory (time t at timestamp T), two adjacent coordinate points are found from the point cloud coordinate trajectory, and these two coordinate points are respectively recorded as coordinates P1 (X1, Y1, Z1) and coordinates P2 (X2, Y2, Z2), and the times are t1 and t2 respectively. i By (X i ,Y i ,Z i ), i is a natural number representing the cyclic sequence, and the formula (X t ,Y t ,Z t )=(X1+S x ×Δt,Y1+S y ×Δt,Z1+S z ×Δt) calculate the point cloud coordinates P at time t t (X t ,Yt ,Z t ), where Δt = t-t1, which indicates the time from t1 to t; the speed at t

[0076] In step S303 of this example, the point cloud coordinate trajectory and the GNSS coordinate trajectory are aligned using Procrustes Analysis to obtain a transformation matrix from the point cloud coordinate system to the GNSS coordinate system, and to establish a mutual conversion relationship between the GNSS positioning data and the point cloud pose data. In step S303, the process of aligning the point cloud coordinate trajectory and the GNSS coordinate trajectory using Procrustes Analysis is as follows: by G = {G1, G2, ... G n} represents the GNSS coordinate trajectory, P = {P1, P2, ... P n} represents the point cloud coordinate trajectory, and the transformation parameters are represented by R as the rotation matrix and T′ as the translation matrix. Perform Protts analysis, where D represents the Protts distance between the point cloud coordinate track and the GNSS coordinate track, and n represents the number of coordinate points.

[0077] In this example, step S303 adopts Procrustes Analysis to find the transformation from the point cloud positioning trajectory to the GNSS positioning trajectory using the least squares method, so that the Procrustes distance between the transformed point cloud positioning trajectory and the GNSS positioning trajectory is minimized; the result of the Procrustes analysis is the transformation matrix from the point cloud coordinate system to the GNSS coordinate system. The transformation matrix can be used to realize the mutual conversion between the GNSS pose data and the point cloud pose data.

[0078] In this example, step S4 preferably includes sub-steps S401 to S403. Specifically, in step S401 of this example, for each key frame point cloud, the point cloud coordinate system pose is converted into the GNSS coordinate system pose using the transformation matrix between the point cloud coordinate system and the GNSS coordinate system (referred to as the point cloud-GNSS transformation matrix); in this example, the point cloud-GNSS transformation matrix is ​​preferably represented by the rotation matrix R and the translation matrix T′, and the rotation matrix R is used to represent the point cloud-GNSS transformation matrix. p and the translation matrix T p Represents the position of the point cloud in the point cloud coordinate system. In practical applications, the complete transformation generally includes rotation, translation and scaling. In this example, no scaling is required, only rotation and translation are required. They can be combined as a 4th-order matrix, or separated into a 3rd-order matrix and a 3D vector. In this example, the separation method is preferred, which is more convenient for calculation. In this example, the formula (R g ,T g )=(R×R p ,T′+T p ) calculates its position and pose in the GNSS coordinate system (R g,T g ):

[0079] In step S402 of this example, for each key frame point cloud P, the voxel size is preferably 1×1×1 meter, and is added to the point cloud map after being converted into the GNSS coordinate system posture; for each point P in the current frame i , using formula X i =R g ×P i +T g Make the transformation.

[0080] Step S403, performing voxel down-sampling on the point cloud map to form a final GNSS coordinate system point cloud map.

[0081] In this example, step S5 preferably includes sub-steps S501 to S504. Specifically, in step S501, according to the GNSS coordinates of the underground hidden danger, the corresponding timestamp is searched and calculated from the GNSS positioning data;

[0082] In step S501 of this example, since the GNSS coordinates of the underground hidden danger must be located on the GNSS coordinate track, the two GNSS coordinates closest to the coordinates of the underground hidden danger are found. The GNSS coordinates of the underground hidden danger are located between the two GNSS coordinates, and the timestamps corresponding to the underground hidden danger are also located between the timestamps of the two GNSS coordinates. The UTM coordinates of the hidden danger are expressed as (E, N, H), and the UTM coordinates of the GNSS are expressed as (E i ,N i ,H i ), using the formula Calculate the distance D between the coordinates of the underground hidden danger and the GNSS coordinates. Then, the process of finding and calculating the corresponding timestamp from the GNSS positioning data is as follows: Find the two GNSS coordinates closest to the underground hidden danger, the timestamps of these two GNSS coordinates are T1′ and T2′, and the distances to the underground hidden danger are D1 and D2, respectively. Use the formula Calculate the timestamp T′ corresponding to the underground hidden danger.

[0083] In step S502 of this example, according to the timestamp, the point cloud with the closest timestamp is found from the Lidar point cloud data; in step S503, the point cloud corresponding to the underground hidden danger is aligned with the GNSS coordinate system point cloud map through the normal distribution transformation algorithm or the iterative closest point algorithm. The voxel size is preferably 1×1×1 meter, the convergence value is 0.01, the step size is 0.1, and the maximum number of iterations is 32 to obtain the position and posture of the underground hidden danger in the GNSS coordinate system point cloud map; the point cloud alignment algorithm includes the use of the NDT (Normal Distribution Transform) algorithm or the ICP (Iterative Closest Point) algorithm.

[0084] In step S504 of this example, the relative coordinates of the underground hidden danger in the GNSS coordinate system point cloud map are added to the UTM coordinates of the origin of the GNSS coordinate system to obtain the global UTM coordinates of the hidden danger, and then converted into longitude and latitude coordinates.

[0085] like Figure 3 As shown, this example also provides a ground-penetrating radar underground hidden danger positioning system based on fusion positioning, which adopts the ground-penetrating radar underground hidden danger positioning method based on fusion positioning as described above, and includes a ground-penetrating radar 1, a global satellite navigation system 2, a laser radar 3 and a collection host 4, wherein the ground-penetrating radar 1 includes a ground-penetrating radar host and a ground-penetrating radar antenna connected to the collection host 4 of the ground-penetrating radar, and the global satellite navigation system 2 (GNSS) includes a GPS global satellite navigation positioning system or a Beidou satellite navigation positioning system, generally including a GNSS host and a GNSS antenna, and the GNSS host 3 has an optional RTK real-time differential function and an INS inertial navigation function. The laser radar 3 (Lidar) can be a multi-line laser radar such as 16 lines, 32 lines, 64 lines and 128 lines. The ground-penetrating radar 1, the global satellite navigation system 2 and the laser radar 3 are respectively connected to the collection host 4 directly or through a network; the ground-penetrating radar 1 is mounted on the tail of a ground-penetrating vehicle 5, and the global satellite navigation system 2 and the laser radar 3 are installed on the top of the ground-penetrating vehicle 5.

[0086] In summary, if Figure 2 As shown, this example targets the special use environment of underground hidden dangers and integrates ground penetrating radar data, GNSS positioning data of the global satellite navigation system, and Lidar point cloud data to form a combination of underground hidden danger data and surface environment data. It can more accurately detect and locate underground hidden dangers, greatly improve the reliability and effectiveness of positioning underground hidden danger data, and has strong anti-interference ability.

[0087] The above contents are further detailed descriptions of the present invention in combination with specific preferred embodiments, and it cannot be determined that the specific implementation of the present invention is limited to these descriptions. For ordinary technicians in the technical field to which the present invention belongs, several simple deductions or substitutions can be made without departing from the concept of the present invention, which should be regarded as falling within the protection scope of the present invention.

Claims

1. A method for locating underground hidden dangers using ground penetrating radar based on fusion positioning, characterized in that: The following steps are involved: Step S1, synchronously collecting and recording ground penetrating radar data, GNSS positioning data of the global satellite navigation system, and Lidar point cloud data of the laser radar, and using a timestamp generated by a unified clock source when recording the GNSS positioning data and the Lidar point cloud data; Step S2, using a point cloud registration algorithm to register the Lidar point cloud data frame by frame, calculating the pose of each frame of the point cloud in the point cloud coordinate system, and obtaining point cloud pose data; Step S3, fusing the GNSS positioning data and the point cloud pose data to establish a mutual conversion relationship; Step S4, performing posture conversion on each frame of point cloud according to the fused posture to construct a point cloud map; Step S5, for the underground hidden dangers interpreted and output by the ground penetrating radar data, the corresponding point cloud is searched according to its GNSS coordinates, the relative coordinates of the point cloud corresponding to the underground hidden danger in the point cloud map are calculated by the point cloud registration algorithm, and the relative coordinates are converted into global coordinates to realize the positioning of the underground hidden dangers; In step S1, the clock of the GNSS positioning data is used as a unified clock source to generate a timestamp mark, the GNSS positioning data and the Lidar point cloud data are recorded, and the timestamp of the Lidar point cloud data is obtained through the following sub-steps: Step S101, obtaining clock data of the first GNSS positioning data, wherein the corresponding GNSS clock is recorded as T1, and the corresponding local clock is recorded as C1; Step S102, obtaining clock data of the second GNSS positioning data, wherein the corresponding GNSS clock is recorded as T2, and the corresponding local clock is recorded as C2; Step S103, obtaining the Lidar point cloud data of the current frame, and the corresponding local clock is recorded as C; Step S104, by formula Record the timestamp T of the Lidar point cloud data of the current frame; The step S2 comprises the following sub-steps: Step S201, preprocessing each frame of point cloud, filtering points whose distance is less than a first threshold or greater than a second threshold, and performing voxel downsampling; Step S202, taking the first frame point cloud as an initial reference point cloud, and taking the initial reference point cloud as the origin of the point cloud coordinate system; Step S203, performing point cloud registration on the current frame point cloud and the reference point cloud using a normal distribution transformation algorithm or an iterative closest point algorithm to obtain the position and posture of the current frame point cloud relative to the reference point cloud; Step S204, adding the point cloud that is successfully registered and whose movement posture is greater than the third threshold to the reference point cloud after posture conversion, and marking it as a key frame; Step S205, discarding point clouds whose poses with the current point cloud are greater than a fourth threshold from the reference point cloud; Step S206, repeating steps S203 to S205 until the last frame; The step S3 comprises the following sub-steps: Step S301: The UTM coordinates of the first GNSS positioning data are used as the origin of the GNSS coordinate system. The origin of the GNSS coordinate system is represented by (E o ,N o ,H o ), the UTM coordinates of the current GNSS are expressed as (E i ,N i ,H i ), H i Refers to the absolute height from the ground, through the formula (X i ,Y i ,Z i )=(E i -E o ,N i -N o ,H i ) Calculate the offset of the UTM coordinates of each GNSS positioning data relative to the UTM coordinates of the origin; Step S302, sampling the point cloud coordinate trajectory so that the point cloud coordinate number is consistent with the GNSS coordinate number, and obtaining the point cloud coordinate trajectory and its coordinates corresponding to the same timestamp for each GNSS positioning data; Step S303, aligning the point cloud coordinate trajectory and the GNSS coordinate trajectory using Proctor analysis to obtain a transformation matrix from the point cloud coordinate system to the GNSS coordinate system, and establishing a mutual conversion relationship between the GNSS positioning data and the point cloud pose data; The step S4 comprises the following sub-steps: Step S401, for each key frame point cloud, using the transformation matrix between the point cloud coordinate system and the GNSS coordinate system, converting its point cloud coordinate system pose into the GNSS coordinate system pose; Step S402, for each key frame point cloud, adding it to the point cloud map after being converted into a GNSS coordinate system pose; Step S403, performing voxel downsampling on the point cloud map to form a final GNSS coordinate system point cloud map; The step S5 comprises the following sub-steps: Step S501, searching and calculating the corresponding timestamp from the GNSS positioning data according to the GNSS coordinates of the underground hidden danger; Step S502, searching the Lidar point cloud data for a point cloud with a timestamp closest to the timestamp according to the timestamp; Step S503, performing point cloud registration on the point cloud corresponding to the underground hidden danger and the GNSS coordinate system point cloud map by using a normal distribution transformation algorithm or an iterative closest point algorithm, and obtaining the position and posture of the underground hidden danger in the GNSS coordinate system point cloud map; Step S504, adding the relative coordinates of the underground hidden danger in the GNSS coordinate system point cloud map to the UTM coordinates of the origin of the GNSS coordinate system to obtain the global UTM coordinates of the hidden danger, and converting them into longitude and latitude coordinates.

2. The method for locating underground hidden dangers by ground penetrating radar based on fusion positioning according to claim 1 is characterized in that: In step S302, the process of sampling the point cloud coordinate trajectory is as follows: for each coordinate point in the GNSS coordinate trajectory, two adjacent coordinate points are found from the point cloud coordinate trajectory, and the two coordinate points are respectively recorded as coordinates P1 (X1, Y1, Z1) and coordinates P2 (X2, Y2, Z2), and the time is t1 and t2 respectively. i By (X i ,Y i ,Z i ), i is a natural number representing the cyclic sequence, and the formula (X t ,Y t ,Z t )=(X1+S x ×Δt,Y1+S y ×Δt,Z1+S z ×Δt) calculate the point cloud coordinates P at time t t (X t ,Y t ,Z t ), where Δt = t-t1, which indicates the time from t1 to t; the speed at t 3. The method for locating underground hidden dangers by ground penetrating radar based on fusion positioning according to claim 2 is characterized in that: In step S303, the process of aligning the point cloud coordinate trajectory and the GNSS coordinate trajectory using Prospectus analysis is as follows: n } represents the GNSS coordinate trajectory, P = {P1, P2, ... P n } represents the point cloud coordinate trajectory, and the transformation parameters are represented by R as the rotation matrix and T′ as the translation matrix. Perform Protts analysis, where D represents the Protts distance between the point cloud coordinate track and the GNSS coordinate track, and n represents the number of coordinate points.

4. The method for locating underground hidden dangers by ground penetrating radar based on fusion positioning according to claim 1 is characterized in that: In step S501, the process of searching and calculating the corresponding timestamp from the GNSS positioning data is as follows: search for the two GNSS coordinates closest to the underground hidden danger, the timestamps of these two GNSS coordinates are T1′ and T2′, and the distances to the underground hidden danger are D1 and D2, respectively, and the formula is used. Calculate the timestamp T′ corresponding to the underground hidden danger.

5. A ground penetrating radar underground hidden danger positioning system based on fusion positioning, characterized in that: A ground-penetrating radar underground hidden danger positioning method based on fusion positioning as described in any one of claims 1 to 4 is adopted, and includes a ground-penetrating radar, a global satellite navigation system, a laser radar and a collection host, wherein the ground-penetrating radar, the global satellite navigation system and the laser radar are respectively connected to the collection host directly or through a network; the ground-penetrating radar is mounted on the rear of a ground-penetrating vehicle, and the global satellite navigation system and the laser radar are installed on the top of the ground-penetrating vehicle.

Images