This module provides a framework for reverse-time migration (RTM) and least-squares RTM with sparsity promotion (SP-LSRTM). We formulate our LSRTM objective function using an elastic net, which is a linear combination of \(\ell_1\) and \(\ell_2\) penalty terms: \[ \begin{equation} \begin{aligned} &\underset{\mathbf{\delta m}}{\operatorname{minimize}} \hspace{.3cm} \lambda||\mathbf{C} \mathbf{\delta m}||_1 + \frac{1}{2} ||\mathbf{C} \mathbf{\delta m}||_2^2\\ &\text{subject to: } ||\hspace{.1cm} \mathbf{J} \delta \mathbf{m} - \delta \mathbf{d}||_2 \leq \sigma. \end{aligned} \label{objective} \end{equation} \] The operator \(\mathbf{J}\) is the linearized Born modeling operator, \(\mathbf{\delta d}\) is the observed reflection data and \(\mathbf{\delta m}\) is the model perturbation (seismic image). \(\mathbf{C}\) is the forward Curvelet transform, \(\sigma\) is size of the noise of the seismic data and \(\lambda\) is a parameter that controls the trade-off between the \(\ell_1\) and the \(\ell_2\)term in the objective function. By formulating our objective this way, we can use the linearized Bregman method to solve this problem, which allows us to work with subsets of shots and rows of \(\mathbf{J}\) in each iteration, thus significantly reducing the cost of LSRTM in terms of PDE solves.
To run sparsity-promoting LSRTM with the linearized Bregman method, the user needs to set up the Jacobian \(\mathbf{J}\) and the observed data \(\mathbf{\delta d}\) as a data container. Refer to the documentation in /Modeling/TimeDomain for a detailed description of how to set up these operators and structures. The linearized Bregman method for SP-LSRTM can then be called by passing `\(\mathbf{J}\) and \(\mathbf{\delta a}\) to the solver:
# Setup operators
F = joModeling(info,model,srcGeometry,recGeometry) # modeling operator
q = joData(srcGeometry,wavelet) # source container
J = joJacobian(F,q) # Setup Born modeling operator
dD = joData(recGeometry,data) # observed reflection data
# SP-LSRTM w/ linearized Bregman
image = linbreg_lsrtm(J,dD;iterations=10,shots_per_iteration=2)The function linbreg_lsrtm requires only J and the observed data dD as input arguments. The number of iterations and the number of shots per iterations can be passed as optional arguments. Further optional arguments are the thresholding parameter lambda and a file name for saving snap shots of each iteration:
image = linbreg_lsrtm(J,dD;iterations=10,shots_per_iteration=2,lambda=1,snapshots="/path/to/snapshots/filename")