Method, medium and equipment for enhancing fracture characterization based on boundary protection filtering

Abstract

The invention relates to the technical field of oil exploration and development, in particular to a method, a medium and equipment for enhancing fracture characterization based on boundary protection filtering. Firstly, signal-to-noise ratio improvement processing is carried out on a three-dimensional seismic data volume to obtain a seismic data volume with a high signal-to-noise ratio; calculating a maximum likelihood data volume based on the seismic data volume to obtain a Likelihood attribute volume for representing fracture; carrying out calibration and three-dimensional space engraving on a Likelihood threshold value representing fracture, and carrying out spatial smoothing filtering on the seismic data volume by using the processed Likelihood attribute volume as a constraint condition to obtain the seismic data volume with clearer fracture boundary; and performing coherent body calculation based on the seismic data volume to obtain a fracture-enhanced coherent data volume. According to the invention, spatial smoothing filtering is carried out on the seismic data volume by using fracture shielding, the signal-to-noise ratio of seismic data is improved, the characteristics of fracture boundaries are protected, the fracture identification precision is improved, and the fracture display effect is enhanced.

Description

Technical Field

[0001] The present invention relates to the technical field of oil exploration and development, and specifically provides a method, medium and device for enhancing fracture characterization based on boundary protection filtering. Background Art

[0002] Fractures are not only channels for oil and gas migration but also containers for oil and gas storage. Fractures are also the main controlling factors for the formation of reservoirs. Therefore, the identification and characterization of fractures play a crucial role in the process of oil and gas exploration. In some areas, there are tall sand dunes on the surface, and the seismic waves are severely absorbed and attenuated. The depth of the Ordovician target layer exceeds 7,000 meters, and the reflection energy of seismic waves is greatly lost, making it difficult to image small deep geological targets. The signal-to-noise ratio of the seismic data inside the Ordovician is low, the reflections are chaotic, the event axes are discontinuous, and the fault points are not clear. Using it for fracture detection has poor effects, the primary and secondary fractures are not distinct, and false fractures are easily formed in seismic coherence detection. Therefore, in order to improve the fracture detection accuracy, it is necessary to perform fracture enhancement processing on the original seismic data.

[0003] In a Chinese invention patent with the authorization announcement number CN 110749924 B and the authorization announcement date December 24, 2021, a fracture zone identification method is disclosed, including step 1: preparing the original seismic data volume; step 2: performing fine preprocessing; step 3: calculating the structure-oriented volume; step 4: calculating the edge detection data volume; step 5: obtaining the fault enhancement data volume; step 6: calculating the seismic gradient structure tensor eigenvalue data volume and obtaining the eigenvalue data volume representing the chaotic reflection characteristics of the fracture zone; step 7: obtaining the fracture zone structure data volume; step 8: performing data volume fusion processing on the fault enhancement data volume and the fracture zone structure data volume to obtain the fracture zone prediction result.

[0004] This solution performs structure-oriented filtering to improve the signal-to-noise ratio and the clarity of fault points of seismic data. On this basis, an edge detection method is used to detect fractures to obtain a data volume representing fractures. This method does not specifically protect the fracture boundary, and the effect of depicting the fracture boundary is weak, not meeting the requirements of coherence detection for basic data.

[0005] Therefore, there is an urgent need for a fracture characterization method that can protect the fracture boundary while improving the clarity of fault points, thereby improving the coherence detection accuracy, enhancing the fracture display effect, and highlighting the primary and secondary relationships of fractures. Summary of the Invention

[0006] To avoid the above problems existing in the prior art, the purpose of the present invention is to provide a method, medium and device for enhancing fracture characterization based on boundary protection filtering.

[0007] To achieve the above purpose, the present invention provides the following technical solution: A method for enhancing fracture characterization based on boundary protection filtering, including the following steps:

[0008] S1: Perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1 to obtain a high signal-to-noise ratio seismic data volume M2;

[0009] S2: Calculate the maximum likelihood data volume based on the seismic data volume M2 to obtain the Likelihood attribute volume C1 for characterizing fractures;

[0010] S3: Calibrate and perform three-dimensional spatial carving on the Likelihood threshold for characterizing fractures to obtain a new Likelihood attribute volume C2;

