A method for dispersion attribute inversion based on improved time-frequency analysis technology

Abstract

This invention discloses a method for dispersion attribute inversion based on improved time-frequency analysis technology, comprising: 1. Obtaining the amplitude spectrum from each seismic record in the pre-stack seismic corner gather; 2. Performing smoothing, energy normalization, and scaling on the amplitude spectrum sequentially; 3. Constructing optimal window parameters and a new window function; 4. Replacing the window function in the conventional time-frequency transform with the new window function and obtaining the time-frequency analysis results; 5. Repeating steps 1-4 to obtain a dataset composed of multiple matrices based on all obtained time-frequency analysis results; 6. Applying the results obtained in step 5 to the P-wave dispersion attribute inversion to obtain the dispersion inversion results; 7. Repeating steps 1-6 to complete the processing of all pre-stack seismic corner gathers and obtain the P-wave dispersion attribute data volume for the entire work area. This invention achieves high-precision time-frequency analysis by adaptively constructing optimal window parameters and applying them to dispersion attribute inversion, resulting in higher accuracy of the obtained dispersion attributes.

Description

Technical Field

[0001] This invention belongs to the field of oil and gas seismic exploration data processing technology, and particularly relates to a method for inverting dispersion attributes based on improved time-frequency analysis technology. Background Technology

[0002] Traditional AVO analysis techniques are based on the Zoeppritz equation, discussing the relationship between the reflection coefficient and the P-wave and S-wave velocities and densities of the layers above and below the interface. By analyzing the influence of changes in rock elastic parameters on seismic amplitude, favorable hydrocarbon reservoir prediction can be achieved. However, fluid-related dispersion attenuation effects lead to frequency-dependent AVO (Amplitude Variation with Offset) responses in seismic waves, meaning the seismic reflection coefficient is not only related to the incident angle but also varies with frequency. To address this, Wilson et al. (2009) extended the conventional two-part AVO formula to a frequency-dependent form and proposed a frequency-varying AVO inversion theory by combining it with time-frequency analysis methods. In recent years, many scholars both domestically and internationally have conducted comprehensive and in-depth research on frequency-varying AVO inversion, and the results have demonstrated the significant importance of this technique for reservoir description and fluid identification. Conventional frequency-varying AVO inversion mainly uses the three-term AVO linear approximation proposed by Smith & Gidlow (1987) to obtain the P-wave and S-wave dispersion factors for fluid identification. However, since the dispersion parameters obtained by the AVO formula with different combinations of elastic parameters are different, frequency-varying AVO inversion based on different combinations of elastic parameters has been further developed. Furthermore, considering the sensitivity of dispersion properties to fluids, dispersion parameter inversion based on sensitivity analysis can obtain dispersion factors that are sensitive to reservoir fluids, thereby improving the accuracy of reservoir prediction.

[0003] Frequency-varying AVO inversion methods also heavily rely on time-frequency analysis techniques. The aforementioned studies mainly utilized commonly used short-time Fourier transform, S-transform or generalized S-transform, wavelet transform, and other analysis techniques. These time-frequency analysis techniques directly affect the accuracy of dispersion attribute inversion results. However, today's oil and gas exploration targets are facing more complex geological conditions. The complex underground structural conditions have brought huge challenges to reservoir prediction. The low accuracy of the aforementioned time-frequency analysis techniques makes it difficult to meet the prediction requirements of complex and thin reservoirs. Therefore, developing higher-precision reservoir prediction techniques has become a key issue in the field of oil and gas exploration.

[0004] In addition, the following techniques have been proposed in existing technologies for dispersive property inversion:

[0005] For example, document CN106291689A discloses a processing method, apparatus, and prediction system for extracting dispersion attributes from seismic data. The method includes: performing a time-frequency transform on pre-stack seismic data to obtain the time-frequency amplitude spectrum of the pre-stack seismic data; extracting seismic wavelets from the pre-stack seismic data and constructing a wavelet window function using the extracted wavelet; using the wavelet window function to perform weighted time-frequency spectrum processing on the time-frequency amplitude spectrum to obtain a processed dispersion time-frequency spectrum; and inverting the processed dispersion time-frequency spectrum to extract the dispersion attribute parameters of the pre-stack seismic data. However, this technique utilizes the conventional generalized S-transform for time-frequency analysis, and the resolution of the obtained dispersion parameters remains low, making high-precision reservoir prediction difficult.

