Ground Penetrating Radar (GPR) plays a crucial role in subsurface imaging, but its utility is often constrained by issues such as noise and data incompleteness. To overcome these limitations, this study explored the application of Denoising Diffusion Probabilistic Models (DDPMs) for generating high-resolution GPR images and proposed a comprehensive strategy integrating dataset preparation, model development, and hyperparameter optimization. Dataset preparation involved diverse data acquisition and sample augmentation to improve model generalizability. DDPM employed a TwoResUnet-based noise prediction network and was trained with 1,000 diffusion time steps and 200,000 training iterations. The results demonstrated that the model could generate GPR images with substantial diversity and feature realism, thereby closely reflecting real-world conditions. It provided a robust framework for synthetic GPR image generation and served as a valuable resource for AI-driven downstream applications. The high quality and variability of the outputs render them suitable for data augmentation, model training, and validation. Moreover, the model’s architecture and training parameters were readily adjustable, enabling customization for specific tasks or broader applications. This study contributes a novel strategy for expanding GPR image datasets and enhancing synthetic image fidelity, thereby improving the interpretability and practical utility of GPR data in archaeology, civil engineering, and environmental monitoring.

Insights into stress-dependent reflection and transmission (R/T) of elastic waves are crucial for geophysical applications involving stress- and angle-dependent seismic data. Stress-induced deformations in rocks typically exhibit an elastic-inelastic progressive transformation due to stress accumulation, crack closure, and crack growth prior to mechanical failure. Based on the truncated Taylor expansion of strain-energy function for small prestress deformations, the classical theory of acoustoelasticity with second- and third-order elastic constants (2oeCs and 3oeCs) is insufficient to handle inelastic deformations caused by the closure or opening of compliant microstructures. Referred to the statically deformed state, Fu and Fu [1] formulate the Padé Piola-Kirchhoff stress equations of motion for wave-induced small dynamic fields superposed on severely deformed rocks under large-magnitude prestress. The resulting Padé acoustoelastic equations of motion describe wave propagation in rocks with elastic-inelastic deformations in terms of 2oeCs, 3oeCs, Padé coefficients, and static stresses. The stress- and angle-dependent R/T of waves at the stressed interface with a large static bias remains largely unaddressed. We first compare conventional and Padé acoustoelastic plane-wave analyses for P- and S-wave velocities. Theoretical predictions by Padé acoustoelasticity agree well with laboratory measurements made on Portland and Massilon sandstones under confining and uniaxial prestresses, respectively. We formulate the Padé acoustoelastic R/T coefficients at an interface between rocks with locally uniform elastic-inelastic deformations, assuming negligible boundary dislocations, which holds under moderate inelastic deformation but may not apply to large plastic flow. We examine three typical loading modes (confining, uniaxial, and pure shear), with each considering the incidence of upgoing and downgoing P- and S-waves. The relevant Knott equations are derived for the estimation of stress-dependent R/T coefficients. We investigate the energy conservation at the interface to confirm the accuracy of stress- and angle-dependent R/T coefficients under hydrostatic stress conditions. Comparisons between different levels of prestresses show the important effect of inelastic strains due to the closure of compliant pores on the energy distribution of R/T waves. Stress-dependent intercepts, slopes, critical angles, and R/T energy ratios are closely related to the sensitivity of upper and lower media to pressure and the type and magnitude of prestress. Particularly, the stress-induced closure of compliant pores must be considered for amplitude-variation-with-offset analyses in geophysical applications.

Fractures are widely existing in rocks and are one of the controlling factors for the formation and distribution of hydrocarbon reservoirs. Electrical surveys are effective means for the quantitative detection and characterization of fractures. However, although all rocks are experiencing pore pressure, the deformation of fractures with pore pressure and its effects on the anisotropic electrical properties of fractured rocks, which can help for the improved interpretation of electrical survey data, remain poorly understood. We bridge this knowledge gap in this work. We first implement dedicated laboratory measurements of the anisotropic electrical conductivity of artificial sandstone samples with and without aligned penny-shaped fractures as a function of pore pressure. We then invert from the experimental data for the fracture parameters that characterize the deformation of the fractures. We finally theoretically model the effects of pore pressure induced fracture deformation on the anisotropic electrical properties of fractured rocks. The results show that the inverted porosity and aspect ratio of the fractures increase exponentially with pore pressure and exhibit linear correlation with each other as an implicit function of varying pore pressure. We also demonstrate that the variations in the fracture porosity caused by the varying pore pressure are the first order parameter affecting the pore pressure dependent electrical properties of fractured rocks. The results not only reveal the deformation of aligned fractures with pore pressure and its effects on the anisotropic electrical properties of fractured rocks, but also provide new insights for the detection and prediction of overpressure in fractured rocks.

