地震波阻抗反演方法综述

地震波阻抗反演方法综述

一、地震反演技术研究现状

地震反演方法是一门综合运用数学、物理、计算机科学等学科发展起来的新技术新方法，每当数学方法、物理理论有了新的认识和发展时，就会有新的地震反演技术、方法的提出。随着计算机技术的不断发展、硬件设施的不断升级，这些方法技术得到了实践验证和提升，反过来地震反演技术运用中出现的新问题、新思路又不断促使数学方法、地球物理学理论的再次发展。时至今日，地震反演技术仍然是一个不断发展、不断成熟、不断丰富着的领域。

反演是正演的逆过程，在地震勘探中正演是已知地下的地质构造情况、岩性物性分布情况，根据地震波传播规律和适当的数学计算方法模拟地震波在地下传播以及接收地震波传输到地表信息的过程。地球物理反演就是使用已知的地震波传播规律和计算方法，将地表接收到的地震数据通过逆向运算，预测地下构造情况、岩性物性分布情况的过程。地震波阻抗正演是对反演的理论基础和实现手段。

1959年美国人Edwin Laurentine Drake在宾夕法尼亚州开凿的第一口钻井揭开了世界石油工业的序幕。从刚开始的查看地质露头、寻找构造高点寻找石油，到通过地震剖面的亮点技术寻找石油，再到现在运用多种科学技术手段进行油气资源的预测，石油勘探经历了一个飞速的发展历程。

声波阻抗（AI）是介质密度和波在介质中传播速度的乘积，它能够反映地下地质的岩性信息。声波阻抗反演技术是20世纪70年代加拿大Roy Lindseth博士提出的，通过反演能够将反映地层界面信息的地震数据变为反映岩性变化的波阻抗（或速度）信息。由于波阻抗与地下岩石的密度、速度等信息紧密联系，又可以直接与已知地质、钻井测井信息对比，因此广泛应用于储层的预测和油藏描述中，深受石油工作者的喜爱。70年代后期，从地震道提取声波资料的合成声波技术得到了快速发展，以此为基础发展的基于模型的一维有井波阻抗反演技术，提高了反演结果的可靠性。进入80年代，Cooke等人将数学中的广义线性方法运用于地震资料反演，提出了广义线性地震反演。此后Seymour等人又提出了测井声波资料和地震数据正反演相结合求取地下声波阻抗的测井约束反演，大大拓宽了反演结果的纵向分辨能力。

90年代，在基于前人对地质统计学研究的基础上Bortoli和Haas提出了地质统计学反演，Dubrule等人对该方法进行了改进和推广。在国内随着油田对地震反演技术的广泛应用，以周竹生为主提出的地震、地质和测井资料联合反演方法，将地质信息引入地震反演中，提高的反演结果与地质认识的联系，克服了线性反演存在的缺陷。1996年，李宏兵等人将宽频带约束方法应用于递推反演并对其进行改进，减弱了噪音对反演结果的影响。

1999年，任职于英国石油公司的Connolly在《弹性波阻抗》一文中介绍了弹性波阻抗（EI）的概念和计算方法，阐述了不同入射角度（偏移距）地震道集部分叠加反演波阻抗随入射角之间的关系，但是该方法求取的弹性阻抗随入射角变化很大，无法与常规叠后反演波阻抗直接比较，因此推广应用较为困难。2002年，Whitcombe通过修正Patrick Connolly的计算公式，得到了弹性波阻抗的归一化求取方法，消除了弹性阻抗随入射角变化大的难题。2003年，西北大学马劲风教授从Zoeppritz方程简化出发提出了广义弹性波阻抗的概念，克服了以往波阻抗反演要求地震波垂直入射到地表的假设条件，推导出了任意入射角下纵波反射系数的递推公式，提高了中等入射角度下弹性波阻抗反演的精度。

二、地震反演基本原理

