Anisotropic seismic migration imaging method, device, electronic equipment and medium

Abstract

The present application discloses an anisotropic seismic migration imaging method, device, electronic device and medium. The method may include: establishing an observation system, reading velocity model and migration parameters; using variable space step size grid partitioning according to velocity model and migration parameters, calculating forward simulation differential coefficients; calculating the source wave field and detection point wave field of each shot at each moment based on the anisotropic wave equation, grid distribution information and forward simulation differential coefficients; performing cross-correlation operation on the source wave field and detection point wave field of the same moment of the shot to obtain a single shot migration result; superimposing all single shot migration results and outputting the final migration result. The grid boundary partitioned by the present invention is highly consistent with the real irregular interface, effectively eliminating the step effect of regular grid seismic wave numerical simulation, with high computational efficiency, effectively reducing the computational complexity of conventional anisotropic migration imaging, and obtaining a high-precision seismic wave field.

Description

Technical Field

[0001] The present invention relates to the field of seismology, and more specifically, to an anisotropic seismic migration imaging method, device, electronic equipment and medium. Background Art

[0002] In the process of complex oil and gas reservoir exploration, seismic wave migration imaging is an important technical means in seismic exploration. The important role of migration imaging is that it can restore the spatial structure and distribution of the underground. Through seismic wave migration imaging, geological areas rich in oil and gas can be efficiently analyzed. As people's understanding of the kinematic and dynamic propagation laws of seismic waves gradually increases, the imaging of anisotropic media has gradually become the research focus of seismic imaging. Anisotropic seismic waves refer to the property that a certain attribute of a medium changes with direction. The anisotropy problem in seismic exploration mainly refers to the change in the velocity of seismic waves propagating in the stratum as the direction of its propagation changes. In three-dimensional exploration with large offset distances, if it is not taken into account, it will seriously affect the accurate imaging and effective use of seismic data.

[0003] Finite difference method is usually used to achieve anisotropic migration imaging. This technology uses a computer to simulate the propagation of seismic waves in underground media using theoretical calculation methods. There are many wave field simulation methods. Although the finite difference method has solved the anisotropic seismic data migration imaging of large geological structures, the global grid division used cannot meet the needs of industrial exploration. For low-speed surfaces or interlayers, the global grid is not fine enough, which will lead to problems such as numerical dispersion and reduced resolution. If a smaller global grid is used, it will lead to memory requirements that are difficult for current computing clusters to accept, and at the same time reduce computing efficiency. When a rectangular grid is used to grid the model, the undulating surface and the underground curved interface can only be replaced by a series of stepped grids, which will cause artificial false scattering and diffraction waves, directly affecting the accuracy of migration imaging.

[0004] Therefore, it is necessary to develop an anisotropic seismic migration imaging method, device, electronic equipment and medium based on variable space grid division with high accuracy and strong stability.

[0005] The information disclosed in the background technology section of the present invention is only intended to deepen the understanding of the general background technology of the present invention, and should not be regarded as acknowledging or suggesting in any form that the information constitutes the prior art already known to those skilled in the art. Summary of the invention

[0006] The present invention proposes an anisotropic seismic migration imaging method, device, electronic equipment and medium, wherein the subdivided grid boundary is highly consistent with the real irregular interface, effectively eliminating the step effect of regular grid seismic wave numerical simulation, having high computational efficiency, effectively reducing the computational complexity of conventional anisotropic migration imaging, and obtaining a high-precision seismic wave field.

[0007] In a first aspect, an embodiment of the present disclosure provides an anisotropic seismic migration imaging method, comprising:

[0008] According to the work area location information, establish an observation system, determine the multi-shot seismic data used for migration, and read the velocity model and migration parameters;

[0009] According to the velocity model and the migration parameters, variable space step size grid generation is adopted to output three-dimensional space grid distribution information, and then forward simulation differential coefficients are calculated;

[0010] According to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficients, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method;

[0011] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single shot migration result;

[0012] All single-shot offset results are superimposed and the final offset result is output.

[0013] Preferably, according to the velocity model and the offset parameter, variable space step size gridding is adopted to output three-dimensional space grid distribution information, including:

[0014] Determining the minimum value of the grid space step length using the dispersion relation of the anisotropic wave equation according to the velocity model and the migration parameters;

[0015] Determining the maximum value of the grid space step length using the stability condition of the anisotropic wave equation according to the velocity model and the migration parameters;

[0016] Determining a range of grid space step lengths to which the seismic wave velocity model is adapted according to the minimum value of the grid space step length and the maximum value of the grid space step length;

[0017] Determining, according to the grid space step range, a mean value of the grid space step adapted to each speed value;

[0018] Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted by each speed value;

[0019] According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

[0020] Preferably, calculating the forward modeling differential coefficients includes:

[0021] Establishing a relationship between the spatial partial derivative and the forward modeling differential coefficient;

[0022] Establishing a simplified expression of the forward modeling differential coefficient;

