Ground penetrating radar (GPR) inversion method based on time-frequency joint physical constraint

Abstract

The invention discloses a ground penetrating radar (GPR) inversion method based on time-frequency joint physical constraint, and relates to the technical field of ground penetrating radar detection. The method comprises the following steps: firstly, acquiring GPR common offset B-scan echo data, and constructing a data set through anti-noise processing and forward simulation; decomposing data through two-dimensional discrete wavelet transform, and outputting a fusion feature map through a time-frequency joint feature fusion module; introducing a Maxwell equation set, an electromagnetic wave propagation equation and an electromagnetic field boundary equation to construct multiple constraints; and finally, outputting a relative dielectric constant distribution diagram, and carrying out visualization and precision verification on the relative dielectric constant distribution diagram for detection. According to the method, the noise immunity and the small target recognition precision are improved through time-frequency deep fusion, the result reasonability is guaranteed through physical constraints, the efficiency is optimized based on the DLFUNet network, the method is suitable for underground pipe network distinguishing, roadbed disease diagnosis and other scenes, and the problems that a traditional method is insufficient in precision and poor in physical interpretability under the low signal-to-noise ratio are solved.

Description

Technical Field

[0001] This invention relates to the field of ground-penetrating radar (GPR) detection technology, and specifically to a ground-penetrating radar (GPR) inversion method based on time-frequency joint physical constraints. Background Technology

[0002] Ground-penetrating radar (GPR) inversion is a core non-destructive detection technology in modern engineering geology and municipal maintenance. By transmitting and receiving electromagnetic waves, it can reconstruct the three-dimensional structure of underground space and accurately obtain physical parameters such as dielectric constant, playing a crucial role in scenarios such as avoiding construction risks, investigating underground hazards, and managing municipal pipeline networks. However, existing GPR inversion technologies still face problems in practical applications, including insufficient inversion accuracy at low signal-to-noise ratios and poor physical interpretability. Insufficient inversion accuracy under low signal-to-noise ratio. The inhomogeneity of the underground medium causes severe scattering and attenuation of electromagnetic waves during propagation, resulting in very weak effective signals that are easily masked by environmental noise and ground clutter. Traditional inversion methods, or time-frequency joint methods relying on Fast Fourier Transform (FFT), have limitations in achieving both time and frequency domain resolution, making it difficult to effectively extract target signals from noise. This leads to low inversion accuracy and weak noise resistance for small targets, making them prone to missed or false detections.

[0003] Poor physical interpretability. Existing GPR inversion methods mostly adopt a purely data-driven approach, lacking effective constraints on the physical laws of electromagnetic wave propagation. Although such methods can utilize deep learning to process massive amounts of data, the inversion results often deviate from the actual electromagnetic propagation laws, making it difficult to explain their physical rationality and limiting the application of the technology in engineering scenarios where the reliability of detection results is highly critical. Summary of the Invention

[0004] This invention addresses the problems of insufficient accuracy and poor physical interpretability of existing GPR inversion techniques at low signal-to-noise ratios. It provides a GPR inversion method based on joint time-frequency physical constraints, which improves noise resistance and accuracy through deep time-frequency fusion, and ensures the physical rationality of the results by combining physical constraints and boundary constraints, thereby achieving efficient and accurate underground target detection.

[0005] To solve the above problems, the technical solution of the present invention is as follows: The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints includes the following steps: S1: Acquire B-scan echo data in GPR co-offset detection mode, construct forward simulation dataset using simulation software based on finite-difference time-domain method, and perform noise reduction and robustness processing on the data. S2: The preprocessed time-domain B-scan data is decomposed using two-dimensional discrete wavelet transform to obtain high-frequency and low-frequency components. The time-frequency joint feature fusion module is used to deeply fuse the time-domain spatial features, high-frequency components and low-frequency components to output a time-frequency fusion feature map. S3: Introduce Maxwell's equations as physical information constraints to build the foundation for data-driven inversion; S4: Construct a physical residual loss function based on the electromagnetic wave propagation equation to constrain the inversion results with physical laws; S5: Construct a boundary loss function based on the electromagnetic field boundary equation to optimize the boundary clarity and structural fidelity of the inversion results; S6: Output relative permittivity distribution map corresponding to the spatial dimensions of the input B-scan data, which is used for underground target detection after visualization processing and accuracy verification.

[0006] Furthermore, the noise immunity and robustness processing in S1 specifically involves introducing Gaussian white noise with a signal-to-noise ratio of -5dB to 5dB into some B-scan samples, and superimposing 1% to 15% of the dielectric constant random perturbation to simulate environmental electromagnetic interference and medium inhomogeneity in real detection.