正问题和反问题就像物理学中的作用力和反作用力一样是同时存在的。地球物理学中的正问题一般是已知我们容易直接测量或者获取的物理量，通过两者内在关系求解未知的不容易获取物理量问题，反问题是已知我们不容易获取的物理量求解问题应该满足的条件或物理量。求解正问题的过程称为正演，求解反问题的过程称为反演。一般情况下正演相对简单，求解结果唯一，而反演的求解大多是欠定的，求解过程较为复杂，其结果也可能出现多解性，因此减少反演结果多解性是反演问题的一个重要研究内容。

在油气地震勘探中，当地下地质体的岩性、物性发生变化时，经常会引起岩石密度和地震波传播速度的变化，导致波阻抗产生变化，波阻抗在岩性、物性界面处的变化会生成一个反射系数脉冲，当地震波传播到波阻抗界面时，就会在其界面处产生反射波。因此，若已知地下介质的波阻抗，通过求取反射系数与并地震子波褶积运算得到地震波通过地下介质后观测记录的过程就是地震波阻抗正演。当我们根据观测到的地震记录，通过地震波反射的褶积理论等求取反映地下介质的岩性、物性信息的波阻抗等信息的过程就是地震波阻抗反演。广义的波阻抗反演包括地震资料的处理、解释以及波阻抗求取，一般讨论的波阻抗反演是在地震资料处理和解释的基础上通过反演方法对波阻抗的求取。

地震资料是反射界面的表现，是一种界面型信息，只能间接表达地下地层的地质特征，不能直接进行储层的描述。而作为储层预测重要工具的地震反演技术可以将地震数据转换成波阻抗数据，波阻抗数据由于是地层型剖面，可以和测井地质信息直接对比，可以直观的进行储层的识别和描述。

地震反演方法有多种分类方法，主要分类方法有按使用的地震资料分类、按反演结果分类、按测井在地震反演中的作用大小来分和按反演方法的实现方式来分等。表2-1是四种分类方法对应的主要反演方法。

2.1 直接反演

直接反演方法是在反褶积理论的基础上，通过对地震道进行数学计算求取地下波阻抗体的反演方法。假设地震记录褶积模型为

()()*()()S t W t R t N t =+ (2-1)

式中，()S t 为地震记录，()W t 为地下介质中传播的地震子波，()R t 为地下波阻抗界面的反射系数，()N t 为采集等原因产生的噪音，一般认为是白化的。这就是地震记录的褶积模型，由模型可以看出，地震记录由地下波阻抗界面反射系数与地震子波的褶积加上一定的白化噪音构成。

地震道直接反演是根据地震波传播理论，用数学计算工具消除地震记录中的子波()W t 、弱化白化噪音()N t ，反射界面的反射系数序列()R t ，再通过积分或递推等方法求取地下波阻抗。道积分反演和递推反演是最常见的直接反演方法。

2.1.1 道积分

道积分反演是最简单的波阻抗反演方法，该方法在地下地层波阻抗是深度变量连续且可微函数的假设条件下，通过地震道自相关方法求取地层反射系数，再对反射系数进行深度积分计算出地层波阻抗。由于道积分反演实际上只是对反射系数的积分，反演波阻抗数据体是相对值。

地震波垂直入射到波阻抗界面时，反射系数的表达式为

111111i i i i i i i i i i i i i

v v AI AI R v v AI AI ρρρρ++++++--==++ (2-2) 式中i R 为第i 、1i +层波阻抗界面的反射系数，i ρ、1i ρ+和i v 、1i v +分别是第i 、1i +层的密度和速度，i AI 、1i AI +为第i 、1i +层的波阻抗。根据地层波阻抗是深度连续可微函数的假设，可以将地层看成是很多很薄的地层，因此相邻地层的波阻抗差异不会太大1i i AI AI +≈，所以1i i AI AI AI +-≈?、12i i AI AI AI ++≈，公式（2-2）可以改写成

2AI R AI

?≈ (2-3) 对反射系数在深度t 上进行积分，可得地层相对波阻抗()AI t

001[ln ()ln (0)]22t

t

AI Rdt dt AI t AI AI ?≈=-?? (2-4) 02()(0)t

Rdt AI t AI e ?= (2-5)

式中R 为地层反射系数，t 为地层深度，()AI t 为地层波阻抗。

