A method for identifying oil and gas-bearing patterns in reservoirs based on data reconstruction

Abstract

The present invention provides a method for identifying reservoir oil and gas patterns based on data reconstruction. This method studies the abnormal characteristics of reservoir oil and gas from the perspectives of seismic signal amplitude, frequency, attenuation, and other aspects to explore new oil and gas identification methods. The method includes: using the Levenberg-Marquardt algorithm to iteratively reconstruct each signal of the seismic data separately; using the instantaneous phase and instantaneous frequency obtained through the Hilbert transform, and the frequency and attenuation factor obtained through the wavelet transform as initial values; performing iterative operations on the original seismic signal data; and finally, using the reconstructed signal data results and the abnormal responses of the reconstructed attenuation factor, frequency, and other parameters, combined with experience and professional knowledge, to determine the oil and gas content of the reservoir. This method helps to determine the oil and gas content of the reservoir and has practical significance for directly finding oil and gas, and reducing the cost of oil and gas extraction.

Description

Technical Field

[0001] The present invention relates to the technical field of oil and gas exploration, and in particular to a method for identifying reservoir oil and gas patterns based on data reconstruction. Background Art

[0002] Today, the contradiction between ever-increasing oil and gas demand and dwindling remaining resources is becoming increasingly prominent. The complex conditions under which oil and gas reservoirs are formed, their deep burial depths, the variable properties of the filling fluids, the limitations of research methods, and the complexity of the geological landscape all contribute to the significant difficulty of exploration and development.

[0003] The oil and gas content of different reservoirs varies greatly. Due to the high cost of exploration and mining, oil and gas content testing needs to be carried out in advance. Only when the test results show that the oil and gas content of the reservoir reaches a certain standard and it is determined to have mining value can it be put into mining.

[0004] Seismic records contain multiple seismic attributes, such as amplitude, frequency, and phase, which can reflect the characteristics of underground reservoirs from different perspectives. However, current identification methods based on seismic attributes are inefficient and unintuitive. Therefore, a new method for effectively identifying oil and gas properties is urgently needed. Summary of the Invention

[0005] In view of the deficiencies in the prior art, the present invention provides a reservoir oil and gas pattern recognition method based on data reconstruction, which solves the problems of low recognition efficiency and non-intuitive recognition results in the prior art.

[0006] The above technical objectives of the present invention are achieved through the following technical solutions:

[0007] A method for identifying reservoir oil and gas patterns based on data reconstruction, characterized by:

[0008] Step 1: Read the earthquake data in sgy format and set the time corresponding to the data;

[0009] Step 2: Read the shear wave velocity v of the corresponding formation of the seismic data s , longitudinal wave velocity v p and density ρ;

[0010] Step 3, calculate the instantaneous frequency f of each sampling point of each seismic trace through wavelet transform;

[0011] Step 4: Use the shear wave velocity v s , longitudinal wave velocity v p and density ρ to find the attenuation factor α of each sampling point of each seismic trace;

[0012] Step 5: Calculate the instantaneous amplitude A and instantaneous phase of each sampling point of each seismic channel

[0013] Step 6: Reconstruct the seismic data, save the reconstructed signal, reconstructed attenuation factor and reconstructed frequency, and set them as data in sgy format;

[0014] Step 7: Use Savitzky-Golay convolution smoothing algorithm to smooth the reconstruction result;

[0015] Step 8: Determine the gas content based on the abnormal information in the calculation results.

[0016] The present invention is further configured as follows: in step 6, the single channel of the original seismic data is processed, and the Levenberg-Marquardt algorithm is used for iterative calculation, and the instantaneous amplitude A and the instantaneous phase The algorithm keeps the attenuation factor α and the instantaneous frequency f unchanged, and uses the two new parameters obtained each time, namely the attenuation factor α and the instantaneous frequency f. The above four parameters are used as initial values. According to the reconstruction signal formula in the Prony algorithm, an initial reconstructed signal value is obtained. The reconstructed signal value is continuously iterated to make it closer and closer to the original signal value. The operation is stopped after the termination condition is reached, and then the next data is calculated. After all the data are calculated, the reconstructed signal, the reconstructed attenuation factor, and the reconstructed frequency are saved.

[0017] The present invention is further configured as follows: In step 3, the method for obtaining the instantaneous frequency f of each sampling point of each seismic trace by wavelet transform is as follows: first determine the sampling frequency f of the seismic data s , sampling period T, number of sampling points N, and then use complex morlet wavelet to perform wavelet transform; suppose the length of the scale sequence is equal to the number of sampling points N, obtain the scale sequence scals, and convert the scale sequence scals into a frequency sequence F, perform wavelet transform on each seismic data with the complex morlet wavelet respectively, and determine the value in the frequency sequence F corresponding to the maximum value of the absolute value of the wavelet coefficient as the frequency size at that location, and the instantaneous frequency f corresponding to all seismic data can be calculated.

