Gas-bearing prediction method and system based on Hilbert-Huang transform time-frequency spectrum

Abstract

The invention discloses a Hilbert-Huang transform time-frequency spectrum-based gas-bearing prediction method and system, and belongs to the technical field of oil and gas field exploration, and the method comprises the steps: carrying out the empirical mode decomposition of a to-be-analyzed segment of each seismic channel, and solving a linear steady-state signal IMF; performing Hilbert transform on each IMF component, and solving an instantaneous frequency and an envelope; synthesizing the frequency spectrums of all the IMF components to obtain a Hilbert-Huang transform two-dimensional time-frequency spectrum; extracting an instantaneous frequency spectrum at a certain seismic channel sampling point, and recording frequencies with the highest amplitude, frequencies with the second highest amplitude and corresponding amplitudes; solving the spectral area of a trapezoid formed by the two spectral lines, and dividing the spectral area by the frequency of the highest amplitude to construct a gas-containing factor; and calculating the gas-bearing factors of all sample points of the to-be-analyzed section to obtain an evaluation data body of the gas-bearing property of the reservoir. According to the method disclosed by the invention, the reliability of gas-bearing prediction is improved through Hilbert transform and construction of gas-bearing factors.

Description

Technical Field

[0001] The invention relates to the technical field of seismic interpretation for oil and gas field exploration and development, and in particular to a gas-bearing property prediction method and system based on Hilbert-Huang transform time-frequency spectrum. Background Art

[0002] In the seismic exploration and interpretation work, for continental clastic sedimentary strata, in addition to conventional structure, reservoir interpretation and prediction, the prediction of reservoir gas content is also an important task. In recent years, natural gas exploration in a certain basin has strengthened the research and review of shallow tight sandstone gas reservoirs, and further characterized the gas-bearing reservoirs in these strata to implement the exploration potential of each block. Among them, the evaluation of reservoir gas content is a difficult point and a key point. The existing seismic technology can already better characterize the plane shape, thickness and scale of the river channel, but the reliability of the prediction of gas content is unsatisfactory, which affects the formulation of well location deployment decisions. Therefore, clarifying the gas content of the reservoir in areas without well control is an important decision-making prerequisite for the deployment of subsequent well location targets.

[0003] Post-stack seismic data is the most widely used, relatively small, and fastest data type in production applications. A large number of reservoir seismic response characteristic analyses, river channel characterization, and reservoir quantitative evaluation work are based on post-stack seismic data. Existing post-stack seismic amplitude and waveform attribute technologies can already better characterize the planar morphology, thickness, and scale of the river channel, but the evaluation of the gas content of the river channel reservoir still requires the help of frequency attributes. According to a large number of theoretical and experimental studies, gas-bearing reservoirs have a strong attenuation effect on the high-frequency components of seismic waves. In the received seismic signals, reservoirs with good gas content show the characteristics of relatively strong low frequencies and weak high frequencies. Therefore, the existing gas content evaluation mainly relies on seismic signal analysis technology, as well as the analysis and transformation of frequency components, such as instantaneous frequency attributes, instantaneous Q values, frequency division attributes, two-dimensional time-frequency transformation and other technologies. The basis of instantaneous frequency technology is Hilbert transform, and then the instantaneous frequency is obtained by taking the derivative of the instantaneous phase with respect to time. However, if it is directly obtained according to the analytical signal formula, negative frequencies without physical meaning will appear, which is not conducive to frequency evaluation. The instantaneous Q value can evaluate the attenuation characteristics of the formation, but it is also based on the calculation of instantaneous frequency and cannot analyze the proportion of different frequency components at the same time. The frequency division attribute is to filter out the signals of different frequency bands of the seismic signal through filtering technology, and analyze the differences between the signals in different frequency bands, so as to evaluate the formation tuning frequency, gas content, etc., but because it is still fundamentally analyzing the amplitude and waveform characteristics, no factors directly related to gas content are found. Two-dimensional time-frequency analysis technologies such as wavelet transform and S transform can decompose non-stationary time domain signals and analyze the signal composition from the time domain and frequency domain at the same time. They are powerful tools for non-stationary signal analysis and are also commonly used basic technologies for seismic gas content prediction. However, due to the serious energy dissipation of the two-dimensional time-frequency spectrum formed by these time-frequency decompositions, the energy cluster is large, and the time and frequency resolution cannot be taken into account at the same time, it is difficult to accurately estimate the frequency composition at a certain time point, and usually cannot produce a good gas content evaluation effect. Therefore, under the existing technical means, the prediction of gas content by conventional post-stack seismic is often ineffective and can only be used as an auxiliary reference. Summary of the invention

