CN108303736A - 各向异性ti介质最短路径射线追踪正演方法 - Google Patents

Abstract

本发明提供一种各向异性TI介质最短路径射线追踪正演方法，其包括如下步骤：S1、建立各向异性介质的群速度模型；S2、通过各向异性介质最短路径射线追踪算法进行各向异性介质最短路径射线追踪。

Description

各向异性TI介质最短路径射线追踪正演方法

技术领域

本发明涉及各向异性介质射线追踪正演技术领域，特别涉及一种各向异性TI介质最短路径射线追踪正演方法。

背景技术

现有技术对各向异性介质射线追踪正演进行研究时，大多数采用各向异性介质的简化形式，即对称轴倾角垂直的介质，也称为VTI介质。选择VTI介质进行研究有两大理由：首先，自然界很多岩石矿物的结构特点类似于VTI介质的组成结构，例如沉积岩中的泥岩、成层性较好的岩石和具有水平裂隙的岩石。其次，VTI介质包含五个弹性模量参数使各向异性介质在一定的程度得到了简化，各向异性介质射线追踪算法里可以避免复杂的几何关系式。Qin与Schuster(1993)根据惠更斯原理提出了二维VTI介质的初至走时计算方法；Ruger与Alkhalifah(1996)针对二维VTI介质提出了一种有效的射线追踪方法；Qian与Symes(2002)针对二维、三维VTI介质提出了利用有限差分法求解体波走时的计算方法。自然界的岩石沉积之后往往会伴随一系列地壳运动，这会导致地层的对称轴倾角存在不同情况，而不仅仅只有垂直对称轴倾角(VTI介质情形)。这里简介的射线追踪方法就是针对具有任意对称轴倾向的地层，这种方法是各向同性介质最短路径射线追踪方法向各向异性介质情形的一种拓展应用。

实际上，将地层假设为弱各向异性介质进行射线追踪具有不合理性，因为自然界的岩石和地层具有任意各向异性强度。Thomsen(1986)指出Thomsen参数既可以描述弱各向异性介质，也可以描述强各向异性介质，差别在于定义的五个参数不一样，表征任意强度各向异性介质的参数为(α₀,β₀,ε,δ^*,γ),表征弱各向异性介质的参数为(α₀,β₀,ε,δ,γ)。虽然这两组参数大部分相同，但是利用它们表达的群速度和相速度表达式差别比较大。弱各向异性介质地震射线追踪的参数模型为v_p0,ε,δ。然而，任意强度各向异性介质地震地震射线追踪的参数模型为α₀,β₀,ε,δ^*,γ，或者对应的弹性模量参数。

发明内容

有鉴于此，本发明提出一种各向异性TI介质最短路径射线追踪正演方法。

一种各向异性TI介质最短路径射线追踪正演方法，其包括如下步骤：

S1、建立各向异性介质的群速度模型；

S2、通过各向异性介质最短路径射线追踪算法进行各向异性介质最短路径射线追踪。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

二维情形的群速度表达式{U_m＝(x,θ₀,n),m＝1,2,3}与三维情形的群速度表达式其中m＝1,2,3分别对应三种波的模式qP,qSV和qSH，三种波模式的群速度都是空间坐标x、岩石对称轴倾角θ₀与方位角以及相慢度方向n的函数。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S2包括：

S21、模型参数步骤；

S22、对走时场进行初始化；

S23、计算一个晶格内两个节点之间的最小走时；

S24、求取所有炮点的走时和射线；

S25、求取反射波与衍射波的射线路径。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S21包括：

以一定的网格单元大小，将二维各向异性介质模型m(x)划分成由主节点组成的网格模型

