High-resolution wave impedance inversion method and device based on spectrum fusion

Abstract

The invention discloses a high-resolution wave impedance inversion method and device based on spectrum fusion, and the method comprises the steps: carrying out the Kriging spatial interpolation of a low-frequency band of logging data under the constraint of horizon data and a spatial variation function, and obtaining a low-frequency-band wave impedance data volume zL; according to the seismic wavelets and the Cauchy norm constraint sparse pulse inversion intermediate frequency model, Cauchy norm constraint sparse pulse inversion is carried out on the seismic data to obtain a reflection coefficient, and channel integration is carried out to obtain an intermediate frequency band relative wave impedance data volume zM; under the constraint of the horizon data and the spatial variation function, performing sequential Gaussian stochastic simulation on the high frequency band of the logging data to obtain a high frequency band wave impedance data volume zH; and performing Fourier transform on zL, zM and zH to obtain frequency domain data volumes ZL, ZM and ZH, performing spectrum fusion to obtain a full-band wave impedance data volume ZMGI, performing Fourier inversion to obtain a high-resolution wave impedance inversion result zMGI, and effectively performing high-resolution wave impedance inversion.

Description

Technical Field

[0001] The present invention relates to the field of geophysical exploration technology, and in particular to a high-resolution wave impedance inversion method and device based on spectrum fusion. Background Art

[0002] This section is intended to provide a background or context to the embodiments of the invention recited in the claims. No admission is made that the description herein is prior art by inclusion in this section.

[0003] Seismic inversion is a process that uses seismic, well logging, geological and other data to reveal the spatial structure and physical properties of underground rock formations. It is currently one of the most effective methods for finding oil and gas reservoirs. Wave impedance (the product of velocity and density) is an important formation attribute that has the ability to link geology and geophysics. The wave impedance obtained through seismic inversion can eliminate the influence of seismic wavelets and convert the interface information in seismic data into formation information, thereby characterizing the spatial distribution characteristics of the reservoir and evaluating the reservoir properties.

[0004] Seismic impedance inversion can generally be divided into two categories: deterministic inversion and random inversion. Deterministic inversion refers to obtaining a single optimal solution estimate under specific solution constraints. This type of method has high computational efficiency, but is limited by the seismic frequency band range. The inversion results cannot identify thin reservoir information outside the seismic frequency band range, and cannot evaluate the uncertainty of the inversion results; random inversion uses a sampling strategy to fully explore the posterior solution space and can output high-resolution inversion results. It has great application potential in thin reservoir prediction and has attracted much attention in recent years.

[0005] Depending on whether the inverse problem is linear or nonlinear, random inversion can be divided into two categories: iterative optimization nonlinear random inversion and Bayesian linearized random inversion. In terms of iterative optimization nonlinear random inversion, Haas and Dubrule (1994) first proposed geostatistical inversion. This method uses seismic data as a constraint to control the selection of random simulation results of logging data in the full frequency band. The calculation efficiency is low, and the channel-by-channel calculation method makes the inversion results have poor lateral continuity and lack of low-frequency information in the vertical direction, which is inconsistent with the actual geological conditions; Sambridge and Mosegaard (2002) proposed a seismic inversion method based on Monte Carlo random sampling. This method also needs to go through a large number of "rejection-acceptance" processes and has low calculation efficiency; Azevedo et al. (2020) combined co-simulation with direct sequential Gaussian random simulation to perform random realization of the entire volume, which improved the lateral continuity of the inversion results, but faced the problem of difficulty in convergence. In the Bayesian linearized stochastic inversion, it is usually assumed that the model parameters obey a multidimensional Gaussian distribution, and that the model parameters and seismic data have a linear forward relationship. Then, the analytical expression of the posterior probability density function of the model parameters can be obtained under the Bayesian theory, and then the posterior probability density function is sampled to obtain the random inversion result (e.g., Buland and Omre, 2003; Hansen et al., 2006; Grana et al., 2021). This method has the advantage of high computational efficiency, but when the entire data is calculated, if the spatial correlation in the horizontal direction is considered, the memory overhead of the model covariance matrix will become very large, and it is difficult to implement on a conventional performance computer. Therefore, this method can generally only consider the longitudinal variogram during inversion, and the random realization result has poor lateral continuity, and the assumptions of the Gaussian prior model and the linear forward operator are usually inconsistent with the actual situation.

[0006] Judging from the published literature, in the field of seismic exploration, there is currently no effective method to integrate seismic inversion and geostatistics in the frequency domain to perform high-resolution wave impedance inversion. Summary of the invention

[0007] The embodiment of the present invention provides a high-resolution wave impedance inversion method based on spectrum fusion, which is used to improve the accuracy of traditional deterministic inversion, solve the problem that traditional random inversion requires a large number of iterations and has low computational efficiency, and effectively perform high-resolution wave impedance inversion. The method includes:

[0008] Acquire seismic data, seismic wavelets, well logging data and horizon data;

[0009] Under the constraints of the layer data and the spatial variation function, the low-frequency band of the logging data is interpolated by Kriging to obtain the low-frequency band wave impedance data volume z L ;

