A method for identifying bioherms

Abstract

The present invention provides a method for identifying bioherms, which includes the following steps: establishing a velocity model; using the velocity model to convert a CRP gather in the time domain into a CRP gather in the depth domain; converting the offset gather of the CRP gather in the depth domain into an angle gather; using the amplitude stacking imaging of the angle gather to determine the target layer segment where the bioherms develop; calculating the seismic morphological index attribute data body of the target layer segment, thereby identifying the bioherms. The bioherm identification method provided by the present invention has low dependence on well drilling information, can improve the accuracy of identification and description of subsalt bioherms, is more suitable for processing and imaging of seismic geometric convexity parameter attributes in large work areas, can provide high-quality drilling targets, reduce drilling risks, and has great application value.

Description

Technical Field

[0001] The invention relates to the field of petroleum exploration, and in particular to a method for identifying bioherms. Background Art

[0002] Bioherms are independent carbonate sediments that have a certain number of in situ reef-building biological frameworks, are able to withstand strong winds and waves, and are often raised in the terrain. They have carbonate rock geometry with (paleo) geomorphic uplifts, and their reservoir space types, diagenetic changes, and porosity and permeability distribution characteristics are quite different from other carbonate reservoirs. Bioherms occupy a very important position in carbonate oil and gas fields due to their good reservoir properties. They contain huge amounts of oil and natural gas resources, and therefore have attracted widespread attention. Seismic identification and description technologies for bioherm reservoirs have developed rapidly at home and abroad.

[0003] At present, the industry mainly studies the evolution of sea level rise and fall and sedimentary cycles based on the theory of sedimentology and sequence stratigraphy, and explores how to use seismic technology to identify and finely describe carbonate bioherm reservoirs based on the analysis of prior well characteristics. For example, patent (CN106353806B) mainly uses the spectrum width and main frequency value of the target area to detect the connectivity of bioherm reservoirs; patent (CN107831554A) mainly uses the input of lithology combination data, logging data, and seismic data of various parts of the bioherm to establish a complex bioherm microfacies identification and determination method; patent (CN104977611B) mainly establishes a set of microfacies identification patterns including logging curve characteristics and seismic reflection structure characteristics, so as to determine the development interval of the bioherm reservoir, and then characterize the internal structure of the bioherm reservoir, thereby realizing a method for carving the bioherm reservoir. Patent (CN103488971B) reports an effective method for geometric morphology identification of bioherm reservoirs. Patent (CN108931813A) reported a method for identifying reef caps, reef cores, and reef bases using seismic information and well logging information. Zhao Jian et al. (Contemporary Chemical Industry 2016, 45(4)) proposed gravity and magnetic, well logging, seismic identification and other methods as common techniques for identifying reefs based on the development characteristics of reefs; Wu Yong et al. (Petroleum Exploration and Development 2017, 44(6)) established a platform-margin reef boundary identification method based on the slope theory in geomorphology and combined with paleo-geomorphological analysis methods to use seismic structural fine interpretation stratigraphic data for high-precision slope attribute calculation. Wang Jingchuan et al. (Computational Technology for Geophysical and Geochemical Exploration 2008, 30(4)) mainly used the seismic reflection characteristics of reefs, conventional seismic amplitude attribute analysis, coherence body, and reservoir spectrum imaging technology to identify the distribution characteristics of reefs.

[0004] After investigation, no relevant research results were found on the identification of bioherms by using the correlation between the developmental morphological characteristics of bioherms and the seismic geometric convexity characteristics. At present, the technology used in the research of carbonate bioherm reservoir identification at home and abroad is mainly based on seismic data with relatively good signal-to-noise ratio and resolution. The research on seismic identification and prediction of subsalt bioherms that are blocked by thick rock salt and have relatively poor seismic data quality is very weak. The main problems are as follows: (1) Due to the influence of the absorption and attenuation of seismic reflection waves by the overlying thick rock salt, the seismic data has weak energy and low frequency, so the identification ability of subsalt bioherms is relatively poor; (2) Subsalt bioherms have special morphology and large thickness variations in the vertical and horizontal directions. At present, effective technical means for accurately describing the characteristics of subsalt bioherms are still very weak; (3) The technical methods for the identification and description of bioherms are usually based on time domain seismic data. The influence of salt body velocity causes a large difference between the time domain imaging and depth domain imaging characteristics, which affects the identification and description accuracy of subsalt bioherms containing oil and gas. Summary of the invention

