技术领域
[0001] 本
发明涉及
地震勘探技术领域，特别涉及一种多波时间域匹配方法及装置。
背景技术
[0002] 由于多波地震数据中含有更为丰富的地下信息，在反演储层物性和含油气预测等方面显示出了良好的应用前景，多波
地震勘探技术于上世纪90年代以来得到广泛应用。相对于单一的纵波勘探无法准确识别的储层信息，多波的地震资料可以保存和挖掘纵波和转换波反射信息，在构造成像、岩性估计、
各向异性分析、
流体检测和储存描述与监测等方面有着广泛的应用。其中，多波地震资料中的纵横波匹配方法研究，是目前困扰多波地震资料处理以及后续反演的一个关键问题。因为纵波和横波传播速度的差异，转换波传播时间要比其对应的纵波传播时间要长。而若想充分利用多波地震资料就必须将两种波进行匹配。因此，良好的纵横波高
精度匹配方法，是充分利用多波进行联合反演、解释的重要
基础和保证，具有非常重要的意义。
[0003] 在现有的常规的纵横波时间匹配方法中，由于地下不均匀介质会对入射的纵波产生不同的响应，产生反射纵波和反射横波。所以早期的纵、横波匹配依赖于纵横波剖面中相似的成像特点，通过处理或者解释人员在纵横波地震剖面中标定出同一
地层的反射轴，并压缩转换波剖面至对应的纵波时间域。这种方法适用于反射特征比较明显的地区。而层位的标定过程可以通过两种途径完成：第一种是通过处理和解释人员的观察找出相似的反射轴；第二种方法是根据
测井资料制作纵横波合成地震记录对纵横波地震剖面进行标定。但是这种匹配方法依赖于人工的地质经验，匹配结果因人而异，具有较大的随意性和误差。
[0004] 假设纵横波具有相似反射系数和相同的子波
频率和
波长，Gaiser扫描如下公式(1)所示的纵横波相似系数对多波资料进行时间匹配，可以得到平均纵横
波速度比，并进而得到层间速度比，
[0005]
[0006] 其中PP(t1)和PS(t2)是纵波和横
波数据。Fomel和Backus采用最小二乘法进行PP波和PS波时间域自动同相轴匹配，减少了人工拾取层位的干扰。这种方法需要一个较好的初始模型，避免收敛结果陷入局部极小。基于纵横波的最大相似性目标函数，J.X.Yuan等利用模拟
退火算法对纵横波进行匹配。这种方法可以避免局部极小的问题，并且对初始模型误差有一个较好的适应性。但是基于纵横波归一化互相关的这种方法对于平均纵横波速度比有一定的要求。Gaiser指出，这种基于相关法的最优化方法要求其平均速度比有一个较为平滑的变化趋势，对纵横波速度比较大的变化适应性较差。所以，这种基于互相关的最优化方法仍具有一系列的问题，一般都需要一个先验信息，并对数据资料的速度比变化率、
信噪比等有较高的要求。
发明内容
[0007] 本发明
实施例提供了一种多波时间域匹配方法，以解决
现有技术中纵横波时间匹配方法对纵横波速度比中较大的变化适应性较差的技术问题。该方法包括：获取纵波和转换波的偏移剖面数据，所述偏移剖面数据为多维地震数据；在所述偏移剖面数据的每一个维向上确定预设数值个
种子点，拉平纵波同相轴，计算纵波
地震道的
叠加道包络，在纵向的时间维度上，纵向种子点是
能量值达到预设值且表征地下地层信息的反射点；在横向的地震道的空间维度上，相邻横向种子点之间保持预设空间间隔，使得所有横向种子点组成的连线完整贯穿在地震数据中；针对纵波的偏移剖面数据，针对每个
采样点，扫描获得该采样点与自身对应的横波的扫描范围内的采样点之间的纵波和横波的幅值平方差，作为该采样点的扫描误差，所述扫描范围为纵波地震记录时间长度；针对每个纵向种子点，按照纵波地震记录时间由小到大的方向，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，并按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累积结果，将该纵向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该纵向种子点处采样点本身对应的扫描误差，获得该纵向种子点处纵波和横波的平滑误差累积结果，所有纵向种子点的平滑误差累积结果形成第一平滑地震数据；将第一平滑地震数据转置90度，针对每个横向种子点，按照地震测线上地震道由小到大的方向，对该横向种子点与上一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的正向误差累积结果，并按照该地震道由大到小的方向，对该横向种子点与下一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的反向误差累积结果，将该横向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该横向种子点处采样点本身对应的扫描误差，获得该横向种子点处纵波和横波的平滑误差累积结果，所有横向种子点的平滑误差累积结果形成第二平滑地震数据；将第二平滑地震数据转置90度，针对除了纵向和横向之外的其他维向的每个种子点，获得每个种子点处纵波和横波的平滑误差累积结果，直至获得所述偏移剖面数据的所有维向的每个种子点处纵波和横波的平滑误差累积结果，并将最终得到的平滑地震数据转置得到与所述偏移剖面数据的排列形式一致的地震数据；对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场；根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场。
[0008] 在一个实施例中，针对每个纵向种子点，按照纵波地震记录时间由小到大的方向，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，包括：按照线性插值的方法得到该纵向种子点与上一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由小到大的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的正向误差累积结果。
[0009] 在一个实施例中，针对每个纵向种子点，按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累计结果，包括：按照线性插值的方法得到该纵向种子点与下一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由大到小的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的反向误差累积结果。
[0010] 在一个实施例中，对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场，包括：通过以下公式回追得到各个纵向种子点处纵波和横波的时间误差：
[0011]
[0012]
[0013] i＝N-1,N-2,...1
[0014] 其中，argmin是指在所有的搜索路径中寻找使d[i-1,l]最小的u[i-1]的路径；u[i-1]是在纵向种子点i-1处平滑误差累积结果最小点对应的纵波和横波的时间误差，是按照时间由大到小的方向逐步回追得到的；d[i-1,l]是在纵向种子点i-1处按照所述扫描范围对扫描误差进行累计得到的平滑误差累积结果；N是纵向种子点的个数；u(0:N-1)是从最后一个纵向种子点到第一个纵向种子点进行反向求取得到的结果，代表从最后一个纵向种子点到第一个纵向种子点处纵波时间域的纵波和横波的时间误差；根据在纵向种子点处纵波和横波的时间误差进行线性插值，得到每个采样点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场。
[0015] 在一个实施例中，通过以下公式根据多维偏移剖面数据的纵波和横波的时差场计算多维偏移剖面数据的速度比场：
[0016]
[0017] 其中，misfit(Tpp)是纵波时间域的纵波和横波的时差场，γi(tpp)是纵波时间域的纵波和横波的速度比场。
[0018] 本发明实施例还提供了一种多波时间域匹配装置，以解决现有技术中纵横波时间匹配方法对纵横波速度比中较大的变化适应性较差的技术问题。该装置包括：数据获取模
块，用于获取纵波和转换波的偏移剖面数据，所述偏移剖面数据为多维地震数据；种子点确定模块，用于在所述偏移剖面数据的每一个维向上确定预设数值个种子点，拉平纵波同相轴，计算纵波地震道的叠加道包络，在纵向的时间维度上，纵向种子点是能量值达到预设值且表征地下地层信息的反射点；在横向的地震道的空间维度上，相邻横向种子点之间保持预设空间间隔，使得所有横向种子点组成的连线完整贯穿在地震数据中；第一误差平滑处理模块，用于针对纵波的偏移剖面数据，针对每个采样点，扫描获得该采样点与自身对应的横波的扫描范围内的采样点之间的纵波和横波的幅值平方差，作为该采样点的扫描误差，所述扫描范围为纵波地震记录时间长度；针对每个纵向种子点，按照纵波地震记录时间由小到大的方向，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，并按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累积结果，将该纵向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该纵向种子点处采样点本身对应的扫描误差，获得该纵向种子点处纵波和横波的平滑误差累积结果，所有纵向种子点的平滑误差累积结果形成第一平滑地震数据；第二误差平滑处理模块，用于将第一平滑地震数据转置90度，针对每个横向种子点，按照地震测线上地震道由小到大的方向，对该横向种子点与上一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的正向误差累积结果，并按照该地震道由大到小的方向，对该横向种子点与下一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的反向误差累积结果，将该横向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该横向种子点处采样点本身对应的扫描误差，获得该横向种子点处纵波和横波的平滑误差累积结果，所有横向种子点的平滑误差累积结果形成第二平滑地震数据；
数据处理模块，用于将第二平滑地震数据转置90度，针对除了纵向和横向之外的其他维向的每个种子点，获得每个种子点处纵波和横波的平滑误差累积结果，直至获得所述偏移剖面数据的所有维向的每个种子点处纵波和横波的平滑误差累积结果，并将最终得到的平滑地震数据转置得到与所述偏移剖面数据的排列形式一致的地震数据；时差数据处理模块，用于对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场；速度比数据处理模块，用于根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场。
[0019] 在一个实施例中，所述第一误差平滑处理模块，包括：正向误差平滑处理单元，用于按照线性插值的方法得到该纵向种子点与上一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由小到大的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的正向误差累积结果。
[0020] 在一个实施例中，所述第一误差平滑处理模块，还包括：反向误差平滑处理单元，用于按照线性插值的方法得到该纵向种子点与下一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由大到小的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的反向误差累积结果。
[0021] 在一个实施例中，所述时差数据处理模块，包括：种子点处时差计算单元，用于通过以下公式回追得到各个纵向种子点处纵波和横波的时间误差：
[0022]
[0023]
[0024] i＝N-1,N-2,...1
[0025] 其中，argmin是指在所有的搜索路径中寻找使d[i-1,l]最小的u[i-1]的路径；u[i-1]是在纵向种子点i-1处平滑误差累积结果最小点对应的纵波和横波的时间误差，是按照时间由大到小的方向逐步回追得到的；d[i-1,l]是在纵向种子点i-1处按照所述扫描范围对扫描误差进行累计得到的平滑误差累积结果；N是纵向种子点的个数；u(0:N-1)是从最后一个纵向种子点到第一个纵向种子点进行反向求取得到的结果，代表从最后一个纵向种子点到第一个纵向种子点处纵波时间域的纵波和横波的时间误差；采样点处时差计算单元，用于根据在纵向种子点处纵波和横波的时间误差进行线性插值，得到每个采样点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场。
[0026] 在一个实施例中，所述速度比数据处理模块通过以下公式根据多维偏移剖面数据的纵波和横波的时差场计算多维偏移剖面数据的速度比场：
[0027]
[0028] 其中，misfit(Tpp)是纵波时间域的纵波和横波的时差场，γi(tpp)是纵波时间域的纵波和横波的速度比场。
[0029] 在本发明实施例中，针对多维偏移剖面数据每一个维向上的种子点，通过扫描获得每个种子点处纵波和横波的正向误差累积结果，同时，通过扫描获得每个种子点处纵波和横波的反向误差累积结果，将每个种子点处的正向误差累积结果和反向误差累积结果累加，再减去该种子点处采样点本身对应的扫描误差，获得每个种子点处纵波和横波的平滑误差累积结果，由于种子点的选择是保持预设空间间隔的，且由于各种子点处纵波和横波的平滑误差累积结果是通过对轨迹进行误差平方的累积得到的，可以增加反演时差的
稳定性；同时，根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场，可以得到变化率较小的纵横和横波的速度比。
附图说明
[0030] 此处所说明的附图用来提供对本发明的进一步理解，构成本
申请的一部分，并不构成对本发明的限定。在附图中：
[0031] 图1是本发明实施例提供的一种多波时间域匹配方法的
流程图；
[0032] 图2是本发明实施例提供的一种一维纵波(PP wave)和横波(PS wave)模型数据；
[0033] 图3是本发明实施例提供的一种一维模型纵波时间域纵横波时差的示意图；
[0034] 图4是本发明实施例提供的一种一维模型纵波时间域纵横波速度比的示意图；
[0035] 图5是本发明实施例提供的一种平滑动态时间调整法的原理示意图；
[0036] 图6是本发明实施例提供的一种合成纵波和转换波叠后剖面的示意图；
[0037] 图7是本发明实施例提供的一种纵波时间域的纵横波时差和速度比场的示意图；
[0038] 图8是本发明实施例提供的一种含噪音的纵横波合成地震数据的示意图；
[0039] 图9是本发明实施例提供的一种采用本申请中的动态平滑时差纵横波时间匹配方法得到的纵横波时差场和速度比场的示意图；
[0040] 图10是本发明实施例提供的一种实际纵波叠后资料的示意图；
[0041] 图11是本发明实施例提供的一种实际转换波叠后资料的示意图；
[0042] 图12是本发明实施例提供的一种纵波资料中种子点分布情况的示意图；
[0043] 图13是本发明实施例提供的一种采用本申请的动态平滑时差纵横波时间匹配方法得到的纵横波时差场和速度比场的示意图；
[0044] 图14是本发明实施例提供的一种时间匹配后的转换波地震资料的示意图。
[0045] 图15是本发明实施例提供的一种多波时间域匹配装置的结构
框图。
具体实施方式
[0046] 为使本发明的目的、技术方案和优点更加清楚明白，下面结合实施方式和附图，对本发明做进一步详细说明。在此，本发明的示意性实施方式及其说明用于解释本发明，但并不作为对本发明的限定。
[0047] 在本发明实施例中，提供了一种多波时间域匹配方法，如图1所示，该方法包括：
[0048] 步骤101：获取纵波和转换波的偏移剖面数据，所述偏移剖面数据为多维地震数据；
[0049] 步骤102：在所述偏移剖面数据的每一个维向上确定预设数值个种子点，拉平纵波同相轴，计算纵波地震道的叠加道包络，在纵向的时间维度上，纵向种子点是能量值达到预设值且表征地下地层信息的反射点；在横向的地震道的空间维度上，相邻横向种子点之间保持预设空间间隔，使得所有横向种子点组成的连线完整贯穿在地震数据中；
[0050] 步骤103：针对纵波的偏移剖面数据，针对每个采样点，扫描获得该采样点与自身对应的横波的扫描范围内的采样点之间的纵波和横波的幅值平方差，作为该采样点的扫描误差，所述扫描范围为纵波地震记录时间长度；针对每个纵向种子点，按照纵波地震记录时间由小到大的方向，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，并按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累积结果，将该纵向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该纵向种子点处采样点本身对应的扫描误差，获得该纵向种子点处纵波和横波的平滑误差累积结果，所有纵向种子点的平滑误差累积结果形成第一平滑地震数据；
[0051] 步骤104：将第一平滑地震数据转置90度，针对每个横向种子点，按照地震测线上地震道由小到大的方向，对该横向种子点与上一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的正向误差累积结果，并按照该地震道由大到小的方向，对该横向种子点与下一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的反向误差累积结果，将该横向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该横向种子点处采样点本身对应的扫描误差，获得该横向种子点处纵波和横波的平滑误差累积结果，所有横向种子点的平滑误差累积结果形成第二平滑地震数据；
[0052] 步骤105：将第二平滑地震数据转置90度，针对除了纵向和横向之外的其他维向的每个种子点，获得每个种子点处纵波和横波的平滑误差累积结果，直至获得所述偏移剖面数据的所有维向的每个种子点处纵波和横波的平滑误差累积结果，并将最终得到的平滑地震数据转置得到与所述偏移剖面数据的排列形式一致的地震数据；
[0053] 步骤106：对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场；
[0054] 步骤107：根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场。
[0055] 由图1所示的流程可知，在本发明实施例中，通过扫描获得每个采样点处纵波和横波的扫描误差，针对多维偏移剖面数据，在每一个维向上确定预设数值个种子点，通过种子点线性插值确定的轨迹，按照纵波地震记录时间由小到大的方向累积采样点对应的扫描误差，获得种子点的正向误差累积结果，同时，按照纵波地震记录时间由大到小的方向累积采样点对应的扫描误差，获得种子点的反向误差累积结果，将每个种子点处的正向误差累积结果和反向误差累积结果累加并减去种子点本身的扫描误差，获得每个种子点处纵波和横波的平滑误差累积结果，这个过程称之为双向平滑过程，由于种子点的选择是保持预设空间间隔的，且由于各种子点处纵波和横波的平滑误差累积结果是通过对轨迹进行误差平方的累积得到的，可以增加反演时差的稳定性；同时，根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场，可以得到变化率较小也更稳定的纵横和横波的速度比。
[0056] 具体实施时，为了得到平滑的时差，本申请
发明人在现有技术中的动态时间校正法(DTW)的基础上进行改进，得到平滑动态时间调整法(Smoothing Dynamic Time Warping，SDTW)。
[0057] 具体的，动态时间校正法原理如下：
[0058] 对于两个
波形上存在相似性而时间上有误差的地震
信号f(i)和g(i)，我们可以使用动态时间校正法计算地震信号f(i)和g(i)的时间误差u(0:n-1)。
[0059]
[0060] 其中，u(0:N-1)代表了一个时间匹配时差的最优解。N是地震信号f(i)和g(i)的长度。当u[j]-u[j-1]＝1，f(i)中两个相邻点对应着g(i)中两个不相邻的点。而当u[j]-u[j-1]＝-1，f(i)中两个相邻点对应着g(i)中的一个点。因为两个相似的
地震波不会变化的很剧烈，所以公式(2)的限制是具有实际的地质意义的。下面进行的过程都是在公式(2)的限制下进行的。
[0061] DTW的第一步是累积。我们定义扫描误差为：
[0062]
[0063] 其中，l[i]表示在第i个采样点处扫描两个地震信号的时间误差的量，扫描范围为L，e[i,l]是纵波的第i个采样点，对应于转换波的第i+l个采样点，时差为l的误差扫描结果；
[0064]
[0065] d[i,l]是从0到第i个采样点的扫描误差累积求和的量。
[0066] 而将一维的DTW方法应用到多维的数据时，如果将公式(4)称之为从第一个采样点至第N-1个采样点的正向累计过程，则对应的从最后一个采样点至第一个采样点的扫描误差的反向累积过程为：
[0067]
[0068] DTW的第二步是回追，
[0069]
[0070]
[0071] i＝N-1,N-2,...1 (6)
[0072] 其中，“argmin”代表着在所有的搜索路径中寻找使d[i-1,l]最小的u[i-1]的路径。这个回追过程与累积相反，从N-1到0进行。求得的u(0:N-1)代表两个地震信号的时间误差的最优解。
[0073] 本申请发明人提出的平滑动态时间调整法原理如下：
[0074] 我们设计一个模型，将深度域一致的随机反射系数序列转换到纵波和转换波时间域，并用Ricker子波褶积得到如图2所示的合成纵波(PP wave)和转换波(PS wave)数据，其时间域的纵横波时差和对应的速度比如图3和图4所示。可以看到，纵横波的时差在纵波时间域是一个平滑的增大过程。采用“以直代曲”思想，我们可以用分段直线来拟合这种平滑的时间差。利用DTW算法中公式(3)，得到其扫描误差如图5所示。如果我们扩大DTW方法中的每一步的累积的时差范围，从一个采样点拓展到dx个采样点，按照新的间隔dx采样点作为分段直线的长度，其两个端点作为地震道中的新的“采样点”，这里我们可以称其为种子点。与DTW方法一致，我们可以假设时差变化幅度为1，从每个种子点的时差为k处，到上一个种子点的时差k至k-dx，由种子点线性插值定义的不同轨迹(如图5中黑色实线所示)对k至k-dx范围内采样点的扫描误差进行累积，可以得到时差为k处的种子点误差累积结果。由于种子点处的时差为整数倍采样间隔，而线性轨迹在非种子点处的采样点大都位于网格中间，但是实践证明，简单的相邻网格的线性插值或者直接取整计算也可以满足绝大多数的需求。最终，如将种子点看作重采样的稀疏的采样点，则可利用公式(6)求得对应种子点的时间误差值，并最终使用与累积过程轨迹相同的线性插值得到每个采样点的时间误差值。由于时间误差的平滑性，这种直线轨迹可以作为一种规则化方法来增加反演时差的稳定性。
这种方法称为平滑动态时间调整法(Smoothing Dynamic time warping)，简称SDTW方法。
[0075] 对于多维的地震数据，我们需要先对一维的纵波地震道和对应的转换波地震道做一个双向的平滑处理，针对每个纵向种子点，按照纵波地震记录时间由小到大的方向，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，包括：按照线性插值的方法得到该纵向种子点与上一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由小到大的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的正向误差累积结果，即实现正向累计过程；
[0076] 针对每个纵向种子点，按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累计结果，包括：按照线性插值的方法得到该纵向种子点与下一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由大到小的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的反向误差累积结果，即实现反向累计过程。
[0077] 例如，我们采用公式(4)对从第一个采样点至第N-1个采样点的扫描误差进行正向累计，得到种子点的正向误差累积结果 则我们可以将地震数据看做一个传播时间从大到小的逆过程，进行一个反向的、相同种子点的反向平滑累积过程，即采用公式(5)对从最后一个采样点至第一个采样点的扫描误差进行反向累计，得到种子点的反向误差累积结果 并通过以下公式(7)将种子点处的正向误差累积结果和反向误差累积结果相加，再减去该种子点处采样点本身对应的扫描误差e[i,l]，获得该种子点处纵波和横波的平滑误差累计结果，即完成一个双向平滑处理：
[0078]
[0079] 公式(7)中从累积结果中减掉种子点自身对应的扫描误差e[i,l]，是对应于每一个采样点的扫描误差的一次求和过程，避免了在每一个采样点处的两次累积。需要注意的是，虽然公式(7)进行了一个双向的平滑滤波，但是对于一维数据，如果我们将平滑结果直接进行如公式(7)所示的回追处理，我们将会发现并没有对结果有任何改善。这是因为公式(7)所进行的只是一个一维的平滑，对于SDTW这种全局寻优方法来说，本质上与一个单方向的累积结果没有区别。而在多维数据中，这种双向的平滑方法需要在不同的维度交替进行，这种方法对于多维数据的良好的抗噪性以及求解的精确性。
[0080] 具体的，由于本申请在进行双向平滑处理过程中，在进行正向累计和反向累计过程中，是先按照线性插值的方法得到该纵向种子点与下一个(或上一个)纵向种子点之间的不同时差的直线轨迹，进而对不同直线轨迹对应的采样点的扫描误差进行累计，不是对两个种子点之间的采样点逐个进行累计，即本申请中的正向累计和反向累计过程添加了时差的轨迹的限定，对纵向种子点与下一个(或上一个)纵向种子点之间的不同直线轨迹对应的采样点的扫描误差进行累计，可以约束采样点直接的变化趋势保持一致。
[0081] 具体实施时，获得偏移剖面数据的所有维向的每个种子点处纵波和横波的平滑误差累积结果后，将最终得到的平滑地震数据转置得到与所述偏移剖面数据的排列形式一致的地震数据，并采用上述公式(6)对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时差，再根据在纵向种子点处纵波和横波的时间差进行线性插值，得到每个采样点处纵波和横波的时间差，得到纵波时间域的纵波和横波的时差场。此时，公式(6)中的argmin是指在所有的搜索路径中寻找使d[i-1,l]最小的u[i-1]的路径；u[i-1]是在纵向种子点i-1处平滑误差累积结果最小点对应的纵波和横波的时间误差，是按照时间由大到小的方向逐步回追得到的；d[i-1,l]是在纵向种子点i-1处按照所述扫描范围对扫描误差进行累计得到的平滑误差累积结果；N是纵向种子点的个数；u(0:N-1)从最后一个纵向种子点到第一个纵向种子点进行反向求取得到的结果，代表从最后一个纵向种子点到第一个纵向种子点处纵波时间域的纵波和横波的时间误差。
[0082] 具体实施时，对于二维和三维地震数据，本方法并无本质区别，因而此处我们以二维数据为例进行说明。对于具有M道和N个采样点二维数据的纵波数据，我们首先对其进行层拉平处理(该层拉平处理方法可以参照现有技术中的方法，例如，Lomark等，2006年提出的方法；Fomel，2010年提出的方法)，并进行叠加、求取叠加道包络，并按选取种子点的条件选择叠加道的纵向种子点
位置及其个数KZ。利层拉平的逆过程，得到纵向种子点在整个剖面地震数据中的分布。然后，针对每个纵向种子点，对纵横波地震数据中的每一对地震道数据进行一个一维的双向平滑处理(只进行双向累积，并不反向回追)，获得每个纵向种子点处纵波和横波的平滑误差累积结果，所有纵向种子点的平滑误差累积结果形成三维的大小为M*KZ*L平滑结果数据。将平滑结果数据进行90度转置处理，看作一个新的KZ道和M个采样点的、大小为KZ*M*L的地震数据误差扫描结果，对于这个新的地震数据误差扫描结果，参考原始数据纵向的双向平滑处理过程，针对每个横向种子点，进行对应的横向的双向平滑处理过程，获得每个横向种子点处纵波和横波的平滑误差累计结果。这里需要注意的是由于叠后地震剖面构造特征多呈横向展布，而在垂向上连续性较差，为减少内存需求，一般可以选取等间隔(即等地震道间隔)的横向种子点分布，无需另外计算其包络大小。我们可以将M个采样点等距离的划分为不同的线段，选取对应的横向种子点，个数为KX。最后将计算得到的平滑结果重新进行一个-90度的转置处理，得到一个数据大小为KX*KZ*L三维平滑结果。
[0083] 最后，对最终的大小为KX*KZ*L平滑结果数据，我们对KX个横向种子点所在的每一地震道的大小为Kz*L的平滑误差累计结果，采用一维的回追方法求取其对应的纵向种子点的时差值，得到二维剖面的纵向种子点处的时差场数据，数据大小为KX*KZ。利用符合累积过程的双线性插值，我们可以得到原始数据的每个采样点的时差，数据大小为M*N。还可以根据二维剖面的纵向种子点处的时差场数据得到对应的纵横波速度比场分布。具体实施时，通过以下公式根据多维偏移剖面数据的纵波和横波的时差场计算多维偏移剖面数据的速度比场：
[0084]
[0085] 其中，misfit(Tpp)是纵波时间域的纵波和横波的时差场，γi(tpp)是纵波时间域的纵波和横波的速度比场。对于地下介质来说，其纵向的纵横波速度比变化率较小，纵横波的时差misfit(Tpp)(纵波时间域)只与其对应的纵横波速度比γi(tpp)有关，因而其时差具有缓慢变化的趋势。
[0086] 以下结合具体的示例来说明上述多波时间匹配方法。
[0087] 示例1.我们建立深度域模型并根据纵波和横波速度转换到对应的纵波和转换波时间域，得到对应的纵波和转换波数据如图6所示。对应的纵波时间域的纵横波时差和速度比分布如图7所示。为验证本申请多波时间匹配方法的稳定性，对原始数据添加随机噪音，使信噪比为2：1，如图8所示。应用本申请的多波时间匹配方法，我们得到含噪音数据对应的纵横波时差和反演速度比如图9所示。从图中可以看到，本方法对于含有较强噪音的数据具有良好的稳定性，求得到的时差符合模型数据，并可以得到符合地下层位分布的良好的速度比场结果。
[0088] 实例2.我们使用实际的经过处理后的纵波数据(图10)和转换波数据测试本方法(图11)。对纵波数据进行层拉平处理，得到其种子点在剖面中的分布情况如图12中白色线所示。使用本申请的多波时间匹配方法，可以得到其纵横波时差和速度比场如图13所示。根据纵横波时差，对转换波进行时间校正，得到结果如图14，从中可以看到，纵波和转换波可以得到良好的自动时间匹配，避免了人工干扰。
[0089] 基于同一发明构思，本发明实施例中还提供了一种多波时间匹配装置，如下面的实施例所述。由于多波时间匹配装置解决问题的原理与多波时间匹配方法相似，因此多波时间匹配装置的实施可以参见多波时间匹配方法的实施，重复之处不再赘述。以下所使用的，术语“单元”或者“模块”可以实现预定功能的
软件和/或
硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现，但是硬件，或者软件和硬件的组合的实现也是可能并被构想的。
[0090] 图15是本发明实施例的多波时间域匹配装置的一种结构框图，如图15所示，该多波时间匹配装置包括：
[0091] 数据获取模块1501，用于获取纵波和转换波的偏移剖面数据，所述偏移剖面数据为多维地震数据；
[0092] 种子点确定模块1502，用于在所述偏移剖面数据的每一个维向上确定预设数值个种子点，拉平纵波同相轴，计算纵波地震道的叠加道包络，在纵向的时间维度上，纵向种子点是能量值达到预设值且表征地下地层信息的反射点；在横向的地震道的空间维度上，相邻横向种子点之间保持预设空间间隔，使得所有横向种子点组成的连线完整贯穿在地震数据中；
[0093] 第一误差平滑处理模块1503，用于针对纵波的偏移剖面数据，针对每个采样点，扫描获得该采样点与自身对应的横波的扫描范围内的采样点之间的纵波和横波的幅值平方差，作为该采样点的扫描误差，所述扫描范围为纵波地震记录时间长度；针对每个纵向种子点，对该纵向种子点与上一个纵向种子点之间的采样点的扫描误差进行累计，获得该纵向种子点的正向误差累积结果，并按照纵波地震记录时间由大到小的方向，对该纵向种子点与下一个纵向种子点之间的采样点的扫描误差进行累计，得到该纵向种子点的反向误差累积结果，将该纵向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该纵向种子点处采样点本身对应的扫描误差，获得该纵向种子点处纵波和横波的平滑误差累积结果，所有纵向种子点的平滑误差累积结果形成第一平滑地震数据；
[0094] 第二误差平滑处理模块1504，用于将第一平滑地震数据转置90度，针对每个横向种子点，按照地震测线上地震道由小到大的方向，对该横向种子点与上一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的正向误差累积结果，并按照该地震道由大到小的方向，对该横向种子点与下一个横向种子点之间的采样点的扫描误差进行累计，获得该横向种子点的反向误差累积结果，将该横向种子点处的正向误差累积结果和反向误差累积结果相加，再减去该横向种子点处采样点本身对应的扫描误差，获得该横向种子点处纵波和横波的平滑误差累积结果，所有横向种子点的平滑误差累积结果形成第二平滑地震数据；
[0095] 数据处理模块1505，用于将第二平滑地震数据转置90度，针对除了纵向和横向之外的其他维向的每个种子点，获得每个种子点处纵波和横波的平滑误差累积结果，直至获得所述偏移剖面数据的所有维向的每个种子点处纵波和横波的平滑误差累积结果，并将最终得到的平滑地震数据转置得到与所述偏移剖面数据的排列形式一致的地震数据；
[0096] 时差数据处理模块1506，用于对转置得到的地震数据回追，得到各个纵向种子点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场；
[0097] 速度比数据处理模块1507，用于根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场。
[0098] 在一个实施例中，所述第一误差平滑处理模块，包括：正向误差平滑处理单元，用于按照线性插值的方法得到该纵向种子点与上一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由小到大的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的正向误差累积结果。
[0099] 在一个实施例中，所述第一误差平滑处理模块，还包括：反向误差平滑处理单元，用于按照线性插值的方法得到该纵向种子点与下一个纵向种子点之间的不同时差的直线轨迹，对不同直线轨迹对应的采样点，按照纵波地震记录时间由大到小的方向，从第一个采样点到最后一个采样点进行扫描误差的累计，得到该纵向种子点的反向误差累积结果。
[0100] 在一个实施例中，所述时差数据处理模块，包括：种子点处时差计算单元，用于通过以下公式回追得到各个纵向种子点处纵波和横波的时间误差：
[0101]
[0102]
[0103] i＝N-1,N-2,...1
[0104] 其中，argmin是指在所有的搜索路径中寻找使d[i-1,l]最小的u[i-1]的路径；u[i-1]是在纵向种子点i-1处平滑误差累积结果最小点对应的纵波和横波的时间误差；d[i-1,l]是在纵向种子点i-1处按照所述扫描范围对扫描误差进行累计得到的平滑误差累积结果；N是纵向种子点的个数；u(0:N-1)是从最后一个纵向种子点到第一个纵向种子点进行反向求取得到的结果，代表从最后一个纵向种子点到第一个纵向种子点处纵波时间域的纵波和横波的时间误差；采样点处时差计算单元，用于根据在纵向种子点处纵波和横波的时间误差进行线性插值，得到每个采样点处纵波和横波的时间误差，得到纵波时间域的纵波和横波的时差场。
[0105] 在一个实施例中，所述速度比数据处理模块通过以下公式根据多维偏移剖面数据的纵波和横波的时差场计算多维偏移剖面数据的速度比场：
[0106]
[0107] 其中，misfit(Tpp)是纵波时间域的纵波和横波的时差场，γi(tpp)是纵波时间域的纵波和横波的速度比场。
[0108] 在本发明实施例中，通过扫描获得每个采样点处纵波和横波的扫描误差，针对多维偏移剖面数据，在每一个维向上确定预设数值个种子点，通过种子点线性插值确定的轨迹，按照纵波地震记录时间由小到大的方向累积采样点对应的扫描误差，获得种子点的正向误差累积结果，同时，按照纵波地震记录时间由大到小的方向累积采样点对应的扫描误差，获得种子点的反向误差累积结果，将每个种子点处的正向误差累积结果和反向误差累积结果累加并减去种子点本身的扫描误差，获得每个种子点处纵波和横波的平滑误差累积结果，这个过程称之为双向平滑过程，由于种子点的选择是保持预设空间间隔的，且由于各种子点处纵波和横波的平滑误差累积结果是通过对轨迹进行误差平方的累积得到的，可以增加反演时差的稳定性；同时，根据纵波时间域的纵波和横波的时差场计算纵波时间域的纵波和横波的速度比场，可以得到变化率较小也更稳定的纵横和横波的速度比。
[0109] 显然，本领域的技术人员应该明白，上述的本发明实施例的各模块或各步骤可以用通用的计算装置来实现，它们可以集中在单个的计算装置上，或者分布在多个计算装置所组成的网络上，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，并且在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤，或者将它们分别制作成各个集成
电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明实施例不限制于任何特定的硬件和软件结合。
[0110] 以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何
修改、等同替换、改进等，均应包含在本发明的保护范围之内。
