ALBERT

Notes: Based on an acoustic assumption (that the shear-wave velocity is zero) and a dispersion relationship, we derive an acoustic wave equation for P-waves in tilted transversely isotropic (TTI) media (transversely isotropic media with a tilted symmetry axis). This equation has fewer parameters than an elastic wave equation in TTI media and yields an accurate description of P-wave traveltimes and spreading-related attenuation. Our TTI acoustic wave equation is a fourth-order equation in time and space. We demonstrate that the acoustic approximation allows the presence of shear waves in the solution. The substantial differences in traveltime and amplitude between data created using vertical transversely isotropic (VTI) and TTI assumptions is illustrated in examples.

Notes: A high-resolution P-wave seismic reflection survey was conducted in the area of Lambro Park within the city of Milan (northern Italy). Several high permeability channel sequences were identified and the stratigraphy of the infill was accurately mapped. The seismic signature imaged a previously undefined unit and also established the lateral correlation of some depositional units within the three major aquifer groups along the survey line.The complexity of the acoustic framework and the reduced depth of the source location limited effective elastic-wave generation. The geological setting also caused propagation of strong coherent noise patterns. Further interference, observed in the recorded data, was due to the traffic noise from the nearby highway. The attenuation of the undesired events required the design of specific filters and their multistep implementation.The results of forward modelling based on borehole information and of noise tests were crucial factors in the design of the processing parameters and in the stratigraphic interpretation of the final stacked section.

Notes: The subsurface geometry of the Sebastián Vizcaíno Basin is obtained from the 2D inversion of gravity data, constrained by a density-versus-depth relationship derived from an oil exploration deep hole. The basin accumulated a thick pile of marine sediments that evolved in the fore-arc region of the compressive margin prevalent along western North America during Mesozoic and Tertiary times. Our interpretation indicates that the sedimentary infill in the Sebastián Vizcaíno Basin reaches a maximum thickness of about 4 km at the centre of a relatively symmetric basin. At the location of the Suaro-1 hole, the depth to the basement derived from this work agrees with the drilled interface between calcareous and volcaniclastic members of the Alisitos Formation. A sensitivity analysis strongly suggests that the assumed density function leads to a nearly unique solution of the inverse problem.

Notes: We describe a method to invert a walkaway vertical seismic profile (VSP) and predict elastic properties (P-wave velocity, S-wave velocity and density) in a layered model looking ahead of the deepest receiver. Starting from Bayes's rule, we define a posterior distribution of layered models that combines prior information (on the overall variability of and correlations among the elastic properties observed in well logs) with information provided by the VSP data. This posterior distribution of layered models is sampled by a Monte-Carlo method. The sampled layered models agree with prior information and fit the VSP data, and their overall variability defines the uncertainty in the predicted elastic properties. We apply this technique first to a zero-offset VSP data set, and show that uncertainty in the long-wavelength P-wave velocity structure results in a sizable uncertainty in the predicted elastic properties. We then use walkaway VSP data, which contain information on the long-wavelength P-wave velocity (in the reflection moveout) and on S-wave velocity and density contrasts (in the change of reflectivity with offset). The uncertainty of the look-ahead prediction is considerably decreased compared with the zero-offset VSP, and the predicted elastic properties are in good agreement with well-log measurements.

Notes: The existence of non-zero reflectivity for zero-offset PS waves in horizontally layered media is discussed. Field measurements have suggested the occurrence of this phenomenon. We show that, in some cases, if anisotropy is considered in the problem, this reflectivity can be predicted. By using an approximated formulation to determine the coefficient of reflection in arbitrarily anisotropic media, it is clear that, for some elastic arrangements of the subsurface, reflected energy is associated with the converted wave resulting from normal P-wave incidence.

Notes: A simple scaling relationship is shown to facilitate comparison, correlation and integration of data recorded using the common experimental configurations in electrical and electromagnetic depth sounding. Applications of the scheme to field data from typical geological and landfill environments show that it is robust and, where transient electromagnetic (TEM) data are available, enables easy identification and quantification of electrical static shift (galvanic distortion) in magnetotelluric and direct current (DC) sounding curves. TEM-based procedures are suggested for both the direct removal of static shift in DC sounding curves and effective joint data inversion with the most-squares criterion in the presence of static shift. A case study of aquifer characterization using sounding data from borehole sites in the Vale of York in England shows that static shift is a common problem in this glacial-covered terrain and demonstrates the effectiveness of the proposed joint DC-TEM inversion strategy in handling distorted soundings.

Notes: A three-dimensional (3D) electrical resistivity modelling code is developed to interpret surface and subsurface data. Based on the integral equation, it calculates the charge density caused by conductivity gradients at each interface of the mesh, allowing the estimation of the potential everywhere without the need to interpolate between nodes. Modelling generates a huge matrix, made up of Green's functions, which is stored by using the method of pyramidal compression. The potential is compared with the analytical and the numerical solutions obtained by finite-difference codes for two models: the two-layer case and the vertical contact case. The integral method is more accurate around the source point and at the limits of the domain for the potential calculation using a pole-pole array. A technique is proposed to calculate the sensitivity (Jacobian) and Hessian matrices in 3D. The sensitivity is based on the derivative with respect to the block conductivity of the potential computed using the integral equation; it is only necessary to compute the electrical field at the source location. A direct extension of this technique allows the determination of the second derivatives. The technique is compared with the analytical solutions and with the calculation of the sensitivity according to the method using the inner product of the current densities calculated at the source and receiver points. Results are very accurate when the Green's function that includes the source image is used. The calculation of the three components of the electric field on the interfaces of the mesh is carried out simultaneously and quickly, using matrix compression.

Notes: In a sedimentary environment, layered models are often capable of representing the actual geology more accurately than smooth minimum structure models. Furthermore, interval thicknesses and resistivities are often the parameters to which non-geophysicist experts can relate and base decisions on when using them in waste site remediation, groundwater modelling and physical planning.We present a laterally constrained inversion scheme for continuous resistivity data based on a layered earth model (1D). All 1D data sets and models are inverted as one system, producing layered sections with lateral smooth transitions. The models are regularized through laterally equal constraints that tie interface depths and resistivities of adjacent layers. Prior information, e.g. originating from electric logs, migrates through the lateral constraints to the adjacent models, making resolution of equivalences possible to some extent. Information from areas with well-resolved parameters will migrate through the constraints in a similar way to help resolve the poorly constrained parameters. The estimated model is complemented by a full sensitivity analysis of the model parameters, supporting quantitative evaluation of the inversion result.Examples from synthetic 2D models show that the model recognition of a sublayered 2D wedge model is improved using the laterally constrained inversion approach when compared with a section of combined 1D models and when compared with a 2D minimum structure inversion. Case histories with data from two different continuous DC systems support the conclusions drawn from the synthetic example.

Notes: We analyse the problem of radiation of seismic waves by a vibroseis source when the baseplate is subject to flexure. A theoretical model is proposed to account for baseplate flexure, generalizing the well-known model of the vibroseis source of Sallas and Weber, which was developed for a rigid plate. Using the model proposed, we analyse the effect of flexure on the properties of seismic waves. We show that the flexure does not contribute to the far-field and mainly affects the readings of the reference accelerometer that is used to measure the force applied to the ground; these readings generally become dependent on the location of the sensor on the plate. For muddy and sandy soils, the effect of flexure on baseplate-acceleration measurements is nonetheless pronounced at the high end of the vibroseis frequency band only (∼100 Hz), and is negligible at all frequencies for stiffer soils. The corresponding phase lags introduced by the flexural vibrations at high frequencies lead to errors in the traveltime measurements (through the cross-correlation function) of up to 0.6 ms for muddy soils and less for denser soils. We show the existence of an optimal position of the reference sensor on the baseplate and also propose a general method of eliminating the phase lag due to the baseplate flexure in acceleration measurements.

Notes: Scalar radiomagnetotelluric measurements were carried out on a contaminated test area close to the Brazi Refinery in Romania in order to detect and to monitor a 1 m thick oil layer expected at 5 m depth. Radio transmitters broadcasting in a frequency range from 10 kHz to 300 kHz were selected to observe the apparent resistivity and the phase data associated with the E- and B-polarizations. They were located parallel and perpendicular to the assumed strike direction of the contamination plume. The data were interpreted by a 2D inversion technique from which the conductivity structure of the area was derived. The 2D inversion models of all profiles on the contaminated area show a poor-conductivity zone above the groundwater table which could be associated with the oil contamination.A first attempt was also made to monitor the contaminated layer: the radiomagnetotelluric measurements were repeated on the same profiles a year later, but this time in a dry period, not in a rainy one. The 2D inversion results of the measurements in the dry period indicate that the high-resistivity layer moved closer to the surface. Additional reference measurements were then carried out on a non-contaminated area situated at a distance from the refinery, in the opposite direction to the flow of the groundwater. These reference measurements were used for the derivation of the unperturbed geology and they were also compared with the measurements of the contaminated test area. There is a significant difference in the frequency dependences of the apparent resistivities of the reference and contaminated areas, which could indicate a contamination at shallow depth. The 2D inversion results show the increase of resistivity at a depth of about 5 m beneath the contaminated area where the oil contamination is expected according to the information from the boreholes.

Notes: Computer simulations are used to calculate the elastic properties of model cemented sandstones composed of two or more mineral phases. Two idealized models are considered – a grain-overlap clay/quartz mix and a pore-lining clay/quartz mix. Unlike experimental data, the numerical data exhibit little noise yet cover a wide range of quartz/cement ratios and porosities. The results of the computations are in good agreement with experimental data for clay-bearing consolidated sandstones.The effective modulus of solid mineral mixtures is found to be relatively insensitive to microstructural detail. It is shown that the Hashin–Shtrikman average is a good estimate for the modulus of the solid mineral mixtures. The distribution of the cement phase is found to have little effect on the computed modulus–porosity relationships. Numerical data for dry and saturated states confirm that Gassmann's equations remain valid for porous materials composed of multiple solid constituents. As noted previously, the Krief relationship successfully describes the porosity dependence of the dry shear modulus, and a recent empirical relationship provides a good estimate for the dry-rock Poisson's ratio.From the numerical computations, a new empirical model, which requires only a knowledge of system mineralogy, is proposed for the modulus–porosity relationship of isotropic dry or fluid-saturated porous materials composed of multiple solid constituents. Comparisons with experimental data for clean and shaly sandstones and computations for more complex, three-mineral (quartz/dolomite/clay) systems show good agreement with the proposed model over a very wide range of porosities.

Notes: With the increasing use of permanently installed seismic installations, many of the issues in time-lapse seismic caused by the lack of repeatability can be reduced. However, a number of parameters still influence the degree of reliability of 4D seismic data. In this paper, the specific impact of seawater velocity variations on time-lapse repeatability is investigated in a synthetic study. A zero-lag time-lapse seabed experiment with no change in the subsurface but with velocity changes in the water column is simulated. The velocity model in the water column is constant for the baseline survey while the model for the repeat survey is heterogeneous, designed from sea salinity and temperature measurements in the West of Shetlands. The difference section shows up to 80% of residual amplitude, which highlights the poor repeatability. A new dynamic correction which removes the effect of seawater velocity variations specifically for permanent installations is developed. When applied to the synthetic data, it reduces the difference residual amplitude to about 3%. This technique shows substantial improvement in repeatability beyond conventional time-lapse cross-equalization.

Notes: A new method for time-lapse signal separation and enhancement using singular-value decomposition is presented. Singular-value decomposition is used to separate a 4D signal into its constituent parts: common geology, time-lapse response and noise. Synthetic tests which demonstrate the advantages of the singular-value decomposition technique over traditional differencing methods are also presented. This signal separation and enhancement technique is used to map out both the original and moved oil–water contacts across the Nelson Field. The singular-value decomposition technique allows the oil–water contact to be mapped across regions which would have been missed using traditional differencing methods. In particular, areas toward the edges of the field are highlighted by the technique. The oil–water contact is observed to move upwards across the field, with the largest movements being associated, as anticipated, with natural production. The results obtained are broadly consistent with those predicted by the reservoir simulator model. Singular-value decomposition is demonstrated to be a useful tool for enhancing the time-lapse signal and for gaining confidence in areas where traditional differencing fails.

