Near-surface anisotropic first-motion wave tomography inversion modeling method and related equipment

Abstract

The invention discloses a near-surface anisotropic first-motion wave tomography inversion modeling method and related equipment, and provides a near-surface model which better conforms to actual geological conditions for pre-stack depth migration from theoretical assumption and actual inversion results which are closer to the actual geological conditions and the propagation law of seismic waves. At present, isotropic static correction and undulating surface pre-stack depth migration are promoted to anisotropic undulating surface pre-stack depth migration, the seismic data imaging quality of a double-complex exploration area is improved, the method is also a huge technical progress, the method is used for undulating surface pre-stack depth migration, a better imaging effect is achieved compared with the prior art, and the seismic data imaging quality of the double-complex exploration area is improved. According to the method, the characteristics that the signal-to-noise ratio of the first-motion waves is high and all the first-motion waves are propagated in the near-surface are utilized, the near-surface anisotropy is considered, and the accuracy and the imaging quality of the imaging position of depth migration are improved.

Description

Technical Field

[0001] The invention belongs to the field of oil and gas seismic exploration, and specifically relates to a near-surface anisotropic first-arrival wave tomography inversion modeling method and related equipment. Background Art

[0002] Figure 1 It is a near-surface topographic map of a certain exploration area. The lateral changes of strata and lithology are dramatic. There are relatively gentle exposed strata, and there are also a large number of nearly vertical or even reversed strata with large dip angles. The stratified structure with directional arrangement, especially when there are fracture fillings between the strata, will show anisotropy. The lithology parallel to the stratification direction is consistent, and the wave propagation speed is high. However, the lithology in the vertical direction of the stratification changes (such as fracture fillings or thin interlayers), and the wave propagation speed is low. The overall performance is anisotropy. Near-surface anisotropy distorts the ray propagation path, affecting the imaging effect and accuracy. Studying the near-surface velocity modeling method of anisotropy is an important way to improve the imaging quality of seismic data in complex exploration areas.

[0003] Figure 2 It shows that when the stratum is tilted (black oblique line in the figure) and anisotropy exists, the propagation of the seismic wave front is shown as the blue arc in the figure, which shows that the propagation speed is higher along the stratum direction. When anisotropy does not exist, the wave front is symmetrical, as shown by the red line in the figure. The imaging points under the two media are different, as shown by the red and blue arrows in the figure respectively.

[0004] The anisotropic characteristics of seismic wave propagation make the data show non-hyperbolic time difference on the CDP gather, which cannot be corrected by standard normal time difference analysis. Whether it is pre-stack time migration or pre-stack depth migration, if the anisotropic properties of the floating layer are ignored, the imaging position will have obvious lateral deviation, especially for complex exploration areas. In order to focus the stacked imaging, we always try to flatten the common imaging point gather. In fact, under the premise of ignoring the anisotropic effect, this flattening will cause the velocity to deviate from the correct value. When the overburden stratum dips to 45 degrees, the lateral deviation distance is the product of the thickness of the floating stratum and the ratio of normal velocity to horizontal velocity. The lateral deviation error is related to the shot offset, that is, the larger the shot offset, the larger the deviation position, which leads to blurred stack imaging. Ray tracing studies have shown that when the stratum is close to 20 degrees or 65 degrees, the blurring of the imaging is greater than that of steeper strata. If the angle of the floating formation is between 10-80 degrees, anisotropy must be considered, otherwise it will cause lateral deviation of the position and blurred imaging. Similarly, if anisotropy is ignored, the propagation speed will be greater than the actual speed, which will also cause the depth of the imaging point to be greater than the drilling result.

[0005] Figure 3Figures a and b are the pre-stack depth migration results of isotropic (left) and anisotropic (right) pre-stack depth migration in a certain work area, respectively. The green dashed line is the horizon that matches the target layer of the drilling, and the red line is the depth obtained by drilling. The target layer is about 2,500 meters. The isotropic migration is 175 meters deeper than the drilling, and the error between the anisotropic migration and the drilling is only 15 meters, within a wavelength range, and the error of the relative burial depth is less than 1%. This shows that the anisotropic characteristics of the floating layer are taken into account, which significantly reduces the depth error.

