2017, Vol. 60
2. 中国石油勘探开发研究院, 北京 100083;
3. 中国石油塔里木油田公司, 新疆库尔勒 841000
2. Research Institute of Petroleum Exploration and Development, PetroChina, Beijing 100083, China;
3. Tarim Oilfield Company, PetroChina, Xinjiang Korla 841000, China
地震记录不是稳态序列，而是随传播时间变大，振幅和相位均不断变化的非稳态序列.为更准确地刻画地下地质情况，需要对地震信号的时变信息进行研究.同时，随着勘探目标埋深和复杂度加大，地震流体识别难度也不断加大，这对地震剖面的分辨率及反演解释的精确度提出更高的要求.地震子波估计的准确程度会直接影响这两者的结果.目前已有的地震子波估计法得到的结果大多是不随时间变化的子波.因此我们需要研究对时变子波的估计方法.
子波估计的误差对后续解释工作带来的影响已有相应的研究.子波相位估计的微小误差会使反演结果与真实值相差甚远(Yuan and Wang，2011; Luo et al.，2014).White和Simm(2003)也给出例子表明，一个微小的错误子波时移估计会显著影响振幅属性切片，从而给解释工作带来很大难度.时变子波估计的结果，不仅能提高井震标定的精确度，修正地震剖面的层位标定；更重要的是将时变子波用于反褶积和反演后续工作，可以为地震解释提供更准确的依据.高静怀等(2009)给出了基于时变子波的地震记录褶积模型的近似数学表达式，并通过自适应时频分析改进谱白化方法，提高地震记录分辨率.Han等(2015)根据合成数据和实际例子，研究时变子波的相位变化与流体预测的关系，指出子波的相位变化是由于流体的吸收衰减性质.对时变子波的研究还有助于研究地震信号的吸收衰减补偿，提高地震记录的分辨率.子波的旁瓣会影响地震记录的横向连续性(Yuan and Wang，2013).通过反褶积消除子波的影响，可以拓宽地震记录的频带，校正相位畸变(Wang et al.，2015).Van der Baan(2012)对比研究反Q滤波与时变维纳反褶积的应用效果，指出纯相位反Q滤波更适合校正频散的影响，而时变维纳反褶积则更适合补偿吸收衰减的影响，在实际应用中两者应先后使用.
现有的地震子波估计方法可以分为两大类(Tygel and Bleistein，2000)：一类是确定型方法，在测井数据丰富的区域，利用测井数据与地震数据进行匹配拟合估算地震子波.估计的子波可用于修正地震剖面的层位标定.实际的合成地震记录层位标定工作中常用到所谓的“拉伸”或“压缩”方法，然而这样的手段既不够科学且缺乏质量控制标准(李国发等，2008).Edgar和Van der Baan(2011)提出以下思路，首先通过地震数据估计子波振幅谱并附以零相位得到初始零相位子波，用其与测井得到的反射系数制作合成记录，通过合成记录与地震记录进行匹配并不断调整子波的相位和振幅，得到最佳匹配结果时对应的子波即为估计子波.Wang等(2015)基于这样的思路，利用蚁群算法搜索匹配最佳子波相位，从而迅速有效地对子波相位进行估计.
Walden和White(1998)提出一种谱分析法.该方法无需对地震子波进行任何假设，通过测井数据和地震记录的功率谱与交互功率谱在一定时窗内估算地震子波.该方法虽然对子波无需假设，但认为在分析时窗内子波是不随时间变化的.White和Naeini(2014)又将这一方法推广到宽频带地震数据子波估计应用中，特别研究了子波在低频端的估计方法.Ma等(2010)根据这一方法提出了测井、垂直地震剖面(Vertical Seismic Profiling，VSP)和地震数据相结合的子波估计方法.
第二类地震子波估计方法是统计型方法，该方法仅利用地震数据，需要对地震子波进行假设，分析反射系数序列的统计分布特征从而估算地震子波.这类方法通常包含反射系数序列白噪及地震子波为最小相位或零相位的假设条件.Wiggins(1978)提出的峰态最大化方法就是一种典型的统计型方法.其中峰态(kurtosis)是统计学中衡量序列非高斯性的参数.除峰态最大准则外，前人提出许多用于衡量反褶积结果的判别准则(Yu et al.，2012)，例如简约准则(Claerbout，1977)、指数变换准则(Ooe and Ulrych，1979)等.Wang等(2014)指出峰态最大准则存在的缺陷，如应用于低频数据时相位估计不准确，并提出一种新的基于Rényi散度的准则，其对于埋藏较深、主频较低的地震记录有更好的效果.刘玉金等(2014)将峰态最大准则和包络最大相似度准则结合，利用测井数据和井旁地震道数据对子波的时变相位进行估计.Ma等(2015)对峰态和偏度(skewness)两种判别准则在子波相位估计中的应用进行了研究，指出两种准则都会受反射系数分布特性的影响，但相比于峰态属性，在同样的相位变化下偏度属性具有更大的动态范围.同时，偏度属性在含噪的情况下表现更稳定.
Van der Baan(2008)，Van der Baan和Fomel(2009)，Van der Baan等(2010)根据Levy和Oldenburg(1987)、White(1988)和Longbottom等(1988)的研究成果，对Wiggins的方法进行改进，提出了非稳态混合相位子波估计方法.假设反射系数序列具有白噪非高斯性，通过对地震数据进行相位旋转至达到非高斯性最强，则可得到子波的相位.再与通过自相关求取的振幅谱结合得到时变子波估计结果.Fomel和Van der Baan(2014)进一步改进，提出一种更稳定的时变子波相位估计方法：局部偏度属性最大化方法.换用偏度统计量判别非高斯性，且通过引入局部属性的概念(Fomel，2007a，2007b)，使用正则化最小二乘优化使得该方法可以更稳定地估计时变信息.Fomel等人使用这一方法对局部观测子波的相位进行估计.本文使用偏度属性最大化方法对传播子波的相位进行估计.
根据以上学者的研究，本文所述的时变子波混合估计法是结合谱分析法和偏度属性最大化法两者的优点，通过谱分析法利用测井数据与地震记录对时变子波的振幅进行估计，同时通过偏度属性最大化法直接利用地震记录对时变子波的相位进行估计.与传统的维纳-莱文森滤波法(Robinson and Treitel，2000)或常相位旋转法(Levy and Oldenburg，1987; Longbottom et al.，1988; White，1988; Gao and Zhang，2010; Wang et al.，2014)相比，这种方法不需要对子波进行任何时不变或相位等的假设，采用对子波的振幅和相位分别求取，实现对时变相位的估计.通过模型测试表明该方法能正确提取时变子波的振幅和相位，从而准确地估计时变子波.随后，将该方法用于实际数据，根据井震拟合的效果对子波估计结果做定性判断，并由拟合度和精确度两个质量控制参数定量判断子波估计结果.
2 方法原理 2.1 时变子波振幅估计实际地震记录可以用以下褶积模型来表示，假设模型中包含的噪声都随机且平稳.对于真实的反射信号st，rt为真实的反射系数序列，wt为地震子波，则
|
(1) |
对于测井数据，xt为测井合成的反射系数序列，ut为随机噪声，则
|
(2) |
对于地震数据，yt为实际接收的地震数据，vt为随机噪声，则
|
(3) |
我们的目的是求取地震子波wt的振幅.针对这一问题，我们基于谱分析法(Walden and White，1998)，利用平稳随机过程的谱分析理论求解这一含噪声系统.在频率域内，根据(1)—(3)式，得到xt与yt的功率谱和交互谱的3个表达式，并联立成方程组(4)—(6).
|
(4) |
|
(5) |
|
(6) |
其中，Φxx(f)、Φyy(f)、Φuu(f)、Φvv(f)、Φrr(f)、Φss(f)分别是xt、yt、ut、vt、rt、st的功率谱，Φxy(f)是xt和yt的交互谱.
(4)—(6)中可以直接由已知的地震和测井数据得到的量共有三个：Φxx(f)、Φyy(f)和Φxy(f)，而未知量却有四个：Φuu(f)、Φvv(f)、Φrr(f)和W(f)，所以此方程组为欠定方程组.这里引入一个额外的表达式v(f)来求解，v(f)定义为输入输出的相对噪声功率比.
|
(7) |
v(f)通过多道地震数据做互相关计算即可求得(Walden and White，1998).将这个量与上述方程组联立求解即可解出估计子波的振幅.
对有限长度序列做功率谱估计时，需使用平滑窗函数减小随机噪声及异常值的影响.使用平滑窗一方面可以减小功率谱估计的误差，但不可避免会对估计谱造成一定畸变.当平滑窗函数在时域内长度过大，则平滑不足从而误差消除不够；如果平滑窗长过小，则会导致平滑过度而引起严重的变形和偏差.因此，平滑窗长的选择需经过实验选择一个平衡点.
2.2 时变子波相位估计高阶统计量是估计地震子波相位的有效工具(Yu et al.，2011).常用的相位校正的算法为峰度最大化算法，是由Wiggins(1978)首先在基于最大峰度准则的盲反褶积算法中提出的.Wiggins算法随后由Levy和Oldenburg(1987)、White(1988)等进行改进，通过常相位假设使计算过程中的自由度减少到一个，从而保证反演的稳定性.这一方法能对子波相位进行校正的原理为：满足白噪假设的时间序列与任意序列褶积，所得结果都比原白噪序列更加高斯化.因此白噪的反射系数序列与子波褶积后，所得的地震记录都比原反射系数序列的高斯性更强(Donoho，1981).基于这一理论，通过对地震记录进行相位旋转并衡量结果的高斯性即可求得子波的相位.
地震记录在时间域内进行相位旋转的表达式为：
|
(8) |
其中，φ为旋转角度，H[.]表示希尔伯特变化.
序列的三阶统计量偏度可以衡量序列的高斯性，通过衡量序列的对称性从而判断高斯性强弱.与序列的峰态属性相比，序列的偏度属性具有动态范围更大、需要样点数更少的特点(Fomel and Van der Baan，2014).对于地震序列y(t)，其偏度可以表示为(Fomel and Van der Baan，2014)：
|
(9) |
其中skew(y)为地震序列的偏度属性，N为地震道的采样点数.
当偏度的绝对值越大，则序列的非高斯性越强.在应用中，通过对原地震记录从-180°到180°进行相位旋转，并提取每次旋转后地震道的偏度属性，当某结果对应的偏度绝对值最大时，此时的旋转角度的负数即对应子波此时的相位角度，即－φ(Feng et al.，2015).
根据任意两个序列an和bn的相关系数定义：
|
(10) |
可以推导skew2(y)的表达式如下：
|
(11) |
根据式(11)计算偏度属性时，根据Fomel(2007a)提出的局部属性估计法，将问题转化为正则化最小二乘反演问题.这一方法对非稳态信息的估计比常用的滑动时窗法具有更强的鲁棒性.在正则化最小二乘反演中，使用的正则化因子长度的大小决定估计的时变相位是局部相位或全局相位.正则化因子长度小则估计的结果为局部观测子波的时变相位信息；反之若正则化因子长度大，局部的地质特征经过平均后无法体现，则估计的结果为传播子波的时变相位信息，即本文的估计结果.
2.3 时变子波估计为完成时变子波的估计，只需将前两步估计得到的时变振幅和相位信息组合即可.第一步，在地震道上给定时窗，通过时变子波振幅估计得到该时窗内的估计子波振幅.通过向下滑动时窗并保证一定的重叠，得到估计子波的时变振幅.第二步，对整个地震道做时变子波相位估计.最后，通过(12)式完成时变子波的估计：
|
(12) |
其中，定义中间时刻为ti的时窗为第i个时窗，|Wi(f)|为第i个时窗内估计的振幅，φskew，i为时变相位在第i个时窗内的均值，sgn(f)为sign函数.
2.4 时变子波估计的准确度评价为评价估计的时变子波是否准确，可以从定性和定量两个角度判断.一方面，可以利用估计子波与从测井数据中提取的反射系数进行褶积，
得到合成记录$\hat s\left( t \right)$将此合成记录与实际地震道y(t)进行匹配，定性地观察匹配效果从而判断估计子波的准确程度.另一方面，可以利用以下两个统计量对匹配效果进行定量评估.
预测能量比例PEP，表示合成记录$\hat s\left( t \right)$与实际地震道y(t)相匹配的能量占实际地震道总能量的比例.PEP越接近1，说明地震数据和测井数据拟合得越好.
|
(13) |
归一化均方误差NMSE，表示合成记录$\hat s\left( t \right)$和真实反射信号s(t)之间的误差.NMSE越接近0，说明预测值越接近真实.
|
(14) |
但真实反射信号无从得知，因此根据预测理论对NMSE进行必要的推导，得到近似NMSE(Walden and White，1984)：
|
(15) |
其中，T为参与匹配的地震道长，b为时变振幅估计中选用的平滑窗带宽.
可见，时变振幅估计中的平滑窗长选择，可以经过多次试验，通过对估计结果定量评价，选择PEP较大且NMSE较小的即最合适的值.
2.5 时变反褶积求得子波w(t)之后可以对地震道进行反褶积.根据最小平方反滤波算法，反滤波因子可以由子波w(t)得到，在频率域内(Berkhout，1977)：
|
(16) |
其中，F为反滤波因子，c2为稳定因子，表示输入地震道的噪声水平；*符号表示复共轭.将反滤波因子作用于地震道即可达到消除子波的滤波效应.
时变反褶积的处理过程类似于常规最小平方反滤波，不同之处在于反滤波因子不是固定不变的(van der Bann，2008).根据式(12)估计的时变子波和式(16)联合可以得到时变的反滤波因子，将这些反滤波因子依次作用于整个地震道可以得到每个反滤波因子对应的反褶积结果.将这些反褶积结果根据(17)式进行整合，即可得到最终的时变反褶积结果.
|
(17) |
其中，d(t)为时变反褶积结果，di(t)为利用时变子波的第i个估计结果的反褶积结果.
3 模型测试为了验证时变子波混合估计法，现通过理论模型试算说明该方法的有效性.这里将子波的振幅和相位的时变性分开讨论.
3.1 固定相位、主频时变的模型首先，利用固定相位、主频时变的雷克子波与随机时间序列褶积得到合成地震记录.在零时刻子波主频为60 Hz，随传播时间主频线性变化，在最终时刻达到15 Hz.这里验证了子波为零相位和90°相位的两种情况.合成记录长为2.1 s，采样率为2 ms.
图 1a为随机时间序列，图 1b和1c分别为零相位和90°相位雷克子波褶积随机时间序列得到的合成记录.合成记录由浅至深表现出由于子波主频降低导致的分辨率降低.同时，对比两个合成记录也可以看出，零相位记录具有很好的对称性，例如0.4 s处.经过时变子波估计后，估计的相位与真实相位的对比如图 1d和1e.估计值与真实值的误差基本保持在20°以内.为进一步验证估计子波的准确性，将估计子波与真实反射系数褶积，与真实合成记录进行匹配，计算PEP和NMSE值，从而评价时变子波估计结果的准确性.这里展示90°相位雷克子波合成模型在0.33~0.45 s、0.93~1.15 s和1.65~1.77 s三个时窗的估计子波结果(图 2a、2b、2c)，可见子波随时间传播主频降低.图 2d、2e和2f分别为三个子波的匹配评价，PEP分别为86.6%、99%和92%，匹配拟合程度较高.与此类似，零相位估计子波的匹配拟合程度也很高，表明子波估计结果的准确性.
|
图 1 0°相位和90°相位变主频雷克子波合成算例 (a)反射系数；(b)0°相位子波的合成地震记录；(c)90°相位子波的合成地震记录(d)0°相位子波的时变相位估计结果；(e)90°相位子波的时变相位估计结果. Fig. 1 Synthetic example for time-varying dominant frequency wavelets with fixed phases 0° and 90° |
|
图 2 90°相位子波合成模型估计时变子波匹配评价 (a)0.33~0.45 s，(b)0.93~1.15 s，(c)1.65~1.77 s三个时窗的估计子波结果；(d)、(e)、(f)为对应子波的匹配结果. Fig. 2 Matching evaluation of 90° phase wavelet estimation on synthetic model for three time-windows (a)0.33~0.45 s;(b)0.93~1.15 s;(c)1.65~1.77 s;(d)，(e)and(f)are matching results of corresponding wavelets. |
其次，验证主频和相位均时变的雷克子波合成模型.在零时刻子波主频和相位分别为(60 Hz，0°)随传播时间主频和相位线性变化，在最终时刻达到(15 Hz，90°).合成记录主频降低，子波的相位也随时间增加而变化，从而模拟实际情况中地层吸收衰减特性对地震记录的影响.合成记录长2.1 s，采样率为2 ms.
图 3a为随机时间序列，图 3b为褶积得到的合成记录.合成记录由浅至深表现出由于子波主频降低导致的分辨率降低.图 4a为0.3~2 s内任选的9个瞬时时刻对应的真实子波波形图.对于时变子波振幅估计，估计功率谱与交互谱时使用120 ms长的Papoulis时窗.估计的子波长为60 ms.对于时变子波相位估计，使用的正则化因子长度为200个采样点，即0.4 s.图 4b为9个瞬时时刻对应的估计子波波形图.图 4c中绿色线表示真实子波的相位随时间变化曲线，蓝色线表示估计子波的相位随时间变化曲线.从图可见估计子波的波形与真实子波的波形符合较好，同时估计相位与真实值的趋势也相符.图 4d为使用White等人提出的谱分析法对整个2 s记录估计的稳态子波，由于模型中真实子波是随时间变化的，因此基于稳态子波的估计结果必然会存在偏差，尤其是估计子波的相位与真实的时变相位存在偏差.
|
图 3 相位、主频均时变的子波合成算例 (a)反射系数；(b)合成地震记录. Fig. 3 Example of synthetic wavelets with time-varying phase and time-varying main frequency (a)Reflectivity coefficients;(b)Synthetic seismic record. |
|
图 4 相位、主频均时变的子波合成算例估计结果 (a)从0.3~2 s间的9个瞬时时刻的真实子波波形；(b)对应时刻的估计子波波形； (c)真实子波相位与估计值的对比图；(d)估计的全局稳态子波波形. Fig. 4 Estimation of wavelets with time-varying phase and time-varying main frequency (a)Instantaneous wavelets waveforms at nine different times between 0.3~2 s;(b)Estimated wavelet waveforms at corresponding times; (c)Comparison of estimated and true phases of wavelets;(d)Estimated waveforms of global stable wavelets. |
图 5为合成记录时变反褶积的结果对比.经过时变反褶积后，记录的分辨率有明显的提高，有助于识别薄层，例如0.3 s和0.6 s附近的薄层经过时变反褶积后更易识别.同时，时变反褶积也使记录旋转至零相位，使结果有了更好的对称性，例如在1.1 s处的反射信号.
|
图 5 相位、主频均时变的子波合成算例时变反褶积结果 (a)合成记录的反射系数；(b)合成地震记录；(c)时变反褶积结果. Fig. 5 Time-varying deconvolution result of synthetic wavelets with time-varying phase and time-varying main frequency (a)Reflectivity coefficients;(b)Synthetic seismic record;(c)Time-varying deconvolution result. |
为了验证效果，我们将该方法用于实际资料，资料来源于塔中低凸起内的塔中72号井附近的三维实际地震数据和测井数据，地震数据取过塔中72号井(inline371，crossline810)、边长1 km的正方形区域内的数据(inline351-391，crossline790-830)，经过静校正、去噪、叠加和时间偏移等流程得到叠后偏移剖面.地震数据采样间隔为4 ms.
在进行时变子波估计前，需要对数据进行以下预处理：测井数据编辑筛选及校正.声波测井数据经过校正可以显著提高拟合质量.通过使用VSP数据得到的CheckShot曲线对声波测井数据进行校正，消除时间误差.此外，该井密度测井数据的浅层部分缺失，可由声波测井数据进行拟合.利用经过校正的声波和密度测井数据可以生成反射系数序列.必要时需要将生成的反射系数序列重采样使之与地震记录的采样间隔相同.此处对生成的反射系数序列重采样至4 ms(如图 6).
|
图 6 实际测井数据的预处理结果 (a)校正后的声波阻抗曲线；(b)合成反射系数；(c)重采样的合成反射系数. Fig. 6 Preprocessing of real log data (a)Corrected impedance log;(b)Synthetic reflectivity;(c)Resampled reflectivity. |
利用测井数据与地震数据，通过谱分析方法求取时变子波的振幅.同时，通过对地震数据使用偏度最大化方法求取时变子波的相位.图 7a为时变相位扫描的结果.使用的正则化因子长度为0.5 s、10道.图 7b为子波相位随时间变化的趋势曲线.可见该地区的地震数据在4 s内相位变化较大.因此，通过时变子波估计及反褶积处理后的地震记录能够更好地满足精细解释的需要.图 8为最终估计得到的时变子波，此处选取0.94~1.44 s，1.83~2.33 s和2.56~3.06 s三个时窗的子波估计结果进行展示.图 8d为利用2.56~3.06 s时窗内的估计子波做井震匹配的结果.地震数据中72.9%的能量可以被估计的子波与由测井数据合成的反射系数褶积得到，归一化均方误差为0.11.考虑到资料来自复杂地质背景的地区及原始资料信噪比较低，这一拟合效果较为满意.
|
图 7 实际地震数据的时变子波相位信息估计 (a)偏度属性极值提取；(b)子波相位随时间变化的曲线. Fig. 7 Estimation of time-varying phase on real seismic data (a)Extracting skewness extremum;(b)Curve of wavelet phase versus time. |
|
图 8 实际数据时变子波估计 (a)0.94~1.44 s，(b)1.83~2.33 s和(c)2.56~3.06 s三个时窗的估计结果；(d)用2.56~3.06 s时窗内的估计子波作井震匹配的结果. Fig. 8 Estimation of time-varying wavelets on real data for three time windows (a)0.94~1.44 s;(b)1.83~2.33 s;(c)2.56~3.06 s;(d)Well matching results to estimated wavelets for 2.56~3.06 s. 1: acoustic wave impedance，2: synthetic reflectivity，3: estimated wavelet，4 and 6: real seismic data，5: synthetic seismic data，7: residuals. |
本文所述的时变子波混合估计方法有以下两个假设条件：地下介质的反射系数序列是白噪的，是非高斯分布的.根据测井数据的分析，地下介质的反射系数序列基本满足非高斯的假设条件(Walden and Hosken，1986).由于该方法中使用偏度属性判别非高斯性，具体来说反射系数序列需要具有不对称的分布性质.实际情况为地下介质的反射系数随深度的增加呈增大的趋势，因此反射系数为正数的概率更大，则反射系数序列总体可以呈现不对称的分布性质，从而保证实际地震数据使用偏度属性最大化法估计时变子波相位行之有效.然而，实际情况中也可能在某一深度段对应介质的反射系数序列是高斯性的或是对称的.
然而，地下介质的反射系数序列并不满足白噪假设条件，而是更趋向为蓝噪声(Walden and Hosken，1985)，即在有限带宽内，功率谱随频率的增加而变大.反射系数序列的非白噪问题可以通过区块内的测井数据来解决(Saggaf and Robinson，2000).Van der Baan(2008)对这一假设开展了更为详细的分析.
Fomel和Van der Baan(2014)提出局部偏度属性最大化方法解决局部观测子波时变相位估计问题.本文使用这种方法，用于估计传播子波时变相位.Van der Baan等(2010)指出传播子波与局部观测子波的差别.传播子波为实际穿过地下介质的子波，是地下介质与几何扩散、固有衰减以及散射的综合产物.而局部观测子波为空间某一点或某一时刻对应的瞬时子波，它反映的是局部的地质特征，与传播子波此时的波形以及该点或该时刻对应的地层反射系数相关.传播子波可用于反褶积或后续反演的输入，局部观测子波则可用作解释过程中的辅助工具，因为它是对目标地质构造或目标层位的反映.
对于本文的时变子波估计方法，要求测井数据和地震数据进行必要的预处理.实际的测井数据既包含噪声，也包含测量误差，其中测井数据合成反射序列时的时间误差会对子波估计结果带来较大的影响(White，1997; 杨培杰和印兴耀，2008)，因此必须对声波测井数据进行仔细校正.对于地震数据，需要进行前期的静校正、去噪、偏移和叠加等处理流程.
6 结论在含有随机噪声的情况下，该方法无需任何先验信息，即可通过谱分析法估计时变子波振幅，通过偏度属性最大化方法估计时变子波相位信息.将两者结合即可得到估计的时变传播子波.为判断估计子波的准确性，使用子波与测井合成反射系数褶积得到合成记录与实际地震记录进行对比，同时使用拟合度(PEP)和准确度(NMSE)两个参数可以对井震匹配的结果进行定量分析.随后可以利用估计时变子波进行时变反褶积，提高剖面保真度.我们将该方法分别用于模型算例和实际数据中，结果证明该方法可较好地处理非平稳地震信号，并从中提取时变子波信息.
| Berkhout A J. 1977. Least-squares inverse filtering and wavelet deconvolution. Geophysics, 42(7): 1369-1383. DOI:10.1190/1.1440798 | |
| Claerbout J F. 1977. Parsimonious deconvolution. Stanford Exploration Project, 13: 1-9. | |
| Donoho D. 1981. On minimum entropy deconvolution.//Applied Time Series Analysis Ⅱ. New York:Academic Press, 565-608. | |
| Edgar J A, Van der Baan M. 2011. How reliable is statistical wavelet estimation?. Geophysics, 76(4): V59-V68. DOI:10.1190/1.3587220 | |
| Feng W, Hu T Y, Zhang Y, et al. 2015. Time-varying wavelet estimation by hybrid method.//77th EAGE Conference and Exhibition. EAGE, doi:10.3997/2214-4609.201413421. | |
| Fomel S. 2007a. Local seismic attributes. Geophysics, 72(3): A29-A33. DOI:10.1190/1.2437573 | |
| Fomel S. 2007b. Shaping regularization in geophysical-estimation problems. Geophysics, 72(2): R29-R36. DOI:10.1190/1.2433716 | |
| Fomel S, Van der Baan M. 2014. Local skewness attribute as a seismic phase detector. Interpretation, 2(1): SA49-SA56. DOI:10.1190/INT-2013-0080.1 | |
| Gao J H, Wang L L, Zhao W. 2009. Enhancing resolution of seismic traces based on the changing wavelet model of the seismogram. Chinese J. Geophys. (in Chinese), 52(5): 1289-1300. DOI:10.3969/j.issn.0001-5733.2009.05.018 | |
| Gao J H, Zhang B. 2010. Estimation of seismic wavelets based on the multivariate scale mixture of Gaussians model. Entropy, 12(1): 14-33. DOI:10.3390/e12010014 | |
| Han L, Liu C C, Yuan S Y. 2015. Can we use wavelet phase change due to attenuation for hydrocarbon detection?//85th Annual International Meeting, SEG, Expanded Abstracts. SEG, 2962-2966, doi:10.1190/segam2015-5890451.1. | |
| Levy S, Oldenburg D W. 1987. Automatic phase correction of common-midpoint stacked data. Geophysics, 52(1): 51-59. DOI:10.1190/1.1442240 | |
| Li G F, Liao Q J, Wang S X, et al. 2008. Discussions about horizon calibration based on well-log synthetic seismogram. Geophysical Prospecting for Petroleum (in Chinese), 47(2): 145-149. | |
| Liu Y J, Li Z C, Wu D, et al. 2014. Well constrained non-stationary phase correction. Chinese J. Geophys. (in Chinese), 57(1): 310-319. DOI:10.6038/cjg20140127 | |
| Longbottom J, Walden A T, White R E. 1988. Principles and application of maximum kurtosis phase estimation. Geophysical Prospecting, 36(2): 115-138. DOI:10.1111/j.1365-2478.1988.tb02155.x | |
| Luo C M, Wang S X, Yuan S Y. 2014. Effect of inaccurate wavelet phase on prestack waveform inversion. Applied Geophysics, 11(4): 479-488. DOI:10.1007/s11770-014-0453-1 | |
| Ma H D, White R E, Hu T Y. 2010. Wavelet estimation by matching well-log, VSP, and surface-seismic data. Applied Geophysics, 7(4): 384-391. DOI:10.1007/s11770-010-0264-y | |
| Ma M, Wang S X, Yuan S Y, et al. 2015. The comparison of skewness and kurtosis criteria for wavelet phase estimation.//85th Annual International Meeting, SEG, Expanded Abstracts. SEG, 5164-5168, doi:10.1190/segam2015-5826646.1. | |
| Ooe M, Ulrych T J. 1979. Minimum entropy deconvolution with an exponential transformation. Geophysical Prospecting, 27(2): 458-473. DOI:10.1111/j.1365-2478.1979.tb00979.x | |
| Robinson E A, Treitel S. 2000. Geophysical Signal Analysis. USA:Society of Exploration Geophysicists. | |
| Saggaf M M, Robinson E A. 2000. A unified framework for the deconvolution of traces of nonwhite reflectivity. Geophysics, 65(5): 1660-1676. DOI:10.1190/1.1444854 | |
| Tygel M, Bleistein N. 2000. An introduction to this special section:wavelet estimation. The Leading Edge, 19(1): 37-37. DOI:10.1190/1.1438448 | |
| Van der Baan M. 2008. Time-varying wavelet estimation and deconvolution by kurtosis maximization. Geophysics, 73(2): V11-V18. DOI:10.1190/1.2831936 | |
| Van der Baan M, Fomel S. 2009. Nonstationary phase estimation using regularized local kurtosis maximization. Geophysics, 74(6): A75-A80. DOI:10.1190/1.3213533 | |
| Van der Baan M, Fomel S, Perz M. 2010. Nonstationary phase estimation:A tool for seismic interpretation?. The Leading Edge, 29(9): 1020-1026. DOI:10.1190/1.3485762 | |
| Van der Baan M. 2012. Bandwidth enhancement:Inverse Q filtering or time-varying Wiener deconvolution?. Geophysics, 77(4): V133-V142. DOI:10.1190/geo2011-0500.1 | |
| Walden A T, Hosken J W J. 1985. An investigation of the spectral properties of primary reflection coefficients. Geophysical Prospecting, 33(3): 400-435. DOI:10.1111/j.1365-2478.1985.tb00443.x | |
| Walden A T, Hosken J W J. 1986. The nature of the non-Gaussianity of primary reflection coefficients and its significance for deconvolution. Geophysical Prospecting, 34(7): 1038-1066. DOI:10.1111/j.1365-2478.1986.tb00512.x | |
| Walden A T, White R E. 1984. On errors of fit and accuracy in matching synthetic seismograms and seismic traces. Geophysical Prospecting, 32(5): 871-891. DOI:10.1111/j.1365-2478.1984.tb00744.x | |
| Walden A T, White R E. 1998. Seismic wavelet estimation:a frequency domain solution to a geophysical noisy input-output problem. IEEE Transactions on Geoscience and Remote Sensing, 36(1): 287-297. DOI:10.1109/36.655337 | |
| Wang S X, Yuan S Y, Ma M, et al. 2015. Wavelet phase estimation using ant colony optimization algorithm. Journal of Applied Geophysics, 122: 159-166. DOI:10.1016/j.jappgeo.2015.09.013 | |
| Wang Z G, Zhang B, Gao J H. 2014. The residual phase estimation of a seismic wavelet using a Rényi divergence-based criterion. Journal of Applied Geophysics, 106: 96-105. DOI:10.1016/j.jappgeo.2014.04.008 | |
| White R E. 1988. Maximum kurtosis phase correction. Geophysical Journal International, 95(2): 371-389. DOI:10.1111/j.1365-246X.1988.tb00475.x | |
| White R E. 1997. The accuracy of well ties:Practical procedures and examples.//69th Annual International Meeting, SEG, Expanded Abstracts. SEG, 816-819, doi:10.1190/1.1886137. | |
| White R E, Simm R. 2003. Tutorial:Good practice in well ties. First Break, 21(10): 75-83. | |
| White R E, Naeini E Z. 2014. Broad-band well tie.//76th EAGE Conference and Exhibition. EAGE, doi:10.3997/2214-4609.20140749. | |
| Wiggins R A. 1978. Minimum entropy deconvolution. Geoexploration, 16(1-2): 21-35. DOI:10.1016/0016-7142(78)90005-4 | |
| Yang P J, Yin X Y. 2008. Summary of seismic wavelet pick-up. Oil Geophysical Prospecting (in Chinese), 43(1): 123-128. DOI:10.13810/j.cnki.issn.1000-7210.2008.01.020 | |
| Yu Y C, Wang S X, Yuan S Y, et al. 2011. Phase estimation in bispectral domain based on conformal mapping and applications in seismic wavelet estimation. Applied Geophysics, 8(1): 36-47. DOI:10.1007/s11770-011-0269-1 | |
| Yu Y C, Wang S X, Yuan S Y, et al. 2012. Phase spectrum estimation of the seismic wavelet based on a criterion function. Petroleum Science, 9(2): 170-181. DOI:10.1007/s12182-012-0197-6 | |
| Yuan S Y, Wang S X. 2011. Influence of inaccurate wavelet phase estimation on seismic inversion. Applied Geophysics, 8(1): 48-59. DOI:10.1007/s11770-011-0273-5 | |
| Yuan S Y, Wang S X. 2013. Spectral sparse Bayesian learning reflectivity inversion. Geophysical Prospecting, 61(4): 735-746. DOI:10.1111/1365-2478.12000 | |
| 高静怀, 汪玲玲, 赵伟. 2009. 基于反射地震记录变子波模型提高地震记录分辨率. 地球物理学报, 52(5): 1289–1300. DOI:10.3969/j.issn.0001-5733.2009.05.018 | |
| 李国发, 廖前进, 王尚旭, 等. 2008. 合成地震记录层位标定若干问题的探讨. 石油物探, 47(2): 145–149. | |
| 刘玉金, 李振春, 吴丹, 等. 2014. 井约束非稳态相位校正方法. 地球物理学报, 57(1): 310–319. DOI:10.6038/cjg20140127 | |
| 杨培杰, 印兴耀. 2008. 地震子波提取方法综述. 石油地球物理勘探, 43(1): 123–128. DOI:10.13810/j.cnki.issn.1000-7210.2008.01.020 | |