[0010] A Cauchy norm constrained sparse pulse inversion intermediate frequency model is established. According to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, the seismic data is subjected to Cauchy norm constrained sparse pulse inversion to obtain the reflection coefficient. The reflection coefficient is integrated to obtain the intermediate frequency band relative wave impedance data volume z. M ;

[0011] Under the constraints of the layer data and the spatial variation function, a sequential Gaussian random simulation is performed on the high-frequency band of the logging data to obtain the high-frequency band wave impedance data volume z H ;

[0012] The low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ;

[0013] According to the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ;

[0014] For the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0015] The embodiment of the present invention further provides a high-resolution wave impedance inversion device based on spectrum fusion, which is used to improve the accuracy of traditional deterministic inversion, solve the problem that traditional random inversion requires a large number of iterations and has low computational efficiency, and effectively perform high-resolution wave impedance inversion. The device includes:

[0016] A data acquisition module, used to acquire seismic data, seismic wavelets, well logging data and layer data;

[0017] The Kriging spatial interpolation module is used to perform Kriging spatial interpolation on the low-frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the low-frequency band wave impedance data volume z L ;

[0018] The Cauchy norm constrained sparse pulse inversion module is used to establish the Cauchy norm constrained sparse pulse inversion intermediate frequency model. According to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, the seismic data is subjected to Cauchy norm constrained sparse pulse inversion to obtain the reflection coefficient; the reflection coefficient is integrated to obtain the intermediate frequency band relative wave impedance data volume z M ;

[0019] The sequential Gaussian random simulation module is used to perform sequential Gaussian random simulation on the high frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the high frequency band wave impedance data volume z H ;

[0020] Fourier transform module is used to transform the low frequency impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ;

[0021] The spectrum fusion module is used to fusion the frequency domain data volume Z according to the selected frequency parameters and scale parameters. L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ;

[0022] The inverse Fourier transform module is used to transform the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0023] An embodiment of the present invention also provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor implements the above-mentioned high-resolution wave impedance inversion method based on spectrum fusion when executing the computer program.

[0024] An embodiment of the present invention further provides a computer-readable storage medium, wherein the computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the high-resolution wave impedance inversion method based on spectrum fusion is implemented.

[0025] An embodiment of the present invention further provides a computer program product, which includes a computer program. When the computer program is executed by a processor, the high-resolution wave impedance inversion method based on spectrum fusion is implemented.

[0026] In the embodiment of the present invention, seismic data, seismic wavelets, well logging data and stratum data are obtained; under the constraints of the stratum data and the spatial variation function, Kriging spatial interpolation is performed on the low frequency band of the well logging data to obtain a low frequency band wave impedance data volume z L ; Establish a Cauchy norm constrained sparse pulse inversion intermediate frequency model, perform Cauchy norm constrained sparse pulse inversion on seismic data according to the selected seismic wavelet and Cauchy norm constrained sparse pulse inversion intermediate frequency model, and obtain the reflection coefficient; perform trace integration on the reflection coefficient to obtain the intermediate frequency band relative wave impedance data volume z M Under the constraints of the layer data and the spatial variation function, a sequential Gaussian random simulation is performed on the high-frequency band of the logging data to obtain the high-frequency band wave impedance data volume z H ; respectively for the low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ; According to the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ; For the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0027] By combining seismic inversion with geostatistical theory and adopting the idea of ​​"divide and conquer", the broadband inversion results are decoupled to obtain low-frequency, medium-frequency and high-frequency results respectively, and finally integrated to obtain high-precision wave impedance inversion results, which improves the accuracy of traditional deterministic inversion and solves the problem that traditional random inversion requires a large number of iterations and has low computational efficiency. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the prior art descriptions. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work. In the drawings:

[0029] Figure 1 is a processing flow chart of a high-resolution wave impedance inversion method based on spectrum fusion in an embodiment of the present invention;

[0030] Figure 2It is a flow chart of the method for establishing a Cauchy norm constrained sparse pulse inversion intermediate frequency model in an embodiment of the present invention;

[0031] Figure 3 This is a schematic diagram of the fusion of wave impedance data volumes and spectra in different frequency bands in an embodiment of the present invention;

[0032] Figure 4 is a logic flow chart of a high-resolution wave impedance inversion method based on spectrum fusion in an embodiment of the present invention;

[0033] Figure 5 This is a schematic diagram of an actual seismic profile of an oil field in an embodiment of the present invention;

[0034] Figure 6 (a) is the result of the traditional sparse pulse inversion in an embodiment of the present invention; Figure 6 (b) is the high-resolution wave impedance inversion result in an embodiment of the present invention;

[0035] Figure 7 (a) is a schematic diagram of the RMS amplitude slice of the inversion result with a high cutoff frequency of 120 Hz in an embodiment of the present invention, Figure 7 (b) is a schematic diagram of the RMS amplitude slice of the inversion result with a high cutoff frequency of 150 Hz in an embodiment of the present invention;

[0036] Figure 8 Schematic diagram of the structure of a high-resolution wave impedance inversion device based on spectrum fusion in an embodiment of the present invention;

[0037] Fig. 9 The figure is a schematic diagram of the structure of a computer device according to an embodiment of the present invention. DETAILED DESCRIPTION