[0006] Figure 4 This is a forward modeling study of anisotropy in a complex exploration area. The figure above is a near-surface velocity model based on first-arrival wave tomography inversion. The black arrows represent the angle field determined based on the near-surface topography map of the area. Here, anisotropy is set only near the surface (above the gray solid line in the figure), with Epsilon = 0.16 and Delta = 0.04. The deep velocity of the model comes from the geological interpretation results, and the anisotropy parameters are all set to zero. Therefore, the model has anisotropy only near the surface and no anisotropy at depth. Figure 5 The depth migration results of considering and ignoring the near-surface anisotropy are shown respectively. It is not difficult to find that if the near-surface anisotropy is ignored, the imaging effect will be significantly worse for complex structures such as overthrust nappes.

[0007] Near-surface anisotropy is widely present in complex exploration areas, and the shallow surface reflection waves have a short effective arrangement and are often suppressed by surface waves and oblique interference, with a low signal-to-noise ratio. It is difficult to establish the anisotropic properties of the near-surface based on the reflection waves. Therefore, the influence of near-surface anisotropy is not considered in the actual data processing. Previous actual data and model forward modeling studies have shown that ignoring near-surface anisotropy will affect the imaging position and imaging quality of depth migration. Summary of the invention

[0008] The present invention provides a near-surface anisotropic first-arrival wave tomography inversion modeling method, which solves the problem in the prior art that near-surface anisotropy is ignored, thus affecting the imaging position and imaging quality of depth migration.

[0009] To achieve the above object, the present invention provides the following technical solutions:

[0010] A near-surface anisotropic first-arrival wave tomographic inversion modeling method, comprising:

[0011] According to the field observation system SPS file and seismic data, the observation system is set up, the first arrival time is picked, and the initial velocity model is established;

[0012] The main shallow surface reflection layer is picked up based on the isotropic deep deflection profile, and the near-surface stratum azimuth and dip field model is established based on the main shallow surface reflection layer picked up based on the isotropic deep deflection profile;

[0013] The first arrival time is subjected to azimuthal isotropic first arrival wave tomographic inversion, and the near-surface Epsilon and Delta initial models are established based on the azimuthal isotropic first arrival wave tomographic inversion results.

[0014] According to the near-surface stratum azimuth and dip field model, initial velocity model, and near-surface Epsilon and Delta initial model, the anisotropic first arrival travel time is forward modeled to obtain the theoretical first arrival time and ray path.

[0015] The anisotropic first-arrival wave tomography inversion equation is established according to the theoretical first-arrival time and ray path, and the anisotropic first-arrival wave tomography inversion equation is subjected to sparse inversion in the wavelet domain to obtain the inversion result;

[0016] The inversion results are evaluated. If the convergence conditions are met, the inversion results are output. If the convergence conditions are not met, the initial model is updated, including the velocity, Epsilon and Delta, the first arrival time and ray path are recalculated, and the inversion is performed until the convergence conditions are met and the inversion results are output.

[0017] Preferably, the azimuthally isotropic first-arrival wave tomography inversion results are specifically the velocity, main energy ray azimuth and inclination of each azimuth inversion.

[0018] Preferably, the initial velocity model is obtained based on a refracted wave delay time method.

[0019] Preferably, the theoretical first arrival time and the ray path are obtained by first obtaining the anisotropic medium wavefront time equation, which is:

[0020]

[0021] Where τ is the wavefront time, α is the phase angle, and v(α) is the phase velocity;

[0022] Thomsen obtained an approximate expression for the phase velocity under weak anisotropy conditions:

[0023] v(α)=v n (1+εsin 4 α+δsin 2 αcos 2 α)

[0024] v n is the velocity in the normal direction of the formation, δ is Delta, ε is Epsilon, and both δ and ε are Thomsen parameters;

[0025] Based on the approximate expression of phase velocity obtained by Thomsen under weak anisotropy conditions, the wavefront time equation of anisotropic media is solved to obtain the theoretical first arrival time and ray path.

[0026] Preferably, the anisotropic first-arrival wave tomographic inversion equation is established according to the theoretical first-arrival time and the ray path. Specifically, according to the previously determined ray path and theoretical first-arrival time, the length of the ray in the grid and the theoretical picking time and the calculated travel time difference are obtained to establish the anisotropic first-arrival wave tomographic inversion equation, which includes three unknowns: the formation normal velocity, Epsilon and Delta.

[0027] Preferably, the wavelet domain sparse equations are solved by an LSQR algorithm.

[0028] Preferably, the output inversion results include a normal velocity field, an Epsilon field and a Delta field.