Notes: This work looks at the application of neural networks in geophysical well-logging problems and specifically their utilization for inversion of nuclear downhole data. Simulated neutron and γ-ray fluxes at a given detector location within a neutron logging tool were inverted to obtain formation properties such as porosity, salinity and oil/water saturation. To achieve this, the forward particle-radiation transport problem was first solved for different energy groups (47 neutron groups and 20 γ-ray groups) using the multigroup code EVENT. A neural network for each of the neutron and γ-ray energy groups was trained to re-produce the detector fluxes using the forward modelling results from 504 scenarios. The networks were subsequently tested on unseen data sets and the unseen input parameters (formation properties) were then predicted using a global search procedure. The results obtained are very encouraging with formation properties being predicted to within 10% average relative error. The examples presented show that neural networks can be applied successfully to nuclear well-logging problems. This enables the implementation of a fast inversion procedure, yielding quick and reliable values for unknown subsurface properties such as porosity, salinity and oil saturation.

Notes: Marine gravimeters mounted on stabilized platforms are commonly used in aircraft to perform airborne gravity measurements. The role of the stabilized platform is to level the sensor mechanically, whatever the aircraft attitude. However, this compensation is generally insufficient due to the sensitivity of modern gravity sensors. Correcting the offlevel error requires that an offlevel correction calculated from positioning data be added to gravimeter measurements, which complicates not only the processing, but also the assessment of precision and resolution.This paper is a feasibility study describing the levelling of a completely strapped-down LaCoste and Romberg gravimeter for airborne gravimetry operation, by means of GPS positioning data. It focuses on the calculation of the sensor offlevel correction needed for the complete gravity data processing. The precision of the offlevel correction that can be achieved, in terms of GPS data precision and gravity wavelengths, is theoretically studied and estimated using the gravity and GPS data acquired during the Alpine Swiss French airborne gravity survey carried out in 1998 over the French Western Alps. While a 1 cm precision of GPS-determined baseline coordinates is sufficient to achieve a 5 mGal precision of the offlevel correction, we maintain that this precision has to reach 1 mm to ensure a 1 mGal precision of the offlevel correction at any wavelength.Without a stabilized platform, the onboard instrumentation becomes significantly lighter. Furthermore, the correction for the offlevel error is straightforward and calculated only from GPS data. Thus, the precision and the resolution of airborne gravity surveys should be estimated with a better accuracy.

Notes: Dip-moveout (DMO) correction is often applied to common-offset sections of seismic data using a homogeneous isotropic medium assumption, which results in a fast execution. Velocity-residual DMO is developed to correct for the medium-treatment limitation of the fast DMO. For reasonable-sized velocity perturbations, the residual DMO operator is small, and thus is an efficient means of applying a conventional Kirchhoff approach. However, the shape of the residual DMO operator is complicated and may form caustics. We use the Fourier domain for the operator development part of the residual DMO, while performing the convolution with common-offset data in the space–time domain. Since the application is based on an integral (Kirchhoff) method, this residual DMO preserves all the flexibility features of an integral DMO. An application to synthetic and real data demonstrates effectiveness of the velocity-residual DMO in data processing and velocity analysis.

Notes: The finite-difference method applied to the full 3D wave equation is a rather time-consuming process. However, in the 2.5D case, we can take advantage of the medium symmetry. By taking the Fourier transform with respect to the out-of-plane direction (the symmetry axis), the 3D problem can be reduced to a repeated 2D problem. The third dimension is taken into account by a sum over the corresponding wave-vector component. A criterion for where to end this theoretically infinite sum derives from the stability conditions of the finite-difference schemes employed. In this way, the computation time of the finite-difference calculations can be considerably reduced. The quality of the modelling results obtained with this 2.5D finite-difference scheme is comparable to that obtained using a standard 3D finite-difference scheme.

Notes: For non-linear kinematic inversion of elastic anisotropy parameters and related investigations of the sensitivity of seismic data, the derivatives of the wavespeed (phase velocity and group velocity) with respect to the individual elastic moduli are required. This paper presents two analytic methods, called the eigenvalue and eigenvector methods, to compute the derivatives of the wavespeeds for wave propagation in a general anisotropic medium, which may be defined by up to 21 density-normalized elastic moduli. The first method employs a simple and compact form of the eigenvalue (phase velocity) and a general form of the group velocity, and directly yields general expressions of the derivatives for the three wave modes (qP, qS1, qS2). The second method applies simple eigenvector solutions of the three wave modes and leads to other general forms of the derivatives. These analytic formulae show that the derivatives are, in general, functions of the 21 elastic moduli as well as the wave propagation direction, and they reflect the sensitivity of the wavespeeds to the individual elastic moduli. Meanwhile, we give results of numerical investigations with some examples for particular simplified forms of anisotropy. They show that the eigenvalue method is suitable for the qP-, qS1- and qS2-wave computations and mitigates the singularity problem for the two quasi-shear waves. The eigenvector method is preferable to the eigenvalue method for the group velocity and the derivative of the phase velocity because it involves simpler expressions and independent computations, but for the derivative of the group velocity the derivative of the eigenvector is required. Both methods tackle the singularity problem and are applicable to any degree of seismic anisotropy for all three wave modes.

Notes: The phase error between the real phase shift and the Gazdag background phase shift, due to lateral velocity variations about a reference velocity, can be decomposed into axial and paraxial phase errors. The axial phase error depends only on velocity perturbations and hence can be completely removed by the split-step Fourier method. The paraxial phase error is a cross function of velocity perturbations and propagation angles. The cross function can be approximated with various differential operators by allowing the coefficients to vary with velocity perturbations and propagation angles. These variable-coefficient operators require finite-difference numerical implementation. Broadband constant-coefficient operators may provide an efficient alternative that approximates the cross function within the split-step framework and allows implementation using Fourier transforms alone. The resulting migration accuracy depends on the localization of the constant-coefficient operators. A simple broadband constant-coefficient operator has been designed and is tested with the SEG/EAEG salt model. Compared with the split-step Fourier method that applies to either weak-contrast media or at small propagation angles, this operator improves wavefield extrapolation for large to strong lateral heterogeneities, except within the weak-contrast region. Incorporating the split-step Fourier operator into a hybrid implementation can eliminate the poor performance of the broadband constant-coefficient operator in the weak-contrast region. This study may indicate a direction of improving the split-step Fourier method, with little loss of efficiency, while allowing it to remain faster than more precise methods such as the Fourier finite-difference method.

Notes: Multi-refractor imaging is a technique for constructing a single two-dimensional image of a number of refractors by stacking multiple convolved and cross-correlated reversed shot records. The method is most effective with high-fold data that have been obtained with roll-along acquisition programs because the stacking process significantly improves the signal-to-noise ratios. The major advantage of the multi-refractor imaging method is that all the data can be stacked to maximize the signal-to-noise ratios before the measurement of any traveltimes. However, the signal-to-noise ratios can be further increased if only those traces that have arrivals from the same refractor are used, and if the correct reciprocal times or traces are employed.A field case study shows that multi-refractor imaging can produce a cross-section similar to the familiar reflection cross-section with substantially higher signal-to-noise ratios for the equivalent interfaces.

Notes: Combined time-lapse reservoir simulation and seismic modelling has been performed on both 1D and 3D models of a channelized turbidite reservoir. The models have been built using core, log, laboratory and seismic data from the Nelson Field (central North Sea) as a template. Oil and water movement in the main channels, channel margins and interchannel regions is investigated, with a particular focus being the effect of poor net-to-gross. The analysis confirms that saturation effects dominate the response whilst stress-sensitivity effects play a minor role. The trough–peak signature in the seismic difference volumes formed by the sweep of the water can be continued and mapped slightly further than the channel margins. This characteristic 4D signature remains roughly intact, despite the complicated depositional architecture, and accurately delineates the area of moved fluid, but it needs additional calibration to combat the detrimental influence of the low net-to-gross. Signal strength is largely dependent on reservoir quality, but is also compromised by the net-to-gross, fluid distribution and, more critically, by the exact timing of the seismic survey. For example, a region of bypassed oil zone remains undetected as it forms early during the production. This work demonstrates that to understand fully the 4D signature at a quantitative level requires adequate knowledge of the fluid properties, but also, more critically, the geology.

Notes: In the central North Sea ‘Gannet-A’ field, a 50 ft oil rim is overlain by a gas cap of variable thickness. Oil is produced from horizontal wells which initially produced dry oil, but as the field became more mature, a significant water cut was seen in several wells. A dedicated 4D seismic monitor survey was acquired in order to assess the remaining distribution of oil reserves. By forward modelling the synthetic seismic response to parameters such as contact movement and residual saturations (using 2D and 3D wedge models), and comparing the results with real seismic data, we are able to decipher the contact movements across the field. It is shown that, in one part of the field, the increased water cut is caused primarily by the vertical displacement of the entire oil rim into the initial gas cap. This oil-rim displacement produces a very different 4D seismic response from the case of a static gas–oil contact and rising oil–water contact (normal production). As a result of these observations, we are able to optimize field production by both re-perforation of existing wells and by drilling sidetracks into the displaced rim: a brown-field development opportunity that might otherwise be missed.

Notes: The cross-calibration of different vintage data is an important prerequisite in attempting to determine the time-lapse seismic effects induced by hydrocarbon production in a reservoir. This paper reports the preprocessing and cross-calibration procedures adopted to modify the data of four seismic vintages (1982, 1989, 1992 and 1999) from the Oseberg field in the North Sea, for optimal conditions for a time-lapse seismic amplitude analysis. The final results, in terms of time-lapse variations, of acoustic impedance and of amplitude-versus-offset, are illustrated for selected data sets. The application of preprocessing to each individual vintage data set reduces the effects of the different acquisition and noise conditions, and leads to consistency in the amplitude response of the four vintages. This consistency facilitates the final amplitude cross-calibration that is carried out using, as reference, the Cretaceous horizon reflections above the Brent reservoir. Such cross-calibration can be considered as vintage-consistent residual amplitude correction.Acoustic impedance sections, intercept and gradient amplitude-versus-offset attributes and coherent amplitude-versus-offset estimates are computed on the final cross-calibrated data. The results, shown for three spatially coincident 2D lines selected from the 1982, 1989 and 1999 data sets, clearly indicate gas-cap expansion resulting from oil production. Such expansion is manifested as a decrease in acoustic impedance and a modification of the amplitude-versus-offset trends in the apical part of the reservoir.

Notes: This paper presents a new algorithm for estimating non-minimum-phase seismic wavelets by using the second- and higher-order statistics (HOS) of the wavelets. In contrast to many, if not most, of the HOS-based methods, the proposed method does not need to assume that subsurface seismic reflectivity is a non-Gaussian, statistically independent and identically distributed random process. The amplitude and phase spectra of the wavelets are estimated, respectively, using the second-order statistics (SOS) and third-order moment (TOM) of the wavelets, which will, in turn, be derived from the HOS of the seismic traces. In our approach, the wavelets can be ‘calculated’ from seismic traces efficiently; no optimization or inversion is necessarily required. Very good results have been obtained by applying this method to both synthetic and real-field data sets.

Notes: We study the propagation of elastic waves that are generated in a fluid-filled borehole surrounded by a cracked transversely isotropic medium. In the model studied the anisotropy and borehole axes coincide. To obtain the effective elastic moduli of a cracked medium we have applied Hudson's theory that enables the determination of the overall properties as a function of the crack orientation in relation to the symmetry axis of the anisotropic medium. This theory takes into account the hydrodynamic mechanism of the elastic-wave attenuation caused by fluid filtration from the cracks into a porous matrix.We have simulated the full waveforms generated by an impulse source of finite length placed on the borehole axis. The kinematic and dynamic parameters of the compressional, shear and Stoneley waves as functions of the matrix permeability, crack orientation and porosity were studied. The modelling results demonstrated the influence of the crack-system parameters (orientation and porosity) on the velocities and amplitudes of all wave types. The horizontally orientated cracks result in maximal decrease of the elastic-wave parameters (velocities and amplitudes).Based on the fact that the shear- and Stoneley-wave velocities in a transversely isotropic medium are determined by different shear moduli, we demonstrate the feasibility of the acoustic log to identify formations with close to horizontal crack orientations.

Notes: The development of an electromagnetic numerical modelling scheme for a magnetic dipole in an arbitrary casing segment in an inhomogeneous conductivity background has been difficult, due to the very high electrical conductivity and magnetic permeability contrasts between the steel casing and the background medium. To investigate the effect of steel casing efficiently, we have developed an accurate but simple finite-element modelling scheme to simulate electromagnetic fields in a medium of cylindrically symmetric conductivity structures. In order to preserve the cylindrical symmetry in the resulting electromagnetic fields, a horizontal loop current source is used throughout. One of the main advantages of the approach is that the problem is scalar when formulated using the azimuthal electric field, even if the casing is both electrically conductive and magnetically permeable. Field calculations have been made inside the cased borehole as well as in another borehole which is not cased. Careful analyses of the numerical modelling results indicate that the anomaly observed in a cross-borehole configuration is sensitive enough to be used for tomographic imaging.