[0038] To make the purpose, technical solution and advantages of the embodiments of the present invention more clear, the embodiments of the present invention are further described in detail below in conjunction with the accompanying drawings. Here, the exemplary embodiments of the present invention and their descriptions are used to explain the present invention, but are not intended to limit the present invention.

[0039] The applicant found that, regardless of deterministic inversion or random inversion, under the action of seismic wavelet convolution (filtering), the inversion results in the mid-frequency band controlled by the seismic data are consistent, that is, the high-frequency band of the inversion results will not affect the degree of matching between the synthetic record and the actual seismic record, which is also the embodiment of multi-solution in geophysical inverse problems.

[0040] Accordingly, in order to solve the problems of low calculation efficiency of traditional geostatistical inversion, poor lateral continuity of inversion results, missing low-frequency information in the vertical direction, and inconsistency with actual geological conditions, the applicant combines deterministic seismic inversion with geostatistical Kriging estimation and random simulation technology, and comprehensively utilizes logging, seismic and geological spatial structure information to provide a high-resolution wave impedance inversion method and device based on spectral fusion, aiming to improve the resolution of wave impedance inversion results, thereby improving the prediction accuracy of thin reservoirs in oilfield development.

[0041] The specific technical solution of this application is introduced below.

[0042] Figure 1 FIG. 4 is a processing flow chart of a high-resolution wave impedance inversion method based on spectrum fusion in an embodiment of the present invention. Figure 1 As shown, the high-resolution wave impedance inversion method based on spectrum fusion in the embodiment of the present invention may include:

[0043] Step 101, obtaining seismic data, seismic wavelets, well logging data and layer data;

[0044] Step 102: Under the constraints of the layer data and the spatial variation function, the low-frequency band of the logging data is subjected to Kriging spatial interpolation to obtain the low-frequency band wave impedance data volume z L ;

[0045] Step 103: Establish a Cauchy norm constrained sparse pulse inversion intermediate frequency model, perform Cauchy norm constrained sparse pulse inversion on the seismic data according to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, and obtain the reflection coefficient; perform trace integration on the reflection coefficient to obtain the intermediate frequency band relative wave impedance data volume z M ;

[0046] Step 104: Under the constraints of the layer data and the spatial variation function, sequential Gaussian random simulation is performed on the high frequency band of the logging data to obtain the high frequency band wave impedance data volume z H ;

[0047] Step 105: respectively analyze the low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ;

[0048] Step 106: According to the selected frequency parameter and scale parameter, the frequency domain data volume Z is L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume ZMGI ;

[0049] Step 107: Full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0050] After obtaining the seismic data, seismic wavelets, logging data and stratigraphic data of the specified underground rock formation, the Kriging low-frequency background model can be established, the low-frequency soft-constrained sparse pulse inversion intermediate-frequency volume can be calculated, the sequential Gaussian random simulation high-frequency model can be established, and the frequency domain integration can be performed.

[0051] Kriging spatial interpolation technology is a generalized regression method that obtains an estimation result with minimum error variance according to the least squares criterion. It not only considers the distance between the known point and the point to be estimated, but also considers the geological correlation between the points. It is widely used in the field of earth science.

[0052] Under different data assumptions, there are different types of Kriging, such as simple Kriging, ordinary Kriging, trend Kriging, etc. Simple Kriging is suitable for estimating stationary variables with known and constant regional means, while ordinary Kriging is suitable for estimating non-stationary variables with unknown regional means. During the calculation process, the local mean is automatically calculated from the neighborhood. Since the wave impedance is non-stationary in space, the present invention adopts the ordinary Kriging method to perform spatial interpolation on the wave impedance logging curve under the constraint of the layer to obtain the low-frequency wave impedance data volume z L .

[0053] For step 102, in one embodiment, Kriging spatial interpolation is performed on the low frequency band of the logging data to obtain a low frequency band wave impedance data volume z L , which may include: using the ordinary Kriging method to perform Kriging spatial interpolation on the low-frequency band of the logging data to obtain the low-frequency band wave impedance data volume z L .

[0054] For step 103, the process of establishing the intermediate frequency model of the Cauchy norm constrained sparse pulse inversion is described below:

[0055] Figure 2 Flow chart of the method for establishing the intermediate frequency model of the Cauchy norm constrained sparse pulse inversion in an embodiment of the present invention. Figure 2 As shown, in one embodiment, the process of establishing the Cauchy norm constrained sparse pulse inversion intermediate frequency model includes:

[0056] Step 201, constructing an objective function with the goal of minimizing the error between seismic data and synthetic seismic records, wherein the seismic data is the actually observed seismic data, and the synthetic seismic records are obtained by convolution of seismic wavelets and reflection coefficients;

[0057] Step 202: Based on the objective function, add the Cauchy norm constraint of the reflection coefficient to obtain a new objective function;

[0058] Step 203: Obtain the derivative of the new objective function with respect to the reflection coefficient, and set the new objective function to 0, to obtain the Cauchy norm constrained sparse pulse inversion intermediate frequency model.

[0059] In specific implementation, assuming that the post-stack seismic record can be described by the Robinson convolution model, we have

