SEP-95 -- TABLE OF CONTENTS

Helix Transform

Multidimensional recursive filters via a helix (ps 96K) (src 148K)
Claerbout J.
A wire is coiled into a helix whose surface is a cylinder. I show that a filter on the 1-D space of the wire mimics a 2-D filter on the cylindrical surface. Thus 2-D convolution can be done with a 1-D convolution program. Likewise, I show some curious examples of 2-D recursive filtering (also called 2-D deconvolution or 2-D polynomial division). In 2-D as in 1-D, the computational advantage of recursive filters is the speed with which they propagate information over long distances. We can estimate 2-D prediction-error filters (PEFs), that are assured of being stable for 2-D recursion. Such 2-D and 3-D recursions vastly speed the solution of a wide class of geophysical estimation problems. Simple tasks that are vastly speeded are (1) estimating missing values on a multidimensional cartesian mesh, and (2) distributing irregularly positioned data on a regular mesh.
Missing data interpolation by recursive filter preconditioning (ps 257K) (src 1450K)
Fomel S., Clapp R., and Claerbout J.
Missing data interpolation problems can be conveniently preconditioned by recursive inverse filtering. A helix transform allows us to implement this idea in the multidimensional case. We show with examples that helix preconditioning can give a magnitude-order speedup in comparison with the older methods.
Solution steering with space-variant filters (ps 579K) (src 3326K)
Clapp R. G., Fomel S., and Claerbout J.
Most geophysical problem require some type of regularization. Unfortunately most regularization schemes produce ``smeared'' results that are often undesirable when applying other criteria (such as geologic feasibility). By forming regularization operators in terms of recursive steering filters, built from a priori information sources, we can efficiently guide the solution towards a more appealing form. The steering methodology proves effective in interpolating low frequency functions, such as velocity, but performs poorly when encountering multiple dips and high frequency data. Preliminary results using steering filters for regularization in tomography problems are encouraging.
Exploring three-dimensional implicit wavefield extrapolation with the helix transform (ps 366K) (src 3296K)
Fomel S. and Claerbout J. F.
Implicit extrapolation is an efficient and unconditionally stable method of wavefield continuation. Unfortunately, implicit wave extrapolation in three dimensions requires an expensive solution of a large system of linear equations. However, by mapping the computational domain into one dimension via the helix transform, we show that the matrix inversion problem can be recast in terms of an efficient recursive filtering. Apart from the boundary conditions, the solution is exact in the case of constant coefficients (that is, a laterally homogeneous velocity.) We illustrate this fact with an example of three-dimensional velocity continuation and discuss possible ways of attacking the problem of lateral variations.

Traveltime Modeling