[0023] The spatial partial derivatives are solved to calculate the forward modeling differential coefficients.

[0024] Preferably, the relationship between the spatial partial derivative and the forward modeling differential coefficient is:

[0025]

[0026] Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, and a(x) is the differential coefficient.

[0027] Preferably, the simplified expression of the forward modeling differential coefficient is:

[0028] a(x)=A -1 (x)B(x)U (2)

[0029] in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

[0030] Preferably, according to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, calculating the source wave field and the receiver wave field of each shot at each time by the finite difference method includes:

[0031] Solving the anisotropic source wave equation using the finite difference method according to the grid distribution information and the differential coefficients to obtain the source wave field of each shot at each moment;

[0032] The anisotropic detection point wave equation is solved by using the finite difference method according to the grid distribution information and the differential coefficient to obtain the detection point wave field of each shot at each moment.

[0033] Preferably, applying imaging conditions to the two wave fields of the shot to obtain a single shot migration result comprises:

[0034] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

[0035] As a specific implementation method of the embodiment of the present disclosure,

[0036] In a second aspect, the present disclosure also provides an anisotropic seismic migration imaging device, including:

[0037] The reading module establishes an observation system based on the work area location information, determines the multi-shot seismic data used for migration, and reads the velocity model and migration parameters;

[0038] A subdivision calculation module, which uses variable space step size grid subdivision according to the velocity model and the migration parameters, outputs three-dimensional space grid distribution information, and then calculates the forward simulation differential coefficient;

[0039] A finite difference module calculates the source wave field and the receiver wave field of each shot at each time by a finite difference method according to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficients;

[0040] The cross-correlation module performs cross-correlation calculations on the source wave field and the detection point wave field at the same time of the shot, and sums the cross-correlation results to obtain the single-shot migration result;

[0041] The superposition module superimposes all single-shot offset results and outputs the final offset result.

[0042] Preferably, according to the velocity model and the offset parameter, variable space step size gridding is adopted to output three-dimensional space grid distribution information, including:

[0043] Determining the minimum value of the grid space step length using the dispersion relation of the anisotropic wave equation according to the velocity model and the migration parameters;

[0044] Determining the maximum value of the grid space step length using the stability condition of the anisotropic wave equation according to the velocity model and the migration parameters;

[0045] Determining a range of grid space step lengths to which the seismic wave velocity model is adapted according to the minimum value of the grid space step length and the maximum value of the grid space step length;

[0046] Determining, according to the grid space step range, a mean value of the grid space step adapted to each speed value;

[0047] Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted by each speed value;

[0048] According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

[0049] Preferably, calculating the forward modeling differential coefficients includes:

[0050] Establishing a relationship between the spatial partial derivative and the forward modeling differential coefficient;

[0051] Establishing a simplified expression of the forward modeling differential coefficient;

[0052] The spatial partial derivatives are solved to calculate the forward modeling differential coefficients.

[0053] Preferably, the relationship between the spatial partial derivative and the forward modeling differential coefficient is:

[0054]

[0055] Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, and a(x) is the differential coefficient.

[0056] Preferably, the simplified expression of the forward modeling differential coefficient is:

[0057] a(x)=A -1 (x)B(x)U (2)

[0058] in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

[0059] Preferably, according to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, calculating the source wave field and the receiver wave field of each shot at each time by the finite difference method includes:

[0060] Solving the anisotropic source wave equation using the finite difference method according to the grid distribution information and the differential coefficients to obtain the source wave field of each shot at each moment;

[0061] The anisotropic detection point wave equation is solved by using the finite difference method according to the grid distribution information and the differential coefficient to obtain the detection point wave field of each shot at each moment.

[0062] Preferably, applying imaging conditions to the two wave fields of the shot to obtain a single shot migration result comprises:

[0063] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

[0064] In a third aspect, an embodiment of the present disclosure further provides an electronic device, the electronic device comprising:

[0065] A memory storing executable instructions;

[0066] A processor runs the executable instructions in the memory to implement the anisotropic seismic migration imaging method.

[0067] In a fourth aspect, an embodiment of the present disclosure further provides a computer-readable storage medium, which stores a computer program, and when the computer program is executed by a processor, the anisotropic seismic migration imaging method is implemented.

[0068] Its beneficial effects are:

[0069] 1) The mesh boundary is highly consistent with the real irregular interface, which can effectively suppress the false wave field generated by the step discretization of the rectangular grid;

[0070] 2) By constraining the variable space step size grid generation based on the dispersion relation and stability condition of the anisotropic wave equation, and comprehensively considering the vertical sampling and lateral oversampling problems, the computational efficiency is improved and the memory usage is reduced compared with the conventional forward modeling algorithm;

[0071] 3) Based on the robust variable spatial step size gridding, high-precision source wave field and receiver wave field information can be obtained. The cross-correlation imaging condition is used to effectively improve the amplitude preservation characteristics of anisotropic migration imaging, providing reliable technical support for energy exploration and development.

