专利名称：用于从单色波场计算有限频率地震迁移传播时间的方法
技术领域：
本发明涉及地球物理探测领域，更具体的涉及用于测定传播时间和用于迁移地震数据的传播时间图的方法。
在石油和天然气工业中，通常使用地球物理探测技术以帮助搜索和评估地下碳氢化合物沉积。通常，使用地震能量源产生传播进入地球的地震信号，并通过地下地震反射体(具有不同声阻抗的地下形成之间的分界面)至少部分对其反射。通过位于地球表面处或在地球表面附近、在水体中、或在钻孔中的已知深度处的地震探测器纪录反射，并且可以处理作为结果的地震数据以收益涉及地下反射物的位置和地下形成的物理特性的信息。
在用于自然资源的地球物理勘测的一种形式中，波场从在地球表面处或在地球的表面附近的源点径向辐射。初始时作为球形展开波阵面传播，辐射声穿透各种地层，由于在岩石密度和声速度中的变化，各种地层常常提供了在地层边界处的声阻抗失谐。由于通过变换器将反射的波场转换成电信号，从各自的层边界对波场进行反射以返回至存在机械地动的表面。将包括地震数据的信号存储于用于将来处理的档案存储媒介。
地震研究的目的是产生在地质和/或经济利益的区域中地下地球形成的量的图像。对于具有恒定速度的各向同性的水平地层，通过半平均速率的放大在源附近的接收器处波场发射和波场接收之间流逝的时间是地层上反射波场的入射点的深度，该地层直接位于源和接收器之间的中点之下。
如果反射器是倾斜的或速度是空间变化的，就不再保持那种简单关系；横向移动入射点以相对于源/接收器中点上升。合适的构图或倾斜的反射器需要从那些倾斜地层起源的波场的迁移。
复合地下环境排除了地震记录的简单堆栈。复合地下环境在多个方向上散射地震波。简单堆栈可以定位错误位置中的地下特征。地震记录的浅迁移能够提供地下特征的更精确的位置。迁移包括返回信号的地质复定位，以给出在它的合适位置中的事件(层边界或其它结构)。传播时间图是预叠加迁移处理的积分部分。一个熟知的迁移技术是基尔霍夫深度迁移。
对于预叠加深度成像，由于它相对于其它方法低的计算成本，常常选择基尔霍夫迁移。基尔霍夫迁移需要使用基于射线跟踪的几种熟知类型中的任何一种波阵面传播时间发生装置。射线跟踪方法在复合地质结构中是有用的。跟踪射线传播时间允许第一波至数据和最大能量波至，其中后面的数据提供了出众的图像。可是那些方法是非常慢和耗散计算处理时间的。
在选择的输出点处波场传播时间的计算是地震数据的成功的基尔霍夫深度迁移中的关键因素。对于复合地质概况，通过来自于源的射线跟踪产生至输出点的波场传播时间。射线跟踪产生方法产生了用于在输出点处所有波阵面波至的传播时间。
Shein-Shen Wang et al于1993年7月20日授权的美国专利No.5,229,938中给出了一种方法。其给出的是用于源和接收器对获得双向传播时间的方法，该方法包括测定用于每个源到多个成像点的一组单向传播时间和测定用于每个接收器到多个成像点的一组单向传播时间。产生射线装置用于源和接收器。通过利用射线装置的两点插值法计算从源位置到成像点的传播时间，通过求和两个装置、每个用于源和接收器位置的装置计算双向传播时间。获得双向传播时间装置以用于具体的源和用于所有成像点的接收器组合。
Eikonal(有限差分)传播时间发生装置是非常快速的并且不会产生影区。有限差分传播时间发生装置常常拾取第一波至传播时间，然而对于射线跟踪传播时间发生装置，必须具体选择波阵面的期望部分。
在由John Vidale，Bulletin of the Seismological Society ofAmerica，v.78，n.6，1998年12月，pp2062-2076出版的论文中给出了熟知的有线差分传播时间发生装置，这里通过有限差分逐点外推法在多向栅格上可以快速计算经过任何速度结构的第一波至地震波的传播时间。取代射线跟踪波阵面。合适的处理折射波并使影区充满合适的波阵面段。该方案是非常快速的并且在线断层反演中是有用的，并通过平滑横向速度梯度表征地质结构中的基尔霍夫迁移。
通过经过具体速率模型的射线跟踪可以直接计算从传播时间图或传播时间表得到的传播时间。通过熟知的方法可以提供速度模型。通常认为速度模型是地球的地震波场速度结构的表示。速度模型是描述了区域或利益的量内部的栅格上速度分布的2D或3D阵列。栅格可以是笛卡尔(x-y-z)栅格，当然，也可以使用球形、多面体或其它栅格。测定用于给出的地震数据量或地球结构的合适的速度模型在本领域中是熟知的，这里就不再详细讨论。
在表面处从源点形成的一束射线可以向下传播进入地球速度模型并经过地下对其跟踪，同时解释速度梯度的变化和具有速度对比的层边界处的折射引起的射线弯曲。识别沿着每个光程的反射点作为与层边界的射线的交叉点。然后通过积分沿着通过与要素联系的速度分开的光程的距离的要素，计算从表面处的源位置到地下中的反射点的传播时间。通过施加倒易性，以相同的方式可以计算从表面处的接收器位置到地下中的反射点的传播时间。最后，对于在表面处的给出源-接收器对和地下中的反射点，通过对从源到反射点的传播时间增加从反射点到接收器的传播时间计算总传播时间。


