叠前深度偏移方法和装置与流程

本发明涉及地震勘探技术领域，尤其是涉及一种叠前深度偏移方法和装置。

背景技术

地震波在传播过程中，一部分能量转化为热能，相应的地震波的幅值会产生衰减效应，这个过程称为吸收。幅值的衰减对地震波的不同频率成份是不同的，频率越高，衰减就越强，这是因为高频成份波长较低频成份波长短，对于一个固定的传播距离来讲，相当于低频成份较少个波长，但相当于高频成份很多个波长，而地震波每传播一个波长的距离，能量损失的程度是固定的。这也导致接收到的反射地震资料的有效频带随反射深度逐渐变窄；而不同频率成份以不同的速度传播，也导致了地震子波的频散，这一频散现象也是反射构造越深，频散越严重。

由于常规叠前深度偏移方法没有补偿粘性吸收导致的幅值衰减，也没有校正频散，因而偏移成像结果的分辨率较低，其特点是，构造越深，分辨率就越低。

技术实现要素：

有鉴于此，本发明的目的在于提供一种叠前深度偏移方法和装置，可以提高地震反射数据对地下构造的分辨能力。

第一方面，本发明实施例提供了一种叠前深度偏移方法，包括：获取目标工区的叠前地震数据；剔除该叠前地震数据中的奇异道得到待处理地震数据；根据该待处理地震数据计算目标工区各激发点的慢度表和等效q值表；根据该慢度表和该等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合；读取预先设定的每个成像点的偏移孔径，并根据该走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值。

结合第一方面，本发明实施例提供了第一方面的第一种可能的实施方式，其中，上述根据该待处理地震数据计算慢度表和等效q值表的步骤包括：根据该待处理地震数据，利用快速波前法计算慢度表和等效q值表。

结合第一方面，本发明实施例提供了第一方面的第二种可能的实施方式，其中，上述根据该走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值的步骤，包括：计算每一个成像点的深度域倾角，并读取该成像点的偏移孔径；识别该深度域倾角是否在该偏移孔径的边界范围内；若是，则计算该成像点的偏移幅值。

结合第一方面的第二种可能的实施方式，本发明实施例提供了第一方面的第三种可能的实施方式，其中，上述计算该成像点的偏移幅值的步骤，包括：分别计算该成像点在预先设定的低、高截止频率内的频率成份对应的第一偏移幅值，以及该成像点在高频端衰减带中的频率成份对应的第二偏移幅值；叠加第一偏移幅值和第二偏移幅值得到该成像点的偏移幅值。

结合第一方面的第三种可能的实施方式，本发明实施例提供了第一方面的第四种可能的实施方式，其中，上述分别计算该成像点在预先设定的低、高截止频率内的频率成份对应的第一偏移幅值，以及该成像点在高频端衰减带中的频率成份对应的第二偏移幅值的公式，包括：

式中，i(h,x,y,z)为第一偏移幅值，im(h,x,y,z)为第二偏移幅值，φ(iδη)为补偿因子光滑性阈值控制函数，m为衰减带厚度，取值范围10～15；g为振幅补偿阈值，取值范围800～1800；h为偏移距，x,y,z为空间位置坐标，表示全部成像点的偏移孔径约束，j为虚部单位，τs为炮点到成像点的走时，τg为检波点到成像点的走时，为炮点走时q比，为检波点走时q比，δω为频率间隔，ω0为地震数据主频，f表示地震数据变换至频率域，表示取实部，l1、l2分别为低、高截止频率对应的频率采样点编号。

结合第一方面，本发明实施例提供了第一方面的第五种可能的实施方式，其中，在上述剔除叠前地震数据中的奇异道得到待处理地震数据的步骤之前，还包括：调用叠前地震数据分配算法将该叠前地震数据按负载平衡规则存放到不同计算节点；该叠前地震数据分配算法为对该叠前地震数据进行虚拟分配和实际分配相结合的算法；该负载平衡规则为使该计算节点上gpu卡所分配的不同偏移距组地震数据包含的地震道数目尽可能均衡。

