A method for viscous medium target imaging based on generalized staining algorithm

Abstract

The present invention belongs to the field of data processing technology, and specifically relates to a viscous medium target area imaging method based on a generalized staining algorithm. 1. Set the staining position; 2. Set the earthquake source, the detection point, and the seismic record information; 3. Obtain the earthquake source seismic wave field value with a time step of t; 4. Calculate the dyeing wave field value with a time step of t; 5. Calculate the real wave field and the dyeing wave field with a time step of t+1; 6. Determine whether the time loop reaches the maximum time step T, otherwise repeat 3 to 5, and end the time loop; 7. Load the detection point wave field data with a time step of T‑t; 8. Solve the detection point seismic wave field value with a time step of T‑t; 9. Calculate the detection point seismic wave field value with a time step of T‑t‑1, and determine whether it reaches the time step 0; otherwise repeat 7 to 9, and end the loop; 10. Obtain the underground structure imaging result after processing using imaging conditions.

Description

Technical Field

[0001] The invention belongs to the technical field of data processing, and in particular relates to a viscous medium target area imaging method based on a generalized staining algorithm. Background Art

[0002] Since 1929, when the U.S. Gulf Coast region successfully detected oil and gas fields using seismic exploration methods, seismic exploration has become the most widely used and indispensable method in oil and gas exploration. After decades of exploration practice, people gradually discovered that not all oil and gas reservoirs are related to anticline structures. The theory of hidden oil and gas reservoirs (Levorsen, 1964) has provided a broader prospect for oil and gas exploration and development. The focus of exploration has also shifted from shallow to deep, from conventional areas to extreme areas, and from structural oil and gas reservoirs to hidden lithologic oil and gas reservoirs. The most eye-catching of these is the subsalt area, which is prone to forming giant oil and gas reservoirs. However, the limited observation system makes it difficult to receive enough effective reflected waves, and the complex overlying structure makes the deep energy complex and weak, resulting in insufficient illumination and low signal-to-noise ratio in the subsalt area during imaging. Therefore, the development of precise imaging technology for subsalt areas has become a consensus in today's industry.

[0003] In response to this problem, the solutions proposed by domestic and foreign scholars can be roughly divided into two categories: one is mainly based on changing the observation system, receiving seismic signals in multiple directions and scales, such as three-dimensional seismic logging technology, which can greatly enrich the illumination of underground structures and is a unique method for imaging in subsalt areas (Chen, 1992). The wide azimuth three-dimensional acquisition system can collect rich azimuth information, which is conducive to imaging deep, high-steep structures, and anisotropic rock masses (Tang Yun, 2003). The other is mainly based on illumination compensation for subsalt area signals, such as using the method of weighting weak illumination areas while retaining strong illumination areas, which improves the imaging effect of subsalt areas (Gherasim, 2014). The reflection path of multiple reflection waves is different from that of primary reflection waves, so multiple waves are also used to enhance the energy signal of weakly illuminated subsalt areas (Li Peng, 2006). Multiple scattered waves are also used to image areas that lack primary reflection wave illumination (Guitton, 2002). The results of multiple wave imaging are also used to correct the results of traditional reverse time migration (Liu, 2011). It is worth noting that a wave equation migration velocity analysis method oriented to the target area is used to synthesize a data set for a specific area by using generalized Born wave field modeling. This data set contains the velocity information necessary for velocity modeling, making this data set suitable for target-oriented and fast velocity model construction (Tang, 2011). Similarly, during the propagation of seismic waves, by marking and tracking the wave field passing through the target area and applying it to migration imaging, high-quality local imaging results can also be obtained. Chen et al. (2014) introduced the idea of ​​"staining" from biology and proposed an elastic medium complex domain staining algorithm. In this algorithm, when the seismic wave contacts the staining area, the seismic wave passing through the staining area will be marked and tracked during the subsequent propagation process. Applying this algorithm can obtain a staining wave field related to the target area, and this wave field is synchronized with the real wave field during the propagation process. However, the staining wave field obtained by the complex domain staining algorithm has the disadvantage that the amplitude is much smaller than the real wave field amplitude and the waveform is distorted due to the influence of complex velocity. Li et al. (2017) proposed the control equations of the generalized coloring algorithm for elastic media. The matching degree between the constructed colored wave field and the real wave field is much higher than that of the complex domain coloring algorithm, which can achieve high-resolution imaging of the target area.

