基于折射波分析地质体结构的方法及存储介质

技术领域

本发明属于地震资料处理领域，具体涉及一种基于折射波分析地质体结构的方法及内部存储有基于折射波分析地质体结构的方法的存储介质。

背景技术

目前，在地球物理勘探中，地震资料大多使用地面地震反射波、井间透射波、VSP(Vertical Seismic Profiling，即垂直地震剖面)上、下行波等进行处理解释。

然而，使用地面地震反射波、井间透射波、VSP上、下行波等进行处理解释存在以下问题：

第一，由于浅层反射的地震信息分布于近道，且受附近面波的严重影响，难以提取到精确的地震反射子波；VSP上、下行波不易分离；井间透射波通常是地面地震反射波的辅助参考。

第二，利用地震反射波对深层基底构造形态进行勘探会因深层反射能量变得极其微弱，有效信息的可靠性也随之变差，并且其无法很好地识别一些地下薄层、裂缝等小尺度地质结构。

第三，地震反射波能量大多仅集中于近偏移距，无法进行远偏移距处信息的采集与成像。

品质因子Q是体现地下介质对地震波能量吸收、衰减的一个重要参数，能够表征地层的特征属性，可用于岩石性质、断层和流体分布的识别与解释，因此精确估计Q值对研究地下介质的岩性、物性、结构和含流体性质等具有重要的意义。地质品质因子可以通过Q层析成像计算。

现在亟须一种基于折射波分析地质体结构的方法及内部存储有基于折射波分析地质体结构的方法的存储介质。

发明内容

本发明所要解决的技术问题是无法精确获取地震子波而导致地质体结构分析结果不准的问题。

为了解决上述技术问题，本发明提供了一种基于折射波分析地质体结构的方法，包括：

S100，将目标储层的地质体剖面图区隔成多个地质网格，以形成多个地质网格点；

S200，设置各地质网格点处地震折射波的地下传播速度的初始值，以及构建地震折射波从震源到各地质网格点的射线路径，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值；

S300，基于地质品质因子反演模型，根据各地质网格点处地震折射波的地下传播速度的修正值和各射线路径的地质吸收特征时间，确定各地质网格点的地质品质因子；

S400，根据各地质网格点的地质品质因子的分布分析目标储层的地质体结构。

优选地，所述步骤S200包括以下步骤：

S210，设置各地质网格点处地震折射波的地下传播速度的初始值；

S220，构建地震折射波从震源到各地质网格点的射线路径，根据地下传播速度的初始值确定地震折射波沿射线路径传播的理论初至旅行时；

S230，基于地下速度反演模型，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值。

优选地，所述步骤S200包括以下步骤：

S210，设置各地质网格点处地震折射波的地下传播速度的初始值；

S220，构建地震折射波从震源到各地质网格点的射线路径，根据地下传播速度的初始值确定地震折射波沿射线路径传播的理论初至旅行时；

S230，基于地下速度反演模型，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值；

S240，判断对应于当前地下传播速度的修正值的理论初至旅行时与实际初至旅行时之间的差值是否小于预设的阈值：

若是，执行步骤S300；

若否，令地下传播速度的初始值等于当前地下传播速度的修正值，返回步骤S220，以对各地质网格点处地震折射波的地下传播速度再次修正。

优选地，所述步骤S200中，设置各地质网格点处地震折射波的地下传播速度的初始值，使其符合均匀变化梯度的二维速度分布；通过射线追踪法构建地震折射波从震源到各地质网格点的射线路径。

优选地，所述步骤S200中，地震折射波从震源到各地质网格点的实际初至旅行时为在未叠加的单炮地震记录剖面上拾取地质网格点处地震折射波的实际初至旅行时。

优选地，所述步骤S400包括以下步骤：

S410，将各地质网格点的地质品质因子对应填入地质体剖面图中，以得到地质品质因子的结构分布；

S420，根据地质体剖面图中地质品质因子的分布分析目标储层的地质体结构。

优选地，所述地下速度反演模型如下：

Bδv＝δt

其中，B为地下速度反演系数矩阵，δv为各地质网格点地震折射波的地下传播速度的修正值与初始值之间的差值，δt为各地质网格点的地震折射波的实际初至旅行时与理论初至旅行时之间的差值。

优选地，所述步骤S300中，通过以下表达式获得各射线路径的地质吸收特征时间：

其中，t

优选地，所述步骤S300中，地质品质因子反演模型为：

其中，Q(x,z)为地质网格点的地质品质因子，v

根据本发明的另一个方面，提供了一种存储介质，其上存储有计算机程序，其特征在于，该程序被处理器执行时实现如上述方法的步骤。

与现有技术相比，上述方案中的一个或多个实施例可以具有如下优点或有益效果：

