基于投影菲涅尔带的三维起伏地表高斯束正演模拟方法
【专利摘要】本发明提供一种基于投影菲涅尔带的三维起伏地表高斯束正演模拟方法，该基于投影菲涅尔带的三维起伏地表高斯束正演模拟方法包括：步骤1，建立起伏地表的三维地质模型；步骤2，进行运动学射线追踪；步骤3，进行动力学射线追踪；步骤4，在三维介质中用投影菲涅尔带椭圆来约束高斯束在射线终点处能量的分布；以及步骤5，采用高斯波包法合成地震记录。本方法中的基于投影菲涅尔带的三维起伏地表高斯束正演模拟方法在中深部地层具有更好的保幅性，对起伏地表和复杂模型具有更好的适应性。
【专利说明】
基于投影菲涅尔带的H维起伏地表高斯束正演模拟方法
技术领域
[0001] 本发明设及地震资料处理，特别是设及到一种基于投影菲涅尔带的=维起伏地表 高斯束演模拟方法。
【背景技术】
[0002] 近年来，国内外许多学者对起伏地表下的正演模拟和高斯束正演模拟进行了很多 研究。首先，在不同介质中应用方面：Cerveny(1982)在引入高斯束理论后，将其应用于单 层二维非均匀介质、侧向变速的层状介质、=维弹性非均匀介质W及=维弹性侧向变速层 状介质的正演模拟中，Muler (1988)将其应用于二维非均匀吸收介质的正演模拟中，徐盛岩 等（1988)将高斯束理论用于二维粘弹性介质地震波场的模拟。其次，在高斯束参数选择方 面：Cerveny(1982)给出了最初的动力学追踪初始参数；Muler(1984)给出了二维介质中几 种初始参数的选择，并讨论了其对高斯束性质的影响；George(1987)对束参数的选择给出 了一些改进，确保了合成记录的稳定性。化UZ (2012)首先提出在二维介质中利用投影菲涅 尔带半径来约束高斯束的有效半宽度，并给出了相应的初始参数，但没有给出投影菲涅尔 带半径的计算方法。近年来，基于起伏地表的=维高斯束正演方法还处于起步阶段。为此 我们发明了一种新的基于投影菲涅尔带的=维起伏地表高斯束演模拟方法，解决了 W上技 术问题。
【发明内容】
[0003] 本发明的目的是提供一种针对起伏地表条件下的=维高斯束正演模拟的基于投 影菲涅尔带的=维起伏地表高斯束演模拟方法。
[0004] 本发明的目的可通过如下技术措施来实现：基于投影菲涅尔带的=维起伏地表高 斯束演模拟方法，该基于投影菲涅尔带的=维起伏地表高斯束演模拟方法包括：步骤1，建 立起伏地表的=维地质模型；步骤2,进行运动学射线追踪；步骤3,进行动力学射线追踪； 步骤4,在=维介质中用投影菲涅尔带楠圆来约束高斯束在射线终点处能量的分布；W及 步骤5,采用高斯波包法合成地震记录。 阳0化]本发明的目的还可通过如下技术措施来实现：
[0006] 在步骤2中，运动学射线追踪在S维坐标系下满足如下方程：
[0007] Cl) 阳00引其中：T为旅行时，V为中屯、射线的P波速度，Vi (i = X，y，Z)为中屯、射线处速度 函数的偏导数，a和0分别为其坐标系下的倾角与方位角；在层内采用四阶龙格-库塔 法进行求解各时间步的坐标，在界面上利用二分法求出交点的近似坐标，并利用矢量Snell 定律计算反射或透射方向。
[0009] 在步骤3中，将动力学追踪方程：
[0012] 化为线性动力学追踪方程：阳01引 dXi/ds = HX-. (i = 1，2) (3)
[0010]
[0011]
[0014] 其中
Xi为X的列向量，
，0为2X2零矩阵，I为2X2 单位矩阵；（S，m，n)为=维中屯、射线坐标，V为速度二阶偏导数矩阵，Q，P为动力射线追踪 方程的参数；
[0015] 在给定四阶单位矩阵的初始条件下，得到X的通解：
[0016] X=n(s).C (4)
[0017] 其中：n(s)为4X4传播矩阵，其列向量为方程（3)四组线性独立的解，C为一个 2X4复值初始参数矩阵； 阳〇1引把n (S)和C分块如下：
[0019]
(5)
[0020] 其中：町郝Ci(i，j = 1，。为2X2子矩阵；
[0021] 运样，根据x = 得到高斯束的动力学参数矩阵： k.尸.J
[0022]
(技)
[0023] 由于射线经过界面时传播矩阵n (S)发生突变，需要在界面处重新计算n (S);
[0024] 在射线与界面交点0点处边界条件为：
阳0巧]
[0026]
[0027]
[0028]
[0029] Eu (0) =-sinisV 2 (0) [ (1+cos2;Ls) v,zi_ e cosisSinisV,z3]
[0030] Ei2 (0) = E21 (0) =-sinisV 2 (0) v,z2 (9)
[0031]
[0032]
[0033]
[0034] 其中：S。为炮点，0为射线与界面交点，资为下一个界面交点，GT为G的转置， 在=*扭/如.1^)，为方向指数山为入射角，1历反射角或透射角；0为界面的曲率矩阵；zl、z2 和z3为局部笛卡尔坐标系中的=个分量；K为射线中屯、坐标系中e。和局部笛卡尔坐标系 中fa的夹角；当发生透射时取等式（9)上面的符号，发生反射时取等式巧）下面的符号。
[0035] 征在于，在步骤4中，在=维介质中用投影菲涅尔带楠圆来约束高斯束在射线终 点处能量的分布；
[0036] 当中屯、射线终点处高斯束的有效半宽度楠圆与射线的投影菲涅尔带楠圆一致时， 可得：
[0037]
( ]-0)
[00測其中：eig为特征值，Gr为中屯、射线终点，化为投影菲涅尔带矩阵，其值可W通过 经典射线传播矩阵n (S)、面面传播矩阵T (S) W及两者转换关系求出；
[0039] 根据传播矩阵n (S)的辛属性及数学推导可得到：
[0040]
( 11) 阳〇4U 其中：入1、入2(入入2)为口11的特征值，n 1、n2(nn 2)为n。的特征值， (1、^ 2)为投影菲涅尔带矩阵化的特征值。 阳0创在步骤5中，记W初始角（a 1，P j)出发在检波点R处t时刻高斯波包为 g(R，t，a 1，0 j)，其近似表达式为： W 创
（12); W44] 其中：fm、丫、V为高斯子波参数，R'为R在中屯、射线上的投影，
qT= (q 1,92), Qi, 92为射线在R点的坐标，A为振幅，在层状介质中
;P为密 度，V为速度，Ri为经过第i个界面反射系数，N表示射线经过界面的个数，"~"表示生成射 线一侧的量值；O = ? /2 JT . I det怕巧'））r2 ? {-det [M很）-Re (M巧'））]}，为能量叠加的 权系数；
[0045] 那么，检波点R处离散的能量叠加表达式为： W46]
( 1恩)
[0047] 其中：1，J分别表示射线在倾角和方位角上离散的个数，A a，A 0为射线间隔。
[0048] 在步骤1中，结合地震和测井资料及其解释结果包括层位、断层等，建立模型框 架，之后加上各层的速度参数，建立起伏地表的=维地质模型。
[0049] 本发明中的基于投影菲涅尔带的S维起伏地表高斯束演模拟方法，引入了S维投 影菲涅尔带的思想，发展了基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，在中 深部地层具有更好的保幅性，对起伏地表和复杂模型具有更好的适应性。相对于常规方法， 在投影法菲涅尔带约束下高斯束的传播更稳定，能量聚集性更好；本方法使高斯束能量分 布在投影菲涅尔带内，赋予了能量分布范围明确的物理意义，使得基于渐进射线理论的高 斯束更符合波动理论；与常规初始参数正演结果相比，本方法在中深部地层具有更好的保 幅性，对起伏地表和复杂模型具有更好的适应性。
【附图说明】
[0050] 图I为本发明的基于投影菲涅尔带的=维起伏地表高斯束演模拟方法的一具体 实施例的流程图；
[0051] 图2为本发明的一具体实施例中用于测试的简单=维水平地表层状模型；
[0052] 图3为本发明的一具体实施例中利用两种初始参数对=维水平地表层状模型进 行=维高斯束正演模拟得到的结果的示意图；
[0053] 图4为本发明的一具体实施例中用于测试的经典起伏地表台地模型；
[0054] 图5为本发明的一具体实施例中起伏地表台地模型盲区示意图；
[0055] 图6为本发明的一具体实施例中利用两种初始参数对起伏地表台地模型进行= 维高斯束正演模拟得到的结果的示意图；
[0056] 图7为本发明的一具体实施例中用于测试的起伏地表复杂模型；
[0057] 图8为本发明的一具体实施例中利用两种初始参数对起伏地表复杂模型进行= 维高斯束正演模拟得到的结果的示意图。
【具体实施方式】
[0058] 为使本发明的上述和其他目的、特征和优点能更明显易懂，下文特举出较佳实施 例，并配合所附图式，作详细说明如下：
[0059] 如图1所示，图1为本发明的基于投影菲涅尔带的=维起伏地表高斯束演模拟方 法的流程图。
[0060] 在步骤101，建立起伏地表的=维地质模型。结合地震和测井资料及其解释结果包 括层位、断层等，建立模型框架，之后加上各层的速度参数，建立起伏地表的=维地质模型。 流程进入到步骤102。
[0061] 在步骤102,进行运动学射线追踪。
[0062] 运动学射线追踪在=维坐标系下满足如下方程：
[0063]
( 1) W64] 其中：T为旅行时，V为中屯、射线的P波速度，Vi (i = X, y, Z)为中屯、射线处速度 函数的偏导数，a和0分别为其坐标系下的倾角与方位角；在层内采用四阶龙格-库塔 法进行求解各时间步的坐标，在界面上利用二分法求出交点的近似坐标，并利用矢量Snell 定律计算反射或透射方向。流程进入到步骤103。
[00化]在步骤103,进行动力学射线追踪。
[0066] 将动力学追踪方程：
[0067]
[0069] 化为线性动力学追踪方程：[0070] dX,/ds = HX, (i = 1, 2) (3)
[0068] 阳07U 其中；
Xi为X的列向量，
：，0为2X2零矩阵，I为2X2 单位矩阵；（s，m，n)为=维中屯、射线坐标，V为速度二阶偏导数矩阵，Q，P为动力射线追踪方 程的参数；
[0072] 在给定四阶单位矩阵的初始条件下，得到X的通解：
[0073] X = n (S) . C (4)
[0074] 其中：n(s)为4X4传播矩阵，其列向量为方程做四组线性独立的解，C为一个 2X4复值初始参数矩阵； 阳〇7引 把n (S)和C分块如下：
[0076]
(百) 阳077] 其中：町郝Ci(i，j = 1，。为2乂2子矩阵；
[0078] 运样，根据X = ，得到高斯束的动力学参数矩阵：
[0079]
( 6)
[0080] 由于射线经过界面时传播矩阵n (S)发生突变，需要在界面处重新计算n (S);
[0081] 在射线与界面交点0点处边界条件为：
[0082]
[0083]
[0084] 阳O化]
[0086]
[00川其中：S。为炮点，0为射线与界面交点，贫为下一个界面交点，护为G的转置，
[0087]
[0088]
[0089]
[0090] 忠=該药據為，为方向指数；is为入射角，i K为反射角或透射角；D为界面的曲率矩阵；zl、z2 和z3为局部笛卡尔坐标系中的=个分量；K为射线中屯、坐标系中e。和局部笛卡尔坐标系 中？2的夹角；
[0092] 当发生透射时取等式（9)上面的符号，发生反射时取等式巧）下面的符号。流程 进入到步骤104。
[0093] 在步骤104,在=维介质中用投影菲涅尔带楠圆来约束高斯束在射线终点处能量 的分布。
[0094] 当中屯、射线终点处高斯束的有效半宽度楠圆与射线的投影菲涅尔带楠圆一致时， 可得： 阳0巧]
(…
[0096] 其中：eig为特征值，G为中屯、射线终点，化为投影菲涅尔带矩阵，其值可W通过经 典射线传播矩阵n (S)、面面传播矩阵T (S) W及两者转换关系求出。
[0097] 根据传播矩阵n (S)的辛属性及数学推导可得到：
[009引
.'11)
[0099] 其中：入1、入2(入入2)为n。的特征值，n 1、n2(ni> n 2)为口12的特征值， (1、C 2)为投影菲涅尔带矩阵化的特征值。流程进入到步骤105。
[0100] 在步骤105,高斯波包法合成地震记录。
[0101] 记W初始角（Qi，0,)出发在检波点R处t时刻高斯波包为g(R，t，Qi, 0,)，其近 似表达式为：
[0102]
(12) 阳10引其中：fm、丫、V为高斯子波参数，R'为R在中屯、射线上的投影，
qT= (q 1,92)，Qi, Q2为射线在R点的坐标，A为振幅，在层状介质中
P为密 度，V为速度，Ri为经过第i个界面反射系数，N表示射线经过界面的个数，"~"表示生成射 线一侧的量值；O = ? /2 JT . I det怕巧'））r2 ? {-det [M很）-Re (M巧'））]}，为能量叠加的 权系数；
[0104] 那么，检波点R处离散的能量叠加表达式为：
[0105]
( 13) 阳106] 其中：1，J分别表示射线在倾角和方位角上离散的个数，A a，A 0为射线间隔。 阳107] 在应用本发明的一具体实施例中，包括：
[0108] 1)选取如图2所示的简单S维水平地表层状模型，模型尺寸为 4000m X 4000m X 4000m，各层速度如图 2 所示，Aa = 2°，AP= 10。；
[0109] 2)观测系统为：S炮四线制（如图2,实线为炮线，虚线为检波点线），道间距25m， 线间距500m，炮间距50m，炮线距500m，采样间隔4ms，记录时间4s，全排列接收。经试算得 到的射线路径如图1中红色线条所示；
[0110] 3)利用两种初始参数进行=维高斯束正演模拟得到的单炮记录如图3所示； 阳111] 在应用本发明的另一具体实施例中，包括：
[0112] 1)选取如图4所示的经典起伏地表台地模型，模型尺寸为4000mX4000mX2500m， 上层速度为2000m/s，下层速度为3000m/s，A a = 2°，A 0 = 2° ;
[011引。观测系统与图2所示模型相同，试算得到的射线路径如图4所示，得到射线追踪 的盲区如图5 ;
[0114] 3)利用两种初始参数进行=维高斯束正演模拟得到的单炮记录如图6所示； 阳115] 在应用本发明的一具体实施例中，包括：
[0116] 1)选取如图6所示的经典起伏地表台地模型，模型尺寸为4000mX4000mX4500m， 接收时间5s，速度分布如图7 ;
[0117] 2)观测系统与图2所示模型相同，试算得到的射线路径如图7所示；
[0118] 3)利用两种初始参数进行=维高斯束正演模拟得到的单炮记录如图8所示；
[0119] 对比图3中（a)、化）两者炮记录，（a)为常规初始参数模拟的炮记录，化）为投影 菲涅尔带约束下模拟的炮记录，可知，（1)两者同相轴的形态和位置是一致的，说明两者的 运动学特征一致；（2)在纵向上，a和b同相轴主要能量分布范围（图3中矩形框）从浅层 到深层都逐渐增大，但b比a增加的更慢些，运是因为投影菲涅尔带楠圆和常规有效半宽度 楠圆都随射线路径增大而增大，但前者增大的速度没有后者速度快；（3)在横向上，同一 深度b中同相轴主要能量分布范围比a小，并且深度越大，两者差异越大；经过此简单模型 的试算和分析，验证了基于投影菲涅尔带的高斯束正演方法的正确性，尤其是本方法还具 有一定的保幅性。
[0120] 图6(a)为常规初始参数模拟的炮记录，化）为投影菲涅尔带约束下模拟的炮记 录，由图5和图6可知，（1)对比两种初始参数合成的地震记录，可知本方法对起伏地表角点 模型是适应的；（2)在普通射线方法盲区内（如图5中的空白区域），两种初始参数试算的 炮记录都有一定的绕射能量（图6中箭头），但b中绕射能量范围（约25道）比a中（约 40道）的小，并且没有a中的能量干设现象（图6a矩形框），究其原因可能是在焦点处投 影菲涅尔带范围小于常规初始参数计算的有效范围，运样使得绕射能量分布在合理的范围 内。经过起伏地表台地模型的计算结果证实了本方法的对起伏地表和焦点绕射一定的适应 性和保幅性。 阳121] 图8(a)为常规初始参数模拟的炮记录，化）为投影菲涅尔带约束下模拟的炮记 录，对比图8中（a)、化）两者炮记录可知，本方法对复杂构造也具有较好的适应性，对绕射 能量具有很好的控制作用，因而模拟结果较常规=维高斯束正演模拟方法具有更好的保幅 性。
【主权项】
1. 基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，其特征在于，该基于投影 菲涅尔带的=维起伏地表高斯束正演模拟方法包括： 步骤1，建立起伏地表的=维地质模型； 步骤2,进行运动学射线追踪； 步骤3,进行动力学射线追踪； 步骤4,在=维介质中用投影菲涅尔带楠圆来约束高斯束在射线终点处能量的分布； W及 步骤5,采用高斯波包法合成地震记录。2. 根据权利要求1所述的基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，其 特征在于，在步骤2中，运动学射线追踪在=维坐标系下满足如下方程：(1) 其中：T为旅行时，V为中屯、射线的P波速度，Vi (i = X，y，Z)为中屯、射线处速度函数 的偏导数，a和0分别为其坐标系下的倾角与方位角；在层内采用四阶龙格-库塔法进行 求解各时间步的坐标，在界面上利用二分法求出交点的近似坐标，并利用矢量Snell定律 计算反射或透射方向。3. 根据权利要求2所述的基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，其 特征在于，在步骤3中，将动力学追踪方程：化为线性动力学追踪方程： dX,/ds = HX, (i = 1, 2) (3) (2) 其C 其中Xi为X的列向量，O为2X2零矩阵，！为2X2单 位矩阵；（S，m，n)为=维中屯、射线坐标，V为速度二阶偏导数矩阵，Q，P为动力射线追踪方程 的参数； 在给定四阶单位矩阵的初始条件下，得到X的通解： X = n (S) ? C (4) 其中：n(s)为4X4传播矩阵，其列向量为方程（3)四组线性独立的解，C为一个2X4 复值初始参数矩阵； 把n (S)和C分块如下：(5) 其中：n。和Ci(i，j = 1，。为2乂2子矩阵； 运样，根袁，得到高斯束的动力学参数矩阵：(邑） 由于射线经过界面时传播矩阵n (S)发生突变，需要在界面处重新计算n (S); 在射线与界面交点O点处边界条件为：其中：S。为炮点，O为射线与界面交点，6为下一个界面交点，护为G的转置， C =喊"(r'n)，为方向指数；is为入射角，i K为反射角或透射角；D为界面的曲率矩阵；zl、z2 和z3为局部笛卡尔坐标系中的=个分量；K为射线中屯、坐标系中e。和局部笛卡尔坐标系 中打的夹角；当发生透射时取等式（9)上面的符号，发生反射时取等式巧）下面的符号。4. 根据权利要求3所述的基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，其 特征在于，在步骤4中，在=维介质中用投影菲涅尔带楠圆来约束高斯束在射线终点处能 量的分布； 当中屯、射线终点处高斯束的有效半宽度楠圆与射线的投影菲涅尔带楠圆一致时，可 得：(10) 其中：eig为特征值，Gf为中屯、射线终点，化为投影菲涅尔带矩阵，其值可W通过经典 射线传播矩阵n (S)、面面传播矩阵T (S) W及两者转换关系求出； 根据传播矩阵n (S)的辛属性及数学推导可得到：(11) 其中：入1、入2(入入2)为口11的特征值，n 1、n2(ni> n 2)为口12的特征值，C 1、 C 2( C 自2)为投影菲涅尔带矩阵化的特征值。5. 根据权利要求4所述的基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法， 其特征在于，在步骤5中，记W初始角（a 1，Pj)出发在检波点R处t时刻高斯波包为 g(R，t，Qi, 0 ,)，其近似表达式为： 邑(R，t 曰 1，0J) = 23lfml 巫 A|exp{-巧 3lfm(t-0) 丫]2+(2 3lfmG/Y)2-23TfmG} (12) X COS (2 JT f*(t-白）+ V -arg (巫 A) + JT /? 其中：fm、丫、V为高斯子波参数，R'为R在中屯、射线上的投影，，，qT= (Qi, 92) ,Qi, 92为射线在 R 点的坐标，A为振幅，在层状介质中为密度，V为速度，Ri为经过第i个界面反射系数，N表示射线经过界面的个数，"~"表示生 成射线一侧的量值；O = ? /2 JT . I det怕巧'））r2 ? {-det [M很）-Re (M巧'））]}，为能量叠 加的权系数； 那么，检波点R处离散的能量叠加表达式为： (13) 其中：1，J分别表示射线在倾角和方位角上离散的个数，A a，A e为射线间隔。6.根据权利要求1所述的基于投影菲涅尔带的=维起伏地表高斯束正演模拟方法，其 特征在于，在步骤1中，结合地震和测井资料及其解释结果包括层位、断层等，建立模型框 架，之后加上各层的速度参数，建立起伏地表的=维地质模型。
【文档编号】G01V1/28GK105988135SQ201510073617
【公开日】2016年10月5日
【申请日】2015年2月11日
【发明人】杨晶, 魏福吉, 黄建平, 段卫星, 徐维秀, 宋建国, 崔汝国, 杨继东, 袁茂林
【申请人】中石化石油工程地球物理有限公司胜利分公司