CN103869366A - 一种确定裂隙裂缝走向的方法及装置 - Google Patents

Abstract

本发明提供了一种确定裂隙裂缝走向的方法及装置，包括：获取裂隙裂缝走向的多波地震测量数据中包括两水平分量数据和垂直分量数据的三分量地震数据；对两水平分量数据进行极化分析确定极化角和极化角方差；确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；确定垂直分量数据旅行时差方位角和方位角方差；对垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差；对极化角、最大特征值对应的特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合确定裂隙裂缝方位角。降低地震裂隙裂缝测量结果的不确定性，提高测量结果的稳定性和精度。

Description

一种确定裂隙裂缝走向的方法及装置

技术领域

本发明涉及地震勘探技术领域，具体的讲是一种确定裂隙裂缝走向的方法及装置。

背景技术

裂隙裂缝是油气的储集空间和渗滤通道，也是水、煤层气富集、储存、运移的场所。探明地下裂隙裂缝走向、提高地震探测精度对于油气预测、提高油气产能以及预防瓦斯爆炸、顶/底板突水等都至关重要。

利用地震波来研究岩石的内部结构在理论上、实验室内和野外观测上都取得了显著的进步，特别是地震横波。定向裂隙裂缝对横波的影响在多年前就载入文献。介质模型由单一垂直的定向裂隙组成时，裂隙的方位各向异性具有两个主方向，分别平行和垂直裂隙走向。横波垂直入射到单一垂直定向裂隙介质时将发生波的偏振，分裂成两个正交偏振的分量快横波和慢横波，分别平行和垂直裂隙方向。相对于横波来说，纵波对介质的各向异性反应灵敏度较低，但由于横波勘探的难度大，上世纪90年代，人们又开始了纵波裂隙裂缝研究，结果表明，在裂隙裂缝引起的方位各向异性介质中，纵波振幅特征、正常时差速度及反射时差均随方位的变化而变化。

利用横波信息测量裂隙裂缝的方法有极化图法、旋转分析和协方差矩阵法等。

极化图法当横波沿z轴（垂直向下）的方向传播时，在x和y（水平面上的两个轴）两方向的质点位移一般为椭圆，椭圆形状与测线方位和裂缝走向的夹角及分裂横波间的延迟时间有关。通常我们通过检测极化图的初始线性起跳的方向确定快慢横波的极化方向和裂隙裂缝走向。极化图解释时主观性比较大，但干涉比较严重的横波分裂来说，判别快横波极化方向比较可靠。

旋转法旋转分析技术可用于二分量或四分量水平记录，最常用的是Alford旋转法。该方法从四分量的地震数据中获得只包含一个波的主时序列，即将在四分量记录剖面按矩阵形式排放，寻找一个旋转角，对四分量数据进行旋转，使非对角线上两个元素达最小，从而求得快横波和慢横波剖面。分析快横波和慢横波剖面，得到快、慢横波到时差或振幅变化，确定各向异性程度及裂隙裂缝的密度。旋转法要求双源激发，震源匹配，难度比较大。

协方差矩阵法用三分量地震记录形成协方差矩阵，协方差矩阵特征值的大小反映质点的平均运动性质，特征向量的方向决定快、慢横波的极化方向。协方差矩阵法即利用协方差矩阵最大特征值对应的特征向量确定裂隙裂缝走向。

胀缩源是目前普遍采用的多分量地震震源激发形式，难度比双源激发小的多，且数据质量好。极化图法和协方差矩阵法距均可用于胀缩源激发数据。

利用纵波信息测量裂隙裂缝的方法有方位振幅分析法、方位NMO速度椭圆法和方位旅行时差法等。

方位振幅分析法纵波在裂隙裂缝性介质中传播时，其振幅不但是偏移距的函数，也与裂隙裂缝方位与测线之间的夹角有关，当纵波在平行裂隙裂缝走向平面内传播时振幅最大，垂直裂隙裂缝走向平面内传播时振幅最小。利用不同方位上的振幅数据，确定最大振幅对应的方位即裂隙裂缝走向。

方位旅行时差在裂隙裂缝性介质中，纵波在不同方位测线上的旅行时的差或层间旅行时差是裂隙裂缝方向的函数。对于两对正交测线，其旅行时的差或层间旅行时差与两倍的方位角呈余弦关系。利用纵波在两对正交测线方向上的旅行时差的比值可以直接计算裂隙裂缝走向的方位角。

