CN111638551A - 地震初至波走时层析方法及装置 - Google Patents

地震初至波走时层析方法及装置 Download PDF

Info

Publication number
CN111638551A
CN111638551A CN201910154540.8A CN201910154540A CN111638551A CN 111638551 A CN111638551 A CN 111638551A CN 201910154540 A CN201910154540 A CN 201910154540A CN 111638551 A CN111638551 A CN 111638551A
Authority
CN
China
Prior art keywords
slowness
updating
travel time
model
ray
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN201910154540.8A
Other languages
English (en)
Inventor
马青坡
黄兴贵
穆志平
王彦铎
潘英杰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China National Petroleum Corp
BGP Inc
Original Assignee
China National Petroleum Corp
BGP Inc
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China National Petroleum Corp, BGP Inc filed Critical China National Petroleum Corp
Priority to CN201910154540.8A priority Critical patent/CN111638551A/zh
Publication of CN111638551A publication Critical patent/CN111638551A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/303—Analysis for determining velocity profiles or travel times
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/303—Analysis for determining velocity profiles or travel times
    • G01V1/305—Travel times
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • G01V1/362—Effecting static or dynamic corrections; Stacking

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

本发明提供了一种地震初至波走时层析方法及装置，其中，该方法包括：获取地震初至波的初至时间；建立初始速度模型；执行第一迭代步骤：按炮检对进行任务分组，并行进行射线追踪、走时残差和慢度更新量的计算；将多个射线密度汇总，得到射线密度汇总值；将多个走时残差汇总，得到走时残差汇总值；将多个慢度更新量汇总，得到慢度更新量汇总值；利用射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新；在满足第一预设条件时，执行更新速度模型的步骤；执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，第二迭代步骤包括第一迭代步骤以及更新速度模型的步骤。上述技术方案提高了初至波走时层析的效率。

Description