道积分反演通过积分方法把地震道记录直接转换为波阻抗数据，反演过程具有积分误差小、计算简单、对计算机要求较低的优点，但是道积分反演过程仅仅依赖于地震资料，地震资料的品质和带宽对反演结果影响很大。因为地震资料由于采集、处理等原因，缺失低频和高频信息，频带较窄，所以该方法求取的波阻抗纵向分辨率低、精度小。道积分反演方法无法使用测井等资料约束，因此反演结果是相对值、纵向分辨率不高，仅适用于地震勘探初期以及没有井数据的研究区。由于反演结果不能反映地层的真实波阻抗，所以道积分反演不能定量计算地层的岩性、物性参数。

2.1.2 递推反演

递推反演是在地下波阻抗界面的反射系数为稀疏分布的假设条件下，首先利用实际地震记录估算波阻抗界面的反射系数序列，然后用递推方法求取地层波阻抗的反演方法。地下介质波阻抗的递推公式推导如下：

公式(2-2)中当第0层反射系数0R 和波阻抗(0)AI 已知时，可以导出第1层的波阻抗(1)AI 为

00

1(1)(0)1R AI AI R +=- (2-6) 当第0层反射系数0R 、第1层波阻抗(1)AI 已知时，可推导出第0层(0)AI 为

00

1(0)(1)1R AI AI R -=+ (2-7) 当已知第(0)m m n ≤<层波阻抗()AI m 和m 到n 层间波阻抗界面的反射系数序列(0,1,2,...,1)i R i n =-时，第n 层地层波阻抗()AI n 为

()(1)AI n AI n =-·1111n n R R --+-=(2)AI n -·1111n n R R --+-·22

11n n R R --+- ...()AI m ==·

1111n n R R --+-·2211n n R R --+-·...·1111m m R R --+-·11m m R R +- 1

1()1n i i m i R AI m R -=+=∏- (0m n ≤<) (2-8) 当已知第n 层波阻抗()AI n 和m 到n 层间波阻抗界面的反射系数(0,1,2,...,1)i R i n =-时，第(0)m m n ≤<层地层波阻抗()AI m 为

()(1)AI m AI m =+·11m m R R -+=(2)AI m +·11m m R R -+·11

11m m R R ++-+ ...()AI n ==·

11m m R R -+·1111m m R R ++-+·...·2211n n R R ---+·1111n n R R ---+ 1

1()1n i i m i R AI n R -=-=∏+ (0m n ≤<) (2-9) 通过公式（2-8）和（2-9）可知，当已知第k 层地层波阻抗()AI k 和该层与待求层第n 层的间各层界面的反射系数序列时，可以递推出k 和n 层间各层的波阻抗值。实际应用中一般选取大套泥岩层作为标准层求取其波阻抗，然后根据反射系数序列递推出所有地层的波阻抗值。

递推反演中最重要的部分是反射系数序列的求取，反射系数的稀疏程度对反演结果影响很大，求取合适的反射系数序列关乎反演的成败。

稀疏脉冲反褶积是实际生产中常用的求取递推反演反射系数的方法。图2-1是稀疏脉冲反褶积求取反射系数序列流程图。首先采用最大似然反褶积等方法估算波阻抗界面的反射系数序列，并经褶积求取合成地震记录与实际地震资料的差异反馈给反射系数序列求取过程，根据反馈结果适当修正反射系数序列，当合成地震记录与实际地震资料的残差满足预设条件时，反射系数序列既是波阻抗界面的稀疏反射系数，然后利用公式（2-8）和（2-9）就可以计算出波阻抗。

图2-1 稀疏脉冲反褶积反射系数求取流程图

递推反演结果能够反映地下波阻抗的分布规律，而且递推算法使得反演结果不存在多解现象；递推反演采用简单的递推计算方法，计算速度快、对计算机要求较低；递推反演纵向分辨率高于道积分反演，反演结果是绝对波阻抗，能够用于储层预测与和砂体雕刻。但是递推反演采用的递推计算可以导致波阻抗误差的层层累计，反演结果受地震资料固有频带宽度的影响，对薄储层的识别能力有限，不能满足薄储层的识别和分析，因此主要用于勘探开发前期。