[0007] Furthermore, the high-frequency components in S2 are used to extract edge, texture, and small target reflection signals from the image, while the low-frequency components are used to extract smooth background and geological layering information from the image.

[0008] Furthermore, the processing flow of the time-frequency joint feature fusion module includes: window segmentation, cutting the input feature map into non-overlapping rectangular local windows according to a preset size; cross-domain attention operation, defining Query as a time-domain spatial feature, Key as a high-frequency component, and Value as a low-frequency component, to achieve weighted fusion of high and low frequency information guided by spatial features; spatial relationship optimization, adding a relative position bias term to the attention score calculation, and applying a shift window mechanism to promote the flow of features between adjacent windows; feature restoration and output, inversely splicing the output results of each window, and outputting a time-frequency fused feature map after dimensionality reduction by convolutional layer channels and fusion of residual structures.

[0009] Furthermore, the Maxwell's equations in S3 are as follows: ; ; in Here, H is the curl operator, and H is the magnetic field strength (unit: A / m). Let E be the dielectric constant, E be the electric field strength (unit: V / m), t be time, and J be the current density. is the magnetic permeability.

[0010] Furthermore, the expression for the physical residual loss function in S4 is: ; in For conductivity, Using the Laplace operator, the residual of the equation is calculated as a physical loss term to ensure that the inversion result conforms to the laws of electromagnetic wave propagation.

[0011] Furthermore, the electromagnetic field boundary equations in S5 include: the continuity equation for the tangential component of the electric field: ; Equation for the tangential component of the magnetic field: Equation for the normal component of the electric displacement vector: ; and the continuity equation of the normal component of magnetic induction intensity ;in , E1t and E2t are the tangential components of the electric field on both sides of the interface, H1t and H2t are the tangential components of the magnetic field on both sides of the interface, and J s D is the free surface current density (unit: A / m). 1n D 2n These are the normal components of the electric displacement vectors on both sides of the interface. Free surface charge density (unit: C / m) 2 ), B 1n B 2n The normal components of the magnetic induction intensity on both sides of the interface are denoted as .

[0012] Furthermore, the accuracy verification in S6 uses root mean square error (RMSE) and structural similarity index (SSIM) as quantitative indicators. The inversion results need to be processed by morphological filtering and edge sharpening to eliminate residual background noise and enhance the clarity of the target body boundary.

[0013] Furthermore, the DLFUNet network includes an encoding module, a decoding module, a DBottleneck module, and an output module. The DBottleneck module is constructed using convolutional layers with dilation rates of 1, 6, 12, and 18 to achieve multi-scale feature extraction.

[0014] Furthermore, the application scenarios of the method include underground pipeline differentiation and quantitative diagnosis of roadbed defects; wherein underground pipeline differentiation is used to identify metal pipelines and non-metal pipelines, and quantitative diagnosis of defects is used to identify roadbed voids and water-bearing subsidence areas in the early stage and to estimate their volume.

[0015] The beneficial effects of this invention are as follows: 1. Significantly Improved Inversion Accuracy and Noise Resistance: This invention achieves accurate decomposition of time-frequency signals through two-dimensional discrete wavelet transform (DWT). Combined with the windowed cross-attention and residual fusion mechanism of the joint time-frequency feature fusion (IFTE) module, it breaks through the limitations of traditional methods that separate time-frequency features. Simultaneously, it captures information on target details (high-frequency components) and environmental background (low-frequency components), significantly improving the fine imaging effect of small-sized targets (such as termite nests and small pipes). Furthermore, the noise simulation and dielectric constant perturbation design in the data preprocessing stage enable the inversion algorithm to effectively extract valid signals even in low signal-to-noise ratio environments ranging from -5dB to 5dB. This avoids the problems of missed and false judgments caused by noise interference in traditional methods, significantly enhancing robustness in complex detection scenarios.

[0016] 2. Significantly Enhanced Physical Interpretability and Result Reliability: An innovative physical residual loss function derived from Maxwell's equations is introduced, forcing the inversion results to adhere to the fundamental laws of electromagnetic wave propagation. This addresses the core deficiency of purely data-driven inversion, which lacks physical meaning. Simultaneously, boundary constraints are constructed based on the four boundary equations of the electromagnetic field, ensuring the conservation and consistency of energy and momentum at the medium interface. This effectively improves the issues of blurred boundaries and smooth images in the inversion results, leading to more accurate target contour localization. These dual constraints ensure that the inversion results not only match the input signal at the data level but also conform to objective laws at the physical level, significantly enhancing their credibility and persuasiveness in engineering applications.