结合第一方面，本发明实施例提供了第一方面的第六种可能的实施方式，其中，上述实际分配的步骤，包括：遍历该叠前地震数据；根据该叠前地震数据的地震道的偏移距确定计算单元；根据该计算单元确定gpu卡；根据该计算单元和该gpu卡将该地震道发送至相应的计算机节点。

结合第一方面，本发明实施例提供了第一方面的第七种可能的实施方式，其中，上述根据该慢度表和该等效q值表计算走时表和走时q比表的步骤包括：选取坐标为yi的成像线，针对所有激发点，结合预先设定的偏移孔径，从该慢度表和该等效q值表中抽取与该成像线临近的慢度面ej(x,y1,z)和ej(x,y2,z)以及等效q值面和y1≤yi≤y2；调用一维插值算法计算该偏移孔径内所有激发点对该成像线有贡献的慢度面集合ej(x,yi,z)与等效q值面集合j＝1,2,…j，j为成像线y方向偏移孔径内所有激发点的个数；逐一选取各地震道，分别读取该地震道的炮点坐标和检波点坐标，从该慢度面集合ej(x,yi,z)与等效q值面集合中选取4个与该炮点临近的激发点的慢度面和等效q值面，并从上述慢度面集合ej(x,yi,z)与等效q值面集合中选取4个与该检波点临近的激发点的慢度面和等效q值面；调用双线性插值算法计算该炮点和该检波点在该成像线处和该成像线的成像点网格相同的慢度面ek(xi,yi,zi)与等效q值面k为s或g，s代表炮点，g代表检波点；根据该慢度面ek(xi,yi,zi)和该等效q值面计算该地震道在该成像线的走时表和走时q比表。

结合第一方面的第七种可能的实施方式，本发明实施例提供了第一方面的第八种可能的实施方式，其中，上述根据该慢度面ek(xi,yi,zi)和该等效q值面计算该地震道在该成像线的走时表和走时q比表的计算公式为：τk(xi,yi,zi)＝ek(xi,yi,zi)/rk(xi,yi,zi)，式中，k为s或g，s代表炮点，g代表检波点，τk(xi,yi,zi)为炮点或检波点沿射线路径到成像点的走时，rk(xi,yi,zi)为炮点或检波点到成像点的距离。

第二方面，本发明实施例还提供了一种叠前深度偏移装置，包括：数据获取模块，用于获取目标工区的叠前地震数据；剔除模块，用于剔除该叠前地震数据中的奇异道得到待处理地震数据；第一计算模块，用于根据该待处理地震数据计算目标工区各激发点的慢度表和等效q值表；第二计算模块，用于根据该慢度表和该等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合；偏移幅值计算模块，用于读取预先设定的每个成像点的偏移孔径，并根据该走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值。

本发明实施例带来了以下有益效果：

本发明实施例提供的一种叠前深度偏移方法和装置，该方法包括获取目标工区的叠前地震数据；剔除该叠前地震数据中的奇异道得到待处理地震数据；根据该待处理地震数据计算目标工区各激发点的慢度表和等效q值表；根据该慢度表和该等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合；读取预先设定的每个成像点的偏移孔径，并根据该走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值；可以提高地震反射数据对地下构造的分辨能力。

本公开的其他特征和优点将在随后的说明书中阐述，或者，部分特征和优点可以从说明书推知或毫无疑义地确定，或者通过实施本公开的上述技术即可得知。

为使本公开的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例提供的一种叠前深度偏移方法的流程图；

图2为本发明实施例提供的一种慢度表在成像线位置处的等值线示意图；

图3为本发明实施例提供的一种等效q值表在成像线位置处的等值线示意图；

图4为本发明实施例提供的一种三维补偿吸收衰减的叠前时间偏移剖面示意图；

图5为本发明实施例提供的一种常规叠前深度偏移结果转换至时间域的剖面示意图；

图6为本发明实施例提供的一种补偿吸收衰减的叠前深度偏移结果转换至时间域的剖面示意图；

图7为本发明实施例提供的一种将图5与图6局部偏移剖面进行频谱分析对比的示意图；