Insights into wave propagation in prestressed porous rocks have great interesting in geophysical applications, such as remote monitoring in-situ stresses. Wave-induced small dynamic fields superposed onto statically deformed objects can be addressed traditionally by acoustoporoelastic theory that extends the classical acoustoelasticity of solids to porous media by incorporating Biot’s theory. Stress-induced deformations in porous rocks are of a progressively scaling feature with increasing prestress, undergoing linear elastic, hyperelastic (nonlinearly elastic), and inelastic deformations prior to mechanical failure. Conventional acoustoporoelastic theory is based on the Taylor expansion for the cubic strain-energy function with linear strains under finite-magnitude prestress. The theory with third-order elastic constants only accounts for stress-induced hyperelasticity, insufficient to handle inelastic deformations with nonlinear strains of compliant microstructures. We replace the Taylor expansion by the Padé approximation to the strain energy function, leading to Padé acoustoporoelastic equations for inelastic deformations under large-magnitude prestress. Theoretical results from plane-wave analyses agree well with the laboratory measurements of fluid-saturated Portland sandstones under confining and uniaxial prestresses. Finite-difference simulations are implemented to solve the first-order velocity-stress formulation of Padé acoustoporoelastic equations for elastic wave propagation in prestressed porous media under isotropic (confining) and anisotropic (uniaxial and pure-shear) prestresses. The resulting wavefield snapshots show the propagation of fast-P and slow-P and S waves in acoustoporoelastic media, illustrating stress-induced velocity orthotropies, strongly related to the direction of prestress. Comparisons with conventional acoustoporoelastic simulations provide a framework to estimate stress-induced inelastic strains from seismic responses in velocity and anisotropy.

The envelope inversion (EI) is an effective method to recover low-wavenumber components, which helps to produce a good initial model for full-waveform inversion (FWI). However, when the initial model is not capable of generating reflections, it brings enormous challenges for EI to invert deep low-wavenumber components, especially for short offset data. In contrast, reflection waveform inversion (RWI) uses demigration data to fit the observed reflection, which focuses on the transmission information with short offset seismic data. However, the cycle skipping and high nonlinearity of the RWI misfit still exist when the low-frequency information is absent. In this letter, we develop a joint traditional and reflection envelope inversion (JREI) that utilizes both reflection and transmission waves with envelope low-frequencies to recover low-wavenumber components in the shallow and deep regions simultaneously. We then use the FWI with high-frequency seismic data to obtain the high-wavenumber components. Applications to the modified Marmousi and Overthrust models demonstrate that the JREI can invert a better starting velocity model for the FWI to achieve a high-resolution inversion result.

Despite a record of great success, seismic imaging is regarded to be far from solved because of a big challenge in velocity building where potential failures are particularly salient for salt-related structures and fault-karst belts that pose strong velocity contrasts and dipping angles. The challenge results from the traveltime equivalence of velocity distortion and structure anamorphose during the imaging. It is closely related to two different wavenumber contents implicit in velocity models and seismic structures, respectively. For salt-related structures, the initial velocity model obtained by current industry techniques is characteristic of a low-wavenumber (LW) distortion, which partially destroys the topological invariance of subsurface structures that are implicitly present in surface seismic data and yields erroneous high-wavenumber (HW) reflectivities correspondingly. For fault-karst belts, the estimated initial velocity is too low in resolution to satisfactorily image the high-frequency part in seismic data. The fundamental question regarding how to couple distorted LW and HW components in velocity building remains largely unaddressed. We propose a general scheme by integrating LW velocity models, HW seismic structures (achieved from automatic interpretation of migrated volumes), and logging salt velocities to achieve interpretative velocity building. For assistance monitoring analysis, we introduce a slowness accumulation method to measure the wavenumber content of velocity fields. It provides a quick and easy way to monitor the change of wavenumber contents during the velocity updating. We also address a tentative scheme to assess the complexity of seismic imaging by a statistical description of velocity contrasts and dipping angles that dominate the wavenumber content of velocity models. Numerical experiments and case studies with overthrust salt structures from western China and fault-karst belts from northwest China demonstrate the applicability of interpretative seismic imaging.

Least-squares reverse time migration (LSRTM) attempts to produce a high-quality image for complicated subsurface structures. However, large amplitude discrepancies between the synthetic and observed seismic data are problematic for high-resolution imaging. Alternatively, correlative LSRTM (CLSRTM) misfit has been proposed to improve the imaging quality of complicated structures. However, the CLSRTM ignores the local characteristics of the 2-D seismic data. Thus, we developed a 2-D local correlative misfit for LSRTM (2-D-LCLSRTM) to improve the imaging resolution. In this case, a 2-D sliding window was used to obtain local-scale seismic data. A 2-D correlation method was then used to measure the similarity between the local-scale synthetic and observed data. Consequently, the 2-D-LCLSRTM misfit could reduce amplitude constraints and emphasize phase similarity, which has a potential for improving deep structure as it can boost weak seismic signals. To suppress the migration artifacts, we incorporated the sparsity promotion method with the 2-D-LCLSRTM misfit and used the fast iterative shrinkage-thresholding algorithm (FISTA) to solve it iteratively. In the numerical examples, a Marmousi model, a Salt model, and a marine field seismic dataset were used to test the effectiveness of the 2-D-LCLSRTM method. Compared with the commonly used RTM and sparsity promotion-based CLSRTM methods, the 2-D-LCLSRTM with sparsity promotion can better image deep reflectors and obtain high-resolution imaging results.