[0017] 3. Significantly Improved Computational Efficiency and Practicality: Leveraging the deep learning architecture of the DLFUNet network, integrating encoder-decoder modules, multi-scale DBottleneck modules, and parallel computing mechanisms, it efficiently processes massive GPR probe data, overcoming the time-consuming bottleneck of traditional iterative optimization algorithms and meeting the needs of modern engineering real-time imaging and rapid response. Simultaneously, the technology is adaptable to most GPR probe scenarios, requiring no manual parameter adjustment; it automatically completes the entire process of time-frequency decomposition, fusion, and constrained inversion simply by inputting raw B-scan data, significantly reducing the operational threshold and improving the convenience and efficiency of engineering applications.

[0018] 4. Deep Expansion of Application Scenarios: Leveraging high-precision quantitative mapping of dielectric constant and pseudo-color visualization technology, this invention successfully overcomes the application limitations of traditional GPR technology. It can effectively distinguish between underground metal and non-metal (PE / PVC) pipelines, solving the pain point of traditional methods' difficulty in detecting non-metallic pipelines. Furthermore, through dielectric constant anomaly gradient analysis, it can achieve early identification and volume estimation of hidden defects such as roadbed voids and water-bearing subsidence areas. Especially in the detection of termite nests in dams, by accurately capturing the difference in electromagnetic properties between low dielectric constant targets and the surrounding soil, it achieves precise location and extent of termite nests, providing reliable technical support for engineering hazard investigation, municipal maintenance, and other fields, with wide application value. Attached Figure Description

[0019] The invention will be further described below with reference to the accompanying drawings: Figure 1 This is a flowchart of the process of the present invention. Figure 2 This is a diagram of the DLFUNet network structure of the present invention. Figure 3 This is a comparison chart of the inversion results of the present invention. Detailed Implementation

