Finite-difference calculation of traveltimes based on rectangular grid

To the most of velocity fields, the traveltimes of the first break that seismic waves propagate along rays can be computed on a 2-D or 3-D numerical grid by finite-difference extrapolation. Under ensuring accuracy, to improve calculating efficiency and adaptability, the calculation method of first-arrival traveltime of finite-difference is de- rived based on any rectangular grid and a local plane wavefront approximation. In addition, head waves and scat- tering waves are properly treated and shadow and caustic zones cannot be encountered, which appear in traditional ray-tracing. The testes of two simple models and the complex Marmousi model show that the method has higher accuracy and adaptability to complex structure with strong vertical and lateral velocity variation, and Kirchhoff prestack depth migration based on this method can basically achieve the position imaging effects of wave equation prestack depth migration in major structures and targets. Because of not taking account of the later arrivals energy, the effect of its amplitude preservation is worse than that by wave equation method, but its computing efficiency is higher than that by total Green′s function method and wave equation method.

Joint tomographic inversion of first-arrival and reflection traveltimes for recovering 2-D seismic velocity structure with an irregular free surface

Irregular surface flattening,which is based on a boundary conforming grid and the transformation between curvilinear and Cartesian coordinate systems,is a mathematical method that can elegantly handle irregular surfaces,but has been limited to obtaining first arrivals only.By combining a multistage scheme with the fast-sweeping method(FSM,the method to obtain first-arrival traveltime in curvilinear coordinates),the reflected waves from a crustal interface can be traced in a topographic model,in which the reflected wavefront is obtained by reinitializing traveltimes in the interface for upwind branches.A local triangulation is applied to make a connection between velocity and interface nodes.Then a joint inversion of first-arrival and reflection traveltimes for imaging seismic velocity structures in complex terrains is presented.Numerical examples all perform well with different seismic velocity models.The increasing topographic complexity and even use of a high curvature reflector in these models demonstrate the reliability,accuracy and robustness of the new working scheme;checkerboard testing illustrates the method's high resolution.Noise tolerance testing indicates the method's ability to yield practical traveltime tomography.Further development of the multistage scheme will allow other later arrivals to be traced and used in the traveltime inversion.

A new scheme for fast calculation of seismic traveltimes and ray paths in heterogeneous media

A new method based on Huygens' principle and Fermat's principle is presented to calculate seismic traveltime and ray-paths in this paper. The algorithm can be easily understood and programmed, and can be applied to heteroge neous media in order to eliminate its disadvantage of slowness, we have improved the basic algorithm to speed its calculation to practical stage without the loss of its accuracy.

An adaptive finite-difference method for seismic traveltime modeling based on 3D eikonal equation

3D eikonal equation is a partial differential equation for the calculation of first-arrival traveltimes and has been widely applied in many scopes such as ray tracing,source localization,reflection migration,seismic monitoring and tomographic imaging.In recent years,many advanced methods have been developed to solve the 3D eikonal equation in heterogeneous media.However,there are still challenges for the stable and accurate calculation of first-arrival traveltimes in 3D strongly inhomogeneous media.In this paper,we propose an adaptive finite-difference(AFD)method to numerically solve the 3D eikonal equation.The novel method makes full use of the advantages of different local operators characterizing different seismic wave types to calculate factors and traveltimes,and then the most accurate factor and traveltime are adaptively selected for the convergent updating based on the Fermat principle.Combined with global fast sweeping describing seismic waves propagating along eight directions in 3D media,our novel method can achieve the robust calculation of first-arrival traveltimes with high precision at grid points either near source point or far away from source point even in a velocity model with large and sharp contrasts.Several numerical examples show the good performance of the AFD method,which will be beneficial to many scientific applications.

基于多声源波速结构成像的岩体异常区域超前辨识方法

异常区域超前辨识对于预防地下岩土工程灾害具有重要作用。为了满足地下工程高精度探测的需求,本文提出一种层析成像方法以辨识复杂岩体结构中的异常区域,结合了走时层析、阻尼最小二乘和高斯滤波等技术。该方法克服了空洞区域辨识中速度差限制,减轻了迭代中孤立速度突变所带来的影响。我们开展了数值和室内实验量化评估最短路径法(Shortest-Path Method,SPM)、动态最短路径法(Dynamic Shortest-Path Method,DSPM)和快速扫描法(Fast Sweeping Method,FSM)等正演模拟的识别精度和计算效率。结果表明,在数值和室内实验中DSPM和FSM均能清晰地辨识出异常区域。陕西震奥矿山现场应用结果证明了该方法可利用矿山开采中爆破、微震等多类声源对矿山内部未知结构进行波速场成像。本研究不仅实现了走时层析成像方法在异常区域识别中的应用,而且为地下岩土工程中潜在风险源的探测提供了新的思路。