1)本发明所分析的地震折射子波易与其他波区分且避免了近道干扰，使得能够精确地提取地震折射子波，便捷可靠地拾取初至，提高了后续地质体结构反演的可靠度；

2)本发明通过对地震折射子波在多个地质网格中的传播过程进行分析，由于地震折射子波同轴能量较强且可追踪性好，能够对一些地下薄层、裂缝等小尺度地质结构进行速度成像，便捷可靠地识别小尺度地质结构，进一步为地质体结构的细节提供解释；

3)本发明的地震折射子波在多个地质网格中的传播过程具备广角特点，能够进行远偏移距处信息的采集与成像，为地质体结构的反演提供更完备的解释。

本发明的其它特征和优点将在随后的说明书中阐述，并且部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点可通过在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

附图说明

附图用来提供对本发明的进一步理解，并且构成说明书的一部分，与本发明的实施例共同用于解释本发明，并不构成对本发明的限制。在附图中：

图1示出了本发明实施例的基于折射波分析地质体结构的方法的流程图；

图2示出了本发明实施例一的基于折射波分析地质体结构的方法的具体步骤流程图；

图3示出了本发明实施例二的基于折射波分析地质体结构的方法的另一具体步骤流程图；

图4示出了本发明实施例四中的试验工区的数据形态；

图5示出了本发明实施例四中的走时层析成像建立的层状初始速度模型，

附图中的标号为：1-折射波，2-反射波，3-水波。

具体实施方式

以下将结合附图及实施例来详细说明本发明的实施方式，借此对本发明如何应用技术手段来解决技术问题，并达成技术效果的实现过程能充分理解并据以实施。需要说明的是，只要不构成冲突，本发明中的各个实施例以及各实施例中的各个特征可以相互结合，所形成的技术方案均在本发明的保护范围之内。

本发明根据走时层析成像原理获取各地质网格点的地下传播速度，走时层析成像原理是以走时为出发点来建立地下介质速度结构的过程。速度的反演过程属于非线性反演问题，采用迭代的方法进行求解，可以将非线性问题线性化。转换为线性化的方法后，首先要给出一个初始速度模型，由最短路径射线追踪法可计算出路径长度及理论上的初至走时。其次，按照选定的正则化准则由走时的差值求得速度模型的残差。最后，由速度的残差，逐步更新速度结构模型直至残差在可接受的足够范围内，即得到了真实地下速度模型，也就得到各地质网格点的地下传播速度。并根据各地质网格点的地下传播速度和Q层析成像原理得到各地质网格点的Q值。Q层析成像原理需要进行二维或三维的射线追踪。Q层析成像原理分为两步：第一步，利用有关波的波形及振幅谱比值求出与Q相关的旅行时t*。第二步，通过第一步求得的t*值进行Q值反演。反演需要将模型进行网格化，在取得速度模型后，每一个网格点处v的值均对应一个待求取的Q值。

具体地，本发明采用频率域的谱比法对折射波射线传播区域进行Q值估计和层析成像。

为解决现有技术中存在的上述技术问题，本发明实施例提供了基于折射波分析地质体结构的方法。

参照图2，本实施例的基于折射波分析地质体结构的方法包括如下步骤：

S100，将目标储层的地质体剖面图区隔成多个地质网格，以形成多个地质网格点；

S200，设置各地质网格点处地震折射波的地下传播速度的初始值，以及构建地震折射波从震源到各地质网格点的射线路径，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值；

具体地，所述步骤S200包括以下步骤：

S210，设置各地质网格点处地震折射波的地下传播速度的初始值；

S220，构建地震折射波从震源到各地质网格点的射线路径，根据地下传播速度的初始值确定地震折射波沿射线路径传播的理论初至旅行时；

S230，基于地下速度反演模型，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值。

S300，基于地质品质因子反演模型，根据各地质网格点处地震折射波的地下传播速度的修正值和各射线路径的地质吸收特征时间，确定各地质网格点的地质品质因子；

S400，根据各地质网格点的地质品质因子的分布分析目标储层的地质体结构。

为解决现有技术中存在的上述技术问题，本发明实施例还提供了基于折射波分析地质体结构的方法的另一具体实施方式。

参照图3，本实施例的基于折射波分析地质体结构的方法包括如下步骤：

S100，将目标储层的地质体剖面图区隔成多个地质网格，以形成多个地质网格点；

在这里，地质体剖面图可以根据目标储层绘制，目标储层为待分析的地质体。地质体剖面图是通过均匀的网格进行参数化的，网格是各个方向上恒定大小的。

S200，设置各地质网格点处地震折射波的地下传播速度的初始值，以及构建地震折射波从震源到各地质网格点的射线路径，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值；

具体地，所述步骤S200包括以下步骤：