[0011] S4: Use the Likelihood attribute volume C2 as a constraint condition to perform spatial smoothing filtering on the seismic data volume M2 to obtain a seismic data volume M3 with clearer fracture boundaries;

[0012] S5: Perform coherence volume calculation based on the seismic data volume M3 to obtain a coherence data volume F1 with enhanced fractures.

[0013] The present invention is further configured such that in step S1, dip-guided filtering is used to perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1.

[0014] The present invention is further configured such that step S1 is specifically as follows: First, calculate the instantaneous frequency:

[0015]

[0016] where ω is the instantaneous frequency, which is a function of time t, space x, and space y, and the unit is radians per second (rad / s); represents the instantaneous phase; μ is the three-dimensional seismic data volume, and μ H is the Hilbert transform of μ corresponding to time t; tan -1 represents the arctangent function, and the output range of variation is (-π, π);

[0017] Calculate the instantaneous wave number k x and k y :

[0018]

[0019]

[0020] where k x and k y are the instantaneous wave numbers in the x and y directions respectively, is the instantaneous phase, μ represents the input data, and μ H represents its Hilbert transform corresponding to time t;

[0021] Calculate the instantaneous time dip:

[0022] θ x = k x / ω

[0023] θ y = k y / ω

[0024] where θ x , θ y are the instantaneous time dips in the x and y directions respectively, and k x , k y are the instantaneous wave numbers in the x and y directions respectively;

[0025] Calculate the formation dip and azimuth:

[0026]

[0027]

[0028] where θ is the formation dip and φ is the azimuth;

[0029] Under the constraints of the formation dip and azimuth parameters, perform spatial smoothing filtering on the 3D seismic data volume M1 along the reflection event axis to obtain the 3D seismic data volume M2.

[0030] The present invention is further configured such that step S2 is specifically to set a known sample set as:

[0031] D = {x 1 , x 2 , …, x N}

[0032] where D is the 3D seismic data set;

[0033] The probability density function P(D|σ) is called the likelihood function of σ with respect to D.

[0034]

[0035] If Θ is the value of σ in the parameter space that maximizes the likelihood function P(D|σ), then Θ is the maximum likelihood estimator of σ, which is a function of the sample set and is denoted as:

[0036] Θ = d(x 1 , x 2 , …, x N ) = d(D)

[0037] Θ(x 1 , x 2 , …, x N ) is the maximum likelihood function estimate value;

[0038] Estimate the maximum probability value Θ of the occurrence of fractures in the 3D seismic data volume, which can be used to characterize the fractures and is denoted as the Likelihood attribute, and obtain the Likelihood attribute volume C1.

[0039] The present invention is further configured such that in step S3, according to the logging reservoir information, the logging reservoir interpretation results, and the fracture manual interpretation results, fine calibration of the reservoir well-seismic is carried out; after the Likelihood attribute is calibrated, the threshold M is obtained.

[0040] The present invention is further configured such that the logging reservoir information includes the blowout, mud loss, and hydrocarbon shows obtained from logging.

[0041] The present invention is further configured such that the three-dimensional spatial carving in step S3 is specifically: using the threshold M to perform three-dimensional spatial carving on the Likelihood attribute volume; if the corresponding Likelihood attribute value is less than the threshold M, the Likelihood attribute value is reset to 0; if the corresponding Likelihood attribute value is greater than or equal to the threshold M, the Likelihood attribute value remains unchanged.

[0042] The present invention is further configured such that in step S4, under the control of the formation dip, spatial smoothing filtering is performed on the 3D seismic data volume, and the specific filtering process is as follows:

[0043] Set m seismic data in the inline direction to participate in the calculation, and n seismic data in the xline direction to participate in the calculation, forming a data bin of m * n channels;

[0044] Assume that the initial values of the smoothing bins m and n are both 3, perform the calculation of the specific smoothing processing formula for a sampling point in the 3D seismic data space, and compare and analyze the signal-to-noise ratio of the horizontal continuous reflection event in the profile of the 3D seismic data volume;

[0045] During the process of spatial smoothing filtering, the fracture attribute volume C2 is used as an occlusion condition. After the spatial smoothing filtering process with fracture occlusion, high signal-to-noise ratio seismic data protected by the fracture boundary is obtained, that is, the seismic data volume M3.