[0072] The methods and apparatus of the present invention have other features and advantages that will be apparent from, or will be described in detail in, the accompanying drawings and subsequent detailed descriptions incorporated herein, which together serve to explain the specific principles of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0073] The above and other objects, features and advantages of the present invention will become more apparent through a more detailed description of exemplary embodiments of the present invention in conjunction with the accompanying drawings, wherein like reference numerals generally represent like components throughout the exemplary embodiments of the present invention.

[0074] Figure 1 A flow chart showing the steps of an anisotropic seismic migration imaging method according to an embodiment of the present invention.

[0075] Figure 2a , Figure 2b , Figure 2c , Figure 2d The v and v of the standard BP velocity model according to an embodiment of the present invention are shown respectively. p Schematic diagram of the model, ε model, δ model, and θ model.

[0076] Figure 3a , Figure 3bSchematic diagrams of superimposed migration sections of a standard BP model obtained by a conventional method and the present method according to an embodiment of the present invention are respectively shown.

[0077] Figure 4a , Figure 4b , Figure 4c , Figure 4d The typical carbonate reservoir model according to one embodiment of the present invention is shown in FIG. p Schematic diagram of the model, ε model, δ model, and θ model.

[0078] Figure 5a , Figure 5b Schematic diagrams showing snapshots of the seismic wave field obtained by a conventional method at 1.5s and 2.25s according to an embodiment of the present invention.

[0079] Figure 6a , Figure 6b Schematic diagrams of seismic wave field snapshots obtained by the method at 1.5s and 2.25s according to an embodiment of the present invention are respectively shown.

[0080] Figure 7a , Figure 7b Schematic diagrams of the migration profiles of a typical carbonate reservoir model obtained by a conventional method and the present method according to an embodiment of the present invention are respectively shown.

[0081] Figure 8 A block diagram of an anisotropic seismic migration imaging device according to an embodiment of the present invention is shown.

[0082] Description of reference numerals:

[0083] 201, reading module; 202, segmentation and calculation module; 203, finite difference module; 204, cross-correlation module; 205, superposition module. DETAILED DESCRIPTION

[0084] The preferred embodiments of the present invention will be described in more detail below. Although the preferred embodiments of the present invention are described below, it should be understood that the present invention can be implemented in various forms and should not be limited to the embodiments set forth herein.

[0085] The present invention provides an anisotropic seismic migration imaging method, comprising:

[0086] According to the location information of the work area, an observation system is established to determine the multi-shot seismic data used for migration, and the velocity model and migration parameters are read; the velocity model and migration parameters include: the longitudinal wave velocity model v p , shear wave velocity model v s , anisotropy parameters ε, δ, inclination model θ, azimuth Offset aperture, time sampling rate dt.

[0087] According to the velocity model and the migration parameters, variable space step size grid generation is adopted to output the three-dimensional space grid distribution information, and then the forward simulation differential coefficient is calculated; in one example, according to the velocity model and the migration parameters, variable space step size grid generation is adopted to output the three-dimensional space grid distribution information including:

[0088] The minimum value of the grid space step is determined based on the velocity model and migration parameters using the dispersion relation of the anisotropic wave equation;

[0089] The maximum value of the grid space step size is determined based on the velocity model and migration parameters using the stability condition of the anisotropic wave equation;

[0090] According to the minimum value of the grid space step length and the maximum value of the grid space step length, the grid space step length range adapted to the seismic wave velocity model is determined;

[0091] According to the range of the grid space step length, determine the average value of the grid space step length that each speed value is adapted to;

[0092] Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted to each speed value;

[0093] According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

[0094] In one example, calculating the forward modeling differential coefficients includes:

[0095] Establish the relationship between the spatial partial derivative and the differential coefficient of forward modeling;

[0096] Establish simplified expressions for the differential coefficients of forward modeling;

[0097] The spatial partial derivatives are solved and then the forward modeling differential coefficients are calculated.

[0098] In one example, the spatial partial derivatives are related to the forward modeling differential coefficients as follows:

[0099]

[0100] Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, and a(x) is the differential coefficient.

[0101] In one example, the simplified expression for the forward modeling differential coefficient is:

[0102] a(x)=A -1(x)B(x)U (2)

[0103] in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

[0104] Specifically, the dispersion condition of the anisotropic wave equation is:

[0105]

[0106] Where v(θ) is the phase velocity of longitudinal waves in anisotropic media, According to formula (3), the minimum value of the grid step size can be obtained.

[0107] The stability condition for the anisotropic wave equation is:

[0108]

[0109] Among them, dx, dy and dz are the grid space steps, v(x, y, z) is the velocity of seismic wave propagation, Δt is the time step, and a m is the difference coefficient, N is the difference order, j is the number of grids from the difference center point, and the maximum value of the grid step size can be obtained according to formula (4).

