pylops.avo.poststack.PoststackInversion#
- pylops.avo.poststack.PoststackInversion(data, wav, m0=None, explicit=False, simultaneous=False, epsI=None, epsR=None, dottest=False, epsRL1=None, **kwargs_solver)[source]#
Post-stack linearized seismic inversion.
Invert post-stack seismic operator to retrieve an elastic parameter of choice from band-limited seismic post-stack data. Depending on the choice of input parameters, inversion can be trace-by-trace with explicit operator or global with either explicit or linear operator.
- Parameters
- data
np.ndarray
Band-limited seismic post-stack data of size \([n_{t_0}\,(\times n_x \times n_y)]\)
- wav
np.ndarray
Wavelet in time domain (must have odd number of elements and centered to zero). If 1d, assume stationary wavelet for the entire time axis. If 2d of size \([n_{t_0} \times n_h]\) use as non-stationary wavelet
- m0
np.ndarray
, optional Background model of size \([n_{t_0}\,(\times n_x \times n_y)]\)
- explicit
bool
, optional Create a chained linear operator (
False
, preferred for large data) or aMatrixMult
linear operator with dense matrix (True
, preferred for small data)- simultaneous
bool
, optional Simultaneously invert entire data (
True
) or invert trace-by-trace (False
) when usingexplicit
operator (note that the entire data is always inverted when working with linear operator)- epsI
float
, optional Damping factor for Tikhonov regularization term
- epsR
float
, optional Damping factor for additional Laplacian regularization term
- dottest
bool
, optional Apply dot-test
- epsRL1
float
, optional Damping factor for additional blockiness regularization term
- **kwargs_solver
Arbitrary keyword arguments for
scipy.linalg.lstsq
solver (ifexplicit=True
andepsR=None
) orscipy.sparse.linalg.lsqr
solver (ifexplicit=False
and/orepsR
is notNone
)
- data
- Returns
- minv
np.ndarray
Inverted model of size \([n_{t_0}\,(\times n_x \times n_y)]\)
- datar
np.ndarray
Residual data (i.e., data - background data) of size \([n_{t_0}\,(\times n_x \times n_y)]\)
- minv
Notes
The cost function and solver used in the seismic post-stack inversion module depends on the choice of
explicit
,simultaneous
,epsI
, andepsR
parameters:explicit=True
,epsI=None
andepsR=None
: the explicit solverscipy.linalg.lstsq
is used ifsimultaneous=False
(or the iterative solverscipy.sparse.linalg.lsqr
is used ifsimultaneous=True
)explicit=True
withepsI
andepsR=None
: the regularized normal equations \(\mathbf{W}^T\mathbf{d} = (\mathbf{W}^T \mathbf{W} + \epsilon_\mathbf{I}^2 \mathbf{I}) \mathbf{AI}\) are instead fed into thescipy.linalg.lstsq
solver ifsimultaneous=False
(or the iterative solverscipy.sparse.linalg.lsqr
ifsimultaneous=True
)explicit=False
andepsR=None
: the iterative solverscipy.sparse.linalg.lsqr
is usedexplicit=False
withepsR
andepsRL1=None
: the iterative solverpylops.optimization.leastsquares.regularized_inversion
is used to solve the spatially regularized problem.explicit=False
withepsR
andepsRL1
: the iterative solverpylops.optimization.sparsity.SplitBregman
is used to solve the blockiness-promoting (in vertical direction) and spatially regularized (in additional horizontal directions) problem.
Note that the convergence of iterative solvers such as
scipy.sparse.linalg.lsqr
can be very slow for this type of operator. It is suggested to take a two steps approach with first a trace-by-trace inversion using the explicit operator, followed by a regularized global inversion using the outcome of the previous inversion as initial guess.