常用叠前深度偏移方法特点分析与实例对比

目前常用叠前深度偏移主要有四种:基于射线理论的Kirchhoff偏移和高斯束偏移以及基于波动方程理论的单程波偏移和逆时偏移.要想在实际应用中取得好的成像效果,必须选取合理、有效的叠前深度偏移算法.本文首先简单阐述了这四种常用叠前深度偏移方法的基本原理,对比分析了四种方法在陡倾角、速度依赖性、起伏地表、非规则数据、低信噪比数据等方面的优劣性以及计算效率和计算精度的差异,然后利用模型试算和实际资料处理分析说明这四种叠前深度偏移方法的成像效果,最后给出不同情况下叠前深度偏移方法的选取建议.     

苏北大陆科学钻探靶区深反射地震的叠前深度偏移

由于深反射地震数据具有信噪比低和记录长度长等特点,叠前深度偏移方法的应用有许多困难．为此,我们研究了一种适合于深反射地震的叠前深度偏移方法;包括：用逆风有限差分方法计算程函方程;在常规速度扫描的基础上,用协方差控制提高速度分析精度;用联合反演算法计算层速度,再插值后得到初始速度模型;用Ｋｉｒｃｈｈｏｆｆ法作为偏移速度分析工具,求得最终的速度模型;最终的速度模型作为有限差分深度偏移的输入,求得最终的偏移结果．用该方法对“中国大陆科学深钻工程”东海二维深反射地震数据ＤＨ－４线进行了叠前深度偏移,取得了良好的效果。     

TI介质局部角度域高斯束叠前深度偏移成像

各向异性射线理论基础上的局部角度域叠前深度偏移方法能够为深度域构造成像与基于角道集的层析反演提供有力支撑,但是对于复杂地质构造而言,高斯度叠前深度偏移在不失高效、灵活等特点的情况下,具有明显的精度优势.为此,本文研究局部角度域理论框架下的高斯束叠前深度偏移方法.为提高算法效率与实用性,文中讨论了一种从经典弹性参数表征的各向异性介质运动学和动力学射线方程演变而来的由相速度表征的简便形式,并提出了一种比较经济的各向异性高斯束近似合成方案.结合地震波局部角度域成像原理,讨论一种适合高斯束偏移的角度参数计算方法.国际上通用的理论模型合成数据试验表明：相比局部角度域Kirchhoff叠前深度偏移成像方法,本文方法具有更高的成像精度与抗噪能力,既适用于复杂构造成像,也可为TI介质深度域偏移速度分析与模型建立提供高效的偏移引擎.     

基于矩形网格的有限差分走时计算方法

对于大多数速度场，地震波沿射线传播的初至波走时，可以用有限差分外推的方法在二维或三维数值网格上计算出来. 在保证精度的条件下，为提高计算效率和适应性，本文推导了基于任意矩形网格和局部平面波前近似的有限差分初至波走时计算方法. 另外，该方法对首波和散射波做了合适的处理，而且不会碰到传统射线法存在的阴影区和焦散区等问题. 简单模型和复杂的Marmousi模型试算的结果表明，该方法精度较高并适用于强纵、横向变速的复杂介质. 基于该方法的Kirchhoff叠前深度偏移, 在主要构造和目的层位置的成像效果上基本达到了波动方程法叠前深度偏移的位置成像效果. 由于未考虑续至波等有效能量，在成像的保幅性上不如波动方程法叠前深度偏移的效果，但其计算效率则明显高于全格林函数法和波动方程法.      

二维各向异性多波高斯束叠前深度偏移

基于多分量地震资料的深度偏移能够更加充分的利用地震记录上的多波信息,而且地下介质中广泛的存在各向异性,本文提出了一种有效的适用于二维各向异性介质的多波高斯束叠前深度偏移方法.首先基于各向异性射线追踪理论,导出了适用于二维各向异性介质的射线追踪方程组,实现了二维各向异性介质中P波和S波的射线追踪.其次将各向异性射线追踪理论引入到高斯束偏移方法中,分别给出适用于二维各向异性介质的PP波和PS波高斯束叠前深度偏移成像方法.本文成像方法充分考虑各向异性因素对地震波场的影响,能够对存在各向异性介质的地下构造准确偏移归位.通过对不同各向异性介质的数值模型进行PP波和PS波高斯束叠前深度偏移成像测试分析表明,本文方法是一种准确有效的适应于各向异性介质条件的多波叠前深度偏移成像算法.     

两种类型共炮点域相关型叠前深度偏移方法分析