"Focusing" eikonal equation and global tomography (ps 594K) (src 1233K)
Biondi B., Fomel S., and Alkhalifah T.
The transformation of the eikonal equation from depth coordinates (z,x) into vertical-traveltime coordinates ($,$)enables the computation of reflections traveltimes independent of depth-mapping. This separation allows the focusing and mapping steps to be performed sequentially even in the presence of complex velocity functions, that otherwise would ``require'' depth migration. The traveltimes satisfying the transformed eikonal equation can be numerically evaluated by solving the associated ray tracing equations. The application of Fermat's principle leads to the expression of linear relationships between perturbations in traveltimes and perturbations in focusing velocity. This linearization, in conjunction with ray tracing, can be used for a tomographic estimation of focusing velocity.
Time-domain anisotropic processing in arbitrarily inhomogeneous media (ps 233K) (src 254K)
Alkhalifah T., Fomel S., and Biondi B.
In transversely isotropic media with a vertical axis of symmetry (VTI media), we can represent the image in vertical time, as opposed to depth, thus eliminating the inherent ambiguity of resolving the vertical P-wave velocity from surface seismic data. In this new $(x-)$-domain, the raytracing and eikonal equations are completely independent of the vertical P-wave velocity, on the condition that the ratio of the vertical to normal-moveout (NMO) P-wave velocity ($$) is laterally invariant. Practical size departures of $$ from lateral homogeneity affect traveltimes only slightly. As a result, for all practical purposes, the VTI equations in the $(x-)$-domain become dependent on only two parameters in laterally inhomogeneous media: the NMO velocity for a horizontal reflector, and an anisotropy parameter, $$. An acoustic wave equation in the $(x-)$-domain is also independent of the vertical velocity. It includes an unsymmetric Laplacian operator to accommodate the unbalanced axis units in this new domain. In summary, we have established the basis for a full inhomogeneous time-processing scheme in VTI media that is dependent on only v and $$, and independent of the vertical P-wave velocity.
Huygens wavefront tracing: A robust alternative to conventional ray tracing (ps 592K) (src 2097K)
Sava P. and Fomel S.
We present a method of ray tracing that is based on a system of differential equations equivalent to the eikonal equation, but formulated in the ray coordinate system. We use a first-order discretization scheme that is interpreted very simply in terms of the Huygens' principle. The method has proved to be a robust alternative to conventional ray tracing, while being faster and having a better ability to penetrate the shadow zones.
Multivalued traveltime interpolation (ps 322K) (src 811K)
Sava P. and Biondi B.
We present a method of interpolating multiple arrival traveltimes and amplitudes. It is based on an assumption of the physical continuity of the traveltimes and provides optimal interpolation between adjacent rays using constrained Delaunay triangulation. We use a graphical, interactive program to visualize the data and for quality control. We have developed the algorithms for both the 2-D and 3-D cases, though implemented only the 2-D so far. The method has been tested on several synthetic velocity models and on the SEG-EAGE salt model.
A variational formulation of the fast marching eikonal solver (ps 594K) (src 3683K)
Fomel S.
I exploit the theoretical link between the eikonal equation and Fermat's principle to derive a variational interpretation of the recently developed method for fast traveltime computations. This method, known as fast marching, possesses remarkable computational properties. Based originally on the eikonal equation, it can be derived equally well from Fermat's principle. The new variational formulation has two important applications: First, the method can be extended naturally for traveltime computation on unstructured (triangulated) grids. Second, it can be generalized to handle other Hamilton-type equations through their correspondence with variational principles.
Implementing the fast marching eikonal solver: Spherical versus Cartesian coordinates (ps 1704K) (src 2009K)
Alkhalifah T. and Fomel S.
Spherical coordinates are a natural orthogonal system to describe wavefronts emanating from a point source. While a regular grid distribution in the Cartesian coordinate system tends to undersample the wavefront description near the source (the highest wavefront curvature) and oversample it away from the source, spherical coordinates, in general, provide a more balanced grid distribution to characterize such wavefronts. Our numerical implementation confirms that the recently introduced fast marching algorithm is both a highly efficient and an unconditionally stable eikonal solver. However, its first-order approximation of traveltime derivatives can induce relatively large traveltime errors for waves propagating in a diagonal direction with respect to the coordinate system. Examples, including the infamous Marmousi model, show that a spherical coordinate implementation of the method results in far fewer errors in traveltime calculation than the conventional Cartesian coordinate implementation, and with practically no loss in computational advantages.

Reservoir Characterization and Interpretation

