基于转换波地震数据求取高精度反射系数的方法
【专利摘要】本发明提供了一种基于转换波地震数据求取高精度反射系数的方法，所述方法包括：采集转换波地震数据，进行预处理和叠前偏移后得到转换波偏移剖面；在偏移剖面上建立稀疏约束和波阻抗约束的多薄层目标函数；对多薄层目标函数利用频域映射与L1范数联合优化算法求取高精度反射系数。本发明的方法能够基于转换波地震数据构建更加合理的多薄层目标函数并获得高精度反射系数，提高了目标函数求解精度和抗噪性，而且所获得的高精度反射系数能够为后续高分辨率剖面重构奠定坚实基础。
【专利说明】
基于转换波地震数据求取高精度反射系数的方法
技术领域
[0001 ]本发明属于转换波地震资料处理技术领域，具体来讲，涉及一种基于转换波地震 数据求取高精度反射系数的方法。
【背景技术】
[0002] 目前，相干体裂缝特征描述是探寻裂缝性油气藏的一项重要手段，而获取裂缝预 测所需高分辨率数据体、裂缝参数的求取一直是裂缝油气藏勘探中重点研究的问题，尤其 转换波低信噪比和分辨率的特点，需要提高其分辨率，便于后续解释的薄层分辨。但是，最 小相位和白噪假设为前提的常规高分辨率处理方法，在主频带范围限制下，提高转换波分 辨率的能力有限，无法达到后续裂缝预测分辨多薄层最小厚度的要求，所以需要采取新的 高分辨率处理思路，g卩：转换波高精度反射系数求取与高频子波高分辨率重构的新思路，该 思路不受任何假设和频带限制，能完成常规手段无法分辨多薄层最小厚度的要求，获取裂 缝预测满意的高品质剖面。
[0003] 通常，反射系数的求取需要构建多薄层目标函数，然而，由于常规方法构建的目标 函数存在诸多不合理之处，故而，在求取过程中的反射系数存在极值和位置、大小上的误 差，这种误差导致后续高分辨率剖面的错像。
[0004] 此外，目标函数的求解方法有很多，比如模拟退火法、共辄梯度法、匹配追踪法等 这些算法，虽然都能求取目标函数的反射系数，但精度和稳定性较差，且对新构建的目标函 数适应性较差。
[0005] 综上，亟需一种能够基于转换波地震数据求取高精度反射系数的方法。
【发明内容】
[0006] 本发明的目的在于解决现有技术存在的上述不足中的至少一项。例如，本发明的 目的之一在于改善现有转换波地震数据领域中为求取反射系数所构建的多薄层目标函数 的合理性，并求解出高精度反射系数。
[0007] 为了实现上述目的，本发明提供一种基于转换波地震数据求取高精度反射系数的 方法。所述方法包括:采集转换波地震数据，进行预处理和叠前偏移后得到转换波偏移剖 面，在转换波偏移剖面上取11道地震记录奴1，〇，1=1，2，一，11;在偏移剖面上建立稀疏约束 和波阻抗约束的多薄层目标函数;对多薄层目标函数利用频域映射与L1范数联合优化算法 求取高精度反射系数；
[0008] 其中，
[0009] 所述目标函数为下式（1):
[0011]其中，r表不反射系数，曰^代表反射系数的偶分量比例，a。代表反射系数的奇分量 比例，&为反射系数的偶分量，r。为反射系数的奇分量，A为变换矩阵，b为地震数据和子波相 关矩阵，A为稀疏因子，#为阻抗因子，C为积分算子矩阵，％为纵波波阻抗矩阵。
[0012] 与现有技术相比，本发明的有益效果包括:能够基于转换波地震数据构建更加合 理的多薄层目标函数，并能有效提高目标函数求解精度和抗噪性，能够获得高精度反射系 数，进而能够为后续高分辨率剖面重构奠定坚实基础。
【具体实施方式】
[0013] 在下文中，将结合示例性实施例来详细说明本发明的基于转换波地震数据求取高 精度反射系数的方法。
[0014] 在本发明的一个示例性实施例中，基于转换波地震数据求取高精度反射系数的方 法可通过以下步骤来实现：
[0015] (1)对采集到的转换波地震数据，经过常规预处理和叠前偏移后得到的转换波偏 移剖面，取n道地震记录x(x，t)，x = l，2,…，n;其中，x表示道数，t表示时间。
[0016] (2)在偏移剖面上建立由稀疏约束和波阻抗约束双重约束的多薄层目标函数； [0017]步骤(2)可按下列方式完成：
[0018] (a)对于不含噪声的褶积模型，地震合成记录在频率域中可表示为：
[0019] S(f)=ff(f) ? R(f) (1)
[0020] 其中，S(f)表示频率域下地震记录（即，上面的奴1，〇)，1(〇表示频率域下地震子 波，R(f)表示频率域下反射系数;f为频率，时间t经过傅里叶变换后变换到频率f中。
[0021] 利用频率域下地震记录与地震子波的丰富信息，利用合适的常规反演算法(如，共 辄梯度法、匹配追踪法、模拟退火等)消除子波（即，地震子波)影响(如，子波主频、旁瓣干扰 等会影响反演求取反射系数精度的内容），从记录中得到所需的反射系数(这里，为了得到 反射系数R(f))。利用反射系数奇偶分解原理，可得到频率域下目标函数表达式：
[0022]
[0023] 其中，ae代表反射系数R( f)的偶分量比例，a。代表反射系数R( f)的奇分量比例;Re 表示R(f)的实部，Im表示R(f)的虚部。
[0024] (b)为建立更符合实际地层情况的反射系数模型，对步骤(a)中的多薄层反射系数 模型下的目标函数进行推导或改进。
[0025] 由奇偶反射系数分解原理，得多薄层反射系数奇偶分量可表示为：
[0026]
[0027] re( t)代表反射系数r (t)的偶分量，r。(t)代表反射系数r (t)的奇分量，t表示时间， T表示时间厚度;r(t)为时间域反射系数，r(t)经过傅里叶变换后到上面的频率域反射系数 R(f);g(t)当前分析位置反射系数，ge(t)代表g(t)的偶分量，g〇(t)代表g(t)的奇分量。
[0028] 对式(3)进行傅里叶变换，得到其频率域实部与虚部表达式：
[0029]
[0030] 其中，31为圆周率。
[0031] 利用式(2)和式(4)得到多薄层反射模型的目标函数表达式为：
[0032]
[0033]其中，re代表反射系数r(t)的偶分量，r。代表反射系数r(t)的奇分量，A t表不时移 量;i无物理意义，是指数与正弦函数变化标识符。
[0034]为了便于求解，将式(5)可写成矩阵形式：
[0038] 其中，t为时间，t、t2、……、tN分别表示第1、2、……川个时间，f为频率，h、 f2、……、fM分别表示第1、2、……、M个分析频率范围内的频率，Ti、T2……、Tn/2分别表示第1、 2、……、N/2个时间厚度，At表示时移量
[0039] (c)由于该目标函数(6)(即，式(6))有多个全局最优解，为了降低求解问题时的多 解性和不确定性，需要加上稀疏化约束，将目标函数改为范数求解形式为：
[0041] 其中，A为稀疏因子，A-般在〇-1范围内取值，| | | U表示L1范数，| |表示欧几里得 (Euclidean)范数的平方。
[0042] (d)由于波阻抗可看作是反射系数对时间的积分，因此发明人把它引入到目标函 数中作为反射系数r的先验约束。
[0043] 由于转换波相对于纵波来说其信噪比和分辨率较低，为了提高反射系数精度，本 发明使用分辨率更高的纵波波阻抗作为转换波的约束条件。
[0044] 若I (t)为纵波波阻抗值，则纵波相对波阻抗e *为：
[0046] 其中，I(to)表示纵波初始波阻抗值，r为反射系数，n表示在时间to与时间t中取值 符合，to表示初始时间，tM-l表示结束时间。
[0047] 把上式(8)写成矩阵形式：
[0048] |P = Cr (9)
[0049] 其中g = [gv，...，，C为积分算子矩阵，其离散形式可表示为：
[0051]则由最小平方定义波阻抗约束得到新的目标函数为：
[0053] 其中，r表示反射系数，ae和a。分别代表反射系数r(t)的偶奇分量比例，是反射系 数的偶分量，r。是反射系数的奇分量;A为变换矩阵，b是地震数据和子波相关矩阵，A是稀疏 因子，3是阻抗因子(其大于零即可），C为积分算子矩阵，％是纵波波阻抗矩阵）
[0054] (3)对于步骤(2)中的式（11)，利用频域映射与L1范数联合优化算法求取高精度反 射系数。
[0055] 本步骤(3)可按下列方式完成：
[0056] (a)对式（11)进行求解中，由于波阻抗约束项是由测井信号得到的已知先验约束 条件，所以式(11)的求解类似于求矩阵方程Ax = b在含有噪声情况下的稀疏解。
[0057] 可将目标函数(11)的求解问题转化为下式：
[0059]其中，u为极小权重值，G为复数集合;本式（12)中，A表不Ae时b表不b e，X表不re; A表 不A。时b表不b。，X表不r。。
[0060] 如果引入变量qGGm，那么式（12)等价于式（13):
[0062]式(13)对应的增广拉格朗日子问题可表示为：
[0064]其中，yT表示拉格朗日乘子的共辄转置运算，0为罚参数。
[0065] 利用式(14)进行频域映射与L1范数联合优化算法实现步骤如下：
[0066] 1)令k = 0对qk，xk，yk赋初始值，并给定常数值(这里，为人为给定，可为大于 零的常值），然后进行以下步骤2)至5)的算法运算，如果满足终止准则(例如，终止准则可以 为:最优解值或者迭代次数），则算法完成运算;否则，进行步骤2);
[0067 ] 2)令x = xk，y = yk，求解r的子问题得到qk+1:
[0069] 3)令q = qk+1，y = yk贝lj关于x的极小化问题式(14)等价于式(16):
[0071] 那么可以通过式(17)近似求解来完成对式(16)的精确求解：
[0073]其中，t为大于零的邻近参数(这里，t为人工给定，其大于零即可），gk为x = xk时二 次项的梯度，那么式(17)通过下式(18)求解得到3+1:
[0075] 其中，。表示逐元素相乘，Shrink(，）表示一维收缩算子，sign表示符号函数。
[0076] 4)令x = xk+1，q = qk+1，求解y 的子问题得到 yk+1:
[0077] yk+1 = yk-yP(Axk+1 +qk+1-b) (19)
[0078] 其中，y为大于零的常数(这里，y人工给定，其大于零即可）。
[0079] 5)令k = k+l，重复步骤1)~步骤4)进行迭代运算。
[0080] 本发明的方法能够基于转换波地震数据构建更加合理的多薄层目标函数并获得 高精度反射系数，提高了目标函数求解精度和抗噪性，而且所获得的高精度反射系数能够 为后续高分辨率剖面重构奠定坚实基础。
[0081] 尽管上面已经结合示例性实施例描述了本发明，但是本领域普通技术人员应该清 楚，在不脱离权利要求的精神和范围的情况下，可以对上述实施例进行各种修改。
【主权项】
1. 一种基于转换波地震数据求取高精度反射系数的方法，所述方法包括:采集转换波 地震数据，进行预处理和叠前偏移后得到转换波偏移剖面，在转换波偏移剖面上取n道地震 记录x(x，t)，x = l，2,…，n; 在偏移剖面上建立稀疏约束和波阻抗约束的多薄层目标函数； 对多薄层目标函数利用频域映射与L1范数联合优化算法求取高精度反射系数； 其中， 所述目标函数为下式(1):其中，r表示反射系数，ae代表反射系数的偶分量比例，a。代表反射系数的奇分量比例， ^为反射系数的偶分量，r。为反射系数的奇分量，A为变换矩阵，b为地震数据和子波相关矩 阵，A为稀疏因子，3为阻抗因子，C为积分算子矩阵，％为纵波波阻抗矩阵。2. 根据权利要求1所述的基于转换波地震数据求取高精度反射系数的方法，其中，其中，t为时间，t、t2、……、tN分别表示第1、2、……、N个时间，f为频率，f i、f2、……、fM 分别表示第1、2、……、M个分析频率范围内的频率，Ti、T2……、TN/2分别表示第1、2、……、N/2 个时间厚度，△ t表示时移量，Re表示频率域下反射系数R(f)的实部，Im表示频率域下反射 系数R(f)的虚部，S表示频率域下地震记录x(x，t)，W表示频率域下地震子波； E = ￡, \，其中，to表示初始时间，tM-i表示结束时间； ~P L. t.o ，j.J， C的离散形式表示为3.如权利要求1所述的基于转换波地震数据求取高精度反射系数的方法，其中，所述求 取高精度反射系数的步骤包括： A、 将多薄层目标函数式(1)的求解问题转化为式(2):其中，u为极小权重值，G为复数集合； B、 引入变量qGGm，则式(2)等价于式(3):C、 式(3)对应的增广拉格朗日子问题表示为式(4):，其中，yT表示乘子的 共辄转置运算，0为罚参数； D、利用式(4)进行频域映射与L1范数联合优化算法求取高精度反射系数。4.如权利要求3所述的基于转换波地震数据求取高精度反射系数的方法，其中，所述步 骤D包括步骤： 1) 令k = 0对rk，xk，yk赋初始值，并给定常数值，然后进行下面步骤2)至5)的算法运 算，若满足终止准则，则完成运算，否则，进行步骤2); 2) 令x = xk，y = yk，求解r的子问题得到是(5):3) 令r = rk+1，y = yk则关于x的极小化问题式⑷等价于式(6):那么可以通过式(7)近似求解来完成对式(6)的精确求解，式(7)为：其中，t为大于零的邻近参数，gk为1 =，时二次项的梯度，那么式(7)通过式(8)求解得 至ljxk+1，式⑶为：其中，。表示逐元素相乘，Shrink(，）表示一维收缩算子； 4) 令x = xk+1，r = rk+1，求解y的子问题得到式(9): yk+1=yk-yP(Axk+1+rk+1-b), 其中，Y为大于零的常数。 5) 令k = k+l，重复步骤1)~步骤4)进行迭代运算。
【文档编号】G01V1/30GK106054247SQ201610357140
【公开日】2016年10月26日
【申请日】2016年5月25日
【发明人】张华 , 何光明, 周强, 赵尧, 金德刚, 罗红明, 谢小翠
【申请人】中国石油集团川庆钻探工程有限公司地球物理勘探公司