Dual complex structure area pre-stack depth migration velocity modeling method and device

Abstract

The present specification provides a pre-stack depth migration velocity modeling method and device for a double complex structure area, for a medium-deep layer, seed point interval velocity value sets are determined based on well logging data of a target straight well and seismic processing velocity spectrum data, and then the interval velocity values in the seed point interval velocity value sets are polynomial fitted to determine a velocity trend surface of a deeper layer; multiple iterations are performed, and at each iteration, interpolation is performed between adjacent two interval velocity trend surfaces, small layers are divided, and tomographic inversion is used to continuously correct a new velocity trend surface. As can be seen, the method can control macro large amplitude structure forms by controlling polynomial fitting parameters and adjusting the fitting degree of the velocity boundary surface and the actual geological conditions, can control the macro large amplitude structure forms, and in the subsequent velocity model correction process through multiple iterations, continuously adjusts and increases the velocity trend surface, gradually depicts the micro small amplitude structure forms, and finally obtains a high-precision pre-stack depth migration velocity model, meeting the needs of more complex migration algorithms for migration imaging.

Description

Technical Field

[0001] This application relates to the field of seismic exploration technology, and in particular to a method and apparatus for modeling pre-stack depth migration velocity in complex tectonic zones. Background Technology

[0002] Dual-complex structures refer to regions with both complex surface and subsurface structures, primarily distributed in mountainous seismic acquisition areas. For example, foreland thrust belts develop inverse-thrust structures, accompanied by imbricate thrust patterns, forming fault blocks, fault noses, and compressional anticlines. The seismic wavefield in dual-complex structural regions is highly complex, and horizontal stacking and post-stack migration techniques based on the assumption of horizontally layered media present certain challenges. Pre-stack time migration, based on the non-zero offset imaging theory of the double square root equation, sums the amplitude along the diffraction curve travel trajectory of the non-zero shot-receiver offset, images each common-shot-receiver offset profile separately, and transfers the reflected wave energy to the true subsurface location. This ensures that the velocity at the analysis point represents the true root-mean-square velocity at the subsurface imaging point location, solving to some extent the problem of multiple solutions for velocities in formations with different dip angles, but still cannot meet the imaging accuracy requirements for dual-complex structural regions. Pre-stack depth migration (PSM), especially wave equation-based reverse time migration based on two-way wave travel time, can adapt to problems such as drastic lateral variations in subsurface seismic wavefields, difficulties in imaging steep structures, and vertical velocity reversal caused by high-velocity geological bodies. It effectively improves the migration imaging quality in complex tectonic zones, providing high-quality and accurate basic data for exploration and development. However, PSM requires extremely high accuracy in velocity modeling. Therefore, how to improve the accuracy of velocity modeling in PSM in complex tectonic zones has always been a hot topic and a challenge in geophysical research and industry. This invention aims to solve this problem. Summary of the Invention

[0003] The purpose of this application is to provide a method and apparatus for pre-stack depth migration velocity modeling in complex structural regions, so as to improve the accuracy of pre-stack depth migration velocity modeling in complex structural regions.

[0004] To address the aforementioned technical problems, this specification provides a method for pre-stack depth migration velocity modeling in a dual-complex structural zone, comprising: determining the interface between a first velocity layer and a second velocity layer on each micrologging well based on micrologging data of the dual-complex structural zone; the depth of the second velocity layer being greater than the depth of the first velocity layer; determining multiple layer velocity values ​​on each target vertical well based on logging data of each target vertical well, and using the layer velocity values ​​as elements in a set of layer velocity values; performing polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer; repeating the following method until a predetermined condition is met to obtain a second target velocity volume of the second velocity layer; interpolating between two adjacent velocity trend surfaces to obtain a new velocity trend surface; and iteratively correcting the new velocity trend surface using tomographic inversion.

[0005] In some embodiments, before performing polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer, the method includes: determining an initial velocity volume for pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; determining the intersection points of the longitudinal and transverse target lines of the second target velocity volume with the layer boundary, and determining the layer velocity value at the intersection points based on the initial velocity volume; and using the layer velocity value at the intersection points as an element of the set of layer velocity values.

[0006] In some embodiments, multiple layer velocity values ​​on each target vertical well are determined based on logging data of each target vertical well. This includes determining multiple layer velocity values ​​on each target vertical well using the following methods: distinguishing sandstone and mudstone marker layers based on the GR curve value of the current target vertical well; performing core repositioning and lithological calibration on the current target vertical well and determining the interfaces of each marker layer; determining the intersection points of the current target vertical well and the interfaces of each marker layer; and determining the layer velocity values ​​at each intersection point based on the sonic logging curve data of the current target vertical well.

[0007] In some embodiments, the predetermined condition includes the accuracy of the velocity trend surface of the second velocity layer reaching a predetermined accuracy threshold; wherein the accuracy of the velocity trend surface of the second velocity layer is calculated by the following method: determining the initial CIP gather in the depth domain based on the logging data of the target vertical well; and determining the accuracy of the velocity trend surface of the second velocity layer by the flattening degree of the in-phase axis after dynamic correction of the CIP gather.

[0008] In some embodiments, after determining the second target velocity body of the second velocity layer, the method further includes: calculating the sum of squares of the differences between the actual values ​​of each layer velocity value in the set of layer velocity values ​​and the fitted values ​​corresponding to the actual values, as a residual value; wherein the fitted value refers to the layer velocity value corresponding to each velocity trend surface in the second target velocity body; and using the residual value as part of the second target velocity body of the second velocity layer.

[0009] In some embodiments, the method further includes: establishing a first target velocity body for the first velocity layer using micrologging data and first-arrival data from a borehole.

[0010] In some embodiments, before determining the velocity trend surface of the second velocity layer based on the set of layer velocity values, the method includes: determining an initial velocity volume for pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; replacing a portion of the initial velocity volume belonging to the first velocity layer with a first target velocity volume of the first velocity layer; determining the intersection points of the longitudinal and transverse target lines of the replaced second target velocity volume with the layer boundary, and determining the layer velocity value of the intersection points based on the initial velocity volume; and using the layer velocity value of the intersection points as an element of the set of layer velocity values.

[0011] The second aspect of this specification provides a pre-stack depth migration velocity modeling device for a dual complex structural zone, comprising: a first determining unit, configured to determine the interface between a first velocity layer and a second velocity layer on each micrologging well based on micrologging data of the dual complex structural zone; wherein the depth of the second velocity layer is greater than the depth of the first velocity layer; a second determining unit, configured to determine multiple layer velocity values ​​on each target vertical well based on logging data of each target vertical well, and to use the layer velocity values ​​as elements in a set of layer velocity values; a fitting unit, configured to perform polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer; and an iterative unit, configured to repeatedly execute the following method until a predetermined condition is met to obtain a second target velocity body of the second velocity layer: interpolating between two adjacent velocity trend surfaces to obtain a new velocity trend surface; and iteratively correcting the new velocity trend surface using tomographic inversion.