[0004] However, underground media usually have viscosity, and seismic waves will produce amplitude attenuation, phase distortion, and band narrowing during their propagation in the strata. If the viscosity of the strata is ignored, the seismic wave energy in the subsalt area will be further weakened and the waveform will be distorted, which will seriously affect the imaging accuracy of the underlying structure.

[0005] At present, there are two widely used viscous acoustic wave equations: one is the viscous wave equation based on the standard linear body model, which was first introduced by Carcione et al. (1988) into the theory of seismic wave propagation by introducing the generalized standard linear body model (GSLS), using memory variables and viscoelastic parameters to characterize the wave field characteristics; the other is the viscous wave equation based on the constant Q model, which is derived from the constant Q theory proposed by Kjartansson (1979). According to the definition of this model, within the seismic exploration frequency band (<150HZ), the quality factor Q hardly changes with frequency, which is more in line with physical laws. In addition, compared with the generalized standard linear body model, the constant Q model has fewer parameters and is more concise in calculation. Combining the above advantages, the constant Q decoupled viscoacoustic wave equation can better meet the needs of seismic exploration.

[0006] Zhu (2014) proposed a viscoacoustic equation of fractional Laplace operator based on constant Q model. This equation contains two Laplace operators corresponding to phase dispersion term and amplitude attenuation term respectively, which can be calculated by fast Fourier transform, simplifying numerical calculation and being a good method to realize viscoacoustic reverse time migration. However, this method still has many problems in its implementation. First, due to the continuous accumulation of the viscosity effect of the medium, the amplitude energy of the imaging result of the deep phase axis is extremely weak, so that the deep underground structure cannot be accurately identified. To solve this problem, Zhu (2014) directly changed the sign of the amplitude attenuation term in the equation for energy compensation, but if this compensation method is used, a large number of high-frequency outliers will be generated during the propagation of seismic waves, reducing the signal-to-noise ratio of the imaging section, so that effective geological structures cannot be identified. In order to solve the problem of high-frequency outliers in compensation, a low-pass filter needs to be introduced. However, when imaging the subsalt area, there are still problems such as insufficient illumination and low signal-to-noise ratio during imaging. Summary of the invention

[0007] In order to solve the problems described in the background technology, the present invention provides a viscous medium target area imaging method based on a generalized staining algorithm, which can further improve the recognition accuracy of complex subsalt oil and gas reservoirs and provide core technical support for the rapid and efficient development of oil and gas resources.

[0008] The technical solution adopted by the present invention is: a viscous medium target area imaging method based on a generalized staining algorithm, and the viscous medium target area imaging method based on a generalized staining algorithm comprises the following steps:

[0009] Step 1: Set the coloring position and choose the same method as the complex field coloring algorithm to mark the coloring position.

[0010] Complex field coloring algorithm equation:

[0011]

[0012]

[0013] In the equation: The original wave field represented by The dyed wave field represented by is the real velocity field, is the velocity field of the dyed region, if and only if This formula is established when When the marked scattering point coincides with the target position, the wavefront of the dyed wavefield will propagate synchronously with the original wavefield; the value of the marked area can only be 1 or 0, 1 represents the marked position, and 0 represents the unmarked position;

[0014] Step 2: Set the earthquake source, detection point, and earthquake record information;

[0015] Step 3: Use the fractional-order Laplace operator constant Q to decouple the viscoacoustic wave equation and boundary conditions to obtain the source seismic wave field value with a time step of t, and save the wave field value;

[0016] Fractional Laplace operator constant Q decoupling viscoacoustic wave equation:

[0017]

[0018]

[0019] In the equation: c is velocity, t is time, p is viscoacoustic wave field; Q is quality factor; c0 is the phase velocity defined at the reference frequency ω0; s is the earthquake source;

[0020] Step 4: Apply the wave field value calculated by the fractional-order Laplace operator constant-Q decoupled viscoacoustic wave equation with a time step of t to the generalized coloring algorithm equation of viscous media to calculate the coloring wave field value with a time step of t, and save the wave field value;

[0021] Generalized coloring algorithm equation for viscous media:

[0022]

[0023]

[0024] In the equation

[0025] Step 5: Use equation (2) and equation (3) to calculate the real wave field and the colored wave field with a time step of t+1;

