Showing papers in "Geophysics in 1987"

Occam's inversion; a practical algorithm for generating smooth models from electromagnetic sounding data

Abstract: The inversion of electromagnetic sounding data does not yield a unique solution, but inevitably a single model to interpret the observations is sought. We recommend that this model be as simple, or smooth, as possible, in order to reduce the temptation to overinterpret the data and to eliminate arbitrary discontinuities in simple layered models.To obtain smooth models, the nonlinear forward problem is linearized about a starting model in the usual way, but it is then solved explicitly for the desired model rather than for a model correction. By parameterizing the model in terms of its first or second derivative with depth, the minimum norm solution yields the smoothest possible model.Rather than fitting the experimental data as well as possible (which maximizes the roughness of the model), the smoothest model which fits the data to within an expected tolerance is sought. A practical scheme is developed which optimizes the step size at each iteration and retains the computational efficiency of layered models, resulting in a stable and rapidly convergent algorithm. The inversion of both magnetotelluric and Schlumberger sounding field data, and a joint magnetotelluric-resistivity inversion, demonstrate the method and show it to have practical application.

Nonlinear two-dimensional elastic inversion of multioffset seismic data

Abstract: The treatment of multioffset seismic data as an acoustic wave field is becoming increasingly disturbing to many geophysicists who see a multitude of wave phenomena, such as amplitude-offset variations and shear-wave events, which can only be explained by using the more correct elastic wave equation. Not only are such phenomena ignored by acoustic theory, but they are also treated as undesirable noise when they should be used to provide extra information, such as S-wave velocity, about the subsurface. The problems of using the conventional acoustic wave equation approach can be eliminated via an elastic approach. In this paper, equations have been derived to perform an inversion for P-wave velocity, S-wave velocity, and density as well as the P-wave impedance, S-wave impedance, and density. These are better resolved than the Lame parameters. The inversion is based on nonlinear least squares and proceeds by iteratively updating the earth parameters until a good fit is achieved between the observed data and the modeled data corresponding to these earth parameters. The iterations are based on the preconditioned conjugate gradient algorithm. The fundamental requirement of such a least-squares algorithm is the gradient direction which tells how to update the model parameters. The gradient direction can be derived directly from the wave equation and it may be computed by several wave propagations. Although in principle any scheme could be chosen to perform the wave propagations, the elastic finite-difference method is used because it directly simulates the elastic wave equation and can handle complex, and thus realistic, distributions of elastic parameters. This method of inversion is costly since it is similar to an iterative prestack shot-profile migration. However, it has greater power than any migration since it solves for the P-wave velocity, S-wave velocity, and density and can handle very general situations including transmission problems. Three main weaknesses of this technique are that it requires fairly accurate a priori knowledge of the low-wavenumber velocity model, it assumes Gaussian model statistics, and it is very computer-intensive. All these problems seem surmountable. The low-wavenumber information can be obtained either by a prior tomographic step, by the conventional normal-moveout method, by a priori knowledge and empirical relationships, or by adding an additional inversion step for low wavenumbers to each iteration. The Gaussian statistics can be altered by preconditioning the gradient direction, perhaps to make the solution blocky in appearance like well logs, or by using large model variances in the inversion to reduce the effect of the Gaussian model constraints. Moreover, with some improvements to the algorithm and more parallel computers, it is hoped the technique will soon become routinely feasible.

On the imaging of reflectors in the earth

Abstract: In this paper, I present a modification of the Beylkin inversion operator. This modification accounts for the band-limited nature of the data and makes the role of discontinuities in the sound speed more precise. The inversion presented here partially dispenses with the small-parameter constraint of the Born approximation. This is shown by applying the proposed inversion operator to upward scattered data represented by the Kirchhoff approximation, using the angularly dependent geometrical-optics reflection coefficient. A fully nonlinear estimate of the jump in sound speed may be extracted from the output of this algorithm interpreted in the context of these Kirchhoff-approximate data for the forward problem. The inversion of these data involves integration over the source-receiver surface, the reflecting surface, and frequency. The spatial integrals are computed by the method of stationary phase. The output is asymptotically a scaled singular function of the reflecting surface. The singular function of a surface is a Dirac delta function whose support is on the surface. Thus, knowledge of the singular functions is equivalent to mathematical imaging of the reflector. The scale factor multiplying the singular function is proportional to the geometrical-optics reflection coefficient. In addition to its dependence on the variations in sound speed, this reflection coefficient depends on an opening angle between rays from a source and receiver pair to the reflector. I show how to determine this unknown angle. With the angle determined, the reflection coefficient contains only the sound speed below the reflector as an unknown, and it can be determined. A recursive application of the inversion formalism is possible. That is, starting from the upper surface, each time a major reflector is imaged, the background sound speed is updated to account for the new information and data are processed deeper into the section until a new major reflector is imaged. Hence, the present inversion formalism lends itself to this type of recursive implementation. The inversion proposed here takes the form of a Kirchhoff migration of filtered data traces, with the space-domain amplitude and frequency-domain filter deduced from the inversion theory. Thus, one could view this type of inversion and parameter estimation as a Kirchhoff migration with careful attention to amplitude.