[0006] Another document, CN110988990A, discloses a high-precision seismic attribute inversion method, which includes the following steps: S1. Converting the pre-stack common center point gather into an angle domain gather; S2. Extracting seismic wavelets from the angle domain gather; S3. Decomposing the angle domain gather using a higher-order synchronous compression transform to convert the seismic amplitude into spectral amplitude; S4. Processing the spectral amplitude in S3 using a weighting function; S5. Introducing the frequency into the approximate equation of the rewritten AVO based on the dispersive medium theory; S6. Obtaining the reflectivity dispersion result through inversion based on the equation obtained in S5. This technique is based on the short-time Fourier transform and performs dispersion inversion through time-frequency analysis using a higher-order compression transform. However, the window function of the short-time Fourier transform is fixed, making it difficult to maintain high time-frequency resolution at both low and high frequencies. Therefore, it is difficult to obtain high-precision dispersion results based on this method, which can lead to interpretation errors in practical applications. Summary of the Invention

[0007] The purpose of this invention is to overcome the above-mentioned technical problems existing in the prior art and to provide a method for inverting dispersion attributes based on improved time-frequency analysis technology. This invention performs high-precision time-frequency analysis by adaptively constructing optimal window parameters and applies them to dispersion attribute inversion, so that the obtained dispersion attributes have higher accuracy, which is beneficial to the accurate identification of reservoirs.

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

[0009] A method for inverting dispersion properties based on improved time-frequency analysis technology includes the following steps:

[0010] Step 1: Input pre-stack seismic angle gathers, transform each seismic record in the pre-stack seismic angle gathers, and obtain the amplitude spectrum;

[0011] Step 2: Perform smoothing, energy normalization, and scaling on the obtained amplitude spectrum in sequence to obtain an amplitude spectrum that meets the set range;

[0012] Step 3: Introduce parameters and construct the optimal window parameters that vary with frequency based on the results obtained in Step 2. Then, use the optimal window parameters to construct a new window function.

[0013] Step 4: Replace the window function in the conventional time-frequency transform with the newly constructed window function, and obtain the time-frequency analysis results;

[0014] Step 5: Repeat steps 1-4 until a complete pre-stack seismic angle gather has been processed, and obtain a dataset consisting of multiple matrices based on all the obtained time-frequency analysis results;

[0015] Step 6: Apply the results obtained in Step 5 to the longitudinal wave dispersion property inversion, select multiple frequencies and determine the reference frequency, and then use the least squares method to obtain high-precision dispersion inversion results;

[0016] Step 7: Repeat steps 1-6 until all pre-stack seismic angle gathers for the entire work area are processed, obtain the P-wave dispersion attribute data volume for the entire work area, and complete the inversion of dispersion attributes.

[0017] In step 1, the amplitude spectrum is obtained by performing a Fourier transform on each seismic record in the pre-stack seismic angle gather.

[0018] In step 2, the smoothed amplitude spectrum results are all greater than 0.

[0019] In step 2, the smoothing process is as follows:

[0020] Y s (f) = S[X(f)]

[0021] In the formula, X(f) is the amplitude spectrum before smoothing, and Y... s (f) is the smoothed amplitude spectrum, where f is the frequency in Hertz (Hz), and S[·] represents the smoothing process.

[0022] In step 2, the method for energy normalization is as follows:

[0023]

[0024] In the formula, This is the amplitude spectrum after energy normalization.

[0025] In step 2, the scaling method is as follows: the amplitude spectrum after energy normalization is scaled by introducing a parameter k, and the scaling formula is:

[0026]

[0027] In the formula, The amplitude spectrum is the scaled-down value, and k is an adjustable parameter, typically set to 3 or 4 for seismic signal analysis.

[0028] In step 3, the optimal window parameters that vary with frequency are constructed as follows:

[0029]

[0030] In the formula, δ(f) is the optimal window parameter, and λ and p are introduced parameters.

[0031] In step 3, the new window function is constructed as follows:

[0032]

[0033] In the formula, t represents time, w represents the constructed window function, and exp(·) represents calculating the exponent.

[0034] In step 5, the time-frequency analysis result is set to X(t,f), and the number of multi-channel seismic records in the pre-stack seismic angle gather is n. The final dataset consisting of multiple matrices is then obtained as follows:

[0035] R(t,n×f)=[X1(t,f),…,X n (t,f)].

[0036] In step 6, the method for obtaining the dispersion inversion result is as follows:

[0037]

