# An introduction to geophysical exploration

An Introduction to Geophysical Exploration Philip Kearey Department of Earth Sciences University of Bristol

Michael Brooks Ty Newydd, City Near Cowbridge Vale of Glamorgan

An Introduction to Geophysical Exploration Philip Kearey Department of Earth Sciences University of Bristol

Michael Brooks Ty Newydd, City Near Cowbridge Vale of Glamorgan

© 2002 by Blackwell Science Ltd Editorial Of?ces: Osney Mead, Oxford OX2 0EL 25 John Street, LondonWC1N 2BS 23 Ainslie Place, Edinburgh EH3 6AJ 350 Main Street, Malden MA 02148-5018, USA 54 University Street, Carlton Victoria 3053, Australia 10, rue Casimir Delavigne 75006 Paris, France

Other Editorial Of?ces: BlackwellWissenschafts-Verlag GmbH Kurfürstendamm 57 10707 Berlin, Germany Blackwell Science KK MG Kodenmacho Building 7–10 Kodenmacho Nihombashi Chuo-ku,Tokyo 104, Japan Iowa State University Press A Blackwell Science Company 2121 S. State Avenue Ames, Iowa 50014-8300, USA

First published 1984 Reprinted 1987, 1989 Second edition 1991 Reprinted 1992, 1993, 1994, 1995, 1996 1998, 1999, 2000 Third edition 2002

Set by SNP Best-setTypesetter Ltd., Hong Kong Printed and bound in Great Britain by TJ International, Padstow, Cornwall.

The right of the Authors to be identi?edastheAuthorsofthisWork has been asserted in accordance with the Copyright, Designs and Patents Act 1988.

All rights reserved. No part of this publication may be reproduced, stored in a retrieval system, or transmitted, in any form or by any means, electronic, mechanical, photocopying, recording or otherwise, except as permitted by the UK Copyright, Designs and Patents Act 1988, without the prior permission of the copyright owner.

A catalogue record for this title is available from the British Library ISBN 0-632-04929-4

Library of Congress Cataloging-in-Publication Data has been applied for

The Blackwell Science logo is a trade mark of Blackwell Science Ltd, registered at the United KingdomTrade Marks Registry DISTRIBUTORS Marston Book Services Ltd PO Box 269 Abingdon, Oxon OX14 4YN (Orders: Tel: 01235 465500 Fax: 01235 465555) The Americas Blackwell Publishing c/o AIDC PO Box 20 50Winter Sport Lane Williston,VT 05495-0020 (Orders: Tel: 800 216 2522 Fax: 802 864 7626) Australia Blackwell Science Pty Ltd 54 University Street Carlton,Victoria 3053 (Orders: Tel: 3 9347 0300 Fax: 3 9347 5001)

Contents Preface, ix 1 The principles and limitations of geophysical exploration methods, 1 1.1 Introduction, 1 3.7 1.2 The survey methods, 1 3.8 1.3 The problem of ambiguity in geophysical interpretation, 6 1.4 The structure of the book, 7

7.8 Aeromagnetic and marine surveys, 164 9.1 7.9 Reduction of magnetic observations, 165 9.2 7.9.1 Diurnal variation correction, 165 7.9.2 Geomagnetic correction, 166 9.3 7.9.3 Elevation and terrain corrections, 166 9.4 7.10 Interpretation of magnetic anomalies, 166 7.10.1 Introduction, 166 7.10.2 Direct interpretation, 168 7.10.3 Indirect interpretation, 170 7.11 Potential ?eld transformations, 172 9.5 7.12 Applications of magnetic surveying, 173 9.6 Problems, 180 9.7 Further reading, 181 9.8 8 Electrical surveying, 183 8.1 Introduction, 183 9.9 8.2 Resistivity method, 183 9.10 8.2.1 Introduction, 183 9.11 8.2.2 Resistivities of rocks and minerals, 183 8.2.3 Current ?ow in the ground, 184 8.2.4 Electrode spreads, 186 8.2.5 Resistivity surveying equipment, 186 9.12 8.2.6 Interpretation of resistivity data, 187 9.13 8.2.7 Vertical electrical sounding interpretation, 188 8.2.8 Constant separation traversing interpretation, 193 10 8.2.9 Limitations of the resistivity method, 196 10.1 8.2.10 Applications of resistivity surveying, 196 10.2 8.3 Induced polarization (IP) method, 199 10.3 8.3.1 Principles, 199 10.4 8.3.2 Mechanisms of induced polarization, 199 8.3.3 Induced polarization measurements, 200 8.3.4 Field operations, 201 8.3.5 Interpretation of induced polarization data, 201 10.5 8.3.6 Applications of induced polarization 10.6 surveying, 202 8.4 Self-potential (SP) method, 203 8.4.1 Introduction, 203 11 8.4.2 Mechanism of self-potential, 203 11.1 8.4.3 Self-potential equipment and survey 11.2 procedure, 203 11.3 8.4.4 Interpretation of self-potential 11.4 anomalies, 204 Problems, 205 Further reading, 207 Contents vii

Preface This book provides a general introduction to the most important methods of geophysical exploration. These methods represent a primary tool for investigation of the subsurface and are applicable to a very wide range of problems. Although their main application is in prospecting for natural resources, the methods are also used, for example, as an aid to geological surveying, as a means of deriving information on the Earth’s internal physical properties, and in engineering or archaeological site investigations. Consequently, geophysical explo- ration is of importance not only to geophysicists but also The book covers the physical principles, methodology, interpretational procedures and ?elds of application of the various survey methods.The main emphasis has been placed on seismic methods because these represent the most extensively used techniques, being routinely and widely employed by the oil industry in prospecting for hydrocarbons. Since this is an introductory text we have not attempted to be completely comprehensive in our coverage of the subject. Readers seeking further infor- mation on any of the survey methods described should refer to the more advanced texts listed at the end of each We hope that the book will serve as an introductory course text for students in the above-mentioned disci- plines and also as a useful guide for specialists who wish to be aware of the value of geophysical surveying to their own disciplines. In preparing a book for such a wide possible readership it is inevitable that problems arise concerning the level of mathematical treatment to be adopted. Geophysics is a highly mathematical subject and, although we have attempted to show that no great mathematical expertise is necessary for a broad under- standing of geophysical surveying, a full appreciation of the more advanced data processing and interpretational techniques does require a reasonable mathematical abil- ity. Our approach to this problem has been to keep the mathematics as simple as possible and to restrict full mathematical analysis to relatively simple cases.We con- sider it important, however, that any user of geophysical surveying should be aware of the more advanced tech- niques of analysing and interpreting geophysical data since these can greatly increase the amount of useful information obtained from the data. In discussing such techniques we have adopted a semiquantitative or quali- tative approach which allows the reader to assess their scope and importance, without going into the details of Earlier editions of this book have come to be accepted as the standard geophysical exploration textbook by nu- merous higher educational institutions in Britain, North America, and many other countries. In the third edition, we have brought the content up to date by taking account of recent developments in all the main areas of geophysical exploration.We have extended the scope of the seismic chapters by including new material on three-component and 4D re?ection seismology, and by providing a new section on seismic tomography. We have also widened the range of applications of refrac- tion seismology considered, to include an account of engineering site investigation.

geophysical exploration methods 1.1 Introduction This chapter is provided for readers with no prior knowledge of geophysical exploration methods and is pitched at an elementary level. It may be passed over by readers already familiar with the basic principles and The science of geophysics applies the principles of physics to the study of the Earth. Geophysical investiga- tions of the interior of the Earth involve taking measure- ments at or near the Earth’s surface that are in?uenced by the internal distribution of physical properties. Analysis of these measurements can reveal how the physical properties of the Earth’s interior vary vertically and By working at different scales, geophysical methods may be applied to a wide range of investigations from studies of the entire Earth (global geophysics; e.g. Kearey & Vine 1996) to exploration of a localized region of Vogelsang 1995, McCann et al. 1997). In the geophysical exploration methods (also referred to as geophysical sur- veying) discussed in this book, measurements within geographically restricted areas are used to determine the distributions of physical properties at depths that re?ect An alternative method of investigating subsurface geology is, of course, by drilling boreholes, but these are expensive and provide information only at discrete locations. Geophysical surveying, although sometimes prone to major ambiguities or uncertainties of interpre- tation, provides a relatively rapid and cost-effective means of deriving areally distributed information on subsurface geology. In the exploration for subsurface resources the methods are capable of detecting and delineating local features of potential interest that could Geophysical surveying does not dispense with the need for drilling but, properly applied, it can optimize explo- ration programmes by maximizing the rate of ground coverage and minimizing the drilling requirement. The importance of geophysical exploration as a means of deriving subsurface geological information is so great that the basic principles and scope of the methods and their main ?elds of application should be appreciated by any practising Earth scientist. This book provides a gen- eral introduction to the main geophysical methods in widespread use.

Method Measured parameter Operative physical property Seismic

Gravity Magnetic Electrical Resistivity Induced polarization

Self-potential Electromagnetic Radar Travel times of re?ected/refracted seismic waves

Spatial variations in the strength of the gravitational ?eld of the Earth Spatial variations in the strength of the geomagnetic ?eld Earth resistance Polarization voltages or frequency- dependent ground resistance Electrical potentials Response to electromagnetic radiation Travel times of re?ected radar pulses Density and elastic moduli, which determine the propagation velocity of seismic waves Density

Magnetic susceptibility and remanence Electrical conductivity Electrical capacitance

Electrical conductivity Electrical conductivity and inductance Dielectric constant