图8为本发明实施例提供的一种叠前深度偏移装置的结构示意图。

图标：

81-数据获取模块；82-剔除模块；83-第一计算模块；84-第二计算模块；85-偏移幅值计算模块。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

在已发展的提高地震勘探分辨率方法中，补偿介质吸收的叠前时间偏移方法通过引入等效q值参数，将介质的吸收补偿与叠前时间偏移结合在一起，实现了沿地震波传播路径恢复地震波被衰减的高频成份，是一种高保真的提高地震勘探分辨率方法。然而，这类方法仅适用于速度横向变化缓慢的地质目标，通常应用于目的层埋深较浅的勘探工区。

对于深层-超深层勘探目标，介质速度存在强的横向变化；构造也较为复杂，地层陡倾。此种情况下，补偿介质吸收的叠前时间偏移方法中，基于均方根速度的走时计算算法难以实现正确的反射波偏移成像，所得到的断层、断点的横向位置也与真实构造存在一些差异。而补偿介质吸收的叠前深度偏移方法基于地层的层速度与层q值考虑地震波传播和吸收衰减过程，它既能在偏移成像过程中恢复地震波传播被衰减的高频成份，又能较好地考虑地震波在复杂构造中的实际传播路径，实现准确的偏移归位，从而在深层-超深层构造存在较强的速度横向变化时也能提高地震成像的分辨率。

上述技术尚存在两个方面的难点：第一，在补偿介质吸收的叠前深度偏移流程中，除了与常规叠前深度偏移类似，需基于工区层速度场计算沿射线传播路径的走时以外，还需考虑基于层速度场和层q场计算沿射线传播路径的走时与层q值的比值的累积效应，其计算量是巨大的。常规叠前深度偏移中，一般通过引入走时表来应对沿射线传播路径的走时计算这一问题，相应的，走时表的存取与插值也成为制约其计算效率与精度的关键。同样，沿射线传播路径的走时表以及走时与层q值的比值的累积效应表的存取与插值也是制约补偿介质吸收叠前深度偏移运算的计算效率和精度的关键环节。第二，补偿介质吸收叠前深度偏移计算的稳定性问题。补偿介质吸收的叠前深度偏移过程是一个对地震波传播过程中发生的吸收衰减进行的一个分频补偿计算过程，包含了振幅校正和相位校正，尤其是要对振幅项做一个e指数的补偿校正，若处理不当，极易在偏移过程中出现失稳现象，从而造成对地下构造成像失败。

目前，利用常规叠前深度偏移方法处理地震反射数据，得到的偏移成像结果仍存在构造越深分辨率就越低的问题，基于此，本发明实施例提供的一种叠前深度偏移方法和装置，可以提高地震反射数据对地下构造的分辨能力。

为便于对本实施例进行理解，首先对本发明实施例所公开的一种叠前深度偏移方法进行详细介绍。

实施例一

如图1所示，为本发明实施例提供的一种叠前深度偏移方法的流程图，由图1可见，该方法包括以下步骤：

步骤s101：获取目标工区的叠前地震数据。

步骤s102：剔除该叠前地震数据中的奇异道得到待处理地震数据。

步骤s103：根据该待处理地震数据计算目标工区各激发点的慢度表和等效q值表。

步骤s104：根据上述慢度表和等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合。

步骤s105：读取预先设定的每个成像点的偏移孔径，并根据上述走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值。

这样，本发明实施例提供的叠前深度偏移方法，通过将补偿介质吸收与叠前深度偏移结合在一起，实现了按照地震波传播路径进行介质吸收补偿，进而提高了地震反射数据对地下构造的分辨能力，对我国深层-超深层复杂构造油气、矿产资源勘探有重要应用价值。

实施例二

本发明实施例二针对叠前深度偏移方法的各个步骤及其实施措施进行更加具体详实的叙述。具体内容如下：

步骤1：获取目标工区的叠前地震数据。

在获取叠前地震数据之后，将该地震数据存放到集群计算机的不同计算节点，以待进行进一步处理。

步骤2：剔除该叠前地震数据中的奇异道得到待处理地震数据。

