Method and device for predicting thin reservoir based on seismic data

Abstract

The invention provides a method and device for predicting a thin reservoir based on seismic data, and relates to the technical field of oil-gas exploration and development, and the method comprises the steps: S1, obtaining actual drilling data, seismic data and geological comprehensive analysis data of a target interval; s2, extracting a wavelet dominant frequency of seismic data of a target layer section, establishing a cross-correlation function of two signals of a seismic record waveform and a synthetic record waveform, obtaining synthetic record drift time, calculating an attenuation coefficient Qf by using the record drift time and the seismic record time, and completing high-precision synthetic record calibration based on a time-frequency attenuation relationship rule; and S3, drawing a cross plot of the sand body thickness and the impedance value, determining a wave impedance identification threshold, and outputting a seismic wave impedance inversion profile of the target layer and a sand body thickness map through threshold value constraint. According to the method, reservoir plane distribution and spatial distribution rules can be depicted finely, thin reservoir sandstone prediction precision and accuracy are improved, and risks caused by drilling are reduced.

Description

Technical Field

[0001] The present invention belongs to the technical field of oil and gas exploration and development, and particularly relates to a method and device for predicting thin reservoirs based on seismic data. Background Art

[0002] Seismic attributes refer to the geometric, kinematic, dynamic or statistical characteristics of seismic waves obtained through mathematical transformation of pre-stack or post-stack seismic data. Seismic attributes mainly include several categories such as time, amplitude, frequency, phase, coherence and attenuation. In the technical field of oil and gas exploration, seismic attributes are mainly applied to the prediction of reservoir distribution characteristics such as reservoir thickness, porosity, permeability, etc. Among them, the application of single seismic attributes such as amplitude and frequency of seismic data is relatively extensive.

[0003] In the prior art, in order to accurately reflect the distribution characteristics of reservoirs, it is necessary to amplify the characteristics of reservoir amplitude and frequency; the currently widely used sweet body attribute is the quotient of amplitude and the square root of frequency, but the sweet body attribute can only be adapted to reservoirs with strong amplitude and low-frequency characteristics. In the case of thin reservoirs (thickness less than one-quarter wavelength), due to interference, the reservoir shows characteristics of strong amplitude and high frequency. The higher the frequency, the smaller the quotient of amplitude and the square root of frequency, resulting in the sweet body attribute being unable to effectively reflect the distribution characteristics of thin reservoirs.

[0004] Chinese Patent 202010911856.X discloses a method and an electronic device for seismic identification of thin reservoirs and seismic weak reflection layers. The method includes: obtaining the geological characteristics and reservoir rock physical characteristics of all reservoirs in the target work area; determining the target reservoir based on the geological characteristics and reservoir rock physical characteristics; obtaining the three-dimensional post-stack seismic data of the target reservoir; performing a feasibility analysis on the three-dimensional post-stack seismic data; obtaining the resampled three-dimensional post-stack seismic data of the target reservoir based on the feasibility analysis result; obtaining the seismic data after differential derivation based on the resampled three-dimensional post-stack seismic data; respectively performing standardization processing on the resampled three-dimensional post-stack seismic data and the seismic data after differential derivation to obtain a first standardized data volume and a second standardized data volume, and performing seismic identification of thin reservoirs and seismic weak reflection layers. The invention can improve the high-frequency amplitude information while preserving the low-frequency components of the original seismic information, and can improve the accuracy of reservoir description.

[0005] In addition, Chinese Patent 201611240442.9 discloses a method for accurate seismic description of thin reservoirs, including: forward modeling to analyze the uncertainty characteristics of seismic reflection waves of thin reservoirs; performing trace integration on the forward modeled seismic data to analyze the uncertainty characteristics of seismic reflection trace integration of thin reservoirs; describing the geometric characteristics of seismic thin reservoirs based on the above analysis results. The invention effectively utilizes the uncertainty characteristics of seismic reflection wave interference of thin reservoirs and can accurately describe the geometric characteristics of reservoirs.

[0006] In the research process of the method for predicting thin reservoirs, there are the following problems in achieving accurate prediction of thin reservoirs based on conventional post-stack seismic data: (1) The layer changes rapidly, the stacking relationship between individual sand bodies is complex, and there are coexisting thick and thin reservoirs with thin interbeds and interlayers; (2) The main frequency is low, the frequency band is narrow, and the resolution restricts the prediction of thin reservoirs. The use of seismic data for reservoir prediction in thin interbeds cannot be fully accepted by people, mainly due to the restriction of seismic vertical resolution. Thin layers are difficult to directly image and identify on seismic reflection profiles, making it difficult to extract thin layer information. Therefore, there is an urgent need for a new method for predicting thin reservoirs based on seismic data to improve the prediction accuracy and accuracy of thin reservoir sandstones and reduce the risks brought by drilling. Summary of the Invention