[0012] In some embodiments, the apparatus further includes: a third determining unit, configured to determine an initial velocity volume for pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; a fourth determining unit, configured to determine the intersection point of the longitudinal and transverse target lines of the second target velocity volume with the stratigraphic interface, and to determine the layer velocity value of the intersection point based on the initial velocity volume; and a fifth determining unit, configured to use the layer velocity value of the intersection point as an element of the set of layer velocity values.

[0013] In some embodiments, the second determining unit includes: a first determining subunit, used to distinguish sandstone and mudstone marker layers based on the GR curve value of the current target vertical well; a second determining subunit, used to perform core repositioning and lithological calibration of the current target vertical well, and determine the interface of each marker layer; a third determining subunit, used to determine the intersection point between the current target vertical well and the interface of each marker layer; and a fourth determining subunit, used to determine the layer velocity value at each intersection point based on the sonic logging curve data of the current target vertical well.

[0014] In some embodiments, the predetermined condition includes the accuracy of the velocity trend surface of the second velocity layer reaching a predetermined accuracy threshold; wherein the accuracy of the velocity trend surface of the second velocity layer is calculated by the following means: determining the initial CIP gather in the depth domain based on the logging data of the target vertical well; and determining the accuracy of the velocity trend surface of the second velocity layer by the flattening degree of the in-phase axis after dynamic correction of the CIP gather.

[0015] In some embodiments, the apparatus further includes: a calculation unit, configured to calculate the sum of squares of the differences between the actual values ​​of each layer velocity value in the set of layer velocity values ​​and the fitted values ​​corresponding to the actual values, as a residual value; wherein the fitted value refers to the layer velocity value corresponding to each velocity trend surface in the second target velocity body; and a sixth determination unit, configured to use the residual value as part of the second target velocity body of the second velocity layer.

[0016] In some embodiments, the apparatus further includes a seventh determining unit, configured to establish a first target velocity body of the first velocity layer using micrologging data and first arrival data from the cannon.

[0017] In some embodiments, before performing polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer, the method includes: an eighth determining unit, configured to determine an initial velocity volume of pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; a replacement unit, configured to replace a portion of the first velocity layer in the initial velocity volume with a first target velocity volume of the first velocity layer; a ninth determining unit, configured to determine the intersection points of the longitudinal and transverse target lines of the replaced second target velocity volume with the layer boundary, and determine the layer velocity value of the intersection points based on the initial velocity volume; and a tenth determining unit, configured to use the layer velocity value of the intersection points as an element of the set of layer velocity values.

[0018] A third aspect of this specification provides an electronic device, comprising: a memory and a processor, wherein the processor and the memory are communicatively connected to each other, the memory stores computer instructions, and the processor executes the computer instructions to implement the steps of the method described in any of the first aspects.

[0019] A fourth aspect of this specification provides a computer storage medium storing computer program instructions that, when executed, implement the steps of the method described in any of the first aspects.

[0020] This specification provides a method and apparatus for pre-stack depth migration velocity modeling in complex structural zones. For shallow near-surface strata, this scheme uses micrologging and first-arrival data from boreholes, employing the interlayer correlation coefficient method to establish a surface velocity model. For intermediate-deep strata, a seed point layer velocity value set is jointly determined based on logging data from the target vertical well and seismically processed velocity spectrum data. Then, polynomial fitting is performed on the layer velocity values ​​in the seed point set to determine the velocity trend surface of deeper strata. Multiple iterations are performed, with interpolation and sub-layering between adjacent velocity trend surfaces in each iteration, and tomographic inversion is used to continuously correct the new velocity trend surface. Therefore, this method, by controlling the polynomial fitting parameters and adjusting the degree of fit between the velocity interface and the actual geological conditions, can effectively control the macroscopic structural morphology. Simultaneously, by replacing the precise velocity with the surface model and adding tomographic inversion modeling constrained by the intermediate-deep velocity trend surface, a more accurate initial depth migration velocity model is obtained. In subsequent iterations of the velocity model, the number and shape of velocity trend surfaces were continuously adjusted and increased, and the constraint scale and range of velocity trend surfaces were refined to gradually characterize the microscopic small-amplitude structure until a high-precision pre-stack depth migration velocity model was finally obtained to meet the migration imaging requirements of more complex migration algorithms. Attached Figure Description

[0021] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this application. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0022] Figure 1 This document presents a flowchart illustrating a pre-stack depth migration velocity modeling method for dual complex structural regions.

[0023] Figure 2 A schematic diagram of the model determined according to the pre-stack depth migration velocity modeling method for dual complex structural regions provided in this specification is shown.

[0024] Figure 3 A flowchart is shown for another pre-stack depth migration velocity modeling method for dual complex structural regions provided in this specification;

[0025] Figure 4 This document presents a flowchart of another pre-stack depth migration velocity modeling method for dual complex structural regions provided in this specification.

[0026] Figure 5 A schematic diagram of micro-logging is shown;

[0027] Figure 6A schematic diagram of micrologging data is shown;

[0028] Figure 7 A schematic diagram is shown showing the interpolation operation using spatially variable interlayer morphological correlation coefficients;

[0029] Figure 8 This specification shows a block diagram illustrating the principle of a pre-stack depth migration velocity modeling device for dual complex structural regions.

[0030] Figure 9 A schematic diagram of the structure of an electronic device provided in this specification is shown. Detailed Implementation

[0031] To enable those skilled in the art to better understand the technical solutions in this application, the technical solutions in the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of this application, and not all of the embodiments. Based on the embodiments in this application, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of this application.

[0032] To overcome the technical challenges of poor seismic data quality and low migration imaging accuracy in the dual-complex tectonic zones of oil and gas basins in central and western China, and to maximize the imaging advantages of pre-stack depth migration technology in these zones, it is urgent to establish a systematic, scientific, and effective pre-stack depth migration velocity modeling method. This will improve the migration imaging accuracy of seismic data in the region, further enhance the quality of seismic data, meet the needs of exploration and development, and contribute to increasing oil and gas reserves and production.

[0033] This specification provides a method for pre-stack depth migration velocity modeling in dual complex structural regions, such as... Figure 1 As shown, the method includes the following steps:

[0034] S10: Based on the micrologging data of the dual complex structural zone, determine the interface between the first velocity layer and the second velocity layer on each micrologging well; the depth of the second velocity layer is greater than the depth of the first velocity layer.

[0035] Micrologging refers to the method of acquiring seismic wave information through methods such as well-hole excitation and surface reception, or surface excitation and well-hole reception, or well-hole excitation and well-hole reception, and further determining near-surface geophysical parameters. The data involved in the micrologging method is called micrologging data. For example... Figure 2 As shown, it includes micrologging well 1 and micrologging well 2.

[0036] like Figure 2 As shown, the first velocity layer may include Figure 2 The low-speed layer and deceleration layer, the second speed layer can be Figure 2 The high-velocity layer is located within the micrologging system. The depth of the micrologging can reach the interface between the first and second velocity layers, meaning it can also reach the interface between the decreasing and high-velocity layers. Therefore, based on micrologging data, the interface between the first and second velocity layers can be determined relatively accurately. Figure 2 The interface between the low-speed layer and the deceleration layer, and the interface between the deceleration layer and the high-speed layer.