Incorporation of attenuation into time-domain computations of seismic wave fields

Abstract: The only numerically tractable way yet found to incorporate attenuation into numerical time-domain computations of seismic wave fields is to approximate the viscoelastic modulus by a low-order rational function of frequency. The coefficients of this function can be determined by the Pade approximation. Our test computations show, however, that this approximation generally is of poor quality. Therefore, we suggest a new approach which is based on the rheological model of the generalized Maxwell body, which has a modulus of the desired rational form. We choose the relaxation frequencies logarithmically equidistant in the frequency band of interest, and determine the weight factors by simple numerical curve fitting to an arbitrary Q law. This approach is superior to the method above both in accuracy and in computational efficiency. For most practical applications, approximations of orders 2 or 3 are sufficient. The computing time and memory requirements for a finite-difference calculation are then approximately twice those of a purely elastic calculation. As a first application of the method, we compute SH channel waves in discontinuous coal seams with Q = 50 within the coal. The results show that the high-frequency Airy phase is strongly attenuated. This indicates that care has to be taken in comparing the results of purely elastic model calculations of the propagation of seam waves with experimental data.

Diffraction tomography and multisource holography applied to seismic imaging

Abstract: Seismic tomography is emerging as an imaging method for determining subsurface structure. When the view-angle coverage is limited and the scale of the medium inhomogeneities is comparable with the wavelength, as is often true in geophysical applications, the performance of ordinary ray tomography becomes poor. Other tomographic methods are needed to improve the imaging process. Here we study diffraction tomography and multisource holography and evaluate their performances for surface reflection profiling (SRP), vertical seismic profiling (VSP), and cross-hole measurements. Theoretical formulations are derived for two-dimensional geometry in terms of line sources along a source line and line receivers along a receiver line. The theory for diffraction tomography is based on the Born or Rytov approximation.The performances of diffraction tomography and multisource holography are evaluated by examining the information coverage in the spatial frequency domain and by numerical examples. Multisource holography, which is similar to Kirchhoff-type migration, often gives distorted images of the object. This distortion causes long tails of the image in the case of SRP and a strong noise belt in the case of VSP and is due to incomplete and nonuniform coverage of the object spectrum. The filtering operation of diffraction tomography helps in correcting the nonuniform coverage (including duplication) of the object spectrum in the reconstruction process and therefore reduces the distortions. On the other hand, multisource holography is better suited for imaging sharp boundaries with large acoustic impedance contrasts since diffraction tomography is restricted, as presently formulated, to weak inhomogeneities. In addition, multisource holography has the flexibility to be used with an arbitrary number of sources (including a single source). Its sampling interval is not restricted by the Nyquist frequency.Numerical examples show that combined data sets (such as surface reflection data combined with VSP data, or cross-hole data combined with surface data, etc.) improve the image quality.

A new slant on seismic imaging: Migration and integral geometry