[0110] make is the grid distribution information, x = (x, y, z). Given the spatial wave field u = u (x, t) of the anisotropic medium wave equation, it represents the seismic wave field in the time-space domain, which can also be expressed as u = u (x, y, z, t). Given the positive definite function Then define the spatial partial derivative The approximate solution is formula (1), which is the relationship between the spatial partial derivative and the differential coefficient of the forward simulation. a(x) is the differential coefficient, which is related to the scope of the influence domain, the number of nodes in the domain, the weight coefficient, etc. Assume that the grid step size of the influence domain around the center point x is R, N is the number of nodes in the whole domain, and the simplified expression of the differential coefficient a(x) is formula (2).

[0111] Spatial partial derivatives The approximate solution function is expressed as:

[0112]

[0113] in, Through the above steps, the differential coefficient data a(x) of each grid can be obtained.

[0114] According to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficient, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method; In one example, according to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method, including:

[0115] According to the grid distribution information and differential coefficients, the finite difference method is used to solve the anisotropic source wave equation to obtain the source wave field of each shot at each time.

[0116] The finite difference method is used to solve the wave equation of anisotropic receiver points according to the grid distribution information and the differential coefficients to obtain the wave field of the receiver points at each shot and at each moment.

[0117] Specifically, according to the anisotropic wave equation:

[0118]

[0119] Among them, the partial differential operator H 1 and H 2 They are

[0120]

[0121] S is a scalar operator,

[0122] The anisotropic source wave equation is:

[0123]

[0124] Among them, u s =u s (x, y, z, t) is the source wave field at the position (x, y, z) at time t, f s is the earthquake source function.

[0125] The wave equation of anisotropic detection point is:

[0126]

[0127] Among them, u r =u r (x, y, z, t) is the wave field of the detection point at the position (x, y, z) at time t, D(x r ,t;x s ) is the multi-shot seismic data used for migration, x r =(x r ,y r ,z r ) is the detection point position, xs =(x s ,y s ,z s ) is the location of the earthquake source.

[0128] The finite difference method is used to solve the anisotropic source wave equation and the receiver point wave equation to obtain the source wave field and the receiver point wave field at each shot and each moment.

[0129] The source wave field and the receiver wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result. In one example, imaging conditions are applied to the two wave fields of the shot to obtain the single-shot migration result, including:

[0130] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

[0131] Specifically, the cross-correlation imaging condition is:

[0132] I(x,y,z)=∫u s (x,y,z,t)u r (x,y,z,t)dt (9)

[0133] Among them, I(x,y,z) is the single shot offset result.

[0134] All single-shot offset results are superimposed and the final offset result is output.

[0135] The present invention also provides an anisotropic seismic migration imaging device, comprising:

[0136] The reading module establishes an observation system based on the work area location information, determines the multi-shot seismic data used for migration, and reads the velocity model and migration parameters; the velocity model and migration parameters include: P-wave velocity model v p , shear wave velocity model v s , anisotropy parameters ε, δ, inclination model θ, azimuth Offset aperture, time sampling rate dt.

[0137] The subdivision calculation module uses variable space step size grid subdivision according to the velocity model and the migration parameters, outputs the three-dimensional space grid distribution information, and then calculates the forward simulation differential coefficient; in one example, according to the velocity model and the migration parameters, the variable space step size grid subdivision is used to output the three-dimensional space grid distribution information including:

[0138] The minimum value of the grid space step is determined based on the velocity model and migration parameters using the dispersion relation of the anisotropic wave equation;

[0139] The maximum value of the grid space step size is determined based on the velocity model and migration parameters using the stability condition of the anisotropic wave equation;

[0140] According to the minimum value of the grid space step length and the maximum value of the grid space step length, the grid space step length range adapted to the seismic wave velocity model is determined;

[0141] According to the range of the grid space step length, determine the average value of the grid space step length that each speed value is adapted to;

[0142] Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted to each speed value;

[0143] According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

[0144] In one example, calculating the forward modeling differential coefficients includes:

[0145] Establish the relationship between the spatial partial derivative and the differential coefficient of forward modeling;

[0146] Establish simplified expressions for the differential coefficients of forward modeling;

[0147] The spatial partial derivatives are solved and then the forward modeling differential coefficients are calculated.

[0148] In one example, the spatial partial derivatives are related to the forward modeling differential coefficients as follows:

[0149]

[0150] Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, and a(x) is the differential coefficient.

[0151] In one example, the simplified expression for the forward modeling differential coefficient is:

[0152] a(x)=A -1 (x)B(x)U (2)

[0153] in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

[0154] Specifically, the dispersion condition of the anisotropic wave equation is given by formula (3), and the minimum value of the grid step size can be obtained according to formula (3).

[0155] The stability condition of the anisotropic wave equation is given by formula (4), and the maximum value of the grid step size can be obtained according to formula (4).

