A thin interbed longitudinal wave impedance prediction method based on deep learning

Abstract

The application discloses a thin interbed longitudinal wave impedance prediction method based on deep learning, and the method comprises the following steps: collecting and processing artificial data to obtain seismic data; well seismic calibration and depth-time conversion are performed to obtain the longitudinal wave impedance of the thin interbed section at the well point and the seismic data; time-frequency analysis and sliding time window are performed to obtain the time-frequency spectrum image of the local time window of the seismic trace beside the well, which is used as the input sample of the deep learning network; the longitudinal wave impedance value of the thin interbed in the local time window is obtained as the sample label; the training of the thin interbed longitudinal wave impedance prediction network model is completed to obtain the best network parameters; the best network parameters are used to traverse all seismic traces to obtain the thin interbed longitudinal wave impedance prediction results in all time windows; and the results in the overlapping time windows are processed to determine the final spatial distribution of the longitudinal wave impedance. Compared with the existing method, the method combines the logging and seismic data, can break through the seismic resolution, does not need the wavelet information, has high prediction accuracy, and is convenient to popularize.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of reservoir prediction in artificial seismic exploration, and further relates to the field of thin interbedded reservoir prediction, in particular to a thin interbedded P-wave impedance prediction method based on deep learning. BACKGROUND

[0002] The artificial seismic exploration method refers to artificially exciting seismic waves on the ground, and when the seismic waves propagate to the underground, they will be reflected and transmitted when encountering the interface of rock layers with different medium properties. The seismic waves are received by geophones on the ground to obtain seismic records. The characteristics on the seismic records are related to the properties and structures of the underground rock layers. Through processing and interpretation of the seismic records, the properties and shapes of the underground rock layers can be inferred. Artificial seismic exploration is superior to other geophysical exploration methods in terms of accuracy and precision, and thus plays an important role in oil and gas reservoir prediction.

[0003] In recent decades, with the continuous deepening of oil and gas exploration, large-scale and easily found structural traps have been almost exhausted, and the targets of oil and gas exploration have gradually shifted to small-thickness and fast-lithology-changing lithology and complex traps. Thin interbedded reservoirs are widely developed in oil and gas basins in China, and large-area thin interbedded lithologic oil and gas reservoirs are the focus of increasing reserves and production in the current conventional oil and gas exploration field. The single-layer thickness of thin interbedded strata is small, the vertical interbedding is frequent, the lateral lithology and thickness change rapidly, and the accurate prediction of P-wave impedance is difficult. Therefore, it is a very meaningful work to carry out the prediction of P-wave impedance of thin interbedded strata in oil and gas exploration.

[0004] In actual artificial seismic exploration, a stratum with a single-layer thickness less than 1 / 4 of the main wavelength of the seismic wave is considered as a thin layer, and multiple thin layers vertically stacked form a thin interbedded layer. Due to the limitation of the longitudinal resolution of seismic data, it is difficult to accurately predict the P-wave impedance of thin interbedded layers. The P-wave impedance prediction using seismic data usually uses seismic inversion method, which can integrate well data to effectively improve the resolution of the results and provide helpful assistance for the identification of thin interbedded layers. However, the inversion method requires seismic wavelet information, and the accuracy of the wavelet directly affects the accuracy of the inversion. The present application provides a thin interbedded P-wave impedance prediction method based on deep learning, which uses the two-dimensional time-frequency spectrum of seismic data as input, automatically extracts the features related to the P-wave impedance of thin interbedded layers from the time-frequency spectrum with the help of deep learning, establishes the mapping relationship between the local time-frequency spectrum and the P-wave impedance, and realizes the accurate prediction of the P-wave impedance of thin interbedded layers. Compared with the existing methods, this method combines logging and seismic data, breaks through the resolution of traditional post-stack seismic inversion method, and does not require seismic wavelet information, with high prediction accuracy and easy to use. SUMMARY

[0005] The purpose of the present application is to provide a thin interbed longitudinal wave impedance prediction method based on deep learning, so as to improve the identification and prediction ability of thin interbeds.

[0006] To achieve the above-mentioned purpose, the present application provides a thin interbed longitudinal wave impedance prediction method based on deep learning, and the content of the present application includes:

