LSJIMP Nonlinear Iterations  

Proponents of so-called ``full-waveform inversion'' Tarantola (1984) seek to solve an ambitious nonlinear inverse problem. In the context of reflection seismology, full-waveform techniques seek to invert a relation of the following qualitative form:

$$\cal{L} \left( \mbox{density, velocity} \right)$$ recorded data =

Most generally, the operator $\cal{L}$ is elastic two-way wave equation modeling, although making the acoustic assumption simplifies matters. In both cases, the modeling operator is non-linearly dependent on the earth model to be estimated. For many reasons (including, but not limited to: computational expense, model non-uniqueness/nullspace, and sensitivity to starting model), full-waveform techniques are rarely applied successfully in today's conventional seismic processing environment. The LSJIMP method can be abstracted in a similar qualitative fashion, using the notation of the previous sections:

$$\bf L \Leftarrow$$ velocity, reflection coefficient, crosstalk model

$$\bf d \bf Lm$$ =

Quantities like imaging velocity, the measured reflection coefficient of the multiple generators, and the crosstalk model are assumed to be fixed. Some LSJIMP implementations might depend only implicitly on velocity or reflection coefficient. For the sake of argument, however, let us assume that the LSJIMP operator, $\bf L$, is a non-linear function of these parameters, which the basic LSJIMP inversion makes no attempt to optimize. A multiple-free estimate of the primaries obviously enhances our ability to estimate imaging velocity, regardless of the method, and also permits us to model crosstalk noise from primaries that are below the onset of the seabed pure multiple. Thus the simplest nonlinear iteration of the LSJIMP method would proceed as follows, where the superscript k denotes that an operator or model vector is attached to the $k^{\rm th}$ nonlinear iteration:

$$\bold L^{[k]} \Leftarrow$$ velocity, reflection coefficients, crosstalk model

$$\bf d \bold L^{[k]} \bold m^{[k]}$$ =

$$\bold m^{[k]} \Rightarrow$$ updated velocity, crosstalk model

$$\}$$

Nonlinear updating of the reflection coefficients of the multiple generators is in general a more difficult, and potentially more valuable, problem. If a multiple generator's pure multiple is obscured by other events, the reflection coefficient estimation scheme outlined in section may produce inaccurate estimates, which the spatial regularization may not account for. For example, in the 2-D field data example shown in Chapter , the pure multiple for one multiple generator happens to be overlapped over almost the entire line by a prominent pegleg from another multiple generator.

I propose a nonlinear reflection coefficient updating scheme which obtains perturbations by fitting unmodeled events in the LSJIMP data residual, $\bold r_d$. If the initial reflection coefficient is perfect, then after convergence, $\bold r_d$ will contain only uncorrelated noise. If it is imperfect, then we also expect to see correlated events left over in the residual. Because LSJIMP separates each multiple mode indpendently in the model space ($\bold m_{i,k,m}$), we can simply apply the forward model for that mode ($\bold L_{i,k,m}$) to obtain an estimate of the particular multiple in data space, $\bold d_{i,k,m}$.

The main idea of my updating scheme is to compute a scalar update to the reflection coefficient of the $m^{\rm th}$ multiple generator, $\Delta \alpha_m$, such that  

$$\Vert \bold r_d - \Delta \alpha_m \; \bold d_{i,k,m} \Vert^2$$ (13)

is minimized. Since $\bold r_d$ may contain correlated events from primaries and other multiple modes, I suggest using $\bold d_{i,k,m}$ to compute a residual weight which is large where $\bold d_{i,k,m}$ is energetic, small where it is not. Using the absolute value of $\bold d_{i,k,m}$ is one possibility:  

$$\left\Vert \left\vert\bold d_{i,k,m}\right\vert \left( \bold r_d - \Delta \alpha_m \; \bold d_{i,k,m} \right) \right\Vert^2.$$ (14)

In the implementation of minimization () used in this thesis, $\Delta \alpha_m$ is computed independently at each midpoint location, y. The updated reflection coefficient, $\alpha_m^{[k]}(y)$, is related simply to $\Delta \alpha_m$ and the previous reflection coefficient $\alpha_m^{[k-1]}(y)$: 

$$\alpha_m^{[k]}(y) = \alpha_m^{[k-1]}(y) \left( 1 + \Delta \alpha_m \right),$$ (15)

where the superscripts again represent nonlinear iteration index.