[0046] The present invention is further configured such that the specific smoothing processing formula for a sampling point in the 3D seismic data space is as follows:

[0047]

[0048] In the formula, Namely, it is the spatially smoothed seismic data; m is the number of seismic traces in the inline direction, and m is an integer greater than 3; n is the number of seismic traces in the xline direction, and n is an integer greater than 3; i is the serial number in the inline direction of the 3D seismic data, j is the serial number in the xline direction of the 3D seismic data, t is the sampling point serial number in the time direction of the 3D seismic data, u is the 3D seismic data, θ is the 3D data of formation dip angle, is the sampling point serial number in the time direction of the 3D seismic data under the control of formation dip angle.

[0049] The present invention is further configured such that when performing the smoothing processing calculation, the signal-to-noise ratio is controlled within the range of 10 - 20; if the signal-to-noise ratio is less than 10, the values of the bin m and n need to be increased; if the signal-to-noise ratio is greater than 20, the values of the bin m and n are decreased; that is, the signal-to-noise ratio of the horizontal continuous reflection event is used as the discrimination criterion for the values of m and n.

[0050] The present invention is further configured such that the occlusion condition, that is, when C 2 (i, j) = 0 at the sample point serial numbers i or j, the smoothing processing formula is calculated, and when C 2 (i, j) ≠ 0, the calculation of the smoothing processing formula is interrupted.

[0051] The present invention is further configured such that step S5 is specifically as follows: First, extract multi-trace seismic data from the given analysis time window to generate sample point vectors, and form a matrix from these sample point vectors; calculate the eigenvalues and eigenvectors of the covariance matrix of this matrix, then the eigenstructure coherence estimation can be defined as:

[0052]

[0053]

[0054] C 3 is the eigen coherence value; Tr(C) is the trace of the covariance matrix, representing the energy of the covariance matrix; λ max is the maximum eigenvalue of the covariance matrix, representing the dominant energy; λ j is the eigenvalue of the covariance matrix.

[0055] Using the seismic data volume with fracture protection filtering for eigenvalue coherence calculation, a coherence data volume F1 with enhanced fractures can be obtained.

[0056] The present invention is further configured such that the matrix formed by the sample point vectors is as follows:

[0057]

[0058] D is the seismic data volume after spatial smoothing filtering, and d njIt represents the nth sampling point of the jth trace within the analysis time window, where n = 1, 2, 3... N; j = 1, 2, 3... J; J is the total number of seismic data traces within the analysis time window, and N is the total number of sampling points for each trace within the analysis time window. If the seismic trace is a non-zero vector, the covariance of the seismic trace is a symmetric, positive semi-definite matrix with a rank of 1 and a non-zero eigenvalue. Therefore, the covariance matrix of the entire data can be expressed as:

[0059]

[0060] From the above symmetric positive semi-definite covariance matrix, J non-negative eigenvalues can be calculated, and then the eigen coherence value is defined.

[0061] The present invention also includes an electronic device, which includes:

[0062] A memory storing executable instructions;

[0063] A processor that runs the executable instructions in the memory to implement the above method for enhancing fracture characterization based on boundary protection filtering.

[0064] A computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, it implements the above method for enhancing fracture characterization based on boundary protection filtering.

[0065] In summary, the beneficial effects of the above technical solutions of the present invention are as follows:

[0066] 1. The present invention uses fracture occlusion to perform spatial smoothing filtering on the seismic data volume, improving the signal-to-noise ratio of seismic data while protecting the characteristics of fracture boundaries. Compared with the prior art, it improves the accuracy of fracture identification and enhances the fracture display effect. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the following will briefly introduce the drawings required for the description of the embodiments. Obviously, the drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, other drawings can be obtained based on these drawings without creative efforts.

[0068] Figure 1 It is a flowchart of the method for enhancing fracture characterization based on boundary protection filtering according to the present invention.

[0069] Figure 2 It is an original seismic profile.

[0070] Figure 3 It is a seismic profile with improved signal-to-noise ratio using the method described in the embodiments of the present invention.