[0037] Layer velocity refers to the speed at which seismic waves propagate within layered strata. Layer velocity directly reflects the lithology of the strata and can be used to classify them. It is generally measured using seismic logging or sonic logging, and specifically refers to the velocity of the P-wave; it can also be calculated using reflection records. In seismic exploration, strata with a layer velocity below 1400 m / s are generally called low-velocity strata, those above 3500 m / s are considered high-velocity strata, and strata with a layer velocity between 1400 and 3500 m / s are considered decelerating strata.

[0038] The micrologging data obtained from seismic acquisitions is lithologically calibrated and reinterpreted to finely delineate the boundaries between low-velocity, decreasing-velocity, and high-velocity layers at each micrologging well. Specifically, different layers have different acoustic propagation velocities, meaning different logging curve slopes; therefore, different layers can be distinguished by the slope of the logging curves. Due to the limited number of micrologging wells, previous modeling methods often used linear thickness interpolation, which is highly susceptible to human factors, resulting in uncertain interpretations. The modeling method provided in this manual, within the spatial range where layer velocities cannot be determined by micrologging, uses first-arrival picking with a borehole to calculate layer velocities, thereby delineating low-velocity, decreasing-velocity, and high-velocity layers.

[0039] S20: Based on the logging data of each target vertical well, determine multiple layer velocity values ​​on each target vertical well, and use the layer velocity values ​​as elements in the layer velocity value set.

[0040] In this specification, the elements in the set of layer velocity values ​​are also called seed point layer velocities. A seed point is a spatial location point corresponding to each layer velocity value in the set of layer velocity values.

[0041] A vertical well, as opposed to a directional well, is a well whose designed trajectory is a vertical line. In the design trajectory of a vertical well, theoretically, the inclination angle of all points on the trajectory is zero. However, in actual drilling, there are errors. As long as the error is within the allowable range of a qualified well, it is still called a vertical well regardless of the actual wellbore inclination angle.

[0042] To ensure the accuracy of pre-stack depth migration velocity modeling, the vertical distance from the well point coordinates of the target vertical well to the longitudinal target line of the target velocity body must, in principle, be less than 100 meters.

[0043] Micro-logging wells and vertical wells can be located at surface protrusions and depressions in complex geological structures.

[0044] like Figure 2 As shown, the depth of a vertical well can reach the high-velocity layer; therefore, velocity modeling of the high-velocity layer can be performed using logging data from a vertical well.

[0045] Log data from vertical wells are typically densely distributed along the depth dimension. To improve modeling efficiency, such dense data is unnecessary when modeling the velocity of high-velocity layers. Therefore, it is necessary to first select the dataset required for modeling from the dense vertical well log data. Specifically, the vertical well log data at the depth where the lateral target line of the target velocity body is located can be selected to form the dataset required for modeling.

[0046] In some embodiments, S20 includes determining multiple layer velocity values ​​on each target vertical well using the following method:

[0047] S21: Distinguish sandstone and mudstone marker layers based on the GR curve value of the current target vertical well.

[0048] S22: Perform core repositioning and lithological calibration for the current target vertical well, and determine the boundaries of each marker layer.

[0049] S23: Determine the intersection points between the current target vertical well and the interfaces of each marker layer.

[0050] S24: Determine the layer velocity values ​​at each intersection point based on the sonic logging curve data of the current target vertical well.

[0051] A marker layer is a rock layer that is known in geological age and has obvious characteristics (such as the properties, color, thickness of the rock and the fossils and impurities it contains) and is easy to identify. It can be used as a stratigraphic correlation marker.

[0052] like Figure 2 As shown, through the above steps S21 to S24, the layer velocity values ​​at points A1, A2, A3, and A4 on target vertical well A can be determined, the layer velocity values ​​at points B1, B2, B3, B4, and B5 on target vertical well B can be determined, and the layer velocity values ​​at points C1, C2, C3, and C4 on target vertical well C can be determined.

[0053] Due to the complex surface morphology of the dual-structure zone, it is not convenient to install many vertical wells. Therefore, the accuracy of the velocity trend surface established based solely on the layer velocity values ​​at multiple points on the target vertical well is relatively low. To improve the accuracy of the established velocity trend surface, the initial velocity volume of the pre-stack depth migration can be determined by using the pre-stack time migration velocity spectrum. Based on this initial velocity volume, the number of layer velocity values ​​in the layer velocity value set can be expanded. The velocity trend surface of the second velocity layer can then be determined using the expanded set of layer velocity values. While the expanded layer velocity values ​​(i.e., the layer velocity values ​​in the initial velocity volume) may not be perfectly precise, they contain information on the changing trends, which can help improve the accuracy of the velocity trend surface. Accordingly, such as... Figure 3 As shown, before step S30, the following steps S70, S80 and S90 are also included.

[0054] S70: Determine the initial velocity volume for pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region.

[0055] The pre-stack time migration velocity volume is converted into the pre-stack depth migration velocity volume using the DIX formula.

[0056] In this specification, "velocity body" refers to a collection of velocity values ​​distributed in space. That is, a velocity body contains multiple velocity values, each with a corresponding spatial location.

[0057] S80: Determine the intersection point of the longitudinal and transverse target lines of the second target velocity body with the layer interface, and determine the layer velocity value at the intersection point based on the initial velocity body.

[0058] Intelligent stratigraphic tracking and picking technology can quickly obtain multiple sets of geological stratigraphic layers with high interpretation accuracy and strong seismic reflection characteristics, such as... Figure 2 As shown, the main types include: 1. High-velocity layer top boundaries (Type I / Type II marker beds); 2. Regionally traceable layers with gentle stratigraphic undulations (Type I marker beds), such as horizontal sediments and monoclinic structures; 3. Regionally traceable layers with dramatic stratigraphic undulations (Type II marker beds), such as progradational structures and thrust structures; 4. Locally traceable geological body reflection interfaces with varying morphologies (Type III marker beds), such as buried hill structures and salt dome structures.

[0059] The target velocity body's transverse and longitudinal target lines (i.e., target lines almost parallel to the vertical well) are the grid lines of the reduced model, equivalent to virtual grid lines in space. The transverse target line lies on a plane perpendicular to the longitudinal target line.

[0060] Step S80 can be as follows: after determining the intersection coordinates of the longitudinal and transverse target lines of the second target velocity body with the layer boundary line, the velocity value corresponding to the intersection coordinates in the initial velocity body can be determined, and the velocity value of the intersection coordinates in the initial velocity body can be used as the velocity value of the intersection coordinates in the second target velocity body.

[0061] S90: Use the layer velocity value at the intersection point as an element of the layer velocity value set.

[0062] Combining steps S20 and S80-S90, the elements in the set of layer velocity values ​​include both the layer velocity values ​​on the target vertical well and the layer velocity values ​​obtained from the initial velocity volume. Therefore, the velocity trend surface determined in step S30 can utilize the key points of the velocity trend surface determination on the target vertical well to establish the velocity trend surface framework, and can also utilize the trend information contained in the initial velocity volume to determine the changing trend of the velocity trend surface.