在每个计算节点，对全部地震道循环，应用快速傅立叶变换算法，对地震道的时间序列离散信号做傅立叶变换。

假设满足快速傅立叶算法的最佳时间采样点数为m，若地震道采样个数小于m则需补充零值，δτ是地震道的时间采样，则角频率采样δω为：

δω＝2π/(mδτ)(1)

令地震资料的主频为ω0，分别计算以ω00＝ω0，ω10＝2ω0，ω20＝2ω0+10，ω30＝2ω0+20为中心，向高频及低频方向分别扩展10hz后得到的频率窗内的能量ek0，其计算公式为：

式中，f表示地震数据变换至频率域，并且，

比较能量值e00，e01，e02，e03的大小，若e01>e00或e02>0.5e00或e03>0.1e00，则认为该地震道是奇异道，可予以剔除；否则予以保留。

步骤3：根据该待处理地震数据计算慢度表和等效q值表；这里，根据待处理地震数据，利用快速波前法计算慢度表和等效q值表。

基于三维深度域层速度模型，沿工区平面x方向和y方向间隔设定间距选取一系列的射线激发点pi，i＝1，2，3…n。对每一个激发点，按照不同的方位角和极角发射射线，形成初始的波前面。将该波前面分解为多个小三角形区域，每个小三角形的中心点作为该三角形区域的等时点，其走时均为τ，其对应的层q值可依据该中心点的位置信息自三维深度域层q模型读取。

对所有小三角形区域，利用射线追踪方法沿射线方向进行推进形成走时为τ+δτ的新波前面。检查构成新的波前面的三角形性质，当三角形面积大于预先设定值或三角形存在钝角时，将该三角形分割为多个小的三角形，同时插值得到新的走时为τ+δτ的等时点。依此类推，直到射线覆盖整个三维深度域层速度模型区域。

在波前面的推进过程中，记录沿射线路径的慢度信息es形成激发点慢度表，激发点慢度计算公式为：

式中，ek为射线激发点慢度，l为波前面序号，为该序号波前面形成用时，τk为激发点沿射线路径到等时点的走时，rk为射线激发点到等时点的距离。同时记录等效q值形成沿射线路径的激发点等效q值表，激发点等效q值的计算公式为：

式中，激发点等效q值，l为波前面序号，为该序号波前面形成用时，为该序号波前面等时点处的层q值，为方便起见，用走时q比符号表示

步骤4：根据上述慢度表和等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合。

其中，根据慢度表和等效q值表计算走时表和走时q比表的步骤包括：

(a1)选取坐标为yi的成像线，针对所有激发点，结合预先设定的偏移孔径，从上述慢度表和等效q值表中抽取与成像线临近的慢度面ej(x,y1,z)和ej(x,y2,z)以及等效q值面和y1≤yi≤y2。

(a2)调用一维插值算法计算该偏移孔径内所有激发点对该成像线有贡献的慢度面集合ej(x,yi,z)与等效q值面集合j＝1,2,…j，j为成像线y方向偏移孔径内所有激发点的个数。

(a3)逐一选取各地震道，分别读取该地震道的炮点坐标和检波点坐标，从该慢度面集合ej(x,yi,z)与等效q值面集合中选取4个与该炮点临近的激发点的慢度面和等效q值面，并从上述慢度面集合ej(x,yi,z)与等效q值面集合中选取4个与该检波点临近的激发点的慢度面和等效q值面。

(a4)调用双线性插值算法计算该炮点和该检波点在成像线处和该成像线的成像点网格相同的慢度面ek(xi,yi,zi)与等效q值面k为s或g，s代表炮点，g代表检波点。

(a5)根据上述慢度面ek(xi,yi,zi)和等效q值面计算各地震道在该成像线的走时表和走时q比表。

这里，在步骤(a5)中走时表和走时q比表的计算公式为：

τk(xi,yi,zi)＝ek(xi,yi,zi)/rk(xi,yi,zi)(7)