Notes: Seismic wavefield reconstruction is posed as an inversion problem where, from inadequate and incomplete data, we attempt to recover the data we would have acquired with a denser distribution of sources and receivers. A minimum weighted norm interpolation method is proposed to interpolate prestack volumes before wave-equation amplitude versus angle imaging.Synthetic and real data were used to investigate the effectiveness of our wavefield reconstruction scheme when preconditioning seismic data for wave-equation amplitude versus angle imaging.

Notes: We estimate the concentration of gas hydrate and free gas at an area located to the north of the Knipovich Ridge (western Svalbard margin). The method is based on P-wave velocities computed by reflection tomography applied to multicomponent ocean-bottom seismometer data. The tomographic velocity field is fitted to theoretical velocities obtained from a poro-elastic model based on a Biot-type approach (the interaction between the rock frame, gas hydrate and fluid is modelled from first physical principles). We obtain average hydrate concentrations of 7% and maximum free-gas saturations of 0.4% and 9%, depending on the saturation model.

Notes: Fluid-flow simulators used in the oil industry model the movement of fluids through a porous reservoir rock. These simulators either ignore coupling between the flow and concurring deformation of the solid rock frame or take it into account approximately, in the so-called loose or staggered-in-time mode. In contrast to existing simulators, the one we describe here fully couples two-phase (oil and water) flow to subsurface deformation and simultaneously accounts for all relevant physical phenomena. As such, our flow simulator inherently links time-dependent fluid pressures, saturations, permeabilities and flow velocities to stresses in the whole subsurface. These stresses relate to strains through the non-linear theory of elasticity, allowing us to model time-lapse changes in seismic velocities and anisotropy. The velocity variations manifest themselves in time shifts and reflection amplitudes that are conventionally measured from 4D seismic data. Changes in anisotropy produce time-dependent shear-wave splitting that can be used for monitoring the horizontal stresses.

Notes: Many joint inversion schemes use 1D forward modelling in the integrated interpretation of various geophysical data. In extending the joint inversion approach to the investigation of 2D structures, the discretization of the model parameters and the appropriate choice of the forward-modelling procedure play a very important role. In this paper, a hybrid seismic–geoelectric joint inversion method is proposed for the investigation of 2D near-surface geological structures. The electric and seismic models are coupled together through the use of common boundaries between the adjacent layers. Assuming a 2D model composed of homogeneous layers with curved boundaries, a fast ray-tracing algorithm is used for the calculation of refraction seismic traveltime data. In the geoelectric forward modelling, a locally 1D approximation is used. The boundary surfaces are written in the form of series expansion; the inversion algorithms are formulated for the expansion coefficients and the petrophysical parameters as unknowns. Two versions of the inversion method are proposed: in versions A and B, interval-wise constant functions and Chebyshev polynomials are, respectively, used as basis functions of the series expansion. The versions are tested by means of synthetic and in situ measured data. The tests show that both methods are stable and accurate.

Notes: Ghawar, the largest oilfield in the world, produces oil from the Upper Jurassic Arab-D carbonate reservoir. The high rigidity of the limestone–dolomite reservoir rock matrix and the small contrast between the elastic properties of the pore fluids, i.e. oil and water, are responsible for the weak 4D seismic effect due to oil production. A feasibility study was recently completed to quantify the 4D seismic response of reservoir saturation changes as brine replaced oil. The study consisted of analysing reservoir rock physics, petro-acoustic data and seismic modelling. A seismic model of flow simulation using fluid substitution concluded that time-lapse surface seismic or conventional 4D seismic is unlikely to detect the floodfront within the repeatability of surface seismic measurements. Thus, an alternative approach to 4D seismic for reservoir fluid monitoring is proposed. Permanent seismic sensors could be installed in a borehole and on the surface for passive monitoring of microseismic activity from reservoir pore-pressure perturbations. Reservoir production and injection operations create these pressure or stress perturbations. Reservoir heterogeneities affecting the fluid flow could be mapped by recording the distribution of epicentre locations of these microseisms or small earthquakes. The permanent borehole sensors could also record repeated offset vertical seismic profiling surveys using a surface source at a fixed location to ensure repeatability. The repeated vertical seismic profiling could image the change in reservoir properties with production.

Notes: For multivalued traveltime computation on dense grids, we propose a wavefront-orientated ray-tracing (WRT) technique. At the source, we start with a few rays which are propagated stepwise through a smooth two-dimensional (2D) velocity model. The ray field is examined at wavefronts and a new ray might be inserted between two adjacent rays if one of the following criteria is satisfied: (1) the distance between the two rays is larger than a predefined threshold; (2) the difference in wavefront curvature between the rays is larger than a predefined threshold; (3) the adjacent rays intersect. The last two criteria may lead to oversampling by rays in caustic regions. To avoid this oversampling, we do not insert a ray if the distance between adjacent rays is smaller than a predefined threshold. We insert the new ray by tracing it from the source. This approach leads to an improved accuracy compared with the insertion of a new ray by interpolation, which is the method usually applied in wavefront construction. The traveltimes computed along the rays are used for the estimation of traveltimes on a rectangular grid. This estimation is carried out within a region bounded by adjacent wavefronts and rays. As for the insertion criterion, we consider the wavefront curvature and extrapolate the traveltimes, up to the second order, from the intersection points between rays and wavefronts to a gridpoint. The extrapolated values are weighted with respect to the distances to wavefronts and rays. Because dynamic ray tracing is not applied, we approximate the wavefront curvature at a given point using the slowness vector at this point and an adjacent point on the same wavefront. The efficiency of the WRT technique is strongly dependent on the input parameters which control the wavefront and ray densities. On the basis of traveltimes computed in a smoothed Marmousi model, we analyse these dependences and suggest some rules for a correct choice of input parameters. With suitable input parameters, the WRT technique allows an accurate traveltime computation using a small number of rays and wavefronts.

Notes: We present a new least-squares migration method called least-squares wave-path migration. The proposed method combines an iterative conjugate-gradient solver with a stationary-phase wave-path migration operator. Numerical tests demonstrate that (i) least-squares wave-path migration is computationally more efficient than and almost as accurate as Kirchhoff least-squares migration, and (ii) many of the artefacts seen in wave-path migration images are suppressed after several conjugate-gradient iterations. Previous results have shown that 3D wave-path migration is up to 100 times faster than a standard 3D Kirchhoff migration, but sometimes at the cost of reduced quality. With the proposed least-squares wave-path migration method, the image quality in wave-path migration can be improved at an acceptable increase in computational cost.

Notes: Imaging pre-salt reflections for data acquired from the coastal region of the Red Sea is a task that requires prestack migration velocity analysis. Conventional post-stack time processing lacks the lateral inhomogeneity capability, necessary for such a problem. Prestack migration velocity analysis in the vertical time domain reduces the velocity–depth ambiguity that usually hampers the performance of prestack depth-migration velocity analysis. In prestack τ-migration velocity analysis, the interval velocity model and the output images are defined in τ (i.e. vertical time). As a result, we avoid placing reflectors at erroneous depths during the velocity analysis process and thus avoid inaccurately altering the shape of the velocity model, which in turn speeds up the convergence to the true model. Using a 1D velocity update scheme, the prestack τ-migration velocity analysis produces good images of data from the Midyan region of the Red Sea. For the first seismic line from this region, only three prestack τ-migration velocity analysis iterations were required to focus pre-salt reflections in τ. However, the second line, which crosses the first line, is slightly more complicated and thus required five iterations to reach the final, reasonably focused, τ-image. After mapping the images for the two crossing lines to depth, using the final velocity models, the placements of reflectors in the two 2D lines were consistent at their crossing point. Some errors occurred due to the influence of out-of-plane reflections on 2D imaging. However, such errors are identifiable and are generally small.

Notes: We have developed a method for imaging magnetic data collected for mineral exploration to yield the following structural information: depth, model type (structural index) and susceptibility. The active nature of mineral exploration data requires we derive the structural information from a robust quantity: we propose that the first- or second-order analytic-signal amplitude is suitably stable. The procedure is to normalize the analytic-signal amplitude by the peak value and then use non-linear inversion to estimate the depth and the structural index for each anomaly. In our field example, different results are obtained depending on whether we inverted for the first- or second-order analytic-signal amplitude. This is probably because the two-dimensional contact, thin sheet or horizontal cylinder models we have assumed are not appropriate. In cases such as these, when our model assumptions are not correct, the results should not be interpreted quantitatively, but they might be useful for giving a qualitative indication of how the structure might vary.With a priori information, it is possible to assume a model type (i.e. set the structural index) and generate estimates of the depth and susceptibility. These data can then be gridded and imaged. If a contact is assumed, the susceptibility contrast is estimated; for the dike model, the susceptibility-thickness is estimated; for the horizontal cylinder, the susceptibility-area is estimated. To emphasize that the results are dependent on our assumed model, we advocate prefixing any derived quantity by the term ‘apparent’.

Notes: An approach to determining the effective elastic moduli of rocks with double porosity is presented. The double-porosity medium is considered to be a heterogeneous material composed of a homogeneous matrix with primary pores and inclusions that represent secondary pores. Fluid flows in the primary-pore system and between primary and secondary pores are neglected because of the low permeability of the primary porosity. The prediction of the effective elastic moduli consists of two steps. Firstly, we calculate the effective elastic properties of the matrix with the primary small-scale pores (matrix homogenization). The porous matrix is then treated as a homogeneous isotropic host in which the large-scale secondary pores are embedded. To calculate the effective elastic moduli at each step, we use the differential effective medium (DEM) approach. The constituents of this composite medium – primary pores and secondary pores – are approximated by ellipsoidal or spheroidal inclusions with corresponding aspect ratios.We have applied this technique in order to compute the effective elastic properties for a model with randomly orientated inclusions (an isotropic medium) and aligned inclusions (a transversely isotropic medium). Using the special tensor basis, the solution of the one-particle problem with transversely isotropic host was obtained in explicit form.The direct application of the DEM method for fluid-saturated pores does not account for fluid displacement in pore systems, and corresponds to a model with isolated pores or the high-frequency range of acoustic waves. For the interconnected secondary pores, we have calculated the elastic moduli for the dry inclusions and then applied Gassmann's tensor relationships. The simulation of the effective elastic characteristic demonstrated that the fluid flow between the connected secondary pores has a significant influence only in porous rocks containing cracks (flattened ellipsoids). For pore shapes that are close to spherical, the relative difference between the elastic velocities determined by the DEM method and by the DEM method with Gassmann's corrections does not exceed 2%. Examples of the calculation of elastic moduli for water-saturated dolomite with both isolated and interconnected secondary pores are presented. The simulations were verified by comparison with published experimental data.

Notes: Shear-wave statics in marine seismic exploration data are routinely too large to be estimated using conventional techniques. Near-surface unconsolidated sediments are often characterized by low values of Vs and steep velocity gradients. Minor variations in sediment properties at these depths correspond to variations in the shear-wave velocity and will produce significant static shifts. It is suggested that a significant proportion of the shear-wave statics solution can be estimated by performing a separate high-resolution survey to target near-surface unconsolidated sediments. Love-wave, shear-wave refraction and geotechnical measurements were individually used to form high-resolution near-surface shear-wave velocity models to estimate the shear-wave statics for a designated survey line. Comparisons with predicted statics revealed that shear-wave statics could not be estimated using a velocity model predicted by substituting geotechnical measurements into empirical relationships. Empirical relationships represent a vast simplification of the factors that control Vs and are therefore not sufficiently sensitive to estimate shear-wave statics. Refraction measurements are potentially sensitive to short-wavelength variations in sediment properties when combined with accurate navigational data. Statics estimated from Love-wave data are less sensitive, and sometimes smoothed in appearance, since interpreted velocity values represent an average both laterally and vertically over the receiver array and the frequency–depth sensitivity range, respectively.For the survey site, statics estimated from near-surface irregularities using shear-wave refraction measurements represent almost half the total statics solution. More often, this proportion will be greater when bedrock relief is less.