[0063] The layer velocity values ​​determined in steps S70, S80, and S90 above can be used only to determine the second velocity volume of the second velocity layer. The first velocity layer is a near-surface layer, which facilitates the installation of multiple micrologging wells. Therefore, a relatively accurate first velocity volume can be established for the first velocity layer using micrologging data and first-reach data from borehole drilling.

[0064] In some embodiments, such as Figure 4 As shown, before step S80, the following step S100 may also be included: replacing a portion of the first velocity layer in the initial velocity body with a first target velocity body of the first velocity layer. Since the first target velocity body is determined based on micrologging data and first arrival data from borehole drilling, it is more accurate than the portion of the first velocity layer in the initial velocity body. Replacing a portion of the first velocity layer in the initial velocity body with the first target velocity body allows step S80 to directly determine the set of layer velocity values ​​based on the replaced initial velocity body, without needing to further filter out the second velocity layer portion in the initial velocity body to determine the velocity trend surface of the second velocity body based on the second velocity layer portion.

[0065] In some embodiments, when determining the first velocity body, the depth of each observation point and the time of receiving the seismic wave can be determined based on micrologging data, a time-distance map can be generated based on the depth and the time, and then interpretation verification and reinterpretation can be performed. Figure 5 A schematic diagram of micro-logging is shown, where multiple seismic sources are set in the well corresponding to S, and multiple observation points are set in the well corresponding to R. x0 represents the distance between wells S and R, and Q1, Q2...Q N Indicates the observation point. Figure 6This diagram illustrates micrologging data. Time represents the time it takes for the seismic wave to be received at each observation point, and Trace represents the observation point number. The distance between each observation point and the excitation point can be determined based on the observation point number and depth. The time-distance plot shows the relationship between the propagation time of the seismic wave from the epicenter to each observation point on the observation line and the distance between the observation point and the excitation point (i.e., the origin). The time-distance curve for a direct wave is a straight line.

[0066] Based on the depth and time corresponding to each observation point, a time-depth curve is plotted. The time corresponding to each depth point can be expressed as T. j (h)=T i +h / V i T j H represents the time corresponding to the target depth point. i V represents the total depth of the layer above the layer containing the target depth point, h represents the depth of the target depth point within its own layer, and V represents the depth of the target depth point within its own layer. i This represents the layer velocity of the layer containing the target depth point, j represents the target depth point number, and i represents the layer number. For example, if the low-velocity layer is 10m deep, and the target depth point is located in a decreasing velocity layer at a depth of 6m, then the time corresponding to the target depth point is T. j (h)=10+6 / V 降速层 In other words, the time corresponding to the target depth point can be recalculated based on the layer velocity and thickness of the associated layers.

[0067] The depth-time curve is divided into multiple slope segments. Strata are then identified based on these slope segments. Linear regression is performed on the depth points of each layer to obtain the regression coefficients a, b, and the mean standard deviation s. When the number of depth points is greater than 3, points whose deviation from the fitted line is greater than 2.5 times the mean standard deviation can be removed. Finally, a final linear regression is performed to obtain the final values ​​of a and b.

[0068] The coefficient 'a' represents the slope of the line, and the coefficient 'b' represents the intercept of the line. This corresponds to the formula: T j (h)=T i +h / V i , 1 / V i It is the slope of the straight line, corresponding to the coefficient 'a', and the unit of speed is m / s. If the speed unit is converted to m / ms, the speed value is 1 / a × 1000; H i It is the intercept of the line, corresponding to the coefficient b.

[0069] The time-depth curve of the i-th layer is: T i+1 (h)=T i +h / V i Based on the fact that the depths at the intersection of the two time-depth curves are equal, we can obtain: T i +h / Vi =T i+1 +h / V i+1 Let the left side represent the depth of the i-th layer, and the right side represent the depth of the (i+1)-th layer. Therefore, we can obtain:

[0070] Thickness of the i-th layer (Note: T is in milliseconds). This method can be used to determine the thickness of each layer.

[0071] The number of observation points set up to obtain micrologging data is usually limited. Therefore, micrologging data can only determine the outline of formation layers. Although this outline is more accurate than that of the second velocity volume, some layer details still cannot be determined. This specification proposes to use the interlayer morphology correlation coefficient with spatial variation for interpolation to further optimize the modeling results.

[0072] like Figure 7 As shown, for a layer interface ( Figure 7 The bold black curve in the image represents the layer interface. The two control points H on the layer interface can be obtained using the method described above. m H m+1 The points between these two control points can be determined using interpolation methods.

[0073] Specifically, for each interface layer, the elevation WH of the control point interface layer between control points m and m+1 is... j The formula is:

[0074] WH 0j =FH 0j +DH 0j ×K0 (Bottom boundary of the first layer)

[0075] WH 1j =FH 1j +DH 1j ×K1 (Bottom boundary of the second layer)

[0076] WH 2j =FH 2j +DH 2j ×K2 (Bottom boundary of the third layer)

[0077] ...

[0078] Where WH represents the physical point measured in the well log, FH represents the point obtained by interpolation, and DH is the formation depth of the previous layer corresponding to the interpolation point. The final depth at this point needs to be determined by adjusting the K value.

[0079] The surface elevation of each seismic channel. Each seismic channel has an elevation, and of course, there is a certain seismic channel with the highest elevation, which should theoretically be the location closest to the highest point of the work area.

[0080] K0 is the interlayer correlation coefficient, that is, the correlation coefficient between the low-velocity layer bottom interface and the surface.

[0081] The value of K is directly proportional to the inter-layer correlation; k=0 indicates no correlation, and k=1 indicates a strong correlation.

[0082] The method for determining the value of K is as follows:

[0083] a) Based on experience, define before modeling and apply in batches;

[0084] b) Calculate the correlation coefficient of spatial variation based on the inter-layer correlation of morphology at surface control points in a set or local area;

[0085] c) For special sections, interactive modeling is used, that is, the model is dynamically generated by interactively changing the correlation coefficient until it is satisfactory. Specifically, it can be: calculate once in the direction of the main survey line, and then calculate again in the direction of the connecting survey line. If the depth error between the two calculations is small enough, or less than a certain threshold value, such as half a sampling interval, it proves that the modeling accuracy requirement has been met.

[0086] S30: Perform polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer.

[0087] like Figure 2 As shown, based on the layer velocity values ​​at each point on each target vertical well, the velocity trend surface of the second velocity layer can be determined. Figure 2 The thick black line represents the speed trend.

[0088] Based on the set of layer velocity values, and according to the grid spacing of the target line for velocity modeling in the longitudinal and transverse depth domains, a grid-like velocity profile framework for the layer velocity of the work area and the main marker layers is constructed, and all current discrete point data are output.