{m(x_k),k＝1,2,...,N₁,N₁＝(N_x+1)(N_z+1)}，其中N_x、N_z分别代表x和z方向的晶格数目；二维各向异性介质模型由多个矩形晶格组成，每个矩形晶格的四个角被定义为主节点；采用较大的晶格单元，同时在矩形晶格的每条边增加次一级节点来提高计算的精度与效率；将网格化的各向异性介质模型转化成为群速度模型U_m(x_k,θ₀,n),m＝1,2,3分别对应qP、qSV和qSH三种波模式，而k∈N₁；群速度写成射线方向r⁰的函数U_m(x_k,θ₀,r⁰)，选择群速度的最小值从而确保群速度曲线的一致连续性；群速度模型晶格内部的群速度根据晶格四个顶点所在节点的群速度值，利用双线性插值函数或者拉格朗日插值函数进行求取

其中矩形晶格四个节点的位置与网格模型的节点位置{x_k,k∈N₁}相重合，这些节点被称为初始节点或者主节点；

在模型参数化过程中增加次节点密度；如果次一级的节点数为N₂，那么整个模型的节点数为N₁+N₂；此外，如果炮点、检波点位置和模型网格节点的位置不重合，在炮点和检波点位置增加额外的节点，因此整个节点网络包含三类节点：主节点、次节点、炮点和炮检点位置的节点，整个节点网络用节点序号l和空间坐标{x_l＝(x_l,z_l),l＝1,2,...,N}来加以识别。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S22包括：

对炮点所在节点的走时赋初值"0"，其它网格节点的走时值赋为无穷大，从炮点所在节点开始，计算炮点周围节点的走时，再从相邻节点计算其它节点的走时，并记录相邻节点的走时最小值。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S23包括：

将一个晶格内两个节点之间的走时通过如下公式近似求取，即利用节点i和节点j之间的距离除以两个节点之间的平均群速度

其中N_j是未计算节点j周围的已知走时节点数目，群速度和均由步骤S21中公式求得。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S24包括：

当所有检波点的最小走时计算完成之后，拾取检点所在节点的走时，并根据入射节点的序号以反向追踪的形式，求取连接检波点与炮点之间的射线路径，得到共炮点道集的走时和射线路径；继续求取下一个炮点对应检波点的走时和射线路径，反复利用这一技术路线，以完成所有炮点的走时和射线求取。

在本发明所述的各向异性TI介质最短路径射线追踪正演方法中，

所述步骤S25包括：

同构一系列密集的节点来代表界面，把界面节点与所述的主节点、次节点、炮点和炮检点位置的节点组合起来，得到一个新的网络，包括主节点、次节点、界面节点和炮、检点节点；每个界面对应的反射波或者衍射波走时τ的计算方法如下：

τ＝min{τ_S(x_int)+τ_R(x_int),x_int∈Ω_int}，其中，τ_S(x_int)和τ_R(x_int)分别是炮点S与检波点R到界面节点位置x_int的初至走时值，Ω_int代表界面节点集合；分别假设炮点S与检波点R为震源点位置，进行两次正演走时计算；τ_S(x_int)在第一次最短路径射线追踪计算时已经有值了，剩下步骤将检波点R作为震源点位置再做一次最短路径射线追踪计算求取走时τ_R(x_int)，然后根据步骤S25中公式寻找相应的界面反射节点，以及反射波的走时与射线路径。

实施本发明提供的各向异性TI介质最短路径射线追踪正演方法与现有技术相比具有以下有益效果：(1)针对的研究对象各向异性介质可以是任意强度各向异性介质情形，即包含强各向异性性质的介质和弱各向异性性质的介质。(2)各向异性介质的倾角可以包含VTI、HTI和TTI三种情形。(3)地震射线追踪算法是基于非均匀节点网格速度模型，这种网格模型的方案是采用较大的晶格单元，同时在矩形晶格的每条边增加次一级节点，这样可以提高射线追踪运算的精度与效率。

附图说明

图1是本发明实施例的各向异性TI介质最短路径射线追踪正演方法流程图；

图2是qP-qP反射波在倾斜界面模型中的射线路径(a)与走时(b)分布图；

