High-order spectrum-driven homomorphic filtering high-resolution seismic data processing method

Abstract

The present invention discloses a high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method. First, a seismic record is input, and a spectral simulation method is used to estimate the amplitude spectrum of the seismic wavelet from the seismic record. The third-order cumulants of the seismic record are calculated and Fourier transformed to obtain its phase spectrum. The phase spectrum of the mixed-phase wavelet at the separation point is calculated. A particle swarm algorithm is used to determine the optimal separation point between the seismic wavelet and the reflection coefficient sequence from the reciprocal spectrum sequence of the seismic record. The mixed-phase seismic wavelet is separated from the seismic record, and the seismic record is subjected to mixed-phase wavelet deconvolution to obtain an output record of the mixed-phase wavelet deconvolution. The method of the present invention reduces the complex high-dimensional functional inversion problem to a single-parameter inversion problem through the combined application of homomorphic filtering and high-order spectral transformation, thereby improving the stability and calculation accuracy of the mixed-phase wavelet estimation and the resolution of the seismic record, and enhancing the ability of the mixed-phase wavelet deconvolution method to reveal thin underground layer structures.

Description

Technical Field

[0001] The present invention belongs to the technical field of seismic data processing in oil and gas geophysical exploration, and particularly relates to a high-resolution seismic data processing method of homomorphic filtering driven by higher-order spectra. Background Art

[0002] Seismic exploration is vividly called CT scanning of the earth, and its goal is to obtain a clear image of the underground structure. The resolution of seismic exploration depends on the waveform characteristics and frequency characteristics of seismic wavelets. Deconvolution is the most commonly used technical method to improve the resolution of seismic data.

[0003] Impulse deconvolution and predictive deconvolution are the most commonly used methods to improve the resolution of seismic data. However, such methods require the assumption that the seismic wavelet is minimum-phase. When the seismic wavelet is mixed-phase, there is a strong residual phase in the seismic data after deconvolution, which seriously affects the resolution ability of the deconvolution method for underground structures.

[0004] To achieve the deconvolution processing of mixed-phase wavelets, geophysicists have proposed the homomorphic deconvolution method. This method uses a homomorphic transform to transform seismic data into a complex spectrum sequence, determines the separation point between the seismic wavelet and the reflection coefficient on the complex spectrum sequence, and performs low-pass filtering on the complex spectrum domain sequence according to the separation point, thereby obtaining the mixed-phase wavelet and realizing the deconvolution processing of the mixed-phase wavelet. However, this technology does not give an implementation method for determining the optimal separation point on the complex spectrum sequence, and the selection of the separation point position directly affects the estimation accuracy of the mixed-phase wavelet and the resolution after deconvolution.

[0005] Higher-order spectrum deconvolution is another type of mixed-phase wavelet deconvolution method. This method calculates the phase spectrum of the seismic wavelet from the higher-order spectrum of the seismic data through the coupling relationship between the higher-order spectrum of the seismic data and the phase spectrum of the seismic wavelet. However, this method involves complex high-dimensional inversion operations and phase wrapping problems, and its stability and estimation accuracy affect the actual application effect of this technology. Summary of the Invention

[0006] To solve the above technical problems, the present invention proposes a high-resolution seismic data processing method of homomorphic filtering driven by higher-order spectra, which reduces the problem of estimating mixed-phase wavelets to an inversion problem of the optimal separation point on the complex spectrum sequence, globally inverses the separation point of the mixed-phase wavelet based on the coupling relationship between the higher-order spectrum of the seismic record and the phase of the seismic wavelet, and determines the optimal separation point of the seismic wavelet on the complex spectrum sequence; through the joint application of homomorphic filtering and higher-order spectrum transformation, the complex high-dimensional functional inversion problem is reduced to a single-parameter inversion problem, improving the stability and calculation accuracy of estimating mixed-phase wavelets, and enhancing the ability of the mixed-phase wavelet deconvolution method to reveal underground thin-layer structures.

[0007] The technical solution adopted by the present invention is as follows: A high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method, and the specific steps are as follows:

[0008] Step 1: Input the seismic record, and use the spectral simulation method to estimate the amplitude spectrum of the seismic wavelet from the seismic record, and calculate the third-order cumulant of the seismic record;

[0009] Step 2: Perform Fourier transform on the third-order cumulant obtained in Step 1 to obtain its phase spectrum, and given the separation range of the seismic wavelet and the reflection coefficient in the complex spectrum domain, calculate the phase spectrum φ α (ω k ) of the mixed-phase wavelet with the separation point being α;