研究不同偏移方法的成像机理是实现复杂构造高精度成像的前提,研究了两类共炮点域的相关型叠前深度偏移成像方法：基于单程波动方程的叠前深度成像方法和基于双程波动方程的叠前深度成像方法,同时对比了它们在计算效率、数据存储量、成像精度、成像机理、速度敏感性等方面的差异及其共性。以复杂构造模型为例,采用了傅里叶有限差分法（FFD）和逆时成像法（RTM）,这两种方法实现了121个共炮点道集的叠前深度偏移成像处理。计算结果表明,当速度准确时,两种深度域波动方程成像方法均可以恢复出各个地质反射界面,其中逆时偏移对陡倾角成像效果显著,当速度存在百分比误差和随机扰动情况时,逆时成像结果要差于单程波方法,因此,逆时偏移方法对速度的敏感性较大,且低频噪声较为严重。     

波动方程深度偏移波场延拓算子的快速重建

基于波动方程的叠前深度偏移对计算机的速度和存储空间都有较高要求.随着并行集群的出现,此类叠前深度偏移问题已经开始应用于工业生产,但与传统偏移方法相比仍然耗时较长.本文应用三次样条函数对波场延拓算子进行光滑处理,然后用抽样函数进行算子重建,既可以保证计算精度,又能减少叠前深度偏移过程所需的计算存储,从而提高效率,缩短整个处理流程的时间.     

最优可分表示法三维波场延拓

适于复杂介质的高精度波场延拓算子是叠前深度偏移研究的重要内容。本文采用最优可分表示方法,运用正反傅立叶变换构造了三维单程波场延拓算子,算子实现了波数域变量与空间(速度)域变量分离。波数域内进行相移计算,在空间域对因介质横向变速引起的时移作修正。脉冲响应显示在区域内各速度的脉冲计算值与理论值基本一致,说明最优可分表示法叠前深度偏移可适用于强变速条件下复杂介质的成像需求。SEG／EAGE模型和实测数据的成像结果验证了本文方法对复杂构造的成像能力。     

波动方程法共成像点道集偏移速度建模

叠前深度偏移的成像效果对偏移速度场相当敏感，建立正确的偏移速度场是实现高质量叠前深度偏移成像的关键。首先应用成像精度高的波动方程法叠前深度偏移抽取共成像点道集；然后基于摄动法通过参数化速度函数和改进的剩余曲率分析建立偏移速度误差和成像深度误差的定量关系；最后采用单参数／多参数联合迭代反演实现偏移速度建模。对Marmousi模型的试算结果表明：该方法对复杂地质体具有较强的适应性和较好的建模和成像效果，一般只需分析和控制主要反射层，通过3－4次近代就可以满足精度要求。     

适于窄线三维地震资料的面炮方法

用小孔径波动方程叠前深度偏移完成大面积数据的成像，不仅总的计算量较大，而且成像畸变较严重。为了进行大孔径波动方程叠前深度偏移，必须将实际的小孔径资料合成大孔径资料，并解决大孔径资料的成像方法。AMO等合成孔径方法，数据合成方式简单，但成像方法是近似的且较复杂。面炮方法从理论上讲数据合成方法较简单，成像方法是精确的且较简单，但其应用适用于远道振幅比近道弱很多的情况。在窄线三维地震资料情况，其联络测线方向不易满足面炮应用条件。针对这种情况，本文将最小二乘原理应用于波场成像，提出适用窄线三维地震资料的面炮方法。为了解三维波动方程叠前深度偏移的面炮方法解决复杂油气藏的可行性、所需计算资源和成像效果，我们利用大庆SW地区三维地震资料进行了实验研究。通过研究，我们获得了面炮方法的实际效果，获得了所需计算资源的测试数据。对实际效果的分析对比表明，面炮方法是一种三维波动方程叠前深度偏移的快速方法，有较好应用前景。     

SEISMIC RAY TRACING USING LINEAR TRAVELTIME INTERPOLATION1

A new ray-tracing method called linear traveltime interpolation (LTI) is proposed. This method computes traveltimes and raypaths in a 2D velocity structure more rapidly and accurately than other conventional methods. The LTI method is formulated for a 2D cell model, and calculations of traveltimes and raypaths are carried out only on cell boundaries. Therefore a raypath is considered to be always straight in a cell with uniform velocity. This approach is suitable to tomography analysis. The algorithm of LTI consists of two separate steps: step 1 calculates traveltimes on all cell boundaries; step 2 traces raypaths for all pairs of receivers and the shot. A traveltime at an arbitrary point on a cell boundary is assumed to be linearly interpolated between traveltimes at the adjacent discrete points at which we calculate traveltimes. Fermat's principle is used as the criterion for choosing the correct traveltimes and raypaths from several candidates routinely. The LTI method has been compared numerically with the shooting method and the finite-difference method (FDM) of the eikonal equation. The results show that the LTI method has great advantages of high speed and high accuracy in the calculation of both traveltimes and raypaths. The LTI method can be regarded as an advanced version of the conventional FDM of the eikonal equation because the formulae of FDM are independently derived from LTI. In the process of derivation, it is shown theoretically that LTI is more accurate than FDM. Moreover in the LTI method, we can avoid the numerical instability that occurs in Vidale's method where the velocity changes abruptly.     