[0018] The present invention is further configured as follows: the calculation method of the attenuation factor α is as follows:

[0019]

[0020] Where, angular frequency ω=2πf,v p is the longitudinal wave velocity and ρ is the density.

[0021] The present invention is further configured as follows: the algorithm needs to calculate the Jacobian matrix for the currently processed seismic trace data, and then calculate the μ value of the Levenberg-Marquardt algorithm; the order of the matrix is N×2, where N is the data length of the seismic trace; the reconstructed signal formula in the Prony algorithm is:

[0022] x p (t) = Ae a·t cos(2πft+θ),

[0023] Obtain the partial derivative with respect to the instantaneous frequency f and the partial derivative with respect to the absorption coefficient α:

[0024]

[0025]

[0026] Finally, the Jacobian matrix J is obtained as:

[0027]

[0028] The present invention is further configured as follows: J is the Jacobian matrix, and N=J T J, μ is x times the maximum diagonal element of matrix N, that is, μ=x·max(N ii ), where i = 1, 2…N, and x is a positive floating point number in the range [0.0001, 0.001].

[0029] The present invention is further configured as follows: the conditions for terminating the iteration only need to meet one of the following three conditions: 1. reaching the set number of iterations; 2. the modulus of the parameter increment is less than the set threshold; 3. the residual between the original signal and the reconstructed signal is less than the set threshold; and at the same time, each time the relationship between the two residual sizes before and after the increment is increased is judged, the size of the μ value used next in the Levenberg-Marquardt algorithm is adjusted: if the residual after the increment is increased is less than the residual before the increment is increased, the μ value is reduced to half of the original value, that is, μ=0.5μ; otherwise, the increment is removed, the μ value is changed to twice the original value, that is, μ=2μ, and then the next iteration is entered.

[0030] The present invention has the following advantages:

[0031] This invention can help technicians identify the oil and gas content of reservoirs, thereby determining whether a specific reservoir is worth developing. By skillfully applying seismic identification and prediction technology, it helps improve the accuracy and reliability of reservoir oil and gas identification, significantly reducing oil and gas extraction costs and the risks of actual oil and gas exploration and development, providing essential support for the development and utilization of oil and gas in my country. BRIEF DESCRIPTION OF THE DRAWINGS

[0032] Figure 1This is a general flow chart of the seismic oil and gas pattern recognition technology of the present invention;

[0033] Figure 2 This is an example diagram of a certain original seismic signal used in the specific embodiment provided by the present invention;

[0034] Figure 3 Flowchart of using wavelet transform to solve frequency in the present invention;

[0035] Figure 4 A flow chart for calculating the attenuation factor using frequency for the present invention;

[0036] Figure 5 Flowchart of the method for solving instantaneous amplitude and phase using Hilbert transform in the present invention;

[0037] Figure 6 is a diagram of the reconstructed signal result of the original seismic signal;

[0038] Figure 7 is a graph of the reconstructed attenuation factor of the original seismic signal;

[0039] Figure 8 is the reconstructed frequency diagram of the original seismic signal;

[0040] Figure 9 This is a flow chart for determining the oil and gas content of the reconstruction results in the present invention. DETAILED DESCRIPTION

[0041] The technical solution of the present invention is further described below with reference to the accompanying drawings and embodiments.

[0042] like Figure 1 As shown, the present invention comprises the following steps:

[0043] Step 1: Read the earthquake data in sgy format and set the time corresponding to the data;

[0044] Step 2: Read the shear wave velocity v of the corresponding formation of the seismic data s , longitudinal wave velocity v p and density ρ, such as Figure 2 shown.

[0045] Step 3: Calculate the instantaneous frequency f of each sampling point of each seismic trace through wavelet transform, as Figure 3 shown.

[0046] Step 4: Use the shear wave velocity v s , longitudinal wave velocity v p And density ρ to find the attenuation factor α of each sampling point of each seismic channel, such as Figure 4 As shown;

[0047] Step 5: Calculate the instantaneous amplitude A and instantaneous phase of each sampling point of each seismic channel like Figure 5 As shown;

[0048] Step 6: Reconstruct the seismic data, save the reconstructed signal, reconstructed attenuation factor and reconstructed frequency, and set them as data in sgy format;

[0049] Step 7: Use Savitzky-Golay convolution smoothing algorithm to smooth the reconstruction result;

[0050] Step 8: Determine the gas content based on the abnormal information in the calculation results.