2.2 模型反演

基于模型的反演是充分利用地震、地质、测井等资料建立研究区的宽频带波阻抗模型，根据地震波传播原理在该模型上进行正演模拟求取合成地震记录，对比其与实际地震资料差异并将两者残差反馈到模型上，根据残差不断修正波阻抗模型再次正演，直到正演模拟结果与实际地震资料匹配较好、残差满足预定要求时，正演模型就是我们要求的反演结果。基于模型的反演是一个不断正演、不断修改模型再次正演的过程，图2-2是模型反演的流程图。

图2-2 模型反演流程图

根据模型反演的思路，适当的初始波阻抗模型和子波以及合理的残差是模型反演的关键，它们决定着反演结果的模型化程度和模型修改次数。因此，基于模型的反演中要注意以下几个方面：

（1）波阻抗模型的建立，模型反演是在初始波阻抗模型的基础上进行正演模拟计算，模型对反演结果影响很大，因此建立一个好的波阻抗模型是进行模型反演的基础。在建立波阻抗模型之前，一般需要对测井数据进行标准化处理，对地质信息（构造发育史、地层断层接触关系、地层岩性信息）有一个充分的理解，搭建合适的地层框架，在建模时需要选择合适的内插函数。

（2）地震子波的求取，进行高精度的地震记录标定是求取子波的关键，只有选取合适的地震子波（包括频率和相位），正演结果和实际地震资料才能匹配。

（3）残差的大小，合成记录与地震资料的残差决定了模型修改次数，只有设定一个合理的残差范围才能使反演结果真实反映地下波阻抗情况，又不会进行过多的迭代次数。

基于模型的反演相比于地震道直接反演方法有很多优点：首先，它不需要假设波阻抗反射截面的反射系数为深度的连续可微函数或者稀疏分布，使得反演结果上更能反映地下真实信息；其次，地震子波可以通过精细层位标定求取，使得子波更接近真实情况；该反演方法充分利用了地震、地质和测井等资料建立地质模型，一定程度上约束了反演结果，使反演误差不随深度累计；最后，模型反演使用测井资料建模，一定程度上拓宽了反演结果频带宽度，提高了纵向分辨率。模型反演方法也有其致命的缺点：由于反演是在模型基础上进行的正演计算，反演结果的模型化严重；不同的地下地质情况可能产生相似的地震记录，反演过程对地质模型有很高的依赖性，使得反演结果具有多解性。

2.3 约束稀疏脉冲反演

约束稀疏脉冲反演是波阻抗模型约束下的地震道递推反演方法，它与稀疏脉冲反演的区别就是初始波阻抗模型和测井资料的约束，约束稀疏脉冲反演既能像递推反演一样快速的求解出反演结果，又不会使递推的误差不随深度累积。

约束稀疏脉冲反演的实现方法是将测井资料中声波和密度曲线生成的波阻抗通过地质框架内插得到一个初始波阻抗模型，然后用波阻抗模型以及测井曲线的趋势来约束反演结果，求取地下地层反射系数序列。约束稀疏脉冲反演得到的波阻抗数据体既与地震资料匹配又与测井数据吻合，更能反映地下地质体的真实分布。

约束稀疏脉冲反演的关键是求取稀疏反射系数序列，在满足合成地震记录与实际地震资料精度的前提下，最少脉冲数目反射系数序列的目标函数为

22()()p q q i i i i i J r d s t z λα=+-+-∑∑∑ (2-10)

式中J 是目标函数，i r 是估算的反射系数序列，i d 为原始地震资料，i s 为合成地震记录，λ是合成地震记录与实际地震资料残差的权重系数，i t 是测井波阻抗变化趋势，i z 是反演选取的以井上波阻抗变化趋势为中心的硬约束范围内的反演值，α为测井波阻抗趋势对反演结果的约束程度。一般情况下取1α=、1p =、2q =。

