一种裂缝性储层的多波联合AVO反演方法、装置及电子设备

一种裂缝性储层的多波联合avo反演方法、装置及电子设备
技术领域
1.本文件涉及地震勘探技术领域，尤其涉及一种裂缝性储层的多波联合avo 反演方法、装置及电子设备。


背景技术：

2.由于褶皱挤压、深部地层抬升等构造运动，地层可以发育出竖直裂缝。这 种裂缝不但可以作为储层的储集空间，也可以改善储层的渗透率，从而增加储 层的经济可采性。在地球物理中，将发育一组直立裂缝的介质近似等效为具有 水平对称轴的横向各向同性(hti)介质，成为致密砂岩、碳酸盐岩、页岩、水 合物等各类储层裂缝预测所使用的基本模型。
3.随着多分量地震勘探的发展，多种反射波的联合反演方法发展迅速，大大 提高了储层预测的精度。在hti介质中，横波可分为快横波(s1波)和慢横波 (s2波)。由于s1波的偏振方向平行于各向同性平面，因此转换快横波(ps1 波)的反射系数仅受各向同性背景参数的影响。然而，转换慢横波(ps2波) 偏振垂直于裂缝面，包含了各向同性和各向异性参数的信息。因此，在振幅随 入射角变化(ava)响应中，分裂横波可以提供更多的反演参数信息。此外， 传播速度差异导致的s1波和s2波之间的时延也表征了裂缝发育的程度。
4.现阶段，由于多波数据在匹配处理上相对困难，尚没有针对pp波、ps1波 和ps2波的联合反演应用。


技术实现要素：

