STANFORD EXPLORATION PROJECT

Gabriel Alvarez, Brad Artman, James Berryman, Biondo Biondi, Morgan Brown, Weitian Chen, Jon Claerbout, Marie Clapp, Robert Clapp, John Etgen, Christopher Liner, Jesse Lomask, Steven Pride, Daniel Rosales, Paul Sava, William Symes, Leon Thomsen, Alejandro Valenciano, and Ioan Vlad

Number 112, September 2002

2002

Preface

The electronic version of this report makes the included programs and applications available to the reader. The markings [ER], [CR], and [NR] are promises by the author about the reproducibility of each figure result. Reproducibility is a way of organizing computational research that allows both the author and the reader of a publication to verify the reported results. Reproducibility facilitates the transfer of knowledge within SEP and between SEP and its sponsors.

ER

denotes Easily Reproducible and are the results of processing described in the paper. The author claims that you can reproduce such a figure from the programs, parameters, and makefiles included in the electronic document. The data must either be included in the electronic distribution, be easily available to all researchers (e.g. SEG-EAGE data sets), or be available in the SEP data library. We assume you have a UNIX workstation with Fortran, Fortran90, C, X-Windows system and the software downloadable from our website (SEP makerules, SEPlib, and the SEP latex package), or other free software such as SU. Before the publication of the electronic document, someone other than the author tests the author's claim by destroying and rebuilding all ER figures. Some ER figures may not be reproducible by outsiders because they depend on data sets that are too large to distribute, or data that we do not have permission to redistribute but are in the SEP data library.

CR

denotes Conditional Reproducibility. The author certifies that the commands are in place to reproduce the figure if certain resources are available. SEP staff have only attempted to make sure that the makefile rules exist and the source codes referenced are provided. The primary reasons for the CR designation is that the processing requires 20 minutes or more, or commercial packages such as Matlab or Mathematica.

${\rm \bf M}\;\;\;\;$

denotes a figure that may be viewed as a movie in the web version of the report. A movie may be either ER or CR.

NR

denotes Non-Reproducible figures. SEP discourages authors from flagging their figures as NR except for figures that are used solely for motivation, comparison, or illustration of the theory, such as: artist drawings, scannings, or figures taken from SEP reports not by the authors or from non-SEP publications.

Our testing is currently limited to LINUX 2.4 (using the Portland Group Fortran90 compiler), but the code should be portable to other architectures. Reader's suggestions are welcome. For more information on reproducing SEP's electronic documents, please visit
<http://sepwww.stanford.edu/research/redoc/>.

SEP-112 -- TABLE OF CONTENTS

Velocity