[0026] Step 6: Determine whether the time loop has reached the maximum time step T, otherwise repeat steps 3 to 5, and end the time loop if it has;

[0027] Step 7: Load the wave field data of the detection point with a time step of Tt;

[0028] Step 8, using the wave equation and boundary conditions of equation (2) to solve the seismic wave field value of the detection point with a time step of Tt;

[0029] Step 9, calculate the seismic wave field value of the detection point with a time step of Tt-1, and determine whether it reaches the time step 0; otherwise, repeat steps 7 to 9, and end the loop if yes;

[0030] Step 10: Use imaging conditions to process the dyed wave fields of all shots, the back-transmitted wave field information of the detection points, and obtain the underground structure imaging results after denoising;

[0031] Imaging conditions:

[0032]

[0033] Where R is the back propagation field of the detection point, is the dyed wave field.

[0034] Furthermore, the finite difference method is used to solve the time derivative on the left side of the equal sign, and the pseudo-spectral method is used to obtain the spatial operator on the right side of the equal sign; therefore, the time accuracy is second-order accuracy, and the spatial accuracy is spectral accuracy. The calculation equation is:

[0035]

[0036]

[0037] Among them, F and F- 1 are the one-dimensional forward and negative Fourier transforms respectively, and k is the discrete wave number.

[0038] Furthermore, the fractional-order viscoacoustic equation realizes the decoupling of the phase dispersion term and the amplitude attenuation term, and compensation is achieved by changing the sign of the amplitude attenuation term during wave field propagation. However, during the compensation process, the compensation operator increases exponentially with the increase of wave number, and the high-frequency components are severely amplified, resulting in numerical instability. In order to suppress the numerical instability phenomenon, a low-pass filter is applied to achieve stable propagation.

[0039] The beneficial effects of the present invention are as follows: a method for imaging target areas of viscous media based on a generalized staining algorithm is provided, which can further improve the recognition accuracy of complex oil and gas reservoirs under salt, and provide core technical support for the rapid and efficient development of oil and gas resources. Its main advantages are as follows:

[0040] (1) Conventional complex domain coloring method has extremely low speed Due to the existence of , the amplitude of the dyed wave field is much smaller than the amplitude of the real wave field. The generalized dyeing algorithm used in this method can realize the propagation of the dyed wave field with amplitude preservation.

[0041] (2) The conventional elastic wave generalized coloring algorithm ignores the absorption and attenuation of the strata and cannot accurately describe the propagation law of seismic waves in the underground. This method uses the fractional-order Laplace operator constant Q decoupling viscoacoustic wave equation to extend the coloring algorithm to viscous media, which can more accurately describe the attenuation mechanism of seismic waves in the underground propagation process. The decoupling characteristics of the equation and the decoupling of amplitude and phase are beneficial to attenuation compensation reverse time migration.

[0042] (3) Introducing a low-pass filter into the wave equation alleviates the problem of high-frequency outliers caused by direct compensation;

[0043] (4) The reverse time migration imaging process of viscoacoustic medium target areas based on the dyeing algorithm improves the signal-to-noise ratio and resolution of imaging sections in complex structural target areas, which is conducive to accurately identifying underground structures and predicting reservoir intervals in seismic data interpretation. BRIEF DESCRIPTION OF THE DRAWINGS

[0044] Figure 1 It is a flow chart of target imaging of viscoacoustic media based on the generalized staining algorithm;

[0045] Figure 2 is the horizontal layer velocity model diagram;

[0046] Figure 3 It is the Q value model diagram of the horizontal layer;

[0047] Figure 4 It is the map of marked positions of the complex domain coloring algorithm;

[0048] Figure 5 It is the generalized coloring algorithm marking position map;

[0049] Figure 6 It is the original wavefield diagram;

[0050] Figure 7 It is the colored wave field diagram required by the viscoacoustic complex domain coloring algorithm;

[0051] Figure 8 This is a comparison chart of the original wave field and the colored wave field obtained by the viscoacoustic complex domain coloring algorithm;

[0052] Fig. 9 It is the dyed wave field diagram required by the viscoacoustic generalized dyeing algorithm;

[0053] Fig.10 It is the original wave field diagram and the dyed wave field diagram required by the viscoacoustic generalized dyeing algorithm;

[0054] Fig.11 It is a salt dome velocity model diagram;