5.本发明实施例目的是提供一种基于裂缝性储层的多波联合avo反演方法、 装置及电子设备，能够联合pp波、ps1波和ps2波的裂缝性vsp储层数据进 行avo反演。
6.为了实现上述目的，本说明书实施例是这样实现的：
7.第一方面，提供一种裂缝性储层的多波联合avo反演方法，包括：
8.将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压缩至零 偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记录；
9.构建初始模型和所述初始模型对应的反演目标函数，其中，所述初始模型 用于根据ava道集数据中包含各向异性参数的测井数据计算裂缝参数，所述 反演目标函数为贝叶斯理论函数；
10.构建所述初始模型对应的反演方程，所述反演方程用于基于所述初始模型 按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导以及 所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初始模型 进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向异性参 数；
11.将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性储层的 各向异性参数；
12.基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝状态。
13.第二方面，提供一种裂缝性储层的多波联合avo反演装置，包括：
14.多波匹配模块，用于将裂缝性储层对应ps波、ps1波和ps2波的vsp地 震观测记录压缩至零偏vsp纵波初至的时间域中进行ava道集转换，得到 ava地震观测记录；
15.模型函数构建模块，用于构建初始模型和所述初始模型对应的反演目标函 数，其中，所述初始模型用于根据ava道集数据中包含各向异性参数的测井 数据计算裂缝参数，所述反演目标函数为贝叶斯理论函数；
16.反演方程构建模型，用于构建所述初始模型对应的反演方程，所述反演方 程用于基于所述初始模型按照精确反射系数正演模拟ava地震合成记录，并 根据正演算子的偏导以及所述反演目标函数的二阶近似解，计算所述初始模型 的梯度以对所述初始模型进行迭代，以根据所述初始模型的迭代结果反演所述 裂缝性储层的各向异性参数；
17.反演执行模块，用于将所述初始模型代入至所述反演方程中进行反演，得 到所述裂缝性储层的各向异性参数；
18.裂缝状态确定模块，用于基于所述裂缝性储层的各向异性参数，确定裂缝 性储层的裂缝状态。
19.第三方面，提供一种电子设备包括：存储器、处理器及存储在所述存储器 上并可在所述处理器上运行的计算机程序，所述计算机程序被所述处理器执行：
20.将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压缩至零 偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记录；
21.构建初始模型和所述初始模型对应的反演目标函数，其中，所述初始模型 用于根据ava道集数据中包含各向异性参数的测井数据计算裂缝参数，所述 反演目标函数为贝叶斯理论函数；
22.构建所述初始模型对应的反演方程，所述反演方程用于基于所述初始模型 按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导以及 所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初始模型 进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向异性参 数；
23.将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性储层的 各向异性参数；
24.基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝状态。
25.第四方面，提供一种算机可读存储介质，所述计算机可读存储介质上存储 有计算机程序，所述计算机程序被处理器执行时实现如下步骤：
26.将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压缩至零 偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记录；
27.构建初始模型和所述初始模型对应的反演目标函数，其中，所述初始模型 用于根据ava道集数据中包含各向异性参数的测井数据计算裂缝参数，所述 反演目标函数为贝叶斯理论函数；
28.构建所述初始模型对应的反演方程，所述反演方程用于基于所述初始模型 按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导以及 所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初始模型 进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向异性参 数；
29.将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性储层的 各向异性参数；
30.基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝状态。
31.本发明能够将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记 录压缩至零偏vsp纵波初至的时间域中进行ava道集转换，从而在不丢失裂 缝信息的情况下，得到ps波、ps1波和ps2波时间匹配的ava地震观测记 录；之后，构建推导裂缝参数的初始模型，并基于初始模型对应的反演方程， 以反演得到裂缝性储层的各向异性参数。其中，反演方程能够基于初始模型按 照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导和反演 目标函数的二阶近似解，计算初始模型的梯度以对初始模型迭代，这种梯度计 算方式能够有效解决精确反射系数与目标函数之间不具有线性关系的问题，可 提高基于精确反射系数的反演性能，最终得到更加准确的各向异性参数以用于 成像裂缝性储层的结构。
附图说明
32.为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施 例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述 中的附图仅仅是本发明实施例中记载的一些实施例，对于本领域普通技术人员 来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
33.图1为本发明实施例提供的裂缝性储层各向异性参数的多波联合反演方法 的第一种流程示意图。
34.图2为本发明实施例提供的裂缝性储层各向异性参数的多波联合反演方法 的第二种流程示意图。
35.图3为零偏vsp数据的p波初至(a)，以及由零偏vsp数据得到的p波 和s波层速度(b)的示意图。
36.图4为预处理后的z、r、t分量的炮记录。
37.图5为分离后的z、r、t分量的炮记录。
38.图6为匹配至t
pfb
域的z、r、t分量的炮记录。
39.图7为匹配至t0
pfb
域的z、r、t分量的炮记录。
40.图8为alford旋转计算的裂缝方位和快慢横波时延τ的地震记录示意 图。
41.图9为匹配至t0
pfb
域的pp波、ps1波和ps2波的地震记录示意图。
42.图10为pp波和ps波的入射角道集(a、b)、反射点x坐标道集(c、d)、 反射点z坐标道集(e)匹配至t0pfb域的地震记录示意图。
43.图11为反演参数λ,μ,和ρ的初始模型的示意图。
44.图12为井旁的裂缝方位和快满波时差τ的示意图。
45.图13为井旁的裂缝密度e的示意图。
46.图14为实际数据(a、c、e)与初始模型(b、d、f)生成的pp、ps1和ps2 波深度域ava道集的示意图。
47.图15为pp、ps1和ps2波深度域叠加记录的示意图。
48.图16为角度独立的深度域pp波、ps1波、ps2波的子波示意图。
49.图17为联合反演的结果示意图。
50.图18为本发明实施例提供的裂缝性储层各向异性参数的多波联合反演装 置的结构示意图。
51.图19为本发明实施例提供的电子设备的结构示意图。
具体实施方式
52.为了使本技术领域的人员更好地理解本说明书中的技术方案，下面将结合 本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述， 显然，所描述的实施例仅仅是本说明书一部分实施例，而不是全部的实施例。 基于本说明书中的实施例，本领域普通技术人员在没有作出创造性劳动前提下 所获得的所有其他实施例，都应当属于本说明书保护的范围。
53.如前文所述，多分量地震勘探能够大大提高了储层预测的精度，是当前地 震勘探领域重要的发展方向。由于多波分量在数据匹配上还存在一定难度，目 前尚没有实用的联合pp波、ps1波和ps2波的反演应用。
54.为此，本发明针对裂缝性储层，提出一种联合pp波、ps1波和ps2波进 行avo反演的技术方案，以用于对裂缝性储层提供更加精确的裂缝参数，比 如岩性和流体的裂缝密度。
55.图1是本发明实施例提供的一种裂缝性储层的多波联合avo反演方法的 流程图，包括如下步骤：
56.s102，将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压 缩至零偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记 录。
57.具体地，本步骤基于纵波初至与检波器深度的关系，将pp波、ps1波和 ps2波的vsp地震观测记录匹配压缩至零偏vsp纵波初至的时间域中，并通 过alford旋转，得到快慢横波时延和裂缝方位；之后，将压缩至零偏vsp纵 波初至的时间域中的vsp地震观测记录进行r分量和t分量的层剥离，得到 ps1波和ps2波各自的vsp地震观测记录；使用快慢横波时延τ校正vsp地 震观测记录中ps2波相对ps1波的时差；基于vsp地震观测记录中p波和s 波的层速度对pp波、ps1波和ps2波进行射线追踪(输入为由零偏vsp地震 观测记录得到的各向同性层速度)，得到每个炮的pp波、ps1波和ps2波的入 射角和反射点的坐标道集；将每个炮的pp波、ps1波和ps2波的入射角和反 射点的坐标道集进行ava道集转换，得到ava地震观测记录。
58.s104，构建初始模型和所述初始模型对应的反演目标函数，其中，所述初 始模型用于根据ava道集数据中包含各向异性参数的测井数据计算裂缝参数， 所述反演目标函数为贝叶斯理论函数。
59.应理解，初始模型还包含有各向同性参数，也就是基于各向同性参数和各 向异性参数的测井数据计算裂缝参数，这里本文不作具限定。
60.可选地，所述初始模型用于根据ava道集数据中的层速度和快慢横波时 延，计算裂缝密度。
61.所述初始模型针对所述裂缝密度的推导公式为：
62.其中，e表示裂缝密度，l表示alford旋转时扫描的时窗利用初至速度转 换成的深
度窗,v
s1
表示s1波的横波速度，τ表示快慢横波时延,g＝μ/(λ 2μ)，μ 和λ为lam
é
反演参数。
63.对应地，所述目标函数的公式为：j(m)＝[a
prs-a
syn
(m
ite
,i)]
t
[a
prs
‑ꢀasyn
(m
ite
,i)]
t
k(m
ite
)；其中，j(m)表示目标函数，i和k分别表示入射角 和采样点；m表示采样点k对应的反演参数向量m＝[λk,μk,ρk,ek]
t
，ρ表示岩石 密度，a
prs
表示反演中输入的ava记录，m
ite
表示迭代更新的反演参数， k(m
ite
)表示正则化柯西分布项，a
syn
表示用于更新初始模型的ava记录， a
syn
(m
ite
,i)＝w*r(w
ite
,i)，w
pp
、w
ps1
、 w
ps2
分别表示pp波、ps1和ps2波的子波矩阵；
[0064]
所述初始模型按照精确反射系数正演模拟ava地震合成记录的公式为 r(m
ite
,i)＝[r
pp
(m
ite
,i),r
ps1
(m
ite
,i),r
ps2
(m
ite
,i)]
t
，r
pp(mite,i)
、 r
ps1
(m
ite
,i)和r
ps2
(m
ite
,i)分别表示pp波、ps1和ps2波对应的使用m
ite
计算的在入射角i时的反射系数矩阵。
[0065]
s106，构建所述初始模型对应的反演方程，所述反演方程用于基于所述初 始模型按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏 导以及所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初 始模型进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向 异性参数。
[0066]
可选地，所述反演方程中对初始模型进行迭代公式为：m
ite
＝m0 δm；其 中，m0表示初始模型，δm表示迭代增量，δm＝-(h λi)-1gt
κ，h表示 hessian矩阵，i表示对角矩阵，κ表示目标函数的jacobian矩阵示目标函数的jacobian矩阵g表示目标函数的jacobian矩阵的二阶形式,的二阶形式,分别表示入射角i下的pp波、ps1波和ps2波的jacobian矩阵。
[0067]
可选地，所述反演方程在所述目标函数j(m)的l2-范数小于预设阈值时， 结束对所述初始模型的迭代。
[0068]
s108，将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性 储层的各向异性参数。
[0069]
应理解，本步骤反演参数还包括各向同性参数，这里本文不对各向同性参 数反演进行举例赘述。
[0070]
s110，基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝发育 状态。
[0071]
其中，裂缝发育状态可以但不限于包括裂缝密度等信息，这里本文不作具 体限定。
[0072]
本发明能够将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记 录压缩至零偏vsp纵波初至的时间域中进行ava道集转换，从而在不丢失裂 缝信息的情况下，得到ps波、ps1波和ps2波时间匹配的ava地震观测记 录；之后，构建推导裂缝参数的初始模型，并基于初始模型对应的反演方程， 以反演得到裂缝性储层的各向异性参数。其中，反演方程能够基于初始模型按 照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导
和反演 目标函数的二阶近似解，计算初始模型的梯度以对初始模型迭代，这种梯度计 算方式能够有效解决精确反射系数与目标函数之间不具有线性关系的问题，可 提高基于精确反射系数的反演性能，最终得到更加准确的各向异性参数以用于 成像裂缝性储层的结构。
[0073]
下面对本发明实施例的方法进行详细介绍。
[0074]
本发明实施例根据bakulin的建模理论，将反演参数简化为了lam
é
参数 λ、μ、岩石密度ρ和裂缝密度e。采用精确反射系数适用于强对比和强各向异 性地层。为了解决精确反射系数带来的非线性问题，提出了一种用于贝叶斯目 标函数的二阶levenberg
–
marquardt(lm)算法，来进行hti介质深度域联合 pp波、ps1波和ps2波的ava反演。
[0075]
为了保护波场中包含的裂缝信息，采用保护数据矢量特征的处理方法对现 场vsp地震观测记录进行预处理。之后利用纵波初至与检波器深度的关系将 三个分量的数据匹配并压缩至零偏vsp纵波初至时间域，并通过alford旋转， 进行r分量和t分量的层剥离，得到ps1波和ps2波的vsp地震观测记录， 并计算出了裂缝方位角和快慢横波时差τ。在使用τ消除ps2波的vsp地震 观测记录的时差后，得到了匹配至零偏vsp纵波初至时间域的pp波、ps1波 和ps2波。
[0076]
在消除了各向异性影响后，通过各向同性的射线追踪，求取了入射角、反 射点坐标的道集，来进行深度域ava道集的转换。并采用零偏vsp地震观测 记录得到的层速度和alford旋转得到的τ提出了一个裂缝密度初始模型计算公 式，实现反演。
[0077]
其中，贝叶斯约束下的非线性最小二乘的反演方程的原理如下：
[0078]
确反射系数的计算需要有刚度系数c
ij
。通过bakulin的建模方式使用反演 参数λ、μ、ρ和e来计算c
ij
。这里，定义裂缝法向和切向弱度分别为δn和δ
t
， 则在干裂缝或充满流体裂缝中有：
[0079][0080][0081]
其中g＝μ/(λ 2μ)。在充满流体的裂缝中有：
[0082]
δn＝0
ꢀꢀꢀ
(3)
[0083]
则根据线性滑移理论，裂缝方位为的hti介质刚度系数矩阵可以写为：
[0084]c11
＝(λ 2μ)(1-δn),
[0085]c13
＝λ(1-δn),
[0086][0087][0088]c44
＝μ,
[0089]c55
＝μ(1-δ
t
),
[0090]c66
＝c
55
.
ꢀꢀꢀ
(4)
[0091]
建立基于柯西分布的后验分布的最大概率得到的目标函数:j(m)＝ [a
prs-a
syn
(m
ite
,i)]
t
[a
prs-a
syn
(m
ite
,i)]
t
k(m
ite
)
ꢀꢀꢀ
(5)
[0092]
其中，j(m)表示目标函数，i和k分别表示入射角和采样点； m表示采样点k对应的反演参数向量m＝[λk,μk,ρk,ek]
t
，ρ表示岩石密度，a
syn
表示反演中输入的ava记录，m
ite
表示迭代更新的反演参数，k(m
ite
)表示正 则化柯西分布项，a
syn
表示用于更新初始模型的ava记录。
[0093]asyn
可以表示为：
[0094]asyn
(m
ite
,i)＝w*r(m
ite
,i)
ꢀꢀꢀ
(6)
[0095][0096]
其中，w
pp
、w
ps1
、w
ps2
分别表示pp波、ps1和ps2波的子波矩阵；
[0097]
公式(6)定义所述初始模型按照精确反射系数正演模拟ava地震合成记 录的公式为：
[0098]
r(m
ite
,i)＝[r
pp
(m
ite
,i),r
ps1
(m
ite
,i),r
ps2
(m
ite
,i)]
t
ꢀꢀꢀ
(8)
[0099]
其中，r
pp
(m
ite
,i)、r
ps1
(m
ite
,i)和r
ps2
(m
ite
,i)分别表示pp波、ps1和 ps2波对应的使用m
ite
计算的在入射角i时的反射系数矩阵。
[0100]
对于hti介质，反射系数一般采用的是线性近似方程。然而，为了使反演 方程适用于更广泛的条件，包括强对比度和强各向异性的情况，公式(8)中采 用了精确反射系数m
ite
由每次迭代中的迭代增量δm和初始模型m0计算得到：
[0101]mite
＝m0 δm
ꢀꢀꢀꢀ
(9)
[0102]
其中δm为在采样点k处的所有反演参数的校正量。
[0103]
δm＝[δλk,δμk,δρk,δ ek]
t
ꢀꢀꢀ
(10)
[0104]
由于反演参数向量m和反射系数矩阵r之间存在非线性关系，这里可以 使用lm算法求解公式(10)。通常，通过求取目标函数的一阶泰勒展开式来 求解迭代参数δm。在弱对比度和弱各向异性条件下，这种近似下的反演结果 是可靠的。然而，在强各向异性或强对比度的情况下，各向异性参数和反射系 数之间的线性关系较差。
[0105]
因此，提出了目标函数的二阶近似：
[0106][0107]
关于δm的偏导为：
[0108][0109]
对于公式(12)中的偏导项，这里定义目标函数的jacobian矩阵为：
[0110][0111]
以及二阶形式：
[0112][0113]
在精确反射系数的计算中，由于矩阵嵌套，很难得到显式导数表达式。因 此，使用中心差分法获得目标函数关于反演参数的一阶和二阶导数的二阶近似 解，
[0114][0115][0116]
其中p表示四个反演参数，表示差分步长。对于通常会根据经验给 出一个小值。但是，如果的给定不在合理的数量级内，则会导致多参数反演 的不收敛。因此，为了提高差分计算的稳定性，使用参数的标准偏差p/σ
p
来限 制差分步长：
[0117][0118][0119]
将公式(13)和(14)带入公式(12)，对于第l次迭代有：
[0120]
κ(m0 δm
l
) g(m0 δm
l
)δm(i)＝0,
ꢀꢀꢀ
(19)
[0121]
则上式可写为：
[0122]
hδm(i)＝-g
t
κ,
ꢀꢀꢀ
(20)
[0123]
其中h为hessian矩阵，有h＝g
t
g。在公式(20)中加入阻尼因子，得 到了利用二阶lm迭代算法的模型更新量：
[0124]
δm＝-(h λi)-1gt
κ
ꢀꢀꢀ
(21)
[0125]
在每次迭代中，将δm加到初始模型m0中以更新m
ite
，然后使用公式(6) 更新a
syn
，直到j(m)的l2-范数足够小以停止反演。
[0126]
对应地，pp波、ps1波和ps2波的层匹配原理下：
[0127]
步骤1，如图3所示，利用零偏vsp地震观测记录的初至波，拾取零偏的 纵波初至时间t0
pfb
和检波器号对应的关系(图3a)，以及p波和s波的层速度 模型(图3b)。这对于后期的层匹配、ava道集转换和初始模型的建立至关重 要。因为零偏vsp地震观测记录的入射角近似为0，地震波在垂向传播时主要 表现介质的各向同性性质。因此获得的层速度模型可以近似为各向同性背景的 速度。
[0128]
步骤2，如图4所示，这里以1012m偏移距的某一炮为例，对于walkaround 三分量(ps波、ps1波和ps2波)vsp地震观测记录，进行水平旋转、噪声 压制和球面扩散补偿等保幅方法在保存三分量之间的振幅特征关系的同时进 行了预处理。
[0129]
步骤3，如图5所示，对步骤2预处理后数据进行上下行波分离，得到分 离后的直达波(下行波)和反射波(上行波)。
[0130]
步骤4在步骤3中得到的直达波记录中可以拾取不同检波器号对应的纵 波初至时间t
pfb
，通过此获得了检波器号-t
pfb
关系。需要注意的是，对于偏移 距不同的炮得到的检波器号-t
pfb
关系不同。
[0131]
步骤5，如图6所示，利用步骤4得到的检波器号-t
pfb
关系，将步骤3中 得到的反射波波场根据检波器号压缩至t
pfb
域中。3个分量上匹配良好的强反 射轴代表的地质标志层位分别是t1j1(示例为：三叠系湘西组-嘉陵江组界面)、 目标层t1j1和t1f1(二叠系-三叠系界面)的底部界面。
[0132]
步骤6此时已完成了z、r、t分量之间的层匹配，但是因为不同炮的检 波器号-t
pfb
关系不同，不同炮之间还没有进行层匹配。如图7所示，将步骤5 中t
pfb
域的z、r、t分量记录根据步骤1中的检波器号-t0
pfb
关系，将所有炮 的z、r、t分量匹配至t0
pfb
域中。此时完成了不同炮的z、r、t分量的层匹 配。经发现，由于裂缝发育导致分裂的s1波和s2波以不同的速度混合，在曲 线包围的范围内，反射同相轴的数量在r分量中比在z分量中更多。此外，由 于分裂横波的偏振方向受裂缝方位角的影响，且不平行于z-r平面，因此图7 呈现t分量记录了分裂横波的投影。所有这些现象都分布在目标地层t1j1附 近，即实际的裂缝储层附近。
[0133]
步骤7,z分量数据可以认为是pp波，但r、t分量是ps1波和ps2波的 混合，需要进行波场分离。如图8所示，将步骤6中t0
pfb
域的r、t分量进行 alford旋转扫描，得到快慢横波时延τ和裂缝方位为了消除上覆层裂缝发 育对下伏层alford旋转结果产生的影响，这里使用对r分量和t分量进行 层剥离，得到ps1和ps2波记录，并通过在快慢横波时延τ消除ps2波时差， 使之与ps1波记录相匹配。如图9所示，在裂缝发育的范围内，分裂横波的波 场混合现象被消除，pp波和分离的ps1波的反射同相轴比z和r分量的匹配 性更好。
[0134]
至此，得到了匹配t0
pfb
域的pp、ps1、ps2波炮记录。
[0135]
对应地，参考图3(b)所示，ava道集生成原理如下：
[0136]
步骤8,在传统的ava道集生成中，入射角是根据层速度进行射线跟踪计 算的。然而，在叠前数据处理中，s2波层速度很难得到，对ava道集生成造 成困难。通过上述步骤7已经将pp波、ps1波、ps2波波的波场拉平，并消除 了各向异性引起的ps2波时差，将其匹配至
t0
pfb
域。对于ps2波，在经过时 差校正后，已经可以使用各向同性背景的s波速度进行模拟。因此，如图3b 所示，采用步骤1获得的各向同性的p波和s波层速度进行射线追踪。通过射 线追踪获得的每个炮的pp和ps波的入射角和反射点坐标道集，临近自然插值 后结果如图10所示。需要注意的是，ps波表示分裂前的上行转换横波。因为 pp波和ps波的具有相同的地质反射点，pp波和ps波反射点的z坐标相同。
[0137]
步骤9，根据步骤8获得的入射角和反射点坐标道集，对步骤7得到的所 有匹配的pp波、ps1波、ps2波的炮记录进行深度域ava道集转换。
[0138]
对应地，参考图2c所示，初始模型建立的原理：
[0139]
步骤10，因为零偏移距vsp数据中p波和s波在各向同性平面中传播， 并且在入射角为零时不受垂直裂缝的影响。因此采用步骤1的层速度(图3b) 作为各向同性背景的纵波速度和横波速度v
p
和vs。对应的反演参数λ,μ,和ρ 的初始模型可以由下公式计算：
[0140]
ρ
gardner
＝-0.0296v
p2
0.461v
p
0.963.
ꢀꢀꢀ
(22)
[0141][0142][0143]
其中，图11为反演参数λ,μ,和ρ的初始模型的示意图。
[0144]
步骤11，将步骤7中的快慢横波时延τ和裂缝方位(图8)通过类似步 骤9的方法利用入射角和反射点坐标道集移位至图12成像范围的空间坐标中。 井旁结果如图12所示。裂缝方位在后续会用于计算反射系数，但是对于初 始模型e，需要用τ换算。
[0145]
步骤12，裂缝密度e的初始模型公式推导。
[0146]
对于thomsen定义的各向异性参数γ有(ruger,2002)：
[0147][0148]
将(24)式与公式(2)联立，得到：
[0149][0150]
已知ps1波旅行时tps1：
[0151][0152]
其中l为alford旋转时扫描的时窗利用初至速度转换成的深度窗。则快慢 波时差：
[0153][0154]
利用公式(2、25、27)可以计算得到裂缝密度公式：
[0155]
[0156]
将步骤1得到的层速度和步骤11得到的井旁τ带入公式(28)计算的e如 图13所示。
[0157]
实际数据反演结果
[0158]
对四川泸州附近某vsp数据进行测试(前述图使得用也是此数据)。由步 骤9得到图14所示的深度域ava道集(图14)，生成的实际数据ava道集与 利用初始模型计算的合成ava道集匹配性较好。ava道集叠加后结果如图16 所示，作为反演范围，曲线为速度比。反演采用由ava道集中提取的角度独立 的深度域统计子波。利用反演发成，得到如图17所示的反演结果。
[0159]
图17代表岩性特征的各向同性背景参数反演结果与地质界面一致。高裂 缝密度范围为2420
–
2490m，这与地质背景和生产数据验证的裂缝性储层的深 度一致。
[0160]
图18是本发明实施例提供的一种裂缝性储层的多波联合avo反演装置的 结构示意图，包括：
[0161]
多波匹配模块1810，用于将裂缝性储层对应ps波、ps1波和ps2波的 vsp地震观测记录压缩至零偏vsp纵波初至的时间域中进行ava道集转换， 得到ava地震观测记录。
[0162]
模型函数构建模块1820，用于基于ava地震观测记录，构建初始模型和 所述初始模型对应的反演目标函数，其中，所述初始模型用于根据ava道集 数据中包含各向异性参数的测井数据计算裂缝参数，所述反演目标函数为贝叶 斯理论函数。
[0163]
反演方程构建模型1830，用于构建所述初始模型对应的反演方程，所述反 演方程用于基于所述初始模型按照精确反射系数正演模拟ava地震合成记录， 并根据正演算子的偏导以及所述反演目标函数的二阶近似解，计算所述初始模 型的梯度以对所述初始模型进行迭代，以根据所述初始模型的迭代结果反演所 述裂缝性储层的各向异性参数。
[0164]
反演执行模块1840，用于将所述初始模型代入至所述反演方程中进行反 演，得到所述裂缝性储层的各向异性参数。
[0165]
裂缝状态确定模块1850，用于基于所述裂缝性储层的各向异性参数，确定 裂缝性储层的裂缝状态。
[0166]
可选地，多波匹配模块1810将裂缝性储层对应ps波、ps1波和ps2波的 vsp地震观测记录压缩至零偏vsp纵波初至的时间域中进行ava道集转换， 得到ava地震观测记录，包括：基于纵波初至与检波器深度的关系，将pp波、 ps1波和ps2波的vsp地震观测记录匹配压缩至零偏vsp纵波初至的时间域 中，并通过alford旋转，得到快慢横波时延和裂缝方位；将压缩至零偏vsp纵 波初至的时间域中的vsp地震观测记录进行r分量和t分量的层剥离，得到 ps1波和ps2波各自的vsp地震观测记录；使用快慢横波时延校正vsp地震 观测记录中ps2波相对ps1波的时差；基于vsp地震观测记录中p波和s波 的层速度对pp波、ps1波和ps2波进行射线追踪，得到每个炮的pp波、ps1 波和ps2波的入射角和反射点的坐标道集；将每个炮的pp波、ps1波和ps2 波的入射角和反射点的坐标道集进行ava道集转换，得到ava地震观测记录。
[0167]
可选地，所述初始模型用于根据ava道集数据中的层速度和快慢横波时 延，计算裂缝密度。
[0168]
可选地，所述初始模型针对所述裂缝密度的推导公式为：：
[0169]
其中，e表示裂缝密度，l表示alford旋转时扫描的时窗利用初至速度转 换成的深度窗,v
s1
表示s1波的横波速度，τ表示快慢横波时延,g＝μ/(λ 2μ)，μ 和λ为lam
é
反演参数。
[0170]
可选地，所述目标函数的公式为：j(m)＝[a
prs-a
syn
(m
ite
,i)]
t
[a
prs
‑ꢀasyn
(m
ite
,i)]
t
k(m
ite
)；其中，j(m)表示目标函数，i和k分别表示入射角 和采样点；m表示采样点k对应的反演参数向量m＝[λk,μk,ρk,ek]
t
，ρ表示岩石 密度，a
prs
表示反演中输入的ava记录，m
ite
表示迭代更新的反演参数， k(m
ite
)表示正则化柯西分布项，a
syn
表示用于更新初始模型的ava记录， a
syn
(m
ite
,i)＝w*r(m
ite
,i)，w
pp
、w
ps1
、 w
ps2
分别表示pp波、ps1和ps2波的子波矩阵；
[0171]
所述初始模型按照精确反射系数正演模拟ava地震合成记录的公式为 r(m
ite
,i)＝[r
pp
(m
ite
,i),r
ps1
(m
ite
,i),r
ps2
(m
ite
,i)]
t
，r
pp
(m
ite
,i)、 r
ps1
(m
ite
,i)和r
ps2
(m
ite
,i)分别表示pp波、ps1和ps2波对应的使用m
ite
计算的在入射角i时的反射系数矩阵。
[0172]
可选地，所述反演方程中对初始模型进行迭代公式为：m
ite
＝m0 δm；其 中，m0表示初始模型，δm表示迭代增量，δm＝-(h λi)-1gt
κ，h表示 hessian矩阵，i表示对角矩阵，κ表示目标函数的jacobian矩阵示目标函数的jacobian矩阵g表示目标函数的jacobian矩阵的二阶形式,的二阶形式,分别表示入射角i下的pp波、ps1波和ps2波的jacobian矩阵。
[0173]
可选地，所述反演方程在所述目标函数的j(m)的l2-范数小于预设阈值时， 结束对所述初始模型的迭代。
[0174]
图18是本说明书的一个实施例电子设备的结构示意图。请参考图18，在 硬件层面，该电子设备包括处理器，可选地还包括内部总线、网络接口、存储 器。其中，存储器可能包含内存，例如高速随机存取存储器(random-accessmemory，ram)，也可能还包括非易失性存储器(non-volatile memory)，例如 至少1个磁盘存储器等。当然，该电子设备还可能包括其他业务所需要的硬件。
[0175]
处理器、网络接口和存储器可以通过内部总线相互连接，该内部总线可以 是isa(industry standard architecture，工业标准体系结构)总线、pci(peripheralcomponent interconnect，外设部件互连标准)总线或eisa(extended industrystandard architecture，扩展工业标准结构)总线等。所述总线可以分为地址总 线、数据总线、控制总线等。为便于表示，图18中仅用一个双向箭头表示，但 并不表示仅有一根总线或一种类型的总线。
[0176]
存储器，用于存放程序。具体地，程序可以包括程序代码，所述程序代码 包括计算机操作指令。存储器可以包括内存和非易失性存储器，并向处理器提 供指令和数据。
[0177]
其中，处理器从非易失性存储器中读取对应的计算机程序到内存中然后运 行，在逻辑层面上形成上述多波联合avo反演装置。对应地，处理器，执行 存储器所存放的程序，
并具体用于执行以下操作：
[0178]
将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压缩至零 偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记录；
[0179]
构建初始模型和所述初始模型对应的反演目标函数，其中，所述初始模型 用于根据ava道集数据中包含各向异性参数的测井数据计算裂缝参数，所述 反演目标函数为贝叶斯理论函数；
[0180]
构建所述初始模型对应的反演方程，所述反演方程用于基于所述初始模型 按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导以及 所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初始模型 进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向异性参 数；
[0181]
将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性储层的 各向异性参数；
[0182]
基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝状态。
[0183]
上述如本说明书图1所示实施例揭示的方法可以应用于处理器中，由处理 器实现。处理器可能是一种集成电路芯片，具有信号的处理能力。在实现过程 中，上述方法的各步骤可以通过处理器中的硬件的集成逻辑电路或者软件形式 的指令完成。上述的处理器可以是通用处理器，包括中央处理器(centralprocessing unit，cpu)、网络处理器(network processor，np)等；还可以是 数字信号处理器(digital signal processor，dsp)、专用集成电路(applicationspecific integrated circuit，asic)、现场可编程门阵列(field－programmablegate array，fpga)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、 分立硬件组件。可以实现或者执行本发明实施例中的公开的各方法、步骤及逻 辑框图。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器 等。结合本发明实施例所公开的方法的步骤可以直接体现为硬件译码处理器执 行完成，或者用译码处理器中的硬件及软件模块组合执行完成。软件模块可以 位于随机存储器，闪存、只读存储器，可编程只读存储器或者电可擦写可编程 存储器、寄存器等本领域成熟的存储介质中。该存储介质位于存储器，处理器 读取存储器中的信息，结合其硬件完成上述方法的步骤。
[0184]
应理解，本发明实施例的电子设备可以使业务处理装置实现对应于图1所 示方法中的步骤和功能。由于原理相同，本文不再赘述。
[0185]
当然，除了软件实现方式之外，本说明书的电子设备并不排除其他实现方 式，比如逻辑器件抑或软硬件结合的方式等等，也就是说以下处理流程的执行 主体并不限定于各个逻辑单元，也可以是硬件或逻辑器件。
[0186]
此外，本发明实施例还提出了一种计算机可读存储介质，该计算机可读存 储介质存储一个或多个程序，该一个或多个程序包括指令。
[0187]
其中，上述指令当被包括多个应用程序的便携式电子设备执行时，能够使 该便携式电子设备执行图1所示方法的步骤，包括：
[0188]
将裂缝性储层对应ps波、ps1波和ps2波的vsp地震观测记录压缩至零 偏vsp纵波初至的时间域中进行ava道集转换，得到ava地震观测记录。
[0189]
构建初始模型和所述初始模型对应的反演目标函数，其中，所述初始模型 用于根
据ava道集数据中包含各向异性参数的测井数据计算裂缝参数，所述 反演目标函数为贝叶斯理论函数。
[0190]
构建所述初始模型对应的反演方程，所述反演方程用于基于所述初始模型 按照精确反射系数正演模拟ava地震合成记录，并根据正演算子的偏导以及 所述反演目标函数的二阶近似解，计算所述初始模型的梯度以对所述初始模型 进行迭代，以根据所述初始模型的迭代结果反演所述裂缝性储层的各向异性参 数。
[0191]
将所述初始模型代入至所述反演方程中进行反演，得到所述裂缝性储层的 各向异性参数。
[0192]
基于所述裂缝性储层的各向异性参数，确定裂缝性储层的裂缝状态。
[0193]
本领域技术人员应明白，本说明书的实施例可提供为方法、系统或计算机 程序产品。因此，本说明书可采用完全硬件实施例、完全软件实施例或结合软 件和硬件方面的实施例的形式。而且，本说明书可采用在一个或多个其中包含 有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd
‑ꢀ
rom、光学存储器等)上实施的计算机程序产品的形式。
[0194]
上述对本说明书特定实施例进行了描述。其它实施例在所附权利要求书的 范围内。在一些情况下，在权利要求书中记载的动作或步骤可以按照不同于实 施例中的顺序来执行并且仍然可以实现期望的结果。另外，在附图中描绘的过 程不一定要求示出的特定顺序或者连续顺序才能实现期望的结果。在某些实施 方式中，多任务处理和并行处理也是可以的或者可能是有利的。
[0195]
以上仅为本说明书的实施例而已，并不用于限制本说明书。对于本领域技 术人员来说，本说明书可以有各种更改和变化。凡在本说明书的精神和原理之 内所作的任何修改、等同替换、改进等，均应包含在本说明书的权利要求范围 之内。此外，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有 其他实施例，都应当属于本文件的保护范围。