Evaluating the Stolt stretch parameter

Sergey Fomel

sergey@sep.stanford.edu

ABSTRACT

The Stolt migration extension to a varying velocity case (Stolt stretch) implies describing a vertical heterogeneity by a constant parameter (W). This paper exploits the connection between modified dispersion relations and traveltime approximations to derive an explicit expression for W. The expression provides theoretically the highest possible accuracy within the Stolt stretch framework. Applications considered include optimal partitioning of the velocity distribution for the cascaded migrations and extension of the Stolt stretch method to transversally isotropic models.

Stolt migration is regarded as the fastest post-stack migration method of all the known algorithms. A known price for that speed is the constant velocity assumption. The time-stretching trick proposed in Stolt's classic paper provides an approximate extension of the method to a variable velocity case. Stolt stretch implicitly transforms reflection traveltime curves to fit an approximate constant velocity pattern (, , ). In other words, the wave equation with variable velocity is transformed by a particular stretch of the time axis to an approximate differential equation with constant coefficients. The two constant coefficients are an arbitrarily chosen frame velocity and a specific nondimensional parameter (W in Stolt's original notation). In the constant velocity case W is equal to 1, and the transformed equation coincides with the exact constant velocity wave equation. In variable velocity media, W is generally assumed to lie between and 1. As shown by Beasley et al. , the cascaded f-k migration approach can move the value of W for each migration in a cascade closer to 1, thus increasing the accuracy of the Stolt stretch approximation.

The W factor was defined by Stolt as an approximate average of a complicated function. Stolt's definition cannot be used directly for computation because it includes a combined dependence on both time and space coordinates. Therefore, in practice, the estimation of this factor is always replaced by a heuristic guess. That's why Levin called the W parameter ``infamous'' (joking, of course), and Beasley et al. called it it ``esoteric.''

This paper develops a method to evaluate the Stolt stretch parameter explicitly. The main idea is to constrain the parameter by fitting the exact and approximated traveltime functions. In the case of isotropic interpretation, the W parameter is connected to the ``parameter of heterogeneity'' (, , ). In the case of anisotropic (transversally isotropic) interpretation, it can be related to the ``parameter of anellipticity'' (, ).

STOLT STRETCH THEORY In order to simplify the references, I will begin with the textbook definitions of the Stolt migration method. The reader familiar with Stolt stretch theory can skip this section and go on to a new piece of theory in the next one.

Post-stack seismic migration is theoretically a two-stage process consisting of wavefield downward continuation in depth z based on the wave equation

$${\partial^2 P \over \partial x^2} + {\partial^2 P \over \par... ...^2(x,z)}}\, {\partial^2 P \over \partial t^2} \EQNLABEL{wave}$$ (1)

and the imaging condition t=0 (here the velocity v is twice as small as the actual wave velocity). Stolt time migration performs both stages in one step, applying the frequency-domain operator

$$\tilde{P_0}\left(k_x, \omega_0\right)= \tilde{P_v}\left(k_x,... ...mega_0\right)}\over{d\omega_0}}\right\vert\;, \EQNLABEL{stolt}$$ (2)

where $P_0\left(x, t_0\right)$ stands for the initial zero-offset (stacked) seismic section defined on the surface z=0, $P_v\left(x, t_v\right)$ is the time-migrated section, and t_v is the vertical traveltime

$$t_v=\int_{0}^{z}{dz \over {v(x,z)}}\;\;. \EQNLABEL{tv}$$ (3)

The function $\omega_v \left(k,\omega_0 \right)$ in stolt corresponds to the dispersion relation of the wave equation wave and in the constant velocity case has the explicit expression

$$\omega_v \left(k,\omega_0 \right)=\mbox{sign}\left(\omega_0\right) \sqrt{\omega_0^2 - v^2 k^2}\;\;. \EQNLABEL{dispersion}$$ (4)

The choice of the sign in dispersion is essential to distinguish between upgoing and downgoing waves. It is the upgoing part of the wave field that is used in migration.