[0007] In view of the problems existing in the prior art, the present invention provides a method and device for predicting thin reservoirs based on seismic data, which can finely depict the planar distribution and spatial distribution law of thin reservoirs, improve the prediction accuracy and accuracy of thin reservoir sandstones, and reduce the risks brought by drilling.

[0008] To achieve the above object, the technical solution adopted by the present invention is as follows:

[0009] First, the present invention provides a method for predicting thin reservoirs based on seismic data, including the steps of:

[0010] S1. Obtain the actual drilling data, seismic data, and geological comprehensive analysis data of the target interval;

[0011] S2. Extract the sub-wave main frequency of the seismic data of the target interval, establish the cross-correlation function of the two signals of the seismic record waveform and the synthetic record waveform, obtain the synthetic record drift time, and calculate the attenuation coefficient Q using the record drift time and the seismic record time; f ; Complete the high-precision synthetic record calibration based on the time-frequency attenuation relationship law;

[0012] S3. Determine the seismic wave impedance identification threshold for sand and mudstone through the intersection diagram of sand body thickness and seismic impedance value, draw the intersection diagram of sand body thickness and impedance value according to the seismic wave impedance inversion result, determine the wave impedance identification threshold, and output the seismic wave impedance inversion profile and sand body thickness map of the target layer through threshold value constraint.

[0013] Preferably, in step S1, the actual drilling data includes drilling data, logging data, and well logging data; the seismic data is three-dimensional post-stack seismic data; the geological comprehensive analysis data is the data for reservoir characteristic analysis, including the fluid properties of the reservoir, single-well sedimentary microfacies characteristics, and oil and gas reservoir and dynamic characteristic parameters.

[0014] Further preferably, the drilling data includes the lithology and physical properties of the target interval; the logging data includes gas logging shows, lithology, fluorescence display levels, etc. of the target interval; the logging data includes acoustic travel time, density curve, well diameter curve, and logging interpretation results, etc.

[0015] Further preferably, the 3D post-stack seismic data is a high-resolution amplitude-preserved data volume for lithology interpretation.

[0016] Preferably, in step S2, the extraction of the dominant frequency of the wavelet of the seismic data of the target interval is specifically: applying the acoustic travel time and density curve of the logging data to extract the dominant frequency of the wavelet of the 3D post-stack seismic data of the target interval.

[0017] Preferably, in step S2, the attenuation coefficient Q f , is obtained by using the cross-correlation method to complete the high-precision synthetic seismogram calibration based on time-frequency attenuation of the seismic data of the target interval.

[0018] Preferably, step S2 is specifically: by establishing the cross-correlation function of two signals, namely the seismic record waveform and the synthetic seismogram waveform, calculating the time delay corresponding to the maximum correlation value, that is, the synthetic seismogram drift time; then, using the synthetic seismogram drift time and the seismic record time to calculate the attenuation coefficient Q f ; after the above steps, the error value between the seismic interval velocity and the interval velocity of the acoustic logging data is approximately the velocity attenuation coefficient, and the accurate seismic synthetic seismogram and the accurate time-depth relationship can be obtained by using the velocity attenuation coefficient.

[0019] Preferably, in step S2, the attenuation coefficient Q f , is calculated by using formulas (1)-(5) as follows:

[0020]

[0021] In formula (1), R xy (m) represents the time series; the cross-correlation function value of x(n - m) and y(n); x(n - m) is the seismic record waveform; y(n) is the random signal; the sampling period is T, t = nT; N represents the number of sampling points within the finite integration time; m is the serial number of the time delay, and n takes an infinite value; the delay time τ = mT;

[0022] Assume that the seismic signal sequence is T seismic (with a total of R sampling points, R being infinite), and the synthetic seismogram signal sequence is T sonic (with a total of S sampling points, S taking about 1500), and the sampling rates of the two signals are the same; symmetrically take values before T sonic and after T sonic to obtain a pre-reference sequence T Ksonic (K sampling points, where K is approximately 1500), formula (2) is obtained using formula (1):

[0023]

[0024] In formula (2), T K is the cross-correlation value of the pre-reference sequence; n is the starting sampling point number involved in the operation in the received signal; Ts is the time interval of the seismic signal sampling points; K is the number of sampling points; i is a natural number (1, 2, 3,..., n, where n takes an infinite value);

[0025] If the calculated cross-correlation value is valid, it is adopted; otherwise, it is cleared; if when n traverses all the previous (S - K) received points and the maximum correlation value T Kmax and the correlation point number n Kmax are obtained, the transit time n KmaxTs of the synthetic record drift is thus obtained; calculate formula (3):