[0020] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0021] The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints includes the following steps: S1: Data preparation and preprocessing, specifically including: Data Acquisition: Using the GPR common offset detection mode, the transmitting antenna (TX) and receiving antenna (RX) maintain a fixed distance and move synchronously along a preset straight survey line trajectory to record the two-way travel time and amplitude information of electromagnetic waves propagating and returning in the underground medium, thereby obtaining B-scan echo data reflecting the underground profile structure. The synchronous movement of the antennas in the common offset mode ensures the consistency of the electromagnetic wave propagation path, accurately captures the reflected signals of the underground medium, and directly obtains the basic data characterizing the reflection characteristics and location information of the underground target, providing reliable input for subsequent inversion. Simulation Dataset Construction: A large-scale forward modeling simulation dataset was constructed using the open-source electromagnetic simulation software GprMax, based on the Finite-Difference Time-Domain (FDTD) method. Diverse underground medium models were built within the simulation environment, covering various soil types with relative permittivity ranging from 3.0 to 15.0. Target objects such as low-dielectric-constant PVC pipes (dielectric constant 2.5), medium-dielectric-constant concrete structures, and high-dielectric-constant water bodies (dielectric constant 81) were embedded in the models. The relative permittivity distribution matrix of the underground space was extracted as training labels, establishing a precise mapping relationship between B-scan data and the distribution of physical parameters of the underground medium. The FDTD method was used to accurately simulate the propagation process of electromagnetic waves in complex media. By covering the media types of real-world detection scenarios with diverse models, the problem of scarce high-quality labeled data in real-world scenarios was addressed, providing sufficient and realistic training samples for deep learning models. Noise and robustness processing: Gaussian white noise with a signal-to-noise ratio of -5dB to 5dB is introduced into some B-scan samples, and 1% to 15% of the dielectric constant random perturbation is superimposed to simulate environmental electromagnetic interference and medium non-uniformity in real detection. By artificially simulating the noise and medium uncertainty in real detection, the model can adapt to complex environments during the training phase, significantly improving the robustness of the inversion algorithm in actual engineering scenarios and avoiding model overfitting caused by the "pure" simulation data. S2: Time-frequency feature decomposition and fusion, specifically including: Time-frequency decomposition: The preprocessed time-domain B-scan data is processed using the two-dimensional discrete wavelet transform (DWT) operator to separate high-frequency components (HFdata) and low-frequency components (LFdata). HFdata is used to extract edges, textures, and small target reflection signals from the image, while LFdata is used to extract smooth background and geological layering information. By utilizing the characteristics of DWT in decomposing signals at different scales, the time-domain and frequency-domain information can be separated and extracted, overcoming the limitation of traditional FFT in achieving both time-frequency resolution and capturing target details and environmental background information, thus laying the foundation for multi-dimensional feature fusion. Feature Fusion: The Inter-Time-Frequency Feature Fusion (IFTE) module performs deep fusion of temporal spatial features, HF data, and LF data. The IFTE module executes processing logic such as window segmentation, cross-domain attention operation, spatial relationship optimization, feature restoration, and output, ultimately outputting a time-frequency fused feature map. The attention mechanism is used to guide the fusion of high and low frequency information through spatial features, and the shift window mechanism is combined to optimize global feature association, breaking through the limitations of single-dimensional processing, realizing deep association of time and frequency information, and improving the comprehensiveness of feature representation. S3: Data-constrained inversion introduces Maxwell's equations as physical information constraints to ensure that the inversion results conform to the basic laws of electromagnetic wave propagation. By using Maxwell's equations, which describe the core physical laws of electromagnetic wave propagation, as constraints, the inversion results are forced to follow the characteristics of electromagnetic propagation, avoiding results that violate physical laws in purely data-driven inversion, and providing a foundation for subsequent constrained inversion. S4: Physically constrained inversion. A physical residual loss function is constructed based on the electromagnetic wave propagation equation. The equation residual is calculated as the physical loss term, and the constrained inversion result conforms to the electromagnetic wave propagation law. This loss function is constructed based on the predicted results of electric field strength and dielectric constant. The equation residual is solved by combining the time derivative and spatial derivative of the electric field with the model's predicted dielectric constant as the physical loss term. This ensures that the dielectric constant output by the model is consistent with the electromagnetic wave propagation law, and improves the physical interpretability of the inversion result. S5: Boundary Constraint Inversion. A boundary loss function is constructed based on the electromagnetic field boundary equation to ensure the clarity of geological boundaries and buried object boundaries in the inversion results, and to ensure that physical quantities such as energy and momentum are conserved and consistent at the boundaries. The constraint function constructed by the electromagnetic field boundary equation solves ill-posed problems and image smoothing defects on the inversion boundary, and balances the numerical accuracy and structural fidelity of the inversion results, making the target body boundary contour clearer. S6: Output and application of inversion results, specifically including: Output results: The output is a relative permittivity distribution map that strictly corresponds to the spatial size of the input B-scan data, realizing the mapping from qualitative radar echo signals to quantitative physical parameter fields; by accurately mapping the electromagnetic properties of underground media, it provides a quantitative basis for target identification and attribute judgment. Visualization: Pseudo-color mapping technology is used to convert the dielectric constant matrix into a visual image, with different color levels corresponding to different media types; the distribution of underground media is presented intuitively through color level differences, clearly outlining the geometric contours, burial depth, and material properties of the target body, reducing the difficulty of manual identification; Accuracy verification: Root mean square error (RMSE) and structural similarity index (SSIM) are introduced as quantitative indicators to evaluate the consistency between the inverted image and the real underground structure model. Morphological filtering and edge sharpening are performed on the inversion results to eliminate residual background noise. The reliability of the inversion is objectively verified through quantitative indicators. The ability to present the details of the target body is further improved through post-processing to ensure that small features are not masked by noise. The method integrates the aforementioned time-frequency fusion, data constraint, physical constraint, and boundary constraint functions through the DLFUNet network. By leveraging the multi-dimensional feature fusion and high-efficiency processing capabilities of deep learning, it solves the technical problems of low GPR inversion accuracy and poor physical interpretability under low signal-to-noise ratio, and achieves high-precision detection of underground targets in complex environments.

[0022] The noise resistance and robustness processing in S1 specifically involves: introducing Gaussian white noise with a signal-to-noise ratio (SNR) of -5dB to 5dB into a portion of the B-scan samples through artificial addition. The intensity of this Gaussian white noise covers the entire spectrum of real-world detection scenarios, including weak interference (SNR 3dB to 5dB), moderate interference (SNR -1dB to 2dB), and strong interference (SNR -5dB to 0dB). Simultaneously, a random perturbation of 1% to 15% is superimposed on the dielectric constant parameter corresponding to this portion of the B-scan samples. The perturbation amplitude is uniformly distributed with a gradient of 1%, accurately simulating the non-uniformity of underground media caused by differences in particle composition and water content. Through this combination of noise introduction and perturbation superposition, the system comprehensively reproduces the environmental electromagnetic interference and medium non-uniformity uncertainty of different intensities in real-world detection. This allows the deep learning model to fully adapt to complex working conditions during the training phase, significantly improving the robustness of the inversion algorithm in actual engineering scenarios. It perfectly adapts to the training requirements of subsequent deep learning algorithms and enhances the model's generalization ability, preventing performance degradation in real-world data applications.

