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.

"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 -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 -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 -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.

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.

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.

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.

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.

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
...

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
...

10/14/1997