[0026]

[0027] In formula (3), T sonic is the synthetic record signal sequence of S sampling points; T seismic is the seismic signal sequence of R sampling points; n Kmax T s is the transit time of the synthetic record drift;

[0028]

[0029] In formula (4), V i Seismic is the seismic layer velocity, and V i Sonic is the acoustic logging layer velocity, in m / s; i is a natural number (1, 2, 3,..., n, where n takes an infinite value);

[0030] Let the expression of the average reflection coefficient (Ri) of the acoustic logging formation thickness be formula (5):

[0031]

[0032] In formula (5), ρ i is the layer density, in g / cm 3 ; V i is the layer velocity, in m / s; i is a natural number (1, 2, 3,..., n, where n takes an infinite value).

[0033] As can be seen from formula (5), calibrating the seismic synthetic record using the velocity attenuation coefficient will not change the calculated high-frequency acoustic wave reflection coefficient and can ensure the stability of the reflection coefficient and the seismic wave group relationship. Based on this, the calculation of formula (5) is effective.

[0034] Preferably, in step S3, the drawing of the crossplot of sand body thickness and impedance value is specifically as follows: Combining the lithology data of the well, a crossplot is drawn through the sand body thickness and the seismic impedance value, and a comprehensive analysis of the crossplot is carried out.

[0035] Further preferably, step S3 is specifically as follows: Combining the lithology data of the actual drilled well, through the analysis of the crossplot of the sand body thickness and the seismic impedance value, the seismic wave impedance identification thresholds of sand and mudstone are determined, and through the constraint of the threshold value, seismic wave impedance inversion is carried out, so as to extract the seismic inversion profile of the target layer and the thickness map of the thin sand body based on time-depth conversion.

[0036] Then, the present invention provides an inversion profile of the target layer and a sand body thickness map, which are established by the above method.

[0037] Finally, the present invention provides a device for seismic prediction of thin reservoirs, including a data acquisition module, a calculation module, and an output module.

[0038] Preferably, the data acquisition module: is used to acquire the drilling, logging, well logging data, 3D post-stack seismic data, and geological comprehensive relevant analysis data of the target layer section;

[0039] Preferably, the calculation module: is used for log curve sensitivity analysis and multi-parameter crossplot analysis, including log curve editing and environmental correction; high-precision synthetic record calibration based on time-frequency attenuation; calculating the attenuation coefficient Q using the drift time and the seismic record time f 。

[0040] Preferably, the output module: is used for broadband wave impedance data volume and high- and low-frequency seismic trace merging, determining the sand-mudstone wave impedance identification threshold, and outputting the thin sandstone seismic wave impedance inversion data volume.

[0041] In the method of the present invention, due to the existence of fractures in the post-stack 3D seismic data volume, it will cause the anisotropy (HTI) of the formation, and the anisotropy of the formation will cause changes in the seismic response characteristics on the wide-azimuth seismic data, that is, the velocity varies with the azimuth angle (VVAz) and the amplitude varies with the azimuth angle (AVAz). These changes in the seismic response characteristics can be used to identify thin reservoirs. Therefore, the wide-azimuth seismic data provides the possibility for thin reservoir prediction.

[0042] The method described in the present invention mainly uses strong wave group reflection boundaries to ensure the calibration accuracy of reservoir groups, changing the magnitude and energy of the reflection coefficient; while the local cumulative error is evenly distributed, which also changes the magnitude of the reflection coefficient.

[0043] In the method described in the present invention, in order to improve the calibration accuracy of thin layers and clarify the internal reflections and seismic wave group characteristics, it is necessary to perform environmental correction on well logging data to eliminate the positive and negative baseline errors introduced by various factors, so as to ensure that the acoustic wave curve can reflect the true formation velocity.

[0044] In the method described in the present invention, regarding the calculation of the attenuation coefficient, it is necessary to eliminate the cumulative time difference caused by the formation absorption effect, obtain an accurate time-depth relationship, and improve the calibration of thin layer reflections within the seismic wave group and the accuracy of seismic inversion while ensuring the stability of the relationship between the thin layer reflection coefficient and the seismic wave group.

[0045] Compared with the prior art, the present invention has the following beneficial effects:

[0046] 1. The method of the present invention can effectively identify glutenite reservoirs with a buried depth greater than 4500 m in the deep layer, efficiently identify thin sand layers, reduce the multi-solution problem existing in the prediction of deep glutenite thin sand layers; optimize the algorithm and workflow of seismic inversion; improve the success rate of deploying and demonstrating development wells to more than 90%; improve the accuracy of the main controlling factors and hydrocarbon-bearing properties of thin reservoirs in deep sandstone and mudstone, laying a theoretical foundation for the efficient development of thin reservoir oil and gas reservoirs with interbedded sandstone and mudstone.