通过式（2-10）可知，约束稀疏脉冲反演最优化目标函数J 由三部分组成，第一项p i r ∑是反射系数序列i r 的p 次方之和，第二项()q q i i d s λ

-∑是反演合成地震记录i s 与实际地质资料i d 残差q 次方之和，第三项22()i i t z α-∑是反演波阻抗与测井趋势均方差之和。式中前两项是互相制约、互相影响的，两者对反演结果的影响依靠λ值大小调控：λ值太大时，过分强调合成地震记录与实际地震资料的残差小，反演结果会引入地震噪声、忽略了反射系数序列的稀疏性；λ值太小时，又过分强调反射系数序列的稀疏性、反射系数个数太少，使得反演结果缺失细节、降低反演结果分辨率，此时合成地震记录与地震资料的吻合性太低，残差太大。因此为了使反演结果分辨率较高，又不引入太多噪音，需要正确选择λ值的大小。 测井曲线以及初始波阻抗模型的约束使得约束稀疏脉冲反演结果的频带宽度较地震资料有所提高，但是地震资料在反演中还是起着主要作用，反演结果的频带没有太大提高，而且容易产生窗帘效应。为了提高反演结果的可靠性及分辨率，需要将地质、钻井测井和地震资料建立的波阻抗模型的低频、高频成分通过滤波补偿到反演结果中，这样才能获得相对宽频带的反演波阻抗体。

约束稀疏脉冲反演结合了直接反演和模型反演的优点，有很广泛的应用，在勘探期可用于储层预测，指导井位部署；在开发期可预测储层的横向变化规律，指导地质统计学反演参数的选择、反演结果的优选。

2.4 地质统计学反演

不管是地震道直接反演还是模型反演，反演过程都依赖于实际地震资料，反演结果受地

震资料品质影响很大。由于地震资料频带宽度大部分都只有十几~几十赫兹，因此反演结果的纵向分辨率较低，对薄储层、薄互层识别能力有限，不能直接应用到精细储层识别的工作中。

以地质统计学和随机模拟为基础的地质统计学反演，将研究区井上资料分析的地质统计学特征应用到随机模拟过程生成井间随机波阻抗数据体，分析合成地震记录与实际地震资料的差异，反过来约束随机模拟过程再次模拟，直至模拟结果满足预设条件，模拟结果就是反演结果。地质统计学反演的优点是利用了地震数据连续性和井上资料纵向高分辨率，因此反演结果能用于薄储层的储层预测。

地质统计学反演由于地质统计学的统计规律和随机模拟的随机性，使得反演结果具有不确定性和多解性，因此一般地质统计学反演会求出多个等概率的数据体，根据前期地质认识和确定性反演结果分析，优选最符合实际情况的结果。

2.4.1 地质统计学

地质统计学（Geostatistics ）是法国G . Matheron 教授通过大量理论实验研究提出的一种统计学方法。地质统计学通过分析样本数据的变差函数和概率密度函数，研究区域化变量的空间分布规律，地质统计学考虑了传统统计学中样本值，同时又考虑了样本间的空间关系。

2.4.1.1 区域化变量

以往我们认为空间中点x 是坐标(,,)u v w x x x 的函数(,,)()u v w x x x z x =，是一个确定的值，只是空间坐标(,,)u v w x x x 的函数，但是在现在我们认为区域化变量是随机变化的，只有进行一次观测我们才能得到它的一个实现()z x ，在进行观测之前我们不能通过空间函数来确定。因此，区域化变量既是随机函数又是空间场的函数。区域化变量的随机性在资源勘查中能够有效的反映地质体的空间局限性、连续性和各向异性。

区域化变量的特征可以用方差函数、协方差函数以及变差函数等进行描述。

2.4.1.2 协方差函数

在地质统计学中，可以用协方差来描述区域化变量的差异。随机向量X 与X 的协方差为

(,)[()()]Cov X Y E X EX Y EY =-- (2-11)

区域化变量()(,,)u v w z x x x x =的（自）协方差可是变量()z x 在空间x 和x h +两点处()z x 和()z x h +的二阶混合中心矩，即