[0005] In view of the problems of difficulty in identifying subsalt reef morphology and low description accuracy in the prior art, the present invention proposes a subsalt reef identification method that enhances the characteristics of seismic geometric uplift parameters. First, based on velocity field modeling, the pre-stack data set in the time domain is converted into a data set in the depth domain, and then the CRP data set is optimized to remove wild values, and denoising and residual time-varying static correction are performed, and then the subsalt offset dominant amplitude stacking imaging is performed. On the basis of dominant offset amplitude imaging, the second-order spatial differential derivative method in the depth domain is used to obtain the geometric convexity parameter data body of the target seismic imaging data body. Using the correlation between the uplift sedimentary body morphology of carbonate reefs and the seismic geometric uplift morphological attributes, the subsalt reef characteristics in the depth domain are identified and described, providing the necessary conditions for optimizing the reef drilling targets.

[0006] A first aspect of the present invention provides a method for identifying a bioherm, comprising the following steps:

[0007] S1: Establish velocity model;

[0008] S2: using the velocity model to convert the CRP gathers in the time domain into CRP gathers in the depth domain;

[0009] S3: converting the offset gather of the CRP gather in the depth domain into an angle gather;

[0010] S4: using the angle gather amplitude stack imaging to determine the target layer section for reef development;

[0011] S5: Calculate the seismic morphological index attribute data volume of the target layer segment to identify the bioherm.

[0012] According to some embodiments of the present invention, S2 further includes using a dominant offset range of the target layer segment to remove amplitude wild values ​​of the offset and NMO distortion amplitudes of far offsets to obtain CRP offset gathers.

[0013] According to some embodiments of the present invention, S2 performs denoising and time-varying static correction on the CRP offset gather.

[0014] According to some embodiments of the present invention, S4 further includes performing smoothing filtering on the imaged seismic trace data.

[0015] According to some embodiments of the present invention, in S5, the seismic morphological index attribute data volume of the target layer segment is calculated by using the spatial second-order difference derivative method in the depth domain.

[0016] According to some embodiments of the present invention, the seismic morphological index attribute data body includes a convexity characteristic parameter K pos , K neg , morphological index S shape .

[0017] According to some embodiments of the present invention, in S5, the method for calculating the seismic morphological index attribute data volume of the target layer segment includes the following steps:

[0018] Step A: Assume that the curvature expression of the two-dimensional function Z(x) is

[0019]

[0020] Step B: Assume that the quadratic surface function expression (6) of a point Z in space is as follows:

[0021] Z(x, y) = ax 2 +by 2 +cxy+dx+ey+f (6)

[0022] Step C: Calculate the curvature of a point on the surface expressed by equation (6), and use the values ​​of the surrounding grid points to perform least squares fitting on the local quadratic surface. The first three terms in equation (6) represent quadratic ellipse fitting, and the last three terms represent plane fitting. Use 3×3 grid units for approximation. According to the definition of first-order and second-order derivatives, the coefficients in equation (6) are obtained as equation (7).

[0023]

[0024] Where: z1,…,z9 are the values ​​of the surface mesh nodes 1…9, Δx, Δy are the distances between the mesh nodes;

[0025] Step D: Calculate the convexity characteristic parameter K of the surface pos , K neg , morphological index Sshape , the calculation formulas are as follows (8), (9), (10):

[0026]

[0027]

[0028]

[0029] Step E: Based on the seismic data volume, calculate the convexity characteristic parameter K of the data volume pos , K neg , morphological index S shape Attribute data body.

[0030] According to some embodiments of the present invention, in S1, the root mean square velocity body is converted into a layer velocity body using Dix formula (1), and the target formation test velocity is used for correction;

[0031]

[0032] In the formula, v n is the interval velocity of the nth layer, τ n and τ n-1 V is the zero-offset two-way travel time of the corresponding formation seismic wave, n and V n-1 is the root mean square velocity of the seismic wave in the corresponding formation.