[0047] 2. The reservoir sweet spots predicted by the method of the present invention can greatly promote the production capacity construction of the reserves of thin reservoir oil and gas reservoirs with interbedded sandstone and mudstone, and then form an effective comprehensive prediction method and technical system for thin sandstone reservoirs, realizing the accurate prediction of thin sand layers with a buried depth of more than 4500 m and a thickness of 3 - 8 m.

[0048] 3. The method of the present invention forms an effective identification method for thin sand layers in the Shushanhe Formation through fine comparison of sand bodies and fine interpretation of structures, combined with means such as attribute, pre-stack, and post-stack inversion. The application of this method effectively promotes the evaluation and development process of the Shushanhe Formation in this area and improves the accuracy of oil and gas discovery and exploitation. BRIEF DESCRIPTION OF THE DRAWINGS

[0049] Figure 1 is a flowchart of the method for predicting thin reservoirs based on seismic data in the present invention.

[0050] Figure 2 is a schematic diagram of the device for predicting thin reservoirs based on seismic data in the present invention.

[0051] Figure 3 is a crossplot of the thickness of the No. 1 sand body in the middle section of the Cretaceous Shushanhe Formation in Well S area and the seismic inversion wave impedance value in the application of the method of the present invention.

[0052] Figure 4 It is a schematic diagram of a seismic time migration profile and a seismic wave impedance inversion profile in the application of the method of the present invention.

[0053] Figure 5 It is a planar distribution map of seismic inversion wave impedance in the S well area in the application of the method of the present invention. Detailed implementation manners

[0054] The following non-limiting embodiments can enable those of ordinary skill in the art to more comprehensively understand the present invention, but do not limit the present invention in any way. The following content is only an exemplary illustration of the scope claimed by the present invention. Those skilled in the art can make various changes and modifications to the present invention based on the disclosed content, and it should also fall within the scope claimed by the present invention. When an embodiment gives a numerical range, it should be understood that unless otherwise specified in the present invention, any value at both ends of each numerical range and any value between the two ends can be selected. Unless otherwise defined, all technical and scientific terms used in the present invention have the same meaning as commonly understood by those of ordinary skill in the technical field to which the present invention belongs. The present invention will be further described below by way of specific embodiments.

[0055] Embodiment 1

[0056] The flow chart of the method for predicting thin reservoirs based on seismic data according to the present invention is shown in Figure 1 .

[0057] S1. Obtain the actual well drilling data, seismic data, and geological comprehensive analysis data of the target interval;

[0058] The actual well drilling data includes drilling data (including the lithology and physical properties of the target interval), logging data (including gas logging shows, lithology, and fluorescence show levels of the target interval), and well logging data (including acoustic travel time, density curve, well diameter curve, and well logging interpretation results);

[0059] The seismic data is 3D post-stack seismic data, which is a high-resolution amplitude-preserved data volume for lithology interpretation;

[0060] The geological comprehensive analysis data is data for reservoir characteristic analysis, including the fluid properties of the reservoir, sedimentary microfacies characteristics, and reservoir and dynamic characteristic parameters.

[0061] S2. Extract the dominant frequency of the wavelet of the 3D post-stack seismic data of the target interval by using the acoustic travel time and density curve of the well logging data. By establishing the cross-correlation function of the two signals of the seismic record waveform and the synthetic record waveform, calculate the time delay corresponding to the maximum correlation value, that is, the synthetic record drift time; then, calculate the attenuation coefficient Q by using the synthetic record drift time and the seismic record time f; After the above steps, the error value between the seismic interval velocity and the interval velocity of the acoustic logging data is approximately the velocity attenuation coefficient; the accurate seismic synthetic record and the accurate time-depth relationship can be obtained by using the velocity attenuation coefficient;

[0062] The attenuation coefficient Q f is confirmed by the following formulas (1)-(5):

[0063]

[0064] In formula (1), R xy (m) represents the time series; the cross-correlation function value of x(n - m) and y(n); x(n - m) is the seismic record waveform; y(n) is the random signal; the sampling period is T, t = nT; N represents the number of sampling points within the finite integration time; m is the serial number of the time delay, and n takes an infinite value; the delay time τ = mT;

