An electromagnetic imaging method, system, medium, device and terminal

Abstract

The application belongs to the technical field of electromagnetic imaging, and discloses an electromagnetic imaging method, system, medium, equipment and terminal, wherein a cross-variation function fitting term of resistivity-velocity is added in a traditional magnetotelluric inversion objective function to represent the relationship between velocity and resistivity; a fuzzy clustering idea is introduced in an NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, clustering input is resistivity and velocity, and the inversion is constrained by structural information by exerting constraint on membership in the clustering process, so that the electromagnetic-seismic combined modeling inversion based on the rock physical relationship and structural similarity double constraints is realized. The application takes the rock physical model as a link, adds the seismic data constraint, suppresses the magnetotelluric inversion multi-solution, and improves the imaging effect of deep exploration. The electromagnetic imaging can reduce the multi-solution of inversion, improve the magnetotelluric imaging precision, exert the complementary characteristics between resistivity and velocity, and highlight the underground physical property distribution and structural characteristics.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of electromagnetic imaging, and particularly relates to an electromagnetic imaging method, system, medium, device and terminal. BACKGROUND

[0002] At present, the magnetotelluric method has low cost, large detection depth, can reach tens of kilometers or even hundreds of kilometers, is not affected by high-resistance layers, is sensitive to low-resistance layers, and can divide structural zones, analyze fracture systems, and identify special rock bodies, and can play an important role in deep exploration. However, with the increase of exploration target depth, the limited deep data information and the non-uniqueness of inversion and interpretation seriously challenge the reliability of exploration. How to reduce the problem of inversion multiplicity has always been one of the key problems to be solved in geophysical inversion.

[0003] Joint inversion is an effective way to reduce the problem of inversion multiplicity and improve imaging quality, and has become a hot issue in the field of geophysical inversion. The premise of reliable joint inversion of electricity and earthquake is reliable rock physics relationship and single inversion result. In the traditional rock physics relationship coupling method, the relationship function between the models is defined first, and the function is applied to the entire region. Since this relationship may change with the geological environment, this method is not reliable and is not suitable for complex geological structure areas with multiple rock physics relationships. The traditional magnetotelluric inversion method has insufficient imaging accuracy of electrical interfaces and is not easy to perform geological interpretation.

[0004] Through the above analysis, the problems and defects of the prior art are that the magnetotelluric inversion method has insufficient imaging accuracy of electrical interfaces and is not easy to perform geological interpretation. The traditional joint inversion rock physics relationship coupling method uses a fixed rock physics relationship, which has low reliability and is not suitable for complex geological structure areas with multiple rock physics relationships. SUMMARY

[0005] In view of the problems existing in the prior art, the present application provides an electromagnetic imaging method, system, medium, device and terminal, and particularly relates to an electromagnetic imaging method, system, medium, device and terminal based on rock physics modeling.

[0006] The present application is implemented as follows: an electromagnetic imaging method comprises:

[0007] A cross-variation function fitting term of resistivity-velocity is added to the traditional magnetotelluric inversion objective function to represent the relationship between velocity and resistivity; a fuzzy clustering idea is introduced into the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, the clustering input is resistivity and velocity, and the inversion is constrained by structural information by applying constraints to the membership in the clustering process, so as to realize joint modeling inversion of electricity and earthquake based on rock physics relationship and structural similarity double constraints.

[0008] Further, the electromagnetic imaging method comprises the following steps:

[0009] Step one, constructing a magnetotelluric fuzzy constraint inversion to obtain a resistivity model;

[0010] Step two, using the obtained resistivity model to establish an initial velocity model with a shallow seismic framework;

[0011] Step three, constructing a seismic tomography inversion to obtain a velocity model;

[0012] Step four, using the obtained velocity model to establish an initial resistivity model with a resistivity-velocity model;

[0013] Step five, based on the rock physical relationship and structural similarity double constraint electromagnetic seismic joint modeling inversion to obtain electromagnetic imaging results.

[0014] Further, the objective function of the double constraint electromagnetic seismic joint modeling inversion in the electromagnetic imaging method comprises:

[0015] Integrate the cross-variogram in multivariate geostatistics into the traditional magnetotelluric inversion objective function equation to represent the relationship between velocity and resistivity, and the cross-variogram fitting term is as follows:

[0016]

[0017] Wherein, γ ρv (m) represents the cross-variogram estimated from the inversion model, γ ρv obs is the prior cross-variogram obtained from prior petrophysical information; the matrix W γ represents the uncertainty in the prior petrophysical relationship, W γ is a diagonal matrix, and the diagonal elements are the inverse of the standard deviation of the cross-variogram mismatch.

[0018] The objective function of the double constraint electromagnetic seismic joint modeling inversion is:

[0019]

[0020] Wherein, Φ(m) is the total objective function, Φ d (m) represents the magnetotelluric data fitting term, Φ m (m) represents the model matching term; Φ γ (m) represents the cross-variogram fitting term, which is the resistivity-velocity relationship fitting term; m is the model parameter, d is the measured MT data, C d is the covariance of data error, f(m) is the forward MT data, m0 is the initial model; λ is the Lagrange operator, which is the regularization parameter; Cm is the model covariance, is the smoothing operator; The parameter a controls the closeness of the cross-variogram obtained by inversion to the prior cross-variogram.

[0021] By adding the cross-variogram fitting term, the magnetotelluric inversion result is fitted to the observed magnetotelluric data and simultaneously satisfies the prior resistivity-velocity relationship.

[0022] Further, the improved nonlinear conjugate gradient method in the electromagnetic imaging method comprises:

[0023] Compared with the traditional magnetotelluric inversion algorithm process, the algorithm process of the double-constrained joint modeling inversion of electric and seismic includes the following three improvements:

[0024] (1) Since the cross-variogram fitting term representing the resistivity-velocity relationship is added in the objective function, the gradient of the cross-variogram fitting term is added when calculating the gradient of the objective function;

[0025] (2) In the NLCG algorithm, the fuzzy clustering idea is introduced in the inversion process to obtain an improved nonlinear conjugate gradient algorithm, and the clustering input is resistivity and velocity, so that the model update is affected by resistivity and velocity at the same time;

[0026] (3) The membership degree is spatially constrained in the clustering process, so that the model update is affected by the resistivity-velocity relationship and the spatial structure constraint.