[0004] In order to solve the estimation problem of gas content prediction of post-stack seismic signals, the present invention proposes a gas content prediction method and system based on Hilbert-Huang transform time-frequency spectrum, uses Hilbert-Huang transform to clearly separate frequency components, and constructs a gas content factor that characterizes the differences in reservoir lithology and gas content, thereby improving the persuasiveness of reservoir gas content prediction.

[0005] In order to achieve the above-mentioned object of the invention, the technical solution of the present invention is as follows:

[0006] A method for predicting gas content based on Hilbert-Huang transform time-frequency spectrum comprises the following steps:

[0007] Step 1: For the segment to be analyzed of the seismic data, perform empirical mode decomposition (EMD) on the segment to be analyzed of each seismic trace to obtain a linear steady-state signal (Intrinsic Mode Function, IMF);

[0008] Step 2: Perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope;

[0009] Step 3: synthesize the instantaneous frequencies and envelopes of all IMF components in the frequency range to be analyzed to obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f);

[0010] Step 4: Extract the instantaneous spectrum P(t i , f), record the frequency with the highest and second highest amplitude and the corresponding amplitude; calculate the spectral area of ​​the trapezoid formed by the two spectral lines and divide it by the frequency with the highest amplitude to construct the gas factor T;

[0011] Step 5: After calculating the gas-bearing factors of all sample points in the section to be analyzed, an evaluation data body of the gas-bearing property of the reservoir is obtained.

[0012] Furthermore, in step 4, the frequency range to be analyzed is the main frequency range of the seismic data, which is 1-100 Hz.

[0013] Furthermore, the method for obtaining the two-dimensional time-frequency spectrum of the Hilbert-Huang transform includes: constructing a time-frequency spectrum matrix with the horizontal axis as frequency, the vertical axis as time, and the Z axis as amplitude; locating each IMF component obtained in the time-frequency spectrum matrix according to its time and instantaneous frequency, and then summing the envelope amplitudes of all IMF components at corresponding positions.

[0014] Furthermore, for each seismic trace sample point, the extracted instantaneous spectrum P(t i , f) is the Hilbert-Huang transform spectrum P(t, f) at the earthquake sample point corresponding to t i The spectrum of the moment.

[0015] Furthermore, in step 4, according to the amplitude, find the highest amplitude F1 and the corresponding frequency f1, the second highest amplitude E2 and the corresponding frequency f2; these two frequencies form a trapezoid on the spectrum, and the area calculation method is S = (E1 + E2) * (f1-f2) / 2.

[0016] Furthermore, the seismic data in step 1 are pre-stack time migration seismic waveform data of the work area.

[0017] Furthermore, the empirical mode decomposition in step 1 is an important part of the Hilbert-Huang transform proposed by Dr. Norden E. Huang. The specific steps of performing the empirical mode decomposition on the segment to be analyzed of each seismic trace include:

[0018] Find all extreme points of the signal x(t);

[0019] Use cubic spline curve to fit the envelope lines s1(t) and s2(t) of the upper and lower extreme points, and calculate the average value of the upper and lower envelope lines m(t) = [s1(t) - s2(t)] / 2, and subtract it from x(t): h(t) = x(t) - m(t);