[0065] Assume that the seismic signal sequence is T seismic (with a total of R sampling points, R being an infinite value), and the synthetic record signal sequence is T sonic (with a total of S sampling points, S taking a value of approximately 1500), and the sampling rates of the two signals are the same. Symmetrically take values before T sonic and after T sonic to obtain a pre-reference sequence T K sonic (K sampling points, K taking a value of approximately 1500) centered around the center point S / 2, and use formula (1) to obtain formula (2):

[0066]

[0067] In formula (2), T K is the cross-correlation value of the pre-reference sequence; n is the starting sampling point serial number of the received signal participating in the operation; Ts is the time interval of the seismic signal sampling points; K is the number of sampling points; i is a natural number (1, 2, 3,..., n, n taking an infinite value);

[0068] If the calculated cross-correlation value is valid, it is adopted, otherwise it is cleared; if when n traverses all the first (S - K) received points and the maximum correlation value T Kmax and the correlation point serial number n Kmax are obtained, thus obtaining the transit time n KmaxTs of the synthetic record drift; calculate formula (3):

[0069]

[0070] In formula (3), T sonic is the synthetic record signal sequence of S sampling points; T seismicis the seismic signal sequence of R sampling points; n Kmax T s is the transit time of the synthetic record drift;

[0071]

[0072] In formula (4), V i Seismic is the seismic interval velocity, V i Sonic is the acoustic logging interval velocity, m / s; i is a natural number (1, 2, 3,..., n, where n takes an infinite value).

[0073] Let the average reflection coefficient of the acoustic logging formation thickness (R i ) have the expression as formula (5):

[0074]

[0075] In formula (5), ρ i is the layer density, g / cm 3 ; V i is the layer velocity, m / s; i is a natural number (1, 2, 3,..., n, where n takes an infinite value).

[0076] It can be seen from formula (5) that calibrating the seismic synthetic record using the velocity attenuation coefficient will not change the calculated high-frequency acoustic wave reflection coefficient and can ensure the stability of the reflection coefficient and the seismic wave group relationship. Based on this, formula (5) is effectively calculated.

[0077] S3. Combine the lithology data of the actual drilled wells, draw a crossplot through the sand body thickness and the seismic impedance value, and conduct a comprehensive analysis of the crossplot; through the analysis of the crossplot of the sand body thickness and the seismic impedance value, determine the seismic wave impedance identification threshold for sand and mudstone, and through the constraint of the threshold value, conduct seismic wave impedance inversion, so as to extract the seismic inversion profile of the target layer and the thickness map of the thin sand body based on the time-depth conversion.

[0078] Figure 2 The device for predicting thin reservoirs based on seismic data according to the method of Embodiment 1 of the present invention is given. The device includes a data acquisition module 201, a calculation module 202, and an output module 203.

[0079] Data acquisition module 201: Acquire the drilling, logging, mud logging data, 3D post-stack seismic data, and geological comprehensive analysis data of the target layer section.

[0080] Calculation module 202: Conduct sensitivity analysis of logging curves for the target layer section, including logging curve editing and environmental correction; high-precision synthetic record calibration based on time-frequency attenuation; seismic wave impedance inversion; extraction of seismic attributes for amplitude values and multi-parameter crossplot analysis of sand body thickness, etc.

[0081] Output module 203: Merge the broadband wave impedance data volume and the high- and low-frequency seismic traces, and output the data volume of thin sandstone seismic wave impedance inversion.

[0082] Application Example 1

[0083] Apply the method described in Example 1 (Steps S1 - S3) to thin reservoir identification. Specifically, select the Cretaceous Shushanhe Formation in the S Well Area of the Tarim Basin as an application example.

[0084] Step (1): Based on forward modeling of the model, select seismic attributes related to lithology and reservoir fluids, and extract seismic attributes for each sub-layer; that is, obtain drilling, logging, well logging data, and 3D post-stack seismic data (high-resolution amplitude-preserved data volume for lithology interpretation) and geological comprehensive analysis data of the target interval. The drilling data includes, for example, the lithology and physical properties of the target interval; the logging data includes, for example, gas logging shows, lithology, and fluorescence display levels of the target interval; the geological comprehensive analysis data includes: fluid properties of the reservoir, sedimentary microfacies characteristics, and reservoir and dynamic characteristic parameters. The fluid properties of the reservoir, sedimentary facies characteristics of a single well, reservoir parameters, etc. are used for reservoir characteristic analysis; the well logging data includes acoustic travel time, density curve, well diameter curve, and well logging interpretation results, etc.

[0085] Step (2): First, it is necessary to perform consistency standardization processing on the well logging curves, including borehole collapse correction. Overall, borehole collapse has a greater impact on density logging. By comparing the target interval, determine the offset of the data in the collapsed interval, and then correct the offset. After correction, select a standard well that penetrates the target layer, has relatively complete well logging items, and generally good well logging curve quality. Based on the characteristics of the same formation thickness and type, correct the curves of multiple wells based on the same set of strata.

