RTM Imaging Research Articles

A comparative analysis between viscoacoustic forward and adjoint wave equations based on Maxwell, Kelvin-Voigt, and SLS rheological models

Most seismic imaging applications still consider the subsurface as a perfect acoustic or elastic medium, ignoring the viscosity of rocky materials. That viscosity converts part of the mechanical energy of seismic waves into heat, affecting the amplitude and phase of the wavefield when propagating through the physical medium. In this work, we derive, implement and analyze a set of forward and adjoint viscoacoustic equations based on the Maxwell, Kelvin-Voigt (KV), and standard linear solid (SLS) rheological models. These equations are described in the time–space domain allowing us to solve with the finite difference method. For this implementation, we used Devito - a domain-specific language (DSL) and code generation framework to design highly optimized finite difference kernels for use in inversion methods. We analyzed the attenuating behavior, comparing each equation’s dissipation and dispersion effects in these models, both forward modeling equations through comparisons between snapshots, seismograms, and traces, and adjoint equations through RTM images.

Multi-parameter full waveform inversion using only the streamer data based on the acoustic-elastic coupled wave equation

The ocean seismic pressure data observed by hydrophones contain the elastic information of the subsurface medium, although few studies of full waveform inversion (FWI) regard the pressure data as generated from elastic media. Based on the acoustic-elastic coupled equation, which can directly obtain the elastic pressure component, we propose an acoustic-elastic FWI (AEFWI) method to use elastic pressure data. First, we illustrated the similarities and differences in wave propagation in acoustic and elastic media using ocean-bottom models. Furthermore, we compared the effect of pressure data with different degrees of elastic effects on acoustic FWI (AFWI) and AEFWI. It was found that the S-wave velocity at the seabed interface had a significant effect on the AFWI. When the seabed is soft, the results of AFWI are acceptable, whereas when the seabed is hard, the results of AFWI are nearly wrong, and the presence of elastic effects has to be considered in the inversion at this time. Contrastingly, the proposed method can obtain more accurate P- and S-wave velocities using only pressure data. Finally, the proposed method was applied to a real streamer dataset of the East China Sea, and high wavenumber P- and S-wave velocities were obtained. The RTM images illustrated the accuracy of the inversion results.

A GPU implementation of least-squares reverse time migration

Least-squares reverse time migration (LSRTM) is a seismic imaging method that can provide higher-resolution image of the subsurface structures compared to other methods. However, LSRTM is computationally expensive. To reduce the computational time of LSRTM, GPU can be utilized. This leads to the objective of this work which is to develop a GPU implementation of LSRTM. In this work, the two-dimensional first-order acoustic wave equations are solved using the second-order finite difference on a staggered grid and a perfectly matched layer is used as an absorbing boundary condition. The adjoint-state method is used to compute the gradient of the objective function. A linear conjugate gradient method is used to minimize the objective function. Both forward- and backward-propagation of wavefields using the finite-difference method are performed on a single GPU using the NVIDIA CUDA library. For a verification purpose, the GPU program of LSRTM was applied to a synthetic data set generated from the Marmousi model. Numerical results show that LSRTM can provide an image with a higher resolution of subsurface structure compared to a conventional RTM image. For a computational cost issue, the GPU-version of LSRTM is significantly faster than the serial CPU-version of LSRTM.

Full waveform inversion of simultaneous long-offset data

SummaryThe standard workflow for velocity model building (VMB) in complex regimes is an interpretive process that requires time-consuming manual intervention, and remains an error-prone process that can produce suboptimal results. We discuss an application of Full Waveform Inversion (FWI) to automate the refinement of legacy velocity models generated by conventional workflows. We demonstrate our solution on a full-azimuth (FAZ) survey acquired in the Gulf of Mexico using dual-sensor streamers and blended sources in the form of simultaneous long-offsets (SLO). The dual-sensor acquisition provides low-frequency data, while the SLO configuration enables the recording of long offsets in excess of 16 km; thereby providing a unique dataset to meet the traditional data dependencies for robust FWI stability. The long offsets and low frequencies were key to using both refractions and reflections to update the deeper parts of the velocity model. Our solution directly inverts the simultaneous data as acquired in the field, and is particularly relevant to the industry growth in blended seismic data acquisition. We also employed an FWI velocity gradient that eliminates the migration isochrones; thereby removing the reflectivity imprint from the model updates. The FWI application to the field survey successfully refined the geometry of the salt bodies, including the base salt and the intra-salt enclosures. It also improved the RTM image; particularly the salt flanks and the subsalt reflectors.

