Traveltime interpolation

Now suppose I know a set of n Green's function traveltime fields $\tau_i({\bf x};{\bf x}_i)$ for n source positions ${\bf x}_i$. Estimating the traveltime field $\tau_k({\bf x};{\bf x}_k)$, due to a single intermediate surface source position ${\bf x}_k$, from the n given fields $\tau_i$, is defined as traveltime interpolation.

For the set of $\tau_i$, (7) can be generalized to a linear matrix equation of the form:

 

$$\nabla\tau_i {\bf \cdot}\nabla\tau_{k} - s^2 \cos\theta_{ik} = 0 \; ; \; i=1,n \;,$$ (12)

or into a general weighted least-squares system ${\bf \cal W}{\bf \cal A}{\bf x}={\bf \cal W}{\bf y}$:

 

$$\pmatrix{ w_1 \cr & w_2 \cr & & w_3 \cr & & & \ddots \c... ...\cos\theta_{3k} \cr \vdots \cr w_n \cos\theta_{nk} \cr } \;.$$ (13)

System (13) can be solved by standard damped least-squares,

 

$${\bf x}\sim ( {\bf \cal A}^T {\bf \cal W}^T{\bf \cal W}{\bf ... ...cal I})^{-1} {\bf \cal A}^T{\bf \cal W}^T{\bf \cal W}{\bf y}\;,$$ (14)

or by the optimal (but slow) Singular Value Decomposition (SVD) method, for the components of the unknown vector $\nabla\tau_k$.I refer the reader to Strang (1980) for an excellent review of damped least squares and SVD linear algebra techniques. The diagonal weighting matrix ${\bf \cal W}$ is of dimensions (nxn), and the weights w_i could be inverse-distance: $w_i \sim \vert{\bf x}_i -{\bf x}_k\vert^{-1}$, for example. The ${\bf \cal A}$ matrix is (nx3) in 3-D and (nx2) in 2-D (ignoring the $\d_y$ terms), and contains the vector components of the known traveltime gradient fields. The data vector ${\bf y}$ is (nx1) and contains the $\cos\theta_{ik}$ values. For a 2-D geometry, at least two traveltime gradient fields must be known, and for a 3-D interpolation at least three traveltime gradient fields must be available.

In the case of interpolation, a refinement to the estimate of the $\cos\theta_{ik}$ values is possible. In this paper, I use the Law of Cosines given by (8) as an initial estimate of, say, $\cos\theta_{12}({\bf x})$ and $\cos\theta_{23}({\bf x})$, where the traveltime field at ${\bf x}_2$ is to be interpolated from two bounding source locations ${\bf x}_1$ and ${\bf x}_3$. Please refer to Figure to consider the appropriate interpolation geometry. However, the angle $\theta_{13}$ is known precisely from the given traveltime (gradient) fields at ${\bf x}_1$ and ${\bf x}_3$:

 

$$\cos\theta_{13} \equiv v^2 \nabla\tau_1 {\bf \cdot}\nabla\tau_3 \;.$$ (15)

Therefore, given the angle estimates $\theta_{12}'$, $\theta_{23}'$ and $\theta_{13}'$, and the true angle $\theta_{13}$, I can refine my initial (primed) estimates such that the total angle $\theta_{13}$ is conserved as a sum of $\theta_{12}$, $\theta_{23}$:

 

$$\theta_{13} = f \theta_{13}' = f (\theta_{12}' + \theta_{23}') = \theta_{12} + \theta_{23} \;.$$ (16)

Then, the correction factor f is given as

$$f = \theta_{13}/ \theta_{13}' \;,$$

and  

$$\theta_{12} \sim f \theta_{12}' \;\;\; ; \;\;\; \theta_{23} \sim f \theta_{23}' \;.$$ (17)

As an example, for a 2-D interpolation of $\tau_2$ given $\tau_1$ and $\tau_3$, I first set up a 2x2 matrix system per (13). I evaluate the required data vector values $\cos\theta_{12}$ and $\cos\theta_{23}$ using (8), and refine the cosine estimates using (17). I then invert system (13) by Cramer's Rule if the determinant is not too small, or by SVD if otherwise. The entire procedure is repeated in a point-by-point independent manner, making it an obvious candidate for massively parallel implementation, for all subsurface grid locations ${\bf x}$.

Finally, the traveltime gradients are spatially integrated to traveltime such that, for a source at (x₂,z₂) on a flat surface z = z₂,

 

$$\tau_2(x,z) = \gamma(x,z_2) + \int_{z_2}^{z}\tau_z(x,z')\,dz' \;,$$ (18)

where  

$$\gamma(x,z_2) = H(x-x_2)\int_{x_2}^{x}\tau_x(x',z_1)\,dx' - H(x_2-x)\int^{x_2}_{x}\tau_x(x',z_1)\,dx' \;,$$ (19)

where H is the Heaviside function. This process completes the traveltime interpolation of $\tau_2$ from $\tau_1$ and $\tau_3$.

 

geometry
Figure 1 The traveltime interpolation geometry. It is assumed that the traveltime field $\tau({\bf x};{\bf x}_2)$ due to a source position at ${\bf x}_2$ can be interpolated from two known traveltime fields $\tau({\bf x};{\bf x}_1)$ and $\tau({\bf x};{\bf x}_3)$ due to sources at positions ${\bf x}_1$ and ${\bf x}_3$ respectively.