式中，k为s或g，s代表炮点，g代表检波点，τk(xi,yi,zi)为炮点或检波点沿射线路径到成像点(xi,yi,zi)的走时，rk(xi,yi,zi)为炮点或检波点到成像点(xi,yi,zi)的距离。

这样，通过步骤4即得到了炮点和检波点的走时表和走时q比表。

步骤5：读取预先设定的每个成像点的偏移孔径，并根据上述走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值。

这里，各个计算节点对于所分配的计算单元地震数据进行补偿介质吸收叠前深度偏移计算。其主要步骤包括：

(b1)计算每一个成像点的深度域倾角，并读取该成像点的偏移孔径。

(b2)识别该深度域倾角是否在该偏移孔径的边界范围内。

(b3)若是，则计算该成像点的偏移幅值。

对于步骤(b3)中成像点的偏移幅值的计算，首先，分别计算该成像点在预先设定的低、高截止频率内的频率成份对应的第一偏移幅值，以及该成像点在高频端衰减带中的频率成份对应的第二偏移幅值；其次，叠加第一偏移幅值和第二偏移幅值得到该成像点的偏移幅值。

其中，第一偏移幅值的计算公式为：

第二偏移幅值的计算公式为：

并且，在上述公式中：

在上述公式中，i(h,x,y,z)为第一偏移幅值，im(h,x,y,z)为第二偏移幅值，φ(iδη)为补偿因子光滑性阈值控制函数，m为衰减带厚度，取值范围10～15；g为振幅补偿阈值，取值范围800～1800；h为偏移距，x,y,z为空间位置坐标，表示全部成像点的偏移孔径约束，j为虚部单位，τs(x,y,z)为炮点到成像点的走时，τg(x,y,z)为检波点到成像点的走时，为炮点走时q比，为检波点走时q比，δω为频率间隔，ω0为地震数据主频，f表示地震数据变换至频率域，表示取实部，l1、l2分别为低、高截止频率对应的频率采样点编号。

通过上述步骤，最终得到成像结果为：

qmig(h,x,y,z)＝i(h,x,y,z)+im(h,x,y,z)(14)

本发明实施例提供的叠前深度偏移方法，通过存储基于工区层速度场计算的沿射线传播路径的慢度表代替走时表，同时以基于层速度场和层q场计算的沿射线路径的等效q值表代替走时q比的累积效应，避免了偏移过程中利用预先存储的走时表和走时q比表插值得到炮检点走时表及走时q比表时所需应用的复杂插值算法，大幅提高了偏移计算的效率。

实施例三

本发明实施例提供了另一种叠前深度偏移方法，该方法在上述实施例一、实施例二的基础上，在剔除叠前地震数据中的奇异道的步骤之前，还调用叠前地震数据分配算法将该叠前地震数据按负载平衡规则存放到不同计算节点。下面对该叠前地震数据分配算法的具体实施进行说明。

首先，上述叠前地震数据分配算法为对叠前地震数据进行虚拟分配和实际分配相结合的算法。并且，负载平衡规则为使计算节点上gpu卡所分配的不同偏移距组地震数据包含的地震道数目尽可能均衡。

假设地震资料的最大偏移距为hmax，最小偏移距为hmin，指定偏移距采样间距δh后，每个成像点处具有相同偏移距区间的地震道偏移结果叠加生成crp道集的一道偏移结果，因此，偏移距个数即成像点crp道集包含道数目为偏移距为h的地震道所处的偏移距索引可用整数表示。

假设集群计算机的计算节点上所有gpu卡的个数总和为ngpu。通常nh一般不等于ngpu，且在三维情况下，地震道数目随crp道集偏移距区间变化也是较大的，一般中间偏移距包含地震道数目多，近、远偏移距包含地震道数目少。需应用地震叠前数据分配算法，使得集群计算机的计算节点上各gpu卡对应的地震道数目均衡。遍历所有地震叠前数据，得到以ih为数组索引，该偏移距地震道数目为数组元素的工区偏移距-地震道数组。

对该数组进行虚拟分配处理的步骤如下：