[0086] Step (3): According to the characteristics of sedimentary cycles and the constraint of local marker beds, combine well and seismic data to complete the small-layer closure comparison of the skeleton profile. On this basis, complete the small-layer division and comparison of the Cretaceous Shushanhe Formation in 19 wells in the S Well Area. Determine the time-depth relationship through synthetic seismogram calibration, calculate and extract seismic inversion wavelets. Then, take the 3D post-stack seismic data as the object, make synthetic seismograms for each well in the seismic inversion area one by one, determine the time-depth relationship of each well, obtain the attenuation coefficient using the cross-correlation method, and extract the seismic wavelet beside the well. Finally, obtain a comprehensive wavelet from these wavelets. After calibrating multiple wells in the S Well Area, it is found that the No. 4 gas-bearing sand body in the middle section of the Cretaceous Shushanhe Formation is calibrated at the trough position, showing low-frequency strong reflection characteristics. This Step (2) - (3) corresponds to Steps S2 and S3 in the method of Example 1.

[0087] Calibration analysis of the reflection position and seismic response characteristics of the No. 1 thin sandstone (with a thickness of 3 - 8 m) within the λ / 4 wavelength in the middle section of the thin-bed synthetic seismic record. Based on the same time-depth relationship, different dominant frequencies, such as 40 Hz, 60 Hz, and 80 Hz, are used to calibrate the original seismic data. Among them, the seismic wave groups of the 40 Hz and 60 Hz synthetic records are all troughs, and there is no response from the sand body. The 80 Hz forward synthetic record can identify the reflection interface of the No. 1 single sand body in the middle section (80 Hz thin-bed trough response), thus determining that the minimum value of the dominant frequency of the wavelet theoretically resolvable for the No. 1 sandstone in the middle section is 80 Hz.

[0088] The thin-bed fully utilizes the identifiability of the wide-band seismic signals of the thin-bed, and uses the identifiable waveforms and well logging high-frequency information for two-way constrained quantitative solution of the thin-bed information. This technology can effectively solve the problem of quantitative prediction of thin reservoirs with a thickness of 3 - 8 m in this area.

[0089] Step (4): Under the control of sedimentary microfacies, through the comparison of the north-south-east-west cross-well profiles in the S well area, it can be seen that the No. 4 sand body in the middle section is encountered by multiple wells and thins out at the tip of Well S6, with a sand body thickness of 5.5 - 10 m. The No. 4 sand body in the middle section develops two sets of upward-graded sand bodies. The oil and gas layers are developed in the upper-graded sand body, with an oil-bearing sand body thickness of 5 - 9 m in the upper part, and no oil and gas in the lower-graded sand body. Through the comparison of the east-west cross-well profiles in the S well area, it can be seen that the No. 4 sand body in the middle section is continuously developed and encountered by multiple wells. Among them, Well S4 (oil and gas layer) and Well S5 (oil and gas layer) encounter oil and gas showing layers. The reservoir characteristics are characterized by determining the reservoir boundary based on "attributes + actual drilling", and the overall shape of the sand body is channel-like.

[0090] The seismic attributes in the S well area are relatively high, and the sand bodies are relatively well developed. Through the statistics of the amplitude values of 19 wells in the S well area, the current amplitude values of the existing wells are > 1000, and reservoirs are developed. Through statistical analysis, the correlation between the amplitude and the sand body thickness is good, and the distribution characteristics of the reservoir and fluid are qualitatively analyzed.

[0091] According to the reservoir thickness and lateral variation characteristics, select appropriate effective samples and high cut-off frequencies to carry out waveform indication inversion. Step (5): Through the cross-plot of the thickness of the No. 1 sand body in the middle section of the Cretaceous Shushanhe Formation in the S well area and the seismic inversion wave impedance value ( Figure 3 ) analysis, determine the impedance identification threshold for sandstone and mudstone. There is a good negative correlation between the wave impedance and the sand body thickness. The lower the impedance value, the greater the sand thickness. The current impedance values of the existing wells are < 7750, and reservoirs are developed. Through statistical analysis, the correlation between the impedance value and the sand body thickness is good, and the coincidence rate reaches 87% ( Figure 4 , Figure 5 ). This step (5) corresponds to step S3 in the method of Example 1. Figure 4The upper figure is the seismic time migration profile, and the lower figure is the seismic wave impedance inversion profile. Among them, in the seismic time migration profile, on the top surface of the No. 4 sand body in the middle section of the target interval, the triangles represent the sand body boundaries; on the seismic wave impedance inversion profile, there are the No. 4 sand body in the middle section and mudstone respectively. The darker the color, the more sandy content, and the lighter the color, the more muddy content; through the above profiles, the corresponding relationship between the sand body boundaries predicted by seismic wave impedance and the amplitude changes in the seismic time migration profile is shown. The amplitude weakens in the southwestern and northeastern directions, and at the same time, the inversion also shows the sand body pinch-out.

