A method for identifying, locating and imaging underground targets

Abstract

The application discloses a kind of underground target identification positioning and imaging method, utilize original B scanning echo data, adopt mean method to remove surface reflection wave and antenna coupling wave;Image processing method is used to convert echo image into binary image;According to the pulse width of transmission waveform, column segment set is constructed on binary image, and according to the width of ricker wave main lobe, column segment set is denoised, and denoised column segment set is obtained;From column segment set, the column segment cluster containing different target corresponding hyperbolic curve label is clustered;Hyperbolic curve fitting is carried out to different target hyperbolic curve column segment cluster using restricted fitting algorithm;According to hyperbolic curve parameter, target position, radius and underground dielectric constant are inverted, and different target characteristic vector is clustered according to similarity;According to binary image and underground dielectric constant, post-projection imaging is carried out, and false image is removed to obtain underground pipeline image.The application combines hyperbolic curve trend characteristics, introduces the concept of column segment trend, and improves the efficiency of clustering target hyperbolic curve.

Description

TECHNICAL FIELD

[0001] The present application belongs to the field of underground pipeline detection, and relates to a method for identifying and positioning underground targets and imaging, in particular, a method for identifying and positioning and imaging multiple underground targets in the case of multiple target echo information interweaving. BACKGROUND

[0002] With the rapid economic development, the modernization level of cities is increasingly high, and the construction of underground pipe network is increasingly important. At present, the underground pipe network includes natural gas pipelines, sewers, cables, water supply pipelines, communication cables, heating pipelines, oil pipelines, etc. The structure of underground pipe network is increasingly complex, and the types are also increasingly diverse. However, the relative lag in the management of underground pipe network, the continuous change of urban environment and focus lead to the need for continuous maintenance and updating of underground pipe network. Underground pipelines have become a hot issue in urban research and a key content of urban scientific management and urban ecological protection.

[0003] In order to accurately detect underground pipelines, thereby providing effective protection for various constructions in cities and the healthy and safe operation of cities, some trenchless technologies are gradually being researched and applied at present. Ground penetrating radar has become one of the mainstream underground pipeline detection methods due to its fast detection speed, high precision and no damage to the ground. Ground penetrating radar transmits electromagnetic waves and receives reflected waves of underground pipelines, and through processing and analyzing waveform, amplitude and time data of the echoes, the position, shape, structure, dielectric constant and burial depth of the underground target are inferred. The most widely used algorithm at present is to generate a B-Scan image (a two-dimensional profile graph under the ground) by moving the receiving antenna, and through a related algorithm, the point set of the hyperbolic curve in the echo image is extracted and clustered and fitted, and finally the position and radius of the underground pipeline are inverted.

[0004] “Automatic GPR B-scan image interpretation model” published by Zhou Xirén et al. in Earth Science and Remote Sensing (2018, 56(99): 3398-3412.) designs a hyperbolic curve opening scanning clustering algorithm based on a line segment set, aiming to quickly and accurately locate the position of the hyperbolic curve from the binary image, and improves the clustering speed and accuracy. However, this method depends on the recognition rate of the hyperbolic curve opening, and cannot adapt to the environment of multiple target echo signals interweaving; at the same time, the preprocessing based on the opening and closing operation will lose a large amount of echo information, resulting in a large deviation of the final fitting method, which cannot adapt to the increasingly complex underground pipe network environment. Therefore, it is very important for future pipe network processing and urban construction to study how to accurately and efficiently extract target information from multiple target echo information and successfully invert underground pipeline information. SUMMARY

[0005] In view of the prior art, the present application aims to provide a method for identifying, positioning and imaging a multi-target underground, which can accurately and quickly position the target and the dielectric constant of the underground medium in the case of multi-target echo information interleaving, and finally image the target.

[0006] To solve the above technical problems, the present application provides a method for identifying, positioning and imaging a target underground, which comprises the following steps:

[0007] Step 1: using the original B-scan echo data, the mean method is used to remove the ground reflection and antenna coupling wave, and the reflection information of the underground target is highlighted;

[0008] Step 2: using an image processing method, the echo image obtained in step 1 is converted into a binary image;

[0009] Step 3: according to the pulse width of the transmitted waveform, a column segment set is constructed from the binary image, and the column segment set is denoised according to the main lobe width of the ricker wave, to obtain a denoised column segment set;

[0010] Step 4: according to the trend characteristics of the target hyperbolic curve, the column segment cluster containing the hyperbolic curve label corresponding to different targets is clustered from the column segment set obtained in step 3;

[0011] Step 5: according to the characteristics of the echo hyperbolic curve of the ground coupling ground penetrating radar detecting the underground target, a restricted fitting algorithm is used to fit the hyperbolic curve of the hyperbolic curve column segment cluster of different targets obtained in step 4;

[0012] Step 6: according to the hyperbolic curve parameters obtained in step 5, the position, radius and underground dielectric constant of the target are inversed, and further clustering is performed according to the similarity of the feature vectors of different targets;