[0089] Based on this, all discrete point data, including current velocity values, layer data, and well point coordinates, are loaded into the velocity modeling software. After setting the gridding parameters (azimuth, smoothing radius, gridding method, grid side length, and grid boundary indentation distance), the discrete points are transformed according to the azimuth and grid boundary indentation distance. The entire region is divided into square grids to obtain the X and Y coordinates of the grid points. Discrete points falling within the smoothing radius of the grid point are used to calculate the Z value of that grid point using the gridding method (distance-weighted average method). If the Z value of a grid point is not allowed to be null, the smoothing radius needs to be gradually expanded until a discrete point falls within the smoothing radius. If smoothing and densification of the grid data are required, smoothing filtering (Box smoothing) and densification interpolation can be performed to ensure that the Z values ​​of the grid points are sufficiently smoothed. The grid data is checked in the form of a color block diagram to see if it meets the requirements. If it does, it is saved to a grid file and associated with the seed point layer velocity database. The coordinates of the nearest adjacent layer velocity seed points are checked in the form of a color block diagram to see if they can be projected into different grids. If not, it indicates that the grid side length is too large, and a denser grid needs to be set and the gridding process repeated. Seed points refer to the spatial location points corresponding to each velocity value in the set of layer velocity values.

[0090] Because previous seismic interpretation schemes were inherently subjective and prone to multiple interpretations, the depth migration velocity modeling process for newly acquired seismic data requires a critical approach and selective use of the data. Adhering strictly to previously interpreted structural frameworks could mislead the next round of seismic data processing, leading to flawed imaging of subsurface structures. Therefore, this specification provides a method for establishing velocity trend surfaces based on polynomial fitting. This method uses a relatively reliable seismic interpretation scheme at the macroscopic scale. By controlling the polynomial fitting parameters and adjusting the fit between the velocity interface and actual geological conditions, it effectively controls large-scale macroscopic structural features. Simultaneously, by replacing the precise velocity with a surface model and adding tomographic inversion modeling constrained by mid-to-deep velocity trend surfaces, a more accurate initial depth migration velocity model is obtained. In subsequent iterations to refine the velocity model, the number and shape of velocity trend surfaces are continuously adjusted and increased, and the scale and range of velocity trend surface constraints are refined, gradually characterizing microscopic, small-scale structural features until a high-precision pre-stack depth migration velocity model is finally obtained, meeting the requirements of more complex migration algorithms for migration imaging. As the accuracy of the velocity model improves with multiple iterations, each iteration modifies the original velocity trend surface to better fit the macroscopic subsurface trend, while adding new, smaller-scale velocity trend surfaces to characterize microstructures one level smaller than the macroscopic structure. Then, in the next velocity model update iteration, more velocity trend surfaces are added to characterize even smaller-scale microstructures, until the modeling accuracy meets the requirements of pre-stack depth migration. This method returns the initiative of migration imaging to more advanced migration algorithms, rather than artificially setting boundaries and obstacles. Based on the aforementioned meshing, velocity modeling is constrained by polynomial fitting of velocity trend surfaces, and a more accurate final layer velocity model is obtained through multiple iterations.

[0091] The principle of polynomial fitting velocity trend surface: decompose "trend value" and "residual value" from actual data to reveal the trend and pattern of spatial distribution of geographical elements.

[0092] a) Establish a velocity trend surface model

[0093] Taking a bivariate function as an example, let Z represent the initial layer velocity value at a certain point, that is, let Z... i (x i ,y i (1, 2, ..., n) represents, using The fitted value representing the layer velocity trend at each point is then... It is a function of the coordinates (x, y) of a point in the plane, that is:

[0094]

[0095]

[0096] Where: ε i This is the velocity residual value (residual value).

[0097] Obviously, when (x i ,y i When the velocity element changes in space, the above equation (2) describes the interaction between the actual distribution surface, trend surface and residual surface of the velocity element.

[0098] The core of velocity trend surface analysis is to deduce the trend surface from actual underground observations. This is generally done using regression analysis to minimize the sum of squared residuals, i.e.:

[0099]

[0100] This is trend surface fitting in the sense of least squares.

[0101] The mathematical equations used to calculate trend surfaces include polynomial functions and Fourier series, with polynomial functions being the most commonly used. This is because any function can be approximated by a polynomial within a suitable range, and adjusting the degree of the polynomial allows the obtained regression equation to suit the needs of the practical problem.

[0102] The form of a polynomial trend surface:

[0103] (1) First-order trend surface model: (Plane: suitable for near-horizontal or monoclinic marker layers)

[0104] z = a0 + a1x + a2y (4)

[0105] (2) Quadratic trend surface model: (Hyperboloid: suitable for marker layers with small fluctuations)

[0106] z = a0 + a1x + a2y + a3x 2 +a4xy+a5y 2 (5)

[0107] (3) Cubic trend surface model: (Cubic surface: suitable for marker layers with drastic fluctuations)

[0108] z = a0 + a1x + a2y + a3x 2 +a4xy+a5y 2 +a6x 3 +a7x 2 y+a8xy 2 +a9y 3 (6)

[0109] b) Extreme value locations of the trend surface

[0110] The necessary condition for the fitted function f(x,y) to reach an extreme value at the point (x0,y0) is:

[0111]

[0112] Taking a quadratic case as an example, let f(x,y) = 2a0 + a1x + a2y + a3x 2 +a4 xy+a5 y 2 Then the geographical coordinates of the extreme point should satisfy the following system of equations:

[0113]

[0114] The above system of equations (8) is a system of two linear equations in two variables, and has one solution:

[0115]

[0116] That is, the two-dimensional quadric surface has an extreme point, which is determined by the above equation (9).

[0117] The regression equation for a two-dimensional cubic surface is:

[0118] z = a0 + a1x + a2y + a3x 2 +a4xy+a5y 2 +a6x 3 +a7x 2 y+a8xy 2 +a9y 3 (6)

[0119] The location of its extreme points should satisfy the following system of equations:

[0120]

[0121] This is a system of two quadratic equations in two variables, which should have four sets of solutions, therefore, a maximum of four extreme points. Similarly, a two-dimensional quartic trend surface should have a maximum of nine extreme points. For a two-dimensional trend surface, let the highest degree of the trend surface be I, then the maximum number of extreme points is (I-1). 2 .

[0122] c) Parameter estimation of trend surface model

[0123] If the trend surface function is an arbitrary polynomial, that is:

[0124]

[0125] In the above formula:

[0126] —This represents the speed trend value;

[0127] L represents the number of trend lines;

[0128]

[0129] The so-called parameter estimation of the trend surface model is to determine the polynomial coefficients a0, a1, ..., a i (x i , y i ), (1, 2, ..., n) to make the sum of squared residuals Q the smallest. Therefore, the polynomial regression problem can be transformed into a multiple linear regression problem to solve. p

[0130] Denote:

[0131] x = x1, y = x2, x 2 = x3, xy = x4, y 2 = x5, … (12)

[0132] Then (Equation 1-12) can be written as

[0133]

[0134] Here, the residual can be written as:

[0135]

[0136] It can be seen from Equation (1-14) that Q is a function of a0, a1, ..., a p .

