VTI medium two-way wave migration imaging method based on De Wolf approximation

Abstract

The invention discloses a VTI medium two-way wave migration imaging method based on De Wolf approximation, and the method comprises the steps: dividing a model into a plurality of thin plates, and sequentially carrying out the calculation of the thin plates from the minimum depth to the deepest depth based on a seismic source, and obtaining a forward transmission downlink migration wave field and a forward transmission primary uplink migration wave field of each thin plate; calculating forward-transmission uplink offset wave fields of the thin plates at different positions from the deepest position to the minimum depth position of the model, and adding the reserved forward-transmission primary uplink offset wave field to finally obtain forward-transmission uplink offset wave fields of all the thin plates; acquiring back propagation downlink and uplink offset wave fields according to the calculation process based on seismic records; and finally, obtaining all imaging results by using cross-correlation imaging conditions. Not only can downlink offset wave field information be utilized, but also uplink offset wave field information can be utilized; a two-way wave migration operator is formed by two one-way migration operators, the influence of anisotropic parameters on the migration imaging result can be corrected, and finally the migration imaging precision of the steep dip angle structure is effectively improved.

Description

Technical Field

[0001] The present invention relates to a geophysical exploration method, in particular to a VTI medium two-way wave migration imaging method based on DeWolf approximation. Background Art

[0002] As oil and gas exploration gradually shifts toward deep, deepwater, and unconventional areas, the anisotropy of media has attracted widespread attention from both academia and industry. Research on anisotropic media parameters is crucial for improving the accuracy of detecting ultra-deep and deep-sea water, small-scale targets, and reservoir (fluid) changes. Research has shown that most sedimentary rocks exhibit anisotropic properties. VTI (Transversely isotropic media with a vertical axis of symmetry) is a model of longitudinal anisotropy that approximates the anisotropic characteristics of actual Earth media (subsurface sedimentary strata are mostly layered and exhibit transverse isotropy). It is characterized by a vertical axis of symmetry. Compared with other anisotropic media models, wavefield propagation in VTI media is easier to describe and more realistic. Therefore, it is often used in seismic migration imaging studies in VTI media.

[0003] The DeWolf approximation splits the scattering potential into forward scattering potential and backscattering potential, and renormalizes the incident wave field and Green's function into forward scattering field and forward propagation Green's function. The forward scattering field is the sum of an infinite subsequence of all multiple forward scattering fields, and the forward propagation Green's function is the sum of similar subsequences of all multiple forward scattering corrections to the Green's function. This method is intended to solve the problem of large calculation errors in the Born approximation when the forward accumulation effect is very strong. In 2025, a scholar published a journal article titled "AVTI medium prestack migration method based on the DeWolf approximation" in the journal "Computers and Geosciences". The article used the one-way wave DeWolf approximation to realize prestack migration imaging of VTI media, effectively correcting the influence of anisotropic parameters on migration imaging. However, due to the dip angle limitation of the one-way wave operator, the imaging accuracy of steep-angle structures is relatively low.

[0004] In order to overcome this problem, how to provide a new migration imaging method that can effectively improve the imaging accuracy of steep-angle structures while ensuring computational efficiency is the research direction required by the present invention. Summary of the Invention

[0005] To address the above-mentioned problems in the prior art, the present invention provides a two-way wave migration imaging method for VTI media based on the DeWolf approximation. The method combines the downlink migration field information with the uplink migration field information, constructs a two-way wave migration operator through two one-way migration operators, and performs cross-correlation, thereby effectively improving the imaging accuracy of steep-angle structures.

[0006] To achieve the above object, the present invention adopts a technical solution: a viscoacoustic medium preprocessing iterative method one-way wave migration method based on thin plate approximation, comprising the following steps: Step 1: Arrange multiple detectors in a row at equal intervals in the required detection area to form an observation system, and arrange multiple seismic sources in a row at equal intervals and number them.

[0007] Step 2: Stimulate each earthquake source in sequence according to the number, and the observation system obtains the earthquake record when each earthquake source is excited.

[0008] Step 3: Set the velocity parameters and anisotropy parameters of the VTI medium model for the propagation of seismic scattered waves in the detection area.

[0009] Step 4: Using the idea of ​​thin plate partitioning, the VTI medium model in step 3 is divided into multiple thin plates along the depth direction, and the velocity parameters and anisotropy parameters in each thin plate are determined.

[0010] Step 5: Load the upper Ricker wavelet at the first source position and convert it into the frequency wavenumber domain using Fourier transform. Use it as the initial value for calculation to obtain the forward down-going migration wavefield and forward up-going migration wavefield of all thin plates.