[0013] Step 7: according to the binary image and the underground dielectric constant obtained in steps 2 and 6, post-projection imaging is performed, and artifact removal processing is performed to obtain an underground target image.

[0014] Further, the image processing method for converting the ground penetrating radar B-scan image into a binary image is as follows:

[0015] The edge of the target hyperbolic curve in the image is extracted using the canny algorithm, the sidelobe of the echo waveform ricker wave is taken as a standard, the mean value of the echo energy greater than the standard on the echo edge is calculated, and the mean value is taken as a threshold for binarization, and the echo image is binarized:

[0016]

[0017] In the formula, S(i,j) represents the amplitude of the i-th row and j-th column in the echo image, R represents the side lobe peak value of the Ricker wave, Edge represents the edge point matrix extracted by the Canny algorithm, num represents the number of sampling points meeting the condition, sum(.) represents the summation operation, and Th represents the threshold value of binaryzation.

[0018] Further, the denoised column segment set is specifically:

[0019] Col seg = {Col, Row up , Row down , Row mid , Trend, Sort}, Col down - Col up ≥ N

[0020] In the formula, Col Seg represents a column segment in the column segment set, Col represents the column number where the column segment is located, Row up represents the upper boundary of the row number where the column segment is located, Row down represents the lower boundary of the row number where the column segment is located, Roe mid represents the row number where the center of the column segment is located, Trend represents the hyperbolic trend of adjacent column segments, the initial state is 0, Sort represents the clustering label of the hyperbolic curve of the column segment, and N represents the threshold value of the received waveform width.

[0021] Further, the clustering of column segment clusters of different hyperbolic curves from the column segment set obtained in step 3 according to the trend characteristics of the target hyperbolic curve in step 4 is specifically:

[0022] Starting from left to right to scan the column segments in the column segment set, if the current column segment is in the initial state, a new hyperbolic clustering label is created, and it is searched whether there is a column segment meeting the trend change condition in the next column of the column segment in the column segment set, the trend change condition including: rising to rising, falling to falling, rising to keeping and rising to falling, if there is a column segment meeting the condition, the clustering label of the current column segment is assigned to the next column segment, and Trend is updated according to the change trend; if not, the current clustering label expansion is terminated, and the clustering label of the hyperbolic curve is re-established, and the number of column segments of the current clustering label is judged, if it is lower than the threshold value Num, the current clustering label is deleted.

[0023] If the next column segment of the current column segment has both rising and falling trend column segments, a new clustering label is established, and the old label is copied to the clustering set of the new label, forming two clustering label sets.

[0024] Further, the restricted fitting algorithm is used to fit the hyperbolic curves of the different target hyperbolic column segment clusters obtained in step 4 in step 5, which is specifically:

[0025] The hyperbola is:

[0026] Ay 2 +Bt+Cx 2 +Dx+E=0

[0027] According to the different target column segment set obtained in step 3, Row mid The hyperbolic parameter fitting is performed with Col;

[0028] The target hyperbola of the ground coupling mode of the ground penetrating radar detection satisfies:

[0029]

[0030] In the formula, (x, y) represents the row and column numbers of the point on the hyperbola, (x0, y0) represents the row and column numbers of the point where the target is located, dx represents the ground penetrating radar detection interval, dt represents the radar sampling point interval time, c represents the light speed, ε represents the dielectric constant of the underground medium, and r represents the radius of the underground pipeline. The fitting parameters are limited by the following conditions:

[0031]

[0032] In the formula, min y represents the minimum y value in the clustering data set, and a limited fitting algorithm is used.

[0033] Further, the position, radius and dielectric constant of the target in step 6 according to the hyperbolic parameters obtained in step 5 are specifically:

[0034]

[0035] In the formula, A, B, D, E are hyperbolic fitting parameters, T x is the target abscissa, T y is the target ordinate, i.e. the target depth, ∈ represents the dielectric constant of the underground medium, dt represents the radar sampling point interval time, c represents the light speed, and r represents the radius of the underground pipeline.

[0036] Further, the similarity of the feature vectors of different targets in step 6 is specifically:

[0037]

[0038] In the formula, (T x1 , T y1 ) represents the abscissa and ordinate of the first target, (T x2 , T y2 ) represents the abscissa and ordinate of the second target, r1 and r2 represent the radii of the detected targets, and Similarity less than a given threshold is classified as the same target.

[0039] Further, the binary image obtained in step 2 and step 6 is post-projected according to the underground dielectric constant in step 7, and the artifact removal process is performed to obtain the underground target image, specifically:

[0040] Step 7.1: First, the detection area is divided into a plurality of pixel points to construct an imaging domain matrix U, and the echo time delay of each pixel point from the antenna is calculated:

[0041]

[0042] In the formula, d is the two-way travel distance of radar detection point i to pixel point (x j ,y j ), c represents the speed of light, ε is the dielectric constant of the medium, and t(i,x j ,y j ) represents the propagation time of radar detection point i to pixel point (x j ,y j ).

