First Arrival Seismic Tomography (FAST) vs. PStomo_eq applied to crooked line seismic data from the Siljan ring area
Introduction
Seismic traveltime tomography is presently one of the most used methods for calculating velocity fields. Several computer programs are available for both 2D and 3D datasets (e.g. Rayinvr, tomo2D, TTT, Jive3D, among others). Seismic data acquired along crooked lines have an inherent 3D component resulting in erroneous solutions from 2D traveltime inversion. The source of the errors is both physical and technical; Zelt and Zelt (1998) showed that there is no 2D model consistent with the data in the case of crooked line geometry and out of the plane velocity variations. This inaccuracy will result in the inverted velocity field having erroneous small-scale structure for 2D and even for 2.5D inversion, i.e. 3D ray tracing on a 2D model (Zelt and Zelt, 1998). Errors also result when 2D programs are used to project the source and receiver locations from the crooked line to a straight line. To illustrate this point we took a 1D constant gradient velocity field and calculated analytically (Shearer, 1999) the traveltimes for seven shots distributed on 2D and crooked line geometries (Fig. 1). We compared the theoretical traveltimes to those calculated by a 2D commercial program for seismic processing (Fig. 2). This program projects the source and receiver locations to a straight line based on CDP (Common Deep Point) numbers. Note the obvious problems encountered by the program when the crooked line geometry is used (Fig. 2B). Profiles that differ from straight lines are common in seismic experiments, mainly due to logistical reasons. Thus, there is a need for algorithms that can handle crooked lines properly.
Furthermore, it is well known that reflection seismic sections can be greatly degraded by shallow structure e.g. weathered layers with variable thickness and very low velocity. Determining a detailed near-surface velocity field is particularly useful for advanced processing, e.g. wave equation datuming and pre-stack migration, to correct for these effects. Taking into account the features due to crooked line geometry will be important for calculating the velocity field.
In this paper we compare results from two publicly available programs that handle crooked line geometry. These algorithms are First Arrival Seismic Tomography (FAST) (Zelt and Barton, 1998) and PStomo_eq (Benz et al., 1996; Tryggvason et al., 2002). The programs use basically the same forward and inverse algorithms, but differ in how smoothing is applied. Both programs have been used to study the upper crustal velocity structure in crystalline rock areas (i.e. Marti et al., 2002; Zelt et al., 2001). We compare the accuracy of first arrival time calculations, inversion of synthetic data and, finally, inversion of real data. In particular, we show that the difference in how smoothing is applied is critical in determining the final velocity model.
Section snippets
Forward modelling
To propagate the traveltime field through a 3D gridded velocity model, both programs use a finite difference approximation to the eikonal equation of ray tracingwhere t is traveltime and s is the slowness (inverse of velocity) of the media. Traveltimes are calculated progressively away from a source on the sides of an expanding cube, completing one side at a time. From the traveltime field, the ray paths may be calculated.
A number of studies (e.g. Vidale, 1988,
Inversion method
Both FAST and PStomo_eq solve the inverse problem iteratively by LSQR, which is a conjugate gradient method, based on the Lanczos bi-diagonalization process (Paige and Saunders, 1982). Regularization is a method to solve mixed-determined problems by applying constraints on the model, in addition to using the data. These constraints require that the final model satisfy some property or conditions. Commonly, we look for models containing “minimum structures”, i.e. the simplest model that fits the
Forward modelling
FAST uses the HZ improvement of Vidale's scheme for forward modelling. It also includes a straight ray calculation in a small cube around the source. This allows the user to assign velocities for nodes above and below the shot inside a 5×5×5 node box around the shot. The user can also specify (at run time) the maximum number of reverse propagations and the side of the cube from which they would start if head wave operators are used in the primary propagation. The user provides the cell size at
Practical and technical features
Table 1 presents the CPU times for each program on a Sun UltraSparc workstation. Times are per non-linear iteration for the tests presented here, i.e., a forward model parameterized using 617 397 nodes and an inversion model parameterized using 151 008 cells. Times shown for FAST are for one-lambda test. In practice, the typical number of lambda values to test is between 4 and 8 (FAST documentation).
The uncompressed package of FAST has a size of about 20 MB. All programs are written in Fortran 77
The Siljan dataset
The real dataset used for this study was acquired in the Siljan Ring impact area, in central Sweden (Juhlin and Pedersen, 1987). It consists of first break traveltimes from two of the deep seismic profiles in the area. These profiles have a general north–south direction and cross lithological boundaries, fracture systems and dike intrusions (Fig. 3). The more northern profile was shot as “end-on” whereas the more southern profile was shot as “split-spread”. There is also a difference in the
Forward calculation of the traveltimes
To test the accuracy of the forward traveltime calculation method, we calculated analytically (Shearer, 1999) the traveltimes on a 1D constant gradient (1/s) velocity field with values of 4500 m/s at the top and 5300 m/s at the bottom of the model (Table 2). The analytical values, from now on referred to as theoretical, were compared to those calculated by FAST and PStomo_eq. Both programs were tested using the HZ algorithm for the forward traveltime calculation; a difference is that in FAST it
Discussion and conclusions
As stated earlier, the smoothness constraint in both programs is handled differently. In FAST, the Laplacian operator is normalized by the prior slowness of the centre cell. FAST is programmed to take the smoothest model with the better-normalized RMS misfit, avoiding overfitting the data. Residual values are normalized by the accuracy of the traveltime picks. Consequently, a good estimation of uncertainties is required. Use of values that are too high results in models lacking information that
Acknowledgments
We are truly thankful to Colin Zelt for his enthusiastic help and quick answers. Ari Tryggvasson also helped us considerably and provided valuable input to this work.
References (23)
Interpretation of reflections in the Siljan Ring area based on results from the Gravberg-1 borehole
Tectonophysics
(1990)- et al.
Study of out-of-plane effects in the inversion of refraction/wide-angle reflection traveltimes
Tectonophysics
(1998) - et al.
Three-dimensional P- and S-wave velocity structure of Redoubt Volcano, Alaska
Journal of Geophysical Research
(1996) - et al.
High Resolution seismic traveltime tomography incorporating static corrections applied in a till covered bedrock environment
Geophysics
(2004) - et al.
3-D finite-difference reflection traveltimes
Geophysical Journal International
(1995) - et al.
Adaptation of back projection tomography to seismic travel time problems
Journal of Geophysical Research
(1988) - et al.
Reflection seismic investigations of the Siljan impact structure
Journal of Geophysical Research
(1987) - Juhlin, C., Palm, H., Bergman, B., 2001. Reflection Seismic Imaging of the Upper Crystalline Crust for Characterization...
- et al.
Mapping brittle fracture zones in three dimensions: high-resolution traveltime seismic tomography in a granitic pluton
Geophysical Journal International
(2002) Geophysical Data Analysis: Discrete Inverse Theory
(1984)
LSQR: an algorithm for sparse linear equations and sparse least squares
ACM Transactions on Mathematical Software
Cited by (3)
Improving horizontal resolution of high-frequency surface-wave methods using travel-time tomography
2016, Journal of Applied GeophysicsCitation Excerpt :The multi-offset phase analysis (MOPA) has been used to detect the lateral S-wave velocity variations by the maximum number of adjacent traces with no detriment to the lateral resolution, however, it demands knowing the exact location of the anomalous body (Vignoli et al., 2011). In recent years, travel-time tomography has become a very popular method in revealing earth structures from regional to global scales (Rodríguez-Tablante et al., 2006; Saygin and Kennett, 2012; Tondi et al., 2012; Bou'e et al., 2014). Considering the high-frequency approximation, surface waves propagate along a straight ray path (Um and Thurber, 1987).
Seismic first-arrival tomography with functional description of traveltimes
2009, Journal of Geophysics and Engineering