Process for characterising the evolution of an oil or gas reservoir over time

Abstract

Disclosed is a method for characterizing the evolution of a reservoir by determining a seismic wavelet which links observed seismic data to a sequence of reflectivities. The method comprises obtaining seismic data (200) representing seismic changes which have occurred between a first time and a second time, said seismic data comprising a plurality of seismic traces; and performing an optimization operation simultaneously (230) on the seismic traces so as to optimize for said seismic wavelet. The optimization operation may be performed without using known reflectivity data as an input.

Claims

1. A method for characterising evolution of a reservoir by determining a seismic wavelet which links observed seismic data to a sequence of reflectivities, said method comprising: obtaining seismic data representing seismic changes which have occurred between a first time and a second time defining a production period over which hydrocarbons have been extracted, said seismic data comprising a plurality of seismic traces; and performing an optimisation operation simultaneously on said plurality of seismic traces so as to optimise for said seismic wavelet, said optimization operation comprising simultaneously optimising for said seismic wavelet and reflectivity change data occurring between said first time and said second time; extracting the time-lapsed changes in the seismic traces collected over the production period; integrating the time-lapsed changes in the seismic traces with production data; using the integrated time-lapsed changes in the seismic traces and the production data to manage at least one of extraction of oil and gas from the reservoir and injection of other fluids into the reservoir; and wherein said optimisation operation is performed without using known reflectivity data as an input.

2. The method as claimed in claim 1 wherein said optimisation comprises minimising a cost function comprising a measurement of a difference between the seismic data and a convolution of the seismic wavelet and reflectivity change data.

3. The method as claimed in claim 1 wherein each reflectivity change described by said reflectivity change data is defined by two parameters describing positions of two reflectors defining the change in reflectivity and a single amplitude parameter representing a magnitude of a measured amplitude change at the positions of said two reflectors, the amplitude change being equal in magnitude and of opposite sign at these positions.

4. The method as claimed in claim 1 comprising the steps of: providing a base survey of the reservoir with a set of seismic traces at said first time; providing a monitor survey of the reservoir, taken at said second time, with a set of seismic traces associated to the same positions as in the base survey; and determining said seismic data from said base survey and monitor survey.

5. The method as claimed in claim 4 wherein said step of determining said seismic data comprises aligning said monitor and base surveys; and subtracting said base survey from said monitor survey.

6. The method as claimed in claim 4 wherein said base and monitor surveys have been performed at a location remote from a well, or adjacent to a horizontal or sub-vertical well.

7. The method as claimed in claim 4, wherein said base and monitor surveys have been performed at a location of interest for the wavelet being determined.

8. The method as claimed in claim 1 comprising a scaling step so as to determine an absolute value for amplitude of the seismic wavelet.

9. The method as claimed in claim 8 wherein said scaling step comprises defining a scaling factor such that ?S=??*1/??R, where ? is the scaling factor, ?S is the seismic data, w is the seismic wavelet and ?R is a reflectivity change data.

10. The method as claimed in claim 8 wherein said scaling step comprises convolving an unscaled wavelet with actual reflectivity measurements taken from another location of the reservoir.

11. The method as claimed in claim 8 wherein said scaling step comprises using travel-time information to estimate the actual reflectivity measurements.

12. The method as claimed in claim 11 wherein said actual reflectivity measurements are determined from a product of measured relative velocity changes and a factor based on an estimated ratio between relative velocity changes and relative density changes.

13. The method as claimed in claim 8 wherein said scaling step comprises directly inferring maximum reflectivity changes from production history data.

14. The method as claimed in claim 1 further comprising the step of using results of said method to aid hydrocarbon recovery from a reservoir.

15. A computer program comprising computer readable instructions which, when run on suitable computer apparatus, cause the computer apparatus to perform the method of claim 1.

16. A computer program carrier comprising the computer program of claim 15.

17. An apparatus specifically adapted to carry out the steps of the method as claimed claim 1.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) Embodiments of the invention will now be described, by way of example only, by reference to the accompanying drawings, in which:

(2) FIGS. 1a and 1b are schematic illustrations of a (1a) base survey and (1b) monitor survey being performed; and

(3) FIG. 2 is a flowchart describing a method according to an embodiment of the invention.

DETAILED DESCRIPTION OF THE EMBODIMENTS

(4) Referring initially to FIGS. 1a and 1b there is illustrated a reservoir, generally indicated by reference numeral 10, containing hydrocarbons 12 in the sub-surface 14. A survey vessel 16, upon which is located a sonar transmitter 18, being an acoustic source, and an array of receivers 20, performs a survey by travelling over the reservoir 10. The first or initial survey, FIG. 1a, may be referred to as a base survey and is typically performed in the exploration phase before production begins.