[0156] make is the grid distribution information, x = (x, y, z). Given the spatial wave field u = u (x, t) of the anisotropic medium wave equation, it represents the seismic wave field in the time-space domain, which can also be expressed as u = u (x, y, z, t). Given the positive definite function Then define the spatial partial derivative The approximate solution is formula (1), which is the relationship between the spatial partial derivative and the differential coefficient of the forward simulation. a(x) is the differential coefficient, which is related to the scope of the influence domain, the number of nodes in the domain, the weight coefficient, etc. Assume that the grid step size of the influence domain around the center point x is R, N is the number of nodes in the whole domain, and the simplified expression of the differential coefficient a(x) is formula (2).

[0157] Spatial partial derivatives The approximate solution function is expressed as formula (5). Through the above steps, the differential coefficient data a(x) of each grid can be obtained.

[0158] The finite difference module calculates the source wave field and the receiver wave field of each shot at each time by the finite difference method according to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficient; In one example, according to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, the source wave field and the receiver wave field of each shot at each time by the finite difference method are calculated, including:

[0159] According to the grid distribution information and differential coefficients, the finite difference method is used to solve the anisotropic source wave equation to obtain the source wave field of each shot at each time.

[0160] The finite difference method is used to solve the wave equation of anisotropic receiver points according to the grid distribution information and the differential coefficients to obtain the wave field of the receiver points at each shot and at each moment.

[0161] Specifically, the anisotropic wave equation is formula (6). The anisotropic source wave equation is formula (7), and the anisotropic receiver wave equation is formula (8). The finite difference method is used to solve the anisotropic source wave equation and the receiver wave equation to obtain the source wave field and the receiver wave field at each shot and each time.

[0162] The cross-correlation module performs cross-correlation calculation on the source wave field and the receiver wave field at the same time of the shot, and sums the cross-correlation results to obtain the single-shot migration result. In one example, the imaging conditions are applied to the two wave fields of the shot to obtain the single-shot migration result, including:

[0163] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

[0164] Specifically, the cross-correlation imaging condition is formula (9).

[0165] The superposition module superimposes all single-shot offset results and outputs the final offset result.

[0166] The present invention also provides an electronic device, which includes: a memory storing executable instructions; and a processor, which runs the executable instructions in the memory to implement the above-mentioned anisotropic seismic migration imaging method.

[0167] The present invention also provides a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, the above-mentioned anisotropic seismic migration imaging method is implemented.

[0168] To facilitate understanding of the solutions and effects of the embodiments of the present invention, four specific application examples are given below. Those skilled in the art should understand that the examples are only for facilitating understanding of the present invention, and any specific details thereof are not intended to limit the present invention in any way.

[0169] Example 1

[0170] Figure 1 A flow chart showing the steps of an anisotropic seismic migration imaging method according to an embodiment of the present invention.

[0171] like Figure 1 As shown, the anisotropic seismic migration imaging method comprises: step 101, establishing an observation system according to the work area location information, determining the multi-shot seismic data used for migration, and reading the velocity model and migration parameters; step 102, using variable space step size grid division according to the velocity model and the migration parameters, outputting three-dimensional space grid distribution information, and then calculating the forward simulation differential coefficient; step 103, according to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficient, calculating the source wave field and the receiver wave field of each shot at each time by the finite difference method; step 104, performing cross-correlation operation on the source wave field and the receiver wave field of the same shot at the same time, and summing the cross-correlation results to obtain the single-shot migration result; step 105, superimposing all the single-shot migration results, and outputting the final migration result.

[0172] Figure 2a , Figure 2b , Figure 2c , Figure 2d The v and v of the standard BP velocity model according to an embodiment of the present invention are shown respectively. pSchematic diagram of the model, ε model, δ model, and θ model. The model has a depth of 11.25km and a lateral width of 94.05km. The spatial grid size used for the offset is 12.5m, with a total of 1641 shots. The starting shot point is located at the starting point of the model, the shot interval is 50m, and the single-sided shot is fired. Each shot is received by 800 channels, the minimum offset distance is 37.5m, the maximum offset distance is 10025m, the channel spacing is 12.5m, the recording time length is 9.2s, the time step is 2ms, and the Ricker wavelet with a main frequency of 20Hz is used as the source time function.

[0173] The mesh boundary of the present invention is highly consistent with the real irregular interface, so this meshing method can effectively suppress the false wave field generated by the step discretization of the rectangular mesh. In addition, the present invention uses the dispersion relationship and stability conditions of the anisotropic wave equation to constrain the variable space step size meshing, comprehensively considers the vertical sampling and lateral oversampling problems, and improves the calculation efficiency and reduces the memory usage compared with the conventional forward simulation algorithm.

[0174] Figure 3a , Figure 3b Schematic diagrams of superimposed migration sections of a standard BP model obtained by a conventional method and the present method according to an embodiment of the present invention are respectively shown.