Abstract: A new approach to seismic migration formalizes the classical diffraction (or common-tangent) stack by relating it to linearized seismic inversion and the generalized Radon transform. This approach recasts migration as the problem of reconstructing the earth’s acoustic scattering potential from its integrals over isochron surfaces. The theory rests on a solution of the wave equation with the geometrical-optics Green function and an approximate inversion formula for the generalized Radon transform. The method can handle both complex velocity models and (nearly) arbitrary configurations of sources and receivers. In this general case, the method can be implemented as a weighted diffraction stack, with the weights determined by tracing rays from image points to the experiment’s sources and receivers. When tested on a finite-difference simulation of a deviated-well vertical seismic profile (a hybrid experiment which is difficult to treat with conventional wave-equation methods), the algorithm accurately reconstructed faulted-earth models. Analytical reconstruction formulas are derived from the general formula for zero-offset and fixed-offset surface experiments in which the background velocity is constant. The zero-offset inversion formula resembles standard Kirchhoff migration. Our analysis provides a direct connection between the experimental setup (source and receiver positions, source wavelet, background velocity) and the spatial resolution of the reconstruction. Synthetic examples illustrate that the lateral resolution in seismic images is described well by the theory and is improved greatly by combining surface data and borehole data. The best resolution is obtained from a zero-offset experiment that surrounds the region to be imaged.

The pseudospectral method: Comparisons with finite differences for the elastic wave equation

Abstract: The pseudospectral (or Fourier) method has been used recently by several investigators for forward seismic modeling. The method is introduced here in two different ways: as a limit of finite differences of increasing orders, and by trigonometric interpolation. An argument based on spectral analysis of a model equation shows that the pseudospectral method (for the accuracies and integration times typical of forward elastic seismic modeling) may require, in each space dimension, as little as a quarter the number of grid points compared to a fourth‐order finite‐difference scheme and one‐sixteenth the number of points as a second‐order finite‐difference scheme. For the total number of points in two dimensions, these factors become 1/16 and 1/256, respectively; in three dimensions, they become 1/64 and 1/4 096, respectively. In a series of test calculations on the two‐dimensional elastic wave equation, only minor degradations are found in cases with variable coefficients and discontinuous interfaces.

A case for upward continuation as a standard separation filter for potential-field maps

Abstract: Separation filtering is incomplete even under the ideal synthetic condition of known power spectra of the regional and residual fields. I have designed some Wiener filters, which minimize the inevitable separation error, from previous statistical source models of Naidu, and Spector and Grant. This formulation includes the classic separation filters of Strakhov and of Elkins as Wiener filters. A proposed generalization of Wiener filters, denoted as uniformly suboptimum filters, quantitatively supports the statement that a wide span of separation problems may be solved adequately using some convenient, small standard filter family. A uniform random‐source model without assumed vertical correlations invokes upward continuation filters. In addition to this role as a Wiener filter, the upward continuation operator is given by elementary functions in both space and wavenumber domains, is numerically stable, and is also physically comprehensible when applied to real, nonrandom anomalies. In view of these disting...

Computing the gravitational and magnetic anomalies due to a polygon; algorithms and Fortran subroutines

Abstract: We present two algorithms for computing the gravitational and magnetic anomalies due to an n‐sided polygon in a two‐dimensional space. Both algorithms have been implemented as subroutines coded in Fortran-77, and listings are provided. Because references to trigonometric functions have been almost completely eliminated, these codes run substantially faster than mosts codes now in existence. Furthermore, anomalies can be computed at any point outside, on, or inside the polygon. Unlike other codes, these algorithms can be used to model subsurface observations.

The dielectric constant of sandstones, 60 kHz to 4 MHz

Abstract: Complex impedance data were collected for eight sandstones at various levels of water saturation (Sw) in the frequency range of 5 Hz to 4 MHz. The measurements were made using a two‐electrode technique with platinum electrodes sputtered onto the flat faces of disk‐ shaped samples. Presentation of the data in the complex impedance plane shows clear separation of the response due to polarization at the sample‐electrode interface from the bulk sample response. Electrode polarization effects were limited to frequencies of less than 60 kHz, allowing us to study the dielectric constant κ′ of the sandstones in the frequency range of 60 kHz to 4 MHz. κ′ of all samples at all levels of saturation shows a clear power‐law dependence upon frequency. Comparing the data from the eight sandstones at Sw=0.36, the magnitude of the frequency dependence was found to be proportional to the surface area‐to‐volume ratio of the pore space of the sandstones. The surface area‐to‐volume ratio of the pore space of each sandstone wa...

Limitations of determining density or magnetic boundaries from the horizontal gradient of gravity or pseudogravity data