[0027] Further, the cross-variogram forward formula and Jacobian matrix derivation in the electromagnetic imaging method comprise:

[0028] The cross-variogram is rewritten as follows:

[0029]

[0030] Where i is the group number of different data groups divided according to different lags h; k and j represent the serial numbers of samples in each data group, and j≠k.

[0031] The cross-variogram differential form with respect to the resistivity model of each serial number j is:

[0032]

[0033] If the resistivity model of serial number j is not in the ith data pair, then And the ith data pair has no effect on The resistivity model differential form is rewritten as:

[0034]

[0035] where N(h) j is the log of the data pair containing resistivity model of index j.

[0036] The differential form of the resistivity model is then rewritten as:

[0037]

[0038] The first term on the right is expanded and combined to give:

[0039]

[0040] The partial derivative of the cross-variogram with respect to the resistivity model of index j is:

[0041]

[0042] The value of the hth row and jth column of the Jacobian matrix of the cross-variogram is found, denoted as J γ:h,j .

[0043] After the forward formula of the cross-variogram and the Jacobian matrix of the cross-variogram with respect to the resistivity model are obtained, the inversion is performed using an iterative solution.

[0044] Further, the cluster center and membership degree calculation in the electromagnetic imaging method comprises:

[0045] The output of the FCM clustering is the cluster center value and the membership degree of the petrophysical information. The cluster center value is the representative physical value of the model, and the membership degree contains the spatial information of the model. When the clustering input is resistivity and velocity, the petrophysical information and structural information of the magnetotelluric and seismic models are exchanged in the clustering process.

[0046] When the cluster center value and the membership degree are known, the model is constructed, which is different from constraining the entire model parameter, and the cluster center value and the membership degree variable are independently constrained respectively.

[0047] (1) Cluster center

[0048] The seismic model S and the magnetotelluric model M are used together to form the fuzzy clustering input data of Z = [M, S], and the seismic model is included in the magnetotelluric inversion constrained by fuzzy clustering. The membership degree and the cluster center value of the magnetotelluric model are affected by the seismic model, and the seismic information indirectly affects the magnetotelluric model through the FCM clustering process.

[0049] (2) Boundary constraint membership degree

[0050] The membership degree is a measure of the degree of belonging of an element to a certain class, and in geological classification or geophysical inversion, it is regarded as the possibility of a certain model unit belonging to a certain lithology or lithostratigraphic unit. The membership degree matrix of the elements in the inversion model contains spatial information, and the prior geometric information is constrained by using the membership degree in the inversion process.

[0051] The spatial information b is simply regarded as the classification feature data in the fuzzy clustering process. The function and the model parameters are combined to form the input data [m, b] for fuzzy clustering, so the distance between the jth model parameter and the kth rock cluster center is represented as:

[0052]

[0053] In the formula, when the model parameters represent resistivity and velocity, c k is the cluster center of the kth rock, w k is a weight coefficient; the first term defines the distance between the model parameter value and the rock cluster center value, and the second term contains boundary information, which is defined as follows:

[0054]

[0055] Another object of the present application is to provide an electromagnetic imaging system applying the electromagnetic imaging method, which comprises:

[0056] A resistivity model construction module is configured to construct a resistivity model obtained by magnetotelluric fuzzy constraint inversion.

[0057] An initial velocity model construction module is configured to establish an initial velocity model by using the obtained resistivity model and a shallow seismic framework.

[0058] A velocity model construction module is configured to construct a velocity model obtained by seismic tomography inversion.

[0059] An initial resistivity model construction module is configured to establish an initial resistivity model by using the obtained velocity model and a resistivity-velocity model.

[0060] An electromagnetic imaging module is configured to obtain electromagnetic imaging results by electroseismic joint modeling inversion based on rock physical relationship and structural similarity double constraints.

[0061] Another object of the present application is to provide a computer device, which comprises a memory and a processor, wherein the memory stores a computer program, and the computer program is executed by the processor to enable the processor to perform the following steps:

[0062] The cross-variogram fitting term of resistivity-velocity is added in the traditional magnetotelluric inversion objective function to represent the relationship between velocity and resistivity; the fuzzy clustering thought is introduced in the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, the clustering input is resistivity and velocity, the clustering process is constrained by the membership to constrain the inversion by the structure information, and the electrical and seismic joint modeling inversion based on the rock physical relationship and structure similarity double constraints is realized.

[0063] Another object of the present application is to provide a computer readable storage medium storing a computer program, which, when executed by a processor, causes the processor to perform the following steps:

[0064] The cross-variogram fitting term of resistivity-velocity is added in the traditional magnetotelluric inversion objective function to represent the relationship between velocity and resistivity; the fuzzy clustering thought is introduced in the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, the clustering input is resistivity and velocity, the clustering process is constrained by the membership to constrain the inversion by the structure information, and the electrical and seismic joint modeling inversion based on the rock physical relationship and structure similarity double constraints is realized.

[0065] Another object of the present application is to provide an information data processing terminal for realizing the electromagnetic imaging system.

[0066] In combination with the above technical solutions and the technical problems solved, the advantages and positive effects of the technical solutions to be protected by the present application are analyzed from the following aspects:

[0067] First, in view of the technical problems existing in the prior art and the difficulty in solving the problems, the technical solutions to be protected by the present application and the results and data in the research and development process are closely combined to analyze in detail and profoundly how the technical solutions solve the technical problems and bring some creative technical effects after solving the problems. The specific description is as follows:

[0068] The electrical and seismic joint modeling inversion based on the rock physical relationship and structure similarity double constraints provided by the present application mainly makes the following improvements compared with the traditional magnetotelluric inversion: 1. The cross-variogram fitting term of resistivity-velocity is added in the traditional magnetotelluric inversion objective function to represent the relationship between velocity and resistivity; 2. The fuzzy clustering thought is introduced in the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, the clustering input is resistivity and velocity, and the clustering process is constrained by the membership to constrain the inversion by the structure information, and the electrical and seismic joint modeling inversion based on the rock physical relationship and structure similarity double constraints is realized.

[0069] The application provides an electromagnetic imaging method based on rock physical modeling, which takes rock physical model as a link, adds the constraint of seismic data, suppresses the multi-solution of magnetotelluric inversion, and improves the imaging effect of deep exploration.