[0092] Figure 5 It is the plane distribution map of seismic inversion wave impedance in the S well area, with low wave impedance, high wave impedance, and sand body pinch-out lines respectively. Through the statistics of the drilled wells in the S well area, the wave impedance value is negatively correlated with the sand body thickness. The lower the wave impedance, the thicker the sand body. The sand body in the No. 4 gas-bearing trap in the middle section is 5 - 12 m thick, with an average porosity of 23%.

[0093] Step (6): Extract the inversion profile of the target interval and the sand body thickness map through threshold constraint. This step corresponds to step S3 in the method.

[0094] It should be noted that the above steps (2) and (4) are preferred steps, and the method further improves the prediction accuracy of the present invention. Practice has proved that in the Cretaceous Shushanhe Formation in the S well area, the sand and mud interbeds are frequent, and the thickness of a single sand body is thin (<10 m). In view of this characteristic, the researchers have carried out a large number of reservoir prediction and identification works. The research results show that the forward simulation experiment confirms that the amplitude attribute is most effective in identifying reservoirs with a thickness of about 20 m (the tuning resolution ability of λ / 4). In addition, reservoir prediction has been carried out on reservoirs with a thickness of about 10 m using seismic amplitude attributes, and good results have been obtained. Using the natural energy of the edge water for development, the recoverable reserves are over 100 million cubic meters of gas and more than 20,000 tons of condensate oil. The designed production capacity is 30,000 cubic meters. And the gas-bearing reservoirs with a thickness of about 10 m in the S well area have basically been put into development, and the reservoirs with a thickness of less than 5 m need to be further characterized. The reservoir thickness is 2.5 - 4 m. Through horizon calibration and extraction of seismic attributes, the traditional method has a low coincidence degree between seismic attributes and sand bodies, and the amplitude attribute basically does not reflect the reservoir distribution. The reservoir identification method needs to be further studied and explored. The method of the present invention forms an effective identification method for thin sand layers in this area through fine comparison of sand bodies and fine interpretation of structures, combined with means such as attribute, pre-stack, and post-stack inversion, laying a foundation for the evaluation and research of thin reservoirs. Through the evaluation and research of key areas in the area using this method, favorable targets are identified, providing a theoretical basis for well placement and effectively promoting the process of evaluation and development.

[0095] Finally, it should be noted that the above content is only used to illustrate the technical solution of the present invention, rather than a limitation on the protection scope of the present invention. Any simple modification or equivalent replacement of the technical solution of the present invention by those of ordinary skill in the art does not depart from the essence and scope of the technical solution of the present invention.

Claims

1. A method for predicting thin reservoirs based on seismic data, characterized in that, Including the steps: S1. Obtain the actual drilling data, seismic data, and geological comprehensive analysis data of the target interval; S2. Extract the dominant frequency of the wavelet of the seismic data in the target interval, establish the cross-correlation function of the two signals of the seismic record waveform and the synthetic seismogram waveform, obtain the drift time of the synthetic seismogram, and calculate the attenuation coefficient Q using the drift time and the seismic record time, f thus completing the high-precision calibration of the synthetic seismogram based on the time-frequency attenuation relationship law; S3. Determine the seismic wave impedance identification thresholds for sand and mudstone through the cross-plot of sand body thickness and seismic impedance value. Draw the cross-plot of sand body thickness and impedance value based on the seismic wave impedance inversion results, determine the wave impedance identification thresholds, and through threshold value constraints, output the seismic wave impedance inversion profile and sand body thickness map of the target layer.

2. The method according to claim 1, characterized in that, In step S1, the actual drilling data includes drilling data, logging data, and well logging data; the seismic data is three-dimensional post-stack seismic data; the geological comprehensive analysis data is the data for reservoir characteristic analysis, including the fluid properties, sedimentary microfacies characteristics, oil and gas reservoir, and dynamic characteristic parameters of the reservoir.

3. The method according to claim 2, wherein The actual drilling data includes the lithology and physical properties of the target interval; the logging data includes gas logging shows, lithology, and fluorescence display levels of the target interval; the well logging data includes acoustic travel time, density curve, well diameter curve, and well logging interpretation results; The three-dimensional post-stack seismic data is a high-resolution amplitude-preserved data volume for lithology interpretation.