Estimating the amount of gas hydrate and free gas from surface seismic (ps 1414K) (src 3271K)
Ecker C., Dvorkin J., and Nur A.
In this study we provide a theoretical tool for quantifying the amount of gas hydrate and gas near a bottom simulating reflector (BSR) at the Blake Outer Ridge from surface seismic. We develop rock-physics models that link the elastic wave velocity in high-porosity marine sediments to density, porosity, effective pressure, mineralogy, and water/gas and hydrate saturation of the pore space. Three models of hydrate deposition are examined: (1) hydrate is part of the pore fluid; (2) hydrate becomes part of the solid frame, thus reducing porosity and weakly affecting the stiffness of the sediment; and (3) hydrate cements grain contacts and therefore strongly reinforces the sediments. Using interval velocities obtained from velocity analysis together with the rock-physics models, we obtain maps of lateral hydrate and gas saturation. Model (1) predicts maximum hydrate saturations between 19% and 33%, model (2) saturations between 16% and 25% and model (3) saturations less than 1%. Maximum gas saturation is about 2% of the pore space. These results are consistent with those that can be obtained using known well-log velocities and porosities from this region. Subsequently, in order to evaluate the effect of the different models on the actual seismic amplitudes, we generated synthetic seismograms using Kirchhoff modeling. The AVO responses showed that models (1) and (2) cannot be differentiated by surface seismic. Comparison with real AVO data from the Blake Outer Ridge suggests that only model (1) or (2) can reproduce the actual observed amplitude trends. Therefore, we conclude that the hydrated sediments at the Blake Outer Ridge are only weakly, if at all, stiffened by the presence of hydrate, which can occupy up to 33% of the pore space.
An amplitude bias correction for 4D seismic cross-equalization (ps 1365K) (src 3915K)
Rickett J., Lumley D., and Martin H.
Amplitude balancing is an important part of the 4D seismic cross-equalization process, which aims to suppress processing and acquisition differences between time-lapse 3D seismic surveys. The matched-filter approach to cross-equalization estimates an amplitude correction that minimizes the norm of the difference section. However, this correction will be biased by the presence of noise, and the amplitude of coherent events will not be equalized correctly. Similarly, a normalization scheme that is based on equalizing the energy in the two surveys, implicitly assumes that their random noise levels are equal. We illustrate the problem with a synthetic example, and present a method that correctly scales the amplitudes based on the relative signal-to-noise levels of each dataset.
Seismic Anisotropy in Trinidad: More parameter estimation   (ps 2705K) (src 236K)
Alkhalifah T. and Rampton D.
New estimates of anisotropy have revealed more details of the subsurface in Trinidad. These estimates are obtained by including more realistic constraints on the anisotropic inversion, which eventually helped boost the stability of the process. The anisotropic parameter $$, which, if not zero, implies the existence of anisotropy, is used to discriminate conservatively between shales and sands. The underlying theory is that shales induce anisotropy, positive $$ in particular, and sands do not. The $$ estimates, through have nice lateral correlation, react to the presence of faults. Correlation of these $$ results with gamma-ray well-log measurements used as a shale estimate proves the credibility of the results. This finding confirms the hypothesis that anisotropy is caused by shales in the subsurface, and, consequently, we can use the inversion for interval $$ to estimate lithology.
Least squares dip and coherency attributes (ps 815K) (src 916K)
Bednar J. B.
This paper demonstates the effectiveness of local least-squares-dip estimates as coherency attributes in 3-D seismic data volumes. The process is based on the recognition that the magnitude of local slowness estimates, calculated from spatial and time derivatives, are good edge detectors and as a consequence are correlated with events which might also be detected by more sophisticated statistical techniques.

Imaging

