非线性射线走时层析成像方法、电子设备及介质与流程

本发明涉及地球物理勘探中的地震走时层析成像领域，更具体地，涉及一种非线性射线走时层析成像方法、电子设备及介质。

背景技术：

1、地震走时层析的目的是找到一个使理论合成走时与观测走时之差在某种意义下达到极小的“最佳”模型。目前基于线性反演思想的方法主要有两种实现方式。第一种是传统的求解层析方程组的方法。这种方法求解的方程组系数矩阵规模庞大，且非常稀疏，病态严重，需要采用正则化来保证反演的正确性与稳定性，因此反演精度严重依赖于方程组的求解方法。第二种是基于局部最优化思想，通过对目标函数的局部线性化近似获得速度的修改方向，然后通过迭代实现该非线性反演问题的线性化求解。但无论哪种实现方式，线性反演方法的主要缺点是对初始表层速度模型依赖性强，当初始表层速度模型选择不合适时，很容易陷入局部极值。尽管如此，这些方法在实际应用中仍可以取得不错的效果。

2、在第二种方法中，每一次迭代中速度修改方向的计算是关键。利用伴随状态法或散射积分法可以获得速度的修改方向。前者利用地震波的一次正传和数据残差的一次反传即可获得目标函数的梯度，无需存储大型矩阵，但其预条件在射线走时层析成像中不易实现。后者通过显式地计算核函数，并与走时残差向量相乘实现梯度的计算。但与传统的层析成像方法一样，散射积分法面临内存占用大的问题，特别是在利用hessian矩阵时。

3、hessian矩阵的利用在局部最优化反演中对提高反演精度与多参数解耦有非常重要的作用。核函数表达的hessian矩阵包含两项，利用准确的hessian矩阵构造的方向称为牛顿方向，只保留hessian第一项的矩阵称为近似hessian矩阵，相对应的方向为高斯-牛顿方向。实际应用中，无论是牛顿方向还是高斯-牛顿方向都涉及到hessian矩阵求逆的操作，因此往往难以准确构造。更通常的做法是利用近似hessian矩阵的对角阵来构造hessian矩阵的逆，所对应的方向是拟牛顿方向的一种。由于hessian矩阵的逆是直接作用在梯度上面，因此这种做法又称为梯度方向的预条件。预条件的实施与否对层析有很大的影响。

4、现有技术在波形反演研究中提出了一种改进的散射积分算法。该方法将大规模的核函数-向量乘运算表示为具有明确物理含义的向量-标量乘的累加运算，因此避免了大规模矩阵的存储，使得传统散射积分法波形反演更加实用化。同时，该方法可以方便地实现hessian矩阵主对角线元素的计算，进而可以方便实现梯度的预条件。但是，与传统的层析成像方法一样，散射积分法也面临内存占用大的问题，特别是在利用hessian矩阵时。

5、因此，有必要开发一种非线性射线走时层析成像方法、电子设备及介质。

6、公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

技术实现思路

1、本发明提出了一种非线性射线走时层析成像方法、电子设备及介质，其能够避免传统层析成像中的大规模矩阵的存储和预条件实现不方便的问题，进而实现带预条件的最速下降法地震走时层析成像。

2、第一方面，本公开实施例提供了一种非线性射线走时层析成像方法，包括：

3、步骤1：对原始地震数据进行去噪和反褶积的预处理，获得预处理地震数据；

4、步骤2：针对所述预处理地震数据进行初至拾取，获得初至波到达时；

5、步骤3：根据地下先验信息建立初始表层速度模型，设定反演参数；

6、步骤4：根据初至波到达时与所述初始表层速度模型通过散射积分法进行非线性射线走时层析成像反演，获得最终表层速度模型并成像。

7、优选地，所述步骤4包括：

8、步骤41：在所述初始表层速度模型上进行射线追踪，获取理论合成初至走时与射线路径，保存所述射线路径；

9、步骤42：判断所述理论合成初至走时与所述初至波到达时是否匹配，若是，则以所述初始表层速度模型作为所述最终表层速度模型，若否，则执行步骤43；

10、步骤43：逐个读取每一个炮检对对应的射线路径作为该炮检对对应的核函数；

11、步骤44：通过核函数-标量乘获得每一个炮检对所对应的梯度，进而计算整个观测系统的全局梯度，并计算所述全局梯度的预条件算子；

12、步骤45：计算预条件最速下降方向、方向更新步长与模型更新量；

13、步骤46：更新表层速度模型，以更新后的表层速度模型为初始表层速度模型，重复步骤41-46，直到获得所述最终表层速度模型。

14、优选地，所述步骤41包括：

15、步骤411：在所述初始表层速度模型网格剖分的网格边界上设定多个关键点；

16、步骤412：自激发点开始以矩形为波前，按惠更斯原理逐层向外扩展，直至扫描完成整个模型空间，再由外向里收缩，直至扫描到激发点，如此重复直至模型空间的最小旅行时场不变为止；

17、步骤413：检波点处的初至波到达时即为所述理论合成初至走时。

18、优选地，所述步骤44包括：

19、步骤441：获取每一个炮检对所对应的射线路径，进而获得所述射线路径与该炮检对对应的走时残差乘，得到该炮检对对应的梯度以及射线路径向量的元素平方向量；

20、步骤442：将所有炮检对的梯度累加求和，得到整个观测系统的全局梯度；

21、步骤443：将所有的炮检对的核函数元素取平方并进行累加，计算全局梯度的预条件算子。

22、优选地，通过公式(1)将所有炮检对的梯度累加求和，得到整个观测系统的全局梯度：

23、

24、其中，ktδt为全局梯度，kij为第i对炮检对对应的核函数在第j个模型网格内的核函数值，δti为旅行时差，kijδti为第i对炮检对的梯度，kt为影响走时的核函数矩阵、δt为旅行时差。

25、优选地，通过公式(2)计算全局梯度的预条件算子：

26、

27、其中，h0为预条件算子。

28、优选地，所述步骤45中，通过公式(3)计算预条件最速下降方向：

29、

30、其中，p为预条件最速下降方向。

31、优选地，所述步骤45中方向更新步长的计算方法为：

32、根据先验信息获取所述初始表层速度模型与真实地下模型的速度差，将速度差除以预先设定的迭代次数，得到每次迭代可更新的最大速度值δvmax，利用待定系数法求取所述方向更新步长t，使得max{|p·t|}＝δvmax。

33、作为本公开实施例的一种具体实现方式，

34、第二方面，本公开实施例还提供了一种电子设备，该电子设备包括：

35、存储器，存储有可执行指令；

36、处理器，所述处理器运行所述存储器中的所述可执行指令，以实现所述的非线性射线走时层析成像方法。

37、第三方面，本公开实施例还提供了一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现所述的非线性射线走时层析成像方法。

38、其有益效果在于：

39、(1)无需存储核函数矩阵即可方便地实现目标函数梯度的计算，在任何时刻只需存储单个核函数—射线路径；

40、(2)计算效率高，易于并行：将大规模的核函数-向量乘运算表示为具有明确物理含义的向量-标量乘的累加运算，无需利用svd或者lsqr方法解大规模矩阵方程组，计算时间大大减少；

41、(3)方便预条件，使得计算效率进一步提高，反演精度也大大提高：该方法无需存储hessian矩阵即可方便实现预条件，预条件算子的应用可以加快收敛速度，并能够获得地下深部介质的信息；

42、(4)反演过程更稳定：优化的梯度导引算法和预条件算子的使用使得反演过程更稳定。

43、本发明的方法和装置具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。