[0029] A near-surface anisotropic first-arrival wave tomography inversion modeling system, comprising:

[0030] The first model building module: according to the field observation system SPS file and seismic data, set up the observation system, pick up the first arrival time, and build the initial velocity model;

[0031] The second model building module: picking up the main shallow surface reflection layer based on the isotropic deep deflection profile, and establishing the near-surface stratum azimuth and dip field model according to the main shallow surface reflection layer picked up by the isotropic deep deflection profile;

[0032] The third model building module: perform azimuthal isotropic first-arrival wave tomographic inversion on the first-arrival time, and establish the near-surface Epsilon and Delta initial models based on the azimuthal isotropic first-arrival wave tomographic inversion results;

[0033] Forward modeling module: forward modeling of anisotropic first arrival travel time based on near-surface formation azimuth and dip field model, initial velocity model, and near-surface Epsilon and Delta initial model to obtain theoretical first arrival time and ray path;

[0034] Solution module: establish anisotropic first-arrival wave tomography inversion equations based on theoretical first-arrival time and ray path, and perform sparse inversion in wavelet domain;

[0035] Evaluation module: Evaluate the inversion results. If the convergence conditions are met, the inversion results are output. If the convergence conditions are not met, the initial model is updated, including velocity, Epsilon and Delta, the first arrival time and ray path are recalculated, and the inversion is performed until the convergence conditions are met and the inversion results are output.

[0036] A computer device comprises a memory, a processor and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, the steps of the near-surface anisotropic first-arrival wave tomography inversion modeling method are implemented.

[0037] A computer-readable storage medium stores a computer program, which, when executed by a processor, implements the steps of a near-surface anisotropic first-arrival wave tomography inversion modeling method.

[0038] Compared with the prior art, the present invention has the following beneficial effects: the present invention provides a near-surface anisotropic first-arrival wave tomographic inversion modeling method, which is closer to the actual geological conditions and the propagation laws of seismic waves in terms of theoretical assumptions and actual inversion results, and provides a near-surface model that is more in line with the actual geological conditions for pre-stack depth migration, and advances the current "static correction" based on isotropy and pre-stack depth migration on undulating surfaces to anisotropic pre-stack depth migration on undulating surfaces, thereby improving the imaging quality of seismic data in dual-complex exploration areas. This is also a huge technological advancement, and by using it for pre-stack depth migration on undulating surfaces, better imaging effects than ever before have been achieved. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Figure 1 It is a typical near-surface topography map of the frontal zone of a certain area;

[0040] Figure 2 Schematic diagram of seismic wave front in anisotropic medium;

[0041] Figure 3 This is the second example of the effect of anisotropy on imaging effect;

[0042] Figure 4 It is a theoretical model of near-surface anisotropy;

[0043] Figure 5 Depth migration results for considering (left) and ignoring (right) near-surface anisotropy;

[0044] Figure 6 The first arrival time curve (above) and elevation (below);

[0045] Figure 7 It is a dip field established based on previous isotropic deep deflection profiles. The shallow surface structure is high and steep, with dramatic lateral changes.

[0046] Figure 8 The velocity field is obtained by tomographic inversion of isotropic first-arrival waves using actual first-arrival times.

[0047] Fig. 9 The velocity field is obtained by using the actual first arrival time to perform anisotropic first arrival wave tomographic inversion.

[0048] Fig.10 Delta field from tomographic inversion of anisotropic first-arrival waves.

[0049] Fig.11 Epsilon field from anisotropic first-arrival tomographic inversion.

[0050] Fig.12 Based on the depth migration results of the near-surface model using isotropic and anisotropic first-arrival wave tomography inversion, the continuity of the steeply dipping strata on the left wing is strengthened, and the imaging quality of the depression between the two mountains is significantly improved.

[0051] Fig.13 The flowchart of near-surface modeling of anisotropic first-arrival wave tomographic inversion includes the steps of setting up the view and picking the first-arrival, picking the main shallow surface layers of the isotropic deep-deviant profile and establishing the near-surface dip and azimuth fields, obtaining the initial values ​​of anisotropic parameters Epsilon and Delta through azimuth inversion, obtaining the first-arrival travel time and ray path through forward modeling of the anisotropic model, establishing the tomographic inversion equation group, and establishing and solving the sparse equation group in the wavelet domain.

[0052] Fig.14 The observation system and quality control examples are shown. The displayed linear correction box (red line in the figure) is consistent with the trend of the first arrival data, indicating that the observation is correct.