being able to survey areas where ground access is dif?cult A wide range of geophysical surveying methods exists, for each of which there is an `operative’ physical property to which the method is sensitive.The methods The type of physical property to which a method Thus, for example, the magnetic method is very suitable for locating buried magnetite ore bodies because of their high magnetic susceptibility. Similarly, seismic or elec- trical methods are suitable for the location of a buried water table because saturated rock may be distinguished from dry rock by its higher seismic velocity and higher Other considerations also determine the type of methods employed in a geophysical exploration pro- gramme. For example, reconnaissance surveys are often carried out from the air because of the high speed of operation. In such cases the electrical or seismic methods are not applicable, since these require physical contact Thus, the initial search for metalliferous mineral deposits often utilizes airborne magnetic and electromagnetic surveying. Similarly, routine reconnaissance of conti- nental shelf areas often includes simultaneous gravity, magnetic and seismic surveying. At the interpretation stage, ambiguity arising from the results of one survey method may often be removed by consideration of results from a second survey method.

Principles of Exploration Methods 3 Application Appropriate survey methods* Exploration for fossil fuels (oil, gas, coal) S, G, M, (EM) Exploration for metalliferous mineral deposits M, EM, E, SP, IP, R Exploration for bulk mineral deposits (sand and gravel) S, (E), (G) Exploration for underground water supplies E, S, (G), (Rd) Engineering/construction site investigation E, S, Rd. (G), (M) Archaeological investigations Rd, E, EM, M, (S) * G, gravity; M, magnetic; S, seismic; E, electrical resistivity; SP, self-potential; IP, induced polarization; EM, electromagnetic; R, radiometric; Rd, ground-penetrating radar. Subsidiary methods in brackets.

+10 0 km 5 0 –40 –30 –20 –10

0 N 0 Fig. 1.1 The gravity anomaly over the Grand Saline Salt Dome,Texas, USA (contours in gravity units — see Chapter 6).The stippled area represents the subcrop of the dome. (Redrawn from Peters & Dugan 1945.)

20 80 60 80 N 100 0 km 5 40 40 Fig. 1.2 Magnetic anomalies over the Grand Saline Salt Dome,Texas, USA (contours in nT — see Chapter 7).The stippled area represents the subcrop of the dome. (Redrawn from Peters & Dugan 1945.)

Principles of Exploration Methods 5 Fig. 1.3 (a) Seismic re?ection section across a buried salt dome (courtesy Prakla-Seismos GmbH). (b) Simple structural interpretation of the seismic section, illustrating some possible ray paths for re?ected rays.

50 35 35 50 35 70 140 100 35 50 35 50 35 50 0 km 2 50 Fig. 1.4 Perturbation of telluric currents over the Haynesville The stippled area represents the subcrop of the dome. (Redrawn from Boissonas & Leonardon 1948.)

con?guration of closely-spaced shots and detectors is moved systematically along a pro?le line and the travel times of rays re?ected back from any subsurface geologi- cal interfaces are measured. If a salt dome is encountered, rays re?ected off its top surface will delineate the shape of 4. Earth materials with anomalous electrical resistivity may be located using either electrical or electromagnetic geophysical techniques. Shallow features are normally investigated using arti?cial ?eld methods in which an electrical current is introduced into the ground and potential differences between points on the surface are measured to reveal anomalous material in the subsurface (Chapter 8). However, this method is restricted in its depth of penetration by the limited power that can be introduced into the ground. Much greater penetration can be achieved by making use of the natural Earth cur- rents (telluric currents) generated by the motions of charged particles in the ionosphere. These currents ex- tend to great depths within the Earth and, in the absence of any electrically anomalous material, ?ow parallel to the surface. A salt dome, however, possesses an anom- alously high electrical resistivity and electric currents preferentially ?ow around and over the top of such a structure rather than through it. This pattern of ?ow causes distortion of the constant potential gradient at the surface that would be associated with a homogeneous subsurface and indicates the presence of the high- resistivity salt. Figure 1.4 presents the results of a telluric current survey of the Haynesville Salt Dome, Texas, USA.The contour values represent quantities describing the extent to which the telluric currents are distorted by subsurface phenomena and their con?guration re?ects the shape of the subsurface salt dome with some accuracy.

conversion of this travel time into a depth requires knowledge of the velocity with which the pulse travelled along the re?ection path and, unlike the velocity of If a velocity is assumed, a depth estimate can be derived And since rocks differ signi?cantly in the velocity with which they propagate seismic waves, it is by no means a straightforward matter to translate the travel time of a seismic pulse into an accurate depth to the geological in- The solution to this particular problem, as discussed in Chapter 4, is to measure the travel times of re?ected pulses at several offset distances from a seismic source because the variation of travel time as a function of range provides information on the velocity distribution with depth. However, although the degree of uncertainty in geophysical interpretation can often be reduced to an acceptable level by the general expedient of taking additional (and in some cases different kinds of ) ?eld measurements, the problem of inherent ambiguity The general problem is that signi?cant differences from an actual subsurface geological situation may give rise to insigni?cant, or immeasurably small, differences in the quantities actually measured during a geophysical survey. Thus, ambiguity arises because many different geological con?gurations could reproduce the observed measurements. This basic limitation results from the unavoidable fact that geophysical surveying attempts to solve a dif?cult inverse problem. It should also be noted that experimentally-derived quantities are never exactly determined and experimental error adds a Principles of Exploration Methods 7

further degree of indeterminacy to that caused by the incompleteness of the ?eld data and the ambiguity associated with the inverse problem. Since a unique solution cannot, in general, be recovered from a set of ?eld measurements, geophysical interpretation is concerned either to determine properties of the subsurface that all possible solutions share, or to introduce assumptions to restrict the number of admissible solutions (Parker 1977). In spite of these inherent problems, however, geophysical surveying is an invaluable tool for the investigation of subsurface geology and occupies a key role in exploration programmes for geological resources.

1.4 The structure of the book The above introductory sections illustrate in a simple way the very wide range of approaches to the geophysical investigation of the subsurface and warn Chapter 2 provides a short account of the more important data processing techniques of general applicability to geophysics. In Chapters 3 to 10 the individual survey methods are treated systematically in terms of their basic principles, survey procedures, Chapter 11 describes the application of these methods to specialized surveys undertaken in boreholes. All these chapters contain suggestions for further reading which provide a more extensive treatment of the material covered in this book. A set of problems is given for all the major geophysical methods.

2.1 Introduction Geophysical surveys measure the variation of some physical quantity, with respect either to position or to time. The quantity may, for example, be the strength of the Earth’s magnetic ?eld along a pro?le across an igneous intrusion. It may be the motion of the ground surface as a function of time associated with the passage of seismic waves. In either case, the simplest way to pre- sent the data is to plot a graph (Fig. 2.1) showing the vari- ation of the measured quantity with respect to distance or time as appropriate. The graph will show some more or less complex waveform shape, which will re?ect physical variations in the underlying geology, superim- posed on unwanted variations from non-geological fea- tures (such as the effect of electrical power cables in the magnetic example, or vibration from passing traf?c for the seismic case), instrumental inaccuracy and data col- lection errors. The detailed shape of the waveform may be uncertain due to the dif?culty in interpolating the curve between widely spaced stations.The geophysicist’s task is to separate the`signal’from the`noise’and interpret Analysis of waveforms such as these represents an es- sential aspect of geophysical data processing and inter- pretation. The fundamental physics and mathematics of such analysis is not novel, most having been discovered in the 19th or early 20th centuries.The use of these ideas is also widespread in other technological areas such as radio, television, sound and video recording, radio- astronomy, meteorology and medical imaging, as well as military applications such as radar, sonar and satellite imaging. Before the general availability of digital com- puting, the quantity of data and the complexity of the processing severely restricted the use of the known tech- niques. This no longer applies and nearly all the tech- niques described in this chapter may be implemented in The fundamental principles on which the various methods of data analysis are based are brought together in this chapter.These are accompanied by a discussion of the techniques of digital data processing by computer that are routinely used by geophysicists.Throughout this chapter, waveforms are referred to as functions of time, but all the principles discussed are equally applicable to functions of distance. In the latter case, frequency (num- ber of waveform cycles per unit time) is replaced by spatial frequency or wavenumber (number of waveform cycles per unit distance).

Geophysical Data Processing 9 600 500 400 300 (a) 200 0 10 Magnetic field (nT) 20 30 40 50 Distance (m) Fig. 2.1 (a) A graph showing a typical magnetic ?eld strength variation which may be measured along a pro?le. (b) A graph of a typical seismogram, showing variation of particle velocities in the ground as a function of time during the passage of a seismic wave.

(a) f (t ) Ground velocity (10-6 m/s) 15 10 5 0 –5 (b) –10

1.0 0 10 20 30 40 50 60 70 80 Time (milliseconds)

?ne electrical power ratios: the ratio of two power values P and P is given by 10 log (P /P ) dB. Since power is 1 2 10 1 2 proportional to the square of signal amplitude A

2 10 log (P P ) = 10 log (A A ) 10 1 2 10 1 2 = 20 log (A A ) (2.1) 10 1 2 t

Thus, if a digital sampling scheme measures ampli- tudes over the range from 1 to 1024 units of amplitude, the dynamic range is given by –1.0 20 log (A A ) = 20 log 1024 ª 60 dB 10 max min 10

0.9 0.9 0.5 0.5 ? 2? 3? 0.0 t (b) g(t ) 1.0

signi?cant loss of information content as long as the fre- quency of sampling is much higher than the highest frequency component in the sampled function. Mathe- matically, it can be proved that, if the waveform is a sine curve, this can always be reconstructed provided that there are a minimum of two samples per period of the Thus, if a waveform is sampled every two milliseconds (sampling interval), the sampling frequency is 500 sam- ples per second (or 500 Hz). Sampling at this rate will preserve all frequencies up to 250 Hz in the sampled function. This frequency of half the sampling frequency is known as the Nyquist frequency ( f ) and the Nyquist N interval is the frequency range from zero up to f N

f = 1 (2Dt) N (2.2)

If frequencies above the Nyquist frequency are pre- sent in the sampled function, a serious form of distortion results known as aliasing, in which the higher frequency Consider the example illustrated in Fig. 2.3 in which sine waves at different frequencies are sampled. The lower frequency wave (Fig. 2.3(a)) is accurately repro- duced, but the higher frequency wave (Fig. 2.3(b), solid line) is rendered as a ?ctitious frequency, shown by the dashed line, within the Nyquist interval. The relation- ship between input and output frequencies in the case of a sampling frequency of 500 Hz is shown in Fig. 2.3(c). It is apparent that an input frequency of 125 Hz, for exam- ple, is retained in the output but that an input frequency To overcome the problem of aliasing, the sampling frequency must be at least twice as high as the highest fre- quency component present in the sampled function. If the function does contain frequencies above the Nyquist frequency determined by the sampling, it must be passed through an antialias ?lter prior to digitization. The antialias ?lter is a low-pass frequency ?lter with a sharp cut-off that removes frequency components above the Nyquist frequency, or attenuates them to an insigni?cant amplitude level.

2.3 Spectral analysis An important mathematical distinction exists between periodic waveforms (Fig. 2.4(a)), that repeat themselves at a ?xed time period T, and transient waveforms (Fig. 2.4(b)), (a)

(b) (c) f 2f 3f 4f NNNN 250 Output frequency (Hz) 125

0 125 250 500 625 750 1000 Input frequency (Hz) (b) Sine wave frequency greater than Nyquist frequency (solid line) showing the ?ctitious frequency that is generated by aliasing (c) Relationship between input and output frequencies for a N

Geophysical Data Processing 11 T (a) 8 8 (b) (a) 2 Amplitude 1

0 (a) (b) Fig. 2.5 Complex waveforms resulting from the summation of two sine wave components of frequency f and 2f. (a)The two sine wave components are of equal amplitude and in phase. (b)The higher frequency component has twice the amplitude of the lower frequency component and is p/2 out of phase. (After Anstey 1965.)

(b) 2 Amplitude 1 0 f 2f Frequency ?/2?/2

0 0 2f Frequency Phase Fig. 2.6 Representation in the frequency domain of the waveforms illustrated in Fig. 2.5, showing their – ?/2 amplitude and phase spectra.

Phase f produce the quite different waveform illustrated in From the above it follows that a periodic waveform can be expressed in two different ways: in the familiar time domain, expressing wave amplitude as a function of time, or in the frequency domain, expressing the amplitude and phase of its constituent sine waves as a function of frequency. The waveforms shown in Fig. 2.5(a) and (b) are represented in Fig. 2.6(a) and (b) in terms of their am- plitude and phase spectra. These spectra, known as line spectra, are composed of a series of discrete values of the amplitude and phase components of the waveform at set frequency values distributed between 0 Hz and the Nyquist frequency.

– ?/2 f 2f Frequency f 2f Frequency

Amplitude density Frequency Phase Frequency Fig. 2.7 Digital representation of the continuous amplitude and phase spectra associated with a transient waveform.

thin frequency slices, with each slice having a frequency equal to the mean frequency of the slice and an ampli- tude and phase proportional to the area of the slice of the appropriate spectrum (Fig. 2.7). This digital expression of a continuous spectrum in terms of a ?nite number of discrete frequency components provides an approximate representation in the frequency domain of a transient waveform in the time domain. Increasing the sampling frequency in the time domain not only improves the time-domain representation of the waveform, but also increases the number of frequency slices in the frequen- cy domain and improves the accuracy of the approxima- Fourier transformation may be used to convert a time function g(t) into its equivalent amplitude and phase spectra A( f ) and f( f ), or into a complex function of frequency G( f ) known as the frequency spectrum, where

G ( f ) = A( f )e if ( f ) (2.3) The time- and frequency-domain representations of a waveform, g(t) and G( f ), are known as a Fourier pair, represented by the notation g (t) ´ G( f ) (2.4) Components of a Fourier pair are interchangeable, such that, if G( f ) is the Fourier transform of g(t ), then g(t) is the Fourier transform of G( f ). Figure 2.8 illus- trates Fourier pairs for various waveforms of geophysical signi?cance. All the examples illustrated have zero phase spectra; that is, the individual sine wave components of the waveforms are in phase at zero time. In this case f( f ) = 0 for all values of f. Figure 2.8(a) shows a spike func- tion (also known as a Dirac function), which is the shortest possible transient waveform. Fourier transformation shows that the spike function has a continuous frequency thus, a spike function contains all frequencies from zero to in?nity at equal amplitude.The`DC bias’waveform of Fig. 2.8(b) has, as would be expected, a line spectrum comprising a single component at zero frequency. Note that Fig. 2.8(a) and (b) demonstrate the principle of interchangeability of Fourier pairs stated above (equa- tion (2.4)). Figures 2.8(c) and (d) illustrate transient waveforms approximating the shape of seismic pulses, together with their amplitude spectra. Both have a band- limited amplitude spectrum, the spectrum of narrower bandwidth being associated with the longer transient waveform. In general, the shorter a time pulse the wider is its frequency bandwidth and in the limiting case a spike Waveforms with zero phase spectra such as those illus- trated in Fig. 2.8 are symmetrical about the time axis and, for any given amplitude spectrum, produce the maximum peak amplitude in the resultant waveform. If phase varies linearly with frequency, the waveform re- mains unchanged in shape but is displaced in time; if the phase variation with frequency is non-linear the shape of the waveform is altered. A particularly important case in seismic data processing is the phase spectrum associated with minimum delay in which there is a maximum con- Analysis of seismic pulses sometimes assumes that they Fourier transformation of digitized waveforms is readily programmed for computers, using a `fast Fourier transform’ (FFT) algorithm as in the Cooley–Tukey method (Brigham 1974). FFT subroutines can thus be routinely built into data processing programs in order to Fourier transformation is supplied as a function to standard spreadsheets such as Microsoft Excel. Fourier transformation can be extended into two dimensions (Rayner 1971), and can thus be applied to areal distribu- tions of data such as gravity and magnetic contour maps.

Geophysical Data Processing 13 Time domain Frequency domain

(a) (b) (c) (d) Fig. 2.8 Fourier transform pairs for (b)A`DCbias’.(c)and(d)Transient Time Frequency waveforms approximating seismic pulses. t = 0

In this case, the time variable is replaced by horizontal distance and the frequency variable by wavenumber (number of waveform cycles per unit distance). The application of two-dimensional Fourier techniques to the interpretation of potential ?eld data is discussed in Chapters 6 and 7.

2.4 Waveform processing The principles of convolution, deconvolution and cor- relation form the common basis for many methods of geophysical data processing, especially in the ?eld of seis- mic re?ection surveying. They are introduced here in general terms and are referred to extensively in later chapters.Their importance is that they quantitatively de- scribe how a waveform is affected by a ?lter. Filtering modi?es a waveform by discriminating between its con- stituent sine wave components to alter their relative am- plitudes or phase relations, or both. Most audio systems are provided with simple ?lters to cut down on high- frequency `hiss’, or to emphasize the low-frequency `bass’. Filtering is an inherent characteristic of any system through which a signal is transmitted.

Amplitude Input displacement W Output displacement Input

Output Time Fig. 2.9 The principle of ?ltering illustrated by the perturbation of a suspended weight system.

response which is de?ned as the output of the ?lter when the input is a spike function (Fig. 2.10). The impulse re- sponse is a waveform in the time domain, but may be transformed into the frequency domain as for any other waveform. The Fourier transform of the impulse re- sponse is known as the transfer function and this speci?es the amplitude and phase response of the ?lter, thus de?ning its operation completely.The effect of a ?lter is described mathematically by a convolution operation such that, if the input signal g(t) to the ?lter is convolved with the impulse response f(t) of the ?lter, known as the con- volution operator, the ?ltered output y(t) is obtained:

y(t) = g (t) * f (t) (2.5)

Figure 2.11(a) shows a spike function input to a ?lter whose impulse response is given in Fig. 2.11(b). Clearly the latter is also the ?ltered output since, by de?nition, the impulse response represents the output for a spike input. Figure 2.11(c) shows an input comprising two 2.11(d)) is now the superposition of the two impulse re- sponse functions offset in time by the separation of the

Spike input input spikes and scaled according to the individual spike amplitudes. Since any transient wave can be represented as a series of spike functions (Fig. 2.11(e)), the general form of a ?ltered output (Fig. 2.11(f )) can be regarded as the summation of a set of impulse responses related to a succession of spikes simulating the overall shape of the The mathematical implementation of convolution involves time inversion (or folding) of one of the func- tions and its progressive sliding past the other function, the individual terms in the convolved output being de- rived by summation of the cross-multiplication products over the overlapping parts of the two functions. In gen- eral, if g (i = 1, 2, . . . , m) is an input function and f ( j = ij 1, 2, . . . , n) is a convolution operator, then the convolu- tion output function y is given by k

m y = Â g f (k = 1, 2, . . . , m + n – 1) (2.6) k i k-i i =1

In Fig. 2.12 the individual steps in the convolution process are shown for two digital functions, a double spike function given by g = g , g , g = 2, 0, 1 and an im- i123 pulse response function given by f = f , f , f , f = 4, 3, 2, i1234 1, where the numbers refer to discrete amplitude values 2.11 it can be seen that the convolved output y = y , y , i12 y , y , y , y = 8, 6, 8, 5, 2, 1. Note that the convolved 3456 output is longer than the input waveforms; if the func- tions to be convolved have lengths of m and n, the con- The convolution of two functions in the time domain becomes increasingly laborious as the functions become longer. Typical geophysical applications may have func- tions which are each from 250 to a few thousand samples long. The same mathematical result may be obtained by transforming the functions to the frequency domain, then multiplying together equivalent frequency terms of their amplitude spectra and adding terms of their phase spectra. The resulting output amplitude and phase spec- Thus, digital ?ltering can be enacted in either the time

Geophysical Data Processing 15 (a) (c) (e) (b) (d) (f) Fig. 2.11 Examples of ?ltering. (a) A spike input. (b) Filtered output equivalent to impulse response of ?lter. (c) An input comprising two spikes. (d) Filtered output given by summation of two impulse response functions offset in time. (e) A complex input represented by a series of contiguous spike functions. (f) Filtered output given by the summation of a set of impulse responses.

4 3 2 1 Cross-products Sum 1 0 2 4×2 = 1 0 2 1 0 2 4×0+3×2 = 4×1+3×0+2×2 = 1 0 2 3×1+2×0+2×1 = 1 0 2 2×1+1×0 = Fig. 2.12 A method of calculating the convolution of two digital functions.

domain or the frequency domain. With large data sets, ?ltering by computer is more ef?ciently carried out in the frequency domain since fewer mathematical opera- Convolution, or its equivalent in the frequency 1 0 2 1×1 = 8 6 8 5 2 1

domain, ?nds very wide application in geophysical data processing, notably in the digital ?ltering of seismic and potential ?eld data and the construction of synthetic seismograms for comparison with ?eld seismograms (see Chapters 4 and 6).