[0023] The high-frequency component (HFdata) in S2 focuses on extracting detailed features of the target, accurately capturing the instantaneous reflection signals of small-sized anomalies (such as tiny ant nests or narrow-diameter pipes). The low-frequency component (LFdata) focuses on extracting the background and macroscopic structure of the underground environment, clearly presenting the stratigraphic interfaces and the distribution characteristics of a large range of media. Relying on the multi-scale decomposition characteristics of DWT, the high- and low-frequency components complement each other in multiple dimensions. The detailed information captured by the high-frequency component provides support for the accurate positioning and identification of the target, while the background structure presented by the low-frequency component provides a basis for judging the environment in which the target is located. The two work together to provide comprehensive and three-dimensional information support for subsequent fusion and inversion, effectively making up for the shortcomings of insufficient information representation by a single component.

[0024] The specific processing flow of the Time-Frequency Joint Feature Fusion (IFTE) module is as follows: Window segmentation: The input temporal spatial feature map, HFdata feature map, and LFdata feature map are each segmented into several non-overlapping rectangular local windows according to a preset size of 32×32 pixels. The local windows are used as the basic units for parallel computing. Parallelization of feature processing is achieved by segmenting windows of fixed size, which greatly improves the efficiency of feature fusion and avoids local feature blurring caused by excessively large window sizes. Cross-domain attention computation: An attention mechanism is executed within each local window. The query vector is explicitly defined as a temporal-spatial feature, the key vector as a high-frequency component (HFdata), and the value vector as a low-frequency component (LFdata). Attention weights are obtained by calculating the similarity between the query and the key, and then the value is weighted and summed based on these weights. Guided by temporal-spatial features, accurate weighted fusion of high- and low-frequency information is achieved, so that the fused features retain spatial location correlation while also containing high-frequency details and low-frequency background information, thereby improving the relevance of feature representation. Spatial relationship optimization: When calculating the attention score, a relative position bias term based on pixel coordinate difference is introduced to correct the spatial association error of features within a local window; at the same time, a shifted window mechanism is applied to translate and overlap adjacent windows by a preset step size; the accuracy of spatial association of local features is improved by the relative position bias term, and the shifted window mechanism breaks the feature isolation of local windows, realizes feature flow and information interaction between adjacent windows, enhances global feature consistency, and avoids local feature fragmentation; Feature Restoration and Output: The output results of each local window are concatenated in reverse order according to the original cutting order to restore a full-size feature map with the same size as the input feature map. Then, a 1×1 convolutional layer is used for channel dimensionality reduction to fuse the multi-channel features corresponding to the three inputs into a single-channel feature map. Then, the residual structure is added to the original input features. The reverse concatenation ensures the integrity of the feature space, the 1×1 convolutional layer achieves efficient channel fusion, and the residual structure effectively preserves the original feature information and enhances the fused features. The final output time-frequency fusion feature map has both comprehensiveness and correlation with the original features, providing high-quality feature input for subsequent inversion.

[0025] The Maxwell's equations in S3 are as follows: (1) (2) in Here, H is the curl operator, and H is the magnetic field strength (unit: A / m). E is the product of the relative permittivity of the medium and the permittivity of vacuum, t is the electric field strength (unit: V / m), t is the time (unit: s), and J is the current density (unit: A / m²). The magnetic permeability of the medium (unit: H / m) is used as a physical information constraint to embed the equations into the model training process. By incorporating the derived electromagnetic propagation law into the loss function calculation, the model learning process is forced to follow the propagation characteristics of electromagnetic waves in the underground medium. This fundamentally avoids deviations from physical laws caused by overfitting the training data in pure data-driven inversion, ensuring that the inversion results are both consistent with the data characteristics and physically reasonable.

[0026] The physical residual loss function expression in S4 is: (3) in The dielectric conductivity is expressed in S / m. Here, E is the Laplace operator, and E is the electric field source term (unit: V / m). This loss function is constructed based on the predicted results of the electric field strength and dielectric constant. The specific calculation process is as follows: first, the dielectric constant at various underground locations is predicted through the model. And the electric field strength E, then calculate the second time derivative of E respectively. Second-order spatial derivative and first-order time derivative Compare the above derivative with the known Predictions and Substituting into equation (3), the difference between the left and right sides of the equation is taken as the physical residual, which is the physical loss term. By minimizing the physical loss term, the dielectric constant of the model output is forced to satisfy the basic equation of electromagnetic wave propagation, so that the inversion result not only matches the input signal at the data level, but also conforms to objective laws at the physical level, significantly improving the physical interpretability of the inversion result and providing a guarantee for the reliability of the result in engineering applications.