[0038] In the formula, D p For the dispersion result of the longitudinal wave, D s The dispersion result is for shear waves. Since shear waves are not sensitive to reservoir fluids, the final result is only the dispersion result of longitudinal waves.

[0039] The expressions for vectors r and e are as follows:

[0040]

[0041]

[0042] In the formula, f0 is the reference frequency, f1, ..., f m This indicates the selection of m frequencies, v p and v s Let A and B be the velocities of the longitudinal and transverse waves, respectively, and A and B be the coefficients in the reflection coefficient formula, expressed as:

[0043]

[0044]

[0045] The advantages of using this invention are:

[0046] 1. This invention uses adaptively constructed optimal window parameters for high-precision time-frequency analysis and applies them to dispersion attribute inversion, resulting in higher accuracy of the obtained dispersion attributes, which is beneficial for accurate reservoir identification. Specifically, step 2 involves smoothing, energy normalizing, and scaling the obtained amplitude spectrum sequentially to obtain an amplitude spectrum that meets the set range. Its advantage lies in ensuring the accuracy of the adaptively calculated window parameters. Step 3 constructs optimal window parameters that vary with frequency by introducing parameters and the results obtained in step 2. This optimal window function is adaptively calculated based on actual seismic data, making the obtained window parameters more reasonable and improving accuracy. Step 4 uses the new window function to obtain time-frequency analysis results that maintain high resolution at both low and high frequencies. Step 5 repeats the aforementioned steps to obtain high-precision post-stack dispersion results corresponding to the pre-stack seismic gather. Step 6 further improves the resolution of the obtained dispersion results. Step 7 utilizes the obtained dispersion results for high-precision reservoir prediction and spatial distribution range characterization.

[0047] 2. The present invention utilizes a time-frequency analysis method with higher time-frequency focusing and resolution to perform dispersion parameter inversion, and the obtained dispersion results have higher resolution, enabling the prediction of complex reservoirs and thin reservoirs.

[0048] 3. The improved window parameter optimization time-frequency analysis technology proposed in this invention can adaptively obtain the optimal window parameter according to the characteristics of the actual seismic signal. The obtained time-frequency analysis results have a more reasonable time-frequency resolution, and the obtained dispersion parameters have higher accuracy, enabling fine characterization of the reservoir. Attached Figure Description

[0049] Figure 1 This is a flowchart of the present invention;

[0050] Figure 2 This is a comparison diagram of the synthesized signal of this invention and its different time-frequency analysis methods; wherein,

[0051] Figure (a) shows the synthesized seismic signal, which is synthesized from two 10Hz signals and two 50Hz signals, respectively;

[0052] Figure (b) shows the time spectrum obtained using the general S-transform time-frequency analysis method;

[0053] Figure (c) shows the time spectrum obtained using the sparse S-transform;

[0054] Figure (d) shows the time spectrum obtained using the generalized S-transform;

[0055] Figure (e) shows the time spectrum obtained using the method of the present invention;

[0056] The horizontal axis of Figure (a) represents amplitude, while the horizontal axis of Figures (b)-(e) represents frequency in Hertz (Hz), and the vertical axis represents time in seconds (s).

[0057] Figure 3 This is a post-stack seismic profile diagram of the present invention. The horizontal axis represents the trace number, ranging from 1 to 131, containing a total of 131 post-stack seismic traces; the vertical axis represents time, in seconds (s), ranging from 0.3 to 0.8 s; and the positions of gas-bearing reservoirs and logging curves are marked.

[0058] Figure 4 This invention uses four different time-frequency analysis methods to obtain longitudinal wave dispersion attribute profiles, where the horizontal axis represents the trace number and the vertical axis represents time, with the unit being seconds (s); and the positions of gas-bearing reservoirs and logging curves are marked.

[0059] Figure (a) shows the dispersion result obtained using the general S-transform;

[0060] Figure (b) shows the dispersion results obtained using the sparse S-transform method;

[0061] Figure (c) shows the dispersion results obtained using the generalized S-transform;

[0062] Figure (d) shows the dispersion results obtained using the method of the present invention. Detailed Implementation

[0063] Example 1

[0064] This invention discloses a method for inverting dispersion properties based on improved time-frequency analysis techniques, such as... Figure 1 As shown, it includes the following steps:

[0065] Step 1: Input pre-stack seismic angle gathers, transform each seismic record x(t) in the pre-stack seismic angle gathers, and obtain the amplitude spectrum X(f). x(t) contains N sampling points, t is time in seconds, and f is frequency in Hertz (Hz). Preferably, the amplitude spectrum is obtained by performing a Fourier transform on each seismic record in the pre-stack seismic angle gathers.

