Spectral analysis for seismic interpretation
Abstract
processing of seismic data. SUBSTANCE: given approach makes it feasible to improve quantitative evaluation and visualization of minor effects of tuning of thin seismic stratum and other side disturbances of continuity of rock. Reflection from thin stratum is characteristically expressed in frequency region and can generate information on thickness of stratum. This reflection has periodic sequence of marks in its amplitude spectrum. These marks are spread by distance which is inversely proportional to temporary thickness of thin stratum. Apart from it this characteristic expression can be employed to follow reflections of thin stratum within limits of 3-D volume and to evaluate thicknesses of strata and their side extent. Salient feature of invention lies in generation of set of seismic paths distributed over specified volume of ground and in selection of specified volume. At least one section of seismic paths is converted with the aid of orthonormal conversion which gives collection of conversion factors. Obtained factors are organized in tuning cube that can be employed jointly with computer for identification of localization of hydrocarbons. It can also result in plotting of map for prospecting for oil and gas. EFFECT: improved quantitative evaluation and visualization of minor effects of tuning of seismic thin stratum. 31 cl, 16 dwg
Term
Term ended
Expired 12 November 2017, 8.9 years ago.
- Priority
- Filed
- Granted
- Expired
- Today
31 claims: 10 independent, 21 dependent
- 1Способ обработки сейсмических данных, отличающийся тем, что он включает в себя следующие операции:(а) получение отображения набора сейсмических трасс, распределенных в пределах заданного объема толщи земли, причем указанные сейсмические трассы содержат цифровые выборки, которые являются характерными по меньшей мере по времени, положению и амплитуде;(б) выбор части указанного объема и содержащихся в нем сейсмических трасс для ограничения зоны интереса в пределах указанного объема;(в) преобразование по меньшей мере одного из участков указанных сейсмических трасс в пределах указанной зоны интереса с использованием дискретного ортонормального преобразования, причем указанное дискретное ортонормальное преобразование дает множество коэффициентов преобразования;(г) организация указанных коэффициентов преобразования в куб настройки, используемый при сейсмической разведке, причем указанный куб настройки может быть использован совместно с компьютером для идентификации потенциальной локализации углеводородов.
- 2Способ по п.1, отличающийся тем, что он дополнительно включает в себя следующую операцию:(д) вывод на индикацию по меньшей мере одного из участков указанного куба настройки.
- 3Способ по п.1, отличающийся тем, что он дополнительно включает в себя следующую операцию:(е) запоминание указанных коэффициентов преобразования в форме, удобной для индикации в виде куба настройки.
- 4Способ по п.1, отличающийся тем, что указанное дискретное ортонормальное преобразование в операции (в) характеризуется наличием множества ортонормальных базисных функций и применено к окну, которое содержит указанные цифровые выборки, для получения множества коэффициентов преобразования, объединенных с указанными ортонормальными базисными функциями.
- 5Способ по п.1, отличающийся тем, что указанное дискретное ортонормальное преобразование представляет собой преобразование Фурье.
- 6Способ по п. 4, отличающийся тем, что весовую функцию применяют в пределах указанного окна, содержащего цифровые выборки, ранее применения преобразования за счет указанного дискретного ортонормального преобразования.
- 7Способ по п.6, отличающийся тем, что указанная весовая функция представляет собой Гауссовскую весовую функцию.
- 8Способ по п.4, отличающийся тем, что указанное дискретное ортонормальное преобразование в операции (в) характеризуется наличием множества ортонормальных базисных функций и применено к окну, которое содержит указанные цифровые выборки, для получения множества коэффициентов преобразования, объединенных с указанными ортонормальными базисными функциями, а операция (г) включает в себя операцию масштабирования указанных коэффициентов преобразования в пределах указанного куба настройки.
- 9Способ по п.8, отличающийся тем, что операция масштабирования указанных коэффициентов преобразования в пределах указанного куба настройки включает в себя следующие операции:(i) выбор базисной функции из указанных в операции (в);(ii) выбор по меньшей мере двух коэффициентов преобразования, соответствующих базисной функции операции (i);(iii) вычисление комплексной величины всех коэффициентов преобразования, выбранных в операции (ii);(iv) вычисление статистической величины из всех величин коэффициентов преобразования, вычисленных в операции (iii);(v) вычисление величины масштабирования из статистической величины;(vi) умножение указанных коэффициентов преобразования операции (ii) на указанную величину масштабирования.
- 10Способ по п.9, отличающийся тем, что указанная статистическая величина представляет собой арифметическое среднее всех вычисленных указанным образом величин коэффициентов преобразования.
- 11Способ по п.2, отличающийся тем, что операция (д) включает в себя операцию записи визуально различимых изображений, отображающих указанный куб настройки на главным образом плоской среде.
- 12Способ по п.11, отличающийся тем, что он дополнительно включает в себя операцию использования указанных визуально различимых изображений для идентификации структурных и осадочных характеристик нижних горизонтов, обычно связанных с захватом и накоплением углеводородов.
- 13Способ по п. 12, отличающийся тем, что он дополнительно включает операцию вывода на индикацию по меньшей мере одного из участков указанного куба настройки, при которой осуществляют запись визуально различаемых изображений, отображающих указанный куб настройки на главным образом плоской среде и операцию использования этих изображений для идентификации структурных и осадочных характеристик нижних горизонтов, обычно связанных с захватом и накоплением углеводородов, которые содержатся по меньшей мере в одном тонком пласте.
- 14Способ по п.3, отличающийся тем, что операция (д) включает в себя операцию запоминания по меньшей мере одного из участков указанного куба настройки в запоминающем устройстве с произвольной выборкой компьютера.
- 15Способ по п.2, отличающийся тем, что операция (д) включает в себя следующие операции:(i) выбор ортонормальной базисной функции;(ii) выбор по меньшей мере двух коэффициентов преобразования, соответствующих указанной ортонормальной базисной функции: (iii) определение положения для каждого выбранного коэффициента преобразования;(iv) вывод на индикацию указанных коэффициентов преобразования в местоположениях, представительных для указанных положений.
- 16Способ по п.15, отличающийся тем, что операции (i) - (iv) повторяют для множества ортонормальных базисных функций.
- 17Способ обработки сейсмических данных с использованием цифрового компьютера, при котором производится обработка пространственно связанных сейсмических трасс в цифровой форме, причем трассы характеризуются по меньшей мере временем, положением и амплитудой, причем цифровой компьютер запрограммирован для выбора части указанного набора пространственно связанных сейсмических трасс для выделения зоны интереса, отличающийся тем, что он включает в себя следующие операции:(а) преобразование по меньшей мере одного из участков указанных сейсмических трасс в пределах зоны интереса с использованием преобразования Фурье, которое характеризуется множеством ортонормальных базисных функций и которое применяют к окну, содержащему указанные цифровые выборки, для получения множества коэффициентов преобразования, объединенных с указанными ортонормальными базисными функциями;(б) организация указанных коэффициентов преобразования в объем спектрального разложения;(в) применение масштабирующего значения к указанным коэффициентам преобразования для получения масштабированного объема разложения, причем этот масштабированный объем разложения является представительным для среднего объема, полученного при помощи следующих операций: (i) выбор базисной функции из функций, полученных в операции (а);(ii) выбор по меньшей мере двух коэффициентов преобразования, соответствующих базисной функции по п.(в) (i);(iii) вычисление комплексных величин для всех коэффициентов преобразования, выбранных по п.(в) (ii), и (iv) вычисление среднего значения из всех величин коэффициентов преобразования по п.(в) (iii);(г) вывод на индикацию указанного масштабированного объема разложения.
- 18Способ обработки сейсмических данных, отличающийся тем, что он включает в себя следующие операции:(а) получение отображения набора пространственно связанных сейсмических трасс, распределенных в пределах заданного объема толщи земли, причем указанные сейсмические трассы содержат цифровые выборки, которые характеризуются по меньшей мере временем, положением и амплитудой;(б) выбор части указанного объема и указанных пространственно связанных сейсмических трасс, содержащихся в нем, для ограничения зоны интереса в пределах указанного объема;(в) очерчивание окна в пределах указанной зоны интереса, причем это окно имеет номер стартовой выборки и включает в себя множество цифровых выборок;(г) преобразование по меньшей мере одного из участков указанных пространственно связанных сейсмических трасс в пределах зоны интереса с использованием дискретного ортонормального преобразования, которое характеризуется множеством ортонормальных базисных функций и которое применяют к указанному окну, содержащему указанные цифровые выборки по п.(в), для получения множества коэффициентов преобразования, объединенных с указанными ортонормальными базисными функциями;(д) организация указанных коэффициентов преобразования в куб настройки, причем указанный куб настройки и коэффициенты преобразования в нем объединены с указанным номером стартовой выборки указанного окна по п.(в);и (е) вывод на индикацию указанного куба настройки, причем операции (в), (г) и (е) повторяют по меньшей мере для одного или нескольких заданных окон, для получения множества кубов настройки.
- 19Способ по п.18, отличающийся тем, что операция (е) включает в себя следующие операции:(i) выбор ортонормальной базисной функции;(ii) выбор куба настройки из указанного множества кубов настройки;(iii) экстрагирование из указанного выбранного куба настройки множества коэффициентов преобразования, объединенных с указанной выбранной ортонормальной базисной функцией;(iv) повтор операций (ii) и (iii) по меньшей мере еще для одного выбранного куба настройки;(v) организация указанных экстрагированных коэффициентов преобразования в куб настройки единственной ортонормальной базисной функции;(vi) вывод на индикацию указанного куба настройки единственной ортонормальной базисной функции.
- 20Способ по п. 19, отличающийся тем, что операцию (v) осуществляют путем упорядочения указанных экстрагированных коэффициентов преобразования при помощи указанного номера стартовой выборки, объединенного с указанным заданным окном.
- 21Способ по п.19, отличающийся тем, что указанное дискретное ортонормальное преобразование представляет собой преобразование Фурье, причем операцию (v) осуществляют путем организации указанных экстрагированных коэффициентов преобразования в куб настройки единственной частоты.
- 22Способ по п. 18, отличающийся тем, что в пределах указанного окна применяют весовую функцию, содержащую цифровые выборки, ранее осуществления преобразования при помощи указанного дискретного ортонормального преобразования.
- 23Способ по п.22, отличающийся тем, что указанная весовая функция представляет собой Гауссовскую весовую функцию.
- 24Способ по п. 18, отличающийся тем, что операция (е) дополнительно включает в себя операцию записи визуально различимых изображений, отображающих один или несколько из указанных кубов настройки на главным образом плоской среде.
- 25Способ по п.24, отличающийся тем, что он дополнительно включает в себя операцию использования указанных визуально различимых изображений для идентификации структурных и осадочных характеристик нижних горизонтов, обычно связанных с захватом и накоплением углеводородов.
- 26Способ обработки сейсмических данных, отличающийся тем, что он включает в себя следующие операции:(а) обращение при помощи компьютера к набору данных, который включает в себя сейсмические трассы, распределенные в пределах заданного объема толщи земли, причем указанные сейсмические трассы содержат цифровые выборки, которые характеризуются по меньшей мере временем, положением и амплитудой;(б) выбор множества пространственно связанных сейсмических трасс;(в) выбор зоны интереса в пределах указанного выбранного множества пространственно связанных сейсмических трасс;(г) преобразование по меньшей мере одного из участков указанных пространственно связанных сейсмических трасс в пределах зоны интереса с использованием дискретного ортонормального преобразования и получение множества коэффициентов преобразования;(д) организация указанных коэффициентов преобразования в куб настройки;(е) вычисление множества значений сейсмических атрибутов из указанных коэффициентов преобразования куба настройки;и (ж) вывод на индикацию указанных значений сейсмических атрибутов в местоположениях, отображающих положения в пределax указанного куба настройки.
- 27Карта для разведки нефти и газа, полученная способом по п.26, предназначенная для использования при разведке углеводородов, когда сейсмические данные, содержащие отраженную сейсмическую энергию, записывают как функцию времени в пределах заданного объема толщи земли, для получения множества пространственно связанных сейсмических трасс, содержащих выборки, которые характеризуются по меньшей мере временем, положением и амплитудой, отличающаяся тем, что указанная карта содержит:(а) главным образом плоскую среду для записи на ней визуально различимых изображений;и (б) по меньшей мере одно визуально различимое изображение на указанной главным образом плоской среде, причем указанное визуально различимое изображение отображает указанные вычисленные значения сейсмического атрибута.
- 28Способ выработки сейсмических атрибутов для использования при разведке углеводородов, отличающийся тем, что он включает в себя следующие операции:(а) получение отображения набора сейсмических трасс, распределенных в пределах заданного объема толщи земли, причем указанные сейсмические трассы содержат выборки, которые характеризуются по меньшей мере временем, положением и амплитудой;(б) выбор части указанного объема и указанных сейсмических трасс, содержащихся в нем, для ограничения зоны интереса в пределах указанного объема;(в) преобразование по меньшей мере одного из участков указанных сейсмических трасс в пределах указанной зоны интереса с использованием очень короткой временной функции дискретного ортонормального преобразования и получение коэффициентов преобразования;(г) организация указанных коэффициентов преобразования в куб настройки;(д) вычисление множества сейсмических атрибутов из указанных коэффициентов преобразования в указанном кубе настройки.
- 29Способ по п.28, отличающийся тем, что он дополнительно включает в себя следующие операции:(е) масштабирование указанных коэффициентов преобразования в пределах указанного куба настройки;и (ж) инвертирование указанного куба настройки с использованием инверсной функции дискретного ортонормального преобразования для получения отфильтрованной версии указанного преобразованного участка указанных пространственно связанных сейсмических трасс.
- 30Способ по п.29, отличающийся тем, что указанная функция преобразования по п. (в) содержит множество ортонормальных базисных функций, причем операция (е) включает в себя следующие операции:(i) выбор одной из базисных функции и по меньшей мере двух соответствующих ей коэффициентов преобразования;(ii) вычисление комплексной величины указанных выбранных коэффициентов преобразования;(iii) вычисление статистического значения из указанной комплексной величины;(iv) вычисление значения масштабирования из указанного статистического значения;(v) применение указанного значения масштабирования к указанным выбранным коэффициентам преобразования.
- 31Способ по п.30, отличающийся тем, что указанное статистическое значение является арифметическим средним указанных комплексных величин коэффициентов преобразования.
Independent claims31
99 paragraphs, as filed
The present invention mainly relates to a method of quantifying and visualizing subtle seismic effects settings for the thin (low-power) of the formation. Disclosed in the present invention, the method comprises decomposition of the reflected seismic signal at its frequency components by using discrete (or other orthonormal) Fourier transform of length (duration), depending on the thickness of the formation to be determined. After decomposition by said discrete transform of use, resulting coefficients are organized and displayed in such a manner that the detected and amplified expression characteristic frequency range reflection events from a thin bed. This allows you to visualize the variation in the thickness of the underlying layer, which otherwise can not be detected. The present invention allows the seismic interpreter to analyze and display a map (map) and geological subsurface stratigraphic features as a function of spatial position, a travel time and frequency, to an extent that was previously unattainable.
By all standards of geophysics exploration is a relatively young science, as the earliest work in this area have appeared in the 20s, and the updated CMP approach only in the 50s. However, since its inception geophysics exploration in the oil industry became the dominant approach in finding oil deposits. Although exploration geophysics that encompasses the three broad regions, namely, a gravitational, magnetic and seismic exploration, seismic currently dominant method is to such an extent that other forms of intelligence are not practically used. In practice, a simple count of the number of seismic exploration in the party became accepted measure of the state of the entire oil industry.
Seismic survey shows try mapping the subsurface of the earth by sending sound energy into the earth and registration "echoes" that come back from the underlying rock layers. The source of the down sonic energy can be burst or seismic vibrators on the surface of the earth or sea surface airgun. During the seismic energy source is moved across the surface of the geologic structure of interest. Each detonation source it is recorded in a large number of points on the surface. Many combinations of explosion / recording is then combined to create an almost continuous profile of the subsurface that can extend for many kilometers. When a two-dimensional (2-D) seismic survey point registration signal usually lie on a straight line, whereas in a three-dimensional (3-D) survey the registration point distributed on the surface in the form of a lattice. At the simplest explanation we can assume that the 2-D seismic line gives a cross-section of each layer, existing under the reception area, and a 3-D survey gives a "cube" of data or the amount of which, at least conceptually, it is a 3-D picture of the lower horizons underlying the shooting area. It should be borne in mind that it is possible to extract individual 2-D line recording of 3-D data volume.
Seismic survey formed a very large number of individual seismic recordings or traces. In a typical 2-D survey is usually a few tens of thousands of traces, whereas in a 3-D survey the number of individual traces may be many millions. Seismic trace is a digital recording of the sound energy reflected from inhomogeneities in the subsurface, and the partial reflection occurs whenever a change of acoustic impedance of the subsurface materials. The digital samples are usually obtained at intervals of 0.004 seconds (4 ms), but often also used and intervals of 2 ms and 1 ms. It should be borne in mind that each sample on the seismic trace is associated with a travel time of wave, and in the case of reflected energy - with a double wave travel time. In addition, the position of each seismic trace on the surface of the earth carefully registers and is usually part of the track (part of the header information of the route). This allows further correlate the seismic information contained in the routes, with specific subsurface locations, thereby obtaining means contouring seismic data and extraction of attributes (characteristic features) to display a map (mapping). The signal that is sent into the earth, is called seismic vibrations or seismic wavelet. Kind of seismic vibrations is dependent on whether the source of air gun, dynamite or a vibrator. The term "signature source" or "source pulse" is generally used to describe the characteristics of the seismic recording of a particular seismic vibrations.
Generated on the surface seismic source pulse starts immediately spread from the point of origin, including into the earth, encountering blocks of rock in the subsurface and passing through them. Each of the interface between two different rocks there is a possibility of seismic reflection. The magnitude of the seismic energy reflected from the specific surface section, which depends on the acoustic impedance contrast between the units of the rock and the reflection coefficient is a known measure of contrast. It is believed that the reflection coefficient is the ratio of the amplitude of the reflected wave to the amplitude of the incident wave. In terms of rock properties: wherein the acoustic impedance of a rock unit is defined as the mathematical product of the rock density (ρ1 and ρ2 are respectively and densities of the upper and lower rock units) on the rate of the signal in said rock. V1 and V2 correspond to velocities in the upper and lower blocks of rock. (Strictly speaking, this expression is completely true only if the wavelet section intersects the surface of the rock at vertical incidence. However, in practice it is usually assumed that the requirement of verticality is satisfied if the wavelet crosses the interface with a deviation approximately 20o from the vertical).
<IMG>
The reflected energy that is recorded on the surface of the earth, can be shown conceptually as a convolution of the seismic imlulsoida reflectivity of the subsurface feature, in the form of so-called "folding model". In short, the model is used to explain the convolution of the seismic signal recorded at the surface as the mathematical convolution downward (into the earth) source wavelet reflectivity function that displays the reflection coefficients at the interfaces between different rock layers in the subsurface. This can be expressed by the equation: x (t) = w (t) • f (t) + n (t) where x (t) is a registered seismogram, w (t) is the seismic source wavelet, f (t) is the reflectivity function of the earth, and n (t) is random ambient noise. "*" Is the mathematical convolution. Additionally convolutional model requires, inter alia, (1) that the source wavelet remains invariant as it passes through the lower layers (i.e., that it be fixed and unchanging), and (2) to record on the surface seismic trace can be represented as the arithmetic sum of the individual parcel source wavelet from each interface in the lower levels (the principle of "superposition" in which the reflection wavelet and signal propagation is a linear system). Despite the fact that few believe that the convolution model completely describes the mechanics of wave propagation, this model is accurate enough for many tasks, which makes it very useful in practice. Some details of construction of the folding model discussed in Section 2.2 the publication "Seismic data processing", Ozdogan Yilmaz, Society of Exploration Geophysicists, 1987, the disclosure of which is incorporated herein by reference.
Seismic data are received and processed properly, provide a wealth of information scout - a specialist oil company whose work is the localization of possible drilling sites. For example, seismic profile gives the Scout a look at the structure of the rock layers of the subsurface and often reveals important characteristics involved in the capture and storage of hydrocarbons, among many other such features such as faults, folds, anticlines, inconsistencies and underlying salt domes and reefs. During the computer processing of seismic data is usually obtained evaluation signal velocity in the subsurface, as well as detect and indicate close to the surface heterogeneity. In some cases, seismic data can be used to directly estimate rock porosity, water saturation and its evaluation for hydrocarbon content. Less obvious is that the attributes of the seismic vibrations, such as phase, peak amplitude, peak to trough, and many other, often can be empirically correlated with known hydrocarbon content, and this correlation is applicable to seismic data collected for new exploration areas. In brief, seismic data provides the best of the existing structural and stratigraphic information regarding the subsurface without drilling.
Seismic data are linked but one fundamental limitation: blocks of rock that are relatively "thin", often do not have sufficient (accurate) resolution. More specifically, while the seismic reflection data can provide maps close to "geological cross-section" subsurface if lithographic layers are relatively "thick" the seismic image, which is obtained from the "thin" layers is much less clear. This phenomenon is known in the art as the seismic resolution problem.
Seismic resolution in the present context refers to vertical resolution within a single seismic trace, and can be defined as the minimum distance between two seismic reflectors in the subsurface that can be identified on the seismic recording as a single interface, rather than as a single composite reflection. As an example, indicate that the block subsurface can be identified ideally on a seismic section as a combination of two quantities: the individual reflection, an outgoing (received) from the top of the block, and the second separate reflection possible with opposite polarity emanating from the base unit. Ideally, the vertex and the base block seen on the recorded seismogram as separate and isolated reflectors that can be individually "are spaced apart in time" (i.e., are marked and identified) in the seismic section, and the seismic data in the interval between two spaced time peaks contain information relative to the intermediate block of rock. On the other hand, if a seismic unit is not sufficiently thick (strong), the return reflection from the top and bottom of the unit overlap, thus creating interference between the two reflection events and blurring occurs subsurface. Blurred image of the subsurface is one example of a phenomenon known in the art as the problem of "thin bed".
FIG. 1 shows in a general form as a problem by using a thin layer of axioms folding model. Consider first reflection from a "thick" layer, shown in Figure 1a. In the left part 1a shows a wavelet source generated at the surface. Source wavelet extends unchanged through the earth ground path P1, so long as it meets the surface of the section "A" block of rock. (It should be borne in mind that this wave trajectory 1a actually are vertical, but a better understanding of the inclined shown. It is usually used in practice). When the downward seismic vibration meets the interface "A" part of its energy is reflected back toward the earth's surface along the path P2 and is recorded on the surface as a reflection of the event R1. It should be borne in mind that the wavelet R1 has reversed polarity compared to the polarity of the source wavelet that displays a negative reflection coefficient of the surface of "A". This polarity reversal is shown as an example only, as experts know that the reflection coefficients may have any polarity.
The remainder of the downward energy (after the partial reflection at the interface "A") continues through the thick layer until it collides with the surface of section "B" at the base of the thick lithographic unit. Upon reaching the surface of the section "B" of the wavelet energy goes deeper into the earth along the path P5 and the rest of the energy is reflected along a path P4 to the ground, where it is registered as a reflection of R2. Note that the reflection from interface "B" occurs later in time than the reflection from interface "A". Current separation in time between the two events depends on the layer thickness between two surfaces and the section from the sound propagation velocity in the layer, the thicker layers or slower rate provides a greater time separation between the reflections from the top and bottom layer. Assessment of the thickness of the layer is the time required for the passage of seismic wavelet this thickness.
At the earth's surface is actually recorded composite (compound) reflected from the thick layer, which is an arithmetic sum (superposition) of the two return reflections given separation in time of these two events. Since the two reflected wavelet have temporal overlap, then the resulting seismic record clearly indicated both events, the information bearing on two discrete (individual) horizons. (It will be appreciated that the separation in time between the two reflected events displayed in Figure 1A is not to scale. The experts know that the time separation is actually twice the characteristic time layer thickness).
Referring now to Figure 1b, which shows the reflection from a thin bed. In this case, the source wavelet is generated on the surface of the earth and is sent along the path P6 until it encounters an interface "C" block of rock. (As with the previous case, the trajectory of the waves in the figure are actually vertical). As shown in Figure 1b, when the downward oscillation meets surface seismic section "C", a part of its energy is reflected back along the path P7 in the direction to the surface where it is recorded as reflection R3. The remainder of the energy continues to move down through the thin bed until the face surface of the section "D". Upon reaching the interface "D" part of the wavelet energy goes deeper into the earth along the path P10, while the remainder of its energy is reflected back along the path P9 to the surface where it is recorded as reflection R4.
In this case, the reflection from interface "D" occurs later in time than the reflection from interface "C", but the separation in time between the two reflections in the case of a thin bed is less because there is less distance that should run the wave front of its reflection from interface "D". In fact, the time separation between the two reflections is so small that the return (reaching to the ground) wavelet overlap. Since in this case the composite reflected from the thin layer is an arithmetic sum of the two return reflection, the actually recorded signal is an event that is not reflected clearly shows from the top and the base unit, so it is quite difficult interpretation. Said composite reflected an uncertain event is an example of a typical problem of a thin layer.
Needless to say that the thickness of the subsurface layer explored has significant economic value to scout the oil company because, other things being equal, the greater the thickness of the lithography unit, the greater the volume of hydrocarbons that it could potentially contain. Taking into account the importance of precise determination of layer thickness, it is not surprising that suggested many approaches to solve the problem of a thin layer.
The first approach, which is used almost universally, is the shortening of the length of the seismic wavelet as wavelet longer generally provide poorer resolution than short. During the data processing phase the recorded seismic wavelet can often be shortened dramatically by the application of well known signal processing techniques. For example, the specialists are well aware that the usual projected inverse convolution can be used for purification of the wavelet spectrum. Similarly, the known technique of wavelet processing, including inverse convolution source signature, as well as a number of other approaches that may alternatively be used in an attempt to achieve the same end result, namely, a more compact shape fluctuations. Although any of these treatments may lead to dramatic changes in the character of the seismic section and may shorten the length of the wavelet significantly, often need to use additional signal processing operation.
Even the best signal processing is ultimately only postpone the inevitable: no matter how compact wavelet, will always be of economic interest rock layers that are too thin for that wavelet so that it can identify them properly. Therefore, it uses a different approach, which is aimed more at the analysis of the nature of the composite reflection. This approach is based on the observation that, even when there is only a single composite reflection and the thickness of the layer can not be directly determined, yet can obtain information by means of the recorded seismic data that may indirectly be used to estimate the actual thickness of the lithographic unit .
For example, in Figure 4a shows the familiar "squeezed" seismic model, in which the thickness of the stratigraphic unit of interest (and the thickness is measured from the time-wavelet) is reduced to extinction (ie "shrinking") on the left end of the drawing. 4b is a set of mathematically derived synthetic seismograms calculated for this model, which illustrate the noise released from the convolution of the seismic wavelet, with interfaces which cover the layer. Note that the right edge 4b, the composite signal recorded on the first trace shows that the reflector is clearly delimited by a negative reflection at the apex unit and a positive reflection at its base. When you move to the left in Figure 4b individual reflections from the top and bottom begin to merge into a single composite reflected vanish as the thickness of the interval to zero. Please note however, that the character of the composite reflection still continues to change even after degeneration of events in a single reflection. Thus, although there is little direct visual evidence that there was a reflection of the two interfaces, changes in reflections with decreasing thickness suggest that they contain information that is related to the thickness of a thin layer.
The pioneering work Uidessa in 1973 (Widess. "How thin thin layer?", Geophysics, Vol. 38, pp. 1176 - 1180) established a popular approach to the analysis of a thin layer, which uses calibration curves, which are obtained using both amplitude peak -vpadina composite reflected by the thin bed events and time division peak-to-trough to provide an indicative estimate of the thickness of the "thin" layer. (See. The publication Neidell and Poggiagliomi, "Stratigraphic models and their interpretation - geophysical principles and technologies" in the book "Application of seismic stratigraphy for the exploration of hydrocarbons", AAPG Memoir 26, 1977). The need for surgery in the calibration process is to set the amplitude of the "settings" for considered reflection from a thin layer, and the amplitude adjustment is obtained by a layer thickness at which a maximum constructive interference between the reflections from the top and bottom of the unit. At least in theory, the width of the setting depends on the dominant wavelength λ wavelet and equal λ / 2, where the coefficients of reflection from the top and bottom of the unit have the same sign, and is equal to λ / 4, when the reflection coefficients have opposite signs.
For the reason that the approaches of the calibration types are flexible, they successfully used in many different directions exploration. However, these are based on the amplitude and time of the calibration methods are highly dependent on the thoroughness of the seismic processing for setting an exact phase wavelet and to control the relative amplitudes of the transition from one to another seismic trace. However, those skilled in the seismic processing know how difficult it is to obtain seismic section, wherein all supported relative amplitudes. Furthermore, as described above, and based on the calibration method is not suitable for the study of thin bed reflections in wide 3-D survey; This method works best if it is applied to an isolated reflector on a single seismic line. Quite a challenge is to develop a calibration curve for a single line; but much more difficult to find a calibration curve that is appropriate for all the lattice 3-D seismic data.
In connection with the above, as is well known in the seismic processing and seismic interpretation, a need exists for a method devoid of the above problems, allowing to extract useful information from a relatively thin layer obtained in the usual way seismic data. Furthermore, this method should also preferably provide attributes for subsequent stratigraphic and structural analysis. It should be recognized, as was done by the present inventors that there is a real need for a method of processing seismic data that takes into account and allows to resolve the problems described above.
Before proceeding to describe the present invention should be however noted that the following description given with reference to the accompanying drawings, given only as examples (or preferred embodiments of the invention), and not limitative. Those skilled in the art, which concerns the present invention, the present invention can be applied in other forms without departing from but outside the scope of the appended claims. Finally, despite the fact that the invention disclosed herein is illustrated with reference to various aspects of the convolutional model, outlined hereinafter methods are not based on any particular model of the recorded seismic trace and work well in the same manner with significant deviations from the standard convolution model.
In accordance with the present invention discloses a new means of using the discrete Fourier transform and mapping to mapping of thin beds and other lateral gaps (discontinuities) of the rock for the conventional 2-D and 3-D seismic data. More particularly, the present invention is motivated by the observation that the reflection from a thin bed has a characteristic expression in the frequency domain, which carries information about the thickness of the reservoir; a homogeneous thin layer introduces a periodic sequence of marks (notches) in the amplitude spectrum of the composite reflection, these markers are moved apart at a distance that is inversely proportional to the time (measured by measuring the travel time of wavelet) the thickness of a thin layer. In addition, if the Fourier transform coefficients obtained properly, this characteristic expression may be used by the interpreter to track thin bed reflections in the 3-D volume and estimate the thickness of the thin bed and extended to such a degree that was previously unattainable. More generally, the methods disclosed herein can be used to detect and identify vertical and lateral discontinuities in the local rock mass. Furthermore, the usefulness of the method in accordance with the present invention is enhanced by the use of a new method of cleaning the frequency domain, which emphasizes the geologic information present in the spectrum. Finally, the present invention is also directed to the detection of seismic attributes that can be correlated. of interest with structural and stratigraphic characteristics of the subsurface, making it possible to obtain quantitative values that can be mapped intelligence and used to predict the presence of hydrocarbons in the subsurface or other accumulation of minerals.
As a general base in accordance with the present invention preferably uses a relatively short discrete Fourier transform to determine frequency components of a seismic trace. As is known, the computation of the Fourier transform of time series even though they have only real values, resulting in a complex Fourier transform coefficients in the form "A + Bi", where "i" indicates "imaginary" number or the square root of minus one. Furthermore, it is well known that the expression "A + Bi" can be written as: A + Bi = re-θ wherein θ = tan-1 (B / A). The expression "θ" is known as the phase angle (or simply as "phase") of the complex quantity A + Bi, the term "g" as amplitude, and the expression | A + Bi | module, also referred to as an absolute value. The frequency spectrum is obtained from the Fourier transform coefficients by calculating the complex value of each transform coefficient. Further, a digital value ('size') of each coefficient in the frequency spectrum of the output frequency is proportional to the initial data. Finally, after applying the Fourier transform to some particular time series, the resulting series of complex obtained coefficients in the frequency domain, while the data is not converted into the time domain.
<IMG>
The present invention is based on the general observation that a frequency spectrum calculated using Fourier transform for the entire route, tends to resemble the spectrum of the source wavelet, whereas shorter window spectra tend to display underlying geological information. This is because long analysis windows encompass a large geological variations, which form with the passage of time "white" (or random and uncorrelated) reflectivity function which has a "flat" amplitude spectrum. Therefore, the shape of the frequency spectrum calculated from the total seismic trace is largely dependent on the frequency content of the source wavelet. (See, e.g., Chapter 2. 1. 2. publication "Seismic data processing", Ozdogan Yilmaz, Society of Exploration Geophysicists, 1987 , incorporated herein by reference). On the other hand, in the case where the analysis window is so short that the earth reflectivity function is not white, then the resulting spectrum is the Fourier components which are dependent on both the wavelet and local geology of. We can say that such a small window geologically acts as a filter, reducing the spectrum of source wavelet and creating a not stationary spectra of small windows.
These ideas are displayed in a general form in Figure 2, which shows a typical seismic trace and some frequency spectra calculated therefrom. The upper part of Figure 2 shows the frequency spectrum of the Fourier transform for the entire seismic trace. This spectrum has the form of a typical wavelet field. However, the spectra calculated over shorter windows for and shown in Figure 2 below, are fixed and do not tend to display the underlying geology which may potentially change dramatically in a very short intervals.
The importance of this observation for the present invention is illustrated in FIG. 3, which generally shows two representative spectrum. The left shows the frequency spectrum of a typical broadband spectrum source wavelet. However right frequency spectrum displays generally expressed in the frequency domain for the composite thin bed reflection. In this latter case, the geology of the thin bed tend to act as a filter in the frequency domain and makes its own contribution to the frequency content of the reflected wavelet. As shown generally in Figure 3, the present inventors found that a homogeneous thin layer affects the amplitude spectrum of the reflected events due to input thereto the "marks" (nicks) or narrow attenuation bands of frequencies having a distinctive appearance. The homogeneous layer is a layer with a constant velocity of propagation and density along its entire length. Furthermore, the distance between marks so introduced is equal to the inverted "time thickness" thin bed, wherein a thickness time understand the time interval that is required for the passage of the wavelet layer in one direction (which is equal to the thickness of the layer divided by speed signal). Thus, the attenuated frequencies in the amplitude spectrum may be used to identify a thin bed reflection and to measure its thickness.
Referring to Figure 4, which received in the previous section results are extended to the analysis of a simplified 2-D geological model for which investigated the frequency domain expression for a thin layer. FIG. 4a shows a typical "compressed" reflectance function (geological model). 4c shows a grayscale image of the amplitudes of the frequency spectrum of the Fourier transform calculated from the model. This image is obtained by creating a series of time locations 50 arranged at equal intervals in the screen pattern, each of which has only two non-zero values: one corresponding to the reflection coefficient at the top layer, and the other - the reflection coefficient at the base layer. Then, calculation was made of the standard discrete Fourier transform for the time series and then calculating the complex magnitude of each coefficient.
FIG. 4c lighter portions correspond to larger amounts of the amplitude spectra, while the darker areas indicate lower values. This "mark" in the amplitude spectra shown darker values on the chart. This 4c shows, in a very real sense, the Fourier transform of the geology, and more particularly, the characteristic signature that is superimposed on the wavelet by this event. Most significant in this graph, in connection with the present invention is that by decreasing the thickness of the pattern interval between the marks increases. In addition, this model for the mark is periodic with a period equal to the time the layer thickness. Thus, if the signature is given, namely, the presence of a periodic frequency markers can be localized in the seismic survey, then this is a clear evidence of a thin bed.
According to a first aspect of the present invention provides a system for interpreting seismic data containing related with the existence of thin bed events, wherein the data are decomposed into a series of Fourier transform 2-D lines or 3-D volume, resulting in increased display layers of said thin formation. In this embodiment, a single Fourier transform window which is used for the isolation of the portion of seismic trace that intersects a zone of interest. This embodiment is shown mainly in Figure 5 as applied to 3-D seismic data, but those skilled will understand that the same method can be advantageously applied to 2-D sets of seismic traces to obtain a gain mapping thin bed reflections which It contains.
The first phase is obtained a set of spatially related to each other, the seismic traces. These traces may comprise, for illustrative purposes only, one or more short records, assembling of seismic traces with constant displacement, mounting assembly seismic traces, VSP survey, a two dimensional seismic line, two-dimensional arranged one above the other seismic line extracted from a 3-D seismic shooting, or, preferably, 3-D plot of 3-D seismic survey. Moreover, the present invention can also be applied to a 2-D or 3-D survey in which data is transposed, i.e. have "offset" or spatial axis (axis "X" or "Y" 3-D data) that oriented in such a manner that replaces the vertical axis or an axis "Time". More generally, any 3-D volume of digital data may be processed using the disclosed according to the present invention methods. In view of the above, for simplicity, the vertical axis is hereinafter referred to as a time axis, despite the fact that experts understand that the digital samples might not be separated by time blocks. Regardless of the selection, the present invention is most effective when applied to a group of seismic traces that have a spatial relationship to some subsurface geological characteristics. Again, for illustration purposes only, the following discussion will be conducted in terms of routes, which are contained in the stack traces when the 3-D shooting, although it can also be used either recruited a group of spatially related seismic traces.
As shown generally in Figure 5, then the selected area of interest within a particular 3-D volume. The zone of interest might be, for example, undulating region bounded by two selected reflectors, as shown in Figure 5. In this case, the reflector is preferably flattened or converted into a reference level (i.e., made flat by shifting the time of individual traces up or down) until analysis and probably Palinspastic reconstructed. Generally it may be set specifically associated time interval (e.g., from 2200 up to 2400 ms), resulting in a "cube", or more precisely, "box" of seismic data within the 3-D volume, namely subobem. Further, the lateral extent of the zone of interest may be limited by specifying the track "in-line" and "cross-line" (in the longitudinal and transverse directions). Also other methods of specifying the zone of interest is known to the authors of the present invention.
The selection and extraction of data, which correspond to the area of interest, known as a step of finding a subset of data (Figure 5). Criterion, which is used in selecting the zone of interest is the desire to have as short (in time) area. This is due to the previously described general philosophy under which the spectra of the Fourier transform long window tend to resemble the source wavelet and Fourier spectra short window tend to contain more information related to geology. It should be borne in mind that there is a "hidden" expansion box, which is often automatically and unexpectedly applied to the windows of the Fourier transform, namely the extension of the window size to a length of degree two. This lengthening of the window is made to improve computational efficiency, since the window with the duration of degree two are candidates for the application of the algorithm of fast Fourier transform (FFT). However, in accordance with the present invention is not used, this widespread practice, and a more general algorithm applies a discrete Fourier transform (although it is less efficient computational capabilities) that allows to keep the minimum possible value of the duration of the analysis window. Taking into account the processing power of computers today, there is little reason to not use the conversion data only within the region of interest.
This operation COMPUTE in Figure 5 ("calculate"), as applied to the present invention comprises at least one operation, namely, the computation of the discrete Fourier transform of the region of interest. The resulting spectral expansion coefficients of the zone of interest is then stored as part of the output volume of the spectral decomposition (such as "tuning cube") for consideration. It should be borne in mind that there is one path (i.e., a set of Fourier transform coefficients) in the output tuning cube volume for each seismic trace processed as part of the input data. It should also be borne in mind that in this preferred output constructing horizontal slices through the volume contain coefficients corresponding to a single common Fourier frequency.
Optionally, the operation COMPUTE may contain additional operations which have the potential to improve the quality of the output volume and subsequent analysis. First of all, the weighting function can be applied to the seismic data within the zone of interest previously computed transformation. The object of the weighting function is a compression or smoothing the data within the Fourier analysis window, thereby reducing the distortion in the frequency domain, which may arise in the analysis window, such as "car". Using the weighting function before the transformation is well known in the art. Advantageously the weighting function in accordance with the present invention is Gaussian in shape, which in many ways optimal for this application. In view of the foregoing it will be appreciated that potentially can be used, and other weighting functions.
Furthermore, since it is usually the amplitude spectrum of greatest interest to the scout, the amplitude spectrum can be calculated from the transform coefficients after they are moved into an auxiliary storage area. Alternatively, the phase spectrum or any other derived attribute can be calculated from the transformation coefficients before being sent to storage, and these calculations made by the present inventors.
Finally, as part of a calculation operation, individual frequency scaling may be applied to each plane (i.e. frequency) on the output screen. As shown generally in Figure 10, the present inventors have found a preference for a separate scale for each frequency slice in the output volume to obtain the same average value before viewing. This type of scalability is one type that can be used, but the inventors prefer this method because it allows to select geologic content stored frequency spectra at the expense of common wavelet information.
After calculating and storing the spectra, they are ready for use in the geophysical exploration for thin beds. It should be borne in mind that in the subsequent display data is important that each spectrum was arranged and examined in the same spatial relationship with the other spectra as the track from which they were calculated. This means that not present in the transformed data space ratio should be maintained in the ratios of transform coefficients. A preferred method of viewing the transform coefficients starts with their formation within the 3-D "volume" (tuning cube), naturally, provided that the input data were originally taken from the 3-D volume. However, it should be understood that the vertical ("z") axis is no longer a "temporary" as it was before the transformation, but rather are now conventionally unit displays the frequency in memorizing the coefficients of the Fourier transform.
FIG. 5 shows that the last step is to consider the tuning cube similarly consider any type of conventional 3-D volume of seismic data. In view of the foregoing, the present inventors have found that the consideration of successive horizontal slices within the scope of the coefficients is preferred to localize and visualize thin bed effects. Note that in the tuning cube using a horizontal section displays all of the coefficients which correspond to a single Fourier frequency, and therefore is a cross section of constant frequency. Furthermore, as an additional tool in the analysis of the data contained within this volume, the present invention advantageously used animation ("recovery") of a series of horizontal species within the volume. In the case where the area of interest is most portion of individual seismic line rather than a volume, the resultant image, which is a set of spectra of the Fourier transform of spatially related seismic traces are displayed in their original spatial relationship, can be considered as a cube configuration, although technically it can not be a "cube" of data.
Animation successive horizontal slices within the spectral volume is the preferred method of reviewing and analyzing the transform coefficients, said animation preferably carried out on a computer monitor having a high speed workstation. As is well known in the art, the animation in the form of interactive panning within the volume is a fast and effective means of treating large volumes of data. The amount of data can be viewed in horizontal, vertical or oblique slices, each of which is a single view of data. However, more importantly in the context of the present invention, a quick review successive horizontal slices one after another, wherein the diagnostic agent is prepared for the recording of a large volume of data and to identify a thin bed reflections therein, as discussed below. It should be noted that the disclosed method is preferable to have the order of the sections in terms of frequency (when it is clear ascending or descending) when producing animations and consideration of these sections.
According to a second aspect of the present invention provides a system for processing seismic data allows amplify seismic records associated with the existence of thin bed events, wherein the data are decomposed into a series of Fourier transform 2-D lines or 3-D volumes by using a series of overlapping Fourier transforms for short windows, whereby the display layers reinforce said thin layer. This embodiment is shown primarily in Figure 6 as applied to 3-D seismic data, but those skilled will understand that the same method can be advantageously applied to 2-D sets of seismic traces to enhance the display of thin bed reflections that therein contained. As shown in Figure 6 and as previously discussed, the first step of this embodiment involves the interpreter mapping the temporal bounds seismic zone of interest. As previously described, resulting in mapping of available seismic data cube or rectangular element inividualnoy seismic line.
In accordance with the present invention does not use Fourier transform of a single window for each trace, and is used instead of Fourier transforms series of overlapping short window. The duration and degree of overlapping windows is altered to a particular application, however in this case the window length need not be equal to "2", but rather should be chosen so as to obtain the best image the underlying geology. It will be appreciated that optional data before the conversion may be applied to a weight function within each short window, and, as before, is preferred Gaussian weighting function.
As shown in Figure 6, after calculating the Fourier transform for each short-window coefficients derived from it are stored in isolation within the individual tuning cube that maintains communication with the source short window. It should be understood that there may be as many tuning cubes however has overlapping windows analysis. Scaling, if applicable, is applied in isolation for each frequency plane in each tuning cube.
Each short-window tuning cube produced by using a sliding window may now be individually analyzed in the same way as has been suggested previously for the first embodiment. In this case, each cube advantageously treated by means of horizontal slices or constant frequency images, thereby obtaining renderer geological changes with changing frequency. Furthermore, since there is now a collection of tuning cubes calculated for different time trace points, in fact get a set of tuning cubes that overlap the depth range of the subsurface.
Finally, in accordance with a third aspect of the present invention provides a system for processing seismic data allows amplify seismic records associated with the existence of thin bed events, wherein the data are decomposed into a series of Fourier transform 2-D lines or 3-D volumes by using Fourier transform for short windows, with subsequent reorganization into cubes with a single tuning frequency, whereby the display layers reinforce said thin layer.
As shown generally in Figure 7, the first operation in the present embodiment, the operation is repeated two preceding embodiments, namely, the data is first interpreted, then partitioned into subsets. Thereafter, from the seismic data within the zone of interest is calculated for a series of Fourier transformations of short overlapping windows prior to transformation using optional weighting function or application of compression within each window. As in the previous embodiment, accumulated coefficients from each short window transform. However, in this case, the calculated coefficients of the Fourier transform is not considered as a tuning cubes and reorganized into single frequency energy cubes which are then analyzed in a horizontal or vertical plane to separate thin bed effects.
More specifically, the reorganization of the present invention conceptually includes the extraction of all the cubes customize each horizontal slice, which corresponds to a specific frequency. Thereafter, the individual sections with the same frequency "are stacked" so that the topmost slice containing coefficients calculated from the topmost sliding window, the next slice containing coefficients calculated from following the upper sliding window, etc. . It should be borne in mind that after the reorganization of the volume coefficient is organized into blocks of "xy" and time block. This is because the vertical axis is the "time" of the sliding window, which gives a special coefficient.
To use information by tuning cubes with a single frequency generated in the preceding operation, the seismic interpreter must select their corresponding frequency and the seismic volume (for example, he can choose the volume ratio corresponding to 10 Hz and / or volume for 11Hz and so on. D .). Each cube constant frequency can be considered as a top or horizontal or in any other manner resulting in a renderer geological changes in the lateral direction for a specific frequency.
Importantly, for all the previously mentioned options, the fact that the original route is not converted spatially related provides additional benefit. More specifically, it is well known that the Fourier coefficients of a short window itself sufficiently noisy and have a poor frequency resolution in comparison with the conversion of a long window. One approach which isp Use This Criterion in accordance with the present invention to improve the reliability of the transformed data is to apply a Gaussian weight function to the data before transformation. However, another important approach in accordance with the present invention is to display the coefficients within a volume in the same spatial relationship as the input data. Since in this manner the displayed route contain spatially correlated information, they display next to each other allows the observe to visually "aliasing" noise and to identify the underlying coherent signal information.
Finally, despite the fact that the present invention is discussed herein in terms of a discrete Fourier transform, in reality the Fourier transform it is just one of a number of discrete time data transformations that could be used similarly. Basic operations, namely, (1) computing a short window transformation (2) combining the resulting coefficient in volume, and (3) the volume of investigation for finding the thin bed effects may be implemented in a variety of discrete data transformations other than the Fourier transform. In these transformations the volume setting is formed by grouping coefficients corresponding to the same base function. Thus, when used as the "single frequency tuning cube" it is not only the tuning cube formed by conventional Fourier transform coefficients, but also any other tuning cube formed by the coefficients of the functions on the same basis.
Experts know that the discrete Fourier transform is only one of a plurality of discrete linear unitary transformations that satisfy the following requirements: (1) they are linear operators that (2) is exactly inverted, and (3) their basic functions form an orthonormal set. In terms of equations, if x (k), k = 1, L, displays the time series, and X (n) represents its "n-ing" transformed value, n = 1, L, then the forward transform of the time series for this class of transformations can be written as: where A (k: n) displays a core direct transformation or set of basis functions. Furthermore, there is an inverse transformation which allows of transformed values back to the original data: where B (k: n) indicates an inverse transform kernel. In accordance with the requirements of orthonormality domestic product between the two basic functions should be zero, and the value of each basic function must be equal to one. This requirement is in compressed form can be expressed by the following equations: where A * (k: n) represents the complex conjugate of A (k: n). For the discrete Fourier transform basis functions that correspond to the direct transformation of length L, typically selected as a set of complex exponentials: A (k; n) = {e-2πikn / L, k = 0, L-1} (5) In this case there L basis functions (or basis vectors), one basis function for each value of "n": n = - (L / 2), ..., 0, ..., (L / 2-1) To summarize, it can be said that each transform coefficient X (n), calculated from a data window corresponds to a particular basis function, and the volume settings established set of all transform coefficients corresponding to a particular area of interest, which are stored in the auxiliary storage area in the same spatial relationship as line, which was calculated from each window.
<IMG>
<IMG>
<IMG>
<IMG>
<IMG>
As another specific example, you can specify that the discrete transformation Welsh (Walsh) can be used instead of the Fourier transform, the coefficients of Welsh in the same way are grouped, displayed and analyzed. As previously described, the conversion can be calculated Welsh within overlapping series of sliding windows and the coefficients obtained may be accumulated and organized in cubes settings. In this case, the calculated transformation coefficients are not displayed frequency and magnitude, referred to as "sequent". In this configuration cubes "sequence identity" may be formed from the transform coefficients Welsh exactly the same way as it happens in the formation of the Fourier tuning cubes. In the future, these cubes setting "the same frequency" (or, more generally, of the same basic functions) will be referred to only configure cubes orthonormal basis function.
Finally, despite the fact that the discrete Fourier transform is a transformation that is characterized by a set of orthonormal basis functions, application of a non-trivial weight function to the basis functions previously calculating transformation destroys their orthonormality. In accordance with the theory, a weight function that is applied within a window that should be regarded as applying to the basis functions rather than the data, whereby the integrity of the underlying data is stored. However, basis functions that were orthogonal before application of the weighting functions are usually not those after its application. In view of the above, regardless of whether the weight function is applied to data or to the basis functions, the final calculation result after the transformation will be the same.
One possible exception of minor theoretical dilemma that arises when using a weighting function with a discrete orthonormal transformation, is the choice of a combination of orthonormal transformation / weight, which is not subject to the specified action. For example, the local cosine (and local sine) transform is a discrete orthonormal transformation, for which the weighting function is selected as a smooth contraction of the special form, wherein orthonormality of the basis functions is preserved by some loss in frequency resolution. Moreover, the rationale of the local cosine / sine transformation serves as a natural bridge to the theoretical field of general wavelet transformation.
1 is a diagram that illustrates in general form a thin bed problem.
FIG. 2 shows a typical seismic trace and a comparison of the spectra of long and short windows, computed from it.
Figure 3 illustrates in general form as expressed in the frequency domain impact seismic wavelet to a thin layer.
FIG. 4 shows a simple model seismic compressed, and convolution of the response mapping in the frequency domain response of said convolution.
5 is a diagram which shows the general approach according to a preferred embodiment of the present invention.
FIG. 6 schematically shows the preferred embodiment of the present invention may be used in exploration.
FIG. 7 schematically illustrates another preferred embodiment of the present invention.
FIG. 8 is a diagram which illustrates a preferred embodiment of the present invention.
FIG. 9 is a diagram that describes the appearance of a thin bed during animation of constant frequency slices.
Figure 10 illustrates a general approach that is used to scale the constant frequency slices to identify the geologic content of the transformed data.
FIG. 11 is a diagram which illustrates another preferred embodiment of the present invention.
In accordance with the present invention there is provided a method of processing seismic data using a discrete Fourier transformation, by means of which increases the usefulness of this transformation for detecting thin beds.
According to a first aspect of the present invention there is provided a method for enhancing and consideration thin bed effects using a discrete Fourier transform wherein a single Fourier transform is calculated for a window spanning the zone of interest and the coefficients obtained by further indicate otherwise. Suppose, as shown in FIG. 5 that x (kjnt) displays a 3-D volume of seismic data, where k = 1, k and j = 1, j represent indices that identify a specific trace within a given 3-D volume. As an example, indicate that these may be the number of the longitudinal and transverse position, although possible other layouts. The variable "nt" is used to display the time (or depth) position within each seismic trace, nt = 0, NTOT - 1 - total number of samples in each individual track. Temporary separation between successive values of x (kjnt) (ie, sampling frequency) is designated as Δt, where Δt is typically measured in milliseconds. Therefore, each track in a 3-D volume contains a record (NTOT) * Δt milliseconds data, the first sample is usually taken at a "zero" time. In view of the foregoing one can understand that seismic data suitable for analysis according to the present invention are not ordered in terms of "time". For example, a sample of the seismic data processed by the program depth migration, stored (preserved) within a seismic trace in order of increasing depth Δz. However, the present invention is fully applicable to such data sample. Therefore, in the following description, the term Δt (and "time") will be used in a broad sense, as referring to the separation between successive digital samples, regardless of the possible forms of measurement division.
As a first step, a scout or seismic interpreter selects a zone of interest within the 3-D volume. This can be accomplished, for example, by digitizing the peak time of seismic events or digitizing table or, more commonly, at a seismic workstation computer. When an event is detected scout trying to find the same characteristics of the reflector (eg, peak, trough, the moment of crossing the zero line, etc.) for each seismic trace; the ultimate goal is to develop a computer file, which contains information about the time and location surface that tracks an event within the 2-D section or a 3-D volume. As shown in Figure 11, if such information can be a computer program to read the peak and finding the zone of interest for any trace within the data volume and / or for implementing the method according to the present invention. This program may be entered into a computer, such as a magnetic disk, tape, optical disk or CD-ROM.
Alternatively, the interpreter might simply specify constant starting and ending that limit event of interest throughout the entire volume, resulting in creating of interest "cube," wherein the term "cube" is used herein in a generic sense to display the 3-D sub-volume of the original a 3-D volume of intelligence. As an example, the following discussion assumes that was removed the 3-D sub-volume, though you can see that the same technique that is discussed below, can be adapted for non-permanent time-window. Again, as an example used equipment is assumed that the time zone of interest after extraction ranges from the first sample in the 3-D sub-volume of sample to the last number "N". Similarly, it is also assumed that the area of interest is present on every track in the sub-volume, although it is known that the area of interest often takes only a part of the 3-D volume.
After selecting a zone of interest is the selection of the next operation window length "L" of the Fourier transform. In the most general form of the length of the conversion must not be longer than is absolutely necessary for the environment (capture) the zone of interest. Typically, the length of the Fourier transform is selected for reasons of computational efficiency and is usually restricted integer power of 2 (e.g., 32, 64, 128, etc.) that enables efficient algorithm for computing the FFT, and not less efficient Fourier transform mixed roots or much less efficient Discrete Fourier Transform. However, in the present context, it is not recommended to increase the length of the selected window, as is usually done, to an integer power of 2; Instead, use the discrete Fourier transform. In view of the foregoing, it should be understood that in the following discussion, for any mention of a discrete Fourier transform of opportunity may have been manufactured by FFT calculation. Otherwise, it may be selected by the general discrete Fourier transform or some variant thereof with mixed roots, if the selected window length is not an integer power of 2.
Before the start of the Fourier transform should create an auxiliary storage volume for calculation of Fourier coefficients. L storage capacity of computer words must be set to store the calculated transformation coefficients for each track, with even more storage capacity is required in the case of storing the values of the seismic data and converted the results of a double (or higher) accuracy. To explain it can be said that the Fourier series in real time requires a memory length L L / 2 complex data values, each of which normally requires two computer words of storage. (In fact only [(L / 2) -1] complex data values, rather than L, because a series of real Fourier transform coefficients which correspond to positive and negative frequencies are directly related to each other and form a complex conjugate pairs. In addition addition, there are two real values: the coefficient at zero ("dc") hertz and the coefficient at the Nyquist frequency, both of which may be stored in a single complex value. Finally, if L is an odd integer, the number of data values equal to [ (L + 1) / 2]. If there are a total (multiplied by K J) seismic traces in the zone (cube) of interest, the minimum required amount of auxiliary storage overall measured in computer words equal to the product of L, J and K. Let that the lattice A (kjnt) displays an auxiliary storage area for the present embodiment.
In a first step of calculating, as shown in Figure 8, the data values within the zone of interest are extracted from an input trace x (j, k, nt), taken in the sub-volume: y (nl) = x (j, k, nt ), nl = 0, L-1 and is used optional weighting function: y (nl) = y (nl) w (nl), nl = 0, L = 1, in which the lattice y (nl) is an area for temporary storage. (It should be borne in mind that in this embodiment, the length of the analysis window equals the length of the storage area). The weight function w (t), or data window as it is referred to some may take various forms. Some of the most popular data windows are windows Hamming, Hanning, Partsa, Bartlett and Blackman. Each of the window functions has certain advantages and disadvantages. However, the inventors have found that the optimum for a given application, for many reasons is to use a Gaussian window. The Gaussian weighting function can be determined by the following equations: w (nl) = στe- (nl-1/2) 2 / σ, nl = 0, L-1 (6) where Generally, the weight function should be a real function of its range not equal to zero.
<IMG>
After applying a weighting function to perform calculations of the discrete Fourier transform according to the following standard expression n = - (L / 2), ..., 0, ..., (L / 2-1), where X (n) indicates a complex factor Fourier transform at frequency f (n), which depend on the window length L. It is known that by means of Fourier transform coefficients are obtained, which give an estimate of the spectral amplitude at the following Fourier frequencies: n = - (L / 2), ..., 0 , ..., (L / 2-1) should be borne in mind that the nominal sampling rate Δt may not correspond to the sampling frequency at which the data were obtained in the field. For example, a common practice is resampled seismic trace at a higher sample rate to save storage when there is little useful information at the highest recording frequency. On the other hand, one of the seismic traces can be produced at a low sampling frequency samples, in combination with other lines, certain higher sampling rate. In any case, the nominal frequency of the data samples may not accurately reflect its true spectral bandwidth. With a simple modification of the previous expression may be obtained contingency: fn = (n / L) Fmax, n = - (L / 2), ..., 0, ..., (L / 2-1) which represents Fmah the highest frequency that is available in the data.
<IMG>
<IMG>
Since the seismic trace is a "real" function (that is, not imaginary), then its Fourier transform is symmetric and Fourier coefficients correspond to positive and negative frequencies, related the following expressions: RE [X (fn)] = RE [X (fn) ] and IM [X (fn)] = - IM [X (fn)] where RE [Z] is a function that extracts the real part of the complex value z, a IM [Z] is a function that extracts its imaginary part. In accordance with the above expressions in each window of the Fourier transform is obtained only the value (L / 2-1). Therefore, for reasons of specificity, the following discussion will consider only the positive frequencies, although it is understood that similar results can be obtained by using only the negative frequencies.
The next step of the process involves placing the calculated values of the complex frequency in the array of auxiliary storage. These tracks are filled calculated complex frequency values as follows: A (j, k, i) = X (i), i = 0, L / 2 where "j" and "k" with the same parameters, the respective initial route data. In practice the lattice constant a (j, k, i) can not be stored entirely at a time in the RAM (a memory with random access RAM), and can be located in whole or in part, on tape or disk, optical disk or other medium Storage. Furthermore, since this embodiment preferred indication thin bed requires more of the frequency spectrum rather than the complex values, it is convenient to simultaneously calculating the integrated value at room each coefficient in the supporting grid storage: A (j, k, i) = | X ( i) |, i = 0, L / 2. However, in many cases, the complex coefficients are necessary and useful, as shown in Figure 8; then mostly the complex coefficients stored in the auxiliary storage grid.
This procedure is repeated for each trace in a given sub-volume, filling the auxiliary grid storage conversion factors prepared for consideration by the scout. However, pending the results of the data again advantageously scaled differently, thereby isolate geologic information into transform coefficients with respect to the contribution of the wavelet. General Procedure, which is used in said frequency domain scaling is illustrated in Figure 10. Scaling method disclosed herein is designed to equalize the average spectral amplitude in each frequency slice, thereby producing a purified tends wavelet spectrum. As shown in FIG. 8, let T (j, k, i) is an auxiliary storage grid, which can be stored the entire tuning cube. For a given frequency slice i may be calculated from the average spectral amplitude: The spectral values are calculated so that T (j, k, i) is potentially a complex value. In the next step the values in this particular frequency slice are set so that their average is equal to some user specified constant value represented by the variable AVG, in accordance with the expression: T (j, k, i) - [AVG / TAVG] (j, k, i), j = 1, J, k = 1, K wherein the primary record is used to indicate that the lattice T (j, k, i) has been modified. In practice, AVG is set to some specific numerical value, such as 100. This scaling operation is repeated in isolation for each frequency slice i = 0, L / 2 in the tuning cube volume. Upon completion of this operation, each slice has the same average amplitude, while achieving some spectral balance. It will be appreciated that this form of single frequency scaling is just one of many scaling algorithms, which, in accordance with the present invention can be applied to the tuning cube data. As an example, indicate that instead of calculating the arithmetic mean of the elements in the slice can be used as other assessment of central tendency or any other statistic (e.g., median, mode, geometric mean, variance, and the like). As another example, that instead of setting the average value in each frequency slice equal to some of the same constant, each slice can have a different constant average value, thereby amplify certain frequencies in the spectrum and suppresses others.
<IMG>
If you invert the data cube after setting the zoom back to the time domain by using inverse Fourier transform of the standard, it is possible to obtain a spectrally balanced version of the original input seismic traces. Let X (k) represents a set of scaled transform coefficients obtained by means of the previously disclosed process, taken from location (j, k) within the scaled tuning cube lattice. Then it may be spectrally purified version of the input data in accordance with the following expression: n1 = 0; L-1 (11) wherein x '(j, k, nl) displays are now modified (spectrally balanced) version of the input data x (j, k, nl). Divider w (nl) used to remove the effects of the weight function, which was used prior to conversion. This element may be omitted if no weight was applied in the direct Fourier transform.
<IMG>
However, rather than inverting the scaled tuning cube, it can be advantageously used as an exploration tool for detecting thin beds. After processing all the runs and placing them into an auxiliary storage area, horizontal (with the constant frequency) amplitude slices Si (j, k), corresponding to the i-th frequency may be extracted from A (j, k, i) for later review and / or animation: Si (j, k) = | A (j, k, i) | When producing animation (quick successive consideration) specified routes, the thin layers can be identified as such events, which are alternating between high and low amplitude values. In addition, for many types of thin layers exist characteristic pattern of moving marks that clearly indicate that the event is caused by a thin layer. It will be appreciated that the method disclosed herein, preferably to produce, in case of considering animations and slicing them ordering frequency (if clear ascending or descending).
FIG. 9 shows the source of the diagnostic mobile pictures. FIG. 9a shows geologic thin bed model lens type, and FIG. 9b shows a stylized Fourier transform of this model, which applied only to the mark. As mentioned previously, marks are periodic with a period that is equal to the inverse temporal thickness of the model at that point. Now imagine a model 9a as displaying a 2-D cross section of the 3-D (disc type) radially symmetric model and 9b - how similarly radially symmetrical set of one-dimensional Fourier transforms said the 3-D model. If the constant frequency plane, called Plane 1, pass through the volume as indicated, in a plan view in this plane can detect low amplitude circular region corresponding to the first mark. Plane Plane 2 passes through two notches and two annular region gives a low amplitude. Finally, the plane Plane 3 passes through three marks and three ring region gives low amplitude. If we now consider these sections in their rapid alternation in ascending order of frequency, the visual impression of the picture appears "bulls-eye" in which the rings are moved outward from the center. This pattern is characteristic of moving markers (diagnostic) for thin layers.
If a thin layer is not a ring, then watch another picture. Instead, there is a series of concentric circles marks that go from thicker zones in the direction of thin zones. For example, let the model 9 is a cross section of the flow path of the lens mold. When considering sequential frequency slice will be observed over the entire length of the channel moving picture out of marks that are moving from the center to the periphery of the channel.
It will be appreciated that if the thin bed is not homogeneous, for example if it contains a stepwise increase or decrease in speed, it can not give a characteristic pattern of "markers" of a homogeneous thin bed, but rather create some expression in other frequency domain. In such cases, the preferred method of identifying the characteristic response is to create a model of the event and calculating the Fourier transform as was illustrated previously in Figure 4b. With this information, the scout can then analyze an animated cube setting for moments predictable response.
Painting markers is not only a qualitative indication of a homogeneous thin bed, but also provides a quantitative assessment of thin bed events. Note that the mark 9 are bounded laterally outermost edges of the pattern. Therefore, when panning the stack of slices of frequency and the determination of the outermost limits of movement marks can be obtained by quantitative evaluation of thin bed events.
Described is an amazing visual effect that can easily be observed in the volume of existing seismic data. Since the event is not a typical thin bed has some definite and slowly changing amplitude spectrum, the thin bed response is significantly different from it and can be easily detected. It should be borne in mind that in the case of a single window for calculating the total area of interest, the actual time position (i.e., depth) of the thin bed within the zone of interest is not important. At any location of the thin bed within the temporal zone of interest, the spectrum for that window will contain a characteristic pattern of moving marks. Experts will easily understand that moving the location of events in time does not change the amplitude spectrum. This movement soon enters a phase change, which is not visible when considering and calculating the amplitude spectrum. Instead of indicating the amplitude spectra in animated plan view in this embodiment can also use any number of other attributes calculated from the complex values stored in the tuning cube. For example, the phase of the complex transform coefficients provides another means of identifying thin bed events and, more generally, lateral discontinuities rock mass. The phase tuning cube may be calculated as follows: where P (j, k, i) comprises a part of the phase of complex Fourier transform coefficients for every point in the original tuning cube. Specialists have long used the phase section to highlight fuzzy reflectors, and the phase section to help highlight the continuity of the seismic data. However, in accordance with the present invention, lateral discontinuities in the spectral phase response carrying information on the lateral variability of the local rock mass, which is the first example of the truncation of thin beds. When viewed in animated form on top of the phase values near the lateral edges are relatively "unstable", with distinct from the normal first derivative. Therefore, the boundary layers thin and, more generally, lateral discontinuities in the rock mass (e.g. faults, cracks, discrepancies, etc.) will have a phase that contrasts with the surrounding phase values and therefore relatively easy to be identified. This behavior can be used in isolation to identify the lateral limits or together with a cube setting the amplitude spectrum for finding variability in local rock mass.
<IMG>
Finally, the present invention contemplates that the disclosed therein tuning cube technology can provide a further understanding of seismic reflection data. It can be produced on the output display and examination of the tuning cube (which contains the phase or amplitude data) to obtain empirical correlation with the underlying rock contents, rock properties from, and the underlying structure or layer stratigraphy. Alternatively, it may be made more manipulation of Fourier transform values stored in the tuning cube, to generate new seismic attributes that can be useful in exploration settings position. As an example, indicate that the attributes that can be calculated from the tuning cube values include the average spectral magnitude or phase, as well as many other attributes. The importance of this aspect of the present invention can best be shown as follows. In the seismic interpretation known that spatial variations in a seismic reflector characteristic can often be empirically correlated with changes in reservoir lithology or fluid content. Since the exact physical mechanism that gives rise to this variation accurately reflect the nature of the unknown, the accepted practice is to calculate the interpreters of the set of seismic attributes, and their subsequent application to the graph or map, and search for an attribute that has a predictable value. Attributes obtained by calculation from the bottom of the settings, display the localized analysis of the properties of the reflector (if the calculation of a short window) and as such potentially have a significant importance in the promotion of interpretation.
According to a second aspect of the present invention, a method of processing seismic data, allowing to strengthen the effects associated with the existence of a thin layer, using the discrete Fourier transform, according to which the calculations of Fourier transforms of a series of sliding short windows arranged in a window overlapping zone of interest and then a new output to the display. This method is shown generally in Figure 6 and in more detail in Figure 8. Conceptually, this embodiment of the present invention may be represented as a series of tuning cubes creating the type described above, one cube for each configuration of the Fourier transform window position specified by the user.
In this case, x (k, j, n) displays a 3-D volume of seismic data, a "L" - the length of the selected short window Fourier transform. In accordance with this embodiment, "L" is typically selected considerably shorter than the length of the zone of interest N. As before, the Fourier transform window length should be selected not on the basis of computational efficiency, but rather so as to obtain better display of certain classes of thin bed events in the lower horizons. As an example, you can specify that a reasonable length of conversion is starting such a length that is sufficient for covering "very thick" thin layer in the area of interest. It will be appreciated that it may be necessary to increase this minimum length in circumstances such as when the waveform is not particularly compact. In this latter case, the minimum length of the windows may be extended to the length of the wavelet measured in samples.
Every citation, both ways
| Document | Relation | Office | Cited during |
|---|---|---|---|
| RU2751088C2 | Cited by | Russian Federation | Search report |
| RU2718135C1 | Cited by | Russian Federation | Search report |
| WO2021049970A1 | Cited by | World Intellectual Property Organization (WIPO) | International search |
23 members in 9 offices
Priority claims3
| Document | Office | Kind | Date |
|---|---|---|---|
| 75965596 | United States of America | A | |
| 08759655 | – | – | – |
| US19960759655 | – | – | – |
Members23
| Document | Office | Kind | |
|---|---|---|---|
| CA2244714A1 | Canada | A1 | |
| WO9825161A2 | World Intellectual Property Organization (WIPO) | A2 | |
| WO9825161A3 | World Intellectual Property Organization (WIPO) | A3 | |
| NO983523D0 | Norway | D0 | |
| NO983523L | Norway | L | |
| EP0882245A2 | European Patent Office (EPO) | A2 | |
| US5870691A | United States of America | A | |
| CN1210591A | China | A | |
| CA2359579A1 | Canada | A1 | |
| WO0043811A2 | World Intellectual Property Organization (WIPO) | A2 | |
| AR011288A1 | Argentina | A1 | |
| US6131071A | United States of America | A | |
| EG21105A | Egypt | A | |
| WO0043811A3 | World Intellectual Property Organization (WIPO) | A3 | |
| NO20013539D0 | Norway | D0 | |
| NO20013539L | Norway | L | |
| EP1149311A2 | European Patent Office (EPO) | A2 | |
| RU2187828C2This record | Russian Federation | C2 | |
| CN1170170C | China | C | |
| CA2244714C | Canada | C | |
| NO328897B1 | Norway | B1 | |
| NO332387B1 | Norway | B1 | |
| CA2359579C | Canada | C |
1 legal event, as the office reported them to INPADOC
Events
| Event | Code | |
|---|---|---|
| The patent is invalid due to non-payment of feesMM4A | MM4A |
Numbers
- Publication, DOCDB
- 2187828
- Publication, EPODOC
- RU2187828
- Application
- 9811695228
- Application, DOCDB
- 98116952
- Application, EPODOC
- RU19980116952
Titles
- English
- SPECTRAL ANALYSIS FOR SEISMIC INTERPRETATION
Classification
- CPC, 3
- G01V1/32
- G01V1/301
- G01V2210/48
- IPC, 2
- G01V1 30
- G01V1 32