[0043] The sample point sequence number of the echo data received at the detection point with a time delay of t(i,x j ,y j ) is calculated:

[0044] M i =round(t(i,x j ,y j ) / dt)

[0045] In the formula, round(.) represents the approximation operation;

[0046] Step 7.2: Repeat step 7.1 to calculate the echo amplitude of a given pixel at all detection points, and add all echo values to obtain the assumed reflection coefficient of the pixel point:

[0047]

[0048] Where N represents the total number of detection points, that is, the column number of the echo matrix, and S(i,M i ) represents the data in the i-th column and the M i -th row of the echo data matrix;

[0049] Step 7.3: Repeat steps 7.1 and 7.2 to calculate the assumed reflection coefficient of all pixel points in the imaging domain to obtain a high-resolution two-dimensional image matrix U;

[0050] Step 7.4: Set the artifact removal threshold to eliminate artifacts in the two-dimensional image matrix U, specifically:

[0051]

[0052] In the formula, n t The number of column segments is the minimum number of cluster column segments, and U(i,j) represents the pixel energy of the image generated by the back projection algorithm at the (i,j) position.

[0053] Further, the cross section of the target is circular or annular.

[0054] The present application has the beneficial effects: the purpose of the present application is to solve the problem that the existing B-Scan image interpretation model cannot effectively cluster different echo hyperbolic curves from the environment where multiple target echo information is intertwined, resulting in missed detection and low positioning accuracy, and the problem that the existing recognition algorithm separates echo information extraction and hyperbolic curve recognition into two independent steps, resulting in slow convergence speed, and a new B-Scan image interpretation model is designed. The model can automatically process the sidelobe interference caused by ricker wave, avoid the problem of too low binary threshold, and at the same time can process the problem of poor shape caused by the destruction of hyperbolic curve opening and the problem of hyperbolic curve fracture caused by curve intersection, identify point clusters with hyperbolic curve characteristics, and obtain accurate underground target position and image. The method can accurately locate the target position and the dielectric constant of the underground medium under the condition that multiple target echo information is intertwined, and finally image the target, enhancing the adaptability to complex scenes. BRIEF DESCRIPTION OF DRAWINGS

[0055] Figure 1 It is a B-Scan image interpretation method process schematic diagram of the present application.

[0056] Figure 2 It is a column segment structure diagram.

[0057] Figure 3 It is a column segment trend structure diagram.

[0058] Figure 4 It is a structure diagram for analyzing hyperbolic curve trend clustering at the pixel level.

[0059] Figure 5 It is a structure diagram for analyzing hyperbolic curve trend clustering at the pixel level.

[0060] Figure 6 It is a structure diagram for analyzing hyperbolic curve trend clustering at the pixel level.

[0061] Figure 7 It is a hyperbolic curve formation schematic diagram under the ground coupling ground penetrating radar detection mode.

[0062] Figure 8 It is a raw echo B-Scan image.

[0063] Figure 9 It is a B-Scan image after binarization.

[0064] Figure 10 is a clustering result based on hyperbolic characteristics.

[0065] Figure 11 is a fitting hyperbolic result.

[0066] Figure 12 is a post-projection algorithm imaging result.

[0067] Figure 13 is an imaging result after removing artifacts. DETAILED DESCRIPTION

[0068] The application will be further described below in conjunction with the accompanying drawings and examples.

[0069] The purpose of the present application is to solve the problems existing in the current underground detection technology field for circular and circular target detection.

[0070] The purpose of the present application is to design a new B-Scan image interpretation method, mainly to solve the problem that the existing B-Scan image interpretation model cannot effectively cluster different echo hyperbolas from the environment where multiple target echo information is intertwined, resulting in missed detection and low positioning accuracy, and the existing recognition algorithm separates echo information extraction and hyperbola recognition into two independent steps, resulting in slow convergence speed. This method combines the trend characteristics of hyperbola, introduces the concept of column segment trend, and improves the efficiency of clustering target hyperbola. At the same time, considering that in the actual detection process, local maximum values will be generated at the intersection of hyperbolas, causing the connection of target hyperbola to be disconnected, this method combines the similarity of clustering centers, fully utilizes the information of hyperbola, and improves the accuracy and accuracy of target positioning. Finally, by using the characteristics of the artifacts generated by the post-projection algorithm, a quantitative artifact concept is introduced to image the underground target, further improving the accuracy of target positioning.

[0071] Example one:

[0072] The technical scheme of the present application includes the following steps:

[0073] The technical scheme of the present application includes the following steps:

[0074] Step 1: Use the original B-Scan echo data, adopt the mean method, remove the ground reflection wave and antenna coupling wave, so as to highlight the reflection information of underground detection target.

[0075] Step 2: According to the echo image obtained in step 1, use image processing method, take the sidelobe of target echo as standard, and convert the echo image into binary image.

