3-D imaging and velocity analysis methods are computationally intensive.
One reason is that they are iterative processes, which need a
forward modeling algorithm as a kernel component. For example, 3-D
Kirchhoff-type imaging usually needs a 3-D Green's function as input.
There are many different ways of estimating Green's function. They can be
roughly categorized into two families according to both accuracy and
efficiency: *limited bandwidth methods* and *infinite frequency
methods* Audebert et al. (1994).

Limited bandwidth methods usually deal with frequency-dependent, complex-valued Green's function. These methods are typically the extrapolation or modeling algorithms working in the time domain or the temporal frequency domain, e.g., phase-shift or finite-difference downward continuation. Nichols 1994 presented a method of calculating band-limited traveltime explicitly, which is essentially a compromise between full wavefield modeling and ray tracing. In terms of accuracy, it is a reduced version of a full wavefield modeling; in terms of efficiency, the cost of this new method is lower than that of full wavefield modeling, but higher than ray tracing method. One main weakness with this method is that no explicit raypaths are calculated and stored. In many applications, e.g., tomographic velocity analysis, raypaths are indeed necessary.

Ray tracing method is the main sub-family of the infinite frequency methods.
The simplest and also the fastest ray tracing algorithm is *geometric ray
tracing*. It transfers the eikonal equation to the ray equation using the
characteristic method and then estimates the raypaths by solving the ray
equation.
Geometric ray tracing can only provide traveltime and
cannot provide any amplitude and phase-shift information. It has been proved
that first arrivals cannot be used to image complex structures correctly
Rekdal and Biondi (1994). *Dynamic ray tracing* Cervený (1985), is
an extension to geometric ray tracing.
Dynamic ray tracing (including paraxial ray approximation and Gaussian Beam
approximation) can provide not only multiple traveltimes but also amplitude
and phase-shift information. Therefore it can also be regarded as a
compromise between full wavefield modeling and ray tracing modeling.

In this paper, on the basis of Cervený's work, we develop an accurate and robust 3-D dynamic ray tracing system. First, our implementation is independent of the model gridding scheme, which can be layer, cube, or tetrahedron, etc. Therefore, it is easy to incorporate our implementation into different 3-D velocity model representations. Second, since dynamic ray tracing is a high-frequency approximation, it has some requirements over the smoothness of the model. In order to make the model more ``ray valid'', we apply a triangle smoothing filter to the model along the x, y, z axis. The smoothed model makes not only our numerical solver more stable, but also the results from dynamic ray tracing more like band-limited seismic data. We use a Fortran90 version Runge-Kutta solver in our application Brankin and Gladwell (1995). This solver has the ability to adapt its step-length in accordance with the local gradient of the slowness, which makes it more accurate than constant step-length solver and also enables us to deal with more complex structures without smoothing it dramatically. Another advantage is that Fortran90 is easily parallelized and we expect our implementation would be more efficient when implemented on a parallel computer.

In the next several sections, we first discuss our representation of the model and then present a brief version of the dynamic ray tracing theory. After that, we show three results from applying our dynamic ray tracing program to synthetic models. The linear model verifies that our implementation is correct when there is a theoretical solution. The layer model shows the importance of smoothing the model. The saltdome model demonstrates our method's robustness.

11/11/1997