[0033] According to some embodiments of the present invention, it is characterized in that

[0034] In S2, the time-depth conversion formula (2) for converting the CRP gathers in the time domain into the CRP gathers in the depth domain is:

[0035]

[0036] where Z is the thickness of the reflecting layer, is the average velocity of the overlying strata, τ 0 is the double-layer reflection time of the formation, z n is the thickness of the nth layer, Vn is the interval velocity of the nth layer, τ n is the double-layer reflection time of the nth layer.

[0037] According to some embodiments of the present invention, in S3, the offset gather of the CRP gather in the depth domain is converted into an angle gather, and the formula (3) for calculating the reflection angle θ of the formation is:

[0038]

[0039] Among them, X represents the size of the offset, and Z is the thickness of the underground reflective layer.

[0040] According to some embodiments of the present invention, in S4, in the preferred amplitude stack imaging using the angle gather, the stacking formula is (4)

[0041] A i =∑a ij (4)

[0042] Among them, A i represents the i-th sample point of a seismic record, a ij Represents the amplitude value of the jth earthquake record on the CRP gather.

[0043] According to some embodiments of the present invention, the angle gather is a dominant angle gather. In the present invention, the dominant angle gather refers to an effective prestack amplitude gather after removing the near angles rich in interlayer multiple waves and amplitude anomalies and the far angles of NMO distortion.

[0044] According to some embodiments of the present invention, in S5, based on the seismic morphological index attribute data body and in combination with the convexity characteristics of carbonate bioherms, bioherm targets are identified and drilling well locations are provided.

[0045] According to a preferred embodiment of the present invention, the bioherm is a carbonate bioherm.

[0046] According to a preferred embodiment of the present invention, the bioherm is a subsalt bioherm.

[0047] The second aspect of the present invention provides an application of the method according to the first aspect in identifying bioherms.

[0048] According to a preferred embodiment of the present invention, the bioherm is a carbonate bioherm.

[0049] According to a preferred embodiment of the present invention, the bioherm is a subsalt bioherm.

[0050] According to some specific implementations of the present invention, the technical solution for deep-domain subsalt carbonate reef identification mainly includes the following 11 parts:

[0051] (1) Establishing the velocity model of the work area: using Dix formula (1), the root mean square velocity volume is converted into a layer velocity volume, and the target formation test velocity obtained by drilling is used for correction;

[0052] The Dix conversion formula is:

[0053]

[0054] In the formula, v n is the interval velocity of the nth layer. τ n and τ n-1V is the zero-offset two-way travel time of the formation seismic wave, n and V n-1 is the root mean square velocity of the seismic wave in the corresponding formation.

[0055] (2) Using the established velocity model, the CRP gathers in the time domain are converted into CRP gathers in the depth domain. The time-depth conversion formula (2) is:

[0056]

[0057] where Z is the thickness of the reflecting layer, The formation is the average velocity, τ 0 is the double-layer reflection time of the formation, z n is the thickness of the nth sublayer, Vn is the layer velocity of the nth sublayer, τ n is the double-layer reflection time of the nth layer.

[0058] (3) Using the advantageous offset range of the target layer segment, the amplitude outliers of the offset distance and the amplitude of the NMO distortion at the far offset distance are removed.

[0059] (4) De-noising and time-varying static correction of CRP offset gathers.

[0060] (5) The offset gather of the CRP gather is converted into an angle gather, and the reflection angle θ of the formation is calculated by formula (3):

[0061]

[0062] Among them, X represents the size of the offset, and Z is the thickness of the underground reflective layer.

[0063] (6) Amplitude stacking imaging of dominant angle gathers, the stacking formula is (4)

[0064] A i =∑a ij (4)

[0065] Among them, A i represents the i-th sample point of a seismic record, a ij Represents the amplitude value of the jth earthquake record on the CRP gather.

[0066] (7) Smoothing and filtering the imaged seismic trace data.

[0067] (8) Determine the target interval for the development of subsalt reefs.

[0068] (9) The second-order spatial difference method in the depth domain is implemented for the target body to calculate the convexity characteristic parameters of the function curvature. We assume that the curvature expression of the two-dimensional function Z(x) is