Ray based tomography using residual Stolt migration (ps.gz 10576K) (pdf 2660K) (src 46677K)
Clapp R. G.
In complex areas, residual vertical movement is not an effective method to calculate traveltime errors for image domain tomography. By scanning over velocity ratios using residual Stolt migration, a different criteria for a coherent image can be defined, and a traveltime error approximated. The resulting traveltime errors are more accurate, therefore the tomography procedure is more robust than the more traditional methodology. Results are shown on a complex 2-D dataset.
Wave-equation MVA applied to 4-D seismic monitoring (ps.gz 441K) (pdf 201K) (src 1270K)
Sava P., Etgen J., and Thomsen L.
4-D seismic processing is gradually maturing as a technique able to aid time lapse monitoring of seismic reservoirs Biondi et al. (1996); Lumley (1995). However, many limitations hamper the ability of 4-D seismic monitoring to produce reliable results in complicated reservoir situations. One such example is that of multi-layer reservoirs where changes at deeper levels are masked by those that occur at the top reservoir. In these cases, ...
Wave-equation MVA using diffracted data (ps.gz 7332K) (pdf 736K) (src 1575K)
Sava P. and Etgen J.
Migration velocity analysis (MVA) using diffracted data is not a new concept. Harlan (1986) addressed this problem and proposed a method to isolate diffraction events around faults. He also proposed a MVA technique applicable to simple geology, constant velocity or v(z), and quantifies the focusing quality using statistical tools. de Vries and Berkhout (1984) use the concept of minimum entropy to evaluate diffraction focusing, and apply this ...
A generalization of wave-equation migration velocity analysis (ps.gz 190K) (pdf 130K) (src 328K)
Sava P. and Symes W. W.
Wave-equation migration velocity analysis is derived from wavefield-continuation migration techniques. The velocity model is updated by optimizing certain properties of the migrated images. Different migration velocity analysis optimization criteria exist, of which two commonly used are fitting a target image and minimizing the differential semblance of migrated images. Both techniques are special cases of a general family of optimization functions. Fitting a target image is an attractive technique because we can guide the solution in the desired direction. However, we can only progress in small steps with the target image being kept within the Born approximation with respect to the reference image. Minimizing differential semblance is an attractive technique, too, because we are operating with small differences of nearby offsets which are likely within the Born approximation. However, this method is not directly guided toward the solution and aliasing or any remnants of coherent noise, like multiples or converted waves, can cause it to diverge.
Dealing with errors in automatic velocity analysis (ps.gz 6328K) (pdf 1179K) (src 23167K)
Clapp R. G.
The lack of human interaction in automatic reflection tomography leads to a larger percentage of ``bad'' data points. The number of data points associated with events with spurious moveouts (such as multiples and converted waves) can be minimized by intelligently controlling the semblance scanning range. The effect of the bad data points can be limited by replacing the standard L₂ norm solution with a norm closer to L₁ by reweighted least-squares. By replacing the standard constant $$ parameter with a diagonal operator, areas with large errors in moveout can be highly regularized with minimal effect on areas with more reliable moveout information. This methodology is applied to a complex 2-D dataset.
Velocity estimation for seismic data exhibiting focusing-effect AVO (part 2) (ps.gz 3347K) (pdf 1375K) (src 8201K)
Vlad I.
Vlad and Biondi (2002) have shown that focusing-effect AVO (FEAVO) exists and is visible in the angle domain. They have conjectured that wave-equation migration velocity analysis (WEMVA) might solve the FEAVO problem. I continue that line of work by showing evidence that WEMVA is the right tool, and by redoing the preprocessing and velocity analysis of the dataset on which FEAVO was defined.

Imaging

Transformation to dip-dependent Common Image Gathers (ps.gz 559K) (pdf 466K) (src 1713K)
Biondi B. and Symes W.
We introduce a new transform of offset-domain Common Image Gathers (CIGs) obtained by wavefield-continuation migration methods. This transformation can be applied to either horizontal-offset CIGs or vertical-offset CIGs. It overcomes the limitations that both kinds of CIGs suffer in the presence of a wide range of reflectors' dips. The result of our transformation is an image cube that is equivalent to the image cube that would have been computed if the offset direction were aligned along the apparent geological dip of each event. The proposed transformation applies a non-uniform dip-dependent stretching of the offset axis and can be efficiently performed in the Fourier domain. Because it is dependent on the image's apparent dip, the offset stretching automatically corrects for the image-point dispersal. Tests on a synthetic data set confirm the potential advantages of the transformation for migration velocity analysis of data containing steeply dipping reflectors.
Deconvolution imaging condition for reverse-time migration (ps.gz 637K) (pdf 473K) (src 14691K)
Valenciano A. A. and Biondi B.
The reverse-time migration imaging condition can be improved by computing the reflection strength at each subsurface point as the zero lag value of the deconvolution of the receiver wavefield by the source wavefield. I show that by using this approach it is possible to eliminate image artifacts due to wavefield multipathing through velocity anomalies. I also show that it has the advantage of handling better amplitudes during imaging.
Multicomponent data regularization (ps.gz 2954K) (pdf 2006K) (src 19689K)
Rosales D. and Biondi B.
Geometry regularization is a key process for obtaining reliable subsurface images with 3D seismic data. 3D regularization is, so far, a technique mainly used on PP land data. Multicomponent ocean bottom cable (OBC) technology simulates 3D land acquisition for multicomponent geophones at the ocean bottom. Reliable subsurface PS images, which provide amplitude information, have to go through a regularization process. Converted wave Azimuth Moveout (PSAMO) acts as a regularization operator in the formulation of the geometry regularization process in the least-squares sense.
Equivalence of source-receiver migration and shot-profile migration (ps.gz 4851K) (pdf 607K) (src 1321K)
Biondi B.
I first review the basic principles of shot-profile migration and source-receiver migration. Then I will show their equivalence.

 

