CN108072897A - 一种混合二维地震波走时计算方法 - Google Patents
一种混合二维地震波走时计算方法 Download PDFInfo
- Publication number
- CN108072897A CN108072897A CN201810077621.8A CN201810077621A CN108072897A CN 108072897 A CN108072897 A CN 108072897A CN 201810077621 A CN201810077621 A CN 201810077621A CN 108072897 A CN108072897 A CN 108072897A
- Authority
- CN
- China
- Prior art keywords
- mrow
- point
- msup
- walked
- mfrac
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000205 computational method Methods 0.000 title claims abstract description 18
- 238000010276 construction Methods 0.000 claims abstract description 13
- 238000000034 method Methods 0.000 claims description 12
- 235000013399 edible fruits Nutrition 0.000 claims description 3
- 238000005516 engineering process Methods 0.000 abstract description 4
- 230000005855 radiation Effects 0.000 abstract 1
- 230000000694 effects Effects 0.000 description 2
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
- 238000013316 zoning Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/30—Analysis
- G01V1/303—Analysis for determining velocity profiles or travel times
- G01V1/305—Travel times
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/282—Application of seismic models, synthetic seismograms
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/67—Wave propagation modeling
- G01V2210/671—Raytracing
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种混合二维地震波走时计算方法,包括以下步骤:读入相关参数、速度模型;沿着震源向不同方向发出射线并计算中心射线信息;使用波前构建法计算射线范围内网格点的地震波走时;对全空间网格点走时属性进行分类并以此为依据建立初始窄带;使用快速推进法计算剩余网格点地震波走时。本发明通过使用窄带技术连接了快速推进法与波前构建法,通过波前构建法提高了震源附近较小区域的地震波走时计算精度,从而提高了快速推进法在剩余区域网格节点的计算精度,实现了一种兼顾计算效率与计算精度的二维地震波走时计算。
Description
技术领域
本发明涉及地震波走时计算领域,特别是一种混合二维地震波走时计算方法。
背景技术
《工程地球物理学报》2009年第3期公开了张双杰等“快速推进法计算精度分析及改进”,介绍了两种提高快速推进法计算精度的方法,一是采用高阶差分格式进行计算;二是对震源点附近网格节点进行局部精细处理,在震源附近采用高阶差分格式进行计算,在剩余区域采用低阶差分格式进行计算。并且通过这两种方法对均匀介质模型进行地震波走时计算,实验结果得到了比较好的效果。
《世界地质》2016年第3期公开了王乾龙等“基于完全三叉树的快速推进法地震波走时计算”,介绍了一种改进的快速推进地震波走时计算方法,将完全三叉树排序方法引入到地震波走时计算中,减少了窄带扩展中寻找最小走时值的时间,提高了整个算法的计算效率。并且通过基于完全三叉树的快速推进地震波走时计算方法对层状模型、Marmousi模型、Sigsbee 2a模型进行了计算,实验结果得到了比较好的效果。
通过以上例子可以看出,现有快速推进地震波走时计算方法在一定程度上能够提升计算精度,但实现过程复杂,提升的计算精度也有限。
发明内容
本发明所要解决的技术问题是提供一种混合二维地震波走时计算方法,通过灵活运用窄带技术将波前构建地震波走时计算方法和快速推进地震波走时计算方法无缝连接起来,即通过在震源附近小范围内使用计算精度较高的波前构建法,在剩余区域使用快速推进法计算走时,在改善原有快速推进法计算精度的同时,仍保留了其高效的特点。
为解决上述技术问题,本发明采用的技术方案是:
一种混合二维地震波走时计算方法,包括以下步骤:
步骤1:读入相关参数文件、速度模型,其中,所述参数文件包含速度模型的网格点数、网格间距和震源位置;
步骤2:从炮点沿着不同方向发射射线并计算中心射线信息;
步骤3:使用波前构建法计算射线覆盖范围内网格节点走时;
步骤4:对所有网格点进行分类,将它们分为接受点、窄带点或远离点,即如果一个网格点的地震波走时已经计算过,并且它周围点的走时全部计算过,则这个点为接受点;如果一个网格点的地震波走时计算过,而它周围点的走时并非全部计算过,则这个点为窄带点;如果一个点的走时没有被计算过,则这个点为远离点;
步骤5:将所有窄带点移入窄带内,构建初始窄带;
步骤6:扩展窄带,直至窄带为空;在这一过程中,网格点走时是通过迎风差分法求解二维程函方程获得的,所述二维程函方程为:
|▽τ|=s
其中,τ为地震波走时,s为慢度,▽为梯度符号,上述公式中梯度项的迎风差分表达式为:
其中,分别为走时在网格节点(i,j)处x或z方向上的向前或向后差分表达式,具体形式分别为:
进一步的,在步骤2中,使用龙哥库塔法求解运动学射线追踪方程组得到中心射线信息,如下式所示:
其中,xi表示空间位置,pi表示慢度分量,τ表示地震波走时,v表示离散点处的速度值。
与现有技术相比,本发明的有益效果是:由于采用了波前构建法提高了震源附近区域网格节点地震波走时计算精度,从而提高了快速推进法计算区域网格点的计算精度;以降低较少的计算效率为代价大幅提高了整个地震波走时算法的计算精度。
附图说明
图1为本发明混合二维地震波走时计算方法流程图。
图2波前构建法区域分割示意图。
图3为均匀介质中快速推进法相对误差。
图4为均匀介质中混合地震波走时计算方法相对误差。
具体实施方式
下面结合附图和具体实施方式对本发明作进一步详细的说明。
图1是混合二维地震波走时计算方法流程图,图中显示了本发明方法的实现流程,具体如下:
1)读入相关参数文件、速度模型,其中,所述参数文件包含速度模型的网格点数、网格间距和震源位置。
2)从炮点沿着不同方向发射射线,射线的发射角度范围为:-90°到+90°,射线之间的角度间隔为3°到10°,使用龙哥库塔法求解运动学射线追踪方程组得到中心射线信息,如下式所示:
其中,xi表示空间位置,pi表示慢度分量,τ表示地震波走时,v表示离散点处的速度值。
3)使用波前构建法计算射线覆盖范围内网格节点走时。在计算过程中这一区域被相邻中心射线上的相邻波前面点划分为多个不规则四边形(如图2所示),每个四边形内的网格节点是通过四边形顶点信息插值获取的。
4)完成波前构建法计算后,对全空间网格节点进行属性分类,规定:如果一个网格点的地震波走时已经计算过,并且它周围点的走时全部计算过,则这个点为接受点;如果一个网格点的地震波走时计算过,而它周围点的走时并非全部计算过,则这个点为窄带点;如果一个点的走时没有被计算过,则这个点为远离点。
5)将所有窄带点移入窄带内,构建初始窄带。
6)扩展窄带,直至窄带为空。在这一过程中,网格点走时是通过迎风差分法求解二维程函方程获得的,所述二维程函方程为:
|▽τ|=s
其中,τ为地震波走时,s为慢度,▽为梯度符号,上述公式中梯度项的迎风差分表达式为:
其中,分别为走时在网格节点(i,j)处x或z方向上的向前或向后差分表达式,具体形式分别为:
下面通过均匀介质模型对本发明方法的计算精度以及稳定性进行分析验证。
图3、图4分别为快速推进法与混合地震波走时计算方法在均匀介质模型中的相对误差,均匀模型的大小为601×601,网格间距为10m×10m,速度为1000m/s。联合走时方法中追出单根射线的最大长度为500m。从图中可以看出混合地震波走时计算方法相对于快速推进法计算精度有了很大的提高。
本发明通过使用窄带技术连接了快速推进法与波前构建法,通过波前构建法提高了震源附近较小区域的地震波走时计算精度,从而提高了快速推进法在剩余区域网格节点的计算精度,实现了一种兼顾计算效率与计算精度的二维地震波走时计算。
Claims (2)
1.一种混合二维地震波走时计算方法,其特征在于,包括以下步骤:
步骤1:读入相关参数文件、速度模型,其中,所述参数文件包含速度模型的网格点数、网格间距和震源位置;
步骤2:从炮点沿着不同方向发射射线并计算中心射线信息;
步骤3:使用波前构建法计算射线覆盖范围内网格节点走时;
步骤4:对所有网格点进行分类,将它们分为接受点、窄带点或远离点,即如果一个网格点的地震波走时已经计算过,并且它周围点的走时全部计算过,则这个点为接受点;如果一个网格点的地震波走时计算过,而它周围点的走时并非全部计算过,则这个点为窄带点;如果一个点的走时没有被计算过,则这个点为远离点;
步骤5:将所有窄带点移入窄带内,构建初始窄带;
步骤6:扩展窄带,直至窄带为空;在这一过程中,网格点走时是通过迎风差分法求解二维程函方程获得的,所述二维程函方程为:
|▽τ|=s
其中,τ为地震波走时,s为慢度,▽为梯度符号,上述公式中梯度项的迎风差分表达式为:
<mrow>
<mo>|</mo>
<mo>&dtri;</mo>
<mi>&tau;</mi>
<mo>|</mo>
<mo>=</mo>
<msup>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mrow>
<mi>max</mi>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>D</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mrow>
<mo>-</mo>
<mi>x</mi>
</mrow>
</msubsup>
<mi>&tau;</mi>
<mo>,</mo>
<mn>0</mn>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>+</mo>
<mi>min</mi>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>D</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mrow>
<mo>+</mo>
<mi>x</mi>
</mrow>
</msubsup>
<mi>&tau;</mi>
<mo>,</mo>
<mn>0</mn>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>+</mo>
<mi>max</mi>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>D</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mrow>
<mo>-</mo>
<mi>z</mi>
</mrow>
</msubsup>
<mi>&tau;</mi>
<mo>,</mo>
<mn>0</mn>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>+</mo>
<mi>min</mi>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>D</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mrow>
<mo>+</mo>
<mi>z</mi>
</mrow>
</msubsup>
<mi>&tau;</mi>
<mo>,</mo>
<mn>0</mn>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mrow>
<mn>1</mn>
<mo>/</mo>
<mn>2</mn>
</mrow>
</msup>
</mrow>
其中,分别为走时在网格节点(i,j)处x或z方向上的向前或向后差分表达式,具体形式分别为:
2.如权利要求1所述的一种混合二维地震波走时计算方法,其特征在于,在步骤2中,使用龙哥库塔法求解运动学射线追踪方程组得到中心射线信息,如下式所示:
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<msub>
<mi>dx</mi>
<mi>i</mi>
</msub>
</mrow>
<mrow>
<mi>d</mi>
<mi>&tau;</mi>
</mrow>
</mfrac>
<mo>=</mo>
<msup>
<mi>v</mi>
<mn>2</mn>
</msup>
<mo>&CenterDot;</mo>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<msub>
<mi>dp</mi>
<mi>i</mi>
</msub>
</mrow>
<mrow>
<mi>d</mi>
<mi>&tau;</mi>
</mrow>
</mfrac>
<mo>=</mo>
<mfrac>
<mo>&part;</mo>
<mrow>
<mo>&part;</mo>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
<mrow>
<mo>(</mo>
<mfrac>
<mn>1</mn>
<mi>v</mi>
</mfrac>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mo>-</mo>
<msup>
<mi>v</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mfrac>
<mrow>
<mo>&part;</mo>
<mi>v</mi>
</mrow>
<mrow>
<mo>&part;</mo>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
其中,xi表示空间位置,pi表示慢度分量,τ表示地震波走时,v表示离散点处的速度值。
Priority Applications (4)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810077621.8A CN108072897A (zh) | 2018-01-23 | 2018-01-23 | 一种混合二维地震波走时计算方法 |
NL2020684A NL2020684B1 (en) | 2018-01-23 | 2018-03-29 | Mixed 2D seismic wave travel time calculation method |
LU100751A LU100751B1 (fr) | 2018-01-23 | 2018-03-30 | Méthode de calcul de temps de propagation sismique mixte en 2D |
BE2018/0040A BE1025488B1 (fr) | 2018-01-23 | 2018-03-30 | Méthode de calcul de temps de propagation sismique mixte en 2D |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810077621.8A CN108072897A (zh) | 2018-01-23 | 2018-01-23 | 一种混合二维地震波走时计算方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN108072897A true CN108072897A (zh) | 2018-05-25 |
Family
ID=61978378
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810077621.8A Pending CN108072897A (zh) | 2018-01-23 | 2018-01-23 | 一种混合二维地震波走时计算方法 |
Country Status (4)
Country | Link |
---|---|
CN (1) | CN108072897A (zh) |
BE (1) | BE1025488B1 (zh) |
LU (1) | LU100751B1 (zh) |
NL (1) | NL2020684B1 (zh) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108957538A (zh) * | 2018-06-21 | 2018-12-07 | 成都启泰智联信息科技有限公司 | 一种虚拟震源二维波前构建地震波走时计算方法 |
CN115201901A (zh) * | 2022-06-30 | 2022-10-18 | 中铁第四勘察设计院集团有限公司 | 隧道波前走时的确定方法、装置、设备及可读存储介质 |
CN117194855A (zh) * | 2023-11-06 | 2023-12-08 | 南方科技大学 | 一种弱各向异性走时的拟合解析分析方法及相关设备 |
CN117607957A (zh) * | 2024-01-24 | 2024-02-27 | 南方科技大学 | 基于等效慢度的快速推进法的地震波走时求解方法及*** |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2002065159A1 (en) * | 2001-02-14 | 2002-08-22 | Hae-Ryong Lim | Method of seismic imaging using direct travel time computing |
CN105425286A (zh) * | 2015-10-30 | 2016-03-23 | 中国石油天然气集团公司 | 地震走时获取方法及基于其的井间地震走时层析成像方法 |
CN106353793A (zh) * | 2015-07-17 | 2017-01-25 | 中国石油化工股份有限公司 | 一种基于走时增量双线性插值射线追踪的井间地震层析反演方法 |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5724310A (en) * | 1996-10-15 | 1998-03-03 | Western Atlas International, Inc. | Traveltime generation in the presence of velocity discontinuities |
US6018499A (en) * | 1997-11-04 | 2000-01-25 | 3Dgeo Development, Inc. | Three-dimensional seismic imaging of complex velocity structures |
US8094513B2 (en) * | 2008-06-03 | 2012-01-10 | Westerngeco L.L.C. | Determining positioning of survey equipment using a model |
-
2018
- 2018-01-23 CN CN201810077621.8A patent/CN108072897A/zh active Pending
- 2018-03-29 NL NL2020684A patent/NL2020684B1/en not_active IP Right Cessation
- 2018-03-30 LU LU100751A patent/LU100751B1/fr active IP Right Grant
- 2018-03-30 BE BE2018/0040A patent/BE1025488B1/fr not_active IP Right Cessation
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2002065159A1 (en) * | 2001-02-14 | 2002-08-22 | Hae-Ryong Lim | Method of seismic imaging using direct travel time computing |
CN106353793A (zh) * | 2015-07-17 | 2017-01-25 | 中国石油化工股份有限公司 | 一种基于走时增量双线性插值射线追踪的井间地震层析反演方法 |
CN105425286A (zh) * | 2015-10-30 | 2016-03-23 | 中国石油天然气集团公司 | 地震走时获取方法及基于其的井间地震走时层析成像方法 |
Non-Patent Citations (3)
Title |
---|
丁鹏程等: "《基于三维多模板快速推进算法的复杂近地表射线追踪》", 《石油物探》 * |
孙章庆等: "《三维复杂山地多级次群推进迎风混合法多波型走时计算》", 《地球物理学报》 * |
李永博: "《快速行进法射线追踪提高旅行时计算精度和效率的改进措施》", 《石油地球物理勘探》 * |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108957538A (zh) * | 2018-06-21 | 2018-12-07 | 成都启泰智联信息科技有限公司 | 一种虚拟震源二维波前构建地震波走时计算方法 |
CN115201901A (zh) * | 2022-06-30 | 2022-10-18 | 中铁第四勘察设计院集团有限公司 | 隧道波前走时的确定方法、装置、设备及可读存储介质 |
CN117194855A (zh) * | 2023-11-06 | 2023-12-08 | 南方科技大学 | 一种弱各向异性走时的拟合解析分析方法及相关设备 |
CN117194855B (zh) * | 2023-11-06 | 2024-03-19 | 南方科技大学 | 一种弱各向异性走时的拟合解析分析方法及相关设备 |
CN117607957A (zh) * | 2024-01-24 | 2024-02-27 | 南方科技大学 | 基于等效慢度的快速推进法的地震波走时求解方法及*** |
CN117607957B (zh) * | 2024-01-24 | 2024-04-02 | 南方科技大学 | 基于等效慢度的快速推进法的地震波走时求解方法及*** |
Also Published As
Publication number | Publication date |
---|---|
BE1025488A1 (fr) | 2019-03-18 |
NL2020684A (en) | 2018-04-20 |
LU100751B1 (fr) | 2018-10-01 |
NL2020684B1 (en) | 2018-11-14 |
BE1025488B1 (fr) | 2019-03-25 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN108072897A (zh) | 一种混合二维地震波走时计算方法 | |
CN107392875A (zh) | 一种基于k近邻域划分的点云数据去噪方法 | |
Verma et al. | 3D building detection and modeling from aerial LIDAR data | |
CN104318622B (zh) | 一种室内场景非均匀三维点云数据的三角网格建模方法 | |
Ben‐Chen et al. | Conformal flattening by curvature prescription and metric scaling | |
CN106528740B (zh) | 基于Delaunay三角网的道路中心线提取方法 | |
CN104200212B (zh) | 一种基于机载LiDAR数据的建筑物外边界线提取方法 | |
CN105574929A (zh) | 一种基于地面LiDAR点云数据的单株植被三维建模方法 | |
CN108734728A (zh) | 一种基于高分辨序列图像的空间目标三维重构方法 | |
CN107886528A (zh) | 基于点云的配电线路作业场景三维重建方法 | |
Alharthy et al. | Detailed building reconstruction from airborne laser data using a moving surface method | |
CN104484868B (zh) | 一种结合模板匹配和图像轮廓的运动目标航拍跟踪方法 | |
CN104700451A (zh) | 基于迭代就近点算法的点云配准方法 | |
CN107843922A (zh) | 一种基于地震初至波和反射波走时联合的层析成像方法 | |
CN108711194A (zh) | 一种基于三次b样条插值的三维网格模型拼接方法 | |
CN105303612B (zh) | 一种基于不规则三角网模型的数字河网提取方法 | |
CN107767453A (zh) | 一种基于规则约束的建筑物lidar点云重构优化方法 | |
CN107424166A (zh) | 点云分割方法及装置 | |
CN105005995A (zh) | 一种计算三维点云模型骨骼的方法 | |
CN108957538A (zh) | 一种虚拟震源二维波前构建地震波走时计算方法 | |
CN109165476A (zh) | 一种模块化风场模型的建模方法及风场模拟方法 | |
CN114332291A (zh) | 一种倾斜摄影模型建筑物外轮廓规则提取方法 | |
CN113327276B (zh) | 一种面向移动测量大体量点云数据配准的方法 | |
CN106960470A (zh) | 三维点云曲面重建方法及装置 | |
CN107886569A (zh) | 一种基于离散李导数的测度可控的曲面参数化方法及*** |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
RJ01 | Rejection of invention patent application after publication |
Application publication date: 20180525 |
|
RJ01 | Rejection of invention patent application after publication |