[0051] In step 6, the single channel of the original seismic data is processed and the Levenberg-Marquardt algorithm is used for iterative calculation. The instantaneous amplitude A and instantaneous phase The algorithm keeps the attenuation factor α and the instantaneous frequency f unchanged, and uses the two new parameters obtained each time, namely the attenuation factor α and the instantaneous frequency f. The above four parameters are used as initial values. According to the reconstruction signal formula in the Prony algorithm, an initial reconstructed signal value is obtained. The reconstructed signal value is continuously iterated to make it closer and closer to the original signal value. The operation is stopped after the termination condition is reached, and then the next data is calculated. After all the data are calculated, the reconstructed signal, the reconstructed attenuation factor, and the reconstructed frequency are saved.

[0052] In step 3, the method of obtaining the instantaneous frequency f of each sampling point of each seismic trace by wavelet transform is as follows: first determine the sampling frequency f of the seismic data s , sampling period T, number of sampling points N, and then use complex morlet wavelet to perform wavelet transform; suppose the length of the scale sequence is equal to the number of sampling points N, obtain the scale sequence scals, and convert the scale sequence scals into a frequency sequence F, perform wavelet transform on each seismic data with the complex morlet wavelet respectively, and determine the value in the frequency sequence F corresponding to the maximum value of the absolute value of the wavelet coefficient as the frequency size at that location, and the instantaneous frequency f corresponding to all seismic data can be calculated.

[0053] The calculation method of the attenuation factor α is as follows:

[0054]

[0055] Where, angular frequency ω=2πf,v p is the longitudinal wave velocity, and ρ is the density. First, in the second formula, the longitudinal wave velocity v p, density ρ are both known data, and the angular frequency ω can be derived from the instantaneous frequency f obtained by wavelet transform using the formula ω=2πf, so η can be derived according to the second formula; then, in the first formula, similarly, the density ρ and angular frequency ω are both known, and based on the calculated η, the attenuation factor α can be calculated.

[0056] In step 5, use Hilbert transform to find the instantaneous amplitude A and instantaneous phase Method:

[0057] Let the original signal be x(t), and its discrete Hilbert transform y(t) is defined as:

[0058]

[0059] The symbol * represents the convolution operation. Then the instantaneous amplitude is:

[0060] A(t)=(x 2 (t)+y 2 (t)) 1 / 2 ,

[0061] Instantaneous phase:

[0062]

[0063] This algorithm needs to calculate the Jacobian matrix of the currently processed seismic trace data and then calculate the μ value of the Levenberg-Marquardt algorithm; the order of this matrix is N×2, where N is the length of the seismic trace data; the reconstructed signal formula in the Prony algorithm is:

[0064] x p (t) = Ae a·t cos(2πft+θ),

[0065] Obtain the partial derivative with respect to the instantaneous frequency f and the partial derivative with respect to the absorption coefficient α:

[0066]

[0067]

[0068] Finally, the Jacobian matrix J is obtained as:

[0069]

[0070] J is the Jacobian matrix, let: N = J T J, μ is x times the maximum diagonal element of matrix N, that is, μ=x·max(N ii ), where i = 1, 2…N, and x is a positive floating point number in the range [0.0001, 0.001].

[0071] The conditions for iterative termination only need to meet one of the following three conditions: 1. The set number of iterations is reached; 2. The modulus of the parameter increment is less than the set threshold; 3. The residual between the original signal and the reconstructed signal is less than the set threshold; at the same time, the size of the μ value used next in the Levenberg-Marquardt algorithm is adjusted by judging the relationship between the two residuals before and after the increment is increased: if the residual after the increment is less than the residual before the increment is increased, the μ value is reduced to half of the original value, that is, μ = 0.5μ; otherwise, the increment is removed and the μ value is doubled, that is, μ = 2μ, and then the next iteration is entered.

[0072] like Figure 6 As shown in the figure, this figure is the reconstructed signal result diagram of the original seismic signal.

[0073] like Figure 7 As shown in the figure, this figure is the reconstructed attenuation factor diagram of the original seismic signal.

[0074] like Figure 8 As shown in the figure, this figure is the reconstructed frequency diagram of the original seismic signal.

[0075] like Figure 9 As shown in the figure, the detailed process for predicting hydrocarbon potential is to save the reconstructed signal, reconstructed attenuation factor, and reconstructed frequency as data in sgy format. This data is then compared with the original signal to determine hydrocarbon potential and the presence, location, and depth of gas. A good prediction of hydrocarbon potential can be made based on comprehensive information such as low dominant frequency, weak amplitude, chaotic reflections, low-frequency anomalies, and strong absorption anomalies.