2.4.2 Deconvolution Deconvolution or inverse ?ltering (Kanasewich 1981) is a process that counteracts a previous convolution (or ?ltering) action. Consider the convolution operation given in equation (2.5) y(t) = g (t) * f (t)

y(t) is the ?ltered output derived by passing the input Knowing y(t) and f(t), the recovery of g(t) represents a de- convolution operation. Suppose that f ¢(t) is the function that must be convolved with y(t) to recover g(t) g(t ) = y(t ) * f ¢(t ) (2.7)

Substituting for y(t) as given by equation (2.5) g(t ) = g(t ) * f (t ) * f ¢(t ) (2.8) Now recall also that

g(t ) = g(t ) * d (t ) (2.9) where d(t) is a spike function (a unit amplitude spike at zero time); that is, a time function g(t) convolved with a spike function produces an unchanged convolution output function g(t). From equations (2.8) and (2.9) it follows that

f (t ) * f ¢(t ) = d (t ) (2.10)

Thus, provided the impulse response f(t) is known, f¢(t) can be derived for application in equation (2.7) to re- cover the input signal g(t). The function f ¢(t) represents Deconvolution is an essential aspect of seismic data processing, being used to improve seismic records by re- moving the adverse ?ltering effects encountered by seis- mic waves during their passage through the ground. In the seismic case, referring to equation (2.5), y(t) is the seismic record resulting from the passage of a seismic wave g(t) through a portion of the Earth, which acts as a ?lter with an impulse response f(t). The particular prob- lem with deconvolving a seismic record is that the input waveform g(t) and the impulse response f(t) of the Earth ?lter are in general unknown. Thus the `deterministic’ approach to deconvolution outlined above cannot be employed and the deconvolution operator has to be designed using statistical methods.This special approach to the deconvolution of seismic records, known as pre- dictive deconvolution, is discussed further in Chapter 4.

2.4.3 Correlation Cross-correlation of two digital waveforms involves cross- multiplication of the individual waveform elements and summation of the cross-multiplication products over the common time interval of the waveforms.The cross- correlation function involves progressively sliding one waveform past the other and, for each time shift, or lag, summing the cross-multiplication products to derive the cross-correlation as a function of lag value. The cross- correlation operation is similar to convolution but does not involve folding of one of the waveforms. Given two digital waveforms of ?nite length, x and y (i = 1, 2, . . . , ii n), the cross-correlation function is given by n -t

f (t ) = Â x y (-m < t < +m) (2.11) xy i+t i i =1

Geophysical Data Processing 17 Waveform 1 Waveform 2 lag Cross-correlation function

– ve lag + ve lag Fig. 2.13 Cross-correlation of two Zero identical waveforms. lag

centred on the time value at which the signal function and its concealed equivalent in the waveform are in A special case of correlation is that in which a wave- form is cross-correlated with itself, to give the autocorre- lation function f (t).This function is symmetrical about xx a zero lag position, so that

fxx (t ) = fxx (-t ) (2.12)

The autocorrelation function of a periodic waveform is also periodic, with a frequency equal to the repetition frequency of the waveform.Thus, for example, the auto- correlation function of a cosine wave is also a cosine wave. For a transient waveform, the autocorrelation These differing properties of the autocorrelation func- tion of periodic and transient waveforms determine one of its main uses in geophysical data processing, namely, the detection of hidden periodicities in any given wave- 2.15) are an indication of the existence of periodicities in the original waveform, and the spacing of the side lobes de?nes the repetition period.This property is particular- ly useful in the detection and suppression of multiple re?ections in seismic records (see Chapter 4).

The autocorrelation function contains all the am- plitude information of the original waveform but none of the phase information, the original phase relation- ships being replaced by a zero phase spectrum. In fact, the autocorrelation function and the square of the am- plitude spectrum A( f ) can be shown to form a Fourier pair

2 f (t ) ´ A( f ) (2.13) xx

Since the square of the amplitude represents the power term (energy contained in the frequency component) the autocorrelation function can be used to compute the power spectrum of a waveform.

Waveform Signal function Cross-correlation function (a) ? (b) S SS 1 23 Fig. 2.14 Cross-correlation to detect Signal positions occurrences of a known signal concealed in waveform in noise. (After Sheriff 1973.)

? (?) xx Fig. 2.15 Autocorrelation of the waveform exhibiting periodicity shown in (a) produces the autocorrelation function with side lobes shown in (b).The spacing of the side lobes de?nes the repetition period of the original ? waveform.

random, and usually due to effects unconnected with the geophysical survey. Coherent noise is, on the other hand, components of the waveform which are generated by the geophysical experiment, but are of no direct interest for the geological interpretation. For example, in a seismic survey the signal might be the seismic pulse arriving at a detector after being re?ected by a geological boundary at depth. Random noise would be back- Coherent noise would be the surface waves generated by the seismic source, which also travel to the detector In favourable circumstances the signal-to-noise ratio (SNR) is high, so that the signal is readily identi?ed and extracted for subsequent analysis. Often the SNR is low and special processing is necessary to enhance the infor- mation content of the waveforms. Different approaches are needed to remove the effect of different types of noise. Random noise can often be suppressed by re- peated measurement and averaging. Coherent noise may be ?ltered out by identifying the particular charac- teristics of that noise and designing a special ?lter to re- move it.The remaining signal itself may be distorted due to the effects of the recording system, and again, if the nature of the recording system is accurately known, suit- able ?ltering can be designed. Digital ?ltering is widely employed in geophysical data processing to improve SNR or otherwise improve the signal characteristics. A very wide range of digital ?lters is in routine use in geo- physical, and especially seismic, data processing (Robin- son & Treitel 2000). The two main types of digital ?lter are frequency ?lters and inverse (deconvolution) ?lters.

Geophysical Data Processing 19 Frequency domain Time domain (a) (b)

Sinc function ff ct (c) (d)

Filter operator Fig. 2.16 Design of a digital low-pass ff ?lter. ct

(BR) in terms of their frequency response. Frequency ?lters are employed when the signal and noise compo- nents of a waveform have different frequency character- Analogue frequency ?ltering is still in widespread use and analogue antialias (LP) ?lters are an essential compo- nent of analogue-to-digital conversion systems (see Sec- tion 2.2). Nevertheless, digital frequency ?ltering by computer offers much greater ?exibility of ?lter design and facilitates ?ltering of much higher performance than can be obtained with analogue ?lters. To illustrate the design of a digital frequency ?lter, consider the case of a LP ?lter whose cut-off frequency is f . The desired out- c put characteristics of the ideal LP ?lter are represented by the amplitude spectrum shown in Fig. 2.16(a).The spec- trum has a constant unit amplitude between 0 and f and c zero amplitude outside this range: the ?lter would there- fore pass all frequencies between 0 and f without atten- c c This amplitude spectrum represents the transfer func- Inverse Fourier transformation of the transfer func- tion into the time domain yields the impulse response of the ideal LP ?lter (see Fig. 2.16(b)). However, this im- pulse response (a sinc function) is in?nitely long and must therefore be truncated for practical use as a convo- lution operator in a digital ?lter. Figure 2.16(c) repre- sents the frequency response of a practically realizable LP ?lter operator of ?nite length (Fig. 2.16(d)). Convolu- tion of the input waveform with the latter will result in LP ?ltering with a ramped cut-off (Fig. 2.16(c)) rather HP, BP and BR time-domain ?lters can be designed in a similar way by specifying a particular transfer func- tion in the frequency domain and using this to design a ?nite-length impulse response function in the time do- main.As with analogue ?ltering, digital frequency ?lter- ing generally alters the phase spectrum of the waveform and this effect may be undesirable. However, zero phase ?lters can be designed that facilitate digital ?ltering with- out altering the phase spectrum of the ?ltered signal.

2.5.2 Inverse (deconvolution) ?lters The main applications of inverse ?ltering to remove the adverse effects of a previous ?ltering operation lie in the ?eld of seismic data processing. A discussion of inverse ?ltering in the context of deconvolving seismic records is given in Chapter 4.

forms themselves are presented in a form in which they simulate an image of the subsurface structure. The most obvious examples of this are in seismic re?ection (Chapter 4) and ground-penetrating radar (Chapter 9) sections, where the waveform of the variation of re?ect- ed energy with time is used to derive an image related to the occurrence of geological boundaries at depth. Often magnetic surveys for shallow engineering or archaeo- logical investigations are processed to produce shaded, coloured, or contoured maps where the shading or colour correlates with variations of magnetic ?eld which are expected to correlate with the structures being sought. Imaging is a very powerful tool, as it provides a way of summarizing huge volumes of data in a format which can be readily comprehended, that is, the visual image. A disadvantage of imaging is that often it can be dif?cult or impossible to extract quantitative informa- In modelling, the geophysicist chooses a particular type of structural model of the subsurface, and uses this The model is then adjusted to give the closest match be- tween the predicted (modelled) and observed wave- forms. The goodness of the match obtained depends on both the signal-to-noise ratio of the waveforms and the initial choice of the model used.The results of modelling are usually displayed as cross-sections through the struc- ture under investigation. Modelling is an essential part of most geophysical methods and is well exempli?ed in gravity and magnetic interpretation (see Chapters 6 and 7).

Problems 1. Over the distance between two seismic recording sites at different ranges from a seismic source, seismic waves have been attenuated by 5 dB. What is the ratio of the wave amplitudes ob- 2. In a geophysical survey, time-series data are (a) What is the Nyquist frequency? (b) In the absence of antialias ?ltering, at what frequency would noise at 200 Hz be aliased back into the 3. If a digital recording of a geophysical time series is required to have a dynamic range of 120 dB, what number of bits is required in each binary word?

4. If the digital signal (-1, 3, -2, -1) is convolved with the ?lter operator (2, 3, 1), what is the con- 5. Cross-correlate the signal function (-1, 3, -1) with the waveform (-2, -4, -4, -3, 3, 1, 2, 2) con- taining signal and noise, and indicate the likely position of the signal in the waveform on the 6. A waveform is composed of two in-phase components of equal amplitude at frequencies f and 3f. Draw graphs to represent the waveform in the time domain and the frequency domain.

Further reading Brigham, E.O. (1974) The Fast Fourier Transform. Prentice-Hall, Camina, A.R. & Janacek, G.J. (1984) Mathematics for Seismic Data Dobrin, M.B. & Savit, C.H. (1988) Introduction to Geophysical Kanasewich, E.R. (1981) Time Sequence Analysis in Geophysics (3rd edn). University of Alberta Press.

Rayner, J.N. (1971) An Introduction to Spectral Analysis. Pion, Sheriff, R.E. & Geldart, L.P. (1983) Exploration Seismology Vol 2: Data-Processing and Interpretation. Cambridge University Press, Cambridge.

3.1 Introduction In seismic surveying, seismic waves are created by a con- Some waves will return to the surface after refraction or re?ection at geological boundaries within the subsur- face. Instruments distributed along the surface detect the ground motion caused by these returning waves and hence measure the arrival times of the waves at different ranges from the source. These travel times may be con- verted into depth values and, hence, the distribution of subsurface geological interfaces may be systematically Seismic surveying was ?rst carried out in the early 1920s. It represented a natural development of the already long-established methods of earthquake seis- mology in which the travel times of earthquake waves recorded at seismological observatories are used to de- Earthquake seismology provides information on the gross internal layering of the Earth, and measurement of the velocity of earthquake waves through the various Earth layers provides information about their physical properties and composition. In the same way, but on a smaller scale, seismic surveying can provide a clear and detailed picture of subsurface geology. It undoubtedly represents the single most important geophysical survey- ing method in terms of the amount of survey activity and the very wide range of its applications. Many of the principles of earthquake seismology are applicable to seismic surveying. However, the latter is concerned solely with the structure of the Earth down to tens of kilometres at most and uses arti?cial seismic sources, such as explosions, whose location, timing and source characteristics are, unlike earthquakes, under the direct control of the geophysicist. Seismic surveying also uses specialized recording systems and associated data pro- Seismic methods are widely applied to exploration problems involving the detection and mapping of sub- surface boundaries of, normally, simple geometry.They also identify signi?cant physical properties of each sub- surface unit. The methods are particularly well suited to the mapping of layered sedimentary sequences and are therefore widely used in the search for oil and gas. The methods are also used, on a smaller scale, for the mapping of near-surface sediment layers, the location of the water table and, in an engineering context, site investigation of foundation conditions including the determination of depth to bedrock. Seismic surveying can be carried out on land or at sea and is used extensively in offshore geological surveys and the exploration for offshore In this chapter the fundamental physical principles on which seismic methods are based are reviewed, starting with a discussion of the nature of seismic waves and going on to consider their mode of propagation through the ground, with particular reference to re?ection and refraction at interfaces between different rock types. To understand the different types of seismic wave that propagate through the ground away from a seismic source, some elementary concepts of stress and strain need to be considered.

Stress Chapter 3 Elastic field Ductile field Yield point Fracture point

Strain known as the principal axes of stress, and the normal stresses acting in these directions are known as the princi- pal stresses. Each principal stress represents a balance of equal-magnitude but oppositely-directed force compo- nents.The stress is said to be compressive if the forces are directed towards each other and tensile if they are If the principal stresses are all of equal magnitude within a body the condition of stress is said to be hydro- static, since this is the state of stress throughout a ?uid body at rest. A ?uid body cannot sustain shearing stresses (since a ?uid has no shear strength), hence there cannot be shear stresses in a body under hydrostatic stress. If the principal stresses are unequal, shearing stresses exist along all surfaces within the stressed body, except for the three orthogonal planes intersecting in the principal A body subjected to stress undergoes a change of shape and/or size known as strain. Up to a certain limit- ing value of stress, known as the yield strength of a ma- terial, the strain is directly proportional to the applied stress (Hooke's Law). This elastic strain is reversible so that removal of stress leads to a removal of strain. If the yield strength is exceeded the strain becomes non-linear and partly irreversible (i.e. permanent strain results), and this is known as plastic or ductile strain. If the stress is in- creased still further the body fails by fracture. A typical The linear relationship between stress and strain in the elastic ?eld is speci?ed for any material by its various elas- tic moduli, each of which expresses the ratio of a particu- lar type of stress to the resultant strain. Consider a rod of original length l and cross-sectional area A which is ex- tended by an increment Dl through the application of a stretching force F to its end faces (Fig. 3.2(a)). The rele- vant elastic modulus isYoung' s modulus E, de?ned by

longitudinal stress F A E= longitudinal strain Dl l Note that extension of such a rod will be accompanied by a reduction in its diameter; that is, the rod will suffer lateral as well as longitudinal strain.The ratio of the lateral The bulk modulus K expresses the stress–strain ratio in the case of a simple hydrostatic pressure P applied to a cubic element (Fig. 3.2(b)), the resultant volume strain being the change of volume Dv divided by the original volume v

volume stress P K= volume strain Dv v In a similar manner the shear modulus (m) is de?ned as the ratio of shearing stress (t) to the resultant shear strain tan q (Fig. 3.2(c))

shear stress t m= shear strain tan q Finally, the axial modulus y de?nes the ratio of longi- tudinal stress to longitudinal strain in the case when there is no lateral strain; that is, when the material is con- strained to deform uniaxially (Fig. 3.2(d))

longitudinal stress F A y= longitudinal strain (uniaxial) Dl l

Elements of Seismic Surveying 23 P (a) l (b)

F F P P l + ?l longitudinal stress F/A volume stress P E= K= longitudinal strain ?l/l volume strain ?v/v

? (c) (d) l ? F F

l + ?l Fig. 3.2 The elastic moduli. (a)Young' s shear stress ? longitudinal stress F/A µ= ?= modulus E. (b) Bulk modulus K. (c) Shear shear strain tan ? longitudinal strain ?l/l modulus m. (d) Axial modulus y. (no lateral strain)

which they pass. There are two groups of seismic waves, 3.3.1 Body waves Body waves can propagate through the internal volume of an elastic solid and may be of two types. Compressional waves (the longitudinal, primary or P-waves of earth- quake seismology) propagate by compressional and dila- Particle motion associated with the passage of a com- pressional wave involves oscillation, about a ?xed point, in the direction of wave propagation (Fig. 3.3(a)). Shear waves (the transverse, secondary or S-waves of earth- quake seismology) propagate by a pure shear strain in a Individual particle motions involve oscillation, about a ?xed point, in a plane at right angles to the direction of wave propagation (Fig. 3.3(b)). If all the particle oscilla- tions are con?ned to a plane, the shear wave is said to be The velocity of propagation of any body wave in any homogeneous, isotropic material is given by:

12 È appropriate elastic modulus of material ? v= Hence the velocity v of a compressional body wave, p which involves a uniaxial compressional strain, is given by 12 v = Èy ? r

or, since y = K + 43 m, by 12 È K + 43 m ? v= r and the velocity v of a shear body wave, which involves a s pure shear strain, is given by

12 v =Èm? r It will be seen from these equations that compres- sional waves always travel faster than shear waves in the same medium. The ratio v /v in any material is deter- ps mined solely by the value of Poisson's ratio (s ) for that material

(a) P-wave Compressions Undisturbed medium Dilatations (b) S-wave

Fig. 3.3 Elastic deformations and ground particle motions associated with the passage (From Bolt 1982.)

independent of density and can be used to derive Pois- son's ratio, which is a much more diagnostic lithological indicator. If this information is required, then both v p s These fundamental relationships between the veloc- ity of the wave propagation and the physical properties of the materials through which the waves pass are inde- pendent of the frequency of the waves. Body waves are non-dispersive; that is, all frequency components in a wave train or pulse travel through any material at the same velocity, determined only by the elastic moduli and Historically, most seismic surveying has used only compressional waves, since this simpli?es the survey technique in two ways. Firstly, seismic detectors which record only the vertical ground motion can be used, and these are insensitive to the horizontal motion of S- waves. Secondly, the higher velocity of P-waves ensures that they always reach a detector before any related S-waves, and hence are easier to recognize. Recording S-waves, and to a lesser extent surface waves, gives greater information about the subsurface, but at a cost of greater data acquisition (three-component recording) and consequent processing effort. As technology ad- vances multicomponent surveys are becoming more One application of shear wave seismology is in engi- neering site investigation where the separate measure- ment of v and v for near-surface layers allows direct ps calculation of Poisson's ratio and estimation of the elas- tic moduli, which provide valuable information on the in situ geotechnical properties of the ground.These may be of great practical importance, such as the value of rip- pability (see Section 5.11.1).

Elements of Seismic Surveying (a) (b) 25 Fig. 3.4 Elastic deformations and ground particle motions associated with the passage of surface waves. (a) Rayleigh wave. (b) Love wave. (From Bolt 1982.)

interior.Analysisoftheobservedpatternofdispersionof earthquake waves is a powerful method of studying the velocity structure of the lithosphere and asthenosphere (Knopoff 1983). The same methodology, applied to the surface waves generated by a sledgehammer, can be used to examine the strength of near-surface materials for If the surface is layered and the surface layer shear wave velocity is lower than that of the underlying layer, a second set of surface waves is generated. Love waves are polarized shear waves with a particle motion parallel to the free surface and perpendicular to the direction of wave propagation (Fig. 3.4(b)). The velocity of Love waves is intermediate between the shear wave velocity of the surface layer and that of deeper layers, and Love waves are inherently dispersive. The observed pattern of Love wave dispersion can be used in a similar way to Rayleigh wave dispersion to study the subsurface structure.

3.3.3 Waves and rays A seismic pulse propagates outwards from a seismic source at a velocity determined by the physical proper- ties of the surrounding rocks. If the pulse travels through a homogeneous rock it will travel at the same velocity in all directions away from the source so that at any subse- quent time the wavefront, de?ned as the locus of all points which the pulse has reached at a particular time, will be a sphere. Seismic rays are de?ned as thin pencils of seismic energy travelling along ray paths that, in isotropic media, are everywhere perpendicular to wavefronts (Fig. 3.5).

Wavefront Source Ray path Fig. 3.5 The relationship of a ray path to the associated wavefront.

of only about 10-10 m.The detection of seismic waves in- 3.4 Seismic wave velocities of rocks grain shape and degree of sorting), porosities and con- tained pore ?uids, rocks differ in their elastic moduli and densities and, hence, in their seismic velocities. Informa- tion on the compressional and shear wave velocities, v p and v , of rock layers encountered by seismic surveys is s important for two main reasons: ?rstly, it is necessary for secondly, it provides an indication of the lithology of a rock or, in some cases, the nature of the pore ?uids To relate rock velocities to lithology, the assumption that rocks are uniform and isotropic in structure must be reviewed. A typical rock texture can be regarded as having mineral grains making up most of the rock (the matrix), with the remaining volume being occupied by void space (the pores). The fractional volume of pore space is the porosity (f) of the rock. For simplicity it may be assumed that all the matrix grains have the same physical properties. This is a surprisingly good approxi- mation since the major rock-forming minerals, quartz, feldspar and calcite, have quite similar physical proper- ties. In this case, the properties of the bulk rock will be an average of the properties of the matrix minerals and the pore ?uid, weighted according to the porosity.The sim- plest case is for the density of a rock, where the bulk density r can be related to the matrix and pore ?uid b densities (r , r ): mf

r = r f + (1 - f ) r bf m For P-wave velocity a similar relationship exists, but the velocity weighting is proportional to the percentage of travel-time spent in each component of the system, which is inversely proportional to velocity, giving the relationship: 1 f (1 - f ) =+ vvv bfm

From the above equations it is possible to produce cross-plot graphs (Fig. 3.6) which allow the estimation of the matrix grain type and the porosity of a rock, purely from the seismic P-wave velocity and density.

6500 6000 5500 5000 -1 4500 4000 3500 Seismic P-wave velocity in m s 3000

2500 2000 1500 1000 Density–Velocity Cross-plot Sandstone Limestone

0% 50% 100% 1000 1500 2000 2500 3000 Density in kg m-3

Fig. 3.6 The relationship of seismic velocity and density to porosity, calculated for mono-mineralic granular solids: open circles – sandstone, calculated for a quartz matrix; solid circles – limestone, calculated for a calcite matrix. Points annotated with the corresponding porosity value 0–100%. Such relationships are useful in borehole log interpretation (see Chapter 11).

and 5. If boreholes exist in the vicinity of a seismic survey, it may be possible to correlate velocity values so derived with individual rock units encountered within borehole sequences. As discussed in Chapter 11, veloc- ity may also be measured directly in boreholes using a sonic probe, which emits high-frequency pulses and measures the travel time of the pulses through a small vertical interval of wall rock. Drawing the probe up through the borehole yields a sonic log, or continuous velocity log (CVL), which is a record of velocity varia- In the laboratory, velocities are determined by meas- uring the travel-time of high-frequency (about 1 MHz) acoustic pulses transmitted through cylindrical rock specimens. By this means, the effect on velocity of vary- ing temperature, con?ning pressure, pore ?uid pressure or composition may be quantitatively assessed. It is im- portant to note that laboratory measurements at low con?ning pressures are of doubtful validity.The intrinsic velocity of a rock is not normally attained in the labora- tory below a con?ning pressure of about 100 MPa (megapascals), or 1 kbar, at which pressure the original solid contact between grains characteristic of the pristine The following empirical ?ndings of velocity studies are noteworthy: 1. Compressional wave velocity increases with con?n- 2. Sandstone and shale velocities show a systematic increase with depth of burial and with age, due to the combined effects of progressive compaction and 3. For a wide range of sedimentary rocks the compres- sional wave velocity is related to density, and well- established velocity–density curves have been published Hence, the densities of inaccessible subsurface layers may be predicted if their velocity is known from seismic 4. The presence of gas in sedimentary rocks reduces the elastic moduli, Poisson's ratio and the v /v ratio. v /v ra- ps ps tios greater than 2.0 are characteristic of unconsolidated sand, whilst values less than 2.0 may indicate either a consolidated sandstone or a gas-?lled unconsolidated sand. The potential value of v in detecting gas-?lled s sediments accounts for the current interest in shear wave Typical compressional wave velocity values and ranges forawidevarietyof EarthmaterialsaregiveninTable3.1.

Elements of Seismic Surveying 27 v (km s-1) p Unconsolidated materials Sand (dry) 0.2–1.0 Sand (water-saturated) 1.5–2.0 Clay 1.0–2.5 Glacial till (water-saturated) 1.5–2.5 Permafrost 3.5–4.0 Sedimentary rocks Sandstones 2.0–6.0 Tertiary sandstone 2.0–2.5 Pennant sandstone (Carboniferous) 4.0–4.5 Cambrian quartzite 5.5–6.0 Limestones 2.0–6.0 Cretaceous chalk 2.0–2.5 Jurassic oolites and bioclastic limestones 3.0–4.0 Carboniferous limestone 5.0–5.5 Dolomites 2.5–6.5 Salt 4.5–5.0 Anhydrite 4.5–6.5 Gypsum 2.0–3.5 Igneous/Metamorphic rocks Granite 5.5–6.0 Gabbro 6.5–7.0 Ultrama?c rocks 7.5–8.5 Serpentinite 5.5–6.5 Pore ?uids Air 0.3 Water 1.4–1.5 Ice 3.4 Petroleum 1.3–1.4 Other materials Steel 6.1 Iron 5.8 Aluminium 6.6 Concrete 3.6

3.5 Attenuation of seismic energy along ray paths As a seismic pulse propagates in a homogeneous ma- terial, the original energy E transmitted outwards from the source becomes distributed over a spherical shell, the wavefront, of expanding radius. If the radius of the wave- front is r, the amount of energy contained within a unit area of the shell is E/4pr2.With increasing distance along a ray path, the energy contained in the ray falls off as r-2 due to the effect of the geometrical spreading of the energy.

Wave amplitude, which is proportional to the square A further cause of energy loss along a ray path arises because, even at the low strains involved, the ground is imperfectly elastic in its response to the passage of seismic waves. Elastic energy is gradually absorbed into the medium by internal frictional losses, leading eventually to the total disappearance of the seismic disturbance.The mechanisms for the absorption of energy are complex, but the loss of energy is usually regarded as being a ?xed proportion of the total energy, for each oscillation of the rock particles involved, during which time the wave- front will have moved forward one wavelength. The ab- sorption coef?cient a expresses the proportion of energy lost during transmission through a distance equivalent to a complete wavelength l.Values of a for common Earth materials range from 0.25 to 0.75 dB l-1 (for a de?nition Over the range of frequencies used in seismic survey- ing the absorption coef?cient is normally assumed to be independent of frequency. If the amount of absorp- tion per wavelength is constant, it follows that higher frequency waves attenuate more rapidly than lower frequency waves as a function of time or distance. To illustrate this point, consider two waves with frequencies of 10 Hz and 100 Hz to propagate through a rock in which v = 2.0 km s-1 and a = 0.5 dB l-1. The 100 Hz p wave (l = 20 m) will be attenuated due to absorption by 5 dB over a distance of 200 m, whereas the 10 Hz wave (l = 200 m) will be attenuated by only 0.5 dB over the same distance. The shape of a seismic pulse with a broad fre- quency content therefore changes continuously during propagation due to the progressive loss of the higher fre- quencies. In general, the effect of absorption is to pro- 3.7). This effect of absorption is a familiar experience as it applies to P-waves in air, sound. The sharp crack of a nearby lightning ?ash is heard far away as the distant `rumble' of thunder.