[0137] According to the principle of the least squares method, select the coefficients a0, a1, ..., a p that make Q reach the minimum, (p < n). For this purpose, it is necessary to find the partial derivatives of Q with respect to a0, a1, ..., a p :

[0138]

[0139] Let That is, we get:

[0140] <00…In the formulas of this specification, the capital "Z" and the lowercase "z" have the same meaning. (Equation 1-16) includes p + 1 linear equations. Solving this system of linear algebraic equations gives p + 1 coefficients a0, a1, ..., a p . (Equation 1-16) is called the normal equation system and can be written in matrix form, denoted as:

[0142]

[0143] ​Then the normal system of equations (16) can be written as:

[0144] X'Xa=X'z (18)

[0145] Where X' is the transpose of X:

[0146]

[0147] Theoretically, the remaining work is merely a computational problem, that is, simply a problem of how to solve the normal equation system. It can be proven that the coefficient matrix of 'a'...

[0148] A = X'X (20)

[0149] It is neither singular nor singular, therefore the solution is...

[0150] a = [X'X] -1 X'z (21)

[0151] It is the only certainty.

[0152] To solve for a0, a1, ..., a p First, we change the way the normal system of equations (16) is written, and denote it as:

[0153]

[0154] And introduce symbols:

[0155]

[0156]

[0157]

[0158] in, —The average value of the extracted layer velocity at a certain point;

[0159] —The j-th variable x j The average value;

[0160] —The total variance of z.

[0161] In the system of equations (16), from the first equation, we can obtain

[0162]

[0163] Using the symbols above, it can be written as

[0164]

[0165] Substituting (23) into equation (14) yields:

[0166]

[0167] Starting from equation (24), find the partial derivative of Q with respect to aj, and make The system of equations can be obtained as follows:

[0168]

[0169] Equations (13) and (16) are equivalent.

[0170] The system of equations (13) consists of p equations, from which a0, a1, ..., a p Substituting this into equation (23) yields a0, thus obtaining all coefficients a0, a1, ..., a p .

[0171] d) Fitting accuracy of the trend surface model

[0172] Based on the extracted value z i and its coordinates (x) i y i The trend surface fitted to a certain time indicates, to some extent, the coordinates (x) i y i ) and observed value z i The relationship is clear: trend surfaces of different orders approximate the original data to varying degrees. Therefore, we need to study how well the trend surfaces approximate the original data.

[0173] The fitting accuracy of a trend surface is a crucial parameter, referring to how closely the mathematical surface represented by the trend surface equation approximates the actual observed data surface. It varies with the degree of the trend surface; generally, a higher degree results in greater fitting accuracy. However, the fitting accuracy of a trend surface has its own characteristics. It does not require extremely high accuracy; on the contrary, excessively high fitting accuracy can cause the mathematical surface to overly approximate the distribution of the measured data, making it difficult to reflect the main patterns of the distribution and thus failing to describe the overall spatial trend. In polynomial trend surface analysis, the fitting accuracy can be controlled by changing the degree of the polynomial to achieve satisfactory analytical results. In practical applications, the fitting accuracy C of the trend surface is generally expressed by the following formula:

[0174]

[0175] in,

[0176] z i : Extracted value of layer velocity;

[0177] Trend values ​​of layer velocity;

[0178] The average layer velocity;

[0179] For z i The sum of squared residuals; that is, S 偏 ;

[0180] For layer velocity data z i Its average The sum of squares of the differences, i.e., S 总 .

[0181] C represents the percentage of the total sum of squared fluctuations in the velocity sampling data reflected by the velocity trend surface. A higher percentage of C indicates a better fit, but excessively high fitting accuracy will cause the disappearance of trend residue, thus missing valuable outliers.

[0182] e) Residual analysis of trend surface model

[0183] The purpose of trend surface analysis includes two aspects: firstly, to determine the regional variation patterns of certain elements, thus requiring the plotting of trend surface graphs; secondly, to identify local anomaly zones and find subtle, hard-to-detect anomalies, which necessitates solving for the residuals and performing residual analysis. According to the definition of residuals, they are the difference between the observed value at a point and the trend value at that point. After solving for all residuals, a residual contour map can be plotted, yielding the residual (error) values:

[0184]

[0185] With the introduction of residual values, the original data is decomposed into two parts: trend values ​​and residual values. Generally, trend values ​​reflect the overall changes over a large regional area, while residual values ​​reflect the characteristics of localized, small-scale changes. Combining these two aspects can help to conduct more in-depth research on layer velocities. First-order trend surface residual values ​​treat complex changes above quadratic curvature as local anomalies, while second-order trend surface residual values ​​treat changes of cubic curvature and above as local anomalies. Therefore, trend surface residual values ​​of different orders from low to high reflect local anomaly characteristics of different scales from large to small. This method can be used to discover meaningful low-gradient anomaly zones.

[0186] Based on this, in some embodiments, after determining the second target velocity body of the second velocity layer, the following steps SA1 and SA2 may also be included.

[0187] SA1: Calculate the sum of squares of the differences between the actual value of each layer velocity value in the set of layer velocity values ​​and the fitted value corresponding to the actual value, and use it as the residual value; wherein the fitted value refers to the layer velocity value corresponding to each velocity trend surface in the second target velocity body.

[0188] SA2: The remaining value is used as part of the second target velocity body of the second velocity layer.

[0189] S40: Interpolate between two adjacent velocity trend surfaces to obtain a new velocity trend surface.

[0190] Specifically, step S40 may involve performing the following operations for each longitudinal target line: obtaining the first intersection point of the current longitudinal target line with the first velocity trend surface and the second intersection point with the second velocity trend surface; determining the average of the velocity values ​​at the first and second intersection points; determining the midpoint of the line connecting the first and second intersection points; and using this average value as the velocity value at the midpoint. This is just one example of interpolation; other interpolation methods may also be used in actual implementation.

[0191] S50: The new velocity trend surface is continuously corrected by using tomographic inversion.

[0192] Tomography and inversion methods are commonly used in the field of seismic exploration. Tomography and inversion can be used to correct for errors in the data.

[0193] S60: Determine whether the predetermined conditions have been met. If yes, obtain the second target velocity body of the second velocity layer; otherwise, jump to step S40 to continue execution.

[0194] First, the linear difference of the velocity between adjacent velocity trend surfaces is calculated. The rate of velocity change is controlled by dividing the velocity into smaller layers. Finally, iterative corrections are made through tomographic inversion. If the predetermined conditions are not met, steps S40 and S50 are repeated, which means continuously dividing the adjacent velocity trend surfaces into smaller layers.

[0195] In some embodiments, the predetermined condition includes the accuracy of the velocity trend surface of the second velocity layer reaching a predetermined accuracy threshold. The accuracy of the velocity trend surface of the second velocity layer can be calculated by the following method: determining the initial CIP gather for the depth domain based on the logging data of the target vertical well; determining the accuracy of the velocity trend surface of the second velocity layer by the flattening degree of the in-phase axis after dynamic correction of the CIP gather. The flattening degree of the in-phase axis can be directly used as the accuracy of the velocity trend surface of the second velocity layer, or the flattening degree of the in-phase axis can be processed to obtain the accuracy of the velocity trend surface of the second velocity layer. Alternatively, the accuracy of the velocity trend surface of the second velocity layer can be calculated using the above formula (26), that is, the fitting accuracy C is used as the accuracy of the velocity trend surface of the second velocity layer.