[0027] The electromagnetic field boundary equations in S5 include: Continuity equation for the tangential component of the electric field: (4) Where E1t is the tangential component of the electric field intensity on the medium 1 side at the interface, and E2t is the tangential component of the electric field intensity on the medium 2 side at the interface, both in volts per meter (V / m). This equation is derived from Faraday's law of electromagnetic induction, ensuring the continuity of the tangential propagation characteristics of the electric field at the interface of the medium. Magnetic field tangential component equation: (5) Where H1t represents the tangential component of the magnetic field intensity on the medium 1 side of the interface, and H2t represents the tangential component of the magnetic field intensity on the medium 2 side of the interface, both in amperes per meter (A / m) and j. s Let J be the free surface current density vector at the interface, in amperes per meter (A / m), and J s The direction is perpendicular to the direction of change of the tangential component of the magnetic field. This equation reflects that the discontinuity of the tangential component of the magnetic field is caused by the free surface current at the interface. The equation for the normal component of the electric displacement vector is: (6) in (D is the electric displacement vector), D 1nD represents the normal component of the electric displacement vector on the medium 1 side at the interface. 2n The electric displacement vector normal component on the medium 2 side at the interface is represented by the unit coulombs per square meter (C / m²). Let C be the free surface charge density at the interface, expressed in coulombs per square meter (C / m²). This equation is a direct manifestation of Gauss's law of electric flux, indicating that the difference in the normal components of the electric displacement is equal to the free surface charge density at the interface. Continuity equation for the normal component of magnetic induction: (7) in (B is the magnetic flux density), B 1n B represents the normal component of the magnetic induction intensity on the medium 1 side at the interface. 2n The normal component of the magnetic flux density on the medium 2 side at the interface is represented by Tesla (T) or Weber per square meter (Wb / m²). This equation is based on the physical fact that there are no magnetic monopoles in nature, ensuring the closure of magnetic field lines so that the magnetic flux entering the interface is equal to the magnetic flux leaving the interface. The boundary constraint function, composed of the above four equations, is transformed into a boundary loss term and incorporated into the model training. By following the laws of energy and momentum conservation at the electromagnetic boundary, it solves the ill-posed problem and image smoothing defects at the inversion boundary, while taking into account the numerical accuracy and structural fidelity of the inversion results. This makes the boundary contour of the underground target and the surrounding medium clearer and the positioning more accurate, avoiding misjudgment or omission of the target due to boundary ambiguity.

[0028] The visualization process in S6 specifically involves using pseudo-color mapping technology to divide the relative permittivity distribution matrix into multiple intervals based on its numerical range. The interval with a permittivity of 2.5 to 5.0 is mapped to blue (corresponding to low-permittivity materials such as cavities, PVC pipes, or ant nests); the interval with a permittivity of 5.1 to 30.0 is mapped to green (corresponding to medium-permittivity media such as soil and concrete); and the interval with a permittivity of 30.1 to 81.0 is mapped to red (corresponding to high-permittivity water bodies or aquifers). The distribution range of different media is visually presented through the differences in color gradients, clearly outlining... The system extracts the geometric contours, burial depth, and material properties of underground targets, reducing the complexity of manual identification. Morphological filtering employs a combination of opening and closing operations. First, opening operations remove small noise spots, and then closing operations fill in tiny cavities inside the target body, preventing noise from interfering with target identification and internal structure assessment. Edge sharpening uses the Laplacian operator to enhance the grayscale differences of boundary pixels, further eliminating residual background noise and strengthening the boundary distinction between the target body and the surrounding medium. This ensures the accurate presentation of detailed features such as tiny cracks and small ant nests, providing support for refined detection needs.