[0060] d obs =Wr+n (1)

[0061] Among them, d obs represents the observed seismic data, W represents the wavelet convolution matrix, r represents the reflection coefficient, and n is the noise term.

[0062] The purpose of seismic inversion is to obtain the observed seismic data from obs The reflection coefficient r is estimated so that the error between the synthetic record and the actual observation data is minimized. At this time, the objective function can be written as:

[0063] min J(r)=(Wr-d obs ) T (Wr-d obs ) (2)

[0064] Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, that is, the minimum value of the objective function, d obs represents seismic data, W represents the seismic wavelet convolution matrix, and r represents the reflection coefficient.

[0065] However, due to the influence of noise in the observed data, directly solving the above objective function will cause strong instability in the inversion result, and it is necessary to further add constraints on the reflection coefficient in the objective function. Compared with the L2 norm constraint, the Cauchy norm constraint can produce sparse reflection coefficients with full frequency band and better noise resistance. The present invention adds the Cauchy norm constraint of the reflection coefficient on the basis of the traditional inversion objective function. The new objective function is

[0066]

[0067] Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, J s (r) is the seismic data constraint term, which controls the similarity between the synthetic seismic record and the seismic data. Cauc h y(r) is the Cauchy norm constraint term of the reflection coefficient, which controls the sparsity of the reflection coefficient, μ is the weight coefficient of the sparse constraint, n is the number of sampling points of the seismic trace, i represents the sampling point sequence number of the seismic trace, and σ represents the variance of the reflection coefficient.

[0068] Take the derivative of the new objective function (Equation 3) with respect to the reflection coefficient r and set it to 0. Then we can get

[0069] (W T W+μQ)r=W T d obs (4)

[0070] Among them, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, μ is the weight coefficient of the sparse constraint, Q is a diagonal weighted matrix, and the diagonal elements i represents the sampling point number of the seismic trace, and σ represents the variance of the reflection coefficient. Since the matrix Q contains the reflection coefficient r to be calculated, Formula 4 is an implicit nonlinear equation group and cannot be directly calculated. The present invention uses iterative reweighted least squares (IRLS) to solve it.

[0071] Finally, the relative wave impedance can be obtained by performing trace integral calculation on the reflection coefficient inversion result. The generated mid-frequency relative impedance data volume is denoted as z M Among them, z(t) represents the relative wave impedance, r represents the reflection coefficient, t represents the time corresponding to a certain sampling point of the seismic trace, and τ represents the seismic trace sampling point before time point t.

[0072] For step 104, the specific process of establishing the sequential Gaussian random simulation high-frequency model is described below:

[0073] Geostatistical random simulation can generate multiple high-resolution geological models that conform to known data and underground spatial structural characteristics. Each model has equal probability and can be used to establish static and dynamic reservoir models and conduct uncertainty assessment. It is one of the key core technologies for fine reservoir description. For continuous variables (such as wave impedance), sequential Gaussian simulation is the most commonly used geostatistical random simulation algorithm, which perfectly combines the Gaussian random function model with the sequential simulation algorithm. The advantage of the Gaussian random function is that it can be fully described only by the mean vector and the covariance matrix (calculated by the variogram); the sequential simulation algorithm can simplify the problem of obtaining a complex high-dimensional space joint distribution to the problem of obtaining a low-dimensional space conditional distribution, thereby realizing random sampling of high-dimensional data.

[0074] The present invention uses sequential Gaussian simulation to establish a high-frequency model for logging data under the control of the horizon. The specific implementation steps are:

[0075] ①Perform normal transformation on logging data;

[0076] Since the basic assumption of the sequential Gaussian simulation algorithm is that the data obeys the Gaussian distribution, it is necessary to first perform a normal (Gaussian) transformation on the logging data to make it meet the basic assumption of the sequential Gaussian simulation. This step also obtains the logging data, but the data distribution has changed and will be transformed back in the end. This step is used to prepare the logging data for sequential Gaussian simulation, thereby obtaining a high-frequency body. This is the step of establishing a high-frequency model. In the specific embodiment of the present invention, only the high-frequency body of the sequential Gaussian model is embodied, and the step of establishing a high-frequency model is not further explored.

[0077] ②Determine a random path and visit all nodes;

[0078] Each sampling point in the three-dimensional data volume can be called a "node", and the path is the order in which all samples are calculated. This order is random; for example, the 100th sampling point is calculated first, and then the 1000th sampling point may be calculated next; this random path is the calculation order of all sampling points, which is the same as the Gaussian distribution and is also required by the sequential Gaussian simulation algorithm.

[0079] ③ At each node, conditional data consisting of the wave impedance of adjacent hard data points (logging data) and simulated points are obtained, and the conditional cumulative probability distribution function (CCDF) of the current node is calculated through Kriging and variogram;

[0080] The hard data points here refer to the points actually measured by logging data, which are equivalent to known points.

[0081] ④ Randomly obtain a value from CCDF as the value of the node and add the simulated value to the data set;

[0082] The data set here refers to the collection of hard data points and previously calculated points. These points are considered as known points when calculating the next point. The data set is the collection of known points. When all points (wave impedance) are calculated, the wave impedance data volume of the high frequency band is obtained.