A review of the influencing factors on teleseismic traveltime tomography

Teleseismic traveltime tomography is an important tool for investigating the crust and mantle structure of the Earth.The imaging quality of teleseismic traveltime tomography is affected by many factors,such as mantle heterogeneities,source uncertainties and random noise.Many previous studies have investigated these factors separately.An integral study of these factors is absent.To provide some guidelines for teleseismic traveltime tomography,we discussed four main influencing factors:the method for measuring relative traveltime differences,the presence of mantle heterogeneities outside the imaging domain,station spacing and uncertainties in teleseismic event hypocenters.Four conclusions can be drawn based on our analysis.(1)Comparing two methods,i.e.,measuring the traveltime difference between two adjacent stations(M1)and subtracting the average traveltime of all stations from the traveltime of one station(M2),reveals that both M1 and M2 can well image the main structures;while M1 is able to achieve a slightly higher resolution than M2;M2 has the advantage of imaging long wavelength structures.In practical teleseismic traveltime tomography,better tomography results can be achieved by a two-step inversion method.(2)Global mantle heterogeneities can cause large traveltime residuals(up to about 0.55 s),which leads to evident imaging artifacts.(3)The tomographic accuracy and resolution of M1 decrease with increasing station spacing when measuring the relative traveltime difference between two adjacent stations.(4)The traveltime anomalies caused by the source uncertainties are generally less than 0.2 s,and the impact of source uncertainties is negligible.

Seismic traveltime inversion of 3D velocity model with triangulated interfaces

Seismic traveltime tomographic inversion has played an important role in detecting the internal structure of the solid earth. We use a set of blocks to approximate geologically complex media that cannot be well described by layered models or cells. The geological body is described as an aggregate of arbitrarily shaped blocks,which are separated by triangulated interfaces. We can describe the media as homogenous or heterogeneous in each block. We define the velocities at the given rectangle grid points for each block,and the heterogeneous velocities in each block can be calculated by a linear interpolation algorithm. The parameters of the velocity grid positions are independent of the model parameterization,which is advantageous in the joint inversion of the velocities and the node depths of an interface. We implement a segmentally iterative ray tracer to calculate traveltimes in the 3D heterogeneous block models.The damped least squares method is employed in seismic traveltime inversion,which includes the partial derivatives of traveltime with respect to the depths of nodes in the triangulated interfaces and velocities defined in rectangular grids. The numerical tests indicate that the node depths of a triangulated interface and homogeneous velocity distributions can be well inverted in a stratified model.

Reflection-based traveltime and waveform inversion with second-order optimization

Reflection-based inversion that aims to reconstruct the low-to-intermediate wavenumbers of the subsurface model, can be a complementary to refraction-data-driven full-waveform inversion(FWI), especially for the deep target area where diving waves cannot be acquired at the surface. Nevertheless, as a typical nonlinear inverse problem, reflection waveform inversion may easily suffer from the cycleskipping issue and have a slow convergence rate, if gradient-based first-order optimization methods are used. To improve the accuracy and convergence rate, we introduce the Hessian operator into reflection traveltime inversion(RTI) and reflection waveform inversion(RWI) in the framework of second-order optimization. A practical two-stage workflow is proposed to build the velocity model, in which Gauss-Newton RTI is first applied to mitigate the cycle-skipping problem and then Gauss-Newton RWI is employed to enhance the model resolution. To make the Gauss-Newton iterations more efficiently and robustly for large-scale applications, we introduce proper preconditioning for the Hessian matrix and design appropriate strategies to reduce the computational costs. The example of a real dataset from East China Sea demonstrates that the cascaded Hessian-based RTI and RWI have good potential to improve velocity model building and seismic imaging, especially for the deep targets.

2-D crustal structure of Changbaishan-Tianchi volcanic region determined by seismic traveltime inversion

2-D velocity structure and tectonics of the crust and upper mantle is revealed by inversion of seismic refraction and wide-angle reflection traveltimes acquired along the profile L1 in the Changbaishan-Tianchi volcanic region. It is used in this study that seismic traveltime inversion for simultaneous determination of 2-D velocity and interface structure of the crust and upper mantle. The result shows that, under Changbaishan-Tianchi crater, there exists a low-velocity body in the shape of an inverted triangle, and the crustal reflecting boundaries and Moho all become lower by a varying margin of 2-6 km, forming a crustal root which is assumed to be the Changbaishan-Tianchi volcanic system. Finally, we make a comparison between our 2-D velocity model and the result from the studies by using trial-and-error forward modeling with SEIS83.

Upper crustal structure under Jingtai–Hezuo profile in Northeastern Tibet from topography-dependent eikonal traveltime tomography