S210，设置各地质网格点处地震折射波的地下传播速度的初始值；

具体地，获取地震折射波在目标储层中传播的初始速度模型并从初始速度模型中提取各地质网格点处地震折射波的地下传播速度的初始值；该初始速度模型为均匀变化梯度的二维速度分布，以设置各地质网格点处地震折射波的地下传播速度的初始值，使其符合均匀变化梯度的二维速度分布。初始速度模型是根据经验数据得到的，与地质体的属性有关。例如，可以根据各地质网格点处地震折射波的实际初至旅行时初步估计，在本发明中，根据测得的海底地形深度、采集范围等已知参数给出均匀变化梯度的二维速度分布。

S220，构建地震折射波从震源到各地质网格点的射线路径，根据地下传播速度的初始值确定地震折射波沿射线路径传播的理论初至旅行时；

具体地，通过射线追踪法构建地震折射波从震源到各地质网格点的射线路径。特别地，旅行时和射线路径的正演过程采用的是Vidale(1988,1990)的射线追踪正演方法，在此不做赘述。

S230，基于地下速度反演模型，根据地震折射波从震源到各地质网格点的实际初至旅行时与地震折射波沿射线路径传播的理论初至旅行时之间的差值对地下传播速度的初始值进行修正，以获得各地质网格点处地震折射波的地下传播速度的修正值。

具体地，地下速度反演模型如下：

Bδv＝δt

其中，B为地下速度反演系数矩阵，δv为各地质网格点的地下传播速度与初始速度的差值，δt为各地质网格点的地震折射波的实际初至旅行时与理论初至旅行时的差值。在这里，B作为地下速度反演系数矩阵，适用于通过δt列向量求取δv的具体矩阵，当δt和δv均为n×1的列向量时，则B为n×n的矩阵。

具体地，在未叠加的单炮地震记录剖面上拾取地质网格点处地震折射波的实际初至旅行时。在这里，基于苏门答腊地区的尼亚斯岛采集的人工源广角地震数据，进行初至折射波旅行时的拾取，选择合适宽度的时间窗将需要的波形完整的截取出来，可以通过人工手动进行，这一步骤截取数据的质量决定了后续反演步骤的可靠度。

S240，判断对应于当前地下传播速度的修正值的理论初至旅行时与实际初至旅行时之间的差值是否小于预设的阈值：

若是，执行步骤S300；

若否，令地下传播速度的初始值等于当前地下传播速度的修正值，返回步骤S220，以对各地质网格点处地震折射波的地下传播速度再次修正。

在这里，当前地下传播速度的修正值为步骤S230中的各地质网格点处地震折射波的地下传播速度的修正值。

以下举例说明各地质网格点处地震折射波的地下传播速度的多次修正过程：

第一步，设置地震折射波从震源到各地质网格点的实际初至旅行时为T，首次构建的射线路径为s0，地下传播速度的初始值为v0；

第二步，根据首次构建的射线路径为s0和地下传播速度的初始值v0正演地震折射波沿s0传播的理论初至旅行时t0，再根据地下传播速度的初始值v0、理论初至旅行时t0、实际初至旅行时T获得地下传播速度的第一次修正值v1；

第三步，判断根据射线路径s0和地下传播速度的第一次修正值v1正演的地震折射波沿s0传播的理论初至旅行时与实际初至旅行时T之间的差值是否满足精度；如果满足，则地下传播速度的第一次修正值v1为地下传播速度的最终修正值，执行步骤S300；如果不满足，则返回步骤S200，重建射线路径s1，根据射线路径s1和地下传播速度的第一次修正值v1正演地震折射波沿s1传播的理论初至旅行时t1，再根据地下传播速度的初始值v1、理论初至旅行时t1、实际初至旅行时为T获得地下传播速度的第二次修正值v2；

第四步，根据射线路径s1和地下传播速度的第二次修正值v2正演的地震折射波沿s1传播的理论初至旅行时与实际初至旅行时T之间的差值是否满足精度；如果满足，则地下传播速度的第二次修正值v2为地下传播速度的最终修正值；如果不满足，则重建射线路径s2，根据射线路径s2和地下传播速度的第二次修正值v2正演地震折射波沿s2传播的理论初至旅行时t2，再根据地下传播速度的初始值v2、理论初至旅行时t2、实际初至旅行时为T获得地下传播速度的第三次修正值v3；

依次类推，直到当前地下传播速度的修正值对应的理论初至旅行时与实际初至旅行时之间的差值满足精度为止。

其中，随着多次迭代δt值不断减小，地下传播速度得到更新，直至δt到达规定的精度范围内，例如δt与实际初至旅行时的比值小于0.1％，得到的地下传播速度即最后的反演结果。