[0175] Depend on Figure 3a , Figure 3b It can be seen that the anisotropic migration method based on robust variable space step size grid subdivision of the present invention is correct. The present invention can obtain high-precision source wave field and detection point wave field information, and adopts cross-correlation imaging conditions to effectively improve the amplitude preservation characteristics of anisotropic migration imaging. In addition, compared with the migration obtained by conventional methods, the migration method obtained by the method of the present invention can better image the complex underground structure, the diffraction wave energy converges more fully, and there are fewer illusions.

[0176] Figure 4a , Figure 4b , Figure 4c , Figure 4d The typical carbonate reservoir model according to one embodiment of the present invention is shown in FIG. p Schematic diagram of the model, ε model, δ model, and θ model.

[0177] The model contains approximately horizontal layered media in the shallow layer to simulate sedimentary strata; the model contains thrust faults in the middle part to depict complex underground structures; the model below 5000m contains multiple groups of fractures of different scales and inclinations, as well as multiple irregular pores and caves filled with oil, gas, water, and rock debris, etc., to depict the geology and geophysics of typical carbonate reservoirs. This model test verifies the effectiveness of the robust variable space step gridding anisotropic migration method in complex carbonate reservoir imaging. The model depth is 9km and the lateral width is 20.8km. The spatial grid size used for the offset is 10m, with a total of 139 shots. The starting shot point is located 50m from the model. The shot points are placed on the surface with a shot interval of 150m. Shots are fired in the middle and received on both sides. Each shot is received with 901 channels. The minimum offset distance is 0m, the maximum offset distance is 9000m, the channel spacing is 10m, the recording time length is 6.0s, the time step is 1ms, and the Ricker wavelet with a main frequency of 30Hz is used as the source time function.

[0178] Figure 5a , Figure 5b Schematic diagrams showing snapshots of the seismic wave field obtained by a conventional method at 1.5s and 2.25s according to an embodiment of the present invention.

[0179] Figure 6a , Figure 6b Schematic diagrams of seismic wave field snapshots obtained by the method at 1.5s and 2.25s according to an embodiment of the present invention are respectively shown.

[0180] Depend on Figure 5a , Figure 5b , Figure 6a , Figure 6b It can be seen from the comparison of multiple groups of seismic wave field snapshots that the anisotropic wave equation of the present invention can simulate the propagation characteristics of seismic waves in anisotropic media more accurately than the wave equation of the conventional method, and its amplitude and phase properties are consistent with the propagation law of seismic waves. In addition, it can be seen from the figure that in the cracks with strong anisotropy, the diffraction waves are obviously developed and the wave field is disordered. These phenomena are also consistent with the actual situation, indicating the correctness and effectiveness of the anisotropic wave equation based on variable space step size grid partitioning of the present invention.

[0181] Figure 7a , Figure 7b Schematic diagrams of the migration profiles of a typical carbonate reservoir model obtained by a conventional method and the present method according to an embodiment of the present invention are respectively shown.

[0182] Depend on Figure 7a , Figure 7bIt can be seen that the migration obtained by the method of the present invention can better accurately image complex structures such as underground fractures and caves. For fractures and caves of different scales, their top and bottom interfaces can be clearly portrayed, and the cross-section of the fault is clear, the breakpoint position is accurate, the diffraction energy converges well, and there are fewer migration artifacts.

[0183] Example 2

[0184] Figure 8 A block diagram of an anisotropic seismic migration imaging device according to an embodiment of the present invention is shown.

[0185] like Figure 8 As shown, the anisotropic seismic migration imaging device comprises:

[0186] The reading module 201 establishes an observation system according to the work area location information, determines the multi-shot seismic data used for migration, and reads the velocity model and migration parameters;

[0187] The meshing calculation module 202 uses variable space step size meshing according to the velocity model and the migration parameters, outputs the three-dimensional space mesh distribution information, and then calculates the forward simulation differential coefficients;

[0188] The finite difference module 203 calculates the source wave field and the receiver wave field at each shot at each time by the finite difference method according to the anisotropic wave equation, based on the grid distribution information and the forward simulation difference coefficients;

[0189] Cross-correlation module 204 performs cross-correlation operation on the source wave field and the detection point wave field at the same time of the shot, and sums the cross-correlation results to obtain the single shot migration result;

[0190] The superposition module 205 superimposes all the single shot migration results and outputs the final migration result.

[0191] As an optional solution, variable space step size gridding is used according to the velocity model and migration parameters, and the output of three-dimensional space grid distribution information includes:

[0192] The minimum value of the grid space step is determined based on the velocity model and migration parameters using the dispersion relation of the anisotropic wave equation;

[0193] The maximum value of the grid space step size is determined based on the velocity model and migration parameters using the stability condition of the anisotropic wave equation;

[0194] According to the minimum value of the grid space step length and the maximum value of the grid space step length, the grid space step length range adapted to the seismic wave velocity model is determined;

[0195] According to the range of the grid space step length, determine the average value of the grid space step length that each speed value is adapted to;