图3是qP-qSV转换波在倾斜界面模型中的射线路径(a)与走时(b)分布图；

图4是qP-qSH转换波在倾斜界面模型中的射线路径(a)与走时(b)分布图；

图5是qP-qP反射波在弯曲界面模型中的射线路径(a)与走时(b)分布图；

图6是qP-qSV转换波在弯曲界面模型中的射线路径(a)与走时(b)分布图；

图7是qP-qSH转换波在弯曲界面模型中的射线路径(a)与走时(b)分布图。

具体实施方式

如图1所示，本发明提供一种各向异性TI介质最短路径射线追踪正演方法，其包括如下步骤：

S1、建立各向异性介质的群速度模型；

S2、通过各向异性介质最短路径射线追踪算法进行各向异性介质最短路径射线追踪。

如图2-7所示，以下对本发明实施例的原理进行进一步说明:

1.各向异性介质的群速度模型

如果将最短路径射线追踪方法从各向同性介质拓展至各向异性介质，首先需要组成各向异性介质的五个弹性模量参数{a₁₁,a₁₃,a₃₃,a₄₄,a₆₆}(量纲为N/m²,牛顿每平方米)与对称轴倾角θ₀(x)，并且进一步研究它们和三个群速度公式{U_m(x,θ₀,r⁰),m＝1,2,3}之间的关系，三个群速度分别对应qP波、qSV波和qSH波的群速度，可以看出群速度表达式不仅是空间位置x的函数，而且还是射线方向r⁰与介质对称轴倾角θ₀(x)的函数。根据群速度模型，可以利用最短路径射线追踪方法计算各向异性介质qP波、qSV波和qSH波的初至走时和射线路径，通过进一步的修改还可以计算反射波、衍射波和转换波的走时和射线路径。VTI介质的相速度表达式，即克里斯托弗矩阵的特征值如下：

式中P和Q的表达式为

式中Q₁、Q₂和Q₃的表达式为

上述关系式中，角度为相慢度向量和各向异性介质对称轴倾角之间的夹角，{a₁₁,a₁₃,a₃₃,a₄₄,a₆₆}是各向异性介质的五个弹性模量参数，方程(4.1)的下标1、2和3分别代表qP波、qSV波和qSH波。计算群速度向量大小的具体表达式如下：

式中m＝1,2,3分别对应qP波、qSV波和qSH波。群速度向量的水平和垂直分量和可以写成如下形式

为了避免m与水平、垂直分量标记h,z相混淆，式中波的模式m写成上标形式。

公式(4.4)和(4.5)的偏导数项可以利用公式(4.1)-(4.3)对角度进行求偏导数得到，它们的具体形式如下

式中和的表达式为

由公式(4.1)-(4.7)可以看出VTI介质的群速度与五个弹性模量参数{a₁₁,a₁₃,a₃₃,a₄₄,a₆₆}以及相慢度向量有关。

为了使公式(4.1)-(4.7)能够应用于一般TI介质，可以利用标准三角关系(坐标旋转)得到入射角的余弦表达式。入射角是指相慢度方向与TI介质对称轴的夹角，而相慢度方向通常是由计算坐标系给出，不是由对称轴坐标系给出。如果知道相慢度方向和各向异性介质的对称轴倾角那么入射角可以表示为然后，才能应用公式(4.1)计算三个波的相速度，所以要由和来求出

将岩石的对称轴从垂直状态旋转至可以得到一般TI介质群速度向量的三个分量形式，它的具体计算方法为：已知岩石对称轴的方向向量为由它可以得到与e_z垂直的x轴方向向量e_x的表达式

然后，由e_y＝(e_z×e_x)/|(e_z×e_x)|得到程序里就是根据这个原理来得到岩石坐标系的。因此我们可以得到群速度的矢量形式为

代入e_x,e_z，可以得到U^(m)的表达式：

这样，就得到了在计算机坐标系的表达式