[0071] Figure 4 It is the seismic profile after fracture boundary protection filtering processing based on the method described in the embodiments of the present invention.

[0072] Figure 5 It is the coherence profile of the original seismic data.

[0073] Figure 6 It is the coherence profile of the fracture-enhanced seismic data by using the method described in the embodiments of the present invention. Specific embodiments

[0074] In order to enable those skilled in the art to better understand the technical solutions of the present invention, the technical solutions of the present invention will be clearly and completely described below in conjunction with the accompanying drawings of the present invention. Based on the embodiments of the present invention, other similar embodiments obtained by those of ordinary skill in the art without creative efforts shall fall within the scope of protection of the present invention.

[0075] In addition, the directional terms mentioned in the following embodiments, such as "up", "down", "left", "right", etc. are only the directions with reference to the accompanying drawings. Therefore, the directional terms used are for illustration rather than limiting the present invention.

[0076] The present invention will be further described below in conjunction with the accompanying drawings and preferred embodiments.

[0077] Embodiment 1:

[0078] As Figure 1 shown, it is a preferred embodiment of the present invention. A method for enhancing fracture characterization based on boundary protection filtering includes the following steps:

[0079] S1: Perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1 to obtain a high signal-to-noise ratio seismic data volume M2.

[0080] The dip-oriented filtering process removes random noise along the direction of the seismic event axis to achieve the purpose of increasing the continuity of the event axis. The key to dip-oriented filtering is to distinguish the dip and azimuth of the reflector while removing the noise superimposed on the reflector. After estimating the dip and azimuth of the reflector, a filter is used to enhance the signal along the reflector. The most important link is the estimation of the dip and azimuth.

[0081] Use dip-oriented filtering to perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1; first calculate the instantaneous frequency:

[0082]

[0083] In the formula, ω is the instantaneous frequency, which is a function of time t, space x, and space y, and the unit is radians per second (rad / s); represents the instantaneous phase; μ is a three-dimensional seismic data volume, and μ H is the Hilbert transform of μ corresponding to the time t; tan -1 represents the arctangent function, and the output range of variation is (-π, π);

[0084] Calculate the instantaneous wave number k x and k y :

[0085]

[0086]

[0087] In the formula, k x , k y are the instantaneous wave numbers in the x and y directions respectively, is the instantaneous phase, μ represents the input data, and μ H represents its Hilbert transform corresponding to the time t;

[0088] Calculate the instantaneous time dip angle:

[0089] θ x = k x / ω

[0090] θ y = k y / ω

[0091] In the formula, θ x , θ y are the instantaneous time dip angles in the x and y directions respectively, and k x , k y are the instantaneous wave numbers in the x and y directions respectively;

[0092] Calculate the formation dip angle and azimuth angle:

[0093]

[0094]

[0095] In the formula, θ is the formation dip angle and φ is the azimuth angle;

[0096] Under the constraints of the formation dip angle and azimuth angle parameters, perform spatial smoothing filtering on the three-dimensional seismic data volume M1 along the reflection event axis to obtain the three-dimensional seismic data volume M2.

[0097] As Figures 2 - 3 shown, Figure 2 is the original seismic data profile, provided by seismic processing and imaging personnel, and this kind of data usually has some random noise.

[0098] Figure 3is the seismic data profile after dip-oriented filtering, and the function of dip-oriented filtering is to eliminate random noise. From Figure 2 , Figure 3 in comparison, the signal-to-noise ratio of the seismic data after dip-oriented filtering is significantly enhanced, providing more reliable data for fault identification.

[0099] S2: Calculate the maximum likelihood data volume based on the seismic data volume M2 to obtain the Likelihood attribute volume C1 for characterizing faults.

[0100] The main purpose of this step is to calculate the coherence volume as the constraint condition for fault boundary protection filtering. The reason for choosing the maximum likelihood attribute is that under the condition of the same seismic signal-to-noise ratio, the signal-to-noise ratio and accuracy of the maximum likelihood attribute are relatively high.

[0101] The maximum likelihood estimation is a statistical method based on the maximum likelihood principle. It provides a method to evaluate model parameters given the observed data, that is: "The model is determined, but the parameters are unknown". Through several experiments, observing the results, if a certain parameter value can make the probability of the sample appearance the largest, it is called the maximum likelihood estimation.