[0020] Determine whether h(t) is an IMF; if the number of extreme points of h(t) is not equal to the number of zero-crossing points, it is not an IMF;

[0021] If h(t) is not an IMF, replace x(t) with h(t) and repeat the above steps until h(t) meets the criterion. Then h(t) is the IMF to be extracted.

[0022] Each time a first-order IMF is obtained, it is subtracted from the original signal, and the above steps are repeated until the remaining part of the signal is r n is a monotonic sequence or a constant sequence;

[0023] In this way, the original signal x(t) is decomposed into a series of IMFs and the linear superposition of the remaining parts, where C i represents the i-th IMF component, r n Represents the residual signal after EMD decomposition:

[0024]

[0025] Furthermore, the steps of performing Hilbert transform on each IMF component include:

[0026] For a signal x(t), its Hilbert transform in the time domain can be expressed as:

[0027]

[0028] in represents the Hilbert transform of the signal x(t), that is, x(t) and the signal Convolution;

[0029] The envelope, instantaneous phase, and instantaneous frequency can be calculated according to the formula:

[0030]

[0031]

[0032]

[0033] Among them, a(t) is the signal envelope, that is, the instantaneous amplitude; φ(t) is the instantaneous phase. By taking the derivative of φ(t), we can get the instantaneous frequency f(t) of the signal x(t).

[0034] The present invention also proposes a gas content prediction system based on Hilbert-Huang transform time spectrum, comprising:

[0035] The data processing module 1 is used for performing empirical mode decomposition on the segment to be analyzed of each seismic trace for the segment to be analyzed, and obtaining the linear steady-state signal IMF;

[0036] Data processing module 2 is used to perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope;

[0037] Data processing module three is used to synthesize the instantaneous frequency and envelope of all IMF components in the frequency range to be analyzed to obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f);

[0038] Data processing module 4 is used to extract the instantaneous spectrum P(t i , f), record the frequency with the highest and second highest amplitude and the corresponding amplitude; calculate the spectral area of ​​the trapezoid formed by the two spectral lines and divide it by the frequency with the highest amplitude to construct the gas factor T;

[0039] The data processing module five is used to calculate the gas-bearing factors of all sample points in the section to be analyzed and obtain an evaluation of the gas-bearing property of the reservoir.

[0040] In summary, the present invention has the following advantages:

[0041] 1. In the present invention, the Hilbert-Huang transform is used to process the IMF components, thereby decomposing the main frequency components of the signal, clearly showing the frequency composition of the signal at a certain moment, providing a good basis for the frequency analysis related to gas content, and facilitating the analysis of frequency characteristics related to gas content;

[0042] 2. The present invention proposes and constructs a gas content factor for characterizing reservoir lithology differences and gas content differences. The gas content factor has a strong correlation with the test production conditions of actual wells, and is therefore more convincing in predicting reservoir gas content. The obtained gas content factor is composed of the spectral area of ​​the numerator and the frequency of the highest amplitude of the denominator, and can simultaneously represent the intensity of lithology differences in the reservoir section, as well as the high-frequency attenuation conditions. It has a high degree of conformity with the actual gas production conditions, and has good reliability in measuring gas content. It is a more reliable means for evaluating gas content in low-velocity sandstone reservoirs. BRIEF DESCRIPTION OF THE DRAWINGS

[0043] Figure 1 is the signal Hilbert-Huang transform t i Schematic diagram of the spectrum area at the time;

[0044] Figure 2 This is a flow chart of the gas content prediction method of the present invention;

[0045] Figure 3 The relationship between the Hilbert transform spectrum of the wellside seismic trace and the test production;

[0046] Figure 4 Comparison of gas factor profiles of high and low production wells;

[0047] Figure 5 This is the gas factor profile of high-yield wells in the same river channel. DETAILED DESCRIPTION