[0076] That is, this scheme uses the Levenberg-Marquardt algorithm to iteratively reconstruct each signal of the seismic data, uses the instantaneous phase and instantaneous frequency obtained by Hilbert transform and the frequency and attenuation factor obtained by wavelet transform as initial values, and performs iterative operations on the original seismic signal data. Finally, the reconstructed signal data results and the reconstructed attenuation factor and frequency abnormal performance are used, combined with experience and professional knowledge, to judge the oil and gas content of the rock formation corresponding to the seismic data, and determine the existence, location, depth and other information of gas.

[0077] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not limiting. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present invention may be modified or replaced by equivalents without departing from the purpose and scope of the technical solutions of the present invention, which should all be included in the scope of the claims of the present invention.

Claims

1. A method for identifying reservoir oil and gas patterns based on data reconstruction, characterized by: The following steps are involved: Step 1: Read the earthquake data in sgy format and set the time corresponding to the data; Step 2: Read the shear wave velocity v of the corresponding formation of the seismic data s , longitudinal wave velocity v p and density ρ; Step 3, calculate the instantaneous frequency f of each sampling point of each seismic trace through wavelet transform; Step 4: Use the shear wave velocity v s , longitudinal wave velocity v p and density ρ to find the attenuation factor α of each sampling point of each seismic trace; Step 5: Calculate the instantaneous amplitude A and instantaneous phase φ of each sampling point of each seismic trace; Step 6: Reconstruct the seismic data, save the reconstructed signal, reconstructed attenuation factor and reconstructed frequency, and set them as data in sgy format; Step 7: Use Savitzky-Golay convolution smoothing algorithm to smooth the reconstruction result; Step 8: Determine the gas content based on the abnormal information in the calculation results; In step 6, the single channel of the original seismic data is processed and iteratively calculated using the Levenberg-Marquardt algorithm. The instantaneous amplitude A and instantaneous phase φ remain unchanged. At the same time, two new parameters, namely the attenuation factor α and the instantaneous frequency f, are obtained each time. The above four parameters are used as initial values. According to the reconstruction signal formula in the Prony algorithm, an initial reconstructed signal value is obtained. Continuous iteration makes the reconstructed signal value closer and closer to the original signal value. The calculation is stopped after the termination condition is reached, and then the next data is calculated. After all the data are calculated, the reconstructed signal, reconstructed attenuation factor, and reconstructed frequency are saved.

2. The method for identifying reservoir oil and gas patterns based on data reconstruction according to claim 1, characterized in that: In step 3, the method of obtaining the instantaneous frequency f of each sampling point of each seismic trace by wavelet transform is as follows: first determine the sampling frequency f of the seismic data s , sampling period T, number of sampling points N, and then use complex morlet wavelet to perform wavelet transform; suppose the length of the scale sequence is equal to the number of sampling points N, obtain the scale sequence scals, and convert the scale sequence scals into a frequency sequence F, perform wavelet transform on each seismic data with the complex morlet wavelet respectively, and determine the value in the frequency sequence F corresponding to the maximum value of the absolute value of the wavelet coefficient as the frequency size at that location, and the instantaneous frequency f corresponding to all seismic data can be calculated.

3. The method for identifying reservoir oil and gas patterns based on data reconstruction according to claim 1, characterized in that: The calculation method of the attenuation factor α is as follows: ， Among them, the angular frequency , v p is the longitudinal wave velocity and ρ is the density.

4. The method for identifying reservoir oil and gas patterns based on data reconstruction according to claim 1, characterized in that: This algorithm needs to calculate the Jacobian matrix for the currently processed seismic trace data, and then calculate the μ value of the Levenberg-Marquardt algorithm; the order of this matrix is ,in N is the data length of the seismic trace; the reconstructed signal formula in the Prony algorithm is: ， Obtain the partial derivative with respect to the instantaneous frequency f and the attenuation factor α Partial derivatives of : ， Finally, the Jacobian matrix J is obtained as: 。 5. The method for identifying reservoir oil and gas patterns based on data reconstruction according to claim 4, characterized in that: J is the Jacobian matrix, let: , μ is x times the largest diagonal element of matrix N, that is, ,in, , x is a positive floating point number in the interval [0.0001, 0.001].

6. The method for identifying reservoir oil and gas patterns based on data reconstruction according to claim 5, characterized in that: The conditions for iterative termination are as follows:

1. The set number of iterations is reached; 2. The modulus of the parameter increment is less than the set threshold; 3. The residual between the original signal and the reconstructed signal is less than the set threshold; at the same time, the size of the μ value used next in the Levenberg-Marquardt algorithm is adjusted by judging the relationship between the two residuals before and after the increment is increased: if the residual after the increment is less than the residual before the increment, the μ value is reduced to half of the original value, that is, Otherwise, remove the increment and change the μ value to 2 times of the original value, that is, , and then enter the next iteration.

Images