Passive multiple reverse time migration imaging based on wave decomposition and normalized imaging conditions

With the development of seismic exploration, passive-source seismic data has attracted increasing attention. Ambient noise passive seismic sources exists widely in nature and industrial production. Passive seismic data is important in logging while drilling (LWD), large-scale structural exploration, etc. In this paper, we proposed a passive multiple reverse time migration imaging (PMRTMI) method based on wavefield decomposition and normalized imaging conditions method. This method differs from seismic interferometry in that it can use raw passive seismic data directly in RTM imaging without reconstruction of virtual active gather, and we use the wavefield decomposition method to eliminate the low frequency noise in RTM. Further, the energy normalized imaging condition is used in full wavefield decomposition, which can not only enhance the image quality of both edge and deep information but also overcome the wrong energy problem caused by uneven distribution of passive sources; furthermore, this method exhibits high efficiency. Finally, numerical examples with the Marmousi model show the effectiveness of the method.

Analytical reverse time migration with new imaging conditions for one-sided nondestructive evaluation of concrete elements using shear waves

Cross-correlation imaging condition of reverse time migration generates high-amplitude artifacts in the reconstructed images. In addition, due to geometrical spreading and scattering attenuations, this imaging condition assigns amplitudes to the points of the reconstructed image that are not a true representative of the reflection coefficient of the scanned medium at those points. These all can lead to a reduction in quality and misinterpretation of the reconstructed images. In this paper, we proposed new imaging conditions to mitigate high-amplitude artifacts and to assign more accurate amplitudes to an RTM image by considering geometrical spreading and scattering attenuation in concrete members when horizontal shear waves are used for imaging. We used data obtained by transmitting horizontal shear waves to 3D synthetic homogeneous and concrete specimens and demonstrated the effectiveness of the new imaging conditions.

Single-Step Data-Domain Least-Squares Reverse-Time Migration Using Gabor Deconvolution for Subsalt Imaging

Least-squares reverse-time migration (LSRTM) distinctly improves seismic imaging quality, but at an expensive computation overhead involving tens of iterations. We herein take a computationally cheaper single-step LSRTM solution, which intrinsically performs deblurring through a data-domain Wiener deconvolution. However, the Wiener filter mainly solves the signal estimation problems for stationary signals. Subsalt imaging often suffers from strong salt-related reflections and artifacts. The former and the latter give rise to strong amplitude variance and changed source wavelets in the demigrated data, increasing its nonstationarity and hindering the data-deblurring operation in the single-step LSRTM. To alleviate the nonstationarity during the data-domain deblurring, we consider a Gabor deconvolution method. Testing on the Sigsbee data sets shows that the Gabor deconvolution method is effective, producing subsalt images of more balanced events and fewer artifacts than the raw RTM image. The Gabor deconvolution-related result also outperforms the standard single-step LSRTM result with more robust behavior and better subsalt imaging quality.

Effective Q-compensated reverse time migration using new decoupled fractional Laplacian viscoacoustic wave equation

Seismic waves are attenuated and distorted during propagation because of the conversion of acoustic energy to heat energy. We focus on intrinsic attenuation, which is caused by [Formula: see text], which is the portion of energy lost during each cycle or wavelength. Amplitude attenuation can decrease the energy of the wavefields, and dispersion effects distort the phase of seismic waves. Attenuation and dispersion effects can reduce the resolution of image, and they can especially distort the real position of interfaces. On the basis of the viscoacoustic wave equation consisting of a single standard linear solid, we have derived a new viscoacoustic wave equation with decoupled amplitude attenuation and phase dispersion. Subsequently, we adopt a theoretical framework of viscoacoustic reverse time migration that can compensate the amplitude loss and the phase dispersion. Compared with the other variable fractional Laplacian viscoacoustic wave equations with decoupled amplitude attenuation and phase dispersion terms, the order of the Laplacian operator in our equation is a constant. The amplitude attenuation term is solved by pseudospectral method, and only one fast Fourier transform is required in each time step. The phase dispersion term can be computed using a finite-difference method. Numerical examples prove that our equation can accurately simulate the attenuation effects very well. Simulation of the new viscoacoustic equation indicates high efficiency because only one constant fractional Laplacian operator exists in this new viscoacoustic wave equation, which can reduce the number of inverse Fourier transforms to improve the computation efficiency of forward modeling and [Formula: see text]-compensated reverse time migration ([Formula: see text]-RTM). We tested the [Formula: see text]-RTM by using Marmousi and BP gas models and compared the [Formula: see text]-RTM images with those without compensation and attenuation (the reference image). [Formula: see text]-RTM results match well with the reference images. We also compared the field data migration images with and without compensation. Results demonstrate the accuracy and efficiency of the presented new viscoacoustic wave equation.

Application of the continuous wavelet transform in the extraction of directional data on RTM imaging condition wavefields

Low-frequency artifacts in reverse time migration result from unwanted cross-correlation of the source and receiver wavefields at non-reflecting points along ray-paths. These artifacts can hide important details in migrated models and increase poor interpretation risk.
 Some methods have been proposed to avoid or reduce the number of these artifacts, preserving reflections, and improving model quality, implementing other strategies such as modification of the wave equation, proposing other imaging conditions, and using image filtering techniques. One of these methods uses wavefield decomposition, correlating components of the wavefields that propagate in opposite directions.
 We propose a method for extracting directional information from the RTM imaging condition wavefields to obtain characteristics allowing for better, more refined imaging. The method works by separating directional information about the wavefields based on the continuous wavelet transform (CWT), and the analysis of the main changes on the frequency content revealed within the scalogram obtained by a Gaussian wavelet family.
 Through numerical applications, we demonstrate that this method can effectively remove undesired artifacts in migrated images. In addition, we use the Laguerre-Gauss filtering to improve the results obtained with the proposed method.

3D acoustic least-squares reverse time migration using the energy norm

We have developed a least-squares reverse time migration (LSRTM) method that uses an energy-based imaging condition to obtain faster convergence rates when compared with similar methods based on conventional imaging conditions. To achieve our goal, we also define a linearized modeling operator that is the proper adjoint of the energy migration operator. Our modeling and migration operators use spatial and temporal derivatives that attenuate imaging artifacts and deliver a better representation of the reflectivity and scattered wavefields. We applied the method to two Gulf of Mexico field data sets: a 2D towed-streamer benchmark data set and a 3D ocean-bottom node data set. We found LSRTM resolution improvement relative to RTM images, as well as the superior convergence rate obtained by the linearized modeling and migration operators based on the energy norm, coupled with inversion preconditioning using image-domain nonstationary matching filters.

Multidirectional-vector-based elastic reverse time migration and angle-domain common-image gathers with approximate wavefield decomposition of P- and S-waves

Elastic reverse time migration (E-RTM) has limitations when the migration velocities contain strong contrasts. First, the traditional scheme of P/S-wave mode separation is based on Helmholtz’s equations, which ignore the conversion between P- and S-waves at the current separation time. Thus, it contains an implicit assumption of the constant shear modulus and requires smoothing the heterogeneous model to approximately satisfy a locally constant condition. Second, the vector-based imaging condition needs to use the reflection-image normal, and it also cannot give the correct polarity of the PP image in all possible conditions. Third, the angle-domain common-image gathers (ADCIGs) calculated using the Poynting vectors (PVs) do not consider the wave interferences that happen at each reflector. Therefore, smooth models are often used for E-RTM. We relax this condition by proposing an improved data flow that involves three new contributions. The first contribution is an improved system of P/S-wave mode separation that considers the converted wave generated at the current time, and thus it does not require the constant-shear-modulus assumption. The second contribution is the new elastic imaging conditions based on multidirectional vectors; they can give the correct image polarity in all possible conditions without knowledge of the reflection-image normal. The third contribution is two methods to calculate multidirectional propagation vectors (PRVs) for RTM images and ADCIGs: One is the elastic multidirectional PV, and the other uses the sign of wavenumber-over-frequency ([Formula: see text]) ratio obtained from an amplitude-preserved approximate-propagation-angle-based wavefield decomposition to convert the particle velocities into multidirectional PRVs. The robustness of the improved data flow is determined by several 2D numerical examples. Extension of the schemes into 3D and amplitude-preserved imaging conditions is also possible.

Q-least-squares reverse time migration with viscoacoustic deblurring filters

