NO347341B1 - Method for time and position corrections for Ocean Bottom Seismic data - Google Patents

Description

METHOD FOR TIME AND POSITION CORRECTIONS FOR OCEAN BOTTOM SEISMIC DATA

Terminology

The word ‘node’ refers to submersible recording device, which may also be described as a data logger, which registers signal input from a hydrophone and/or motion sensor(s) and/or accelerometer(s). A node records seismic pressure wave and/or shear waves. A node is commonly referred to in the industry as an Ocean Bottom Node (OBN) or seafloor node.

The word ‘span’ refers to a limited period of time, usually a period of a few hours.

Scope

The invention aims to improve and make more economical the methods for acquiring geophysical data in water covered areas. Such geophysical data is used for academic and commercial investigations into the earth’s geological structures. The experimental method used is often referred to as ‘The Seismic Method’ and is well known and used in the oil and gas exploration industry and also for academic studies.

The Seismic Method is a generic experimental procedure which aims to record the earth’s response to excitation by a seismic source. The earth’s response consists of a complicated series of acoustic waves and shear waves which contain information about the earth’s sub-surface. The seismic excitation may be generated by an airgun or array of arrays, but it may also be generated by a vibrator, or a source in a borehole, or a natural source such as an earthquake or any other method. Generally the seismic excitation will produce both acoustic and shear waves, which may collectively be termed seismic waves or seismic energy. Water has a shear modulus of zero and therefore seismic excitation in the water column will only produce acoustic waves. From the point of excitation acoustic and shear waves propagate in patterns dictated by physical properties of the materials that the waves travel through. As waves travel into the earth, they are reflected, refracted, transmitted, attenuated and mode-converted according to the laws of physics and the physical properties of the rock types encountered. Some of the acoustic and shear wave energy after an excitation event can be recorded at the seafloor by sensitive recording devices. Analysis of the recorded data provides information about the physical properties of the materials through which the seismic energy propagated. The acoustic waves may be recorded by pressure sensors, but may also be recorded by particle motion sensors or accelerometers. Shear waves may be recorded by particle motion sensors or accelerometers. As shear waves cannot propagate through water, it is important that the sensors are located at the sea-floor and that they are coupled such that shear wave energy will propagate into the sensor unit.

Description of the invention

In the seismic method, the seismic signals are recorded for subsequent analysis. For these recordings to be of use, the timing of the recordings relative to source event timing as well as the location of the recordings is required. This invention concerns a method for correcting for timing errors and establishing the locations for such sensor recordings. Both time relative to source event timing and the position need to be established repeatedly as the local time in each node may deviate from the source event time and the node position may change under external influence such as currents. Source event timing is commonly synchronised with an external reference time such as GPS time or World Time.

The local time within a node is controlled by a local oscillator chip or ‘clock’. During a survey the nodes may be deployed for extended periods of time ranging from several days to 120 days or more. At the time of deployment, the local clock time relative to source event controller time may be established, however, the local clock time may run faster or slower than the source event controller clock. The accuracy of the seismic method is compromised if the local node clock time is more than 1ms different from the source event reference time. The 1ms limit corresponds to the Nyquist sampling rate of a 500Hz signal, which is considered here as the upper physical limit of useful seismic signals. This accuracy may in practice be relaxed to 2ms or 4ms or 8ms depending on the geological and geophysical survey requirements. For the purpose of this patent, the 1ms requirement is used to illustrate positioning and clock accuracy requirements.

Taking a 30 day deployment as an example; at the end of the 30 day period the 1ms accuracy is still required. This corresponds to an accuracy requirement in the order of 0.4ppb or 4x10<-10 >[1ms cumulative error per 30days]. This accuracy can currently only be achieved by the use of specialised clocks, which are expensive and have relatively high electric power consumption. By using the data driven method as described in this invention, the recorder clock time variations are calculated and corrected for in post processing, enabling the local node clock accuracy to be relaxed to as much as 10ppm or 1x10<-5>.

A possible implementation of the method is for an Ocean Bottom Seismic Recording Node system. A possible node design contains a sensor for recording pressure, a digital clock, electronics to record the signal data, a memory storage system for storing data, a battery or rechargeable power source or other power source. Optionally the node are fitted with a wireless communications system and a wireless charging system.