On the inversion of 3D multichannel data (ps 106K) (src 149K)
Chemingui N. and Biondi B.
Spatial sampling is an important consideration in the design of seismic surveys. It affects acquisition costs and is directly tied to processing and imaging requirements. During the acquisition stage, economic constraints, obstructions, cable feathering, environmental objectives, and many other factors cause seismic data to be sampled irregularly. Therefore, 3D surveys typically have sparse and irregular geometry that often results in spatial aliasing. ...
Results in depth focusing analysis for 3-D migration velocity estimation (ps 1292K) (src 19834K)
Malcotti H. and Biondi B.
This paper is the continuation of a previous paper (Malcotti and Biondi, 1997) in which we discussed the theoretical and practical bases of depth focusing analysis for 2-D and 3-D. In this work, we present some results of the depth focusing analysis methodology based on a 3-D prestack common-azimuth migration (Biondi et. al., 1996). We apply this methodology to two differents data sets with complex velocities, a synthetic data set and a real 3-D data set. In both data sets, we show one iteration in the methodology loop for interval velocity estimation using depth focusing analysis. For the real seismic data set, our goal is to update the velocity model where the image has focusing problems. Therefore, we used the best velocity field available. For the interpretation of the depth error gathers, we use a three dimensional graphic tool that help to simplify the picking in complex error gathers.
A theoretical comparison of equivalent offset migration and dip moveout prestack imaging (ps 66K) (src 66K)
Bednar J. B.
Equivalent-offset migration is a methodology for prestack-Kirchhoff time migration that partially reverses the order of velocity analysis, normal moveout correction, stack, and migration. Although claimed to be computationally and analytically superior to earlier time-domain approaches, it is not independent of velocity. Velocity-independent dip-moveout followed by prestack imaging is similar to equivalent-offset migration in that it also postpones normal moveout correction to the post-migration stage, but it is independent of velocity. In this paper, I investigate the theoretical relationship between these two processes, showing that equivalent-offset migration and prestack imaging are asymptotically equivalent. Moreover, equivalent-offset migration has little, if any, computational or analytic advantage. The fact that the combination of dip-moveout followed by prestack imaging is independent of velocity is a major advantage and suggests that this latter method is better suited to problems of velocity estimation.

Anisotropic Modeling

An anisotropic Marmousi model (ps 2842K) (src 2715K)
Alkhalifah T.
I use the acoustic wave equation for transversely isotropic media with vertical symmetry axis (VTI media), to generate synthetic VTI data for an anisotropic version of the Marmousi model. This acoustic equation, though it represents a physically impossible medium, provides an extremely accurate approximation of the widely used elastic wavefield. The anisotropic Marmousi model has the same NMO velocity as the original Marmousi model and an anisotropy $$ distribution that possesses the same layering characteristics as the velocity model. Interval $$ spans the range of 0 to 0.27, which are values that are commonly observed in practice. The traveltime and amplitude differences between the synthetic seismograms of this new anisotropic model and those produced by the original isotropic Marmousi data set are quite apparent. As a result, a prestack isotropic migration failed to properly image the anisotropic data when using the exact original velocity model.
An acoustic wave equation for anisotropic media (ps 513K) (src 863K)
Alkhalifah T.
A wave equation, derived under the acoustic medium assumption for P-waves in transversely isotropic media with a vertical symmetry axis (VTI media), though physically impossible, yields good kinematic approximation to the familiar elastic wave equation for VTI media. The VTI acoustic wave equation is fourth-order and has two sets of complex conjugate solutions. One set of solutions is just perturbations of the familiar acoustic wavefield solutions in isotropic media for incoming and outgoing waves. The second set describes a wave type that propagates at speeds slower than the P-wave for the positive anisotropy parameter, $$, and grows exponentially, becoming unstable, for negative values of $$. Luckily, most $$values corresponding to anisotropies in the subsurface have positive values. Placing the source in an isotropic layer, a common occurrence in marine surveys where the water layer is isotropic, eliminates most of the energy of this additional wave type. Numerical examples prove the usefulness of this acoustic equation in simulating wave propagation in complex models. From this acoustic wave equation, the eikonal and transport equations that describe the ray theoretical aspects of wave propagation are derived. These equations, based on the acoustic assumption (shear wave velocity equals zero), are much simpler than their elastic counterparts, and yet yield exceptionally accurate description of traveltime and geometrical amplitude, or wavefront spreading.

Inversion