二维各向异性介质模型观测方位角与介质的对称轴方位角均为零，而且介质的对称轴倾角θ₀、相慢度方向倾角θ与射线入射角之间的关系可以进一步简化为因此二维各向异性介质群速度的表达式可以写成

可以采用如下表达式，对二维与三维各向异性介质模型进行描述

m(x)＝{a₁₁(x),a₁₃(x),a₃₃(x),a₄₄(x),a₆₆(x),θ₀(x)}，2D case (4.14)

模型向量的所有分量都随着空间坐标(x)的变化而变化，将公式(4.14)和(4.15)代入公式(4.1)和(4.8)可以得到二维情形的群速度表达式{U_m＝(x,θ₀,n),m＝1,2,3}与三维情形的群速度表达式(m＝1,2,3分别对应三种波的模式qP,qSV和qSH)，三种波模式的群速度都是空间坐标x、岩石对称轴倾角θ₀与方位角以及相慢度方向n的函数。将公式(4.14)和(4.15)定义的弹性参数模型代入二维与三维群速度表达式就可以得到弹性参数与群速度之间的关系式，这种关系式是后续进行射线追踪算法的基础。群速度可以用相慢度方向n或者射线方向r⁰表示，利用公式(4.13)可以计算出二维情形的射线角

因此射线方向可以写成r⁰＝(sinθ′,cosθ′)，群速度{U_m＝(x,θ₀,r⁰),m＝1,2,3}随射线方向r⁰的变化而变化。

2各向异性介质最短路径射线追踪算法

2.1各向异性介质最短路径射线追踪的原理

如果已知射线路径，那么可以沿着射线路径利用线性积分计算走时，计算表达式如下

式中ds＝|dx|是射线段的长度，x是局部射线向量；由公式(4.17)可知走时计算需要确定弹性参数模型m(x)的群速度U(x,θ₀,r⁰)与射线路径R(x)。上一节介绍了怎样建立弹性参数模型m(x)与群速度U(x,θ₀,r⁰)之间的关系，这一节将描述如何确定射线路径段R(x)，并进一步利用群速度U(x,θ₀,r⁰)计算走时τ。

根据费马原理，射线路径R(x)是使积分公式(4.17)不变(δτ＝0)的一条轨迹，这意味着在空间位置x_B处的走时计算表达式如下

式中Ω_B是模型内x_B点的邻域，r⁰是射线方向从x_A点到x_B点的单位向量，公式(4.18)为各向同性介质的射线追踪拓展至各向异性介质提供了可能性。可以明显地看出这种扩展过程需要考虑到基于射线方向r⁰的群速度与搜索从x_A点到x_B点最小走时的有效方法。这使得我们可以将Nakanishi和Yamaguchi提出的最短路径方法从各向同性介质扩展至各向异性介质。

2.2各向异性介质最短路径射线追踪的步骤

最短路径射线追踪方法遵循惠更斯原理，根据这种原理各向异性介质的最短路径方法需要包括以下几个步骤。

(1)模型的参数化：以一定的网格单元大小，将二维各向异性介质模型m(x)划分成由主节点组成的网格模型

{m(x_k),k＝1,2,...,N₁,N₁＝(N_x+1)(N_z+1)} (4.19)

式中N_x、N_z分别代表x和z方向的晶格数目。二维各向异性介质模型由许多矩形晶格组成，每个矩形晶格的四个角被定义为主节点，根据这种规则网格方法，在实际走时计算时须将晶格划分的非常小，但是我们可以采用非规则网格方法：采用较大的晶格单元，同时在矩形晶格的每条边增加次一级节点来提高计算的精度与效率(Bai^[82]，2007)。将网格化的各向异性介质模型转化成为群速度模型U_m(x_k,θ₀,n),m＝1,2,3分别对应qP、qSV和qSH三种波模式，而k∈N₁；群速度可以写成射线方向r⁰的函数U_m(x_k,θ₀,r⁰)，但是qSV波的群速度存在多值现象，在处理qSV波的多值问题时，我们选择群速度的最小值从而确保群速度曲线的一致连续性。群速度模型晶格内部的群速度可以根据晶格四个顶点所在节点的群速度值，利用双线性插值函数或者拉格朗日插值函数进行求取