3.6 Ray paths in layered media At an interface between two rock layers there is gen- erally a change of propagation velocity resulting from At such an interface, the energy within an incident seis- mic pulse is partitioned into transmitted and re?ected pulses. The relative amplitudes of the transmitted and re?ected pulses depend on the velocities and densities Input spike

20 ms After 1 s After 2 s After 3 s After 4 s After 5 s Fig. 3.7 The progressive change of shape of an original spike pulse during its propagation through the ground due to the effects of absorption. (After Anstey 1977.)

of the two layers, and the angle of incidence on the 3.6.1 Re?ection and transmission of normally incident seismic rays Consider a compressional ray of amplitude A normally 0 incident on an interface between two media of differ- ing velocity and density (Fig. 3.8). A transmitted ray of amplitude A travels on through the interface in the 2 same direction as the incident ray and a re?ected ray of amplitude A returns back along the path of the 1 The total energy of the transmitted and re?ected rays must equal the energy of the incident ray. The relative proportions of energy transmitted and re?ected are de- termined by the contrast in acoustic impedance Z across the interface. The acoustic impedance of a rock is the prod- uct of its density (r) and its wave velocity(v); that is, Z = rv

Incident ray, amplitude A0 Elements of Seismic Surveying 29

Reflected ray, v1, ?1 amplitude A1 Transmitted ray, v2, ?2 amplitude A2 ?2v2 =/ ?1v1 Fig. 3.8 Re?ected and transmitted rays associated with a ray normally incident on an interface of acoustic impedance contrast.

interface, and more energy is re?ected the greater the contrast. From common experience with sound, the best echoes come from rock or brick walls. In terms of physical theory, acoustic impedance is closely analogous to elec- trical impedance and, just as the maximum transmission of electrical energy requires a matching of electrical im- pedances, so the maximum transmission of seismic ener- The re?ection coef?cient R is a numerical measure of the effect of an interface on wave propagation, and is calcu- lated as the ratio of the amplitude A of the re?ected ray 1 to the amplitude A of the incident ray 0

R=A A 10 To relate this simple measure to the physical properties of the materials at the interface is a complex problem. As we have already seen, the propagation of a P-wave de- pends on the bulk and shear elastic moduli, as well as the density of the material. At the boundary the stress and strain in the two materials must be considered. Since the materials are different, the relations between stress and strain in each will be different. The orientation of stress and strain to the interface also becomes important. The formal solution of this physical problem was derived early in the 20th century, and the resulting equations are named the Zoeppritz equations (Zoeppritz 1919; and for explanation of derivations see Sheriff & Geldart 1982). Here, the solutions of these equations will be accepted. For a normally incident ray the relationships are fairly simple, giving:

rv-rv Z-Z R= 22 11= 2 1 rv+rv Z+Z 22 11 2 1 where r , v , Z and r , v , Z are the density, P-wave 111 222 velocity and acoustic impedance values in the ?rst and second layers, respectively. From this equation it follows that -1 £ R £ +1. A negative value of R signi?es a phase The transmission coef?cient T is the ratio of the ampli- tude A of the transmitted ray to the amplitude A of the 20 incident ray

