Description of the method

Aki and Richards (1980) stated that, up to $40^{\circ}$, the p-wave reflection coefficients, $R(\theta)$, as a function of the incidence angle, $\theta$, can be approximated by  

$$R(\theta) = A + B \sin^{2}{\theta} + C \tan^2{\theta},$$ (1)

where A, B and C are parameters related to the average and contrasts of the elastic properties of the limiting layers. Therefore, if it is possible to determine these parameters, one can simulate the amplitudes of primary events at every depth or time step. I refer to the simulation of the primary amplitudes as AVA modeling. The aim of this work is not to invert the AVA parameters for rock and fluid properties, but just to obtain a reasonable estimate of the primaries to be further used in a residual-multiple-attenuation sequence.

The AVA modeling problem can be addressed by solving locally, for every different depth (or time) step, the normal equations $\bold {(G^{^T}G)} \bold {\hat m} = \bold {G^{^T}} \bold d$,where $\bold {\hat m}$ are the estimated AVA parameters, $\bold G$ has the column vector form $[\bold 1$, $\bold{sin^2 \theta_i}$, $\bold{tan^2 \theta_i}]$ and $\bold d$ is the ADCIG (or a CMP gather after NMO). The simulated primaries are obtained by using the estimated AVA parameters in equation (1). One problem that arises from this local solution is that the AVA parameters are not constrained to be smooth along depth. This can lead to anomalous amplitudes in the simulated primaries for a certain depth.

In the present approach, the AVA modeling problem consists of determining the three parameters A, B and C of equation (1) by solving the following data fitting problem:  

$$\bold L \bold m = \bold d ,$$ (2)

where $\bold L$ in matrix form is  

$$\left[ \begin{array} {ccc} {\bold 1}_1&{\bold S}_1&{\bold T... ...1}&{\bold S}_{n_1}&{\bold T}_{n_1}\\ \end{array} \right] \ \\ end{displaymath}$$ (3)

which has dimension $n_{_1}n_{_2} \times 3n_{_1}$, $\bold m$ is the model-parameter vector (3n₁ elements), and $\bold d$ is the data vector (n₁n₂ elements), where n₁ and n₂ are the number of samples in depth (or time) and the number of traces in the ADCIG, respectively. ${\bold 1}_j$, ${\bold S}_j$ and ${\bold T}_j$ are $n_{_2} \times n_{_1}$ matrices consisting of 1, $(\sin^2 \theta_i)$ and $(\tan^2 \theta_i)$ along their j-th column entries, respectively. i stands for the reflection-angle indexes.

Generally, in ADCIGs, residual multiples are more persistent at near angles, because of insufficient moveout difference between them and the flattened primaries. Additionally, at the farthest angles, stretch occurs. To avoid these imperfections in the input data, the fitting goal is to minimize the residual,  

$$\bold 0 \approx \bold M (\bold L \bold m - \bold d),$$ (4)

where $\bold M$ is a selector operator which applies the appropriate internal and external mutes.

If unrealistic variations of the AVA parameters with depth are an issue, the following regularization goal can be introduced:  

$$\bold 0 \approx \epsilon \bold D \bold m ,$$ (5)

where $\bold D$ is the derivative operator along depth and $\epsilon$ is the regularization parameter. As usual, care must be taken when choosing $\epsilon$ not to destroy any recoverable residual moveout information and not to spread simulated primaries to angles at which they originally do not occur.

To accelerate the solution, the problem can be solved with preconditioning Claerbout and Fomel (2001), using the transformation $\bold m = \bold C \bold p$, where $\bold C = D^{-1}$ and $\bold p$ is the preconditioned variable. Finally, the fitting goals reduce to  

$$\begin{array} {c} \bold 0 \approx \bold M (\bold {LC} \bold p - \bold d)\\ \bold 0 \approx \epsilon \bold p .\\ \end{array}$$ (6)

The final model is obtained with $\bold m = \bold C \bold p$. I use conjugate-gradients to solve the inverse problem.

After the determination of the three AVA parameters for every depth step, primaries are simulated by computing the reflection coefficient for all reflection angles using equation (1). Of course, the method relies on the flatness of the reflectors in the CIG to correctly extract the 3 parameters of the AVA curve.

To get the residual-multiple-attenuated data, the adaptive subtraction must be applied in two steps. The first one aims to obtain an estimate of the residual multiples. This is done by subtracting the adjusted version of the simulated primaries from the original data. Amplitude and phase adjustments are achieved by the convolution of the simulated primaries with a prediction-error Wiener filter.

As the estimated residual multiples may contain some primary information (mainly at small reflection angles), I use the strategy proposed by Guitton et al. (2001), the so-called `subtraction method,' in which two different prediction-error filters (PEFs), which model primaries and multiples, are computed. Their method allows regularization, to decrease the crosstalk between multiples and primaries. The corresponding fitting goals are

$\bold 0 \approx \bold m_n$
$\bold 0 \approx \epsilon \bold m_s$
subjected to $\leftrightarrow \bold d = \bold S^{-1} \bold m_s + \bold N^{-1}\bold m_n$

after preconditioning, with $\bold s = \bold {S^{-1}m_s}$ and $\bold n = \bold {N^{-1}m_n}$. In the preconditioning equations, $\bold s$ represents the primaries, $\bold n$ the multiples; $\bold N^{-1}$ and $\bold S^{-1}$ represent deconvolution with PEFs for multiples and primaries, respectively; $\bold {m_n}$ is the multiple model content and $\bold {m_s}$ is the primary model content. Finally, the estimated primaries are obtained by $\bold {\hat s} = \bold d - \bold N^{-1} \bold {m_n}$. All PEFs are computed with conjugate-gradients.