4. The method according to claim 1, characterized in that, In step S2, the specific method for extracting the dominant frequency of the seismic wavelet of the target interval is: apply the acoustic travel time and density curve of the well logging data to extract the dominant frequency of the three-dimensional post-stack seismic data of the target interval.

5. The method according to claim 1, wherein In step S2, the attenuation coefficient Q f , is obtained by using the cross-correlation method to complete the high-precision synthetic record calibration of the seismic data in the target interval based on time-frequency attenuation.

6. The method according to claim 1, wherein Step S2 specifically includes: by establishing the cross-correlation function of two signals, namely the seismic record waveform and the synthetic record waveform, calculating the time delay corresponding to the maximum correlation value, which is the synthetic record drift time; then, using the synthetic record drift time and the seismic record time to calculate the attenuation coefficient Q f ; After the above steps, the error value between the seismic layer velocity and the layer velocity of the acoustic well logging data is approximately the velocity attenuation coefficient. Using the velocity attenuation coefficient, an accurate seismic synthetic record and an accurate time-depth relationship can be obtained.

7. The method according to claim 1, wherein In step S2, the attenuation coefficient Q f , is calculated using formulas (1)-(5) as follows: In formula (1), R xy (m) represents a time series; the cross-correlation function value of x(n - m) and y(n); x(n - m) is the seismic record waveform; y(n) is a random signal; the sampling period is T, t = nT; N represents the number of sampling points within a finite integration time; m is the sequence number of the time delay, n takes an infinite value; the delay time τ = mT; Assume that the seismic signal sequence is T seismic (with a total of R sampling points, where R is an infinite value), and the synthetic record signal sequence is T sonic (with a total of S sampling points), and the sampling rates of the two signals are the same; for T sonic before and T sonic after symmetrically take values, and obtain the pre-reference sequence T K sonic (with a total of K sampling points), and use formula (1) to obtain formula (2): In formula (2), T K is the cross-correlation value of the pre-reference sequence; n is the starting sampling point number involved in the operation in the received signal; Ts is the time interval of the seismic signal sampling points; K is the number of sampling points; i is a natural number (1, 2, 3,..., n, where n takes an infinite value); If the calculated cross-correlation value is valid, it is adopted; otherwise, it is cleared. If, when n takes values that traverse all the previous (S - K) receiving points, the maximum correlation value T Kmax and the correlation point sequence number n Kmax are obtained, thereby obtaining the transit time n of the synthetic record drift KmaxTs ; Calculation formula (3): In formula (3), T sonic is the synthetic record signal sequence of S sampling points; T seismic is the seismic signal sequence of R sampling points; n Kmax T s is the transit time of the synthetic record drift; V i Seismic = V i Sonic Q f , formula (4), In formula (4), V i Seismic is the seismic layer velocity, and V i Sonic is the sonic logging layer velocity, in m / s; i is a natural number (1, 2, 3,…, n, where n takes an infinite value); Assume that the expression for the average reflection coefficient (Ri) of the acoustic well logging formation thickness is formula (5): In formula (5), ρ i is the layer density, g / cm 3 ; V i is the layer velocity, m / s; i is a natural number (1, 2, 3,…, n, where n takes an infinite value); It can be seen from formula (5) that using the velocity attenuation coefficient to calibrate the seismic synthetic record will not change the calculated high-frequency acoustic wave reflection coefficient and can ensure the stability of the reflection coefficient and seismic wave group relationship. Based on this, the calculation of formula (5) is effective.

8. The method according to claim 1, characterized in that, In step S3, the specific method for drawing the cross-plot of sand body thickness and impedance value is: combine the lithology data of the drilling, draw the cross-plot through the sand body thickness and seismic impedance value, and conduct comprehensive analysis on the cross-plot.

9. The method according to claim 1, wherein Step S3 is specifically: combine the lithology data of the actual drilling, analyze through the cross-plot of sand body thickness and seismic impedance value, determine the seismic wave impedance identification thresholds for sand and mudstone, and through threshold value constraints, conduct seismic wave impedance inversion, so as to extract the seismic inversion profile of the target layer and the thickness map of thin sand bodies based on time-depth conversion.

10. An inversion profile of the target layer and a sand body thickness map, characterized in that, Established by the method described in any one of claims 1-9.

11. An apparatus for seismic prediction of thin reservoirs, characterized in that, Including a data acquisition module, a calculation module, and an output module; The calculation module: for well logging curve sensitivity analysis and multi-parameter cross-plot analysis, including well logging curve editing and environmental correction; high-precision synthetic record calibration based on time-frequency attenuation; calculating attenuation coefficient Q using drift time and seismic record time f .

Images