2015 Vol. 26, No. 4
With the development of computational power, there has been an increased focus on data-fitting related seismic inversion techniques for high fidelity seismic velocity model and image, such as full-waveform inversion and least squares migration. However, though more advanced than conventional methods, these data fitting methods can be very expensive in terms of computational cost. Recently, various techniques to optimize these data-fitting seismic inversion problems have been implemented to cater for the industrial need for much improved efficiency. In this study, we propose a general stochastic conjugate gradient method for these data-fitting related inverse problems. We first prescribe the basic theory of our method and then give synthetic examples. Our numerical experiments illustrate the potential of this method for large-size seismic inversion application.
We present a method based on least-squares reverse time migration with plane-wave encoding (P-LSRTM) for rugged topography. Instead of modifying the wave field before migration, we modify the plane-wave encoding function and fill constant velocity to the area above rugged topography in the model so that P-LSRTM can be directly performed from rugged surface in the way same to shot domain reverse time migration. In order to improve efficiency and reduce I/O (input/output) cost, the dynamic encoding strategy and hybrid encoding strategy are implemented. Numerical test on SEG rugged topography model show that P-LSRTM can suppress migration artifacts in the migration image, and compensate amplitude in the middle-deep part efficiently. Without data correction, P-LSRTM can produce a satisfying image of near-surface if we could get an accurate near-surface velocity model. Moreover, the pre-stack P-LSRTM is more robust than conventional RTM in the presence of migration velocity errors.
Seismic illumination plays an important role in subsurface imaging. A better image can be expected either through optimizing acquisition geometry or introducing more advanced seismic migration and/or tomographic inversion methods involving illumination compensation. Vertical cable survey is a potential replacement of traditional marine seismic survey for its flexibility and data quality. Conventional vertical cable data processing requires separation of primaries and multiples before migration. We proposed to use multi-scale full waveform inversion (FWI) to improve illumination coverage of vertical cable survey. A deep water velocity model is built to test the capability of multi-scale FWI in detecting low velocity anomalies below seabed. Synthetic results show that multi-scale FWI is an effective model building tool in deep-water exploration. Geometry optimization through target oriented illumination analysis and multi-scale FWI may help to mitigate the risks of vertical cable survey. The combination of multi-scale FWI, low-frequency data and multi-vertical-cable acquisition system may provide both high resolution and high fidelity subsurface models.
Estimation of an accurate macro velocity model plays an important role in seismic imaging and model parameter inversion. Full waveform inversion (FWI) is the classical data-domain inversion method. However, the misfit function of FWI is highly nonlinear, and the local optimization cannot prevent convergence of the misfit function toward local minima. To converge to the global minimum, FWI needs a good initial model or reliable low frequency component and long offset data. In this article, we present a wave-equation-based reflection traveltime tomography (WERTT) method, which can provide a good background model (initial model) for FWI and (least-square) pre-stack depth migration (LS-PSDM). First, the velocity model is decomposed into a low-wavenumber component (background velocity) and a high-wavenumber component (reflectivity). Second, the primary reflection wave is predicted by wave-equation demigration, and the reflection traveltime is calculated by an automatic picking method. Finally, the misfit function of the l2-norm of the reflection traveltime residuals is minimized by a gradient-based method. Numerical tests show that the proposed method can invert a good background model, which can be used as an initial model for LS-PSDM or FWI.
A velocity model is an important factor influencing microseismic event locations. We review the velocity modeling and inversion techniques for locating microseismic events in exploration for unconventional oil and gas reservoirs. We first describe the geological and geophysical characteristics of reservoir formations related to hydraulic fracturing in heterogeneity, anisotropy, and variability, then discuss the influences of velocity estimation, anisotropy model, and their time-lapse changes on the accuracy in determining microseismic event locations, and then survey some typical methods for building velocity models in locating event locations. We conclude that the three tangled physical attributes of reservoirs make microseismic monitoring very challenging. The uncertainties in velocity model and ignoring its anisotropies and its variations in hydraulic fracturing can cause systematic mislocations of microseismic events which are unacceptable in microseismic monitoring. So, we propose some potential ways for building accurate velocity models.
In order to improve the efficiency of 3D near-surface velocity model building, we develop a layer-stripping method using seismic first-arrival times. The velocity model within a Common Mid-Point (CMP) gather is assumed to be stratified into thin layers, and the velocity of each layer varies linearly with depth. The thickness and velocity of the top layer are estimated using minimum-offset first-arrival data in a CMP gather. Then the top layer is stripped and the second layer becomes a new top layer. After removing the effect of the top layer from the former first-arrival data, the new first-arrival data are obtained and then used to estimate the parameters of the second layer. In this manner, the velocity model, being regarded as that at a CMP location, is built layer-by-layer from the top to the bottom. A 3D near-surface velocity model is then formed using the velocity models at all CMP locations. The tests on synthetic and observed seismic data show that the layer-stripping method can be used to build good near-surface velocity models for static correction, and its computation speed is approximately hundred times faster than that of grid tomography.
Reverberation is significant in shallow water and produces obvious notches in OBC spectra. It also degrades the quality of sections and increases the difficulty of processing and interpretation. This article presents the relationship between notch, shooting depth, and seabed depth based on the seismic convolution model. Forward modelling based on wave equation theory is used to verify this relationship. Dual-sensor summation is applied to suppress receiver-side multiples and remove notches according to the opposite response of geophones and hydrophones to down-going wave fields based on a detailed analysis of the OBC technique. The good results obtained in practical applications reveal the effectiveness of this method.
As the key technique of improved Hilbert-Huang transform (HHT), ensemble empirical mode decomposition (EEMD) has a good performance of eliminating mode mixing phenomenon, which has a strong impact on the observation of seismic information. However, the intrinsic mode functions (IMF) obtained from EEMD contain noises, so that it is required to find a more robust frequency estimation method to calculate the instantaneous frequency (IF) of IMF. For this reason, the improved HHT algorithm based on the damped instantaneous frequency (DIF) is proposed to overcome the shortage of EEMD. Compared with other IF estimation methods, the DIF has strong antinoise ability and high estimation accuracy. The test results of synthetic and real seismic data show that the proposed algorithm is feasible and effective for extracting seismic instantaneous attributes.
Common-reflection-point (CRP) gather is a bridge that connects seismic data and petrophysical parameters. Pre-stack attributes extraction and pre-stack inversion, both of them are important means of reservoir prediction. Quality of CRP gather usually has great impact on the accuracy of seismic exploration. Therefore, pre-stack CRP gathers noise suppression technology becomes a major research direction. Based on the vector decomposition principle, here we propose a method to suppress noise. This method estimates optimal unit vectors by searching in various directions and then suppresses noise through vector angle smoothing and restriction. Model tests indicate that the proposed method can separate effective signal from noise very well and suppress random noise effectively in single wavenumber case. Application of our method to real data shows that the method can recover effective signal with good amplitude preserved from pre-stack noisy seismic data even in the case of low signal to noise ratio (SNR).
The present study aims to reveal the recovering period of the postseismic fluid pressure in fault zone, offering an insight into earthquake recurrence. Numerical modeling is performed based on a 2D simple layered fault-valve model to simulate the fluid activities within the earthquake fault. In order to demonstrate the features of postseismic fluid pressure in natural state, the interference of tectonic movements is not considered. The recovering period of postseismic fluid pressure includes a suddenchanging period and a much longer fluctuating period. Modeling results show that fault permeability and porosity are sensitive parameters and reversely proportional to the recovering period of the fluid pressure in earthquake fault zone. When the permeability reduces from 10-15 to 10-18 m2, the recovering period increases from 400 to 2 000 yrs, correspondently. The upper and lower fluid pressures are separated by the valve seal, causing their fluctuations in opposite tendencies.
Lushan Earthquake (~Mw 6.6) occurred in Sichuan Province of China on 20 April 2013, was the largest earthquake in Longmenshan fault belt since 2008 Wenchuan Earthquake. To better understand its rupture pattern, we focused on the influences of fault parameters on fault slips and performed fault slip inversion using Akaike's Bayesian Information Criterion (ABIC) method. Based on GPS coseismic data, our inverted results showed that the fault slip was mainly confined at depths. The maximum slip amplitude is about 0.7 m, and the scalar seismic moment is about 9.47×1018 N·m. Slip pattern reveals that the earthquake occurred on the thrust fault with large dip-slip and small strike-slip, such a simple fault slip represents no second sub-event occurred. The Coulomb stress changes (ΔCFF) matched the most aftershocks with negative anomalies. The inverted results demonstrated that the source parameters have significant impacts on fault slip distribution, especially on the slip direction and maximum displacement.
Clustering earthquakes refer to the seismic events that occur closely with each other in time and space. Because their overlapping waveform records make it difficult to pick the first arrivals, the hypocenters of clustering earthquakes cannot be determined accurately by traveltime location methods. Here we apply a reverse-time imaging (RTI) method to map clustering earthquakes. Taking the advantage of directly using waveforms, the RTI method is capable to map either a single small earthquake or some densely distributed clustering earthquakes beneath a 2-D seismic array. In 3-D case the RTI method is successfully applied to locate the long-offset doublet earthquakes using the data from a set of sparsely distributed surface stations. However, for the same acquisition geometry, the RTI encounters challenges in mapping densely distributed clustering earthquakes. While it is obvious that improving the mapping of clustering earthquakes requires a denser receiver network with wider range of illumination angles, it is necessary to verify the actual resolution of the RTI method with synthetic data. In our study area in the Three Gorges region, our tests in 3-D case suggest that some events beneath the linear aligned sub-arrays have reasonable resolution.
Elastic impedance (EI) inversion has been widely used in industry to estimate kinds of elastic parameters to distinguish lithology or even fluid. However, it is found that conventional three-term elastic impedance formula is unstable even with slight random noise in seismic data, due to the ill-conditioned coefficient matrix of elastic parameters. We presented two-term Fatti elastic impedance inversion method, which is more robust and accurate than conventional three-term elastic impedance inversion. In our method, density is ignored to increase the robustness of inversion matrix. Besides, P-impedance and S-impedance, which are less sensitive to random noise, are inverted instead of VP and VS in conventional three-term elastic impedance. To make the inversion more stable, we defined the range of K value as a constraint. Synthetic tests claim that this method can obtain promising results with low SNR (signal noise ratio) seismic data. With the application of the method in a 2D line data, we achieved λρ, μρ and VP/VS sections, which matched the drilled well perfectly, indicating the potential of the method in reservoir prediction.
Seismic while drilling (SWD) is an emerging borehole seismic imaging technique that uses the downhole drill-bit vibrations as seismic source. Without interrupting drilling, SWD technique can make near-real-time images of the rock formations ahead of the bit and optimize drilling operation, with reduction of costs and the risk of drilling. However, the signal to noise ratio (SNR) of surface SWD-data is severely low for the surface acquisition of SWD data. Here, we propose a new method to retrieve the drill-bit signal from the surface data recorded by an array of broadband seismometers. Taking advantages of wavefield analysis, different types of noises are identified and removed from the surface SWD-data, resulting in the significant improvement of SNR. We also optimally synthesize seismic response of the bit source, using a statistical cross-coherence analysis to further improve the SNR and retrieve both the drill-bit direct arrivals and reflections which are then used to establish a reverse vertical seismic profile (RVSP) data set for the continuous drilling depth. The subsurface images derived from these data compare well with the corresponding images of the three-dimension surface seismic survey cross the well.
Multiscale strategies are very important in the successful application of waveform-based velocity inversion. The strategy that sequentially preceeds from long to short scale of velocity model, has been well developed in full waveform inversion (FWI) to solve the local mininum problem. In contrast, it's not well understood in the image-domain waveform tomography (IWT), which back-projects incoherent waveform components of the common image gather into velocity updates. IWT is less prone to local minimum problem but tends to build long-scale model with low resolution. In order to build both long- and short-scale model by IWT, we discuss several multiscale strategies restricted in the image domain. The strategies include model reparameterization, objective function switching and gradient rescaling. Numerical tests on Marmsousi model and real data demonstrate that our proposed multiscale IWT is effective in buidling velocity model with wide wavenumber spectrum.
Conventional time imaging techniques are not capable of producing accurate seismic imaging of the subsurface in the mountain front of the Tarim Basin, China. Their imaged structures have led to some major drilling failures before, bearing a disrepute that "their structural closures have wheels and their structural highs have springs". This article first lists the imaging challenges, and explains in a schematic why the time imaging techniques fail in this area. Then through a series of real data examples, it demonstrates that when there exist lateral velocity variations, depth imaging is the only solution to tackle the imaging challenges in this area. Depth imaging accounts for the complexity of the wavefield, therefore produces superior and geological plausible images. The core task in properly performing depth imaging is building the velocity model. This article stresses some the main aspects in this regard.
Marine seismic reflection surveys are often masked by strong water-bottom multiples that limit the use of data beyond the first multiple waves. In this study, we have successfully suppressed much of the multiple artifacts in the depth images of two of the marine seismic reflection profiles from the Los Angeles regional seismic experiment (LARSE) by applying reverse time migration (RTM). In contrast to most seismic reflection methods that use only primary reflections and diffractions, the two-way RTM migrates both primaries and multiple reflections to their places of origination: seabed multiples to the sea bottom and primaries to the reflecting interfaces. Based on the RTM depth sections of LARSE lines 1 and 2, we recognize five stratigraphic units from the sea bottom to a depth of 6 km. These units are Pliocene and younger strata, probably Miocene syntectonic strata, two deeper sequences of unknown age and lithology as well as Miocene volcanic layers on Catalina ridge. Several inferred igneous intrusions in the upper crust comprise a sixth unit. The existence of a thick sedimentary section in the Catalina Basin, which might include Paleogene and Cretaceous fore-arc strata, has important geologic significance. If borne out by further studies, significant revisions of current structural and stratigraphic interpretations of the California borderland would be warranted.
Reservoir architecture of meandering river deposition is complex and traditional seismic facies interpretation method cannot characterize it when layer thickness is under seismic vertical resolution. In this study, a seismic sedimentology interpretation method and workflow for point bar characterization is built. Firstly, the influences of seismic frequency and sandstone thickness on seismic reflection are analyzed by outcrop detection with ground penetrating radar (GPR) and seismic forward modeling. It is found that (1) sandstone thickness can influence seismic reflection of point bar architecture. With the increasing of sandstone thickness from 1/4 wavelength (λ) to λ/2, seismic reflection geometries various from ambiguous reflection, "V" type reflection to "X" type reflection; (2) seismic frequency can influence reservoirs' seismic reflection geometry. Seismic events follow inclined lateral aggradation surfaces, which is isochronic depositional boundaries, in high frequency seismic data while the events extend along lithologic surfaces, which are level, in low frequency data. Secondly, strata slice interpretation method for thin layer depositional characterization is discussed with seismic forward modeling. Lastly, a method and workflow based on the above study is built which includes seismic frequency analysis, 90º phasing, stratal slicing and integrated interpretation of slice and seismic profile. This method is used in real data study in Tiger shoal, the Gulf of Mexico. Two episodes of meandering fluvial deposition is recognized in the study layer. Sandstone of the lower unit, which is formed in low base level stage, distributes limited. Sandstone distribution dimension and channel sinuosity become larger in the upper layer, which is high base level deposition.
According to Vening Meinesz-Moritz (VMM) global inverse isostatic problem, either the Moho density contrast (crust-mantle density contrast) or the Moho geometry can be estimated by solving a non-linear Fredholm integral equation of the first kind. Here solutions to the two Moho parameters are presented by combining the global geopotential model (GOCO-03S), topography (DTM 2006) and a seismic crust model, the latter being the recent digital global crustal model (CRUST1.0) with a resolution of 1º×1º. The numerical results show that the estimated Moho density contrast varies from 21 to 637 kg/m3, with a global average of 321 kg/m3, and the estimated Moho depth varies from 6 to 86 km with a global average of 24 km. Comparing the Moho density contrasts estimated using our leastsquares method and those derived by the CRUST1.0, CRUST2.0, and PREM models shows that our estimate agrees fairly well with CRUST1.0 model and rather poor with other models. The estimated Moho depths by our least-squares method and the CRUST1.0 model agree to 4.8 km in RMS and with the GEMMA1.0 based model to 6.3 km.