[0053] Fig.15 An example of first arrival picking is shown. The red dot in the figure is the picked first arrival time, which coincides with the first arrival wavelet, indicating that the picking is correct.

[0054] Fig.16 It shows how to pick reflection horizons and build shallow surface dip and azimuth fields on an isotropic deep deflection profile. Fig.16 The solid lines of different colors on the section in (a) represent the different picked horizons; Fig.16 (b) Displaying the layers picked up from different sections in three-dimensional space, the spatial distribution of each geological layer can be clearly seen; Fig.16 (c) is the calculation of the dip field based on the picked horizon, with the dip range being 0° to 90°; Fig.16 (d) is the azimuth field calculated based on the picked horizon.

[0055] Fig.17 To illustrate the azimuth, the seismic traces are intercepted according to the angle between the shot-detection line and the detection point direction. Here, they are divided into 6 different data sets: (-30°-30°), (0°-60°), (30°-90°), (60°-120°), (90°-150°) and (120°-180°);

[0056] Fig.18is the velocity field obtained by inversion in six directions;

[0057] Fig.19 It is the ray density distribution map obtained by inversion in six directions;

[0058] Fig. 20 It is the azimuth distribution diagram of the main energy rays obtained by inversion in six azimuths;

[0059] Fig.21 It is the main energy ray inclination distribution diagram obtained by inversion in six azimuths. The main energy is the main part of all rays in each grid;

[0060] Fig. 22 is the anisotropic normal velocity distribution diagram calculated based on the azimuth inversion results;

[0061] Fig.23 is the Delta field distribution map obtained from the azimuth inversion results;

[0062] Fig.24 The Epsilon field distribution map obtained from the azimuth inversion results;

[0063] Fig.25 This is an example of anisotropic first arrival travel time forward modeling. A single-layer near-surface anisotropic model is designed, where the inclination is 18°, the azimuth is 20°, Epsilon is equal to 0.25, and Delta is equal to 0.1. (a) is the gradually varying normal velocity field, (b) is the anisotropic wavefront time field, and (c) is a schematic diagram of ray tracing time shift;

[0064] Fig.26 The invention discloses a near-surface anisotropic first-arrival wave tomography inversion modeling system. DETAILED DESCRIPTION

[0065] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Generally, the components of the embodiments of the present invention described and shown in the drawings here can be arranged and designed in various different configurations.

[0066] Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention claimed for protection, but merely represents selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0067] It should be noted that similar reference numerals and letters denote similar items in the following drawings, and therefore, once an item is defined in one drawing, further definition and explanation thereof is not required in subsequent drawings.

[0068] In the description of the embodiments of the present invention, it should be noted that if the terms "upper", "lower", "horizontal", "inner", etc. indicate an orientation or positional relationship based on the orientation or positional relationship shown in the drawings, or the orientation or positional relationship in which the product of the invention is usually placed when in use, it is only for the convenience of describing the present invention and simplifying the description, and does not indicate or imply that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present invention. In addition, the terms "first", "second", etc. are only used to distinguish the description, and cannot be understood as indicating or implying relative importance.

[0069] In addition, if the term "horizontal" appears, it does not mean that the component must be absolutely horizontal, but can be slightly tilted. For example, "horizontal" only means that its direction is more horizontal than "vertical", which does not mean that the structure must be completely horizontal, but can be slightly tilted.

[0070] In the description of the embodiments of the present invention, it is also necessary to explain that, unless otherwise clearly specified and limited, the terms "set", "install", "connect", and "connect" should be understood in a broad sense, for example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be a direct connection, or it can be indirectly connected through an intermediate medium, or it can be the internal connection of two components. For ordinary technicians in this field, the specific meanings of the above terms in the present invention can be understood according to specific circumstances.

[0071] In order to enable those skilled in the art to better understand the technical solution of the present invention, the present invention will be further described in detail below with reference to the accompanying drawings.

[0072] like Fig.13 As shown, a near-surface anisotropic first-arrival wave tomography inversion modeling system comprises the following steps:

[0073] S1: Set up the observation system and pick the first arrival of seismic data. According to the given field observation system SPS file and seismic data, set up the observation system, perform quality control and inspection, and then accurately pick the first arrival time.