(c1)给定一个分组标准，即每一计算单位内包含的地震道数目nidcal；对每一偏移距区间，将偏移距区间内包含的偏移距个数δh记为数组data的长度，偏移距区间内每个偏移距的序号记为数据的索引，其对应的地震道数目记为数组data的元素值，数组元素的和即偏移距区间内包含的地震道总数目为ntotal，利用非连续分组算法将该偏移距区间分为ntotal/nidcal组或ntotal/nidcal+1组，将每一个分组视为一个计算单元，使得各计算单元地震道数目平均值与nidcal相差最小，同时不同计算单元之间的地震道数目相差最小；记录达到这一要求的各偏移距区间的计算单元个数及每一计算单元的地震道数目。

其中，上述非连续分组算法介绍如下：

设有数组，其元素为data[i]，i＝1,2,3…n；

将数组data的元素分成m(m<<n)组，每组元素各自组合成一个新的数组，记为数组集合kr，r＝1,2,3…m。分组后，数组集合kr中每个数组的元素的和为sumr，r＝1,2,3…m，要求sumr,r＝1,2,3…m之间的差最小；

对数组data的元素按照降序排序，得到数组datasort[j]，j＝1,2,3…n；

将datasort[j]进行(int)(n/m+1)次分配；

对于第k次分配，k＝1,2,3…(int)(n/m+1)，

若k为奇数即k％2＝1，则将datasort本次分配的元素按照顺序分配给数组集合kr，r＝1,2,3…m，即：

datasort[km]分配给k1，datasort[km+1]分配给k2，…datasort[km+m-1]分配给km；

若k为偶数即k％2＝0，则将datasort本次分配的元素按照逆序分配给数组集合kr，r＝1,2,3…m；

datasort[km+m-1]分配给k1，datasort[km+m-2]分配给k2，…datasort[km]分配给km；

所得到的数组集合kr，r＝1,2,3…m中，每个数组的元素仍然是降序排列的；

计算分组后每个数组的元素的和sumr，r＝1,2,3…m，并搜索数值sumr，r＝1,2,3…m中两两之间差的绝对值的最大值maxdiff；

对于数组集合kr，r＝1,2,3…m，试验两两交换数组最后一个元素，并计算新的2个数组元素的和及其二者之间的差值的绝对值，看能否使得maxdiff变小，若能，执行交换；否则不执行交换。重复该操作，直到满足所有交换均无法使得maxdiff变小结束。

(c2)将所有计算单元的总个数npart记为数组data的大小，每个计算单元对应的地震道数目记为数组data的元素，再次应用非连续分组算法将所有的计算单元分为ngpu组，将每一个分组包含的计算单元分配给集群计算机的1块gpu卡，使得不同gpu卡分配的地震道数目相差最小。记录此时集群式计算机地震道数目最大的gpu卡与最小gpu卡之间的差值百分比d及集群式计算机分配的计算单元总个数npart。

(c3)改变计算量nidcal，重复上述步骤即可得到序列di，nparti，i＝1,2,3…100。

(c4)搜寻di及nparti，先比较di，将达到一定标准(如小于1％)的d组成新的序列dj及npartj；再比较npartj，选取更利于偏移结果回收的最小的npart，这一方案即为实现了各集群计算节点的偏移任务负载平衡的较佳方案，同时记录这一方案对应的数据偏移距分组及计算节点分配信息表。

在进行实际分配时，遍历所有地震叠前数据，首先依据地震道的偏移距确定计算单元，再依据所属计算单元确定gpu卡，然后，根据确定的gpu卡及计算单元信息即可把当前地震道发送至相应的集群计算机节点上。

本发明实施例提供的叠前深度偏移方法，通过使用叠前地震数据分配算法，可使集群计算机各计算节点分配的地震道数目基本一致，进而可使各计算节点偏移处理完成时间基本一致，提高了利用集群计算机进行补偿介质吸收叠前深度偏移运算的并行计算加速比。

实施例四

本发明实施例四介绍了利用叠前深度偏移方法对三维反射地震工业数据的处理情况，以进一步说明该方法的实际应用效果。