地震初至波走时层析方法及装置
技术领域
本发明涉及地球物理地震勘探技术领域，特别涉及一种地震初至波走时层析方法及装置。
背景技术
在地球物理地震勘探领域中，初至波走时层析是建立近地表速度模型的主要方法之一，其建立的近地表速度模型可以用于层析静校正量的计算，也可以作为叠前深度偏移和全波形反演的近地表速度模型。在数据量较小的情况下，传统的保存射线路径方法(将射线追踪出的射线路径保存为文件)或全内存操作方法可以满足计算效率的要求。但随着地震采集密度的日益增大，有效的地震道数成倍增加，这样一方面数据量成几何倍数增加，导致了射线路径文件非常庞大，有时达到TB级，全内存操作已成为不可能的事情。而频繁的读写这一庞大的射线路径文件，给计算效率带来了很大的挑战，导致计算效率显著降低，无法满足常规生产对计算效率的要求。
针对上述问题，目前尚未提出有效的解决方案。
发明内容
本发明实施例提供了一种地震初至波走时层析方法，用以提高初至波走时层析的效率，该方法包括：
获取地震初至波的初至时间；
建立地震初至波走时层析反演的初始速度模型；
根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
本发明实施例还提供了一种地震初至波走时层析装置，用以提高初至波走时层析的效率，该装置包括：
获取单元，用于获取地震初至波的初至时间；
建立单元，用于建立地震初至波走时层析反演的初始速度模型；
第一迭代单元，用于根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
速度模型更新单元，用于在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
第二迭代单元，用于执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
本发明实施例还提供了一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述地震初至波走时层析方法。
本发明实施例还提供了一种计算机可读存储介质，所述计算机可读存储介质存储有执行地震初至波走时层析方法的计算机程序。
本发明实施例提供的技术方案通过：获取地震初至波的初至时间；建立地震初至波走时层析反演的初始速度模型；根据初至时间、初始速度模型，执行如下第一迭代步骤：按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；利用射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，第二迭代步骤包括第一迭代步骤以及更新速度模型的步骤；近地表速度模型用于地震初至波走时层析静校正量的计算，实现了：以地震波的初至时间为基础，射线追踪后直接进行慢度更新量的求取，不需要保存射线路径，其计算量仅与有效的炮检对数量和迭代次数有关，具有计算效率高、易于并行，对大数据和大模型的适应性强等特点，提高了初至波走时层析的效率。
附图说明
此处所说明的附图用来提供对本发明的进一步理解，构成本申请的一部分，并不构成对本发明的限定。在附图中：
图1是本发明实施例中地震初至波走时层析方法的流程示意图；
图2是本发明另一实施例中地震初至波走时层析方法的流程示意图；
图3是本发明实施例中一个炮检对，炮点到检波点的射线路径以及穿过每个网格的射线路径长度示意图；
图4是利用本发明实施例提供的地震初至波走时层析方案效果与传统串并行效果的对比示意图；
图5是本发明实施例中地震初至波走时层析装置的结构示意图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明白，下面结合实施方式和附图，对本发明做进一步详细说明。在此，本发明的示意性实施方式及其说明用于解释本发明，但并不作为对本发明的限定。
在介绍本发明实施例之前，首先对本发明实施例涉及的专业技术名词进行介绍：
1、炮检对：一个炮点和对应的一个检波点组成一个炮检对。
2、射线及射线追踪：按照几何地震学理论，地震波从炮点传播检波点的路径可用一条射线表示。射线追踪就是根据地震学理论得到炮点到检波点的射线路径，从而构建层析反演方程，射线追踪的精度和效率决定初至波走时层析的精度和效率。
3、地震道：单个检波器记录地面震动形成的数据。
4、网格：对所求区域进行剖分的单元，一般采用规则网格进行剖分，如二维采用矩形，三维采用长方体等。
5、初至波走时层析：一种按照射线理论，结合层析反演技术，采用初至波的走时进行近地表速度反演的方法。
发明人发现：面对背景技术存在的技术问题，现有技术中有些层析算法采用分块计算的方法，但存在着块间耦合的问题，导致了整体模型反演不合理；与此同时，采集密度的增加为反演更为精细的近地表速度模型提供了可能，这就需要采用更小的网格尺寸进行层析反演，而网格尺寸的减小会显著增大射线路径文件的大小，如在其他参数一致的情况下，网格尺寸在三个方向上减小一倍，会导致射线路径文件增大八倍，这会更加影响计算效率。
由于考虑到如上技术问题，发明提出了一种新的地震初至波走时层析方案，该方案结合初至波走时层析反演算法的特点以及目前多核多线程计算机的硬件特点，提出了一种适用于海量数据的初至波走时层析方法和装置，其计算量仅仅与实际的有效炮检对个数和迭代次数成正比，其内存需求与模型大小成正比，且易于并行化，计算效率得到显著提升，非常适合大数据量和大模型初至波走时层析反演的计算。
具体地，本发明涉及地球物理勘探近地表模型建立方法，以地震记录初至时间为基础，在初始速度模型的基础上，通过内外两层循环实现，内层循环通过射线追踪、走时残差和慢度模型的更新来实现，多次迭代之后进行速度模型的更新，完成一次外部迭代，即外层循环，这样反复迭代，直到满足给定的条件后，便得到最终的近地表速度模型。本发明在内层循环中，射线追踪后可以直接利用反投影算法(反投影算法是一种经典的层析反演算法，在层析理论方面的文献中有详细介绍)进行慢度更新量的求取，不需要保存射线路径，具有计算效率高、易于并行，对大数据和大模型的适应性强等特点，同时可以充分地利用高性能计算机的硬件资源，其计算效率得到显著提升。下面对该地震初至波走时层析方案进行详细介绍。
图1是本发明实施例中地震初至波走时层析方法的流程示意图，如图1所示，该方法包括如下步骤：
步骤101：获取地震初至波的初至时间；
步骤102：建立地震初至波走时层析反演的初始速度模型；
步骤103：根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
步骤104：在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
步骤105：执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
本发明实施例提供的技术方案通过：获取地震初至波的初至时间；建立地震初至波走时层析反演的初始速度模型；根据初至时间、初始速度模型，执行如下第一迭代步骤：按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；利用射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，第二迭代步骤包括第一迭代步骤以及更新速度模型的步骤；近地表速度模型用于地震初至波走时层析静校正量的计算，实现了：以地震波的初至时间为基础，射线追踪后直接进行慢度更新量的求取，不需要保存射线路径，其计算量仅与有效的炮检对数量和迭代次数有关，具有计算效率高、易于并行，对大数据和大模型的适应性强等特点，提高了初至波走时层析的计算效率。
本发明实施例提供的技术方案为适用于大数据量和大模型的初至波走时层析方案。具体实施时，近地表速度模型除了可以用于地震初至波走时层析静校正量的计算，还可以作为叠前深度偏移和全波形反演的近地表速度模型。
下面结合附图2至图4，对本发明实施例涉及的各个步骤进行介绍如下。
一、首先，介绍上述步骤101。
具体实施时，采集地震数据，进行初至拾取：采集地震数据，采用高精度的自动拾取或交互拾取初至的方法拾取初至时间，获取的初至时间作为初至波走时层析的基础数据之一。
二、其次，介绍上述步骤102。
具体实施时，可以按照简单梯度模型法或其他方法(如基于初至的层状梯度化模型法等)建立层析反演的初始速度模型，其中，简单梯度模型法的计算公式可以为：
Vi＝V0+G*(H0-Hi)；
其中，Vi为地下某一网格的速度，V0为给定的地表速度，如1500m/s，G为梯度因子，即为深度上每增加一个单位长度速度的变化量，如可取3s-1，H0为地表所在网格中心点的高程值，Hi为对应某一深度网格中心点的高程值(高程沿向下方向其值逐渐减小)。
三、再次，介绍上述步骤103。
具体实施时，上述步骤103的步骤即为上述提到的“内层循环”，内层循环通过射线追踪、走时残差和慢度模型的更新来实现，如图2所示，该步骤103可以包括①并行计算的步骤、②汇总的步骤、③慢度模型更新的步骤，下面对该步骤103进行详细介绍。
1、首先介绍并行计算的步骤：
在一个实施例中，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量，可以包括：
根据CPU核数，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，每一CPU负责一组的计算任务；每一CPU根据组内每一炮检对对应的初至时间，进行射线追踪、走时残差和慢度更新量的计算，该计算的步骤包括：
在当前速度模型下进行射线追踪，得到当前炮检对之间的射线穿过每个网格的射线路径长度；
根据所述射线路径长度，以及当前炮检对对应的地震道初至时间，计算当前慢度模型下的所述地震道的走时残差；
根据所述走时残差和射线路径长度，计算当前射线对每个网格的慢度更新量；
对每一任务组内的网格慢度模型的慢度更新量进行汇总，得到任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量；对每个网格穿过的射线条数进行汇总，得到任务组内对应的射线条数；对每个网格的走时残差的平方和信息进行汇总，得到任务组内对应的走时残差的平方和信息；多个任务组内的总贡献量为所述多个慢度更新量，多个任务组内对应的射线条数为所述多个射线密度，多个任务组内对应的走时残差的平方和信息为所述多个走时残差。
具体实施时，网格慢度是对于每一个网格来说的，相当于一个点的值，而慢度(速度的倒数)模型对于整体来说的，二维是一个面，三维是一个体。
具体实施时，将所有参与层析反演计算的有效初至，按照所使用的CPU核数进行任务分解，每个CPU核只负责对应任务内的数据，对这些数据进行射线追踪、计算走时残差和慢度更新量。具体包括：
首先，读取初至信息，包括总炮数、每炮的有效接收道数、炮检点的坐标和初至信息，具体地，可以利用总炮数、每炮的有效接收道数信息进行如下工作：1)进行分组使用；2)得到总的有效道数，相当于总的射线条数已知；炮检点坐标主要是计算炮检点在那个网格上，并计算炮检点之间的射线穿过那些网格使用；初至信息主要是后续中计算走时残差(射线追踪的走时减去初至时间)使用，按照CPU核数进行任务分组，分组的个数和CPU的核数一致，各个CPU核只负责对应任务组内的数据计算，例如CPU的核数为8，炮检对的数量为24，那么就分成8组，每一CPU负责3个炮检对对应数据的计算，具体包括以下几个方面(请参见图3)：
a)在当前速度模型下进行射线追踪，得到某一炮检对之间的射线穿过每个网格的射线路径长度li(例如图3中的l24)；
b)计算当前慢度模型下的走时残差：
Figure BDA0001982486900000081
其中，Δt为某一道(地震道)的走时残差，t0为对应道的初至时间，即为步骤101拾取的结果(初至时间)，li为步骤a)中射线追踪得到的射线穿过每个网格的射线路径长度，si为射线穿过每个网格所对应的慢度值(该信息已知，速度模型已知后，该信息通过速度模型求倒数即可得到)。
c)计算当前射线对每个网格慢度模型更新量Δsi的贡献(慢度更新量)，其计算公式为：
Figure BDA0001982486900000082
d)重复步骤a)、步骤b)和步骤c)完成当前任务组内的所有炮检对对应数据的计算，同时对每个网格的慢度更新量Δsi，每个网格穿过的射线条数，以及走时残差的平方和等信息进行汇总(累加求和)，得到任务组内所有数据对每个网格慢度更新量的总贡献，以及对应的射线条数，走时残差的平方和等信息(走时残差的平方和随后再除以总射线数，再求平方根作为迭代收敛的一个参数，相当于均方根误差，可以作为后续迭代终止的判断条件数据使用)。
2、其次，介绍汇总的步骤：
具体实施时，将慢度更新量、射线条数和走时残差等信息进行汇总，完成慢度模型的更新：每个CPU得到的每个网格的慢度更新量、射线条数和走时残差进行累加，然后将每个网格慢度更新量的总和除以穿过该网格的射线条数即可得到每个网格的慢度更新量，即为最终的慢度更新量，然后进行慢度模型的更新。
3、接着介绍更新慢度模型的步骤：
具体实施时，对于某一个网格的慢度模型更新公式为：
Snew＝Sold+ΔS；
其中，Snew为更新后的慢度模型，Sold为原始的慢度模型(速度模型的倒数即为慢度模型)，ΔS为慢度模型的更新量。
四、接着，介绍上述步骤104。
具体实施时，重复上述步骤103进行多次迭代，更新速度模型，并行实现射线追踪、慢度模型的更新，这样反复迭代多次，当满足给定的条件(该给定条件请参见下述实施例的介绍)后，终止迭代，然后再更新速度模型，即将最后生成的慢度模型取倒数后赋值给速度模型即可。
在一个实施例中，在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤，可以包括：
在第一预设的迭代次数大于第一预设值，或前后两次迭代的走时残差的均方根误差小于第二预设值时，根据更新后的慢度模型，执行更新速度模型的步骤。
具体实施时，走时残差的均方根误差可以根据上述提到的走时残差的平方和信息确定，上述第二预设值可以是1.0e-3。上述第一预设值可以是：默认值可以是8，用户可以根据实际需要，进行修改。
五、最后，介绍上述步骤105。
在一个实施例中，执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，可以包括：
执行第二迭代步骤，直至得到第二预设的迭代次数大于第三预设值(默认值一般选10，用户可以根据实际需要，进行修改)，或前后两次迭代的走时残差的均方根误差小于第四预设值时对应的近地表速度模型。
具体实施时，重复步骤103和104进行多次迭代，生成最终的近地表速度模型，在步骤104生成的速度模型的基础上，重复步骤103和104，这样进行多次迭代，直到满足给定的条件(该给定的条件请参见下述实施例的介绍)，如设定的迭代次数(内部循环的迭代次数与外部循环的迭代次数不一致)或前后两次迭代的均方根误差小于给定的值(该给定值与第一次迭代时的给定值可以不一样也可以一样，有一个默认值，用户也可以修改)后，终止迭代，得到的速度模型(图2中的最终速度模型)即为最终的近地表速度模型。
具体实施时，走时残差的均方根误差可以根据上述提到的走时残差的平方和信息确定，上述第三预设值可以是：10。
经过发明人反复的实验，如图4所示，与传统的串行(传统单线程算法)方案和分块计算(传统算法并行)的方案相比较，本发明实施例提供的技术方案的提高了初至波走时层析的计算效率。
基于同一发明构思，本发明实施例中还提供了一种地震初至波走时层析装置，如下面的实施例。由于地震初至波走时层析装置解决问题的原理与上述地震初至波走时层析方法相似，因此地震初至波走时层析装置的实施可以参考上述地震初至波走时层析方法的实施，重复之处不再赘述。以下所使用的，术语“模块”或者“单元”可以实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现，但是硬件，或者软件和硬件的组合的实现也是可能并被构想的。
图5是本发明实施例中地震初至波走时层析装置的结构示意图，如图5所示，该装置包括：
获取单元01，用于获取地震初至波的初至时间；
建立单元02，用于建立地震初至波走时层析反演的初始速度模型；
第一迭代单元03，用于根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
速度模型更新单元04，用于在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
第二迭代单元05，用于执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
在一个实施例中，所述第一迭代单元具体用于：
根据CPU核数，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，每一CPU负责一组的计算任务；每一CPU根据组内每一炮检对对应的初至时间，进行射线追踪、走时残差和慢度更新量的计算，该计算的步骤包括：
在当前速度模型下进行射线追踪，得到当前炮检对之间的射线穿过每个网格的射线路径长度；
根据所述射线路径长度，以及当前炮检对对应的地震道初至时间，计算当前慢度模型下的所述地震道的走时残差；
根据所述走时残差和射线路径长度，计算当前射线对每个网格的慢度更新量；
对每一任务组内的网格慢度模型的慢度更新量进行汇总，得到任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量；对每个网格穿过的射线条数进行汇总，得到任务组内对应的射线条数；对每个网格的走时残差的平方和信息进行汇总，得到任务组内对应的走时残差的平方和信息。
在一个实施例中，所述第一迭代单元具体用于：
将每个CPU得到的任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量进行累加，得到慢度更新量汇总值；
将每个CPU得到的任务组内对应的射线条数进行累加，得到射线密度汇总值；
将慢度更新量汇总值除以射线密度汇总值，得到每个网格的慢度更新量；
根据所述每个网格的慢度更新量对慢度模型进行更新，得到更新后的慢度模型；根据走时残差的平方和信息或迭代次数，确定是否继续更新慢度模型。
在一个实施例中，所述速度模型更新单元具体用于：在第一预设的迭代次数大于第一预设值，或前后两次迭代的走时残差的均方根误差小于第二预设值时，根据更新后的慢度模型，执行更新速度模型的步骤。
在一个实施例中，所述第二迭代单元具体用于：执行第二迭代步骤，直至得到第二预设的迭代次数大于第三预设值，或前后两次迭代的走时残差的均方根误差小于第四预设值时对应的近地表速度模型。
本发明实施例还提供了一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述地震初至波走时层析方法。
本发明实施例还提供了一种计算机可读存储介质，所述计算机可读存储介质存储有执行地震初至波走时层析方法的计算机程序。
本发明实施提供的技术方案的有益技术效果为：本发明目的是提供一种适用于海量数据的初至波走时层析方法和装置，解决常规走时层析无法有效适应性大数量、大模型以及计算效率低等问题，有效地利用目前多核多线程高性能计算机性能，减小编程复杂度，极大地提高了初至波走时层析的计算效率。
显然，本领域的技术人员应该明白，上述的本发明实施例的各模块或各步骤可以用通用的计算装置来实现，它们可以集中在单个的计算装置上，或者分布在多个计算装置所组成的网络上，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，并且在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明实施例不限制于任何特定的硬件和软件结合。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种地震初至波走时层析方法，其特征在于，包括：
获取地震初至波的初至时间；
建立地震初至波走时层析反演的初始速度模型；
根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
2.如权利要求1所述的地震初至波走时层析方法，其特征在于，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量，包括：
根据CPU核数，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，每一CPU负责一组的计算任务；每一CPU根据组内每一炮检对对应的初至时间，进行射线追踪、走时残差和慢度更新量的计算，该计算的步骤包括：
在当前速度模型下进行射线追踪，得到当前炮检对之间的射线穿过每个网格的射线路径长度；
根据所述射线路径长度，以及当前炮检对对应的地震道初至时间，计算当前慢度模型下的所述地震道的走时残差；
根据所述走时残差和射线路径长度，计算当前射线对每个网格的慢度更新量；
对每一任务组内的网格慢度模型的慢度更新量进行汇总，得到任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量；对每个网格穿过的射线条数进行汇总，得到任务组内对应的射线条数；对每个网格的走时残差的平方和信息进行汇总，得到任务组内对应的走时残差的平方和信息。
3.如权利要求2所述的地震初至波走时层析方法，其特征在于，将多个慢度更新量进行汇总，得到慢度更新量汇总值，包括：将每个CPU得到的任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量进行累加，得到慢度更新量汇总值；
将多个射线密度进行汇总，得到射线密度汇总值，包括：将每个CPU得到的任务组内对应的射线条数进行累加，得到射线密度汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型，包括：
将慢度更新量汇总值除以射线密度汇总值，得到每个网格的慢度更新量；
根据所述每个网格的慢度更新量对慢度模型进行更新，得到更新后的慢度模型；根据走时残差的平方和信息或迭代次数，确定是否继续更新慢度模型。
4.如权利要求1所述的地震初至波走时层析方法，其特征在于，在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤，包括：
在第一预设的迭代次数大于第一预设值，或前后两次迭代的走时残差的均方根误差小于第二预设值时，根据更新后的慢度模型，执行更新速度模型的步骤。
5.如权利要求1所述的地震初至波走时层析方法，其特征在于，执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，包括：
执行第二迭代步骤，直至得到第二预设的迭代次数大于第三预设值，或前后两次迭代的走时残差的均方根误差小于第四预设值时对应的近地表速度模型。
6.一种地震初至波走时层析装置，其特征在于，包括：
获取单元，用于获取地震初至波的初至时间；
建立单元，用于建立地震初至波走时层析反演的初始速度模型；
第一迭代单元，用于根据所述初至时间、初始速度模型，执行如下第一迭代步骤：
按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，并行进行射线追踪、走时残差和慢度更新量的计算，得到多个射线密度、多个走时残差和多个慢度更新量；
将多个射线密度进行汇总，得到射线密度汇总值；将多个走时残差进行汇总，得到走时残差汇总值；将多个慢度更新量进行汇总，得到慢度更新量汇总值；
利用所述射线密度汇总值、走时残差汇总值和慢度更新量汇总值，对慢度模型进行更新，得到更新后的慢度模型；
速度模型更新单元，用于在满足第一预设条件时，根据更新后的慢度模型，执行更新速度模型的步骤；
第二迭代单元，用于执行第二迭代步骤，直至得到满足第二预设条件的近地表速度模型，所述第二迭代步骤包括所述第一迭代步骤以及更新速度模型的步骤；所述近地表速度模型用于地震初至波走时层析静校正量的计算。
7.如权利要求6所述的地震初至波走时层析装置，其特征在于，所述第一迭代单元具体用于：
根据CPU核数，按炮检对，将射线追踪、走时残差和慢度更新量的计算任务进行分组，每一CPU负责一组的计算任务；每一CPU根据组内每一炮检对对应的初至时间，进行射线追踪、走时残差和慢度更新量的计算，该计算的步骤包括：
在当前速度模型下进行射线追踪，得到当前炮检对之间的射线穿过每个网格的射线路径长度；
根据所述射线路径长度，以及当前炮检对对应的地震道初至时间，计算当前慢度模型下的所述地震道的走时残差；
根据所述走时残差和射线路径长度，计算当前射线对每个网格慢度模型的慢度更新量；
对每一任务组内的网格慢度模型的慢度更新量进行汇总，得到任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量；对每个网格穿过的射线条数进行汇总，得到任务组内对应的射线条数；对每个网格的走时残差的平方和信息进行汇总，得到任务组内对应的走时残差的平方和信息。
8.如权利要求7所述的地震初至波走时层析装置，其特征在于，所述第一迭代单元具体用于：
将每个CPU得到的任务组内所有数据对每个网格慢度模型慢度更新量的总贡献量进行累加，得到慢度更新量汇总值；
将每个CPU得到的任务组内对应的射线条数进行累加，得到射线密度汇总值；
将慢度更新量汇总值除以射线密度汇总值，得到每个网格的慢度更新量；
根据所述每个网格的慢度更新量对慢度模型进行更新，得到更新后的慢度模型；根据走时残差的平方和信息或迭代次数，确定是否继续更新慢度模型。
9.一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至5任一所述方法。
10.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有执行权利要求1至5任一所述方法的计算机程序。
CN201910154540.8A 2019-03-01 2019-03-01 地震初至波走时层析方法及装置 Pending CN111638551A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910154540.8A CN111638551A (zh) 2019-03-01 2019-03-01 地震初至波走时层析方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910154540.8A CN111638551A (zh) 2019-03-01 2019-03-01 地震初至波走时层析方法及装置