[0074] S2: Pick up the main shallow surface layers of the isotropic depth migration section. The anisotropic processing is based on the isotropic processing results. The early isotropic deep deviation section can be used to more accurately obtain the angle field information of the formation interface. First, pick up the layers on each Inline section, and each layer forms a surface in space.

[0075] S3: Based on the multiple geological layer "surfaces" picked up in the previous step, the normal direction of each grid point is calculated according to the following formula, and then the inclination and azimuth are calculated. The angle between the normal direction and the Z axis of the model is defined as the inclination, and the angle between the normal direction and the north direction of the geodetic coordinate (i.e., the Y coordinate) is defined as the azimuth.

[0076] Suppose P0(x0,y0,z0) is a known point on the stratum. To obtain the dip and azimuth of the stratum at point P0, take another n known points P1(x1,y1,z1), P2(x2,y2,z2), ..., P0 at different directions near P0. n (x n ,y n ,z n ), find the value of P0, P1, ..., P n The plane where the sum of the squares of the distances is the smallest is found, and the inclination and azimuth of this plane are taken as the inclination and azimuth of the formation at point P0. The calculation steps are as follows:

[0077] (1) Input P0, P1, ..., P n The coordinates of the point are x0, y0, z0; x1, y1, z1; ...; x n ,y n ,z n .

[0078] (2) Calculate a according to the following formula ij (i=1,2,3;j=1,2,3),d j (j=1,2,3)

[0079]

[0080]

[0081]

[0082] (3) Calculate D, D1, and D2 according to the following formulas.

[0083]

[0084]

[0085]

[0086] (4) Calculate c1 and c2 according to the following formula.

[0087] c1=D1 / D, c2=D2 / D.

[0088] (5) Calculate the azimuth of the formation at point P0 by the following formula:

[0089]

[0090]

[0091] (6) Press to calculate the inclination angle θ of the formation at point P0.

[0092]

[0093] S4: The first arrival time is divided into azimuths, and azimuthally isotropic first arrival wave tomographic inversion is performed to obtain the velocity, azimuth angle and inclination angle of the main energy ray for each azimuth inversion result.

[0094] S5: Due to the influence of anisotropy, the results of different azimuthal inversions are different. The present invention backcalculates Epsilon and Delta based on these differences. Theoretically, the isotropic tomographic inversion still performed in the azimuthal inversion here is inconsistent with the actual geological conditions. Therefore, it is used as the initial field of the anisotropic attribute parameters Epsilon and Delta.

[0095] (1) Assume that there is a ray I passing through a near-surface imaging unit, and the inclination angle θ of this ray is known: i (i=1,2,L,I), azimuth and the group velocity V on each ray i (i=1,2,L,I), the corresponding group slowness is Assume that the imaging unit is a TTI medium, and the inclination and azimuth of its symmetry axis are θ0 and The angles between each ray and the axis of symmetry (i.e., group angles) are φ i (i=1,2,L,I). Establish the following system of equations about the parameters a1, a2 and a3.

[0096]

[0097] in,

[0098]

[0099] (2) Use the least squares method to solve the above equations and obtain the parameters a1, a2 and a3.

[0100] (3) Solving the following equations, we can obtain the velocity V in the direction of the symmetry axis. P⊥ And Thomsen parameters Epsilon and Delta.

[0101]

[0102] S6: Under the anisotropic model, the propagation velocity of seismic waves varies with the angle between the velocity and the normal direction of the stratum. The velocity here refers to the propagation velocity perpendicular to the interface, that is, the velocity in the direction of the normal direction of the interface. Here, the initial velocity model is still established based on the refraction wave delay time method.

[0103] S7: Forward modeling of the anisotropic first arrival travel time to obtain the theoretical first arrival time.

[0104] The wavefront time equation for anisotropic media is:

[0105]

[0106] Where τ is the wavefront time, α is the phase angle, and v(α) is the phase velocity. Thomsen obtained an approximate expression for the phase velocity under weak anisotropy conditions:

[0107] v(α)=v n (1+εsin 4 α+δsin 2 αcos 2 α)

[0108] v n is the velocity in the normal direction of the formation, δ and ε are Thomsen parameters.

[0109] Solve equation (1) to obtain the first arrival time and ray path.

[0110] S8: According to the previously determined ray path and first arrival time, the length of the ray in the grid Δl and the picking time and calculated travel time difference Δt are obtained, and the anisotropic first arrival wave tomography inversion equation is established, which includes three unknowns: formation normal velocity, Epsilon and Delta.

[0111]

[0112] in:

[0113]

[0114] Where subscript i is the ray number and j is the grid number.

[0115] S9: The sparse inversion equation group in wavelet domain is established, and the tomographic inversion equation is subjected to wavelet transformation to form a sparse equation group.

[0116] Let I, J, K be the number of grids in the z, y and x directions. Any ray satisfies the following equation:

[0117]

[0118] For example, if the grid slowness Δs is transformed by orthogonal wavelet in the x and y directions, it can be decomposed into the smoothness coefficient x and the detail coefficient d. Since the detail coefficient is approximately 0, Δs can be approximately expressed by x as:

[0119]

[0120] where h x and h y is the reconstruction filter in the x and y directions. Substituting Δs in the equation with the above formula, we obtain the following equation with the wavelet coefficient x as the unknown number:

[0121]

[0122] In the above formula, the dimension of the wavelet coefficient x is N=I*J*K / 4, which is one-fourth of the dimension I*J*K of Δs. The coefficient Δl of the original equation is decomposed by wavelet to obtain the coefficient of the new equation, and the unknown number of the new equation is one-tenth of the original equation.

[0123] S10: LSQR least squares QR decomposition to solve the inverse equation.

[0124] Since the equations of tomographic inversion are all ill-conditioned and have no exact solution, we can only find its optimal solution. Different optimization methods will result in different solutions. The LSQR algorithm takes up less memory, has moderate speed, fast convergence and stable results. This method is often used in conventional tomographic inversion.

[0125] The LSQR algorithm is an algorithm for calculating large sparse linear equations, proposed by Paige and Saunders in 1982. The method is mainly used to solve the following linear equations (A is an m*n matrix, m>n) and keep the second-order residual norm as small as possible:

[0126] Ax=b

[0127]

[0128] in

[0129] A={c ij} M×N

[0130] b=[Δt1,Δt2,...,Δt M ] T

[0131] x=[x1,x2,...,x N ] T

[0132] For traditional least squares problems, the solution to the equation can be directly found, but for large sparse matrices, it is often an ill-conditioned problem (overdetermined or underdetermined system of equations). Therefore, the LSQR method does not directly solve x, but adheres to the principle of least squares and seeks the best solution in the space defined by the normal equation.

[0133] S11: There are two evaluation methods for inversion result evaluation: one is that the iterative convergence difference is less than a given threshold, and the other is that the convergence difference of the last few iterations tends to be stable.

[0134] Based on the velocity field, Epsilon field and Delta field generated by each iterative inversion, anisotropic first arrival travel time forward modeling is performed to obtain the theoretical synthetic time. The mean square error between this set of time and the actual first arrival time is called the iterative convergence difference. When the iterative convergence difference is less than the given threshold, the inversion result can be considered reliable and the inversion result is output. Due to different iterative convergence differences in different work areas, it is difficult to give a reasonable threshold in actual production. Another evaluation method is to consider that the iterative convergence difference is stable, that is, as the number of iterations increases, the iterative convergence difference no longer decreases, and the inversion result is considered reasonable and the inversion result is output.

[0135] S12: If the iterative convergence conditions are met, the inversion results are output, including the normal velocity field, Epsilon field and Delta field; if the iterative convergence conditions are not met, the model is updated and S6 to S10 are repeated until the convergence conditions are met.

[0136] The present invention is closer to the actual geological conditions and the propagation laws of seismic waves in terms of theoretical assumptions and actual inversion results, providing a near-surface model that is more in line with the actual geological conditions for pre-stack depth migration, and advancing the current "static correction" based on isotropy and pre-stack depth migration on undulating surfaces to pre-stack depth migration on anisotropic undulating surfaces, thereby improving the imaging quality of seismic data in dual-complex exploration areas. This is also a huge technological advancement. By using it for pre-stack depth migration on undulating surfaces, better imaging effects have been achieved than ever before. The application effect of the two-dimensional survey line in a mountainous area in Xiqiu is demonstrated here.

[0137] Figure 6 It is the first arrival time and topography of the two-dimensional survey line, with an elevation difference of 700 meters; Figure 7 It is a dip field established based on previous isotropic deep deflection profiles, with high and steep shallow surface structures and dramatic lateral variations; Figure 8 and Fig. 9 These are the velocity models obtained by isotropic and anisotropic tomographic inversion. It is not difficult to find that the shallow surface velocity law of anisotropic inversion is more reasonable. Fig.10 and 11 They are the anisotropic attribute parameters Epsilon and Delta values ​​obtained by inversion; Fig.12The results of prestack depth migration obtained by using isotropic inversion and anisotropic inversion velocity models respectively were demonstrated. The continuity of the steeply dipped strata on the left wing was strengthened, and the imaging quality of the depression between the two mountains was significantly improved.