For the case of a varying velocity Stolt suggested the following change of the time variable (referred to in the literature as Stolt stretch):

$$s(t)={\left({{2 \over v_0^2}\,\int_0^t\eta d \tau}\right)}^{1/2}\;, \EQNLABEL{ss}$$ (5)

where v₀ is an arbitrarily chosen constant velocity, and $\eta$ is a function defined by the parametric expressions

$$\eta(\zeta)=\int_0^{\zeta} v(x,z) \,dz \;,\; \tau(\zeta)=\int_0^{\zeta} { {dz} \over {v(x,z)}}\;\;.$$ (6)

With the stretch ss, seismic time migration can be related to the transformed wave equation

$${\partial^2 P \over \partial x^2} + W\,{\partial^2 P \over \... ..., {\partial^2 P \over \partial \hat{t}^2}\;. \EQNLABEL{swave}$$ (7)

Here $\hat{z}$ and $\hat{t}$ are the transformed depth and time coordinates that possess the following property: if $\hat{z}=0$, $\hat{t}=s\left(t_0\right)$, and if $\hat{t}=0$, $\hat{z}=v_0 s\left(t_v\right)$. W is a varying coefficient defined as

$$W=a^2+2b\,(1-a^2)\;, \EQNLABEL{wstolt}$$ (8)

where Stolt's idea was to replace the slowly varying parameter W with its average value. Thus equation swave is approximated by an equation with constant coefficients, which has the dispersion relation

$$\widehat{\omega}_v\left(k,\widehat{\omega}_0 \right)= \left(... ...{\widehat{\omega}_0^2 - W v_0^2 k^2}\;. \EQNLABEL{sdispersion}$$ (9)

Stolt's approximate method for migration in heterogeneous media consists of the following steps:

1.

stretching the time variable according to ss,

2.

interpolating the stretched time to a regular grid,

3.

double Fourier transform,

4.

f-k time migration by operator stolt with the dispersion relation sdispersion,

5.

inverse Fourier transform,

6.

inverse stretching (shrinking) the vertical time variable on the migrated section.

The value of W must be chosen prior to migration. According to Stolt's original definition wstolt, the depth variable z gradually changes in the migration process from zero to $\zeta$, causing the coefficient b in wstolt to change monotonically from 0 to 1. If the velocity v monotonically increases with depth, then $\eta''(z)={\partial v \over \partial z}\geq 0$, and the average value of b is

$$\bar{b}={1 \over {\zeta \eta(\zeta)}}\, {\int_0^\zeta \eta(z... ...\zeta) {z \over \zeta} }dz}= {1 \over 2}\;. \EQNLABEL{bishalf}$$ (10)

As follows from wstolt and bishalf, in the case of monotonically increasing velocity, the average value of W has to be less than 1 (W equals 1 in a constant velocity case). Analogously, in the case of a monotonically decreasing velocity, W is always greater than 1. In practice, W is included in migration routines as a user-defined parameter, and its value is usually chosen to be somewhere in the range of 1/2 to 1.

In this paper I will describe a straightforward way to determine the most appropriate value of W for a given velocity distribution.

A useful tool for that purpose is Stewart Levin's formula for the traveltime curve. Levin applied the stationary phase technique to the dispersion relation sdispersion to obtain an explicit formula for the summation curve of the integral migration operator analogous to the Stolt stretch migration. The formula evaluates the summation path in the stretched coordinate system, as follows:

$$s\left(t_0\right)= \left(1-{1\over W}\right) s\left(t_v\righ... ...\, {{\left(x-x_0\right)^2} \over v_0^2}}}\;. \EQNLABEL{levin}$$ (11)