The Northeastern Tibetan plateau records Caledonian Qilian orogeny and Cenozoic reactivation by continental collision between the Indian and Asian plates. In order to provide the constraint on the Qilian orogenic mechanism and the expansion of the plateau,wide-angle seismic data was acquired along a 430 km-long profile between Jingtai and Hezuo. There is strong height variation along the profile,which is dealt by topography flattening scheme in our crustal velocity structure reconstruction. We herein present the upper crustal P-wave velocity structure model resulting from the interpretation of first arrival dataset from topography-dependent eikonal traveltime tomography. With topography flattening scheme to process real topography along the profile,the evenness of ray coverage times of the image area（upper crust）is improved,which provides upper crustal velocity model comparable to the classic traveltime tomography（with model expansion scheme to process irregular surface）. The upper crustal velocity model shows zoning character which matcheswith the tectonic division of the Qaidam-Kunlun-West Qinling belt,the Central and Northern Qilian,and the Alax blocks along the profile. The resultant upper crustal P-wave velocity model is expected to provide important base for linkage between the mapped surface geology and deep structure or geodynamics in Northeastern Tibet.

Teleseismic P-Wave Tomography of the New Guinea-Solomon Arc System

A P-wave tomographic traveltime inversion was applied to obtain a new model of seismic velocity anomalies beneath the New Guinea-Solomon arc system(PN-SL).The P-wave traveltime data,obtained from the revised International Seismological Center catalog,were recorded by 82 seismic stations in the PN-SL.Under the constraints of the epicenter distance,magnitude,and the number of stations recorded,15009 effective P-wave traveltime data were selected from 2011 teleseismic events.The obtained model showed that the Solomon Sea Plate subducted beneath the New Britain Island along the New Britain Trench at an angle of>70°and that the slab can be traced down to a depth of approximately 800 km.Conversely,we cannot observe a high-velocity anomaly exhibited by the subducted Solomon Sea Plate in the deep mantle at the Trobriand Trench,and the slab stopped at a depth of<200 km.The double subduction of the Solomon Sea Plate strongly modified the subduction patterns of the early subducted Pacific and Australian plates in the mantle along the West Melanesian Trench and the Pocklington Trough,respectively.In addition,the subducted Solomon Sea Plate induced the melting of the upper mantle to form a low-velocity anomaly,which provided the deep dynamic source for the expansion of the Bismarck Sea.Based on the joint consideration of the tomography results and a petrological analysis,the low-velocity anomalies beneath the Solomon Sea and Woodlark Basin are closely related to the early subduction of the Pacific and Australian Plates,respectively.

Application of teleseismic tomography to the study of shallow structure beneath Shizigou in the western Qaidam basin

Teleseismic body wave traveltime tomography is used to inverse the three-dimensional seismic velocity structure beneath Shizigou in the western Qaidam basin. The travel time are picked from the continuous observation data on a small seismic array of stations deployed during 2004-2007. The tomographic results obtained indicate that a NW-trending low velocity anomaly just beneath the target region insert northeastwards with a high dip angle. In the north, northeast and east of the low velocity anomaly, some high-velocity anomalies distribute with the same strike and coverage as those of Shizigou anticline.

Finite difference calculation of traveltime on non-orthogonal grid

Finite difference methods have been widely employed in solving the eikonal equation so as to calculate traveltime of seismic phase. Most previous studies used regular orthogonal grid. However, much denser grid is required to sample the interfaces that are undulating in depth direction, such as the Moho and the 660 km discontinuity.Here we propose a new finite difference algorithm to solve the eikonal equation on non-orthogonal grid（irregular grid）.To demonstrate its efficiency and accuracy, a test was conducted with a two-layer model. The test result suggests that the similar accuracy of a regular grid with ten times grids could achieve with our new algorithm, but the time cost is only about 0.1 times. A spherical earth model with an undulant660 km discontinuity was constructed to demonstrate the potential application of our new method. In that case, the traveltime curve fluctuation corresponds to topography. Our new algorithm is efficient in solving the first arrival times of waves associated with undulant interfaces.

Source-independent wave-equation based microseismic source location using traveltime inversion

A new source location method using wave-equation based traveltime inversion is proposed to locate microseismic events accurately. With a sourceindependent strategy, microseismic events can be located independently regardless of the accuracy of the source signature and the origin time. The traveltime-residuals-based misfit function has robust performance when the velocity model is inaccurate. The new Fréchet derivatives of the misfit function with respect to source location are derived directly based on the acoustic wave equation, accounting for the influence of geometrical perturbation and spatial velocity variation. Unlike the mostly used traveltime inversion methods, no traveltime picking or ray tracing is needed.Additionally, the improved scattering-integral method is applied to reduce the computational cost. Numerical tests show the validity of the proposed method.