[0138] like Fig.26 The present invention also provides a near-surface anisotropic first-arrival wave tomography inversion modeling system, comprising:

[0139] The first model building module: according to the field observation system SPS file and seismic data, set up the observation system, pick up the first arrival time, and build the initial velocity model;

[0140] The second model building module: picking up the main shallow surface reflection layer based on the isotropic deep deflection profile, and establishing the near-surface stratum azimuth and dip field model according to the main shallow surface reflection layer picked up by the isotropic deep deflection profile;

[0141] The third model building module: perform azimuthal isotropic first-arrival wave tomographic inversion on the first-arrival time, and establish the near-surface Epsilon and Delta initial models based on the azimuthal isotropic first-arrival wave tomographic inversion results;

[0142] Forward modeling module: forward modeling of anisotropic first arrival travel time based on near-surface formation azimuth and dip field model, initial velocity model, and near-surface Epsilon and Delta initial model to obtain theoretical first arrival time and ray path;

[0143] Solution module: establish anisotropic first-arrival wave tomography inversion equation according to theoretical first-arrival time and ray path, perform wavelet transform on the anisotropic first-arrival wave tomography inversion equation to obtain a sparse equation group in the wavelet domain, and then solve the equation group;

[0144] Evaluation module: Evaluate the inversion results. If the convergence conditions are met, the inversion results are output. If the convergence conditions are not met, the initial model (including velocity, Epsilon and Delta) is updated, the first arrival time and ray path are recalculated, and the inversion is performed until the convergence conditions are met and the inversion results are output.

[0145] A terminal device is provided in one embodiment of the present invention. The terminal device of this embodiment includes: a processor, a memory, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, the steps in the above-mentioned method embodiments are implemented. Alternatively, when the processor executes the computer program, the functions of the modules / units in the above-mentioned device embodiments are implemented.

[0146] The computer program may be divided into one or more modules / units, and the one or more modules / units are stored in the memory and executed by the processor to accomplish the present invention.

[0147] The terminal device may be a computing device such as a desktop computer, a notebook, a PDA, a cloud server, etc. The terminal device may include, but is not limited to, a processor and a memory.

[0148] The processor can be a central processing unit (CPU), or other general-purpose processors, digital signal processors (DSP), application-specific integrated circuits (ASIC), field-programmable gate arrays (FPGA) or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc.

[0149] The memory may be used to store the computer programs and / or modules, and the processor implements various functions of the terminal device by running or executing the computer programs and / or modules stored in the memory and calling the data stored in the memory.

[0150] If the module / unit integrated in the terminal device is implemented in the form of a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the present invention implements all or part of the processes in the above-mentioned embodiment method, and can also be completed by instructing the relevant hardware through a computer program. The computer program can be stored in a computer-readable storage medium, and the computer program can implement the steps of the above-mentioned various method embodiments when executed by the processor. Among them, the computer program includes computer program code, and the computer program code can be in source code form, object code form, executable file or some intermediate form. The computer-readable medium may include: any entity or device capable of carrying the computer program code, recording medium, U disk, mobile hard disk, disk, optical disk, computer memory, read-only memory (ROM, Read-Only Memory), random access memory (RAM, Random Access Memory), electric carrier signal, telecommunication signal and software distribution medium. It should be noted that the content contained in the computer-readable medium can be appropriately increased or decreased according to the requirements of legislation and patent practice in the jurisdiction. For example, in some jurisdictions, according to legislation and patent practice, computer-readable media do not include electric carrier signals and telecommunication signals.

[0151] Although the embodiments of the present invention are described above in conjunction with the accompanying drawings, the present invention is not limited to the above-mentioned specific embodiments and application fields, and the above-mentioned specific embodiments are only illustrative and instructive, rather than restrictive. Under the guidance of the specification, a person skilled in the art can also make many forms without departing from the scope of protection of the claims of the present invention, all of which belong to the scope of protection of the present invention.

Claims

