Seismic tuning of dispersive thin layers

The wedge model is commonly used to study the limits of seismic resolution where conventionally the thickness of sub‐resolution seismic layers can be determined from thin layer tuning curves. Tuning curves relate the layer temporal thickness to the wavelet frequency, but, to our knowledge, no systematic study has been done to date of the effects of velocity dispersion and attenuation on tuning. In this work, we study the tuning properties of a thin layer dispersing according to the standard linear solid model. We show that the first tuning curve is sensitive to the attenuation, relative polarity and magnitude of the two reflection coefficients at the top and base of the layer. We provide an analytic formula for the upper bound in the mismatch between the elastic and dispersive tuning curves. We conclude that highly attenuative thin layers can appear thicker or thinner than they actually are depending on polarity and relative magnitude of reflection coefficients at the layer interfaces. Our results are particularly relevant in the quantification of CO2 where knowledge of the temporal thickness of CO2 layers is reliant on their tuning curves.


I N T RO D U C T I O N
One of the goals of seismic interpretation is to relate seismic amplitude variations to rock and fluid properties from variations in acoustic impedance. However, when geological layers are thinner than the seismic wavelength, interference from their top and base distorts the relation between seismic amplitude and acoustic impedance. Starting from the systematic approach of Widess (1973), much research has been done in the past to understand this interference by studying the tuning properties of waves encountering thin layers.
Notably properties of thin layer tuning are used in carbon monitoring for the quantification of CO 2 injected in depleted reservoirs using spectral characteristics (Williams and Chadwick, 2012;Huang et al., 2016;Ivandic et al., 2018;White et al., 2018) or time delays (Arts et al., 2004;Furre et al., 2015). These are based on earlier modelling works (De Voogd * E-mail: giorgos.papageorgiou@ed.ac.uk and Den Rooijen, 1983;Chung andLawton, 1995b, 1995a) all of which assume thin layers to be elastic and the reflection coefficients real-valued.
There is increasing evidence, however, that a range of dispersion mechanisms affect seismic wave propagation all of which are relevant in carbon monitoring. For example, patchy fluid saturation as well as wave-induced flow from large-scale fractures (see Müller et al., 2010;Caspari et al., 2011;Germán Rubino and Velis, 2011;Papageorgiou and Chapman, 2017;Falcon-Suarez et al., 2018) have been shown to affect seismic wave propagation. But the tuning properties of such dispersive thin layers have not been investigated.
In this paper, we revisit the wedge model and use it to understand tuning properties of dispersive thin layers. We assume a thin layer that disperses according to the standard linear solid model (Ursin and Toverud, 2002) Figure 1 A wedge model of linearly increasing temporal thickness and acoustic impedance Z 2 is encased within two half-spaces of acoustic impedance Z 1 , Z 3 , respectively. Varying amplitude, apparent temporal thickness and spectral characteristics of such a configuration help in the study of seismic resolution limits and tuning.
to solid-solid friction or contact-line movement (Winkler and Nur, 1982;Rozhko, 2020) in which dispersion depends on wave amplitude rather than frequency, but postpone such considerations for future work. Likewise, the wedge model we discuss here is assumed homogeneous although this is not always the case, particularly in the context of CO 2 monitoring (Glubokovskikh et al., 2020).

