Lithology identification method and system based on independent component and amplitude offset

Abstract

The invention discloses a lithology identification method and system based on an independent component and an amplitude offset. The lithology identification method comprises the following steps: establishing a pseudo well model containing noise according to rock physical data; determining a seismic reflection coefficient sequence based on the noise-containing pseudo well model; determining sparse pulse reflection and a weight matrix from the seismic reflection coefficient sequence; and comparing the weight matrix with a preset rock physical template, and establishing a longitudinal lithology sequence of the input result. Through application of ICA, sparse pulse reflection can be extracted from the AVO gather more efficiently, the signal weight result is processed, and through comparison and classification with the preset rock physical template, the purpose of rapid lithology identification is achieved.

Description

Technical Field

[0001] The invention belongs to the technical field of lithology identification, and in particular relates to a lithology identification method and system based on independent components and amplitude with offset. Background Art

[0002] In oil and gas exploration tasks, lithology identification is an important basis for reservoir division, a prerequisite for accurate determination of porosity and oil saturation, the basis for reservoir characteristic research, reserve calculation and geological modeling, and a strong evidence for subsequent lithology spatial distribution target determination and complex three-dimensional modeling. Therefore, lithology prediction is of great significance in the exploration and development of oil and gas.

[0003] As exploration in complex areas becomes increasingly difficult, using only post-stack seismic information can no longer fully meet the needs of single-factor prediction of lithology, physical properties and gas content in refined exploration. Although pre-stack seismic data has more dimensional information than post-stack data, the most commonly used method in the industry is to obtain parameters with strong correlation with lithology through pre-stack simultaneous inversion to achieve the purpose of lithology identification, but the data volume is usually large and the computational efficiency is low. There are few mature algorithms for identifying lithology from pre-stack seismic data.

[0004] Independent component analysis (ICA) originated in the 1980s and is an unsupervised machine learning algorithm for feature analysis. During this period, amplitude versus offset analysis (AVO) was also more widely used by geophysicists in the field of oil and gas exploration. Although some scholars have applied the ICA algorithm in the field of geophysics, there is a lack of research linking ICA and AVO. Summary of the invention

[0005] In view of the above problems, the present invention discloses a lithology identification method based on independent components and amplitude with offset, comprising:

[0006] A pseudo-well model with noise is established based on petrophysical data;

[0007] determining a seismic reflection coefficient sequence based on the pseudo-well model containing noise;

[0008] determining a sparse pulse reflection and a weight matrix from the seismic reflection coefficient sequence;

[0009] The weight matrix is ​​compared with a preset rock physics template to establish a vertical lithology sequence of the input results.

[0010] Furthermore, the rock physical data includes porosity, longitudinal wave velocity, shear wave velocity and density.

[0011] Furthermore, the establishment of a pseudo-well model containing noise according to petrophysical data comprises the following steps:

[0012] Select rocks of different lithologies and measure rock physics data;

[0013] Set the number of lithology interfaces;

[0014] Assign rock physical data to the upper and lower layers of each lithology interface to obtain a pseudo-well model without noise;

[0015] Noise is added to the pseudo-well model without noise to obtain the pseudo-well model with noise.

[0016] Furthermore, the lithology includes normal porosity water-bearing sandstone, high porosity water-bearing sandstone, tight gas-bearing sandstone, normal porosity water-bearing limestone and normal porosity water-bearing dolomite.

[0017] Furthermore, the specific steps of determining the seismic reflection coefficient sequence based on the pseudo well model containing noise are as follows:

[0018] Based on the pseudo-well model containing noise, a seismic reflection coefficient sequence with pseudo-seismic incident angles ranging from 0° to 40° is determined.

[0019] Furthermore, the determining of sparse pulse reflections and weight matrices from the seismic reflection coefficient sequence comprises the following steps:

[0020] Perform centralization and whitening preprocessing on the seismic reflection coefficient sequence to generate a random orthogonalized initial matrix