Optionally the node is completely sealed without external connectors to ensure that no water penetrates the node. Optionally the node may contain sensors to record motion and/or acceleration and/or an inclinometer and/or a compass and/or a temperature sensor. Sensor data, and optionally data from the pressure sensor and/or the inclinometer and/or the compass and/or the temperature sensor will be recorded and stored while the device is in operation. Optionally the node functionality may be extended to recording of electromagnetic waves.

The method of this patent is equally applicable to ocean bottom cable systems.

Discussion of existing Prior Art

Designs and methods known in this field rely on expensive chip scale atomic clocks, known as CSAC or similar to ensure the local clock time does not deviate more allowed while the device is in operation. The node position may be obtained using an external device such as an acoustic location device (USBL, acoustic beacon system or similar), or may be derived from the deployment vehicle location, or may be derived from the seismic data if the clock is sufficiently accurate.

An attempt has been made to solve the timing and positioning simultaneously by post-processing the seismic data, but in practice it is hard to implement in a stable fashion because some of the variables have very similar error expression. Specifically, Prior Art as described in US Patent US 8,995,222 B2 titled “System and method for accurate determination of ocean bottom seismometer positioning and timing” describes an instant method utilizing linearized inversion in conjunction with a conventionally accurate clock to provide both time and positioning for each OBS unit with high accuracy as compared with the prior art approach. The suggested simultaneous determination of the node clock timing and node positioning problematic due to relationships of the variables. It is very difficult to ensure convergence to a plausible solution as the clock error, the velocities and the node depth are not well resolved as independent variables due to their similar expression in terms of travel times. The use of a single step in this Prior Art is clear from the use of the word simultaneous in claim 1f and from the mathematical formulation as specified in claim 4. As the position and time update are updated in a single solver step the earlier method is fundamentally different to the method proposed in the current disclosure.

The present invention solves this serious shortcoming by solving for time error and position in autonomous steps, aided by a priori knowledge of the bathymetry. This method has been shown to work in a stable manner on real world data sets.

The main novelty lies in the time regression, which is an autonomous part of the method. This step is based on the fact that within a period of a few hours, which shall be referred to as a ‘span’, the time and temperature variant clock error, hence referred to as ‘static error’, is relatively small. Hence when data is analysed within a span, the clock variation within this period may be largely ignored and a reasonable estimate of clock static error correction within this period may be established based on approximate distance from source events to the node and observed arrival times of seismic events such as direct waves, refractions and reflections. By extension the static error correction may then be established for other spans and a relationship may be established between the static error and the known time of the source event, usually GPS time. This relationship may be made more accurate by taking into account a priori knowledge that crystal oscillators used in the clocks have a known logarithmic aging profile, and the clock frequency varies with temperature. The established relationship may then be used to correct the observed seismic arrival times for node clock deviations.

US 8,995,222 B2 includes the step of determining a length of time since a clock associated with said selected receiver was last synchronized [claim 1 e]. The method considered in the present invention does not require any such determination of time since last synchronisation as it relies on correlating source events with observed seismic events registered by the node to obtain a first estimate of clock error as disclosed in WO2019224354 A1.

Within a span, the number of observations is necessarily limited. This has a negative impact on the accuracy of the calculated position. However by using an approach similar to height aiding as for example described in patent US8736487B2 for improving GPS satellite positioning solution, bathymetry data may be used to improve the position solution within a span. This is also implicitly referred to in WO2019224354 A1.

The use of various common mathematical techniques, such as correlation, regression, least squares positioning algorithm need to be combined in a specific way to give a sable result.

In summary, the present method combined the insight that the problem may be simplified by considering spans of time, with elements from GPS positioning technology, seismic recording and processing technology and mathematics to allow the use of economical, low power crystal clocks in nodes.

Nature and use of the invention

The method described in the present patent may be used to determine the clock corrections as well as the location of a node or a plurality of nodes. These node(s) may or not be connected by any means and may or may not be in contact by whatever means with each other.

In the preferred embodiment, multiple nodes are deployed on the seafloor and left there for extended periods of time in order to record seismic data. Multiple traverses are made across or near to the deployed nodes with one or more seismic sources which fire regularly creating groupings of data called sail line passes. Such traverses, often referred to as sail line passes, are usually several hours in duration and are useful spans for the purpose of grouping recorded data.

