叠前各向异性系数反演及在裂缝预测中的应用

图1 各向同性和HTI介质模型纵波反射系数随方位角变化关系

Fig.1 Variation of P-wave reflection coefficient with azimuth in isotropic and HTI media models

为进一步研究各向异性系数与地震波场特征的关系,采用入射角及其他参数不变,改变各向异性参数取值范围和方位取值范围的方式进行模拟计算,并对各个方位下3个各向异性参数所引起的反射系数的变化规律进行分析。研究结果表明,不同的各向异性参数对纵波反射系数影响不同,其中,各向异性参数ε对反射系数的影响最大,即引起的反射系数的变化最大(图2)。

图2

图2 HTI介质模型纵波反射系数随各向异性参数变化关系

Fig.2 Variation of P-wave reflection coefficient of HTI medium model with anisotropic parameters

图2给出了方位角为45°时,在固定入射角及其他参数不变条件下,单一改变各向异性参数取值范围所得到反射系数随各向异性系数变化的结果。当研究各向异性参数中δ与反射系数的变化关系时,单一改变各向异性参数δ取值范围,使其从-1 ~1每隔0.2取值一次,γ、ε为0.25并保持不变。同理,研究其他各向异性参数与反射系数之间的关系时,固定另外两个各向异性参数值为0.25不变,使用Rüger公式进行反射系数计算。其他各参数取值见表1。

表1 HTI介质纵波反射系数随各向异性参数变化正演模型参数

Table 1 Variation of P-wave reflection coefficient of HTI medium with anisotropic parameters

参数	v_p/(m·s^-1)	v_s/(m·s^-1)	ρ	φ	θ
Layer 1	2000	1253	2.2	45	30
Layer 2	2337	1437	2.3	45	30

新窗口打开| 下载CSV

根据以上研究结果,可以得到几点认识:①由裂缝诱导的各向异性可以引起地震振幅的方位各向异性,基于宽方位地震资料可以反映这种各向异性;②用于表征高角度裂缝模型的HTI介质的各向异性参数(ε、δ、γ)中,γ对反射系数的影响最大,各向异性越强,表明裂缝越发育。因此,可以利用宽方位资料反演各向异性参数γ,从而得到裂缝发育程度的反演结果。

为了得到反射系数与方位角之间的定量关系,在忽略式(1)中的高阶项sin²θtan²θ后,可以将式(1)化简为截距项和梯度项,其中梯度项可以表示为各向同性和方位各向异性两项^[8]:

(2)

R (ϕ, θ) = F + G (ϕ) s i n^{2} θ,

(3)

F = \frac{1}{2} \frac{Δ Z}{\bar{Z}},

(4)

G = G_{i s o} + G_{a n i} c o s^{2} ϕ,

(5)

G_{i s o} = \frac{1}{2} [\frac{Δ α}{\bar{α}} - {(\frac{2 Δ β}{\bar{α}})}^{2} \frac{Δ G}{\bar{G}}],

(6)

G_{a n i} = \frac{1}{2} [Δ δ + 2 {(\frac{2 \bar{β}}{\bar{α}})}^{2} Δ γ] 。

式中:F为截距,与纵波速度有关;G(ϕ)为与方位角有关的梯度项,且可分为各向同性梯度项G_iso和各向异性梯度项G_ani。其他参数含义与式(1)中的一致。

2 各向异性系数反演的方法原理

由弱各向异性近似理论可知,各向同性背景项和各向异性项之和可以用来描述各向异性介质中的反射系数^[9]。因而,反射系数公式可以进一步化简为:

(7)

\begin{array}{l} R_{p} (θ, φ) = R_{p}^{i s o} (θ) + R_{p}^{a n i} (θ, φ), \\ R_{p}^{i s o} (θ, ϕ) = \frac{1}{2} \frac{Δ α}{α} + \frac{1}{2} [\frac{Δ α}{\bar{α}} - {(2 \frac{\bar{β}}{\bar{α}})}^{2} \frac{Δ G}{\bar{G}}] s i n^{2} θ + \\ \frac{1}{2} (\frac{Δ α}{\bar{α}}) s i n^{2} θ t a n^{2} θ, \\ R_{p}^{a n i} (θ, ϕ) = \frac{1}{2} \frac{Δ α}{α} + \\ \frac{1}{2} \{[Δ δ^{(v)} + 2 (2 \frac{Δ \bar{β}}{\bar{α}}) Δ γ^{(v)}] c o s^{2} ϕ\} s i n^{2} θ + \\ \frac{1}{2} [Δ ε^{(v)} c o s^{4} ϕ + Δ δ^{(v)} s i n^{2} ϕ c o s^{2} ϕ] s i n^{2} θ t a n^{2} θ 。 \end{array}