[0021] For the initial matrix Perform orthogonalization iterations;

[0022] For the initial matrix Perform decorrelation processing to obtain the new process matrix W n+1 ;

[0023] Determine the initial matrix and the new process matrix W n+1 Are they mutually replaced? If yes, the iteration process ends and the new process matrix W n+1 Transpose to get the weight matrix; if not, continue orthogonalization iteration;

[0024] Sparse pulse reflections are determined based on a weight matrix and an input noisy signal.

[0025] Furthermore, the sparse pulse reflection is obtained by the following formula:

[0026] X=T - 1*W T *s e +μ

[0027] Where X is the input noisy target signal; T-1 is the whitening preprocessing operator; W T is the weight matrix; s e is the sparse pulse reflection; μ is the central preprocessing residual.

[0028] Furthermore, the specific steps of comparing the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input result are as follows:

[0029] Combine different lithologies in pairs to determine the seismic reflection coefficient sequence for seismic incident angles ranging from 0 to 40°;

[0030] By intersecting the incident angle and the seismic reflection coefficient intensity, the preset AVO change curve corresponding to the lithology combination is obtained, that is, the preset rock physics template;

[0031] Determine the similarity index between the AVO curve of the n-layer lithology interface in the weight matrix and all the preset AVO change curves;

[0032] Arrange all convergent lithology combinations of n-layer lithology interfaces;

[0033] Determine the sum of similarity indices of all convergent lithology combinations of n-layer lithology interfaces, and select the maximum similarity index and the corresponding lithology combination as the output result;

[0034] The redundant information of the upper and lower layers of the lithology interface in the output results is merged to obtain a multi-layer lithology medium model, that is, a vertical lithology sequence.

[0035] The present invention also discloses a lithology identification system based on independent components and amplitude with offset, comprising:

[0036] A building unit is used to build a pseudo-well model containing noise according to rock physical data;

[0037] A determination unit, configured to determine a seismic reflection coefficient sequence based on the pseudo-well model containing noise;

[0038] A weight matrix unit for determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence;

[0039] The comparison unit is used to compare the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input results.

[0040] Furthermore, the establishing unit is specifically used for:

[0041] Select rocks of different lithologies and measure rock physics data;

[0042] Set the number of lithology interfaces;

[0043] Assign rock physical data to the upper and lower layers of each lithology interface to obtain a pseudo-well model without noise;

[0044] Noise is added to the pseudo-well model without noise to obtain the pseudo-well model with noise.

[0045] Compared with the prior art, the embodiments of the present invention have at least the following advantages: through the use of ICA, sparse pulse reflections can be extracted from AVO gathers more efficiently, and the signal weight results can be processed, and the purpose of rapid lithology identification can be achieved by comparing and classifying with preset rock physics templates.

[0046] Other features and advantages of the present invention will be described in the following description, and partly become apparent from the description, or understood by practicing the present invention. The purpose and other advantages of the present invention can be realized and obtained by the structures pointed out in the description and the drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] 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 description of the prior art. Obviously, the drawings described below are some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0048] Figure 1 A flow chart of a lithology identification method according to an embodiment of the present invention is shown;

[0049] Figure 2 A schematic diagram of a Gaussian function with an intercept window and intercept outlets with different sample point lengths according to an embodiment of the present invention is shown;

[0050] Figure 3 (a) shows a pseudo well model without noise generated based on a three-layer horizontal lithologic medium model according to an embodiment of the present invention;

[0051] Figure 3 (b) shows a pseudo well model containing noise generated based on a three-layer horizontal lithologic medium model according to an embodiment of the present invention;

[0052] Figure 4 (a) shows a seismic reflection coefficient sequence according to an embodiment of the present invention;

[0053] Figure 4 (b) shows a 40 Hz Ricker wavelet according to an embodiment of the present invention;