[0048] In order to explain the present invention more clearly, the present invention is further described below in conjunction with preferred embodiments and drawings. It should be understood by those skilled in the art that the content described below is illustrative rather than restrictive, and should not be used to limit the scope of protection of the present invention.

[0049] The technical solution for evaluating the gas content of reservoirs by using instantaneous frequency and instantaneous Q value calculates the Hilbert transform of the signal, analyzes the attenuation of the formation by estimating the instantaneous frequency, and thus characterizes the characteristics of the seismic signal related to the fluid in the formation. However, this type of method only uses the Hilbert transform, and cannot parse the frequency components of complex and non-stationary signals such as actual seismic signals. Using only a single instantaneous frequency cannot accurately characterize the true composition of the signal at a certain moment. In order to solve this problem, the present invention provides a gas content prediction method based on the Hilbert-Huang transform time-frequency spectrum, including:

[0050] Step 1: Input the seismic data volume and determine the data range and frequency range of the segment to be analyzed;

[0051] Step 2: Perform empirical mode decomposition on the segment to be analyzed of each seismic trace to obtain the linear steady-state signal IMF;

[0052] Step 3: Perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope;

[0053] Step 4: Summarize the instantaneous frequencies and envelopes of all IMF components, synthesize them within the frequency range to be analyzed, and obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum (t, f);

[0054] Step 5: For each seismic trace sample point, extract the instantaneous spectrum at the sample point, record the highest and second highest frequencies f1, f2 and the corresponding amplitudes E1, E2, the area S of the trapezoid formed by these two spectrum lines is called the spectrum area, calculate the area and divide it by the frequency of the highest amplitude to construct the gas factor: T = S / f1;

[0055] Step 6: After calculating the gas-bearing factors of all sample points in the section to be analyzed, an evaluation data body for the gas-bearing property of the reservoir is obtained.

[0056] The overall idea of ​​this method is based on the two-dimensional time-frequency decomposition method. The Hilbert-Huang transform is selected as the time-frequency decomposition method. It has the characteristic of clearly splitting the instantaneous frequency components of non-stationary signals, which is convenient for calculating the high-frequency attenuation, analyzing the instantaneous frequency components, and evaluating the reservoir lithology differences and high-frequency attenuation, thereby evaluating the gas content of the reservoir. It avoids the problem of unclear frequency composition analysis caused by the energy dissipation phenomenon of methods such as wavelet transform and S transform, and is conducive to the analysis of frequency characteristics related to gas content.

[0057] At the same time, by implementing this technical solution, we finally get the evaluation data body of gas-bearing property of the reservoir, which can directly show the favorable gas-bearing parts in the reservoir. According to these favorable parts, we can provide data basis for geological modeling or provide arguments for well location deployment. This method can be used for seismic data in time domain and depth domain, and is suitable for low-speed and low-impedance river sand body gas-bearing reservoirs.

[0058] Example 1

[0059] The following specific application examples are provided to illustrate a gas content prediction method based on the Hilbert-Huang transform time-frequency spectrum of the present invention. This embodiment is the application of the present invention to the Shaximiao Formation gas reservoir in a certain structure in a certain basin. The gas reservoir is a structural-lithological gas reservoir with river channel deposits. Several main river channels are developed in the work area. The river channel sand body is the main reservoir with low-speed and low-wave impedance characteristics. The morphology of the river can be characterized by conventional seismic amplitude attributes, but the conventional seismic frequency attributes related to the fluid cannot establish a good relationship with the actual test gas production of the drilling, and the prediction effect of the gas content of the river channel is not good. Therefore, in the case of only post-stack data, a gas content prediction method based on the seismic time-frequency spectrum area method was carried out in this work area to characterize the gas content in the main river channel.

[0060] Refer to the instruction manual Figure 2 , and its implementation process is as follows:

[0061] Step 1: Input the seismic data volume, and determine the data range and frequency range of the segment to be analyzed based on the tracking of the target layer; the seismic data volume input in this step is the pre-stack time migration seismic waveform data volume of the work area, and the seismic waveform data volume is collectively referred to as seismic data in the following text.

