CN120652538A - Shale reservoir azimuth anisotropy time difference solving method and device - Google Patents
Shale reservoir azimuth anisotropy time difference solving method and deviceInfo
- Publication number
- CN120652538A CN120652538A CN202411956930.4A CN202411956930A CN120652538A CN 120652538 A CN120652538 A CN 120652538A CN 202411956930 A CN202411956930 A CN 202411956930A CN 120652538 A CN120652538 A CN 120652538A
- Authority
- CN
- China
- Prior art keywords
- point
- reflection
- shot point
- shot
- offset
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
本申请属于地震勘探领域，涉及该领域的地球物理专业的地震资料处理环节，具体为一种页岩储层方位各向异性时差求取方法及装置，该方法包括：基于每个炮点对应的所有地震波反射信号的分布特征，得到偏移道集；分析地震波不同传播路径之间的差异，以及地震波传播过程中的噪声，得到白噪声抑制系数、距离可信度系数，进而获得炮点各方位角的偏移叠加函数；得到每个检波点的成像道函数，获取互相关判断函数；得到炮点各方位角的每个检波点的最佳时间平移大小；对相应检波点的地震波反射信号进行处理，得到每个炮点各方位角在各反射点对应反射位置的地震数据特征、各向异性特征。本申请旨在获取不同方位地震波数据之间的方位时差。
This application belongs to the field of seismic exploration and involves the seismic data processing link of the geophysical profession in this field. Specifically, it is a method and device for obtaining the azimuthal anisotropic time difference of shale reservoirs. The method includes: obtaining an offset track set based on the distribution characteristics of all seismic wave reflection signals corresponding to each shot point; analyzing the differences between different seismic wave propagation paths and the noise during the seismic wave propagation process to obtain a white noise suppression coefficient and a distance credibility coefficient, and then obtaining an offset superposition function for each azimuth angle of the shot point; obtaining an imaging track function for each detection point and obtaining a cross-correlation judgment function; obtaining the optimal time shift size for each detection point at each azimuth angle of the shot point; processing the seismic wave reflection signals of the corresponding detection points to obtain the seismic data characteristics and anisotropy characteristics of each azimuth angle of each shot point at the corresponding reflection point. This application aims to obtain the azimuthal time difference between seismic wave data of different azimuths.
Description
Technical Field
The application relates to the field of seismic exploration, in particular to a shale reservoir azimuth anisotropy time difference solving method and device.
Background
The azimuth anisotropy time difference characteristic plays a crucial role in shale reservoir description, provides information about complex structures and physical characteristics in shale reservoirs, and generally has the purposes of crack identification and evaluation, pore structure analysis, reservoir dessert prediction, ground stress analysis, seismic imaging quality improvement and the like, so that the accuracy of a geological model can be improved to a great extent by deeply researching the azimuth anisotropy time difference characteristic of the shale reservoirs, and effective exploration and development of shale oil and gas resources are promoted.
In shale reservoir description, the current azimuth anisotropy time difference solving method has the problem of strong data dependence in analyzing azimuth anisotropy time difference, only high-quality seismic data can be used as a basis for accurately analyzing azimuth anisotropy parameters, data acquired by actual seismic are often influenced by factors such as noise, accuracy of the azimuth anisotropy time difference is reduced, uncertainty of inversion results is large, and meanwhile, the problem of poor interpretation of the seismic data is caused by neglecting coupling relations of different azimuth data in reflected seismic waves.
Disclosure of Invention
In view of the above, it is necessary to provide a shale reservoir azimuth anisotropy time difference obtaining method and device, which solve the problems of large data error and poor data interpretability in the seismic data obtaining process of the shale layer.
The application provides a shale reservoir azimuth anisotropy time difference solving method, which comprises the following steps:
Acquiring seismic wave reflection signals of all wave detection points in each reflection point in each azimuth of each shot point;
Obtaining offset gathers of each reflection point at each azimuth of each shot point based on the distribution characteristics of all seismic wave reflection signals corresponding to each shot point;
analyzing differences among seismic wave propagation paths corresponding to different elements in the offset gather to obtain propagation path difference values between two elements in the offset gather;
Obtaining white noise suppression coefficients of the elements based on the change and the numerical difference of the propagation path difference values between the elements and the rest elements of each offset gather, comparing travel time of the seismic waves from the shot point to the reflection point and from the reflection point to each detection point, and obtaining the distance reliability coefficient of each element in the offset gather corresponding to the reflection point;
fusing the white noise suppression coefficient and the distance reliability coefficient of each element, and weighting the derivative results of each element of all offset gathers based on the fusion results to obtain offset superposition functions of all azimuth angles of shot points;
Obtaining an imaging channel function of each detector according to the seismic wave signal change characteristics of each detector of each azimuth of each shot point at all reflection points, and obtaining a cross-correlation judging function of the offset superposition function of each azimuth of each shot point and the imaging channel function of each detector under different time delays in a preset time period;
processing the seismic wave reflection signals of the corresponding detection points according to the optimal time translation, and synthesizing the seismic wave reflection signals processed by all the detection points of each azimuth of each shot point to obtain the seismic data characteristics of each azimuth of each shot point at the corresponding reflection position of each reflection point;
and obtaining the anisotropic characteristic of each shot point based on the distribution of the seismic data characteristics of all azimuth angles of each shot point at the corresponding reflection position of each reflection point.
Preferably, the offset gather is specifically a set of seismic wave reflection signal intensities of each reflection point at corresponding moments of all detection points of each azimuth angle of the shot point.
Preferably, the obtained propagation path difference value between two elements in the offset gather is specifically:
based on the distance between the shot point and the wave detection point, combining the propagation speed of the seismic wave to obtain the horizontal offset time from the shot point to the wave detection point;
According to the travel time of the seismic waves from the shot point to the corresponding reflection point of any offset gather, and then to the travel time of each element in any offset gather, combining the horizontal offset time to obtain the propagation triangle from the shot point to each element in any offset gather;
the difference of the areas of the propagation triangles between the different elements in each offset gather is taken as the propagation path difference value between the different elements in each offset gather.
Preferably, the horizontal offset time from shot to geophone is obtained by dividing the euclidean distance between shot and geophone by the propagation velocity of the seismic wave.
Preferably, the acquiring process of the propagation triangle specifically includes:
The travel time from the shot point to each offset gather corresponding to the reflection point is used as a first right-angle side of the propagation triangle, the horizontal offset time from the shot point to each detection point on each azimuth angle is used as a second right-angle side of the propagation triangle, and the travel time from each offset gather corresponding to the reflection point to each detection point is used as a hypotenuse of the propagation triangle.
Preferably, the obtaining the white noise suppression coefficient of each element includes:
And obtaining a fusion result of a negative correlation mapping result and the numerical value difference value of the propagation path difference value, normalizing the reciprocal of the accumulation sum of the fusion result of each element of each offset gather and all the remaining elements, and obtaining the white noise suppression coefficient of each element of each offset gather.
Preferably, the distance reliability coefficient of each element in the offset gather corresponding to the reflection point is specifically the square of the ratio of the first right-angle edge to the hypotenuse of the propagation triangle corresponding to each element in the offset gather corresponding to the reflection point.
Preferably, the formula of the offset superposition function of each azimuth of the shot point is as follows:
The kth element of the offset superposition function corresponding to the kth offset gather of each azimuth of the shot is marked as I (x, y, T k, n), and the formula is as follows: Wherein x and y represent the abscissa and ordinate of the shot at the earth surface, T k represents the depth time of the seismic wave from the shot to the kth reflection point, N represents the number of the detection points in each azimuth of the shot, N represents the total number of the detection points in each azimuth of the shot, WF n represents the white noise suppression coefficient of the nth element in the kth offset gather of each azimuth of the shot, DK n represents the distance reliability coefficient of the nth element in the kth offset gather of each azimuth of the shot, f n ′(τs,k+τg,n) represents the first derivative value of the nth element in the kth offset gather of each azimuth of the shot, τ s,k represents the travel time of the seismic wave from the shot to the kth reflection point, and τ r,n represents the travel time of the seismic wave from the kth reflection point to the nth detection point.
Preferably, the imaging channel function of each wave detection point is specifically the first derivative of the seismic reflection signal of each wave detection point at all reflection points.
Preferably, the cross-correlation judgment function of the offset superposition function of each azimuth of each shot point and the imaging channel function of each detector point in a preset time period is obtained, and the specific formula is as follows:
Wherein phi j,n (tau) represents the cross-correlation function values of G j (T) and G n,j (T+tau), wherein G n,j(T+τ)＝pn(T+τ)-gj(T+τ);gj (T) represents the offset channel function of each azimuth of the jth shot point, G j (T+tau) represents the offset channel function of the jth shot point after each azimuth time offset tau, p n (T+tau) represents the imaging channel function of the jth shot point after each azimuth time offset tau and X 1、X2 represents the starting time and the ending time of the seismic wave from the shot point to the depth of the target shale reservoir respectively;
taking the difference value of the square of the cross-correlation function and 1 as the denominator of the cross-correlation judgment function, and taking the cross-correlation function as the numerator of the cross-correlation judgment function.
Preferably, the optimal time shift of each wave detection point of each azimuth of the shot point is specifically an independent variable corresponding to the maximum function value of the cross-correlation judging function.
Preferably, the seismic data of each azimuth of each shot point at the reflection position of the corresponding reflection point is specifically a mean value of the seismic wave reflection correction signals of all the detection points acquired by each shot point at each moment.
Preferably, the anisotropic characteristic of each shot point comprises an anisotropic time difference accumulation result corresponding to the azimuth of the maximum anisotropic time difference;
taking the azimuth corresponding to the maximum value of the seismic data characteristic of all azimuth angles of each shot point at the reflection position corresponding to each reflection point as the azimuth corresponding to the maximum anisotropic time difference of each shot point at the reflection position corresponding to each reflection point;
and taking the mean value of the seismic data characteristics of all azimuth angles of each shot point at the corresponding reflection positions of the reflection points as an anisotropic time difference accumulation result of each shot point at the corresponding reflection positions of the reflection points.
In a second aspect, the embodiment of the application also provides a shale reservoir azimuth anisotropy time difference obtaining device, which comprises a memory, a processor and a computer program stored in the memory and running on the processor, wherein the processor realizes the steps of any one of the methods when executing the computer program.
The application has at least the following beneficial effects:
The method comprises the steps of setting shot points and wave detection points, collecting seismic wave data, converting the seismic wave data into an offset gather, representing system noise difference in the seismic wave propagation process by calculating propagation path difference values aiming at the condition that elements in the offset gather are interfered by noise, calculating white noise suppression coefficients by combining data differences among element values, screening out elements greatly influenced by white noise in the offset gather, calculating distance reliability coefficients by using seismic wave propagation paths corresponding to the elements in the offset gather, screening out elements greatly influenced by systematic noise in the offset gather, and finally carrying out superposition calculation on the offset gather by combining the white noise suppression coefficients and the distance reliability coefficients to obtain an offset gather function serving as a standard function for subsequent time correction. The white noise and systematic noise in the offset channel are respectively eliminated through the white noise suppression coefficient and the distance reliability coefficient, and the accuracy of the standard function is improved.
Furthermore, the azimuth anisotropy time difference information is represented through the offset set function and the cross-correlation function of the seismic wave data, so that the seismic wave data are aligned in time, the description precision of the underground geological structure can be greatly improved, particularly in three-dimensional seismic interpretation, a crack system model which is more in line with the actual situation is formed, the data are formatted uniformly, data noise caused by time delay is eliminated, and a seismic wave data set is obtained.
Compared with the prior art, the method eliminates the time delay interference of noise, utilizes the irregular travel time difference correction of the overlying stratum to carry out time difference correction on the seismic reflection on the target layer, peels off the influence of the stratum near the upper part of the target layer on the azimuth anisotropy time difference of the target layer, can effectively reduce the time error between seismic wave data, and is favorable for carrying out azimuth anisotropy analysis and azimuth inversion on shale reservoirs.
Drawings
FIG. 1 is a flow chart of steps of a shale reservoir azimuth anisotropy moveout calculation method according to an embodiment of the present application;
FIG. 2 is a schematic view of a shot point and surrounding detector placement at a measurement location according to an embodiment of the present application;
FIG. 3 is a schematic view of a seismic reflection signal of a detector at a reflection point according to an embodiment of the present application;
FIG. 4 is a schematic diagram of a seismic reflection correction signal for a geophone at a reflection point, according to an embodiment of the present application.
Detailed Description
In describing embodiments of the present application, words such as "exemplary," "or," "such as," and the like are used to mean serving as an example, instance, or illustration. Any embodiment or design described herein as "exemplary" or "for example" should not be construed as preferred or advantageous over other embodiments or designs. Rather, the use of words such as "exemplary," "or," "such as," and the like are intended to present related concepts in a concrete fashion.
Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application belongs. The terminology used in the description of the application herein is for the purpose of describing particular embodiments only and is not intended to be limiting of the application.
It should be further noted that the terms "first" and "second" in the present disclosure and the accompanying drawings are used to distinguish similar objects from each other, and are not used to describe a specific order or sequence. The method disclosed in the embodiments of the present application or the method shown in the flowchart includes one or more steps for implementing the method, and the execution sequence of the steps may be interchanged with each other, where some steps may be deleted without departing from the scope of the present application.
Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application belongs.
The application provides a shale reservoir azimuth anisotropy time difference solving method and a shale reservoir azimuth anisotropy time difference solving device.
Referring to fig. 1, a flowchart illustrating a shale reservoir azimuth anisotropy time difference calculation method according to an embodiment of the present application is shown, and the method includes the following steps:
The first step is to obtain the seismic wave reflection signals of all the wave detection points in each reflection point in each azimuth of each shot point.
The application aims to deeply understand physical attribute differences of underground media, namely rock stratum, reservoir and the like in different directions, analyze azimuth anisotropy characteristics, set a shot point at a preset measuring position on the ground, arrange detection points of a net structure around the shot point, acquire three-dimensional coordinates of the shot point and the detection points, and acquire preset number of azimuth angles of the shot point by taking the shot point as a circle center and taking the north of the shot point as an initial direction and clockwise, wherein in the embodiment, the preset number is 12, namely each azimuth angle of the shot point is 30 degrees. Wherein, the schematic diagram of the arrangement of the shots at one measuring position and the surrounding detection points is shown in fig. 2, the square in fig. 2 represents the shots at the measuring position, the black dots represent the detection points around the shots, and the broken lines represent the division of azimuth angles. The method comprises the steps of obtaining seismic wave reflection signals of all detection points in each azimuth of each shot point at each reflection point, wherein the reflection points are specifically signal reflection points of the shot points in each shale reservoir, it is understood that each reflection point corresponds to a specific geological interface and can be determined through depth time of seismic waves, sampling frequency of each seismic wave reflection signal is 1000Hz, corresponding time length of each seismic wave reflection signal is 4.5 seconds, each seismic wave reflection signal consists of seismic wave reflection signal intensity at each moment, and the detection points on the common edge of two adjacent azimuth angles are divided into azimuth angles which are clockwise and sequentially behind.
The seismic waves are transmitted through the vibration device at the shot point, the seismic waves downwards propagate, and are reflected at different shale reservoirs to propagate to the periphery, and the wave detection points around the shot point receive the seismic wave reflection signals. The seismic response frequencies of different shale intervals are different, the shallow structure may have better response at the high frequency band, and the deeper structure is more prone to the low frequency band. The application aims at the seismic wave reflection signals of different layer segments to select dominant frequency bands, and the dominant frequency bands are selected by considering the type of a seismic source, the sensitivity of the frequency response of a receiving seismic instrument and the signal-to-noise ratio of materials in different frequency bands. The optimal response frequency band of the geological structure of the target stratum is ensured to be acquired, and the information of the dominant frequency band is highlighted through frequency band fusion and enhancement means.
In addition, the seismic data in the dominant frequency band is subjected to earth surface consistency static correction by adopting the dominant frequency band earth surface consistency static correction principle, so that the problem that the single integral static correction cannot be perfectly suitable for the difference of the propagation speeds of seismic waves with different frequencies is solved, the influence of earth surface topography fluctuation, an underground low-speed band and the irregularity of earth surface and underground interfaces on the propagation time of the seismic waves is compensated, and the time for the seismic waves to reach all detectors under the assumption of a flat earth surface becomes consistent. The problem of seismic wave propagation time deviation caused by a complex earth surface structure can be more accurately solved through dominant frequency band earth surface consistency static correction, so that quality and accuracy of subsequent azimuth time difference solving are improved.
The partial schematic diagram of the seismic wave reflection signal of the wave detection point at one reflection point is shown in fig. 3, the abscissa in fig. 3 shows the serial number value of the wave detection point, each curve in fig. 3 shows the seismic wave reflection signal obtained at each wave detection point, and the ordinate shows the receiving time of the seismic wave reflection signal, and it can be observed from fig. 3 that although the seismic wave reflection signal is the seismic wave signal reflected by the same reflection point, the receiving time of the maximum peak data of the seismic wave reflection signal of each wave detection point swings about 200ms, which indicates that the propagation time of the seismic wave is influenced by a medium, and if the medium has anisotropy, the propagation paths of the wave are different, and the time of the wave reaching the receiving point is different.
And the second step is to obtain offset gathers of each reflection point at each azimuth angle of each shot point based on the distribution characteristics of all seismic wave reflection signals corresponding to each shot point.
Since seismic wave data is susceptible to interference from ambient noise, a mesh structure is typically used to arrange the geophones during the measurement process so that seismic wave reflection signals can be acquired by the geophones at different locations. Typically, similar portions of the seismic reflection signal data from different geophones represent the effective signal, while dissimilar portions are typically caused by noise generated during propagation.
In order to remove the noise, the conventional seismic wave data processing method performs comprehensive analysis on the data of a plurality of detection points so as to realize denoising. However, since the distance between each detector point and the seismic reflection point is different, there is a certain time difference between the data, and thus time difference correction is required. To solve this problem, the present embodiment first converts the seismic reflection signal into an offset gather (offset gather), so that processing is performed under a unified framework, ensuring that the time difference is corrected.
The method comprises the steps of receiving signals at different moments, and corresponding to different reflection points at different moments, so that a set formed by the seismic wave reflection signal intensities of all wave detection points of each reflection point at each azimuth angle of a shot point is used as an offset gather, and in the embodiment, the seismic wave reflection signal intensity of an nth wave detection point at the moment corresponding to a kth reflection point in each azimuth angle of the shot point is recorded asAbbreviated as f n(τs,k+τr,n), wherein x and y represent the abscissa and ordinate of a shot at the ground surface, T k represents the depth time of a seismic wave from the shot to a kth reflection point, n represents the number of a detection point in each azimuth of the shot, τ s,k represents the travel time of the seismic wave from the shot to the kth reflection point, and τ r,n represents the travel time of the seismic wave from the kth reflection point to the nth detection point. It should be understood that one reflection point corresponds to one offset gather at each azimuth of the shot point, and that the travel time from the same reflection point to different detection points is different because the paths from the same reflection point to different reflection points are different.
And thirdly, analyzing differences among seismic wave propagation paths corresponding to different elements in the offset gather to obtain propagation path difference values between two elements in the offset gather.
The horizontal offset time from the shot to the geophone is obtained by combining the propagation velocity of the seismic wave based on the distance between the shot and the geophone, and in this embodiment, the horizontal offset time is obtained by dividing the Euclidean distance between the shot and the geophone by the propagation velocity of the seismic wave.
Further, the propagation path may form a triangle from the shot point to the reflection point and then to the detection point to be detected. The travel time of the seismic wave from the shot point to the reflection point is used for describing the path distance from the shot point to the reflection point, the travel time of the seismic wave from the reflection point to the detection point is used for describing the path distance from the reflection point to the detection point, the travel time of the shot point to each offset gather is used as a first right-angle side of a propagation triangle, the horizontal offset time of the shot point to each detection point on each azimuth angle is used as a second right-angle side of the propagation triangle, the travel time of the reflection point to each detection point is used as a hypotenuse of the propagation triangle, each element in each offset gather corresponds to one propagation triangle based on the travel time, and the area difference of the propagation triangles between different elements in each offset gather is used as the propagation path difference value between different elements in each offset gather. In this embodiment, for each offset gather, the absolute value of the difference of the areas between propagation triangles corresponding to different elements is obtained, and linear normalization is adopted to obtain the propagation path difference value between the elements.
It should be understood that the smaller the propagation path difference value, the smaller the numerical difference of the two elements should be. If the difference value of the propagation paths of the two elements is smaller, the numerical value difference of the two elements is more obvious, which means that one element is more influenced by noise.
And a fourth step of obtaining white noise suppression coefficients of the elements based on the variation and the numerical difference of the propagation path difference values between the elements and the rest elements of each offset gather, and comparing the travel time of the seismic waves from the shot point to the reflection point and from the reflection point to each detection point to obtain the distance reliability coefficient of each element in the offset gather corresponding to the reflection point.
The seismic wave is generally affected by noise in different reflection paths, one is white noise generated by the influence of the environment on the whole system, the noise usually appears randomly without fixed occurrence rules, and the other is systematic noise caused by the influence of the rock stratum structure when the seismic wave propagates in the rock stratum, and the farther the propagation distance of the seismic wave is, the larger the systematic noise is.
And for white noise, based on the change and the numerical difference of the propagation path difference value between each element and the rest elements of each offset gather, obtaining the white noise suppression coefficient of each element, namely obtaining the numerical difference value and the propagation path difference value between each element and the rest elements of each offset gather, obtaining the fusion result of the negative correlation mapping result and the numerical difference value of the propagation path difference value, and normalizing the inverse of the accumulation sum of the fusion result of each element and the rest elements of each offset gather to obtain the white noise suppression coefficient of each element of each offset gather. In this embodiment, the fusion result of the second difference negative correlation mapping result and the first difference is specifically the ratio of the first difference to the second difference, and it should be noted that, when the second difference is equal to 0, the denominator is 0, so that the ratio is meaningless, so that a preset value needs to be added to the denominator, the value is 0.01, and the normalization method adopts linear normalization.
It should be appreciated that the greater the difference in values between each element of the offset gather and the remaining elements, the greater the likelihood of being affected by white noise, the greater the difference in propagation path difference values, indicating that the difference between the element values may be due to propagation path differences rather than white noise, and therefore the greater the white noise suppression coefficient of an element, the greater the likelihood that the element will be able to characterize the geologic features of the shale reservoir in which the corresponding reflection point is located.
For system noise which is not white noise, the noise is related to the propagation distance of the reflected seismic waves in the shale reservoir, and the farther the propagation distance is, the larger the influence of the system noise on data is, and the lower the data reliability is. Therefore, the distance reliability coefficient of each element in the offset channel corresponding to the reflection point is obtained by comparing the travel time of the seismic waves from the shot point to the reflection point and from the reflection point to each detection point, namely, the square of the ratio of the first right-angle side to the hypotenuse of the propagation triangle corresponding to each element in the offset channel corresponding to the reflection point is used as the distance reliability coefficient of each element.
It should be understood that the longer the hypotenuse of the element propagation triangle in each offset gather, i.e. the longer the travel time of the seismic wave from the corresponding reflection point to the detection point of the offset gather, the longer the path thereof, the higher the received systematic noise, the lower the corresponding distance reliability coefficient, the lower the data reliability of the corresponding element, and the larger the value of the distance reliability coefficient indicates that the smaller the received systematic noise of the element, and the more can characterize the geological features of the corresponding shale layer position.
And fifthly, fusing the white noise suppression coefficient and the distance reliability coefficient of each element, and weighting the derivative results of each element of all offset gathers based on the fusion results to obtain offset superposition functions of all azimuth angles of the shot point.
The application focuses on the change condition of the geological structure, usually the more severe the amplitude change of the seismic wave, so the application derives the elements in the offset trace, namely f n′(τs,k+τr,n, and is used for representing the change condition of the seismic wave reflection signal intensity of the nth detection point in each azimuth angle of the shot point at the moment corresponding to the kth reflection point.
It should be appreciated that the formula corresponding to the value of the nth element of the kth offset gather of each azimuth of the shot is still f n(τs,k+τr,n), so that for the kth offset gather of each azimuth of the shot, the kth element of the offset superposition function corresponding to the kth offset gather of each azimuth of the shot is calculated, denoted as I (x, y, T k, n), the formula of which is:
Where x and y represent the abscissa and ordinate of the shot at the earth's surface, T k represents the depth time of the seismic wave from the shot to the kth reflection point, N represents the number of the spots in each azimuth of the shot, N represents the total number of the spots in each azimuth of the shot, WF n represents the white noise suppression coefficient of the nth spot in the kth reflection point in each azimuth of the shot, that is, the white noise suppression coefficient of the nth element in the kth offset gather of each azimuth of the shot, DK n represents the distance reliability coefficient of the nth spot in each azimuth of the shot from the kth reflection point, that is, the distance reliability coefficient of the kth offset gather of each azimuth of the shot from the kth reflection point, f n ′(τs,k+τg,n) represents the first derivative value of the intensity of the seismic wave reflected by the nth spot in each azimuth of the kth reflection point, that is, the first derivative value of the kth element in the kth offset gather of the seismic wave in each azimuth of the shot, s,k represents the distance reliability coefficient of the nth spot in each azimuth from the kth reflection point to the kth reflection point, and f n ′(τs,k+τg,n represents the travel time of the seismic wave from the kth reflection point to the kth reflection point.
It should be noted that, since each element in the offset trace set is expressed in the form of a function, the weighted overlap-add can still be expressed in the form of a function, and one shot point corresponds to one offset overlap-add function.
It should be understood that, for each azimuth of the shot point, the first derivative corresponding to the nth element value of the kth offset gather is taken as the characteristic value of the offset gather, and the kth element value of the offset stacking function is obtained by weighting and summing the white noise suppression coefficient and the distance reliability coefficient. The characteristic value of the offset gather is in a derivative form, so that the offset superposition function can directly pay attention to the change condition of the intensity of the seismic wave reflection signal and is used for representing the change degree of the geological structure of a geological layer, wherein the high-noise data in the offset gather is restrained through a white noise restraining coefficient and a distance credibility coefficient, and the reliability of a final result is higher.
And a sixth step of obtaining an imaging channel function of each detector according to the seismic wave signal change characteristics of each detector at all reflection points of each shot point, obtaining a cross-correlation judging function of the offset superposition function of each shot point and the imaging channel function of each detector under different time delays in a preset time period, and obtaining the optimal time shift of each detector of each shot point according to the numerical distribution of the cross-correlation judging function.
The elements of each offset gather can represent the reflection condition of the same reflection point position, the elements are comprehensively calculated, the distribution condition of the seismic wave reflection signal intensity values in all azimuth angles of the shot points obtained at the reflection point position is represented, and the reflection paths from the reflection points to the detection points of the seismic waves are different due to the fact that the positions of the detection points in all azimuth angles are different, and then the influence condition of environmental factors is also different.
The method comprises the steps of taking an offset superposition function of each azimuth angle of a shot point as a standard signal, performing time translation operation on all seismic wave reflection signals of all detection points at corresponding reflection points, performing time error correction on seismic wave reflection data, eliminating interference of environmental noise, obtaining a first derivative vector of an nth seismic wave reflection signal, and recording the first derivative vector as an imaging channel function p n (T), wherein n represents the nth detection point corresponding to the nth seismic wave reflection signal, and T represents depth time of the seismic wave.
And meanwhile, acquiring an offset superposition function of each azimuth of the jth shot point, and recording the offset superposition function as an offset channel function g j (T), wherein T represents the depth time of the seismic wave.
The embodiment adopts a cross-correlation function as a method for acquiring the translation time, and the specific formula is as follows:
Wherein phi j,n (tau) represents the cross-correlation function values of G j (T) and G n,j (T+tau), wherein G n,j(T+τ)＝pn(t+τ)-gj(x+τ);gj (T) represents the offset channel function of each azimuth of the jth shot, G j (T+tau) represents the offset channel function of the jth shot after each azimuth time offset tau, p n (T+tau) represents the imaging channel function of the nth detector after each azimuth time offset tau, X 1、X2 represents the starting moment and the ending moment of the seismic wave from the shot to the depth of the target shale reservoir, and X 1 takes 1 second and X 2 takes 1.8 seconds in the embodiment.
It should be noted that [ X 1,X2 ] is a time window, and an operator can adjust the time window by himself, when determining the time window, a horizon with strong energy, obvious characteristics and high signal to noise ratio of the reflected wave group should be selected, while the width of the time window is related to the signal to noise ratio of the data, and a narrower time window should be selected under the condition that the signal to noise ratio of the target layer is low, so as to improve the accuracy of the related calculation, but too narrow time window may cause artifacts, and the target layer is particularly a page layer corresponding to the reflection position of the reflection point.
It should be understood that, for the two functions G j (T) and G n,j (t+τ), the correlation between the two functions is larger, the larger the function value is, the larger the G n,j (t+τ) is used to describe the difference between the imaging channel function and the shifting channel function shifted by τ units from the time axis of delay, then the cross correlation function is calculated with the original shifting channel function G j (T), the cross correlation between the shifting channel function G j (T) and the imaging channel function p n (T) is described, the autocorrelation of the shifting channel function G j (T) is subtracted, this quantity reflects the correlation between the shifting channel function G j (x) and the imaging channel function p n (T), and the interference effect of the shifting channel function G j (T) is eliminated, the larger the cross correlation function is described that the imaging channel function is shifted by τ units and then the more similar to the shifting channel function is, and the time error of the imaging channel function is described to be closer to τ.
In this embodiment, τ is changed continuously, which is equivalent to shifting the imaging function at different times, and when τ corresponds to the maximum value of the cross correlation function, it is explained that after shifting according to τ, the offset channel function is most similar to the imaging channel function, and then τ corresponds to the optimal time shift of the nth seismic wave reflection signal.
The ideal cross-correlation function should have a relatively pronounced peak, the main peak being significantly larger than its side lobes. The low signal-to-noise ratio data may cause the peak value of the cross-correlation function to be not obvious, and in order to highlight the difference of the correlation function values, the cross-correlation judgment function W j,n (τ) of G j (x) and G n,j (x+τ) is obtained, where the formula is: Finding the independent variable tau corresponding to the maximum function value in the definition domain of the cross-correlation judging function, and corresponding to the optimal time translation of the nth seismic wave reflection signal at each azimuth of the jth shot point.
And seventhly, processing the seismic wave reflection signals of the corresponding demodulation points according to the optimal time translation, and synthesizing the seismic wave reflection signals processed by all demodulation points of each azimuth angle of each shot point to obtain the seismic data characteristics of each azimuth angle of each shot point at the corresponding reflection position of each reflection point.
And obtaining the optimal time translation tau corresponding to the vector of the nth seismic reflection wave according to the cross-correlation judging function, and translating the nth seismic reflection signal by tau unit length along the coordinate axis corresponding to time to obtain the seismic reflection correction signal. The partial schematic diagram of the seismic reflection correction signal of the geophone at one reflection point is shown in fig. 4, the abscissa in fig. 4 represents the serial number value of the geophone, each curve in fig. 4 represents the seismic reflection correction signal of each geophone, and the ordinate is the receiving time of the seismic reflection correction signal.
Finally, intercepting the data of all seismic wave reflection correction signals from 1.0 to 1.8 seconds, and recording the data as a seismic wave data set with an azimuth angle of 0 to 30 degrees and a target layer depth of 1.0 to 1.8 seconds at a j-th shot point coordinate (x, y, z) j, and recording the first two elements, namely (x, y) j, of which the j-th shot point horizontal coordinate is the shot point coordinate (x, y, z) j. In addition, the mean value of the seismic wave reflection correction signals of all the wave detection points acquired by each shot at each moment is used as the seismic data characteristic of each azimuth of each shot at the corresponding reflection position.
It should be understood that the seismic wave reflection correction signal is a signal after time correction, that is, each moment corresponds to a reflection point, so that the signal at each moment can be averaged to uniformly reflect geological reflection data, that is, seismic data characteristics, of the shot point in the shale reservoir where the reflection point is located.
And eighth step, obtaining the anisotropic characteristic of each shot point based on the distribution of the seismic data characteristics of all azimuth angles of each shot point at the corresponding reflection positions of each reflection point.
The seismic wave data set eliminates the environmental influence by a denoising method, improves the reliability of the seismic wave data, and further obtains anisotropic characteristics by the seismic wave data. The anisotropic characteristic comprises an azimuth corresponding to the maximum anisotropic time difference and an anisotropic time difference accumulation result, and specifically comprises the following steps:
The method comprises the steps of obtaining the maximum value and the average value of seismic data characteristics of all azimuth angles of each shot point at the corresponding reflection position of each reflection point, taking the obtained average value as an anisotropic time difference accumulation result of each shot point at the corresponding reflection position of each reflection point, and taking the azimuth angle corresponding to the obtained maximum value as the azimuth corresponding to the maximum anisotropic time difference of each shot point.
Based on the same inventive concept as the method, the embodiment of the application also provides a shale reservoir azimuth anisotropy time difference solving device, which comprises a memory, a processor and a computer program stored in the memory and running on the processor, wherein the processor realizes the steps of any one of the shale reservoir azimuth anisotropy time difference solving methods when executing the computer program.
In summary, the embodiment of the application collects seismic wave data by setting shot points and detection points and converts the seismic wave data into an offset gather, characterizes system noise difference in the seismic wave propagation process by calculating propagation path difference values aiming at the condition that elements in the offset gather are interfered by noise, calculates white noise suppression coefficients by combining data difference among the element values, screens out elements greatly influenced by white noise in the offset gather, calculates distance reliability coefficients by the seismic wave propagation paths corresponding to the elements in the offset gather, screens out elements greatly influenced by systematic noise in the offset gather, and finally performs superposition calculation processing on the offset gather by combining the white noise suppression coefficients and the distance reliability coefficients to obtain an offset gather function serving as a standard function for subsequent time correction. The white noise and systematic noise in the offset channel are respectively eliminated through the white noise suppression coefficient and the distance reliability coefficient, and the accuracy of the standard function is improved.
Furthermore, the azimuth anisotropy time difference information is represented through the offset set function and the cross-correlation function of the seismic wave data, so that the seismic wave data are aligned in time, the description precision of the underground geological structure can be greatly improved, particularly in three-dimensional seismic interpretation, a crack system model which is more in line with the actual situation is formed, the data are formatted uniformly, data noise caused by time delay is eliminated, and a seismic wave data set is obtained.
Compared with the prior art, the method eliminates the time delay interference of noise, utilizes the irregular travel time difference correction of the overlying stratum to carry out time difference correction on the seismic reflection on the target layer, peels off the influence of the stratum near the upper part of the target layer on the azimuth anisotropy time difference of the target layer, can effectively reduce the time error between seismic wave data, and is favorable for carrying out azimuth anisotropy analysis and azimuth inversion on shale reservoirs.
The flowcharts and block diagrams in the figures illustrate the architecture, functionality, and operation of possible implementations of systems, methods and computer program products according to embodiments of the present application. In this regard, each block in the flowchart or block diagrams may represent a module, segment, or portion of code, which comprises one or more executable instructions for implementing the specified logical function(s). In some alternative implementations, the functions noted in the block may occur out of the order noted in the figures. For example, two blocks shown in succession may, in fact, be executed substantially concurrently, or the blocks may sometimes be executed in the reverse order, depending upon the functionality involved. In the description corresponding to the flowcharts and block diagrams in the figures, operations or steps corresponding to different blocks may also occur in different orders than that disclosed in the description, and sometimes no specific order exists between different operations or steps. For example, two consecutive operations or steps may actually be performed substantially in parallel, they may sometimes be performed in reverse order, which may be dependent on the functions involved. Each block of the block diagrams and/or flowchart illustration, and combinations of blocks in the block diagrams and/or flowchart illustration, can be implemented by special purpose hardware-based systems which perform the specified functions or acts, or combinations of special purpose hardware and computer instructions.
It will be evident to those skilled in the art that the application is not limited to the details of the foregoing illustrative embodiments, and that the present application may be embodied in other specific forms without departing from the essential characteristics thereof. Therefore, the above-mentioned embodiments of the present application should be regarded as exemplary and non-limiting in all respects, and modifications of the technical solutions described in the above-mentioned embodiments or equivalent substitutions of some technical features thereof should be included in the protection scope of the present application without making the essence of the corresponding technical solutions deviate from the scope of the technical solutions of the embodiments of the present application.
Claims (14)
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202411956930.4A CN120652538A (en) | 2024-12-29 | 2024-12-29 | Shale reservoir azimuth anisotropy time difference solving method and device |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202411956930.4A CN120652538A (en) | 2024-12-29 | 2024-12-29 | Shale reservoir azimuth anisotropy time difference solving method and device |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN120652538A true CN120652538A (en) | 2025-09-16 |
Family
ID=97008977
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202411956930.4A Pending CN120652538A (en) | 2024-12-29 | 2024-12-29 | Shale reservoir azimuth anisotropy time difference solving method and device |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN120652538A (en) |
-
2024
- 2024-12-29 CN CN202411956930.4A patent/CN120652538A/en active Pending
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Sabbione et al. | Automatic first-breaks picking: New strategies and algorithms | |
| US8948463B2 (en) | Method for analyzing seismic data | |
| WO2007087435A2 (en) | Method for processing acoustic reflections in array data to image near-borehole geological structure | |
| WO2007021857A2 (en) | Method of accounting for wavelet stretch in seismic data | |
| CA2750200A1 (en) | Processing time series data embedded in high noise | |
| EA032186B1 (en) | Seismic adaptive focusing | |
| CN102939546A (en) | System and method for local attribute matching in seismic processing | |
| CN105116448B (en) | A kind of converted wave azimuthal anisotropy bearing calibration and device | |
| US20230333272A1 (en) | Quality control and preconditioning of seismic data | |
| EA001509B1 (en) | A method of processing seismic data signals | |
| US20160054465A1 (en) | Seismic data spectrum restoring and broadening | |
| US20180017692A1 (en) | Device and method for estimating pre-stack wavelet model from seismic gathers | |
| CN111880223B (en) | Attribute analysis method for well-side sound wave remote detection imaging fracture-cave carving | |
| NO20240186A1 (en) | Method of simulating seismic data | |
| US6502038B1 (en) | Method for non-hyperbolic moveout analysis of seismic data | |
| Zhao et al. | Wavelet-crosscorrelation-based interferometric redatuming in 4D seismic | |
| CN111691876B (en) | Method, device and storage medium for imaging adjacent well by using acoustic logging | |
| CN110579798B (en) | Seismic acquisition observation method and system with equal reflection angle intervals | |
| EP3259618B1 (en) | Black hole boundary conditions | |
| CN120652538A (en) | Shale reservoir azimuth anisotropy time difference solving method and device | |
| US20090052279A1 (en) | Method and system for determining azimuth parameters for seismic data processing | |
| Sun et al. | A comparison of coherency measurement using semblance and multiple signal classification, from a seismic-while-drilling perspective | |
| CN116774293B (en) | Method, system, electronic equipment and medium for automatically picking up same phase shaft | |
| JP2019536050A (en) | A method for improving seismic acquisition utilizing an active ultralight seismic detection system | |
| US20220397690A1 (en) | Detection and Removal of Delayed Seismic Travel Times Produced by Velocity Inversions |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination |