一种高精度的三维弹性波全波形反演方法

本发明涉及一种地下三维弹性介质的成像方法，尤其涉及一种对于石油、煤矿等资源勘探的三维弹性波全波形反演方法，属于地震勘探速度建模。

背景技术：

1、地下资源的开采离不开目标靶区的精准识别，诸如油气勘探中的储层精准识别并布置恰当的钻井，煤炭开采中关键层位的精准识别及其空间内隐蔽致灾因素的圈定等。为了保障地下资源的安全开采，探明目标区域的地质构造分布及其岩性属性变化是一项必须且迫切的任务。地球物理勘探中的地震勘探具备精细构造的探查能力，在资源勘查中扮演着主导地位。作为一种高精度的地震勘探技术，全波形反演通过分析地震波在地下的传播情况，推断出地下介质的结构和性质。

2、由于地质体在空间中是客观的三维分布，二维全波形反演无法准确地刻画其空间分布。三维全波形反演技术已经成为地球物理勘探领域的研究热点之一。相比于二维全波形反演，三维全波形反演需要更为复杂的方程及其离散方程的求解，所需的计算资源和计算效率是二维问题的数十倍乃至上百倍，同时反演的复杂性和多解性增加。

3、目前的三维反演相比于二维反演较为缓慢。其中由于时间域反演存在着主频优势的影响，会存在频率信息使用不充分的情况。进而对三维全波形反演的精度造成较大影响，尤其是全波形反演对初始模型有着强依赖性，会导致计算精度的降低。究其原因是数据中的低频信息无法有效地被用来更新模型中的高波数信息，对于全波形反演这类迭代算法，当低频反演不准确时，高波数的更新会存在很大的偏差；而导致数据中低频信息无法被有效地利用的原因有两点：1.数据本身存在着非常有限且能量较弱的低频信息；2.目前的全波形反演方法无法有效地使用这部分弱低频信号。尽管目前常用的基于卷积类型的多尺度反演策略在一定程度上缓解上述问题，但仍然无法解决上述问题，且还存在反演时参考道质量差导致反演结果不理想的问题，上述这些问题都会降低反演精度。

4、因此，如何实现一种高精度三维全波形反演方法，能将观测数据中低频信息较好的用于更新模型中的高波数信息，并通过多尺度反演的方式，最终有效保证三维反演的精度。

技术实现思路

1、针对上述现有技术存在的问题，本发明提供一种高精度的三维弹性波全波形反演方法，能将观测数据中低频信息较好的用于更新模型中的高波数信息，并通过多尺度反演的方式，最终有效保证三维反演的精度。

2、为了实现上述目的，本发明采用的技术方案是：一种高精度的三维弹性波全波形反演方法，具体步骤为：

3、步骤一、获取观测系统采集的观测数据，先给定一个初始速度场模型，从观测数据中选取数个离散频点，结合常规带限子波震源中多个频点进行振幅归一化处理后形成一个新的震源；

4、步骤二、利用步骤一得出的新震源，采用三维弹性波方法计算正传波场，并结合观测系统记录的波场数据得到合成数据；

5、步骤三、根据步骤二获取的合成数据建立最小二乘目标函数并进行计算，计算后的数据残差作为反传震源；

6、步骤四、根据步骤三获取的反传震源，计算反传震源条件下的反传波场；

7、步骤五、利用步骤二获取的正传波场和步骤四获取的反传波场构建模型参数的梯度；

8、步骤六、根据步骤五得到的梯度，采用具有二阶收敛性的l-bfgs方法计算海森矩阵的逆矩阵并对步骤一的初始速度场模型进行持续更新，直至数据残差达到设定阈值或迭代更新次数达到预设值，停止迭代更新获得最终的速度场模型参数，此时的速度场模型参数即为所需区域的三维高精度成像。

9、进一步，所述步骤一中归一化形成新的震源的具体公式为：

10、

