MONO-PLANE DECONVOLUTION

Consider two-dimensional filters. The dip and deconvolution character of such filters is described in considerable detail in chapter 8 of PVI. The coefficients of these two-dimensional convolution filters are chosen to produce minimum power out of a given input. Here we limit attention to the filter coefficient pattern below:  

$$\begin{array} {cc} c &d \\ b &e \\ a &f \\ 1 &g \\ \cdot &h \\ \cdot &i \\ \cdot &j \end{array}$$ (1)

The vertical axis is time. The horizontal axis is space. The ``''s are zeros. The filter can be lengthened in time but not in space. The choice of exactly two columns is a choice to have an analytic form that can exactly destroy a single plane, but cannot destroy two. (The code shown can be used for more columns and more planes.) With two signals that are statistically independent, the filter reduces to the well-known prediction-error filter in the left column and zeros in the right column. If the filter coefficients were extended in both directions on t and to the right on x, it would flatten the two-dimensional spectrum of the input. (An appendix proves these whitening properties). Imagine all the coefficients in (1) vanished but d=-1 and the given 1. Such a filter would annihilate an appropriately sloping plane wave.

A subroutine for two-dimensional filtering is cinlof() from PVI. That subroutine handles any number of columns so we would limit it to two columns by making the choice a2=2. This routine is overly general for the task at hand, but it is a well tested starting place. With cinlof() the filter stops before it runs off the end of the data so there are no edge transients and as a necessary byproduct, the output space is somewhat smaller than the input space. Also, the output aligns with the input under the coefficient that is constrained (in autoregression) to unity.

# CINLOF -- Convolution INternal with Lags.  Output is conjugate with FILTER.
#
subroutine cinlof( conj, add, lag1,lag2,  n1,n2,xx,  nb1,nb2,bb,  yy)
integer            conj, add, lag1,lag2,  n1,n2,     nb1,nb2
real                                      xx(n1,n2), bb(nb1,nb2), yy(n1,n2)
integer y1,y2, x1,x2, b1, b2
call conjnull( conj, add,                            bb,nb1*nb2,  yy,n1*n2)
if( conj == 0 )
	do b2=1,nb2 { do y2= 1+nb2-lag2, n2-lag2+1 {  x2= y2 - b2 + lag2
	do b1=1,nb1 { do y1= 1+nb1-lag1, n1-lag1+1 {  x1= y1 - b1 + lag1
				yy(y1,y2) = yy(y1,y2) + bb(b1,b2) * xx(x1,x2)
				}} }}
else
	do b2=1,nb2 { do y2= 1+nb2-lag2, n2-lag2+1 {  x2= y2 - b2 + lag2
	do b1=1,nb1 { do y1= 1+nb1-lag1, n1-lag1+1 {  x1= y1 - b1 + lag1
				bb(b1,b2) = bb(b1,b2) + yy(y1,y2) * xx(x1,x2)
				}} }}
return; end