[0070] The electromagnetic imaging method based on rock physical modeling can reduce the multi-solution of inversion, improve the magnetotelluric imaging precision, exert the complementary characteristics between resistivity and velocity, and highlight the underground physical property distribution and structural characteristics.

[0071] Secondly, the technical effect and advantages of the technical solution to be protected by the application are described as follows from the perspective of the product as a whole.

[0072] Theoretical model test shows that the electromagnetic imaging method based on rock physical modeling can reduce the multi-solution of inversion, improve the magnetotelluric imaging precision, exert the complementary characteristics between resistivity and velocity, and highlight the underground physical property distribution and structural characteristics.

[0073] Thirdly, the creativity of the claims of the application is also embodied in the following important aspects.

[0074] (1) The expected income and commercial value of the technical solution of the application after transformation are as follows:

[0075] The technical solution of the application improves the imaging precision of resistivity, exerts the complementary characteristics between resistivity and velocity, highlights the underground physical property distribution and structural characteristics, and is preliminarily applied in a certain area of the piedmont zone of the Tarim Basin. The electromagnetic imaging result based on rock physical modeling has better consistency with the known prior information, and better describes the distribution of the conglomerate layer in the piedmont zone. It is proved that the technical solution of the application after transformation can lay a foundation for subsequent structural analysis and correct interpretation of oil and gas, provide necessary geophysical basis, and support the deep-sea and deep-land oil and gas strategy of China.

[0076] (2) The technical solution of the application fills the technical blank in the industry at home and abroad:

[0077] The technical solution of the application fills the technical blank of constructing an underground resistivity-velocity model which meets the prior rock physical relationship and has similar structure by using the fuzzy clustering and geostatistical methods, and improving the electromagnetic imaging precision.

[0078] (3) The technical scheme of the present application solves the technical problems that people have been eager to solve but have failed to succeed:

[0079] The technical scheme of the present application solves the technical problems of the conventional magnetotelluric inversion method, which has low imaging accuracy of the electrical interface, is not easy to be geologically interpreted, and the joint inversion based on the conventional rock physical relationship coupling method, which uses fixed rock physical relationship, has low reliability and is not suitable for complex geological structure regions with multiple rock physical relationships.

[0080] (4) The technical scheme of the present application overcomes the technical prejudice:

[0081] The technical scheme of the present application overcomes the technical prejudice that electromagnetic imaging is generally considered to have low accuracy and is not easy to be geologically interpreted. BRIEF DESCRIPTION OF DRAWINGS

[0082] In order to more clearly illustrate the technical scheme of the embodiments of the present application, the drawings needed in the embodiments of the present application will be briefly introduced as follows. Obviously, the drawings described below are only some embodiments of the present application, and other drawings can be obtained by those skilled in the art without creative labor on the basis of these drawings.

[0083] Figure 1 is a flowchart of the electromagnetic imaging method provided by the embodiments of the present application;

[0084] Figure 2 is a structural block diagram of the electromagnetic imaging system provided by the embodiments of the present application;

[0085] Figure 3 is a flowchart of the improved nonlinear conjugate gradient algorithm provided by the embodiments of the present application;

[0086] Figure 4 is a schematic diagram of using the FCM method in geophysical inversion provided by the embodiments of the present application;

[0087] Figure 5 is a schematic diagram of the model being divided into two geological units A (dark gray) and B (light gray) provided by the embodiments of the present application (the model units [a1, b1] and [a2, b2] have the same physical property values, but belong to different geological units);

[0088] Fig. 6 is a schematic diagram of a theoretical model and rock physical relationship provided by the embodiments of the present application;

[0089] Fig. 6(a) is a schematic diagram of a three-dimensional resistivity model provided by the embodiments of the present application;

[0090] Fig. 6(b) is a schematic diagram of a three-dimensional velocity model provided by the embodiments of the present application;

[0091] Figure 6(c) is a resistivity-velocity cross plot provided by an embodiment of the present application;

[0092] Figure 6(d) is a resistivity-velocity cross-variability function diagram provided by an embodiment of the present application;

[0093] Figure 7 is a schematic diagram of a theoretical model over y=4000 two-dimensional slice provided by an embodiment of the present application;

[0094] Figure 7(a) is a schematic diagram of a two-dimensional resistivity model provided by an embodiment of the present application;

[0095] Figure 7(b) is a schematic diagram of a two-dimensional velocity model provided by an embodiment of the present application;

[0096] Figure 8 Figure 7(c) is a schematic diagram of an inverted resistivity model over y=4000 two-dimensional slice provided by an embodiment of the present application;

[0097] Figure 9 Figure 7(d) is a schematic diagram of an inverted velocity model over y=4000 two-dimensional slice provided by an embodiment of the present application;

[0098] Figure 10 is a schematic diagram of a final magnetotelluric imaging based on rock physics modeling and rock physics relationship provided by an embodiment of the present application;

[0099] Figure 10(a) is a schematic diagram of a final magnetotelluric imaging result based on rock physics modeling provided by an embodiment of the present application;

[0100] Figure 10(b) is a comparison diagram of an inverted cross-variability function and a true cross-variability function provided by an embodiment of the present application;

[0101] Figure 10(c) is a schematic diagram of a comparison of an inverted result and a theoretical model cross plot provided by an embodiment of the present application;

[0102] Figure 11 Figure 10(d) is a drilling lithology map provided by an embodiment of the present application;

[0103] Figure 12 Figure 10(e) is a lithology velocity-resistivity cross plot of a formation in a study area provided by an embodiment of the present application, which shows a velocity-resistivity value cross plot for different lithologies of the formation;

[0104] Figure 13 Figure 10(f) is a different direction cross-variability function diagram provided by an embodiment of the present application;

[0105] Figure 14 Figure 10(g) is a flowchart of a prior resistivity model construction based on cluster analysis provided by an embodiment of the present application;

[0106] Figure 15The comparison of the electromagnetic imaging results provided by the embodiment of the present invention with the electromagnetic imaging interpretation results based on rock physics modeling;

[0107] Figure 16 This is a resistivity distribution map obtained by electromagnetic imaging of rock physics modeling provided by an embodiment of the present invention;

[0108] Figure 17 This is a comparison diagram of electromagnetic imaging results and prior information for rock physics modeling provided by an embodiment of the present invention;

