A new bidirectional deconvolution method that overcomes the minimum phase assumption

Introduction

In the previous report (Zhang and Claerbout, 2010), we introduced the spiking deconvolution problem using the hybrid norm solver (Claerbout, 2009a). Synthetic examples (Zhang and Claerbout, 2010) showed that given a minimum-phase wavelet, it retrieved the sparse reflectivity model almost perfectly even with a reflection series that is far from white, while conventional L2 deconvolution did a poor job. However, if the assumption of a minimum-phase wavelet was removed, the hybrid norm spiking deconvolution failed quickly and gave a poor result similar to the conventional L2 deconvolution.

In this paper, we still rely on the hybrid norm solver to retrieve the sparse model, but we use a slightly more complex formulation that avoids the minimum-phase wavelet constraint.

We start by realizing that any (mixed-phase) wavelet $C(Z)$ can be decomposed into a minimum-phase part $A(Z)$ and a maximum-phase part $B(1/Z)$ plus a certain time shift:

$$C(Z)=A(Z)B(1/Z)Z^{k},$$ (1)

where $B(Z)$ is also a minimum-phase wavelet (therefore $B(1/Z)Z^k$ is a maximum-phase wavelet) and the exponent $k$ is the order of $B(Z)$ . This $Z^{k}$ term makes the wavelet $C(Z)$ causal. In the time domain, (1) can be written as

$$c = a*b^r*\delta(n-k),$$ (2)

where $b^r$ stands for the time reverse of series $b$ .

Our original spiking deconvolution can find only a minimum-phase wavelet which has the same spectrum of real wavelet $c$ . It can be defined as an inverse problem as follows:

$$[d]f_c = r,$$ (3)

where $[d]$ is the data convolution operator, and $f_c$ is the unknown filter. In this formulation, the filter is the only unknown, the hybrid norm is applied on the residual term $r$ to enforce the sparseness constraint. In theory, the residual $r$ itself is the reflectivity model. Such a method requires the wavelet in the data to be minimum-phase because only a minimum-phase wavelet has a causal stable inverse.

The following bidirectional deconvolution formulation utilizes a pair of conventional deconvolutions, trying to invert components $a$ and $b$ separately:

$$\begin{array}{l} \left[(d*f_b^r) \right] f_a = r_{a}, \\ \left[(d*f_a)^r \right] f_b = r_{b}, \end{array}$$ (4)

in which $f_a$ and $f_b$ are the corresponding filters that corresponds to the inverses of $a$ and $b$ denoted above, the superscript $r$ means time-reverse. The operator in each equation is the convolution operator. Again the hybrid norm is applied to $r_a$ and $r_b$ , and the reflectivity model is simply $r_{a}$ plus a time shift. Notice that this is a non-linear inversion, since the operator itself depends on the unknown $f_a$ and $f_b$ . In practice we have to solve these two inversions alternately and therefore iteratively.

To understand the meaning of (4), let

$$d = m * c = m*a*b^r*\delta(n-k),$$ (5)

where m is the reflectivity model and the $\delta$ term is just a time shift. Assume $f_a$ and $f_b$ are perfectly known in the operators (which is not true in reality), i.e.

$$f_a*a = \delta(n), f_b*b=\delta(n)$$

Substituting (5) into (4), since

$$d*f_b^r = m*\delta(n-k)*a,$$ (6)

$$(d*f_a)^r = (m*b^r*\delta(n-k))^r = m^r*\delta(n+k)*b,$$ (7)

we have

$$\begin{array}{l} \left[(m*\delta(n-k))*a \right] f_a = r_{a}, \\ \left[(m^r*\delta(n+k))*b \right] f_b = r_{b}. \end{array}$$ (8)

From (8) it is easier to see what is behind the bidirectional deconvolution formulation (4): It tries to separate the two parts of the wavelet, turning each one into a traditional deconvolution problem in which the wavelet ($a,b$ ) is always minimum-phase.

As with all non-linear estimation, iteration is required. Convergence is assured if the starting solution is close enough. We expect the traditional PEF for $a$ and an impulse function for $b$ to be a pretty good first guess. The following section shows several examples (complexity varies from low to high) illustrating the effectiveness and limitations of the method.

A new bidirectional deconvolution method that overcomes the minimum phase assumption