[0054] Figure 4 (c) shows a forward synthesis record according to an embodiment of the present invention;

[0055] Figure 5 An ICA algorithm iteration flow chart according to an embodiment of the present invention is shown;

[0056] Figure 6 A schematic diagram showing when the number of predicted independent components is equal to the number of actual sparse pulse reflections according to an embodiment of the present invention;

[0057] Figure 7 A schematic diagram showing a case where the number of predicted independent components is less than the number of actual sparse pulse reflections according to an embodiment of the present invention;

[0058] Figure 8 A schematic diagram showing a case where the number of predicted independent components is greater than the number of actual sparse pulse reflections according to an embodiment of the present invention;

[0059] Fig. 9 A flow chart of estimating the number of independent components in a target signal according to an embodiment of the present invention is shown;

[0060] Fig.10 A schematic diagram of an ICA algorithm result according to an embodiment of the present invention is shown;

[0061] Fig.11 A schematic diagram of similarity calculation of two curves according to an embodiment of the present invention is shown;

[0062] Fig.12 A schematic diagram showing comparison and sorting of a weight matrix result of a noise-free four-layer horizontal lithologic medium model and all preset AVO change curves according to an embodiment of the present invention is shown;

[0063] Fig.13 A flow chart for establishing a multi-layer lithologic medium model according to an embodiment of the present invention is shown. DETAILED DESCRIPTION

[0064] In order to make the purpose, technical solution and advantages of the embodiments of the present invention clearer, the technical solution in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0065] like Figure 1 As shown, the present invention proposes a lithology identification method based on independent components and amplitude with offset, comprising the following steps:

[0066] Step S100: establishing a pseudo well model containing noise according to rock physical data;

[0067] Step S200: determining a seismic reflection coefficient sequence based on the pseudo-well model containing noise;

[0068] Step S300: determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence by an ICA algorithm;

[0069] Step S400: Compare the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input results.

[0070] Through the use of ICA, sparse pulse reflections and reflection coefficients can be extracted from AVO gathers more efficiently, and then the lithology identification can be completed more accurately according to the preset rock physics template based on the reflection coefficient. The results show that the ICA algorithm can effectively extract sparse pulse reflection signals from AVO gathers to achieve the role of lithofacies identification.

[0071] Step S100: establishing a pseudo well model containing noise according to rock physical data;

[0072] In order to better interpret the technical process of this method, a multi-layer medium pseudo-well model with noise is first established as the research object based on the measured rock physical data. Five different lithologies were selected in this study, and their combination covers all AVO types, namely: general porosity water-bearing sandstone, high porosity water-bearing sandstone, tight gas-bearing sandstone, general porosity water-bearing limestone, and general porosity water-bearing dolomite. Their measured average P-wave velocity, S-wave velocity and density data are listed in Table 1. Rock physical data include porosity, P-wave velocity, S-wave velocity and density. Among them, general porosity refers to a porosity of 10% to 15%, high porosity refers to a porosity greater than 15%, and tight refers to a porosity less than 10%.

[0073] Table 1 Lithological elastic parameters

[0074]

[0075] After selecting the lithology and its related data, this study adopts the following flow chart ( Figure 2 ) is a two-step process to establish a pseudo-well model of multi-layer media with noise. Figure 2 (a) is the Gaussian function with intercept window used in the study; Figure 2 (b) is the correlated Gaussian noise with an intercept window length of 10 sample points; Figure 2 (c) is the correlated Gaussian noise with an intercept window length of 20 sample points; Figure 2 (d) is the correlated Gaussian noise with an intercept window length of 40 sample points. First, the number of lithologic interfaces in a pseudo-well model is randomly determined, and then the lithologic elastic parameter data (P-wave velocity, S-wave velocity, and density) in Table 1 are assigned to the upper and lower layers of each lithologic interface, thereby establishing a pseudo-well model without noise.