[0076] Step 3: According to the binary image obtained in step 2, a set of column segments (a segment of a column in the binary image with consecutive 1) is constructed from the binary image, and the main lobe width of the ricker wave is combined to perform denoising work, reduce the set of column segments, and speed up the subsequent operation.

[0077] Step 4: According to the column segment set obtained in step 3, the column segment set is clustered to contain different hyperbolic curve label clusters in combination with the trend characteristics of the detected target hyperbolic curve.

[0078] Step 5: According to the column segment clusters of different detected target hyperbolic curves obtained in step 4, in combination with the characteristics of ground coupling ground penetrating radar detection, a restricted fitting algorithm is used for hyperbolic curve fitting, and the related parameters of the hyperbolic curve are obtained.

[0079] Step 6: According to the hyperbolic curve related parameters obtained in step 5, the position, radius and underground dielectric constant of the target pipeline are inverted, and the feature vectors of different targets are further calculated. Then, according to the similarity of the feature vectors, further clustering is performed.

[0080] Step 7: According to the binary image obtained in step 2 and the underground dielectric constant obtained in step 6, back projection imaging is performed, and artifacts are removed to obtain the final underground pipeline image.

[0081] The application also includes the following technical features:

[0082] Step 2 includes: using the canny algorithm to preliminarily extract the edge of the target hyperbolic curve, taking the sidelobe of the target echo as the standard, calculating the mean value of the echo energy on the echo edge greater than the standard, and taking the mean value as the threshold value of binarization, and binarizing the echo image:

[0083]

[0084] In the formula, S(i,j) represents the amplitude of the i-th row and j-th column in the echo image, R represents the sidelobe peak value of the ricker wave, Edge represents the edge point matrix extracted by the canny algorithm, num represents the number of sampling points meeting the condition, sum(.) represents summation, and Th represents the threshold value of binarization.

[0085] Step 3 includes: constructing a set of column segments from the binary image according to the binary image obtained in step 2 and in combination with the width of the transmitted waveform, and each column segment has 6 parts, as follows:

[0086] Col Seg ={Col,Row up ,Row down ,Row mid ,Trend,Sort},Col down -Col up≥N (2)

[0087] where Col Seg represents a column segment in the set, Col is the column number where the column segment is located, Row up is the row number of the upper boundary of the column segment, Row down is the row number of the lower boundary of the column segment, Row mid is the row number of the center of the column segment, Trend is the hyperbolic trend of the adjacent column segment, the initial state is 0. Sort is the clustering label of the hyperbolic curve, and N represents the threshold of the received waveform width.

[0088] Further, the column segment greater than the threshold of the received waveform width is constructed into a column segment set.

[0089] Step 4 includes that according to the characteristics of the target hyperbolic curve, the column segment on the hyperbolic curve has three trends: rising, maintaining and falling trends. According to this characteristic, the trend change between the column segments of adjacent columns on the hyperbolic curve includes rising-rising, rising-maintaining, maintaining-maintaining, maintaining-falling and falling-falling, and the trend change abnormality caused by interference is excluded.

[0090] The algorithm starts to scan the column segments on the binary image from left to right, calculates the trend relationship between the column segments of adjacent columns according to the positional relationship between the column segments, and updates Trend. The specific updating principle is as follows:

[0091] When the positional relationship between the adjacent column segments meets the rising trend, Trend of the current column segment is updated to 1;

[0092] When the positional relationship between the adjacent column segments meets the maintaining trend, Trend of the current column segment is updated to 0;

[0093] When the positional relationship between the adjacent column segments meets the falling trend, Trend of the current column segment is updated to -1;

[0094] Further, if the current column segment is in the initial state, a new hyperbolic clustering label is given;

[0095] Further, if the current column segment and the column segment of the next column meet the trend change of the hyperbolic curve, the current clustering label is extended to the right; if they do not meet the trend change of the hyperbolic curve, for example, falling trend-rising trend, the extension is terminated, and a new hyperbolic label is created.

[0096] Further, if the next column of the column segment of the current hyperbolic label has two trends meeting the hyperbolic trend at the same time, such as rising trend-falling trend and rising trend-rising trend, a new hyperbolic point cluster is established, and the column segments of the old point cluster are copied to the new point cluster set, so as to form the column segment set of two labels.

[0097] Step 5 includes: according to the different target hyperbolic curve segment set obtained in step 4, combining the characteristics of ground coupling ground penetrating radar, using a restricted fitting algorithm to perform hyperbolic curve fitting. The hyperbolic curve can be expressed as:

[0098] Ay 2 +By+Cx 2 +Dx+E=0 (3)

[0099] Because the target hyperbolic curve detected by the ground coupling mode ground penetrating radar satisfies:

[0100]