[0011] Step 6: Perform Fourier transform on the seismic record corresponding to the first source excitation, convert it into the frequency-wavenumber domain, and use it as the initial value to calculate using the calculation process of step 5 to obtain the back-propagation downward migration wave field and the back-propagation upward migration wave field of all thin plates.

[0012] Step 7: Using cross-correlation imaging conditions, perform migration imaging on the forward and reverse down-going migration wave fields and up-going migration wave fields obtained in steps 5 and 6 to obtain the imaging result of the first earthquake source.

[0013] Step 8: Repeat steps 5 to 7 for each earthquake source in order of earthquake source number, thereby obtaining imaging results for each earthquake source.

[0014] Step 9: Superimpose the imaging results of all earthquake sources to finally obtain the offset imaging results of the detection area.

[0015] Furthermore, in step 1, a coordinate system is established along the arrangement direction of each detector and source as the X-axis and perpendicular to the arrangement direction (i.e., the depth direction of the detection area) as the Z-axis, with the coordinate origin at the position (0.0, 0.0); the number of detectors is nx, the distance between adjacent detectors is dr, and the position of the first detector is (0.0, 0.0); the number of sources is ns, the distance between adjacent sources is ds, and the position of the first source is (0.0, 0.0), the main frequency of the wavelet is fHz, the sampling interval is set to 1ms, and the sampling time length is nt.

[0016] Furthermore, the velocity model parameters in step 3 are ,in is the coordinate parameter, the grid parameter of the X axis is nx, the grid parameter of the Z axis is nz, and the grid spacing is dx and dz respectively; the anisotropy parameter of the VTI medium model is and The expression is: (1) in and are the maximum speed and minimum speed respectively.

[0017] Furthermore, the step five is specifically as follows: A. In the jth thin plate, the velocity parameter and anisotropy parameters and Decomposed into background parameters and disturbance parameters, where the background parameters include background velocity and background anisotropy parameters and , the disturbance parameters include velocity disturbance and anisotropic parameter perturbations and .

[0018] B. Load the Ricker wavelet at the first source location , and transform it into the frequency wavenumber domain using Fourier transform ;make is the initial value, that is .

[0019] C. Perform phase shift calculation in the frequency-wavenumber domain, including the frequency-wavenumber domain 、 、 、 .

[0020] (2) (3) (4) (5) in: (6) in for The corresponding transverse wave number, , .

[0021] D. will 、 、 and Transformed into the frequency space domain, they are 、 、 and Correction is performed and the above values ​​are substituted into formulas (7) and (8), and the forward downward offset field of all thin plates is obtained. and the positive transmission uplink offset wave field .

[0022] (7) (8) in .

[0023] E. After the calculation reaches the bottom of the model, the forward migration wave field of the second thin plate sorted from deep to shallow is As its positive upward offset wave field The initial value of is calculated according to formula (7) (when calculating the upgoing wave, the imaginary part separation i needs to be multiplied by the sign), and the previous positive upgoing migration wave field of the thin plate is added. As the initial value of the adjacent thin plate above it, the calculation continues until the top of the model, and the forward upward migration wave field of all thin plates is obtained. .

[0024] Furthermore, the calculation formula for the cross-correlation imaging in step seven is: (9) Where * is the conjugate and W is the cutoff frequency.

[0025] Compared with the traditional one-way wave migration method for VTI media based on the DeWolf approximation, the present invention first adopts the idea of ​​thin plate partitioning to divide the model parameters into thin plates. This method not only effectively improves the computational efficiency, but also improves the computational accuracy of subsequent calculations by using each thin plate. Then, based on the earthquake source, the forward downlink migration wavefield and the forward first-pass uplink migration wavefield of each thin plate are calculated in sequence from the thin plate with the smallest depth to the deepest depth. Then, the forward uplink migration wavefield is calculated for each thin plate at different positions from the deepest to the smallest depth of the model. The previous forward first-pass uplink migration wavefield is added to each thin plate as the initial value for the calculation of the next thin plate. The above steps are repeated to obtain the forward uplink migration wavefield of all thin plates. Then, based on the seismic records, the reverse downlink migration wavefield and reverse uplink migration wavefield are obtained according to the above calculation process. Finally, the complete imaging results are obtained using the cross-correlation imaging condition. The present invention can utilize not only the wavefield information of multiple downlink migrations but also the wavefield information of one uplink migration. A two-way wave migration operator is formed by two one-way migration operators, which can correct the influence of anisotropic parameters on the migration imaging results, and ultimately effectively improve the migration imaging accuracy of steep-angle structures. BRIEF DESCRIPTION OF THE DRAWINGS

[0026] Figure 1 It is the overall flow chart of the present invention.

[0027] Figure 2 are velocity model parameters and anisotropy parameters of the embodiment of the present invention.