方位NMO速度椭圆法对任意非均匀各向异性介质，单一波型NMO速度的方位变化典型地描述为一个椭圆，对裂隙裂缝介质，水平反射层的纵波NMO速度受对称轴方位、垂直入射的纵波速度及两个各向异性参数控制，其NMO速度椭圆半轴与对称面一致，长半轴指示裂缝走向。

分裂的横波基本上沿同一射线路径旅行，它们的速度、反射振幅之间的差异主要与地层裂隙裂缝特征相关，而与横向非均匀性关系不大。正因为这一点，横波方法是现行地震法测量裂隙裂缝最常用方法，但横波技术通常需要在地面激发S波并用多分量观测，成本高，而且难以采集到适于作可靠的偏振及旅行时分析所需要的高质量的横波数据，于是，人们又开始研究P-S转换波。P-S转换波含有与S-S波相同的信息，但震源是纵波震源，其采集成本低，且劳动强度小，但转换横波的数据质量仍不及纵波数据，这些都直接限制了利用横波数据测量裂隙裂缝的准确度。

纵波裂隙裂缝检测方法应该说是一种有效可行的方法，但由于影响振幅的因素很多，除炮检距和方位之外，较为敏感的还有采集偏差，上覆层非均匀性等，再者，纵波方位NMO速度不但受裂隙裂缝方位的影响，也与地层倾角、倾斜界面的走向及介质的非均匀性有关，这些因素都直接影响纵波裂隙裂缝测量的精度。另一方面，纵波特征对裂隙裂缝的敏感性比较差，只利用纵波数据进行裂隙裂缝测量，结果也存在一定的不确定性。

目前地震裂隙裂缝测量要么只利用纵波数据（即常规地震勘探），要么只利用横波数据，而且基本上只采用一种方法。由于地下环境复杂、地震数据采集的不完备和测量误差等，只采用一种波型数据或一种方法测量地层裂隙裂缝存在两个缺陷：一是抗干扰能力弱，一是不能充分利用多分量数据提供的纵、横波数据裂隙裂缝信息的冗余性，测量结果稳定性差，精度低。

发明内容

为获得高精度、稳定性好的裂隙裂缝测量结果，降低地震裂隙裂缝测量结果的不确定性，本发明实施例提供了一种确定裂隙裂缝走向的方法，方法包括：

获取多波地震测量数据中的三分量地震数据，所述三分量地震数据包括两水平分量数据和垂直分量数据；

对所述的两水平分量数据进行极化分析确定极化角和极化角方差；

根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

根据所述垂直分量数据旅行时确定旅行时差方位角和旅行时差方位角方差；

对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和椭圆拟合方差；

根据所述的极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差对所述的极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合，获得定裂隙裂缝走向。

此外，本发明还公开了一种确定裂隙裂缝走向的装置，装置包括：

数据获取模块，用于获取多波地震测量数据中的三分量地震数据，所述三分量地震数据包括两水平分量数据和垂直分量数据；

极化角确定模块，用于对所述的两水平分量数据进行极化分析确定极化角和极化角方差；

特征向量方位角确定模块，用于根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

旅行时差方位角确定模块，用于根据垂直分量数据确定旅行时差方位角和旅行时差方位角方差；

椭圆长轴方位角确定模块，用于对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差；

裂隙裂缝方向确定模块，用于根据所述的极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差对所述的极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合，获得定裂隙裂缝走向。

本发明利用纵波数据和横波数据提供的裂隙裂缝信息具有冗余性、互补性和协作性，同时利用纵波数据和横波数据测量裂隙裂缝，并通过一种数据融合方法对多种方法获得的裂隙裂缝信息进行融合，降低地震裂隙裂缝测量结果的不确定性，提高测量结果的稳定性和精度。

为让本发明的上述和其他目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附图式，作详细说明如下。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。

图1为本发明提供的一种确定裂隙裂缝走向的方法的流程图；

图2为本发明实施例中确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差的流程图；

图3为本发明实施例中确定所述垂直分量数据的裂隙裂缝方位角和裂隙裂缝方位角方差的流程图；

图4为本发明提供的一种确定裂隙裂缝走向的装置的框图；

图5为本发明实施例中特征向量方位角确定模块的框图；

图6为本发明实施例中垂直分量数据裂隙裂缝方位角确定模块的框图；