初至波相位走时层析

近地表速度结构通常是利用射线走时层析或菲涅尔体走时层析等反演方法得到的,但它们的目标函数仍利用射线走时残差构建,导致反演精度不高.为此,本文提出了基于散射积分算法的初至波相位走时层析成像方法.该方法的核心是:(1)提出了依赖于频率的相位走时概念;(2)利用依赖于频率的相位走时信息,而非单一的无限频率射线走时;(3)发展了一种改进的相位展开方法,即通过监测相位不连续性和2π周期判定来消除相位折叠现象;(4)考虑了地震波传播的有限频特征,即基于波动理论而非传统的射线路径或有限空间的菲涅尔体构建核函数.通过利用Overthrust模型的数值实验及与传统射线走时层析和菲涅尔体走时层析的对比表明:本文提出的方法是一种有效的初至波走时反演方法.同时,基于Overthrust模型的数值试验还证明了下列结论,即通过挖掘更多的走时信息的确可以获得更高的反演精度和分辨率.     

PS转换波界面二次源法射线追踪

界面二次源法是最近提出的一种最小走时射线追踪方法，尤其适合层状介质中走时和射线路径的计算.该方法相对于传统的最小走时树方法（如Moser法），仅在物性界面上设置二次源，射线路径的方向只在层界面处发生改变，该方法最大程度地消除了射线路径的锯齿状现象，同时也避免了低变速区的射线路径多值现象，因此，它具有更高的追踪精度和效率.本文采用界面二次源法在各向同性介质中实现了PS转换波射线追踪，理论模型的计算证实了界面二次源法追踪PS转换波的准确性和高效性，同时该方法在各向异性介质中也很好地追踪出分离的PSV波和PSH波，因此该方法有利于横波分裂在地震勘探中的研究和应用     

Second-order interpolation of later-arrival traveltimes

The performance of a 3D prestack migration of the Kirchhoff type can be significantly enhanced if the computation of the required stacking surface is replaced by an efficient and accurate method for the interpolation of diffraction traveltimes. Thus, input traveltimes need only be computed and stored on coarse grids, leading to considerable savings in CPU time and computer storage. However, interpolation methods based on a local approximation of the traveltime functions fail in the presence of triplications of the wavefront or later arrivals. This paper suggests a strategy to overcome this problem by employing the coefficients of a hyperbolic traveltime expansion to locate triplications and correct for the resulting errors in the interpolated traveltime tables of first and later arrivals.     

Numerical modeling of global seismic phases and its application in seismic phase identification

Earthquake data include informative seismic phases that require identification for imaging the Earth's structural interior. In order to identify the phases, we created a numerical method to calculate the traveltimes and raypaths by a shooting technique based upon the IASP91 Earth model, and it can calculate the traveltimes and raypaths for not only the seismic phases in the traditional traveltime tables such as IASP91, AK135, but also some phases such as pPcP, pPKIKP, and PPPPP. It is not necessary for this method to mesh the Earth model, and the results from the numerical modeling and its application show that the absolute differences between the calculated and theoretical traveltimes from the ISAP91 tables are less than 0.1 s. Thus, it is simple in manipulation and fast in computation, and can provide a reliable theoretical prediction for the identification of a seismic phase within the acquired earthquake data.     

矩形网格三点Fermat射线追踪技术

矩形网格三点Fermat射线追踪法是基于矩形网格三点扰动法的一种提高计算速度的方法．取矩形网格三个点,在Fermat最小旅行时原则下求取扰动中间点的位置,而不象扰动法那样依次扰动．因此,计算速度比扰动法提高2倍多,同时不受扰动摄动量大小选择的困扰．该方法继承了矩形网格三点扰动法的优点,对任意离散的速度场,总能找到最短时间路径,避免了射线盲区和追踪路径并非时间最短路径问题．     

3D adaptive tomography using Delaunay triangles and Voronoi polygons

The solutions of traveltime inversion problems are often not unique because of the poor match between the raypath distribution and the tomographic grid. However, by adapting the local resolution iteratively, by means of a singular value analysis of the tomographic matrix, we can reduce or eliminate the null space influence on our earth image: in this way, we get a much more reliable estimate of the velocity field of seismic waves. We describe an algorithm for an automatic regridding, able to fit the local resolution to the available raypaths, which is based on Delaunay triangulation and Voronoi tessellation. It increases the local pixel density where the null space energy is low or the velocity gradient is large, and reduces it elsewhere. Consequently, the tomographic image can reveal the boundaries of complex objects, but is not affected by the ambiguities that occur when the grid resolution is not adequately supported by the available raypaths.     