[0196] In the velocity modeling process, the target line grid is gradually increased in density. For example, a velocity trend surface is established by using an 800*800m grid to establish the target line grid. This can correct the layer velocity residual value at the velocity line in the 400*400m grid. By gradually increasing the density and iteratively correcting the grid, a high-precision layer velocity model that conforms to the laws of underground geological structure can be obtained.

[0197] The aforementioned modeling method replaces traditional pre-stack depth migration tomography velocity modeling with stratigraphic constraints by establishing a velocity trend surface. This method allows for manual adjustment of the fit between the velocity interface and actual geological conditions by controlling polynomial fitting parameters. It controls large-scale macroscopic structural features, filters out small-scale microscopic structural features, replaces precise velocities with surface models, and incorporates tomographic inversion modeling constrained by mid-to-deep velocity trend surfaces. This results in a final pre-stack depth migration velocity model that meets the migration imaging requirements of more complex migration algorithms.

[0198] The pre-stack depth migration velocity modeling method for complex structural zones presented in this specification, for deeper formations, determines a set of seed point layer velocity values ​​based on logging data from the target vertical well. Then, polynomial fitting is performed on the layer velocity values ​​in the seed point set to determine the velocity trend surface of the deeper formations. This is iterated multiple times, with interpolation and subdivision between adjacent velocity trend surfaces in each iteration. Tomographic inversion is then used to continuously correct the new velocity trend surface. Therefore, this method, by controlling the polynomial fitting parameters and adjusting the degree of fit between the velocity interface and the actual geological conditions, can control large-scale macroscopic structural features, filter out small-scale microscopic structural features, replace precise velocities with surface models, and add tomographic inversion modeling constrained by mid-deep velocity trend surfaces to obtain the final pre-stack depth migration velocity model, meeting the migration imaging requirements of more complex migration algorithms.

[0199] This specification provides a pre-stack depth migration velocity modeling device for dual complex structural regions, such as... Figure 8 As shown, the device includes a first determining unit 10, a second determining unit 20, a fitting unit 30, and an iterating unit 40.

[0200] The first determining unit 10 is used to determine the interface between the first velocity layer and the second velocity layer on each micro-logging well based on the micro-logging data of the dual complex structural zone; the depth of the second velocity layer is greater than the depth of the first velocity layer.

[0201] The second determining unit 20 is used to determine multiple layer velocity values ​​on each target vertical well based on the logging data of each target vertical well, and to use the layer velocity values ​​as elements in the layer velocity value set.

[0202] The fitting unit 30 is used to perform polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer.

[0203] The iteration unit 40 is used to repeatedly execute the following method until a predetermined condition is met to obtain the second target velocity body of the second velocity layer: interpolation is performed between two adjacent velocity trend surfaces to obtain a new velocity trend surface; the new velocity trend surface is continuously iteratively corrected using tomographic inversion.

[0204] In some embodiments, the apparatus further includes: a third determining unit, configured to determine an initial velocity volume for pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; a fourth determining unit, configured to determine the intersection point of the longitudinal and transverse target lines of the second target velocity volume with the stratigraphic interface, and to determine the layer velocity value of the intersection point based on the initial velocity volume; and a fifth determining unit, configured to use the layer velocity value of the intersection point as an element of the set of layer velocity values.

[0205] In some embodiments, the second determining unit includes: a first determining subunit, used to distinguish sandstone and mudstone marker layers based on the GR curve value of the current target vertical well; a second determining subunit, used to perform core repositioning and lithological calibration of the current target vertical well, and determine the interface of each marker layer; a third determining subunit, used to determine the intersection point between the current target vertical well and the interface of each marker layer; and a fourth determining subunit, used to determine the layer velocity value at each intersection point based on the sonic logging curve data of the current target vertical well.

[0206] In some embodiments, the predetermined condition includes the accuracy of the velocity trend surface of the second velocity layer reaching a predetermined accuracy threshold; wherein the accuracy of the velocity trend surface of the second velocity layer is calculated by the following means: determining the initial CIP gather in the depth domain based on the logging data of the target vertical well; and determining the accuracy of the velocity trend surface of the second velocity layer by the flattening degree of the in-phase axis after dynamic correction of the CIP gather.

[0207] In some embodiments, the apparatus further includes: a calculation unit, configured to calculate the sum of squares of the differences between the actual values ​​of each layer velocity value in the set of layer velocity values ​​and the fitted values ​​corresponding to the actual values, as a residual value; wherein the fitted value refers to the layer velocity value corresponding to each velocity trend surface in the second target velocity body; and a sixth determination unit, configured to use the residual value as part of the second target velocity body of the second velocity layer.

[0208] In some embodiments, the apparatus further includes a seventh determining unit, configured to establish a first target velocity body of the first velocity layer using micrologging data and first arrival data from the cannon.

[0209] In some embodiments, before performing polynomial fitting on the layer velocity values ​​in the set of layer velocity values ​​to determine the velocity trend surface of the second velocity layer, the method includes: an eighth determining unit, configured to determine an initial velocity volume of pre-stack depth migration based on the pre-stack time migration velocity spectrum of the dual complex structural region; a replacement unit, configured to replace a portion of the first velocity layer in the initial velocity volume with a first target velocity volume of the first velocity layer; a ninth determining unit, configured to determine the intersection points of the longitudinal and transverse target lines of the replaced second target velocity volume with the layer boundary, and determine the layer velocity value of the intersection points based on the initial velocity volume; and a tenth determining unit, configured to use the layer velocity value of the intersection points as an element of the set of layer velocity values.

[0210] For a detailed description of the above-mentioned device and its beneficial effects, please refer to the corresponding method embodiments, which will not be repeated here.

[0211] This invention also provides an electronic device, such as... Figure 9 As shown, the electronic device may include a processor 901 and a memory 902, wherein the processor 901 and the memory 902 may be connected via a bus or other means. Figure 9 Taking the example of a connection between China and Israel via a bus.

[0212] Processor 901 can be a Central Processing Unit (CPU). Processor 901 can also be other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, or combinations of the above types of chips.