Using Traveltime Tomography to Reconstruct Electrical Conductivity

A new method for reconstructing electrical conductivity distribution from electromagnetic (EM) data by using traveltime tomography is presented in this paper. Diffusive EM fields can be mathematically transformed to wavefields defined in a time like variable. The transform uniquely relates a field satisfying a diffusion equation in time, or in frequency, to an integral of the corresponding wavefield. This paper first transforms numerically calculated transient magnetic fields to wavefields. Traveltime data from a source to the receivers are estimated from the transformed wavefields. Then an iterative reconstruction algorithm is used to obtain the slowness distribution of a medium. This algorithm is an improved ART algorithm taking account of bending ray paths. The slowness distribution is transformed to electrical conductivity distribution according to their relation. The simulation result is presented in this paper.

A Neural Network for Weighted Least－Squares Criteria of Traveltime Tomography

TH network has been applied successfully to linear optimum and some quadratic optimum problems.This paper discusses how to determine the gain function of the net amplifiers and choose appropriate parameters to solve weighted least-squares problems.To test the performance of the net,it is used for seismic traveltime tomography in computer simulation. Two media of different contrast are taken in simulation.The experimental results show that if the parameters are determined appropriately, the performance of the network is good and the results are close to the ideal ones.

An Improved Iterative Reconstruction Algorithm for Traveltime Tomography

Traveltime tomography is a technique to reconstruct acoustic, seismic, or electromagnetic wave-speed distributions from first arrival traveltime data. The ray paths that should be used for tomographic techniques strongly depend on the wave-speed distribution. In this paper, a new method is proposed for finding out the ray paths from Fermat's principle, that means the traveltime of the ray path should be a minimum value. The problem of finding out the ray path is actually an optimum problem. Our new method uses the idea to find out the shortest path in a weighted directed graph to solve the problem. The ray paths found out by this method are used in the iterative reconstruction algorithm. Computer simulation result produced by this reconstruction algorithm is better than that by the conventional ones. It also shows that the new algorithm is effective with good convergency and stability.

mFAST:A MATLAB toolbox for ocean bottom seismometer refraction first-arrival traveltime tomography

First-arrival seismic traveltime tomography(FAST)is a well-established technique to estimate subsurface velocity structures.Although several existing open-source packages are available for first-arrival traveltime tomography,most were written in compiled languages and lack sufficient extendibility for new algorithms and functionalities.In this work,we develop an open-source,selfcontained FAST package based on MATLAB,one of the most popular interpreted scientific programming languages,with a focus on ocean bottom seismometer refraction traveltime tomography.Our package contains a complete traveltime tomography workflow,including ray-tracing-based first-arrival traveltime computation,linearized inversion,quality control,and high-quality visualization.We design the package as a modular toolbox,making it convenient to integrate new algorithms and functionalities as needed.At the current stage,our package is most efficient for performing FAST for two-dimensional ocean bottom seismometer surveys.We demonstrate the efficacy and accuracy of our package by using a synthetic data example based on a modified Marmousi model.

Application of Snell's law in reflection raytracing using the multistage fast marching method

Accurate calculations of travel times and raypaths of reflection waves are important for reflection travel time tomography.The multistage shortest path method(MSPM)and multistage fast marching method(MFMM)have been widely used in reflection wave raytracing,and both of them are characterized by high efficiency and ac-curacy.However,the MSPM does not strictly follow Snell's law at the interface because it treats the interface point as a sub-source,resulting in a decrease in accuracy.The MFMM achieves high accuracy by solving the Eikonal equation in local triangular mesh.However,the implementation process is complex.Here we propose a new method which uses linear interpolation to compute the incident travel time of interface points and then using Snell's law to compute the reflection travel time of grid points just above the interface.Our new method is much simpler than the MFMM;furthermore,numerical simulations show that the accuracy of the MFMM and our new method are basically the same,thus the reflection tomography algorithms which use our new method are easier to implement without decreasing accuracy.Besides,our new method can be extended easily to other grid-based raytracing methods.

Influences of coarse grid selection on Kirchhoff beam migration

Kirchhoff beam migration is a beam migration method, which focuses on rapid imaging of geological structures. Although this imaging method ignores the amplitude information in the calculation process, it can calculate multi-arrival traveltime. This migration method takes into account both imaging accuracy and computational efficiency. Kirchhoff beam migration employs coarse grid techniques in several key steps such as traveltime calculation, weight function calculation, and imaging calculation. The selection of the coarse mesh size has an important influence on the computational efficiency and imaging accuracy of the migration imaging method. This paper will analyze this influence and illustrate the analysis results by the Marmousi data sets.