[0101] In the formula, (x, y) represents the row and column number of the point on the hyperbolic curve, (x0, y0) represents the row and column number of the point where the target is located, dx represents the ground penetrating radar detection interval, dt represents the radar sampling point interval time, c represents the speed of light, ε represents the dielectric constant of the underground medium, and r represents the radius of the underground pipeline. Because dx, dt and c are known, the fitting parameters can be restricted by the following conditions:

[0102]

[0103] In the formula, min y represents the minimum y value in the clustering data set.

[0104] Step 6 includes: according to the hyperbolic curve related parameters obtained in step 5, combining the characteristics of ground coupling ground penetrating radar, inverting the position, radius and underground dielectric constant of the target pipeline, and establishing a feature vector according to the position and radius of different targets, and further clustering according to the similarity of the vector.

[0105]

[0106] In the formula, T x is the horizontal coordinate of the detected target, T y is the vertical coordinate of the detected target, that is, the depth of the detected target.

[0107] Further calculate the similarity of the feature vector composed of the positions and radii of different targets:

[0108]

[0109] Example two:

[0110] In combination Figure 1 , the application includes the following steps:

[0111] Step 1: using the original B-Scan echo data, using the mean method to remove the ground reflection wave and antenna coupling wave, so as to highlight the reflection information of the underground target.

[0112]

[0113] In the formula, S(i,j) represents the echo amplitude of the ground penetrating radar B-scan image at (i,j), i and j represent the row and column numbers of the echo data matrix respectively, i∈[1,M], j∈[1,N], M represents the row number, and N represents the column number.

[0114] Step 2: According to the echo image obtained in step 1, an image processing method is used to convert the echo image into a binary image.

[0115] Step 2.1: For the echo image processed in step 1, first, the canny algorithm is used to extract the edges of the target in the image to obtain a set of edge points Edge, and the ratio of the main lobe to the side lobe of one of the received echoes is obtained.

[0116]

[0117] Step 2.2: Find the maximum value of the points in the edge point set on the echo image, and multiply the maximum value by the ratio of the main lobe to the side lobe to obtain a preliminary threshold value. Calculate the mean value of the edge points greater than the preliminary threshold value, and use this value as the threshold value for binarization.

[0118]

[0119] In the formula, S(i,j) represents the amplitude of the i-th row and j-th column in the echo image, R represents the side lobe peak value of the ricker wave, Edge represents the edge point matrix extracted by the canny algorithm, num represents the number of sampling points that meet the condition, sum(.) represents summation, and Th represents the threshold value for binarization.

[0120] Step 3: According to the binary image obtained in step 2, since the echo of the underground target belongs to the reflection information of static targets and does not exist Doppler effect, the pulse width of the echo is equal to the pulse width of the transmitted signal.

[0121] Step 3.1: According to the pulse width of the transmitted waveform, a set of column segments is constructed from the binary image, as shown in the accompanying drawings. Figure 2 Each column segment has 6 parts, as follows:

[0122] Col Seg ={Col,Row up ,Row down ,Row mid ,Trend,Sort},Col down -Col up ≥N (11)

[0123] In the formula, Col SegCol is the column number of the column segment, Row up Row is the row number of the upper boundary of the column segment down Row is the row number of the lower boundary of the column segment mid Row is the row number of the center of the column segment, Trend is the hyperbolic trend of the adjacent column segment, the initial state is 0. Sort is the cluster label of the hyperbolic curve, and N represents the threshold of the received waveform width.

[0124] Step 3.2: Construct the column segment set greater than the threshold of the received waveform width into a column segment.

[0125] Step 4: According to the characteristics of the target hyperbolic curve, the column segment on the hyperbolic curve has three trends: rising, maintaining and falling. According to this characteristic, the trend change between the column segments of adjacent columns on the hyperbolic curve includes: rising-rising, rising-maintaining, maintaining-maintaining, maintaining-falling, falling-falling, and any other trend change is abnormal due to interference. Set the threshold of the number of column segments Num.

[0126] Step 4.1: The algorithm scans the column segments on the binary image from left to right, calculates the trend relationship between them according to the position relationship between the column segments of adjacent columns, as shown in the accompanying Figure 3 , and updates Trend. The specific updating principle is as follows:

[0127] When the position relationship between adjacent column segments meets the rising trend, update the Trend of the current column segment to 1;

[0128] When the position relationship between adjacent column segments meets the maintaining trend, update the Trend of the current column segment to 0;

[0129] When the position relationship between adjacent column segments meets the falling trend, update the Trend of the current column segment to -1;

[0130] Step 4.2: If the current column segment is in the initial state, assign a new hyperbolic cluster label, as shown in the accompanying Figure 4 ;

[0131] Step 4.3: If the current column segment and the column segment of the next column meet the trend change of the hyperbolic curve, the current cluster label is extended to the right; if they do not meet the trend change of the hyperbolic curve, for example, falling trend-rising trend, terminate the extension, and create a new hyperbolic label, as shown in the accompanying Figure 5 .