(5) The base survey of the reservoir 10 provides a set of seismic traces at a first time T. For a given trace, the base survey provides amplitudes that are a function of time. With digital recording and processing the trace is sampled at a set of values; typical trace lengths correspond to around 1000 samples. The trace is then handled as a set of values.

(6) One or more wells 22 may be drilled in order to extract the hydrocarbons 12. As the reservoir 10 is produced, hydrocarbons will be substituted by other fluids and the fluid pressure will change. Additionally, enhanced oil recovery techniques may be applied wherein a fluid is injected into the reservoir at one or more locations giving changes in fluid pressure and saturation. Changes within the reservoir may also change the stress and strain state of the surrounding rocks. Thus when a further survey is carried out, FIG. 2b, these changes will be seen due to a consequential change in the velocity field. These changes to velocity will produce time shifts in the seismic expression of underlying reflectors and associated changes in reflectivity, causing a change in the local wavefield.

(7) Thus reservoir monitoring performs a monitor survey of the reservoir 10, taken at a second time T+? T, with a set of seismic traces. In the simplest assumption, ? T is a positive quantity, and the monitor survey is taken at a time later than the base survey; however, the order in which the surveys are taken is irrelevant to the operation of the process of the invention and, in principle, the time lapse ?T could as well be negativewhich amounts to comparing the earlier survey to the later one. As for the base survey, a sampled trace in the monitor survey is represented as a set of values.

(8) Ideally, the traces in the monitor survey are associated to the same positions as in the base survey. This is carried out by using, inasmuch as possible, the same equipment, acquisition geometry and processes for running the base and monitor surveys. Techniques such as interpolation may be used where traces in the monitor survey and in the base survey do not fulfil this condition.

(9) In order to link the seismic data (traces) on which interpretations are based and the geology (reflection coefficients) being interpreted, the seismic wavelet (an impulse response which travels through the surface) is defined, based on the known convolutional model of a seismic trace, where the seismic trace is a convolution of the seismic wavelet with the subsurface reflectivity (plus noise). The seismic wavelet is the waveform which would be recorded by a seismic system for a reflection from a single plane reflecting boundary in the subsurface. This wavelet can be thought of as acting like a filter through which the geology is viewed when interpreting the image provided by seismic data. However, this wavelet needs to be estimated accurately.

(10) The problem of wavelet estimation has been a key problem in geophysics for some time. Present methods of wavelet estimation comprise deconvolving a seismic trace with a sequence of reflection coefficients from a seismogram. A common way of estimating wavelets is to estimate the wavelet that would best match the reflectivities measured at the well using a convolution process on the 3D seismic data. This technique has many drawbacks which include:

(11) 1. Wavelets are not stable neither laterally nor vertically and if the wells are distant to the place of interest, the wavelets may not be well adapted;

(12) 2. Wavelets estimated in this way are adapted to 3D and may not be adequate for 4D processing;

(13) 3. This technique requires a large time window and so the estimated wavelet may be an average wavelet;

(14) 4. With better drilling technology the majority of wells (producer and injectors) are now drilled sub-vertically or in some cases horizontally through the reservoir. This makes conventional wavelet estimation difficult and sometimes impossible.

(15) Other techniques attempt to estimate the wavelets directly from the seismic data: the spectrum can be estimated from the square root of the autorcorrelation (with the assumption that reflectivity has a white spectrum) and phase through Kurtotsis analysis (with assumption of uniform distribution); however none of these approaches are reliable and they are hardly used.

(16) The problem of an inaccurately defined wavelet is more acute when dealing with 4D seismic data since only tiny anomalies are being looked for, and an inaccurate wavelet can be a huge source of error. Despite this, wavelets estimated as above from 3D seismic data at the wellhead are typically used when interpreting 4D seismic data.

(17) Therefore, instead of using the reflectivities measured at the well, it is proposed to use the reflectivities provided by the 4D signal. Advantages of this approach include:

(18) 1. 4D reflectivities are limited to only a small number of layers where dynamic changes have occurred;

(19) 2. the reflectivities can be measured away from the wells where the 4D signal actually exists.

(20) FIG. 2 is a flowchart describing a method using the 4D signal to estimate the wavelet. At step 200, the Survey data (Base and Monitor) is obtained. This data may have been gathered as described in relation to FIGS. 1a and 1b.

(21) At step 210, the monitor survey data is aligned to the base survey data. While the base and monitor surveys each show a large number of reflections (for every boundary), once aligned, most of the reflections overlap. Only in regions where changes have occurred over the time between the surveys will there be non-overlapping reflections. Consequently the change in the seismic data ?S between base and monitor is sparse, with few dynamic reflections, and can be defined as:
?S=?+R.sub.m??*R.sub.B=?*?R
where ? is the wavelet, R.sub.B are the base reflectivities, R.sub.m are the monitor reflectivities and ?R is the change in reflectivities between base and monitor.