[0109] In the figure: 1. Resistivity model construction module; 2. Initial velocity model construction module; 3. Velocity model construction module; 4. Initial resistivity model construction module; 5. Electromagnetic imaging module. DETAILED DESCRIPTION

[0110] In order to make the purpose, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with the embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not intended to limit the present invention.

[0111] In order to solve the problems existing in the prior art, the present invention provides an electromagnetic imaging method, system, medium, device and terminal. The present invention is described in detail below with reference to the accompanying drawings.

[0112] 1. Explanatory Examples In order to enable those skilled in the art to fully understand how to implement the present invention, this section provides an illustrative example that expands upon the technical solutions of the claims.

[0113] like Figure 1 As shown, the electromagnetic imaging method provided by an embodiment of the present invention includes: adding a resistivity-velocity cross-variogram fitting term to the traditional magnetotelluric inversion objective function to characterize the relationship between velocity and resistivity; introducing the fuzzy clustering concept into the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, in which the clustering inputs are resistivity and velocity, and by imposing constraints on the membership in the clustering process, the inversion is constrained by structural information, thereby realizing electro-seismic joint modeling inversion based on the dual constraints of rock physical relationship and structural similarity.

[0114] The specific steps include:

[0115] S101, constructing the resistivity model by using magnetotelluric fuzzy constraint inversion;

[0116] S102, establishing an initial velocity model using the obtained resistivity model and the mid-shallow seismic framework;

[0117] S103, constructing a seismic tomography inversion to obtain a velocity model;

[0118] S104, establishing an initial resistivity model by using the obtained velocity model and a resistivity-velocity model;

[0119] S105, obtaining electromagnetic imaging results by electromagnetic seismic joint modeling inversion based on the dual constraints of rock physical relationship and structural similarity.

[0120] As shown in Figure 2 The electromagnetic imaging system provided by the embodiment of the present application comprises:

[0121] A resistivity model construction module 1 is configured to construct a resistivity model by magnetotelluric fuzzy constraint inversion.

[0122] An initial velocity model construction module 2 is configured to establish an initial velocity model by using the obtained resistivity model and a shallow seismic structure.

[0123] A velocity model construction module 3 is configured to construct a velocity model by seismic tomography inversion.

[0124] An initial resistivity model construction module 4 is configured to establish an initial resistivity model by using the obtained velocity model and a resistivity-velocity model.

[0125] An electromagnetic imaging module 5 is configured to obtain electromagnetic imaging results by electromagnetic seismic joint modeling inversion based on the dual constraints of rock physical relationship and structural similarity.

[0126] The technical solutions of the present application will be further described in combination with specific embodiments.

[0127] In order to combine magnetotelluric and seismic data, make the best of the advantages of different geophysical data in resolution and sensitivity, and obtain electromagnetic imaging results with higher vertical resolution, the present application finally adopts the following process for electromagnetic imaging:

[0128] 1. A resistivity model is obtained by magnetotelluric fuzzy constraint inversion.

[0129] 2. An initial velocity model is established by using the resistivity model obtained in 1 and a shallow seismic structure.

[0130] 3. A velocity model is obtained by seismic tomography inversion.

[0131] 4. An initial resistivity model is established by using the velocity model obtained in 3 and a resistivity-velocity model.

[0132] 5. More reliable electromagnetic imaging results are obtained by electromagnetic seismic joint modeling inversion based on the dual constraints of rock physical relationship and structural similarity.

[0133] 1. Objective function of dual constraint electromagnetic seismic joint modeling inversion

[0134] Joint inversion of different geophysical data has attracted much attention in order to reduce the non-uniqueness of inversion. For joint inversion based on rock physics relations, it is critical to establish reasonable rock physics relations. The cross-variogram in multivariate geostatistics is integrated into the objective function of magnetotelluric inversion in the present invention to characterize the relationship between velocity and resistivity.

[0135] The variogram is a good geostatistical method for quantitatively characterizing a model, and can reflect the vertical and horizontal variation of physical properties. When the statistical characteristics of the physical properties of the subsurface geological body cannot be accurately known, the variogram can be obtained from indirect sources, such as geophysical data, outcrops, well logs and prior geological information.

[0136] The cross-variogram in multivariate geostatistics is integrated into the commonly used objective function equation of magnetotelluric inversion in the present invention to characterize the relationship between velocity and resistivity, and the cross-variogram fitting term is as follows:

[0137]

[0138] wherein γ ρv (m) represents the cross-variogram estimated from the inversion model, γ ρv obs is the prior cross-variogram obtained from prior rock physics information. The matrix W γ represents the uncertainty in the prior rock physics relation, W γ is a diagonal matrix whose diagonal elements are the inverse of the standard deviation of the cross-variogram misfit.

[0139] The objective function of the dual-constrained joint electrical and seismic modeling inversion is as follows:

[0140]

[0141] wherein Φ(m) is the total objective function, Φ d (m) represents the magnetotelluric data fitting term, Φ m (m) represents the model matching term, Φ γ (m) represents the cross-variogram fitting term, i.e. the resistivity-velocity relation fitting term; m is the model parameter, d is the measured MT data, C d is the covariance of data error, f(m) is the forward MT data, m0 is the initial model, and λ is the Lagrange operator (regularization parameter). C m is the model covariance, also known as the smoothing operator, The parameter α controls the closeness of the cross-variogram obtained by the present invention to the prior cross-variogram.

[0142] By adding the cross-variogram fitting term, the magnetotelluric inversion results can not only fit the observed magnetotelluric data, but also satisfy the priori resistivity-velocity relationship.

[0143] 2. Improved nonlinear conjugate gradient method

[0144] The inversion algorithm further improves the nonlinear conjugate gradient method, and the improved algorithm process is as shown in Figure 3

[0145] Compared with the algorithm process of the magnetotelluric ambiguity constraint inversion algorithm, the algorithm process of the double constraint joint modeling inversion of the electric and seismic has three improvements:

[0146] 1. Since the cross-variogram fitting term representing the resistivity-velocity relationship is added to the objective function, the gradient of the cross-variogram fitting term is added when calculating the gradient of the objective function;

[0147] 2. In the inversion process, the input of clustering is changed from resistivity to resistivity and velocity, so that the update of the model is affected by resistivity and velocity at the same time;