图7为本发明实施例中椭圆长轴方位角确定模块的框图；

图8本实施例的裂隙裂缝地质模型图；

图9为方位测线分布图；

图10为两对正交测线分布图；

图11为15°方位测线方向上的三分量地震数据；

图12为极化分析和特征值分析结果图

图13为方位旅行时差交会分析图；

图14为叠加速度分析速度谱和解释结果图；

图15为NMO速度椭圆图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明提供了一种确定裂隙裂缝走向的方法，如图1所示，该方法包括：

步骤S101，获取多波地震测量数据中的三分量地震数据，其中，三分量地震数据包括两水平分量数据和垂直分量数据；

步骤S102，对两水平分量数据进行极化分析确定极化角和极化角方差；

本发明中，采用两分量极化图分析法确定极化角和极化角方差。根据步骤S101获得的两个水平分量数据进行极化分析，获得极化角β₁(i)，并计算方差δ₁；

极化分析采用常规两分量极化图分析方法；

方差δ₁按下式计算，

${\mathrm{\δ}}_1=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{({\mathrm{\β}}_1\left(i\right)-{\stackrel{\‾}{\mathrm{\β}}}_1)}^2},{\stackrel{\‾}{\mathrm{\β}}}_1=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\β}}_1\left(i\right)---\left(1\right)$

其中i为道号（i=1,2,……,N），N为道数。

步骤S103，根据两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

步骤S104，根据垂直分量数据确定旅行时差方位角和旅行时差方位角方差；

步骤S105，对垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差；具体为对每条方位测线方向上数据进行速度分析，获得不同方位测线上的NMO速度，用常规椭圆拟合法进行拟合，获得椭圆长轴的方位角β₄和拟合方差δ₄；

步骤S106，根据极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差对所述的极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合，获得定裂隙裂缝走向。对获得的数据按式（5）进行融合，获得裂隙裂缝走向。

$\mathrm{\α}=\frac{\frac{{\stackrel{\‾}{\mathrm{\β}}}_1}{{\mathrm{\δ}}_1^2}+\frac{{\stackrel{\‾}{\mathrm{\β}}}_2}{{\mathrm{\δ}}_{\text{2}}^2}\text{+}\frac{{\mathrm{\β}}_3}{{\mathrm{\δ}}_3^2}+\frac{{\mathrm{\β}}_4}{{\mathrm{\δ}}_4^2}+\frac{{\mathrm{\β}}_0}{{\mathrm{\δ}}_0}}{\frac{1}{{\mathrm{\δ}}_1^2}+\frac{1}{{\mathrm{\δ}}_2^2}+\frac{1}{{\mathrm{\δ}}_3^2}+\frac{1}{{\mathrm{\δ}}_4^2}+\frac{1}{{\mathrm{\δ}}_0}}$

β₀和δ₀按下式计算，

${\mathrm{\β}}_0=\frac{1}{4}({\stackrel{\‾}{\mathrm{\β}}}_1+{\stackrel{\‾}{\mathrm{\β}}}_2+{\mathrm{\β}}_3+{\mathrm{\β}}_4)$

${\mathrm{\δ}}_0=\sqrt{\frac{1}{4}[{({\stackrel{\‾}{\mathrm{\β}}}_1-{\mathrm{\β}}_0)}^2+{{({\stackrel{\‾}{\mathrm{\β}}}_2-{\mathrm{\β}}_0)}^2({\stackrel{\‾}{\mathrm{\β}}}_3-{\mathrm{\β}}_0)}^2+{({\mathrm{\β}}_4+{\mathrm{\β}}_0)}^2]}---\left(5b\right)$

上两式中，α为裂隙裂缝走向，

和δ₁为步骤S102获得的数据，和δ₂为步骤S103获得的数据，

和δ₃为步骤S104中确定的垂直分量数据的裂隙裂缝方位角和裂隙裂缝方位角方差，

和δ₄为步骤S105中的椭圆长轴的方位角和拟合方差。

其中，如图2所示，步骤S103中根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差包括：

步骤S1031，计算选取的目标层时窗内数据的方差和协方差；

步骤S1032，根据所述目标层时窗内数据的方差和协方差确定协方差矩阵；

步骤S1033，根据所述协方差矩阵确定最大的特征值对应的特征向量方位角和特征向量方位角方差。

具体确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差的步骤为：