[0076] Secondly, in order to add noise to the model built in the first step, the present invention defines a correlation function equation f(x, w). After studying different intercept windows w, the present invention selects the number of sample points w=20 as the standard value for establishing the noise model ( Figure 2 (c)).

[0077]

[0078] S′=S*f (2)

[0079] Where S is a random number sequence generated using Gaussian distribution; S′ is correlated Gaussian noise.

[0080] Then, the random numbers that conform to the normal distribution are convolved with the equation, and the band energy of the noise-free model and the noise model is adjusted according to the defined signal-to-noise ratio. The signal-to-noise ratio is the ratio of signal energy to noise energy. The higher the signal-to-noise ratio, the more obvious the signal energy. What is known is the signal energy. Given the signal-to-noise ratio value, the band energy of the noise can be obtained.

[0081] Before adding noise to the noise-free pseudo-well model, its energy is proportionally adjusted according to a given SNR (signal to noise ratio, ie, the bandwidth energy ratio (band power) of the noise channel and the pseudo-well model channel). Figure 3 (a) is a pseudo-well model generated based on a three-layer horizontal lithologic medium model (water-bearing sandstone-gas-bearing sandstone-water-bearing sandstone); Figure 3 (b) is the noisy pseudo-well model generated after adding correlated Gaussian noise. Figure 3 The addition of noise channels to the noise-free pseudo-well model to generate the noisy pseudo-well model used in this study is explained.

[0082] Step S200: determining a seismic reflection coefficient sequence based on the pseudo-well model containing noise;

[0083] Finally, based on the established multi-layer medium pseudo-well model with noise, the Zoeppritz formula (Formula 3) is used to calculate the seismic reflection coefficient sequence (e.g. Figure 4 (a)), and then use the 40Hz Ricker wavelet (as shown in Figure 4 (b)), by convolving the reflection coefficient and the wavelet, a forward synthetic record based on the pseudo-well model can be produced (as shown in Figure 4 (c)).

[0084]

[0085] in, x = sinθ;

[0086] ρ1 is the density of the upper medium in the double-layer lithology model; ρ2 is the density of the lower medium in the double-layer lithology model; Vp1 is the P-wave velocity of the upper medium in the double-layer lithology model; Vp2 is the P-wave velocity of the lower medium in the double-layer lithology model; Vs1 is the S-wave velocity of the upper medium in the double-layer lithology model; Vs2 is the S-wave velocity of the lower medium in the double-layer lithology model; θ is the incident angle; R pp is the longitudinal wave reflection coefficient; R ps is the converted shear wave reflection coefficient; T pp is the transmission longitudinal wave coefficient; T ps is the transmitted shear wave coefficient.

[0087] Step S300: determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence by an ICA algorithm;

[0088] In the blind source separation problem processed by ICA, the target signal matrix x m×n The matrix s is considered to be a series of independent components j×n After the transformation matrix α m×j In this study, the seismic pre-stack angle gather is considered to be the target signal x m×n (θ, t), which is the sparse pulse reflection s from the underlying strata j×n (t) and the characteristic curve matrix α containing the amplitude variation of each layer's emission as the angle changes m×j (θ) is mixed and superimposed:

[0089] x m×n (θ, t) = α m×j (θ)*s j×n (t) (4)

[0090] Where t is time.

[0091] An important theoretical basis of ICA is the central limit theory: when a series of independent components are superimposed, the result will tend to be more normally distributed. In view of this, the above formula is changed as follows:

[0092]

[0093] Where X is the input noisy target signal; W T is the weight matrix; W is W T The transpose of .

[0094] Since the matrix s contains a series of independent components, any other transformation matrix W T α will make the final result X closer to the normal distribution than s, unless X is completely equal to s, so if the non-Gaussianity of X can be maximized, then the target signal x m×nThe independent component s in (θ, t) j×n (t) can be calculated.