[0213] Memory 902, as a non-transitory computer-readable storage medium, can be used to store non-transitory software programs, non-transitory computer-executable programs, and modules, such as the program instructions / modules corresponding to the pre-stack depth migration velocity modeling method for dual complex structural regions in this embodiment of the invention (e.g., Figure 8 (Iteration unit 40 shown). The processor 901 executes various functional applications and data classification by running non-transient software programs, instructions and modules stored in the memory 902, that is, it implements the pre-stack depth migration velocity modeling method for dual complex structural regions in the above method embodiment.

[0214] The memory 902 may include a program storage area and a data storage area. The program storage area may store the operating system and applications required for at least one function; the data storage area may store data created by the processor 901, etc. Furthermore, the memory 902 may include high-speed random access memory and may also include non-transitory memory, such as at least one disk storage device, flash memory device, or other non-transitory solid-state storage device. In some embodiments, the memory 902 may optionally include memory remotely located relative to the processor 901, and these remote memories may be connected to the processor 901 via a network. Examples of such networks include, but are not limited to, the Internet, corporate intranets, local area networks, mobile communication networks, and combinations thereof.

[0215] The one or more modules are stored in the memory 902, and when executed by the processor 901, the above-described pre-stack depth migration velocity modeling method for dual complex structural regions is performed.

[0216] For details regarding the specific electronic devices described above, please refer to the relevant descriptions and effects in the corresponding embodiments for further understanding; they will not be repeated here.

[0217] The systems, devices, modules, or units described in the above embodiments can be implemented by computer chips or entities, or by products with certain functions.

[0218] For ease of description, the above devices are described separately by function as various units. Of course, in implementing this application, the functions of each unit can be implemented in one or more software and / or hardware.

[0219] While this specification provides method operation steps as described in the embodiments or flowcharts, more or fewer operation steps may be included based on conventional or non-inventive means. The order of steps listed in the embodiments is merely one possible execution order among many and does not represent the only execution order. In actual device or end product execution, the methods shown in the embodiments or drawings may be executed sequentially or in parallel (e.g., in a parallel processor or multi-threaded processing environment, or even a distributed data processing environment). The terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, product, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, product, or apparatus. Without further limitations, the presence of other identical or equivalent elements in the process, method, product, or apparatus that includes said elements is not excluded.

[0220] The embodiments described in this specification can be described in the general context of computer-executable instructions, such as program modules, that are executed by a computer. Generally, program modules include routines, programs, objects, components, data structures, etc., that perform a specific task or implement a specific abstract data type. The embodiments of this specification can also be practiced in distributed computing environments where tasks are performed by remote processing devices connected via a communication network. In distributed computing environments, program modules can reside in local and remote computer storage media, including storage devices.

[0221] The above description is merely an embodiment of the present specification and is not intended to limit the embodiments of the present specification. For those skilled in the art, various modifications and variations can be made to the embodiments of the present specification. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principle of the embodiments of the present specification should be included within the scope of the claims of the embodiments of the present specification.

Claims

1. A method for pre-stack depth migration velocity modeling in a dual-complicated structure area, characterized in that, The method comprises the following steps: determining the interface between the first velocity layer and the second velocity layer on each microlog according to the microlog data of the double complex structure area; the depth of the second velocity layer is greater than the depth of the first velocity layer; determining a plurality of interval velocity values on each target straight well according to the logging data of each target straight well, and taking the interval velocity values as elements in an interval velocity value set; polynomial fitting is performed on the interval velocity values in the interval velocity value set to determine the velocity trend surface of the second velocity layer; the following method is repeatedly performed until the second target velocity body of the second velocity layer is obtained after a predetermined condition is reached: interpolation is performed between two adjacent interval velocity trend surfaces to obtain a new velocity trend surface; the new velocity trend surface is iteratively corrected by tomographic inversion; the method further comprises the following steps: establishing a first target velocity body of the first velocity layer through microlog data and big-gun first-arrival data; before determining the velocity trend surface of the second velocity layer based on the interval velocity value set, the following steps are included: determining an initial velocity body for pre-stack depth migration according to the pre-stack time migration velocity spectrum of the double complex structure area; replacing part of the first velocity layer in the initial velocity body with the first target velocity body of the first velocity layer; determining the intersection of the longitudinal and horizontal target lines of the second target velocity body after replacement and the horizon interface, and determining the interval velocity value of the intersection point according to the initial velocity body; 2. The method of claim 1, wherein, taking the interval velocity value of the intersection point as an element of the interval velocity value set. determining a plurality of interval velocity values on each target straight well according to the logging data of each target straight well comprises the following steps: distinguishing sandstone and mudstone marker layers according to the GR curve value of the current target straight well; core homing and lithology calibration are performed on the current target straight well, and the interface of each marker layer is determined; determining the intersection of the current target straight well and the interface of each marker layer; 3. The method of claim 1, wherein, determining the interval velocity value at each intersection point according to the acoustic logging curve data of the current target straight well. The predetermined condition includes that the accuracy of the velocity trend surface of the second velocity layer reaches a predetermined accuracy threshold; wherein the accuracy of the velocity trend surface of the second velocity layer is calculated by the following method: determining the initial CIP gather in the depth domain according to the logging data of the target straight well; 4. The method of claim 1, wherein, determining the accuracy of the velocity trend surface of the second velocity layer according to the flattening degree of the same-phase axis after CIP gather dynamic correction. After determining the second target velocity body of the second velocity layer, the following steps are further included: calculating the sum of squares of the difference between the actual value of each interval velocity value in the interval velocity value set and the fitting value corresponding to the actual value as a residual value; wherein the fitting value refers to the corresponding interval velocity value on each velocity trend surface in the second target velocity body; 5. A device for pre-stack depth migration velocity modeling in a dual complex structure area, characterized in that, taking the residual value as part of the second target velocity body of the second velocity layer. The method comprises the following steps: a first determining unit is configured to determine the interface between the first velocity layer and the second velocity layer on each microlog according to the microlog data of the double complex structure area; the depth of the second velocity layer is greater than the depth of the first velocity layer; The second determining unit is configured to determine a plurality of interval velocity values on each target straight well according to logging data of the target straight well, and take the interval velocity values as elements in an interval velocity value set; The fitting unit is configured to perform polynomial fitting on the interval velocity values in the interval velocity value set to determine a velocity trend surface of the second velocity layer; The iteration unit is configured to repeatedly perform the following method until a second target velocity body of the second velocity layer is obtained after a predetermined condition is reached: interpolating between two adjacent interval velocity trend surfaces to obtain a new velocity trend surface; and iteratively correcting the new velocity trend surface by tomographic inversion; The device further comprises a seventh determining unit configured to establish a first target velocity body of the first velocity layer by using micro-logging data and big-gun first-arrival data; The device further comprises: An eighth determining unit configured to determine an initial velocity body of pre-stack depth migration according to a pre-stack time migration velocity spectrum of a double-complex structure area; A replacing unit configured to replace a part of the first velocity layer in the initial velocity body with the first target velocity body of the first velocity layer; A ninth determining unit configured to determine an intersection of a vertical-horizontal target line of the second target velocity body after replacement and a horizon interface, and determine an interval velocity value of the intersection according to the initial velocity body; A tenth determining unit configured to take the interval velocity value of the intersection as an element in the interval velocity value set.

6. An electronic device, comprising: The device further comprises: A memory and a processor, which are communicatively connected, and the memory stores computer instructions, and the processor implements the steps of the method according to any one of claims 1 to 4 by executing the computer instructions.

7. A computer storage medium, characterized in that The computer storage medium stores computer program instructions, and the computer program instructions are executed to implement the steps of the method according to any one of claims 1 to 4.

Images