[0102] Set a known sample set as:

[0103] D = {x 1 , x 2 , …, x N}

[0104] In the formula, D is the three-dimensional seismic data set;

[0105] The probability density function P(D|σ) is called the likelihood function of σ with respect to D.

[0106]

[0107] If Θ is the value of σ in the parameter space that can make the likelihood function P(D|σ) the largest, then Θ is the maximum likelihood estimator of σ. It is a function of the sample set and is denoted as:

[0108] Θ = d(x 1 , x 2 , …, x N ) = d(D)

[0109] Θ(x 1 , x 2 , …, x N ) is the estimated value of the maximum likelihood function;

[0110] Estimate the maximum probability value Θ of the occurrence of faults in the three-dimensional seismic data volume, which can be used to characterize faults, denoted as the Likelihood attribute, and obtain the Likelihood attribute volume C1.

[0111] S3: Calibrate the Likelihood threshold representing fractures and perform three-dimensional spatial carving to obtain a new Likelihood attribute volume C2.

[0112] Based on the blowout, mud loss, and hydrocarbon shows obtained from mud logging, the three types of reservoirs interpreted from well logging curves, and the fractures interpreted manually; carry out fine reservoir calibration between well and seismic data by integrating this information. After calibrating the Likelihood attribute, obtain the threshold M. Use the threshold M to perform three-dimensional spatial carving on the Likelihood attribute volume; if the corresponding Likelihood attribute value is less than the threshold M, reset the Likelihood attribute value to 0; if the corresponding Likelihood attribute value is greater than or equal to the threshold M, the Likelihood attribute value remains unchanged. After well-seismic calibration and three-dimensional carving, obtain a new Likelihood attribute volume C2.

[0113] S4: Use the Likelihood attribute volume C2 as a constraint condition to perform spatial smoothing filtering on the seismic data volume M2 to obtain a seismic data volume M3 with clearer fracture boundaries.

[0114] In this step, the Likelihood attribute volume is used to represent the fracture boundary. Under the constraint of the formation dip angle θ parameter, the seismic data volume is subjected to spatial smoothing filtering along the reflection event axis direction. The fracture acts as an occlusion condition during the spatial smoothing filtering process, that is, the spatial smoothing filtering is interrupted when encountering a fracture, and the samples at the breakpoints do not participate in the spatial smoothing filtering, thus protecting the characteristics of the fracture boundary, and the signal-to-noise ratio of the seismic data far from the fracture is significantly improved.

[0115] Perform spatial smoothing filtering on the three-dimensional seismic data volume under the control of the formation dip angle. The specific filtering process is as follows:

[0116] Set m seismic data traces in the inline direction and n seismic data traces in the xline direction to participate in the calculation, forming a data bin of m * n traces. The specific smoothing processing formula for a sampling point in the three-dimensional seismic data space is as follows:

[0117]

[0118] In the formula, is the seismic data after spatial smoothing; m is the number of seismic traces in the inline direction, m is an integer greater than 3; n is the number of seismic traces in the xline direction, n is an integer greater than 3; i is the serial number of the inline direction of the three-dimensional seismic data, j is the serial number of the xline direction of the three-dimensional seismic data, t is the sampling point serial number in the time direction of the three-dimensional seismic data, u is the three-dimensional seismic data, θ is the three-dimensional data of the formation dip angle, is the sampling point serial number in the time direction of the three-dimensional seismic data under the control of the formation dip angle.

[0119] Assume that the initial values of the smooth bins m and n are both 3, and complete the calculation of, and compare and analyze the three-dimensional seismic data volume whether the signal-to-noise ratio of the horizontal continuous reflection event in the middle section has increased. The signal-to-noise ratio is appropriately controlled between 10 and 20. If the signal-to-noise ratio is less than 10, the values of bins m and n need to be appropriately increased. If the signal-to-noise ratio is greater than 20, the values of bins m and n can be appropriately decreased. That is, the signal-to-noise ratio of the horizontal continuous reflection event is used as the discrimination criterion for the values of m and n.

