May 2007 | Volume 32 Issue 05 | View IssueClick Here for a PDF of this Article
Introduction
Horizontal transverse isotropy (HTI) is a widely used method to describe reservoirs with vertical fracture s (Thomsen, 1998). Anisotropic effects on moveout yield abnormal NMO velocities unless the seismic line is parallel to the fractures. Here we consider a subsurface velocity model with horizontal layers. Each of the layers may include a set of aligned, vertical, circular flat cracks. We use HTI anisotropy with a horizontal axis of rotational symmetry to describe azimuthal anisotropy. These horizontal axes may have different orientations for different layers. The anisotropic co-efficients, determined from P-wave moveout, can be combined to estimate crack parameters. The HTI model provides not just anisotropic coefficients but also true vertical velocities (Tsvankin, 1997, Grechka et al. 1999, Contreras et al., 1999). It implies that azimuthally dependent normal moveout (NMO) velocities (three parameters of NMO velocity ellipse) acquired over a layered HTI subsurface can be used to build a depth velocity model and perform anisotropic depth imaging without using any additional information. For a vertically transverse isotropic layered model (VTI), to estimate vertical velocities, we need additional information. It might be joint P and S inversion (Blias, 1983) or additional well information about vertical velocity (Tsvankin and Thomsen, 1995). Blias (1989) derived explicit approximate formulas for the zero-offset time and azimuthally-dependent NMO velocity for a slightly tilted symmetry axis in VTI layered media. To describe fractured reservoirs, orthorhombic symmetry is believed to fit the best. The main reason to use orthorhombic symmetry in sedimentary basins is that this anisotropy describes the verticallytransverse layer (composed of the finally layered sediments) with embedded parallel vertical cracks (Wild and Crampin, 1991, Schoenberg and Helbig, 1997).
For a horizontally stratified anisotropic subsurface above the dipping reflector, the formula for NMO velocity was derived by Grechka et al. (1999b, formula (28)). Their formula includes interval 2×2 matrices Wl-1 (l is the interval number) in terms of the component of the slowness vector. These matrices describe the NMO ellipses corresponding to the individual layers. To compute the azimuthally-dependent NMO velocity, we have to calculate an inverse matrix for each matrix Wl and sum them up with the interval times as scalars.
For a model composed of a stack of homogeneous HTI layers, I will derive an explicit formula for the azimuthally-dependant NMO velocity. This formula does not include matrices. This provides more convenient visual analysis and enables us to analytically obtain uncertainties in the estimation of anisotropic parameters. To derive such a formula we use an approach designed for a vertically transverse layered medium (Blias, 1983, 1987) and based on a Taylor series expansion for phase and interval velocities with respect to phase and group angle. These formulae lead to an explicit Dix-type analytical NMO inversion for a horizontally-layered HTI subsurface (no raytracing needed) that can be easily implemented into any program code. This approach is extended to the orthorhombic layered subsurface with a horizontal plane of symmetry.
VTI horizontally layered medium
For VTI horizontally-stratified subsurface, an explicit presentation of reflection traveltime for P and SV waves was derived by Blias (1983):
EQUATION 01
where t0 is a vertical time, VNMO is NMO velocity and the coefficient c2 describes the non-hyperbolic part of the NMO function:
EQUATION 02
Here hk is the thickness of the k-th layer; vk is a vertical velocity in the layer. Coefficients qk and Bk describe the influence of VTI anisotropy on reflected traveltimes of P or SV waves. For SH wave, formulas (2) and (3) are the same with B = 1 and q = (VSH,hor/VSH,vert)2, where Vhor and Vvert are horizontal and vertical velocities. In terms of elastic modulus Cαβ, the formula for NMO velocity for P and SV waves was derived by Lyakhovitskiy and Nevskiy (1971). A short English review of their results can be found in Lyakhovitskiy (1981). They called coefficients
EQUATION
“apparent coefficients of anisotropy”, and defined them as the ratio of the horizontal and vertical axes of the ellipses which approximate the group velocity for small angles (< 20 – 250). They derived formulae in terms of the elastic modulus and showed that for both P and SV NMO velocities
EQUATION
where VDix is the interval velocity estimation from Dix formula.
Using the Thomsen parameters ε and δ (Thomsen, 1986) we can write for the apparent coefficients of anisotropy:
EQUATION
where γ = Vs,Vert/VP,Vert , Vvert is a vertical velocity. For the SH wave, k̃SH = VSH,Hor/VSH,Vert where Vhor is a horizontal velocity, qP and qSV are coefficients in formulas (2) for P and SV waves correspondingly. Analyzing thin-layered two-component media, Lyakhovitskiy and Nevskiy suggested four types of P anisotropy according to the difference in velocities VS and VP in these components, (Fig. 1). Possible values for the coefficients qk and Bk from formula (2) for different types of anisotropy were determined in Blias (1989, table 1). For a VTI medium, for both P and SV waves, coefficients Bk may differ significantly from unity, which leads to a significant difference between the NMO function and its hyperbolic approximation even for a single-layer model. For SH waves, Bk = 1 and traveltimes behave similarly to the isotropic subsurface.
EQUATION 03
After substituting (3) into (2), we find coefficients t0, VNMO and c2 for a new model :
EQUATION 04
Formulas (4) show that the three-term NMO function depends on three velocity parameters: velocity uk, coefficient Bk and thickness dk. For each layer, we have three unknown parameters: dk, uk and Bk; equations (4) do not include the coefficient qk. Consequently, we cannot determine coefficients qk from equations (2) or, by implication, from three-term velocity analysis. Formulas (3) and (4) show that we cannot determine true depth (as a sum of the thicknesses hk) and vertical interval velocities vk but only the products
EQUATION
This implies that we need additional information to determine interval vertical velocities. It may be well information about vertical velocities or boundaries. Using the Thomsen parameter δ and formula q = 1 + 2δ, we conclude that parameter δ (the same as parameter q) cannot be restored from three-term P wave velocity analysis. This fact was discussed in Tsvankin and Thomsen (1995). The same result holds for the VTI model with gently curvilinear boundaries (Blias, 1987).
Coefficients q and B can be expressed in terms of the group velocity ν(ψ) as a function of group angle ψ and its derivatives in the vertical direction (Blias, 1983, 1987):
EQUATION
To calculate the derivatives of group velocity v with respect to group angle ψ, I used the differentiation rule of the implicit function. Formulas for the coefficients Bk were derived in Blias (1983) in terms of effective anisotropy coefficients. Unfortunately, there were some errata in this paper. The same formulas were provided later (Blias 1987, 1989):
EQUATION
where γ = VS,vert/VP,vert, k̃p and k̃SV are apparent anisotropy coefficients introduced by Lyakhovitskiy and Nevskiy (1971). Using Thomsen’s parameters ε and δ (Thomsen, 1986), formulas for the coefficients qk and Bk for P waves can be written as:
EQUATION 05
where
EQUATION
and η is a parameter introduced by Alkhalifah and Tsvankin (1995):
EQUATION
Formulas (5) show that B = 1 if η = 0. Since for elliptical anisotropy in the k-th, layer Bk = 1 (Blias, 1983, 1987, 1989), parameter B (the same as parameter η) is responsible for the nonelliptical part of anisotropy. For weak anisotropy, when we take into account only linear terms with respect to ε and δ, formulas for BP and BSV are much simpler:
EQUATION
Formulas for the higher order terms in Taylor series (1) were derived by Ursin and Stovas (2006). An approximation to the group velocity as a function of the group angle ψ was suggested in Blias (1989) in a form:
EQUATION
The first term is an ellipse; the second term is a non-elliptical correction. This approximation usually provides an accurate description of the group velocity in the interval [00 to 450].
In Blias (1983), formulas (2) were used to solve the inverse problem in a layered VTI subsurface, using three coefficients from the P reflected traveltime Taylor series (t0, VNMO and c2) and the two first coefficients for SV waves (t0 and VNMO). In that paper, a Dix type inversion was designed for a multilayered VTI media to find the anisotropic parameters for each layer. The inversion problem was reduced to a quadratic equation with respect to
EQUATION
A similar approach for joint traveltime inversion of P and SV waves, which also uses three coefficients of the P wave Taylor series and two coefficients from the S wave expansion, and leads to a quadratic equation with respect to V2P ,vert' was described in Ursin and Stovas (2005).
Formulas (1) and (2) provide different NMO approximations, similar to the isotropic model. For example, we can write two fractional approximations to the NMO function (Blias, 2007, formulas (10) and (13)):
EQUATION
where VNMO is determined by the second formula (2); S is a dimensionless parameter, defined by the formula:
EQUATION
For the isotropic model, formulas (t) show higher accuracy compared to formula (1) or to a shifted hyperbola approximation (Malovichko, 1978, Blias, 1982):
EQUATION
We can expect that for anisotropic media, formulas (t) also have higher accuracy. Numerical investigation of these formulas can be performed similarly to Blias (2007) but it is beyond the scope of this paper.
Comparing formulas (4) with formulas for an isotropic medium (Taner and Koehler, 1969), we see that the only difference between them is the coefficient Bk. If the phase velocity (and, consequently, group velocity) is elliptical, then Bk = 1 (Blias, 1983, 1987, 1989). This implies that coefficient Bk is responsible for nonelliptical anisotropy, similarly to parameter η. Formulas (3) and (4) also demonstrate that for the horizontally-layered model with elliptical anisotropy, reflected traveltimes are the same as for the isotropic model (Bk = 1) with thicknesses dk and interval velocities uk. This leads to the fact that we cannot distinguish a VTI subsurface from the isotropic velocity model using only conventional hyperbolic velocity analysis for P waves and SH waves. For a conventional offset to depth ratio, close to 1, the hyperbolic approximation of reflected traveltimes for the isotropic model is quite accurate. The final of the set of formula (4) shows that anisotropy-induced deviation in the moveout from hyperbolic (additionally to the influence of vertical heterogeneity) depends on the magnitudes of the coefficients Bk, that is, on the nonhyperbolic part of the velocity and the non-elliptic part of the phase/ group velocity. For the isotropic (or elliptically anisotropic) model (Bk = 1) the coefficient c2 is always negative, and the RMS velocity is larger than the average velocity (Al- Chalabi, 1974). In a VTI horizontally-stratified medium, the sign of c2 depends on the values of Bk. If Bk >1 then c2 is negative (ηk > 0) and the deviation from hyperbolic moveout increases. In this case (corresponding to the most common positive values of the difference εk – δk, (Tsvankin and Thomsen, 1994)) vertical heterogeneity and anisotropy work together to increase the non-hyperbolic component of the moveout. If Bk <1 (ηk < 0) then the anisotropy-induced part of the coefficient c2 is positive. Vertical heterogeneity and anisotropy lead c2 to the opposite direction and the non-hyperbolic moveout deviation becomes smaller. Contrary to isotropic vertically heterogeneous subsurface, in a VTI medium, the coefficient c2 can be positive when coefficients Bk < 1. Consequently, for a VTI horizontally-layered media, the RMS velocity may be less than the average velocity.
NMO velocity ellipse
Author Biography
Emil Blias received his M. Sc. (Mathematics) from Samara University in 1974, Ph.D. in Geophysics from Novosibirsk Research Geophysical Institute (Russian Academy of Science) in 1984 and Ph.D. in Mathematics and Physics from St.-Petersburg University in 1990 (all Russia). He began his geophysical career in 1977 in Samara Geological Research Institute as a re s e a rch scientist. In 1986 he moved to Murmansk and became a chief of a Laboratory of Marine Research Geophysical Institute, where he was leading R&D for VSP processing and interpretation. In 1994 – 1998 he served as a full professor of Mathematics in Murmansk Technical University. In 1999 he joined Scott Pickford (later Core Laboratories) as a Research Scientist. In Core Labs, he has developed 3-D VSP software and invented a new method for enhancing surface seismic data using downgoing VSP wavefield, now knowing as HFR (High Frequency Restoration). From August 2005, E. Blias is with VSFusion, doing R&D in 3-D VSP data processing and interpretation.
Emil derived an analytical solution for the joint (P and S) Dix type inversion in VTI layered medium in 1983. He established the analytical connection between NMO velocities and second lateral derivatives of the boundaries and interval velocities, which explains anomalous fluctuations of stacking velocities. This analytical result forms the base for the processing non-first-break technology in the presence of unresolved shallow velocity anomalies. He designed a new wave-by-wave method for 3-C VSP wavefield separation, as well as concepts, algorithms and software for 3-D VSP data processing.
His field of interest includes seismic velocities, velocity analysis and depth velocity model building, shallow velocity anomaly problem, anisotropy, AVOInversion and VSP (including 3-D VSP) data processing and interpretation. Emil Blias has over 100 publications and is an author of three university text books in mathematics and signal processing. He is a member of SEG, CSEG and EAGE.