[0007] Field acquisition of artificial seismic data, indoor processing to obtain post-stack seismic data;

[0008] Well seismic calibration and depth-time conversion are carried out to obtain thin interbed longitudinal wave impedance data and seismic data at well points;

[0009] Time-frequency analysis is performed on the seismic data of the target thin interbed section to obtain a two-dimensional time-frequency spectrum of the seismic trace;

[0010] The two-dimensional time-frequency spectrum image in the local time window is obtained from the well seismic trace, and after preprocessing, it is used as the input sample of the deep learning network;

[0011] The continuous value of the thin interbed longitudinal wave impedance in the local time window is obtained from the logging data, and after normalization, it is used as the label of the sample;

[0012] The deep learning network architecture is built, the loss function is defined, the network model evaluation method is selected, the thin interbed longitudinal wave impedance prediction network model training is completed, and the best network parameter is obtained;

[0013] All two-dimensional time-frequency spectrum images in the local time window are obtained by traversing all seismic traces;

[0014] The preprocessed two-dimensional time-frequency spectrum image set is converted into thin interbed longitudinal wave impedance prediction results by using the best network parameter;

[0015] The longitudinal wave impedance prediction results in the overlapping time window are processed to determine the final spatial distribution of the thin interbed longitudinal wave impedance.

[0016] The method of the present application has the advantages that, compared with the existing methods, the method combines seismic and well data, breaks through the resolution of traditional post-stack seismic inversion methods, and does not require seismic wavelet information, so that the absolute longitudinal wave impedance signal is obtained, the prediction accuracy is high, and the method is easy to use and promote. BRIEF DESCRIPTION OF DRAWINGS

[0017] In order to more clearly illustrate the technical solutions in the embodiments of the present application or the prior art, the following will briefly introduce the drawings needed to be used in the embodiments or prior art description. Obviously, the drawings in the following description are only some embodiments of the present application, and for those skilled in the art, other drawings can also be obtained without creative labor.

[0018] Figure 1 A flowchart of a deep learning-based thin interbed P-wave impedance inversion method of an embodiment of the present application;

[0019] Figure 2 Seismic data containing a thin interbed section of an embodiment of the present application;

[0020] Figure 3 A seismic trace and its two-dimensional time-frequency spectrum image of an embodiment of the present application;

[0021] Figure 4 A sample and label making method of an embodiment of the present application;

[0022] Figure 5 A deep learning network architecture of an embodiment of the present application;

[0023] Figure 6 A data set division proportion diagram of an embodiment of the present application;

[0024] Figure 7 A network hyperparameter table of an embodiment of the present application;

[0025] Figure 8 A loss curve of an embodiment of the present application;

[0026] Figure 9 A thin interbed P-wave impedance inversion result obtained by using commercial software of the present embodiment;

[0027] Figure 10 A deep learning-based thin interbed P-wave impedance prediction result of an embodiment of the present application. DETAILED DESCRIPTION

[0028] The technical solutions in the embodiments of the present application will be described clearly and completely below with reference to the accompanying drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those of ordinary skill in the art without creative work fall within the protection scope of the present application.

[0029] Figure 1 A flowchart of a deep learning-based thin interbed P-wave impedance prediction method of an embodiment of the present application. As shown in Figure 1 , the method comprises the following steps:

[0030] S101, artificial earthquake data acquisition on the ground, indoor seismic data processing, obtain poststack seismic data. Specifically, in step S101, artificial earthquake refers to the purpose of exploration, artificially excited earthquake. In the field of artificial earthquake acquisition, it can be three-dimensional acquisition along a surface, or two-dimensional acquisition along a line, and the obtained poststack seismic data can be three-dimensional data or two-dimensional data. Select a three-dimensional seismic data of a two-dimensional line of a certain area in the underground exploration target, and carry out the application display of the embodiment of the application. The exploration target in the area includes a set of thin interbedded section, which is delta front sediment, sand and mudstone lithology, and develops multiple layers of sand and mudstone interbedded structure from bottom to top, the spatial relationship of reservoir sand body is complex, the lateral variation is fast, and the heterogeneity is strong. The seismic data of the target thin interbedded section is intercepted, Figure 2 For the seismic data of the embodiment of the application containing the thin interbedded section, the upper arrow indicates the position of the well, that is, there are five drilling wells in the area, and the well data is complete.