[0148] 3. The membership in the clustering process is spatially constrained, so that the model update is not only affected by the resistivity-velocity relationship, but also constrained by the spatial structure.

[0149] The specific improvements of the present application will be described in more detail below.

[0150] 2.1 Cross-variogram forward formula and Jacobian matrix derivation

[0151] Calculation of the gradient of the objective function is an important content of applying NLCG algorithm to solve magnetotelluric inversion. Because the cross-variogram fitting term is introduced into the objective function of the traditional magnetotelluric inversion in the present application, the Jacobian matrix of the cross-variogram needs to be derived in order to perform NLCG inversion.

[0152] In order to calculate the subsequent cross-variogram Jacobian matrix, the cross-variogram is rewritten as follows:

[0153]

[0154] Where i is the group number of different data groups divided according to different lags h; k and j represent the serial numbers of samples in each data group, and j≠k. For example, a group of sample numbers [m 1,j ,m 2,j ,m 3,j ,…] all have the same serial number j, but they belong to different data pairs.

[0155] ​The partial derivative of the cross-variogram with respect to the resistivity model of index j is

[0156]

[0157] Note that if the resistivity model of index j is not in the ith data pair, then and the ith data pair has no effect on Therefore, the equation (4) is rewritten as

[0158]

[0159] where N(h) j is the number of data pairs containing the resistivity model of index j.

[0160] The equation (5) can be written as

[0161]

[0162] The first term on the right side is expanded and combined, and we have

[0163]

[0164] Therefore, the partial derivative of the cross-variogram with respect to the resistivity model of index j is

[0165]

[0166] The equation (8) gives the value of the hth row and jth column of the Jacobian matrix of the cross-variogram, denoted as J γ:h,j

[0167] After obtaining the forward equation of the cross-variogram and the Jacobian matrix of the cross-variogram with respect to the resistivity model, the inversion can be performed by iteration.

[0168] 2.2 Cluster center and membership calculation

[0169] The output of the FCM clustering is the cluster center value and the membership of the petrophysical information. The cluster center value is the representative physical value of the model, and the membership contains the spatial information of the model. Therefore, when the clustering input is resistivity and velocity, the petrophysical information and structural information of the magnetotelluric and seismic models can be exchanged during the clustering process.

[0170] When the cluster center value and the membership are known, the model can be constructed as shown in Figure 4 Therefore, unlike constraining the entire model parameter, the cluster center value and the membership can be independently constrained.

[0171] Figure 4In the present application, the model parameters of the geophysical model can be decomposed into membership and cluster center value, the membership represents the possibility of the model element belonging to the cluster, and also contains the spatial information of the model, the cluster center value represents the representative physical value of the model.

[0172] 2.2.1 Cluster center

[0173] The present application has shown that introducing fuzzy clustering idea in traditional magnetotelluric inversion can well improve the inversion effect. In order to include the seismic model in the magnetotelluric inversion constrained by fuzzy clustering, the seismic model (S) is used together with the magnetotelluric model (M) to form the fuzzy clustering input data of Z = [M, S]. Therefore, the membership and cluster center value of the magnetotelluric model will also be affected by the seismic model, that is, the seismic information will indirectly affect the magnetotelluric model through the FCM clustering process.

[0174] 2.2.2 Boundary-constrained membership

[0175] The membership is a measure of the degree of belonging of an element to a certain class, which can be regarded as the possibility of a certain model unit belonging to a certain lithology or rock formation unit in geological classification or geophysical inversion. It has been pointed out by predecessors that the membership matrix of the elements in the inversion model contains spatial information. Therefore, the prior geometric information can be constrained by using the membership in the inversion process.

[0176] In the specific algorithm, the spatial information b is simply regarded as the classification feature data in the fuzzy clustering process. This function is combined with the model parameters to form the input data [m, b] for fuzzy clustering. Therefore, the distance between the jth model parameter and the kth rock cluster center can be represented as:

[0177]

[0178] In formula (9), when the model parameters represent resistivity and velocity, c k is the cluster center of the kth rock, w k is the weight coefficient; the first term defines the distance between the model parameter value and the rock cluster center value, and the second term contains the boundary information, which is defined as follows:

[0179]

[0180] Next, an example is given to illustrate how the spatial information constrains the clustering process. It is assumed that the model is divided into multiple geological units by the boundary, and the physical characteristics of these geological units are different. Figure 5The 16 model cells of the geophysical model during the inversion process are shown. The model cells belonging to two geological units A and B should have different physical property values. The physical property values of almost all model cells in group A are smaller than those in group B. However, two cells al and bl have the same physical property value of 2.6, while another two cells a2 and b2 have the same physical property value of 2.3. al and bl, a2 and b2 are located in different geological units, which leads to the fact that they cannot be classified into the correct geological units according to their physical property characteristics. Therefore, other information besides the physical property information is needed to help the inversion process to classify these model cells into different correct geological units. In this case, the spatial information can be used as additional constraint information for the clustering process, and the clustering of the model parameters is not only based on the physical property characteristics, but also constrained by the spatial information.

[0181] Table 1 The membership values of the FCM clustering results do not contain boundary information (membership A and membership B), while the color labels represent the model cells with the same physical property values but belonging to different geological units

[0182]

[0183] Table 1 shows the comparison between the fuzzy clustering results with and without the boundary information constraint. If the spatial information is not constrained during the clustering process, according to the membership, model cells al and bl are classified into geological unit B, and model cells a2 and b2 are classified into geological unit A. On the other hand, with the help of the further boundary information constraint, according to the membership, model cells al and a2 are classified into geological unit A, and model cells bl and b2 are classified into geological unit B. The clustering results with the boundary information constraint classify the model cells into the correct geological units, which shows that the model update is indeed constrained by the prior boundary information.

[0184] In addition, the clustering center values of the two clustering results are also not the same due to the influence of the boundary information on the membership change. In this example, the clustering centers are 2.08 and 2.90 without the boundary information constraint, and the clustering centers are 2.05 and 2.82 with the boundary information constraint. This shows that by constraining the membership with spatial information, not only the spatial distribution of the physical property parameters can be changed, but also the values of the physical property parameters can be affected. Both the geometric characteristics and the physical characteristics help to update the model in the inversion process, which means that the prior geological information can be incorporated into the inversion process. These information can guide the update process of the inversion to follow the geological information and establish a model that can be geologically interpreted.

