A complex structural velocity modeling method based on structural constraints using well-seismic combined methods

Abstract

The present invention discloses a well-seismic combined complex structure velocity modeling method based on structural constraints, which is applied to the field of seismic data processing and interpretation. In view of the fact that the velocity model in the medium with complex structure still does not conform to the structure and the sedimentation pattern in the prior art, the present invention first extracts the low-frequency information of the seismic velocity model as the background field, then performs structural modeling according to the structural interpretation data, analyzes the sedimentation pattern, and adds this information as constraints to the process of velocity correction. Finally, a more refined initial velocity model that is more consistent with the geological structure and sedimentation characteristics is obtained.

Description

Technical Field

[0001] The invention belongs to the field of seismic data processing and interpretation, and particularly relates to a velocity inversion technology. Background Art

[0002] Velocity analysis and modeling are key issues in seismic data processing and interpretation. With the development and maturity of VSP technology, the accuracy and credibility of the one-dimensional velocity results of VSP velocity analysis are relatively high. Therefore, the velocity modeling of well-seismic joint has always been the focus of the industry. Since the traditional stacking velocity analysis method based on horizontal layered media cannot adapt to the situation of lateral velocity changes and complex structures, it has become increasingly difficult to apply to the increasingly high-precision imaging requirements of complex structures.

[0003] In recent years, velocity analysis methods based on tomographic inversion are important tools for realizing velocity modeling. This is due to the rapid development of computer performance on the one hand, and the increasingly high-precision acquisition systems on the other hand, which provide high-quality data support for calculations. VSP zero-bias data can provide more accurate one-dimensional velocity, but the problem of lateral velocity changes caused by complex structures requires a reliable method to solve. The traditional approach is to use seismic offset profiles as a guide for velocity inversion to obtain more accurate velocity modeling. However, this method still has the problem that the velocity model does not match the structure and the sedimentary pattern in complex structural media. The reason is that the multi-solution of velocity inversion will lead to singular values ​​in the inversion process of some areas.

[0004] Because velocity inversion using first arrival travel time is an underdetermined problem, there must be multiple solutions. In the case of complex structures, if no constraints are added, the calculated velocity model will mostly have singular values. The conventional method is to carry out a velocity modeling method combining well and seismic data, make full use of ground seismic data to calculate seismic velocity, and then use VSP first arrival information to correct the seismic velocity. This solves the problem of insufficient seismic velocity accuracy and the problem of singular values ​​appearing in the velocity model using only VSP information for iterative inversion. However, this method still has the following problems:

[0005] 1) The near-offset VSP direct wave only passes through the near-well medium and cannot directly correct the medium velocity at medium and far well source distances.

[0006] 2) Directly use the ground seismic velocity as the initial model for iteration. The high-frequency information in the seismic velocity will seriously affect the inversion effect. If the error of the high-frequency information is large, the accuracy of the velocity model will be reduced.

[0007] 3) Although the initial velocity model contains structural information, no structural constraints are added during the inversion process, and the structural information contained in the ground seismic velocity is not accurate enough, which leads to a low degree of match between the velocity model and the structural model, and the velocity model does not conform to the sedimentary pattern. Summary of the invention

[0008] In order to solve the above technical problems, the present invention proposes a well-seismic combined complex structure velocity modeling method based on structural constraints.

[0009] The technical solution adopted by the present invention is: a well-seismic combined complex structure velocity modeling method based on structural constraints, comprising:

[0010] S1. Establish a unified observation system;

[0011] S2, based on the observation system established in step S1, obtaining a ground seismic velocity model, performing low-pass filtering on the ground seismic velocity model, and obtaining an initial velocity model;

[0012] S3. Pick up the actual first arrival time according to the VSP data provided by the work area;

[0013] S4, using the shortest path ray tracing algorithm to perform forward modeling on the initial velocity model, and obtain a set of first arrival travel time information and ray path information;

[0014] S5, calculating the difference between the first arrival travel time information obtained in step S3 and the actual first arrival time obtained in step S3, thereby obtaining the first arrival difference;

[0015] S6. Calculate the velocity correction on the ray path based on the ray path information and the first arrival difference, introduce the sedimentary structure model as a constraint, correct the initial velocity from shallow to deep, and obtain the final velocity model.

[0016] The sedimentary structure model is specifically:

[0017] Based on the acquired seismic data, the interpreters combined the ground seismic post-stack data to interpret the stratigraphic position of each well layer and obtain the original seismic interpretation data; based on the original seismic interpretation data, they carried out structural modeling to form a block geological model;

[0018] On the basis of structural modeling, the sedimentation pattern of each geological block is analyzed.