[0095] In this study, a variation based on kurtosis (Formula 6) is used to calculate the non-Gaussianity of a matrix. The technical process of solving the non-Gaussianity maximization matrix mainly includes two steps: data preprocessing and ICA iteration process.

[0096] G=E[(X-μ) 4 ]-3(E[(X-μ) 2 ]) 2 (6)

[0097] Among them, G is non-Gaussian; μ is the centered preprocessing residual.

[0098] Before starting the algorithm, the input signal needs to be preprocessed in two steps: centering and whitening. The centering process changes the matrix X into a zero-mean matrix Xc (Formula 7); the whitening preprocessing process further changes the matrix after centering into a unit-variance matrix Xcw (Formula 7). Through these two steps of preprocessing, Formula 6 can be simplified to Formula 8. Similarly, the target signal X calculated at the end is the same as the independent component matrix s calculated. e The relationship is as described in Formula 9:

[0099] X=Xc+μ;Xc=T -1 *Xcw (7)

[0100] G * =E[(Xcw) 4 ]-3 (8)

[0101] X=T -1 *W T *s e +μ (9)

[0102] Among them, T -1 is the whitening preprocessing operator; G * is the simplified non-Gaussian property; X is the input noisy target signal; W T is the weight matrix; s e It is a sparse pulse reflection.

[0103] After preprocessing, the following ICA algorithm can be understood as Figure 5 It is an iterative process consisting of four main steps as shown in .

[0104] Generate a random orthogonalized initial matrix

[0105] Using formula 10, the initial matrix Perform orthogonalization iterations;

[0106] Use formula 11 to decorrelate the previous step matrix to obtain the new process matrix W n+1 ;

[0107] Check the initial matrix And the new process matrix W after iteration n+1 Are they mutually replaced? If yes, the iteration process ends and the new process matrix W n+1 Transpose to get the weight matrix; if not, continue the orthogonalization iteration:

[0108] The ultimate goal of this process is to make W T The non-Gaussian property of x is maximized, that is, to maximize the G value in formula 6. First, generate a random orthogonalized initial matrix Then iterate it according to formula 10:

[0109]

[0110] in, is the initial matrix after iteration.

[0111] At the same time, in order to prevent different vectors in the weight matrix from converging to the same maximum value, they need to be decorrelated after each iteration according to formula 11:

[0112]

[0113] Among them, W n+1 is the new process matrix after decorrelation.

[0114] If the iterative process is convergent, the new process matrix W obtained by the iterative algorithm is n+1 With the original initial matrix They should be mutually replaced, that is, the dot product of the two should be 1. If this condition is met, the loop is exited and the weight matrix is ​​obtained, which is regarded as the target signal with the maximum non-quotient. Then, the target signal x is obtained according to formula 9. m×n The sparse pulse reflection s in (θ, t) j×n (t).

