还剩3页未读,继续阅读
文本内容:
旅行时线性插值射线追踪算法的改进0按列/按行扫描算法地震波传播理论的重要路径之一是,地震波传播理论也被用作研究不均匀介质和复杂结构的地震波传播问题利用旅行时信息进行层析成像最常用的正演方法就是射线追踪,该方法还被广泛用于偏移、地震波反演及其他正演模拟中,计算地震波的旅行时和波前Asawaka在文献[4]中指出,LTI算法向前处理过程中计算节点旅行时,应先按列扫描再进行按行扫描,但其并没有对按行扫描进行进一步的分析,也没有给出计算实例说明若不这样做会存在什么样的问题国内学者在对该算法进行研究时,也没有提到在按列扫描之后还要进行按行扫描,相反让人误以为LTI算法向前处理只需按列扫描即可但实际上,在向前处理过程中仅按列或按行扫描均没有考虑到逆向传播的射线,在处理速度变化剧烈的介质时,会导致所追踪出来的射线并不满足最小旅行时因此,笔者通过实例分析,说明LTI算法在向前处理过程中应采用全方位循环扫描的方法,并给出了其具体实现步骤;在保证计算精度的同时,在边界加入次生节点的方法可以大大提高该方法的计算效率1C点旅行时计算LTI方法存在一个线性假设前提,如图1所示:设A(x其中t若射线从A,B之间的D点穿过到达C点,则C点的旅行时t其中:v为介质的波传播速度;@为线段AB与AC之间的夹角式
(2)中C点的旅行时是关于r的函数根据最短旅行时原理,射线从D点传到C点应满足联立公式
(2)、
(3)可以解出将式
(4)代入式
(2)可得C点最小旅行时为由式
(4)可得D点的坐标为该方法可以分为向前处理和向后处理两步向前处理是从炮点开始计算炮点到各网格边界节点的最短旅行时,可由式
(5)求得;向后处理是从接收点开始,确定所有炮点一接收点的射线路径与各单元网格边界的交点,可由式
(6)求得2Iti算法的问题和改进
2.1旋转后过程中的射线路径对比LTI算法在向前处理过程中若只按列或按行扫描来计算节点旅行时没有考虑到逆向传播的射线以按列扫描为例,如图2a所示,S点表示炮点所在位置,R点是其中的一个计算节点,位于炮点右方LTI算法在向前处理过程中对R点最小旅行时的计算只考虑了炮点S,计算点R所在左半平面的所有射线;而在复杂的速度结构中,射线极有可能先绕到R点的右半平面再传到R点,如图2b所示同理,当计算点位于炮点左边,也存在同样的问题这样会导致在横向速度变化剧烈的介质中(即出现逆向传播),计算出来的节点旅行时不是最小,向后处理中追踪出来的射线路径也不满足最短旅行时原理,降低了该方法对复杂模型的适应能力图3为Asakawa和Kawanaka的一个经典模型为说明射线逆向传播的问题,将图3a的模型及观测系统整体顺时针旋转90°,如图4a所示,重新追踪出的初至波射线路径见图4b对比图3b,4b中的射线路径可以看出,旋转后追踪的
1、2两条射线路径有误由表1也可以看出,旋转后仅用按列扫描的LTI算法计算出来的
1、2两条射线的旅行时不是最小由图3b中的射线路径可知,正确情况下图4b中标示
1、2的2条射线应该是由左侧的高速层折射回来究其原因是只按列扫描的LTI法在向前处理过程中,计算震源左侧单元节点最小旅行时只考虑了来自右半平面的射线,并没有考虑来自左半平面的射线同理,只按行扫描也存在射线考虑不全的问题因此,笔者认为只按列或按行扫描的LTI法在向前处理过程并没有考虑来自各个方向的射线,只有采用全方位循环的方法才能确保计算得到的各节点的旅行时满足最小在具体实施中,先按LTI算法按列扫描计算各节点的旅行时,再按行扫描的方法来重新计算各网格节点旅行时;将按行扫描计算的节点旅行时与按列扫描计算的节点旅行时相比较,取其中旅行时最小的值作为网格节点的旅行时这样就可以确保在向前处理过程中考虑来自各个方向的射线,即所计算出来的节点旅行时是最小的全方位循环法向前处理实现步骤如下1)先用LTI算法按列扫描计算出整个区域各网格节点的旅行时2)除震源所在网格外,重新计算震源所在行各节点的旅行时(按行扫描),如图5a所示,并与原先按列扫描计算的结果对比,取其中最小值作为节点的旅行时3)由下到上逐个计算震源所在行上面所有节点的旅行时
(1)计算各单元水平边界及垂直边界上各节点的旅行时,这时只考虑来自下边界的射线(图5b),并与原先按列计算的比较,取其中旅行时最小的;
(2)从左往右,利用单元左边界上的已知点计算右边界和上边界节点的旅行时(图5c),并与节点已有值比较,取其中旅行时最小的;
(3)从右往左,利用单元右边界上的已知点计算左边界和上边界节点的旅行时(图5d),并与节点已有值比较,取其中旅行时最小的通过
(1)、
(2)、
(3)步计算,就可求得震源所在行上面所有节点的最小旅行时4)同理,从上到下可以计算出震源所在行下面各节点的旅行时,并与按列计算的旅行时比较取其中最小的最后可得到整个模型区域内所有节点的最小旅行时图6为LTI算法全方位循环对图4a追踪的结果与图3b相比,追踪的射线其路径和理论值是一样的,尤其是标示为
1、2的两条射线,其旅行时与理论值一致(表1)由试算结果可知,LTI算法全方位循环法在向前处理过程中可以确保计算出来的节点旅行时最小,从而能够避免仅按列或按行扫描未考虑逆向传播的射线导致在横向变化剧烈的模型中所追踪出来的初至射线路径不满足最小旅行时原理的问题
2.2相同节点间距的算法LTI法是基于线性假设的,模型网格剖分的精细程度直接影响射线旅行时的计算精度:在一定范围内,网格剖分越精细,射线旅行时精度越高;但网格剖分得越精细,在计算时,所耗费的时间也就越长因此在保证计算精度同时,如何提高计算效率也是目前需解决的问题为了保证LTI法的计算精度,同时提高该方法的计算效率,笔者提出在网格单元边界加入次生节点的技巧图7a是改进前LTI算法网格节点示意图,图7b为改进后在网格边界插入次生节点的示意图可以看出,在剖分精度相同的情况下,传统LTI计算点数多改进后只在主节点上给定速度值,边界上插入的次生节点速度可由单元边界上相邻两主节点的速度线性插值得到;而原来的做法(图7a)则需在每个节点都赋速度值因此在将该方法用于层析成像时,插入次生节点的LTI法在不影响计算精度的前提下所需反演的参数更少,更有利于反演为了验证加入次生节点后LTI法的计算精度,建立一均匀速度模型模型参数及观测系统如图8a所示图8b显示的是网格剖分个数为10X10,单元边界次生节点数为99X99,即节点间距为1m情况下的初至波射线路径;表2列出了传统的LTI法与网格边界插入次生节点的LTI法在节点间距为5m和1m两种情况下各道所接收到的初至波旅行时;图9给出了不同节点间距下插入次生节点前后LTI法的计算耗时由表2可以看出在节点间距为5m的情况下,传统的LTI法绝对误差最大值为
0.391ms,在单元边界插入次生节点后其绝对误差最大值降为
0.006ms;在节点间距为1m的情况下,传统的LTI法绝对误差最大值为
0.084ms,单元边界插入次生节点后绝对误差最大值为0mso从图9中可以看到在节点间距相同的情况下,传统的LTI法比插入次生节点的LTI法耗时少,而且随着节点间距逐渐减小,耗时的下降更为显著从试算结果可以看到:在节点间距相同的情况下,网格边界插入次生节点的LTI法较传统的LTI法计算精度至少可以提高一个数量级同时,计算速度也更快,随着节点间距剖分的越精细,计算耗时下降也越明显,计算速度较传统的方法可提高n倍3预测结果及旅行时残差为了综合检验改进后LTI法对复杂模型的适应能力,笔者引用了经典的Marmousi模型对改进前后的LTI法进行试算及误差分析图10为Marmousi模型的试算结果模型参数为网格数384义122,网格大小lOmXlOm,速度位于网格节点观测系统设置参数:激发点位于模型底部2800,1220,在模型顶部均匀放置77个接收点,接收点分别位于0,0,50,0,100,0,3800,0o图10是改进后的LTI算法追踪出来的初至射线路径,次生节点数为4X4由射线路径可以看出,在模型底部,射线沿着高速层走,符合最小旅行时射线路径的规律图11是传统的LTI法计算的各道初至旅行时减去改进后的LTI法计算的各道初至旅行时得到的旅行时残差曲线可以看到,所计算出来的旅行时残差均处于残差为的参考线上方,也就是说改进后的LTI算法算出来的初至旅行时小于传统的LTI法所计算出来的初至旅行时可以验证改进后的LTI算法较传统的LTI法计算精度更高,且改进后的LTI法对复杂模型具有更强的适应能力及稳定性4网格边界采用cti算法1笔者用实例详细分析了LTI算法仅按列或按行扫描没有考虑到逆向传播的射线,因此建议在向前处理时采用全方位循环的方法计算网格节点走时,这样可以确保追踪出来的射线路径满足最短旅行时2在保证计算精度的前提下,为能有效地提高计算效率,笔者采用在网格边界加入次生节点的措施试算结果表明,较传统的LTI法,本文所采用的LTI算法其适应能力更强、计算精度及计算效率均大幅度提高3本文采用2D模型,由于其高精度、高效率的特点,稍加扩展后即可较好地适用三维射线追踪及三维层析成像。