[0028] Where (a) is the velocity model parameter; (b) is the anisotropic Parameters; (c) anisotropic parameter; Figure 3 This is the migration imaging result of the embodiment of the present invention.

[0029] Among them, (a) is the one-way wave imaging result; (b) is the two-way wave imaging result. DETAILED DESCRIPTION

[0030] The present invention will be further described below.

[0031] like Figure 1 As shown, the present invention includes the following steps: Step 1: Arrange multiple detectors in a row at equal intervals in the required detection area to form an observation system, and arrange multiple seismic sources in a row at equal intervals and number them; establish a coordinate system along the arrangement direction of each detector and seismic source as the X-axis and perpendicular to the arrangement direction as the Z-axis, with the coordinate origin at the position (0.0, 0.0); the number of detectors is nx, the distance between adjacent detectors is dr, and the position of the first detector is (0.0, 0.0); the number of seismic sources is ns, the distance between adjacent seismic sources is ds, and the position of the first seismic source is (0.0, 0.0); the main frequency of the wavelet is fHz, the sampling interval is set to 1ms, and the sampling time length is nt.

[0032] Step 2: Stimulate each source in sequence according to the number. The observation system obtains the earthquake records when each source is excited. .

[0033] Step 3: Set the velocity parameters and anisotropy parameters of the VTI medium model for the propagation of seismic scattered waves in the detection area. Specifically, the velocity model parameters are: ,in is the coordinate parameter, the grid parameter of the X axis is nx, the grid parameter of the Z axis is nz, and the grid spacing is dx and dz respectively; the anisotropy parameter of the VTI medium model is and The expression is: (1) in and are the maximum speed and minimum speed respectively.

[0034] Step 4: Using the thin plate partitioning concept, the VTI medium model in step 3 is divided into multiple thin plates along the depth direction, a total of j thin plates, with a thickness of dp, which is usually the thickness of a grid dz; 、 and Represent the velocity parameter and anisotropy parameter in a thin plate respectively. Figure 2 shown.

[0035] Step 5: Load the Ricker wavelet at the first source location (t is the time parameter), and is converted to the frequency wavenumber domain using Fourier transform ( is the circular frequency), and the forward down-transmission migration wavefield and the forward up-transmission migration wavefield of all thin plates are obtained after calculation using it as the initial value, specifically: A. In the jth thin plate, the velocity parameter and anisotropy parameters and Decomposed into background parameters and disturbance parameters, where the background parameters include background velocity and background anisotropy parameters and , the disturbance parameters include velocity disturbance and anisotropic parameter perturbations and .

[0036] B. Load the Ricker wavelet at the first source location , and transform it into the frequency wavenumber domain using Fourier transform ;make is the initial value, that is .

[0037] C. Perform phase shift calculation in the frequency-wavenumber domain, including the frequency-wavenumber domain 、 、 、 .

[0038] (2) (3) (4) (5) in: (6) in for The corresponding transverse wave number, , .

[0039] D. will 、 、 and Transformed into the frequency space domain, they are 、 、 and Correction is performed and the above values ​​are substituted into formulas (7) and (8), and the forward downward offset field of all thin plates is obtained. and the positive transmission uplink offset wave field .

[0040] (7) (8) in .

[0041] E. After the calculation reaches the bottom of the model, the forward migration wave field of the second thin plate sorted from deep to shallow is As its positive upward offset wave field The initial value of is calculated according to formula (7) (when calculating the upgoing wave, the imaginary part separation i needs to be multiplied by the sign), and the previous positive upgoing migration wave field of the thin plate is added. As the initial value of the adjacent thin plate above it, the calculation continues until the top of the model, and the forward upward migration wave field of all thin plates is obtained. .

[0042] Step 6: Perform Fourier transform on the seismic record corresponding to the first source excitation and convert it into the frequency-wavenumber domain , and use it as the initial value to calculate using the calculation process of step 5 to obtain the corresponding reverse propagation downward migration wave field and reverse uplink wave field .

[0043] Step 7: Use the cross-correlation imaging condition to perform migration imaging on the forward and reverse down-transmitted wave fields and up-transmitted wave fields obtained in steps 5 and 6 to obtain the imaging result of the first earthquake source. The specific calculation formula is: (9) Where * is the complex conjugate and W is the cutoff frequency.

[0044] Step 8: Repeat steps 5 to 7 for each earthquake source in sequence according to the earthquake source number, so as to obtain the imaging results of each earthquake source, which are 、 … .

[0045] Step 9: Superimpose the imaging results of all earthquake sources, such as Figure 3 As shown, the final offset imaging result of the detection area is obtained , the specific formula is: (10).