Notes: Parsimonious post-stack migration is extended to three dimensions. By tracing single rays back along each incident wave direction (as determined by a local slant stack at the receivers), the ray tracing can be embedded in the migration. This approach significantly reduces the computer time and disk space needed because it is not necessary to build and save image time maps; 3D migration can be performed on a workstation or personal computer rather than using a supercomputer or cluster.The location of a reflector in the output image is defined by tracing a zero-offset ray to the one-way traveltime (the image condition); the orientation of the reflector is defined as a surface perpendicular to the raypath. The migration impulse response operator is confined to the first Fresnel zone around the estimated reflection point, which is much smaller than the large isochronic surface in traditional Kirchhoff depth migration. Additional efficiency is obtained by applying an amplitude threshold to reduce the amount of data to be migrated. Tests on synthetic data show that the proposed implementation of parsimonious 3D post-stack Kirchhoff depth migration is at least two orders of magnitude faster than traditional Kirchhoff migration, at the expense of slightly degraded migration image coherence. The proposed migration is expected to be a useful complement to conventional time migrations for fast initial imaging of subsurface structures and for real-time imaging of near-offset sections during data acquisition for quality control.

Notes: The semi-automated detection of objects has been quite successful in detecting various types of seismic object, such as chimneys. The same technique can be applied successfully to detect faults in 3D seismic data. We show that several different attributes – among others, similarity, frequency and curvature, all of which potentially enhance the visibility of faults – can be combined successfully by an artificial neural network. This results in a fault ‘probability’ cube in which faults are more continuous and noise is suppressed compared with single-attribute cubes. It is believed that the fault-cube can be improved further by applying image-processing techniques to enhance the fault prediction.

Notes: Filters for migrated offset substacks are designed by partial coherence analysis to predict ‘normal’ amplitude variation with offset (AVO) in an anomaly free area. The same prediction filters generate localized prediction errors when applied in an AVO-anomalous interval. These prediction errors are quantitatively related to the AVO gradient anomalies in a background that is related to the minimum AVO anomaly detectable from the data. The prediction-error section is thus used to define a reliability threshold for the identification of AVO anomalies. Coherence analysis also enables quality control of AVO analysis and inversion. For example, predictions that are non-localized and/or do not show structural conformity may indicate spatial variations in amplitude–offset scaling, seismic wavelet or signal-to-noise (S/N) ratio content. Scaling and waveform variations can be identified from inspection of the prediction filters and their frequency responses. S/N ratios can be estimated via multiple coherence analysis.AVO inversion of seismic data is unstable if not constrained. However, the use of a constraint on the estimated parameters has the undesirable effect of introducing biases into the inverted results: an additional bias-correction step is then needed to retrieve unbiased results. An alternative form of AVO inversion that avoids additional corrections is proposed. This inversion is also fast as it inverts only AVO anomalies. A spectral coherence matching technique is employed to transform a zero-offset extrapolation or near-offset substack into P-wave impedance. The same technique is applied to the prediction-error section obtained by means of partial coherence, in order to estimate S-wave velocity to P-wave velocity (VS/VP) ratios. Both techniques assume that accurate well ties, reliable density measurements and P-wave and S-wave velocity logs are available, and that impedance contrasts are not too strong. A full Zoeppritz inversion is required when impedance contrasts that are too high are encountered. An added assumption is made for the inversion to the VS/VP ratio, i.e. the Gassmann fluid-substitution theory is valid within the reservoir area. One synthetic example and one real North Sea in-line survey illustrate the application of the two coherence methods.

Notes: A local plane-wave approach of generalized diffraction tomography in heterogeneous backgrounds, equivalent to Kirchhoff summation techniques when applied in seismic reflection, is re-programmed to act as repeated synthetic aperture radar (SAR) imaging for seismic prestack depth migration. Spotlight-mode SAR imaging quickly provides good images of the electromagnetic reflectivity of the ground via fast Fourier transform (FFT)-based signal processing. By calculating only the Green's functions connecting the aircraft to the centre of the illuminated patch, scattering structures around that centre are also recovered. SAR technology requires us to examine seismic imaging from the local point of view, where the quantity and quality of the available information at each image point are what are important, regardless of the survey geometry. When adapted to seismics, a local image of arbitrary size and sampling is obtained by FFT of seismic energy maps in the scattering wavenumber domain around each node of a pre-calculated grid of Green's functions. These local images can be used to generate a classic prestack depth-migrated section by collecting only their centres. However, the local images also provide valuable information around the centre, as in SAR. They can therefore help to pre-analyse prestack depth migration efficiently, and to perform velocity analysis at a very low cost. The FFT-based signal-processing approach allows local, efficient and automatic control of anti-aliasing, noise and resolution, including optimized Jacobian weights. Repeated local imaging could also be used to speed up migration, with interpolation between local images associated with a coarse grid of Green's functions, as an alternative to interpolation of Green's functions. The local images may, however, show distortions due to the local plane-wave approximation, and the velocity variations across their frame. Such effects, which are not necessarily a problem in SAR, should be controlled and corrected to further enhance seismic imaging. Applications to realistic models and to real data show that, despite the distortion effects, the local images can yield similar information to prestack depth migration, including common-image-point gathers for velocity analyses and AVO/AVA effects, at a much lower cost when a small target is considered.

Notes: Vibroseis is a method that imparts coded seismic energy into the ground. The energy is recorded with geophones and then processed using the known (coded) input signal. The resulting time-domain representation of vibroseis data is an impulsive wavetrain with wavelet properties consistent with the coded input signal convolved with the earth's reflectivity series. Historically, vibratory seismic surveys collect data from one source location at a time, summing one or more sources at each location. We present a method of designing orthogonal sweeps using the concept of combisweeps. The orthogonal sweeps allow simultaneous recording and later separation of two or more unique source locations. Orthogonality of sweeps permits separation of the data into unique source-location field records by a conventional correlation procedure. The separation power of the orthogonal sweeps is demonstrated by a comparison between separated data and data acquired with one vibrator. Separation noise was at a negligible level for our demonstration data sets when two vibrators were located 50 m to 200 m apart. Coincident generation and recording of two vibroseis sweeps at different locations would allow almost double the amount of data to be recorded for a given occupation time and requires only half the storage medium.

Notes: Shales are a major component of sedimentary basins, and they play a decisive role in fluid flow and seismic-wave propagation because of their low permeability and anisotropic microstructure. Shale anisotropy needs to be quantified to obtain reliable information on reservoir fluid, lithology and pore pressure from seismic data, and to understand time-to-depth conversion errors and non-hyperbolic moveout. A single anisotropy parameter, Thomsen's δ parameter, is sufficient to explain the difference between the small-offset normal-moveout velocity and vertical velocity, and to interpret the small-offset AVO response. The sign of this parameter is poorly understood, with both positive and negative values having been reported in the literature. δ is sensitive to the compliance of the contact regions between clay particles and to the degree of disorder in the orientation of clay particles. If the ratio of the normal to shear compliance of the contact regions exceeds a critical value, the presence of these regions acts to increase δ, and a change in the sign of δ, from the negative values characteristic of clay minerals to the positive values commonly reported for shales, may occur. Misalignment of the clay particles can also lead to a positive value of δ. For transverse isotropy, the elastic anisotropy parameters can be written in terms of the coefficients W200 and W400 in an expansion of the clay-particle orientation distribution function in generalized Legendre functions. For a given value of W200, decreasing W400 leads to an increase in δ, while for fixed W400, δ increases with increasing W200. Perfect alignment of clay particles with normals along the symmetry axis corresponds to the maximum values of W200 and W400, given by 〈inlineGraphic alt="inline image" href="urn:x-wiley:00168025:GPR495:GPR_495_mu1" location="equation/GPR_495_mu1.gif"/〉 and 〈inlineGraphic alt="inline image" href="urn:x-wiley:00168025:GPR495:GPR_495_mu2" location="equation/GPR_495_mu2.gif"/〉. A comparison of the predictions of the theory with laboratory measurements shows that most shales lie in a region of the (W200, W400)-plane defined by W400/W200≤Wmax400/Wmax200.

Notes: We introduce a method of wavefield separation from multicomponent data sets based on the use of the continuous wavelet transform. Our method is a further generalization of the approach proposed by Morozov and Smithson, in that by using the continuous wavelet transform, we can achieve a better separation of wave types by designing the filter in the time–frequency domain. Furthermore, using the instantaneous polarization attributes defined in the wavelet domain, we show how to construct filters tailored to separate different wave types (elliptically or linearly polarized), followed by an inverse wavelet transform to obtain the desired wave type in the time domain. Using synthetic and experimental data, we show how the present method can be used for wavefield separation.

Notes: This paper describes the measurements of the acoustic and petrophysical properties of two suites of low-shale sandstone samples from North Sea hydrocarbon reservoirs, under simulated reservoir conditions. The acoustic velocities and quality factors of the samples, saturated with different pore fluids (brine, dead oil and kerosene), were measured at a frequency of about 0.8 MHz and over a range of pressures from 5 MPa to 40 MPa.The compressional-wave velocity is strongly correlated with the shear-wave velocity in this suite of rocks. The ratio VP/VS varies significantly with change of both pore-fluid type and differential pressure, confirming the usefulness of this parameter for seismic monitoring of producing reservoirs.The results of quality factor measurements were compared with predictions from Biot-flow and squirt-flow loss mechanisms. The results suggested that the dominating loss in these samples is due to squirt-flow of fluid between the pores of various geometries. The contribution of the Biot-flow loss mechanism to the total loss is negligible. The compressional-wave quality factor was shown to be inversely correlated with rock permeability, suggesting the possibility of using attenuation as a permeability indicator tool in low-shale, high-porosity sandstone reservoirs.

Notes: We have studied 56 unfractured chalk samples of the Upper Cretaceous Tor Formation of the Dan, South Arne and Gorm Fields, Danish North Sea. The samples have porosities of between 14% and 45% and calcite content of over 95%. The ultrasonic compressional- and shear-wave velocities (VP and VS) for dry and water-saturated samples were measured at up to 75 bar confining hydrostatic pressure corresponding to effective stress in the reservoir. The porosity is the main control of the ultrasonic velocities and therefore of the elastic moduli. The elastic moduli are slightly higher for samples from the South Arne Field than from the Dan Field for identical porosities. This difference may be due to textural differences between the chalk at the two locations because we observe that large grains (i.e. filled microfossils and fossil fragments) that occur more frequently in samples from the Dan Field have a porosity-reducing effect and that samples rich in large grains have a relatively low porosity for a given P-wave modulus. The clay content in the samples is low and is mainly represented by either kaolinite or smectite; samples with smectite have a lower P-wave modulus than samples with kaolinite at equal porosity. We find that ultrasonic VP and VS of dry chalk samples can be satisfactorily estimated with Gassmann's relationships from data for water-saturated samples. A pronounced difference between the VP/VS ratios for dry and water-saturated chalk samples indicates promising results for seismic amplitude-versus-offset analyses.

Notes: Seismic wavefield scattering from a statistically randomly rough interface in a multilayered piecewise homogeneous medium is studied in 3D. The influence of the surface roughness on the scattered wavefield is analysed numerically by using a finite-difference operator in the acoustic domain. Since interface scattering in the real practical sense is a 3D physical phenomenon, we show in this work that the scattering response of a randomly rough interface is not the same in 3D situations as in the 2D cases described in some earlier works. For a given interface roughness height in 3D, an interface roughness height at least three times greater is required to produce an equivalent phase scattering effect in 2D situations, for a given correlation length of the interface roughness scale. Based on observations from spectral analysis, we show that scattering results principally in de-phasing and frequency band-limiting of the incident wavefront, the frequency band-limiting properties being comparable to cases reported in the literature for absorption and thin-layer filtering. The interface scattering phenomenon should be critically considered when using amplitude and phase information from seismic signal during inversion processes.

Notes: Seismic waveforms contain valuable information about the media, but the waveform inversion is a non-linear problem. We present a waveform inversion method that combines a local optimization method with a genetic algorithm to determine the anisotropic parameters of a horizontally stratified medium. Synthetic seismograms for a horizontally stratified anisotropic medium are calculated using the reflectivity technique. In the initial stage of the inversion, the global space-sampling properties of the genetic algorithm are used to direct the search to the region close to the global solution. This solution is then further improved using a conjugate-gradient method. The numerical experiments performed with noisy synthetic data show that our hybrid optimization method satisfactorily reconstructs the anisotropic parameters at a reasonable computing cost while the range of slowness is adequate. We found that (i) for small-angle data neither single- nor multiple-component data are sufficient to determine the anisotropic parameters uniquely; (ii) for medium-angle data the multiple-component data are sufficient to determine the anisotropic parameters exactly whereas the single-component data are not sufficient; and (iii) for wide-angle data, either single- or multiple-component data are sufficient to determine the anisotropic parameters accurately.