Abstract: The horizontal‐gradient method has been used since 1982 to locate density or magnetic boundaries from gravity data (Cordell, 1979) or pseudogravity data (Cordell and Grauch, 1985). The method is based on the principle that a near‐vertical, fault‐like boundary produces a gravity anomaly whose horizontal gradient is largest directly over the top edge of the boundary. Magnetic data can be transformed to pseudogravity data using Fourier techniques (e.g., Hildenbrand, 1983) so that they behave like gravity data; thus the horizontal gradient of pseudogravity also has maximum magnitude directly over the boundary. The method normally is applied to gridded data rather than to profiles. The horizontal‐gradient magnitude is contoured and lines are drawn or calculated (Blakely and Simpson, 1986) along the contour ridges. These lines presumably mark the top edges of magnetic or density boundaries. However, horizontal‐gradient magnitude maxima (gradient maxima) can be offset from a position directly over the boundary f...

Elastic reverse-time migration

Abstract: Elastic, prestack, reverse-time, finite-difference migration of two-component seismic surface data requires data extrapolation and application of an imaging condition. Data extrapolation involves synchronous driving of the vertical-component and horizontal-component finite-difference meshes with the time reverse of the recorded vertical and horizontal traces, respectively. Extrapolation uses the coupled elastic wave equation for variable velocity solved with a second-order, explicit finite-difference scheme. The imaging condition at any point in the grid is the one-way traveltime from the source to that point.Elastic migrations of both synthetic test data and real two-component common-source gathers produce simpler images than acoustic migrations because of the coalescing of double reflections (compressional waves and shear waves) into single loci.

Tomographic inversion via the conjugate gradient method

Abstract: Tomographic inversion of seismic traveltime residuals is now an established and widely used technique for imaging the Earth’s interior. This inversion procedure results in large, but sparse, rectangular systems of linear algebraic equations; in practice there may be tens or even hundreds of thousands of simultaneous equations. This paper applies the classic conjugate gradient algorithm of Hestenes and Stiefel to the least‐squares solution of large, sparse systems of traveltime equations. The conjugate gradient method is fast, accurate, and easily adapted to take advantage of the sparsity of the matrix. The techniques necessary for manipulating sparse matrices are outlined in the Appendix. In addition, the results of the conjugate gradient algorithm are compared to results from two of the more widely used tomographic inversion algorithms.

Wave‐equation trace interpolation

Abstract: Spatial aliasing in multichannel seismic data can be overcome by solving an inversion in which the model is the section that would be recorded in a well sampled zero-offset experiment, and the data are seismic data after normal moveout (NMO). The formulation of the (linear) relation between the data and the model is based on the wave equation and on Fourier analysis of aliasing. A processing sequence in which one treats missing data as zero data and performs partial migration before stacking is equivalent to application of the transpose of the operator that actually needs to be inverted. The inverse of that operator cannot be uniquely determined, but it can be estimated using spatial spectral balancing in a conjugate-gradient iterative scheme. The first iteration is conventional processing (including prestack partial migration). As shown in a field data example in which severe spatial aliasing was simulated, a few more iterations are necessary to achieve significantly better results.

Practical aspects of reflectivity modeling

Abstract: This paper is intended to help those not familiar with the “lore” of layered earth modeling to avoid some common problems. In the computation of the reflectivity function, an easily incorporated phase‐integral approximation is used away from turning points when the velocity gradient is smaller than the frequency. Hanning windows, or segments thereof, work well for both the slowness integral and the frequency integral. For the quadrature of the slowness integral the Filon method of Frazer is easily coded and vectorizes well; Levin’s Filon method and the Clenshaw‐Curtis‐Filon method of Xu and Mal are more difficult to vectorize, but more powerful because they require fewer evaluations of the reflectivity function. A modification of Strick’s power law is a convenient way to calculate complex frequency‐dependent seismic velocities. The complex frequency technique for avoiding time aliasing is explained by use of the Poisson sum formula. In writing code for vector computers, such as the CRAY, if frequency‐inde...

Absorbing boundary conditions and surface waves

Abstract: One of the major problems in numerically simulating waves traveling in the Earth is that an artificial boundary must be introduced to produce unique solutions. To eliminate the spurious reflections introduced by this artificial boundary, we use a damping expression based on analogies to shock absorbers. This method can reduce the amplitude of the reflected wave to any pre‐specified value and is successful for waves at any angle of incidence. The method can eliminate unwanted reflections from the surface, reflections at the corners of the model, and waves reflected off an interface that strike the artificial boundary. Many of the boundary conditions currently used in the numerical solution of waves are approximations to perfectly absorbing boundary conditions and depend upon the angle of incidence of the incoming wave at the artificial boundary. Stability problems often occur with these boundary conditions. The method we use at the artificial boundary allows use of stable Dirichlet or von Neumann condition...