[0019] Beneficial effects of the invention: The invention fully utilizes and mines known structural information and sedimentary laws, adds them as important constraints to the inversion process, realizes velocity modeling of complex high and steep structures in mountainous areas, and improves the resolution and accuracy of the velocity model. The effectiveness of the invention is verified by actual data velocity modeling, improves the quality of velocity modeling and provides new ideas. BRIEF DESCRIPTION OF THE DRAWINGS

[0020] Figure 1 Schematic diagram of shortest path wavelet propagation;

[0021] Figure 2 This is a schematic diagram of the travel time of the downlink wave;

[0022] Figure 3 It is a flow chart of the scheme of the present invention;

[0023] Figure 4 is a schematic diagram of the observation system;

[0024] Figure 5 It is a schematic diagram of low-pass filtering;

[0025] Figure 6 Modeling schematics for construction;

[0026] Figure 7 Schematic diagram of the deposition mode;

[0027] Among them, (a) is the deposition mode parallel to the top, (b) is the deposition mode parallel to the bottom, and (c) is the isoproportional deposition mode;

[0028] Figure 8 This is a schematic diagram of the division of small layers within the block;

[0029] Fig. 9 It is a schematic diagram of the velocity grid data structure;

[0030] Fig.10 This is a schematic diagram of speed correction;

[0031] Fig.11 This is a schematic diagram of VSP first arrival picking;

[0032] Fig.12 Provide a schematic diagram for construction interpretation;

[0033] Fig.13 Modeling schematics for geological structures;

[0034] Fig.14 Schematic diagram for modeling sedimentary structures;

[0035] Fig.15 This is a schematic diagram of the ground seismic velocity model;

[0036] Fig.16 Schematic diagram of the effect of modeling the initial velocity. DETAILED DESCRIPTION

[0037] In order to facilitate those skilled in the art to understand the technical content of the present invention, the relevant prior art is first described:

[0038] The relevant prior art is as follows:

[0039] 1. Shortest Path Ray Tracing Algorithm

[0040] Ray tracing technology is an important part of earthquake research. Traditional methods include wavefront expansion methods based on the eikonal equation, Fermat's principle and Huygens' principle, seismic ray symplectic geometry algorithms based on Hamiltonian systems, and shortest path algorithms based on Fermat's principle in graph theory. Among them, the shortest path algorithm is widely used because of its stable and flexible characteristics and is suitable for various complex structural models.

[0041] The shortest path algorithm was first proposed by Nakanishi and further studied by Moser. This method divides the velocity model into a network of nodes, with adjacent nodes connected by lines to indicate the adjacency relationship and record the distance. Starting from the earthquake source point, it expands to the surrounding nodes, and finds the shortest time and path for the seismic wave to reach the surrounding nodes. Then, each node calculates the shortest time and path of its surrounding nodes ( Figure 1 , Figure 2 ). When the calculation is completed, a shortest path tree can be obtained, which records the shortest time and propagation path of all nodes receiving seismic waves. This method can improve the calculation accuracy by encrypting the distribution of nodes, but it will also increase the amount of calculation and reduce the efficiency. Thanks to the improvement of computer computing power and the optimization of this method by predecessors, the availability of the shortest path ray tracing algorithm has been greatly improved.

[0042] 2. First arrival travel time tomographic velocity inversion

[0043] The basic idea of ​​the tomographic velocity inversion method based on first arrival travel time is to pick up the first arrival time of the seismic wave and use the first arrival time to invert the underground velocity information. The basic principle of conventional layer velocity calculation was proposed by Dix, but this method is only applicable to horizontal layered media. Later generations have conducted a lot of research and improvement on this method, making the tomographic velocity inversion method applicable to complex geological conditions. Among them, the inversion iteration method has been widely used and studied. Its principle is to set an initial velocity model and observation system, and perform forward simulation on the model to obtain the forward first arrival travel time information. Then, the difference between the forward first arrival and the actual earthquake first arrival is used to modify the velocity on the seismic wave path, and repeated iterations are performed to achieve a velocity model that meets the accuracy requirements.

[0044] The process of tomographic velocity inversion can be understood as a process of seeking the optimal solution. In order to obtain a more accurate velocity model, on the one hand, it is necessary to improve the accuracy of the forward travel time, and on the other hand, it is necessary to add constraints to deal with the problem of singular values. Under the condition of ensuring the above two points, the first arrival travel time tomographic velocity inversion method can obtain a high-resolution and more accurate velocity model.

[0045] The present invention is further explained below in conjunction with the accompanying drawings. Figure 3 As shown, the method of the present invention comprises the following steps:

[0046] 1. Data preparation

[0047] First, you need to preprocess the existing data. The steps are as follows:

[0048] 11. Establish a unified observation system: In order to ensure the uniformity of the coordinates of the velocity model, sedimentary structure model, etc., it is necessary to normalize them according to the actual observation system. The accuracy of the observation system determines the accuracy of the calculation results and also affects the efficiency of the forward modeling. Figure 4 As shown, it mainly includes model boundary, velocity grid size, excitation point coordinates, receiving point coordinates, etc., among which the grid size affects the accuracy and efficiency of calculation. In this embodiment, the grid size is set to 10m.

[0049] 12. If Figure 5 As shown, the ground seismic velocity model is low-pass filtered: a two-dimensional Butterworth low-pass filter is designed to perform low-pass filtering on the ground seismic velocity, filter out the high-frequency components, and obtain the low-frequency velocity background field as the initial velocity model.

[0050] 13. Based on the VSP data provided by the work area, the interpreters accurately pick up the first arrival time of the near-offset VSP data to obtain the first arrival travel time data.

[0051] 2. Sedimentation pattern analysis and velocity structure analysis

[0052] 21. Based on the seismic data obtained from oil fields and geophysical service companies, the interpreters will interpret each well layer in combination with the ground seismic post-stack data to obtain the original seismic interpretation data. Figure 6 As shown, structural modeling is carried out on this basis to form a block geological model. It is generally believed that the velocity in each block has a similar variation pattern. The data of the block geological model should be completely consistent with the velocity model in scale.

[0053] 22. Based on the structural modeling, analyze the sedimentation pattern of each geological block. There are three main types of sedimentation patterns, such as Figure 7 Different sedimentary patterns reflect the sedimentary process of each stratum. This information is used to divide the small layers within the block, establish the sequence structure within the block, and form the longitudinal velocity change within the block, that is, the deeper the small layer, the higher the velocity. According to this law, the velocity layer structure is obtained as follows Figure 8 As shown, the thickness of the small layer is designed to be 3-5 grid widths.

[0054] 23. After the above steps, we get an initial velocity grid model. Each grid in the model records velocity information, block information, and layer index information. The velocity data structure of a single grid is as follows: Fig. 9As shown in the figure. The block index marks which velocity block the current grid is in, and the control layer index marks which sublayer in the velocity block the current grid is in. The velocity model grid index is also added to the block information and the control layer index information, which makes it easy to reversely search for grids in the same block or layer. The weight information is calculated based on the location of the current block. The closer to the wellhead and the shallower the depth, the higher the weight of the block, which means the credibility of the impact of the VSP data on the velocity here.

[0055] 3. Modify the velocity based on sedimentary characteristics and structural constraints

[0056] 31. The initial velocity model is forward modeled using the shortest path ray tracing algorithm. The number of nodes on each grid edge is designed to be 50. When the grid width is 10m, the forward modeling travel time accuracy is about 1ms. A set of first arrival travel time information and ray path information is obtained through calculation, and then the first arrival obtained by the forward modeling is subtracted from the actual first arrival to obtain the first arrival difference.

[0057] 32. According to the stratification of the sedimentary structure model, first calculate the velocity correction of the shallowest layer based on the ray path length and the first arrival difference and modify all grid velocities of the current layer, and then calculate the velocity of the next layer in the same way. Fig.10 As shown, if the corrected speed of point i is vel i＝ vel i0 +Δvel i , then the velocity of another far-off point j in the layer sequence Q where point A is located is:

[0058]

[0059] In the above formula, vel A0 Andvel B0 are the initial velocities of points A and B respectively. In this way, the velocity points at medium and far offsets in the same layer are constrained by the velocity at near offsets, that is, the velocity changes in the same layer are similar, so the velocities of all points in a layer can be calculated based on the near offset points.

[0060] The calculation formula for the velocity correction at point i on the first floor is as follows:

[0061]

[0062] Among them, t0 is the first arrival travel time of the ray, t is the first arrival travel time obtained by forward modeling, Δt=t-t0 is the first arrival difference, and L is the path length of the ray. Because the velocity is corrected from shallow to deep layers, when calculating the velocity correction amount of the next layer, the velocity of the upper layer has been corrected, so the actual remaining time difference of the ray in this layer is:

[0063]

[0064] The meaning of the above formula is that the sum of the compensated time difference after velocity correction of each layer is equal to the total first arrival difference, vel 0k +Δvel k is the corrected velocity of the point on the kth layer of the ray path, vel 0k is the original velocity of the point, Δvel k is the velocity correction corresponding to the point on the kth layer, L k is the path length of the ray at the kth layer. Then the velocity correction at point i can be calculated as:

[0065]

[0066] For a grid with multiple ray constraints, the weight is introduced into the formula. The rays emitted by the gun with a farther offset have a lower velocity modification weight for the grid. Therefore, the final velocity correction is:

[0067]