[0055] Fig.12 It is the Q value model diagram of salt dome;

[0056] Fig.13 is the attenuation imaging diagram;

[0057] Fig.14 It is a sonic imaging picture;

[0058] Fig.15 is the compensated imaging diagram;

[0059] Fig.16 It is the imaging image of the attenuated complex domain staining algorithm;

[0060] Fig.17 It is the imaging image of the acoustic wave complex domain dyeing algorithm;

[0061] Fig.18 It is the image of compensated complex domain staining algorithm;

[0062] Fig.19 It is the attenuated generalized staining algorithm imaging image;

[0063] Fig. 20 It is the imaging image of the acoustic generalized staining algorithm;

[0064] Fig.21 It is the image of compensated generalized staining algorithm;

[0065] Fig. 22 This is a single-channel information map at a horizontal position of 2.4 km. DETAILED DESCRIPTION

[0066] Embodiment 1

[0067] Refer to the figures.

[0068] A viscous medium target area imaging method based on a generalized staining algorithm, characterized in that: the viscous medium target area imaging method based on a generalized staining algorithm comprises the following steps:

[0069] Step 1: Set the dyeing position and select the same method as the complex domain dyeing algorithm to mark the dyeing position. The complex domain dyeing algorithm equation is:

[0070]

[0071]

[0072] In the equation: The original wave field represented by The dyed wave field represented by is the real velocity field, is the velocity field of the dyed region, if and only if This formula is established when When the marked scattering point coincides with the target position, the wavefront of the dyed wavefield will propagate synchronously with the original wavefield; the value of the marked area can only be 1 or 0, 1 represents the marked position, and 0 represents the unmarked position;

[0073] Step 2: Set the earthquake source, detection point, and earthquake record information;

[0074] Step 3: Use the fractional-order Laplace operator constant Q to decouple the viscoacoustic wave equation and boundary conditions to obtain the source seismic wave field value with a time step of t, and save the wave field value;

[0075] Fractional Laplace operator constant Q decoupling viscoacoustic wave equation:

[0076]

[0077]

[0078] In the equation: c is velocity, t is time, p is viscoacoustic wave field; Q is quality factor; c0 is the phase velocity defined at the reference frequency ω0; s is the earthquake source;

[0079] Step 4: Apply the wave field value calculated by the fractional-order Laplace operator constant-Q decoupled viscoacoustic wave equation with a time step of t to the generalized coloring algorithm equation of viscous media to calculate the coloring wave field value with a time step of t, and save the wave field value;

[0080] Generalized coloring algorithm equation for viscous media:

[0081]

[0082]

[0083] In the equation

[0084] Step 5: Use equation (2) and equation (3) to calculate the real wave field and the colored wave field with a time step of t+1;

[0085] Step 6: Determine whether the time loop has reached the maximum time step T, otherwise repeat steps 3 to 5, and end the time loop if it has;

[0086] Step 7: Load the wave field data of the detection point with a time step of Tt;

[0087] Step 8, using the wave equation and boundary conditions of equation (2) to solve the seismic wave field value of the detection point with a time step of Tt;

[0088] Step 9, calculate the seismic wave field value of the detection point with a time step of Tt-1, and determine whether it reaches the time step 0; otherwise, repeat steps 7 to 9, and end the loop if yes;

[0089] Step 10: Use imaging conditions to process the dyed wave fields of all shots, the back-transmitted wave field information of the detection points, and obtain the underground structure imaging results after denoising;

[0090] Imaging conditions:

[0091]

[0092] Where R is the back propagation field of the detection point, is the dyed wave field.

[0093] The finite difference method is used to solve the time derivative on the left side of the equal sign, and the pseudo-spectral method is used to obtain the spatial operator on the right side of the equal sign; therefore, the time accuracy is second-order accuracy, and the spatial accuracy is spectral accuracy. The calculation equation is:

[0094]

[0095]

[0096] Where F and F-1 are the one-dimensional forward and negative Fourier transforms, respectively, and k is the discrete wave number.

[0097] The fractional-order viscoacoustic equation realizes the decoupling of the phase dispersion term and the amplitude attenuation term. Compensation is achieved by changing the sign of the amplitude attenuation term during wave field propagation. However, the compensation operator increases exponentially with the increase of wave number during the compensation process, and the high-frequency component is severely amplified, resulting in numerical instability. In order to suppress the numerical instability, a low-pass filter is applied to achieve stable propagation.