Notes: The conventional spectral analysis method for interpretation of magnetic data assumes stationary spatial series and a white-noise source distribution. However, long magnetic profiles may not be stationary in nature and source distributions are not white. Long non-stationary magnetic profiles can be divided into stationary subprofiles following Wiener filter theory. A least-squares inverse method is used to calculate the scaling exponents and depth values of magnetic interfaces from the power spectrum. The applicability of this approach is demonstrated on non-stationary synthetic and field magnetic data collected along the Nagaur–Jhalawar transect, western India. The stationarity of the whole profile and the subprofiles of the synthetic and field data is tested. The variation of the mean and standard deviations of the subprofiles is significantly reduced compared with the whole profile. The depth values found from the synthetic model are in close agreement with the assumed depth values, whereas for the field data these are in close agreement with estimates from seismic, magnetotelluric and gravity data.

Notes: A recently developed laboratory method allows for simultaneous imaging of fluid distribution and measurements of acoustic-wave velocities during flooding experiments. Using a specially developed acoustic sample holder that combines high pressure capacity with good transparency for X-rays, it becomes possible to investigate relationships between velocity and fluid saturation at reservoir stress levels. High-resolution 3D images can be constructed from thin slices of cross-sectional computer-tomography scans (CT scans) covering the entire rock-core volume, and from imaging the distribution of fluid at different saturation levels. The X-ray imaging clearly adds a new dimension to rock-physics measurements; it can be used in the explanation of variations in measured velocities from core-scale heterogeneities. Computer tomography gives a detailed visualization of density regimes in reservoir rocks within a core. This allows an examination of the interior of core samples, revealing inhomogeneities, porosity and fluid distribution. This mapping will not only lead to an explanation of acoustic-velocity measurements; it may also contribute to an increased understanding of the fluid-flow process and gas/liquid mixing mechanisms in rock. Immiscible and miscible flow in core plugs can be mapped simultaneously with acoustic measurements. The effects of core heterogeneity and experimentally introduced effects can be separated, to clarify the validity of measured velocity relationships.

Notes: 4D reservoir monitoring is now becoming established as a normal part of business in some companies, in some regions such as the North Sea, and for some applications, such as thermally enhanced oil recovery. Some spectacular technical successes have been reported, mostly in monitoring water floods. We are beginning to understand how to acquire and process surveys to achieve good results. There is still, however, much to improve in how we use 4D technology to capture its potential value in field management and still some way to go in the tool itself, to make it faster, cheaper and more sensitive to small production effects.

Notes: Recent improvements in the local wavenumber approach have made it possible to estimate both the depth and model type of buried bodies from magnetic data. However, these improvements require calculation of third-order derivatives of the magnetic field, which greatly enhances noise. As a result, the improvements are restricted to data of high quality. We present an alternative method to estimate both the depth and model type using the first-order local wavenumber approach without the need for third-order derivatives of the field. Our method is based on normalization of the first-order local wavenumber anomalies and provides a generalized equation to estimate the depth of some 2D magnetic sources regardless of the source structure. Information about the nature of the sources is obtained after the source location has been estimated. The method was tested using synthetic magnetic anomaly data with random noise and using three field examples.

Notes: Three-dimensional residual moveout analysis is the basic step in velocity model refinement. The analysis is generally carried out using horizontal and/or vertical semblances defined on a sparse set of in-lines or cross-lines with densely sampled source–receiver offsets. An alternative approach, which we call dense residual moveout analysis (DRMA), is to use all the bins of a three-dimensional survey but sparsely sampled offsets. The proposed technique is very fast and provides unbiased and statistically efficient estimates of the residual moveout. Indeed, for the sparsest possible offset distribution, when only near- and far-angle stacks are used, the variance of the residual moveout estimate is only 1.4 times larger than the variance of the least-squares estimate obtained using all offsets.The high performance of DRMA makes it a useful tool for many applications, of which azimuthal velocity analysis is considered here. For a horizontal transverse isotropy (HTI) model, a deterministic procedure is proposed to define, at every point of residual moveout estimation, the azimuthal angle of the HTI axis of symmetry, the Thomsen anisotropy coefficients, and the interval (or root-mean-square) velocities in both the HTI isotropy and symmetry planes. The procedure is not restricted by DRMA assumptions; for example, it is also applicable to semblance-based residual moveout estimates.The high resolution of the technique is illustrated by azimuthal velocity analysis over an oilfield in West Siberia.

Notes: A part of the Békés Basin (an extensional sub-basin of the Pannonian Basin, where the basement under thick Pannonian sediments is well known from deep boreholes and from seismic measurements, and where many magnetotelluric (MT) soundings have been carried out for frequencies ranging from 1 to 10−3 Hz) was selected as a test area to assess the imaging performances of various apparent-resistivity definitions computed with rotational invariants of either the real part of the complex impedance tensor, or its imaginary part, or both. A comparison (based on earlier 3D numerical studies) has been made between the magnetotelluric images obtained in this way and the depths to the high-resistivity basement, as known from boreholes and seismic investigations. The correlation coefficient between the series of basement depth values at 39 MT sites and the apparent-resistivity values was found to be stronger and high correlation appeared at a shorter period when it was computed with apparent resistivities based on the real tensor rather than with apparent resistivities based on the imaginary tensor. In the light of our studies, ρReZ and the impedance phase seem to be more informative than any other combination of magnetotelluric interpretation parameters.

Notes: We use the integral equation for a d.c. electric field, published in the literature, to introduce the concept of the electrical reflectivity tensor into d.c. electric field modelling. It is shown that in d.c. electric field modelling, the electric reflectivity tensor can be obtained in exactly the same way as in electromagnetic modelling. As a result, for a d.c. electric field, the quasi-linear and the quasi-analytical approximations, as well as the quasi-analytical series, can be constructed in exactly the same way as in electromagnetic modelling. If the primary field is uniform, and if the anomalous body is a uniform circular cylinder or a uniform sphere, the reflectivity tensor is zero order (constant), relating to the free surface charge density. Thus, for some homogeneous bodies that have simple shapes and are embedded in a uniform primary field, the electrical reflectivity tensor is not only a mathematical mechanism for obtaining approximate solutions, but also a physical reality. Indeed, the free surface charge density is defined as the change of the electric displacement vector across the boundary surface under consideration. If the primary field is caused by a point source, and if the anomalous body is a uniform sphere, the reflectivity tensor is second order, varying slowly within the sphere. The relationship to the free surface charge density can be established only when both the reflectivity tensor and the free surface charge density are approximated by the first terms of their series solutions. If the point source is far from the centre of the sphere, the corresponding reflectivity tensor reduces to zero order, and is independent of the observation position within the sphere, i.e. it is a constant. Therefore, the basic idea of the quasi-analytical approximation, i.e. taking the reflectivity tensor outside the integral operator, is justified in the case considered here.

Notes: Rock physics analysis plays a vital role in time-lapse seismic interpretation because it provides the link between changes in rock and fluid properties and the resulting seismic data response. In this case study of the Schiehallion Field, we discuss a number of issues that commonly arise in rock physics analyses for time-lapse studies. We show that:〈list xml:id="l1" style="custom"〉1Logarithmic fits of dry bulk (Kdry) and shear (Gdry) moduli vs. effective pressure (Peff) are superior to polynomial fits.22D surface fits of Kdry and Gdry over porosity (φ) and effective pressure using all the core data simultaneously are more useful and accurate than separate 1D fits over φ and Peff for each individual core.3One average set (facies) of Kdry(φ, Peff) and Gdry(φ, Peff) can be chosen to represent adequately the entire Schiehallion reservoir.4Saturated velocities and densities modelled by fluid substitution of Kdry(φ, Peff), Gdry(φ, Peff) and the dry bulk density ρdry(φ) compare favourably with well-log velocities and densities.5P- and S-wave impedance values resulting from fluid substitution of Kdry(φ, Peff), Gdry(φ, Peff) and ρdry(φ) show that the largest impedance changes occur for high porosities and low effective pressures.6Uncertainties in Kdry(φ, Peff) and Gdry(φ, Peff) derived for individual cores can be used to generate error surfaces for these moduli that represent bounds for quantifying uncertainties in seismic modelling or pressure–saturation inversion.

Notes: Legacy streamer data and newer 3D ocean-bottom-cable data are cross-matched and analysed for time-lapse analysis of geomechanical changes due to production in the Valhall Field. The issues relating to time-lapse analysis using two such distinctly different data sets are addressed to provide an optimal cross-matching workflow that includes 3D warping. Additionally an assessment of the differences between the imaging using single-azimuth streamer and multi-azimuth ocean-bottom-cable data is provided. The 3D warping utilized in the cross-matching procedure is sensitive to acquisition and processing differences but is also found to provide valuable insight into the geometrical changes that occur in the subsurface due to production. As such, this work also provides a demonstration of the use of high-resolution 3D interpreted warping to resolve the 3D heterogeneity of the compaction and subsidence. This is an important tool for Valhall, and possibly other fields, where compaction and subsidence (and monitoring thereof) are key factors in the reservoir management since the predominant observed production-induced changes are compaction of the soft, high-porosity chalk reservoir, due to pore-pressure reduction, and the resultant overburden subsidence. Such reservoir compaction could have significant implications for production by changing permeabilities and production rates. Furthermore the subsidence effects could impact upon subsea installations and well-bore stability. Geomechanical studies that have previously been used to model such subsidence and compaction are only constrained by observed surface displacements and measured reservoir pressure changes, with the geological overburden being largely neglected. The approaches suggested herein provide the potential for monitoring and assessment in three dimensions, including the probable heterogeneity and shearing, that is needed for full understanding of reservoir compaction and the resultant effects on the overburden to, for example, mitigate well-bore failures.

Notes: Walkaway vertical seismic profile (VSP) data acquired in basalt-covered areas can be used to improve knowledge of the sub-basalt structure. A synthetic example and a case study from the North Atlantic (UK) show that elastic two-way downward-continuation migration combined with the stationary-phase principle is well suited to the processing of VSP data. Vector data are processed using decoupled elastic migration algorithms in both isotropic and anisotropic media. To illustrate the value of decoupled imaging equations, conventional PP imaging is carried out on the enhanced VSP data and compared with the decoupled scheme. Decoupled vector migration operates directly on the displacement vector, and uses various wave modes. Downgoing waves are migrated to image basalt lava flows and measure their anisotropy. Upgoing waves are used for high-resolution sub-basalt imaging.

Notes: This paper develops an automatic method for interpretation of magnetic data using derivatives of the analytic signal. A linear equation is derived to provide source location parameters of a 2D magnetic body without a priori information about the nature of the source. Then using the source location parameters, the nature of the source can be ascertained. The method has been tested using theoretical simulations with random noise for two 2D magnetic models placed at different depths with respect to the observation height. In both cases, the method gave a good estimate for the location and shape of the sources. Good results were obtained on two field data sets.

Notes: Results of an experimental study of shear- and compressional-wave propagation in an orthorhombically anisotropic medium are presented. The experiments were performed on a physical model consisting of two sets of fractures. The first set consisted of orientated rubber inclusions simulating weak material-filled cracks. The second set consisted of a system of closely spaced parallel fractures simulated by thin plates of epoxy resin, superimposed orthogonally on the first set. Three cases of fracture orientations within the model were identified and studied. Case 1 is analogous to a jointed fracture reservoir with one vertical set of fluid-filled cracks or fractures and one non-filled horizontal set. This case is referred to as JFV. Case 2 is analogous to a double fracture reservoir with one horizontal set of fluid-filled fractures or cracks and one non-filled vertical set. This is referred to as DFH. Case 3 is analogous to a double fracture reservoir with two vertical sets of fractures or cracks, with only one fluid-filled. Case 3 is referred to as DFV. A pulse transmission method was performed on all three modelled cases along the three principal axes. A directional variation in the compressional- and shear-wave velocities, as well as distinct shear-wave splitting, was observed. The elastic constants for each case were determined and differences between them were noted and compared with the controlled results of both layered (transverse isotropy, TI) and vertically fractured (azimuthally anisotropic models, VF) media. The differences in elastic moduli and velocities indicate the potential of recognizing the different fracture orientations and suggest an approach to designing a method of drilling to further enhance oil recovery and reservoir exploitation.

Notes: Most methods for velocity macromodel estimation require considerable operator input, mainly concerning the regularization and the picking of events in the data set or in the migrated images. For both these aspects, slope tomography methods offer interesting solutions. They consider locally coherent events characterized by their slopes in the data cube. Picking is then much easier and consequently denser than in standard traveltime tomography. Stereotomography is the latest slope tomography method. In recent years it has been improved significantly, both from an algorithmic point of view and in terms of practical use. Robust and fast procedures are now available for 2D stereotomographic picking and optimization.Concerning the picking, we propose simple criteria for the selection of relevant data among the automatically picked events. This enables an accurate smooth velocity macromodel to be estimated quite rapidly and with very limited operator intervention. We demonstrate the method using a 2D line extracted from the Oseberg NH8906 data set.