选取目标层，确定时窗[t₁,t₂]，根据步骤S101获得的两个水平分量数据计算方差和协方差，构成协方差矩阵，

$\left[\begin{array}{cc}{C}_{\mathrm{xx}}& C_{\mathrm{xy}}\\ C_{\mathrm{yx}}& C_{\mathrm{yy}}\end{array}\right]---\left(2\right)$

该步骤中方差按下式计算，

${\stackrel{\‾}{u}}_x=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}{u}_x\left(t\right),{\stackrel{\‾}{u}}_y=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}{u}_y\left(t\right),---\left(2a\right)$

$C_{\mathrm{xx}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_x\left(t\right)-{\stackrel{\‾}{u}}_x{]}^2,C_{\mathrm{yy}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_y\left(t\right)-{\stackrel{\‾}{u}}_y{]}^2,---\left(2b\right)$

协方差按下式计算，

$C_{\mathrm{xy}}=C_{\mathrm{yx}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_x\left(t\right)-{\stackrel{\‾}{u}}_x][u_y\left(t\right)-{\stackrel{\‾}{u}}_y];---\left(2c\right)$

其中，u_x为x方向水平分量数据，u_y为y方向水平分量数据；

时窗[t₁,t₂]的选取根据目标层底部反射时间和反射子波长度确定，包含整个目标层底部反射波数据；

计算获得的协方差矩阵（2）的特征值和向量，获得最大特征值对应的特征向量的方位角β₂(i)(i为道号，i=1,2,……,N，N为道数)，并计算方差；

协方差矩阵（2）的特征值和向量采用常规矩阵求解方法；

方差按式（3）计算，

${\mathrm{\δ}}_2=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{({\mathrm{\β}}_2\left(i\right)-{\stackrel{\‾}{\mathrm{\β}}}_2)}^2},{\stackrel{\‾}{\mathrm{\β}}}_2=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\β}}_2\left(i\right)---\left(3\right)$

其中i为道号（i=1,2,……,N），N为道数。

如图3所示，步骤S104中根据垂直分量数据旅行时确定旅行时差方位角和方位角方差包括：

步骤S1041，从垂直分量数据中抽取两对正交测线上的数据，所述两对正交测线的方位相差为45度；

步骤S1042，计算所述两对正交测线的不同偏移距的旅行时差，并对所述不同偏移距的旅行时差进行交会分析，确定旅行时差方位角和方位角方差。

从步骤101获得的垂直分量数据中抽取两对正交测线上的数据，用方位旅行时差法获得裂隙裂缝方位角β₃，并计算偏差δ₃；

第i偏移距上方位旅行时差获得的裂隙裂缝方位角按下式计算，

${\mathrm{\θ}}_i=\frac{1}{2}{\mathrm{tg}}^{-1}(-\frac{{\mathrm{\Δt}}_2({\mathrm{\θ}}_i,x_i)}{\mathrm{\Δ}{t}_1({\mathrm{\θ}}_i,x_i)})---\left(4b\right)$

其中，Δt₁(θ_i,x_i)为第一对正交测线方向上在偏移距xi处的旅行时差，Δt₂(θ_i,x_i)为第二对正交测线方向上在偏移距x_i处的旅行时差，x_i为第i偏移距上的偏移距。

偏差δ₃按下式计算，

${\mathrm{\δ}}_3\sqrt{\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{({\mathrm{\θ}}_i-\stackrel{\‾}{\mathrm{\θ}})}^2}---\left(4\right)$

$\stackrel{\‾}{\mathrm{\θ}}=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{\mathrm{\θ}}_i---\left(4a\right)$

上式中，M为偏移距个数，θ_i为第i偏移距上方位旅行时差获得的裂隙裂缝方位角。

此外本发明还公开了一种确定裂隙裂缝走向的装置，如图4所示，装置包括：

数据获取模块401，用于获取多波地震测量数据中的三分量地震数据，所述三分量地震数据包括两水平分量数据和垂直分量数据；

极化角确定模块402，用于对所述的两水平分量数据进行极化分析确定极化角和极化角方差；本发明实施例中，采用两分量极化图分析法确定极化角和极化角方差。

特征向量方位角确定模块403，用于根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

旅行时差方位角确定模块404，用于根据旅行时差法确定所述垂直分量数据的裂隙裂缝方位角和裂隙裂缝方位角方差；

椭圆长轴方位角确定模块405，用于对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差；