式中: $R_{p}^{i s o}$ (θ)为Rüger公式各向同性项^[10]; $R_{p}^{a n i}$ (θ,φ)为由(式(2)~(6))根据Rüger推导HTI介质各向同性项而推导出的各向异性项,其过程为在推导各向异性项时忽略各向同性项。

通过对Rüger 公式中各向同性项和各向异性项的研究,一些研究者提出了适用于不同情况下的 HTI 介质反射系数计算公式,来满足当下对实际反演中弱各向异性近似理论的需求,其中Aki-Richard 和 Fatti 近似式在各向同性介质中的应用最为广泛。Aki-Richard 近似式反映了纵、横波速度反射率和密度反射率二者与反射系数之间的关系,Fatti 近似式则反映了纵、横波阻抗反射率和密度反射率与反射系数的关系:

Aki-Richard 近似式

(8)

\begin{array}{l} R_{p}^{i s o} (θ) = \frac{1}{2} (1 - 4 g s i n^{2} θ) \frac{Δ ρ}{\bar{ρ}} + \\ \frac{1}{2} s e c^{2} θ \frac{Δ α}{\bar{α}} - 4 g s i n^{2} θ \frac{Δ β}{\bar{β}}, \end{array}

Fatti 近似式

(9)

\begin{array}{l} R_{p}^{i s o} (θ) = \frac{1}{2} s e c^{2} θ \frac{Δ I_{p}}{\bar{I_{p}}} - 4 g s i n^{2} θ \frac{Δ I_{s}}{\bar{I_{s}}} - \\ (\frac{1}{2} t a n^{2} θ - 2 g s i n^{2} θ) \frac{Δ ρ}{\bar{ρ}}, \end{array}

式中:g= $\frac{v_{s}^{2}}{v_{p}^{2}}$ ,α为纵波速度,β为横波速度,θ为入射角,ρ为密度, $\frac{Δ I_{p}}{\bar{I_{p}}}$ 为纵波阻抗反射系数, $\frac{Δ I_{s}}{\bar{I_{s}}}$ 为横波阻抗反射系数。

通过对比式(8)、式(9)与式(1)中的各向同性项可以发现,Fatti 各向同性表达式更适合在较大入射角情况下的反射系数计算。本文通过应用 Gardner 提出的密度和纵波速度的经验公式,对式(8)中的第三项进行化简,达到提高反演稳定性和减少误差的目的,弥补Fatti近似式中为达到在反演过程中减少未知数的个数而忽略密度项使得对较大入射角范围的反射系数计算会出现明显误差时的不足。具体推导过程如下:

(10)

ρ = c v_{p}^{k},

式中:c为常数,k为幂指数;可根据研究区测井资料中纵波速度和密度之间的关系进行拟合来确定 c 和 k 的值。

对式(10)两边进行求导,可得密度和纵波速度之间的关系:

(11)

\frac{Δ ρ}{\bar{ρ}} = L \frac{Δ α ρ}{\bar{α}},

式中,L=c×k。假设目的层处密度变化不大,则式(11)可以改写为:

(12)

\frac{Δ ρ}{\bar{ρ}} = L \frac{Δ α ρ}{\bar{α ρ}} = L \frac{Δ I_{p}}{\bar{I_{p}}},

将式(12)代入到Fatti 公式中的密度项进行化简得:

(13)

\begin{array}{l} R_{p}^{i s o} (θ) = (\frac{1}{2} s e c^{2} θ - \frac{L}{2} t a n^{2} θ + 2 g L s i n^{2} θ) \frac{Δ I_{p}}{\bar{I_{p}}} - \\ 4 g s i n^{2} \frac{Δ I_{s}}{\bar{I_{s}}}, \end{array}

上式为消除了密度项影响后,基于Fatti 公式的反射系数近似计算公式,可用于替换Rüger 公式中各向同性项的表达。