Notes: The theme of the 2003 EAGE/SEG imaging workshop concerned the contrast between different philosophies of ‘model building’: whether an explicit, user-determined model should be imposed throughout the processing, with user updates at each step; or alternatively, whether user intervention should be kept to a minimum so as to avoid preconceived bias, and instead to allow the data itself to guide some heuristic process to converge to an optimal solution.Here we consider a North Sea study where our initial approach was to build the subsurface model using interpreted horizons as a guide to the velocity update. This is common practice in the North Sea, where the geology ‘lends itself’ to a layer-based model representation. In other words, we encourage preconceived bias, as we consider it to be a meaningful geological constraint on the solution.However, in this instance we had a thick chalk sequence, wherein the vertical compaction gradient changed subtly, in a way not readily discernible from the seismic reflection data. As a consequence, imposing the explicit top and bottom chalk horizons, with an intervening vertical compaction gradient (of the form v(x, y, z) =v0(x, y) +k(x, y).z), led to a misrepresentation of the subsurface.To address this issue, a gridded model building approach was also tried. This relied on dense continuous automatic picking of residual moveout in common-reflection point gathers at each iteration of the model update, followed by gridded tomography, resulting in a smoothly varying velocity field which was able to reveal the underlying local changes within the chalk.

Notes: In recent years, the advances in velocity model building and depth imaging have provided a better understanding of complex subsalt plays. The tomographic approach to subsurface velocity modelling, using interpretive processes, has led to significant progress in solving subsalt imaging problems, which were once considered to be impenetrable barriers. We show how gravity data, as an alternative data source, can be integrated into iterative velocity–depth model building to constrain the overburden velocity model and delineate the shape of the salt body above the target reflector. In this way, a structurally accurate image of subsalt reflectors is achieved.

Notes: A tomographic inversion method is presented that uses kinematic information in the form of zero-offset traveltimes and kinematic wavefield attributes (first and second spatial traveltime derivatives) to determine smooth, laterally inhomogeneous 3D subsurface velocity models for depth imaging. The kinematic wavefield attributes can be extracted from the seismic prestack data by means of the common reflection surface (CRS) stack. The input for the tomography is then taken from the resulting attribute volumes at a number of pick locations in the CRS stacked zero-offset volume. As a smooth model description based on B-splines is used and reflection points are treated independently of each other, only locally coherent events in the stacked volume are required and very few picks are needed. Thus, picking is considerably simplified.During the iterative inversion process, the required forward-modelled quantities are obtained by dynamic ray tracing along normal rays pertaining to the input data points. Fréchet derivatives for the tomographic matrix are calculated with ray perturbation theory. The inversion algorithm is demonstrated on a 3D synthetic data example, where the kinematic wavefield attributes have directly been obtained by forward modelling.

Notes: The use of the differential semblance misfit function on common-image-point gathers in the angle domain lends itself to an automated tomographic approach through a gradient-based search in the model space. The velocity model is described by a layer-based model with linear velocity trends and a superimposed bicubic B-spline. The interfaces of the layer-based model are computed by map migration of the PP zero-offset traveltimes of key reflectors. The common-image-point gathers are produced by a restricted inverse generalized Radon transform or amplitude-versus-angle-compensated migration. We present a complete description of all 2.5D formulae for isotropic velocity analysis of PP reflections and the results for ocean-bottom seismic data.

Notes: Seismic migration can be formulated in terms of two consecutive downward extrapolation steps: refocusing the receivers and refocusing the sources. Applying only the first focusing step with an estimate of the focusing operators results in a common focal point (CFP) gather for each depth point at a reflecting boundary. The CFP gathers, in combination with the estimates of the focusing operators, can be used in an iterative procedure to obtain the correct operators. However, current 3D seismic data acquisition geometries do not contain the dense spatial sampling required for calculation of full 3D CFP gathers. We report on the construction of full 3D CFP gathers using a non-full 3D acquisition geometry. The proposed method uses a reflector-orientated data infill procedure based on the azimuthal redundancy of the reflection data. The results on 3D numerical data in this paper show that full 3D CFP gathers, which are kinematically and dynamically correct for the target event, can be obtained. These gathers can be used for iterative updating of the 3D focusing operators.

Notes: We present a migration velocity analysis (MVA) method based on wavefield extrapolation. Similarly to conventional MVA, our method aims at iteratively improving the quality of the migrated image, as measured by the flatness of angle-domain common-image gathers (ADCIGs) over the aperture-angle axis. However, instead of inverting the depth errors measured in ADCIGs using ray-based tomography, we invert ‘image perturbations’ using a linearized wave-equation operator. This operator relates perturbations of the migrated image to perturbations of the migration velocity. We use prestack Stolt residual migration to define the image perturbations that maximize the focusing and flatness of ADCIGs.Our linearized operator relates slowness perturbations to image perturbations, based on a truncation of the Born scattering series to the first-order term. To avoid divergence of the inversion procedure when the velocity perturbations are too large for Born linearization of the wave equation, we do not invert directly the image perturbations obtained by residual migration, but a linearized version of the image perturbations. The linearized image perturbations are computed by a linearized prestack residual migration operator applied to the background image. We use numerical examples to illustrate how the backprojection of the linearized image perturbations, i.e. the gradient of our objective function, is well behaved, even in cases when backprojection of the original image perturbations would mislead the inversion and take it in the wrong direction.We demonstrate with simple synthetic examples that our method converges even when the initial velocity model is far from correct. In a companion paper, we illustrate the full potential of our method for estimating velocity anomalies under complex salt bodies.

Notes: A method is presented to estimate the elastic parameters and thickness of media that are locally laterally homogeneous using P-wave and vertically polarized shear-wave (SV-wave) data. This method is a ‘layer-stripping’ technique, and it uses many aspects of common focal point (CFP) technology. For each layer, a focusing operator is computed using a model of the elastic parameters with which a CFP gather can be constructed using the seismic data. Assuming local homogeneity, the resulting differential time shifts (DTSs) represent error in the model due to anisotropy and error in thickness. In the (τ−p) domain, DTSs are traveltimes Δτ that connect error in layer thickness z, vertical slowness q, and ray parameter p. Series expansion is used to linearize Δτ with respect to error in the elastic parameters and thickness, and least-squares inversion is used to update the model.For stability, joint inversion of P and SV data is employed and, as pure SV data are relatively rare, the use of mode-converted (PSV) data to represent SV in the joint inversion is proposed. Analytic and synthetic examples are used to demonstrate the utility and practicality of this inversion.

Notes: An integrated multiscale seismic imaging flow is applied to dense onshore wide-aperture seismic data recorded in a complex geological setting (thrust belt).An initial P-wave velocity macromodel is first developed by first-arrival traveltime tomography. This model is used as an initial guess for subsequent full-waveform tomography, which leads to greatly improved spatial resolution of the P-wave velocity model. However, the application of full-waveform tomography to the high-frequency part of the source bandwidth is difficult, due to the non-linearity of this kind of method. Moreover, it is computationally expensive at high frequencies since a finite-difference method is used to model the wave propagation. Hence, full-waveform tomography was complemented by asymptotic prestack depth migration to process the full-source bandwidth and develop a sharp image of the short wavelengths. The final traveltime tomography model and two smoothed versions of the final full-waveform tomography model were used as a macromodel for the prestack depth migration.In this study, wide-aperture multifold seismic data are used. After specific preprocessing of the data, 16 frequency components ranging from 5.4 Hz to 20 Hz were inverted in cascade by the full-waveform tomography algorithm. The full-waveform tomography successfully imaged SW-dipping structures previously identified as high-resistivity bodies. The relevance of the full-waveform tomography models is demonstrated locally by comparison with a coincident vertical seismic profiling (VSP) log available on the profile. The prestack depth-migrated images, inferred from the traveltime, and the smoothed full-waveform tomography macromodels are shown to be, on the whole, consistent with the final full-waveform tomography model. A more detailed analysis, based on common-image gather computations, and local comparison with the VSP log revealed that the most accurate migrated sections are those obtained from the full-waveform tomography macromodels. A resolution analysis suggests that the asymptotic prestack depth migration successfully migrated the wide-aperture components of the data, allowing medium wavelengths in addition to the short wavelengths of the structure to be imaged.The processing flow that we applied to dense wide-aperture seismic data is shown to provide a promising approach, complementary to more classical seismic reflection data processing, to quantitative imaging of complex geological structures.

Notes: Three models for the dynamics of seismic airgun-generated bubbles and their associated far-field signals are developed and compared with geophysical data. The first model of an airgun-generated bubble uses a spherical approximation, the second is an approximate Lagrangian model which allows for small deformations from a spherical shape, whilst the final model is an axisymmetric boundary-integral method which permits the bubble to evolve into highly non-spherical geometries. The boundary-integral method also allows both geometric interference and strong dynamic interactions in multi-bubble studies.When comparing the spherical model to experimental data there are three apparent, significant differences: the magnitude of the primary pressure peak, which is greater in the model; the subsequent decay of the pressure peaks and motion – the experimental data demonstrating greater decay and a slower rise rate; and the frequency of oscillation, which is slower in the experimental data. It is believed that the first discrepancy is due to the initial stages of expansion where the compressed air is forced to sparge through the airgun ports. The other differences indicate that there is some other energy-loss mechanism which is not accounted for in the spherical bubble model. Non-spherical bubble behaviour is investigated through the use of two different deformable many-bubble codes and their predictions are compared with the spherical model and experimental data.The Lagrangian model predicts the formation of a buoyancy-driven liquid jet on the first collapse of a typical airgun bubble; however, the model breaks down when the bubble becomes significantly deformed, due to a low-order spherical-harmonic approximation for the potential. The axisymmetric boundary-integral code models the jet shape accurately and it is found that these bubbles evolve to toroidal geometries when the jet impacts on the opposite surface of the bubble. This highly non-spherical behaviour is readily observed on high-speed films of airgun bubbles, and is one key source of energy loss; it damps the pulsations of the bubble and slows its rise speed.Inter-bubble interactions are investigated using the two deformable bubble models, and the predictions are compared to field data. It was found that as the bubbles approach each other, their periods of oscillation increase in accordance with observations, and jets are formed in the direction of motion upon collapse.

Notes: Helicopter electromagnetic (HEM) systems are commonly used for conductivity mapping and the data are often interpreted using an isotropic horizontally layered earth model. However, in regions with distinct dipping stratification, it is useful to extend the model to a layered earth with general anisotropy by assigning each layer a symmetrical 3 × 3 resistivity tensor.The electromagnetic (EM) field is represented by two scalar potentials, which describe the poloidal and toroidal parts of the magnetic field. Via a 2D Fourier transform, we obtain two coupled ordinary differential equations in the vertical coordinate. To stabilize the numerical calculation, the wavenumber domain is divided into two parts associated with small and large wavenumbers. The EM field for small wavenumbers is continued from layer to layer with the continuity conditions. For large wavenumbers, the EM field behaves like a DC field and therefore cannot be sensed by airborne EM systems. Thus, the contribution from the large wavenumbers is simply ignored.The magnetic fields are calculated for the vertical coaxial (VCX), horizontal coplanar (HCP) and vertical coplanar (VCP) coil configurations for a helicopter EM system. The apparent resistivities defined from the VCX, VCP and HCP coil responses, when plotted in polar coordinates, clearly identify the principal anisotropic axes of an anisotropic earth. The field example from the Edwards Aquifer recharge area in Texas confirms that the polar plots of the apparent resistivities identify the principal anisotropic axes that coincide well with the direction of the underground structures.