裂隙裂纹方向确定模块406，用于对所述极化角、极化角方差、最大特征值对应的特征向量方位角、特征向量方位角方差、裂隙裂缝方位角、裂隙裂缝方位角方差、椭圆长轴的方位角、拟合方差进行数据融合确定裂隙裂缝方位角。

如图5所示，特征向量方位角确定模块403包括：

窗内数据方差计算单元4031，用于计算选取的目标层时窗内数据的方差和协方差；

协方差矩阵确定单元4032，用于根据所述目标层时窗内数据的方差和协方差确定协方差矩阵；

特征向量方位角确定单元4033，用于根据所述协方差矩阵确定最大的特征值对应的特征向量方位角和特征向量方位角方差。

如图6所示，旅行时差方位角确定模块404包括：

正交测试数据抽取单元4041，用于从垂直分量数据中抽取两对正交测线上的数据，所述两对正交测线的方位相差为45度；

旅行时差确定单元4042，用于计算两对正交测线上不同偏移距的旅行时差；

旅行时差方位角确定单元4043，用于对两对正交测线上不同偏移距的旅行时差进行交会分析，确定旅行时差方位角和方位角方差。

如图7所示，椭圆长轴方位角确定模块405包括：

方位测线数据抽取单元4051，用于从垂直分量数据抽取不同方位测线上的数据；

NMO速度获取单元，用于对所述的不同方位测线上数据进行速度分析，获得不同方位的NMO速度；

椭圆拟合单元，用于对所述的不同方位NMO速度进行椭圆拟合，获得椭圆长轴的方位角和拟合方差。

本发明提供一种提高多波地震裂隙裂缝测量稳定性方法，通过多分量地震测量提供的纵波数据和横波数据，利用不同方法测量裂隙裂缝走向，并通过一种数据融合技术对多方法结果进行融合，获得高精度、稳定性好的裂隙裂缝测量结果。

下面结合具体的实施例对本发明做进一步详细说明。本实施例是一个预测裂隙裂缝走向的实例。该实例有三层地层，中间层发育有裂缝，裂缝走向与x1轴平行，即方位角为0°。如图8所示，为本实施例的裂隙裂缝地质模型图。图9为方位测线分布图，测线间隔15°，图10为两对正交测线分布图。

具体实现步骤为：

1.获取用于裂隙裂缝预测的三分量数据，两个水平分量和一个垂直分量，水平分量数据为共转换点道集数据，垂直分量为共中心点道集数据；图11为15°方位测线方向上的三分量地震数据，x和y为两水平分量，z为垂直分量。

2.根据两个水平分量数据，对每道进行极化分析，获得每一道的极化角β₁(i)；

3.按下式对多道获得的极化角进行平均，并计算方差，

${\mathrm{\δ}}_1=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{({\mathrm{\β}}_1\left(i\right)-{\stackrel{\‾}{\mathrm{\β}}}_1)}^2},{\stackrel{\‾}{\mathrm{\β}}}_1=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\β}}_1\left(i\right)---\left(1\right)$

其中i=1,2,……,N，N=100；

4.选取目标层时窗，根据两个水平分量数据：u_x(t)和u_y(t)，按下式计算时窗内数据方差和协方差，

${\stackrel{\‾}{u}}_x=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}{u}_x\left(t\right),{\stackrel{\‾}{u}}_y=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}{u}_y\left(t\right),---\left(2a\right)$

$C_{\mathrm{xx}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_x\left(t\right)-{\stackrel{\‾}{u}}_x{]}^2,C_{\mathrm{yy}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_y\left(t\right)-{\stackrel{\‾}{u}}_y{]}^2,---\left(2b\right)$

$C_{\mathrm{xy}}=C_{\mathrm{yx}}=\frac{1}{t_2-t_1}\underset{t=t_1}{\overset{t_2}{\mathrm{\Σ}}}[u_x\left(t\right)-{\stackrel{\‾}{u}}_x][u_y\left(t\right)-{\stackrel{\‾}{u}}_y];---\left(2c\right)$

确定协方差矩阵：

$\left[\begin{array}{cc}{C}_{\mathrm{xx}}& C_{\mathrm{xy}}\\ C_{\mathrm{yx}}& C_{\mathrm{yy}}\end{array}\right]---\left(2\right)$