对于Rüger 公式中的各向异性项反射系数计算,常采用的方法是在不考虑高阶项影响的情况下进行小入射角假设,以此实现各向异性的梯度反演。但在现实情况下,裂缝的发育导致了方位各向异性现象更多地出现在中远偏移距(较大入射角处),故此方法将导致 Rüger 公式在进行各向异性项反射系数化简时会产生一定的偏差,损失裂缝预测的精度。为此,本文通过应用对裂缝中含有流体时给出的各向异性参数之间的关系^[10],对 Rüger公式中的各向异性参数项进行化简,得到新的各向异性项表达式。

当裂缝中含有饱和流体时,Bakulin 给出了各向异性参数之间的关系式:

(14)

ξ^{(v)} = 0, γ^{(v)} = \frac{δ^{(v)}}{4 g} 。

将式(14)代入到Rüger 公式中的各向异性项,并保留高阶项的影响,得到

(15)

\begin{array}{l} R_{p}^{a n i} (θ, φ) = (6 g s i n^{2} θ c o s^{2} φ + \\ 2 g s i n^{2} φ c o s^{2} φ s i n^{2} θ t a n^{2} θ) Δ γ^{(v)} 。 \end{array}

利用式(15)可对Δγ⁽^v⁾进行各向异性表征,其中Δγ⁽^v⁾可以直接由地震数据反演得到,进一步减少了各向异性参数反演的未知量个数和计算所需的工作量。

结合式(13)与式(15),可得基于研究区含流体裂缝型储层新的反射系数随入射角和方位角的计算公式:

(16)

\begin{array}{l} R_{p} (θ, φ) = \\ (\frac{1}{2} s e c^{2} θ - \frac{L}{2} t a n^{2} θ + 2 g L s i n^{2} θ) \frac{Δ I_{p}}{\bar{I_{p}}} - 4 g s i n^{2} θ \frac{Δ I_{s}}{\bar{I_{s}}} + \\ (6 g s i n^{2} θ c o s^{2} φ + 2 g s i n^{2} φ c o s^{2} θ s i n^{2} θ t a n^{2} θ), \end{array}

式中,L为根据Gardner 公式对研究区拟合得到的常数。

基于式(16),根据N道OVT域叠前地震道集,第i道叠前地震资料对应的测线方位角为ϕ_i,入射角为θ_i,实际地震道集记录振幅为R(θ_i,ϕ_i)。N道地震记录中反射振幅与纵横波阻抗反射系数和各向异性参数之间的表达式为:

(17)

