陆地地震资料的剩余静校正量的自动获取方法与流程

本发明属于勘探地震信号处理领域，涉及陆地地震资料自动剩余静校正技术，尤其涉及一种基于累积量相干积累法的时延估计自动获取陆地地震资料的初至剩余静校正量的方法。

背景技术：

目前，据文献[1]记载，中国西部地区成为油气勘探的重点区域，该地区地表条件相对复杂，如山前带沉积扇，戈壁和沙漠，近地表的非成岩低速介质区往往存在近地表高程、厚度和速度的空间变化，静校正问题突出，影响叠加成像效果以及地表风化层速度异常的检测。

折射类方法在静校正问题中有重要应用，如用于基准面静校正的折射静校正方法和层析静校正方法，以及基于统计相关实现的剩余静校正方法。这些方法的有效性都依赖初至拾取的准确度，但复杂地表条件往往伴随强背景噪声，降低初至信噪比(signal-to-noiseratio,snr)，难以实现准确的初至自动拾取，手动拾取在面对大量地震数据时人力时间成本极大。

现有方法中，文献[2]～[3]提出使用滤波类方法或地震干涉法提高初至波低信噪比并用于自动拾取。文献[4]提出折射波干涉法，通过相同炮点激发在不同检波点接收的折射波做互相关提取检波点间的折射波响应，利用不同炮对相同检波点间的折射响应具有多次覆盖的特点，实现同相叠加提高折射波信噪比，由此提高时延估计的准确性。文献[5]采用超虚干涉法(super-virtualinterferometry，svi)，在折射波干涉法的基础上再通过褶积叠加的步骤，进一步加强炮检点间的初至折射波的snr，用于提高自动初至拾取的准确性以及自动静校正的效率。折射波干涉法得到的虚折射道集包含的折射波时延信息可用于计算近地表折射层的厚度和速度，文献[6]将这一思路进一步发展用于计算剩余静校正量。

上述干涉类方法都利用互相关实现时延估计，针对低信噪比大数据量地震资料需要的初至拾取工作的人力时间成本可能很高；而且，由于互相关很容易受高斯噪声影响，尤其是相干高斯噪声，因此，这类方法抗噪性不强，实现的时延估计准确性不高，进一步计算得到的陆自动初至剩余静校正量的准确性也不高。

参考文献：

[1]符力耘,肖又军,孙伟家等,库车坳陷复杂高陡构造地震成像研究.地球物理学报,2013,56(6):1985–2001.

[2]bharadwaj,p.,schuster,g.,mallinson,i.,anddai,w.,theoryofsupervirtualrefractioninterferometry.geophysicaljournalinternational.2012,188,263–273.

[3]安圣培,胡天跃,崔永福等,起伏地表条件下复杂路径初至自动拾取.应用地球物理,2015,12(1):92–99.

[4]dong,s.,sheng,j.,andschuster,g.t.,theoryandpracticeofrefractioninterferometry.76thannualinternationalmeeting,segexpandedabstracts.2006,3021–3025.

[5]an,s.,hu,t.,liu,y.,peng,g.,andliang,x.,automaticfirst-arrivalpickingbasedonextendedsuper-virtualinterferometrywithqualitycontrolprocedure.explorationgeophysics.2016,online.doi:10.1071/eg14120.

[6]mikesell,d.,vanwijk,k.,calvert,a.,andhaney,m.,thevirtualrefraction:usefulspuriousenergyinseismicinterferometry.geophysics.2009,74(3),a13–a17.

技术实现要素：

为了克服上述现有技术的不足，本发明提供一种基于累积量的相干积累法(cumulant-basedcoherentintegration，cci)，针对折射波初至，通过高阶累积量求取检波点或炮点间的初至波时延，并进一步求取炮检点的剩余静校正量。

本发明的核心是：本发明对cci方法进行发展和改进，针对折射波初至，通过高阶累积量(higher-ordercumulant，hoc)求取检波点间的初至波时延，并用于求取检波点间的相对剩余静校正量差，由此逐点计算全部炮检点的剩余静校正量。针对实际数据可能存在的复杂噪声类型导致误差积累产生虚假长波长分量，应用局部回归平滑的方法校正剩余静校正量。相比折射波干涉法，本发明方法抗噪性更强，且有效压制相干的高斯噪声；相比传统的基于初至旅行时的静校正方法，本发明方法避免了面对低信噪比大数据量地震资料时可能需要的人力时间成本极高的初至拾取工作。

本发明的技术方案如下：

一种陆地地震资料的自动剩余静校正量的获取方法，基于累积量的时延估计，求取检波点间的初至波时延，将初至波时延用于求取检波点间的相对剩余静校正量差，通过逐点计算全部炮检点的剩余静校正量，实现陆地地震资料的自动剩余静校正；包括如下步骤：