[0196] Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted to each speed value;

[0197] According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

[0198] As an option, calculate the forward modeling differential coefficients including:

[0199] Establish the relationship between the spatial partial derivative and the differential coefficient of forward modeling;

[0200] Establish simplified expressions for the differential coefficients of forward modeling;

[0201] The spatial partial derivatives are solved and then the forward modeling differential coefficients are calculated.

[0202] As an optional solution, the relationship between the spatial partial derivative and the forward modeling differential coefficient is:

[0203]

[0204] Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, and a(x) is the differential coefficient.

[0205] As an alternative, the simplified expression of the forward modeling differential coefficient is:

[0206] a(x)=A -1 (x)B(x)U (2)

[0207] in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

[0208] As an optional solution, according to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method, including:

[0209] According to the grid distribution information and differential coefficients, the finite difference method is used to solve the anisotropic source wave equation to obtain the source wave field of each shot at each time.

[0210] The finite difference method is used to solve the wave equation of anisotropic receiver points according to the grid distribution information and the differential coefficients to obtain the wave field of the receiver points at each shot and at each moment.

[0211] Optionally, imaging conditions are applied to the two wave fields of the shot to obtain single-shot migration results including:

[0212] The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

[0213] Example 3

[0214] The present disclosure provides an electronic device, which includes: a memory storing executable instructions; and a processor, which runs the executable instructions in the memory to implement the above-mentioned anisotropic seismic migration imaging method.

[0215] An electronic device according to an embodiment of the present disclosure includes a memory and a processor.

[0216] The memory is used to store non-temporary computer-readable instructions. Specifically, the memory may include one or more computer program products, which may include various forms of computer-readable storage media, such as volatile memory and / or non-volatile memory. The volatile memory may, for example, include random access memory (RAM) and / or cache memory (cache), etc. The non-volatile memory may, for example, include read-only memory (ROM), hard disk, flash memory, etc.

[0217] The processor may be a central processing unit (CPU) or other forms of processing units having data processing capabilities and / or instruction execution capabilities, and may control other components in the electronic device to perform desired functions. In one embodiment of the present disclosure, the processor is used to run the computer-readable instructions stored in the memory.

[0218] Those skilled in the art should be able to understand that in order to solve the technical problem of how to obtain a good user experience, the present embodiment may also include well-known structures such as a communication bus and an interface, and these well-known structures should also be included in the protection scope of the present disclosure.

[0219] For detailed description of this embodiment, reference may be made to the corresponding descriptions in the aforementioned embodiments, which will not be repeated here.

[0220] Example 4

[0221] An embodiment of the present disclosure provides a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, the anisotropic seismic migration imaging method is implemented.

[0222] According to the computer-readable storage medium of the embodiment of the present disclosure, non-transitory computer-readable instructions are stored thereon. When the non-transitory computer-readable instructions are executed by a processor, all or part of the steps of the above-mentioned methods of each embodiment of the present disclosure are executed.

[0223] The above-mentioned computer-readable storage media include, but are not limited to: optical storage media (e.g., CD-ROM and DVD), magneto-optical storage media (e.g., MO), magnetic storage media (e.g., magnetic tape or mobile hard disk), media with built-in rewritable non-volatile memory (e.g., memory card) and media with built-in ROM (e.g., ROM box).

[0224] Those skilled in the art should understand that the purpose of the above description of the embodiments of the present invention is only to exemplarily illustrate the beneficial effects of the embodiments of the present invention, and is not intended to limit the embodiments of the present invention to any given examples.

[0225] The embodiments of the present invention have been described above, and the above description is exemplary, not exhaustive, and is not limited to the disclosed embodiments. Many modifications and changes will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments.

Claims

1. A method for anisotropic seismic migration imaging, characterized in that: include: According to the work area location information, establish an observation system, determine the multi-shot seismic data used for migration, and read the velocity model and migration parameters; According to the velocity model and the migration parameters, variable space step size grid generation is adopted to output three-dimensional space grid distribution information, and then forward simulation differential coefficients are calculated; According to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficients, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method; The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single shot migration result; All single-shot offset results are superimposed and the final offset result is output; Wherein, according to the velocity model and the offset parameter, variable space step size gridding is adopted to output three-dimensional space grid distribution information, including: The minimum value of the grid space step is determined according to the velocity model and the migration parameters using the dispersion relation of the anisotropic wave equation; wherein the dispersion relation of the anisotropic wave equation is: Where v(θ) is the phase velocity of longitudinal waves in anisotropic media, According to formula (3), the minimum value of the grid space step length, v p is the longitudinal wave velocity, v s is the shear wave velocity, ε and δ are anisotropy parameters, and θ is the inclination angle; The maximum value of the grid space step size is determined according to the velocity model and the migration parameters using the stability condition of the anisotropic wave equation; wherein the stability condition of the anisotropic wave equation is: Among them, dx, dy and dz are the grid space steps, v(x, y, z) is the velocity of seismic wave propagation, Δt is the time step, and a m is the difference coefficient, N is the difference order, j represents the number of grids from the difference center point, and the maximum value of the grid space step length is obtained according to formula (4); Determining a range of grid space step lengths to which the seismic wave velocity model is adapted according to the minimum value of the grid space step length and the maximum value of the grid space step length; Determining, according to the grid space step range, a mean value of the grid space step adapted to each speed value; Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted by each speed value; According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