[0185] This example also illustrates that expert knowledge can be incorporated into the geophysical inversion process simultaneously by fuzzy clustering. Because fuzzy clustering is used, imprecision and uncertainty can be tolerated. This means that "soft constraints" can be introduced in the geophysical inversion, which better resolves the contradiction between the prior information and the inversion model than "hard constraints".

[0186] 3、The present application provides an electromagnetic imaging method based on rock physics modeling. With rock physics model as a link, the method adds the constraint of seismic data, suppresses the multi-solution of magnetotelluric inversion, and improves the imaging effect of deep exploration. The method adds a cross-variability function fitting term to make the inversion result reflect the vertical and horizontal variation law of physical property parameters, and uses a fuzzy clustering method of boundary constraint to adaptively correct the model in the inversion iteration process, to constrain the rock physics relationship while increasing the structural constraint.

[0187] Theoretical model test shows that the electromagnetic imaging based on rock physics modeling is realized according to the process of the present application, compared with the magnetotelluric fuzzy constraint inversion method alone, the multi-solution of inversion can be reduced, the magnetotelluric imaging accuracy can be improved, the complementary characteristics between resistivity and velocity can be played, and the underground physical property distribution and structural characteristics can be highlighted.

[0188] II. Evidence of the effects of the embodiments. The embodiments of the present application have achieved some positive effects in the development or use process, and indeed have great advantages compared with the prior art. The following content is described in combination with the data, charts and the like of the test process.

[0189] The theoretical model test experiment provided by the embodiments of the present application is as follows:

[0190] The electromagnetic imaging method based on rock physics modeling is described in detail above. The present application further illustrates the implementation process of the inversion strategy and verifies its effectiveness through theoretical model test. A layered model is designed, which contains a strong shielding layer, a high and steep structure and a volcanic rock body, and has seven strata. The third stratum and the fourth stratum have wave impedance difference, forming a strong shielding layer. The models shown in Figs. 6(a) and 6(b) show the theoretical resistivity and velocity models respectively, which are the same in structure, and in the crossplot of Fig. 6(c), seven single points are regarded as circular clustering, and Fig. 6(d) shows the prior cross-variability function of the resistivity model and the velocity model.

[0191] Fig. 7 is a theoretical model over y=4000 two-dimensional slice, Fig. 7(a) is a two-dimensional resistivity model, and Fig. 7(b) is a two-dimensional velocity model.

[0192] The present application first performs magnetotelluric fuzzy constraint inversion, only using magnetotelluric data, and the initial model is a homogeneous half-space of 50Ω.m. Figure 8The resistivity model slice for the fuzzy constrained inversion of the magnetotelluric, the inversion result reflects the structural fluctuation of the real model, the high-speed volcanic rock body corresponds to the middle high uplift, but the accuracy of the description of the subsurface is limited, which shows that it is necessary to join the seismic data to further improve the quality of the electromagnetic imaging.

[0193] Due to the existence of the strong shielding layer, the seismic record is difficult to receive the signal from the deep part, in order to improve the reliability of the seismic tomographic inversion, the initial model of the seismic inversion is constructed by using the magnetotelluric inversion result, instead of selecting the progressive model as the initial model. Figure 9 The velocity model slice for the seismic tomographic inversion, the velocity inversion of the stratum above the salt dome top boundary is relatively accurate, and the inversion velocity below the strong shielding layer deviates from the real model velocity.

[0194] On the basis of the single inversion result of the magnetotelluric and the seismic, the electromagnetic-seismic joint modeling inversion is carried out by using the prior petrophysical relationship to improve the electromagnetic imaging effect. For the cross-variogram, only the prior information needs to be introduced into the objective function of the double-constrained electromagnetic-seismic joint modeling inversion (formula 2). For the fuzzy constrained inversion method, it has been explained in the part of the petrophysical modeling that different geometric structures need different distance measurements, as described above, the FCM clustering algorithm is suitable for circular clustering. In the inversion iteration process, the guided FCM clustering is selected, and the spatial information is introduced for constraint by using the membership constraint.

[0195] The resistivity model, the resistivity-velocity cross plot and the cross-variogram obtained by the final inversion are shown in Fig. 10. As shown in Fig. 10 (a), the shape and boundary of the seven strata and the volcanic rock body of the model have higher consistency with the real model in Fig. 7 (a). Compared with the single inversion, the layered structure of the background of the joint modeling inversion result is better, and Fig. 10 (b), (c) shows that the final inversion result conforms to the prior petrophysical relationship, and the inversion result is more reliable.

[0196] III. Application Examples. In order to prove the creativity and technical value of the technical scheme of the application, this part is the application examples of the technical scheme of the claims on the specific products or related technologies.

[0197] The application examples of the application are as follows:

[0198] The single magnetotelluric fuzzy constrained inversion and the electromagnetic imaging based on the petrophysical modeling are carried out based on the data of a certain zone of a certain basin piedmont belt, and the imaging effects are compared and analyzed.

[0199] The research area is located in a zone of the piedmont zone of Tarim Basin, and two sets of thick gravel layers are mainly developed. The gravel layers have the characteristics of developing thickly, and the lithology of the gravel layers has great changes, which leads to great lateral velocity variation and lateral discontinuity, and the deep strata cannot be accurately imaged. Influenced by multiple tectonics, the top and bottom interfaces of the gypsum salt layer are complex, the fracture zone under the salt layer is severely broken, and the thickness of the gypsum salt layer has great changes, which increases the difficulty of accurately identifying the deep strata and has an influence on the determination of drilling wells.

[0200] Through years of exploration, the predecessors have had certain understanding of identifying the gravel layers and the gypsum layers in the area by using the electrical method. Although the electrical data can reflect the gravel layers and the gypsum layers, it is difficult to finely depict the strata distribution by using only the electrical data. Although the seismic data cannot accurately reflect the distribution of the gravel layers and the gypsum layers, the framework of the strata can be correctly depicted. The gravel layers of the Neogene System are mainly coarse-grained, and have little influence on the velocity, so that the fine layer depiction can be carried out. The joint inversion of the electrical data and the seismic data can make the best of both worlds, so as to track the strata and realize the distribution depiction of the gravel layers and the gypsum layers.

[0201] Figure 11 The drilling lithology map shows the main lithology of different strata obtained by statistical drilling data, and the drilling data mainly includes seven stratum information.