T=A A 20 For a normally incident ray this is given, from solution of Zoeppritz's equations, by

2Z T= 1 Z +Z 21 Re?ection and transmission coef?cients are some- times expressed in terms of energy rather than wave am- plitude. If energy intensity I is de?ned as the amount of energy ?owing through a unit area normal to the direc- tion of wave propagation in unit time, so that I , I and I 01 2 are the intensities of the incident, re?ected and trans- mitted rays respectively, then

where R¢ and T¢ are the re?ection and transmission coef- This is the case when there is no contrast of acoustic im- pedance across an interface, even if the density and ve- 12 A good approximation to this situation occurs at the free surface of a water layer: rays travelling upwards from an explosion in a water layer are almost totally re?ected back from the water surface with a phase change of p (R Values of re?ection coef?cient R for interfaces be- tween different rock types rarely exceed ±0.5 and are typically much less than ±0.2.Thus, normally the bulk of seismic energy incident on a rock interface is transmitted and only a small proportion is re?ected. By use of an em- pirical relationship between velocity and density (see also Section 6.9), it is possible to estimate the re?ection coef?cient from velocity information alone (Gardner et al. 1974, Meckel & Nath 1977):

R=0.625ln(v1v2) Such relationships can be useful, but must be applied with caution since rock lithologies are highly variable and laterally heterogeneous as pointed out in Section 3.4.

3.6.2 Re?ection and refraction of obliquely incident rays When a P-wave ray is obliquely incident on an inter- face of acoustic impedance contrast, re?ected and trans- mitted P-wave rays are generated as in the case of normal incidence. Additionally, some of the incident com- pressional energy is converted into re?ected and trans- mitted S-wave rays (Fig. 3.9) that are polarized in a vertical plane. Zoeppritz's equations show that the am- plitudes of the four phases are a function of the angle of incidence q. The converted rays may attain a signi?cant magnitude at large angles of incidence. Detection and identi?cation of converted waves can be dif?cult in seis- mic surveys, but they do have potential to provide more constraints on the physical properties of the media at the interface. Here consideration will be con?ned to the In the case of oblique incidence, the transmitted P- wave ray travels through the lower layer with a changed direction of propagation (Fig. 3.10) and is referred to as a refracted ray.The situation is directly analogous to the be- Reflected S

Incident P ? Reflected P v1 v2 > v1 Refracted P

Refracted S Fig. 3.9 Re?ected and refracted P- and S-wave rays generated by a P-wave ray obliquely incident on an interface of acoustic impedance contrast.

Incident P ?1 ?1 Reflected P

v1 v2 > v1 Refracted P ?2 Fig. 3.10 Re?ected and refracted P-wave rays associated with a P-wave rays obliquely incident on an interface of acoustic impedance contrast.

haviour of a light ray obliquely incident on the boundary between, say, air and water and Snell’s Law of Refraction applies equally to the optical and seismic cases. Snell de- ?ned the ray parameter p = sin i/v, where i is the angle of inclination of the ray in a layer in which it is travelling with a velocity v. The generalized form of Snell’s Law states that, along any one ray, the ray parameter remains a For the refracted P-wave ray shown in Fig. 3.10, therefore

sin q v 1=1 sin q v 22 Note that if v > v the ray is refracted away from the 21 normal to the interface; that is, q > q . Snell’s Law also 21 applies to the re?ected ray, from which it follows that the 3.10).

3.6.3 Critical refraction When the velocity is higher in the underlying layer there is a particular angle of incidence, known as the critical angle q , for which the angle of refraction is 90°. This c gives rise to a critically refracted ray that travels along the interface at the higher velocity v . At any greater angle of 2 incidence there is total internal re?ection of the incident energy (apart from converted S-wave rays over a further range of angles).The critical angle is given by

sin q sin 90? 1 c= = vvv 122 so that q = sin – 1 ( ) c v1 v 2

The passage of the critically refracted ray along the top of the lower layer causes a perturbation in the upper layer that travels forward at the velocity v , which is greater 2 than the seismic velocity v of that upper layer. The 1 situation is analogous to that of a projectile travelling through air at a velocity greater than the velocity of sound in air and the result is the same, the generation of a shock wave.This wave is known as a head wave in the seis- Elements of Seismic Surveying 31

mic case, and it passes up obliquely through the upper layer towards the surface (Fig. 3.11). Any ray associated with the head wave is inclined at the critical angle i . By c means of the head wave, seismic energy is returned to the surface after critical refraction in an underlying layer of higher velocity.

3.6.4 Diffraction In the above discussion of the re?ection and transmission of seismic energy at interfaces of acoustic impedance contrast it was implicitly assumed that the interfaces were continuous and approximately planar. At abrupt discontinuities in interfaces, or structures whose radius of curvature is shorter than the wavelength of incident waves, the laws of re?ection and refraction no longer apply. Such phenomena give rise to a radial scattering of incident seismic energy known as diffraction. Common sources of diffraction in the ground include the edges of faulted layers (Fig. 3.12) and small isolated objects, such as boulders, in an otherwise homogeneous layer.

Ray paths Head wave generated in overlying layer ?c

v1 A B v2 > v1

Wavefront expanding in lower layer Fig. 3.11 Generation of a head wave in the upper layer by a wave propagating through the lower layer.

(a) (b) t 1 v 1 1 v 2 xx x crit cros (c) Fig. 3.13 (a) Seismogram showing the S D output traces of 24 geophones distributed X along the Earth’s surface as a function of time. (b)Travel-time curves for direct, z re?ected and refracted rays in the case of a v1 simple two-layer model. (c) Direct, re?ected and refracted ray paths from a near surface source to a surface detector v2 in the case of a simple two-layer model.

Diffracted phases are commonly observed in seismic recordings and are sometimes dif?cult to discriminate from re?ected and refracted phases, as discussed in Chapter 4.

refracted ray travels obliquely down to the interface at ve- locity v , along a segment of the interface at the higher 1 21 The travel time of a direct ray is given simply by t =xv dir 1

which de?nes a straight line of slope l/v passing through 1 The travel time of a re?ected ray is given by

12 (x2 + 4z2) t= refl v 1 which, as discussed in Chapter 4, is the equation of an The travel time of a refracted ray (for derivation see Chapter 5) is given by

x 2z cosq t=+ c refr vv 21 which is the equation of a straight line having a slope of l/v and an intercept on the time axis of 2

2z cosq c v 1 Travel-time curves, or time–distance curves, for 3.13. By suitable analysis of the travel-time curve for re?ected or refracted rays it is possible to compute the depth to the underlying layer. This provides two inde- pendent seismic surveying methods for locating and mapping subsurface interfaces, re?ection surveying and refraction surveying. These have their own distinctive methodologies and ?elds of application and they are dis- cussed separately in detail in Chapters 4 and 5. However, some general remarks about the two methods may be made here with reference to the travel-time curves and seismogram of Fig. 3.13. The curves are more compli- cated in the case of a multilayered model, but the follow- The ?rst arrival of seismic energy at a surface detector offset from a surface source is always a direct ray or a re- fracted ray. The direct ray is overtaken by a refracted ray at the crossover distance x . Beyond this offset distance cros the ?rst arrival is always a refracted ray. Since critically re- Elements of Seismic Surveying 33

fracted rays travel down to the interface at the critical angle there is a certain distance, known as the critical dis- tance x , within which refracted energy will not be re- cr it turned to the surface. At the critical distance, the travel times of re?ected rays and refracted rays coincide because they follow effectively the same path. Re?ected rays are never ?rst arrivals; they are always preceded by direct rays The above characteristics of the travel-time curves determine the methodology of refraction and re?ection surveying. In refraction surveying, recording ranges are chosen to be suf?ciently large to ensure that the cross- over distance is well exceeded in order that refracted rays Indeed, some types of refraction survey consider only these ?rst arrivals, which can be detected with unsophis- ticated ?eld recording systems. In general, this approach means that the deeper a refractor, the greater is the range over which recordings of refracted arrivals need to be In re?ection surveying, by contrast, re?ected phases are sought that are never ?rst arrivals and are normally of very low amplitude because geological re?ectors tend to have small re?ection coef?cients. Consequently, re?ec- tions are normally concealed in seismic records by higher amplitude events such as direct or refracted body Re?ection surveying methods therefore have to be capable of discriminating between re?ected energy and many types of synchronous noise. Recordings are nor- mally restricted to small offset distances, well within the critical distance for the re?ecting interfaces of main interest. However, in multichannel re?ection surveying recordings are conventionally taken over a signi?cant range of offset distances, for reasons that are discussed fully in Chapter 4.

· record and display the seismic waveforms on a suitable The general methodology of examining hidden structures by studying their effects on arti?cially gener- ated acoustic or seismic waves has an enormously wide range of applications covering a wide range of spatial scales. Perhaps the smallest scale is ultrasound imaging in medicine, which can also be applied industrially to examining engineering structures. Within the more geophysical applications, the scales range from depths of a metre or less in engineering, environmental or archae- ological surveys to tens of kilometres for crustal and For each application there is a limit to the smallest structures that can be detected, known as the resolution of the survey.The resolution is basically determined by the pulse length: for a pulse of any particular length there is a minimum separation below which the pulses will overlap in time in the seismic recording. Although the pulse length may be shortened at the processing stage by deconvolution (see Section 4.8.2), this is only possible if the data are of good quality, and is a complement to, not a substitute for, good survey design. The pulse width is determined by both the maximum frequency Since Earth materials absorb seismic energy in a fre- quency-selective way (Section 3.5), the optimum wave- form will be speci?c to each survey. It is an important characteristic of all geophysical surveys, and particularly seismic ones, that they must be designed individually for each speci?c case. The general aspects of the equipment used for seismic surveys are reviewed here; speci?c variations for re?ection and refraction surveying are described in Chapters 4 and 5.

3.8.1 Seismic sources and the seismic/acoustic spectrum A seismic source is a localized region within which the sudden release of energy leads to a rapid stressing of the surrounding medium. The archetypal seismic source is an explosion. While explosives are still used, there is an increasing number of more sophisticated and ef?cient The main requirements of the seismic source are: · Suf?cient energy across the broadest possible fre- quency range, extending up to the highest recordable · Energy should be concentrated in the type of wave energy which is required for a speci?c survey, either P- wave or S-wave, and generate minimum energy of other wave types. Such other unwanted energy would degrade · The source waveform must be repeatable. Seismic sur- veys almost always involve comparing the seismograms Variations on the seismograms should be diagnostic of the ground structure, not due to random variations of · The source must be safe, ef?cient, and environmen- tally acceptable. Most seismic surveys are commercial operations which are governed by safety and environ- mental legislation.They must be as cost-effective as pos- sible. Sometimes the requirements for ef?ciency lead to higher safety and environmental standards than legally enforced.Whether involving personal injury or not, ac- cidents are referred to as `lost-time incidents’. Safety aids ef?ciency as well as being desirable from many other The complete seismic/acoustic spectrum is shown in Fig. 3.14.There is a very wide variety of seismic sources, characterized by differing energy levels and frequency characteristics. In general, a seismic source contains a wide range of frequency components within the range from 1 Hz to a few hundred hertz, though the energy is Source characteristics can be modi?ed by the use of several similar sources in an array designed, for example, to improve the frequency spectrum of the transmitted pulse. This matter is taken up in Chapter 4 when dis- cussing the design parameters of seismic re?ection surveys.

Elements of Seismic Surveying 35 Echo sounders Pingers Boomers Sparkers Air guns Vibroseis Quarry blasts Earthquake body waves Earthquake surface waves

10–2 10–1 1 101 102 103 104 105 Frequency (Hz) (log scale)

and usually three, of the basic requirements for modern surveys, their use is steadily declining and limited to lo- cations where alternative sources cannot be used.

Vibroseis sweep signal Reflection from base of layer 1 Reflection from base of layer 2

Reflection from base of layer 3 (phase inverted) Field recording (superposition of above reflections)

Output trace resulting from correlation of field recording with sweep signal t=0 Time

Fig. 3.15 Cross-correlation of aVibroseis® seismogram with the input sweep signal to locate the positions of occurrence of re?ected arrivals.

Firing pin pulse-encoded signal from buried interfaces. Peaks in the cross-correlation function reveal the positions of Weight drops and hammers. Perhaps the simplest land seismic source is a large mass dropped on to the ground surface. Weight drops have been manufactured in a wide variety of forms from eight-wheel trucks dropping a weight of several tonnes, to a single person with a sledgehammer. If the source energy required is relatively low, these types of sources can be fast and ef?cient. The horizontal impact of a weight or hammer on to one side of a vertical plate partially embedded in the ground can be used as a source for shear wave Shotguns, buffalo guns and ri?es. One solution to gaining additional energy for small-scale surveys is to use the Ri?es have been used as seismic sources by ?ring the bullet into the ground. While effective as a very high- frequency source, this is banned by legislation in many countries. An alternative is to ?re a blank shotgun car- tridge in a hole using a suitable device, commonly termed a buffalo gun (Fig 3.16).The blank shotgun car- tridge offers an impulsive source giving considerably more energy than a sledgehammer, with few of the safety problems of explosives.

Ground surface Auger flight Firing chamber Fig. 3.16 Schematic cross-section of a typical buffalo gun.The cartridge is ?red by dropping a simple ?ring pin on to the cartridge.

(a) Variable chamber size (b) Air Elements of Seismic Surveying 37

Water Fig. 3.17 Schematic cross-sections through (a) a Bolt air gun and (b) a Sodera water gun to illustrate the principles of operation. (Redrawn with permission of Bolt Associates and Sodera Ltd.)

(a) Single 270 in3 air gun 0 (b) 0.1 0.2 0.3 0.4 0.5 0.6

Seven – gun array (1222 in3 total volume) water.When the piston stops, a vacuum cavity is created behind the advancing water jet and this implodes under the in?uence of the ambient hydrostatic pressure, gener- Since the implosion represents collapse into a vacuum, no gaseous material is compressed to ` bounce back’ as a bubble pulse. The resulting short pulse length offers a potentially higher resolution than is achieved with air guns but at the expense of a more complex initial source Several marine sources utilize explosive mixtures of gases, but these have not achieved the same safety and reliability, and hence industry acceptance, as air guns. In sleeve exploders, propane and oxygen are piped into a sub- merged ?exible rubber sleeve where the gaseous mix- 0.7 0.8 s

Fig. 3.18 Comparison of the source signatures of (a) a single air gun (peak pressure: 4.6 bar metres) and (b) a seven- Note the effective suppression of bubble pulses in the latter case. (Redrawn with permission of Bolt Associates.)