11、其中，f表示傅里叶变换，表示频率为ωi的谱的傅里叶逆变换，xs表示震源的位置，表示信号a的归一化，n表示使用的频率个数，并且ω＝(ω1,ω2,…,ωn)，s(xs,t)为常规带限子波震源，t表示时间，i表示循环下标。

12、进一步，所述步骤二中采用三维弹性波方法计算正传波场的具体计算公式为：

13、

14、其中，ρ表示介质的密度，正传波场中质点沿着三个方向的振动速度为v＝[vx vyvz]，正传波场中质点三个方向的正应力和三个方向的剪应力波场为τ＝[τxx τyy τzz τxy τxzτyz]；表示时间偏导数算子，t表示时间，链接矩阵与弹性参数相关的虎克矩阵为

15、生成合成数据的具体公式为：

16、ψsyn(xr,t；m)＝ψ(x,t；m)δ(xr),

17、其中，x表示计算区域，δ(xr)表示脉冲函数，合成数据ψ(x,t；m)＝[v τ]t是在检波点xr处记录的波场数据。

18、进一步，所述步骤三中最小二乘的目标函数为：

19、

20、其中，模型参数m＝[α β]t，α和β分别为纵横波速度，ψsyn表示合成数据，ψobs表示观测数据，上述目标函数计算后对应的数据残差为：

21、

22、进一步，所述步骤四中反传波场的计算公式为：

23、

24、其中，反传波场反传波场中质点沿着三个方向的振动速度为υ＝[υxυy υz]；反传波场中质点三个方向的正应力和三个方向的剪应力波场为σ＝[σxx σyy σzz σxyσyz σxz]。

25、进一步，所述步骤五中利用正传波场和反传波场构建模型参数的梯度，采用的模型参数为拉梅参数，则拉梅参数的梯度表示为：

26、

27、

28、其中，λ和μ为拉梅参数，nt表示记录的时间步数，和表示分别对x，y和z方向的空间偏导数算子；通过拉梅参数与纵横波速度之间的关系λ+2μ＝ρα2和μ＝ρβ2，α和β分别为纵横波速度，求取目标函数对纵横波速度的梯度：

29、

30、

31、进一步，所述步骤六中采用具有二阶收敛性的l-bfgs方法计算海森矩阵的逆矩阵并对步骤一的初始速度场模型进行持续更新的具体过程为：

32、采用l-bfgs方法求取海森矩阵的逆矩阵dk+1，具体公式为：

33、

34、其中，i为单位对角矩阵，梯度变化量yk＝gk+1-gk，速度场模型的变化量sk＝mk+1-mk，该公式中的梯度由步骤五获得；逆矩阵dk+1公式中初始的dk为单位矩阵，其中下标k表示迭代次数，模型更新如下：

35、mk+1＝mk-dk+1gk+1.

36、采用上述公式，对初始速度场模型进行持续更新，直至数据残差达到设定阈值或迭代更新次数达到预设值，停止迭代更新获得最终的速度场模型参数。

37、与现有技术相比，本发明具有如下优点：

38、1、本发明先从观测数据中选取数个离散频点，结合常规带限子波震源及其位置进行归一化处理后形成一个新的震源，然后进行反演时采用该新型子波震源，从而克服了反演方法中使用常规带限子波震源数据而引起的低波数更新不足的缺陷。

39、2、本发明通过结合新震源获得正传波场和反传波场，进而构建模型参数的梯度，通过这种方式使用由低到高的几组频组反演，可以有效实现多尺度反演，避免了卷积型多尺度反演中的参考道质量差导致反演结果不理想的问题。

40、3、本发明采用l-bfgs方法结合获取的梯度对初始速度场模型进行持续更新，采用更新后的模型再进行反演，这种方式能实现对反演中低频信号的能量进一步提升，增强了低频信号对模型低波数信息的贡献，从而有效地减弱反演对初始速度场模型的依赖性，进而提升了反演精度，实现了高精度反演目标。

41、4、本发明的三维高精度全波形反演方案，相比于二维反演成像方法，可以获得准确地探查地下三维介质的空间分布，能为地下资源的勘探开发提供数据指导作用。