[0132] Step 4.4: If the next column of the current hyperbola label contains both hyperbola trends, such as an upward trend followed by a downward trend and an upward trend followed by an upward trend, create a new hyperbola point cluster and copy the old point cluster column to the new cluster set, thus forming a column segment set with both labels, as shown in the appendix. Figure 5 As shown.

[0133] Step 4.5: Calculate the number of column segments contained in different cluster labels and compare it with Num, discarding cluster labels with fewer than Num.

[0134] Step 5: Based on the characteristics of the hyperbolic echo of underground pipelines detected by ground-coupled ground-penetrating radar, a constrained fitting algorithm is used to fit the hyperbolic segment clusters of different targets obtained in Step 4.

[0135] Step 5.1: From the standard hyperbola equation, we know that a hyperbola can be formulated as:

[0136] Ay 2 +By+Cx 2 +Dx+E=0 (12)

[0137] Step 5.2: As attached Figure 6 As shown, based on the analysis of ground-penetrating radar detection methods coupled to the ground, the formula for the hyperbola generated by an underground pipeline is as follows:

[0138]

[0139] In the formula, (x, y) represents the row and column numbers of the hyperbola, (x0, y0) represents the row and column numbers of the target location, dx represents the ground-penetrating radar detection interval, dt represents the radar sampling point interval time, c represents the speed of light, ∈ represents the dielectric constant of the underground medium, and r represents the radius of the underground pipeline. The fitting parameters can be constrained by the following conditions:

[0140]

[0141] In the formula, min y This represents the smallest y value in the clustered dataset.

[0142] Step 5.3: Based on the above analysis, the Col and Row values ​​in the column segment sets of different targets obtained in Step 4 are... mid A restricted fitting algorithm can be used to achieve fast fitting of hyperbolas.

[0143] Step 6: Based on the hyperbolic parameters obtained in Step 5 and the hyperbolic characteristics analyzed in Step 5.2, the location, radius, and underground dielectric constant of the target pipeline are deduced, and clustering is further performed based on the similarity of the location vectors of different targets.

[0144] Step 6.1: Inversion of the relative permittivity of the subsurface medium and the relevant parameters of the subsurface targets according to the fitting results obtained in Step 5.

[0145]

[0146] where T x is the horizontal coordinate of the target, T y is the vertical coordinate of the target, i.e. the depth of the target.

[0147] Step 6.2: Calculation of the similarity of the feature vectors composed of different target positions and radii:

[0148]

[0149] where (T x1 , T y1 ) represents the horizontal and vertical coordinates of the No. 1 target, (T x2 , T y2 ) represents the horizontal and vertical coordinates of the No. 2 target, and r1 and r2 represent the radii of the detected targets. Those with a similarity less than 0.2 are classified as the same subsurface pipeline.

[0150] Step 7: Post-projection imaging according to the binary image obtained in Step 2 and the subsurface permittivity obtained in Step 6, and removal of artifacts to obtain the final subsurface pipeline image.

[0151] Step 7.1: First, divide the detection area into a number of pixel points to construct the imaging domain matrix U. And calculate the echo time delay of each pixel point from the antenna:

[0152]

[0153] where d is the two-way travel distance from radar detection point i to pixel point (x j , y j ), c is the speed of light, ε is the permittivity of the medium, and t(i, x j , y j ) represents the propagation time from radar detection point i to pixel point (x j , y j ).

[0154] Further calculate the sample point sequence number with a time delay of t(i, x j , y j ) in the echo data received at this detection point:

[0155] M i = round(t(i, x j , y j ) / dt) (18)

[0156] where M i denotes the M i th sample point in the echo signal of the i

[0157] Step 7.2: Repeat step 7.1 to calculate the echo amplitude of a given pixel at all detection points, and add all echo values to obtain the assumed reflection coefficient u of the pixel, where N represents the total number of detection points, that is, the column number of the echo matrix, and S(i, M i ) denotes the data in the i i th column and M

[0158]

[0159] Step 7.3: Repeat steps 7.1 and 7.2 to calculate the assumed reflection coefficient u of all pixel points in the imaging domain to obtain a high-resolution two-dimensional image U.

[0160] Step 7.4: Set an energy threshold to eliminate artifacts in the two-dimensional image U, and the filtering operation is as follows:

[0161]

[0162] where n t denotes the number of hyperbolic curves identified in the echo image, and U(i, j) denotes the pixel energy at the (i, j) position in the image generated by the back-projection algorithm.

[0163] The following embodiments are given in conjunction with specific parameters.

[0164] Table 1 shows the relevant parameters of the detection model, including the position of the buried object, the frequency of the transmitted electromagnetic wave, the number of detection points, the sampling frequency, the transmitted waveform, the radius of the buried object, and the relative permittivity of the underground medium.

[0165] Table 1

