Theory

A time-domain seismic dataset $d(t,\textbf{g},\textbf{s})$ , where $\textbf{s}$ and $\textbf{g}$ are source and receiver vectors, can be digitized into a 3D array $\textbf{d}_{it,ig,is}$ ( $it=1,2,..,n_{t};ig=1,2,..,n_{g};is=1,2,..,n_{s}$ ), assuming there are $n_{t}$ time samples, $n_{g}$ receivers for a shot, and $n_{s}$ shots in total. Given time sampling $dt$ , the time domain array can be transformed to the frequency domain as $\tilde{\textbf{d}}_{i\omega,ig,is}$ ( $i\omega=1,2,..,n_{\omega})$ and the angular frequency sampling is $d\omega=\frac{2\pi}{n_{t}*dt}$ . In the frequency domain, only these samples that fall into the frequency band of the seismic data are kept, so for a dataset with peak frequency $f$ (frequency band $0\sim 2.5f$ ), $n_{\omega}$ can be calculated as $n_{\omega}$ = $\frac{2\pi*2.5f}{d\omega}$ .

With the frequency-selection encoding, the encoding function is defined as

$$N_{s}(i\omega,i\omega_s)=\left\{\begin{matrix} 1 & when~i\omega=i\omega_s\\ 0 & otherwise, \end{matrix}\right.$$ (31)

where $i\omega_s$ is a function of shot index $is$ , and it represents the selected frequency for the shot. Similar to conventional blended source technique, all the shots are encoded with the encoding function and blended together to form a supergather

$$\tilde{\textbf{d}}_{i\omega_s,ig}=\sum_{is=1}^{n_{s}}N_{s}(i\omega,i\omega_s)\tilde{\textbf{d}}_{i\omega,ig,is}.$$ (32)

Now the supergather $\tilde{\textbf{d}}_{i\omega_s,ig}$ becomes a 2D array, and each frequency component corresponds to a different shot. Note that a supergather can only accommondate up to $n_{\omega}$ shots. It is obvious that the frequency-selection encoding method is applicable to seismic data with a marine streamer acquisition geometry, because at each receiver position, the data components from different shots can be distinguished from one another according to their frequency contents.

In least-squares reverse time migration, a reflectivity model vector $\textbf{m}$ is sought to best fit the observed data with a Born modeling operator $\textbf{L}$ by minimizing the misfit functional

$$f(\textbf{m})=\frac{1}{2}\vert\vert\textbf{Lm}-\tilde{\textbf{d}}\vert\vert^{2}+\frac{\gamma}{2}\vert\vert\textbf{m}\vert\vert^2,$$ (33)

where $\tilde{\textbf{d}}$ is a vector representing a supergather $\tilde{\textbf{d}}_{i\omega_s,ig}$ and $\gamma$ is damping coefficient.

The following numerical scheme can be implemented with Born modeling and the reverse time migration method:

$$\textbf{g}^{(k)} = \textbf{L}^{T}[\textbf{L}(\textbf{m}^{(k)})-\tilde{\textbf{d}}]+\gamma \textbf{m}^{(k)},$$

$$\alpha = \frac{(\textbf{g}^{(k)})^{T}\textbf{g}^{(k)}}{(\textbf{Lg}^{(k)})^{T}\textbf{}\textbf{Lg}^{(k)}+\gamma \vert\vert\textbf{g}^{(k)}\vert\vert^2},$$

$$\textbf{m}^{(k+1)} = \textbf{m}^{(k)}-\alpha \textbf{g}^{(k)}.$$ (34)

At each iteration, a new supergather with new encoding functions should be used to sample a different frequency for each shot. Therefore, if $N_{\omega}$ frequencies are needed to avoid wrap-around effects, the LSM procedure should be iterated at least $N_{\omega}$ iterations to ensure that all the frequencies are visited by a shot. In contrast, the iterative stacking method can be applied to those $N_{\omega}$ supergathers to produce an image with less computational cost than conventional RTM, because usually $N_{\omega}$ can be much smaller than $n_{\omega}$ due to data reduntancy in the frequecy domain (Mulder and Plessix, 2004b).

In the following section, I demonstrate the numerical implementation of modeling and migration of a supergather with a time-domain finite-difference method.