式中矩形晶格四个节点的位置与网格模型的节点位置{x_k,k∈N₁}相重合，这些节点被称为初始节点或者主节点。为了提高计算的效率和精度，模型参数化过程中需要增加次节点密度。如果次一级的节点数为N₂，那么整个模型的节点数为N₁+N₂。另外，如果炮点、检波点位置和模型网格节点的位置不重合，那么必须在炮点和检波点位置增加额外的节点，因此整个节点网络包含三类节点：主节点、次节点、炮点和炮检点位置的节点，整个节点网络需要用节点序号l和空间坐标{x_l＝(x_l,z_l),l＝1,2,...,N}来加以识别。

(2)初始化走时场：对炮点所在节点的走时赋初值"0"，其它网格节点的走时值赋为无穷大，例如赋初值"10⁷"，从炮点所在节点开始，计算炮点周围节点的走时，再从相邻节点计算其它节点的走时，并记录相邻节点的走时最小值。

(3)最小走时计算方法：一个晶格内两个节点之间的走时可以利用如下公式近似求取，即利用节点i和节点j之间的距离除以两个节点之间的平均群速度

式中N_j是未计算节点j周围的已知走时节点数目，群速度和均由公式(4.20)进行计算求得；由公式(4.20)可以看出如果模型是各向同性介质，则群速度与射线方向无关，公式(4.21)可以简化为标准的最短射线路径走时计算方程。

(4)射线路径的求取：当所有检波点的最小走时计算完成之后，可以拾取检点所在节点的走时，并根据入射节点的序号以反向追踪的形式，求取连接检波点与炮点之间的射线路径，这样就可以得到共炮点道集的走时和射线路径；根据这种走时和路径的求取方法，可以继续求取下一个炮点对应检波点的走时和射线路径，反复利用这一技术路线，便可以完成所有炮点的走时和射线求取。

(5)反射波与衍射波的射线路径求取：如果复杂的各向异性介质模型存在界面，那么地震波在界面就会产生反射波与衍射波。我们可以用一系列密集的节点来代表界面，把界面节点与前面所述的三类节点组合起来，可以得到一个新的网络，它包括主节点、次节点、界面节点和炮、检点节点。网络中大部分晶格的节点分布是一样的，除了含有界面节点的那些晶格，但是它们只占整个模型晶格数目的一小部分。每个界面对应的反射波或者衍射波走时τ的计算方法如下

τ＝min{τ_S(x_int)+τ_R(x_int),x_int∈Ω_int} (4.22)

式中，τ_S(x_int)和τ_R(x_int)分别是炮点S与检波点R到界面节点位置x_int的初至走时值，Ω_int代表界面节点集合。分别假设炮点S与检波点R为震源点位置，进行两次正演走时计算。τ_S(x_int)在第一次最短路径射线追踪计算时就已经有值了，剩下步骤是将检波点R作为震源点位置再做一次最短路径射线追踪计算求取走时τ_R(x_int)，然后根据公式(4.22)寻找相应的界面反射节点，以及反射波的走时与射线路径。值得注意的问题是：在计算反射波与衍射波走时的时候，群速度模型只与界面以上的部分有关，也就是说在计算τ_S(x_int)和τ_R(x_int)时，只截取了模型界面以上部分的节点，这样可以大大提高计算效率。如果一个晶格包含界面节点，那么需要将界面以下节点的群速度值用界面以上节点的群速度值来替代，再用群速度计算公式(4.20)来求取τ_S(x_int)和τ_R(x_int)。同理，这种计算方法也可以运用于转换波走时与射线路径的求取。

3各向异性介质最短路径射线追踪数值模拟