[0069]

[0070] Let the quadratic surface function expression (6) of a point Z in space be as follows:

[0071] Z(x, y) = ax 2 +by 2 +cxy+dx+ey+f (6)

[0072] Here, the second-order difference method is mainly used to calculate the derivative. To calculate the curvature of a point on the surface expressed by equation (6), the values ​​of the surrounding grid points are used to fit the local quadratic surface using the least squares method. The first three terms in equation (6) represent the quadratic ellipse fitting, and the last three terms represent the plane fitting. We use a 3×3 grid unit for approximation ( Figure 2 ), according to the definition of first-order and second-order derivatives, the coefficients in equation (6) are obtained as expression (7).

[0073]

[0074] Where: z1,…,z9 are Figure 1 The values ​​of the mesh nodes 1…9 at the level shown in the figure, Δx and Δy are the distances between the mesh nodes, and the convexity characteristic parameter K of the surface is calculated on this basis pos , K neg , morphological index S shape , the calculation formulas are as follows (8), (9), (10):

[0075]

[0076]

[0077]

[0078] (10) Based on the seismic data volume, calculate the convexity characteristic parameter K of the data volume for the target layer segment pos , K neg , morphological index S shape Attribute data body.

[0079] (11) Based on the calculation of K of the data body pos , K neg Convexity characteristic parameter, morphology index S shape Characteristics, combined with the convexity characteristics of carbonate bioherms, select subsalt carbonate bioherms targets and provide drilling well locations.

[0080] In view of the low accuracy and difficulty in identifying the above-mentioned subsalt reef morphology description, the present invention proposes a subsalt reef identification method that enhances the characteristics of seismic geometric convexity parameters. This technical method can improve the accuracy of subsalt reef identification description, is more suitable for the processing imaging effect of seismic geometric convexity parameter attributes, and provides an important basis for providing high-quality targets. In the reef identification application of a certain oil field in Brazil, the morphology and targets of subsalt reefs were identified by using the characteristics of seismic convexity parameters. Subsequent drilling confirmed that the effective thickness of the biological reefs reached more than 300 meters, which has great commercial value. It also proves the effectiveness and feasibility of the technology and has great application value. By using the correlation between the convex sedimentary body morphology of carbonate reefs and the seismic geometric convexity morphological attributes, the subsalt reef reservoir is identified and described in the depth domain, and high-quality reef drilling targets are selected, which can reduce the risk of subsalt reef drilling. This subsalt reef identification method has low dependence on well drilling information, avoids the problem of large differences between time domain imaging and depth domain imaging caused by the influence of salt body velocity, and affects the accuracy of seismic attribute analysis of subsalt reservoirs. This technology can improve the accuracy of identification and description of subsalt reefs, is more suitable for processing and imaging of seismic geometric convexity parameter attributes in large work areas, can provide high-quality drilling targets, reduce drilling risks, and has great application value. BRIEF DESCRIPTION OF THE DRAWINGS

[0081] Figure 1 It is a technical flow chart of the bioherm identification method according to the present invention.

[0082] Figure 2 Schematic diagram of a grid unit for second-order difference according to the present invention.

[0083] Figure 3 The velocity model established according to Embodiment 1 of the present invention is shown.

[0084] Figure 4 The CRP gather data in the depth domain according to Example 1 of the present invention is shown.

[0085] Figure 5 The CRP gather data in the depth domain after de-wiring, denoising and time-varying static correction according to Example 1 of the present invention are shown.

[0086] Figure 6 The converted depth-domain CRP angle gather according to Example 1 of the present invention is shown.

[0087] Figure 7 It is a comparison diagram of subsalt depth domain and time domain imaging according to Example 1 of the present invention, wherein (A) is subsalt seismic dominant offset imaging in the depth domain, and (B) is subsalt seismic dominant offset imaging in the time domain.

[0088] Figure 8The depth domain subsalt dominant offset bioherm imaging according to Example 1 of the present invention is shown.

[0089] Fig. 9 The seismic convexity attribute parameter K according to Example 1 of the present invention is shown. pos Identification profile of bioherms. DETAILED DESCRIPTION