[0120] During the spatial smoothing filtering process, the fracture attribute volume C2 is used as an occlusion condition, that is, when C 2 (i, j) = 0 at the sample point serial numbers i or j, calculate C 2 when C(i, j) ≠ 0, the calculation is interrupted

[0121] After the spatial smoothing filtering process with fracture occlusion, high signal-to-noise ratio seismic data with fracture boundary protection is obtained that is, the seismic data volume M3.

[0122] As Figure 4 shown, it is the seismic data section after fracture boundary protection filtering. This step mainly uses the Likelihood attribute characterizing fractures for fracture occlusion and performs spatial smoothing filtering on the seismic data. The purpose is to improve the signal-to-noise ratio of seismic data while protecting the clarity of fracture boundaries. Compare Figure 2 and Figure 3 , Figure 4 the signal-to-noise ratio of the section is higher, the breakpoints are clearer, and the longitudinal fracture extension characteristics are more obvious.

[0123] S5: Perform coherence cube calculation based on the seismic data volume M3 to obtain the fracture-enhanced coherence data volume F1.

[0124] Perform coherence cube calculation on the spatially smoothed seismic data volume M3. This step uses the eigenstructure-based coherence cube technique. This method is obtained by calculating the eigenvalues of the seismic data volume. First, extract multi-trace seismic data from the given analysis time window to generate sample point vectors, and form a matrix from these sample point vectors:

[0125]

[0126] D is the spatially smoothed seismic data volume, d njrepresents the nth sampling point of the jth trace within the analysis time window, where n = 1, 2, 3... N; j = 1, 2, 3... J; J is the total number of seismic traces within the analysis time window, and N is the total number of sampling points for each trace within the analysis time window. If a seismic trace is a non-zero vector, the covariance of the seismic trace is a symmetric, positive semi-definite matrix of rank 1 with a non-zero eigenvalue. Therefore, the covariance matrix of the entire data can be expressed as:

[0127]

[0128] From the above symmetric positive semi-definite covariance matrix, J non-negative eigenvalues can be calculated. The magnitudes of the eigenvalues can quantitatively describe the degree of variation of the data volume. Then, the eigen-coherence value of the coherence cube can be defined as:

[0129]

[0130]

[0131] C 3 is the eigen-coherence value; Tr(C) is the trace of the covariance matrix, representing the energy of the covariance matrix; λ max is the largest eigenvalue of the covariance matrix, representing the dominant energy; λ j is the eigenvalue of the covariance matrix.

[0132] The coherence cube algorithm based on eigen-structure analysis applies the idea of principal component analysis. In principal component analysis, the first principal component always coincides with the major axis of the divergence ellipsoid through the solid angle formed by the eigenvectors of the covariance matrix. The principal axes of the ellipsoid are determined by the eigenvectors of the matrix, and the eigenvalues are equal to the lengths of the semi-axes of the ellipsoid. The largest eigenvalue reflects the common part of the original observation point information, that is, coherence, and C3 is the eigenvalue coherence.

[0133] By performing eigen-coherence calculation on the seismic data volume using fracture protection filtering, a fracture-enhanced coherence data volume F1 can be obtained.

[0134] As Figure 5 、 Figure 6 shown, Figure 5 is the coherence profile of the original seismic data, Figure 6 is the coherence profile after fracture enhancement. Comparing Figure 5 and Figure 6 , the signal-to-noise ratio of the coherence profile after fracture enhancement is higher, the characteristics of the main fractures are more prominent, and the fractures are more continuous longitudinally. From the comparison of the profiles, it can be seen that the spatial smoothing filtering for fracture occlusion has a huge improvement on the coherence calculation accuracy, and the fracture enhancement effect is obvious.

[0135] Example 2:

[0136] An electronic device, the electronic device comprising:

[0137] A memory storing executable instructions;

[0138] A processor that runs the executable instructions in the memory to implement the above-mentioned method for enhancing fracture characterization based on boundary protection filtering.

[0139] Embodiment 3:

[0140] A computer-readable storage medium storing a computer program, which when executed by a processor implements the above-mentioned method for enhancing fracture characterization based on boundary protection filtering.