\{\begin{array}{l} A m = B, \\ A = [\begin{array}{l} \frac{1}{2} s e c^{2} θ_{1} - \frac{L}{2} t a n^{2} θ_{1} + 2 g L s i n^{2} θ_{1} & - 4 g s i n^{2} θ_{1} & 6 g s i n^{2} θ_{1} c o s^{2} φ_{1} + 2 g s i n^{2} φ_{1} c o s^{2} φ_{1} s i n^{2} θ_{1} t a n^{2} θ_{1} \\ \frac{1}{2} s e c^{2} θ_{2} - \frac{L}{2} t a n^{2} θ_{2} + 2 g L s i n^{2} θ_{2} & - 4 g s i n^{2} θ_{2} & 6 g s i n^{2} θ_{2} c o s^{2} φ_{2} + 2 g s i n^{2} φ_{2} c o s^{2} φ_{2} s i n^{2} θ_{2} t a n^{2} θ_{2} \\ \frac{1}{2} s e c^{2} θ_{3} - \frac{L}{2} t a n^{2} θ_{3} + 2 g L s i n^{2} θ_{3} & - 4 g s i n^{2} θ_{3} & 6 g s i n^{2} θ_{3} c o s^{2} φ_{3} + 2 g s i n^{2} φ_{3} c o s^{2} φ_{3} s i n^{2} θ_{3} t a n^{2} θ_{3} \end{array}], \\ m^{T} = [\begin{array}{l} \frac{Δ I_{p}}{\bar{I_{p}}} & \frac{Δ I_{s}}{\bar{I_{s}}} & Δ γ^{(v)} \end{array}], B = [\begin{array}{l} R_{1} \\ R_{2} \\ R_{3} \end{array}] 。 \end{array}

对实际数据进行处理时,会存在一个由裂缝对称面与测线初始方位的夹角φ=ϕ-φ_sym,其中φ_sym能够根据实际测井和地质资料等先验信息进行估计。

从式(17)可以看出,只需要3个不同方位角和入射角的叠前道集,就可以求解此方程,得到每一反射界面上的纵波波阻抗反射率ΔI_p/ $\bar{I_{p}}$ 、横波波阻抗反射率ΔI_s/ $\bar{I_{s}}$ ,以及各向异性系数差值Δγ。其中Δγ值越大,表示各向异性越强,裂缝越发育。但对于实际 OVT域宽方位叠前道集数据来说^[11],每个CMP点覆盖次数非常高,不同的道具有不同的入射角和方位角信息,远远大于实际解方程所需要的方程数。因此方程(17)的求解过程是一个超定方程的求解过程。本文根据实际资料选取适当的超定解,一般为15个左右,然后利用最小二乘法求最优化解作为方程(17)最后的解。图3给出了本文提出的各向异性系数反演的工作流程。

图3

图3 各向异性系数反演流程

Fig.3 Flow chart of anisotropic coefficient inversion

3 模型实验

3.1 新公式精确性分析

在进行模型实验时,尽量采用了研究区的已知信息,主要是已钻井的测井数据和裂缝发育情况等,建立相关HTI模型,其模型参数如表2所示。

表2 根据研究区资料所建模型参数

Table 2 According to the data of the study area, the model parameters are established

v_p/(m·s^-1)	v_s/(m·s^-1)	ρ/(g·cm^-3)	ξ	δ	γ	θ	φ
4000	2300	2.8	0	0	0	0~60	0~180
5500	3100	2.65	0.05	0.25	0.2

新窗口打开| 下载CSV

采用表2中的参数分别研究新公式(式(16))与其他不同反射系数计算方程在各向同性项随入射角的变化情况和各向异性项随入射角和方位角的变化情况,以及新公式(16)相比Rüger公式和其他常用的反射系数计算方程,得到的反射系数对比结果见图4和图5。

图4

图4 不同方程各向同性计算结果随入射角变化对比

Fig.4 Comparison of isotropic calculation results of different equations with incident angle

图5

图5 不同方程各向异性项随入射角(a)和方位角(b)变化对比

Fig.5 Comparison of anisotropic terms of different equations with incident angle(a) and azimuth(b)

图4分别为Aki-Richard、Fatti三项式和常用简化后的Fatti两项式、Rüger方程各向同性项以及基于Fatti-Gardner公式推导的新方程在各向同性项下固定方位角计算得到的反射系数随入射角的变化关系。可以看出,在小入射角情况下各方程计算得到的各向同性项反射系数基本相同,但随着入射角的增大,根据式(16)计算的反射系数精度更高,其结果与Rüger公式的计算结果更加接近,表明了公式的正确性。

图5为考虑了各向异性项的各反射系数方程随入射角和方位角变化的对比。从图中可以看出,不同入射角和方位角情况下,式(16)在中远偏移距(大入射角)情况下较常规使用的Fatti两项式具有更高的精度,说明该方程可用于实际资料求取反射系数。

3.2 基于实际叠前道集资料的模型反演

模型建立基于研究区BZ-B井点实际叠前地震资料,将已有的宽方位叠前道集沿目的层进行不同的方位道集提取,得到不同方位角情况下的合成数据。

1)通过时间域的波阻抗反射系数曲线和各向异性参数值曲线采用新推导的式(16)进行反射系数计算,得到基于研究区BZ-B井处实际资料不同入射角(0~45°)和方位角(30°,60°,90°)利用40 Hz的雷克子波进行褶积得到的合成地震记录(图6)。

图6

图6 研究区对应井目的层处不同方位角合成叠前道集

Fig.6 Prestack gathers are synthesized at different azimuth angles at the target layer of the corresponding well in the study area

2)选取BZ研究区内典型井BZ-B,根据已有的井资料提取其目的层处的纵横波速度和随机产生的各向异性参数,利用波阻抗反射系数计算公式得到对应的波阻抗反射系数曲线,得到基于BZ研究区实际测井资料的纵横波阻抗反射系数和各向异性参数γ⁽^v⁾值曲线。利用以上合成数据,根据反演方程式(17)和反演流程(图3),得到基于BZ-B井处的模型反演结果。图7为反演得到的纵、横波阻抗反射系数和各向异性参数曲线反演结果与理论结果对比。