[0166] Parameter Correlation value Electromagnetic wave form Ricker pulse Electromagnetic wave frequency 1 GHz Relative permittivity of the underground medium 6 Number of probe points 90 Probe point step length 0.02m Sampling frequency 2.1199e+11 Position of the 1st pipe buried object (0.36 m, 0.3 m) Radius of the 1st buried object 0.01m Position of the 2nd pipe buried object (0.7 m, 0.4 m) Radius of the 2nd buried object 0.01m Position of the 3rd pipe buried object (1.04 m, 0.3 m) Radius of the 3rd buried object 0.01m Position of the 4th pipe buried object (1.7 m, 0.4 m) Radius of the 4th buried object 0.01m

[0167] The GPRmax 3.0 electromagnetic simulation software is used to simulate the actual underground scene, and the method in embodiment one of the present application is used for solving. The accuracy of the results is shown in Table 2:

[0168] Table 2 shows the inversion result error

[0169]

[0170] Table 2 shows the simulation results. From the table, it can be seen that the present method can well detect the position information of the underground pipeline and the dielectric constant information of the underground environment.

[0171] Secondly, Figure 7Display the original echo image, Figure 8 Display the binary echo image, under this model, because the opening of the hyperbolic curve of target 2 is destroyed, the existing clustering algorithm cannot focus on the position of target 2, thereby causing the phenomenon of missing target 2. Figure 9 The results of clustering the intersection of hyperbolic curves using the method herein. Figure 10 The results of the clustering algorithm based on hyperbolic trend proposed by the method herein, it can be seen that the echo hyperbolic curve of target 2 has a complex opening structure, but the algorithm can still well cluster the left and right curve information of the target hyperbolic curve. Figure 11 The fitting results of the fitting algorithm using the method herein. Figure 12 The results of inversion imaging using the back projection algorithm, Figure 13 The results after using hierarchical filtering.

[0172] The present application relates to the field of underground pipeline detection, and particularly relates to a method for processing echo information and realizing multi-target identification and positioning imaging under the condition of multi-target echo information interweaving. The present application aims to design a new type of B-Scan image interpretation method, mainly to solve the problems of missing detection and information loss of the existing interpretation model under the condition of multi-target echo information interweaving. The present method combines the hyperbolic curve characteristics and introduces the concept of column segment trend to improve the clustering efficiency of the target hyperbolic curve. At the same time, considering that in the actual detection process, local maximum values will be generated at the intersection of hyperbolic curves, causing the problem of disconnection of target hyperbolic curves, the present method combines the method of calculating similarity, fully utilizes the information of hyperbolic curves, and improves the positioning accuracy and accuracy of targets. Finally, by using the characteristics of artifacts generated by the back projection algorithm, a quantitative artifact concept is introduced to image underground targets, and the accuracy of target positioning is further improved.

[0173] The above only describes the preferred embodiments of the present application and is not used to limit the present application. For those skilled in the art, the present application can have various changes and modifications. Any modification, equivalent replacement, improvement, etc. made within the spirit and principles of the present application shall be included in the protection scope of the present application. The components in the experimental examples can be realized by using the existing technology.

Claims

1. A method for identification, location and imaging of subsurface targets, characterized in that, The method comprises the following steps: Step 1: using original B-scan echo data, using mean method to remove ground reflection and antenna coupling wave, and highlighting underground target reflection information; Step 2: using image processing method, converting the echo image obtained in step 1 into a binary image; Step 3: constructing a column segment set from the binary image according to the pulse width of the transmitted waveform, and performing denoising on the column segment set according to the main lobe width of the ricker wave to obtain a denoised column segment set; Step 4: clustering column segment clusters containing hyperbolic curve labels corresponding to different targets from the column segment set obtained in step 3 according to the trend characteristics of the target hyperbolic curve; Step 5: according to the characteristics of the echo hyperbolic curve of the ground coupling ground penetrating radar detecting the underground target, using a limited fitting algorithm to perform hyperbolic curve fitting on the hyperbolic curve column segment clusters of different targets obtained in step 4; Step 6: according to the hyperbolic curve parameters obtained in step 5, inverting the position, radius and underground dielectric constant of the target, and further clustering according to the similarity of the feature vectors of different targets; Step 7: according to the binary image and the underground dielectric constant obtained in steps 2 and 6, performing back projection imaging and removing artifacts to obtain an underground target image.

2. The method of claim 1, wherein: The method of converting the ground penetrating radar B-scan image into a binary image by using image processing method is specifically: using a canny algorithm to extract the edge of the target hyperbolic curve in the image, taking the sidelobe of the echo waveform ricker wave as a standard, calculating the mean value of the echo energy of the echo edge greater than the standard, and the mean value is used as the threshold value for binaryzation, and the echo image is binaryzation: wherein S(i,j) represents the amplitude of the i-th row and j-th column in the echo image, R represents the sidelobe peak value of the ricker wave, Edge represents the edge point matrix extracted by the canny algorithm, num represents the number of sampling points meeting the condition, sum(.) represents the summation operation, and Th represents the threshold value for binaryzation.