5.求解协方差矩阵（2）的特征值和特征向量，取最大特征值对应的特征向量的方位角β₂(i)；

6.按下式对多道获得的特征向量方位角进行平均，并计算方差，

${\mathrm{\δ}}_2=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{({\mathrm{\β}}_2\left(i\right)-{\stackrel{\‾}{\mathrm{\β}}}_2)}^2},{\stackrel{\‾}{\mathrm{\β}}}_2=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\β}}_2\left(i\right)---\left(3\right)$

其中i=1,2,……,N，N=100；

7.从垂直分量数据中抽取两对正交测线上的数据，分别计算每对方位测线上的方位旅行时差，进行交会分析，获得裂隙裂缝方位角β₃；

8.并按下式计算方位旅行时差法获得裂隙裂缝方位角的偏差δ₃：

${\mathrm{\δ}}_3\sqrt{\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{({\mathrm{\θ}}_i-\stackrel{\‾}{\mathrm{\θ}})}^2}---\left(4\right)$

其中M=100，按下式计算，

$\stackrel{\‾}{\mathrm{\θ}}=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{\mathrm{\θ}}_i---\left(4a\right)$

θ_i按下式计算，

${\mathrm{\θ}}_i=\frac{1}{2}{\mathrm{tg}}^{-1}(-\frac{{\mathrm{\Δt}}_2({\mathrm{\θ}}_i,x_i)}{\mathrm{\Δ}{t}_1({\mathrm{\θ}}_i,x_i)})---\left(4b\right)$

其中，Δt₁(θ_i,x_i)为第一对正交测线方向上在偏移距xi处的方位旅行时差，Δt₂(θ_i,x_i)为第二对正交测线方向上在偏移距x_i处的方位旅行时差；

9.根据垂直分量数据，对每条方位测线方向上数据进行叠加速度分析，获得不同方位测线上的NMO速度，用常规椭圆拟合法进行拟合，获得椭圆长轴的方位角β₄和拟合方差δ₄；

10.按下式计算β₀和δ₀，

${\mathrm{\β}}_0=\frac{1}{4}({\stackrel{\‾}{\mathrm{\β}}}_1+{\stackrel{\‾}{\mathrm{\β}}}_2+{\mathrm{\β}}_3+{\mathrm{\β}}_4),$

${\mathrm{\δ}}_0=\sqrt{\frac{1}{4}[{({\stackrel{\‾}{\mathrm{\β}}}_1-{\mathrm{\β}}_0)}^2+{{({\stackrel{\‾}{\mathrm{\β}}}_2-{\mathrm{\β}}_0)}^2({\stackrel{\‾}{\mathrm{\β}}}_3-{\mathrm{\β}}_0)}^2+{({\mathrm{\β}}_4+{\mathrm{\β}}_0)}^2]}---\left(5b\right)$

11.按下式计算裂隙裂缝方位角α，

$\mathrm{\α}=\frac{\frac{{\stackrel{\‾}{\mathrm{\β}}}_1}{{\mathrm{\δ}}_1^2}+\frac{{\stackrel{\‾}{\mathrm{\β}}}_2}{{\mathrm{\δ}}_{\text{2}}^2}\text{+}\frac{{\mathrm{\β}}_3}{{\mathrm{\δ}}_3^2}+\frac{{\mathrm{\β}}_4}{{\mathrm{\δ}}_4^2}+\frac{{\mathrm{\β}}_0}{{\mathrm{\δ}}_0}}{\frac{1}{{\mathrm{\δ}}_1^2}+\frac{1}{{\mathrm{\δ}}_2^2}+\frac{1}{{\mathrm{\δ}}_3^2}+\frac{1}{{\mathrm{\δ}}_4^2}+\frac{1}{{\mathrm{\δ}}_0}}---\left(5\right)$

图12为极化分析和特征值分析结果图，图13为方位旅行时差交会分析图，图14为叠加速度分析速度谱和解释结果图，图15为NMO速度椭圆图。四种方法所得的裂隙裂缝层方位角列于表1。用本发明方法进行数据融合后的角度值为-0.48度，即裂隙裂缝的优势走向为-0.48度，与每种方法预测结果比较，稳定性好，精度明显提高。

表1.方位角度表

本发明中应用了具体实施例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。

Claims (12)

1.一种确定裂隙裂缝走向的方法，其特征在于，所述的方法包括：