上述内容介绍了各向异性介质最短路径射线追踪算法的基本原理与实现步骤，这里建立倾斜界面模型和弯曲界面模型，分别研究qP-qP、qP-qSV与qP-qSH波在两种模型中的反射波射线路径和走时分布。

3.1倾斜界面模型的射线追踪数值模拟

倾斜界面模型水平方向的长度为800米，垂直方向的深度为500米，第一层介质的弹性模量参数大小为a₁₁＝9.08,a₁₃＝2.98,a₃₃＝7.53,a₄₄＝2.27,a₆₆＝3.84，相应的对称轴倾角和方位角分别为θ₀＝90°,第二层介质的弹性模量参数大小为a₁₁＝20.3,a₁₃＝9.58,a₃₃＝22.3,a₄₄＝8.35,a₆₆＝11.35，相应的对称轴倾角和方位角分别为θ₀＝90°,第三层介质的弹性模量参数大小为a₁₁＝13.86,a₁₃＝4.31,a₃₃＝10.93,a₄₄＝3.31,a₆₆＝4.34，相应的对称轴倾角和方位角分别为θ₀＝0°,炮点位于地表400米处，共有33个检波点，检波点间距为25米。水平与垂直方向的网格间距为5米。下面依次研究qP-qP、qP-qSV与qP-qSH波在倾斜界面模型的反射波射线路径和走时分布。

(1)qP-qP反射波在倾斜界面模型中的射线路径与走时分布，如图2所示。

(2)qP-qSV转换波在倾斜界面模型中的射线路径与走时分布，如图3所示。

(3)qP-qSH转换波在倾斜界面模型中的射线路径与走时分布，如图4所示。

3.2弯曲界面模型的射线追踪数值模拟

弯曲界面模型水平方向的长度为800米，垂直方向的深度为500米，第一层介质的弹性模量参数大小为a₁₁＝9.08,a₁₃＝2.98,a₃₃＝7.53,a₄₄＝2.27,a₆₆＝3.84，相应的对称轴倾角和方位角分别为θ₀＝90°,第二层介质的弹性模量参数大小为a₁₁＝20.3,a₁₃＝9.58,a₃₃＝22.3,a₄₄＝8.35,a₆₆＝11.35，相应的对称轴倾角和方位角分别为θ₀＝90°,第三层介质的弹性模量参数大小为a₁₁＝13.86,a₁₃＝4.31,a₃₃＝10.93,a₄₄＝3.31,a₆₆＝4.34，相应的对称轴倾角和方位角分别为θ₀＝0°,炮点位于地表400米处，共有33个检波点，检波点间距为25米。水平与垂直方向的网格间距为5米。下面依次研究qP-qP、qP-qSV与qP-qSH波在弯曲界面模型的反射波射线路径和走时分布。

(1)qP-qP反射波在弯曲界面模型中的射线路径与走时分布，如图5所示。

(2)qP-qSV转换波在弯曲界面模型中的射线路径与走时分布，如图6所示。

(3)qP-qSH转换波在弯曲界面模型中的射线路径与走时分布，如图7所示。

可以理解的是，对于本领域的普通技术人员来说，可以根据本发明的技术构思做出其它各种相应的改变与变形，而所有这些改变与变形都应属于本发明权利要求的保护范围。

Claims (8)

1.一种各向异性TI介质最短路径射线追踪正演方法，其特征在于，其包括如下步骤：

S1、建立各向异性介质的群速度模型；

S2、通过各向异性介质最短路径射线追踪算法进行各向异性介质最短路径射线追踪。

2.如权利要求1所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

二维情形的群速度表达式{U_m＝(x,θ₀,n),m＝1,2,3}与三维情形的群速度表达式其中m＝1,2,3分别对应三种波的模式qP,qSV和qSH，三种波模式的群速度都是空间坐标x、岩石对称轴倾角θ₀与方位角以及相慢度方向n的函数。