图7

图7 基于最小二乘法的BZ-B井处理论结果与反演结果

a—纵波阻抗反射系数;b—横波阻抗反射系数;c—各向异性参数γ

Fig.7 Anisotropic parameter curve inversion results theoretical results and inversion results at well BZ-B based on least square method

a—longitudinal wave impedance reflection coefficient;b—transverse wave impedance reflection coefficient;c—anisotropic parameter γ

以上结果证明了基于新推导的式(17)进行各向异性系数反演的正确性,为实际资料的应用奠定了基础。

4 应用效果

将本文提出的相关方法应用到渤海某地太古宇潜山裂缝型储层的预测中,以检验方法的适用性。研究区的目的层为深埋(平均埋深5 000 m)潜山储层,断缝系统十分发育并且是优质储层的主控因素。所用实际地震资料为该地区新采集的横纵比为0.8的宽方位地震资料。图8展示了研究区的一条典型的联井地震剖面,图中3 000~3 500 ms间的红色曲线为太古界潜山顶面目的层的层位解释结果,剖面上标注了根据现有钻探成果对地震资料的解释结果。主要勘探结论是:BZ-B、A、C及D井断缝系统较为发育,为优质储层发育区,而BZ-E井断缝系统及储层发育不好。

图8

图8 研究区地震连井剖面

Fig.8 Seismic well connection profile in the study area

图9为用商业软件进行AVAZ反演得到的同一剖面的断缝系统反演结果,图中3 000~3 500 ms间的红色区域代表异常振幅带,表示断缝系统发育区。从图中可以看到反演结果与实际钻探结果吻合较好(图中蓝色区域)。

图9

图9 AVAZ叠前裂缝预测连井线剖面

Fig.9 AVAZ prestack fracture prediction well connection profile

图10为利用本文提出的方法得到的反演结果,图中3 000~3 500 ms的红紫色区域代表各向异性强的区域,说明断缝系统较为发育。对比图8和图10可以看到,本文所提出的方法所得到的反演结果与实际钻探结果较为吻合,并且好于AVAZ反演结果(图中蓝色区域),从而证明了所提方法的正确性和适用性。

图10

图10 叠前各向异性参数反演裂缝预测连井线剖面

Fig.10 Prestack anisotropic parameter inversion fracture prediction well connection profile

需要指出的是,利用本方法按照图1的反演流程进行实际资料反演时,需要考虑如何更好地利用宽方位资料对式(17)的超定方程进行求解。由于实际宽方位叠前地震资料数据量巨大,为了进一步减少计算量,同时实现基于不同方位的实际资料反演,可采取以下方法进行:①在实际进行反演的时候,通过对叠前宽方位道集进行不同方位角和入射角的叠加^[12];②将每个叠加数据体对应的中心入射角和方位角作为输入数据的方位角和入射角,用于实际资料的反演。本次实际资料反演中,选用了4个方位叠加道集和4个入射角叠加道集进行反演,具体反演选用的实际资料入射角和方位角划分情况如表3所示。

表3 实际资料反演选用的入射角和方位角划分情况

Table 3 Division of incident angle and azimuth selected for actual data inversion

数据体编号	入射角范围/(°)	中心入射角/(°)	方位角范围/(°)	中心方位角/(°)
1	10~30	15	30~60	45
2	20~30	25	90~120	105
3	20~40	30	120~150	135
4	30~40	35	150~180	165

新窗口打开| 下载CSV

5 结论

通过研究HTI介质模型下的地震波场随方位的变化特征,总结归纳了各向异性系数与反射系数的关系及其随方位变化的规律。结合Gardner公式和Bakulin等提出的含流体裂缝各向异性参数之间的相互关系,对Rüger公式进行了近似简化,推导了一个新的反射系数与各向异性系数γ及方位角和入射角间的定量公式。基于裂缝诱导的各向异性理论和反射系数随方位变化的相关理论,提出了一种基于各向异性参数反演的裂缝预测方法,通过理论模型和实际资料应用证明了方法的有效性和适用性,为应用宽方位叠前地震资料进行裂缝预测提供了一种新的方法。在研究区变质岩潜山储层裂缝预测中取得了较好的应用效果,为研究区变质岩潜山裂缝型储层的预测提供了一种可行适用的方法技术。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