本次三维反射地震工业数据在进行地震采集时，采用共炮点观测，单点检波器接收，其中，炮线距为160米，接收线距为160米，每炮最大接收道数21504道，覆盖次数672次，最大炮检距5880米，cpp面元inline方向10米，crossline方向20米，时间采样间隔为2毫秒，样点数3001个，共12000余炮数据。

参见图2与图3，分别为自某激发点计算的慢度表和等效q值表在415号成像线位置处的等值线。由图2、图3可见，两种类型的等值线均比较平缓，这有利于提高基于该慢度表和等效q值表进行慢度和等效q值插值的精度。

如图4所示，为415号成像线上利用该数据集进行补偿介质吸收叠前时间偏移处理得到的剖面，考察可见，由于深层-超深层构造复杂，图4中上部矩形所示的火山口及其围岩以及下部矩形所示的致密天然气储层均未得到较好成像。

图5为415号成像线上利用该数据集，忽略地层层q模型，进行常规叠前深度偏移处理得到的剖面。为便于比较，这一结果转到了时间域，图中上部矩形所示的火山口及其围岩以及下部矩形所示的致密天然气储层均得到较好成像，但由于未考虑介质的吸收补偿，剖面的分辨率较低，导致对火山口及其围岩以及致密天然气储层的刻画不够精细。

图6为415号成像线上利用该数据集，考虑地层层q模型，利用本发明提供的叠前深度偏移方法进行处理得到的剖面。同样，这一结果也转到了时间域，观察可见，这一处理结果的分辨率明显提高，精细刻画了火山口及其围岩以及致密天然气储层的构造形态。

图7展示了将图5与图6中的致密天然气储层局部偏移剖面进行频谱分析对比的结果，在图7中，虚线为利用本发明提供的叠前深度偏移方法处理得到的剖面的频谱曲线，实线为利用常规叠前深度偏移方法处理得到的剖面的频谱曲线，由图7可见，利用本发明提供的叠前深度偏移方法比利用常规叠前深度偏移方法处理得到的剖面分辨率提高了20hz以上。

实施例五

本发明实施例还提供了一种叠前深度偏移装置，参见图8，为本发明实施例提供的一种叠前深度偏移装置的结构示意图，由图8可见，该装置包括依次相连的数据获取模块81、剔除模块82、第一计算模块83、第二计算模块84和偏移幅值计算模块85。其中，各个模块的功能如下：

数据获取模块81，用于获取目标工区的叠前地震数据；

剔除模块82，用于剔除该叠前地震数据中的奇异道得到待处理地震数据；

第一计算模块83，用于根据该待处理地震数据计算目标工区各激发点的慢度表和等效q值表；

第二计算模块84，用于根据该慢度表和该等效q值表计算走时表和走时q比表；该走时表包括炮点走时表和检波点走时表，该走时q比表包括炮点走时q比表和检波点走时q比表；该走时q比表为地震波沿射线路径传播的走时与相应射线路径处层q值的比值的累计数值的集合；

偏移幅值计算模块85，用于读取预先设定的每个成像点的偏移孔径，并根据该走时表、走时q比表和偏移孔径计算该目标工区各个成像点对应不同偏移距的偏移幅值。

本发明实施例所提供的装置，其实现原理及产生的技术效果和前述方法实施例相同，为简要描述，装置实施例部分未提及之处，可参考前述方法实施例中相应内容。

除非另外具体说明，否则在这些实施例中阐述的部件和步骤的相对步骤、数字表达式和数值并不限制本发明的范围。

在这里示出和描述的所有示例中，任何具体值应被解释为仅仅是示例性的，而不是作为限制，因此，示例性实施例的其他示例可以具有不同的值。

附图中的流程图和框图显示了根据本发明的多个实施例的系统、方法和计算机程序产品的可能实现的体系架构、功能和操作。在这点上，流程图或框图中的每个方框可以代表一个模块、程序段或代码的一部分，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。也应当注意，在有些作为替换的实现中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用执行规定的功能或动作的专用的基于硬件的系统来实现，或者可以用专用硬件与计算机指令的组合来实现。

另外，在本发明实施例的描述中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

在本发明的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。

最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以所述权利要求的保护范围为准。