获取确定裂隙裂缝走向的多波地震测量数据中的三分量地震数据，所述三分量地震数据包括两水平分量数据和垂直分量数据；

对所述的两水平分量数据进行极化分析确定极化角和极化角方差；

根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

根据所述垂直分量数据旅行时确定旅行时差方位角和旅行时差方位角方差；

对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和椭圆拟合方差；

对所述极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合，获得裂隙裂缝走向。

2.如权利要求1所述的确定裂隙裂缝走向的方法，其特征在于，所述的对所述的两水平分量数据进行极化分析确定极化角和极化角方差包括：

采用两分量极化图分析法确定极化角和极化角方差。

3.如权利要求1所述的确定裂隙裂缝走向的方法，其特征在于，所述的根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差包括：

计算选取的目标层时窗内数据的方差和协方差；

根据所述目标层时窗内数据的方差和协方差确定协方差矩阵；

根据所述协方差矩阵确定最大的特征值对应的特征向量方位角和特征向量方位角方差。

4.如权利要求1所述的确定裂隙裂缝走向的方法，其特征在于，所述的根据垂直分量数据旅行时确定旅行时差方位角和方位角方差包括：

从所述垂直分量数据中抽取两对正交测线上的数据，所述两对正交测线的方位相差为45度；

计算所述两对正交测线的不同偏移距的旅行时差；

对所述的两对正交测线上不同偏移距的旅行时差进行交会分析，确定旅行时差方位角和方位角方差。

5.如权利要求1所述的确定裂隙裂缝走向的方法，其特征在于，所述的对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差包括：

从垂直分量数据抽取不同方位测线上的数据；

对所述的不同方位测线上的数据进行速度分析，获得不同方位的NMO速度；

对所述的不同方位的NMO速度进行椭圆拟合，获得椭圆长轴的方位角和拟合方差。

6.如权利要求1所述的确定裂隙裂缝走向的方法，其特征在于，所述的对极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合包括：

根据所述的极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差对所述的极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角和下式（1）进行数据融合：

$\mathrm{\α}=\frac{\frac{{\stackrel{\‾}{\mathrm{\β}}}_1}{{\mathrm{\δ}}_1^2}+\frac{{\stackrel{\‾}{\mathrm{\β}}}_2}{{\mathrm{\δ}}_{\text{2}}^2}\text{+}\frac{{\mathrm{\β}}_3}{{\mathrm{\δ}}_3^2}+\frac{{\mathrm{\β}}_4}{{\mathrm{\δ}}_4^2}+\frac{{\mathrm{\β}}_0}{{\mathrm{\δ}}_0}}{\frac{1}{{\mathrm{\δ}}_1^2}+\frac{1}{{\mathrm{\δ}}_2^2}+\frac{1}{{\mathrm{\δ}}_3^2}+\frac{1}{{\mathrm{\δ}}_4^2}+\frac{1}{{\mathrm{\δ}}_0}}---\left(1\right)$

${\mathrm{\β}}_0=\frac{1}{4}({\stackrel{\‾}{\mathrm{\β}}}_1+{\stackrel{\‾}{\mathrm{\β}}}_2+{\mathrm{\β}}_3+{\mathrm{\β}}_4)$

${\mathrm{\δ}}_0=\sqrt{\frac{1}{4}[{({\stackrel{\‾}{\mathrm{\β}}}_1-{\mathrm{\β}}_0)}^2+{{({\stackrel{\‾}{\mathrm{\β}}}_2-{\mathrm{\β}}_0)}^2({\stackrel{\‾}{\mathrm{\β}}}_3-{\mathrm{\β}}_0)}^2+{({\mathrm{\β}}_4+{\mathrm{\β}}_0)}^2]}$

其中，

分别为所述极化角、特征向量方位角的平均值，β₃，β₄为所述旅行时差方位角、椭圆长轴的方位角；

δ₁，δ₂，δ₃，δ₄分别为所述极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差。

7.一种确定裂隙裂缝走向的装置，其特征在于，所述的装置包括：

数据获取模块，用于获取确定裂隙裂缝走向的多波地震测量数据中的三分量地震数据，所述三分量地震数据包括两水平分量数据和垂直分量数据；

极化角确定模块，用于对所述的两水平分量数据采用两分量极化图分析法进行极化分析确定极化角和极化角方差；

特征向量方位角确定模块，用于根据所述两水平分量数据确定协方差矩阵最大特征值对应的特征向量方位角和特征向量方位角方差；