[0090] The present invention is described in detail by the following examples, which can enable those skilled in the art to more fully understand the present invention, but these examples do not constitute any limitation to the scope of the present invention.

[0091] The present invention takes a subsalt oil field in the Santos Basin of Brazil as an example. Figure 1 The technical process was used to establish the velocity field model of the work area, and the CRP gathers in the time domain were converted into CRP gathers in the depth domain. After the pre-stack CRP gathers were de-fielded, de-noised, and subjected to time-varying static correction, the CRP offset gathers were converted into angle gathers to complete the imaging of the dominant offset angle gathers. The second-order spatial differential derivative method in the depth domain was used to obtain the geometric morphological attribute parameter data body of the target seismic imaging data body, and the correlation between the convex sedimentary body morphology of carbonate reefs and the seismic geometric convexity attribute parameters was used to identify the subsalt reef development zone in the depth domain, and select high-quality reef drilling targets. These microbial reefs, which are mainly mound-like or imbricate in morphology, are difficult to identify in seismic profiles, but can be clearly imaged and effectively identified in seismic geometric convexity parameter profiles, see. Figure 8 and Fig. 9 ,in, Figure 8 is in accordance with Figure 1 The deep-domain subsalt carbonate reef imaging is completed by the technical process, and the reef morphology is not easy to identify; Fig. 9 is the seismic convexity attribute parameter K of a two-dimensional seismic profile in the study area after processing at the top interface of BVE100. pos The bioherm identification profile shows the top hillock reflection feature and the internal weak or even chaotic reflection feature, which are typical bioherm seismic response features.

[0092] Specifically,

[0093] The area of ​​pre-stack 3D seismic data in the study area is 400 km 2 Due to the absorption and attenuation of the overlying salt rock, the apparent main frequency of the subsalt seismic data in this area is about 17Hz. The main oil layer in this area is a set of microbial reef reservoirs. Figure 1 The deep-domain subsalt carbonate reef seismic convexity parameter attribute imaging processing completed by the technical process has achieved good application results:

[0094] (1) Establishing the velocity model of the work area: using Dix formula (1), the root mean square velocity volume is converted into a layer velocity volume, and the target formation test velocity obtained by drilling is used for correction;

[0095] The Dix conversion formula is:

[0096]

[0097] In the formula, v n is the formation velocity limited by the (n-1)th formation interface above and the nth formation interface below. τ n and τ n-1 is the two-way zero offset time, V n and V n-1 is the corresponding RMS velocity.

[0098] The velocity model established is as follows Figure 3 As shown, the speed changes greatly due to the influence of salt.

[0099] (2) Using the established velocity model, the CRP gathers in the time domain are converted into CRP gathers in the depth domain. The time-depth conversion formula (2) is:

[0100]

[0101] where Z is the thickness of the reflecting layer, The formation is the average velocity, τ 0 is the double-layer reflection time of the formation, z n is the thickness of the nth sublayer, Vn is the layer velocity of the nth sublayer, τ n is the double-layer reflection time of the nth layer.

[0102] The CRP original offset gather in the time domain is converted into the CRP original gather data in the depth domain as follows: Figure 4 As shown in the figure, the amplitude wild values ​​of the acquired signal and the long-channel NMO distortion seriously affect the quality of subsalt reef imaging.

[0103] (3) Using the advantageous offset range of the target layer segment, the amplitude outliers of the offset distance and the amplitude of the NMO distortion at the far offset distance are removed.

[0104] (4) De-noising and time-varying static correction of CRP offset gathers.

[0105] The CRP gathers in the depth domain are de-noised, de-noised and time-varying statically corrected to obtain the CRP gathers in the depth domain. Figure 5 As shown, the quality has been significantly improved.

[0106] (5) The offset gather of the CRP gather is converted into an angle gather, and the reflection angle θ of the formation is calculated by formula (3):

[0107]

[0108] Among them, X represents the size of the offset, and Z is the thickness of the underground reflective layer.

[0109] The converted depth domain CRP angle gathers are as follows: Figure 6 shown.

[0110] (6) Amplitude stacking imaging of dominant angle gathers, the stacking formula is (4)

[0111] A i =∑a ij (4)