Here x₀ is the midpoint location on a zero-offset seismic section, and x is the space coordinate on the migrated section. Formula levin shows that, with the stretch of the time coordinate, the summation curve has the shape of a hyperbola with the apex at $\left\{x,s\left(t_v\right)\right\}$ and the center (the intersection of the asymptotes) at $\left\{x,{\left(1-{1\over W}\right)}\,s\left(t_v\right)\right\}$. In the case of homogeneous media, W=1, $s(t)\equiv t$, and levin reduces to the well-known hyperbolic diffraction traveltime curve. It is interesting to note that inverting formula levin for $s\left(t_v\right)$ determines the impulse response of the migration operator, which can be interpreted as the wavefront from a point source in the $\{x,\hat{z},\hat{t}\}$ domain of equation swave:

$$\hat{z}-\hat{z_0}= \left({1\over Q}-1\right) R \pm {1\over Q... ...\sqrt{R^2 - {Q\,{\left(x-x_0\right)^2}}}\;, \EQNLABEL{sfront}$$ (12)

where $R=v_0 \hat{t}$, and Q=2-W. According to equation sfront, wavefronts from a point source in the stretched coordinates for W<2 have an elliptic shape, with the center of the ellipse at $\{x,\hat{z_0}+\left( {1 \over Q}-1\right)\, R \}$ and the semi-axes $a_x={R \over \sqrt{Q}}$ and $a_z={r \over Q}$. The ellipses stretch differently for W<1 and W>1 (Figure stofro). In the upper part that corresponds to the upgoing waves, they look nearly spherical, since the radius of the front curvature at the top apex equals the distance from the source.

 

stofro
Figure 1 Wavefronts from a point source in the stretched coordinate system. Left: velocity decreases with depth (W=1.5). Right: velocity increases with depth (W=0.5).

EVALUATING THE W PARAMETER AND STOLT STRETCH ACCURACY

Formula levin belongs to the three-parameter class of traveltime approximations. The key result of this paper uses a remarkable formal similarity between levin and Malovichko's approximation for the reflection traveltime curve in vertically inhomogeneous media (, , ) defined by

$$t_0= \left(1-{1\over {S\left(t_v\right)}}\right) \,t_v+ {1\o... ...ht)^2} \over {v_{rms}^2\left(t_v\right)}}}}\;, \EQNLABEL{ussr}$$ (13)

where v_rms stands for the effective (root mean square) velocity along the vertical ray

$$v_{rms}^2\left(t_v\right)={\eta(z)\over t_v}= {1 \over t_v}\,\int_{0}^{t_v} v^2 \,dt_v\;, \EQNLABEL{vrms}$$ (14)

and S is the parameter of heterogeneity:

$$S\left(t_v\right)={1 \over{v_{rms}^4 t_v}}\,\int_{0}^{t_v} v^4 \,dt_v\;. \EQNLABEL{heter}$$ (15)

In terms of the S parameter, the variance of the squared velocity distribution along the vertical ray is

$$\sigma^2={1 \over t_v}\,\int_{0}^{t_v} v^4 \,dt_v - v_{rms}^4=v_{rms}^4 (S-1)\;. \EQNLABEL{sigma}$$ (16)

As follows from equality sigma, $S\geq 1$ for any type of velocity distribution (S equals 1 in a constant velocity case). For most of the distributions occurring in practice, S ranges between 1 and 2.

Malovichko's formula ussr is known as the most accurate three-parameter approximation of the NMO curve in vertically inhomogeneous media. Since reflection from a horizontal reflector in that class of media is kinematically equivalent to diffraction from a point, formula ussr can be similarly regarded as an approximation of the summation path of the post-stack Kirchhoff-type migration operator. In this case, it has the same meaning as formula levin. An important difference between the two formulae is the fact that equation ussr is written in the initial coordinate system and includes coefficients varying with depth, while equation levin applies the transformed coordinate system and constant coefficients. Using this fact, the rest of this section compares the accuracy of the approximations and relates Stolt's W factor to Malovichko's parameter of heterogeneity.

Equations levin and ussr both approximate the traveltime curve in the neighborhood of the vertical ray. Therefore, to compare their accuracy, it is appropriate to consider series expansion of the diffraction traveltime in the vicinity of the vertical ray:

$$t_0(l)={\left.t_0\right\vert _{l=0}}+ {1 \over 2}\,{\left.{d... ...ver {dl^4}\right\vert _{l=0}}l^4+\cdots\;\;, \EQNLABEL{taylor}$$ (17)

where l=x-x₀. Expansion taylor contains only even powers of l because of the obvious symmetry of t₀ as a function of l.

The special choice of parameters t_v, v_rms, and S allows Malovichko's formula ussr to provide correct values for the first three terms of expansion taylor:

$$\begin{eqnarray} \left.t_0\right\vert _{l=0} & = & t_v\;; \EQNLABEL{mal0}\\ \lef... ...ht)} \over {t_v^3 v_{rms}^4\left(t_v\right)}}\;\;. \EQNLABEL{mal4}\end{eqnarray}$$ (18) (19) (20)

Considering Levin's formula levin as an implicit definition of the function $t_0\left(t_v\right)$, we can iteratively differentiate it following the rules of calculus:

$$\left.{d^2s}\over {dl^2}\right\vert _{l=0} = \left.\left(s'... ...l=0} = {1\over {v_0^2 \,s\left(t_v\right)}}\;; \EQNLABEL{lev2}$$ (21)

$$\begin{eqnarray} \left.{d^4s}\over {dl^4}\right\vert _{l=0} & = & \left(6\,s'''\... ...-{{3\,W} \over {v_0^4 \,s^3\left(t_0\right)}}\;\;. \EQNLABEL{lev4}\end{eqnarray}$$ (22)

Substituting the definition of Stolt stretch transform ss into lev2 produces an equality similar to mal2, which means that approximation levin is theoretically accurate in depth-varying velocity media up to the second term in taylor. It is this remarkable property that proves the validity of the Stolt stretch method (, ). Formula levin will be accurate up to the third term if the value of the fourth-order traveltime derivative in lev4 coincides with mal4. Substituting equation mal4 into lev4 transforms the latter to the form

$${{1-W}\over {v_0^2\,s^2\left(t_v\right)}}= {{v^2\left(t_v\ri... ...over {v_{rms}^4\left(t_v\right)\,t_v^2}}\;. \EQNLABEL{malvsst}$$ (23)

It is now easy to derive from equation malvsst the desired explicit expression for the Stolt stretch parameter W, as follows:

$$W=1-{{v_0^2\,s^2\left(t_v\right)} \over{v_{rms}^2\left(t_v\r... ...t(t_v\right)}} -S\left(t_v\right) \right)\;\;. \EQNLABEL{main}$$ (24)

Expression main is derived so as to provide the best possible value of W for a given depth (vertical time t_v). To get a constant value for a range of depths one should take an average of the right hand side of main in that range. The error associated with Stolt stretch can be approximately estimated from taylor as the difference between the fourth-order terms:

$$\delta={{l^4 \over 8}\,{{W\left(t_v\right)-W} \over {t_v s^... ...ght) v_{rms}^2\left(t_v\right) v_0^2}}}\;, \EQNLABEL{maxerror}$$ (25)

where $W\left(t_v\right)$ is the right-hand side of main, and W is the constant value of W chosen for Stolt migration.

To estimate the best possible accuracy that the Stolt stretch method can achieve, we must take into account the sixth-order term in taylor related to the sixth-order derivative of the traveltime curve. For the true traveltime curve, the expression for the sixth-order derivative in the vicinity of the vertical ray is known from the literature (, ) to be

$$\left.{d^6t_0}\over {dl^6}\right\vert _{l=0} = {{45} \over ... ...right)}}\, \int_0^{t_v} v^6 \,dt_v\right)\;\;. \EQNLABEL{true6}$$ (26)

First, let us estimate the error of Malovichko's approximation ussr. Differentiating ussr six times and setting the offset l to zero yields

$$\left.{d^6t_0}\over {dl^6}\right\vert _{l=0} = {{45\,S^2\left(t_v\right)} \over {t_v^5\, v_{rms}^6}}\;. \EQNLABEL{mal6}$$ (27)

The estimated error is proportional to the difference between mal6 and true6:

$$\delta_M={{l^6} \over {6!}}\,\left[ {{45} \over {t_v^5\, v_{... ...dt_v-S^2\left(t_v\right) \right)\right]\;. \EQNLABEL{malerror}$$ (28)

It is interesting to note that replacing the parameter of heterogeneity S by its definition heter changes the expression in the round brackets to the following form:

$${1 \over {t_v \,v_{rms}^6}}\, \int_0^{t_v} v^6 \,dt_v-S^2= {... ...({\int_0^{t_v} v^4 \,dt_v}\right)^2\right)\;. \EQNLABEL{Cauchy}$$ (29)

According to the Schwarz inequality from calculus ( also known as the Cauchy-Bunyakovski inequality), the value of expression Cauchy can never be less than zero; hence $\delta_M\geq 0$ for any velocity distribution. This conclusion indicates that Malovichko's approximation tends to increase the traveltime at large offsets beyond its true value.

Differentiating lev4 twice and eliminating terms that vanish at l=0 produces

$$\begin{eqnarray} \left.{d^6s}\over {dl^6}\right\vert _{l=0} & = & \left.\left(15... ...{{45\,W^2} \over {s\left(t_v\right)^5\, v_0^6}}\;. \EQNLABEL{lev6}\end{eqnarray}$$ (30)

Evaluating the sixth-order traveltime derivative from lev6 and subtracting true6, we get a somewhat lengthy but explicit expression for the error associated with Stolt stretch approximation in the case of the best possible choice of W:

$$+ {l^6 \over 6!}\,\left[ {{45\, \left( 1-W \right) } \over ... ..._v^4\,v_{rms}^8\left(t_v\right)}} \right]\;. \EQNLABEL{leverr}$$ (31)

ISOTROPIC HETEROGENEITY VERSUS ANELLIPTIC ANISOTROPY A controversial issue associated with the topic of this paper is whether the non-hyperbolicity of the traveltime curves is caused mainly by heterogeneity or by anisotropy. To find a connection between the two different descriptions of media, we can consider an alternative three-parameter traveltime approximation (the anelliptic anisotropic moveout formula), proposed by Muir and Dellinger :

$$t_0 = {{t_v^4+(f+1)\,t_v^2\,{\left(x-x_0\right)^2 \over v_{r... ...f\,{\left(x-x_0\right)^2 \over v_{rms}^2}}}\;. \EQNLABEL{muir}$$ (32)

Here f is the parameter of anellipticity. Differentiating muir four times, setting l=x-x₀ to zero, and equating the result with mal4 results in the following formal relationship between f and Malovichko's parameter of heterogeneity:

$$S = 1 + 4\,f - 4\,f^2\;. \EQNLABEL{Svsf}$$ (33)

Equation Svsf clearly demonstrates the uncertainty between the anisotropic and heterogeneous isotropic interpretations. Both of them can explain the cause of the nonhyperbolicity of traveltime curves. An important difference is that the parameter of heterogeneity is uniquely determined by the velocity distribution according to heter, while the f parameter is assumed to be an independent functional. The definition heter, applied in combination with main, is suitable for calculating the Stolt stretch factor in an isotropic model for a given velocity function. If the correction parameter is measured experimentally by a non-hyperbolic velocity analysis in the form of either equation ussr or equation muir, it accumulates both heterogeneous and anisotropic factors and can be used for an explicit determination of W in main independently of the preferred explanation. In the case of the anisotropic moveout velocity analysis, we merely need to substitute the connection formula Svsf into main to find W. An alternative approach to Stolt-type migration in transversally isotropic media was proposed recently by Ecker and Muir . However, Stolt stretch migration is superior to that method in its ability to cope with varying rms velocities.

EXAMPLES A simple analytic example of isotropic heterogeneity is the case of a constant velocity gradient. In this case the velocity distribution can be described by the linear function $v\left(z\right)=v\left(0\right)(1+\alpha z)$. Stolt stretch transform is found from ss as

$$s(t)=\left({e^{2 \alpha v\left(0\right)\,t} -1 - 2 \alpha v\... ...\,t} \over {2 \alpha^2 v_0^2}\right)^{1/2}\;. \EQNLABEL{ssvz}$$ (34)

Let $\kappa$ be the logarithm of the velocity change v(z)/v(0). Then an explicit expression for W factor follows from main:

$$W={{2\,\kappa}\over{e^{2\,\kappa}-1}}={v^2\left(0\right) \over v_{rms}^2(\kappa)}\;. \EQNLABEL{wvz}$$ (35)

For $\kappa \rightarrow 0$ (a small depth or a small velocity gradient), $W \approx 1-\kappa$. For $\kappa \rightarrow \infty$ (a large positive change of velocity) W monotonically approaches zero. Formula wvz can be a useful rule of thumb for a rough estimation of W.

Numerical example of the Stolt stretch parameter computation is illustrated in Figures stovwt and stocvw. The left side of Figure stovwt shows a smoothed interval velocity curve from the Gulf of Mexico. The corresponding optimal values of the W factor as a function of vertical time (in the isotropic model) are shown on the right. Though the velocity function is smooth, substantial changes in W occur, making its mean value for the times $t_v \leq 6$ sec equal to 0.631.

The theory of cascaded migrations (, ) proves that Stolt-type f-k migration for a nonuniform velocity $v\left(t_v\right)$ can be performed as a cascaded process consisting of migrations with the smaller velocities $v_i\left(t_v\right)\,,\;i=1,2,\ldots,n$, such that $v_1^2+v_2^2+\cdots+v_n^2=v^2$ . As shown by Larner and Beasley , it is important to partition the velocity so that for each particular t_v all the velocities in the cascade, except maybe the last one, are constant. The advantage of the cascaded f-k migration method is based on the fact that each small velocity v_i describes a more homogeneous medium than the initial $v\left(t_v\right)$ function. Therefore, the W factor for each migration in a cascade is closer to 1, and the Stolt stretch approximation is more accurate. This fact is illustrated in Figure stocvw, which shows an optimal partitioning of the velocity and the corresponding values of the W factor. In accordance with the empirical conclusions of Beasley et al. , a cascade of only four migrations was sufficient to increase the value of W to more than 0.8. With a further increase of the number of cascaded migrations, the method becomes as accurate with respect to vertical velocity variations as phase-shift migration. Theoretically, this limit corresponds to the velocity continuation concept (). Note that the theory of cascaded f-k migration is strictly valid for isotropic models. The anisotropic interpretation does not support it, since the intrinsic anisotropy factor is not supposed to change with the velocity partitioning.

 

stovwt
Figure 2 Smoothed interval velocity distribution from the Gulf of Mexico (left) and the corresponding W factor as a function of vertical time (right). The mean value of W is 0.631.

 

stocvw
Figure 3 Left: Optimal partitioning of the velocity function for the method of cascaded migrations. Right: corresponding mean values of W. Top: four cascaded migrations. Bottom: seven cascaded migrations.

The main result of this paper is an analytic explicit expression main that allows us to choose the most appropriate value for the Stolt stretch factor. Possible applications include the optimal design of interval velocities partitioning for the method of cascaded f-k migrations and extension of the Stolt stretch method to a transversally isotropic model.

Nowadays the topic of this paper seems to be out of fashion. When everyone is interested in prestack depth migration in the time-space domain, it is difficult to attract any attention to post-stack time migration in the frequency domain. Nevertheless, I believe the art of approximation demonstrated by Robert Stolt in his famous paper to be a good example to follow when working on many different problems, which was the main reason for this research.

[SEP,MISC,GEOTLE,paper,SEGCON]