(22) Considering a layer where the saturation or pressure (or both) have changed, there will be a change of impedance ?Ip at the top of this layer and an opposite change of impedance ??Ip at the bottom of this layer. Therefore the change in reflectivities ?R for any dynamic layer can be defined using only three parameters relating to the positions of the boundaries defining the layer and the magnitude of the reflectivity signal at these boundaries. In one embodiment, each layer is defined by the position of the layer top t (or bottom or any other position identifier), the thickness of the layer ?t and the change in the reflectivity signal amplitude ?A for the layer. The 4D seismic data observed will be the result of the convolution of the wavelet by a dipole of opposite sign at position t, t+?t.

(23) At step 220 initial values are given to the unknowns in ??R (amplitudes, initial reflectivities and the wavelet. A general initial guess of the number of layers is the minimum input. More complex initial guesses can optionally be made by picking horizons from the 4D or 3D data.

(24) At step 230, all the seismic traces (or a subset thereof) are optimised simultaneously. The expected ranges of variation of the initial values can also be specified. This optimisation may be performed by minimising a cost function such as:
cost=?S??*?R
although it should be appreciated that any norm or difference measurement of ?S and ?*?R can be used to calculate the cost.

(25) While it is true that both the wavelet ? and ?R are both unknown, their values relative to each another can be determined. The sparse distribution of the reflection coefficients in the 4D seismic data is the key of this technique. There are theorems in the super resolution literature showing that under certain conditions of sparsity the inversion process is exact. Of course the seismic data are not noise free and noise will perturb the inversion response. Fortunately as the single inversion uses data from several or many seismic traces, the wavelet can be constrained and a unique wavelet can be solved for. If instead of a single layer there are several layers where changes occur, the situation is the same since the convolution process is linear.

(26) Therefore, for a signal composed of N traces, with M layers (composed of 2 reflectors) and a wavelet of length L, there are in total 3?M?N+L parameters to optimise (where the 3 in this total results from the three parameters which define ?R). This total is much smaller than N?S, where S is the number of samples. Therefore the problem is overconditioned.

(27) With the wavelet estimated in this way, the absolute amplitude of the wavelet co or of the change in reflectivity ?R is unknown. In an embodiment, the relative impedance inversion values may be used as they are. In another embodiment, the wavelet may be scaled (step 240) to determine its absolute amplitude. A scaling factor ? (positive or negative) may be defined, such that:
?S=??*1/??R

(28) The scaling step determines the constant ? that determines correctly the true scaled wavelet:
?.sub.true=??,
and the true scaled reflectivity:
?R.sub.true=1/??R,

(29) There are a number of different options for performing scaling step 240, which include: 1. Convolving the unscaled wavelet with some true reflectivity measured (for example) at a well location: S=?*?R.sub.true (local to the estimation or not). The computed synthetic trace (S) can be compared with true seismic (S.sub.true=?.sub.true*?R.sub.true) at the well location using the ratio:
S.sub.true/S=?.sub.true/?=?. 2. Using travel-time information. By using the fact that 4D data has a time shift that is given by the integral of relative velocity changes it is possible to use the base and non-aligned monitor traces to determine the magnitude of ?Vp/Vp. Thus using the commonly known expression: ?R??Ip/Ip??Vp/Vp+??/????VpVp and assuming a value for ?, the true reflectivity can be approximated along with the scaling factor ? (where ?Vp/Vp is the change in p-wave velocity value, ?Ip/Ip the change in impedance, and ??/? is the change in density). The factor ? may represent an estimated ratio between ?Vp/Vp, and ??/? (such an estimation can be made based upon knowledge of the subsurface composition, as is understood by the skilled person) and may in a specific embodiment equate to 1+?Vp/Vp/??/?. 3. It is also possible to use prior information based on production information, reservoir simulation and rock physics modelling to directly infer the maximum reflectivity changes (max(?R.sub.true)) expected, and to scale the data accordingly by max(?R)/max(?R.sub.true)=?.

(30) One or more steps of the methods and concepts described herein may be embodied in the form of computer readable instructions for running on suitable computer apparatus, or in the form of a computer system comprising at least a storage means for storing program instructions embodying the concepts described herein and a processing unit for performing the instructions. As is conventional, the storage means may comprise a computer memory (of any sort), and/or disk drive, optical drive or similar. Such a computer system may also comprise a display unit and one or more input/output devices.

(31) The concepts described herein find utility in all aspects (real time or otherwise) of surveillance, monitoring, optimisation and prediction of hydrocarbon reservoir and well systems, and may aid in, and form part of, methods for extracting hydrocarbons from such hydrocarbon reservoir and well systems.

(32) It should be appreciated that the above description is for illustration only and other embodiments and variations may be envisaged without departing from the spirit and scope of the invention.