Improving modelling and inversion in refraction seismics with a first-order Eikonal solver

A first-order Eikonal solver is applied to modelling and inversion in refraction seismics. The method calculates the traveltime of the fastest wave at any point of a regular grid, including head waves as used in refraction. The efficiency, robustness and flexibility of the method give a very powerful modelling tool to find both traveltimes and raypaths. Comparisons with finite-difference data show the validity of the results. Any arbitrarily complex model can be studied, including the exact topography of the surface, thus avoiding static corrections. Later arrivals are also obtained by applying high-slowness masks over the high-velocity zones. Such an efficient modelling tool may be used interactively to invert for the model, but a better method is to apply the refractor-imaging principle of Hagedoorn to obtain the refractors from the picked traveltime curves. The application of this principle has already been tried successfully by previous authors, but they used a less well-adapted Eikonal solver. Some of their traveltimes were not correct in the presence of strong velocity variations, and the refractor-imaging principle was restricted to receiver lines along a plane surface. With the first-order Eikonal solver chosen, any topography of the receiving surface can be considered and there is no restriction on the velocity contrast. Based on synthetic examples, the Hagedoorn principle appears to be robust even in the case of first arrivals associated with waves diving under the refractor. The velocities below the refractor can also be easily estimated, parallel to the imaging process. In this way, the model can be built up successively layer by layer, the refractor-imaging and velocity-mapping processes being performed for each identified refractor at a time. The inverted model could then be used in tomographic inversions because the calculated traveltimes are very close to the observed traveltimes and the raypaths are available.     

2D inversion of refraction traveltime curves using homogeneous functions

A method using simple inversion of refraction traveltimes for the determination of 2D velocity and interface structure is presented. The method is applicable to data obtained from engineering seismics and from deep seismic investigations. The advantage of simple inversion, as opposed to ray‐tracing methods, is that it enables direct calculation of a 2D velocity distribution, including information about interfaces, thus eliminating the calculation of seismic rays at every step of the iteration process. The inversion method is based on a local approximation of the real velocity cross‐section by homogeneous functions of two coordinates. Homogeneous functions are very useful for the approximation of real geological media. Homogeneous velocity functions can include straight‐line seismic boundaries. The contour lines of homogeneous functions are arbitrary curves that are similar to one another. The traveltime curves recorded at the surface of media with homogeneous velocity functions are also similar to one another. This is true for both refraction and reflection traveltime curves. For two reverse traveltime curves, non‐linear transformations exist which continuously convert the direct traveltime curve to the reverse one and vice versa. This fact has enabled us to develop an automatic procedure for the identification of waves refracted at different seismic boundaries using reverse traveltime curves. Homogeneous functions of two coordinates can describe media where the velocity depends significantly on two coordinates. However, the rays and the traveltime fields corresponding to these velocity functions can be transformed to those for media where the velocity depends on one coordinate. The 2D inverse kinematic problem, i.e. the computation of an approximate homogeneous velocity function using the data from two reverse traveltime curves of the refracted first arrival, is thus resolved. Since the solution algorithm is stable, in the case of complex shooting geometry, the common‐velocity cross‐section can be constructed by applying a local approximation. This method enables the reconstruction of practically any arbitrary velocity function of two coordinates. The computer program, known as godograf , which is based on this theory, is a universal program for the interpretation of any system of refraction traveltime curves for any refraction method for both shallow and deep seismic studies of crust and mantle. Examples using synthetic data demonstrate the accuracy of the algorithm and its sensitivity to realistic noise levels. Inversions of the refraction traveltimes from the Salair ore deposit, the Moscow region and the Kamchatka volcano seismic profiles illustrate the methodology, practical considerations and capability of seismic imaging with the inversion method.     

复杂地表条件下基于线性插值和窄带技术的地震波走时计算

为了在复杂地表条件下实现地震波走时计算,提出了一种基于线性插值和窄带技术的走时计算新方法.其中,线性插值用于局部走时计算,窄带技术用于局部波前捕获和追踪.为了逼近起伏地表,采用三角网和矩形网相结合的方法对速度模型进行剖分.为了得到局部走时计算公式,利用费马（Fermat）原理和关于入射点位置的限定条件.有关编程实践和数值试验表明：新方法不仅可以有效、灵活地处理地表高程的剧烈变化,而且还具有很好的适应性和稳定性,得到的计算结果满足波前传播规律.