图1A是具有一组4个光程段用于简单的3层地球模型的速度模型举例。每个3层，V1、V2和V3具有速率，其中层V1的速率异于V2，并且V2异于V3。来自于地震能量在地球表面处的源位置S作为射线段(通过箭头表示)经过上部模型层V1传播。当它横越下一个地球层V2时，下一个射线段然后获得不同的路径。从源S到反射点M的总传播时间是第一到射线段的和。在第二层V2和第三层V3的分界面处，在反射点M处向着具有速度V2的表面向后反射射线并到达上部层V1，横越该上部层V1，并然后到达接收器位置R。对于这种源-接收器结构，通过沿着4独立的段路径增加射线段传播时间计算全部源到接收器的传播时间。
对具有与速度模型类似的节点结构的栅格可以构图计算的传播时间。这些传播时间图是表示了在利益区域内部的栅格上传播时间分布的2D或3D阵列。节点的这些栅格也与笛卡尔(x-y-z)栅格、球形、多面体或其它栅格中的速度模型相似。
图1B是2D计算栅格的举例，可以使用其对来自于速度模型的传播时间构图。在图1B中包括S和R的所有圆形表示栅格点GP，该栅格点表示用于可以被评估的传播时间的节点。该栅格在以箭头表示的x方向和通过箭头以z方向表示的深度上横向延伸。
从有限差分解法到eikonal方程或通过如先前给出的射线跟踪可以计算传播时间。可是两个方法都是计算用于渐进高频波的传播时间。由于波将以频率相关方式通过空间变化的速度区域，这些高频方法不对任何有限频率进行精确构图。
通过求和图像相素处记录的反射地震数据，析取用于在需要用于从源到图像点到接收器的波的评估时间处求和的数据取样，基尔霍夫迁移产生图像。通过对高频eikonal方程用于超过十进位的一些解法已经评估了这种传播时间。解决eikonal方程的使用第一报告的“逆风”有限差分方法是Wiggins，et，al，1986给出的。在基尔霍夫深度迁移中使用高频传播时间引起的问题是熟知的，并从早期的1990s(Van Trierand色，1991；Vinje，et al，1993；Geoltrain andBrac，1993；Audebert，et，al，1994a).Biondi，1992中做出了深入的讨论，其通过实施速度模型的频率相关平滑设法改进高频方法。
基尔霍夫迁移过程中的振幅保持是计算上昂贵的。对于振幅保持，必须对迁移轨迹加权合适的迁移振幅因子。可以参考1993年12月由Tygel，published in Geophysics，vol.59，No.12给出的文章，其是使用的测量技术。例如，必须计算等于包括在迁移中的取样的数量的权重的数量，并且这是计算上昂贵的。
在1996年Stanford，Dave Nichols的学术论文中给出了计算的单色2D格林函数，并在试图避免高频假设的限制中实施2D迁移。用于给出的差分方程的格林函数是对具有如源的空间德耳塔函数的单色波场的解法。Nichols利用在极坐标中的有限差分方法计算单色2D格林函数。Nichols正在寻找传播时间以计算在时间上是德耳塔函数的格林函数的近似值。通过处理单色格林函数作为正弦波的近似值，可以计算多重频率并将作为结果的正弦波叠加在一起作为预叠加函数。可以选择该函数峰值作为德耳塔函数的位置。
先前技术给出了不精确计算通过地球模型具有变化的速度结构的传播时间的高频方法和系统。先前技术给出了计算迁移振幅因子的方法，其是计算上昂贵的。因此存在一种用于该方法和系统的需要比先前技术方法更有效和精确的计算传播时间。存在一种用于该方法和系统的需要有效计算迁移振幅因子以精确的迁移地震图像和表示地球结构。
一种用于测定在多个地下位置点内部的从表面位置到至少一个计算点(CP)的地震波的传播时间的方法。定义区域的速度模型。利用速度模型，对于该至少一个CP，计算复合单色波场。在最接近所述多个位置点中的至少一个CP的至少一个邻近计算点(NCP)处测定地震波的展开相位。利用展开的NCP相位和用于至少一个CP的计算的单色波场到至少一个CP的地震波的传播时间。从测定的传播时间值可以产生传播时间图或表。
该专利的文本包括至少一副带有颜色的图基于请求和必须成本的支付通过专利和商标局将提供具有颜色图的该专利的副本。通过参考下面的详细描述和涉及的图可以更好的理解本发明和它的优点，其中图1A给出了地球速度模型。
图1B给出了计算栅格的举例，可以使用该计算栅格对从速度模型得到的传播时间进行构图。
图1C给出了相对于图1B的栅格计算点的近邻计算点的选择结构。
图1D给出了传播时间图如何计算与涉及用于本发明的CP的展开相位结合的图1A的简单速度模型的图1B的栅级的示意性举例。
图1E给出了展开处理的结构的举例。
图2给出了本发明优选实施例的流程图。
图3A给出了用于模型1的速度模型。
图3B给出了用于模型2的速度模型。
图4给出了波场和用于模型1的得到的传播时间。
图5给出了波场和用于模型2的传播时间。
图6给出了在Marmousi模型中具有最大振幅传播时间覆盖的30Hz波场的实部。
图7给出了用于Marmousi模型的速度模型。
图8给出了利用从30Hz点源波场获得的传播时间的迁移Marmousi数据预叠加。
结合它的优选实施例给出本发明的同时，可以理解本发明不受该处的限制。相反，通过附件权利要求的限定，其期望覆盖可以包括在本发明的精神和范围中的所有替换、修改以及等价物。
本发明是一种用于测定传播时间和产生用于迁移地震数据的传播时间图和传播时间表的方法。本发明的技术提供了对传播时间有效和精确的构图，同时提供了对迁移振幅因子的有效计算。因此，利用本发明的传播时间图得到的迁移图像在位置和振幅上是更连续和更精确的。基于下面的详细描述本发明的其它优点将更容易清楚地展现于本领域的熟练技术人员面前。就下面的详细描述来说，其对具体实施例和本发明的具体利用是清楚的，其是用于描述而不是进行解释以限定本发明的范围。
下文中描述的本发明提供了用于计算传播时间的改进方法，该方法避免了先前高频技术的缺点。在适用于地震数据的单独、有限频率或少数频率处计算传播时间。这些传播时间比高频方法能给出更好的地震图像。此外，可将地震数据过滤成几个频率范围，通过选择频率个数将每个子范围分别迁移成在它的子范围中心处计算的传播时间图和传播时间表，从而给出连续质量/成本权衡。本发明的方法比对渐进高频传播时间实施基尔霍夫求和法也进行了改进的其它方法在计算上更加有效。
本发明的一个优点是避免了峰化时间域格林函数的需要。传播时间涉及单色格林函数的相位。可从单一单色格林函数直接获得传播时间图或传播时间表。从单一频率或选择出的频率计算传播时间，以给出用于那个频率或选择出的频率的成像问题的精确解法(使用在精确度之内的方法以计算格林函数)。对传播时间进行近似修正以用于邻近频率。
单色波场也给出了可以用以计算迁移振幅因子的振幅。由射线跟踪计算振幅通常是不连续和无规律可变的。在有限频率处获得的振幅给出了平滑的迁移图像和更精确的反射图像。
图2给出了本发明优选实施例的流程图。201定义了以阵列形式给出了在利益区域内部的栅级上速度分布的速度模型。203选择了源位置。205计算用于该源位置的波场。在优选实施例中，利用单向声波方程计算波场。可替代的实施例包括对声波方程和弹性波方程的单向和双向解法。通过有限元素或有限差别方法解答每个方程，并且其可以是隐含或明显的。计算的波场给出了含混、包裹的相位和在每个计算点处的振幅。可以测定初始传播时间和初始展开相位值，并可将其设计成阵列的任何节点。例如，在源位置处的展开相位和传播时间是0。207选择在已经计算的展开相位处毗邻邻近计算点(NCP)的计算点(CP)。初始时，源点是单独的候选NCP。可是，当计算更加展开的相位时，CP的选择从源点向外移动，并且选择的CP将具有几个候选NCP。通过涉及最近的邻近计算点(NCP)的先前展开相位测定到达CP的地震波的展开相位。209选择NCP以提供用于获得CP展开相位的参考展开相位。在优选实施例中，从在一个或多个NCP211处的展开相位获得参考相位。通过任何合适的方法获得或选择用于测定CP展开相位的参考相位下面给出举例。参考相位用于测定CP展开相位213。利用NCP参考相位测定的CP展开相位用于测定CP215处的传播时间。在本发明的可替代的实施例中，利用在215处计算的传播时间可以测定传播时间图和传播时间表。
计算用于在图2中给出的优选实施例的传播时间的本质步骤是测定展开相位211；也就是将在范围-π百到π中的周期相位转换成与传播时间成比例的绝对相位。
在本发明中，从单一单色格林函数或合适的几个选择的频率计算传播时间。从单一单色格林函数计算传播时间，其中s(r，t)是用于声波方程的单色点源激励s(r,t)=δ(r)eiw0t,---(1)]]>其中，r=|x→-x0→|]]>是2维或3维的点源、 和任意表面下的点、 之间的距离。激励的频率是ω0，时间是t，以及σ( )是迪拉克德耳塔函数。
现在， 是用于激励s(r，ω)＝FT(s(r，t))的声波的解法▿2U+ω2v2U=s(r,ω)=δ(r)δ(ω-ω0).---(2)]]>通过用于对声波方程的计算解法的任何熟知方法，可以对于给出的速度、 计算U。该解法可以分成振幅和相位U=A(x→,x→0,ω)eiΦ(x→,x→0,ω),---(3)]]>其中，A(x→,x0→,ω)=|U|]]>是振幅，以及Φ(x→,x0→,ω)=arctan(Im(U)/Re(U))]]>是相位。没有近似值包括在U的除法中其仅仅是用于代表复合函数的替代方法。通过傅里叶逆变换U(r，ω)获得时间域解法。由于U具有唯一的单一偏率，它的傅里叶逆变换简单的是U(x→,x→0,t)=Aei(Φ-ω0t).---(4)]]>当指数是零时，(4)中波场的峰值是φ-ω0t＝0(5)其意思是t=Φω0---(6)]]>可是，就给它增加任何多重的2π时它将给出相同的波场来说，φ是含混的。在t＝φ/ω0变成用于计算传播时间的有用公式之前不得不移除含混性。如果在半波长或更小的波长的间距栅级处计算波场，然后可以使用从一个计算点到邻近点的相位的限定变化移除含混性。
由计算的解释在源点处知道波场的相位和传播时间φsource＝tsource＝0然后由于我们知道邻近源的计算点处的波场需要正传播时间，并且其小于来自于源的波长，我们可以选择它们的绝对或展开相位以包括在范围
中设置的差别φ-φsource的多重的2π。
一旦已将毗邻源的点指派为展开相位，下一个邻近它们的点可以接收展开任务。从源点向外迭代处理该任务直到已将所有的计算点指派为展开相位。在每个迭代处，从与被讨论的点邻近的点的先前展开相位计算参考相位、φref。选择多重的2π用于被讨论的点的相位、φ，以使φ-φref存在于范围
中。可替代的，可使在点处的展开的相位的预测作为φprediction＝φref+dφ，其中计算dφ作为从参考点到预测点的局部传播时间。然后可以选择多重的2π以设置在范围[-π，+π]中的参考相位、φ-φprediction。
使用振幅 以提供在基尔霍夫迁移中需要的振幅因子。通过恒定速度事件提供计算来自于单色解法的传播时间的简单举例。通过解析获得解法为U=A0rei(ωvr-ωt).---(7)]]>传播时间简单的为t=Φω=ωvrw=rv.]]>该时间、t是直线传播时间。
通过施加在上面给出的用于指派对每个计算点的展开相位的迭代方法，通过对在传播时间图中对邻近计算点的分析，就实现了通过相位展开(图2中的211)解决相位含混性。单色波场的计算获得了在传播时间图中的每个计算点(CP)处的复合波函数。在图1B中代表的任何栅级点GP可以是CP。
利用复合函数的实部和虚部可以完成计算。下面给出了了波函数、U的相位、φΦ=arctan(Im(U)Re(U)).---(9)]]>就给它增加任何多重的2π时将不会引起波函数的变化来说，φ是含混的。因此，为了使用相位作为利用方程6的传播时间计算的基础，必须解决含混性。
该含混性的解决需要在CP处在复合值之上的附加信息。该附加信息来自于NCP、波场的空间演化和用于限定传播时间的辅助选择。
从相对于图1C中的CP的结构的下面非详尽举例选择用于图2的209中选择的NCP组的举例。结构101中的NCP在深度上位于101的CP之上。给出的计算栅级节点是显示装置101和107中“垂直结构”的举例。作为需要向左/右颠倒“拐角结构”以设置距离源点最近一测NCP(在深度上该NCP近邻于CP)。当使用拐点结构时，使用垂直结构以计算一行新深度CP上的第一CP。其它的结构对本领域的精通技术人员来说是显然的。这些结构是直接从2D到3D外推。
在离散的数字计算中，在CP处计算的波场是足够的紧密，从而测定从一个CP到它的邻近点的绝对相对相位变化。在有限积分计算中，必须将CP分隔成相隔小于半波长，这样将相位差限定在小于π。在一些傅里叶方法中，CP必须间隔大于半波长，但是计算方法采用在那些CP之间的相位的已知的函数相关。因此，在一组NCP处的展开相位，通过需要其在来自于相位φCP、来自于NCP的预测相位的±π的相位差别可以测定在CP处波的相位。
当从与来自于另一个NCP的预测展开相位不同的一个NCP预测展开相位(图2中的211)时引起了第二含混性。通过在可能性之中的选择可以解决该含混性。从选择的结构之中的最强振幅NCP可以预测展开相位，或形成从振幅加权NCP获得的预测值。可替代的，可以对沿着从源点到CP的线最接近的NCP给与选择权。例如，可在清楚的有限差分向下连续波传播中使用这些方法。在结合测定CP的相位的下面举例中讨论解决NCP含混性。
图1D给出了传播时间图如何计算与涉及用于本发明的CP的展开相位结合的图1A的简单速度模型的图1B的栅级的示意性举例。射线从源位置SP经过具有速度V1的层，并到达具有速度v2的层，在具有速度v2的层它到达计算点CP。在图1D中给出的CP同时在该例子中，其也是一个反射点，但是任何格点可以是CP。轮廓箱103包围了具有3NCP的CP的选择的拐角结构。
图1C给出了几个NCP结构101、103、105和107。在波场的直角栅格、向下延拓计算中，这些结构是有用的。如果相位展开在增加深度的水平行中是连续的，并是朝向来自于行内部的源位置的两侧是连续的计算，然后在图1C中给出的NCP结构常常具有在用于给出的CP的所有NCP处可用的展开相位。
可以使用最大振幅预测值以解决用于利用了NCP的选择设置的图2的213中CP相位含混性，例如来自于图1D的结构103。在CP处预测的相位是ΦCP=ΦMAX+ω·drv,---(10)]]>其中，φMAX是参考相位，在NCP处具有最大振幅的先前展开相位，以及dr是分离CP的距离和NCP具有最大振幅。在CP处的相位是通过具有增加了多重的2π以将其设置于范围φCP±π的方程9给出的展开相位。
同样给出了在图2中用于211在NCP中展开处理的结构的举例，波场的局部方向的特征在于在图1E中的波数向量k，其中k是表现了波方向的局部评估的波数。在用以测定在节点CP处的相位、φCP的振幅预测方法中，使用具有结构103的最大振幅的三个节点中的一个节点作为参考相位、φMAX，并且用于在图2中给出的优选实施例的预测相位213是ΦCP=ΦMAX+k·dx→,---(11)]]>其中，k是来自于NCP的评估，并且 是从参考点到CP的向量。可替代的平面波方法将选择取代所有三个NCP的两个左侧NCP节点的最大振幅，其中k在图1E结构103之间通过。这是对沿着从源点到CP的线最接近的NCP给出的参考的举例。同时这些是2D举例，它是直通的以延伸这些事件至3D。
可以选择加权振幅最小平方平面波预测值以解决相位含混性。举例包括诸如在图1C中的CP最直接之上的五个NCP的垂直结构107的一组NCP，或作为在最接近于源点的CP的一侧的正NCP的CP最直接之上的三个NCP的拐角结构105。在用于参数φ0的2D事件中的加权线性最小平方解法中利用这些NCP相位和振幅，参考相位同样来自于图2，kx，侧向尺寸，以及kz，深度尺寸。给出预测相位φCP＝φ0+kx·dx+kz·dz (12)
这可以延伸至3D，φCP＝φ0+kr·dr，其中r延伸至任何选择的坐标系统。
同样，通过选择来自于图1C的结构可以使用用于解决相位含混性的简单平面波预测值。例如，NCP组的结构103包括直接在CP之上的NCP，在CP一侧的NCP最接近于源点SP，并且通过CP和第一两个NCP(也就是在图1C和图1D和图1E中的103)限定NCP完全直角。这些三个NCP的相位唯一的限定了同样在方程12中的平面波预测值。
通过选择任何组NCP可以使用惠更斯子波预测值。构造ΦCP=ΣNCPs(ΦNCP+ω·drv).---(13)]]>也可以加权子波和，也就是通过振幅ΦCP=ΣNCPsANCP(ΦNCP+ω·drv).---(14)]]>可以使用向下延拓波推断方法用于传播时间计算。首先，计算用于单一频率的波场、U。然后，在每个深度步骤处，计算来自于方程9的局部相位作为φ＝tan-1(Im(U)/Re(U))。
利用比如C语言atan2()函数计算角度，其返回在范围-π到π中的值。如果值是负的，可以增加2π以给出用于局部相位的0到2π的范围。返回的值可以与以在先前深度步骤处的源位置的方向最接近的栅格单元处已测定的相位进行比较。可对局部相位增加多重的2π以给出比在比较的位置处的球形相位大0到2π的球形相位。除了这些非常简单的相位预测值，多种预测值给出了同样或类似的导致了速度模型的传播时间。
这里用于图3A，3B，4和5讨论的例子使用了下面的参数计算的栅格尺寸nx＝nz＝100。
计算的栅格间距dx＝dz＝10m。
频率10Hz。
速度在图3A中的模型1在x和y中总计1000m。速度恒定贯穿模型，1000m/s。
在图3B中的模型2两个500m厚度层。上部层是1000m/s和下部层是2000m/s。
图4给出了U401的实部和用于恒定速度模型、来自于图3A的模型1的展开相位。展开相位如来自于方程8的预期相位随大约恒定的倾斜增加。如锥形面给出了从用于恒定速度模型的波场相位获得的传播时间。从水平平面401向上标绘时间。在有颜色的图401中的等相等高线给出了在缺乏任何速度不连续时在深度上波的稳定传播。
图5给出了用于速度模型、来自于图3B的模型2的波场501和获得的传播时间503，其中在最大深度的一半处加倍速度，在505给出了速度不连续。以在水平平面501之上的有颜色的图示出计算的波场的实部。以具有通过在水平片面之上的高度代表的时间的表面503给出了从波场相位得到的传播时间。
图5给出了具有水平速度不连续505从模型2得到的结果。在不连续505之下的速度是在不连续之上它的值的两倍。如预期的在不连续之下对传播时间表面503的倾斜进行均分。当它遇到速度不连续505时，有颜色的图501中的等相等高线507给出了波褶皱。
图6给出了在Marmousi模型中具有最大振幅传播时间覆盖的30Hz波场的实部。释放Marmousi合成数据集(Bourgeois等，Proceedings of the 1990 EAEG Workshop on Practical Aspects ofData InversionEur.Assoc.Expl.Geophys.，5-16，1991)作为用于速度评估的测试。常常使用数据集用于测试迁移和速度地震数据处理运算法则。它是基于Angola中Cuanza盆地的地质概况的复杂声波2-D数据。在图7中给出了速度模型。通过生长断层支配结构形式，生长断层起因于盐沼蠕变，并引起在模型的上部部分中的复杂速度结构。主要成像目标是在盐沼之下的背斜结构中的水库。
图8给出了具有Marmousi数据使用的本发明的方法。利用从30Hz点源波场获得的传播时间已对Marmousi数据进行迁移预叠加。在模型的主要成像目标处的反射的底部处的反射的清晰度和质量描述了迁移的保真度。
可将单色波场表述成具有相对高振幅的河川，例如在图6中603是高振幅河川，其中相位具有平滑单调增加。在具有相对低振幅代表数据区域的河川之间可以有明显的缝隙。这些缝隙是多重波至交迭的地方，产生了相消干扰，并且这样仅产生了弱反射。相位展开发现在高振幅区域中的合适的展开相位，并且使用这些相位作为在低振幅时缝隙中的参考，可以使用其改进展开处理的稳定性。
用于计算宽缝隙波场的可替代方法是为了计算在展开多面层上的波场，例如，由于对二十面体的每一面的法线大约是径向的，在二十面体上限定的栅格是有用的。利用有限差分方法或其它数字方法在大约球形多面层上可以计算波场。
来自于波方程的单色解法的波传播时间的计算可以采取几种不同的形式。声波方程的解法可以包括声波方程的单向解法和声波方程的双向解法。本发明也可以使用弹性波方程的解法。弹性波方程解法也包括单向解法和双向解法。数字方法解法包括有限差分方法、初始元素方法和傅里叶变换方法、以及其它变换。通过分析先前展开相位和一组邻近计算点(NCP)的展开相位用于展开在给出的计算点(CP)处波的相位的方法包括选择来自于源点和CP的相对位置的NCP，通过分析传播时间或展开相位的局部梯度选择NCP，发现最高的振幅NCP，发现NCP的最小展开相位，对NCP拟合平面波，对NCP拟合加权振幅或另外的加权平面波，并且实施振幅或从NCP获得另外的加权惠更斯子波计算。
从用于单色波场解法计算的振幅可以获得用于基尔霍夫迁移的迁移振幅或振幅加权参数。利用如在这里给出的计算的传播时间可以实施基尔霍夫迁移。利用在选择的频率处或在那个频带中的频率处用于计算的每个频带的传播时间在多重频带中可以实施基尔霍夫迁移。
从单色源点波场或格林函数可以计算用于深度迁移的传播时间图或表。这些传播时间图对在计算的频率处的波方程是精确的解法，并且与以高频渐近方法计算的传播时间相比，是用于邻近频率的更好的近似值。
在地震频带中的少数频率处可以计算传播时间。可将地震数据分解成在这些频率周围的频带，并且可以以它的联合的传播时间对每个频带进行迁移。由于处理是线性的，通过对这些带限图像求和可以恢复完整的地震图像。通过选择为了使用的频带的数量，我们可以期望图像质量中的一个增量的改进具有在计算成本中的增量的增加。
作为本发明的最后概括，定义速度模型并选择源位置。对于每个CP计算复合波场。在每个CP处计算展开相位，例如利用反正切函数。在源点处开始迭代循环。初始时，由于它是零，已知展开相位在源点处是单一的。由于它是具有在迭代处的展开相位的单一NCP，专用迭代以用于邻近源点的CP。由于在它们没有行，专用迭代以用于CP阵列的上部行。如果缺少普通事件CP中的一些，在模型的边缘处可以专用迭代。
在已经处理的专用事件之后使用普通事件。选择已具有展开相位的CP的CP。通过合适的设置中的具有展开CP的必要性支配选择CP的顺序。从源点径向地或多或少的继续该顺序。计算来自于NCP的展开相位的参考相位。从参考相位和传播延迟在CP处计算预测的展开相位。选择多重的2π以增加给CP的展开相位，从而使作为结果的展开相位尽可能的接近预测的展开相位，或者交替地，其尽可能的接近，但是差分＞0。
本领域的熟练技术人员可以理解如给出的，可以实践这里给出的方法，其包括但不受给出的实施例的限制。进一步的，可以理解本发明不受已经给出的用于描述作用的前面所述的不正当限制。如在下面权利要求定义的，只要不脱离本发明的真实范围，各种修改和替代对本领域的熟练技术人员是明显的。
权利要求
1.用于测定在多个地下位置点内部的从表面位置到至少一个计算点(CP)的地震波的传播时间的方法包括(a)定义包括所述多个地下位置点的区域的速度模型；(b)利用速度模型，对于该至少一个CP，计算单色波场；(c)在最接近所述多个位置点中的至少一个CP的至少一个邻近计算点(NCP)处测定地震波的展开相位；以及(d)对于该至少一个CP，从所述测定的展开相位和所述的计算单色波场，测定到该至少一个CP的地震波的传播时间。
2.权利要求1的方法，其中所述表面位置是至少一个i)源位置，和ii)接收器位置。
3.权利要求1的方法，其中计算用于该至少一个CP的单色波场进一步包括对至少一个i)单向弹性波方程，ii)双向弹性波方程，iii)单向声波方程，iv)双向声波方程测定单色解法。
4.权利要求3的方法，其中计算单色波场进一步包括以至少一个i)有限元方法，ii)有限差分方法，iii)显式方法，iv)隐式方法解决波方程解法。
5.权利要求3的方法，其中计算单色波场进一步包括使用至少一个格林函数。
6.权利要求3的方法进一步包括从计算的单色波场获得基尔霍夫迁移振幅因子。
7.权利要求1的方法，其中至少一个CP包括多个CP。
8.权利要求1的方法，其中至少一个NCP包括多个NCP。
9.权利要求1的方法进一步包括从至少一个计算的传播时间测定I)传播时间图，和ii)传播时间表。
10.权利要求1的方法，其中至少一个NCP包括多个NCP，该方法进一步包括通过形成多个NCP的相位的加权和来测定多个NCP的相位。
11.权利要求1的方法，其中测定至少一个CP的展开相位包括将至少一个CP的展开相位限定于用于至少一个NCP的测定相位的±π之内。
12.权利要求1的方法，其中测定传播时间进一步包括计算多面体表面之上的单色波场。
13.权利要求1的方法，其中测定传播时间进一步包括计算是二十面体的多面体表面之上的单色波场。
14.权利要求1的方法，其中测定到至少一个CP的地震波的传播时间进一步包括i)分析传播时间的局部梯度，ii)分析展开相位的局部梯度，iii)测定至少一个NCP的最大振幅，iv)测定至少一个NCP的最小展开相位，v)对至少一个NCP拟合平面波，vi)对至少一个NCP拟合加权平面波，以及vii)计算来自于至少一个NCP的加权惠更斯子波，中的至少之一。
15.用于测定在多个地下位置点内部的从表面位置到至少一个计算点(CP)的地震波的传播时间的方法包括(a)定义包括所述多个地下位置点的区域的速度模型；(b)利用速度模型，对于该至少一个CP，计算单色波场；(c)在最接近所述多个位置点中的至少一个CP的至少一个邻近计算点(NCP)处测定地震波的展开相位；(d)对于该至少一个CP，从所述测定的展开相位和所述的计算单色波场，测定到该至少一个CP的地震波的传播时间，以用于；以及(e)以图的形式存储传播时间值。
16.权利要求15的方法，其中所述表面位置是i)源位置，和ii)接收器位置中的至少之一。
17.权利要求15的方法，其中至少一个CP包括多个CP。
18.权利要求15的方法，其中至少一个NCP包括多个NCP。
19.权利要求15的方法，其中至少一个NCP包括多个NCP，该方法进一步包括通过形成多个NCP的相位的加权和来测定多个NCP的相位。
20.权利要求15的方法，其中测定至少一个NCP的展开相位包括将至少一个CP的展开相位限定于用于至少一个NCP的测定相位的±π之内。
21.权利要求15的方法，其中使用至少一个格林函数计算测定的单色波场。
22.权利要求15的方法进一步包括从计算的单色波场获得基尔霍夫迁移振幅因子。
23.权利要求15的方法，其中测定单色波场进一步包括计算多面体表面之上的波场。
24.权利要求15的方法，其中测定单色波场进一步包括计算是二十面体的多面体表面之上的波场。
25.权利要求15的方法，其中测定单色波场进一步包括计算i)单向弹性波方程，ii)双向弹性波方程，iii)单向声波方程，iv)双向声波方程中至少之一的解。
26.权利要求15的方法，其中测定单色波场进一步包括以i)有限元方法，ii)有限差分方法，iii)显式方法，iv)隐式方法中至少之一解决波方程解法。
27.权利要求15的方法，其中测定到至少一个CP的地震波的传播时间进一步包括i)分析传播时间的局部梯度，ii)分析展开相位的局部梯度，iii)测定至少一个NCP的最大振幅，iv)测定至少一个NCP的最小展开相位，v)对至少一个NCP拟合平面波，vi)对至少一个NCP拟合加权平面波，以及vii)计算来自于至少一个NCP的加权惠更斯子波，中的至少一个。
28.用于测定在多个地下位置点内部的从表面位置到至少一个计算点(CP)的地震波的传播时间的方法包括(a)定义包括所述多个地下位置点的区域的速度模型；(b)利用速度模型，对于至少一个CP，对于至少一个选择频率计算波场，；(c)在最接近所述多个位置点中的至少一个CP的至少一个邻近计算点(NCP)处测定地震波的展开相位；以及(d)对于该至少一个CP，从所述测定的展开相位和所述的计算单色波场，测定到该至少一个CP的地震波的传播时间。
29.权利要求28的方法，其中所述表面位置是i)源位置，和ii)接收器位置中的至少一个。
30.权利要求28的方法，其中对于该至少一个CP计算波场进一步包括对于至少一个选择频率测定i)单向弹性波方程，ii)双向弹性波方程，iii)单向声波方程，iv)双向声波方程中至少之一的波场解法。
31.权利要求30的方法，其中计算波场进一步包括以i)有限元方法，ii)有限差分方法，iii)显式方法，iv)隐式方法中的至少之一解决波方程解法。
32.权利要求30的方法，其中计算波场包括进一步包括使用用于至少一个选择频率的至少一个格林函数。
33.权利要求30的方法进一步包括从计算的波场获得基尔霍夫迁移振幅因子。
34.权利要求28的方法，其中至少一个CP包括多个CP。
35.权利要求28的方法，其中至少一个NCP包括多个NCP。
36.权利要求28的方法，进一步包括从至少一个计算传播时间测定i)传播时间图，和ii)传播时间表中的至少之一。
37.权利要求28的方法，其中至少一个NCP包括多个NCP，该方法进一步包括通过形成多个NCP的相位的加权和来测定多个NCP的相位。
38.权利要求28的方法，其中测定至少一个CP的展开相位包括将至少一个CP的展开相位限定于用于至少一个NCP的测定相位的±π之内。
39.权利要求1的方法，其中测定传播时间进一步包括对于至少一个选择频率计算在多面体表面之上的波场。
40.权利要求1的方法，其中测定传播时间进一步包括对于至少一个选择频率计算在二十面体的多面体表面之上的波场。
41.权利要求15的方法，其中测定到该至少一个CP的地震波的传播时间进一步包括i)分析传播时间的局部梯度，ii)分析展开相位的局部梯度，iii)测定至少一个NCP的最大振幅，iv)测定至少一个NCP的最小展开相位，v)对至少一个NCP拟合平面波，vi)对至少一个NCP拟合加权平面波，以及vii)计算来自于至少一个NCP的加权惠更斯子波，中的至少之一。
全文摘要
一种用于测定在多个地下位置点内部的从表面位置到至少一个计算点(CP)的地震波传播时间的方法。定义区域的速度模型。利用速度模型，计算复合单色波场，以用于至少一个CP。在最接近所述多个位置点中的至少一个CP的至少一个邻近计算点(NCP)处测定地震波的展开相位。利用展开的NCP相位和用于至少一个CP的计算的单色波场测定到至少一个CP的地震波的传播时间。从测定的传播时间值产生转播时间图和表。
文档编号G01V1/00GK1613021SQ02826701
公开日2005年5月4日 申请日期2002年12月30日 优先权日2002年1月4日
发明者詹姆斯·W.·威金斯 申请人:维斯特恩格科有限责任公司