[0062] Step 2: Perform empirical mode decomposition on the segment to be analyzed of each seismic trace to obtain the linear steady-state signal IMF. The specific steps are as follows:

[0063] 1) Find all extreme points of the signal x(t);

[0064] 2) Use cubic spline curve to fit the envelope of the upper and lower extreme points s1(t) and s2(t), and calculate the average value of the upper and lower envelopes m(t) = [s1(t) - s2(t)] / 2, and subtract it from x(t): h(t) = x(t) - m(t);

[0065] 3) Determine whether h(t) is an IMF. If the number of extreme points of h(t) is not equal to the number of zero-crossing points, it is not an IMF;

[0066] 4) If h(t) is not an IMF, replace x(t) with h(t) and repeat the above steps until h(t) meets the criterion, then h(t) is the IMF to be extracted;

[0067] 5) Every time you get a first-order IMF, subtract it from the original signal and repeat the above steps until the remaining part of the signal is r n is a monotonic sequence or a constant sequence.

[0068] 6) In this way, the original signal x(t) is decomposed into a series of IMFs and the linear superposition of the remaining parts, where C i represents the i-th IMF component, r n Represents the residual signal after EMD decomposition:

[0069]

[0070] Step 3: Perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope;

[0071] Hilbert transform is a commonly used method in signal processing. For a signal x(t), its Hilbert transform in the time domain can be expressed as:

[0072]

[0073] in represents the Hilbert transform of the signal x(t), that is, x(t) and the signal convolution.

[0074] The envelope, instantaneous phase, and instantaneous frequency can be calculated according to the formula:

[0075]

[0076]

[0077]

[0078] Among them, a(t) is the signal envelope, that is, the instantaneous amplitude; φ(t) is the instantaneous phase. By taking the derivative of φ(t), we can get the instantaneous frequency f(t) of the signal x(t).

[0079] Step 4: Summarize the instantaneous frequencies and envelopes of all IMF components, synthesize them within the frequency range to be analyzed, and obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum (t, f); this spectrum represents the frequency component composition at each time point; in this step, the frequency range analyzed is the main frequency range of seismic data, generally 1-100 Hz;

[0080] The method for calculating the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f) is as follows: construct a time-frequency spectrum matrix with the horizontal axis as frequency, the vertical axis as time, and the Z axis as amplitude; locate each IMF component obtained in steps two to three in the time-frequency spectrum matrix according to its time and instantaneous frequency, and then sum the envelope amplitudes of all IMF components at the corresponding positions.

[0081] Step 5: For each seismic trace sample point, extract the instantaneous spectrum at the sample point, record the highest and second highest frequencies f1, f2 and the corresponding amplitudes E1, E2, the area S of the trapezoid formed by these two spectrum lines is called the spectrum area, calculate the area and divide it by the frequency of the highest amplitude to construct the gas factor: T = S / f1;

[0082] like Figure 1 As shown, for each earthquake sample point, the extracted instantaneous spectrum P(t i , f) is the Hilbert-Huang transform spectrum P(t, f) at the earthquake sample point corresponding to t i Since the number of IMFs is small, the spectrum P(t i , f) is sparse and easy to determine. According to the amplitude, find the highest amplitude E1 and the corresponding frequency f1, the second highest amplitude E2 and the corresponding frequency f2. These two frequencies can form a trapezoid on the spectrum, and the area calculation method is S = (E1 + E2) * (f1-f2) / 2, which is called the spectrum area.

[0083] like Figure 3 As shown in the figure, the construction basis of the gas factor T is: from the production situation of the actual well and the Hilbert-Huang transform two-dimensional time-frequency spectrum of the seismic trace near the well, a rule is obtained. The spectrum area S of the spectrum in the reservoir part of the high-yield gas well is larger, and the f1 value is lower. On the contrary, the spectrum area of ​​the low-yield well, dry well or non-reservoir is smaller, and the f1 value is higher. According to this rule, the spectrum area S is used as the numerator and f1 is used as the denominator to construct the gas factor.

