地震第3章反褶积ppt课件.ppt

Chapter3 反褶积,本章首先介绍反褶积的概念及其在地震数据处理中的作用，然后讨论各种反褶积方法原理及其实现问题。,3. 1 反褶积及褶积模型,一、反褶积的概念,反褶积是地震数据处理中一个基本的处理环节。反褶积的基本作用是压缩地震记录中的地震子波，同时，可以压制鸣震和多次波，因而反褶积可以明显提高地震的垂直分辨率。反褶积通常是用于叠前地震数据处理，也可以用于叠后数据处理，通常在一个地震数据处理流程中，为了提高地震垂直分辨率，在叠前和叠后不止一次用到反褶积处理。,反褶积处理是褶积处理的反过程，因而称为反褶积。在前面第一章第一节中，我们曾讲过，一个滤波器的滤波过程在时间域的输出是输入信号与滤波器滤波因子的褶积。因此，时间域的褶积处理就相当于一个滤波过程。而反褶积则相当于时间域的一个反滤波过程。地震记录可以看作是地震子波与地层脉冲响应的褶积，即,为地震道记录;为地震子波;为地层脉冲响应，为震源是单位脉冲 时零炮检距自激自收的地震记录。,(3-1),(3一1)式可视为一个滤波过程，如图3-1所示。,图3-1 褶积滤波过程,这个滤波过程的输入为地震子波。 滤波器的滤波因子为地层脉冲响应 ，输出为地震道记录 。 或者输入为地层脉冲响应 ，滤波器滤波因子为子波 ，输出为地震道记录城 。,如果我们设计一个滤波器，其滤波因子 具有与滤波器约恰好相反的性质，即当输入为地震道记录 约时，其输出为地层脉冲响应 。我们称这个反过程为反滤波或反褶积，如图3-2所示。,二、褶积模型,反褶积是以地震褶积模型为基础的。图3一3所示为地震褶积过程图。图3一3(a)为一段声波测井曲线，表明该井处的地层层速度随深度而变化 。图3一3(b)表示图3一3(a)得到的反射系数随深度的变化 ，计算时假定平面波垂直入射，井忽略了地层密度的变化。图3一3(c)表示图3一3（b)中的反射系数随深度变化 。利用图3-3(a）的速度信息进行深度一时间转换后，得到的反射系数随双程旅行时变化 。图3一3(d)表示由图3一3(c)的反射系数序列 得到的地层脉冲响应，其中包括一次反射和各种多次波的响应。最后，图3一3(e)为图3-3（d）的地层脉冲响应与图3一4中的震源子波按(3一1）式褶积得到的人工合成地震记录。,(3-1)式所表示的地震道记录 是不包含随机环境噪声的。合成地震记录，则需要在(3-l)式中加入适量的随机环境噪声 。 为了得到更真实的人工，这时，人工合成地震记录褶积公式为,式中 随机环境噪声。 整个人工合成地震记录褶积过程如图3一4所示。而反褶积过程则是与之相反，试图由所得到的地震记录 恢复地层脉冲响应 。 或反射系数序列。 由(3一2)式可以得到,(3-2),(3-3),式中 分别为 和 的频谱。,由于随机噪声 和地层脉冲响应两者均接近白噪声，它们的振幅谱 和 在接近全频带范围内是近似于相对平坦的，因而地震子波的振幅谱 近似于光滑后的地震记录的振幅谱 ,两者的自相关函数 与 也是近似的这一性质对后面将讨论的一些反褶积方法是重要的。,3. 2 反滤波,一、反滤波的概念,在反射波法地震勘探中，由炸药爆炸等震源产生一个尖锐的脉冲，在地层介质中传播，并经反射界面反射后返回地面，其理想的地震记录应该是如图3一5所示的一系列尖脉冲，其中每个脉冲表明地下存在一个反射界面，整个脉冲系列表明了地下一组反射界面。这种理想的地震记录x(t)可表示为,式中。 震源脉冲值，为一常数; 反射界面的反射系数。 但是，由于地层介质具有滤波作用，这种大地的滤波作用相当于一个滤波器。因此，由震源发出的尖脉冲经过大地滤波器的滤波作用后，变成一个具有一定时间延续的波形 ，通常叫作地震子波(图3一6)。这时，地震记录是许多反射波叠加的结果，即地震记录 是地震子波 与反射系数 的褶积,(3-4),实际的地震记录城 除了(3一4)式所表示的一系列反射波 而外，还存在着干扰波 ，因此，地震记录双 的一般模型可以写为,其结果为一复杂的地震记录形式(图3-7)。,(3-5),在普通的地震记录上，一个界面的反射波一般是一个延续时间为几十毫秒的波形。由于地下反射界面一般是相距为几米至几十米的密集层，他们的到达时间差仅为几毫秒到几十毫秒，因此，在反射地震记录上它们彼此干涉，难于区分开来。为了提高反射地震记录的分辨能力，我们希望在所得到的地震记录上，每个界面的反射波表现为一个窄脉冲，每个脉冲的强弱与界面的反射系数的大小成正比，而脉冲的极性反映界面反射系数的符号。那么，怎样把延续几十毫秒的地震子波 ，压缩成为一个反映反射系数 的窄脉冲呢?这就是反褶积所要解决的问题。 如果，地震记录是(3-4)式所表示的地震子波 与反射系数 的褶积，即地震记录中只有反射波 ，而没有干扰波 。这时反褶积问题很简单。 根据(3一4)式，在频率域相应有,(3-6),式中， 分别是地震记录 、地震子波 和反射系数 的频谱。显然,如果令,则得到,(3-7),(3-8),(3-9),在时间域，得到,其中 是 的时间函数。由(3-10)式得到,叫作反子波或逆子波。由此可知，已知地震子波 。求出反子波 利用(3一10)式，将反子波 。与地震记录 褶积，即可求出反射系数,(3-10),(3-11),(3-12),这个过程叫作反褶积(图3一8)。,因而，所谓的反褶积或反滤波实际就是一个滤波过程，只不过是这种滤波过程其作用恰好与某个滤波过程的作用相反。,二、地震子波的求取,在进行反褶积处理时，通常必须知道地震子波 的形状。地震子波求取得是否准确对反褶积结果的影响很大，求取地震子波的方法较多，这里只讲在反褶积处理中常用的几种求取地震子波的方法。,1.直接观测法 这种方法是用专门布置在震源附近的检波器直接记录地震子波 ，此方法只适用于海上地震勘探。 在某些地区的海上地震勘探中，在地震记录上海底反射波到达之前曾记录到一个地震波。经过分析知道这是由于海水含盐量有分层性所形成的。由于海水的含盐量有分层性使海水明显地分成上下两层。下层的含盐量较上层含盐量高，形成了一个较为清楚的界面。由震源出发的地震波到达这个界面引起反射返回到海面下的检波器，被记录下来。由于这个波没有与其他波干涉，所以可以作为地震子波 。使用这样求取的地震子波进行反褶积，得到了良好的效果。,2.自相关法 对某个地震记录道选记录质量高的一段，取时窗起点为时间起点，时窗长度为T，则该段地震记录,其Z变换为,假设反射系数 为白噪声，其z变换为 ，则 自相关 的Z变换,从(3-4)式，可知地震记录 的z 变换,式中 地震子波 的Z变换。地震记录 自相关 的Z变换为,(3-13),(3-14),(3-15),将(3-13)式代入上式，得到,将 代入_上式，得到,由于 和 的系数 。 和 都是实数，因而,其中 和 分别是 和 的共扼复数。由于,(3-16),因而，由(3一1)式得到,所以,其中， 是未知的。现在要确定出 。假设地震子波 是最小相位的。则地震子波 满足因果关系，即其Z变换,(3-17),(3-19),(3-18),地震子波， 还满足稳定性条件，即,地震子波 的Z变换 在单位圆内没有根，即当 时,在上述条件下，对(3一19)式两边取对数，得到,令,(3-20),因而得到,根据复变函数理论，,其中 c 是常数，和,利用(3-21）式求出 后，代入(3一19）式，可以求出 ，再利用,(3-21),(3-22),(3-23),(3-24),其中，M要求取的足够大。即可求出地震子波 。 如果地震子波 不是最小相位的，而是零相位的。假设反射系数 为白噪声，则其振幅谱,地震记录 的振幅谱,这时地震记录的振幅谱 与地震子波的振幅谱 相同。 据此，可以对地震记录以 求自相关,(3-25),(3-26),(3-27),再计算其频谱,由于当反射系数 为白噪声时,即可得到地震子波 的功率谱,由此得到地震子波的振幅谱 。当地震子波 零相位时 ，其频谱,(3-28),(3-29),(3-30),(3-31),上式也可由(3一1）式中，令 ，直接得出。将所得到的地震子波的频谱 进行傅里叶反变换，得到,因而，当地震子波 为零相位时，在对地震记录 的自相关进行傅里叶变换，求出其振幅谱 后，利用（3-32)式，再对 进行傅里叶反变换，即可求出地震子波 。,3.多项式求根法,对某个地震记录道选择地震记录质量较好的一段,(3-32),假设反射系数 白噪声，即,且,其中 E 表示数学期望。则由地震记录,可知其自相关函数,即地震记录 的自相关函数与地震子波 的自相关函数相同,并且,其Z变换为,(3-33),(3-34),将上式两边同乘以 ，得到,对上式进行因式分解，即可求出多项式(3一35)的ZM个根，在这ZM个根中，共有M对互为倒数的根。设在这ZM个根中有M个模大于1的根 。 如果假设地震子波。 是最小相位的，则它的Z变换 的根都在单位圆外，因而，得到最小相位的地震子波的z变换为,因此，（3-34)式可以写成,(3-35),(3-36),则(3-35)式可以写成,令Z=O，得到,因而,即,将上式所得到的 ，代入(3-36)式，得到,(3-37),由此得出最小相位的地震子波,4、利用测井资料求子波的方法,这种方法要求有良好的声波测井和密度测井资料，并且在井旁有质量较高的地震记录。 其方法是根据声波测井和密度测井资料得到声速曲线 和密度曲线 ，因而求出声阻抗曲线 ，把深度H转换成垂直双程旅行时,其中 ，得到随反射时间变化的声阻抗曲线 。然后利用反射系数公式,计算出反射系数 ，从井旁的地震记录得到 。,利用傅里叶变换，求出反射系数 和地震记录 的频谱 和 ，因为地震记录的频谱,得到地震子波的频谱,(3-39),(3-40),(3-41),最后，对 进行傅里叶反变换，就得到地震子波,这种方法不必假设反射系数是白噪声，也不必预先知道地震子波的相位特性。,5. 对数分解法,这种方法也不需要假设反射系数是白噪声和地震子波是最小相位的。 假设地震记录 是地震子波。 与反射系数 褶积的结果，即,在频率域，则有,(3-42),为了将地震子波 和反射系数 从地震记录 中分离开来，对上式等式两端取对数得到,频谱的对数叫作对数谱。上式表明频谱 和 的乘积在对数谱中变成 相加。因为对数谱是频率 的函数， 可以用傅里叶反变换，求出相应的时间序列,对数谱的时间序列叫作对数谱序列。由(3-43)式得到,(3-43),(3-44),(3-45),上式表明，地震记录 的对数谱序列 是地震子波对数谱序列 和反射系数对数谱序列 之和。但是， 与 两者在时间轴上分布的位置是不同的。地震子波对数谱序列 分布在时间轴靠近原点附近，而反射系数对数谱序列 分布在离时间轴原点较远的区域。,因此，对对数谱序列 在时间轴上进行低通滤波。选择时间域低通滤波器的时间特性,其中， 是滤波门限时间，则滤波后的输出,(3-46),(3-47),如果，地震子波的对数谱序列 与反射系数对数谱序列 在时间轴上分离较好，适当选择时间域低通滤波门限时间 ，使它接近于地震子波对数谱序列 在时间轴分布区间的上限，则滤波后的结果将非常接近地震子波对数谱序列 。即,根据滤波后的结果所得到的地震子波对数谱序列 即可计算出地震子波 。为此，首先对 进行傅里叶变换，求出 的对数谱,再对对数谱 取指数，得到 的频谱,(3-48),(3-49),(3-50),然后，再对频谱 进行傅里叶反变换，就可以得到地震子波,但是，由于难于准确确定地震子波对数谱序列 在时间轴上分布的区域，进行时间域低通滤波的门限时间 不易选择得很适当，再加上随机噪声的干扰对对数分解法求取子波的影响很大，实践证明，即使较低的噪声干扰水平(例如信号与噪声的振幅比为10)也可以对这种方法所求取的子波的结果产生较为严重的影响。因此，在用对数分解法求取地震子波时，要注意选择记录质量好的地震记录进行计算。 为此，发展了一种对数谱序列平均法，用这种方法求取地震子波效果较好。,(3-51),取 n 个地震记录道 ，假设各地震记录道上的地震子波是相同的;各道的反射系数是随机分布的;各道的噪声也是随机分布的。对各地震记录道分别计算对数谱序列，得到 。根据(3-45)式，有,由于假定反射系数在各道是随机分布的，所以对各地震道的对数谱序列进行平均时，其平均值,(3-52),中，各地震道的反射系数对数谱序列平均值,因而，得到,上式表明对各地震道的地震记录对数谱序列 进行平均，所得到的平均值即接近于地震子波对数谱序列 。,在平均后的地震记录对数谱序列 时间轴原点附近进行时间域低频滤波，所得到的结果就是比较准确的地震子波对数谱序列 。对这个 进行前面叙述过的一系列计算，即可获得比较准确的地震子波 。,(3-53),(3-54),三、反滤波的实现,在应用上面各种方法获得地震子波 之后，再根据(3-8)式,所表示的反子波 与地震子波 频谱之间的关系求取反子波。求取反子波 时最便利的方法是利用反子波 与地震子波 的Z变换之间的关系,首先，根据地震子波时间序列,得到其Z变换,(3-55),(3-56),然后，利用(3-55)式，求出反子波 的Z变换,从而得到反子波时间序列,将反子波 作为反滤波的滤波因子，与输入的地震记录 褶积，即可得到反射系数序列,(3-57),当地震子波 是最小相位时，其反子波 也是最小相位的。这时，反滤波的滤波因子系数为一收敛序列，反滤波器是稳定的。否则，如果地震子波 是最大相位或混合相位的，则其反滤波的滤波因子 的系数是发散的。这时，反滤波器是不稳定的。,3. 3 最佳维纳滤波及最小平方反褶积,最佳维纳滤波即最小平方滤波。本节应用最佳维纳滤波原理实现最小平方反褶积的目的。,一、最佳维纳滤波,最佳维纳滤波是数字滤波中的一大类滤波方法。它是在滤波器实际输出与期望输出的误差平方和为最小的情况下，确定滤波器的滤波因子的，因而称为最小平方滤波。其基本原理如下:,已知输入为,现在要求设计一个滤波器，其滤波因子为,使得滤波后的实际输出为,与期望输出,在最小平方意义下最接近。即使滤波器的实际输出 与期望输出 的误差平方和,(3-58),为最小，即使滤波因子 满足,由此得出,令,(3-59),其中 是时间延迟为 的输入 的自相关， 是时间延迟为 s 的输入 与期望输出 的互相关。于是，(3-59）式可以写成,或写成矩阵形式 ，即,解方程组(3-60)或解矩阵方程(3-61)，得到滤波因子,(3-61),(3-62),(3-60),其中,是矩阵 的逆矩阵。,（3-61)式左端的输入 的自相关矩阵为一对称矩阵，称为托布里兹(Toeplitz) 矩阵。求解这种矩阵方程的方法很多，其中莱文森(Levision)递推算法可对其快速求解，从而得到最佳维纳滤波器的滤波因子 。将输入 经滤波器 作用后，得到输出为,此输出 与设计维纳滤波器的期望输出 在最小平方意义下最为接近。,(3-63),设计最佳维纳滤波器的关键在于针对输入 的特性，根据滤波处理的目的设定期望输出 。期望输出可以是:(1）零延迟尖脉冲;(2)具有延迟时间的尖脉冲;(3)向前时移输入预测距 任意形状期望输出等。设定不同的期望输出可以得到不同的反褶积的应用。,二、最小平方反褶积,将上述最佳维纳滤波原理应用于反褶积(反滤波)问题，就可以得到最小平方反褶积(反滤波)方法。最小平方反褶积方法原理叙述如下: 已知地震记录,假设地震记录 为,(3-64),现在，要求设计一个滤波器,使地震子波 变成窄脉冲,同时，使干扰 得到最大限度的压制。 和前述一样，在确定滤波器 时，使地震记录 经过滤波器 作用后的实际输出,与期望输出一列窄脉冲,(3-65),(3-66),的误差平方和,为最小。则,由此得出,令,(3-68),(3-67),其中 是时间延迟为 的地震记录自相关 是时间延迟为 s 的地震记录与期望输出的互相关。,于是，(3-68)式可以写成,现在，来计算地震记录 的自相关函数。假设，反射系数 为白噪声，它满足,随机噪声 约与反射系数 不相关，即满足,(3-69),在上述假设条件下，地震记录 的自相关,其中,(3-70),令,再来计算地震记录 与期望输出 的互相关函数。,(3-73),令(3-70)式中的 ,(3 -73)式中的 ，并代入(3-69)式中，得到方程组,将上面的方程组写成矩阵形式，得到,或,(3-74),(3-75),(3-76),如果输入的地震记录 不包含 ，即,此时，(3-71)式中的 ，或 e=0 。由(3-70)式可知,(3-75)式和(3-76)式，可以写成,(3-77),(3-78),上式即为当期望输出 为任意形状时的最小平方反褶积滤波方程。利用这个矩阵方程即可在已知地震子波 情况下，计算出反滤波器的滤波因子 ，经过反滤波后，得到与所选择的任意期望输出 在最小平方意义下最接近的输出结果。因此，利用(3-78)式可以进行子波整形反褶积处理。,3. 4 脉冲反褶积,一、脉冲反褶积原理 在上节的最小平方反褶积中如果期望输出不是具有一定延续时间的波形 ，而是一个尖脉冲,则地震子波 与期望输出 的互相关函数,只有当 j=0 时， 当j为其他值时 ，(3-75)式和(3-76)式变成,或,(3-79),(3-80),解矩阵方程(3-79)或矩阵方程(3-80)，即可得到期望输出为尖脉冲 的反滤波因子 。 求出反滤波因子 之后，对输入地震记录 进行反褶波，即可得到反滤波后的输出,一般在进行反滤波时，对地震子波 有一定的要求，即要求地震子波 的相位是最小延迟的。但是，在进行最小平方反滤波时，可以避免这个条件，对非最小延迟的地震子波找到近似解，其方法是将反滤波因子 延迟到适当的位置。 假设地震子波,的相位性质是任意的，既可以是最小相位或最小延迟，也可以是最大相位或最大延迟和混合相位或混合延迟的。要得到期望输出 为一尖脉冲,所设计的反滤波器的滤波因子 在理论上应为无穷项，但实际上只能取有限项，并且取其主要部分,其中， 表示滤波因子 主要部分起始的时间位置， m 表示滤波因子 主要部分的长度，这里 都是正整数。反滤波因子 的形状和位置要根据地震子波 的形状来决定。 利用Z变换，由(3-57)式可以得出，反滤波因子 的形状和位置与地震子波 之间的关系，如图3-9表示。,图3-9表明，对于不同相位性质的子波，反滤波因子的主要部分的位置和形状是不同的。按照地震子波的类型可以将反滤波分为三种情况，各种情况的反滤波因子的位置和形状不同，确定反滤波因子的方程的形式也不相同。,(1)当地震子波为最小延迟时,反滤波因子为,由于当 t0 时 w(t)=0 矩阵方程为,上式即前面得到的(3-80)式的结果。(2)当地震子波是最大延迟时,由于 tn 时 w(t)=0 矩阵方程为,当地震子波 为最大延迟时，可用上式求出反滤波因子 (3)当地震子波是混合延迟时,矩阵方程为,(3-81),当地震子波 混合延迟时，可用上式求出反滤波因子 。,二、参数选择,最小平方反褶积程序中的几个主要参数选择如下: (l)反滤波算子长度的选择:反滤波算子 长度 m 可以任意选择。一般在一个地区或一段测线上，通过试验来进行选择。,(3-82),反滤波算子长度可以选择 80ms,120ms,160ms,200ms,240ms 等，计算时要将它们换算为采样点数。,(2)相关时窗长度的选择:相关时窗长度 m+n 的选择，最小不应小于反滤波算子长度 m 的2倍，最长为地震记录的有效长度。 一个地震道可以开几个时窗，计算几个反滤波算子，在相邻两算子之间进行线性内插，对整个地震记录道进行时变反褶积。 (3)解矩阵方程时，稳定常数的选择:解矩阵方程时，表示随机干扰水平的常数。可以根据噪声水平进行选择。一般取为 r(0) 的百分数。 当地震记录h干扰较小时，稳定常数 e 选择为 在干扰较大的地震记录上，稳定常数 e 取作 稳定常数 e 最大可取作 0.100r(0) 。,一、预测滤波原理,预测问题就是已知某个物理量的过去值和现在值，通过对已知的信息进行加工处理获得未来某个时刻的预测值。例如，我们要击毁敌人的导弹，而在发射反导弹时，就不能对准敌弹现在的即时位置。因为当反导弹到达该位置时，敌弹已飞走了，肯定是击不中的。这就要求我们自动跟踪敌弹的运动，发射反导弹朝着射击目标在未来某个时刻的位置，使发射的反导弹与敌弹在该时刻同时到达空间相同的位置。为此就必须解决预测导弹未来位置的方法问题。预测滤波就可以解决这类间题。根据从实践中和理论上总结出来的规律，就可以设计一个算子预测因子，利用导弹现在和以前的信息(包括位置、速度、加速度等)来预测出它在未来某个时刻的位置。因而，预测滤波就是要设计一个预测因子 ，对输入 已知的过去值 和现在值 进行滤波处理，获得未来某个时刻通 时的预测值,3. 5 预测反褶积,(3-83),使它与实际的未来值 的误差预测误差,最小，按照最小平方原理，即使预测误差以 的平方和,为最小。得到,即,(3-84),或,令,分别表示时间延迟为 和 。的输入 的自相关函数。则(3- 85）式，可以写成,(3-85),(3-86),得到一个方程组，将上述方程组写成矩阵形式，得到,由输入 求出自相关函数 ，解矩阵方程(3-87)，即可得到预测滤波的因子 ，然后用预测因子 对输入进行预测滤波，得到未来的预测值,这就是利用 已知的过去值 和现在值 通过预测滤波所得到的未来 时预测 叫做预测步长。,(3-87),二、预测反褶积原理,上面讲的是一般的预测滤波问题。但是，在地震数据处理中所用的预测反褶积也是预测滤波)是用预测的方法，根据地震记录一次反射和干扰的信息预测出纯干扰部分，再由包括一次波和干扰的地震记录中减去纯干扰部分，得到消除干扰后的一次反射信号，以消除一次反射后面的海上鸣震等多次波干扰。在预测反滤波问题中，设计一个预测因子 对输入地震记录 的过去值 和现在值 进行预测，所得到的未来的预测值 是海上鸣震等多次波干扰，将它从包括一次反射和干扰的地震记录 中减去，所得到的预测误差,就是消除干扰后的一次反射信号。 设地震记录的数学模型为,其中，地震子波 满足最小相位条件，反射系数 是白噪声，它满足,其中， E 为数学期望。根据(3-88)式可以写出 时的 表达式,(3-88),(3-89),令 则上式可以写成,在(3-90）式中，要计算等式右边第一项，需要利用 等 t 时刻以后未来的信息。在时刻 t 这些值是未知的。因此，用 t 时刻已有的信息不能估计出这一项的值。而(3- 90)式等式右边第二项则可以用时刻 t 及 t 时刻以前的信息估算出来。这个估算值为,(3-90),(3-91),为了用 t 时刻及 t 时刻以前的值 求得未来某个时刻 的预测值 就必须将(3-91)式表示为值 与 的关系而不是 与 的关系 。,由平稳时间序列的预测分解定理可知，在满足(3-89)式条件下，地震子波 必然存在一个反子波 ，使得,将(3-92)式代入(3-91)式，得到,(3-92),(3-93),令 ，则(3-93)式可以写成,(3-94)式中的 为 t 时刻及t 时刻以前的值 等。 为预测因子。(3-94)式表明未来某个时刻 的预测值 为，时刻 t 及 t 刻以前的输入值 与预测因子 的褶积。,为了利用(3-94）式计算 时的预测值 ，除了要利用。 t 时刻和以前的 值，还必须知道预测因子 。,在确定预测因子。（l）时，仍然按照最小平方原理，使未来的预测值,(3-94),与实际的未来值 之间的预测误差,最小，或使预测误差 的平方和,为最小，得到,(3-95),或,令,分别表示延迟时间为 l-s 和 的地震记录 的自相关函数，则(3-96）式可以写成一个方程组,将上述方程写成矩阵形式，则为前述的矩阵方程(3-87),(3-96),(3-97),由输入地震记录 求出自相关函数 ，解矩阵方程(3-87)，即可求出预测因子 。然后用预测因子 与输入地震记录 进行褶积运算，得到未来 时的预测值,这就是预测的干扰。 将预测的未来值从实际未来值中减去，得到预测误差,(3-98),就得到了消除干扰的地震记录。,三、预测反褶积压缩反射脉冲,预测反褶积对反射脉冲或地震子波有着明显的压缩作用。它可以将原来一长度为 n 的地震子波 ，压缩为长度为 的窄脉冲 特别是预测步长很小时，地震子波 将发生明显的压缩。如果预测步长 ，则地震子波 将被压缩为尖脉冲,下面，进一步讨论这个问题,1. 的情况,在上一节最小平方反滤波中，将地震子波 压缩成为尖脉冲，其期望输出 为尖脉冲 ，根据最小平方原理，得到反滤波因子所满足的矩阵方程(3-79),其等式右边的列向量只有第一个元素为1，其他元素均为零。如果能证明预测因子所满足的矩阵方程(3-87),也能变成(3-79)式的形式，即其右边列向量的第一个元素不为零，其他元素均为零，即可证明预测反褶积在顶测步长 时，也可将地震子波压缩为尖脉冲。将 代入（3-87)式，并将 简写成 ，则得到,写成联立方程组的形式，得到,将上面方程组中各等式右边的项移到左边，并都乘以一1，再在方程组的前面补上一个恒等式,(3-99),(3-100),其中,(3-100)式变成,将上式写成矩阵形式，得到,(3-101),(3-102),将上式与(3一79)式比较，可以看出，两个矩阵方程在形式上是完全相同的。这就是说，如果，用(3-102)式所确定的滤波因子 对输入地震记录 进行最小平方反滤波，其输出为一尖脉冲，即,其中， 是由矩阵方程(3-87)所确定的预测因子， 是地震子波开始出现的时间。,如果，用(3-87)式所确定的预测因子 对输入地震记录 进行预测反褶积，当预测步长 时，所得到的预测值由(3-83)式可知为,(3-103),其输出的预测误差为,令上式中的 ，则得到,比较(3-103)式与(3-105)式，可以看出，两者等式右边完全相同，因而得出,或,(3-104),(3-105),由此可知，当预测步长 时，预测反褶积的输出可以将地震子波 压缩成为一个尖脉冲 。,2. 的情况,当预测步长 时，预测反褶积可以将长度为 的地震子波压缩成为长度为 的短脉冲。,在最小平方反褶积中，要将地震子波 压缩成为长度为 的短脉冲，则其期望输出必须为长度是 的短脉冲。这时，反滤波因子所满足的矩阵方程为,(3-106),其等式右边的列向量的前。个元素不为零，其他各元素均为零。如果能证明预测因子所满足的矩阵方程(3-87),也能写成(3-106)式的形成，即其右边列向量的前 个元素不为零，其他元素都为零，即可证明预测反褶积在预测步长为 时，也可将地震子波压缩为长度是 的短脉冲。将预测因子所满足的矩阵方程(3-87)式写成联立方程组的形式,(3-87),(3-107),将上面方程组中各等式右边的项移到左边，并都乘以一1，得到,再取,在(3-108)式中补充 个恒等式，得到,(3-108),写成矩阵形式，得到,(3-109),(3-110),比较(3-106)式和(3-110)式，可以看出，两者等式右边的形式相同。如果，用（3-110)式所确定的滤波因子,对输入地震记录 进行最小平方反褶积，其输出为一长度是 。 的短脉冲，其输出可以写成,其中， 是由矩阵方程(3-87)所确定的预测因子。,(3-111),如果，用(3-87)式所确定的预测因子 对输入地震记录 进行预测反褶积，当预测步长为 时，所得到的预测值为,其输出的预测误差为,令上式中的 ，则得到,比较(3-111)式与(3-113)式，两者等式右边完全相同，因而得出,由此可知，当预测步长为 时，预测反褶积可以将长度为 的地震子波 压缩为长度 的短脉冲。,(3-112),(3-113),四、参数选择,预测反褶积程序中的几个主要参数选择如下:,预测步长 ；是预测烦褶积的关键参数，预测反褶积的效果与预测步长 的选择有很大的关系。 为了选择最佳预测步长，可以将叠加前地震记录作为输入，进行预测反褶积试验，在其他参数（例如预测因子长度和预白化量）不变的情况下，改变预测步长 ，从不同预测步长的预测反褶积输出记录直接确定最佳预测步长 。 选择最佳预测步长 ，也可以用人工合成的地震记录作为作为输入，进行预测反褶积试验，同样，在其他参数不变的情况下，改变预测步长，从预测反褶积输出结果确定最佳预测步长 。,如图3-10所示，预测因子长度 m 不变为128ms，预白化量e为0%，采样间隔为2ms，取预测步长 分别为2ms，8ms，20ms，32ms，44ms，和88ms。图中（a）为脉冲相应，（b）为合成地震记录，地震子波为已知最小相位子波。从图中可以看出预测步长 = 2ms时，预测反褶积效果最好。,图3-10 用合成地震记录试验预测步长,（2）预测因子长度m；对预测反褶积的结果也有明显的影响。总的来说，预测因子不能太小。长度太小预测反褶积波形尾部出现明显波动，长度增大，这种尾部波动减弱，反褶积效果较好。达到一定长度时，反褶积效果将趋于稳定。预测因子长度也不能太大，长度太大不仅增加计算量，而且会增加相邻反射及噪声的影响。,图3-11所示为用合成地震记录确定预测因子长度的试验，预测步长 不变为2ms，预白化量 为0%，取预测因子长度分别为20ms，44ms，94ms，128ms，192ms，292ms，图中（a）为脉冲响应，（b）为合成地震记录，地震子波为已知最小相位子波。从图中可以看出。当预测因子长度太小（例如20ms和44ms）时，预测反褶积波形尾部波动明显，当预测因子长度达到 =128ms时，预测反褶积效果已经较好。,（3）预白化量 ：在进行预测反褶积计算预测因子时，要根据（3-87）式求解地震记录相关矩阵。为了求解矩阵方程的稳定性，需要根据随机噪声干扰水平确定随机噪声自相关值e或者它与地震子波自相关的百分比。 图3-12所示为合成地震记录确定预白化量得试验，预测步长a为2ms，预测因子长度m为128ms，预白化量分别为0,0.1%，1%，5%，10%和20%，图中（a）为脉冲响应，（b）为合成地震记录，地震子波为最新相位子波。从图中可以看出，当预白化量 =0.1%时，预测反褶积的输出结果基本不变，而当预白化量很大（ 5%），预测反褶积波形尾部波动明显增大，预测反褶积效果明显变差。,3. 6 子波整形反褶积,子波整形反褶积是通过所设计的整形滤波器，将地震子波转换成所期望的形状。,一、子波整形反褶积原理,为了将已知的地震子波 整形成为我们所期望的形状 我们要设计一个整形滤波器，其滤波因子为 ，根据最小平方反褶积原理，我们得到（3-78）式的结果,(3-78),上式是当期望输出 为任意形状时的最小平方反褶积滤波方程。利用这个矩阵方程，在已知地震子波 情况下，计算出整形滤波器的滤波因子 ，经过整形滤波后，即可得到与期望波形 在最小平方意义下最接近的整形后的处理结果。通常选择一个宽频带零相位窄脉冲作为期望输出 ，以获得尽可能高的分辨率。,二、子波整形反褶积,子波整形反褶积采用的期望输出波形 通常视输入地震子波 的相位性质不同，选择不同延迟的尖脉冲。正如第四节脉冲反褶积中讲过的，对于最小相位子波，其期望输出波形为零延迟尖脉冲，即 对子最大相位子波，其期望输出波形为最大延迟尖脉冲，即 而对于混合相位子波，其期望输出波形在零延迟尖脉冲与最大延迟尖脉冲之间的非零延迟尖 (图3-9)。,图3-13显示用不同延迟尖脉冲作为期望输出 ，输入为一混合相位子波 进行子波整形反褶积所得到一系列的结果。图中(0)表示一个混合相位输入子波， 8个框表示8个不同延迟尖脉冲作为期望输出的结果，从 (a) 至 (h) 尖禾冲的延迟时间逐渐加大。每个框中的(1）表示不同延迟的期望输出尖脉冲 （2）示整形滤波因子 （3）表示滤波后的实际输出 。我们的目标是将混合相位子波 经(3- 78）式计算得到的整形滤波因子 滤波后得到的实际输出与不同延迟尖脉冲 的误差平方和最小。,从图3-13中可以看出， (e) 框中具有60ms延迟的尖脉冲期望输出的子波整形反褶积效果最好，这时的实际输出 最接近非零延迟尖脉冲期望输出 。,从图中还可以看出，当实际输出 接近最佳延迟尖脉冲期望输出 时，尖脉冲延迟时间的改变对整形反褶积的结果影响不大。,在子波整形反褶积处理中，经常采用的期望输出除了上面的零延迟或非零延迟尖脉冲期望输出外，还经常采用零相位子波作为期望输出。即将最小相位、混合相位和最大相位的子波，保持其振幅谱不变，仅仅改变子波的相位，使其变为零相位子波。即将原来非零相位的子波，经过整形滤波后，将其整形为零相位子波。这是因为在具有相同振幅谱的各种相位的子波中，零相位子波具有最高的分辨率。,这样，就完成了最小相位子波整形为零相位子波的整形反褶积过程，在这个过程中，保持其振幅谱不变，仅仅改变子波的相位。,图3-14显示一个子波整形反褶积过程，它将一个最小相位子波通过子波整形反褶积处理整形成为一个振幅谱相同的零相位子波。,图3-14中(a)表示由地震道得到的自相关图，(b)是(a）平滑后的结果，（c)是(b）的单边曲线，(d)为计算得到的尖脉冲反褶积因子，(e)是反褶积因子(d）的逆得到的最小相位子波，(f)是与（e）具有相同振幅谱的零相位子波，(g)是根据零相位子波(f）设计的整形滤波因子，通过整形滤波可以将最小相位子波(e）整形为零相位子波(f)，(h)是整形滤波器的实际输出，它与零相位子波期望输出(f)是非常接近的。,3. 7 同态反褶积,同态反褶积不同于前面所讲过的最小平方反褶积或预测反褶积方法，它不需要假设地震子波的最小相位延迟性质，也不需要假设反射系数的白噪声性质。同态反褶积可以对任意相位延迟性质的地震子波进行反褶积。它主要是通过对地震记录频谱取对数，将地震子波和反射系数分离开来，原则上可以同时求取地震子波和确定反射系数，达到反褶积的目的。因而，这种方法又叫对数分解法。,下面讨论同态反褶积的原理。我们已经知道，反射地震记录 是由地震子波 和反射系数 褶积形成，即,由于褶积的原因，使地震子波 与反射系数 混在一起了，因而在地震记录 上，不能明显地显示出地震子波 的形状和反射系数 的位置。能不能设法将地震子波 和反射系数 分开呢?为此，将上式变换到频率域，得到,再对上式两边取对数，得到,其中， 分别表示地震记录 ，地震子波 和反射系数 频谱的对数，它们分别叫做 x 和 的对数谱，并分别用 表示，(3-114)式可以写作,由于,(3-114),(3-115),其中 和 分别表示频谱 的模和幅角，即分别是 的振幅谱和相位。因为对数函数是多值函数，为了使 是单值函数，只要 是 的连续函数就行,因为 是频率 的函数，就可以用傅里叶反变换将它从频率域变换到时间域，它的相应的时间序列,我们将 叫做 的对数谱序列。对(3-115）式进行傅里叶反变换，得到,其中， 和 分别是。 和 的对数谱序列。,(3-116),(3-117),上式表明，地震记录的对数谱序列分 是地震子波对数谱序列 和反射系数对数谱序列 之和。但是，地震子波对数谱序列 与反射系数对数谱序列 两者分布的位置几不同的。为了说明这个问题，举一个例子(图3-15）,假设有一个地震记录,(3-118),上式表明 是地震子波 与反射系数 的褶积。,由(3-118)式可知， 的频谱,其对数谱为,其对数谱序列为,(3-119),(3-120),(3-121),(3一121)式的对数谱序列示于图3一15（d）中，从该图可以看出地震子波的对数谱序列 集中在原点附近，而反射系数的对数谱序列 则是离原点较远的一系列尖脉冲。这就是说地震子波对数谱序列 与反射系数对数谱序列 分布在:轴的不同位置上，它们在一定程度上是彼此分离的。,因此，可以在原点附近沿 t 轴开一时窗，时窗长度为 ，将时窗 里的对数谱序列 保留，而将时窗外 的对数谱序列 去掉，这样，只要时窗的长度 取得合适，使 而又接近 ，就可以使时窗里所保留下来的 接近于地震子波的对数谱序列 图3一16（a）。相反的，如果离开原点 处开一时窗，将时窗 里的 保留，而将时窗外面 的 去掉，就有可能得到反射系数对数谱序列 。,上述开时窗保留或去掉 ，实质上就是对地震记录对数谱序列 进行滤波。 类似于频率域的低通、带通和高通滤波，可以