Viscoacoustic least-squares reverse time migration, also denoted as Q-LSRTM, linearly inverts for the subsurface reflectivity model from lossy data. Compared with conventional migration methods, it can compensate for the amplitude loss in the migrated images due to strong subsurface attenuation and can produce reflectors that are accurately positioned in depth. However, the adjoint [Formula: see text] propagators used for backward propagating the residual data are also attenuative. Thus, the inverted images from [Formula: see text]-LSRTM with a small number of iterations are often observed to have lower resolution when compared with the benchmark acoustic LSRTM images from acoustic data. To increase the resolution and accelerate the convergence of [Formula: see text]-LSRTM, we used viscoacoustic deblurring filters as a preconditioner for [Formula: see text]-LSRTM. These filters can be estimated by matching a simulated migration image to its reference reflectivity model. Numerical tests on synthetic and field data demonstrate that [Formula: see text]-LSRTM combined with viscoacoustic deblurring filters can produce images with higher resolution and more balanced amplitudes than images from acoustic RTM, acoustic LSRTM, and [Formula: see text]-LSRTM when there is strong attenuation in the background medium. Our preconditioning method is also shown to improve the convergence rate of [Formula: see text]-LSRTM by more than 30% in some cases and significantly compensate for the lossy artifacts in RTM images.

Elastic least-squares reverse time migration

Least-squares migration can produce images with improved resolution and reduced migration artifacts, compared with conventional imaging. We have developed a method for elastic least-squares reverse time migration (LSRTM) based on a new perturbation imaging condition that yields scalar images of squared P- and S-velocity perturbations. These perturbation images do not suffer from polarity reversals that are common for more conventional elastic imaging methods. We use 2D synthetic and field-data examples to demonstrate the proposed LSRTM algorithm using the perturbation imaging condition. Our results show that elastic LSRTM improves the energy focusing and illumination of the elastic images and it attenuates artifacts resulting, for instance, from sparseness in the wavefield sampling and crosstalk of the P- and S-modes. Compared with RTM images, the LSRTM images provide more accurate relative amplitude information that is useful for reservoir characterization.

Three-dimensional acoustic wave equation modeling based on the optimal finite-difference scheme

Generally, FD coefficients can be obtained by using Taylor series expansion (TE) or optimization methods to minimize the dispersion error. However, the TE-based FD method only achieves high modeling precision over a limited range of wavenumbers, and produces large numerical dispersion beyond this range. The optimal FD scheme based on least squares (LS) can guarantee high precision over a larger range of wavenumbers and obtain the best optimization solution at small computational cost. We extend the LS-based optimal FD scheme from two-dimensional (2D) forward modeling to three-dimensional (3D) and develop a 3D acoustic optimal FD method with high efficiency, wide range of high accuracy and adaptability to parallel computing. Dispersion analysis and forward modeling demonstrate that the developed FD method suppresses numerical dispersion. Finally, we use the developed FD method to source wavefield extrapolation and receiver wavefield extrapolation in 3D RTM. To decrease the computation time and storage requirements, the 3D RTM is implemented by combining the efficient boundary storage with checkpointing strategies on GPU. 3D RTM imaging results suggest that the 3D optimal FD method has higher precision than conventional methods.

Improved seismic image by Q-compensated reverse time migration: Application to crosswell field data, west Texas

To test the effectiveness of the [Formula: see text]-compensated reverse time migration ([Formula: see text]-RTM) method, we applied it to crosswell seismic data from western Texas. This crosswell field survey was aimed at determining the boundaries and even the internal features of the reservoir. In this area, the reservoir geologic body exhibits strong attenuation that reduces high frequencies more rapidly. Thus, conventional acoustic RTM produces a dimmed image (reduced amplitude and low resolution) of the reservoir area and structures underneath. In contrast, [Formula: see text]-RTM is able to compensate for the attenuation effects during imaging. The [Formula: see text] and [Formula: see text] profiles needed for [Formula: see text]-RTM were produced by joint traveltime and frequency shift tomography. Preprocessing of the data was carried out to reduce noise, remove tube waves, and to separate up- and downgoing wavefields. Along with recovered high wavenumbers, the final [Formula: see text]-RTM image provided many details about geologic layers and structures. The lateral and vertical extent and internal structures within the reservoir unit were clearly determined. These geologic features were also correlated to the velocity profile and sonic logs. We concluded that [Formula: see text]-RTM imaging practically improved the image resolution in attenuating geologic media.