2. The anisotropic seismic migration imaging method according to claim 1, wherein: Calculation of the forward modeling differential coefficients includes: Establishing a relationship between the spatial partial derivative and the forward modeling differential coefficient; Establishing a simplified expression of the forward modeling differential coefficient; The spatial partial derivatives are solved to calculate the forward modeling differential coefficients.

3. The anisotropic seismic migration imaging method according to claim 2, wherein: The relationship between the spatial partial derivative and the forward modeling differential coefficient is: Among them, u = u (x, t) is the seismic wave field in the time-space domain, p T (x) = [1, x, y, z] is the basis function, m is the number of terms of the basis function, a(x) is the differential coefficient, x is the grid distribution information, and x = (x, y, z).

4. The anisotropic seismic migration imaging method according to claim 2, wherein: The simplified expression of the forward modeling differential coefficient is: a(x)=A -1 (x)B(x)U (2) in, A(x) is the weight function correlation matrix, and B(x) is the basis function correlation matrix.

5. The anisotropic seismic migration imaging method according to claim 1, wherein: According to the anisotropic wave equation, based on the grid distribution information and the differential coefficient, the source wave field and the receiver wave field of each shot at each time are calculated by the finite difference method, including: Solving the anisotropic source wave equation using the finite difference method according to the grid distribution information and the differential coefficients to obtain the source wave field of each shot at each moment; The anisotropic detection point wave equation is solved by using the finite difference method according to the grid distribution information and the differential coefficient to obtain the detection point wave field of each shot at each moment.

6. The anisotropic seismic migration imaging method according to claim 1, wherein: Imaging conditions are applied to the two wave fields of the shot to obtain single-shot migration results including: The source wave field and the detection point wave field of the shot at the same time are cross-correlated, and the cross-correlation results are summed to obtain the single-shot migration result.

7. An anisotropic seismic migration imaging device, characterized in that: include: The reading module establishes an observation system based on the work area location information, determines the multi-shot seismic data used for migration, and reads the velocity model and migration parameters; A subdivision calculation module, which uses variable space step size grid subdivision according to the velocity model and the migration parameters, outputs three-dimensional space grid distribution information, and then calculates the forward simulation differential coefficient; A finite difference module calculates the source wave field and the receiver wave field of each shot at each time by a finite difference method according to the anisotropic wave equation, based on the grid distribution information and the forward simulation differential coefficients; The cross-correlation module performs cross-correlation calculations on the source wave field and the detection point wave field at the same time of the shot, and sums the cross-correlation results to obtain the single-shot migration result; The superposition module superimposes all the single-shot offset results and outputs the final offset result; Wherein, according to the velocity model and the offset parameter, variable space step size gridding is adopted to output three-dimensional space grid distribution information, including: The minimum value of the grid space step is determined according to the velocity model and the migration parameters using the dispersion relation of the anisotropic wave equation; wherein the dispersion relation of the anisotropic wave equation is: in, v (θ) is the phase velocity of longitudinal waves in anisotropic media, According to formula (3), the minimum value of the grid space step length, v p is the longitudinal wave velocity, v s is the shear wave velocity, ε and δ are anisotropy parameters, and θ is the inclination angle; The maximum value of the grid space step size is determined according to the velocity model and the migration parameters using the stability condition of the anisotropic wave equation; wherein the stability condition of the anisotropic wave equation is: Among them, dx, dy and dz are the grid space steps, v(x, y, z) is the velocity of seismic wave propagation, Δt is the time step, and a m is the difference coefficient, N is the difference order, j represents the number of grids from the difference center point, and the maximum value of the grid space step length is obtained according to formula (4); Determining a range of grid space steps to which the seismic wave velocity model is adapted according to the minimum value of the grid space step and the maximum value of the grid space step; Determining, according to the grid space step range, a mean value of the grid space step adapted to each speed value; Determine the grid space step size corresponding to each speed value according to the average value of the grid space step size adapted by each speed value; According to the grid space step size corresponding to each velocity value, it is substituted into the Halton grid generation method to output the three-dimensional space grid distribution information.

8. An electronic device, characterized in that: The electronic device comprises: A memory storing executable instructions; A processor, wherein the processor runs the executable instructions in the memory to implement the anisotropic seismic migration imaging method according to any one of claims 1 to 6.

9. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the anisotropic seismic migration imaging method according to any one of claims 1 to 6 is implemented.

Images