[0029] The DLFUNet network includes an encoding module, a decoding module, a DBottleneck (Dilated Bottleneck) module, and an output module. Its specific structure and workflow are as follows: The encoding module consists of a Patch embedding layer, a Conv1x1 layer, a BatchNorm layer, a ReLU activation function, and a Conv3x3 convolutional layer. The input B-scan data is segmented into fixed-size feature blocks by the Patch embedding layer and its dimensions are increased. The number of channels is then adjusted by the Conv1x1 layer, normalized by the BatchNorm layer, and non-linearly transformed by the ReLU activation function. Finally, the Conv3x3 convolutional layer extracts basic features, completing the feature encoding of the data layer by layer. This hierarchical encoding process gradually refines the core features of the data, while normalization and non-linear transformations enhance the feature representation capability, laying the foundation for subsequent multi-scale feature processing. The DBottleneck module uses Conv3x3 convolutional layers with dilation rates of 1, 6, 12, and 18 to construct multi-scale feature extraction channels. Each convolutional layer is followed by a BatchNorm layer and a ReLU activation function. By using convolutional kernels with different dilation rates, it captures feature information of different ranges. Small dilation rates focus on local details, while large dilation rates cover global structural features, realizing the extraction and enhancement of multi-scale features. This significantly improves the model's adaptability to targets of different sizes, effectively identifying both tiny ant nests and large-scale holes. The time-frequency fusion layer includes a DWT decomposition layer, an Inter-Time-Frequency Feature Fusion (IFTE) module, a window attention layer, and a shifted window mechanism layer. Encoded features are decomposed by the DWT layer to obtain high-frequency and low-frequency components, which are then fed into the IFTE module for time-frequency fusion. The window attention layer further strengthens local feature associations, and the shifted window mechanism layer promotes global feature interaction. This progressive processing achieves deep fusion of time-frequency information, preserving the accuracy of local features while ensuring the consistency of global features, thus improving the comprehensiveness of feature representation. The decoding module includes an Upsample layer, a Skip connection layer, a Conv3x3 convolutional layer, a BatchNorm layer, and a ReLU activation function. The Upsample layer restores the size of the fused feature map to match the input data. The Skip connection layer fuses the shallow detail features output from the encoding module with the deep features from the decoding process. After processing by the Conv3x3 convolutional layer, BatchNorm layer, and ReLU activation function, the feature decoding and reconstruction are completed. Upsampling restores the feature space size, and the Skip connection effectively preserves shallow detail information, avoiding detail loss during decoding and ensuring the spatial accuracy of the inversion result. Output module: sequentially passes through a Conv1x1 channel dimensionality reduction layer, a (conv3x3)x3+relu combined layer, and an AdaptiveAvgPolls layer; the Conv1x1 channel dimensionality reduction layer compresses multi-channel features into a single channel, the (conv3x3)x3+relu combined layer further refines the features, and the AdaptiveAvgPolls layer adaptively adjusts the feature map size, finally outputting a relative permittivity distribution map with the same spatial size as the input B-scan data; The total loss function of the network is a weighted sum of data loss, physical loss, and boundary loss. The data loss is the mean squared error (MSE) loss, the physical loss is the physical residual loss as described in claim 6, and the boundary loss is the boundary constraint function derived loss as described in claim 7. The weighting coefficients are set to 0.4, 0.3, and 0.3, respectively. The model is iteratively trained until convergence using the gradient descent algorithm, achieving synergistic optimization of data-driven and physical-driven approaches. This enables the model to efficiently process massive amounts of data while strictly adhering to physical laws, thus balancing inversion efficiency, accuracy, and physical interpretability.

[0030] The specific application scenarios and implementation methods of the method are as follows: Underground Pipeline Network Differentiation: By combining the relative permittivity distribution map obtained from the inversion with pseudo-color visualization results, and based on the difference in permittivity between metal pipelines (whose permittivity is mostly in the upper-middle range) and non-metallic (PE / PVC) pipelines (whose permittivity is in the lower range), automatic differentiation of the two types of pipelines is achieved by setting a permittivity threshold. Accurate classification is achieved by leveraging the differences in permittivity characteristics of pipelines of different materials, solving the problem of difficulty in detecting non-metallic pipelines in traditional methods. This method is suitable for municipal pipeline network inspection scenarios, improving the efficiency and accuracy of pipeline network sorting. Quantitative diagnosis of road subgrade defects: By analyzing the gradient changes in areas with abnormal dielectric constants, an abnormal gradient threshold is set to identify hidden defects such as road subgrade voids and water-bearing subsidence areas. Combining the pixel count, spatial size, and burial depth information of abnormal areas in the dielectric constant distribution map, a volume estimation model is established (volume = horizontal area of ​​abnormal area × average burial depth × medium density correction coefficient). Based on the correlation between abnormal dielectric constants and defects, early defect identification is achieved. The volume of defects is accurately estimated through the quantitative model, providing data support for early defect treatment and reducing the engineering risks caused by defect expansion. Dam termite nest detection: Utilizing the difference in dielectric constant between termite nests (with dielectric constants close to PVC pipes, falling within the low dielectric constant range) and the surrounding soil (with medium dielectric constant range), the location and extent of the termite nest are pinpointed through pseudo-color visualization. Combined with detailed features after morphological filtering, accurate nest identification and boundary delineation are achieved. By leveraging the difference in electromagnetic properties between the termite nest and the soil, this method overcomes the limitations of traditional detection methods in the complex environment of dams, accurately locating the nest's position and extent, providing reliable technical support for termite control in dams, and ensuring the structural safety of dams. In the above application scenarios, no manual parameter adjustment is required. Only the original B-scan data needs to be input, and the model can automatically complete the entire process of time-frequency decomposition, fusion, and constraint inversion. This adapts to the actual needs of different detection scenarios, greatly reduces the cost of manual operation, and improves detection efficiency and reliability.