[0202] Figure 12 The stratum lithology velocity-resistivity crossplot of the research area shows the velocity-resistivity crossplot of different stratum lithology. The crossplot shows that the joint seismic and magnetotelluric data can reflect the spatial distribution of the gravel layers and the gypsum layers, because the resistivity-velocity relationship of the gravel rock and the gypsum layer is obviously distinguished from the surrounding rock. At the same time, it can also be seen that the resistivity-velocity relationship of the area is suitable for the petrophysical modeling based on fuzzy clustering.

[0203] Through the prior information such as logging data, outcrop data and the like, cross-variogram in different directions such as Figure 13 The cross-variogram in different directions is shown in the figure.

[0204] The present application uses the result of single magnetotelluric fuzzy constraint inversion to combine logging information and seismic data to obtain a comprehensive interpretation result, and then constructs a new resistivity model as an initial model of dual constraint electroseismic joint modeling inversion. The steps of constructing the prior resistivity model are illustrated by L1 line of three wells:

[0205] 1. As shown in Figure 14 The flow chart of constructing the prior resistivity model based on clustering analysis shows that the resistivity model ( Figure 14 a) is obtained by magnetotelluric fuzzy constraint inversion;

[0206] 2. PSDM data ( Figure 14b) While normal reflections can be formed for Neogene coarse-grained gravel layers, it is difficult to determine the distribution of gravel and gypsum layers. Calibration is performed using well logging information, and the PSDM data is combined with the resistivity model derived from magnetotelluric fuzzy constraint inversion to obtain interpretation results for key horizons.

[0207] 3. Velocity model obtained by seismic inversion ( Figure 14 c);

[0208] 4. Using rock physics relationships to obtain resistivity values ​​as the initial model for dual-constrained electric-seismic joint modeling inversion ( Figure 14 d).

[0209] After comprehensive analysis of the well logging data, core outcrops, and one-dimensional magnetotelluric block inversion results of the work area, and based on the purpose of identifying gravel layers and gypsum-salt layers, the present invention ultimately determined the number of cluster centers to be 6, the velocity cluster centers to be 4050 m / s, 4500 m / s, 5200 m / s, 5500 m / s, 5700 m / s, and 6000 m / s, and the resistivity cluster centers to be 25 Ω·m, 50 Ω·m, 100 Ω·m, 165 Ω·m, 275 Ω·m, and 676 Ω·m as prior information for the clustering method to optimize the inversion results.

[0210] Using the above prior information, electromagnetic imaging based on rock physics modeling is implemented, and the cross-variogram function is established using the well logging data and the velocity and resistivity profiles. Figure 15 The comparison of electromagnetic imaging results with the electromagnetic imaging interpretation results based on rock physics modeling shows that the obtained resistivity model is as follows Figure 15 a The electromagnetic imaging results based on rock physics modeling are shown, which is different from the single magnetotelluric fuzzy constraint inversion ( Figure 15 b) Compared with the results of magnetotelluric fuzzy constraint inversion imaging, the resistivity profile obtained by the inversion strategy of the present invention has better stratification, especially the division of the target layer gravel layer and gypsum salt layer is clearer, and the inversion result is more consistent with the physical property information. The dual-constrained electro-seismic joint modeling inversion combines the advantages of seismology and electromagnetics, that is, the advantage of seismology is to obtain the stratigraphic framework, and the advantage of electromagnetics is to invert the volume distribution of electrical properties. The two have different focuses, and the comprehensive interpretation results are as follows. Figure 15 c Comparison of electromagnetic imaging results and electromagnetic imaging interpretation results based on rock physics modeling.

[0211] The corresponding histogram of the inverted resistivity model is shown in Figure 16 In the resistivity distribution map obtained from electromagnetic imaging of rock physics modeling, six tight clusters can be easily identified, whose resistivity values ​​in the histogram are close to the target cluster center. Figure 17 As shown in the figure, the electromagnetic imaging results of rock physics modeling are compared with the prior information. Figure 17a shows that the cross-variability function obtained from the joint inversion of the electrical and seismic data is consistent with the prior cross-variability function. Figure 13 Figure 17 b shows that the inversion results are consistent with the resistivity logging data, which confirms the reliability of the inversion strategy proposed in the present application.

[0212] It should be noted that the embodiments of the present application can be realized by hardware, software, or a combination of software and hardware. The hardware part can be realized by special logic; the software part can be stored in a memory and executed by a suitable instruction execution system, such as a microprocessor or a specially designed hardware. Those skilled in the art can understand that the above-mentioned devices and methods can be realized by computer executable instructions and / or included in processor control codes, for example, such codes are provided on a carrier medium, such as a magnetic disk, a CD or a DVD-ROM, a programmable memory, such as a read-only memory (firmware), or a data carrier, such as an optical or electronic signal carrier. The devices of the present application and their modules can be realized by hardware circuits, such as very large scale integrated circuits or gate arrays, semiconductors, such as logic chips, transistors, etc., or programmable hardware devices, such as field programmable gate arrays, programmable logic devices, etc., by software executed by various types of processors, or by a combination of the above-mentioned hardware circuits and software, such as firmware.

[0213] The above description is only a specific implementation of the present application, but the protection scope of the present application is not limited thereto, and any modification, equivalent replacement and improvement within the technical range disclosed in the present application and within the spirit and principle of the present application should be covered within the protection scope of the present application.​

Claims