3. The method of claim 1, wherein: The denoised column segment set is specifically: Col Seg = {Col, Row up , Row down , Row mid , Trend, Sort}, Col down - Col up ≥ N where Col Seg represents one of the set of column segments, Col represents the column number of the column segment, Row up represents the upper boundary of the row number of the column segment, Row down represents the lower boundary of the row number of the column segment, Row mid represents the row number of the center of the column segment, Trend represents the hyperbolic trend of the adjacent column segment, the initial state is 0, Sort represents the clustering label of the hyperbolic curve of the column segment, and N represents the threshold of the received waveform width.

4. The method of claim 3, wherein: The clustering of different hyperbolic curve column segment clusters from the column segment set obtained in step 3 according to the trend characteristics of the target hyperbolic curve in step 4 is specifically: starting to scan the column segments in the column segment set from left to right, if the current column segment is in the initial state, creating a new hyperbolic curve clustering label, and searching whether there is a column segment in the next column of the column in the column segment set that meets the trend change condition, the trend change condition includes: rising to rising, falling to falling, rising to keeping and rising to falling, if there is a column segment that meets the condition, the clustering label of the current column segment is assigned to the next column segment, and the trend is updated according to the change trend; if not, terminate the current clustering label expansion, and re-establish the clustering label of the hyperbolic curve, and judge the number of column segments of the current clustering label, if it is lower than the threshold value Num, delete the current clustering label; if the next column segment of the current column segment simultaneously exists the column segments with rising and falling trends, a new clustering label is established, and the old label is copied to the clustering set of the new label, forming two clustering label sets.

5. The method of claim 1, wherein: The hyperbolic curve fitting on the hyperbolic curve column segment clusters of different targets obtained in step 4 by using a limited fitting algorithm in step 5 is specifically: the hyperbolic curve is: Ay 2 +By+Cx 2 +Dx+E=0 Row in different target column segment set according to step 3 mid Hyperbolic parameter fitting with Col; The target hyperbola of ground penetrating radar detection in ground coupling mode satisfies: In the formula, (x, y) represents the row and column numbers of points on the hyperbola, (x0, y0) represents the row and column numbers of points where the target is located, dx represents the ground penetrating radar detection interval, dt represents the sampling point interval time of the radar, c represents the speed of light, ε represents the dielectric constant of the underground medium, and r represents the radius of the underground pipeline. The fitting parameters are limited by the following conditions: where min y represents the minimum y value in the clustered data set, using the restricted fitting algorithm.

6. The method of claim 1, wherein: The position, radius and dielectric constant of the target are inversely calculated from the hyperbola parameters obtained in step 5 according to step 6, and the specific steps are as follows: where A, B, D, E are hyperbolic fitting parameters, T x is the target horizontal coordinate, T y is the target vertical coordinate, i.e., the target depth, ∈ represents the dielectric constant of the underground medium, dt represents the radar sampling point interval time, c represents the speed of light, and r represents the radius of the underground pipeline.

7. The method of claim 1, wherein: The similarity of the feature vectors of different targets in step 6 is as follows: where (T x1 , T y1 ) denotes the horizontal and vertical coordinates of the first target, (T x2 , T y2 ) denotes the horizontal and vertical coordinates of the second target, r1 and r2 denote the radii of the detected targets, and Similarity is less than a given threshold value.

8. The method of claim 1, wherein: According to the binary image obtained in step 2 and step 6 and the dielectric constant of the underground, the back projection imaging is performed, and the artifact removal processing is performed to obtain the underground target image, and the specific steps are as follows: Step 7.1: First, divide the detection area into a plurality of pixel points, construct an imaging domain matrix U, and calculate the echo time delay of each pixel point from the antenna: where d is the two-way travel distance from radar detection point i to pixel point (x j ,y j ), c represents the speed of light, ε is the dielectric constant of the medium, and t(i, x j ,y j ) represents the propagation time from radar detection point i to pixel point (x j ,y j ). The sample point number with time delay t(i, x j ,y j ) in the echo data received at the detection point is calculated. M i = round(t(i, x j ,y j ) / dt) In the formula, round(.) represents the approximation operation; Step 7.2: Repeat step 7.1 to calculate the echo amplitude of the given pixel at all detection points, and add all echo values to obtain the assumed reflection coefficient of the pixel point: Where N represents the total number of detection points, which is the number of columns in the echo matrix, S(i,M) i ) represents the Mth column of the i-th column in the echo data matrix. i Row data; Step 7.3: Repeat steps 7.1 and 7.2 to calculate the assumed reflection coefficient of all pixel points in the imaging domain, and obtain a high-resolution two-dimensional image matrix U; Step 7.4: Set a threshold to filter out artifacts to eliminate the artifacts of the two-dimensional image matrix U, and the specific steps are as follows: In the formula, n t U(i,j) represents the minimum number of clustered segments, and U(i,j) represents the pixel energy at position (i,j) in the image generated by the post-projection algorithm.

9. The method of any one of claims 1-8, wherein: The cross section of the target is circular or annular.

Images