第1章叠后波阻抗反演技术
1.1地震反演
由地下地质信息得到地震信息的过程，称为正演；反过来，由地震信息得到地下地质信息的过程，称为地震反演。地球物理反演就是从有限频带宽度的地震数据中恢复宽带波阻抗（图1.1），恢复地下地质结构和岩石性质（Sarkheil，2020；撒利明等，2015；Ning，2006）。
通常在地震解释过程中是沿同相轴解释的，同相轴的变化对应着地下地层的变化，也就是说同相轴实际上就是地下岩层的分界面（杨文采，1993）。反演是在解释层位的基础上建立地质模型并反演出波阻抗界面，通过已知的地震资料、测井资料及其他资料，利用反褶积等方法推算出地下地质体的空间分布规律，反演出的波阻抗界面实际上就是地下地质体的反映，其本质是地震资料解释的一部分，可以看成地下地质体的一种地震属性（Brossier et al.，2015；张国栋等，2010）。因此，地震反演技术将地震反射剖面转换为地层波阻抗剖面，可直接进行地质解释和分析，是地震数字处理结果之一，比解释的层位更能反映地下地层的变化，逐渐成为煤层勘探的重要手段（管永伟等，2016）。
煤层预测最重要的是横向和纵向上的分辨率，然而一般地震数据只能提供中频信息（姚逢昌和甘利灯，2000），可得到比较高的横向分辨率，但是对于提高纵向分辨率的低频和高频信息却没有什么手段，而测井资料中含有高低频信息，因此，在反演过程中将测井资料作为约束地震反演的条件，从而达到提高地震分辨率的目的（Li et al.，2011；Li，2004；沈财余和阎向华，1999）。
反演结果的可靠性与地震资料及测井资料有关，一般认为井位越多测井资料越全，对反演越有利。当然井口数量多时，一定程度上能提高反演准确度，但并不是越多越好，关键是能找到测井曲线和地下地质体的正确对应关系（沈财余等，2002）。有些测井资料和地震资料在采集时受到各种影响，并不能真实反映地下地质体的真实特征，因此需要对已知的资料做处理，比如对测井曲线去异常值、标准化、曲线重构等，这样处理后的测井资料更具有横向对比性。最终反演的结果是通过层位建立地质模型，在测井约束的条件下进行反演，然后利用反演出来的波阻抗数据体沿层位提取属性进行综合分析，并进行勘探预测和评价（张超英等，2004；杨立强，2003）。
1.2波阻抗反演原理
地下层位速度的不同对应着波阻抗的不同，速度的差异界面即波阻抗的分界面，而界面的反射系数可以表示为界面上下波阻抗之差与波阻抗之和的比值（Tarantola，1984）。由地震反射波理论可知，地震反射记录是地下反射系数序列与地震子波的褶积。由此，建立了地震反射记录叠加剖面与地震波阻抗剖面的联系（张宏兵等，2005；张永华等，2004）。
在地震反射波勘探中，地下介质的反射界面主要靠地震波的波阻抗界面来反映。地下岩性分界面所产生的岩性差异，会使各地层的速度和密度产生变化，而波阻抗为速度与密度之积，因而岩性差异会导致波阻抗差异。所以利用反演得到的波阻抗数据可以有效识别地下地质体的差异，并进行目的层的预测（Gan and Goulty，2010；Koesoemadinata and Mcmechan，2003）。基于褶积模型的地震波阻抗正演和反演方法的基本原理如图1.2所示。
在正演过程中，设s（t）为最终得到的地震道，vi和vi+1分别是第i层和第i+1层的速度，单位为m/s，ρi和ρi+1分别是第i层和第i+1层的密度，单位为g/cm3，当地震波垂直入射岩性分界面时，反射系数ri为
（1.1）
根据式（1.1），当波阻抗从小变大时，反射系数为正；反之当波阻抗从大变小时，反射系数为负。将反射系数序列r（t）与地震子波w（t）褶积，再加上随机噪声n（t），即为最终得到的地震道s（t）：
（1.2）
反演是正演的逆过程，原理与正演恰好相反。反演通过一系列数学算法，消除随机噪声n（t）和地震子波w（t）对地震道s（t）的影响，得到反射系数序列r（t），再将反射系数序列代入递推公式（1.3）：
（1.3）
式中，Z0是初始波阻抗；Zj+1是j+1层的波阻抗。如果知道初始波阻抗，就可以依据式（1.3）逐层计算出各层的波阻抗，得到波阻抗数据体。
波阻抗反演的基本思路：首先根据研究区域的测井相关信息建立初始模型，然后对建立的初始模型进行地震正演，求得合成地震记录；将合成地震记录与实际地震记录进行比较，根据比较结果，对地震波阻抗模型的密度、速度、子波及深度进行修改，随后正演求取合成地震记录，再与实际地震记录进行比较，继续修改地震波阻抗模型；进行多次反复修改，直至最终的合成地震记录与实际地震记录达到最为接近状态，这样反复修改所得到的结果，就是波阻抗反演结果（Figueiredo et al.，2017；孙振涛等，2002）。
1.3褶积模型
基于褶积模型的反演是地震反演中的一类重要方法，主要包括递推反演、稀疏脉冲反演和基于模型反演等方法，在煤炭资源勘探中发挥着重要作用（Wang et al.，2021；左博新等，2012）。
在地震勘探中，通常把地震记录面貌的形成过程概括为以下的数学模型（韦瑜等，2017）。假设地震道f（t）是由有效波s（t）和干扰波n（t）叠加组成的，即
（1.4）
此处的有效波是指一次反射波。对反射波地震勘探而言，除一次反射波以外的所有波都是干扰波。层状介质的一次反射波通常用线性褶积模型表示：
（1.5）
式中，s（t）是地震道；w（t）是地震子波；r（t）是反射系数函数；符号“最”表示褶积运算；T是地震记录的双程旅行时，s;τ是地震子波的延续时间，s;t是岩层间双程旅行时，s。
式（1.5）称为地震记录道的时间域褶积模型。根据傅里叶展开式中的褶积定理，式（1.5）在频率域中就是乘积关系，即
（1.6）
式中，S（jω）、W（jω）、R（jω）分别是s（t）、w（t）、r（t）的傅里叶变换；ω=2π/T，是圆频率。
它们的复数形式可分解为振幅谱和相位谱两部分，即
（1.7）
三者之间的振幅谱和相位谱分别为
（1.9）
式（1.8）和式（1.9）表明，地震道的振幅谱是地震子波振幅谱及反射系数振幅谱的乘积，地震道的相位谱是地震子波相位谱和反射系数相位谱之和，若以离散形式表示，则式（1.5）可写成
（1.10）
式中，
（1.11）
式中，si，wi，ri分别是有效波记录、地震子波及反射系数的时间序列，长度各为n+1、m+1及p+1，一般有m （1.12）
由于地震子波为物理可实现序列，故w－k=0；又因为地震子波长度为m，故wk>m=0，因而有
（1.13）
矩阵式（1.13）可写成如下的向量形式：
（1.14）
式中，s是n+1阶地震记录道的列向量；w是（m+1）×（p+1）阶的子波矩阵；r是p+1阶反射系数的列向量。
上述表征地震记录数学关系的一系列表达式，如时间域表达式（1.5）、频率域表达式（1.6）、离散表达式（1.10）、向量表达式（1.14），其表示的物理过程都是等效的，统称为地震记录的褶积模型。
大量事实表明，利用声波测井资料和其他资料换算出反射系数函数r（t），并选用合适的地震子波w（t），按上述褶积模型计算出的人工合成地震记录与对应的井旁地震记录大都符合较好。由此可见，这一套地震记录形成的理论（也称地震记录的褶积模型理论，此外还有波动理论）是基本符合客观实际的，并且是正确合理的。
1.4地震子波
在地震数据处理中，特别是高分辨率处理中，常常涉及子波，但在不同场合下子波的含义并不完全相同（高静怀等，2009）。
一种场合是指原始记录中包含子波，它与震源、大地滤波、接收仪器特性等因素有关。在数据处理过程中，常要以很大的努力来估算子波和消除子波的影响，但解释人员基本上与这种子波不产生联系。
另一种场合是指处理的数据及向解释人员提供的数据中包含子波。高分辨率处理虽然以获得反射系数剖面为理想目标，但受原始数据信号频带的限制及信噪比的限制，处理结果只能是反射系数的带限版本。这就是说处理结果仍然是反射系数与带限子波的褶积，这种子波在一定程度上是可以选择的，与解释工作有直接联系（Rui and Zhi，2018）。
在高分辨率数据的解释和反演中，子波扮演着非常重要的角色，在进行目的层标定合成记录过程中以及反演中都要对子波进行分析和提取，因此很有必要在这里对子波进行详细的研究（戴永寿等，2008）。
在地震反演过程中，子波提取是最复杂的问题之一，主要分为统计子波提取和雷克子波提取两大类。统计子波提取无须测井信息就能得到估计子波，但要对地震资料和反射系数序列的分布进行假设，得到的子波精度与假设条件有关；雷克子波是根据规律计算得到的，进行标定时要确定其主频和相位（Sun and Feng，2010；杨培杰和印兴耀，2008）。
1.4.1地震资料提取统计子波
地震反演系统中的统计子波提取，只是使用地震道来提取子波。该方法不计算相位谱，默认为最小相位。振幅谱是用地震道的自相关计算出来的，对每一个用作子波提取的地震道而言计算步骤如下。
（1）提取分析时窗。
（2）将时窗的起始和结束时间用一个斜坡长度创造斜坡，其长度取10个采样点或1/4个时窗长度二者中的较小长度。
（3）对时窗内的数据计算自相关系数，自相关的长度等于期望子波长度的1/2。
（4）计算自相关的振幅谱。
（5）取自相关振幅谱的均方根，近似为子波的振幅谱。
（6）加期望相位。
（7）采用反快速傅里叶变换（fast Fourier transform，FFT）产生子波。
（8）将该道计算的结果与分析时窗内其他道计算的子波求和。
在上述流程中，子波长度是一个关键的参数，它决定了地震道振幅谱的平滑程度。当子波长度增加时，子波谱就接近于分析时窗内地震数据的谱。
在提取统计子波的过程中，仅仅使用了地震资料，没有考虑测井资料，忽略了子波相位的重要性。在反演流程中，这种子波是最初始的子波，往往是为了进一步提取高精度的子波做基础，很少直接用作反演子波。对于信噪比高的资料，统计子波的提取方法是较为可取的。
1.4.2理论子波
根据地震子波的特点用一些具有特殊数学表达式的波形来表示，如雷克子波等，计算时要确定其主频和相位，对目的层段提取的统计子波进行主频分析，该主频就可以作为雷克子波的主频；对于相位而言，实际的地震记录是混合相位，但是通常都对地震资料做了零相位化处理，所以雷克子波的相位近似选用零相位。雷克子波在时间域和频率域的表达式如下：
（1.15）
（1.16）
式中，fp是雷克子波的主频；t是时间；w（t）是子波；w（f）是振幅谱；θ（f）是相位谱。
1.5基于模型的反演
由于使用的原始资料不同，地震反演可分为叠前反演和叠后反演两大类；由于所用的地震信息不同，地震反演还可以分为地震波振幅反演和地震波旅行时反演；还有一些其他的因素可以影响分类，如地质结果，地震反演可分为波阻抗反演、构造反演、参数反演等（潘新朋等，2018；赵小龙等，2016；Gunning and Glinsky，2007）。
随着地震反演技术的日渐成熟，地震反演的种类随着新反演思路的不断提出以及计算方法的不断变化大幅度增加，不过从最根本的层面也就是实现地震反演的方式来分类，叠后波阻抗反演可以分为以下三种：道积分反演、递推反演及基于模型的反演。为满足煤田精细探测的需要，主要采用基于模型的反演方法（梁光河和顾贤明，1993）。
基于模型的反演是在反褶积反演的基础上发展起来的，其一定程度上解决了褶积反演方法频带宽度有限的问题，同时通过实际的测井数据以及多种地质信息的约束，该方法的精确度和分辨率都得到了提升（贺维胜等，2007；孟宪军等，2005）。
基本技术思路：先建立一