[0010] Step 3: Set an objective function for solving the optimal separation point, use the particle swarm algorithm to solve the optimal separation point, and calculate the phase spectrum corresponding to the optimal separation point;

[0011] Step 4: Use the inverse Fourier transform to determine the mixed-phase wavelet w m (ω) of the seismic wavelet and the phase spectrum φ α (ω k ) of the mixed-phase wavelet with the separation point being α; m (t);

[0012] Step 5: Use Wiener filtering to calculate the deconvolution operator c(t) of the mixed-phase wavelet w m (t), and use the deconvolution operator to perform deconvolution operation on the input seismic record to obtain the output record b(t) of the mixed-phase wavelet deconvolution.

[0013] Furthermore, the specific content of Step 1 is as follows:

[0014] Let d(t) represent the seismic record, and |w m (ω)| represent the amplitude spectrum of the seismic wavelet. Then, the third-order cumulant of the seismic record d(t) is calculated as follows:

[0015]

[0016] where t represents the reflection time, ω represents the angular frequency, and τ k (k = 1, 2) represents the time delay.

[0017] Furthermore, the specific content of Step 2 is as follows:

[0018] Perform Fourier transform on the third-order cumulant to obtain its phase spectrum ψ(ω1, ω2). The separation range of the seismic wavelet and the reflection coefficient in the complex spectrum domain is α1 ≤ α ≤ α2.

[0019] where ω k(k = 1, 2) represents the angular frequency, α represents the separation point, and α1, α2 represent the pre-set separation point range parameters.

[0020] Perform a homomorphic transformation on the seismic record d(t) to obtain its complex cepstrum sequence For the complex cepstrum sequence Perform low-pass filtering to obtain the complex cepstrum sequence of the seismic wavelet The calculation formula is as follows:

[0021]

[0022] Then for the complex cepstrum sequence Perform an inverse homomorphic transformation to obtain the seismic wavelet w α (t), and then perform a Fourier transform on the seismic wavelet w α (t) to obtain the phase spectrum φ α (ω k ) of the mixed-phase wavelet at the separation point α.

[0023] Furthermore, the specific steps of step 3 are as follows:

[0024] Set an objective function ε(α) for solving the optimal separation point, and the expression is as follows:

[0025]

[0026] According to the objective function ε(α), use the particle swarm optimization algorithm to solve the optimal separation point α from α1 ≤ α ≤ α2 m , and calculate the corresponding phase spectrum φ m (ω), φ α (ω k ) represents the phase spectrum of the mixed-phase wavelet at the separation point α (k = 1, 2).

[0027] Advantages of the present invention: The method of the present invention first inputs a seismic record, estimates the amplitude spectrum of the seismic wavelet from the seismic record using a spectral simulation method, calculates the third-order cumulant of the seismic record and performs a Fourier transform to obtain its phase spectrum, calculates the phase spectrum of the mixed-phase wavelet at the separation point, uses the particle swarm optimization algorithm to determine the optimal separation point between the seismic wavelet and the reflection coefficient sequence from the complex cepstrum sequence of the seismic record, separates the mixed-phase seismic wavelet from the seismic record, and performs mixed-phase wavelet deconvolution on the seismic record to obtain the output record of the mixed-phase wavelet deconvolution. The method of the present invention reduces the complex high-dimensional functional inversion problem to a single-parameter inversion problem through the combined application of homomorphic filtering and high-order spectral transformation, improves the stability and calculation accuracy of the mixed-phase wavelet estimation and the resolution of the seismic record, enhances the ability of the mixed-phase wavelet deconvolution method to reveal the underground thin-layer structure, and effectively improves the ability of seismic data to detect underground structures. Description of the Drawings

[0028] Figure 1 It is a flowchart of a high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method of the present invention.

[0029] Figure 2 It is the input seismic record diagram after processing the A block of a certain oilfield by using the existing seismic method in the embodiment of the present invention.

[0030] Figure 3 It is the result diagram after processing the seismic data of the A block of a certain oilfield by the method of the present invention in the embodiment of the present invention.

[0031] Figure 4 It is the input seismic data diagram of the B block of a certain oilfield in the embodiment of the present invention.

[0032] Figure 5 It is the result diagram after processing the seismic data of the B block of a certain oilfield by the method of the present invention in the embodiment of the present invention. Specific implementation manners