[0098] Due to the huge amount of calculation of the viscoacoustic equation, in order to reduce the time complexity of the algorithm and meet industrial needs, multi-GPU parallel acceleration of wave field simulation calculations is introduced, and CUFFT in the CUDA command is used to accelerate the Fourier transform of the spatial derivative.

[0099] from Fig.10 It can be concluded that the color wave field value calculated by the generalized coloring algorithm equation of viscous media has the same amplitude and phase as the original wave field value, and the color wave field value with preserved amplitude can be obtained; in practice, if no compensation is performed, Fig.13 As shown in the figure, the deep-seated imaging effect under salt is very weak. Using the acoustic wave data imaging as the reference solution, Fig.15 Although the deep energy is compensated to a certain extent after compensation, the structure is still unclear. Fig.18 The dyed wave field is due to The existence of imaging amplitude, resolution and Fig.21 There is a big gap. Fig.21 Compared to Fig.18The imaging resolution of the target area has been significantly improved, the target area structure is clear, the fault is crisp, and the breakpoints are clear, which can better achieve accurate imaging of the deep structure under the salt. Fig.14 , 15 , 18, 21 single-channel information at 2.4 km in the horizontal direction, where the blue dotted line is the first layer of the target area. Fig. 22 It can also be verified that the new staining algorithm and the generalized staining algorithm proposed in this paper can obtain imaging results with high imaging amplitude, high resolution, and correct structural positioning.

Claims

1. A method for imaging a target area of ​​a viscous medium based on a generalized staining algorithm, characterized in that: The viscous medium target imaging method based on the generalized staining algorithm includes the following steps: Step 1: Set the coloring position and choose the same method as the complex field coloring algorithm to mark the coloring position. Complex field coloring algorithm equation: In the equation: The original wave field represented by The dyed wave field represented by is the real velocity field, is the velocity field of the dyed region, if and only if This formula is established when When the marked scattering point coincides with the target position, the wavefront of the dyed wavefield will propagate synchronously with the original wavefield; the value of the marked area can only be 1 or 0, 1 represents the marked position, and 0 represents the unmarked position; Step 2: Set the earthquake source, detection point, and earthquake record information; Step 3: Use the fractional-order Laplace operator constant Q to decouple the viscoacoustic wave equation and boundary conditions to obtain the source seismic wave field value with a time step of t, and save the wave field value; Fractional Laplace operator constant Q decoupling viscoacoustic wave equation: In the equation: c is velocity, t is time, p is viscoacoustic wave field; Q is quality factor; c0 is the phase velocity defined at the reference frequency ω0; s is the earthquake source; Step 4: Apply the wave field value calculated by the fractional-order Laplace operator constant-Q decoupled viscoacoustic wave equation with a time step of t to the generalized coloring algorithm equation of viscous media to calculate the coloring wave field value with a time step of t, and save the wave field value; Generalized coloring algorithm equation for viscous media: In the equation Step 5: Use equation (2) and equation (3) to calculate the real wave field and the colored wave field with a time step of t+1; Step 6: Determine whether the time loop has reached the maximum time step T, otherwise repeat steps 3 to 5, and end the time loop if it has; Step 7: Load the wave field data of the detection point with a time step of Tt; Step 8, using the wave equation and boundary conditions of equation (2) to solve the seismic wave field value of the detection point with a time step of Tt; Step 9, calculate the seismic wave field value of the detection point with a time step of Tt-1, and determine whether it reaches the time step 0; otherwise, repeat steps 7 to 9, and end the loop if yes; Step 10: Use imaging conditions to process the dyed wave fields of all shots, the back-transmitted wave field information of the detection points, and obtain the underground structure imaging results after denoising; Imaging conditions: Where R is the back propagation field of the detection point, is the dyed wave field.

2. The method for imaging a target area of ​​a viscous medium based on a generalized staining algorithm according to claim 1, characterized in that: The finite difference method is used to solve the time derivative on the left side of the equal sign, and the pseudo-spectral method is used to obtain the spatial operator on the right side of the equal sign; therefore, the time accuracy is second-order accuracy, and the spatial accuracy is spectral accuracy. The calculation equation is: Among them, F and F -1 are the one-dimensional forward and negative Fourier transforms respectively, and k is the discrete wave number.

Images