Marine controlled source electromagnetic inversion method based on image structure tensor guidance

Abstract

The present invention relates to a marine controlled source electromagnetic inversion method guided by image structure tensor, comprising the following steps: A) obtaining a seismic image, establishing an inversion objective function based on the roughness of a model vector m, the degree of deviation from a preferred model vector m′ and the difference between measured data and the model; B) obtaining an inversion structure tensor matrix based on the structure tensor and combining a Gaussian difference operator; C) decomposing the matrix into two eigenvalues ​​and two orthogonal eigenvectors u and v, wherein u and v are the directions with the maximum and minimum intensity value change rates in the seismic image, respectively; D) determining the geological boundary in the seismic image by the eigenvalue, calculating the absolute value of the cosine value of the angle formed by the eigenvector v of a grid adjacent to the geological boundary and the eigenvector v of the corresponding central grid; E) assigning an adjustment coefficient to the absolute value and using it as the coefficient of the roughness of the model vector m. The recognition method of the present invention can improve the model resolution of traditional marine controlled source electromagnetic inversion.

Description

Technical Field

[0001] The invention relates to a geological exploration data processing method, and in particular to an ocean controllable source electromagnetic inversion method based on image structure tensor guidance. Background Art

[0002] Seismic exploration is a widely used method for oil and gas exploration. Different strata have different pore structures and are filled with different pore fluids, thus corresponding to different rock physical parameters. The differences in these parameters are ultimately reflected in the seismic migration images and velocity inversion results, from which we can roughly determine the location of the oil and gas reservoir.

[0003] Since the resistivity of reservoir rocks containing oil and gas depends largely on their porosity and the resistivity and content of the fluid in the pore space, as the content of hydrocarbon fluid filling the rock pores increases, the resistivity of the rock will also increase. Due to its high sensitivity to underground high-resistance thin layers, marine controlled source electromagnetic method has become one of the most effective methods for marine oil and gas resource exploration in the past two decades, and inversion of data obtained by marine controlled source electromagnetic method is an important means to obtain the boundaries of oil and gas fields.

[0004] The inversion methods in the prior art are generally forced to solve a very large set of ill-conditioned linear equations, and the results are non-unique and unstable. In order to make the results stable, the prior art often adds a regularization term to the inversion objective function to obtain a smooth resistivity model. However, this method blurs the edges of the true abnormal area, resulting in low resolution of the imaging results, large errors, and low practicality.

[0005] In view of this, it is necessary to provide an ocean controlled source electromagnetic inversion method based on image structure tensor guidance. Summary of the invention

[0006] The present invention aims to provide an ocean controlled source electromagnetic inversion method based on image structure tensor guidance, which can improve the resolution of boundary determination and adjust the sharpness of the inversion boundary according to actual needs.

[0007] To achieve the above purpose of the invention, the marine controlled source electromagnetic inversion method based on image structure tensor guidance provided by the present invention comprises the following steps:

[0008] A) acquiring a seismic image and establishing an inversion objective function based on the roughness of a model vector m, the deviation between the model vector m and a preferred model vector mʹ in a priori preferred model, and the difference between measured data and the model F(m);

[0009] B) Based on the structure tensor, the inverted structure tensor matrix is ​​obtained by combining the Gaussian difference operator;

[0010] C) decomposing the inversion structure tensor matrix into two eigenvalues ​​and two orthogonal eigenvectors, eigenvector u and eigenvector v, wherein the eigenvector u is the direction with the largest intensity value change rate in the seismic image, and the eigenvector v is the direction with the smallest intensity value change rate in the seismic image;

[0011] D) determining a geological boundary in the seismic image according to the characteristic value, and calculating the absolute value of the cosine value of an angle formed by a characteristic vector v of a grid adjacent to the geological boundary and a characteristic vector v of a central grid on the boundary corresponding to the grid;

[0012] E) taking the absolute value of the cosine value as the coefficient of the roughness of the model vector m, and assigning an adjustment coefficient to the absolute value of the cosine value to adjust the weight of the roughness in the objective function.

