*********************************** Program package ANRAY, version 4.00 *********************************** Dirk Gajewski *), Ivan Psencik **) *) Inst.of Geophysics, University of Hamburg, Bundesstrasse 55, 20146 Hamburg, FRG **) Geophysical Institute, Czechosl.Acad.Sci., Bocni II, 141 31 Praha 4, Czechoslovakia *************************************** Working version, not yet fully debugged *************************************** Main differences from the version 3.02 -------------------------------------- There are several new options offerred in this version. The main news are the following: (a) The distribution of elastic parameters in individual layers can be specified in grid points of a 3-D rectangular grid. Each of the 21 elastic parameters is specified in the grid individually. In case of an isotropic layer, only specification of P and S wave velocities is sufficient. If larger grid than specified in the present version is to be used, changes in the dimensions in the routines of the block MODBS.FOR must be made, see descritpion of input data for the B-spline approximation. (b) Graphical output in all the programs of the package is produced in the form of postscript files. Calcomp-GKS interface was developed by L.Klimes. (c) In order to increase the precision of the determination of points of intersection of rays and dynamic ray tracing with interfaces and boreholes, and thus to improve the convergence of the two-point ray tracing, procedures proposed in [4] are adopted in this package. (d) Tests, making possible a check of the precission of the ray tracing and dynamic ray tracing, are offered. The tests can be used both passively and actively, see the new switch IPREC in the description of input data for the program ANRAY, sub 10. Further investigation of the precission of the ray and dynamic ray tracing and its control is under way. (e) Many errors and problems have been removed and other small scale modifications have been made. Most errors were detected and removed in the routines for transformation of dynamic ray tracing across interfaces. Most of them were spotted by P.Bakker (Shell Int., Rijswijk). I am very grateful to Peter for it and I shall be grateful to anybody else for indication of other problems since with the new options, new bugs were surely introduced. Description of the package -------------------------- Program package ANRAY can be used for computation of rays, travel times, ray amplitudes and ray synthetic seismograms in 3-D laterally varying structures containing isotropic and/or anisotropic layers. Synthetic seismograms can be constructed at receivers distributed regularly or irregularly along surface, interfaces or vertical profiles. Program package consists of nine programs, basic program ANRAY and programs ANRAYPL, SYNTAN, SEISPLOT, POLARPLOT, FRESAN, SYNFAN, BPLOT and VELPL. Program ANRAY is an updated version of program with the same name which has been used in packages ANRAY86 and ANRAY89 written by Gajewski & Psencik [1],[2]. It is designed for ray, travel time and ray amplitude computations. Two ways of approximation of distribution of elastic parameters are available. In the first, "B-spline approximation", elastic parameters in an arbitrary point of a layer are determined by B-spline approximation from values of parameters specified at grid points of a 3-D rectangular network. B-splines with tension by A.K.Cline [3] with additions by L.Klimes, also used in [4], are used. In the second way, "isosurface interpolation", elastic parameters are determined by vertical linear interpolation of values specified at interfaces, which represent surfaces of constant values of elastic parameters. Rays can be computed in two modes. In the first one, rays are specified by the point source location and the initial orientation of the slowness vector at the source: initial-value ray tracing. In the second mode, rays are specified by the point source location and a system of regularly or irregularly distributed receivers situated on surface or on an interface or on a vertical profile: two-point ray tracing. The point source can be situated at any point of the model. Polarization vectors, geometrical spreading and reflection, transmission and conversion coefficients may be evaluated along the rays. Program ANRAY can produce three files, one for plotting ray diagrams, travel times and ray amplitudes, second for computation of synthetic seismograms and third for plotting plane sections of slowness, phase velocity or group velocity surfaces. Program ANRAY consists of 6 blocks of routines, ANRAY.FOR, A2.FOR, A3.FOR, A4.FOR, A5.FOR, CALCOMP.FOR and MODBS.FOR or MODIS.FOR, which must be linked together. If B-spline approximation is to be used, MODBS.FOR must be used, if isosurface interpolation is to be used, MODIS.FOR must be used. The program ANRAYPL can be used for plotting of horizontal and vertical ray diagrams, time-distance and amplitude-distance curves of individual elementary waves computed in the program ANRAY. Several types of point sources, namely explosive source, single force, double couple source can be considered. Program ANRAYPL is a modification of the program RAYPLOT from the package SEIS83 [5]. Program ANRAYPL consists of routines: ANRAYPL.FOR, BORDER.FOR, SERV.FOR, SOURCE.FOR and CALCOMP.FOR. The program SYNTAN can be used for computation of ray synthetic seismograms from the results computed in the program ANRAY. Gabor wavelet is used as a source-time function. Hilbert transform of this wavelet is evaluated by approximate formula, see [6]. Several types of point sources, namely explosive source, single force, double couple source, can be considered. Program SYNTAN is a modification of the program SYNTPL from the program package SEIS83 [5]. Program SYNTAN consists of routines: SYNTAN.FOR, SERV.FOR, SOURCE.FOR and CALCOMP.FOR. Programs FRESAN and SYNFAN can be used for frequency domain computation of ray synthetic seismograms. Frequency response is computed in the program FRESAN. In the program SYNFAN it is multiplied by the spectrum of the considered source-time function. Inverse Fourier transform is then used to evaluate synthetic seismograms. Programs FRESAN and SYNTAN are modifications of the programs GB and SYNTGB from the package BEAM87 by V.Cerveny. Program FRESAN consists of routines FRESAN.FOR, BORDER.FOR, SOURCE.FOR, SERV.FOR and CALCOMP.FOR. Program SYNFAN consists of routines SYNFAN.FOR, SERV.FOR and CALCOMP.FOR. Programs SEISPLOT, BPLOT and POLARPLOT can be used for plotting results computed in progams SYNTAN and SYNFAN. Program SEISPLOT can be used for plotting synthetic seismograms generated by the program SYNTAN, program BPLOT can be used for plotting synthetic seismograms generated by the program SYNFAN. The program POLARPLOT serves for plotting particle motion diagrams and can be used in both cases. Programs SEISPLOT and POLARPLOT are modifications of the programs with the same name from the package SEIS83, program BPLOT is a modification of the program with the same name from the package BEAM87. Program SEISPLOT consists of routines SEISPLOT.FOR, BORDER.FOR, SERV.FOR and CALCOMP.FOR. Program BPLOT consists of routines BPLOT.FOR, BORDER.FOR, SERV.FOR and CALCOMP.FOR. Program POLARPLOT consists of routines POLARPLOT.FOR, BORDER.FOR, SERV.FOR and CALCOMP.FOR. Program VELPLOT can be used for plotting plane sections of slowness, phase velocity and group velocity surfaces from the file generated in the program ANRAY. Program VELPLOT consists of three routines: VELPL.FOR, SERV.FOR and CALCOMP.FOR. FORTRAN77 is used throughout the program package ANRAY (with exception of the routine SERV). Standard CALCOMP routines, namely PLOTS, PLOT, SYMBOL and NUMBER, are used in plotting programs. Calcomp-GKS interface containing the above listed routines and developed by L.Klimes, is a part of the package. This makes possible to generate the graphical output in the form of ASCII files coded in the postscript level 2 language - encapsulated postscript file format version 3.0. The Calcomp-GKS interface is a part of the block of routines CALCOMP.FOR. This block must be linked with each program of the package. The package is supplemented by a complete set of test input data for an "anisotropic" model, in which isosurface approximation is used. At present version the test input data for a model, in which B-splines are used, are given only for the program ANRAY. A simple three-layered model consisting of two homogeneous isotropic layers bounding from above and below an anisotropic layer with horizontal variation of elastic parameters in x- and y-direction is considered. At the receivers situated in a borehole both in isotropic and anisotropic layer all direct waves and one reflected qS wave are calculated. The program package is NOT yet fully debugged. Problems can be encountered especially in singular regions in which two slowness surfaces of the shear waves are close to each other. The ray method gives in such regions distorted results and the two-point ray tracing procedure can fail there. Introduction of criteria checking precission of the ray tracing and dynamic ray tracing indicates that present formulation of these procedures is very sensitive to the above mentioned singular regions. Alternative formulations have to be, therefore, studied. Search for alternative ways of solving the polynomial equation of the sixth order, necessary for the determination of slowness vectors of waves generated by reflection or transmission at an interface, is also desirable. Any report on other errors, problems or inconsistencies, which may be encountered by users, is welcome. In case of publication of results obtained with the program package ANRAY, the users are asked kindly to refer to this report. Good luck! Praha, September 1996 REFERENCES [1] Gajewski,D., Psencik,I., 1986. Numerical modelling of seismic wave fields in 3-D laterally varying layered anisotropic structures - Program ANRAY86, Internal Report Inst.of Earth and Planet.Phys., University of Alberta, Edmonton. [2] Gajewski,D., Psencik,I.,1989. Ray synthetic seismograms in 3-D laterally inhomogeneous anisotropic structures - Program ANRAY89, Internal Report Centre for Computational Seismology, LBL, Berkeley. [3] Cline,A.K.,1981. FITPACK, Dept. of Computer Sciences, Univ. of Texas at Austin. [4] Cerveny,V., Klimes,L., Psencik,I.,1988. Complete seismic ray tracing in three-dimensional structures. In: Seismological Algorithms, D.J.Doornbos edit., Academic Press, New York, 89-168. [5] Cerveny,V., Psencik,I.,1984. Documentation of earthquake algorithms. SEIS83 - Numerical modeling of seismic wave fields in 2-D laterally varying layered structures by the ray method. E.R.Engdahl, edit., Report SE-35, Boulder, 36-40. [6] Cerveny,V., Molotkov,I.A., Psencik,I.,1977. Ray method in seismology, Univerzita Karlova, Praha.