behind the survey vessel. Operating voltages are typi- This electrical discharge leads to the formation and rapid growth of a plasma bubble and the consequent genera- tion of an acoustic pulse. For safety reasons, sparkers are Boomers comprise a rigid aluminium plate attached below a heavy-duty electrical coil by a spring-loaded mounting. A capacitor bank is discharged through the coil and the electromagnetic induction thus generated forces the aluminium plate rapidly downwards, setting up a compressional wave in the water. The device is typically towed behind the survey vessel in a catamaran Sparkers and boomers generate broad-band acoustic pulses and can be operated over a wide range of energy levels so that the source characteristics can to some ex- tent be tailored to the needs of a particular survey. In gen- eral, boomers offer better resolution (down to 0.5 m) but more restricted depth penetration (a few hundred metres Pingers consist of small ceramic piezoelectric transducers, mounted in a towing ?sh, which, when ac- tivated by an electrical impulse, emit a very short, high- frequency acoustic pulse of low energy.They offer a very high resolving power (down to 0.1 m) but limited pene- tration (a few tens of metres in mud, much less in sand or rock). They are useful in offshore engineering applica- tions such as surveys of proposed routes for submarine Chirp systems are electro-mechanical transducers that produce an extended, repeatable, source waveform which allows greater energy output. This longer signal can be compressed in processing to give greater resolu- Further discussion of the use of air guns, sparkers, boomers and pingers in single-channel seismic re?ec- tion pro?ling systems is given in Section 4.15.

3.8.2 Seismic transducers Conversion of the ground motion to an electrical signal requires a transducer which is sensitive to some compo- nent of the ground motion, and can record the required The ?rst issue is which component of the motion to measure. As the ground oscillates, it is possible to mea- sure either the displacement, velocity or acceleration of the ground particles as the wave passes.The ground mo- tion also takes place in three dimensions. To record it Elements of Seismic Surveying 39

Shunt resistor Leaf spring supportCoil output to seismic amplifier

Coil N S N Permanent magnet Fig. 3.19 Schematic cross-section through a moving-coil geophone.

Resonant frequency (a) (b) 180 0.80 150 – 1 ) 0.40h = 0.2 120 s –1 0.20h = 0.5 90 h = 0.7 60 Output phase Output (V cm 0.10 30

0 0 10 20 30 40 Frequency (Hz) Resonant frequency

10 20 30 40 Frequency (Hz) h = 0.7 h = 0.5 h = 0.2 Fig. 3.20 Amplitude and phase responses of a geophone with a resonant frequency of 7 Hz, for different damping factors h. Output phase is expressed relative to input phase. (AfterTelford et al. 1976.)

hydrophones are made up into hydrophone streamers by distributing them along an oil-?lled plastic tube. The tube is arranged to have neutral buoyancy and is manu- factured from materials with an acoustic impedance close to that of water to ensure good transmission of seismic energy to the hydrophone elements. Since piezoelectric elements are also sensitive to accelerations, hydrophones are often composed of two elements mounted back to back and connected in series so that the effects of accelerations of the streamer as it is towed through the water are cancelled out in the hydrophone outputs. The response of each element to pressure change is, however, unaffected and the seismic signal is Arrays of geophones or hydrophones may be con- nected together into linear or areal arrays containing tens or even hundreds of transducers whose individual out- puts are summed. Such arrays provide detectors with a directional response that facilitates the enhancement of signal and the suppression of certain types of noise as dis- cussed further in Chapter 4.

3.8.3 Seismic recording systems Recording a seismogram is a very dif?cult technical operation from at least three key aspects: 1. The recording must be timed accurately relative to 2. Seismograms must be recorded with multiple trans- ducers simultaneously, so that the speed and direction of The least dif?cult of these problems is the timing. For nearly all seismic surveys, times need to be accurate to better than one thousandth of a second (one milli- second). For very small-scale surveys the requirement may be for better than 0.1 ms. In fact, with modern electronics, measuring such short time intervals is not dif?cult. Usually the biggest uncertainty is in deciding how to measure the instant when the seismic source started the wave. Even in a simple case, as for a sledge- hammer hitting the ground, is the correct instant when the hammer ?rst hits the ground, or when it stops compressing the ground and a seismic wave radiates out- ward?The ?rst is easy to measure, the second is probably more important, and they are usually separated by more In order to determine the subsurface path of the seis- mic energy, the direction from which the wave arrives at the surface must be determined.This is achieved by hav- Elements of Seismic Surveying 41

and potentially reveals more information about the Distributed systems In seismic surveying the outputs of several detectors are fed to a multichannel recording system mounted in a recording vehicle. The individual detector outputs may be fed along a multicore cable.The weight and complex- ity of multicore cables becomes prohibitive as the num- ber of channels of data rises into the hundreds. Modern systems distribute the task of ampli?cation, digitization and recording of data from groups of detectors to indi- vidual computer units left unattended in the ?eld.These are connected together to make a ?eld computer net- work using lightweight ?bre-optic cables or telemetry links. The separate units can then be controlled by a central recording station, and upload their digital seis- mograms to it on command.

Problems 1. How does the progressive loss of higher fre- quencies in a propagating seismic pulse lead to 2. A 10 Hz seismic wave travelling at 5 km s-1 propagates for 1000 m through a medium with an absorption coef?cient of 0.2 dB l-1. What is the wave attenuation in decibels due solely to 3. A wave component with a wavelength of 100 m propagates through a homogeneous medium from a seismic source at the bottom of a borehole. Between two detectors, located in boreholes at radial distances of 1 km and 2 km from the source, the wave amplitude is found to be attenuated by 10 dB. Calculate the contribu- tion of geometrical spreading to this value of at- tenuation and, thus, determine the absorption coef?cient of the medium.

4. What is the crossover distance for direct and critically refracted rays in the case of a horizontal interface at a depth of 200 m separating a top layer of velocity 3.0 km s-1 from a lower layer of 5. A seismic pulse generated by a surface source is returned to the surface after re?ection at the tenth of a series of horizontal interfaces, each of which has a re?ection coef?cient R of 0.1. What is the attenuation in amplitude of the pulse caused by energy partitioning at all interfaces 6. At what frequency would a 150 Hz signal be recorded by a digital recording system with a sampling rate of 100 Hz?

Further reading Al-Sadi, H.N. (1980) Seismic Exploration. Birkhauser Verlag, Anstey, N.A. (1981) Seismic Prospecting Instruments. Vol 1: Signal Characteristics and Instrument Speci?cations. Gebruder Born- Dobrin, M.B. & Savit, C.H. (1988) Introduction to Geophysical Gregory, A.R. (1977) Aspects of rock physics from laboratory and log data that are important to seismic interpretation. In: Payton, C.E. (ed.), Seismic Stratigraphy — Applications to Hydrocarbon Exploration. Memoir 26, American Association of Petroleum Geologists,Tulsa.

SEG (1997) Digital Tape Standards (SEG-A, SEG-B, SEG-C, SEG-Y and SEG-D formats plus SEG-D rev 1&2). Compiled by SEG Technical Standards Committee. Society of Exploration Sheriff, R.E. & Geldart. L.P. (1982) Exploration Seismology Vol 1: History.Theory and DataAcquisition. Cambridge University Press, Sheriff, R.E. & Geldart, L.P. (1983) Exploration Seismology Vol 2: Data-Processing and Interpretation. Cambridge University Press, Waters, K.H. (1978) Re?ection Seismology — A Tool For Energy Re- Zoeppritz, K, 1919. Uber re?exion und durchgang seismischer wellen durch Unstetigkerls?aschen. Berlin, Uber Erdbeben- wellen VII B, Nachrichten der Koniglichen Gesellschaft der Wissensschaften zu Gottingen, math-phys. Kl. pp. 57–84.

4.1 Introduction Seismic re?ection surveying is the most widely used and well-known geophysical technique.The current state of sophistication of the technique is largely a result of the enormous investment in its development made by the hydrocarbon industry, coupled with the development of advanced electronic and computing technology. Seismic sections can now be produced to reveal details of geo- logical structures on scales from the top tens of metres of drift to the whole lithosphere. Part of the spectacu- lar success of the method lies in the fact that the raw data are processed to produce a seismic section which is an image of the subsurface structure. This also provides a trap for the unwary, since the seismic section is similar to, but fundamentally different from, a depth section of the geology. Only by understanding how the re?ection method is used and seismic sections are created, can the geologist make informed interpretations. This chapter provides the essential knowledge and understanding to support interpretation of seismic re?ection data. It builds up systematically from the basics of seismic wave re- ?ection from rock layers, and refers back to relevant material in Chapters 2 and 3.

4.2 Geometry of re?ected ray paths In seismic re?ection surveys seismic energy pulses are re?ected from subsurface interfaces and recorded at near-normal incidence at the surface. The travel times are measured and can be converted into estimates of depths to the interfaces. Re?ection surveys are most commonly carried out in areas of shallowly dipping sedimentary sequences. In such situations, velocity varies as a function of depth, due to the differing physical properties of the individual layers. Velocity may also vary horizontally, due to lateral lithological changes within the individual layers. As a ?rst approxi- mation, the horizontal variations of velocity may be Figure 4.1 shows a simple physical model of horizontally-layered ground with vertical re?ected ray paths from the various layer boundaries. This model assumes each layer to be characterized by an interval veloc- ity v , which may correspond to the uniform velocity i within a homogeneous geological unit or the average velocity over a depth interval containing more than one unit. If z is the thickness of such an interval and t is the ii one-way travel time of a ray through it, the interval velocity is given by

z v=i it i The interval velocity may be averaged over several depth intervals to yield a time-average velocity or, simply, average velocityV .Thus the average velocity of the top n layers in Fig. 4.1 is given by

nn Âz Âvt i ii V = i =1 = i =1 nn Ât Ât ii i =1 i =1

or, if Z is the total thickness of the top n layers and T is nn the total one-way travel time through the n layers,

v 1 v 2 v 3 v n–1 v n Fig. 4.1 Vertical re?ected ray paths in a horizontally-layered ground.

?ector lying at a depth z beneath a homogeneous top layer of velocity V. The equation for the travel time t of the re?ected ray from a shot point to a detector at a hori- zontal offset, or shot–detector separation, x is given by the ratio of the travel path length to the velocity 12 t = (x2 + 4z2) V (4.1) In a re?ection survey, re?ection time t is measured at an offset distance x.These values can be applied to equation (4.1), but still leave two unknown values which are re- lated to the subsurface structure, z and V. If many re?ec- tion times t are measured at different offsets x, there will be enough information to solve equation (4.1) for both these unknown values. The graph of travel time of re?ected rays plotted against offset distance (the time– distance curve) is a hyperbola whose axis of symmetry is Substituting x = 0 in equation (4.1), the travel time t 0 of a vertically re?ected ray is obtained:

2z t = (4.2) 0 V (a) xx

z V (b)t t x ?T t 0

–x O + x x Fig. 4.2 (a) Section through a single horizontal layer showing the geometry of re?ected ray paths and (b) time–distance curve for re?ected rays from a horizontal re?ector. DT = normal moveout (NMO).

This is the intercept on the time axis of the time–distance curve (see Fig. 4.2(b)). Equation (4.1) can be written

4z2 x2 t 2 = + (4.3) V2 V2

Thus x2 t 2 = t 2 + (4.4) 0 V2

velocity is by considering the increase of re?ected travel time with offset distance, the moveout, as discussed Equation (4.3) can also be rearranged 12 12 È 2? 2 2z Ê x ^ È Ê x ^ ? Ë ¯ 0 Ë ¯ (4.5) V Î 2z ° Î Vt ° 0

This form of the equation is useful since it indicates clearly that the travel time at any offset x will be the vertical travel time plus an additional amount which 0 This relationship can be reduced to an even simpler form with a little more rearrangement. Using the standard binomial expansion of equation (4.5) gives È24 1Ê x ^ 1Ê x ^ ? 0 Î 2 ËVt0 ¯ 8 ËVt0 ¯ °

Remembering that t = 2z/V, the term x/Vt can be 00 written as x/2z. If x = z, the second term in this series becomes 1/8 of (1/2)4, i.e. 0.0078, which is less than a 1% change in the value of t. For small offset/depth ratios (i.e. x/z << 1), the normal case in re?ection surveying, this equation may be truncated after the ?rst term to obtain the approximation

È 1Ê x ^2? x2 t ª t Í1 + . ª t + (4.6) 0 Î 2 ËVt ¯ ° 2V 2t0 0 0

This is the most convenient form of the time–distance equation for re?ected rays and it is used extensively in the Moveout is de?ned as the difference between the travel times t and t of re?ected-ray arrivals recorded at 12 two offset distances x and x . Substituting t , x and t , x 1 2 11 22 in equation (4.6), and subtracting the resulting equations gives

x2 - x2 t -t ª 2 1 2 1 2V 2t 0 Normal moveout (NMO) at an offset distance x is the difference in travel time DT between re?ected arrivals at x and at zero offset (see Fig. 4.2)

x2 DT = t - t ª (4.7) x 0 2V 2t 0 Seismic Re?ection Surveying 45

Note that NMO is a function of offset, velocity and re- ?ector depth z (since z = Vt /2).The concept of move- 0 out is fundamental to the recognition, correlation and enhancement of re?ection events, and to the calculation of velocities using re?ection data. It is used explicitly or implicitly at many stages in the processing and interpre- As an important example of its use, consider the T –DT method of velocity analysis. Rearranging the terms of equation (4.7) yields

x V ª (4.8) 12 (2t DT ) 0

Using this relationship, the velocity V above the re?ector can be computed from knowledge of the zero-offset re?ection time (t ) and the NMO (DT ) at a particular 0 offset x. In practice, such velocity values are obtained by computer analysis which produces a statistical estimate based upon many such calculations using large numbers of re?ected ray paths (see Section 4.7). Once the velocity has been derived, it can be used in conjunction with t to compute the depth z to the re?ector using 0 0

4.2.2 Sequence of horizontal re?ectors In a multilayered ground, inclined rays re?ected from the nth interface undergo refraction at all higher inter- faces to produce a complex travel path (Fig. 4.3(a)). At offset distances that are small compared to re?ector depths, the travel-time curve is still essentially hyper- bolic but the homogeneous top layer velocity V in equa- tions (4.1) and (4.7) is replaced by the average velocityV or, to a closer approximation (Dix 1955), the root-mean- square velocity V of the layers overlying the re?ector. As r ms the offset increases, the departure of the actual travel- time curve from a hyperbola becomes more marked The root-mean-square velocity of the section of ground down to the nth interface is given by

12 nn È? V = Âv2t Ât Î i =1 i =1 °

where v is the interval velocity of the ith layer and t is the ii one-way travel time of the re?ected ray through the ith layer.

(a) (b)

v 1 v 2 zV v rms, n 3 v n–1 v n t x Actual curve

Hyperbolic curve 2 2 x2 t =t + 02 V rms Fig. 4.3 (a)The complex travel path of a re?ected ray through a multilayered ground, showing refraction at layer boundaries. (b)The time–distance curve for re?ected rays following such a travel path. Note that the divergence from the hyperbolic travel-time curve for a r ms

Thus at small offsets x (x << z), the total travel time t n of the ray re?ected from the nth interface at depth z is given to a close approximation by

12 t =(x2+4z2) V cf.equation(4.1) n rms

and the NMO for the nth re?ector is given by x2 DTn ª cf. equation (4.7) 2V 2 t rms,n 0

The individual NMO value associated with each re?ec- tion event may therefore be used to derive a root-mean- Values of V down to different re?ectors can then r ms be used to compute interval velocities using the Dix formula. To compute the interval velocity v for the nth n interval

12 ÈVr2ms,n tn -V n -1 tn -1 ? 2 rms, n Î t -t ° n n -1

where V 2 , t and V , t are, respectively, the rms,n-1 n-1 rms,n n root-mean-square velocity and re?ected ray travel times to the (n - 1)th and nth re?ectors (Dix 1955).

4.2.3 Dipping re?ector In the case of a dipping re?ector (Fig. 4.4(a)) the value of dip q enters the time–distance equation as an additional unknown. The equation is derived similarly to that for horizontal layers by considering the ray path length divided by the velocity:

12 (x2+4z2+4xzsinq) t = cf. equation (4.1) V

The equation still has the form of a hyperbola, as for the horizontal re?ector, but the axis of symmetry of the Proceeding as in the case of a horizontal re?ector, using a truncated binomial expansion, the following expres- sion is obtained:

(x2+4xzsinq) t ª t + (4.9) 0 2V 2t 0

Seismic Re?ection Surveying 47 (a) (b) xx Fig. 4.4 (a) Geometry of re?ected ray paths and (b) time–distance curve for re?ected rays from a dipping re?ector. DT = dip moveout. –x 0 +x x d z

t ?T d t x t –x V ? difference in travel times t and t of rays re?ected from x -x the dipping interface to receivers at equal and opposite offsets x and -x DT = t - t d x -x

Using the individual travel times de?ned by equation (4.9) DT = 2x sin q V d

Rearranging terms, and for small angles of dip (when sin q ª q ) q ª V DT 2x d

(a) Primary Double-path Near-surface Peg-leg multiple multiples multiple

(b) Short-path multiples extend pulse length Long-path multiples generate discrete pulse

Fig. 4.5 (a) Various types of multiple re?ection in a layered ground. (b)The difference between short-path and long-path multiples.

recorded pulse. Such multiples are known as short-path multiples (or short-period reverberations) and these may be contrasted with long-path multiples whose additional path length is suf?ciently long that the multiple re?ec- tion is a distinct and separate event in the seismic record Misidenti?cation of a long-path multiple as a primary event, for example, would lead to serious interpretation error. The arrival times of multiple re?ections are pre- dictable, however, from the corresponding primary re- ?ection times. Multiples can therefore be suppressed by suitable data processing techniques to be described later (Section 4.8).

4.3 The re?ection seismogram The graphical plot of the output of a single detector in a re?ection spread is a visual representation of the local pattern of vertical ground motion (on land) or pressure variation (at sea) over a short interval of time following the triggering of a nearby seismic source. This seismic trace represents the combined response of the layered ground and the recording system to a seismic pulse. Any display of a collection of one or more seismic traces is termed a seismogram. A collection of such traces repre- senting the responses of a series of detectors to the energy from one shot is termed a shot gather. A collection of the traces relating to the seismic response at one surface mid-point is termed a common mid-point gather (CMP gather). The collection of the seismic traces for each CMP and their transformation to a component of the image presented as a seismic section is the main task of seismic re?ection processing.

4.3.1 The seismic trace At each layer boundary a proportion of the incident en- The proportion is determined by the contrast in acoustic impedances of the two layers, and for a vertically travel- ling ray, the re?ection coef?cient can be simply calcul- ated (see Section 3.6). Figure 4.6 shows the relationship of the geological layering, the variation in acoustic im- pedance and the re?ection coef?cients as a function of depth. The detector receives a series of re?ected pulses, scaled in amplitude according to the distance travelled and the re?ection coef?cients of the various layer boundaries.The pulses arrive at times determined by the depths to the boundaries and the velocities of propaga- Assuming that the pulse shape remains unchanged as it propagates through such a layered ground, the resultant seismic trace may be regarded as the convolution of the input pulse with a time series known as a re?ectivity func- tion composed of a series of spikes. Each spike has an am- plitude related to the re?ection coef?cient of a boundary and a travel time equivalent to the two-way re?ection time for that boundary. This time series represents the impulse response of the layered ground (i.e. the output for a spike input). The convolution model is illustrated schematically in Fig. 4.6. Since the pulse has a ?nite length, individual re?ections from closely-spaced boundaries are seen to overlap in time on the resultant seismogram.

Seismic Re?ection Surveying 49 Depth Fig. 4.6 The convolutional model of the re?ection seismic trace, showing the trace as the convolved output of a re?ectivity function with an input pulse, and the relationship of the re?ectivity function to the physical properties of the geological layers.

(a) Detectors Central shot Geological Acoustic Reflection section impedance coefficient log log Reflectivity function Time vvv vv vvv

(b) Detectors End shot Input Seismic * pulse = trace

Detectors Fig. 4.7 Shot–detector con?gurations used in multichannel seismic re?ection pro?ling. (a) Split spread, or straddle spread. (b) Single- ended or on-end spread.

In practice, as the pulse propagates it lengthens due to the progressive loss of its higher frequency components by absorption. The basic re?ection seismic trace may then be regarded as the convolution of the re?ectivity function with a time-varying seismic pulse. The trace will be further complicated by the superposition of various types of noise such as multiple re?ections, direct and re- fracted body waves, surface waves (ground roll), air waves and coherent and incoherent noise unconnected with the seismic source. In consequence of these several effects, seismic traces generally have a complex appear- ance and re?ection events are often not recogniz- able without the application of suitable processing In seismic re?ection surveying, the seismic traces are recorded, and the purpose of seismic processing can be viewed as an attempt to reconstruct the various columns of Fig. 4.6, moving from right to left.This will involve: · removing noise · determining the input pulse and removing that to give the re?ectivity function · determining the velocity function to allow conversion from time to depth axis · determination of the acoustic impedances (or related properties) of the formations.

particular shot. In shot gathers, the seismic traces are plotted side by side in their correct relative positions and the records are commonly displayed with their time axes arranged vertically in a draped fashion. In these seismic records, recognition of re?ection events and their corre- lation from trace to trace is much assisted if one half of the normal `wiggly-trace' waveform is blocked out. Figure 4.8 shows a draped section with this mode of display, de- rived from a split-spread multichannel survey. A short time after the shot instant the ?rst arrival of seismic ener- gy reaches the innermost geophones (the central traces) and this energy passes out symmetrically through the two arms of the split spread.The ?rst arrivals are followed by a series of re?ection events revealed by their hyper- bolic moveout.

4.3.3 The CMP gather Each seismic trace has three primary geometrical factors which determine its nature.Two of these are the shot po- sition and the receiver position. The third, and perhaps most critical, is the position of the subsurface re?ection point. Before seismic processing this position is un- known, but a good approximation can be made by as- suming this re?ection point lies vertically under the position on the surface mid-way between the shot and Older terminology is to call this point the depth point, but the former term is a description of what the position is, rather than what it is wished to represent, and is hence preferred. Collecting all the traces with a common mid- The seismic industry and the literature use the older term common depth point (CDP) interchangeably for The CMP gather lies at the heart of seismic processing for two main reasons: 1. The simple equations derived in Section 4.2 assume horizontal uniform layers. They can be applied with less error to a set of traces that have passed through the same geological structure. The simplest approximation to such a set of traces is the CMP gather. In the case of horizontal layers, re?ection events on each CMP gather are re?ected from a common depth point (CDP – see Fig. 4.9(a)). For these traces, the variation of travel time with offset, the moveout, will depend only on the veloc- ity of the subsurface layers, and hence the subsurface 2. The re?ected seismic energy is usually very weak. It is imperative to increase the signal-to-noise ratio of most Fig. 4.8 A draped seismic record of a shot gather from a split spread (courtesy Prakla-Seismos GmbH). Sets of re?ected arrivals from individual interfaces are recognizable by the characteristic hyperbolic alignment of seismic pulses.The late-arriving, high- amplitude, low-frequency events, de?ning a triangular-shaped central zone within which re?ected arrivals are masked, represent surface waves (ground roll).These latter waves are a typical type of coherent noise.

Seismic Re?ection Surveying 51 Shot – detector (a) midpoint S S S S CMP D D D D 43211234

CDP (b) S S S S CMP D D D D 43211234

Fig. 4.9 Common mid-point (CMP) re?ection pro?ling. (a) A set of rays from different shots to detectors re?ected off a common depth point (CDP) on a horizontal re?ector. (b)The common depth point is not achieved in the case of a dipping re?ector.

random and coherent noise. Combining all the traces in a CMP together will average out the noise, and increase the signal-to-noise ratio (SNR). This process is termed Strictly, the common mid-point principle breaks down in the presence of dip because the common depth point then no longer directly underlies the shot– detector mid-point and the re?ection point differs for rays travelling to different offsets (see Fig. 4.9(b)). Never- theless, the method is suf?ciently robust that CMP stacks almost invariably result in marked improvements In two-dimensional CMP surveying, known as CMP pro?ling, the re?ection points are all assumed to lie in three-dimensional surveying, the re?ection points are distributed across an area of any subsurface re?ector, and the CMP is de?ned as a limited area on the surface.

across a series of closely-spaced survey lines or around a grid of lines. However, as discussed later in Section 4.10, three-dimensional surveys provide a much better means of mapping three-dimensional structures and, in areas of structural complexity, they may provide the only means Re?ection pro?ling is normally carried out along pro?le lines with the shot point and its associated spread of detectors being moved progressively along the line to build up lateral coverage of the underlying geological section. This progression is carried out in a stepwise fashion on land but continuously, by a ship under way, at The two most common shot–detector con?gurations in multichannel re?ection pro?ling surveys are the split 4.7), where the number of detectors in a spread may be several hundred. In split spreads, the detectors are dis- tributed on either side of a central shot point; in single- ended spreads, the shot point is located at one end of the detector spread. Surveys on land are commonly carried out with a split-spread geometry, but in marine re?ec- tion surveys single-ended spreads are the normal con- ?guration due to the constraint of having to tow equipment behind a ship. The marine source is towed close behind the ship, with the hydrophone streamer (which may be several kilometers long) trailing behind.

4.4.1 Vertical and horizontal resolution Re?ection surveys are normally designed to provide a speci?ed depth of penetration and a particular degree of resolution of the subsurface geology in both the vertical and horizontal dimensions. The vertical resolution is a measure of the ability to recognize individual, closely- spaced re?ectors and is determined by the pulse length on the recorded seismic section. For a re?ected pulse represented by a simple wavelet, the maximum resolu- tion possible is between one-quarter and one-eighth of the dominant wavelength of the pulse (Sheriff & Gel- dart 1983). Thus, for a re?ection survey involving a sig- nal with a dominant frequency of 50 Hz propagating in sedimentary strata with a velocity of 2.0 km s-1, the dominant wavelength would be 40 m and the vertical This ?gure is worth noting since it serves as a reminder that the smallest geological structures imaged on seismic sections tend to be an order of magnitude larger than the Since deeper-travelling seismic waves tend to have a x

0.5x Fig. 4.10 The horizontal sampling of a seismic re?ection survey is half the detector spacing.

Source ? 4 Reflector Fresnel zone Fig. 4.11 Energy is returned to source from all points of a re?ector.The part of the re?ector from which energy is returned within half a wavelength of the initial re?ected arrival is known as the Fresnel zone.

structively to build up the re?ected signal, and the part of the interface from which this energy is returned is known as the ?rst Fresnel zone (Fig. 4.11) or, simply, the Fresnel zone. Around the ?rst Fresnel zone are a series of annular zones from which the overall re?ected energy tends to interfere destructively and cancel out. The width of the Fresnel zone represents an absolute limit on the horizontal resolution of a re?ection survey since re- ?ectors separated by a distance smaller than this cannot be individually distinguished.The width w of the Fresnel zone is related to the dominant wavelength l of the source and the re?ector depth z by

12 w = (2zl ) (for z >> l ) The size of the ?rst Fresnel zone increases as a function of re?ector depth. Also, as noted in Section 3.5, deeper- travelling re?ected energy tends to have a lower domi- nant frequency due to the effects of absorption. The lower dominant frequency is coupled with an increase in interval velocity, and both lead to an increase in the wavelength. For both these reasons the horizontal reso- lution, like the vertical resolution, reduces with increas- As a practical rule of thumb, the Fresnel zone width for the target horizons should be estimated, then the Seismic Re?ection Surveying 53

geophone spacing ?xed at no more than one-quarter of that width. In this case the horizontal resolution will be limited only by the physics of the seismic wave, not by the survey design.

4.4.2 Design of detector arrays Each detector in a conventional re?ection spread con- sists of an array (or group) of several geophones or hy- drophones arranged in a speci?c pattern and connected together in series or parallel to produce a single channel of output. The effective offset of an array is taken to be the distance from the shot to the centre of the array. Ar- rays of geophones provide a directional response and are used to enhance the near-vertically travelling re?ected pulses and to suppress several types of horizontally trav- elling coherent noise. Coherent noise is that which can be correlated from trace to trace as opposed to random noise (Fig. 4.12). To exemplify this, consider a Rayleigh surface wave (a vertically polarized wave travelling along the surface) and a vertically travelling compressional wave re?ected from a deep interface to pass simultane- ously through two geophones connected in series and spaced at half the wavelength of the Rayleigh wave. At any given instant, ground motions associated with the Rayleigh wave will be in opposite directions at the two geophones and the individual outputs of the geophones at any instant will therefore be equal and opposite and be cancelled by summing. However, ground motions asso- ciated with the re?ected compressional wave will be in phase at the two geophones and the summed outputs of the geophones will therefore be twice their individual The directional response of any linear array is gov- erned by the relationship between the apparent wave- length l of a wave in the direction of the array, the a The response is given by a response function R

sin nb R= sin b where b = pDx l a R is a periodic function that is fully de?ned in the inter- aa Typical array response curves are shown in Fig. 4.13.

1.0 Array response function R n=2 n=4 n=8 0 0 0.5 Detector spacing ?x Wavelength ?

Arrays comprising areal rather than linear patterns of geophones may be used to suppress horizontal noise The initial stage of a re?ection survey involves ?eld trials in the survey area to determine the most suitable combination of source, offset recording range, array geometry and detector spacing (the horizontal distance between the centres of adjacent geophone arrays, often referred to as the group interval) to produce good seis- Source trials involve tests of the effect of varying, for example, the shot depth and charge size of an explosive source, or the number, chamber sizes and trigger delay times of individual guns in an air gun array.The de- tector array geometry needs to be designed to suppress the prevalent coherent noise events (mostly source- generated). On land, the local noise is investigated by means of a noise test in which shots are ?red into a spread of closely-spaced detectors (noise spread) consisting of in- dividual geophones, or arrays of geophones clustered to- gether to eliminate their directional response. A series of shots is ?red with the noise spread being moved progres- sively out to large offset distances. For this reason such a test is sometimes called a walk-away spread. The purpose of the noise test is to determine the characteristics of the coherent noise, in particular, the velocity across the spread and dominant frequency of the air waves (shot noise travelling through the air), surface waves (ground roll), direct and shallow refracted arrivals, that together tend to conceal the low-amplitude re?ections. A typical noise section derived from such a test is shown in Fig.

1.0 Fig. 4.13 Response functions for different detector arrays. (After Al-Sadi 1980.)

4.12(a). This clearly reveals a number of coherent noise events that need to be suppressed to enhance the SNR of re?ected arrivals. Such noise sections provide the neces- sary information for the optimal design of detector phone arrays. Figure 4.12(b) shows a time section ob- tained with suitable array geometry designed to suppress the local noise events and reveals the presence of re?ec- tion events that were totally concealed in the noise It is apparent from the above account that the use of suitably designed arrays can markedly improve the SNR of re?ection events on ?eld seismic recordings. Further improvements in SNR and survey resolution are achiev- able by various types of data processing discussed later in the chapter. Unfortunately, the noise characteristics tend to vary along any seismic line, due to near-surface geological variations and cultural effects.With the tech- nical ability of modern instrumentation to record many hundreds of separate channels of data, there is an increas- ing tendency to use smaller arrays in the ?eld, record more separate channels of data, then have the ability to experiment with different array types by combining recorded traces during processing. This allows more so- phisticated noise cancellation, at the cost of some increase in processing time.

Each seismic trace then represents a unique sampling of some point on the re?ector. In common mid-point (CMP) pro?ling, which has become the standard method of two-dimensional multichannel seismic sur- veying, it is arranged that a set of traces recorded at dif- ferent offsets contains re?ections from a common depth The fold of the stacking refers to the number of traces in the CMP gather and may conventionally be 24, 30, 60 or, exceptionally, over 1000.The fold is alternatively ex- pressed as a percentage: single-fold = 100% coverage, six-fold = 600% coverage and so on.The fold of a CMP pro?le is determined by the quantity N/2n, where N is the number of geophone arrays along a spread and n is the number of geophone array spacings by which the spread is moved forward between shots (the move-up rate). Thus with a 96-channel spread (N = 96) and a move-up rate of 8 array spacings per shot interval (n = 8), the coverage would be 96/16 = 6-fold. A ?eld procedure for the routine collection of six-fold CMP coverage using a single-ended 12-channel spread con?guration progressively moved forward along a The theoretical improvement in SNR brought about by stacking n traces containing a mixture of coherent in- Stacking also attenuates long-path multiples. They have travelled in nearer-surface, lower velocity layers and have a signi?cantly different moveout from the primary re- ?ections. When the traces are stacked with the correct velocity function, the multiples are not in phase and do not sum. The stacked trace is the equivalent of a trace recorded with a vertical ray path, and is often referred to as a zero-offset trace.

4.4.4 Display of seismic re?ection data Pro?ling data from two-dimensional surveys are con- ventionally displayed as seismic sections in which the in- dividual stacked zero-offset traces are plotted side by side, in close proximity, with their time axes arranged verti- cally. Re?ection events may then be traced across the section by correlating pulses from trace to trace and in this way the distribution of subsurface re?ectors beneath the survey line may be mapped. However, whilst it is tempting to envisage seismic sections as straightforward images of geological cross-sections it must not be forgot- ten that the vertical dimension of the sections is time, not depth.

Seismic Re?ection Surveying 57 DS 21 Path 1 D 4S 2 Path 2 D 6S 3

Path 3 D 8S 4 Path 4 D 10S 5 Path 5 D 12S 6 Path 6 Common depth (reflection) point

Fig. 4.14 A ?eld procedure for obtaining six-fold CDP coverage with a single-ended 12-channel detector spread moved progressively along the survey line.

4.5 Time corrections applied to seismic traces Two main types of correction need to be applied to re- ?ection times on individual seismic traces in order that the resultant seismic sections give a true representation of geological structure. These are the static and dynamic corrections, so-called because the former is a ?xed time correction applied to an entire trace whereas the latter varies as a function of re?ection time.

(a) (b) D 3 Surface D 1 D 2Datum S

V w r e lya d h B s f at ere ao e we S 3 S 1 S 2 V 1

ttt 12 12 12 1 Seismic Re?ection Surveying 59

D SD SD 1 22 33 Datum

V 1 SD 11 SD 22 SD 33 t SD 11 SD 22 SD 33 ttt 23 23 23

Fig. 4.15 Static corrections. (a) Seismograms showing time differences between re?ection events on adjacent seismograms due to the different elevations of shots and detectors and the presence of a weathered layer. (b)The same seismograms after the application of elevation and weathering corrections, showing good alignment of the re?ection events. (After O’Brien, 1974.)

to the shot hole to measure the vertical time (VT ) or uphole time, from which the velocity of the surface layer above The complex variations in velocity and thickness within the weathered layer can never be precisely de- ?ned. The best estimate of the static correction derived It always contains errors, or residuals, which have the effect of diminishing the SNR of CMP stacks and reducing the coherence of re?ection events on time sections. These residuals can be investigated using so- This purely empirical approach assumes that the weath- ered layer and surface relief are the only cause of irregu- larities in the travel times of rays re?ected from a shallow interface. It then operates by searching through all the data traces for systematic residual effects associated with individual shot and detector locations and applying these as corrections to the individual traces before the CMP stack. Figure 4.16 shows the marked improvement in SNR and re?ection coherence achievable by the appli- cation of these automatically computed residual static In marine re?ection surveys the situation is much simpler since the shot and receivers are situated in a The static correction is commonly restricted to a con- version of travel times to mean sea-level datum, without removing the overall effect of the water layer. Travel times are increased by (d + d )v , where d and d are s hw s h the depths below mean sea-level of the source and hydrophone array and v is the seismic velocity of sea w water. The effect of marine tidal height is often signi?- cant, especially in coastal waters, and demands a time- variant static correction. Tidal height data are usually readily available and the only complexity to the correc- tion is their time-variant nature.

computer analysis of moveout in the groups of traces from a common mid-point (CMP gathers). Prior to this velocity analysis, static corrections must be applied to the individual traces to remove the effect of the low- velocity surface layer and to reduce travel times to a common height datum.The method is exempli?ed with reference to Fig. 4.17 which illustrates a set of statically corrected traces containing a re?ection event with a zero-offset travel time of t . Dynamic corrections are cal- 0 culated for a range of velocity values and the dynamically corrected traces are stacked. The stacking velocity V is st de?ned as that velocity value which produces the maxi- mum amplitude of the re?ection event in the stack of traces.This clearly represents the condition of successful removal of NMO. Since the stacking velocity is that which removes NMO, it is given by the equation x2 t 2 = t 2 + (cf. equation (4.4)) 0 V2 st Fig. 4.16 Major improvement to a seismic section resulting from residual (Courtesy Prakla Seismos GmbH.)

Seismic Re?ection Surveying 61 Offset Velocity DDDDDDx VVV 123456 123

Fig. 4.17 A set of re?ection events in a t 0 CMP gather is corrected for NMO using a range of velocity values.The stacking velocity is that which produces peak cross-power from the stacked events; that is, the velocity that most successfully removes the NMO. In the case illustrated, V represents the stacking velocity. (After 2 Taner & Koehler 1969.) Reflection time Fig. 4.18 The velocity spectrum is used to determine the stacking velocity as a function of re?ection time.The cross-power function (semblance) is calculated over a large number of narrow time windows down the seismic trace, and for a range of possible velocities for each time window.The velocity spectrum is typically displayed alongside the relevant CMP gather as shown. Peaks in the contoured semblance values correspond to appropriate velocities for that travel time, where a re?ection phase occurs in the CMP gather.

t 0 V 3 Cross-power V 2Peak power defines correct stacking velocity V 1

stacked wavelet. A velocity function de?ning the in- crease of velocity with depth for that CMP is derived by picking the location of the peaks on the velocity spec- Velocity functions are derived at regular intervals along a CMP pro?le to provide stacking velocity values for use in the dynamic correction of each individual trace.

Fig. 4.19 Filter panels showing the frequency content of a panel of re?ection records by passing them through a series of narrow-band frequencies.This plot allows the geophysicist to assess the frequency band that maximises the signal-to-noise ratio. Note that this may vary down the traces due to frequency-dependent absorption. (From Hatton et al. 1986, p. 88)

main types of waveform manipulation are frequency ?l- tering and inverse ?ltering (deconvolution). Frequency ?ltering can improve the SNR but potentially damages the vertical resolution, while deconvolution improves the resolution, but at the expense of a decrease in the SNR. As with many aspects of seismic processing, com- promises must be struck in each process to produce the optimum overall result.

4.8.1 Frequency ?ltering Any coherent or incoherent noise event whose domi- nant frequency is different from that of re?ected arrivals may be suppressed by frequency ?ltering (see Chapter 2). Thus, for example, ground roll in land surveys and several types of ship-generated noise in marine seismic surveying can often be signi?cantly attenuated by low- cut ?ltering. Similarly, wind noise may be reduced by high-cut ?ltering. Frequency ?ltering may be carried out at several stages in the processing sequence. Nor- mally, shot records would be ?ltered at a very early stage in the processing to remove obvious noise. Later applica- tions of ?lters are used to remove artefacts produced by other processing stages. The ?nal application of ?lters is to produce the sections to be used by the seismic inter- preters, and here the choice of ?lters is made to produce the optimum visual display.

Since the dominant frequency of re?ected arrivals decreases with increasing length of travel path, due to the selective absorption of the higher frequencies, the char- acteristics of frequency ?lters are normally varied as a function of re?ection time. For example, the ?rst second of a 3 s seismic trace might typically be band-pass ?ltered between limits of 15 and 75 Hz, whereas the frequency limits for the third second might be 10 and 45 Hz. The choice of frequency bands is made by inspection of ?lter panels (Fig. 4.19). As the frequency characteristics of re?ected arrivals are also in?uenced by the prevailing geology, the appropriate time-variant frequency ?lter- ing may also vary as a function of distance along a seismic pro?le. The ?ltering may be carried out by computer in the time domain or the frequency domain (see Chapter 2).

processing, each designed to remove some speci?c ad- verse effect of ?ltering in the ground along the transmis- Deconvolution is the analytical process of removing the effect of some previous ?ltering operation (convolu- tion). Inverse ?lters are designed to deconvolve seismic traces by removing the adverse ?ltering effects associated with the propagation of seismic pulses through a layered ground or through a recording system. In general, such effects lengthen the seismic pulse; for example, by the generation of multiple wave trains and by progressive ab- sorption of the higher frequencies. Mutual interference of extended re?ection wave trains from individual inter- faces seriously degrades seismic records since onsets of re?ections from deeper interfaces are totally or partially concealed by the wave trains of re?ections from shal- Examples of inverse ?ltering to remove particular ?ltering effects include: · dereverberation to remove ringing associated with · deghosting to remove the short-path multiple asso- ciated with energy travelling upwards from the source and re?ected back from the base of the weathered layer or the surface; and · whitening to equalize the amplitude of all frequency components within the recorded frequency band (see All these deconvolution operations have the effect of shortening the pulse length on processed seismic sec- Consider a composite waveform w resulting from an k initial spike source extended by the presence of short- path multiples near source such as, especially, water layer reverberations. The resultant seismic trace x will be k given by the convolution of the re?ectivity function r k with the composite input waveform w as shown k schematically in Fig. 4.6 (neglecting the effects of attenuation and absorption)

xk = rk * wk (plus noise) Re?ected waveforms from closely-spaced re?ectors will overlap in time on the seismic trace and, hence, will in- terfere. Deeper re?ections may thus be concealed by the reverberation wave train associated with re?ections from shallower interfaces, so that only by the elimination of Note that short-path multiples have effectively the same normal moveout as the related primary re?ection and are Seismic Re?ection Surveying 63

therefore not suppressed by CDP stacking, and they have similar frequency content to the primary re?ection so Deconvolution has the general aim, not fully realiz- able, of compressing every occurrence of a composite waveform w on a seismic trace into a spike output, in k order to reproduce the re?ectivity function r that would k fully de?ne the subsurface layering.This is equivalent to the elimination of the multiple wave train.The required deconvolution operator is an inverse ?lter i which, k when convolved with the composite waveform w , k yields a spike function d k

ik * w = d kk Convolution of the same operator with the entire seismic trace yields the re?ectivity function

2. that the composite waveform w for an impulsive k source is minimum delay (i.e. that its contained energy is concentrated at the front end of the pulse; see From assumption (1) it follows that the autocorrela- tion function of the seismic trace represents the autocor- relation function of the composite waveform w . From k assumption (2) it follows that the autocorrelation func- tion can be used to de?ne the shape of the waveform, the

* Input waveform Filter operator

= Filtered output Desired output (a) ? (?) xx necessary phase information coming from the Such an approach allows prediction of the shape of the composite waveform for use in Wiener ?ltering. A par- ticular case of Wiener ?ltering in seismic deconvolution This is the basis of spiking deconvolution, also known as whitening deconvolution because a spike has the amplitude spectrum of white noise (i.e. all frequency components A wide variety of deconvolution operators can be de- signed for inverse ?ltering of real seismic data, facilitat- ing the suppression of multiples (dereverberation and deghosting) and the compression of re?ected pulses.The presence of short-period reverberation in a seismogram is revealed by an autocorrelation function with a series of decaying waveforms (Fig. 4.21(a)). Long-period rever- berations appear in the autocorrelation function as a series of separate side lobes (Fig. 4.21(b)), the lobes occurring at lag values for which the primary re?ection aligns with a multiple re?ection.Thus the spacing of the side lobes represents the periodicity of the reverberation pattern.The ?rst multiple is phase-reversed with respect to the primary re?ection, due to re?ection at the ground surface or the base of the weathered layer. Thus the ?rst side lobe has a negative peak resulting from cross- correlation of the out-of-phase signals.The second mul- tiple undergoes a further phase reversal so that it is in phase with the primary re?ection and therefore gives rise 4.21(b)). Autocorrelation functions such as those shown in Fig. 4.21 form the basis of predictive deconvolution operators for removing reverberation events from seismograms.

? (b) ? (?) xx

? Fig. 4.21 Autocorrelation functions of (a) A gradually decaying function (b) A function with separate side lobes indicative of long-period reverberation.

Fig. 4.22 Removal of reverberations by predictive deconvolution. (a) Seismic record dominated by strong reverberations. (b) Same section after spiking deconvolution. (Courtesy Prakla Seismos GmbH.)

Practically achievable inverse ?lters are always approximations to the ideal ?lter that would produce a re?ectivity function from a seismic trace: ?rstly, the ideal ?lter operator would have to be in?nitely long; second- ly, predictive deconvolution makes assumptions about the statistical nature of the seismic time series that are only approximately true. Nevertheless, dramatic im- provements to seismic sections, in the way of multiple suppression and associated enhancement of vertical resolution, are routinely achieved by predictive decon- volution. An example of the effectiveness of predictive deconvolution in improving the quality of a seismic section is shown in Fig. 4.22. Deconvolution may be carried out on individual seismic traces before stacking (deconvolution before stacking: DBS) or on CMP stacked traces (deconvolution after stacking: DAS), and is com- monly employed at both these stages of data processing.

Seismic Re?ection Surveying 65 4.8.3 Velocity ?ltering The use of velocity ?ltering (also known as fan ?ltering or pie slice ?ltering) is to remove coherent noise events from seismic records on the basis of the particular angles at which the events dip (March & Bailey 1983). The angle of dip of an event is determined from the apparent velocity with which it propagates across a spread of A seismic pulse travelling with velocity v at an angle a to the vertical will propagate across the spread with an apparent velocity v = v/sin a (Fig. 4.23). Along the a spread direction, each individual sinusoidal component of the pulse will have an apparent wavenumber k related a to its individual frequency f, where