[0013] Specifically, the seismic image in step A) is obtained by electromagnetic exploration of underground media.

[0014] More specifically, the objective function in step A) is:

[0015]

[0016] in, is the roughness of the model vector m, is the degree of deviation between the model vector m and the preference model vector mʹ in the prior preference model,

[0017] Represents the mismatch value between the measured data d and the model F(m) and the target X 2 The difference between , W is the diagonal matrix of the data covariance weighting function, μ -1 is the Lagrange multiplier that balances roughness and poor fit.

[0018] Preferably, the roughness It can be discretized using an unstructured grid as:

[0019]

[0020] Where np is the total number of grids, nb is the number of adjacent grids of each grid, which changes with the change of the center grid number i; is the area of ​​the jth adjacent grid of the ith grid, is the area of ​​all adjacent grids of the ith grid the sum of is the distance between the center grid numbered i and the adjacent j-th grid, is the relative difference between the center grid numbered i and the adjacent j-th grid.

[0021] Specifically, the inverted structure tensor matrix in step B) is for:

[0022]

[0023] Among them, * is the convolution operator, is the Gaussian difference operator, is the gradient of the intensity value of the pixel in the seismic image in the horizontal direction; is the gradient of the intensity value of the pixel in the horizontal and vertical directions, = , , , , x represents the position number of the pixel in the horizontal direction, and y represents the position number of the pixel in the vertical direction.

[0024] More specifically, the Gaussian difference operator for:

[0025]

[0026] Where e is the natural logarithm, is the scaling factor.

[0027] Specifically, the inverted structure tensor matrix The decomposition formula is as follows:

[0028]

[0029] in, ,

[0030] ,

[0031] is the transposed matrix of the eigenvector u, is the transposed matrix of the eigenvector v.

[0032] More specifically, the The value and When the values ​​of are all tending to 0, the pixel is in the flat area of ​​the seismic image; The value is much greater than 0 and the When the value of tends to 0, the pixel is in the boundary area of ​​the seismic image; The value is oriented to 0 and the When the value of is much greater than 0, the pixel is in the corner area of ​​the seismic image.

[0033] More specifically, the coefficient of the roughness is expressed as:

[0034]

[0035] in, is the feature vector v of the i-th grid; is the feature vector v of the j-th adjacent grid of the i-th grid.

[0036] More specifically, the adjustment coefficient K is:

[0037]

[0038] Among them, w is a constant whose size can be adjusted.

[0039] The marine controlled source electromagnetic inversion method based on image structure tensor guidance provided by the present invention first introduces an unstructured grid in the roughness calculation of the inversion. The unstructured grid has a stronger adaptability to the terrain and is more suitable for inversion under real terrain; then, a structure tensor matrix is ​​introduced in the roughness calculation of the inversion, so that the geological boundary in the seismic image can be judged according to the eigenvalue of the structure tensor matrix, and further the difference in the characteristic direction of the adjacent grid and the central grid corresponding to the adjacent grid is judged by the angle formed by the characteristic vectors of the adjacent grids adjacent to the boundary and the central grid, and the absolute value of the cosine value of the angle is used as the coefficient of the roughness, so that the boundary obtained after the inversion can be smoother and have higher accuracy and resolution; in addition, an artificially adjustable adjustment coefficient is introduced in the roughness calculation to strengthen the influence of the absolute value of the cosine value of the angle on the roughness, and the sharpness of the inversion boundary can be adjusted by changing the adjustment coefficient according to actual needs.

[0040] Other features and advantages of the embodiments of the present invention will be described in detail in the following specific embodiment section. BRIEF DESCRIPTION OF THE DRAWINGS

[0041] Figure 1 It is a flow chart of the marine controlled source electromagnetic inversion method based on image structure tensor guidance of the present invention;

[0042] Figure 2 It is the original seafloor resistivity model diagram of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention;

[0043] Figure 3 It is the abnormal body boundary screened by eigenvalue in the marine controlled source electromagnetic inversion method guided by image structure tensor of the present invention;