[0141] In summary, the present invention uses fracture occlusion to perform spatial smoothing filtering on seismic data volumes, protecting the characteristics of fracture boundaries while improving the signal-to-noise ratio of seismic data, and improving the accuracy of fracture identification and enhancing the fracture display effect compared with the prior art.

[0142] The above are only the preferred embodiments of the present invention, and the protection scope of the present invention is not limited to the above embodiments. All technical solutions falling within the idea of the present invention belong to the protection scope of the present invention. It should be noted that for those of ordinary skill in the art, several improvements and refinements made without departing from the principle of the present invention should also be regarded as within the protection scope of the present invention.

Claims

1. A method for enhancing fracture characterization based on boundary protection filtering, characterized in that, it includes the following steps: S1: Perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1 to obtain a high signal-to-noise ratio seismic data volume M2; S2: Calculate the maximum likelihood data volume based on the seismic data volume M2 to obtain the Likelihood attribute volume C1 for characterizing fractures; S3: Calibrate the Likelihood threshold for characterizing fractures and perform three-dimensional spatial carving to obtain a new Likelihood attribute volume C2; S4: Use the Likelihood attribute volume C2 as a constraint condition to perform spatial smoothing filtering on the seismic data volume M2 to obtain a seismic data volume M3 with clearer fracture boundaries; S5: Perform coherence volume calculation based on the seismic data volume M3 to obtain a fracture-enhanced coherence data volume F1.

2. The method for enhancing fracture characterization based on boundary protection filtering according to claim 1, characterized in that, in step S1, dip-oriented filtering is used to perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1.

3. The method for enhancing fracture characterization based on boundary protection filtering according to claim 2, characterized in that, using dip-oriented filtering to perform signal-to-noise ratio improvement processing on the three-dimensional seismic data volume M1 specifically means that, under the constraints of formation dip and azimuth parameters, perform spatial smoothing filtering on the three-dimensional seismic data volume M1 along the reflection event axis to obtain the three-dimensional seismic data volume M2; Formation dip and azimuth: θ x = k x / ω θ y = k y / ω where θ x and θ y are the instantaneous time dips in the x and y directions respectively, k x and k y are the instantaneous wave numbers in the x and y directions respectively, ω is the instantaneous frequency; θ is the formation dip, and φ is the azimuth angle.

4. The method for enhancing fracture characterization based on boundary protection filtering according to claim 1, characterized in that, step S2 is specifically as follows: Set a known sample set as: D = {x 1 , x 2 , …, x N} where D is a three-dimensional seismic data set; The probability density function P(D|σ) is called the likelihood function of σ with respect to D; If Θ is the value of σ in the parameter space that can maximize the likelihood function P(D|σ), then Θ is the maximum likelihood estimator of σ, which is a function of the sample set and is denoted as: Θ = d(x 1 , x 2 , …, x N ) = d(x) Θ(x 1 ,x 2 ,…,x N ) is the maximum likelihood function estimate value; Estimate the maximum probability value Θ of the occurrence of fractures in the three-dimensional seismic data volume, which can be used to characterize fractures and is denoted as the Likelihood attribute, to obtain the Likelihood attribute volume C1.

5. The method for enhancing fracture characterization based on boundary protection filtering according to claim 1, characterized in that, in step S3, according to the logging reservoir information, logging reservoir interpretation results, and fracture manual interpretation results, carry out fine calibration of reservoir well-seismic; after calibrating the Likelihood attribute, obtain the threshold M.

6. The method for enhancing fracture characterization based on boundary protection filtering according to claim 5, characterized in that, the logging reservoir information includes blowout, mud loss, and oil and gas shows obtained from logging.

7. The method for enhancing fracture characterization based on boundary protection filtering according to claim 6, characterized in that, The three-dimensional spatial carving described in step S3 specifically involves performing three-dimensional spatial carving on the Likelihood attribute volume using the threshold M; if the corresponding Likelihood attribute value is less than the threshold M, the Likelihood attribute value is reset to 0; if the corresponding Likelihood attribute value is greater than or equal to the threshold M, the Likelihood attribute value remains unchanged.