E L A S T I C L AY E R T U N I N G
Let us consider the properties of the wedge model depicted in Fig. 1, and its normal incidence reflection coefficient time series given by or, in the frequency domain: where the normal incidence reflection coefficients are given by the impedance contrast between consecutive layers: and index E indicates elastic modelling. The wedge model of Fig. 1 is useful for studying tuning which occurs when interference from top and base of layer 2 is maximum. Based on equation (2), the amplitude spectrum takes extreme value (minimum/maximum for opposite/same polarity reflections, respectively) when f t = n/2, n = 1, 2, . . . . The first of these, that is, defines the first tuning curve of the thin bed (Chung and Lawton, 1995b), and it has been used in seismic imaging to resolve layer thickness. For example, a zero phase Ricker wavelet with spectrum W (ω) and central angular frequency ω c = 2π f c defined by (see Wang, 2015) will have an effective frequency corresponding to the wavelength between its two troughs f eff. = π √ 6 f c and tune at temporal thickness t = √ 3/2 π fc according to (4). The spectral and amplitude characteristics of such a configuration are shown in Fig. 2 with the tuning thickness √ 3/2 π fc clearly marked to show the maximum (resp. minimum) amplitude as well as the first tuning curve. It is clear that for the reflectivity time series in equation (2) these characteristics depend solely on relative polarity between top and base of the layer.

D I S P E R S I V E L AY E R T U N I N G
For the remainder of this paper, we investigate how spectral and amplitude characteristics are affected when layer 2 of Fig. 1 is dispersive, where dispersion implies that layer 2 has a Figure 3 Shape of viscoelastic impedance Z(ω) and attenuation 1/Q derived from the standard linear solid model. Different values of the parameter α lead to different attenuation curves. The model has a complex phase that is significant within a bandwidth of roughly two orders of magnitude around the relaxation frequency 1/τ . complex-valued, frequency-dependent acoustic impedance Z 2 given by where ρ is the density of the medium and the frequency dependence is dictated by its phase velocity c(ω). Likewise, the normal incidence reflection coefficientsr 1 ,r 2 become frequency dependent and complex valued. They are derived by matching stress and strain boundary conditions at each interface (see Innanen, 2011): Withr 1 ,r 2 , we obtain the modified reflectivity series: where we use index D to indicate dispersive modelling. It is clear that equations (7) and (8) depend on the choice of the velocity model c(ω) of equation (6), which we choose to be a standard linear solid (SLS) (e.g. Ursin and Toverud, 2002) whose impedance Z 2 (ω) is the product of the elastic impedance Z 2 and a complex function where τ is a characteristic timescale and α is the ratio of unrelaxed to relaxed wave moduli (see Fig. 3). It is straightforward to see that These limits are useful if the timescale of the observation can be assumed small relative to the timescale of the model, that is if the impedance Z 2 is observed at reference frequency ω 0 = 0 and its behaviour is extrapolated to non-zero frequency by means of the SLS. Alternatively, one may choose reference frequency ω 0 = 1/τ , where attenuation peaks to and define some reference, zero-frequency real-valued acoustic impedance Z 0 such that which is used to extrapolate the impedance when ω = ω 0 . In practical terms, the former choice ω 0 = 0, involves contribution of both amplitude and phase from the SLS whereas the latter choice ω 0 = 1/τ of only phase. Here we chose the latter to isolate the potentially neglected effect of phase shift  (16) for the upper bound Q 10 at reference frequency ω 0 = 1/τ . Eight reflection coefficient pairs are marked to represent different types of reflectivity as per Chung and Lawton (1995a) but also distinguishing between positive/negative first reflection as detailed in Table 1. at the wedge interfaces. Using equation (12) and using only first orders in α − 1, the reflection coefficients becomẽ and at ω 0 = 1/τ : Inserting the reflection coefficients of equation (13) in the dispersive reflectivity series in equation (8) we find that the first tuning curve at the extreme of |R D (ω)| has the form: where φ takes the following simple form at ω = ω 0 : A reasonable upper bound for φ of equation (16), is obtained if we constrain the attenuation to 1/Q min ≤ 0.1, (1 ≤ α 1.22) which still retains a good approximation in using the first order expansion in α − 1. In Fig. 4, we show the quantity φ against the magnitude of the reflection coefficients at the top (r 1 ) and base (r 2 ) of the wedge for Q = 10 when ω 0 = 1/τ . It is apparent from Fig. 4 that the tuning of the viscoelastic wedge depends nontrivially on the relative magnitude of reflection coefficient pairs. We  Chung and Lawton (1995b). Spectral amplitude together with the shifted tuning curves for the viscoelastic wedge are shown above (when the first reflection is negative) and below (when the first reflection is positive) the amplitude plots in each figure, with a red point marking the elastic tuning thickness for the effective wavelet frequency. The filled/hollow square, cross, circle and triangle symbols refer to reflectivity pairs marked in Fig. 4. choose eight representative reflection coefficient pairs marked in Fig. 4 based on the classification of Chung and Lawton (1995b) (see also Puryear and Castagna, 2008) where reflectivity pairs belong in one of four types (all permutations of equal/unequal magnitude and same/opposite polarity).
Additionally, we also distinguish between positive/ negative polarity of the first reflection which determines whether the wavelet acquires phase retardation or advancement at the dispersive interface. It is evident from Fig. 4 that cases where the top reflection is smaller (in absolute terms) than the base reflection do not need separate consideration due to the symmetry. Together with their respective Chung and Lawton (1995b) types the following reflectivity pairs are marked on Fig. 4: • Type I: Opposite polarity-equal magnitude (white square/ circle) with r 1 = ∓0.15, r 2 = ±0.15, respectively • Type II: Equal polarity-equal magnitude (white cross/ triangle) with r 1 = ∓0.15, r 2 = ∓0.15, respectively • Type III: Opposite polarity-unequal magnitude (black square/circle) with r 1 = ∓0.15, r 2 = ±0.03, respectively • Type IV: Equal polarity-unequal magnitude (black cross/ triangle) with r 1 = ∓0.15, r 2 = ∓0.03, respectively Explicitly the tuning curves from the viscoelastic reflectivity series, as well as their amplitude characteristics from convolution with a zero-phase Ricker wavelet of equation (5) are shown in Fig. 5 together with the corresponding elastic models for comparison. Apart from Type I, the rest of the reflectivity types show measurable difference illustrating the potential for mis-estimation of layer thickness if viscoelasticity is neglected. We have calculated this error in Table 1 which represent what we consider to be an upper bound, rather than typical values. This is reflected on our choice of low-quality factor (Q 10) as well as the ratio 5:1 of top-base reflectivity for pairs of types III and IV; a higher ratio of reflectivity pairs would yield the smallest of them very difficult to observe.