Stereology as inverse problem (ps 89K) (src 59K)
Berryman J. G.
Stereology is the part of imaging science in which the three-dimensional structure of a body is determined from two-dimensional views. Although it is relatively easy to determine volume information from 2-D slices, it is nontrivial in general to determine other physical properties such as internal surface areas unless the medium is known to have some simple symmetry such as isotropy. For this reason, stereology can be viewed as a type of inverse problem. In earlier work I showed that an anisotropic spatial correlation function of a random porous medium could be used to compute the specific surface area when it is stationary as well as anisotropic by first performing a three-dimensional radial average and then taking the first derivative with respect to lag at the origin. This result generalized the earlier result for isotropic porous media of Debye et al. (1957). Here I provide more detailed information about the use of spatial correlation functions for anisotropic porous media and in particular I show that, for stationary anisotropic media, the specific surface area can be related to the derivative of the two-dimensional radial average of the correlation function measured from cross sections taken through the anisotropic medium. The main concept is first illustrated using a simple pedagogical example for an anisotropic distribution of spherical voids. Then, a general derivation of formulas relating the derivative of the planar correlation functions to surface integrals is presented. When the surface normal is uniformly distributed (as is the case for any distribution of spherical voids), my formulas can be used to relate specific surface area to easily measureable quantities from any single cross section. When the surface normal is not distributed uniformly (as would be the case for an oriented distribution of ellipsoidal voids), my results show how to obtain valid estimates of specific surface area by averaging measurements on three orthogonal cross sections. One important general observation for porous media is that the surface area from nearly flat cracks may be underestimated from measurements on orthogonal cross sections if any of the cross sections happen to lie in the plane of the cracks. This result is illustrated by taking the very small aspect ratio (penny-shaped crack) limit of an oblate spheroid, but holds for other types of flat surfaces as well.
Comparison of autoregressive and multitaper spectral analysis for long time series (ps 685K) (src 2296K)
Fodor I. K., Berryman J. G., and Stark P. B.
The periodogram (i.e., the square of the Fourier transform) of a time series generally provides a poor estimate of the spectrum if that spectrum has a wide dynamic range. So the spectrum of any process that includes either one or many resonant modes (sharp peaks) can be expected to be poorly computed by such elementary means. Multitaper spectral analysis is a nonparametric method designed to provide a rigorous method of resolving the spectrum of such complex processes. There are some practical difficulties with this method, such as deciding what tapers to use and how many, that can make the method some what difficult for the uninitiated user. Another approach to spectral analysis is the parametric method known as autoregressive (AR) analysis (related but not identical to maximum entropy spectral analysis). AR analysis can be used to approximate the dominant modes in the spectrum, remove them from the time series and thus prewhiten the remaining time series, thereby eliminating the main problem with analysis based on the periodogram. Furthermore, if the main purpose of the spectral analysis is to determine the number and location of the modes, the AR method provides a set of parameters that can be used to describe the mode structure if desired. The present paper provides a set of examples comparing the use of both the multitaper method and the autoregressive method for analyzing the same (long) time series data. We find that both methods give very comparable results in these examples, thus providing a type of empirical cross-validation showing that both methods are therefore doing an excellent job of estimating the true spectrum of the time series.

Noise Removal

Multi-source experiment for ground roll removal (ps 415K) (src 2241K)
Crawley S.
Shallow seismic data are typically recorded with multiple source impacts per shot gather, with individual impacts stacked in the recorder to economize storage. Don Steeples has observed that in some cases, separately recorded impacts at the same shot point share similar spectra at ground roll frequencies, but differ at higher signal frequencies. This observation suggests that multiple source gathers at a single station can be subtracted to remove ground roll ...
Upward continuation multiple suppression applied to the SEG multiples datasets (ps 13716K) (src 34144K)
Crawley S.
Two SEG multiple suppression workshop data sets were used as input to a multiple suppression algorithm based on upward continuation, originally formulated by Berryhill and Kim 1986. The recorded data are upward continued a distance equal to twice the water depth, with water velocity, in order to add a lag equal to the multiple period. A primary in the upward continued data should then line up with its first multiple in the original recorded data, a first order multiple in ...

Ratfor90

RATional FORtran == Ratfor90 (ps 54K) (src 57K)
Clapp R. G. and Claerbout J. F.
Fortran is generally accepted as the most universal computer language for computational physics. However, for general programming, it has been surpassed by C. Ratfor is Fortran with C-like syntax. Ratfor was invented by the Kernighan and Plauger1976, the same people who invented C. Ratfor uses C-like syntax, the syntax that is also found in the popular languages ...