1. A near-surface anisotropic first-arrival wave tomographic inversion modeling method, characterized in that: include: According to the field observation system SPS file and seismic data, the observation system is set up, the first arrival time is picked, and the initial velocity model is established; The main shallow surface reflection layer is picked up based on the isotropic deep deflection profile, and the near-surface stratum azimuth and dip field model is established based on the main shallow surface reflection layer picked up based on the isotropic deep deflection profile; The first arrival time is subjected to azimuthal isotropic first arrival wave tomographic inversion, and the near-surface Epsilon and Delta initial models are established based on the azimuthal isotropic first arrival wave tomographic inversion results. According to the near-surface stratum azimuth and dip field model, initial velocity model, and near-surface Epsilon and Delta initial model, the anisotropic first arrival travel time is forward modeled to obtain the theoretical first arrival time and ray path. The anisotropic first-arrival wave tomography inversion equation is established according to the theoretical first-arrival time and ray path, and the anisotropic first-arrival wave tomography inversion equation is subjected to sparse inversion in the wavelet domain to obtain the inversion result; The inversion results are evaluated. If the convergence conditions are met, the inversion results are output. If the convergence conditions are not met, the initial model is updated, including the velocity, Epsilon and Delta, the first arrival time and ray path are recalculated, and the inversion is performed until the convergence conditions are met and the inversion results are output.

2. A near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The results of azimuthal isotropic first-arrival tomographic inversion are specifically the velocity, main energy ray azimuth and dip angle of each azimuthal inversion.

3. A near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The initial velocity model is obtained based on the refraction wave delay time method.

4. A near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The theoretical first arrival time and ray path are obtained by first obtaining the wavefront time equation of the anisotropic medium, which is: Where τ is the wavefront time, α is the phase angle, and v(α) is the phase velocity; Thomsen obtained an approximate expression for the phase velocity under weak anisotropy conditions: v(a)=v n (1+sin 4 α+δsin 2 acos 2 a) v n is the velocity in the normal direction of the formation, δ is Delta, ε is Epsilon, and both δ and ε are Thomsen parameters; Based on the approximate expression of phase velocity obtained by Thomsen under weak anisotropy conditions, the wavefront time equation of anisotropic media is solved to obtain the theoretical first arrival time and ray path.

5. The near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The anisotropic first-arrival wave tomographic inversion equation is established according to the theoretical first-arrival time and the ray path. Specifically, according to the previously determined ray path and theoretical first-arrival time, the length of the ray in the grid and the theoretical picking time and the calculated travel time difference are obtained to establish the anisotropic first-arrival wave tomographic inversion equation, which includes three unknowns: the formation normal velocity, Epsilon and Delta.

6. A near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The wavelet domain sparse equations are solved by the LSQR algorithm.

7. The near-surface anisotropic first-arrival wave tomographic inversion modeling method according to claim 1, characterized in that: The output inversion results include normal velocity field, Epsilon field and Delta field.

8. A near-surface anisotropic first-arrival wave tomographic inversion modeling system, characterized in that: include: The first model building module: according to the field observation system SPS file and seismic data, set up the observation system, pick up the first arrival time, and build the initial velocity model; The second model building module: picking up the main shallow surface reflection layer based on the isotropic deep deflection profile, and establishing the near-surface stratum azimuth and dip field model according to the main shallow surface reflection layer picked up by the isotropic deep deflection profile; The third model building module: perform azimuthal isotropic first-arrival wave tomographic inversion on the first-arrival time, and establish the near-surface Epsilon and Delta initial models based on the azimuthal isotropic first-arrival wave tomographic inversion results; Forward modeling module: forward modeling of anisotropic first arrival travel time based on near-surface formation azimuth and dip field model, initial velocity model, and near-surface Epsilon and Delta initial model to obtain theoretical first arrival time and ray path; Solution module: establish anisotropic first-arrival wave tomography inversion equations based on theoretical first-arrival time and ray path, and perform sparse inversion in wavelet domain; Evaluation module: Evaluate the inversion results. If the convergence conditions are met, the inversion results are output. If the convergence conditions are not met, the initial model is updated, including velocity, Epsilon and Delta, the first arrival time and ray path are recalculated, and the inversion is performed until the convergence conditions are met and the inversion results are output.

9. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the steps of the near-surface anisotropic first-arrival wave tomography inversion modeling method as described in any one of claims 1 to 7 are implemented.

10. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the steps of a near-surface anisotropic first-arrival wave tomography inversion modeling method as described in any one of claims 1 to 7 are implemented.

Images