[0031] The embodiments described in this specification are merely examples of implementations of the inventive concept. The scope of protection of this invention should not be considered as limited to the specific forms stated in the embodiments. The scope of protection of this invention also extends to equivalent technical means that can be conceived by those skilled in the art based on the inventive concept.

Claims

1. A ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints, characterized in that, Includes the following steps: S1: Acquire B-scan echo data in GPR co-offset detection mode, construct forward simulation dataset using simulation software based on finite-difference time-domain method, and perform noise reduction and robustness processing on the data. S2: The preprocessed time-domain B-scan data is decomposed using two-dimensional discrete wavelet transform to obtain high-frequency and low-frequency components. The time-frequency joint feature fusion module is used to deeply fuse the time-domain spatial features, high-frequency components and low-frequency components to output a time-frequency fusion feature map. S3: Introduce Maxwell's equations as physical information constraints to build the foundation for data-driven inversion; S4: Construct a physical residual loss function based on the electromagnetic wave propagation equation to constrain the inversion results with physical laws; S5: Construct a boundary loss function based on the electromagnetic field boundary equation to optimize the boundary clarity and structural fidelity of the inversion results; S6: Output relative permittivity distribution map corresponding to the spatial dimensions of the input B-scan data, which is used for underground target detection after visualization processing and accuracy verification.

2. The ground-penetrating radar (GPR) inversion method based on time-frequency joint physical constraints according to claim 1, characterized in that, The noise immunity and robustness processing in S1 specifically involves introducing Gaussian white noise with a signal-to-noise ratio of -5dB to 5dB into some B-scan samples, and superimposing 1% to 15% of the dielectric constant random perturbation to simulate environmental electromagnetic interference and medium inhomogeneity in real detection.

3. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The high-frequency components in S2 are used to extract the edge, texture, and small target reflection signals in the image, while the low-frequency components are used to extract the smooth background and geological layering information in the image.

4. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The processing flow of the time-frequency joint feature fusion module includes: window segmentation, which cuts the input feature map into non-overlapping rectangular local windows according to a preset size; cross-domain attention operation, which defines Query as temporal spatial feature, Key as high-frequency component, and Value as low-frequency component, to achieve weighted fusion of high and low frequency information guided by spatial features; spatial relationship optimization, which adds a relative position bias term to the attention score calculation and applies a shift window mechanism to promote the flow of features between adjacent windows; and feature restoration and output, which conversely splices the output results of each window, and outputs a time-frequency fused feature map after dimensionality reduction by convolutional layer channels and fusion of residual structures.

5. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The Maxwell's equations in S3 are as follows: ； ； in Here, H is the curl operator, and H is the magnetic field strength (unit: A / m). Let E be the dielectric constant, E be the electric field strength (unit: V / m), t be time, and J be the current density. is the magnetic permeability.

6. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The physical residual loss function expression in S4 is: ； in For conductivity, Using the Laplace operator, the residual of the equation is calculated as a physical loss term to ensure that the inversion result conforms to the laws of electromagnetic wave propagation.

7. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The electromagnetic field boundary equations in S5 include: the continuity equation for the tangential component of the electric field. ; Equation for the tangential component of the magnetic field: Equation for the normal component of the electric displacement vector: ; and the continuity equation of the normal component of magnetic induction intensity ;in , E1t and E2t are the tangential components of the electric field on both sides of the interface, H1t and H2t are the tangential components of the magnetic field on both sides of the interface, and J s D is the free surface current density (unit: A / m). 1n D 2n These are the normal components of the electric displacement vectors on both sides of the interface. Free surface charge density (unit: C / m) 2 ), B 1n B 2n The normal components of the magnetic induction intensity on both sides of the interface are denoted as .

8. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, In S6, the accuracy verification uses root mean square error (RMSE) and structural similarity index (SSIM) as quantitative indicators. The inversion results need to be processed by morphological filtering and edge sharpening to eliminate residual background noise and enhance the clarity of the target body boundary.

9. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The DLFUNet network includes an encoding module, a decoding module, a DBottleneck module, and an output module. The DBottleneck module is constructed using convolutional layers with dilation rates of 1, 6, 12, and 18 to achieve multi-scale feature extraction.

10. The ground-penetrating radar (GPR) inversion method based on joint time-frequency physical constraints according to claim 1, characterized in that, The application scenarios of the method include the differentiation of underground pipe networks and the quantitative diagnosis of roadbed defects. The underground pipe network is used to distinguish between metal and non-metal pipes, while the quantitative diagnosis of defects is used to identify roadbed voids and water-bearing subsidence areas in the early stage and to estimate their volume.

Images