• Shot-profile migration ...
  Damped imaging condition for reverse-time migration. (ps.gz 212K) (pdf 121K) (src 779K)
  Valenciano A. A.
  Reverse-time migration of shot-profiles Etgen (1986) has been proposed Biondi (2002) as an alternative to downward continuation methods to perform imaging in complex subsurface environments (e.g., under complex and rugose salt bodies). In these situations, due to the poor illumination given by near-offset ray-paths, all the events present in the data (such as overturned reflections and prismatic reflections) are needed to generate interpretable images. Not only the kinematic response of the migration is important, seismic data amplitudes also have the potential to provide information on reservoir properties. However, the most common implementation of shot-profile reverse-time migration uses the zero lag of the cross-correlation of the source and the receiver wavefields as imaging condition. This implementation has the advantage of being robust and honoring the kinematics of Claerbout's imaging principle Claerbout (1971) but does not honor the dynamics of the problem, resulting in the loss of amplitude accuracy. ...
  Multicomponent Stolt residual migration: a real data example (ps.gz 1584K) (pdf 387K) (src 20286K)
  Rosales D.
  Residual migration is the process of updating an image without re-running the entire migration process. Its main application is in refining of velocity models. Rosales et al. 2001 first introduced Stolt residual migration for converted waves (PS). Stolt migration is a constant velocity process. Although ...
  Migrating passive seismic data (ps.gz 346K) (pdf 104K) (src 2957K)
  Artman B.
  It is possible to migrate raw passive seismic data with a modified shot-profile migration algorithm to produce a subsurface image. This skips time-intensive and space-consuming pre-processing steps as has heretofore been assumed necessary. Further, output sections are better focused and more accurately imaged using less computer time.
  Deconvolving passive data (ps.gz 185K) (pdf 85K) (src 1247K)
  Artman B. and Claerbout J.
  Deconvolution prior to autocorrelation processing for a passive seismic data set has the potential to ameliorate wave-parameter and azimuthal inconsistency of arriving energy during acquisition. If any particular subset of plane-wave energy dominates the passive recording sequence, full illumination of the model-space may not be achieved. Further limitation of the result could also arise from the fact that the bulk of the ambient energy recorded in the experiment will likely be ground-roll energy that does not probe the subsurface. Thus, damping over-represented energy components by convolving the data with a prediction error filter (PEF) prior to ...
  

  Optimization

  Flattening without picking (ps.gz 1358K) (pdf 615K) (src 2597K)
  Lomask J. and Claerbout J.
  We introduce an analytical method for integrating dip information to flatten uninterpreted seismic data. First, dips are calculated over the entire seismic volume. The dip is then integrated in the Fourier domain, returning for each sample a time shift to a flat datum. Then each sample is shifted in the seismic data to remove all structural folding deformation in a single non-interpretive step. Using the Fourier domain makes it a quick process but requires that the boundaries are periodic. This method does not yet properly handle faults because of their discontinuous nature, but is presently very effective at removing warping and folding.
  Combined inversion: preconditioning with regularization (ps.gz 3118K) (pdf 535K) (src 1490K)
  Clapp M. L.
  Iterative inversion schemes are becoming more common in seismic processing. The high cost of the operators generally used in these inversion schemes makes it very important to minimize the number of iterations needed to obtain a good model. In complex environments, inversion schemes can be improved by styling the model through regularization or preconditioning. At early iterations, regularization provides a result that has a frequency content comparable to that of the ``ideal'' model. Preconditioning defines a solution at every model point at earlier iterations than regularization. An ``improved'' model should combine these two characteristics. This paper examines a scheme that uses the result of preconditioned inversion as an initial model for regularized inversion. I show that this scheme allows us to obtain an improved model in fewer iterations than would be needed for preconditioned inversion or regularized inversion alone.
  Conjugate gradient total least-squares in geophysical optimization problems (ps.gz 3341K) (pdf 667K) (src 1300K)
  Brown M.
  Total least-squares (TLS) optimization is a methodology to solve least-squares optimization problems when the modeling operator has errors. In standard least-squares optimization, errors are assumed to be concentrated in the data only. Golub and Loan (1980) presented a numerically-stable TLS algorithm which utilizes the singular value decomposition (SVD). Subsequent refinements to the method predominantly use SVD, and much of the current literature emphasizes stabilization of the inverse and implicit model regularization by SVD truncation Fierro et al. (1997). ...
  

  Non-linear problems

  Nonlinear pairwise alignment of seismic traces (ps.gz 2857K) (pdf 802K) (src 1513K)
  Liner C. L. and Clapp R. G.
  Alignment of seismic traces is a recurring need in seismic processing and interpretation. For global alignment via static shift there are robust tools available, including cross correlation. However, another kind of alignment problem arises in applications as diverse as associating synthetic seismograms to field data, harmonizing P-wave and mode converted data, final multilevel flattening of common image gathers, and so on. These cases require combinations of trace compression, extension, and shift - all of which are time variant. The difficulty is to find a mapping between the traces which is in some sense optimum. This problem is solved here using a modified form of the Needleman-Wunsch algorithm, which was originally developed for amino acid sequence alignment in proteins. Applied to seismic traces, this global optimization algorithm provides a nonlinear mapping of one seismic trace onto another. The method extends to alignment of any number of traces since that problem can be broken down into a cascade of pairwise alignments. The Needleman-Wunch algorithm is discussed, extended to the seismic case, and applied to field data. The results show a promising new tool for nonlinear alignment or flattening of seismic data.
  Simultaneous estimation of two slopes from seismic data, applied to signal/noise separation (ps.gz 2313K) (pdf 545K) (src 3638K)
  Brown M.
  I present an efficient new approach to simultaneously estimate two slopes from seismic data. I employ a Newton iteration to overcome the problem's nonlinearity. In spite of my method's theoretical inability to handle aliased data, it robustly estimates two independent slopes in many circumstances. I apply my method to the problem of signal/noise separation on synthetic and real data examples. The estimated slopes provide approximate inverse signal and noise covariance operators good enough to obtain an excellent separation, with only a limited amount of prior information required.
  Can we make genetic algorithms work in high-dimensionality problems? (ps.gz 167K) (pdf 513K) (src 2509K)
  Alvarez G.
  In this paper I compare the performance of a standard genetic algorithm versus a micro-genetic algorithm for matching a randomly-generated seismic trace to a reference trace with the same frequency spectrum. A micro genetic algorithm evolves a very small population that must be restarted whenever the population loses its genetic diversity. I show that the micro-genetic algorithm is more efficient in solving this problem in terms of improved rate of converge, especially in the first few generations. This characteristic may make the method useful for locating the most promising valleys in the search space which can then be searched with more traditional gradient-based methods. An additional benefit is a significant reduction in the number of evolution parameters that needs to be adjusted making the method more easy to use.
  Velocity inversion of a seismic trace with a micro-genetic algorithm (ps.gz 71K) (pdf 148K) (src 1958K)
  Alvarez G.
  An ever present goal of seismic processing and inversion is to extract meaningful geologic information from seismic data. The simplest possibility is to invert for seismic impedance, or, considering density constant, invert for velocity. Here I use a micro-genetic algorithm program to achieve this goal. A micro-genetic algorithm is different from a standard genetic algorithm in that it evolves a very small population that must be restarted whenever genetic diversity is lost. I use a real sonic log to compute a synthetic seismic trace and then use that trace as input to the micro-genetic algorithm program to invert for the velocity log. Without further input, the velocity inversion cannot hope to recover the general velocity-depth trend because this information corresponds to the very low frequencies which are absent in the seismic data. In order to achieve a good match between the real and the inverted sonic log, the inversion must be suppled with an estimate of the velocity-depth trend. Here I use a 33-point Savitzky-Golay filter with a sixth-degree polynomial to smooth the velocity log and show that with this trend the velocity log is reasonably well recovered. I used both the L₁ and the L₂ norms of the sample-to-sample difference between the reference and the inverted trace and show that the results with both norms are similar.
  

  Rock Physics

  Modeling high-frequency acoustics velocities in patchy and partially saturated porous rock using differential effective medium theory (ps.gz 49K) (pdf 70K) (src 256K)
  Berryman J. G.
  Differential effective medium (DEM) theory is applied here to the problem of modeling physical properties of poroelastic media that are partially saturated with liquid. Typical fluid saturants are air and water, or gas and oil. If the liquid and gas saturants are homogeneously mixed, then we say the medium is partially saturated. If the liquid and gas saturants are poorly mixed, so each constituent occupies separate, but contiguous, regions of the porous medium, we say the medium has patchy saturation. Some examples are presented to show that a reasonable approach to modeling the effects of patchy saturation at high frequencies (200 kHz and above) is produced by treating the medium as if it is a composite of gas-saturated and liquid-saturated porous inclusions that are homogeneously mixed together. Estimates of the properties are obtained using differential effective medium theory. The results differ dramatically from those predicted by Gassmann's equations for homogeneous mixing of the fluids in individual pores. In particular, the shear modulus depends on the elastic properties of the fluid constituents, unlike the quasi-static behavior predicted by Gassmann.
  Dispersion of waves in porous cylinders with patchy saturation (ps.gz 101K) (pdf 153K) (src 247K)
  Berryman J. G. and Pride S.
  Laboratory experiments on wave propagation through saturated and partially saturated porous media have often been conducted on porous cylinders that were initially fully saturated and then allowed to dry while continuing to acquire data on the wave behavior. Since it is known that drying typically progresses from outside to inside, a sensible physical model of this process is concentric cylinders having different saturation levels - the simplest example being a fully dry outer cylindrical shell together with a fully wet inner cylinder. We use this model to formulate the equations for wave dispersion in porous cylinders for patchy saturation (i.e., drainage) conditions. In addition to multiple modes of propagation obtained numerically from these dispersion relations, we find two distinct analytical expressions for torsional wave modes. We solve the dispersion relation for torsional waves for two examples: Massillon sandstone and Sierra White granite. The drainage analysis appears to give improved agreement with the data for both these materials.
  Exploring the relationship between uncertainty of AVA attributes and rock information (ps.gz 172K) (pdf 283K) (src 7871K)
  Chen W. and Clapp R.
  Amplitude versus Angle (AVA) attributes include information about rock properties. Using a dataset from South America, we performed a multiple realization method to get multiple equal-probable AVA intercepts, gradients, and their products. We generated a 3-D histogram to evaluate the variability of those AVA attributes. In the same area, we chose a 2-D section by matching it to three wells. Then we calculated the shale volume along these three wells and found the well with low shale volume has high AVA uncertainty, which made us guess the low shale/sand ratio may cause high AVA uncertainty. The further work need to be done is to use more real data to exam our conjecture, namely, whether there exist an empirical relationship between AVA uncertainty and rock information, such as shale volume, impedance or velocity.
  

   

  Ray based tomography using residual Stolt migration

  Robert G. Clapp

  ABSTRACT

  In complex areas, residual vertical movement is not an effective method to calculate traveltime errors for image domain tomography. By scanning over velocity ratios using residual Stolt migration, a different criteria for a coherent image can be defined, and a traveltime error approximated. The resulting traveltime errors are more accurate, therefore the tomography procedure is more robust than the more traditional methodology. Results are shown on a complex 2-D dataset.

  Depth migration is often necessary for complex structures, but requires an accurate interval velocity model. Estimating this velocity model is one of the essential problems in reflection seismology. One of the most common methods to estimate an interval velocity model is from ray-based reflection tomography after migration (), which I will refer to as image domain tomography. In previous work (, ), I showed how to use angle gathers (, ) in conjunction with downward continuation based migration (, ) for back-projection. Post-migration reflection tomography is based on the fact that an offset gather (Kirchhoff) or angle gather (wavefield continuation) should be flat after migration. Deviation from flatness indicates a velocity error and can be converted into a travel time error and back-projected. (), among many others, pointed out that, in complex environments, looking at vertical moveout in gathers is not the optimal method to describe moveout errors. () presented one alternative approach, constructing gathers where moveout is normal to an event.

  Another approach is to use residual migration (, , ) to find the best focusing velocity. () and (, ) showed how to do residual migration for wave continuation methods. These methods allow scanning over slowness field ratios. () noticed that when scaling the slowness field by a constant, ray behavior is unchanged. () described a method of updating the velocity model by back projecting along the normal ray.

  In this paper I take these works a step further. After performing downward continuation based migration, I use residual migration to find a smooth field of $\gamma$ values that best focuses the data. I then convert these $\gamma$ values to approximate travel-time errors and back-project.

  I begin by outlining a method of selecting my $\gamma$ field. I then describe the approximations used to convert $\gamma$ to $\bf \Delta t$. Finally, I show the procedure applied to a 2-D North Sea example.