[0115] Finally, a technical problem encountered in the algorithm and its solution are explained: underdetermination and overdetermination in the application of the ICA algorithm. For the ICA algorithm, a preliminary estimate of the number of independent components is required before the algorithm starts. Generally speaking, this number is equal to the number of channels in the target signal X. However, in this study, underdetermination will cause the algorithm to be unable to calculate all the sparse pulse reflections in the objective function, and overdetermination will lead to repeated sparse pulse reflections in the final calculation results, causing changes in the characteristic information corresponding to the sparse pulse reflections in the weight matrix, thereby affecting the final classification results. Taking the same four-layer horizontal lithology medium model as an example, when the number of independent components is correctly selected as 3, the algorithm can accurately classify the sparse pulse reflections ( Figure 6 (b) from the noisy seismic prestack angle gather ( Figure 6 (a)), but if the number of initial selections is small and leads to underdetermination (for example, only 2 are selected), then the algorithm can only convert the noisy prestack angle gathers ( Figure 7 (a) shows the first two sparse pulse reflections with larger G values ​​( Figure 7 (b)). When the predicted number is greater than the actual number, two situations may occur. One of them is that there are repeated sparse pulse reflections in the calculation result (such as Figure 8 (c)), Figure 8 (d) is the corresponding weight matrix. In one case, the calculated result contains noise channels ( Figure 8 (a) Figure 8 (b) is the corresponding weight matrix.

[0116] Therefore, there needs to be a workflow in the algorithm to ensure that the number of independent components input to the algorithm is as close to the actual number as possible. Fig. 9 The processing flow is explained: the number of predictions starts from 1 at the beginning of each run, and increases by 1 at the beginning of each iteration. The result after the increase is checked to see if there is exactly the same sparse pulse reflection signal in the sparse pulse reflection matrix. If the overdetermined feature does not appear, the iteration is returned until the calculated weight result shows an overdetermined feature.

[0117] Fig.10 This paper explains the example of applying the ICA algorithm to a seismic prestack angle gather based on a three-layer horizontal lithologic medium model. The ICA algorithm can accurately extract sparse pulse reflections from the seismic prestack angle gather and calculate the angle-varying information matrix corresponding to each pulse signal. Fig.10 (a) is the input data of the ICA algorithm. Fig.10 (b) and (c) are the calculated sparse pulse reflections and the weight matrix corresponding to each signal, respectively.

[0118] Step S400: Compare the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input results.

[0119] The algorithm summary before this study uses the blind source separation algorithm to decompose the target signal into two matrices, one of which contains the location information of the sparse pulse reflection, and the other contains the corresponding information of the amplitude variation with angle (AVO), which we will call the weight matrix. The following program flow will compare the weight matrix with the preset rock physics template to calculate the formation geological model corresponding to the target pseudo-well function (that is, each formation reflection signal obtained). This problem is also divided into two parts: lithology determination based on the similarity law of AVO curves and the establishment of a multi-layer lithology medium model.

[0120] Lithology determination based on the similarity law of AVO curves:

[0121] First, according to all possible lithologies and their corresponding lithology elastic parameter data provided in Table 1, the seismic reflection coefficient sequence (i.e., seismic prestack angle gather) with seismic incident angles varying from 0° to 40° can be calculated by combining them in pairs and according to the Zoeppritz formula (Formula 3), and the preset AVO change curve of the corresponding lithology combination can be obtained by intersecting the incident angle and the reflection coefficient intensity. The AVO curves of all lithology combinations appearing in Table 1 are integrated to form a preset AVO change curve sample data dictionary library (i.e., preset rock physics template) for target comparison.

[0122] Secondly, by comparing the similarity between the AVO curve in the weight matrix and the preset AVO change curve in the sample data dictionary library, the corresponding lithology combination information is obtained through the preset curve with the strongest similarity (largest similarity index).

[0123] like Fig.11 As shown in the figure, this study uses the minimum square distance to calculate the similarity of two curves: For each preset AVO change curve, the similarity index I is defined, and this study uses formula 12 to calculate the similarity between the relevant element in the weight matrix and it. After calculating the similarity between a relevant element and all preset AVO change curves, they are compared and the possible lithology combinations are arranged. Fig.12 The similarity calculation between the weight matrix results of a noise-free four-layer horizontal lithologic medium model (each row in the weight matrix represents an AVO curve) and all preset AVO variation curves is explained. Fig.12 (a) is the noise-free seismic pre-stack angle gather; Fig.12 (b) is the calculated sparse pulse reflection; Fig.12(c) is the corresponding weight matrix. Table 2 shows the correlation factors of all possible lithology combinations.

[0124] Table 2 Correlation factors for all possible lithology combinations

[0125] Reflection 1 Reflection 2 Reflection 3 <![CDATA[L1&L2]]> 0.0808 0.1151 0.0528 <![CDATA[L1&L3]]> 0.1405 0.1289 0.0932 <![CDATA[L1&L4]]> 0.1417 0.1271 0.1246 <![CDATA[L1&L5]]> 0.0157 0.0810 0 <![CDATA[L2&L3]]> 0.1119 0.0947 0.1363 <![CDATA[L2&L4]]> 0.1478 0.1251 0.1177 <![CDATA[L2&L5]]> 0.1206 0.0866 0.1315 <![CDATA[L3&L4]]> 0 0 0.1091 <![CDATA[L3&L5]]> 0.1049 0.1084 0.1301 <![CDATA[L4&L5]]> 0.1362 0.1330 0.1047

[0126]

[0127] Among them, I is the similarity index; is the vertical distance between the AVO curve and all preset AVO change curve sample points in the weight matrix; i is the i-th one.

[0128] Since the potential model is a multi-layer lithologic medium model, it is necessary to judge and eliminate redundant stratum lithologic information before obtaining the final lithologic medium model.

[0129] like Fig.13 As shown in the figure, the similarity index between the AVO curve of the n-layer lithology interface in the weight matrix and all the preset AVO change curves is determined.

[0130] Arrange all convergent lithology combinations of n-layer lithology interfaces;

[0131] Determine the sum of similarity indices of all convergent lithology combinations of n-layer lithology interfaces, and select the maximum similarity index and the corresponding lithology combination as the output result;

[0132] The redundant information of the upper and lower layers of the lithology interface in the output results is merged to obtain a multi-layer lithology medium model, that is, a vertical lithology sequence.

[0133] First, for a multi-layer lithologic medium model with n+1 layers, there will be a total of n stratum interface information, corresponding to n AVO curves generated by different lithologic combinations. Through the process in step 1, the similarity of the AVO curves of these n interface information and all the preset AVO curves is calculated respectively.

[0134] Secondly, for a multi-layer model, the lower segment of the lithology combination of the previous model should be the upper segment of the lithology combination of the next model. If the similarity of the lithology combination of each layer is maximized, the lithology combination sequence in the vertical direction may not converge. Therefore, firstly, all the converged lithology combinations of the n-layer lithology medium model are arranged, and the sum of their similarity indexes is calculated.

[0135] Finally, the sequence with converged lithology combination and the largest similarity index is taken as the output result, and the final lithology combination is obtained by merging the redundant information of the upper and lower layers.

[0136] Based on the above-mentioned lithology identification method based on independent components and amplitude with offset, this embodiment proposes a lithology identification system based on independent components and amplitude with offset, including:

[0137] A building unit is used to build a pseudo-well model containing noise according to rock physical data;

[0138] A determination unit, configured to determine a seismic reflection coefficient sequence based on the pseudo-well model containing noise;

[0139] A weight matrix unit for determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence;

[0140] The comparison unit is used to compare the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input results.

[0141] The establishing unit is specifically used for:

[0142] Select rocks of different lithologies and measure rock physics data;

[0143] Set the number of lithology interfaces;

[0144] Assign rock physical data to the upper and lower layers of each lithology interface to obtain a pseudo-well model without noise;

[0145] Noise is added to the pseudo-well model without noise to obtain the pseudo-well model with noise.

[0146] Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent substitutions for some of the technical features therein; and these modifications or substitutions do not make the essence of the corresponding technical solutions deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A lithology identification method based on independent components and amplitude versus offset, characterized in that: include: A pseudo-well model with noise is established based on petrophysical data; determining a seismic reflection coefficient sequence based on the pseudo-well model containing noise; determining a sparse pulse reflection and a weight matrix from the seismic reflection coefficient sequence; The weight matrix is ​​compared with a preset rock physics template to establish a vertical lithology sequence of the input results.

2. The lithology identification method based on independent components and amplitude with offset according to claim 1 is characterized in that: The rock physical data include porosity, P-wave velocity, S-wave velocity and density.

3. The lithology identification method based on independent components and amplitude with offset according to claim 1 is characterized in that: The method of establishing a pseudo well model containing noise according to petrophysical data comprises the following steps: Select rocks of different lithologies and measure rock physics data; Set the number of lithology interfaces; Assign rock physical data to the upper and lower layers of each lithology interface to obtain a pseudo-well model without noise; Noise is added to the pseudo-well model without noise to obtain the pseudo-well model with noise.

4. The lithology identification method based on independent components and amplitude with offset according to claim 3 is characterized in that: The lithology includes normal porosity water-bearing sandstone, high porosity water-bearing sandstone, tight gas-bearing sandstone, normal porosity water-bearing limestone and normal porosity water-bearing dolomite.

5. The lithology identification method based on independent components and amplitude with offset according to claim 1, characterized in that: The specific steps of determining the seismic reflection coefficient sequence based on the pseudo well model containing noise are as follows: Based on the pseudo-well model containing noise, a seismic reflection coefficient sequence with pseudo-seismic incident angles ranging from 0° to 40° is determined.

6. The lithology identification method based on independent components and amplitude with offset according to claim 1, characterized in that: Determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence comprises the following steps: Perform centralization and whitening preprocessing on the seismic reflection coefficient sequence to generate a random orthogonalized initial matrix For the initial matrix Perform orthogonalization iterations; For the initial matrix Perform decorrelation processing to obtain the new process matrix W n+1 ; Determine the initial matrix and the new process matrix W n+1 Are they mutually replaced? If yes, the iteration process ends and the new process matrix W n+1 Transpose to get the weight matrix; if not, continue orthogonalization iteration; Sparse pulse reflections are determined based on a weight matrix and an input noisy signal.

7. The lithology identification method based on independent components and amplitude with offset according to claim 6 is characterized in that: The sparse pulse reflection is obtained by the following formula: X=T -1 *W T *s e +μ Where X is the input noisy target signal; T -1 is the whitening preprocessing operator; W T is the weight matrix; s e is the sparse pulse reflection; μ is the central preprocessing residual.

8. The lithology identification method based on independent components and amplitude with offset according to claim 1 is characterized in that: The specific steps of comparing the weight matrix with the preset rock physics template to establish the vertical lithology sequence of the input result are as follows: Different lithologies are combined in pairs to determine the seismic reflection coefficient sequence for seismic incident angles ranging from 0 to 40°; By intersecting the incident angle and the seismic reflection coefficient intensity, the preset AVO change curve corresponding to the lithology combination is obtained, that is, the preset rock physics template; Determine the similarity index between the AVO curve of the n-layer lithology interface in the weight matrix and all the preset AVO change curves; Arrange all convergent lithology combinations of n-layer lithology interfaces; Determine the sum of similarity indices of all convergent lithology combinations of n-layer lithology interfaces, and select the maximum similarity index and the corresponding lithology combination as the output result; The redundant information of the upper and lower layers of the lithology interface in the output results is merged to obtain a multi-layer lithology medium model, that is, a vertical lithology sequence.

9. A lithology identification system based on independent components and amplitude with offset, characterized in that: include: A building unit is used to build a pseudo-well model containing noise according to rock physical data; A determination unit, configured to determine a seismic reflection coefficient sequence based on the pseudo-well model containing noise; A weight matrix unit for determining sparse pulse reflections and a weight matrix from the seismic reflection coefficient sequence; The comparison unit is used to compare the weight matrix with a preset rock physics template to establish a vertical lithology sequence of the input results.

10. The lithology identification system based on independent components and amplitude with offset according to claim 9, characterized in that: The establishing unit is specifically used for: Select rocks of different lithologies and measure rock physical data; Set the number of lithology interfaces; Assign rock physical data to the upper and lower layers of each lithology interface to obtain a pseudo-well model without noise; Noise is added to the pseudo-well model without noise to obtain the pseudo-well model with noise.

Images