[0066] Step 2: The obtained amplitude spectrum is sequentially smoothed, energy normalized, and scaled. After smoothing, the result must be greater than 0. After scaling, an amplitude spectrum that meets the set range is obtained. Furthermore,

[0067] The smoothing process is as follows:

[0068] Y s (f) = S[X(f)]

[0069] In the formula, X(f) is the amplitude spectrum before smoothing, and Y... s (f) is the smoothed amplitude spectrum, where f is the frequency in Hertz (Hz), and S[·] represents the smoothing process.

[0070] The energy normalization process is as follows:

[0071]

[0072] In the formula, This is the amplitude spectrum after energy normalization.

[0073] The scaling method is as follows: the amplitude spectrum after energy normalization is scaled by introducing a parameter k, and the scaling formula is:

[0074]

[0075] In the formula, The amplitude spectrum is the scaled-down value, and k is an adjustable parameter, typically set to 3 or 4 for seismic signal analysis.

[0076] Step 3: Introduce parameters and construct the optimal window parameters that vary with frequency based on the results obtained in Step 2. Then, use the optimal window parameters to construct a new window function.

[0077] Furthermore, this step introduces parameters λ and p and constructs the optimal window parameters that vary with frequency based on the results obtained in step 2, as follows:

[0078]

[0079] In the formula, δ(f) is the optimal window parameter, and λ and p are introduced parameters.

[0080] Furthermore, the new window function constructed using the above formula is:

[0081]

[0082] In the formula, t represents time, w represents the constructed window function, and exp(·) represents calculating the exponent.

[0083] Step 4: Replace the window function in the conventional time-frequency transform with the newly constructed window function, and obtain the time-frequency analysis result X(t,f). The conventional time-frequency transform refers to the short-time Fourier transform, S-transform, and generalized S-transform.

[0084] Step 5: Repeat steps 1-4 until a complete pre-stack seismic angle gather is processed, and a dataset consisting of multiple matrices is obtained based on all the time-frequency analysis results.

[0085] In this step, the time-frequency analysis result is set as X(t,f), and the number of multi-channel seismic records in the pre-stack seismic angle gather is n. The final dataset consisting of multiple matrices is then obtained as follows:

[0086] R(t,n×f)=[X1(t,f),…,X n (t,f)].

[0087] Step 6: Apply the results obtained in Step 5 to the longitudinal wave dispersion property inversion, select multiple frequencies and determine the reference frequency, and then use the least squares method to obtain high-precision dispersion inversion results.

[0088] Specifically, the method for obtaining the dispersion inversion result is as follows:

[0089]

[0090] In the formula, D p For the dispersion result of the longitudinal wave, D s The dispersion result is for shear waves. Since shear waves are not sensitive to reservoir fluids, the final result is only the dispersion result of longitudinal waves.

[0091] The expressions for vectors r and e are as follows:

[0092]

[0093]

[0094] In the formula, f0 is the reference frequency, f1, ..., f m This indicates the selection of m frequencies, v p and v s Let A and B be the velocities of the longitudinal and transverse waves, respectively, and A and B be the coefficients in the reflection coefficient formula, expressed as:

[0095]

[0096]

[0097] Step 7: Repeat steps 1-6 until all pre-stack seismic angle gathers for the entire work area are processed, obtain the P-wave dispersion attribute data volume for the entire work area, and complete the inversion of dispersion attributes.

[0098] In summary, by adopting the specific technical solutions described above, this invention can perform high-precision time-frequency analysis by adaptively constructing optimal window parameters and applying them to dispersion attribute inversion, thereby obtaining more accurate dispersion attributes, which is beneficial for accurate reservoir identification.

[0099] Example 2

[0100] To more clearly illustrate the technical advantages of this invention, the invention will now be described in further detail with reference to the accompanying drawings, as follows:

[0101] Figure 2 As shown in the figure, the time-frequency spectrum obtained by the present invention (e) has higher time-frequency focusing and the highest resolution, which demonstrates the superiority of the present invention.

[0102] Figure 3 The input is a post-stack seismic profile, which contains 131 traces with a time range of 0.3 to 0.8 seconds and a time sampling interval of 0.002 seconds. The black circles in the figure mark the locations of gas-bearing reservoirs (between approximately 0.6 and 0.7 seconds). The inserted curve is the P-wave velocity curve. It can be seen that at the location of the gas-bearing reservoir, the P-wave velocity exhibits a significant low-velocity anomaly.