[0044] Figure 4 This is the inversion result diagram of the traditional method;

[0045] Figure 5 It is a graph showing the relationship between the target fitting error and the number of iterations during the inversion process of the traditional method;

[0046] Figure 6 It is a graph showing the relationship between roughness and the number of iterations during the inversion process of the traditional method;

[0047] Figure 7 This is the inversion result diagram of the marine controlled source electromagnetic inversion method based on image structure tensor guidance when w=0.1;

[0048] Figure 8 It is a relationship diagram between the target fitting error and the number of iterations in the inversion process of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.1;

[0049] Fig. 9 It is a relationship diagram between roughness and the number of iterations in the inversion process of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.1;

[0050] Fig.10 This is the inversion result diagram of the marine controlled source electromagnetic inversion method based on image structure tensor guidance when w=0.05;

[0051] Fig.11 It is a relationship diagram between the target fitting error and the number of iterations in the inversion process of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.05;

[0052] Fig.12 It is a relationship diagram between roughness and the number of iterations in the inversion process of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.05;

[0053] Fig.13 This is the inversion result diagram of the marine controlled source electromagnetic inversion method based on image structure tensor guidance when w=0.01;

[0054] Fig.14 It is a relationship diagram between the target fitting error and the number of iterations in the inversion process of the marine controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.01;

[0055] Fig.15 It is a relationship diagram between roughness and the number of iterations in the inversion process of the ocean controlled source electromagnetic inversion method guided by the image structure tensor of the present invention when w=0.01. DETAILED DESCRIPTION

[0056] The specific implementation of the embodiment of the present invention is described in detail below in conjunction with the accompanying drawings. It should be understood that the specific implementation described here is only used to illustrate and explain the embodiment of the present invention, and is not used to limit the embodiment of the present invention.

[0057] like Figure 1 As shown, in one embodiment of the marine controlled source electromagnetic inversion method based on image structure tensor guidance of the present invention, the marine controlled source electromagnetic inversion method based on image structure tensor guidance includes the following steps:

[0058] Step A): Obtain a seismic image, and establish an inversion objective function based on the roughness of the model vector m, the deviation between the model vector m and the preferred model vector mʹ in the prior preference model, and the difference between the measured data and the model F(m). Taking the detection of oil and gas fields under the seabed as an example, the marine controlled source electromagnetic method can be used, that is, electromagnetic exploration of the medium under the seabed is performed to obtain seismic image data of the seabed, and then the resistivity model of the detection area is obtained by inverting the seismic image data, and the geological boundary, that is, the boundary of the oil and gas field, is determined by the obtained resistivity model. The specific inversion formula can use the inversion objective function of the two-dimensional marine controlled source electromagnetic method:

[0059]

[0060] in, is the roughness of the model vector m, is the degree of deviation between the model vector m and the preference model vector mʹ in the prior preference model, represents the difference between the measured data d and the model F(m), W is the data covariance weighting function diagonal matrix, X 2 is the target mismatch value, μ -1 is the Lagrange multiplier that balances roughness and poor fit.

[0061] In order to make the objective function better suitable for the inversion under the real terrain, the roughness The calculation is extended to the inversion roughness of unstructured grids, i.e. triangular grids. Specifically, the roughness After discretization by unstructured grid, it can be expressed as:

[0062]

[0063] Where np is the total number of grids, nb is the number of adjacent grids of each grid, which changes with the change of the center grid number i; is the area of ​​the jth adjacent grid of the ith grid, is the area of ​​all adjacent grids of the ith grid sum; is the distance between the center grid numbered i and the adjacent j-th grid, is the relative difference between the center grid numbered i and the adjacent j-th grid.