Notes: The diffusion of electromagnetic fields is dependent not only on conductivity, but also on magnetic permeability, dielectric permittivity and polarizability, i.e. dispersive conductivity. The long-offset transient electromagnetic (LOTEM) method is mainly used to determine the spatial distribution of conductivity in the subsurface. However, earlier work on loop-loop TEM suggests that transient EM methods can also be affected by induced polarization (IP).Numerous 1D forward calculations were carried out to study the IP effect on LOTEM data, using the Cole-Cole relaxation model to simulate the polarizability of the ground. Besides the polarizability of each layer, the IP effect depends on the LOTEM field set-up and the spatial distribution of conductivity in the ground. In particular, near-surface layers with high chargeabilities can significantly distort the late time transients of the electric field components in the vicinity of the transmitter. The influence of polarizable layers on the magnetic field components can be neglected under normal circumstances.In 1997 and 1999, LOTEM measurements were carried out at Mt. Vesuvius in Italy to explore the geological structure of the volcano. Sensitivity studies on the effect of polarizable layers suggest that high chargeabilities in connection with conductive layers at greater depths would result in a detectable distortion of the electric field transients. Although the simultaneous IP measurements revealed high chargeabilities in a near-surface layer, no evidence of IP effects could be found in the measured LOTEM data. We conclude that the observed chargeabilities are local and that 3D effects are probably present in the data.Another aspect is the measurement of the system response, which is usually measured by placing a receiver very close to the transmitter. Therefore, large distortions can be expected if near-surface polarizable layers exist. This was verified in practice by field measurements in an area with high chargeabilities in Longerich, Cologne.

Notes: The shallowest few hundred metres of the earth cannot be adequately imaged by conventional seismic when tuned for deeper targets. Adding independent measurements (such as uphole or shallow refraction surveys) reduces this information gap, but in some arid areas (such as Saudi Arabia) the near-surface complexities are not well resolved, even in this way. The joint tomographic inversion of different wave types can contribute to reducing these ambiguities further, by complementing the different penetration ranges and propagation directions of reflected, refracted and diving waves. Here, we demonstrate the weakness of diving waves when used alone, and the value of complementing them by available reflected and refracted arrivals.

Notes: Effective communication between seismic specialists should be facilitated by a shared process model that can be used at different levels of abstraction. In this shared model the seismic work-flow is presented as a value chain, showing the complex interrelationships between the broad range of specialized activities that are needed in today's practice. One of these activities is velocity estimation, providing the relationship between seismic time and geological depth. Excellence in the seismic value chain will depend on the quality of the specialized tools and skills (abilities) involved, as well as on the capability of the organization to combine these abilities in an integrated work-flow to realize maximum value at the end of the chain.

Notes: Seismic data volumes, which require huge transmission capacities and massive storage media, continue to increase rapidly due to acquisition of 3D and 4D multiple streamer surveys, multicomponent data sets, reprocessing of prestack seismic data, calculation of post-stack seismic data attributes, etc. We consider lossy compression as an important tool for efficient handling of large seismic data sets. We present a 2D lossy seismic data compression algorithm, based on sub-band coding, and we focus on adaptation and optimization of the method for common-offset gathers. The sub-band coding algorithm consists of five stages: first, a preprocessing phase using an automatic gain control to decrease the non-stationary behaviour of seismic data; second, a decorrelation stage using a uniform analysis filter bank to concentrate the energy of seismic data into a minimum number of sub-bands; third, an iterative classification algorithm, based on an estimation of variances of blocks of sub-band samples, to classify the sub-band samples into a fixed number of classes with approximately the same statistics; fourth, a quantization step using a uniform scalar quantizer, which gives an approximation of the sub-band samples to allow for high compression ratios; and fifth, an entropy coding stage using a fixed number of arithmetic encoders matched to the corresponding statistics of the classified and quantized sub-band samples to achieve compression. Decompression basically performs the opposite operations in reverse order. We compare the proposed algorithm with three other seismic data compression algorithms. The high performance of our optimized sub-band coding method is supported by objective and subjective results.

Notes: The estimation of velocity models is still crucial in seismic reflection imaging as it controls the quality of the depth-migrated image, which is the basis of geological interpretation. Among the numerous existing methods for velocity determination, tomographic methods are very attractive for their efficiency and ability to retrieve heterogeneities of the medium. We present three tomographic methods in order to estimate heterogeneous velocity models from 2D prestack PP reflection data: a traveltime tomography in the time-migrated domain, a traveltime and slope tomography in the non-migrated time domain, and a slope tomography in the depth-migrated domain.The first method (traveltime tomography in the time domain) is based on continuous picked events, whereas the two slope tomographic methods, one in the time domain and the other in the depth domain, are based on locally coherent events, with no assumptions about reflector geometry or the unknown velocity field. The purpose of this paper is not to describe in detail the theoretical basis and implementation of the methods, but to apply and compare their output using the same marine real data set. Based on the estimated velocity models, the migrated images and the common-image gathers from the three processing routes, the relative strengths and weaknesses of the methods are discussed. Finally, similarities are indicated and potential alternative approaches are proposed.

Notes: Skin depth is an electromagnetic (EM) scale length that provides a measure of the degree of attenuation experienced by a particular frequency of an EM system. As has been discussed in the literature, skin depth is not a complete measure of the depth of investigation, but the two may be related. Frequency-domain airborne EM systems employ pairs of transmitter and receiver coils that use a frequency range from several hundred hertz to over 100 kHz. For elevated dipoles, both geometrical and frequency-dependent attenuation of the induced fields must be considered. For airborne EM systems it is possible to define a skin depth based only on the electric field induced by the transmitter. A vertical decay scale length, here defined from the at-surface position of maximum electric field, enables the same skin-depth estimate to be obtained for both cases of vertical and horizontal dipolar excitation. Such dipolar skin depths associated with towed-bird and fixed-wing airborne systems are studied in relation to frequency, conductivity and sensor elevation. Dipolar skin depths are found to be much smaller than their plane-wave counterparts except at high frequency (〉50 kHz) and in combination with high conductivity. For the majority of airborne systems the influence of altitude on skin depth is highly significant. Dipolar skin depths increase with increasing sensor elevation. Low frequencies display the greatest sensitivity. At low elevation (〈40 m), geometrical attenuation dominates the behaviour of the skin depth. The study indicates that typical low-altitude airborne surveys provide vertically compact assessments of subsurface conductivity, well suited to near-surface, environmental applications.

Notes: The image-wave equation for depth remigration is a partial differential equation that is similar to the acoustic wave equation. In this work, we study its finite-difference solution and possible applications. The conditions for stability, dispersion and dissipation exhibit a strong wavenumber dependence. Where higher horizontal than vertical wavenumbers are present in the data to be remigrated, stability may be difficult to achieve. Grid dispersion and dissipation can only be reduced to acceptable levels by the choice of very small grid intervals. Numerical tests demonstrate that, upon reaching the true medium velocity, remigrated images of curved reflectors propagate to the correct depth and those of diffractions collapse to single points. The latter property points towards the method's potential for use as a tool for migration velocity analysis. A first application to inhomogeneous media shows that in a horizontally layered medium, the reflector images reach their true depth when the remigration velocity equals the inverse of the mean medium slowness.

Notes: Wave-equation redatuming can be a very efficient method of overcoming the overburden imprint on the target area. Owing to the growing amount of 3D data, it is increasingly important to develop a feasible method for the redatuming of 3D prestack data.Common 3D acquisition designs produce relatively sparse data sets, which cannot be redatumed successfully by applying conventional wave-equation redatuming. We propose a redatuming approach that can be used to perform wave-equation redatuming of sparse 3D data. In this new approach, additional information about the medium velocity below the new datum is included, i.e. redatumed root-mean-square (RMS) velocities, which can be extracted from the input data set by conventional velocity analysis, are used. Inclusion of this additional information has the following implications: (i) it becomes possible to simplify the 4D redatuming integral into a 2D integral such that the number of traces needed to calculate one output time sample and the computational effort are both reduced; (ii) the information about the subsurface enables an infill of traces which are needed for the integral calculation but which are missing in the sparse input data set.Two tests applying this new approach to fully sampled 2D data show satisfactory results, implying that this method can certainly be used for the redatuming of sparse 3D data sets.

Notes: Subsalt imaging is strongly dependent on the quality of the velocity model. However, rugose salt bodies complicate wavefield propagation and lead to subsalt multipathing, illumination gaps and shadow zones, which cannot be handled correctly by conventional traveltime-based migration velocity analysis (MVA). We overcome these limitations by the wave-equation MVA technique, introduced in a companion paper, and demonstrate the methodology on a realistic synthetic data set simulating a salt-dome environment and a Gulf of Mexico data set. We model subsalt propagation using wave paths created by one-way wavefield extrapolation. Those wave paths are much more accurate and robust than broadband rays, since they inherit the frequency dependence and multipathing of the underlying wavefield. We formulate an objective function for optimization in the image space by relating an image perturbation to a perturbation of the velocity model. The image perturbations are defined using linearized prestack residual migration, thus ensuring stability, relative to the first-order Born approximation assumptions. Synthetic and real data examples demonstrate that wave-equation MVA is an effective tool for subsalt velocity analysis, even when shadows and illumination gaps are present.