[0083] ⑤Return to step ③ until all nodes are calculated;

[0084] ⑥ The simulation results are transformed into the original space through inverse normal transformation, and finally the high-frequency wave impedance data volume is obtained, which is recorded as z H .

[0085] For step 105, after obtaining the low-frequency band wave impedance data volume z according to the above steps, L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z HAfterwards, the low-frequency impedance data volume z can be L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H Then, in step 106, the frequency domain data volume Z is processed according to the selected frequency parameter and scale parameter. L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI .

[0086] Figure 3 FIG. 1 is a schematic diagram of the fusion of impedance data volumes of different frequency bands in an embodiment of the present invention. Based on this, in one embodiment, the frequency domain data volume Z can be fused according to the selected frequency parameter and scale parameter according to the following formula: L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI :

[0087] Z MGI =Z L H(0,0,f1,f2)+αZ M ·H(f1,f2,f3,f4)+βZ H H(f1,f2,f3,f4)

[0088] Among them, Z MGI It represents the full-band wave impedance data volume, α represents the scale factor of the mid-frequency data volume, β represents the scale factor of the high-frequency data volume, which is used to adjust the amplitude of the components in each frequency band, H represents a filter with an amplitude of 1, and f1, f2, f3, and f4 represent the frequency parameters of the frequency domain filter.

[0089] Finally, in step 107, the full-band impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0090] To simply and clearly present the logical method of the present invention, Figure 4 FIG. 2 shows a logic flow chart of a high-resolution wave impedance inversion method based on spectrum fusion in an embodiment of the present invention. Figure 4 As shown, it is helpful to quickly understand the specific processing logic of the high-resolution wave impedance inversion method based on spectrum fusion in an embodiment of the present invention.

[0091] The present invention can obtain high-resolution wave impedance inversion results and has a good application effect in the prediction of thin reservoirs in oil field development. Figure 5 Schematic diagram of an actual seismic profile of an oil field in an embodiment of the present invention. Figure 5 As shown, the curve inserted in the profile is the wave impedance logging curve. The main frequency of the seismic data is about 30 Hz. The target layer is a thin sandstone reservoir. Due to the limitation of the seismic frequency band, these thin reservoirs cannot be directly identified in the seismic profile.

[0092] Figure 6 (a) is the result of the traditional sparse pulse inversion in the embodiment of the present invention. Figure 6 (b) is the high-resolution wave impedance inversion result in the embodiment of the present invention. From the arrow, it can be seen that the resolution of the wave impedance inversion result obtained by the present invention is significantly higher than the traditional sparse pulse inversion result, and can clearly depict the spatial distribution characteristics of thin reservoirs.

[0093] Figure 7 (a) is a schematic diagram of the RMS amplitude slice of the inversion result with a high cutoff frequency of 120 Hz in an embodiment of the present invention, Figure 7 (b) is a schematic diagram of the RMS amplitude slice of the inversion result with a high cutoff frequency of 150 Hz in an embodiment of the present invention. Figure 7 (a) and Figure 7 (b) shows the comparison of the RMS amplitude slices of the inversion results with different high cutoff frequencies. From the comparison between the two, it can be seen that although the high-frequency components of the inversion results are different, they have good consistency on the slices as a whole, and the lateral distribution is clear and natural, without the "bull's eye" modeling phenomenon, which also reflects that the lateral information of the inversion results of the present invention mainly depends on seismic data and has good predictability.

[0094] The present invention also provides a high-resolution wave impedance inversion device based on spectrum fusion in an embodiment, as described in the following embodiment. Since the principle of solving the problem by the device is similar to that of the high-resolution wave impedance inversion method based on spectrum fusion, the implementation of the device can refer to the implementation of the high-resolution wave impedance inversion method based on spectrum fusion, and the repeated parts will not be repeated.

[0095] Figure 8 FIG. 1 is a schematic diagram of the structure of a high-resolution wave impedance inversion device based on spectrum fusion in an embodiment of the present invention. Figure 8 As shown, the use risk detection device for an electric vehicle in an embodiment of the present invention may specifically include:

[0096] Data acquisition module 801, used to acquire seismic data, seismic wavelets, well logging data and layer data;

[0097] The Kriging spatial interpolation module 802 is used to perform Kriging spatial interpolation on the low frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the low frequency band wave impedance data volume z L ;

[0098] The Cauchy norm constrained sparse pulse inversion module 803 is used to establish a Cauchy norm constrained sparse pulse inversion intermediate frequency model, perform Cauchy norm constrained sparse pulse inversion on the seismic data according to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, and obtain the reflection coefficient; perform trace integration on the reflection coefficient to obtain the intermediate frequency band relative wave impedance data volume z M ;

[0099] The sequential Gaussian random simulation module 804 is used to perform sequential Gaussian random simulation on the high frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the high frequency band wave impedance data volume z H ;

[0100] Fourier transform module 805 is used to transform the low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ;

[0101] The spectrum fusion module 806 is used to fusion the frequency domain data volume Z according to the selected frequency parameters and scale parameters. L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ;

[0102] The inverse Fourier transform module 807 is used to transform the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0103] In one embodiment, the Kriging spatial interpolation module 802 is specifically used for:

[0104] The ordinary Kriging method is used to perform Kriging spatial interpolation on the low-frequency band of the logging data to obtain the low-frequency band wave impedance data volume z L .

[0105] In one embodiment, the Cauchy norm constrained sparse pulse inversion module 803 is specifically used to:

[0106] The objective function is constructed with the goal of minimizing the error between seismic data and synthetic seismic records, where seismic data refers to the actual observed seismic data and synthetic seismic records are obtained by convolution of seismic wavelets and reflection coefficients.

[0107] Based on the objective function, the Cauchy norm constraint of the reflection coefficient is added to obtain a new objective function;

[0108] The derivative of the new objective function with respect to the reflection coefficient is obtained, and the new objective function is set to 0, and the Cauchy norm constrained sparse pulse inversion intermediate frequency model is obtained.

[0109] In one embodiment, the objective function can be expressed as follows:

[0110] minJ(r)=(Wr-d obs ) T (Wr-d obs )

[0111] Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, that is, the minimum value of the objective function, d obs represents seismic data, W represents the seismic wavelet convolution matrix, and r represents the reflection coefficient.

[0112] In one embodiment, the new objective function can be expressed as follows:

[0113]

[0114] Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, J s (r) is the seismic data constraint term, which controls the similarity between the synthetic seismic record and the seismic data. Cauc h y (r) is the Cauchy norm constraint term of the reflection coefficient, which controls the sparsity of the reflection coefficient, μ is the weight coefficient of the sparse constraint, n is the number of sampling points of the seismic trace, i represents the sampling point sequence number of the seismic trace, and σ represents the variance of the reflection coefficient.

[0115] In one embodiment, the Cauchy norm constrained sparse pulse inversion intermediate frequency model is expressed as follows:

[0116] (W T W+μQ)r=W T d obs

[0117] Among them, d obsrepresents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, μ is the weight coefficient of the sparse constraint, Q is a diagonal weighted matrix, and the diagonal elements i represents the sampling point number of the seismic trace, and σ represents the variance of the reflection coefficient.

[0118] In one embodiment, the Cauchy norm constrained sparse pulse inversion module 803 is specifically used for:

[0119] Perform trace integration on the reflection coefficient to obtain the relative wave impedance;

[0120] According to the relative wave impedance, the mid-frequency band relative wave impedance data volume z is obtained. M .

[0121] In one embodiment, the Cauchy norm constrained sparse pulse inversion module 803 is specifically used to perform trace integration on the reflection coefficient according to the following formula to obtain the relative wave impedance:

[0122] z(t)=exp[2∫0 t r(τ)dτ]

[0123] Among them, z(t) represents the relative wave impedance, r represents the reflection coefficient, t represents the time corresponding to a certain sampling point of the seismic trace, and τ represents the seismic trace sampling point before time point t.

[0124] In one embodiment, the spectrum fusion module 806 is specifically used to: according to the following formula, based on the selected frequency parameter and scale parameter, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI :

[0125] Z MGI =Z L H(0,0,f1,f2)+αZ M ·H(f1,f2,f3,f4)+βZ H H(f1,f2,f3,f4)

[0126] Among them, Z MGI It represents the full-band wave impedance data volume, α represents the scale factor of the mid-frequency data volume, β represents the scale factor of the high-frequency data volume, which is used to adjust the amplitude of the components in each frequency band, H represents a filter with an amplitude of 1, and f1, f2, f3, and f4 represent the frequency parameters of the frequency domain filter.

[0127] Based on the above invention concept, Fig. 9As shown, the present invention also proposes a computer device 900, including a memory 910, a processor 920, and a computer program 930 stored in the memory 910 and executable on the processor 920, wherein the processor 920 implements the aforementioned high-resolution wave impedance inversion method based on spectral fusion when executing the computer program 930.

[0128] An embodiment of the present invention further provides a computer-readable storage medium, wherein the computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the high-resolution wave impedance inversion method based on spectrum fusion is implemented.

[0129] An embodiment of the present invention further provides a computer program product, which includes a computer program. When the computer program is executed by a processor, the high-resolution wave impedance inversion method based on spectrum fusion is implemented.

[0130] In summary, in the embodiment of the present invention, seismic data, seismic wavelets, well logging data and stratum data are obtained; under the constraints of the stratum data and the spatial variation function, Kriging spatial interpolation is performed on the low-frequency band of the well logging data to obtain the low-frequency band wave impedance data volume z L ; Establish a Cauchy norm constrained sparse pulse inversion intermediate frequency model, perform Cauchy norm constrained sparse pulse inversion on seismic data according to the selected seismic wavelet and Cauchy norm constrained sparse pulse inversion intermediate frequency model, and obtain the reflection coefficient; perform trace integration on the reflection coefficient to obtain the intermediate frequency band relative wave impedance data volume z M Under the constraints of the layer data and the spatial variation function, a sequential Gaussian random simulation is performed on the high-frequency band of the logging data to obtain the high-frequency band wave impedance data volume z H ; respectively for the low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ; According to the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ; For the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