3.如权利要求2所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S2包括：

S21、模型参数步骤；

S22、对走时场进行初始化；

S23、计算一个晶格内两个节点之间的最小走时；

S24、求取所有炮点的走时和射线；

S25、求取反射波与衍射波的射线路径。

4.如权利要求3所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S21包括：

以一定的网格单元大小，将二维各向异性介质模型m(x)划分成由主节点组成的网格模型

{m(x_k),k＝1,2,...,N₁,N₁＝(N_x+1)(N_z+1)}，其中N_x、N_z分别代表x和z方向的晶格数目；二维各向异性介质模型由多个矩形晶格组成，每个矩形晶格的四个角被定义为主节点；采用较大的晶格单元，同时在矩形晶格的每条边增加次一级节点来提高计算的精度与效率；将网格化的各向异性介质模型转化成为群速度模型U_m(x_k,θ₀,n),m＝1,2,3分别对应qP、qSV和qSH三种波模式，而k∈N₁；群速度写成射线方向r⁰的函数U_m(x_k,θ₀,r⁰)，选择群速度的最小值从而确保群速度曲线的一致连续性；群速度模型晶格内部的群速度根据晶格四个顶点所在节点的群速度值，利用双线性插值函数或者拉格朗日插值函数进行求取

其中矩形晶格四个节点的位置与网格模型的节点位置{x_k,k∈N₁}相重合，这些节点被称为初始节点或者主节点；

在模型参数化过程中增加次节点密度；如果次一级的节点数为N₂，那么整个模型的节点数为N₁+N₂；此外，如果炮点、检波点位置和模型网格节点的位置不重合，在炮点和检波点位置增加额外的节点，因此整个节点网络包含三类节点：主节点、次节点、炮点和炮检点位置的节点，整个节点网络用节点序号l和空间坐标{x_l＝(x_l,z_l),l＝1,2,...,N}来加以识别。

5.如权利要求3所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S22包括：

对炮点所在节点的走时赋初值"0"，其它网格节点的走时值赋为无穷大，从炮点所在节点开始，计算炮点周围节点的走时，再从相邻节点计算其它节点的走时，并记录相邻节点的走时最小值。

6.如权利要求5所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S23包括：

将一个晶格内两个节点之间的走时通过如下公式近似求取，即利用节点i和节点j之间的距离除以两个节点之间的平均群速度

其中N_j是未计算节点j周围的已知走时节点数目，群速度和均由步骤S21中公式求得。

7.如权利要求6所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S24包括：

当所有检波点的最小走时计算完成之后，拾取检点所在节点的走时，并根据入射节点的序号以反向追踪的形式，求取连接检波点与炮点之间的射线路径，得到共炮点道集的走时和射线路径；继续求取下一个炮点对应检波点的走时和射线路径，反复利用这一技术路线，以完成所有炮点的走时和射线求取。

8.如权利要求6所述的各向异性TI介质最短路径射线追踪正演方法，其特征在于，

所述步骤S25包括：

同构一系列密集的节点来代表界面，把界面节点与所述的主节点、次节点、炮点和炮检点位置的节点组合起来，得到一个新的网络，包括主节点、次节点、界面节点和炮、检点节点；每个界面对应的反射波或者衍射波走时τ的计算方法如下：

τ＝min{τ_S(x_int)+τ_R(x_int),x_int∈Ω_int}，其中，τ_S(x_int)和τ_R(x_int)分别是炮点S与检波点R到界面节点位置x_int的初至走时值，Ω_int代表界面节点集合；分别假设炮点S与检波点R为震源点位置，进行两次正演走时计算；τ_S(x_int)在第一次最短路径射线追踪计算时已经有值了，剩下步骤将检波点R作为震源点位置再做一次最短路径射线追踪计算求取走时τ_R(x_int)，然后根据步骤S25中公式寻找相应的界面反射节点，以及反射波的走时与射线路径。