The plane-wave destructor (PWD)

 A plane wave in a wavefield u(t,x)=u(t-px) with stepout p can be extinguished with a partial differential operator that we write as a matrix $\bold A$ where

$$\begin{eqnarray} 0 \quad\approx\quad v(t,x) &=& \left( {\partial \over \partial ... ... \ u(t,x) \\ \bold 0 \quad\approx\quad \bold v &=& \bold A \bold u\end{eqnarray}$$ (13) (14)

The parameter p is called the wave number or Snell parameter and |p| can take on any value less than 1/v where v is the medium velocity. The wave's angle of propagation is $\sin\theta = pv$.

#   vv = (aa Dx - pp Dt) uu
#
subroutine wavekill1( aa, pp, m1,n1,n2 , uu, vv)
integer conj, i1,i2, m1,n1,n2
real aa, pp, s11, s12, s21, s22, uu(m1,n2), vv(m1,n2-1)
s11 = -aa+pp;	s12 = aa+pp
s21 = -aa-pp;	s22 = aa-pp
call zero( m1*(n2-1), vv)
do i1=1, n1-1 {		# vv is one point shorter than uu on both axes.
do i2=1, n2-1 {
	vv(i1,i2) = vv(i1,i2) +
		uu(i1  ,i2) * s11 + uu(i1  ,i2+1) * s12 +
		uu(i1+1,i2) * s21 + uu(i1+1,i2+1) * s22
	}}
return; end

In a finite difference representation, the difference operator $\delta_x - p_i\delta_t$is a two-dimensional filter that fits on a $2\times 2$ differencing star. As a matrix operation, this two dimensional convolution is denoted $\bold A$.(More details can be found in my previous book ``''.) The program wavekill1() applies the operator $a \delta_x -p\delta_t$which can be specialized the operators $\delta_x$, $\delta_t$,$\delta_x - p_i\delta_t$.I carefully arranged the side boundaries so the filter never runs off the sides of the data. Thus the output is shorter than the input by one point both the t-axis and the x-axis. The reason for using these side boundaries is so that large data sets can be chopped into independent sections without the boundaries themselves affecting the result. By chopping a large data set into sections, we can handle curved events as piecewise linear.

When only one wave is present and the data is adequately sampled, then finding the best value of p is a single-parameter, linear least-squares problem. Let $\bold x$ be a vector whose components are values of $\partial u / \partial x$ taken on a mesh in (t,x). Likewise let $\bold t$ contain $\partial u / \partial t$.Since we wish $\bold x - p\, \bold t \approx \bold 0$,we minimize the quadratic function of p

$$Q(p) {=} (\bold x - p\, \bold t) \cdot (\bold x - p\, \bold t)$$ (15)

by setting to zero the derivative getting

$$p {=}{\bold x \cdot \bold t \over \bold t \cdot \bold t}$$ (16)

As data will not always fit the model very well, it may be helpful to have some measure of goodness of fit. I suggest

$$Q^2 {=}1 \ -\ { (\bold x -p\, \bold t) \cdot (\bold x -p\, \bold t) \over \bold x \cdot \bold x}$$ (17)

which on inserting $p=(\bold x \cdot \bold t)/(\bold t \cdot \bold t)$leads to

$$Q {=} { \bold x \cdot \bold t \over \sqrt{ (\bold x \cdot \bold x) (\bold t \cdot \bold t) } }$$ (18)

which is known as the normalized correlation.  The program for this calculation is straightforward. I named the program puck() to denote picking on a continuum.

# measure coherency and dip, and compute residual  res = (Dx - p Dt) uu
#
subroutine puck ( m1, n1, n2, uu, coh, pp, res )
integer i1, i2, m1, n1,n2
real	uu(m1,n2), res(m1,n2), xx, xt, tt, coh, pp
temporary real dx(m1,n2-1), dt(m1,n2-1)
call wavekill1( 1., 0.,  m1,n1,n2 , uu, dx)
call wavekill1( 0.,-1.,  m1,n1,n2 , uu, dt)
xx = 1.e-30; tt = 1.e-30; xt = 0.
do i1=1,n1-1 {
do i2=1,n2-1 {
	xt = xt + dt(i1,i2) * dx(i1,i2)
	tt = tt + dt(i1,i2) * dt(i1,i2)
	xx = xx + dx(i1,i2) * dx(i1,i2)
	}}
coh = sqrt((xt/tt)*(xt/xx))
pp  =  xt/tt 
call wavekill1( 1., pp,  m1,n1,n2 , uu, res)
return;	end

Finally, an undesirable feature of equation (13) is that the derivatives roughen the data so the spectrum of $\bold x$ and $\bold t$ is not that of $\bold u$.What I am doing now is putting the data through leaky integration on the time axis before I start. I am planning to switch this from t to x so the operator can be thought of as $1-p\delta_t^x$.Then the residual should be in phase with the input signal so I can subtract it from the data leaving ``cleaned up'' data.