关宝文, 郭建明, 杨燕, 等.

油气储层裂缝预测方法及发展趋势

[J]. 特种油气藏, 2014, 21(1):12-17.

Guan

B W

, Guo

J M

, Yang

, et al.

Methods of fracture prediction in oil & gas reservoirs and their development trend

[J]. Special Oil & Gas Reservoirs, 2014, 21(1):12-17.

[2]

张昕, 郑晓东.

裂缝发育带地震识别预测技术研究进展

[J]. 石油地球物理勘探, 2005, 40(6):724-731.

Zhang

, Zheng

X D

Research progress of seismic identification and prediction technology in fracture zone

[J]. Oil Geophysical Prospecting, 2005, 40(6):724-731.

[3]

印兴耀, 张洪学, 宗兆云.

OVT数据域五维地震资料解释技术研究现状与进展

[J]. 石油物探, 2018, 57(2):155-178.

Yin

X Y

, Zhang

H X

, Zong

Z Y

Research status and progress of 5D seismic data interpretation in OVT domain

[J]. Geophysical Prospecting for Petroleum, 2018, 57(2):155-178.

DOI:10.1111/j.1365-246X.1978.tb03754.x URL [本文引用: 1]

[4]

Crampin

Seismic-wave propagation through a cracked solid:polarization as a possible dilatancy diagnostic

[J]. Geophysical Journal International, 1978, 53(3):467-496.

[5]

Thomsen

Elastic anisotropy due to aligned cracks in porous rock

[J]. Geophysical Prospecting, 1995, 43(6):805-829.

DOI:10.1111/j.1365-2478.1995.tb00282.x URL [本文引用: 1]

[6]

张琬璐

. 基于Rüger公式的地震各向异性参数反演方法研究[D]. 北京: 中国石油大学, 2016.

Zhang

W L

. Research on seismic anisotropic parameter inversion method based on Rüger’s formula[D]. Beijing: China University of Petroleum, 2016.

[7]

任涛, 秦军, 周阳, 等.

利用叠前方位各向异性预测准噶尔盆地西北缘车476井区火山岩裂缝

[J]. 石油地质与工程, 2021, 35(1):1-7.

Ren

, Qin

, Zhou

, et al.

Prediction of volcanic fractures in Che 476 well block in northwestern margin of Junggar basin by using pre-stack azimuthal anisotropy

[J]. Petroleum Geology and Engineering, 2021, 35(1):1-7.

[8]

Rüger

P-wave reflection coefficients for transversely isotropic models with vertical and horizontal axis of symmetry

[J]. Geophysics, 1997, 62(3):713-721.

DOI:10.1190/1.1444181 URL [本文引用: 2]

[9]

陈怀震, 印兴耀, 杜炳毅, 等.

裂缝型碳酸盐岩储层方位各向异性弹性阻抗反演

[J]. 地球物理学进展, 2013, 28(6):3073-3079.

Cheng

H Z

, Yin

X Y

, Du

B Y

, et al.

Current situations and recent progress in different azimuths seismic exploration

[J]. Progress in Geophysics, 2013, 28(6):3073-3079.

[10]

Bakulin

, Vladimir

Tsvankin Ilya.Estimation of fracture parameters from reflection seismic data—Part I:HTI model due to a single fracture set

[J]. Geophysics, 2000, 65(6):1788-1802.

DOI:10.1190/1.1444863 URL [本文引用: 3]

[11]

龚明平, 张军华, 王延光, 等.

分方位地震勘探研究现状及进展

[J]. 石油地球物理勘探, 2018, 53(3):642-658.

Gong

M P

, Zhang

J H

, Wang

Y G

, et al.

Current situations and recent progress in different azimuths seismic exploration

[J]. Oil Geophysical Prospecting, 2018, 53(3):642-658.

[12]

张保庆, 周辉, 左黄金, 等.

宽方位地震资料处理技术及应用效果

[J]. 石油地球物理勘探, 2011, 46(3):396-400.

Zhang

B Q

, Zhou

, Zuo

H J

, et al.

Wide azimuth data processing techniques and their applications

[J]. Oil Geophysical Prospecting, 2011, 46(3):396-400.