D I S C U S S I O N
We found that when intra-layer multiples are taken into account, the form for φ in equation (15) is similar both in form and in magnitude to Fig. 4 but we could not find a simple analytic expression like equation (16) in this case. Our work implies that highly attenuative thin layers introduce negligible energy loss but significant phase shifts. Therefore, in order to improve layer thickness estimates, we need a measurement independent of seismic amplitude like a well log, or an estimate of the attenuation in the layer which further implies that confident knowledge of one of thickness or attenuation (and the reflection coefficients at the top and base layer) can be used to invert for the other via equation (16).

C O N C L U S I O N S
We demonstrated that tuning of dispersive layers depends on attenuation as well as the polarity and magnitude of the reflection coefficients at the top and base of the layer. We showed that traditional tuning analysis of highly dispersive thin layers (Q 10) can lead to more than 10% mis-estimation in layer thickness if dispersion is ignored. The potential for misestimation is larger if the reflection coefficients from top and base layers have different magnitudes.

AC K N OW L E D G E M E N T S
We acknowledge financial support from the Research Council of Norway (RCN grant number: 267765). We also thank Alexey Stovas and Stanislav Glubokovskikh for their constructive comments on drafts of this work, as well as Alexander Rozhko and an anonymous reviewer for their reviews. GP acknowledges financial support from the UK CCSRC (funded by the EPSRC as part of the RCUK Energy Programme) for an IRC spent at Curtin University. He is grateful for the hospitality of Boris Gurevich and Stanislav Glubokovskikh and the creative and stimulating environment they contributed to.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data sharing is not applicable to this article as no new data were created or analysed in this study.