next up previous print clean
Next: Computational cost savings Up: Lomask and Fomel: Flattening Previous: Introduction

Method

The flattening method described in Lomask et al. (2005) creates a time-shift (or depth-shift) field $\boldsymbol{\tau}(x,y,t)$ such that its gradient approximates the dip $\bf{p}\rm (x,y,\boldsymbol{\tau})$. The dip is a function of $\boldsymbol{\tau}$ because for any given horizon, the appropriate dips to be summed are the dips along the horizon itself. Using the gradient operator ($\boldsymbol\nabla=[ \frac{\partial }{\partial x} \; \; \frac{\partial }{\partial y}]^T$) and the estimated dip (${\bf p} = [{{\bf p}_x} \; \; {{\bf p}_y}]^T$), our regression is
\begin{displaymath}
\boldsymbol{\nabla} {\boldsymbol \tau}(x,y,t)\quad = \quad {\bf p}(x,y,{\boldsymbol \tau}).\end{displaymath} (1)
Note that because the estimated dip $\bf{p}\rm (x,y,\boldsymbol{\tau})$ field is a function of the unknown $\boldsymbol{\tau}(x,y,t)$ field, this problem is non-linear and, therefore, difficult to solve directly. We solve this using a Gauss-Newton approach by iterating over equations (2)-(4), i.e.,


		 iterate {    
         \begin{eqnarray}
\bf r \quad &=& \quad [\boldsymbol{\nabla} \boldsymbol{\tau}_k ...
 ...} \quad &=& \quad \boldsymbol{\tau}_{k} + \Delta \boldsymbol{\tau}\end{eqnarray} (2)
(3)
(4)

		  }  ,  
where the subscript k denotes the iteration number.

To greatly improve efficiency, we solve equation (3) in the Fourier domain. We apply the divergence to the estimated dips and divide by the Z-transform of the discretized Laplacian in the Fourier domain as  
 \begin{displaymath}
\Delta \boldsymbol{\tau} \quad \approx \quad {\rm FFT_{\rm 2...
 ... r} \right]}\ \over { -Z_x^{-1} -Z_y^{-1} +4 -Z_x -Z_y} \right]\end{displaymath} (5)
where $Z_x = e^{i w \Delta x} \ $and$ \ Z_y = e^{i w \Delta y}$. This is extremely efficient compared to iterative methods for solving equation (3).

Unfortunately, mirroring is required to avoid boundary affects in the Fourier domain. Mirroring is concatenating reversed copies of the data along each dimension so that the input to the tranform has periodic boundaries. In 1D the input is increased by a factor of 2, in 2D the input is increased by a factor of four, and so on. In this case, we mirror the divergence of the residual ($\nabla\it{^T}{\bf r}$). Mirroring is described in more detail by Ghiglia and Pritt (1998) for application to 2D phase unwrapping.

Ghiglia and Pritt (1998) also describe an alternative method using discrete cosine transforms that does not need mirroring. This method does not require mirroring because the cosine transform assumes that image is even (periodic) implicitly. Since it makes this assumption, then only the unmirrored input itself must be transformed. The alternative formulation is  
 \begin{displaymath}
\Delta \boldsymbol{\tau} \quad \approx \quad {\rm DCT_{\rm 2...
 ...]}\ \over { -2\cos( w \Delta x) -2\cos(w \Delta y)+4 } \right].\end{displaymath} (6)