[0103] Figure 4 The figure shows the P-wave dispersion profiles obtained using four different time-frequency analysis methods. As can be seen, the dispersion result obtained using the S-transform (Figure a) has the lowest resolution and struggles to identify reservoir locations. The dispersion result obtained using the sparse S-transform (Figure b) shows improved resolution and can delineate reservoir locations to some extent. The dispersion result obtained using the generalized S-transform (Figure c) shows significantly improved resolution and can better identify reservoirs. The dispersion result obtained using this invention has even higher resolution and better reservoir focusing, enabling more accurate reservoir identification. This demonstrates that this invention has higher accuracy than other methods for reservoir prediction. Therefore, this technology can be used for fine-grained detection of fluid-bearing reservoirs, providing strong support for subsequent seismic interpretation and drilling deployment.

[0104] The above description is merely a specific embodiment of the present invention. Any feature disclosed in this specification may be replaced by other equivalent or similar features unless otherwise specified. All features or steps in the disclosed methods or processes may be combined in any way, except for mutually exclusive features and / or steps.

Claims

1. A method for inverting dispersion properties based on improved time-frequency analysis technology, characterized in that... Includes the following steps: Step 1: Input pre-stack seismic angle gathers, transform each seismic record in the pre-stack seismic angle gathers, and obtain the amplitude spectrum; Step 2: Perform smoothing, energy normalization, and scaling on the obtained amplitude spectrum in sequence to obtain an amplitude spectrum that meets the set range; Step 3: Introduce parameters and construct the optimal window parameters that vary with frequency based on the results obtained in Step 2. Then, use the optimal window parameters to construct a new window function. Step 4: Replace the window function in the conventional time-frequency transform with the newly constructed window function, and obtain the time-frequency analysis results; Step 5: Repeat steps 1-4 until a complete pre-stack seismic angle gather has been processed, and obtain a dataset consisting of multiple matrices based on all the obtained time-frequency analysis results; Step 6: Apply the results obtained in Step 5 to the longitudinal wave dispersion property inversion, select multiple frequencies and determine the reference frequency, and then use the least squares method to obtain high-precision dispersion inversion results; Step 7: Repeat steps 1-6 until all pre-stack seismic angle gathers in the entire work area are processed, obtain the P-wave dispersion attribute data volume of the entire work area, and complete the inversion of dispersion attributes. In step 2, the smoothing process is as follows: In the formula, The amplitude spectrum before smoothing. The amplitude spectrum after smoothing. Frequency, measured in Hertz (Hz). Indicates smoothing processing; In step 2, the method for energy normalization is as follows: ； In the formula, The amplitude spectrum after energy normalization; In step 2, the scaling process is performed by introducing parameters. The amplitude spectrum after energy normalization is scaled using the following scaling formula: In the formula, This is the scaled amplitude spectrum. As an adjustable parameter, it is typically set to 3 or 4 for seismic signal analysis; In step 3, the optimal window parameters that vary with frequency are constructed as follows: In the formula, For optimal window parameters, and For the introduced parameters; In step 3, the new window function is constructed as follows: In the formula, For time, For the constructed window function, This indicates the calculation of the exponent; In step 5, the time-frequency analysis result is set as follows: If the number of seismic records in the pre-stack seismic angle gather is n, then the final dataset consisting of multiple matrices is as follows: 。 2. The method for inverting dispersion properties based on improved time-frequency analysis technology according to claim 1, characterized in that: In step 1, the amplitude spectrum is obtained by performing a Fourier transform on each seismic record in the pre-stack seismic angle gather.

3. The method for inverting dispersion properties based on improved time-frequency analysis technology according to claim 1, characterized in that: In step 2, the smoothed amplitude spectrum results are all greater than 0.

4. The method for inverting dispersion properties based on improved time-frequency analysis technology according to claim 1, characterized in that: In step 6, the method for obtaining the dispersion inversion result is as follows: In the formula, This is the dispersion result of the longitudinal wave. The dispersion result is for shear waves. Since shear waves are not sensitive to reservoir fluids, the final result is only the dispersion result of longitudinal waves. Where, vector and The expressions are as follows: In the formula, For reference frequency, This indicates the selection of m frequencies. and These are the longitudinal wave and transverse wave velocities, respectively. and The coefficient of the reflection coefficient formula is expressed as: 。

Images