[0112] Among them, A i represents the i-th sample point of a seismic record, a ij Represents the amplitude value of the jth earthquake record on the CRP gather.

[0113] (7) Smoothing and filtering the imaged seismic trace data.

[0114] The imaging of carbonate reefs under salt in a two-dimensional seismic profile in the depth domain of the study area is as follows: Figure 8 As shown, the reef morphology is not easy to identify.

[0115] (8) Based on the previous drilling research results and seismic horizon interpretation results, determine the target layer for the development of subsalt reefs.

[0116] (9) The second-order spatial difference method in the depth domain is implemented for the target body to calculate the convexity characteristic parameters of the function curvature. We assume that the curvature expression of the two-dimensional function Z(x) is

[0117]

[0118] Let the quadratic surface function expression (6) of a point Z in space be as follows:

[0119] Z(x, y) = ax 2 +by 2 +cxy+dx+ey+f (6)

[0120] Here, the second-order difference method is mainly used to calculate the derivative. To calculate the curvature of a point on the surface expressed by equation (6), the values ​​of the surrounding grid points are used to fit the local quadratic surface using the least squares method. The first three terms in equation (6) represent the quadratic ellipse fitting, and the last three terms represent the plane fitting. We use a 3×3 grid unit for approximation ( Figure 2 ), according to the definition of first-order and second-order derivatives, the coefficients in equation (6) are obtained as expression (7).

[0121]

[0122] Where: z1,…,z9 are Figure 1 The values ​​of the mesh nodes 1…9 at the level shown in the figure, Δx and Δy are the distances between the mesh nodes, and the convexity characteristic parameter K of the surface is calculated on this basis pos , K neg , morphological index S shape , the calculation formulas are as follows (8), (9), (10):

[0123]

[0124]

[0125]

[0126] (10) Based on the seismic data volume, calculate the convexity characteristic parameter K of the data volume for the target layer segment pos , K neg , morphological index S shape Attribute data body.

[0127] (11) Based on the calculation of K of the data body pos , K neg Convexity characteristic parameter, morphology index S shape Characteristics, combined with the convexity characteristics of carbonate bioherms, select subsalt carbonate bioherms targets and provide drilling well locations.

[0128] After processing, the two-dimensional seismic profile is as follows: Fig. 9 As shown in Figure 2, the seismic convexity attribute parameter K under the top interface of BVE100 is pos The bioherm identification profile. The top hillock reflection feature and the weak or even chaotic reflection feature inside can be seen, which are initially considered to be the seismic response characteristics of the bioherm beach. The seismic profile cannot effectively describe the specific morphological characteristics of the bioherm, and its distribution range cannot be determined. However, the morphological characteristics of the bioherm can be clearly shown on the seismic convexity attribute profile after depth domain processing, and subsequent drilling confirmed the existence of the bioherm.

[0129] Comparison of subsalt depth domain and time domain imaging Figure 7 As shown in the figure, (A) is the subsalt seismic advantageous offset imaging in the depth domain, and (B) is the subsalt seismic advantageous offset imaging in the time domain. Comparative analysis shows that the imaging morphologies of subsalt carbonate rocks in different domains are very different. Only the subsalt bioherm imaging in the depth domain can truly reflect the development morphology of bioherms, which illustrates the importance of depth domain imaging and the importance of obtaining seismic convexity index in the depth domain.

[0130] The above examples prove that in the application of bioherm identification in a Brazilian subsalt oil field, the seismic uplift parameter characteristics are used to identify the morphology and target of subsalt bioherms. Subsequent drilling confirms that the effective thickness of the bioherms is more than 300 meters, which has great commercial value. It also proves the effectiveness and feasibility of the technology and has great application value.

[0131] It should be noted that the embodiments described above are only used to explain the present invention and do not constitute any limitation to the present invention. The description of the present invention by reference to typical embodiments should be interpreted as the words used therein are descriptive and explanatory words, rather than restrictive words. The present invention may be modified as specified within the scope of the claims of the present invention, and the present invention may be revised without departing from the spirit and scope of the present invention. Although the present invention described therein relates to specific methods, materials and embodiments, it does not mean that the present invention is limited to the specific examples disclosed therein. On the contrary, the present invention can be extended to all other methods and applications with the same functions.