((),())[()()][()][()]Cov z x z x h E z x z x h E z x E z x h +=+-+ (2-12)

由公式可知，区域化变量()z x 的协方差是空间点x 和偏移量h 的函数。假设()i z x 为()z x 在空间位置i x 处的一个观测值（实现），()i z x h +是()z x 在i x +h 的观测值（实现）

[1,2,...,()]i N h =，由式（2-12）可得：

()

11()[()()][()()]()N h i i i i i c h z x z x z x h z x h N h ?

==-+-+∑ (2-13) 式中()i z x 、()]i z x h +为空间x 和x h +两点的平均期望，假设()()i i z x z x h m =+=（m 为常数），()c h ?

可简化为 ()

211()[()()]()N h i i i c h z x z x h m N h ?

==+-∑ (2-14) 2.4.1.3 变差函数

地质统计学分析中，以直方图和变差函数的形式来反映地质变量空间分布的随机性和相关性。

变差函数是一种能定量表示空间相关关系或地质变量连续性的数学方法。它是指地质变量()z x 在x 和x h +两点处的增量的方差之半，即区域变化量在相距为h 的任意两点处的平方差值的一半。其数学表达式为

21()[()()]2

h E z x z x h γ=-+ (2-15)

对于()z x 的观测值（实现），()i z x 是离散的，表达式（2-15）可改写为 ()

21

1()[()()]2()N h i i i h z x z x h N h γ==-+∑ (2-16) 在地震储层反演中，变差函数中有三个重要参数：块金值0c 、基台值c 和变程a 。图2-3为变差函数示意图，当两个点距离较小h 时，变差函数随距离h 增大而增大，当距离h 过大时，变差函数不再增大。变差函数两点间距离增大临界距离称为变程a ，变差函数到达稳定是的值称为基台值c ，两点重合距离为零时的变差函数为块金值0c 。

图2-3 变差函数示意图

从图上不难理解：变程a 代表了研究区内任意两点观测值是否具有相关性的最大距离，也反映了在某个方向储层参数变化的快慢，当距离h a <时，观测值有相关性，当h a ≥时，观测值没有相关性；基台值c 是变差函数达到平稳时对应的值(,)x h γ，反映了储层参数在某个方向变化的幅度。因此，储层参数在不同方向上的不同变化规律以通过求取变差函数来体现。

测井曲线等数据是区域化变量的观测值，由于测量误差等原因导致观测数据产生混乱，为了较好的反映地下储层参数的真实规律、便于定量研究，实际生产中采用理论变差函数模型拟合观测数据。采用变差函数模型既可以用数学工具研究空间规律，又可以求得连续的变差函数。常用的变差函数类型及其表达式如下：

(1)高斯模型

()0010220>=??

???

???? ??-+=-h h e c c h a h γ

(2-17) (2)球状模型 ()a h a h h c c a h a h c c h >≤≤=???

????+???? ??-+=00212300330γ (2-18)

(3)指数模型 ()00100>=????????? ??-+=-h h e c c h a h γ (2-19)

2.4.2 反演过程

通过地质统计分析测井曲线、地震、地质等数据分析，可以获得地下储层的统计特征。在随机模拟过程中，可以根据储层统计特征进行模拟计算。

地质统计学模拟从储层参数开始，以变差函数和概率密度分布为基础，运用随机模拟算法预测储层的过程。现在最常用的随机模拟算法是序贯模拟算法和马尔科夫链-蒙托卡罗模拟算法。序贯模拟需要对地质体进行网格划分，对各个节点进行随机模拟，而且序贯模拟要求概率密度函数呈正态分布，因此模拟结果对模型依赖较为严重，容易出现局部最优化解；马尔科夫链-蒙托卡罗模拟算法可以根据实际的概率密度分布函数得到统计意义上真正的随机样点分布，通过类似于优化算法的增量调整算法实现全局优化。

地质统计学反演是将随机模拟的波阻抗数据体或者其他属性数据体通过相关变换或云变换转换成波阻抗数据体，进而求取反射系数序列合成地震记录，根据合成地震记录与实际地震资料的残差约束随机过程，求取反演数据体。