[0084] Step 6: After calculating the gas-bearing factors of all sample points in the section to be analyzed, an assessment of the gas-bearing nature of the reservoir is obtained.

[0085] Attached Figure 4 The time domain gas factor profiles of two over-high and under-yield wells are shown. The line marked J2s2 in the profile indicates the bottom boundary of the upper sub-member of the Jurassic Shaxi Temple, and the river reservoir is developed within the range of 15ms-25ms below this interface. From the calculated gas-bearing data body, it can be seen that, firstly, the high-value area of ​​the data body can better reflect the location of the reservoir and will not deviate too much from the main body of the river channel; secondly, the gas-bearing favorable area reflected by the gas factor corresponds to the production of the well: Profile 1 is the gas factor profile of the over-high-yield well A, which has a higher gas factor (greater than 0.2) at the target layer position and a larger lateral extension; Profile 2 is the gas factor profile of the under-yield well C, with a lower gas factor value (less than 0.1).

[0086] Attached Figure 5 The difference in gas content of three wells in the same river channel is shown. The production of Well D is lower than that of Wells A and B. On the gas factor profile, it can be seen that a long section of its trajectory is located in the low gas factor section; the production of Well B is higher than that of Well A. First, it is because its horizontal section is longer, and second, the gas factor of the formation traversed by its trajectory is higher than that of Well A. The gas factor of Well A is about 0.22, and some areas of Well B are greater than 0.24. It can be seen that the gas factor calculated by this technical solution is highly consistent with the actual gas production situation, and the reliability of gas content measurement is better.

[0087] Example 2

[0088] This embodiment also proposes a gas content prediction system based on Hilbert-Huang transform time spectrum, including:

[0089] The data processing module 1 is used for performing empirical mode decomposition on the segment to be analyzed of each seismic trace for the segment to be analyzed, and obtaining the linear steady-state signal IMF;

[0090] Data processing module 2 is used to perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope;

[0091] Data processing module three is used to synthesize the instantaneous frequency and envelope of all IMF components in the frequency range to be analyzed to obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f);

[0092] Data processing module 4 is used to extract the instantaneous spectrum P(t i , f), record the frequency with the highest and second highest amplitude and the corresponding amplitude; calculate the spectral area of ​​the trapezoid formed by the two spectral lines and divide it by the frequency with the highest amplitude to construct the gas factor T;

[0093] The data processing module five is used to calculate the gas-bearing factors of all sample points in the section to be analyzed and obtain an evaluation of the gas-bearing property of the reservoir.

[0094] The functional implementation of each data processing module of the prediction system in this embodiment corresponds one-to-one to each method step in the above embodiment, and will not be repeated here.

[0095] The above description is only a preferred embodiment of the present invention and does not limit the present invention in any form. Any simple modification or equivalent change made to the above embodiment based on the technical essence of the present invention shall fall within the protection scope of the present invention.

Claims

1. A method for predicting gas content based on Hilbert-Huang transform time-frequency spectrum, characterized in that: The steps include: Step 1: For the segment to be analyzed of the seismic data, perform empirical mode decomposition on the segment to be analyzed of each seismic trace to obtain the linear steady-state signal IMF; Step 2: Perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope; Step 3: synthesize the instantaneous frequencies and envelopes of all IMF components in the frequency range to be analyzed to obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f); Step 4: Extract the instantaneous spectrum P(t i , f), record the frequency with the highest and second highest amplitude and the corresponding amplitude; calculate the spectral area of ​​the trapezoid formed by the two spectral lines and divide it by the frequency with the highest amplitude to construct the gas factor T; Step 5: After calculating the gas-bearing factors of all sample points in the section to be analyzed, an evaluation data body of the gas-bearing property of the reservoir is obtained.