Magnetotelluric strike rules

Abstract: The distortion of the magnetotelluric impedance tensor by complex “near‐surface” structure leads to leakage between the elements of the tensor. The magnetotelluric impedance tensor for our principal model, which has both a local and a regional strike, can be written in the long‐period limit as a sum of the regional, undistorted impedance and a perturbed impedance. The latter can be written as a product of a local distortion (which can be regarded as thin‐sheet distortion in the long‐period range) tensor and the regional impedance. Local and regional strikes are found by rotating the impedance tensor into directions in which diagonal elements are proportional and column elements are proportional, respectively. The regional impedance tensor is calculated assuming that the strikes are uniquely defined. An example from a crystalline area with well conducting fracture zones illustrates the model concepts. A weighted least‐squares procedure is used for the estimation of distortion parameters.

On the theory of sea-floor conductivity mapping using transient electromagnetic systems

Abstract: The electrical conductivity of the sea floor is usually much less than that of the seawater above it. A theoretical study of the transient step-on responses of some common controlled-source, electromagnetic systems to adjoining conductive half-spaces shows that two systems, the horizontal, in-line, electric dipole-dipole and horizontal, coaxial, magnetic dipole-dipole, are capable of accurately measuring the relatively low conductivity of the sea floor in the presence of seawater. For these systems, the position in time of the initial transient is indicative of the conductivity of the sea floor, while at distinctly later times, a second characteristic of the transient is a measure of the seawater conductivity. The diagnostic separation in time between the two parts of the transient response does not occur for many other systems, including several systems commonly used for exploration on land. A change in the conductivity of the sea floor produces a minor perturbation in what is essentially a seawater response. Some transient responses which could be observed with a practical, deep-towed coaxial magnetic dipole-dipole system located near the sea floor are those for the half-space, the layer over a conductive or resistive basement, and the half-space with an intermediate resistive zone. The system response to two adjoining half-spaces, representing seawater and sea floor, respectively, is derived analytically. The solution is valid for all time, provided the conductivity ratio is greater than about ten, or less than about one-tenth. The analytic theory confirms the validity of numerical evaluations of closed-form solutions to these layered-earth models. A lateral conductor such as a vertical, infinite, conductive dike outcropping at the sea floor delays the arrival of the initial crustal transient response. The delay varies linearly with the conductance of the dike. This suggests that time delay could be inverted directly to give a measure of the anomalous integrated conductance of the sea floor both between and in the vicinity of the transmitter and the receiver dipoles.

An effective anisotropy parameter in transversely isotropic media

Abstract: An interesting physical meaning is presented for the anisotropy parameter δ, previously introduced by Thomsen to describe weak anisotropy in transversely isotropic media. Roughly, δ is the difference between the P-wave and SV-wave anisotropies of the medium. The observed systematic depth errors in the North Sea are reexamined in view of the new interpretation of the moveout velocity through δ. The changes in δ at an interface adequately describe the effects of transverse isotropy on the P-wave reflection amplitude, The reflection coefficient expression is linearized in terms of changes in elastic parameters. The linearized expression clearly shows that it is the variation of δ at the interface that gives the anisotropic effects at small incidence angles. Thus, δ effectively describes both the moveout velocity and the reflection amplitude variation, two very important pieces of information in reflection seismic prospecting, in the presence of transverse isotropy.

Estimating slowness dispersion from arrays of sonic logging waveforms

Abstract: Acoustic wave propagation in a fluid‐filled borehole is affected by the type of rock which surrounds the hole. More specifically, the slowness dispersion of the various body‐wave and borehole modes depends to some extent on the properties of the rock. We have developed a technique for estimating the dispersion relations from data acquired by full‐waveform digital sonic array well‐logging tools. The technique is an extension of earlier work and is based on a variation of the well‐known Prony method of exponential modeling to estimate the spatial wavenumbers at each temporal frequency. This variation, known as the forward‐backward method of linear prediction, models the spatial propagation by purely real‐valued wavenumbers. The Prony exponential model is derived from the physics of borehole acoustics under the assumption that the formation does not vary in the axial or azimuthal dimensions across the aperture of the receiver array, but can vary arbitrarily in the radial dimension. The exponential model fits...