[0033] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0034] Taking the 3D seismic data of the A block of a certain oilfield as an example in Embodiment 1, the area is a thin sandstone oil and gas reservoir. To detect and characterize the depth and scope of the oil and gas reservoir, high-resolution seismic data is required.

[0035] As Figure 1 shown, a flowchart of a high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method of the present invention is as follows:

[0036] Step 1: Input the seismic record, and use the spectral simulation method to estimate the amplitude spectrum of the seismic wavelet from the seismic record, and calculate the third-order cumulant of the seismic record;

[0037] Step 2: Perform Fourier transform on the third-order cumulant obtained in Step 1 to obtain its phase spectrum, and given the separation range of the seismic wavelet and the reflection coefficient in the complex spectrum domain, calculate the phase spectrum φ α (ω k ) of the mixed-phase wavelet with the separation point α;

[0038] Step 3: Set an objective function for solving the optimal separation point, use the particle swarm optimization algorithm to solve the optimal separation point, and calculate the phase spectrum corresponding to the optimal separation point;

[0039] Step 4: Use the inverse Fourier transform to determine the mixed-phase wavelet w from the amplitude spectrum |w m (ω)| of the seismic wavelet and the phase spectrum φ α (ω k ) of the mixed-phase wavelet with the separation point α;m (t);

[0040] Step 5: Calculate the deconvolution operator c(t) of the mixed-phase wavelet w m (t), and perform deconvolution on the input seismic record using the deconvolution operator to obtain the output record b(t) of the mixed-phase wavelet deconvolution.

[0041] In this embodiment, the specific steps of step 1 are as follows:

[0042] Input the seismic record d(t) as shown in Figure 2 , |w m (ω)| represents the amplitude spectrum of the seismic wavelet. In this embodiment, the spectral simulation method is used to estimate the amplitude spectrum of the seismic wavelet, and the third-order cumulant of the seismic record d(t) The calculation formula is as follows:

[0043]

[0044] where t represents the reflection time, ω represents the angular frequency, and τ k (k = 1, 2) represents the time delay.

[0045] In this embodiment, the specific steps of step 2 are as follows:

[0046] Perform Fourier transform on the third-order cumulant to obtain its phase spectrum ψ(ω1, ω2). The separation range of the seismic wavelet and the reflection coefficient in the bispectrum domain is α1 ≤ α ≤ α2. In this embodiment, α1 = 10 ms and α2 = 50 ms.

[0047] where ω k (k = 1, 2) represents the angular frequency, α represents the separation point, and α1 and α2 represent the pre-set separation point range parameters.

[0048] Perform homomorphic transformation on the seismic record d(t) to obtain its bispectrum sequence Perform low-pass filtering on the bispectrum sequence to obtain the bispectrum sequence of the seismic wavelet The calculation formula is as follows:

[0049]

[0050] Then perform inverse homomorphic transformation on the bispectrum sequence to obtain the seismic wavelet w α (t), and then perform Fourier transform on the seismic wavelet w α (t) to obtain the phase spectrum φ α (ω k ) of the mixed-phase wavelet at the separation point α.

[0051] In this embodiment, step 3 is specifically as follows:

[0052] Set a target function ε(α) for solving the optimal separation point, and the expression is as follows:

[0053]

[0054] According to the target function ε(α), use the particle swarm algorithm to solve the optimal separation point α from α1 ≤ α ≤ α2 m = 26 ms, and calculate the corresponding phase spectrum φ m (ω) based on step 2, where φ α (ω k ) represents the phase spectrum of the mixed-phase wavelet at the separation point α (k = 1, 2).

[0055] Figure 2 represents the seismic record after being processed by the existing seismic method. Although the seismic data has a high signal-to-noise ratio, its resolution cannot well depict the formation details and sequence content. Therefore, the method of the present invention is used for resolution improvement processing. Figure 3 represents the seismic record after being processed by the method of the present invention. The resolution of the seismic data has been significantly improved, clearly revealing the interlayer details and formation interior.

[0056] As Figure 4 , Figure 5 shown, the present invention also provides Embodiment 2.

[0057] This embodiment is an application example in Block B of a certain oilfield. Figure 4 is the input seismic data of this embodiment. The resolution of the seismic data is low. Although the seismic record shows the overall outline and basic framework of the underground structure well, it cannot well distinguish the sequence relationship and formation interior between large formations. Figure 5 represents the seismic data after improving the resolution by using the method of the present invention. The resolution of the seismic record has been significantly improved, and the sequence structure and structural details between large formations are shown more clearly.