[0064] Using unstructured grids to After discretization, the model roughness can be calculated based on the grid area and the distance between grids as the metric of the size of the grid's corresponding weight coefficient. Adjacent grids with large areas and close distances have higher weights, while adjacent grids with small areas and far distances have lower weights. However, the inversion of the two-dimensional marine controlled source electromagnetic method aims to seek a smooth resistivity model. The inversion result may not be very accurate, but in terms of roughness. Under the constraint of , the inversion is stable, and the inversion result also reduces the possibility of false anomalies. However, if the weight information cannot reflect the characteristic boundary information of the target body, the inversion will only obtain a smooth resistivity model, and the resolution of the inversion result is low, which cannot meet the accuracy requirements of exploration. Therefore, it is necessary to improve the resolution of the resistivity model. Therefore, the structural tensor information of the seismic image with high resolution can be further introduced into the inversion to redefine the roughness term in the objective function.

[0065] That is, execute step B): Based on the structure tensor, combined with the Gaussian difference operator The inverted structure tensor matrix is ​​obtained :

[0066]

[0067] Among them, * is the convolution operator, is the Gaussian difference operator, is the horizontal gradient of the intensity value of the pixel in the seismic image; is the gradient of the pixel intensity value in the horizontal and vertical directions, = , , , , x represents the position number of the pixel in the horizontal direction, y represents the position number of the pixel in the vertical direction, and the Gaussian difference operator Specifically:

[0068]

[0069] Where e is the natural logarithm, is the scaling factor.

[0070] Then perform step C): decompose the inverted structure tensor matrix into two eigenvalues ​​and two orthogonal eigenvectors: eigenvector u and eigenvector v, where , then the eigenvector u is the direction with the largest intensity value change rate in the seismic image, and the eigenvector v is the direction with the smallest intensity value change rate in the seismic image. The decomposed inversion structure tensor matrix It is expressed as:

[0071]

[0072] in, ,

[0073] ,

[0074] is the transposed matrix of the eigenvector u,

[0075] is the transposed matrix of the eigenvector v.

[0076] Then, step D is executed: the geological boundary in the seismic image is determined according to the decomposed eigenvalues. The specific determination method is: when the inversion structure tensor matrix of a certain pixel of and When the values ​​of tend to 0, the pixel is in the flat area of ​​the seismic image; when the inversion structure tensor matrix of a pixel of The value is much greater than 0 and When the value of tends to 0, the pixel is in the boundary area of ​​the seismic image; when the inversion structure tensor matrix of a pixel of The value tends to 0 and When the value of is much greater than 0, the pixel is in the corner area of ​​the seismic image. Figure 2 Taking the original seafloor resistivity model diagram shown in the figure as an example, the anomaly boundary (geological boundary) obtained by the above method is as follows Figure 3 As shown by Figure 3 It can be seen that the geological boundary obtained by the above method is a relatively wide range, that is, the resolution of the inversion result is low. In order to highlight the more detailed geological boundary structure and improve the resolution of the inversion result, it is necessary to apply a smaller weight value on the anomaly boundary to relax the roughness. The calculation highlights the resistivity difference near the boundary, and at the same time increases the weight value in the area where the linear feature exists, limits the parameter change in this area during the iteration process, and makes the model parameters approach the same, so that the smoothing direction of the inversion result is guided by higher-resolution image information, which can be more consistent with the linear feature, and the resolution of the inversion result will also be improved accordingly.

[0077] Specifically, the absolute value of the cosine value of the angle formed by the characteristic vector v of the grid adjacent to the geological boundary and the characteristic vector v of the central grid corresponding to the grid on the geological boundary can be selected as the roughness The weight coefficient is expressed as:

[0078]

[0079] in, is the characteristic vector v of the i-th grid; is the feature vector v of the jth adjacent grid of the ith grid. The size of can reflect the difference between the adjacent grid and the central grid in the characteristic direction. The adjacent grid with almost the same characteristic direction as the central grid will have a larger value, otherwise it has a smaller Value, use Roughness Correction can achieve relaxation at the boundary and smoothing of linear features. is always less than or equal to 1. In order to increase the difference between the two grids and give the approximately isotropic grid a larger weight, step E can be further performed): the absolute value of the cosine value Assign adjustment coefficient K:

[0080]

[0081] Among them, w is a constant with adjustable size. A smaller w value produces a larger penalty weight difference and a sharper inversion boundary, while a larger w value makes the inversion result closer to the result of smooth inversion of the two-dimensional marine controlled source electromagnetic method.

[0082] Adjusted roughness The weight coefficient It is expressed as:

[0083]

[0084] Therefore, the adjusted roughness term It can be specifically expressed as:

[0085]

[0086] In order to verify the influence of the w value on the inversion results, we first establish Figure 2 The original seabed resistivity model is shown in FIG. 1 . Then, the original seabed resistivity model is first inverted using the traditional inversion method to obtain the following Figure 4 The inversion result diagram shown is Figure 4It can be seen that traditional methods can only roughly determine the general direction of abnormal geological areas (such as oil and gas fields, etc.), with low resolution and accuracy, and cannot accurately determine their boundaries.

[0087] Then, the inversion objective function provided by the present invention is used, and w is assigned a value of 0.1 to invert the original seabed resistivity model, and the following can be obtained: Figure 7 The inversion result diagram shown is Figure 7 It can be seen that when w is assigned to 0.1, the overall outer contour of the geological anomaly area can be judged more clearly and accurately, but the accuracy still cannot meet the needs of geological exploration; then w is assigned to 0.05 to invert the original seabed resistivity model, and the following can be obtained: Fig.10 The inversion result diagram shown is Fig.10 It can be seen that when w is assigned a value of 0.05, the overall outer contour of the geological anomaly area can be judged more clearly and accurately. Its accuracy can basically meet the needs of geological exploration and can reflect the boundaries between geological anomaly areas such as areas with different natural gas contents in gas fields. Finally, w is assigned a value of 0.01 to invert the original seabed resistivity model, and the following can be obtained: Fig.13 The inversion result diagram shown is Fig.13 It can be seen that when w is assigned to 0.01, the overall outer contour of the geological anomaly area can be basically consistent with the geological boundary in the original seabed resistivity model established, and can accurately and clearly reflect the geological anomaly area, such as the boundary between areas with different natural gas contents in the gas field. Therefore, the larger the w value, the smoother the image-guided inversion result, and the closer it is to the traditional smooth inversion result. When a smaller w value is selected, the model boundary in the inversion result is sharper and clearer; in addition, if Figure 5 , Figure 8 , Fig.11 and Fig.14 As shown in FIG. 1 , in terms of the number of inversion iterations, the inversion method provided by the present invention requires more iterations to reach the target RMS Misfit compared with the traditional inversion method. The larger the value of w, the fewer iterations are required, and vice versa. Figure 6 , Fig. 9 , Fig.12 and Fig.15 As shown in the figure, the model roughness of the inversion method provided by the present invention is greater than that of the traditional inversion, which means that after being processed by the inversion method provided by the present invention, the parameter differences between the grids change from smooth to obvious, so the final image resolution is improved.

[0088] Based on the above test results, it can be seen that the marine controlled source electromagnetic inversion method based on image structure tensor guidance of the present invention can add higher-precision underground structure information detected by earthquakes or other methods as weight constraints in the inversion to the inversion objective function. Compared with the joint inversion, the marine controlled source electromagnetic inversion method based on image structure tensor guidance of the present invention does not introduce other data into the objective function, so the amount of calculation is smaller. At the same time, since the form of the inversion objective function of the two-dimensional marine controlled source electromagnetic method is not changed, the implementation is also more convenient. Under the same number of iterations, the resolution of the image-guided inversion result is significantly better than that of the traditional inversion method. In addition, the degree of inversion smoothness can be adjusted artificially by controlling the adjustment coefficient K (specifically the w value).

[0089] The optional implementation modes of the embodiments of the present invention are described in detail above in conjunction with the accompanying drawings. However, the embodiments of the present invention are not limited to the specific details in the above implementation modes. Within the technical concept of the embodiments of the present invention, various simple modifications can be made to the technical scheme of the embodiments of the present invention, and these simple modifications all belong to the protection scope of the embodiments of the present invention.

