METHOD FOR RECONSTRUCTING AT LEAST ONE TRACE IN A SEISMIC IMAGE

20230066911 · 2023-03-02

    Inventors

    Cpc classification

    International classification

    Abstract

    The present invention is related to a method for reconstructing at least one trace in a seismic image of a common receiver and time domain, the image comprising traces in time domain with seismic data and one or more traces to be reconstructed. A first aspect of the invention is a method that is characterized by a specific use of a convolutional neural network trained under an unsupervised learning approach with a modified receptive field. A second aspect of the invention is a deblending method based on the use of a reconstructing method according to the first aspect of the invention applied to a denoising step of a deblending process allowing a very effective data acquisition while keeping a high quality output data sets after being processed according to the first and/or second aspects of the invention.

    Claims

    1. A method for reconstructing at least one trace in a seismic image, the method comprising the steps of: a) providing a seismic image having seismic traces extracted from a data set acquired via a seismic survey, wherein the seismic image includes at least one seismic trace with data that is unavailable from the acquired data set; b) training a convolutional neural network to predict values for the unavailable data of the at least one trace, the convolutional neural network includes at least one layer having a kernel function configured to evaluate a bounded domain that defines a blind-trace receptive field, wherein the blind-trace receptive field encompasses data corresponding to a plurality of traces in the seismic image that are adjacent to, but exclude, the at least one trace having unavailable data; c) predicting values for the unavailable data of the at least one trace by inputting the seismic image into the trained convolutional neural network; d) reconstructing the at least one trace having unavailable data using the data values predicted in step c); and e) generating a reconstructed seismic image with the at least one trace reconstructed in step d).

    2. The method according to claim 1, wherein the convolutional neural network is a U-net.

    3. The method according to claim 1, wherein the at least one layer is a down-sampling layer.

    4. The method according to claim 1, wherein the blind-trace receptive field of the kernel function is determined by the following processing steps: limiting the receptive field to one adjacent side of the at least one trace having unavailable data, the one adjacent side being determined based on a direction of the at least one trace having unavailable data; generating a first copy of the seismic image and a second copy of the seismic image, the second copy of the seismic image being rotated 180° with respect to the first copy of the seismic image; inputting the first copy and the second copy of the seismic image into the convolutional neural network with the one-adjacent-sided receptive field resulting in two output images that are subsequently combined to form a single image.

    5. The method according to claim 4, wherein the at least one layer with a blind-trace receptive field is a down-sampling layer and wherein any layer of the convolutional neural network includes an output for outputting a feature map and, wherein: before inputting a feature map, output of a previous layer or the image if the current layer is the first layer, into the current layer, the feature map is padded by adding rows of zeros at the end of the feature map located at one side of the trace; carrying out a convolution operation; cropping out the same number of rows previously added wherein cropping the output feature map is carried out on the side opposite to the side on which rows were previously added.

    6. The method according to claim 2, wherein training the convolutional neural network includes using a converging criterion based on an approximation error estimation E.sub.s for measuring an interpolation loss when predicting the reconstructed trace, and the approximation error estimation E.sub.s being determined as a linear combination of a misfit loss and a regularization loss, the regularization loss being determined by the following steps: determining a main energy area of the seismic image; calculating a norm of the reconstructed seismic image in the frequency domain limited to an area outside of the main energy area based on a predetermined norm.

    7. The method according to claim 6, wherein the linear combination determining the approximation error estimation is:
    E.sub.S=∥misfit loss∥+αƒregularization loss∥ being α a positive weighting value and ∥⋅∥ the predetermined norm.

    8. The method according to claim 7, wherein the misfit loss is determined as a difference between the seismic image and the reconstructed seismic image after removing the at least one reconstructed trace.

    9. A method for deblending seismic data in a receiver domain, the method comprising: deploying a plurality of n.sub.s acoustic sources in the upper surface of the reservoir domain wherein the n.sub.s acoustic sources are grouped in B groups σ.sub.i, l=1 . . . B, of acoustic sources, each acoustic source being in only one group σ.sub.i of sources and at a location x.sub.i.sup.s, i=1 . . . n.sub.s and, deploying a plurality of n.sub.r, acoustic receivers in the upper surface of the reservoir domain at a location x.sub.j.sup.r, j=1 . . . n.sub.r; a) for each group σ.sub.i, l=1 . . . B of acoustic sources, each acoustic source is shot with a random delay time τ.sub.li and the response in the acoustic receivers stored in a data structure that may be represented by P.sub.ljk.sup.b=(x.sub.i.sup.s, x.sub.i.sup.r, t.sub.k−τ.sub.li; wherein t.sub.k is the k.sup.th time sample in the time domain; b) calculating, for each group σ.sub.l, l=1 . . . B, the Fourier-transform Π.sup.b(:,:,ω.sub.k)=F{P.sub.ljk.sup.b} wherein ω.sub.k is the k.sup.th frequency and “:” denoting variables depending in index i or index j; c) for each frequency ω.sub.k determining Π.sup.LS(:,:,ω.sub.k)=Γ*Π.sup.b(:,:,ω.sub.k) wherein Γ* is
    Γ*=Γ.sub.k.sup.HD being D a diagonal matrix and Γ.sub.k.sup.H the conjugate transpose of Γ.sub.k the blending matrix that may be calculated from the random delay times τ.sub.li as ( Γ k ) li = { e - - 1 ω k τ l i for x i s ϵσ l 0 else d) calculating an Inverse Fourier Transform of F.sup.−1(Π.sup.LS(:,:,ω.sub.k))=P.sub.ljk; e) the shot-gather ordering in the output P.sub.ljk is sorted to get the image I.sub.j of the trace data in the common-receiver domain or the common-midpoint receiver or the common-offset receiver and, time domain; f) for each trace carrying out a deblending step by reconstructing a coherent signal of the traces using the reconstructing method according to claim 1.

    10. The method for deblending seismic data according to claim 9, wherein D=(Γ.sub.kΓ.sub.k.sup.H).sup.−1.

    11. The method for deblending seismic data according to claim 9, wherein D=I, wherein I denotes the identity matrix.

    12. The method for deblending seismic data according to claim 9, wherein each source is shot once and wherein for all m≠n, σ.sub.m∩σ.sub.n=Ø.

    13. A non-transitory computer program product stored on a computer-readable medium and comprising computer-implementable instructions, which, when executed by a computer, cause the computer to carry out the method according to claim 1.

    14. A computer system having a processor and a non-transitory computer-readable medium storing computer-executable instructions which, when executed by the processor, cause the processor to carry out the method according to claim 1.

    Description

    DESCRIPTION OF THE DRAWINGS

    [0091] These and other features and advantages of the invention will be seen more clearly from the following detailed description of a preferred embodiment provided only by way of illustrative and non-limiting example in reference to the attached drawings.

    [0092] FIG. 1 This figure shows an example of the failed result by the normal convolutional neural network, according to the prior art, which trained with the blending loss. (a) is the input noisy data and (b) is the prediction after the convergence.

    [0093] FIG. 2 This figure shows schematically the concept of a blind-trace method according to the invention. In order to recover the thick black line representing a trace in (b), only information in the dashed patches in (a) are used, ignoring the original noisy trace itself.

    [0094] FIG. 3 This figure shows schematically the structure of the blind-trace network according to an embodiment of the invention. The thick black line marks a sample target trace. The dashed patches (the first patch with inclined lines and the second patch with Cartesian grid) indicates the areas used for predicting the predicted result. In the first layer after the blind-trace U-net, the patches cover the target trace, and with padding and cropping operations, the thick black line is excluded and the model is fully blind-trace.

    [0095] FIG. 4 This figure shows schematically the modified convolutional layer and max-pooling layer in a blind-trace network according to an embodiment of the invention.

    [0096] FIGS. 5 and 5 (continuation)

    [0097] FIG. 5 ranges from (a) to (f) showing examples of blended-unblended pairs using three different blending schemes. Figures (a), (b) and (c) are pseudo-deblended common receiver gathers from the three blending schemes identified as (“alternate”, “half”, and “continuous”) respectively. Figures (d), (e) and (f) are the corresponding unblended data.

    [0098] FIG. 6 This figure shows results of blended SEAM II-unconventional data. Figures (a)-(c) are the blended data (input), unblended data (ground truth), and the deblending result (prediction) respectively in the common receiver domain. Figures (d)-(f) are the blended data unblended data, and the deblending result in the common shot domain.

    [0099] FIG. 7 This figure shows the zoom-in of the boxes shown on FIG. 6.

    [0100] FIG. 8 This figure shows the results of the “half” blended SEAMII-unconventional data. Figures (a)-(c) are the blended data (input), unblended data (ground truth), and the deblending result (prediction) respectively in the common receiver domain. Figures (d)-(f) are the blended data unblended data, and the deblending result in the common shot domain.

    [0101] FIG. 9 This figure shows the zoom-in of the boxes shown on FIG. 8.

    [0102] FIG. 10 This figure shows results of the “continuous” blended BP2004 data. Figures (a)-(c) are the blended data (input), unblended data (ground truth), and the deblending result (prediction) respectively in the common receiver domain. Figures (d)-(f) are the blended data unblended data, and the deblending result in the common shot domain.

    [0103] FIG. 11 This figure shows the zoom-in of the boxes shown on FIG. 10.

    DETAILED DESCRIPTION OF THE INVENTION

    [0104] As will be appreciated by one skilled in the art, aspects of the present invention may be embodied method that may be implemented as a computer program product, at least those parts manipulating the acquired seismic data.

    [0105] In order to show an specific embodiment of the disclosed methods in a detailed manner, let custom-character represent the 1-D Fourier transformation along the time axis, that is, for what has been named a trace, custom-character.sup.−1 its inverse and f.sub.θ(⋅) is a denoising convolutional neural network (CNN) with parameters θ. Further, let C denote the array concatenation operator which takes arrays and stacks them into a single higher-dimensional array, e.g. Π=CΠ.sub.k. With no unblended data available, a model that minimizes the unsupervised blending loss function is trained, being the loss function

    [00004] L θ ( x , f θ ( .Math. ) ) = 1 N .Math. k = 1 N .Math. x ( k ) - - 1 C Γ m H Γ m ( f θ ( x ( k ) ) ) m .Math. 1

    [0106] and ={x.sup.(1), . . . , x.sup.(N)}, where the input x=custom-character.sup.−1CΓ.sub.m.sup.HΠ.sub.k.sup.b, is the pseudo-deblended data in the common receiver domain. f.sub.θ(x) gives the prediction of the deblended data. The prediction of the model is blended again by the known blending matrix Γ. This result is the blended prediction: Γ.sub.m(custom-characterf.sub.θ(x.sup.(k))).sub.m. The term custom-character.sup.−1CΓ.sub.m.sup.HΓ.sub.m(custom-characterf.sub.θ(x.sup.(k))).sub.m is the pseudodeblended output of the model's prediction. The mean absolute error is taken from the pseudodeblended output of the model's prediction and the original pseudo-deblended data, i.e. the input x. The action of the blending followed by the pseudodeblending operator custom-character.sup.−1CΓ.sub.m.sup.HΓ.sub.m(custom-character(⋅)).sub.m can be implemented efficiently and directly in the time domain by first combining the traces according to the blended-acquisition time dithering code, then each bended shot-gather is copied a number of times equal to the corresponding source group size, and finally each of these copies is time-shifted (dithering decoded) to undo the delays introduced in the field. By adopting this blending loss, the deblending task is accomplished in a unsupervised way with the blended data and the blending matrix.

    [0107] Direct minimization on this blending loss with a traditional deep CNN leads to identity mapping. It will produce the pseudo-deblended result which is exactly the least square solution. This phenomenon could be observed in FIG. 1. In this test, the blending factor is 2, i.e. every two shots are blended. Therefore, a simple scaling of ½ on both the noise and signal can perfectly satisfies a blending loss of 0.

    [0108] To avoid this local minimum, a strong constraint for the model by making it trace-wise blind is proposed. According to this proposal, the coherent signal of the traces are essentially reconstructed from their adjacent traces, while not looking at themselves in the input. Thus for a specific trace, there would be no chance for its blending noise to be mapped from input to output. For those noises in its adjacent area, the model automatically neglect them due to their discontinuity and the convergence of the blending noise.

    [0109] The illustration is shown in FIG. 2 where at the left side, two vertical bands with dashed with inclined lines shows two patches used when inputting information to recover the trace marked at the right side, ignoring the original noisy trace itself. This is why we name “blind-trace” the reference to the receptive field of the CNN.

    [0110] According to this embodiment, the blind-trace convolutional neural network is constructed for the deblending task.

    [0111] The network structure is shown in FIG. 3. A deep CNN is constructed based on a modified version of the U-net with all convolutional layers and down-sampling layers changed for the receptive field restriction. In the beginning, the input seismic image is clockwise rotated twice with angles 90° and 270° so that both the left and right sections of the original receiver gather is moved upward, and they are separately inputted to the network. The network is composed of two parts, namely the blind-trace U-net and the merging layers, the merging layers being previously disclosed as the responsible for carrying out the combination of the two feature maps.

    [0112] In the blind-trace U-net, the receptive field of each layer is strictly restricted to the upper half area for each row so that the model can extract the coherent features for each trace based on its left or right adjacent area in the original input respectively.

    [0113] After the blind-trace U-net, the two output feature maps are cropped at the bottom and padded at the top such that the target traces can be excluded from its receptive field to fulfill the “blind-trace” purpose.

    [0114] Then, they are rotated back (270° and 90° accordingly) and concatenated, followed by two consecutive 1×1 convolutional layers that integrate the feature maps and squeeze the channel size to 1.

    [0115] In FIG. 3, a trace example is marked in black, the informative left and right areas of it at each stage are marked using two different patterns, a first pattern with inclined lines for the left patch and a second pattern with a square grid for the right patch respectively.

    [0116] According to a preferred embodiment that may be applied to any of the disclosed examples, before the padding and cropping operations, the trace is buried under the patches, and at the end of the process, the trace will be excluded from the patches. For the edges, there is only one side of the patch that is informative and zeros are padded on the other side.

    [0117] The image is processed only once and the network gives the predictions for all traces at the same time.

    [0118] Specifically, to constrain the receptive field within the upper half for all rows in the blind-trace U-net, the convolutional layers and down-sampling layers are changed as shown in FIG. 4.

    [0119] Convolution: The feature maps are padded with zeros before each convolutional layer. Given the size of the filter k×k, └k/2┘ lines are padded on the top of the feature maps, and crop the └k/2┘ lines at the bottom after the convolution.

    [0120] Max-pooling: The feature maps are padded on the top with one line and crop one at the bottom before max-pooling. There is no need to change the up-sampling layers since it does not affect the receptive filed after the modification for the down-sampling layers.

    [0121] With the aforementioned changes, for a single trace in the feature maps, it takes the information purely from itself and the area above in the input. No information below it can leak into the result.

    [0122] According to this embodiment, the rotated inputs are projected to 32 feature maps in the first layer. There are 4 contracting blocks and 4 expanding blocks in the blind-trace U-net, each of which consists of two consecutive modified convolutional layers and a max-pooling/up-sampling layer. The feature maps are doubled in the last convolutional layer of each contracting block and halved correspondingly in the expanding blocks.

    [0123] Each blind-trace convolutional layer in the encoder is replaced by a blind-trace residual block (i.e. a residual block with all convolutional layers modified). Batch Normalization is adopted before each activation in the decoder.

    [0124] The last two 1×1 convolutional layers project the concatenated 64 feature maps to 32 and 1 respectively. It has been used ReLU activation in the intermediate convolutional layers and leaky ReLU for the last 1×1 convolutional layer.

    [0125] The proposed blind-trace network with blending loss combines the merits of both the conventional filter-based method and the inversion-based method.

    [0126] The large number of the weights in the U-Net endows the network to achieve a complex non-linear filter to extract information from the coherent signal, meanwhile, for the low-SNR area where filters could not obtain any coherent information, the minimization of the blending loss reconstructs the coherent signal underneath through nonlinear inversion. This is a one-stage deblending framework and does not require much exhaustive and meticulous hyperparameter tuning.

    [0127] Regarding the experiment, three different types of blending schemes are used to demonstrate the performance of the proposed method. The synthetic data examples can be seen in FIG. 5.

    [0128] “alternate” blending. In this blending strategy, two consecutive shots are blended using short and random dithering time. After pseudo-deblending, the noise will cover the whole common receiver gather as shown in FIG. 5(a). Since the delay is short, the amplitudes of both signals and noises are at the same level. SEAMII-unconventional data is used for deblending test (SEAM denoting Society of Exploration Geophysicists and “II” denoting Phase II, Phase I refers to Subsalt imaging in tertiary basins and Phase II refers to Land seismic challenges. https://seg.org/News-Resources/Research-Data/Land-Seismic-Data)

    [0129] 200 unblended shots are collected with the size (1600(receiver numbers)×750(time samples)). Random delays in the range of (0, 1]s are added to even shots and blended with the odd shots. The sampling period is 6 ms.

    [0130] “half” blending. The second half of the shot gathers are shifted with relatively long delays and added to the first half. Noises in the pseudo-deblended common receiver gathers concentrate in the two corners of the image as shown in FIG. 5(b). Since the amplitude of blending noise from early arrivals is often much higher than late arrival signals, the SNR (signal/noise rate) is very low. SEAMII-unconventional data is used for deblending test. 200 unblended shots are extracted from the SEAMII-unconventional dataset with size (1600(receiver numbers)×750(time samples)). Shot 101 to 200 are blended accordingly to shot 1 to 100 with added random delays in the range of [2, 3]s. the sampling rate is 6 ms.

    [0131] “continuous” blending. Simultaneous source of streamer acquisition is simulated in this experiment. A BP2004 velocity model is used for generating 200 consecutive unblended shots.

    [0132] (BP2004 denotes the 2004 British Petroleum Velocity Benchmark. https://software.seg.org/datasets/2D/2004_BP_Vel_Benchmark/eage_abstract.pdf)

    [0133] The shots have been recorded every 5 to 6 seconds continuously. There are 959 receivers for each record and 3334 samples along the time axis (dt=6 ms, hence T≈20 s). Approximately 4 shots are blended at different locations in the common receiver gathers. The noise level, in this case, is much higher than the other blending schemes as shown in FIG. 5(c).

    [0134] FIGS. 5(d), 5(e) and 5(f) are the corresponding unblended data of that shown in FIGS. 5(a), 5(b) and 5(c) respectively.

    [0135] Results of One-Stage Deblending

    [0136] In order to avoid massive matrix multiplication in the f-x domain in loss computation, the blending and pseudo-deblending is performed in the time-space by applying dithering codes to the shots, followed by summation and cropping. In the first stage, the weights and bias are randomly initialized and trained with the Adam optimizer.

    [0137] FIGS. 6 to 11 shows the deblending results using the blind-trace network and one-stage deblending for different blending schemes.

    [0138] FIG. 6 presents the result of the “alternate” blending scheme. Short shifting between two consecutive shots prevents the early arrival noises from mixing with the late arrivals, but the noises distributed across the whole image. In FIG. 7 a zoom-in area of FIG. 9 is presented. From the demonstration in both common receiver domain and common shot domain, it is shown that the one-stage deblending can effectively remove the noises for both early arrivals and late arrivals, preserving the details for both low and high-frequency patterns.

    [0139] “Half” blending scheme introduces more challenges for deblending models. Due to the concentration of the noises and long dithering time, the weak amplitudes in late arrivals are severely contaminated by strong noises from the early arrivals. As can be seen in FIG. 8 and the emphasized area of it in FIG. 9, traditional filters tend to fail due to the loss of coherency. Besides, the heteroscedasticity problem would increase the difficulty for the iterative method, in which harsh requirements are imposed for the accuracy of amplitude estimation. As shown in FIGS. 9(c) and (f), the proposed method not only works as a denoising model but also acts as non-linear inversion from the blending loss with the blind-trace regularization. The signals submerged by strong blending noises are reconstructed by optimization.

    [0140] In “continuous” blending scheme, more unblended shots are collected within one blended shot. The noisy area enlarges to the whole image, and the signal is covered by multiple layers of noises. Less coherent information can be attained by the model.

    [0141] However, the one-stage deblending can still give a very good result as demonstrated in FIGS. 10 and 11. The following table shows the blending loss comparison in different stages:

    TABLE-US-00001 Blending scheme Pseudo-deblended Training-stage Tuning-stage “alternate” 0.5178 0.1594 0.1141 “half” 0.0674 0.0169 0.0072 “continuous” 0.3250 0.1176 0.0789