Three‐dimensional seismic monitoring of an enhanced oil recovery process

Abstract: Seismic reflection data were used to monitor the progress of an in‐situ combustion, enhanced oil recovery process. Three sets of three‐dimensional (3-D) data were collected during a one‐year period in order to map the extent and directions of propagation in time. Acquisition and processing parameters were identical for each survey so that direct one‐to‐one comparison of traces could be made. Seismic attributes were calculated for each common‐depth‐point data set, and in a unique application of seismic reflection data, the preburn attributes were subtracted from the midburn and postburn attributes. The resulting “difference volumes” of 3-D seismic data showed anomalies which were the basis for the interpretation shown in this case study. Profiles and horizon slices from the data sets clearly show the initiation and development of a bright spot in the reflection from the top of the reservoir and a dim spot in the reflection from a limestone below it. Interpretation of these anomalies is supported by informa...

Imaging quasi-layered conductive structures by simple processing of transient electromagnetic data

Abstract: An “imaged” conductivity section of a layered earth can be obtained by simple transformation of step‐response electromagnetic data measured in the quasi‐static zone. This method of data transformation is presented as an alternative to conventional apparent conductivity transformations. At each delay time, the variation of the step response as a function of geometry (transmitter and receiver location) is transformed to an equivalent reference depth h, which can be related to the depth of electromagnetic field diffusion. The behavior of h as a function of delay time is nearly independent of the source‐receiver geometry. The slowness dt/dh divided by the magnetic permeability is almost exactly proportional to the cumulative conductance measured from the surface down to a depth h. Thus we can estimate an apparent conductivity, which we call the “imaged conductivity,” at depth to be d2t/μ0dh2. The cost of this transformation is a fraction of the cost of conventional data inversion, and it does not require an a...

Surface emittance, temperature, and thermal inertia derived from Thermal Infrared Multispectral Scanner (TIMS) data for Death Valley, California

Abstract: The NASA airborne Thermal Infrared Multispectral Scanner (TIMS) was flown over Death Valley, California, on both a daytime flight and a nighttime flight within a two‐day period in July 1983. This Daedulus scanner has six channels in the thermal infrared, between 8 and 12 μm. Calibrated digital spectral radiance data from these flights, along with Landsat Thematic Mapper (TM) reflectance data, permit the calculation of both spectral emittance and thermal inertia. Spectral emittance images were derived for the test area for data sets from both the day and night tests, and they show good qualitative agreement. Comparison of the numerical values of emittance derived from these day and night images shows a decrease in spectral contrast at night. This is probably due primarily to an increased atmospheric contribution to the radiance reaching the sensor at night when the ground is cold, rather than to a change in spectral characteristics of the surface at night. These spectral emittance data contribute to an und...

Physical properties of shale at temperature and pressure

Abstract: Deformational properties, P‐wave and S‐wave velocities, and electrical resistivity were measured for three North Sea Malm shales in the laboratory under pressures to 800 bars and temperatures to 100 °C. These data were used to evaluate how factors such as mineralogy, microstructure, compaction, and pore‐fluid conductivity affect a shale’s seismic and electrical responses. Deformation in the shales is dominated by inelastic processes which cause time‐dependent changes in velocity, resistivity, and pore pressure. Overall, shales are less sensitive to pressure changes as compared to sandstones of similar porosity. However, changes in temperature result in large changes in physical properties as compared to sandstones or shaly sands. P‐wave and S‐wave velocities may decrease by as much as 10 percent over the temperature range studied, and calculated activation energies for surface conduction are nearly twice those observed in shaly sands. These comparisons emphasize fundamental differences in fabric among the...

Use of ground-probing radar and resistivity surveys for archaeological investigations

Abstract: In Japan, geophysical methods are normally used to estimate the distribution of cultural relics before digging. Objects of archaeological interest are usually located within a few meters of the surface. Therefore, geophysical methods suitable for archaeological exploration are those which provide high resolution at shallow depths. The most commonly used geophysical methods are ground‐probing radar, resistivity, and magnetometry. Of these methods, we used mainly ground‐probing radar and resistivity surveys in archaeological investigations at four sites. Three of the sites were in Gumma Prefecture (Japan); they were covered with volcanic deposits (loam or pumice). Using ground‐probing radar, we were able to locate ancient dwellings, burial mounds, and a distribution of archaeologically significant “culture layers.” At the other site, in Nara Prefecture, we located part of the remains of an ancient city. In this investigation, the resistivity method and ground‐probing radar were combined to determine the loc...