[0131] By combining seismic inversion with geostatistical theory and adopting the idea of ​​"divide and conquer", the broadband inversion results are decoupled to obtain low-frequency, medium-frequency and high-frequency results respectively, and finally integrated to obtain high-precision wave impedance inversion results, which improves the accuracy of traditional deterministic inversion and solves the problem that traditional random inversion requires a large number of iterations and has low computational efficiency.

[0132] Those skilled in the art will appreciate that embodiments of the present invention may be provided as methods, systems, or computer program products. Therefore, the present invention may take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0133] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to produce a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0134] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0135] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0136] The specific embodiments described above further illustrate the objectives, technical solutions and beneficial effects of the present invention in detail. It should be understood that the above description is only a specific embodiment of the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A high-resolution wave impedance inversion method based on spectrum fusion, characterized in that: include: Acquire seismic data, seismic wavelets, well logging data and horizon data; Under the constraints of the layer data and the spatial variation function, the low-frequency band of the logging data is interpolated by Kriging to obtain the low-frequency band wave impedance data volume z L ; A Cauchy norm constrained sparse pulse inversion intermediate frequency model is established. According to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, the seismic data is subjected to Cauchy norm constrained sparse pulse inversion to obtain the reflection coefficient. The reflection coefficient is integrated to obtain the intermediate frequency band relative wave impedance data volume z. M ; Under the constraints of the layer data and the spatial variation function, a sequential Gaussian random simulation is performed on the high-frequency band of the logging data to obtain the high-frequency band wave impedance data volume z H ; The low frequency band wave impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ; According to the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ; For the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

2. The method according to claim 1, characterized in that Perform Kriging spatial interpolation on the low-frequency band of the logging data to obtain the low-frequency band wave impedance data volume z L ,include: The ordinary Kriging method is used to perform Kriging spatial interpolation on the low-frequency band of the logging data to obtain the low-frequency band wave impedance data volume z L .

3. The method according to claim 1, characterized in that Establish a Cauchy norm constrained sparse pulse inversion intermediate frequency model, including: The objective function is constructed with the goal of minimizing the error between seismic data and synthetic seismic records, where seismic data refers to the actual observed seismic data and synthetic seismic records are obtained by convolution of seismic wavelets and reflection coefficients. Based on the objective function, the Cauchy norm constraint of the reflection coefficient is added to obtain a new objective function; The derivative of the new objective function with respect to the reflection coefficient is obtained, and the new objective function is set to 0, and the Cauchy norm constrained sparse pulse inversion intermediate frequency model is obtained.

4. The method according to claim 3, characterized in that The objective function is expressed as follows: minJ(r)=(Wr-d obs ) T (Wr-d obs ) Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, that is, the minimum value of the objective function, d obs represents seismic data, W represents the seismic wavelet convolution matrix, and r represents the reflection coefficient.

5. The method according to claim 3, characterized in that The new objective function is expressed as follows: Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, J s (r) is the seismic data constraint term, which controls the similarity between the synthetic seismic record and the seismic data. Cauchy (r) is the Cauchy norm constraint term of the reflection coefficient, which controls the sparsity of the reflection coefficient, μ is the weight coefficient of the sparse constraint, n is the number of sampling points of the seismic trace, i represents the sampling point sequence number of the seismic trace, and σ represents the variance of the reflection coefficient.

6. The method according to claim 3, characterized in that The Cauchy norm constrained sparse pulse inversion intermediate frequency model can be expressed as follows: (W T W+μQ)r=W T d obs Among them, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, μ is the weight coefficient of the sparse constraint, Q is a diagonal weighted matrix, and the diagonal elements i represents the sampling point number of the seismic trace, and σ represents the variance of the reflection coefficient.

7. The method according to claim 1, characterized in that Perform channel integration on the reflection coefficient to obtain the mid-frequency relative wave impedance data volume z M ,include: Perform trace integration on the reflection coefficient to obtain the relative wave impedance; According to the relative wave impedance, the mid-frequency band relative wave impedance data volume z is obtained. M .

8. The method according to claim 7, characterized in that This includes integrating the reflection coefficient according to the following formula to obtain the relative wave impedance: z(t)=exp[2∫0 t r(τ)dτ] Among them, z(t) represents the relative wave impedance, r represents the reflection coefficient, t represents the time corresponding to a certain sampling point of the seismic trace, and τ represents the seismic trace sampling point before time point t.

9. The method according to claim 1, characterized in that Including according to the following formula, based on the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI : Z MGI =Z L ·H(0,0,f1,f2)+αZ M ·H(f1,f2,f3,f4)+βZ H ·H(f1,f2,f3,f4) Among them, Z MGI It represents the full-band wave impedance data volume, α represents the scale factor of the mid-frequency data volume, β represents the scale factor of the high-frequency data volume, which is used to adjust the amplitude of the components in each frequency band, H represents a filter with an amplitude of 1, and f1, f2, f3, and f4 represent the frequency parameters of the frequency domain filter.