旅行时差方位角确定模块，用于根据垂直分量数据确定旅行时差方位角和旅行时差方位角方差；

椭圆长轴方位角确定模块，用于对所述的垂直分量数据进行叠加速度分析和椭圆拟合确定椭圆长轴的方位角和拟合方差；

裂隙裂缝方向确定模块，用于对所述极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角进行数据融合，获得裂隙裂缝走向。

8.如权利要求7所述的确定裂隙裂缝走向的装置，其特征在于，所述的极化角确定模块采用两分量极化图分析法确定极化角和极化角方差。

9.如权利要求7所述的确定裂隙裂缝走向的装置，其特征在于，所述的特征向量方位角确定模块包括：

窗内数据方差计算单元，用于计算选取的目标层时窗内数据的方差和协方差；

协方差矩阵确定单元，用于根据所述目标层时窗内数据的方差和协方差确定协方差矩阵；

特征向量方位角确定单元，用于根据所述协方差矩阵确定最大的特征值对应的特征向量方位角和特征向量方位角方差。

10.如权利要求7所述的确定裂隙裂缝走向的装置，其特征在于，所述的旅行时差方位角确定模块包括：

正交测线数据抽取单元，用于从所述垂直分量数据中抽取两对正交测线上的数据，所述两对正交测线的方位相差为45度；

旅行时差确定单元，用于计算所述的两对正交测线上不同偏移距的旅行时差；

旅行时差方位角确定单元，用于对所述的两对正交测线上不同偏移距的旅行时差进行交会分析，确定旅行时差方位角和方位角方差。

11.如权利要求7所述的确定裂隙裂缝走向的装置，其特征在于，所述的椭圆长轴方位角确定模块包括：

方位测线数据抽取单元，用于从垂直分量数据抽取不同方位测线上的数据；

NMO速度获取单元，用于对所述的不同方位测线上数据进行速度分析，获得不同方位的NMO速度；

椭圆拟合单元，用于对所述的不同方位NMO速度进行椭圆拟合，获得椭圆长轴的方位角和拟合方差。

12.如权利要求7所述的确定裂隙裂缝走向的装置，其特征在于，所述的裂隙裂缝方向确定模块根据所述的极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差对所述的极化角、特征向量方位角、旅行时差方位角、椭圆长轴的方位角和下式（1）进行数据融合：

$\mathrm{\α}=\frac{\frac{{\stackrel{\‾}{\mathrm{\β}}}_1}{{\mathrm{\δ}}_1^2}+\frac{{\stackrel{\‾}{\mathrm{\β}}}_2}{{\mathrm{\δ}}_{\text{2}}^2}\text{+}\frac{{\mathrm{\β}}_3}{{\mathrm{\δ}}_3^2}+\frac{{\mathrm{\β}}_4}{{\mathrm{\δ}}_4^2}+\frac{{\mathrm{\β}}_0}{{\mathrm{\δ}}_0}}{\frac{1}{{\mathrm{\δ}}_1^2}+\frac{1}{{\mathrm{\δ}}_2^2}+\frac{1}{{\mathrm{\δ}}_3^2}+\frac{1}{{\mathrm{\δ}}_4^2}+\frac{1}{{\mathrm{\δ}}_0}}---\left(1\right)$

${\mathrm{\β}}_0=\frac{1}{4}({\stackrel{\‾}{\mathrm{\β}}}_1+{\stackrel{\‾}{\mathrm{\β}}}_2+{\mathrm{\β}}_3+{\mathrm{\β}}_4)$

${\mathrm{\δ}}_0=\sqrt{\frac{1}{4}[{({\stackrel{\‾}{\mathrm{\β}}}_1-{\mathrm{\β}}_0)}^2+{{({\stackrel{\‾}{\mathrm{\β}}}_2-{\mathrm{\β}}_0)}^2({\stackrel{\‾}{\mathrm{\β}}}_3-{\mathrm{\β}}_0)}^2+{({\mathrm{\β}}_4+{\mathrm{\β}}_0)}^2]}$

其中，

分别为所述极化角、特征向量方位角的平均值，β₃，β₄为所述旅行时差方位角、椭圆长轴的方位角；

δ₁，δ₂，δ₃，δ₄分别为所述极化角方差、特征向量方位角方差、旅行时差方位角方差、椭圆拟合方差。

Images