[0068] In the above formula, Δvel i，j is the velocity correction of the ray emitted by the jth shot to point i, P j is the weight of the gun, taking the inverse of the offset. The velocity is corrected in order from shallow to deep to get the final velocity.

[0069] 33. The corrected velocity model is forward modeled again using the shortest ray path algorithm to obtain new first arrival information and calculate the first arrival difference again with the actual first arrival. Set the accuracy requirements to the average first arrival difference absolute value less than 50ms and the maximum first arrival difference less than 200ms, and determine whether the first arrival difference meets the accuracy requirements. If it does not meet the accuracy requirements, repeat step 32 and use the new velocity model as the initial model to recalculate the velocity correction amount and perform correction until it meets the accuracy requirements.

[0070] 34. According to the correction results of the layered model, there will be a sudden change in velocity at the boundary of the layer, so the boundary needs to be moderately smoothed. Perform a two-dimensional mean filter on the velocity model results, and design the window size to be equal to the thickness of the small layer. When the thickness of the small layer is 3, the convolution kernel is designed as:

[0071]

[0072] The technical effects of the present invention are described below in conjunction with specific data:

[0073] This example shows the application of the method of the present invention in the WVSP treatment of the Dabei 4 well, including the following steps:

[0074] First, accurately pick the first arrival of the VSP and obtain the first arrival travel time information as follows: Fig.11Then the post-stack seismic data is interpreted structurally. Fig.12 To obtain a detailed structural interpretation of the Dabei 4 well, it is necessary to interpret the control horizon in each velocity layer and well layer depth. In theory, the more detailed the horizon interpretation, the more accurate the velocity control.

[0075] Using the results of fine structural interpretation, a block geological model is constructed. Fig.13 As shown in Figure 2, we believe that each geological block has similar velocity variation patterns. In other words, velocity correction requires constraints from geological structures.

[0076] Secondly, based on the structural modeling, the sedimentation pattern of each block is analyzed. According to the sedimentation pattern, the law of longitudinal velocity change in each block is analyzed. Fig.14 As shown in Figure 2, under the constraints of the structural model, the sequence structure of each block is formed, thereby obtaining the structural constraints of the velocity model.

[0077] Finally, the surface seismic velocity model is corrected using the near-offset VSP first arrival and tectonic constraint velocity structure. Fig.15 The ground seismic velocity and Fig.16 It can be seen from the corrected velocity that the corrected velocity maintains the low-frequency background of the ground seismic velocity and is consistent with the near-offset VSP velocity. The corrected velocity is more precise and reasonable.

[0078] Those skilled in the art will appreciate that the embodiments described herein are intended to help readers understand the principles of the present invention, and should be understood that the scope of protection of the present invention is not limited to such specific statements and embodiments. For those skilled in the art, the present invention may have various changes and variations. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of the claims of the present invention.

Claims

1. A complex structure velocity modeling method based on structural constraints and combined with well-seismic data, characterized in that: include: S1. Establish a unified observation system; S2, based on the observation system established in step S1, obtaining a ground seismic velocity model, performing low-pass filtering on the ground seismic velocity model, and obtaining an initial velocity model; S3. Pick up the actual first arrival time according to the VSP data provided by the work area; S4, using the shortest path ray tracing algorithm to perform forward modeling on the initial velocity model, and obtain a set of first arrival travel time information and ray path information; S5, calculating the difference between the first arrival travel time information obtained in step S3 and the actual first arrival time obtained in step S3, thereby obtaining the first arrival difference; S6, calculating the velocity correction amount on the ray path according to the ray path information and the first arrival difference, and introducing the sedimentary structure model as a constraint, correcting the initial velocity from shallow to deep layers, and obtaining the final velocity model; The sedimentary structure model is specifically: Based on the acquired seismic data, the interpreters combined the ground seismic post-stack data to interpret the stratigraphic position of each well layer and obtain the original seismic interpretation data; based on the original seismic interpretation data, they carried out structural modeling to form a block geological model; On the basis of structural modeling, the sedimentation pattern of each geological block is analyzed; The calculation formula for velocity correction of each layer in each geological block is: Among them, Δvel i represents the velocity correction of the i-th layer, Δt i Indicates the actual remaining time difference of the ray in this layer, L i represents the path length of the ray at layer i, vel 0i represents the initial velocity of the i-th point; When calculating the first floor velocity correction, Δt i The value is Δt, where Δt is the first arrival difference; Δt i The calculation formula is: Among them, Δt is the first arrival difference, vel 0k is the original velocity of the kth layer, L k is the path length of the ray at layer k, Δvel k Represents the velocity correction amount of the kth layer.

2. The method for complex structural velocity modeling based on structural constraints and combined well-seismic data according to claim 1 is characterized in that: After step S6, the method further includes performing two-dimensional mean filtering on the velocity model result obtained in step S6.

Images