此外，在本实施例中，由理论的旅行时和拾取的实际旅行时进行反演迭代，当迭代完成后，可以然后通过GMT(the Generic Mapping tools，通用绘图工具)绘图脚本进行绘制，得到地层速度结构的示意图。

此外，在步骤S100至S240中，优选采用折射波初至走时层析成像法获得本地区的地层速度结构和射线路径。

S300，基于地质品质因子反演模型，根据各地质网格点处地震折射波的地下传播速度的修正值和各射线路径的地质吸收特征时间，确定各地质网格点的地质品质因子；

具体地，各射线路径的地质吸收特征时间t*的求取采用Wilcock(1995)提出的方法。在时间域，地震信号可以以褶积的形式表示出来：

x(t)＝s(t)*c(t)*i(t)

其中，t是时间，s(t)是震源信号，c(t)是地层响应，i(t)是仪器响应。

相应地，在频率域，地层响应的频率谱C可以通过地震信号的频谱X除以震源信号的频谱S和仪器响应的频谱I得到：

其中，C(f)为地层响应频谱，X(f)为接收地震信号频谱，S(f)为震源信号频谱，I(f)为仪器响应频谱，f是频率。

假设用于频谱测量的时间窗只有一个相位并且地震折射子波的传播分量中唯一与频率有关的只有衰减，那么C可以表示为：

C(f,s)＝G(s)exp[-2πft(s)]

其中，s是波的路径，G是波的散度，为一个常数。

t*值可以通过画出ln[C(f)]在频率域的曲线再进行线性化来获得。

其中，t

因此，通过以下表达式获取地质吸收特征时间：

其中，t

具体地，Q层析成像原理即t

其中，Q(x,z)为地质网格点的地质品质因子，v(x,z)为地质网格点的地下传播速度。

在速度v(x，z)已知的情况下，运用走时t就可以进行Q(x，z)的反演。

这个求取Q的过程可以用下面的表达式表示：

其中，A是系数矩阵是和射线路径与网格交点有关的量；λ为约束参量；L是拉普拉斯算子；t*是求取的t*的列向量。

进行离散后的地质品质因子反演模型如下：

其中，Q(x,z)为地质网格点的地质品质因子，v

以上为基于所有射线对应的t*值与Q值模型之间的线性关系建立的线性方程组，利用最小二乘方法求解Q值。

S400，根据各地质网格点的地质品质因子的分布分析目标储层的地质体结构。

本实施例二，通过迭代方式控制当前地下传播速度的修正值对应的理论初至旅行时与实际初至旅行时之间的差值的精度来确保当前地下传播速度的修正值的精确度。

为解决现有技术中存在的上述技术问题，本发明实施例三还提供了一种存储介质，其上存储有计算机程序，该程序被处理器执行时实现如上述方法的步骤。

需要说明的是，存储介质上存储的程序被执行实现的方法与实施例一中基于折射波分析地质体结构的方法相同，因此本实施例存储介质上存储的程序被执行实现的方法在此也不再对其进行赘述。

下面结合具体应用实例来进一步说明利用本发明的基于折射波分析地质体结构的方法对X轴范围为0-120km，Z轴范围为0-20km的地质体进行构造分析，以显示沉积层及其下部地壳存在的Q值的纵向与横向变化。图4示出了本发明实施例四中的试验工区的数据形态，其中可以清楚看到该数据的折射波1、反射波2以及直达水波3的形态。由图4可知，本发明所分析的地震折射子波易与其他波区分且避免了近道干扰，使得能够精确地提取地震折射子波，便捷可靠地拾取初至，提高了后续地质体结构反演的可靠度，并且，本发明的地震折射子波在多个地质网格中的传播过程具备广角特点，能够进行远偏移距处信息的采集与成像，为地质体结构的反演提供更完备的解释。图5示出了本发明实施例四中的走时层析成像建立的层状初始速度模型，图5中的箭头表示层状初始速度由小变大的方向。

首先，采用折射波初至走时层析成像法获得选定试验工区地震折射波的初始速度模型和射线路径；其次，利用折射波初至和水层直达波(即震源信号频谱)的谱比值斜率获得与Q值有关的t*值；最后基于所有射线对应的t*值与Q值模型之间的线性关系，建立线性方程组，利用最小二乘方法求解Q值模型。

Q层析成像结果显示在5-12km深度，沉积层及其下部地壳存在的Q值的纵向与横向变化。从X＝0km到X＝120km的范围内，Q

虽然本发明所公开的实施方式如上，但所述的内容只是为了便于理解本发明而采用的实施方式，并非用以限定本发明。任何本发明所属技术领域内的技术人员，在不脱离本发明所公开的精神和范围的前提下，可以在实施的形式上及细节上作任何的修改与变化，但本发明的保护范围，仍须以所附的权利要求书所界定的范围为准。