[0058] In summary, the method of the present invention reduces the complex high-dimensional functional inversion problem to a single-parameter inversion problem through the combined application of homomorphic filtering and high-order spectral transformation, improves the stability and calculation accuracy of mixed-phase wavelet estimation and the resolution of seismic records, enhances the ability of the mixed-phase wavelet deconvolution method to reveal underground thin-layer structures, and effectively improves the ability of seismic data to detect underground structures.

[0059] Those of ordinary skill in the art will realize that the embodiments described herein are provided to assist the reader in understanding the principles of the present invention, and it should be understood that the scope of protection of the present invention is not limited to such specific statements and embodiments. Those of ordinary skill in the art can make various other specific deformations and combinations that do not depart from the essence of the present invention based on the technical revelations disclosed in the present invention, and these deformations and combinations are still within the scope of protection of the present invention.

Claims

1. A high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method, the specific steps are as follows: Step 1: Input the seismic record, and use the spectral simulation method to estimate the amplitude spectrum of the seismic wavelet from the seismic record, and calculate the third-order cumulant of the seismic record; Step 2: Perform Fourier transform on the third-order cumulant obtained in step 1 to obtain its phase spectrum. Given the separation range of the seismic wavelet and the reflection coefficient in the complex spectrum domain, calculate the phase spectrum φ of the mixed phase wavelet with the separation point α. α (ω k ); Step 3: Use the particle swarm algorithm to solve the optimal separation point, and calculate the phase spectrum corresponding to the optimal separation point; The specific content of step 3 is as follows: Set an objective function ε(α) for solving the optimal separation point, and the expression is as follows: According to the objective function ε(α), the particle swarm algorithm is used to solve the optimal separation point α from α1 ≤ α ≤ α2 m , and calculate the corresponding phase spectrum φ m (ω) based on step 2, where φ α (ω k ) represents the phase spectrum of the mixed-phase wavelet at the separation point α (k = 1, 2); Among them, ψ(ω1,ω2) represents the phase spectrum obtained by the Fourier transform of the third-order cumulant in step 2; ω1,ω2 represent the angular frequencies; α1, α2 represent the pre-set separation point range parameters; Step 4: Use the inverse Fourier transform to determine the mixed-phase wavelet \(w\) m (\(\omega\)) and the phase spectrum \(\varphi\) α (\(\omega\) k ) of the mixed-phase wavelet at the separation point \(\alpha\), and determine the mixed-phase wavelet \(w\) m (t); Step 5: Calculate the deconvolution operator c(t) of the mixed-phase wavelet w(t) using Wiener filtering, and perform deconvolution operation on the input seismic record using the deconvolution operator to obtain the output record b(t) of the mixed-phase wavelet deconvolution. m (t), and perform deconvolution operation on the input seismic record using the deconvolution operator to obtain the output record b(t) of the mixed-phase wavelet deconvolution.

2. The high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method according to claim 1, wherein, The specific content of step 1 is as follows: d(t) represents the seismic record, and |w m (ω)| represents the amplitude spectrum of the seismic wavelet. Then, the third-order cumulant of the seismic record d(t) is calculated as follows: where t represents the reflection time, ω represents the angular frequency, and τ k (k = 1, 2) represents the time delay.

3. A high-order spectrum-driven homomorphic filtering high-resolution seismic data processing method according to claim 1, characterized in that The specific content of step 2 is as follows: The third-order cumulant is Fourier-transformed to obtain its phase spectrum ψ(ω1,ω2), and the separation range of the seismic wavelet and the reflection coefficient in the bispectrum domain is α1 ≤ α ≤ α2; where ω k (k = 1, 2) represents the angular frequency, α represents the separation point, and α1, α2 represent the pre-set separation point range parameters; Perform a homomorphic transformation on the seismic record d(t) to obtain its cepstrum sequence For the cepstrum sequence Perform low-pass filtering to obtain the cepstrum sequence of the seismic wavelet The calculation formula is as follows: Next, perform an anti-homomorphic transformation on the semi-final spectral sequence to obtain the seismic wavelet w α (t). Then, perform a Fourier transform on the seismic wavelet w α (t) to obtain the phase spectrum φ α (ω k ) of the mixed-phase wavelet at the separation point α.

Images