a.针对在炮点sk激发在检波点rj接收的地震记录uk；j(t)(t为地震道的旅行时)可能存在的长波长静校正问题，对uk；j(t)进行高程静校正或野外静校正；

b.通过高通滤波去除地震数据uk；j(t)中的异常高频噪声，采用时窗大致截取远偏移距的初至波信号；

c.利用cci方法求取检波点间的初至波时延；

c1.求取炮点sk在不同检波点rj,rj+1,...,rj+n接收的n+1道地震记录的n+1阶累积量；

c2.利用不同炮点得到的累积量包含相同折射波响应的特点，通过相干叠加提高折射波响应的snr。

c3.对叠加得到的高阶累积量，以最大值位置确定检波点间折射波响应的时延。

d.求取检波点间或炮点间的相对剩余静校正量，再求取全部检波点和炮点的剩余静校正量；

e.通过局部回归方法，校正剩余静校正量。

本发明的有益效果是：

本发明提供一种基于累积量的相干积累法(cci)，针对折射波初至，通过高阶累积量求取检波点或炮点间的初至波时延，并进一步求取炮检点的剩余静校正量。与现有技术相比，本发明具有以下技术优势：

(一)本发明直接估计检波点间的初至时延并计算静校正量，相比传统基于初至旅行时的静校正方法，本发明方法可以避免可能耗费大量人力时间成本的初至拾取工作；

(二)本发明采用cci方法，可以提高折射波初至的信噪比，适用于低snr数据的处理，相比传统折射波干涉法，本发明方法的抗噪性提高；

(三)相比现有传统初至拾取方法拾取初至的起跳点，本发明通过峰值估计检波点间的初至时延，适用于可控震源数据。

附图说明

图1是本发明提供的自动剩余静校正量获取方法的流程框图。

图2是通过n+1阶cci方法加强折射波信号的示意图；

其中，sk是炮点；rj,rj+1,...,rj+n是检波点；v表示虚源位置；灰色实线和虚线分别表示具有正、负的旅行的折射波射线路径；使用不同炮点sk在检波点rj,rj+1,...,rj+n接收的地震记录求取的累积量叠加可以实现折射波信号的相干积累，并加强折射波信号。

图3是cci求取剩余静校正量示意图；

其中，si和sk分别是检波点左侧和右侧的炮点，rj和rj+1是相邻检波点，和是虚源位置；灰色细线表示两检波点共用的折射波路径，灰色粗线表示使用左侧炮点si估计的rj和rj+1间的折射波旅行时，黑色粗线表示使用右侧炮点sk估计的rj+1和rj间的折射波旅行时；实(虚)线表示旅行时为正(负)的部分。

图4是本发明实施例对含噪数据应用4阶cci静校正得到的炮集结果：

其中，(a)为未做静校正的单炮地震道集；(b)为应用炮点剩余静校正量的单炮地震道集；(c)和(d)分别是(a)和(b)中方框标示区域的局部放大图。

图5是实际数据应用本发明方法得到的剩余静校正量；

其中，(a)为检波点的剩余静校正量；(b)为炮点的剩余静校正量；对比本发明方法的结果(黑色细线)和基于手动拾取初至方法的结果(灰色粗线)，两者在大部分区域都吻合，而本发明可以避免手动初至拾取的繁重工作量。

图6是本发明实施例应用cci剩余静校正的炮集结果；

其中，(a)为静校正前的炮集；b)为本发明方法的结果；(c)为基于商业软件中初至方法的结果；(d)～(f)分别是(a)～(c)中方框标示区域的局部放大图。

图7是本发明实施例中的叠加剖面结果；

其中，(a)为未做静校正的叠加剖面；(b)为本发明方法的结果；(c)为基于商业软件中初至方法的结果。

具体实施方式

下面结合附图，通过实施例进一步描述本发明，但不以任何方式限制本发明的范围。

本发明提供一种基于累积量的相干积累法(cci)，针对折射波初至，通过高阶累积量求取检波点或炮点间的初至波时延，并进一步求取炮检点的剩余静校正量。

基于累积量的时延估计实现陆地地震资料的自动剩余静校正的具体实施方式如下：

a.针对在炮点sk激发在检波点rj接收的地震记录uk；j(t)(t为地震道的旅行时)可能存在的长波长静校正问题，对uk；j(t)进行高程静校正或野外静校正；

b.通过高通滤波去除地震数据uk；j(t)中的异常高频噪声，用时窗大致截取远偏移距的初至波信号；

c.利用cci方法求取检波点间的初至波时延；