Claims

1. A method for identifying a bioherm, comprising the following steps: S1: Establish velocity model; S2: using the velocity model to convert the CRP gathers in the time domain into CRP gathers in the depth domain; S3: converting the offset gather of the CRP gather in the depth domain into an angle gather; S4: using the angle gather amplitude stack imaging to determine the target layer section for reef development; S5: Calculate the seismic morphological index attribute data volume of the target layer segment to identify the bioherm; In S1, the root mean square velocity volume is converted into the layer velocity volume using Dix formula (1), and the target formation test velocity is used for correction; In the formula, v n is the interval velocity of the nth layer, τ n and τ n-1 is the two-way zero-offset travel time of the corresponding formation, V n and V n-1 is the root mean square velocity of the corresponding formation; In S2, the time-depth conversion formula (2) for converting the CRP gathers in the time domain into the CRP gathers in the depth domain is: where Z is the thickness of the reflecting layer, is the average velocity of the overlying strata, τ0 is the double-layer reflection time, z n is the thickness of the nth layer, Vn is the interval velocity of the nth layer, τ n is the double-layer travel time of the seismic wave in the nth layer; In S5, the method for calculating the seismic morphological index attribute data volume of the target layer segment includes the following steps: Step A: Let the curvature expression (5) of the two-dimensional function Z(x) be Step B: Assume that the quadratic surface function expression (6) of a point Z in space is as follows: Z(x,y)=ax 2 +by 2 +cxy+dx+ey+f (6) Step C: Calculate the curvature of a point on the surface expressed by equation (6), and use the values ​​of the surrounding grid points to perform least squares fitting on the local quadratic surface. The first three terms in equation (6) represent quadratic ellipse fitting, and the last three terms represent plane fitting. Use 3×3 grid units for approximation. According to the definition of first-order and second-order derivatives, the coefficients in equation (6) are obtained as equation (7). Where: z1, ..., z9 are the values ​​of the surface mesh nodes 1 ... 9, Δx, Δy are the distances between the mesh nodes; Step D: Calculate the convexity characteristic parameter K of the surface pos , K neg , morphological index S shape , the calculation formulas are as follows (8), (9), (10): Step E: Based on the seismic data volume, calculate the convexity characteristic parameter K of the data volume pos , K neg , morphological index S shape Attribute data body.

2. The method according to claim 1, characterized in that: S2 also includes using the dominant offset range of the target layer segment to remove the amplitude wild value of the offset and the amplitude of the NMO distortion of the far offset to obtain a CRP offset gather; and / or denoising and performing time-varying static correction on the CRP offset gather; and / or S4 also includes smoothing and filtering the imaged seismic trace data.

3. The method according to claim 1, characterized in that In S5, the seismic morphological index attribute data volume of the target layer segment is calculated by the spatial second-order difference derivation method in the depth domain.

4. The method according to claim 3, characterized in that: The seismic morphological index attribute data body includes a convexity characteristic parameter K pos , K neg , morphological index S shape .

5. The method according to any one of claims 1 to 4, characterized in that In S3, the offset gather of the CRP gather in the depth domain is converted into an angle gather, and the formula (3) for calculating the reflection angle θ of the formation is: Among them, X represents the size of the offset, and Z is the thickness of the underground reflective layer.

6. The method according to any one of claims 1 to 4, characterized in that In S4, the amplitude stacking imaging using the angle gathers is performed, and the stacking formula is (4) A i =∑a ij (4) Among them, A i represents the i-th sample point of a seismic record, a ij Represents the amplitude value of the jth earthquake record on the CRP gather.

7. The method according to any one of claims 1 to 4, characterized in that The angle gather is a dominant angle gather.

8. The method according to any one of claims 1 to 4, characterized in that In S5, based on the seismic morphological index attribute data body and in combination with the convexity characteristics of carbonate bioherms, bioherm targets are identified and drilling well locations are provided.

9. Use of the method according to any one of claims 1 to 8 in identifying bioherms.

10. The use according to claim 9, characterized in that: The bioherm is a carbonate bioherm.

11. The use according to claim 10, characterized in that: The bioherm is a subsalt bioherm.

Images