[0031] S102, well-to-seismic calibration and depth-time conversion are carried out on the target layer section by using logging data and well seismic trace, and longitudinal wave impedance data and seismic data of the thin interbedded section at the well point are obtained; specifically, in step S102, for the embodiment of the application, the logging curves used are acoustic travel time curve and density curve. In the present example, well-to-seismic calibration adopts the method of making artificial synthetic seismogram, first, the acoustic travel time curve is used to obtain the acoustic velocity curve by inverse operation, the acoustic velocity curve is multiplied by the density curve to obtain the longitudinal wave impedance curve, the longitudinal wave velocity curve is subjected to point-by-point reflection coefficient operation to obtain the discrete longitudinal wave reflection coefficient sequence, the wavelet matched with the main frequency of the seismic data and the reflection coefficient sequence are subjected to convolution operation to obtain the synthetic seismogram, the over-well or well-side seismic trace is extracted from the actual seismic data, the well-side seismic trace is matched and calibrated with the synthetic seismogram, and the time-depth relationship at the well is determined. In well-to-seismic calibration, the matching degree of the target thin interbedded section is focused on. The same calibration operation is carried out on the wells in the work area. Using the obtained time-depth relationship, the longitudinal wave impedance data calculated from each well are converted into longitudinal wave impedance data with time coordinates, the converted longitudinal wave impedance data are consistent with the time coordinates of the seismic data, and have the same time sampling interval. In the present example, the time sampling interval of the seismic data is 1ms.

[0032] S103, time-frequency analysis is performed on the seismic data of the target thin interbedded section, and one-dimensional waveform data of the seismic trace is converted into two-dimensional time-frequency spectrum of the seismic trace. The time-frequency analysis in step S103 can be linear time-frequency analysis, such as short-time Fourier transform, continuous wavelet transform and S transform, or nonlinear time-frequency analysis, such as Wigner-Ville distribution, pseudo Wigner-Ville distribution and Cohen class time-frequency distribution. In the present example, wavelet transform time-frequency analysis method is adopted.