[0046] Figure 3 In the figure, a comparison between the conventional single-way wave imaging and the double-way wave imaging of the embodiment of the present invention shows that the imaging of the present invention is more precise.

[0047] The above is only a preferred embodiment of the present invention. It should be pointed out that for ordinary technicians in this technical field, several improvements and modifications can be made without departing from the principles of the present invention. These improvements and modifications should also be regarded as the scope of protection of the present invention.

Claims

1. A two-way migration imaging method for VTI media based on DeWolf approximation, characterized by: The following steps are involved: Step 1: Arrange multiple geophones in a row at equal intervals in the desired detection area to form an observation system, and arrange multiple seismic sources in a row at equal intervals and number them; Step 2: Stimulate each seismic source in sequence according to the number, and the observation system obtains the earthquake record when each seismic source is excited; Step 3: Set the velocity parameters and anisotropy parameters of the VTI medium model for the propagation of seismic scattered waves in the detection area; Step 4: Using the thin plate partitioning concept, the VTI medium model in step 3 is divided into multiple thin plates along the depth direction, and the velocity parameters and anisotropy parameters in each thin plate are determined; Step 5: Load the upper Ricker wavelet at the first source location and convert it into the frequency wavenumber domain using Fourier transform. Use it as the initial value to calculate and obtain the forward down-transmission migration wavefield and forward up-transmission migration wavefield of all thin plates. Step 6: Perform Fourier transform on the seismic record corresponding to the first source excitation to convert it into the frequency-wavenumber domain, and use it as the initial value to perform calculations using the calculation process of step 5 to obtain the back-propagation down-migration wavefield and the back-propagation up-migration wavefield of all thin plates; Step 7: Using cross-correlation imaging conditions, perform migration imaging on the forward and reverse down-transmitted wave fields and the up-transmitted wave fields obtained in steps 5 and 6 to obtain an imaging result of the first earthquake source; Step 8: Repeat steps 5 to 7 for each earthquake source in order of earthquake source number, thereby obtaining imaging results for each earthquake source; Step 9: Superimpose the imaging results of all earthquake sources to finally obtain the offset imaging results of the detection area.

2. The VTI medium two-way wave migration imaging method based on DeWolf approximation according to claim 1, characterized in that: In the step 1, a coordinate system is established along the X-axis and perpendicular to the Z-axis of the layout direction of each detector and source, with the origin of the coordinate system at the position (0.0, 0.0); the number of detectors is nx, the distance between adjacent detectors is dr, and the position of the first detector is (0.0, 0.0); the number of sources is ns, the distance between adjacent sources is ds, and the position of the first source is (0.0, 0.0); the main frequency of the wavelet is fHz, the sampling interval is set to 1ms, and the sampling time length is nt.

3. The VTI medium two-way wave migration imaging method based on DeWolf approximation according to claim 1, characterized in that: The velocity model parameters in step 3 are: ,in is the coordinate parameter, the grid parameter of the X axis is nx, the grid parameter of the Z axis is nz, and the grid spacing is dx and dz respectively; the anisotropy parameter of the VTI medium model is and The expression is: （1） in and are the maximum speed and minimum speed respectively.

4. The VTI medium two-way wave migration imaging method based on DeWolf approximation according to claim 1, characterized in that: The step five is specifically as follows: A. In the jth thin plate, the velocity parameter and anisotropy parameters and Decomposed into background parameters and disturbance parameters, where the background parameters include background velocity and background anisotropy parameters and , the disturbance parameters include velocity disturbance and anisotropic parameter perturbations and ; B. Load the Ricker wavelet at the first source location , and transform it into the frequency wavenumber domain using Fourier transform ,make is the initial value, that is ; C. Perform phase shift calculation in the frequency-wavenumber domain, including the frequency-wavenumber domain 、 、 、 ; （2） （3） （4） （5） in: （6） in for The corresponding transverse wave number, , ; D. will 、 、 and Transformed into the frequency space domain, they are 、 、 and Correction is performed and the above values ​​are substituted into formulas (7) and (8), and the forward downward offset field of all thin plates is obtained. and the positive transmission uplink offset wave field : （7） （8） in ; E. After the calculation reaches the bottom of the model, the forward migration wave field of the second thin plate sorted from deep to shallow is As its positive upward offset wave field The initial value of is calculated according to formula (7), and the previous positive up-transmission wave field of the thin plate is added. As the initial value of the adjacent thin plate above it, the calculation continues until the top of the model, and the forward upward migration wave field of all thin plates is obtained. .

5. The VTI medium two-way wave migration imaging method based on DeWolf approximation according to claim 1, characterized in that: The calculation formula for the cross-correlation imaging in step seven is: （9） Where * is the conjugate and W is the cutoff frequency.

Images