8. A method for enhancing fracture characterization based on boundary protection filtering according to claim 1, characterized in that, in step S4, spatial smoothing filtering is performed on the three-dimensional seismic data volume under the control of the formation dip angle, and the specific filtering process is as follows: Set m seismic data traces in the inline direction to participate in the calculation, and n seismic data traces in the xline direction to participate in the calculation, forming a data bin of m*n traces; Assume that the initial values of the smoothing bins m and n are both 3, perform specific smoothing processing formula calculations for a sampling point in the three-dimensional seismic data space, and compare and analyze the signal-to-noise ratio of the horizontal continuous reflection event in the section of the three-dimensional seismic data volume; During the process of spatial smoothing filtering, the fracture attribute volume C2 is used as an occlusion condition. After the spatial smoothing filtering process with fracture occlusion, high signal-to-noise ratio seismic data with fracture boundary protection is obtained, that is, the seismic data volume M3.

9. A method for enhancing fracture characterization based on boundary protection filtering according to claim 8, characterized in that, the specific smoothing processing formula for a sampling point in the three-dimensional seismic data space is as follows: In the formula, is the seismic data after spatial smoothing; m is the number of seismic traces in the inline direction, and m is an integer greater than 3; n is the number of seismic traces in the xline direction, and n is an integer greater than 3; i is the serial number in the inline direction of the 3D seismic data, j is the serial number in the xline direction of the 3D seismic data, t is the sampling point serial number in the time direction of the 3D seismic data, u is the 3D seismic data, and θ is the 3D data of the formation dip angle, is the sampling point serial number in the time direction of the 3D seismic data under the control of the formation dip angle.

10. A method for enhancing fracture characterization based on boundary protection filtering according to claim 9, characterized in that, When performing smoothing processing calculations, the signal-to-noise ratio is controlled within the range of 10 - 20; if the signal-to-noise ratio is less than 10, the values of bins m and n need to be increased; if the signal-to-noise ratio is greater than 20, the values of bins m and n are decreased; that is, the signal-to-noise ratio of the horizontal continuous reflection event is used as the discrimination criterion for the values of m and n.

11. A method for enhancing fracture characterization based on boundary protection filtering according to claim 10, characterized in that, The occlusion condition, that is, C at sample point serial number i or j 2 When C 2 (i,j) = 0, perform the smoothing processing formula calculation. 2 When C 2 (i,j) ≠ 0, interrupt the smoothing processing formula calculation.

12. A method for enhancing fracture characterization based on boundary protection filtering according to claim 1, characterized in that, Step S5 is specifically as follows: First, extract multi-trace seismic data from the given analysis time window to generate sample point vectors, and form a matrix from these sample point vectors; calculate the eigenvalues and eigenvectors of the covariance matrix of this matrix, then the eigenstructure coherence estimation is defined as: C 3 is the intrinsic coherence value; Tr(C) is the trace of the covariance matrix, representing the energy of the covariance matrix; λ max is the maximum eigenvalue of the covariance matrix, representing the dominant energy; λ j is the eigenvalue of the covariance matrix; Perform eigenvalue coherence calculation on the seismic data volume using fracture protection filtering to obtain the fracture-enhanced coherence data volume F1.

13. A method for enhancing fracture characterization based on boundary protection filtering according to claim 12, characterized in that, The matrix formed by the sample point vectors is as follows: D is the seismic data volume after spatial smoothing filtering, d nj represents the nth sampling point of the jth trace within the analysis time window, where n = 1, 2, 3... N; j = 1, 2, 3... J; J is the total number of seismic traces within the analysis time window, and N is the total number of sampling points per trace within the analysis time window; if the seismic trace is a non-zero vector, then the covariance of the seismic trace is a symmetric, positive semi-definite matrix of rank one with a non-zero eigenvalue, so the covariance matrix of the entire data can be expressed as: J non-negative eigenvalues can be calculated from the above symmetric positive semi-definite covariance matrix, and then the eigen-coherence value described in claim 12 is defined.

14. An electronic device, characterized in that, the electronic device includes: a memory storing executable instructions; A processor that runs the executable instructions in the memory to implement a method for enhancing fracture characterization based on boundary protection filtering as recited in any one of claims 1-13.

15. A computer-readable storage medium, characterized in that the computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, it implements a method for enhancing fracture characterization based on boundary protection filtering as recited in any one of claims 1-13.

Images