[0033] In step S103, the frequency range of the time-frequency analysis is based on the dominant frequency f of the target thin interbedded section seismic data. m and its effective frequency band range [f1, f h [Confirmed, low frequencies start from 0Hz, high frequencies start from 2Hz.] m +(f h -f i ) / 2, frequency domain sampling interval 1 Hz. In this example, the dominant frequency f of the seismic data m The frequency range is approximately 30 Hz, with an effective frequency band ranging from 8 Hz to 88 Hz. The frequency range for time-frequency analysis is 0 Hz to 100 Hz, and the frequency domain sampling interval is 1 Hz. Figure 3 This paper presents a seismic trace and its two-dimensional time-spectrum image according to an embodiment of the present invention. The left side shows the waveform of the seismic trace, and the right side shows the two-dimensional time-spectrum image obtained by wavelet transform of the seismic trace.

[0034] S104. Obtain a two-dimensional time-spectral image within a local time window from the seismic trace near the well using a sliding time window, and preprocess it as input samples for the deep learning network. Specifically, in step S104, the sliding time window's length is based on the dominant frequency f of the target thin interbedded section seismic data. m It is determined that the minimum value is 1 / f. m Twice as much, with a maximum of 1 / f m The time window's main frequency is four times the time sampling interval and is an integer multiple of the time sampling interval; the sliding step size of the sliding window is at least one time sampling interval and at most one time window length. In this example, the seismic data's dominant frequency f m The frequency is approximately 30 Hz, the length of the sliding window is determined to be 100 ms, and the sliding step size adopts the minimum step size, that is, one time sampling interval of 1 ms.

[0035] The preprocessing in step S104 specifically involves normalizing and resampling the two-dimensional time-spectrum image within the local time window. Normalization is necessary, while resampling is selective, depending on the required recognition accuracy and computational efficiency. Resampling can be performed in both the vertical and horizontal directions of the two-dimensional time-spectrum image within the local time window, or only in one direction. In this example, normalization and resampling are performed during preprocessing. Normalization, based on the pixel mean and standard deviation, maps pixels to the range [-1, 1]. After normalization, the data distribution is centered at 0, which is beneficial for gradient descent optimization. The formula is:

[0036]

[0037] Where x is the pixel value of the image; μ is the mean used for normalization; σ is the standard deviation used for normalization; x norm It is the normalized pixel value.

[0038] In the embodiment, the resampling processing utilizes a linear interpolation manner to improve the resolution of the time-frequency spectrum, and the specific operation is to perform three rounds of linear interpolation on each row of time-frequency data. Assuming that the original data is x = [x1, x2,..., xN], after one round of interpolation, the interpolated data is: n

[0039]

[0040] Repeating the above step three times, the final interpolated data can be obtained.

[0041] S105, using the same sliding time window as in step S104, obtaining the continuous value of the P-wave impedance in the local time window from the well logging data, and performing minimum-maximum normalization processing to obtain the label of the sample. The formula of the minimum-maximum normalization processing is:

[0042]

[0043] wherein I is the original P-wave impedance value; I max is the global maximum value of the P-wave impedance, that is, the maximum value of all P-wave impedance values in the sample set; I min is the global minimum value of the P-wave impedance, that is, the minimum value of all P-wave impedance values in the sample set; and I norm is the normalized P-wave impedance value, which ranges between [0, 1]. The normalized P-wave impedance label is a one-dimensional vector with a numerical value distributed in the interval [0, 1], and the length of the vector is the same as the time sampling point number of the sample. In step S105, the continuous value of the P-wave impedance in the local time window is the P-wave impedance data with the same time sampling interval and sampling number as the seismic data in the local time window. The sliding time window is the same as in step S104, and has the same time window length and the same sliding step length. Figure 4 The sample and label making manner of the embodiment of the application is shown. In the example, the time window length and the sliding step length are 100 ms and 1 ms respectively, I max and I min are 6600 m / s*g / cm 3 and 6380 m / s*g / cm 3 respectively.

[0044] In step S105, the sample label is a one-dimensional vector after normalization processing, and the length of the vector is the same as the time sampling point number of the sample. In the example, the label of each sample is a one-dimensional vector of continuous numerical values normalized to 0-1, and the vector length is 100. In the example, a single label indicates the P-wave impedance value at each time sampling point position in the corresponding 100 ms local time window, and with the sliding of the time window, the numerical value change of the adjacent label indicates the vertical change of the P-wave impedance.​

[0045] S106, match the sample and the label, divide the training set and the verification set, build the deep learning network architecture, define the loss function, select the network model evaluation method, test the network hyperparameters, complete the thin interlayer longitudinal wave impedance prediction network model training, and obtain the best network parameters. The deep learning network refers to a deep network with three or more feature extraction layers, and the feature extraction layer adopts two-dimensional convolution operation. In this example, a deep network with a U-shaped structure is used, Figure 5 The deep learning network architecture diagram of the embodiment of the application is shown in FIG. 6. The network is an encoder-decoder structure, the encoder is composed of four progressive down-sampling modules, each module contains two consecutive 3x3 convolution operations and one 2x2 max pooling operation, the decoder is composed of four progressive up-sampling modules, each module contains one 2x2 deconvolution and two 3x3 convolution operations, and the output of each layer in the encoder is spliced with the up-sampling result of the corresponding layer in the decoder by using a jump connection to realize the fusion of image features of different scales.

[0046] Figure 6 The data set division ratio diagram of the embodiment of the application is shown in FIG. 7. In this example, there are a total of 5 wells, the total number of labeled samples is 3840, the samples of 4 wells are used as the training set, and the samples of the remaining well are used as the verification set. The training set and the verification set account for 80% and 20% respectively. Figure 2

[0047] In step S106, the loss function adopts a mean square error function (MSE):

[0048]

[0049] Wherein, N is the number of samples; Y represents the label data; y represents the predicted data.

[0050] In step S106, the network model evaluation method adopts a loss curve evaluation method, and the loss curve is generated by the loss function. The loss curve is generated by the value of the loss function with the change of the training iteration number, and reflects the optimization process of the model. Figure 7 The network hyperparameter table of the embodiment of the application is shown in FIG. 8. In this example, the Adam optimization method is used, and the ideal batch size obtained by testing is 1, the epoch is 150, and the learning rate is 0.001. Figure 8 The loss curve of the embodiment of the application is shown in FIG. 9. The loss curves of the training set and the verification set can converge normally, which indicates that the network training achieves the expected effect, the thin interlayer longitudinal wave impedance prediction network model training is completed, and the best network parameters are obtained.

[0051] ​S107, traverse the 2D time-frequency spectrum of all seismic traces with the same length of sliding time window as in step S104 to obtain the whole set of 2D time-frequency spectrum images in local time window, and pre-process all images in the set. In step S107, the length of sliding time window is the same as that in step S104; the sliding step can be the same as that in step S104 or different. In the present example, the same sliding time window as in step S104 is used, and the length of time window and the sliding step are 100 ms and 1 ms respectively. In step S107, the same pre-processing method as in step S104 is used, and in the present example, normalization and resampling pre-processing are performed.

[0052] S108, convert the whole set of pre-processed 2D time-frequency spectrum images in local time window into the prediction results of thin interbed P-wave impedance in local time window with the optimal network parameters. Specifically, in step S108, the optimal network is the deep network with optimal weights obtained in step S106, and the conversion process refers to inputting each image in the set of pre-processed 2D time-frequency spectrum images in local time window into the optimal network, and the output is the prediction result of thin interbed P-wave impedance in the corresponding local time window.

[0053] S109, perform averaging and de-normalization processing on the prediction results of P-wave impedance in overlapping time window according to the vertical position of the center point of sliding time window of each seismic trace to determine the final spatial distribution of the target thin interbed P-wave impedance. Specifically, in step S109, the sliding time window is the same sliding time window as in step S107, and has the same length of time window and sliding step; the overlapping time window refers to the overlapping part of time window when the sliding step is less than the length of time window, and the averaging method of the prediction results of P-wave impedance in overlapping time window is:

[0054]

[0055] wherein M is the number of overlapping time windows, which is related to the length of sliding time window and the sliding step; N is the number of overlapping sampling points, which is related to the length of seismic trace and the length of sliding time window; I ij is the prediction value of P-wave impedance of the i-th sampling point in the j-th overlapping time window; I i is the average prediction value of P-wave impedance of the i-th sampling point, which is in the normalized state and has a value range of [0, 1]. The existence of overlapping time window realizes multiple prediction of P-wave impedance at the same position in the subsurface, and the comprehensive use of multiple prediction results of P-wave impedance determines the final P-wave impedance value at the position, avoiding the uncertainty of single prediction result and improving the prediction accuracy and robustness. In the present example, M = 101 and N = 450.

[0056] Further, in step S109, Ii The inverse normalization processing is performed to restore the actual P-wave impedance range, and the final P-wave impedance prediction value is obtained, and the inverse normalization method is:

[0057] I f = I i ×(I max -I min )+I min

[0058] wherein, I f is the final P-wave impedance value after inverse normalization; I i is the normalized average P-wave impedance prediction value of the i-th sampling point; I max and I min are the global maximum value and the global minimum value of the P-wave impedance, which are consistent with the normalization in step S105, and in the present example, I max and I mm are 6600 m / s*g / cm 3 and 6380 m / s*g / cm 3 , respectively. In step S109, the final spatial distribution of the P-wave impedance of the thin interbed is the P-wave impedance prediction result consistent with the P-wave impedance data range at the well point after the P-wave impedance prediction values of the sampling points in the space are averaged and inverse normalized. In the present example, it means that the one-dimensional vector with the value range between 0 and 1 predicted by the network is converted into the actual P-wave impedance value, and the range is between 6600 m / s*g / cm 3 and 6380 m / s*g / cm 3 . In step S109, the obtained P-wave impedance spatial distribution has the same dimension as the post-stack seismic data in step S101, that is, when the post-stack seismic data is three-dimensional data, the final P-wave impedance spatial distribution of the thin interbed obtained is three-dimensional; when the post-stack seismic data is two-dimensional data, the final P-wave impedance spatial distribution of the thin interbed obtained is two-dimensional. In the present example, Figure 10 is the thin interbed P-wave impedance inversion result based on deep learning of the embodiment of the present application, and shows a well line. Figure 9 is the P-wave impedance inversion result obtained by the traditional method using commercial software, and Figure 9 and Figure 10 It can be seen that the method proposed in the present application has obvious advantages for the identification and prediction of the P-wave impedance value of the thin interbed, the thin interbed sand body structure is clear and separable, the vertical superposition relationship and the horizontal extinction point of the sand body are clearly revealed, and the P-wave impedance value prediction accuracy is high.

[0059] Although the present application has been described in detail with general description and specific embodiments above, it is obvious to those skilled in the art that some modifications or improvements can be made on the basis of the present application. Therefore, these modifications or improvements made on the basis of not deviating from the spirit of the present application, all belong to the scope of protection claimed by the present application.

Claims

1. A method for predicting the impedance of thin interlayered longitudinal waves based on deep learning, characterized in that... Includes the following steps: (1) Artificial seismic data acquisition and processing are carried out on the Earth's surface to obtain post-stack seismic data; (2) Use well logging data and wellside seismic traces to perform well-seismic calibration and depth-time conversion on the target section to obtain P-wave impedance data and seismic data of thin interbedded sections at the well point; (3) Perform time-frequency analysis on the seismic data of the target thin interbedded section, converting the one-dimensional waveform data of the seismic trace into a two-dimensional time-frequency spectrum of the seismic trace. The frequency range of the time-frequency analysis is based on the dominant frequency f of the seismic data of the target thin interbedded section. m and its effective frequency band range [f1, f h [Confirmed, low frequencies start from 0Hz, high frequencies start from 2Hz.] m +(f h -f1) / 2 or higher, with a frequency domain sampling interval of no more than 1 Hz; (4) A two-dimensional time-spectrum image within a local time window is obtained from the seismic trace near the well using a sliding time window, and the image is preprocessed and used as input samples for a deep learning network. The time window length of the sliding time window is based on the dominant frequency f of the seismic data of the target thin interbedded section. m It is determined that the minimum value is 1 / f. m Twice as much, with a maximum of 1 / f m 4 times; the minimum sliding step size is one time sampling interval, and the maximum is one time window length; (5) Using the same sliding time window as in step (4), obtain the P-wave impedance value within the local time window from the logging data, and perform minimum-maximum normalization on it as the label of the sample. The formula for minimum-maximum normalization is: Where I is the original longitudinal wave impedance value; I max This represents the global maximum value of the P-wave impedance, i.e., the maximum value of all P-wave impedance values ​​in the sample set; I min This represents the global minimum of the P-wave impedance, i.e., the minimum value of all P-wave impedance values ​​in the sample set; I norm The normalized longitudinal wave impedance value ranges between [0, 1]. The normalized longitudinal wave impedance label is a one-dimensional vector with values ​​distributed in the interval [0, 1], and the length of the vector is the same as the number of time sampling points of the sample. (6) Match samples and labels, divide training and validation sets, select loss function and network model evaluation method, build deep learning network architecture, test network hyperparameters, complete the training of thin interlayer longitudinal wave impedance prediction network model, and obtain the best network parameters. (7) Use a sliding time window of the same length as in step (4) to traverse the two-dimensional time spectrum of all seismic traces, obtain the complete set of two-dimensional time spectrum images within the local time window, and preprocess all images in the set. (8) Using the optimal network parameters, the entire set of two-dimensional time-spectrum images within the preprocessed local time window is converted into the prediction result of the longitudinal wave impedance of thin interlayer within the local time window. (9) Based on the vertical position of the center point of each seismic trace sliding window, the P-wave impedance prediction results within the overlapping windows are averaged and denormalized to determine the final spatial distribution of the target thin interlayer P-wave impedance. The overlapping window refers to the overlapping portion of adjacent sliding windows when the sliding step size is less than the window length. The averaging method for the P-wave impedance prediction results within the overlapping window is as follows: Where M is the number of overlapping time windows, which is related to the length of the sliding time window and the sliding step size; N is the number of overlapping sampling points, which is related to the time length of the seismic trace and the length of the sliding time window; I ij I represents the predicted longitudinal wave impedance of the i-th sampling point in the j-th overlapping time window; i is the predicted average longitudinal wave impedance of the i-th sampling point, and is the value under normalized conditions, which ranges from [0, 1]. to I i After performing inverse normalization to restore the actual P-wave impedance range, the final predicted P-wave impedance value is obtained. The inverse normalization method is as follows: I f =I i ×(I max -I min )+I min Among them, I f This is the final longitudinal wave impedance value after inverse normalization; I i I is the normalized average predicted P-wave impedance at the i-th sampling point; max I min The global maximum and global minimum values ​​of the longitudinal wave impedance are consistent with those in step (5) during normalization.

2. The method according to claim 1, characterized in that the artificial earthquake in step (1) is an earthquake artificially triggered for the purpose of exploration; When conducting artificial seismic data acquisition in the field, three-dimensional acquisition can be performed along a surface or two-dimensional acquisition along a line. The resulting post-stack seismic data can be either three-dimensional or two-dimensional.

3. The method according to claim 1, characterized in that the logging data in step (2) includes sonic transit time curves and density curves; The P-wave impedance data of the thin interbedded section at the well point refers to the P-wave impedance curve obtained by dividing the logging density curve by the sonic transit time curve at the well point, and the time domain P-wave impedance data obtained after depth-time conversion is consistent with the time coordinate of the seismic data and has the same time sampling interval.

4. The method according to claim 1, wherein the time-frequency analysis in step (3) can be linear time-frequency analysis or nonlinear time-frequency analysis.

5. The method according to claim 1, characterized in that the sliding window is the same in steps (4) and (5), having the same window length and the same sliding step size; The preprocessing described in step (4) refers to normalizing and resampling the two-dimensional time spectrum image within the local time window. Resampling can be performed in both the vertical and horizontal directions of the two-dimensional time spectrum image within the local time window, or only in one direction.

6. The method according to claim 1, characterized in that the continuous P-wave impedance value of thin interlayered layers within the local time window in step (5) is P-wave impedance data with the same time sampling interval and number of samples as the seismic data within the local time window.

7. The method according to claim 1, characterized in that the deep learning network in step (6) refers to a deep network with three or more feature extraction layers, wherein the feature extraction layers are performed using two-dimensional convolution operations; The loss function mentioned above uses the mean square error function; The network model evaluation method described above uses a loss curve evaluation method, where the loss curve is generated by the loss function.

8. The method according to claim 1, characterized in that the length of the sliding window in step (7) is the same as the length of the sliding window in step (4), and the sliding step size is the same as or different from the sliding step size in step (4); the preprocessing adopts the same method as the preprocessing in step (4).

9. The method according to claim 1, characterized in that the optimal network in step (8) is a deep network with optimal weights obtained in step (6); The aforementioned transformation refers to inputting each image from the set of two-dimensional time-spectrum images within the preprocessed local time window into the optimal network, and its output is the prediction result of the longitudinal wave impedance of the thin interlayer within the corresponding local time window.

10. The method according to claim 1, characterized in that the existence of the overlapping time window in step (9) enables multiple predictions of P-wave impedance at the same underground location, and the results of multiple P-wave impedance predictions are combined to determine the final P-wave impedance value at that location, avoiding the uncertainty of a single prediction result and improving the accuracy and robustness of the prediction. The final spatial distribution of the longitudinal wave impedance of thin interlayer refers to the result obtained by averaging and inverse normalization of the final predicted values ​​of each sampling point in space according to step (9); The spatial distribution has the same dimension as the post-stack seismic data in step (1). That is, when the post-stack seismic data is three-dimensional, the final spatial distribution of the thin interlayered P-wave impedance is three-dimensional; when the post-stack seismic data is two-dimensional, the final spatial distribution of the thin interlayered P-wave impedance is two-dimensional.

Images