1. An electromagnetic imaging method, characterized in that: The electromagnetic imaging method includes: adding a resistivity-velocity cross-variogram fitting term to the magnetotelluric inversion objective function to characterize the relationship between velocity and resistivity; utilizing fuzzy clustering introduced into the NLCG algorithm to obtain an improved nonlinear conjugate gradient algorithm, wherein the clustering inputs are resistivity and velocity; and imposing constraints on the membership degree in the clustering process so that the inversion is constrained by structural information, thereby achieving electro-seismic joint modeling inversion based on the dual constraints of rock physical relationships and structural similarity; The cluster center and membership degree calculation in the electromagnetic imaging method includes: The output of FCM clustering is the cluster center value and membership degree of rock physical information; the cluster center value is the representative physical property value of the model, and the membership degree contains the spatial information of the model. When the clustering input is resistivity and velocity, the rock physical information and structural information of the magnetotelluric and seismic models are exchanged during the clustering process. When the cluster center value and membership degree are known, the model is constructed. Instead of constraining the entire model parameters, the cluster center value and membership degree variables are constrained independently. (1) Cluster center The seismic model S is used together with the magnetotelluric model M to form fuzzy clustering input data Z = [M, S], and the seismic model is included in the magnetotelluric inversion constrained by fuzzy clustering. The membership degree and cluster center value of the magnetotelluric model are affected by the seismic model, and the seismic information indirectly affects the magnetotelluric model through the FCM clustering process. (2) Boundary constraint membership Membership is a measure of the degree to which an element belongs to a certain class. In geological classification or geophysical inversion, it is considered the probability that a model unit belongs to a certain lithology or lithostratigraphic unit. The membership matrix of the elements in the inversion model contains spatial information, and the membership is used to constrain the prior geometric information during the inversion process. The spatial information b is simply regarded as categorical feature data in the fuzzy clustering process; the features are combined with model parameters to form the input data for fuzzy clustering , so the distance between the j-th model parameter and the k-th rock cluster center is expressed as: ； Where, when the model parameters represent resistivity and velocity, ; is the cluster center of the k-th rock, w k is the weight coefficient; the first term defines the distance between the model parameter value and the rock cluster center value, and the second term contains the boundary information, which is defined as follows: 。 2. The electromagnetic imaging method according to claim 1, wherein: The electromagnetic imaging method comprises the following steps: Step 1: construct the magnetotelluric fuzzy constraint inversion to obtain the resistivity model; Step 2: Use the obtained resistivity model and the mid-shallow seismic framework to establish an initial velocity model; Step 3: construct seismic tomography inversion to obtain a velocity model; Step 4: Using the obtained velocity model and resistivity-velocity model to establish an initial resistivity model; Step 5: Electromagnetic imaging results are obtained by inversion of the joint electro-seismic modeling based on the dual constraints of rock physics relationship and structural similarity.

3. The electromagnetic imaging method according to claim 1, wherein: The objective function of the dual-constrained electric-seismic joint modeling inversion in the electromagnetic imaging method includes: The cross-variogram in multivariate geostatistics is integrated into the traditional magnetotelluric inversion objective function equation to characterize the relationship between velocity and resistivity. The cross-variogram fitting term is as follows: ； in, represents the cross-variogram estimated from the inversion model, is the prior-cross variogram obtained from the prior petrophysical information; the matrix represents the uncertainty in the a priori rock physics relationships, is a diagonal matrix whose diagonal elements are the inverse of the standard deviation of the cross-variogram misfit; The objective function of the dual-constrained electro-seismic joint modeling inversion is: ； in, is the overall objective function, represents the magnetotelluric data fitting term, represents the model matching item; represents the cross-variogram fitting term, which is the resistivity-velocity relationship fitting term; are model parameters, It is the measured MT data. is the covariance of the data errors, , is the forward MT data, is the initial model; is the Lagrangian operator, is the regularization parameter; is the model covariance, is the smoothing operator; , ,parameter Controls how close the desired inverted cross-variogram is to the prior cross-variogram; By adding the cross-variogram fitting term, the MT inversion results are made to fit the observed MT data while satisfying the a priori resistivity-velocity relationship.

4. The electromagnetic imaging method according to claim 1, wherein: The improved nonlinear conjugate gradient method in the electromagnetic imaging method includes: Compared with the MT fuzzy constrained inversion algorithm, the dual-constrained electro-seismic joint modeling inversion algorithm includes the following three improvements: (1) Since the cross-variogram fitting term that characterizes the resistivity-velocity relationship is added to the objective function, the gradient of the cross-variogram fitting term is added when calculating the objective function gradient; (2) During the inversion process, the input of clustering is changed from resistivity to resistivity and velocity, so that the update of the model is affected by both resistivity and velocity; (3) Spatial constraints are imposed on membership during the clustering process, so that the model update is affected by the resistivity-velocity relationship and the spatial structure constraints.

5. The electromagnetic imaging method according to claim 1, wherein: The cross-variogram forward formula and Jacobian matrix derivation in the electromagnetic imaging method include: Rewrite the cross-variogram as follows: ； Where i is the group number of different data groups divided according to different lag distances h; k and j represent the sequence numbers of samples in each data group. ; is the number of pairs of data comprising the resistivity model; The differential form of the cross-variogram function with respect to each resistivity model with serial number j is: ； If the resistivity model with serial number j is not in the i-th data pair, then And the i-th data pair does not have If there is any influence, the resistivity model differential form can be rewritten as: ； in, is the number of pairs of data containing the resistivity model with sequence number j; Then the resistivity model differential form is rewritten as: ； Expand and merge the first term on the right side to get: ； The partial derivative of the cross-variogram with respect to the resistivity model with serial number j is: ； The value obtained is the value of the hth row and jth column of the Jacobian matrix of the cross-variogram function, which is recorded as ; After obtaining the forward formula of the cross-variogram and the Jacobian matrix of the cross-variogram to resistivity model, the inversion is performed using an iterative solution method.

6. An electromagnetic imaging system using the electromagnetic imaging method according to any one of claims 1 to 5, characterized in that: The electromagnetic imaging system comprises: Resistivity model building module, used to build the resistivity model obtained by magnetotelluric fuzzy constraint inversion; An initial velocity model building module is used to build an initial velocity model using the obtained resistivity model and the medium-shallow seismic framework; Velocity model building module, used to build velocity models obtained by seismic tomography inversion; An initial resistivity model building module is used to build an initial resistivity model using the obtained velocity model and resistivity-velocity model; The electromagnetic imaging module is used to obtain electromagnetic imaging results through joint electro-seismic modeling and inversion based on the dual constraints of rock physical relationships and structural similarity.

7. A computer device, characterized in that: The computer device includes a memory and a processor, the memory stores a computer program, and when the computer program is executed by the processor, the processor performs the steps of the electromagnetic imaging method according to any one of claims 1 to 5.

8. A computer-readable storage medium storing a computer program, wherein when the computer program is executed by a processor, the processor is caused to perform the steps of the electromagnetic imaging method according to any one of claims 1 to 5.

9. An information data processing terminal, characterized in that: The information data processing terminal is used to implement the electromagnetic imaging system as claimed in claim 6.

Images