Nonstationarity, moving windows

Seismic reflections are not generally planar though they are locally so. The puck() method was designed to be ignorant of side boundaries so that it can be applied in a small window and the window moved freely about the data. A strength of the method is that the window can be smaller than a wavelength--it can be merely two traces wide. A sample is shown in Figures 2 through 4.

 

puckin Figure 2 Input synthetic data.

 

residual Figure 3 Residuals, i.e. an evaluation of Ux - p Ut.

The residual shows several features: It is small in the central region of the data; it is large where there is spatial aliasing; and it is rough because it is made from derivatives and the synthetic data was made with nearest neighbor interpolation.

 

puckout Figure 4 Output values of p are shown by the slope of short line segments.

The residuals are small. If a weaker dip were present it would be revealed in the residuals. In Figure 3 there is no weaker dip, only artifacts of inadequate sampling of the signals and the transient onset of the signal. Later in chapter 8 several strong dips will be estimated simultaneously.

# slide a window around on a wall of data measuring coherency, dip, residual
#
subroutine slider( n1,n2, w1,w2, k1,k2, data, coh, pp, res)
integer	n1,n2, w1,w2, k1,k2, s1,s2, e1,e2, i1,i2
integer	p1,p2		# math size of pp(p1,p2), coh(p1,p2)
real	data(n1,n2)			# input
real	res(n1,n2)			# outputs. math size (n1-1,n2-1)
real	pp(n1,n2), coh(n1,n2)		# outputs.
temporary real  window(n1,n2), tres(n1,n2), count(n1,n2)
call zero( n1*n2, res)
call zero( n1*n2, count)
p2=1; s2=1; e2=s2+w2-1; while( e2<=n2) {
p1=1; s1=1; e1=s1+w1-1; while( e1<=n1) {
	do i1=1,w1
		do i2=1,w2
			window(i1,i2) = data(s1+i1-1, s2+i2-1)
	call puck ( n1, w1, w2, window, coh(p1,p2), pp(p1,p2), tres)
	do i1=1,w1-1 {
	do i2=1,w2-1 {
		res(  s1+i1-1, s2+i2-1) =    res(s1+i1-1, s2+i2-1) + tres(i1,i2)
		count(s1+i1-1, s2+i2-1) =  count(s1+i1-1, s2+i2-1) + 1.
		}}
	  p1=p1+1; s1=s1+k1; e1=s1+w1-1
	} p1=p1-1
	  p2=p2+1; s2=s2+k2; e2=s2+w2-1
	} p2=p2-1
do i2=1,n2-1
	do i1=1,n1-1
		if( count(i1,i2) > 0. )
			res(i1,i2) = res(i1,i2) / count(i1,i2)
return; end

A disadvantage of the method is that the finite difference operator is susceptible to spatial aliasing as well as distortions at spatial frequencies that are high but not yet aliased. This suggests the next step--estimating missing interlaced traces--that is taken up in chapter 8.