Least-squares diffraction imaging using shaping regularization by anisotropic smoothing

Discussion

The dependency of workflow's ability to produce sharp images of edge diffractions upon the migration velocity accuracy requires further investigation. On one hand, the approach incorporates least-squares migration framework known to be quite sensitive towards velocity model (Nemeth et al., 1999). On the other hand, path-summation integral provides a velocity-model-independent weighting of the misfit, which is expected to increase the method's tolerance towards velocity model errors.

The second field example illustrates the efficiency of the proposed approach in a complex geological environment. Most edge diffractions, and especially those associated with major discontinuities, are extracted and denoised. Some reflection energy remainders are present but are limited to locations characterized by a peculiar reflection pattern, which was not picked up by the dip estimation tuned to perform reflection-diffraction separation over the whole area. This is the dip estimation problem, results of which can be improved, for instance, by subdividing the area into smaller fragments, and, thus, does not question the validity of the approach presented. The latter statement is also supported by the high performance of the developed approach when applied to the first field data example. First field data example is characterized by a low signal-to-noise ratio and is highly contaminated by the acquisition footprint. The approach allows successfully untangling edge diffractions from reflections and noise and provides high resolution images of subsurface discontinuities. We expect that in production-like environment, geological knowledge can be used to further adjust the parameter values. For instance, in the first field data example expected channel sinuosity could guide the smoothing strength for the edges. In this paper, we focus on highlighting the advantages of the method rather than on delivering final results for drilling decisions.

The natural extension of the approach is to include reflection modeling into the inversion. As demonstrated by Merzlikin et al. (2019), both reflections and diffractions can be inverted for by the same forward modeling operator whereas the separation into the components can be done on the regularization level. Regularization of diffractions can stay the same whereas reflections, for instance, can be penalized by a strong isotropic smoothing operator along the dominant local slopes: specular events locally do not exhibit lateral symmetry as opposed to edges. Extension of the model space to include reflections can help to eliminate the reflection remainders in the diffraction image domain.

High complexity of the overburden often leads to the interference between seismic events, which results in the presence of multiple dominant local slopes at a single data sample. While, in this case, the effectiveness of the proposed inversion scheme in general and of AzPWD in particular will be degraded, the performance could be improved by pre-applying migration with approximate velocity model to untangle interfering events, running AzPWD, and then going back to the original data domain.

Anisotropic smoothing is capable of emphasizing one edge direction at once. Edges with conflicting orientations can be a challenge. The ``brute force" way to tackle the challenge can be scanning for various edge diffraction orientations and picking the desired ones (Merzlikin et al., 2016). Then, inversion results with alternative orientations can be compared. At the same time, if coherent noise with a consistent spatial orientation is present in the data, it can be emphasized by anisotropic smoothing. Structure tensor orientation determination will treat this noise as signal. Poor illumination and velocity model errors can also reduce the accuracy of structure-tensor based edge diffraction orientation determination. The latter will degrade the performance of anisotropic-smoothing regularization operator. We expect the problem can be alleviated by utilizing a priori information about geologic discontinuities' orientation, e.g. by using predominant azimuths of the faults in the region extracted from a geomechanical model.

The workflow described in this paper is a 3D extension of the method proposed by Merzlikin and Fomel (2016). In two dimensions one cannot discriminate between point and edge diffractions. Distribution of scatterers is spiky and intermittent, which leads to a natural choice of sparsity constraints on the diffractivity model. The new inversion scheme is based on two regularization operators: thresholding and anisotropic smoothing. While thresholding operator imposing sparsity constraints applicable to both point and edge diffractions remains to be the same as in the 2D counterpart of the workflow, anisotropic smoothing enforces continuity along the directions picked up by the structure tensor and thus is only applicable to edge diffractions. Point diffractions, which are not elongated in space, will be smeared under anisotropic smoothing operator action. Currently, our method is biased towards edge diffractions.

Anisotropic smoothing can have spatially variable diffusion coefficient defining its strength (Weickert, 1998; Hale, 2009). For instance, the coefficient and the direction of smoothing can depend on the linearity, which can be computed as a ratio of eigenvalues of the PWD-based structure tensor and which can distinguish between the edges and regions of continuous amplitude variation (Hale, 2009; Wu, 2017). Spatially variable diffusion coefficient could help to alleviate smearing of point diffractions. Point diffractions in the diffractivity model will exhibit low linearity. For these samples smoothing power can be reduced and thresholding will be a predominant regularization operator.

The proposed approach is equivalent to total variation (TV) regularization (e.g., Strong and Chan (2003)), in which minimizing $l_1$ norm of a second derivative penalizes the model. Hessian of TV is guided by a structure tensor, which forces model smoothing to be applied along the edges with no smearing across them. Thus, TV regularization is similar to the one proposed in this paper and can also be used to penalize edge diffractions. TV implementation challenges are associated withregularization term differentiation during optimization. While this obstacle can be accommodated by using sophisticated optimization methods and representing $l_1$ -norm as a square root with a damper (Chan et al., 1999; Anagaw and Sacchi, 2012; Burstedde and Ghattas, 2009), shaping regularization by anisotropic diffusion appears to be a viable, simple to implement and fast to converge alternative with no approximations required. Anisotropic smoothing can be used to regularize full wavefield images in ``conventional" least squares migration and even in iterative velocity-model building methods.

The workflow can be utilized to extract and denoise diffractions for their subsequent depth imaging. Alternatively, a depth imaging operator can replace Kirchhoff time migration in forward modeling to allow for depth-domain model conditioning, while misfit weighting by path-summation integral is still performed in time migration domain.

The inversion can be extended to pre-stack domain. In this case, pre-stack counterparts of the chained forward modeling operators should be used: pre-stack path-summation integral (Merzlikin and Fomel, 2017), pre-stack migration engine and pre-stack AzPWD. While expressions for the former two exist, AzPWD has not been applied in the pre-stack domain. The corresponding method can be derived based on the approach developed by Taner et al. (2006).

Least-squares diffraction imaging using shaping regularization by anisotropic smoothing