c1.求取炮点sk在不同检波点rj,rj+1,...,rj+n接收的n+1道地震记录的n+1阶累积量；

二阶统计工具，如互相关方法，常被用于时延估计，但互相关方法对高斯噪声敏感。在实际问题中，信号往往是非高斯过程，而加性噪声具有高斯性，因此，使用高阶统计工具，如hoc方法，进行时延估计更加合理，因为高斯噪声的hoc恒为零。文献[7](mendel,j.m.,tutorialonhigher-orderstatistics(spectra)insignalprocessingandsystemtheory:theoreticalresultsandsomeapplications.proceedingsoftheieee,1991,79(3),278–305.)记载了hoc的数学表达和性质。对在炮点sk激发在检波点rj接收的地震记录uk；j(t)，假设它由折射波信号和高斯随机噪声组成。求取炮点sk在不同检波点rj,rj+1,...,rj+n接收的n+1道地震记录的n+1阶累积量[cj；j+1,...,j+n(τ1,...,τn)]k(如图2所示通过n+1阶cci方法加强折射波信号)，其中，下标k表示该累积量是通过使用炮点sk有关的地震记录得到，(n+1)阶累积量包含n个时延τ1到τn，用于表示uk；j(t)与uk；j+1(t+τ1)到uk；j+n(t+τn)所包含的折射波响应间的时延。假设数据包含的噪声是高斯随机噪声，因此噪声的高阶累积量恒为零。图2是通过n+1阶cci方法加强折射波信号的示意图；其中，sk是炮点；rj,rj+1,...,rj+n是检波点；v表示虚源位置；灰色实线和虚线分别表示具有正、负的旅行的折射波射线路径；使用不同炮点sk在检波点rj,rj+1,...,rj+n接收的地震记录求取的累积量叠加可以实现折射波信号的相干积累，并加强折射波信号。

c2.利用不同炮点得到的高阶累积量包含相同折射波响应的特点，通过相干叠加提高折射波响应的snr，获得高阶累计量；

如图2所示，炮点sk得到的累积量包含的τp(p＝1,2,…,n)表示检波点rj和rj+p间的折射波响应的时延，假设v表示rj对应的虚源位置，sk到rj的折射波旅行时为sk到rj+p的折射波旅行时为共用的旅行时在计算累积量的过程中被抵消，因此时延τp不包含与炮点sk有关的旅行时，而v的位置仅与检波点位置有关，因此使用不同炮点得到的累积量[cj；j+1,...,j+n(τ1,...,τn)]k包含的折射波响应的时延相同，可以实现检波点间折射波响应的相干叠加。同时，反射波和非相干的随机噪声由于具有不同的时延，不能相干叠加，因此被相对压制。

c3.对相干叠加得到的高阶累积量，以最大值位置确定检波点间折射波响应的时延；作为检波点间的初至波时延。

d.利用cci得到的检波点间的初至波时延，求取检波点间或炮点间的相对剩余静校正量，再求取全部检波点和炮点的剩余静校正量；

如图3，假设远偏移距初至是来自同一层出射的折射波，在地表一致性条件下，假设检波点rj和rj+1具有的剩余静校正量为和则检波点间的初至波时延为：

其中，δtj,j+1是检波点j和j+1间的不包含剩余静校正量的时延，而δtj+1,j表示j+1和j间的不包含剩余静校正量的时延。δt^left和δt^right分别表示使用左侧炮点si和右侧炮点sk估计的检波点间的时延，其中包含了剩余静校正量。

为求取剩余静校正量，假设地表高程和低速层厚度引起的长波长静校正量已被消除，则δtj,j+1和δtj+1,j包含的上下行射线具有相同的旅行时，同时高速层顶界面的起伏在长波长静校正中被修正，δtj,j+1和δtj+1,j包含的沿折射界面的传播距离都近似等于检波点间距，因此在折射界面速度不发生剧烈变化的情况下，沿折射界面传播的旅行时也相同，因此δtj,j+1和δtj+1,j相等。文献[8](zhang,c.,andzhang,j.,2dseismicresidualstaticsderivedfromrefractioninterferometry.journalofappliedgeophysics,2016,130,145–152.)记载了上述推导过程。

将表达式(式1)和(式2)相减，得到：

其中，δt^left和δt^right可以由cci方法估计得到，因此得到了相邻检波点的剩余静校正量的差值。假定任一检波点的剩余静校正量为零，可以根据相邻检波点的剩余静校正量差，逐一求解所有检波点的剩余静校正量。

e.通过局部回归方法，校正剩余静校正量，去除误差积累导致的长波长分量。

