next up previous [pdf]

Next: MULTIVARIATE LEAST SQUARES Up: UNIVARIATE LEAST SQUARES Previous: Formal path to the

The plane-wave destructor

We address the question of shifting signals into best alignment. The most natural approach might seem to be via cross correlations, which is indeed a good approach when signals are shifted by large amounts. Here, we assume signals are shifted by small amounts, often less than a single pixel. We take an approach closely related to differential equations. Consider this definition of a residual.
\begin{displaymath}
0 \quad \approx \quad \hbox{residual}(t,x) \eq \left( \frac{...
...al}{\partial x} + p \frac{\partial}{\partial t} \right) u(t,x)
\end{displaymath} (14)

By taking derivatives we see the residual vanishes when the two-dimensional observation $u(t,x)$ matches the equation of moving waves $u(t-px)$. The parameter $p$ has units inverse to velocity, the velocity of propagation.

In practice, $u(t,x)$ might not be a perfect wave but an observed field of many waves that we might wish to fit to the idea of a single wave of a single $p$. We seek the parameter $p$. First, we need a method of discretization that allows the mesh for $\partial u/\partial t$ to overlay exactly $\partial u /\partial x$. To this end, I chose to represent the $t$-derivative by averaging a finite difference at $x$ with one at $x+\Delta x$.

\begin{displaymath}
\frac{\partial u}{\partial t} \quad \approx \quad \frac{1}{2...
...u(t+\Delta t,x+\Delta x) - u(t,x+\Delta x) }{\Delta t}
\right)
\end{displaymath} (15)

Likewise, there is an analogous expression for the $x$-derivative with $t$ and $x$ interchanged. The function $u(t,x)$ lies on a grid, and the differencing operator $\delta_x + p\delta_t$ lies atop it and convolves across it. The operator is a $2\times 2$ convolution filter. We may represent equation (14) as a matrix operation,
\begin{displaymath}
{\bf0} \quad \approx\quad {\bf r} = {\bf A}{\bf u}
\end{displaymath} (16)

where the two-dimensional convolution with the difference operator is denoted ${\bf A}$.

Now, let us find the numerical value of $p$ that fits a plane wave $u(t-px)$ to observations $u(t,x)$. Let $\bold x$ be an abstract vector having components with values $\partial u /\partial x$ taken everywhere on a 2-D mesh in $(t,x)$. Likewise, let $\bold t$ contain $\partial u/\partial t$. Because we want ${\bf x} + p {\bf t} \approx {\bf0}$, we minimize the quadratic function of $p$,

\begin{displaymath}
Q( p) = ({\bf x} + p {\bf t}) \cdot ({\bf x} + p {\bf t})
\end{displaymath} (17)

by setting to zero the derivative by $p$. We get:
\begin{displaymath}
p \eq - \ \frac{ {\bf x} \cdot {\bf t} }{ {\bf t} \cdot {\bf t} }
\end{displaymath} (18)

Because data does not always fit the model very well, it may be helpful to have some way to measure how good the fit is. I suggest:
\begin{displaymath}
C^2 \eq 1 \ -\ \frac{ ({\bf x} + p {\bf t}) \cdot ({\bf x} + p {\bf t}) }{ {\bf x} \cdot {\bf x}}
\end{displaymath} (19)

which, on inserting $p=-({\bf x} \cdot {\bf t})/({\bf t} \cdot {\bf t})$, leads to $C$ , where
\begin{displaymath}
C \eq \frac{ {\bf x} \cdot {\bf t} } { \sqrt{ ( {\bf x} \cdot {\bf x} ) ( {\bf t}\cdot{\bf t}) }}
\end{displaymath} (20)

is known as the ``normalized correlation.''

To suppress noise, the quadratic functions $\bold x \cdot \bold x$, $\bold t \cdot \bold t$, and $\bold x \cdot \bold t$ were smoothed over time with a triangle filter.

puckin
Figure 3.
Input synthetic seismic data includes a low level of noise.
puckin
[pdf] [png] [scons]

residual5
Figure 4.
Residuals, i.e., an evaluation of $U_x + pU_t$.
residual5
[pdf] [png] [scons]

puckout
Figure 5.
Output values of $p$ are shown by the slope of short line segments.
puckout
[pdf] [png] [scons]

Subroutine puck2d shows the code that generated Figures 3-5. An example based on synthetic data is shown in Figures 3 through 5. The synthetic data in Figure 3 mimics a reflection seismic field profile, including one trace that is slightly delayed as if recorded on a patch of unconsolidated soil.

Figure 4 shows the residual. The residual is small in the central region of the data; it is large where the signal is not sampled densely enough, and it is large at the transient onset of the signal. The residual is rough because of the noise in the signal, because it is made from derivatives, and the synthetic data was made by nearest-neighbor interpolation. Notice that the residual is not particularly large for the delayed trace.

Figure 5 shows the dips. The most significant feature of this figure is the sharp localization of the dips surrounding the delayed trace. Other methods based on ``beam stacks'' or Fourier concepts might lead us to conclude that the aperture must be large to resolve a wide range of angles. Here, we have a narrow aperture (two traces), but the dip can change rapidly and widely.

Once the stepout $p = d t /d x$ is known between each of the signals, it is a simple matter to integrate to get the total time shift. A real-life example is shown in Figure 6.

twod
Figure 6.
A seismic line before and after flattening.
twod
[pdf] [png] [scons]

In this case the flattening was a function of $x$ only. More interesting (and more complicated) cases arise when the stepout $p = d t /d x$ is a function of both $x$ and $t$. The code shown here should work well in such cases.

A disadvantage, well known to people who routinely work with finite-difference solutions to partial differential equations, is that for short wavelengths a finite difference operator is not the same as a differential operator; therefore, the numerical value of $p$ is biased. This problem can be overcome in the following way. First, estimate the slope $p = d t /d x$ between each trace. Then, shift the traces to flatten arrivals. Now, there may be a residual $p$ because of the bias in the initial estimate of $p$. This process can be iterated until the data is flattened. Everywhere in a plane we have solved a least squares problem for a single value $p$.


next up previous [pdf]

Next: MULTIVARIATE LEAST SQUARES Up: UNIVARIATE LEAST SQUARES Previous: Formal path to the

2014-12-01