Upon retrieval, the data from the nodes are transferred to a computer system, whereupon the data are analysed and processed to determine the relative time within each node compared to source event controller time and the position of each node while it was deployed. This process consists of several steps: getting a starting point [deriving approximate node position and span static error], static function update [F2.2: time regression], lateral position update [F2.3 least squares estimation of node x and y coordinates], velocity model update [F2.4: distance regression], full position update [F2.5: least squares estimation of node x, y and z coordinates].

Table 3 describes a workflow which will work without any a priori knowledge of node location or node timing. As no initial position or synchronisation is needed, the workflow may be thought of as solving a ‘double blind’ problem.

List of figures

Figure 1 STAGE F3.5: without timing correction, the first break arrival times make no sense.

Figure 2 STAGE F3.10 without timing corrections, and with incorrect node position and incorrect velocity model, a near linear trend is observed.

Figure 3 STAGE F3.15 with time corrected, first breaks and calculated offsets have a linear trend for direct arrivals.

Figure 4 STAGE F3.15 with time corrected, there are still significant errors in c-o.

Figure 5 STAGE F3.16 with time and position corrections, the first breaks and offset linearization is optimised.

Figure 6 STAGE F3.16 with time and position corrections, the c-o is minimised.

Detailed explanation

The recordings made by a node on the sea floor may be used to relate timing within a node to source sequence timing, and to derive the position of said node. In a possible embodiment, data from spans are first correlated with source event timing to establish causal relationships between source events and observed recordings for each span. Specifically, the arrival times of seismic events recorded in a node may be correlated to source events of known location and timing. Approximate positions may be obtained from the deployment method (node drop positions, dead reckoning, USBL, deployment vehicle) or they may be deduced from signal strength analysts. A starting velocity model may be a single layer constant velocity model represented the water layer [TABLE 1 STEP 0]. With these initial values the travel time from each source to a node may be calculated for the various expected seismic events: direct wave, refractions, reflections and converted waves. Possible corrections to these calculated arrival times may be based on time and space varying water currents, source array dimensions, individual source element timing, and calculated tidal height [TABLE 1 STEP 2]. The calculated arrival times may be compared to actual arrival times observed in the node recording. The difference between calculated arrival times and observed arrival times is referred to as c-o (calculated – observed).

As local node clock time is not synchronous with the timing of the source events, the calculated arrival time will not match the observed arrival time. We must remove the effect of the local clock error. The crystal clocks types considered are commonly referred to as OCXO (Oven Controlled Crystal Oscillators) or TCXO (Temperature Compensated Crystal Oscillator) or equivalent. These clocks may have a frequency offset relative to design specification, causing a cycle to be shorter or longer than nominally specified. The cumulative effect of this frequency offset expresses itself as clock drift. The oscillating frequency is known to change with the length of time of use of the clock and expresses itself as a slow, logarithmic change of oscillator frequency with oscillator age. The temperature of the oscillator also influences the frequency of oscillation. The net effect is that the local clock error is almost linear due to the frequency offset, but there is also a logarithmic aging and a temperature dependence that should be incorporated in the correction. It is not unusual for clocks to have an internal compensation method to correct for oscillator frequency offset and/or temperature, but there is usually still a residual error.

In a possible embodiment, the c-o from several spans are first used to derive an approximate relationship between clock age, clock temperature on the one side and c-o on the other side in a step referred to as ‘time regression’ [TABLE 1 STEP 3].

In this possible embodiment data from several spans are processed to establish an approximate lateral node position through lateration. This lateration may consist of a least squares iterative approach, a direct triangulation method using source event positions and associated estimated travel distances or another method. This step is called ‘XY lateration’ [TABLE 1 STEP 4]

This may be followed by another instance of time regression, after which the data may be used to update the velocity model, which in the simplest embodiment would entail updating the water velocity by regressing the calculated offset on to the observed first arrival times in a step called ‘distance regression’ [TABLE 1 STEP 5]. A relevant aspect observation is that the water velocity is not static or uniform and may vary temporally and spatially.

After updating the velocity model a further refinement of the node position in three dimensions is normally possible in a step referred to as ‘XYZ solver’. [TABLE 1 STEP 6]

At any stage one may elect to reject observations that are spurious in a step called ‘outlier rejection’.