Analysis of electrical conduction in the grain consolidation model

Abstract: In the grain consolidation model the diagenetic processes of compaction and cementation are represented in terms of the growth of an array of originally spherical grains. Grain growth toward the nodes of the pore space leads to an electrical formation factor F(ϕ) that increases slowly as the porosity ϕ decreases. By contrast, grain growth toward the throats of the pore space leads to a rapidly increasing F(ϕ). In all the cases we have examined, the value of the percolation threshold, ϕc is less than 0.055. Network simulation techniques have been developed to calculate the electrical conductivity of the ordered versions of the grain consolidation model. We find that the minimum‐area approximation employed in our earlier work is generally quite satisfactory. The network techniques can also be used to model the effects of mixed pore‐space fluid saturation, with results that are physically reasonable although not necessarily in agreement with empirical rules regarding saturation.

Two and one‐half dimensional Born inversion with an arbitrary reference

Abstract: Multidimensional inversion algorithms are presented for both prestack and poststack data gathered on a single line. These algorithms both image the subsurface (i.e., give a migrated section) and, given relative true amplitude data, estimate reflection strength or impedance on each reflector. The algorithms are “two and one‐half dimensional” (2.5-D) in that they incorporate three‐dimensional (3-D) wave propagation in a medium which varies in only two dimensions. The use of 3-D sources does not entail any computational penalty, and it avoids the serious degradation of amplitude incurred by using the 2-D wave equation. Our methods are based on the linearized inversion theory associated with the “Born inversion.” Thus, we assume that the sound speed profile is well approximated by a given background velocity, plus a perturbation. It is this perturbation that we seek to reconstruct. We are able to treat the case of an arbitrary continuous background profile. However, the cost of implementation increases as one...

The effects of transverse isotropy on reflection amplitude versus offset

Abstract: Studies have shown that elastic properties of materials such as shale and chalk are anisotropic. With the increasing emphasis on extraction of lithology and fluid content from changes in reflection amplitude with shot‐to‐group offset, one needs to know the effects of anisotropy on reflectivity. Since anisotropy means that velocity depends upon the direction of propagation, this angular dependence of velocity is expected to influence reflectivity changes with offset. These effects might be particularly evident in deltaic sand‐shale sequences since measurements have shown that the P-wave velocity of shales in the horizontal direction can be 20 percent higher than the vertical P-wave velocity. To investigate this behavior, a computer program was written to find the P- and S-wave reflectivities at an interface between two transversely isotropic media with the axis of symmetry perpendicular to the interface. Models for shale‐chalk and shale‐sand P-wave reflectivities were analyzed.

Automatic phase correction of common-midpoint stacked data

Abstract: The residual wavelet on a processed seismic section is often not zero phase despite all efforts to make it so. In this paper we adopt the convolutional model for the processed seismogram, assume that the residual phase shift can be approximated by a frequency-independent constant, and use the varimax norm to generate an algorithm to estimate the residual phase directly. Application of our algorithm to reflectivities from well logs suggests that it should work in the majority of cases so long as the reflectivity is non-Gaussian. An application of our algorithm to stacked data enhances the interpretability of the seismic section and leads to an improved match between the recovered relative acoustic impedance and a measured velocity log.

Wave-field transformations of vertical seismic profiles

Abstract: Vertical seismic profile (VSP) data may be partitioned in a variety of ways by application of wave‐field transformations. These transformations provide insights into the nature of the data and aid in the design of processing operations. Transformations are implemented in a reversible sequence that takes the observed VSP data from the depth‐time (z-t) domain through the slowness‐time intercept (p-τ) domain (by a slant stack), to the slowness‐frequency (p-ω) domain (by a 1-D Fourier transform over τ), to the wavenumber‐frequency (k-ω) domain (by resampling using the Fourier central‐slice theorem), and finally back to the z-t domain (by an inverse 2-D Fourier transform). Multidimensional wave‐field transformations, combined with k-ω, p-ω, and p-τ filtering, can be applied to wave‐field resampling, interpolation, and extrapolation; separation of P-waves and S-waves; separation of upgoing and downgoing waves; and wave‐field decomposition for isolation, identification, and analysis of arrivals.