Publications (1)

Publication Number Publication Date
CN111638551A true CN111638551A (zh) 2020-09-08

Family

ID=72327022

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910154540.8A Pending CN111638551A (zh) 2019-03-01 2019-03-01 地震初至波走时层析方法及装置

Country Status (1)

Country Link
CN (1) CN111638551A (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112630833A (zh) * 2020-11-19 2021-04-09 安徽理工大学 一种基于测井曲线的快速地震初至波走时联合反演方法
CN114428292A (zh) * 2020-09-22 2022-05-03 中国石油化工股份有限公司 近地表速度模型的构建方法和存储介质
CN114624768A (zh) * 2020-12-14 2022-06-14 中国石油化工股份有限公司 训练地震初至拾取模型的方法及装置
CN117724166A (zh) * 2024-02-07 2024-03-19 中国石油大学(华东) 基于大炮初至的近地表三维速度建模方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103135132A (zh) * 2013-01-15 2013-06-05 中国科学院地质与地球物理研究所 Cpu/gpu协同并行计算的混合域全波形反演方法
CN104133240A (zh) * 2014-07-29 2014-11-05 中国石油天然气集团公司 大规模并行的克希霍夫叠前深度偏移方法及装置
CN105277978A (zh) * 2015-10-15 2016-01-27 中国石油天然气集团公司 一种确定近地表速度模型的方法及装置
CN107783185A (zh) * 2017-09-14 2018-03-09 中国石油天然气股份有限公司 一种层析静校正的处理方法及装置
CN107843922A (zh) * 2017-12-25 2018-03-27 中国海洋大学 一种基于地震初至波和反射波走时联合的层析成像方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103135132A (zh) * 2013-01-15 2013-06-05 中国科学院地质与地球物理研究所 Cpu/gpu协同并行计算的混合域全波形反演方法
CN104133240A (zh) * 2014-07-29 2014-11-05 中国石油天然气集团公司 大规模并行的克希霍夫叠前深度偏移方法及装置
CN105277978A (zh) * 2015-10-15 2016-01-27 中国石油天然气集团公司 一种确定近地表速度模型的方法及装置
CN107783185A (zh) * 2017-09-14 2018-03-09 中国石油天然气股份有限公司 一种层析静校正的处理方法及装置
CN107843922A (zh) * 2017-12-25 2018-03-27 中国海洋大学 一种基于地震初至波和反射波走时联合的层析成像方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
潘英杰 等: "一种针对物探应用的多机并行计算框架", 《地球物理学进展》 *
陈楠: "基于D－Bus技术的地震作业管理架构设计", 《计算机技术与发展》 *
马青坡 等: "基于近地表模型约束的初至波走时层析方法", 《CPS/SEG北京2018国际地球物理会议暨展览电子论文集》 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114428292A (zh) * 2020-09-22 2022-05-03 中国石油化工股份有限公司 近地表速度模型的构建方法和存储介质
CN112630833A (zh) * 2020-11-19 2021-04-09 安徽理工大学 一种基于测井曲线的快速地震初至波走时联合反演方法
CN114624768A (zh) * 2020-12-14 2022-06-14 中国石油化工股份有限公司 训练地震初至拾取模型的方法及装置
CN117724166A (zh) * 2024-02-07 2024-03-19 中国石油大学(华东) 基于大炮初至的近地表三维速度建模方法

Similar Documents

Publication Publication Date Title
CN111638551A (zh) 地震初至波走时层析方法及装置
CN106772583B (zh) 一种地震绕射波分离方法和装置
CN105277978B (zh) 一种确定近地表速度模型的方法及装置
CN109917454B (zh) 基于双基准面的真地表叠前深度偏移成像方法及装置
CN106526674A (zh) 一种三维全波形反演能量加权梯度预处理方法
CN105005072B (zh) 利用cuda的pml边界三维地震波传播模拟方法
CN105093319B (zh) 基于三维地震数据的地面微地震静校正方法
CN107843922A (zh) 一种基于地震初至波和反射波走时联合的层析成像方法
CN116774292B (zh) 一种地震波走时确定方法、系统、电子设备及存储介质
CN114839673B (zh) 多震源高效采集波场分离方法、分离系统及计算机设备
CN102692644A (zh) 生成深度域成像道集的方法
CN108828668A (zh) 一种叠前时间偏移数据处理方法及装置
CN107817516A (zh) 基于初至波信息的近地表建模方法及系统
CN108680968B (zh) 复杂构造区地震勘探数据采集观测系统评价方法及装置
CN105204063B (zh) 地震数据速度模型建立方法和装置
CN106353798A (zh) 多分量联合高斯束叠前逆时偏移成像方法
CN105301638B (zh) 一种提取风化层底界面的方法和装置
CN110187386B (zh) 一种自动快速识别地质构造的dtw地震体属性分析方法
CN109738944B (zh) 基于广角反射的地震采集参数确定方法及装置
CN108983291B (zh) 一种近地表散射波获取方法、装置及系统
CN110161565A (zh) 一种地震数据重建方法
CN105572730B (zh) 三维复杂结构声波正演方法
CN108072895A (zh) 一种基于gpu的各向异性叠前逆时偏移优化方法
CN114428342A (zh) 基于MPI的Marchenko格林函数重构方法及系统
CN113970789A (zh) 全波形反演方法、装置、存储介质及电子设备

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
RJ01 Rejection of invention patent application after publication

Application publication date: 20200908