[0090] It should also be noted that the various specific technical features described in the above specific embodiments can be combined in any suitable manner without contradiction. To avoid unnecessary repetition, the embodiments of the present invention will not further describe various possible combinations.

[0091] In addition, various implementation modes of the embodiments of the present invention may be arbitrarily combined, and as long as they do not violate the concept of the embodiments of the present invention, they should also be regarded as the contents disclosed by the embodiments of the present invention.

Claims

1. A method for marine controlled source electromagnetic inversion based on image structure tensor guidance, characterized in that: The steps include: A) acquiring a seismic image and establishing an inversion objective function based on the roughness of a model vector m, the deviation between the model vector m and a preferred model vector mʹ in a priori preferred model, and the difference between measured data and the model F(m); B) Based on the structure tensor, the inverted structure tensor matrix is ​​obtained by combining the Gaussian difference operator; C) decomposing the inversion structure tensor matrix into two eigenvalues ​​and two orthogonal eigenvectors, eigenvector u and eigenvector v, wherein the eigenvector u is the direction with the largest intensity value change rate in the seismic image, and the eigenvector v is the direction with the smallest intensity value change rate in the seismic image; D) determining a geological boundary in the seismic image according to the characteristic value, and calculating the absolute value of the cosine value of an angle formed by a characteristic vector v of a grid adjacent to the geological boundary and a characteristic vector v of a central grid on the boundary corresponding to the grid; E) taking the absolute value of the cosine value as a coefficient of the roughness of the model vector m, and assigning an adjustment coefficient to the absolute value of the cosine value to adjust the weight of the roughness in the objective function; The objective function in step A) is: in, is the roughness of the model vector m, is the model vector m and the preference model vector in the prior preference model The deviation Represents the mismatch value between the measured data d and the model F(m) and the target X 2 The difference between them, W is the diagonal matrix of the data covariance weighting function, μ -1 Lagrange multiplier to balance roughness and poor fit; The roughness It can be discretized using an unstructured grid as: Where np is the total number of grids, nb is the number of adjacent grids of each grid, which changes with the change of the center grid number i; is the area of ​​the jth adjacent grid of the ith grid, is the area of ​​all adjacent grids of the ith grid the sum of is the distance between the center grid numbered i and the adjacent j-th grid, is the relative difference between the center grid numbered i and the adjacent j-th grid.

2. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 1 is characterized in that: The seismic image in step A) is obtained by electromagnetic exploration of underground media.

3. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 1, characterized in that: Step B) The inverted structure tensor matrix for: Among them, * is the convolution operator, is the Gaussian difference operator, is the gradient of the intensity value of the pixel in the seismic image in the horizontal direction; is the gradient of the intensity value of the pixel in the horizontal and vertical directions, = , , , x represents the position number of the pixel in the horizontal direction, and y represents the position number of the pixel in the vertical direction.

4. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 3 is characterized in that: The Gaussian difference operator for: Where e is the natural logarithm, is the scaling factor.

5. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 4 is characterized in that: The inverted structure tensor matrix The decomposition formula is as follows: in, , ， is the transposed matrix of the eigenvector u, is the transposed matrix of the eigenvector v.

6. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 5, characterized in that: Said The value and When the values ​​of are all tending to 0, the pixel is in the flat area of ​​the seismic image; The value is much greater than 0 and the When the value of tends to 0, the pixel is in the boundary area of ​​the seismic image; The value tends to 0 and the When the value of is much greater than 0, the pixel is in the corner area of ​​the seismic image.

7. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 6 is characterized in that: The coefficient of the roughness is expressed as: in, is the feature vector v of the i-th grid; is the feature vector v of the j-th adjacent grid of the i-th grid.

8. The method for marine controlled source electromagnetic inversion based on image structure tensor guidance according to claim 7 is characterized in that: The adjustment coefficient K is: Among them, w is a constant whose size can be adjusted.

Images