一种利用gpu进行逆时偏移提取角度道集的方法
【专利摘要】本发明提供了一种利用GPU进行逆时偏移角度道集提取的方法，属于石油工业的地震勘探与开发领域。本发明方法包括以下步骤：(1)炮点波场正向延拓；(2)并行执行逆时偏移和角度信息统计；(3)角度道集提取。本发明方法提供了一种快速稳定的逆时偏移角度道集提取方法，其具逆时偏移和角度信息统计并行实施的能力，更适合于GPU程序的实现，解决了以往逆时偏移和角度道集提取必须串行执行的问题，解决了计算效率低，存储需求过大，容错能力弱，不易反复试验的缺点。
【专利说明】一种利用GPU进行逆时偏移提取角度道集的方法
【技术领域】
[0001]本发明属于石油工业的地震勘探与开发领域，具体涉及一种利用GPU进行逆时偏移提取角度道集的方法，用于人工地震记录的偏移成像软件系统中的共像点道集提取。
【背景技术】
[0002]GPU的英文全称是Graphic Processing Unit,中文翻译为“图形处理器”,它是从CPU(中央处理器)衍生出来的的一个概念，它是显卡的心脏，和显卡板载内存构成一个子系统，共同决定PC系统的图形处理性能。越来越多的图形应用使得GPU在现代的计算机中变得越来越重要，这种芯片在CPU主宰PC性能走向几十年后的近两年中迅速崛起，其浮点计算能力已经达到同期CPU的十倍甚至几十倍，在高性能计算领域，差不多达到了和CPU平起平坐的地位。
[0003]将图形处理器(GPU)用于解决通用计算(GPGPU)及其相关方面的问题目前已成为一个热门话题。事实上，整个IT产业都已经敏锐地意识到了 GPU通用计算将给PC带来革命性的变化，进而影响到CPU的发展。GPU是专为图像处理设计的，其存储系统实际上是一个二维的分段存储空间，包括一个区段号(从中读取图像)和二维地址(图像中的X、Y坐标)。此外，没有任何间接写指令。输出写地址由光栅处理器确定，而且不能由程序改变。这对于自然分布在存储器之中的算法而言是极大的挑战，存储器的访问效率是提升算法效率的关键问题。以往利用GPU做高性能计算是非常困难的工作，2007年以后CUDA和OpenCL等高级编程语言的问世使编程新手也能够掌握GPU的计算性能优势。
[0004]逆时偏移的应用中一个基础的研究工作就是角度道集的提取。进入21世纪以来，波动方程输出共反射角度像集就成为业界研究的热点。其实现方法可以划分为两大类:一类是成像前的提取方法(如图1所示)，后文称为方法一，该方法需要在逆时偏移的过程中实时计算波场角度，并且实时输出角度道集，三者不可分割实现(如图5所示，图5是图1的细化流程，其中角度道集提取和波场延拓耦合在一个模块中执行)，并且需要大量(TB级)的内存空间作为保障，是目前计算机技术难以满足的条件；另一类是成像后的提取方法(如图2所示)，后文称为方法二，该方法采用平移或时移成效条件首先生成偏移距共像点道集，然后通过Tau-P变换等方法转换生成角度道集。成像后的提取方法适合于以往的偏移模块，只需要做较少的改动和后续的转换处理就可以得到角度道集，也是目前使用较多的方法。但是，由于通常偏移距与角度的转换关系并不能做到一一对应，所以生成的角度道集同样也存在假象。
[0005]如今的地震资料采集已经发展到了三维甚至四维，其数据量更是十分庞大，已经达到TB级甚至IOTB级的叠前单炮数据，如此庞大的数据的计算时间和IO时间也是逆时偏移角度道集提取工作必须面临的难题。
[0006]由于GPU高性能计算的飞速发展，可以大幅度提高逆时偏移和角度道集计算的问题，但是在现有技术的实施过程中仍然没有解决存储问题。
【发明内容】
[0007]本发明的目的在于解决上述现有技术中存在的难题，提供一种利用GPU进行逆时偏移角度道集提取的方法，能够快速稳定地提取角度道集。
[0008]本发明是通过以下技术方案实现的:
[0009]一种利用GPU进行逆时偏移角度道集提取的方法，所述方法包括以下步骤:
[0010](I)炮点波场正向延拓；
[0011](2)并行执行逆时偏移和角度信息统计；
[0012](3)角度道集提取。
[0013]所述步骤(I)包括以下步骤:
[0014](Al)在炮点位置放入人工子波震源；
[0015](A2)利用有限差分方法对步骤(Al)中的人工子波震源激发出的地震波的波场沿时间方向进行正向延拓，得到不同时刻的炮点波场；
[0016]所述步骤(A2)是在GPU上完成的。
[0017]所述步骤(2)中，所述逆时偏移包括以下步骤:
[0018](BI)将步骤(A2)中得到的不同时刻的炮点波场存储在波场数据缓冲池中；
[0019](B2)将从步骤(A2)得到的炮点波场沿时间方向正向延拓至最大时刻，得到炮点最大时刻的波场；正向延拓是指沿时间递增的方向延拓波场；
[0020](B3)从野外采集的地震数据中读取检波点最大时刻的波场，然后从最大时刻开始，将炮点波场和检波点波场沿时间方向反向延拓，得到炮点和检波点不同时刻的波场；反向延拓是指沿时间递减的方向延拓波场；
[0021](B4)将每个时刻的炮点波场与检波点波场相乘得到成像结果；
[0022](B5)每个炮点对应一个单炮，将每个单炮对应的步骤(B4)中的所述成像结果作为单炮偏移成像结果单独存储在磁盘上；
[0023](B6)每个从野外采集的地震数据包含很多单炮，利用每个单炮对应的步骤(B4)中的所述成像结果做多炮叠加得到多炮叠加的偏移成像结果；所述多炮叠加是指把多个单炮成像结果针对地下相同位置处的成像值进行累加。
[0024]所述步骤(2)中，所述角度信息统计包括以下步骤:
[0025](Cl)角度信息统计:将所述波场数据缓冲池中的炮点波场按照时间顺序存储，并且提取每个时刻每个空间位置的振幅和传播角度；
[0026](C2)完成步骤(A2)后，将步骤(Cl)中得到的炮点波场在每个空间位置最大振幅时刻的传播角度作为角度信息数据体存储存储在磁盘上；
[0027](C3)在所述步骤(B6)得到的所述多炮叠加的偏移成像结果上，提取地层倾角作为地层倾角信息数据体存储在磁盘上。
[0028]其中，所述步骤(B2)至(B5)与步骤(Cl)至(C2)并行执行，步骤(B6)执行完毕后再执行步骤(C3)。
[0029]所述步骤(3)包括以下步骤:
[0030](Dl)从磁盘上读取步骤(C2)和(C3)中得到的角度信息数据体和地层倾角信息数据体；
[0031](D2)根据所述角度信息数据体和地层倾角信息数据体计算入射角:[0032]入射角Θ等于炮点波场在每一空间位置处达到最大振幅时刻的传播角度β与地层倾角α —半的差，公式表达为:θ = β -0.5 α ；[0033](D3)从磁盘上读取步骤(Β5)输出的单炮偏移成像结果，然后将单炮偏移成像结果输出到对应于步骤(D2)得到的入射角的道集中，得到最终的角度道集数据体。[0034]所述数据缓冲池中的基本数据结构包括一个浮点型数据、一个整型数据和一个浮点型数组，其中浮点型数据存储当前时刻值、整型数据存储所述浮点数组的长度、浮点型数组存储当前时刻的波场。[0035]与现有技术相比，本发明的有益效果是:[0036](I)本发明中利用数据缓冲池技术，实现了逆时偏移、角度信息统计的并行化实施方案。避免以往二者只能串行实现的方法，提闻计算效率。[0037](2)由于炮点波场的角度信息是在成像前统计生成的，保持了成像前的角度计算精度条件下，实现逆时偏移完成后再进行角度道集的提取工作。这样可以减少方法一中逆时偏移和角度道集提取同时进行的内存开销，同时可以消除方法二中成像后角度道集提取的成像假象。[0038](3)本专利由于只统计波场最大能量时刻的角度信息，所以在不同时刻的角度计算结果可存放在相同的内存或磁盘空间，实现存储空间的重复利用，这样可以大幅度的减少内存和磁盘的需求量。
【专利附图】
【附图说明】[0039]图1是现有技术中的方法一的步骤框图。
[0040]图2是现有技术中的方法二的步骤框图。
[0041]图3是本发明利用GPU逆时偏移提取角度道集的方法中设计的模块结构图。
[0042]图4是本发明利用GPU提取逆时偏移角度道集方法的详细步骤框图。
[0043]图5是采用方法一提取角度道集的详细步骤框图。
[0044]图6是本发明根据炮点波场传播角度和地层倾角计算入射角的示意图。
[0045]图7是本发明中一个单炮对应地下成像点一个入射角的示意图。
[0046]图8 是 sigseeb2a 模型图。
[0047]图9-1是利用本发明方法对图8的模型图进行角度道集提取得到的效果图。
[0048]图9-2是利用方法一对图8的模型图进行角度道集提取得到的效果图。
[0049]图9-3是利用方法二对图8的模型图进行角度道集提取得到的效果图。
【具体实施方式】
[0050]下面结合附图对本发明作进一步详细描述:
[0051]本发明针对GPU逆时偏移提取角度道集的方法包括三个模块:逆时偏移模块(图4中的(一))、角度信息统计模块(图4中的(二))和角度道集提取模块(图4中的(三))，模块的细分如图4，每个步骤都可以用CPU和GPU来实现，但是在GPU上效率会大大提高，只有利用GPU才能满足工业应用的需要。图3是CPU/GPU协同计算实现逆时偏移角度道集提取的模块图。其中逆时偏移模块和角度信息统计模块可以并行实施，为角度道集提取模块提供必要的输入数据；逆时偏移模块中的炮点波场正向延拓步骤需要为数据缓冲池提供波场数据，角度信息统计模块直接从数据缓冲池中读取波场数据。
[0052]图4是模块并行处理的详细流程图(图4中的(一)、(二)、(三)是表示具体步骤分属于哪个模块)，其中白色背景表示逆时偏移流程(属于现有技术)，深色背景表示角度信息统计模块的实施流程(包括波场数据缓冲池、角度信息统计、最大振幅时刻的传播角度计算与输出和成像剖面的地层倾角计算与输出)，这两个流程完成后利用逆时偏移的单炮成像结果、空间位置的最大振幅角度信息数据体和成像剖面的地层倾角信息提取角度道集。具体流程如下:
[0053](I)炮点波场正向延拓
[0054](2)并行执行逆时偏移和角度信息统计
[0055](3)角度道集提取。
[0056]所述步骤⑴中，所述炮点波场正向延拓包括以下步骤:
[0057](Al)在炮点位置放入人工子波震源；
[0058](A2)利用有限差分方法正向延拓，可得到不同时刻的炮点波场。此为现有技术，必须在GPU上执行，因为GPU上执行才能达到满意的计算效率，符合工业应用要求；所述有限差分方法是基本的波场延拓的算法之一；
[0059]所述步骤(2)中，所述逆时偏移包括以下步骤:
[0060](BI)将(A2)步骤中得到的不同时刻的炮点波场存储在波场数据缓冲池中；角度信息统计模块中，为了可以保证与逆时偏移模块的并行化实施，本发明采用了波场数据缓冲池的中间环节，此波场数据缓冲池主要是将炮点波场正演过程中的不同时刻的波场信息暂时存储在内存或者磁盘上，用于角度信息的统计工作。其中并行执行应包含GPU上的波场延拓并行，和CPU上或者另一个GPU上的角度信息统计模块的并行化执行。这样可以在不影响原有逆时偏移效率的基础上实现角度道集的提取工作；
[0061](B2)将从步骤(A2)得到的炮点波场沿时间方向正向延拓(又称为正演)至最大时刻，得到炮点最大时刻的波场；
[0062](B3)从野外采集的地震数据中读取检波点最大时刻的波场，然后从最大时刻开始，将炮点波场和检波点波场反向延拓，得到炮点和检波点不同时刻的波场(因为逆时偏移的成像条件是沿着时间反向积分的，所以需要反向传播过程中成像，若是在正向延拓过程中把每个时刻的波场存储下来，需要很多的磁盘空间，三维工区时可以达到几十基至几百TB的需求，反向计算一次可以避免这个存储，此为已有的技术)；反向延拓就是指沿时间递减的方向延拓波场；
[0063](B4)将每个时刻的炮点波场与检波点波场相乘得到成像结果(即图4中的应用相关成像条件)；波场是空间和时间的函数其物理含义是振幅。
[0064](B5)每个炮点对应一个单炮，将每个单炮对应的步骤(B4)中的所述成像结果作为单炮偏移成像结果单独存储在磁盘上；
[0065](B6)每个从野外采集的地震数据包含很多单炮，利用每个单炮的所述成像结果做多炮叠加得到多炮叠加的偏移成像结果；所述多炮叠加是指把多个单炮成像结果，针对地下相同位置的成像值进行累加。
[0066]所述步骤(2)中，所述角度信息统计包括以下步骤:
[0067](Cl)角度信息统计:将所述波场数据缓冲池中的波场按照时间顺序存储，角度信息提取模块直接访问波场数据缓冲池中的波场，并且提取每个时刻每个空间位置的振幅和传播角度(提取方法可以是Poynting矢量法或者局部波分解法，此为已有技术)；当波场数据缓冲池中波场时刻的个数大于I个时，首先比较不同时刻的波场振幅然后提取最大振幅时刻的传播角度；统计每个空间位置达到最大振幅时刻的传播角度(计算中时间是离散化的，以等间隔的时刻序列来体现，延拓过程中会把波场暂存在缓冲池中，用以统计角度信息，每个时刻的角度信息统计完后会被下一个时刻的波场覆盖，若角度信息统计的速度比延拓速度慢时，缓冲池中会同时存在多个时刻的波场，此时需要比较现存在于缓冲池中不同时刻的波场振幅以选择其中振幅最大的波场来计算角度信息，统计完成后缓冲池清空，等待后面时刻的波场进入缓冲池，再进行相同的统计。)。波场是空间和时间的函数，统计的结果也是空间的函数。
[0068]角度信息统计的实现方法也是本发明的重要保护内容，波场数据缓冲池的基本数据结构包括一个浮点型数据、一个整型数据和一个浮点型数组，其中浮点数据存储当前时刻值、浮点数组存储当前时刻的波场(波场正向延拓指的就是沿时间递增的方向延拓波场，其结果即为不同时刻的波场，那么在波场数据缓冲池中包含的波场是为了计算角度而保留的，每个时刻计算一次即可被下一时刻覆盖掉，所以需要指明所计算的角度是那个时刻的，即所述的当前时刻。)、整型数据存储浮点数组的长度。角度信息统计模块直接从波场数据缓冲池中将不同时刻的波场值读出，然后计算每个时刻每个空间位置的振幅和传播角度。需要指出的是，当波场数据缓冲池中波场时刻较多时首先比较振幅然后提取最大振幅的角度信息(缓冲池中可能同时存在多个时刻的波场，此时需要比较振幅)，用以后续的角度道集提取工作。这样可以保证与逆时偏移中的波场延拓的同步性，而且不需要任何通信，只需访问波场数据缓冲池中的时间和波场即可。
[0069](C2)完成步骤(A2)后(即全部时刻延拓完毕后)，将步骤(Cl)中得到的炮点波场在每个空间位置最大振幅时刻的传播角度作为角度信息数据体存储(因为在Cl执行过程中不能保证所有时刻都已包含，必须全部时刻延拓完毕才能得到该角度信息数据体。)存储在磁盘上；所述角度信息数据体作为后续角度道集提取的输入数据的一部分；
[0070](C3)在所述步骤(B6)得到的多炮叠加的偏移成像结果上，提取地层倾角作为地层倾角信息数据体存储在磁盘上；所述地层倾角信息数据体作为后续角度道集提取的输入数据的一部分；
[0071]其中，步骤(B2)至(B5)可以与(Cl)至(C2)并行执行，(B6)执行完毕后再执行C3 ；
[0072]所述步骤(3)中，所述角度道集提取包括以下步骤:
[0073](Dl)从磁盘上读取步骤(C2)和(C3)中得到的角度信息数据体和地层倾角信息数据体；
[0074](D2)根据所述角度信息数据体和地层倾角信息数据体计算入射角:本发明的角度信息统计模块最终保存下来的是炮点波场在每一空间位置处达到最大振幅时刻的波场传播角度和逆时偏移多炮叠加成像后的地层倾角。利用这两个角度信息，推算地震波的入射角信息来抽取成像结果中的成像值生成角度道集。具体推算方法如图6所示(图6中S表示炮点、R表示检波点)，入射角Θ等于炮点波场在每一空间位置处达到最大振幅时刻的波场传播角度β与地层倾角(又称为反射同相轴倾角)α —半的差，公式表达为:Θ =β _0.5 ct 0
[0075](D3)从磁盘上读取步骤(Β5)输出的单炮偏移成像结果，然后将单炮偏移成像结果输出到对应于步骤(D2)得到的入射角的道集中，得到最终的角度道集数据体。图7描述了不同单炮对应不同的入射角是自然形成的。
[0076]在角度道集抽取模块中，每个成像点需要不同入射角的信息，本发明把不同入射角映射为不同的单炮位置，如图7所示(Xs表示炮点位置坐标)，一个单炮数据的成像结果对应一个入射角，在实施过程中分别保存单炮数据的成像结果(即步骤(Β5)中的单炮偏移结果单独存储的结果)，做为后续角度道集提取模块的输入数据。
[0077]所述步骤(Α2)是在GPU上完成的，其它步骤在GPU或CPU上完成。本发明主要强调此流程的划分方法有利于并行和提高效率，(Α2)步骤的有限差分计算量很大需要使用GPU，但是角度提取针对不同方法可采取GPU或CPU上完成。
[0078]下面以国际SEG标准模型sigseeb2a(如图8所示)作为本发明的实例数据对本发明的效果进行说明，该模型的地震工区信息如表1所示。
[0079]
【权利要求】
1.一种利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述方法包括以下步骤: (1)炮点波场正向延拓； (2)并行执行逆时偏移和角度信息统计； (3)角度道集提取。
2.根据权利要求1所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述步骤(I)包括以下步骤: (Al)在炮点位置放入人工子波震源； (A2)利用有限差分方法对步骤(Al)中的人工子波震源激发出的地震波的波场沿时间方向进行正向延拓，得到不同时刻的炮点波场； 所述步骤(A2)是在GPU上完成的。
3.根据权利要求2所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述步骤(2)中，所述逆时偏移包括以下步骤: (BI)将步骤(A2)中得到的不同时刻的炮点波场存储在波场数据缓冲池中； (B2)将从步骤(A2)得到的炮点波场沿时间方向正向延拓至最大时刻，得到炮点最大时刻的波场；正向延拓是指沿时间递增的方向延拓波场； (B3)从野外采集的地震数据中读取检波点最大时刻的波场，然后从最大时刻开始，将炮点波场和检波点波场沿时间方 向反向延拓，得到炮点和检波点不同时刻的波场；反向延拓是指沿时间递减的方向延拓波场； (B4)将每个时刻的炮点波场与检波点波场相乘得到成像结果； (B5)每个炮点对应一个单炮，将每个单炮对应的步骤(B4)中的所述成像结果作为单炮偏移成像结果单独存储在磁盘上； (B6)每个从野外采集的地震数据包含很多单炮，利用每个单炮对应的步骤(B4)中的所述成像结果做多炮叠加得到多炮叠加的偏移成像结果；所述多炮叠加是指把多个单炮成像结果针对地下相同位置处的成像值进行累加。
4.根据权利要求3所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述步骤(2)中，所述角度信息统计包括以下步骤: (Cl)角度信息统计:将所述波场数据缓冲池中的炮点波场按照时间顺序存储，并且提取每个时刻每个空间位置的振幅和传播角度； (C2)完成步骤(A2)后，将步骤(Cl)中得到的炮点波场在每个空间位置最大振幅时刻的传播角度作为角度信息数据体存储存储在磁盘上； (C3)在所述步骤(B6)得到的所述多炮叠加的偏移成像结果上，提取地层倾角作为地层倾角信息数据体存储在磁盘上。
5.根据权利要求4所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述步骤(B2)至(B5)与步骤(Cl)至(C2)并行执行，步骤(B6)执行完毕后再执行步骤(C3)。
6.根据权利要求5所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述步骤(3)包括: (Dl)从磁盘上读取步骤(C2)和(C3)中得到的角度信息数据体和地层倾角信息数据体； (D2)根据所述角度信息数据体和地层倾角信息数据体计算入射角: 入射角Θ等于炮点波场在每一空间位置处达到最大振幅时刻的传播角度β与地层倾角α —半的差，公式表达为:θ = β -0.5 α ； (D3)从磁盘上读取步骤(Β5)输出的单炮偏移成像结果，然后将单炮偏移成像结果输出到对应于步骤(D2)得到的入射角的道集中，得到最终的角度道集数据体。
7.根据权利要求3所述的利用GPU进行逆时偏移角度道集提取的方法，其特征在于:所述数据缓冲池中的基本数据结构包括一个浮点型数据、一个整型数据和一个浮点型数组，其中浮点型数据存储当前时刻值、整型数据存储所述浮点数组的长度、浮点型数组存储当前时刻 的波场。
【文档编号】G01V1/32GK103472481SQ201210184592
【公开日】2013年12月25日 申请日期:2012年6月6日 优先权日:2012年6月6日
【发明者】李敏, 李振华, 任鹏 申请人:中国石油化工股份有限公司, 中国石油化工股份有限公司石油物探技术研究院