The output of each stage may be verified visually by plotting the observed arrival time minus the calculated arrival time for each observation. Such plots are well known in the geophysical industry and are referred to as LMO corrected gathers or NMO corrected gathers. The considered seismic events should line-up as a flat and centred on zero time if the velocity model, node timing correction and node position are correct.

A possible measure for the quality of the result is the absolute sum of all c-o values.

Having achieved sufficient accurate estimate of the clock error, it is possible to update node X, Y, Z -coordinates with a stable least squares solver – velocity model update loop.

Tables

Table 1 GENERIC STEPS

Table 2 MODEL WORK FLOW

Table 3 EXAMPLE ACTUAL WORKFLOW “XY LATERATION”

BASED ON TAYLOR EXPANSION OF THE DISTANCE EQUATION

Or written in vector and matrix format:

with least squares solution

is the driving vector containing ‘observations - provisional values’ is the output vector containing estimated residuals ‘update - provisional’ A is the design matrix

Table 4 XY LATERATION

“XYZ SOLVER”

BASED ON TAYLOR EXPANSION OF THE DISTANCE EQUATION – IGNORING Z COMPONENT

Or written in vector and matrix format:

with least squares solution

is the driving vector containing ‘observations - provisional values’

is the output vector containing estimated residuals ‘update - provisional’ A is the design matrix

Table 5 XYZ SOLVER

Relevant prior art:

EP0958511A1 titled “Acoustic positioning of seismic ocean bottom cable”

US Patent 6,005,828 titled “Acoustic Positioning of Seismic Ocean Bottom Cable”

US Patent 7,660,189 B2 titled “Apparatus, systems and methods for determining position of marine seismic acoustic receivers”

US 8,995,222 B2 titled “System and method for accurate determination of ocean bottom seismometer positioning and timing”

US 8,736,487 B2 titled “Method and apparatus of using height aiding from a contour table for GNSS positioning”

US Patent 2016/0146958 A1 titled “Method and computer system for determining seismic node position” US Patent US 2017/0090052 A1 titled “Determining node depth and water column transit velocity” WO2017105885A1 titled “System and method for correction of seismic receiver clock drift”

US 5,757,722 “Method for verifying the location of an array of detectors”

US 5,696,733 “Method for verifying the location of an array of sensors”

Textbook (1997) “Geodetic Applications of GPS,” published by the Swedish Land Survey.

WO2019224354 A1. “METHOD FOR ESTABLISHING POSITION AND TIMING OF SEISMIC

RECORDING DEVICES THAT ARE DEPLOYED ON A SEAELOOR”

Claims (9)

Claims

1) A method for establishing node timing and node position by continuously recording underwater sound and/or vibrations with a node with known approximate position, while recording source events within the recording time of the node, at known location and time, the method comprising a post recording correction characterised by the steps:

i. Calculate travel times from the approximate node position, the known source event positions, known source event timings, and a starting velocity model,

ii. Calculate the expected arrival times of the source events at the nodes by using the recorded source event times and the calculated travel times,

iii. Calculate node clock drift parameters by time regression, deriving a relationship between clock time, clock frequency offset and the difference between calculated arrival times and observed arrival times to best reduce the difference between expected arrival time and observed arrival time,

iv. Calculate the expected arrival times with the updated node clock drift parameters,

v. Updating node position and velocity model to best reduce the difference between expected arrival time and observed arrival time.

2) The method of claim 1 whereby an estimate of node position is achieved by assuming that the difference between node time and source event timing is stationary during a limited time span.

3) The method of claim 1 whereby an estimate of node time variation is derived from observed clock error in several spans.

4) The method of claim 1 whereby the derivation of relationship between source event timing versus node timing, and position-velocity updates are repeated.

5) The method in claim 1 whereby data the results are improved by:

- Making use of bathymetric data,

- Using data from multiple nodes is used to improve the accuracy of the velocity model.

6) The method as described in claim 1 whereby the velocity model allows for:

- lateral variation,

- anisotropy,

- temporal and spatially varying water column.

7) The method described in claim 1 whereby the regression of source event timing and observed seismic event timing in the node includes:

- temperature,

- working age of the clock crystal.

8) The method described in claim 1 whereby outlier rejection, currents, source array dimensions, tidal heights, environmental variations and overlapping source events are accounted for.

9) The method described in claim 1 whereby the node functionality is extended to record electromagnetic waves.