2. The method according to claim 1, characterized in that In step 4, the frequency range to be analyzed is the main frequency range of the seismic data, which is 1-100 Hz.

3. The method according to claim 1, characterized in that The method for obtaining a two-dimensional time-frequency spectrum of a Hilbert-Huang transform includes: constructing a time-frequency spectrum matrix with a horizontal axis as frequency, a vertical axis as time, and a Z axis as amplitude; locating each IMF component obtained in the time-frequency spectrum matrix according to its time and instantaneous frequency, and then summing the envelope amplitudes of all IMF components at corresponding positions.

4. The method according to claim 1, characterized in that: For each seismic trace sample point, the extracted instantaneous spectrum P(t i , f) is the Hilbert-Huang transform spectrum P(t, f) at the earthquake sample point corresponding to t i The spectrum of the moment.

5. The method according to claim 1, characterized in that In step 4, according to the amplitude, find the highest amplitude E1 and the corresponding frequency f1, the second highest amplitude E2 and the corresponding frequency f2; these two frequencies form a trapezoid on the spectrum, and the area calculation method is S = (E1 + E2) * (f1-f2) / 2.

6. The method according to claim 1, characterized in that The seismic data in step 1 are pre-stack time migration seismic waveform data of the work area.

7. The method according to claim 1, characterized in that The specific steps of performing empirical mode decomposition on the segment to be analyzed of each seismic trace include: Find all extreme points of the signal x(t); Use cubic spline curve to fit the envelope lines s1(t) and s2(t) of the upper and lower extreme points, and calculate the average value of the upper and lower envelope lines m(t) = [s1(t) - s2(t)] / 2, and subtract it from x(t): h(t) = x(t) - m(t); Determine whether h(t) is an IMF; if the number of extreme points of h(t) is not equal to the number of zero-crossing points, it is not an IMF; If h(t) is not an IMF, replace x(t) with h(t) and repeat the above steps until h(t) meets the criterion. Then h(t) is the IMF to be extracted. Each time a first-order IMF is obtained, it is subtracted from the original signal, and the above steps are repeated until the remaining part of the signal is r n is a monotonic sequence or a constant sequence; In this way, the original signal x(t) is decomposed into a series of IMFs and the linear superposition of the remaining parts, where C i represents the i-th IMF component, r n Represents the residual signal after EMD decomposition:

8. The method according to claim 1, characterized in that The steps of Hilbert transforming each IMF component include: For a signal x(t), its Hilbert transform in the time domain can be expressed as: in represents the Hilbert transform of the signal x(t), that is, x(t) and the signal Convolution; The envelope, instantaneous phase, and instantaneous frequency can be calculated according to the formula: Among them, a(t) is the signal envelope, that is, the instantaneous amplitude; φ(t) is the instantaneous phase. By taking the derivative of φ(t), we can get the instantaneous frequency f(t) of the signal x(t).

9. A gas content prediction system based on Hilbert-Huang transform time-frequency spectrum, characterized in that: include: The data processing module 1 is used for performing empirical mode decomposition on the segment to be analyzed of each seismic trace for the segment to be analyzed, and obtaining the linear steady-state signal IMF; Data processing module 2 is used to perform Hilbert transform on each IMF component and obtain the instantaneous frequency and envelope; Data processing module three is used to synthesize the instantaneous frequency and envelope of all IMF components in the frequency range to be analyzed to obtain the Hilbert-Huang transform two-dimensional time-frequency spectrum P(t, f); Data processing module 4 is used to extract the instantaneous spectrum P(t i , f), record the frequency with the highest and second highest amplitude and the corresponding amplitude; calculate the spectral area of ​​the trapezoid formed by the two spectral lines and divide it by the frequency with the highest amplitude to construct the gas factor T; The data processing module five is used to calculate the gas-bearing factors of all sample points in the section to be analyzed and obtain an evaluation of the gas-bearing property of the reservoir.

Images