Viscous effects cause strong energy decay and waveform changes of seismic waves. These distortions can be corrected usingQ-compensated reverse time migration ($Q$-RTM) algorithms, and high-resolution migration images can be obtained. However, all$Q$-RTM methods require a relatively accurate$Q$model. The traditional wave-equation$Q$tomography can invert the$Q$model by eliminating the difference in peak frequency between the observed and synthetic early arrivals. However, this approach only can be used to invert the$Q$value only for large-scale applications or on the surface. Moreover, the reflected wave can also be applied in the extended domain, but its computational efficiency is low compared to that of the early arrivals. To overcome these problems, this work proposes a new wave-equation-based$Q$inversion methodology to evaluate more accurate underground$Q$values in local domain. The proposed approach is applicable both to the early arrivals and reflected waves. Accordingly, we first transform the seismic data into the local domain using a sliding Gaussian window to alleviate the crosstalk noise in nearby seismic waves. Then, we use an improved cross correlation algorithm between the amplitude spectra of the observed and synthetic data to calculate the peak frequency shift of each seismic event in local domain. Thus, the inversion accuracy of$Q$can be improved by using different kinds of waves. The numerical inversion examples demonstrate the ability of our proposed method to produce satisfactory inversion results, especially in high-attenuation and deep areas. The$Q$-RTM images further illustrate the accuracy of our proposed$Q$tomography method.

Q-compensated reverse-time migration ( Q-RTM) has been proven as an efficient method for seismic imaging with high fidelity. However, the source (forward) and receiver (backward) wavefields propagate along the opposite direction of time, and the recursive computation with the out-of-order access requires that all the wavefields of source propagation should be stored on the hard disk. For massive amounts of seismic data, saving the source wavefield from the central processing unit (CPU) [or graphics processing unit (GPU)] device to the disk and loading these data from the hard disk to the CPU (or GPU) device become extremely intensive in time and storage, which has been a bottleneck of Q-RTM. Several methods have been developed to reduce the huge wavefield storage in acoustic media, but are not applicable in the attenuated media. In this letter, we present a reversible hybrid absorbing boundary condition for Q-RTM, which is implemented by mixing the reversible attenuation and the random boundary conditions. Based on our developed new boundary, we just need to save the wavefield at the last one or two time steps in the forward process and then reconstruct the source wavefield in the time-reversal order. Numerical results demonstrate that the method can avoid the huge seismic data input and output (I/O) requirement and improve the computational efficiency dramatically.

Q-compensated reverse time migration (Q-RTM) can compensate seismic attenuation caused by the anelastic behavior of subsurface media. Although, the traditional Q-RTM has high computational efficiency, it is instable because the high frequency or wavenumber ambient noise is exponentially boosted during forward and backward seismic wavefield propagation. The existing stable Q-RTM method costs twice as much computing time and memory compared to the traditional Q-RTM. In this letter, we propose a new Q-RTM method to address the above issues simultaneously. First, a new viscoacoustic wave equation is derived based on a wavefield decomposition method to obtain the velocity-dispersion-only and viscoacoustic wavefields efficiently. Then, a theoretical framework of stable and high-efficiency Q-RTM method is proposed based on the velocity-dispersion-only and viscoacoustic wavefields. The synthetic example shows that the new stable Q-RTM results match well with the reference images (without attenuation images). Moreover, the field data images also demonstrate the stability and high-efficiency of our proposed Q-RTM method.

Singular value decomposition is an effective way to remove ground clutter in ground-penetrating radar (GPR) applications. The main limitation of this method is the selection of principal components to completely reconstruct the ground clutter or the target. To date, no effective criteria or technology have been developed. To solve this problem, a new method is proposed in this paper. The research and analysis presented herein reveal that the root-mean-square height (RMSH) of the first-arrival curve corresponding to the ground clutter has a welldefined positive relationship with the number of singular values associated with the principal components of the ground clutter. The number of singular values of these principal components (N) can be precisely determined based on the ground clutter by a linear function, N = 0.2634D + 1.3086, where D represents the RMSH value. In addition, an algorithm called developed histogram equalization was developed to improve the contrast to highlight the targets in denoized GPR data sets. The proposed strategy of extracting the principal components of the ground clutter and highlighting the contrast between the target signal and environmental reflections was successfully applied to the field GPR data, thus demonstrating the practicality and validity of the proposed approach.