实际资料中的噪声和其他干扰会影响时延估计的准确性，即求取的检波点间的静校正量差可能存在误差，由此，逐个求取全部检波点的绝对静校正量值会产生误差积累，在静校正量曲线中产生明显的长波长分量，这部分分量与真实的地下构造形态无关，因此要进一步对cci方法得到的剩余静校正量做长短波长分离，去除误差积累导致的长波长分量。将检波点的静校正量视为连续排列的散点，使用局部回归平滑的方法提取长波长分量，即对每一点j周围一定平滑半径内的数据做加权最小二乘拟合，加权过程使距离点j距离近的点具有相对大的权重，由拟合得到的多项式可以得到点j平滑后的数值，由此提取的平滑曲线即长波长分量。文献[9](clevelandw.s.,anddevlin,s.j.,locallyweightedregression:anapproachtoregressionanalysisbylocalfitting.journaloftheamericanstatisticalassociation,1988,83(403),596–610.)和文献[10](fotiadis,d.a.,astaras,a.,kalfas,a.,andbamidis,p.d.,anovelsignalprocessingmethodbasedonthefrequencymodalityforintra-bodymedicalinstrumenttracking.5thinternationalconferenceonmoderncircuitsandsystemstechnologies(mocast).2016,1–4.)分别记载了上述局部回归拟合方法。

由于误差累积得到的长波长分量不具有特定的曲线形态，因此，本发明选择了局部而非整体的回归过程，好处在于不要求被拟合的散点具有某种特定的全局函数形态或理论模型。文献[11](taner,m.t.,wagner,d.e.,baysal,e.,andlu,l.,aunifiedmethodfor2-dand3-drefractionstatics.geophysics.1998,63(1),260–274)记载了常规的长短波长分离选取菲涅尔带半径作为平滑半径。但是，本发明涉及的误差累积导致的长波长分量不具有真实物理意义，即与空间排列长度无关，而是由检波点数量决定，因误差累积是逐点计算产生的。平滑半径应选取适中：过大会导致长波长分量残留，过小会导致过度拟合损伤短波长分量，本发明采用文献[12](sun,y.,tonellot,t.,kamel,b.,andbakulin,a.,atwo-phaseautomaticstaticcorrectionmethod.geophysicalprospecting.2016.online.doi:10.1111/1365-2478.12453.)记载的方法，具体数值可根据实际的剖面叠加的能量最大化标准和对构造的改变程度综合决定。由于检波点和炮点的关系可以互易，因此同样过程可以求解炮点的剩余静校正量。

下面通过实例对本发明做进一步说明。

实施例一：

利用合成数据实例验证cci方法计算剩余静校正量的有效性。通过声波有限差分正演得到合成数据，所用模型参数见表1，炮距为40m，道距为10m，炮点数100，检波点数400。

表1模型参数

剩余静校正量由随机数序列产生，满足均值为0ms和标准差为8ms的正太分布。使用含噪数据测试，加入的噪声仍是相干高斯噪声，含噪数据snr为–7db。将4阶cci方法得到的剩余静校正量应用到含噪炮集中，看到反射波的同相轴连续性明显提高(图4)。该实例验证了cci方法求取剩余静校正量的正确性以及较好的抗噪性。

实施例二：

选取中国西部地区的实际地震资料，检验cci方法的有效性。该工区地表存在较强的近地表速度异常，仅校正高程导致的长波长静校正量后仍存在剩余静校正问题。图5a的黑色细线给出了使用4阶cci方法及局部回归校正得到的检波点剩余静校正量，并对比了某商业软件中的基于拾取初至的初至波剩余静校正模块的结果(灰色粗线)，使用的是广义互换方法[11]，两者基本吻合，只在检波点号为900和1050附近位置存在微弱差异。相比之下，为了保证基于拾取初至的初至波剩余静校正效果，需要拾取全部道集的初至，在远偏移距处，受该工区复杂近地表的影响，初至snr极低，使用手动初至拾取耗费大量的时间人力成本，而本发明可以避免拾取初至的工作量，也取得了相似的结果。图5b为炮点剩余静校正量的结果，同样，与基于拾取初至方法的结果基本吻合。进一步将剩余静校正量应用到全部地震道集，选取了剩余静校正问题比较突出的原始单炮道集(图6a)，可以看到明显的初至波同相轴的不连续，图6b和6c分别是应用本发明和基于拾取初至方法的结果，都明显改善了初至波同相轴的连续性。图6d～f分别是6a～c的局部放大图，反射波同相轴的连续性也明显改善。图7是时间叠加剖面的对比图，未做剩余静校正的剖面的同相轴连续性较差，应用本发明和基于初至方法的效果相当，都明显改善了同相轴的连续性，但本发明避免了拾取初至的工作量。

需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的精神和范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。