Notes: Angle-domain common-image gathers (ADCIGs) are an essential tool for migration velocity analysis (MVA). We present a method for computing ADCIGs in 3D from the results of wavefield-continuation migration. The proposed methodology can be applied before or after the imaging step in a migration procedure. When computed before imaging, 3D ADCIGs are functions of the offset ray parameters (p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉x〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-1.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-1"/〉, p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉y〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-3.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-3"/〉); we derive the geometric relationship that links the offset ray parameters to the aperture angle γ and the reflection azimuth φ. When computed after imaging, 3D ADCIGs are directly produced as functions of γ and φ.The mapping of the offset ray parameters (p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉x〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-5.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-5"/〉, p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉y〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-7.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-7"/〉) into the angles (γ, φ) depends on both the local dips and the local interval velocity; therefore, the transformation of ADCIGs computed before imaging into ADCIGs that are functions of the actual angles is difficult in complex structure. By contrast, the computation of ADCIGs after imaging is efficient and accurate even in the presence of complex structure and a heterogeneous velocity function. On the other hand, the estimation of the offset ray parameters (p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉x〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-9.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-9"/〉, p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉y〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-11.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-11"/〉) is less sensitive to velocity errors than the estimation of the angles (γ, φ). When ADCIGs that are functions of the offset ray parameters (p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉x〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-13.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-13"/〉, p〈mml:math〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mrow〉〈mml:mtext fontstyle="italic"〉y〈/mml:mtext〉〈mml:msub〉〈mml:mtext〉 〈/mml:mtext〉〈mml:mtext〉h〈/mml:mtext〉〈/mml:msub〉〈/mml:mrow〉〈/mml:msub〉〈/mml:math〉〈inlineGraphic alt="math image" location="equation/tex2gif-sub-15.gif" href="urn:x-wiley:00168025:GPR444:tex2gif-sub-15"/〉) are adequate for the application of interest (e.g. ray-based tomography), the computation of ADCIGs before imaging might be preferable.Errors in the migration velocity cause the image point in the angle domain to shift along the normal to the apparent geological dip. By assuming stationary rays (i.e. small velocity errors), we derive a quantitative relationship between this normal shift and the traveltime perturbation caused by velocity errors. This relationship can be directly used in an MVA procedure to invert depth errors measured from ADCIGs into migration velocity updates. In this paper, we use it to derive an approximate 3D residual moveout (RMO) function for measuring inconsistencies between the migrated images at different γ and φ. We tested the accuracy of our kinematic analysis on a 3D synthetic data set with steeply dipping reflectors and a vertically varying propagation velocity. The tests confirm the accuracy of our analysis and illustrate the limitations of the straight-ray approximation underlying our derivation of the 3D RMO function.

Notes: The common focal point (CFP) method and the common reflection surface (CRS) stack method are compared. The CRS method is a fast, highly automated procedure that provides high S/N ratio simulation of zero-offset (ZO) images by combining, per image point, the reflection energy of an arc segment that is tangential to the reflector. It uses smooth parametrized two-way stacking operators, based on a data-driven triplet of attributes in 2D (eight parameters in 3D). As a spin-off, the attributes can be used for several applications, such as the determination of the geometrical spreading factor, multiple prediction, and tomographic inversion into a smooth background velocity model. The CFP method aims at decomposing two-way seismic reflection data into two full-aperture one-way propagation operators. By applying an iterative updating procedure in a half-migrated domain, it provides non-smooth focusing operators for prestack imaging using only the energy from one focal point at the reflector. The data-driven operators inhibit all propagation effects of the overburden. The CFP method provides several spin-offs, amongst which is the CFP matrix related to one focal point, which displays the reflection amplitudes as measured at the surface for each source–receiver pair. The CFP matrix can be used to determine the specular reflection source–receiver pairs and the Fresnel zone at the surface for reflection in one single focal point. Other spin-offs are the prediction of internal multiples, the determination of reflectivity effects, velocity-independent redatuming and tomographic inversion to obtain a velocity–depth model. The CFP method is less fast and less automated than the CRS method. From a pointwise comparison of features it is concluded that one method is not a subset of the other, but that both methods can be regarded as being to some extent complementary.

Notes: Numerical simulations are used to compare the resolution and efficiency of 2D resistivity imaging surveys for 10 electrode arrays. The arrays analysed include pole-pole (PP), pole-dipole (PD), half-Wenner (HW), Wenner-α (WN), Schlumberger (SC), dipole-dipole (DD), Wenner-β (WB), γ-array (GM), multiple or moving gradient array (GD) and midpoint-potential-referred measurement (MPR) arrays. Five synthetic geological models, simulating a buried channel, a narrow conductive dike, a narrow resistive dike, dipping blocks and covered waste ponds, were used to examine the surveying efficiency (anomaly effects, signal-to-noise ratios) and the imaging capabilities of these arrays. The responses to variations in the data density and noise sensitivities of these electrode configurations were also investigated using robust (L1-norm) inversion and smoothness-constrained least-squares (L2-norm) inversion for the five synthetic models.The results show the following. (i) GM and WN are less contaminated by noise than the other electrode arrays. (ii) The relative anomaly effects for the different arrays vary with the geological models. However, the relatively high anomaly effects of PP, GM and WB surveys do not always give a high-resolution image. PD, DD and GD can yield better resolution images than GM, PP, WN and WB, although they are more susceptible to noise contamination. SC is also a strong candidate but is expected to give more edge effects. (iii) The imaging quality of these arrays is relatively robust with respect to reductions in the data density of a multi-electrode layout within the tested ranges. (iv) The robust inversion generally gives better imaging results than the L2-norm inversion, especially with noisy data, except for the dipping block structure presented here. (v) GD and MPR are well suited to multichannel surveying and GD may produce images that are comparable to those obtained with DD and PD. Accordingly, the GD, PD, DD and SC arrays are strongly recommended for 2D resistivity imaging, where the final choice will be determined by the expected geology, the purpose of the survey and logistical considerations.

Notes: Many prospective passive ocean margins are covered by large areas of basalts. These basalts are often extremely heterogeneous and scatter the seismic energy of the conventional seismic reflection system so that it becomes difficult to obtain information on deeper reflectors. Since high frequencies are scattered more than low frequencies, we argue that the acquisition system for sub-basalt targets should be modified to emphasize the low frequencies, using much larger airguns, and towing the source and receivers at about 20 m depth. In the summer of 2001 we obtained seismic reflection data over basalt in the northeast Atlantic using a system modified to enhance the low-frequency energy. These new data show deep reflections that are not visible on lines shot in the same places with a conventional system.

Notes: We present a new method for producing a ‘brute’ velocity image rapidly and automatically from traveltimes picked from densely sampled refraction data. The procedure involves imaging by data transformation from the time–offset domain into the tau–p (intercept–slope) domain, and does not include conventional modelling steps. Differences in apparent velocity and tau along reciprocal paths in the up- and downdip directions allow the estimation of the true velocity and geometrical position of the ray turning points. The tau–velocity–turningpoint (τ–ν–x) map distributes phases automatically on the basis of geometry and velocity to give a two-dimensional representation of subsurface structure. This map may be converted simply to depth and two-way-time images. Such images have potential for direct geological interpretation, for use as a starting model for seismic inversion, for superimposition on to conventional reflection images, or for input into prestack depth migration and other processing routines.

Notes: Reflection and refraction data are normally processed with tools designed to deal specifically with either near- or far-offset data. Furthermore, the refraction data normally require the picking of traveltimes prior to analysis. Here, an automatic processing algorithm has been developed to analyse wide-angle multichannel streamer data without resorting to manual picking or traveltime tomography. Time–offset gathers are transformed to the tau–p domain and the resulting wavefield is downward continued to the depth–p domain from which a velocity model and stacked section are obtained. The algorithm inputs common-depth-point (CDP) gathers and produces a depth-converted stacked section that includes velocity information. The inclusion of long-offset multichannel streamer data within the tau–p transformation enhances the signal from high-velocity refracted basalt arrivals. Downward continuation of the tau–p transformed wavefield to the depth–p domain allows the reflection and refraction components of the wavefield to be treated simultaneously. The high-slowness depth–p wavefield provides the velocity model and the low-slowness depth–p wavefield may be stacked to give structural information. The method is applied to data from the Faeroe Basin from which sub-basalt velocity images are obtained that correlate with an independently derived P-wave model from the line.

Notes: We demonstrate that the use of long-offset seismic data allows wide-angle reflections and diving waves to be recorded, and that these can be used in conjunction with prestack depth migrations to constrain and to image the base of the basalt flows and the underlying structure in the Faeroe-Shetland Basin. Crustal velocity models are built first by inverting the traveltimes of the recorded reflections and diving waves using ray-tracing methods. Finer details of the velocity structure can then be refined by analysis of the amplitudes and waveforms of the arrivals. We show that prestack depth migration of selected wide-angle arrivals of known origin, such as the base-basalt reflection, using the crustal velocity model, allows us to build a composite image of the structure down to the pre-rift basement. This has the advantage that the wide-angle first-arriving energy must be primary, and not from one of the many multiples or mode-converted phases that plague near-offset seismic data. This allows us to ‘tag’ these primary arrivals with confidence and then to identify the same arrivals on higher-resolution prestack migrations that include data from all offsets. Examples are drawn from the Faeroe-Shetland Basin, with a series of regional maps of the entire area showing the basalt depths and the thickness of the basalt flows and underlying sediment down to the top of the pre-rift basement. The maps show how the basalts thin to the southeast away from their presumed source west of the present Faeroe Islands, and also show the extent to which the structure of the pre-rift basement controls the considerable variations in sediment thickness between the basement and the cap formed by the overlying basalt flows.

Notes: The idea of imaging beneath a high-velocity layer using converted waves has been popular since 1990. Because these wave types have their maximum amplitudes at mid- to far-offsets, the search for pure P-waves at the highly multiple-contaminated near-offsets can be avoided. For the Atlantic Margin, with buried thin-layered basalts, our quantitative study shows that the initial single-layered approach is not viable. Even in an unrealistic ideal geological setting, the amplitude of the symmetrical PSP-mode is far too weak to be recognized on towed streamer data. Furthermore, in the far-offset window, where locally converted waves have their strongest amplitudes, there is a multitude of other reflections, refractions and interbedded multiples, which have similar moveouts and, often, higher amplitudes. Without the removal of these events, a reliable image of the subsurface cannot be produced. We show that even if this problem were solved, it would be far easier to use the P-wave reflection from beneath the basalt at near-offsets. Our study shows that this wave type is by far the strongest response. A borehole-derived model using a thin-layered basalt sequence reveals that the strongest locally converted wave has an asymmetrical path and is 10 times weaker. All our results indicate that the pure P-modes provide the best chance of imaging sub-basalt sedimentary interfaces.

Notes: In this study, a set of 20 2D seismic lines, acquired over the Golan Heights basaltic plateau, was processed and analysed. Although the data were acquired and processed by standard techniques, in some cases good-quality seismic images were obtained under several hundred metres of basalts. We describe how the seismic characteristics of the top basalt layer were defined and show the effect of the numerous widespread volcanic sources on the quality of the final images. The new data reveal the first images of the sedimentary sequence under the basalt flows, and indicate that strands of the Dead Sea Transform extend into this area. The entire region was found to be very deformed. Several attractive traps for hydrocarbon exploration were also identified on the output sections.

Notes: Artificial neural networks were used to implement an automatic inversion of frequency-domain airborne electromagnetic (AEM) data that do not require a priori information about the survey area. Two classes of model, i.e. homogeneous half-space models and horizontally layered half-space models with two layers, are used in this 1D inversion, and for each data point the selection of the class of 1D model is performed prior to the inversion, also using an artificial neural network. The proposed inversion method was tested in a survey area situated in Austria, northwest of Vienna in the Bohemian Massif. The results of the inversion were compared with the geological setting, logging results, and seismic and gravimetric measurements. This comparison shows a good correlation between the AEM models and the known geological and geophysical data.

Notes: Wide-angle multicomponent ocean-bottom cable (OBC) data should further enhance sub-basalt imaging by using both compressional and converted shear wavefields. The first step in analysing multicomponent OBC data is to decompose the recorded wavefields into pure P- and pure S-wavefields, and extract the upgoing P- and S-waves. This paper presents a new scheme to separate P- and S-wavefields from wide-angle multicomponent OBC data in the τ–p domain. By considering plane-wave components with a known horizontal slowness, the P- and S-wavefields are separated into the directions of observed P- and S-wave oscillations using the horizontal and vertical components of the data. The upgoing P- and S-waves are then extracted from the separated P- and S-wavefields. The parameters used in the separation are the seismic wave velocities and the density at the receiver location, which can be estimated from the first reflection phase observed on the horizontal and vertical components. Numerical tests on synthetic data for a plane-layered model show good performance and demonstrate the accuracy of the scheme. Separation of wavefields from a basalt model is performed using synthetic wide-angle multicomponent OBC data. The results show that both near-offset and wide-angle reflections and conversions from within and below basalt layers are enhanced and clearly identified on the separated wavefields.

Notes: This paper is concerned with the problem of interpretation of anomalous seismic amplitudes, induced by the amplitude-scattering phenomenon. This phenomenon occurs in the vicinity of a crack distribution at the interface between elastic layers. The purpose of this work is to obtain a better understanding of the physics of this distinctive phenomenon, in order to interpret correctly the amplitudes of the reflected events. By analogy with studies in optics and in acoustics, we suggest that diffraction is widely involved in the amplitude-scattering phenomenon. Analytical evaluation of the amount of energy carried by the reflected and the diffracted waves shows that neglecting diffraction in numerical models leads to local underestimation of the amplitude of waves reflected at interfaces with gas-filled crack distribution.

Notes: Fluid flow in many hydrocarbon reservoirs is controlled by aligned fractures which make the medium anisotropic on the scale of seismic wavelength. Applying the linear-slip theory, we investigate seismic signatures of the effective medium produced by a single set of ‘general’ vertical fractures embedded in a purely isotropic host rock. The generality of our fracture model means the allowance for coupling between the normal (to the fracture plane) stress and the tangential jump in displacement (and vice versa). Despite its low (triclinic) symmetry, the medium is described by just nine independent effective parameters and possesses several distinct features which help to identify the physical model and estimate the fracture compliances and background velocities. For example, the polarization vector of the vertically propagating fast shear wave S1 and the semi-major axis of the S1-wave normal-moveout (NMO) ellipse from a horizontal reflector always point in the direction of the fracture strike. Moreover, for the S1-wave both the vertical velocity and the NMO velocity along the fractures are equal to the shear-wave velocity in the host rock.Analysis of seismic signatures in the limit of small fracture weaknesses allows us to select the input data needed for unambiguous fracture characterization. The fracture and background parameters can be estimated using the NMO ellipses from horizontal reflectors and vertical velocities of P-waves and two split S-waves, combined with a portion of the P-wave slowness surface reconstructed from multi-azimuth walkaway vertical seismic profiling (VSP) data. The stability of the parameter-estimation procedure is verified by performing non-linear inversion based on the exact equations.

Notes: We obtain the wave velocities and quality factors of gas-hydrate-bearing sediments as a function of pore pressure, temperature, frequency and partial saturation. The model is based on a Biot-type three-phase theory that considers the existence of two solids (grains and gas hydrate) and a fluid mixture. Attenuation is described with the constant-Q model and viscodynamic functions to model the high-frequency behaviour. We apply a uniform gas/water mixing law that satisfies Wood's and Voigt's averages at low and high frequencies, respectively. The acoustic model is calibrated to agree with the patchy-saturation theory at high frequencies (White's model). Pressure effects are accounted by using an effective stress law for the dry-rock moduli and permeabilities. The dry-rock moduli of the sediment are calibrated with data from the Cascadia margin. Moreover, we calculate the depth of the bottom simulating reflector (BSR) below the sea floor as a function of sea-floor depth, geothermal gradient below the sea floor, and temperature at the sea floor.

Notes: This paper presents a traveltime inversion approach, using the reflection traveltimes from offset VSP data, to reconstruct the horizontal and vertical velocities for stratified anisotropic media. The inverse problem is reduced to a set of linear equations, and solved by the singular value decomposition (SVD) technique. The validity of this inversion scheme is verified using two sets of synthetic data simulated using a finite-difference method, one for an isotropic model and the other for an elliptically anisotropic model. The inversion result demonstrates that our anisotropic velocity inversion scheme may be applied to both isotropic and anisotropic media. The method is finally applied to a real offset VSP data set, acquired in an oilfield in northwestern China.