10. A high-resolution wave impedance inversion device based on spectrum fusion, characterized in that: include: A data acquisition module, used to acquire seismic data, seismic wavelets, well logging data and layer data; The Kriging spatial interpolation module is used to perform Kriging spatial interpolation on the low-frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the low-frequency band wave impedance data volume z L ; The Cauchy norm constrained sparse pulse inversion module is used to establish the Cauchy norm constrained sparse pulse inversion intermediate frequency model. According to the selected seismic wavelet and the Cauchy norm constrained sparse pulse inversion intermediate frequency model, the seismic data is subjected to Cauchy norm constrained sparse pulse inversion to obtain the reflection coefficient; the reflection coefficient is integrated to obtain the intermediate frequency band relative wave impedance data volume z M ; The sequential Gaussian random simulation module is used to perform sequential Gaussian random simulation on the high frequency band of the logging data under the constraints of the layer data and the spatial variation function to obtain the high frequency band wave impedance data volume z H ; Fourier transform module is used to transform the low frequency impedance data volume z L , mid-frequency relative impedance data volume z M , high frequency band wave impedance data volume z H Perform Fourier transform to obtain the corresponding frequency domain data volume Z L , Z M , Z H ; The spectrum fusion module is used to fusion the frequency domain data volume Z according to the selected frequency parameters and scale parameters. L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI ; The inverse Fourier transform module is used to transform the full-band wave impedance data volume Z MGI Perform inverse Fourier transform to obtain high-resolution wave impedance inversion result z MGI .

11. The device according to claim 10, characterized in that The Kriging spatial interpolation module is specifically used for: The ordinary Kriging method is used to perform Kriging spatial interpolation on the low-frequency band of the logging data to obtain the low-frequency band wave impedance data volume z L .

12. The device according to claim 10, characterized in that The Cauchy norm constrained sparse pulse inversion module is specifically used for: The objective function is constructed with the goal of minimizing the error between seismic data and synthetic seismic records, where seismic data refers to the actual observed seismic data and synthetic seismic records are obtained by convolution of seismic wavelets and reflection coefficients. Based on the objective function, the Cauchy norm constraint of the reflection coefficient is added to obtain a new objective function; The derivative of the new objective function with respect to the reflection coefficient is obtained, and the new objective function is set to 0, and the Cauchy norm constrained sparse pulse inversion intermediate frequency model is obtained.

13. The device according to claim 12, characterized in that The objective function is expressed as follows: minJ(r)=(Wr-d obs ) T (Wr-d obs ) Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, that is, the minimum value of the objective function, d obs represents seismic data, W represents the seismic wavelet convolution matrix, and r represents the reflection coefficient.

14. The device according to claim 12, characterized in that The new objective function is expressed as follows: Where minJ(r) represents the minimum error between seismic data and synthetic seismic records, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, J s (r) is the seismic data constraint term, which controls the similarity between the synthetic seismic record and the seismic data. Cauc h y (r) is the Cauchy norm constraint term of the reflection coefficient, which controls the sparsity of the reflection coefficient, μ is the weight coefficient of the sparse constraint, n is the number of sampling points of the seismic trace, i represents the sampling point sequence number of the seismic trace, and σ represents the variance of the reflection coefficient.

15. The device according to claim 12, characterized in that The Cauchy norm constrained sparse pulse inversion intermediate frequency model is expressed as follows: (W T W+μQ)r=W T d obs Among them, d obs represents seismic data, W represents the seismic wavelet convolution matrix, r represents the reflection coefficient, μ is the weight coefficient of the sparse constraint, Q is a diagonal weighted matrix, and the diagonal elements i represents the sampling point number of the seismic trace, and σ represents the variance of the reflection coefficient.

16. The device according to claim 10, characterized in that The Cauchy norm constrained sparse pulse inversion module is specifically used for: Perform trace integration on the reflection coefficient to obtain the relative wave impedance; According to the relative wave impedance, the mid-frequency band relative wave impedance data volume z is obtained. M .

17. The device according to claim 16, characterized in that The Cauchy norm constrained sparse pulse inversion module is specifically used to perform trace integration on the reflection coefficient according to the following formula to obtain the relative wave impedance: z(t)=exp[2∫0 t r(τ)dτ] Among them, z(t) represents the relative wave impedance, r represents the reflection coefficient, t represents the time corresponding to a certain sampling point of the seismic trace, and τ represents the seismic trace sampling point before time point t.

18. The device according to claim 10, characterized in that The spectrum fusion module is specifically used to: According to the following formula, based on the selected frequency parameters and scale parameters, the frequency domain data volume Z L , Z M , Z H Perform spectrum fusion to obtain the full-band wave impedance data volume Z MGI : Z MGI =Z L ·H(0,0,f1,f2)+αZ M ·H(f1,f2,f3,f4)+βZ H ·H(f1,f2,f3,f4) Among them, Z MGI It represents the full-band wave impedance data volume, α represents the scale factor of the mid-frequency data volume, β represents the scale factor of the high-frequency data volume, which is used to adjust the amplitude of the components in each frequency band, H represents a filter with an amplitude of 1, and f1, f2, f3, and f4 represent the frequency parameters of the frequency domain filter.

19. A computer device comprising a memory, a processor and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the method according to any one of claims 1 to 9 is implemented.

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

21. A computer program product, characterized in that The computer program product comprises a computer program, and when the computer program is executed by a processor, the method according to any one of claims 1 to 9 is implemented.

Images