A Deep Learning-Based Method for Reconstructing Multiple Wave Driven Seismic Data

Abstract

This application relates to the field of seismic signal processing technology and discloses a deep learning-based method for reconstructing multiple wave-driven seismic data. The method first constructs a focusing operator using the acquired data domain seismic shot gather and performs a focusing transformation to obtain the corresponding focusing domain data. Then, it removes focal points and noise regions from the focusing domain data, establishing a mapping relationship between the data domain and the focusing domain. Next, it extracts data blocks using a sliding window to construct training samples and establishes a nonlinear mapping between focusing domain features and data domain signals through a deep neural network. Finally, it inputs the focusing domain data to be reconstructed into the trained model, performs block prediction and Gaussian weighted smoothing of overlapping areas, and obtains the reconstructed seismic data. Finally, it uses the reconstructed data to fill in the missing positions of the original seismic data. This invention introduces focusing transformation into a deep learning reconstruction framework, achieving high-precision reconstruction of seismic data by mining the effective information contained in multiple waves, while maintaining the dynamic characteristics and structural continuity of the wavefield.

Description

Technical Field

[0001] This invention relates to the field of seismic signal processing technology, specifically to a method for reconstructing multi-wave driven seismic data based on deep learning. Background Technology

[0002] During the acquisition of seismic exploration data, the raw seismic data often suffers from missing traces or irregular spatial sampling due to objective factors such as surface obstacles, restricted areas, acquisition costs, and detector malfunctions. This data incompleteness disrupts the continuity of the wavefield, leading to severe spurious frequency or repositioning errors in subsequent processing steps such as multiple wave suppression, velocity analysis, and migration imaging, thus affecting the accuracy of subsurface geological structure interpretation. Therefore, high-precision reconstruction of missing seismic data at the beginning of the data processing workflow is a crucial step.

[0003] Traditional seismic data reconstruction methods primarily rely on signal processing theory, utilizing the sparsity of seismic data in the transform domain for interpolation. These methods typically assume that the seismic wavefield exhibits local linear or parabolic characteristics, mapping the data to the sparse domain using mathematical tools such as Fourier transform, Radon transform, or curvelet transform, and then using inverse transforms to recover missing data. However, when faced with large-scale continuous data gaps or complex subsurface structures leading to non-stationary wavefield variations, these mathematically transform-based methods struggle to accurately describe the wavefield's dynamic characteristics, are prone to spurious frequency interference, and have limited ability to recover effective weak signals.

[0004] With the development of artificial intelligence technology, data-driven methods based on deep learning have been introduced into the field of seismic data reconstruction. Existing deep learning reconstruction methods mostly establish a direct mapping relationship from missing data to complete data. However, in practical processing, conventional methods often treat surface-correlated multiples as coherent noise and suppress them, ignoring the fact that multiples are wavefields that have undergone multiple reflections from subsurface interfaces, and that they often contain richer subsurface illumination information than primary waves. Without fully utilizing the physical coupling relationship between multiples and primary waves, directly training the network can easily lead to it being dominated by strong-energy primary waves. This makes it difficult for the network to extract effective information about missing locations from the complex wavefield through this physical correlation, especially in low signal-to-noise ratio regions. The amplitude fidelity and phase axis continuity of the reconstruction results often fail to meet the requirements of high-precision exploration. Summary of the Invention

[0005] To address the shortcomings of existing technologies, this invention provides a deep learning-based method for reconstructing multi-wave driven seismic data. This method solves the problems that existing seismic data reconstruction techniques struggle to handle large-scale continuous gaps, and that conventional deep learning methods fail to fully utilize the physical information contained in multi-waves, resulting in limited wavefield reconstruction accuracy and fidelity under complex geological conditions.

[0006] To achieve the above objectives, the present invention is implemented through the following technical solution: a deep learning-based multi-wave driven seismic data reconstruction method, which establishes a mapping relationship between the features of the focusing domain and the original seismic data domain by combining physical transformation and data-driven approach.

[0007] This method first performs a focusing transformation. Using the acquired data domain seismic shot gather, the primary wavefield with strong energy characteristics is separated through preprocessing and used as the focusing operator. A sampling operator matrix is ​​constructed, which numerically identifies the spatial existence of receivers relative to shot points in the seismic record. Correlation operations are performed between the focusing operator and the seismic shot gather containing multiples, and the results are combined with the sampling operator matrix to generate focused domain data. Physically, this focused domain data consists of a focal point located at the origin and a background term represented by a reduced-order multiple, where the focal point corresponds to the focusing of primary wave energy, and the background term contains wavefield characteristics that implicitly contain subsurface medium information transformed from multiples.

[0008] Secondly, feature data for network input is constructed. A spatiotemporal cutoff window is defined, which is a parabolic function that moves with the shot point, and the curvature of the parabola is determined by the phase axis of the seismic record. The sample values ​​within the window are set to zero or smoothed and attenuated to remove the focal point and its surrounding noise region, thereby extracting the focal domain background term. Amplitude normalization is performed on this background term to ensure that its numerical distribution meets the input requirements of the deep neural network. At this point, a missing pattern consistent with the data domain is set, thus establishing a physical mapping basis between the focal domain background term (containing missing information) and the missing seismic data in the data domain.

[0009] Then, a deep neural network was constructed and trained. A sliding window was used to simultaneously scan and extract data blocks from the focal domain data and the corresponding data domain seismic shot gathers to construct training samples. The training samples included the focal domain data as input and the raw seismic data as labels. To enhance the model's generalization ability, spatial trace order flipping, phase reversal, and random gain adjustment based on the physical properties of the wave equation were performed on the training samples, and band-limited Gaussian noise was superimposed. A ResUNet network structure containing an encoder, decoder, and skip connections was constructed. The encoder extracts multi-scale energy distribution features from the focal domain data through cascaded convolutions, the decoder restores the spatiotemporal resolution of the seismic data through upsampling, and the skip connections transfer shallow high-resolution features to deeper layers to preserve wavefield details.

[0010] During model training, a hybrid loss function is constructed to constrain the updating of network parameters. This hybrid loss function consists of a weighted sum of an L1 norm loss term and a structural similarity loss term. The L1 norm loss term constrains the absolute error between the predicted and actual data at sample points; the structural similarity loss term, calculated based on the local mean, standard deviation, and covariance of the predicted and labeled data, constrains the continuity of seismic phase axes and texture structure features. An adaptive moment estimation optimizer iteratively updates the network parameters based on the gradient of the hybrid loss function until a stable nonlinear mapping relationship is established between the features of the focused domain and the signal in the data domain.

[0011] Finally, data reconstruction and backfilling are performed. The focal domain data to be reconstructed is input into a trained deep neural network model, and block prediction is performed according to a set overlap rate to obtain predicted data blocks. To eliminate the boundary effect caused by block prediction, a two-dimensional weight matrix with the same size as the predicted data blocks is constructed. The values ​​of this matrix follow a two-dimensional Gaussian distribution, with the largest weight at the center and exponentially decreasing towards the edges. A full-size signal buffer matrix and weight buffer matrix are established. The Hadamard accumulation of each predicted data block and the weight matrix is ​​added to the signal buffer matrix, and the weight matrix is ​​accumulated to the weight buffer matrix. By dividing the elements of the signal buffer matrix by the corresponding elements of the weight buffer matrix and performing inverse normalization, globally consistent reconstructed data is obtained. The dimensionality consistency between the reconstructed data and the original data is verified. The spatial coordinates of the original data are traversed, and according to the identifier of the sampling operator matrix, the values ​​of the reconstructed data are written only at the positions marked as missing, finally obtaining a complete seismic data profile.

[0012] This invention provides a deep learning-based method for reconstructing multi-wave driven seismic data. It offers the following advantages: 1. This invention utilizes a focusing transformation to map seismic data to a focusing domain. The focusing transformation recovers missing seismic traces through a multiple wave reduction process and highlights background features containing multiple wave information by removing focal points. This allows deep neural networks to focus on learning the nonlinear mapping relationship between multiple waves and missing data. Compared with traditional methods that only use single waves or direct interpolation, this invention effectively utilizes multiple wave information as a constraint, significantly improving the accuracy of recovering effective wavefield features in large-area continuous missing scenarios.

[0013] 2. This invention constructs a ResUNet network model and uses a hybrid loss function combining L1 norm and structural similarity for training constraints. The L1 norm constrains the absolute error of the reconstructed data in terms of amplitude values, ensuring amplitude fidelity. The structural similarity loss term strengthens the preservation of the continuity of seismic phase axes and texture features. This dual constraint mechanism avoids waveform blurring or structural distortion that may be caused by a single loss function, ensuring that the reconstructed seismic profile has both accurate dynamic characteristics and maintains a clear geological structure.

[0014] 3. This invention adopts a prediction strategy based on overlapping blocks and combines it with Gaussian weighted smoothing technology. By applying a two-dimensional Gaussian weight matrix with the largest center weight and exponential decay at the edges to the overlapping areas of the predicted data blocks, it effectively eliminates numerical abrupt changes and block artifacts generated at the splicing boundaries of the block prediction. This processing ensures the continuity and consistency of the reconstructed seismic data in the entire spatial range and improves the signal-to-noise ratio and overall quality of the final seismic data volume. Attached Figure Description

[0015] Figure 1 This is a flowchart of the method of the present invention; Figure 2 This is a system architecture diagram of the present invention; Figure 3 The seismic data reconstruction framework proposed for the implementation of this invention; Figure 4 The ResUnet network model used in the embodiments of this invention; Figure 5 This is a schematic diagram of the residual block structure used in an embodiment of the present invention; Figure 6 The original profile and corresponding focused domain profile of the seismic record synthesized by the planar model in this embodiment of the invention; Figure 7 This is a cross-section showing the comparison results between the seismic data reconstructed from the level model in this embodiment of the invention and two other methods; Figure 8 A waveform comparison of seismic data reconstructed from a level model in an embodiment of the present invention with that of two other methods; Figure 9 The speed model of the complex model used in the embodiments of the present invention; Figure 10 This refers to the original profile and corresponding focused domain profile of the complex model synthesized seismic record in an embodiment of the present invention; Figure 11 This is a cross-section showing the comparison results of seismic data reconstructed from a complex model in an embodiment of the present invention with two other methods; Figure 12 A waveform comparison of seismic data reconstructed from a complex model in an embodiment of the present invention with that of two other methods; Figure 13 This is a comparison diagram of single-channel waveforms of the present invention; Figure 14 This is a bar chart comparing the signal-to-noise ratio and quantitative indicators of structural similarity of the reconstruction results from the convex set projection method, the supervised learning method, and the method of the present invention in the embodiments of the present invention. Detailed Implementation

[0016] The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0017] Please see the appendix Figure 1 - Appendix Figure 2 This invention provides a deep learning-based method for reconstructing multi-wave driven seismic data, comprising the following steps: S1. Construct a focusing operator using the acquired seismic shot gathers and perform a focusing transformation on the seismic shot gathers to obtain the focusing domain data corresponding to the data domain. S2, in the focal domain data, remove focal points and noise areas, set a missing mode consistent with the data domain, and construct a mapping relationship between the data domain and the focal domain; S3 uses a sliding window to sample the data profile with a preset overlap rate, extracting data blocks as training samples. The data blocks include focal domain data as input samples and original seismic data as labels. Data augmentation is performed on the training samples, and they are divided into training set and validation set according to a preset ratio. S4. Construct a deep neural network model, use the training set to train the model, learn and establish the nonlinear mapping relationship between the features of the focal domain and the seismic signals in the data domain; S5. Input the test set data to be processed into the trained deep neural network model, divide it into blocks according to the overlap rate at the time of sampling, predict the output, and perform Gaussian weighted smoothing on the overlapping area to obtain the reconstructed missing seismic data. S6 fills the missing seismic data after reconstruction into the missing positions of the original seismic data, realizing seismic data interpolation and obtaining a complete seismic data profile.

[0018] Please see the appendix Figure 2 The seismic data reconstruction framework in this embodiment mainly includes a focus transformation module, a mapping construction module, a sample processing module, a model training module, a prediction smoothing module, and a data imputation module. Each module is connected in sequence according to the data flow direction and works together to complete the reconstruction process from the original seismic data with missing data to the complete seismic data.

[0019] The focusing transformation module is configured to acquire the collected seismic shot gather data and construct focusing operators based on the strong energy phase axes or first waves in the seismic data. The focusing transformation module then performs focusing transformation operations on the seismic shot gather using these focusing operators to obtain focused domain data corresponding to the data domain. This process utilizes the physical order reduction mechanism of multiple waves to refocus the energy of higher-order wavefields, aiming to achieve physical reconstruction of the wavefield at missing sources or receivers.

[0020] The mapping construction module connects to the focusing transformation module and receives the focusing domain data. Since the primary wave energy within the focusing domain is highly concentrated near the origin, forming focal points, the mapping construction module removes these focal points and noise regions from the focusing domain data. This module sets a missing data pattern consistent with the data domain to construct a mapping relationship between the data domain and the focusing domain. The focusing domain data, after removing focal points, serves as a weighted representation of the original seismic wavefield in a physical sense, providing feature input for subsequent network learning.

[0021] The sample processing module is used to construct the training dataset required for the neural network. This module uses a sliding window to sample the data profile at a preset overlap rate, extracting data blocks as training samples. These data blocks include focal region data as input samples and raw seismic data as labels. The sample processing module performs data augmentation operations on the training samples and divides them into training and validation sets according to a preset ratio.

[0022] The model training module constructs a deep neural network model and trains it using a training set. Through an iterative optimization process, this module learns and establishes a nonlinear mapping between the features of the focal domain and the seismic signals in the data domain. During training, the module uses a validation set to evaluate model performance and adjusts network parameters to improve the model's generalization ability.

[0023] The prediction smoothing module inputs the test set data to be processed into the trained deep neural network model. This module performs block prediction output according to the overlap rate at the time of sampling, and applies Gaussian weighted smoothing to the overlapping areas. Through weighted smoothing, the boundary effects that may be caused by block prediction are eliminated, resulting in the reconstructed missing seismic data.

[0024] The data infilling module fills in the missing locations of the original seismic data with the reconstructed missing seismic data. This infilling operation performs seismic data interpolation, ultimately obtaining a complete seismic profile. This profile contains restored near-offset information, which can be used for subsequent seismic imaging and inversion processing.

[0025] In step S1, a focusing operator is constructed using the acquired seismic shot gather, and a focusing transformation is performed on the seismic shot gather to obtain the focusing domain data corresponding to the data domain. This embodiment of the invention utilizes the rich illumination information carried by multiple waves during underground propagation to compensate for the lack of near-offset information through a physical transformation mechanism.

[0026] Specifically, the implementation of step S1 includes the following logical process: In sub-step S11, a focusing operator for wavefield focusing is constructed. The focusing operator is based on strong-energy in-phase axes or primary waves in the seismic data. Its physical essence is the impulse response representing the main subsurface reflection interface. In this embodiment, after performing conventional preprocessing on the original seismic record (such as removing direct waves and suppressing random noise), the primary wavefield with strong energy characteristics is separated using FK filtering or Radon transform, and it is defined as the focusing operator. This operator serves as the benchmark for wavefield propagation and is used for subsequent cross-correlation calculations with the full wavefield data to achieve energy refocusing.

[0027] In sub-step S12, to mathematically describe the insufficient spatial sampling during the acquisition process, a sampling operator matrix is ​​defined and introduced. The sampling operator matrix is ​​used to identify the spatial existence of receivers relative to shot points in seismic records. elements The definition is as shown in the formula: ; In the formula, Represents the first in the matrix Line number The element values ​​of the column, Indicates the spatial coordinate index of the receiver point. Indicates the spatial coordinate index of the shot point. Represents the absolute value of the shot-receiver offset. In this embodiment, the threshold value representing the missing range on one side of the shot point is... The value is determined based on the physical configuration of the field seismic acquisition system, specifically corresponding to the minimum safe offset between the seismic source and the nearest geophone or the blind zone radius in the actual observation system. When the calculated source-receiver offset is greater than this threshold, it is marked as 1, indicating that the data is valid; otherwise, it is marked as 0, indicating that the data is missing.

[0028] In substep S13, a focusing positive transform is performed to generate focused domain data, utilizing the constructed focusing operator. Seismic records containing multiple waves Cross-correlation calculations are performed. Based on the physical principle of surface correlation multiple cancellation, multiples can be considered as the convolution of the primary wave and the surface reflection coefficient. Therefore, by eliminating the propagation path of the primary wave through correlation calculations, the energy of higher-order multiples can be focused at time zero and zero offset. Considering missing data, this focusing forward transformation process is shown in the formula: ; In the formula, This indicates the result of data transformation in the focal region; Indicates the focusing operator; This represents data containing multiples; The sampling operator matrix; It is a surface correlation operator that physically characterizes the combined effect of surface reflectivity and subsurface medium response; It is an identity matrix. It should be noted that although the original records... There are deficiencies, but due to the focusing operator Typically constructed in a transform domain (such as the Radon domain) or after interpolation regularization, it possesses spatial continuity. Therefore, the above operations can physically extract the full-wave field background information implicit in the multiple waves.

[0029] Please see the appendix Figure 6 To visually illustrate the wavefield morphology evolution during the above transformation process, we will take the composite data of the flat-layer model as an example: Figure 6 (a) shows the original complete seismic record under ideal conditions, with its phase axis exhibiting a complete hyperbolic shape; Figure 6 (b) shows the data state in the simulated actual acquisition where there is a missing near offset, and there is no effective signal in the middle area; Figure 6 (d) visually demonstrates the... Figure 6 (b) Focus domain data obtained after performing focus transformation As can be seen, the seismic signal in the spatiotemporal domain is mapped to a specific energy distribution in the focusing domain, which includes a strong energy focusing point located at the top center and a background texture extending into the depth. This background texture physically implies full-wave field information, forming the basis for the input of subsequent deep neural networks.

[0030] During this calculation, taking into account Potential bandwidth limitations or numerical instabilities may exist when calculating its inverse matrix. In the case of deconvolution, as a preferred implementation, a small pre-whitening factor is added to the denominator or a damped least squares method is used for stabilization to prevent numerical divergence.

[0031] In step S2, which involves removing focal points and noise regions from the focal domain data, setting a missing pattern consistent with the data domain, and constructing a mapping relationship between the data domain and the focal domain, this embodiment aims to solve the problem of how to adapt the physically transformed data to the input characteristics of deep learning networks.

[0032] Specifically, the implementation of step S2 includes the following logical process: In sub-step S21, a physical gating process is performed on the energy distribution characteristics of the focused domain data. Given the focused domain data... It is composed of the identity matrix With background items It is formed by linear superposition, where the identity matrix is In the physical spectrum, this manifests as a strong energy pulse located at time zero and zero offset. Although this pulse has enormous energy, it is essentially an autocorrelation main lobe of the focusing transformation and does not contain the scattering information of the subsurface medium to be reconstructed. If directly retained, its extremely large dynamic range would suppress the weak effective signal, making it difficult for network training to converge. Therefore, this embodiment performs a cut-off operation. The specific cut-off range is not arbitrarily set, but determined based on the energy distribution characteristics of the primary wave within the focusing domain. As a preferred approach, a spatiotemporal cut-off window is set, which is a parabolic function that moves with the shot point, and the curvature of the parabola is determined by the phase axis of the seismic record. Data sample values ​​within this window are forced to zero or smoothed by a Hanning window, thereby removing interference from high-energy focusing points while fully preserving the background term containing multiple wave path information and subsurface medium reflection characteristics. .

[0033] In sub-step S22, the background item after removing the focal point is processed. Numerical preprocessing is performed. To accommodate the sensitivity of deep neural networks to the distribution of input data, amplitude normalization is executed. When calculating normalization factors (such as maximum absolute amplitude or root mean square amplitude), a small regularization constant is introduced into the denominator to prevent anomalies of zero denominators caused by silent segments or dead channels. (e.g., 10) -6 The processed data is mapped to a standard numerical range (e.g., [-1, 1]) to ensure that the data energy of different shot sets is on the same order of magnitude, eliminating training bias caused by differences in the excitation energy of the seismic source.

[0034] In sub-step S23, the mapping logic between the background terms in the focus domain and the original data in the data domain is constructed. This is a crucial step connecting the physical model and the data-driven model. Based on the physical principles of the wave equation, the background terms... In physical terms, it is the raw seismic data. Surface-related operators The convolutional response can be viewed as a blurred or encoded representation of the original wavefield. Although It cannot be directly equated with However, there is a strict physical causal relationship between the two. Therefore, this embodiment does not require explicitly solving a complex deconvolution operator, but directly constructs input-output pairs: the background term after cropping and normalization... As input features to the network, the known raw seismic data at the corresponding locations are used. As a supervisory label for the network. Through this construction, the deep neural network is trained to learn from the blurred wave field after convolution. Restore clear original wave field The nonlinear inversion process.

[0035] In step S3, a sliding window is used to sample the data profile with a preset overlap rate, extracting data blocks as training samples, and performing data augmentation and set partitioning. This embodiment constructs the data foundation for driving a deep neural network to learn the mapping relationship between physical transformation features and the real wave field.

[0036] Specifically, the implementation of step S3 includes the following logical process: In sub-step S31, a sliding window sampling strategy based on spatiotemporal overlap is performed. Given that the original seismic profile typically has a large data dimension in both the temporal and spatial directions, directly inputting it in its entirety into a deep learning network would exhaust existing resources and make it difficult to capture fine local texture features. Therefore, this embodiment uses a fixed-size two-dimensional sliding window to simultaneously scan and extract data from the focused domain data profile and the corresponding original data profile. To ensure that the sampling process covers different dip angles of the seismic wavefield and to ensure high continuity of the subsequent reconstruction results at the block stitching points, the sliding window's step size is set to a value smaller than the window size. As a preferred implementation, the step sizes in the temporal and spatial directions are set to 25% to 50% of the corresponding dimensions of the window, thereby forming at least 50% spatial overlap between adjacent samples. The specific size of the window... (in The number of time sampling points, The number of spatial traces is determined based on the dominant frequency and maximum apparent velocity of the seismic data in the work area, and is usually set to a range of 64×64 to 256×256 to ensure that each data block contains complete local linear features of the same phase axis.

[0037] In sub-step S32, a pairwise training sample set under physical constraints is constructed. For each window truncation operation, the sub-block extracted from the preprocessed focal region data is defined as the input sample. Sub-blocks at the same location extracted from the corresponding original seismic data are defined as label samples. The selection of training samples strictly follows the principle of "known supervision of unknowns," meaning that only the effective data area actually acquired by the detector in the original seismic record (i.e., the sampling operator) is included. Samples are extracted from the region. Although the ultimate goal of this invention is to reconstruct missing near-offset data, by training the network across the entire observation area, the model can learn a robust nonlinear mapping between the "focused domain wavefield after removing the focal point" and the "true wavefield". Based on the statistical stationarity and physical self-similarity of seismic wavefield propagation in the space medium, this mapping law is spatially transferable and can therefore be generalized to wavefield prediction in missing regions.

[0038] In sub-step S33, data augmentation based on the physical properties of the wave equation is implemented. To prevent overfitting of the deep neural network on limited seismic samples and to improve the model's adaptability to complex geological conditions and environments with different signal-to-noise ratios, this embodiment expands the basic sample set. Unlike the geometric transformations in conventional image processing, the data augmentation operations used in this embodiment strictly follow the physical laws of seismic wave propagation: including simulating spatial trace order reversal based on the wavefield reciprocity principle, simulating phase reversal based on different excitation polarities, and simulating random gain adjustment based on different source energies. Furthermore, to address potential interference in actual field data, band-limited Gaussian noise of varying intensities is dynamically superimposed on the training samples, forcing the network to learn the ability to separate and extract coherent signals from low signal-to-noise ratio inputs, thereby improving the model's robustness in processing real-world data.

[0039] In sub-step S34, the training set and validation set are divided. The augmented sample set is randomly shuffled and distributed according to a preset ratio. To ensure the objectivity of the model evaluation, the division of the validation set should avoid direct spatial overlap with the training set. That is, the block-based leave-out method is preferred over completely random sampling to prevent data leakage caused by the spatial correlation of seismic data.

[0040] In step S4, a deep neural network model is constructed, and the model is trained using a training set to learn and establish a nonlinear mapping relationship between the features of the focused domain and the seismic signals in the data domain.

[0041] Please see the appendix Figure 3 The area within the blue dashed box at the top represents the missing focus data during the training phase. Feedback to the network as input data The output result is obtained through forward propagation. Construct the loss function and calculate the output results. With the original data missing as labels The differences between them are analyzed, and the network parameters are iteratively updated using the backpropagation algorithm.

[0042] Specifically, the implementation of step S4 includes the following logical process: In sub-step S41, a deep neural network model for feature mapping is constructed. Given the complex wavefield characteristics of seismic data across spatiotemporal scales, this embodiment preferably uses the ResUNet model, which combines residual learning with the U-Net architecture, as the backbone network; please refer to the appendix. Figure 4The network has a U-shaped symmetrical structure. Input data enters the network from the left. The left path consists of multiple cascaded residual blocks and alternating max pooling operations (shown by the orange downward arrow in the figure), used for progressive downsampling to extract features. The right path consists of bicubic upsampling operations (shown by the green upward arrow in the figure) and residual blocks, used to progressively restore the spatial dimension of the data. Skip connections (shown by the blue horizontal arrow in the figure) are set between corresponding layers on the left and right sides, directly concatenating the feature maps from the encoder to the decoder. The network ends with a 3x3 convolutional layer for feature aggregation, finally generating the output data. The network mainly consists of an encoder, a decoder, and skip connections. The encoder is configured to extract the multiple wave energy distribution features at different scales in the focal domain data through multi-layer cascaded convolutional operations and compress the data into a low-dimensional high-level semantic feature vector.

[0043] Please see the appendix Figure 5 The basic unit of the cascaded convolution operation is Figure 5 (a) shows the pre-activated residual module. This module adopts a dual-path structure, including a main path consisting of two cascaded convolutional modules and a side path branch containing a skip module. The outputs of the two paths are fused element-wise using an adder (+). The specific construction of the convolutional module is as follows: Figure 5 As shown in (b), following the pre-activation design principle, the data flow sequentially passes through a two-dimensional group normalization layer, a sigmoid linear unit activation function (SiLU), and a 3×3 two-dimensional convolutional layer. This structural design can effectively prevent the gradient vanishing problem in deep network training and improve the non-linear expressive power of feature extraction.

[0044] The decoder performs the corresponding upsampling operation, aiming to gradually restore the abstract feature vector to the spatiotemporal resolution of the original seismic data. Skip connections are introduced between the corresponding levels of the encoder and decoder. Their physical function is to directly transmit the high-frequency phase information and small amplitude details of the seismic signal retained in the shallow network to the decoding path, effectively alleviating the waveform distortion problem caused by multiple downsampling in the deep network.

[0045] In sub-step S42, the optimization objective for network training is defined. To ensure that the network parameters accurately represent the mapping from the residual components of the focused domain to the original complete wavefield, a global optimization problem needs to be constructed. The constructed training dataset is input into the network for iterative training, and the parameter optimization process is described by the following formula: ; In the formula, This represents the set of optimal network parameters after optimization and convergence. This represents the total number of training samples; This represents the loss function, used to quantify the difference between the predicted and actual values; Indicates the first Label data for each sample (i.e., the original earthquake record); Indicates the first The network input data of each sample (i.e., the preprocessed focal domain data); Indicates the current parameter The network predicts the output. To prevent the model from overfitting to noise on the training data and thus reducing its generalization ability, a regularization term is introduced into the objective function. ,in The L2 norm (k=2) is usually chosen. is the regularization coefficient, used to control the strength of the penalty for model complexity.

[0046] In sub-step S43, a hybrid loss function that conforms to the physical characteristics of seismic data is designed. A single loss function is insufficient to simultaneously ensure the amplitude fidelity and structural continuity of the seismic data. Therefore, this embodiment constructs a hybrid loss function that combines L1 norm loss and structural similarity loss. This hybrid loss function is defined by the formula: ; In the formula, It refers to the prediction data generated by the network; Refers to the corresponding actual label data; This represents the total number of samples within a single data block; Indicates the first The absolute error of each sample point at a given location. (Symbol) It is a non-negative balancing weight factor used to adjust the proportion of structural similarity constraints in the total loss. In actual training, The value is usually set between 0.1 and 1.0, and the specific value is dynamically adjusted according to the convergence curve of the validation set.

[0047] In substep S44, the structural similarity index (SSIM) is calculated to quantify the similarity between the reconstructed seismic profile and the original profile in terms of brightness, contrast, and structure. The calculation is performed based on the following formula: ; In the formula, and These represent the predicted data. and tag data The local mean; and They represent their standard deviations respectively; The covariance of the two variables is represented by a minimum constant. To ensure the stability of numerical calculations and prevent the denominator from becoming zero due to flat data in local regions (zero variance), a minimum constant is introduced. and These two constants are typically based on the dynamic range of the data. Determine, for example, take , .

[0048] In substep S45, gradient-based parameter updates are performed. An adaptive moment estimator optimizer is used to update the network parameters based on the calculated gradient of the mixture loss function. Perform backpropagation updates.

[0049] In step S5, the trained model is used to predict the test data, and weighted smoothing is applied to generate the complete focal region prediction data; please refer to the appendix. Figure 3 The area within the lower red dashed box represents the region where, during the testing phase, network parameters are kept constant, and complete focus domain data is displayed. Input to the trained network In the process, inference operations are performed, and the complete interpolation result is obtained directly through mapping. .

[0050] Specifically, the implementation of step S5 includes the following logical process: In sub-step S51, model inference under the consistent block-based strategy is performed. The seismic test data to be reconstructed, containing missing data (i.e., preprocessed focal domain data), is input into the trained deep neural network model. To ensure strict consistency between the scale space of feature extraction and the training phase, the test data input is not a single, complete profile fed in at once, but rather segmented using the same sliding window size and step parameters as when constructing the training set. This strategy ensures that the local texture features processed by the network during inference match the effective receptive field learned during training, thereby guaranteeing the effectiveness of feature extraction. The neural network performs forward propagation calculations on each input focal domain data block and outputs the corresponding predicted data block, which physically corresponds to the originally missing near-offset wavefield data obtained from the inversion of focal domain features.

[0051] In sub-step S52, a Gaussian kernel-based weighted smoothing process for overlapping regions is implemented. Because convolutional neural networks are limited by zero-fill boundary conditions and edge decay of the effective receptive field when processing image edges, the prediction accuracy of data block edge regions is usually lower than that of the central region. If the predicted blocks are directly hard-stitched, significant numerical jumps, i.e., block effects, will occur at the stitching boundaries, which will severely disrupt the lateral continuity of the seismic phase axis. Therefore, this embodiment introduces a distance-related weighted fusion mechanism during overlapping region recombination. Specifically, a two-dimensional weight matrix with the same size as the predicted data block is constructed, and the numerical distribution of this matrix follows a two-dimensional Gaussian distribution. Its physical characteristics are: the weight value at the center of the matrix is ​​the largest (close to 1), and the weight value decays exponentially (approaching 0) with increasing distance from the center.

[0052] In sub-step S53, stitching and inverse normalization are performed. During the stitching process, for each pixel on the full-size profile, its final value is obtained by weighted averaging of the predicted values ​​at that point from all the prediction blocks covering that point. The calculation logic is as follows: A full-size accumulated signal buffer matrix and an accumulated weight buffer matrix are established respectively; each prediction block is multiplied by its corresponding Gaussian weight matrix and then accumulated into the signal buffer matrix, while the weight matrix is ​​accumulated into the weight buffer matrix. Finally, a normalization operation is performed, dividing each element of the signal buffer matrix by the corresponding element of the weight buffer matrix. To prevent the risk of numerical overflow due to the denominator approaching zero caused by the weights of some edge points being too small, a very small positive number is added to the denominator of the division operation. (e.g., 10) -8 Subsequently, the original data statistical parameters (such as maximum absolute amplitude or root mean square amplitude) recorded in step S2 are used to perform an inverse transformation to restore its original physical dynamic range, thereby obtaining reconstructed data consistent with the energy level of the actual earthquake record.

[0053] In step S6, based on the original observation geometry, the reconstructed seismic data is filled into the missing positions of the original seismic data, realizing full-field data interpolation and obtaining a complete seismic data profile. This embodiment achieves the physical fusion of predicted and measured data.

[0054] Specifically, the implementation of step S6 includes the following logical process: In sub-step S61, precise data backfilling based on the observation mask is performed. This step aims to ensure that the final output data volume retains the original high-fidelity observation information while only filling in unsampled areas. Before performing the merging operation, the dimensionality consistency of the two sets of data to be merged is first checked to ensure that the predicted data volume after inverse transformation is completely aligned with the original acquisition records in terms of the number of temporal sampling points and spatial channels. Subsequently, this embodiment utilizes the sampling operator matrix constructed in the data preprocessing stage. (i.e., a binary mask that records the actual physical location of the detector) is used as a backfill criterion.

[0055] The specific fusion logic is as follows: traverse every spatial coordinate point of the seismic data volume. Query the value at the corresponding position in the sampling operator matrix. If This indicates the presence of actual geophone records at that location. The system strictly preserves the original acquired seismic values ​​and does not use network predictions for overlay, thus avoiding the loss of high-frequency details or amplitude fidelity of the true signal due to the smoothing effect of the deep learning model; if This indicates that the location belongs to a near-offset missing region or a bad channel region. The system reads the reconstructed data generated in step S5, which has been Gaussian-weighted and smoothed, and writes it to this coordinate location. Through this point-by-point discrimination and filling mechanism based on physical coordinates, seamless stitching of virtual source records and real source records in terms of wavefield energy and phase is achieved.

[0056] In sub-step S62, a complete full-wavefield seismic profile is generated and output. After the above backfilling operation, the near-offset information that was missing due to the limitations of the acquisition aperture is effectively supplemented, forming a complete two-dimensional or three-dimensional data volume that is continuously distributed in both the time and spatial domains. This complete data volume contains first arrival waves, reflected waves, and multiple wave information recovered after focusing transformation from shallow layers. As a preferred embodiment, this embodiment encapsulates the reconstructed data volume according to the seismic data exchange format of the petroleum industry standard (such as SEGY format) and writes it to the storage medium. The output data has complete wavefield characteristics and can be directly used as standard input data for subsequent processing procedures (such as surface multiple adaptive subtraction SRME, pre-stack depth migration PSDM, or full waveform inversion FWI).

[0057] In summary, the technical solution of this embodiment establishes an efficient and robust seismic data loss recovery mechanism by combining physical focusing transformation with data-driven deep learning technology. Its technical effects are mainly reflected in the following aspects: First, this invention achieves effective utilization of multiple illumination information. It overcomes the limitation of traditional processing methods that treat multiples as coherent noise and remove them. Based on the coupling relationship between multiples and primary waves along their physical propagation paths, a focusing transformation is used to map multiples into high-energy features in the focused domain. By learning their kinematic and dynamic correlations through a deep neural network, implicit subsurface medium reflection features can be extracted from the background field containing multiples, thereby achieving high-precision signal reconstruction in areas lacking near-offset illumination.

[0058] Secondly, it constructs a feature space with clear physical meaning, enhancing the network's generalization ability. Unlike direct end-to-end image inpainting in the spatiotemporal domain, this invention generates focused domain data rich in underground geological structure information through focused transformation, providing deep neural networks with input features that have clear physical orientation. This preprocessing method significantly reduces the complexity of the network learning nonlinear mappings, enabling the model to maintain high prediction accuracy and structural continuity even when facing complex laterally variable geological structures such as faults and pinch-outs.

[0059] Third, it reduces the reliance on high-precision velocity models. This method falls under the category of data-driven reconstruction. While ensuring the quality of data reconstruction, it replaces the complex iterative solution process in traditional wave equation inversion with convolutional filtering techniques from deep learning. This eliminates the need for high-precision explicit velocity models or complex time-frequency transformation operators in the data reconstruction process, avoiding reconstruction errors caused by difficulties or inaccuracies in velocity model establishment, thereby effectively improving the efficiency and adaptability of data processing in seismic exploration.

[0060] Specific application example: Deepwater towed cable seismic data reconstruction To more clearly illustrate the technical details and practical application effects of the present invention, the following detailed explanation is provided in conjunction with deep-sea towed cable collection scenarios and experimental verification data.

[0061] Application scenarios and data foundation: This embodiment selects two-dimensional towed seismic acquisition data from a deep-water work area (such as a deep-water basin in the South China Sea) as the application object. In this observation system, the seismic source is excited by an air gun, the towed cable is 6000 meters long, and the geophone channel spacing is 12.5 meters.

[0062] During actual data acquisition, due to the strong wake noise generated by the tugboat and the need to maintain a safe physical distance between the seismic source and the first geophone, significant near-offset gaps inevitably exist in the original seismic records. Specifically, the gaps range from 0 to 150 meters in minimum shot-receiver distance. This gap directly leads to insufficient imaging illumination of shallow seabed reflections and key geological interfaces, severely affecting subsequent AVO analysis and high-precision imaging.

[0063] Specific implementation steps: This invention utilizes joint focusing transformation and deep learning methods to reconstruct the aforementioned missing data. The specific steps are as follows: S1: Preprocessing and Focus Transformation The system first reads the single-shot record and uses FK filtering to separate the primary wave with strong energy characteristics as the focusing operator. Meanwhile, to mathematically describe the insufficient spatial sampling during the acquisition process, a sampling operator matrix is ​​constructed. In this embodiment, the grid within the range of 0-150 meters of the shot-receiver distance is marked as 0 (representing missing), and the remaining range is marked as 1.

[0064] Subsequently, according to the aforementioned formula This maps the raw data, containing multiple waves, to the focused domain. Within this focused domain data, (Focus point) corresponds to a high-energy pulse with zero offset at time zero, while The background information carries implicit, near-offset full-wave field information that needs to be extracted.

[0065] S2-S3: Feature Processing and Network Mapping Spatiotemporal stripping is performed on the focused domain data to remove strong-energy focal points. Then, a region with an offset of 200 to 2000 meters (i.e., the data-complete region) from the full wavefield data is selected as the training set to construct a pairwise mapping relationship between the focused domain background term and the original wavefield. Considering the highly developed multiples in deep water, the network will focus on learning the nonlinear correspondence between the morphological characteristics of the multiples in the focused domain and the original primary wave.

[0066] S4-S5: Model Training and Inference A ResUNet deep neural network model was constructed and trained using a hybrid loss function combining L1 norm and structural similarity (SSIM). After training, complete focal domain data containing the 0-150 meter missing region was input into the model. Based on the learned physical laws, the model outputs the predicted wavefield for the missing region.

[0067] S6: Data Backfilling and Output After eliminating the boundary effects that may arise from block prediction by using Gaussian weighted smoothing, the sampling operator matrix is ​​then used. Only the predicted data within the 0-150 meter range is backfilled into the original record to generate a complete SEGY data volume. This data volume completes the near-offset information, providing a complete wavefield basis for subsequent processing.

[0068] Experimental verification and effect comparison: To verify the effectiveness of the method of the present invention, a simulation experiment based on the Pluto sound wave velocity model was designed and compared with the prior art.

[0069] Please see the appendix Figure 9To further test the adaptability of the method under extreme geological conditions, this embodiment also introduces a complex velocity model with strong lateral velocity variation for verification. As shown in the figure, the model has a giant, irregularly shaped, high-velocity rock body at its center (dark blue area in the figure), whose velocity value is significantly higher than that of the surrounding sedimentary strata (red and light-colored layered areas in the figure). This dramatic velocity difference will produce strong refraction and multiple reflections during seismic wave propagation. Steep faults and fold structures are developed on the flanks and below the rock body; these structures are often blind spots for conventional interpolation methods. Based on this velocity model, forward modeling using the acoustic equations was performed to generate complex wavefield data containing rich multiple wave information, serving as a benchmark for verifying the reconstruction accuracy of the deep learning network of this invention.

[0070] Please see the appendix Figure 10 In order to incorporate physical constraints into deep learning networks, this embodiment first processes the simulated raw full-wave field data (such as...) Figure 10 (a) Preprocessing is performed. Based on the actual observation system's geometry, near-offset data is removed to simulate acquisition gaps (such as...). Figure 10 (b) shows). Subsequently, the primary wave in-phase axis with strong energy and clear morphology is separated using FK filtering technology (e.g., Figure 10 (c) shows that it is defined as the focusing operator. Finally, this operator is used to perform a focusing transformation operation on the missing data to obtain the focusing domain data (such as...). Figure 10 (d) shown). Comparison Figure 10 (b) and Figure 10 (d) It can be seen that the complex hyperbolic wave field in the spatiotemporal domain is mapped to specific texture features in the focal domain. The X-shaped energy at the top intersection corresponds to the focal point, while the complex background texture extending below it contains multiple effective information of the wave to be mined, providing high-dimensional physical feature input for the subsequent nonlinear mapping learning of the network.

[0071] Please see the appendix Figure 11 The figure shows the final reconstructed profiles of the three methods described above in areas where complex geological models are missing.

[0072] like Figure 11 As shown in (a), the seismic profile reconstructed by the method of the present invention exhibits extremely high continuity in the near offset region (i.e., near the horizontal coordinate 0 meters in the figure), the curvature of the shallow strong energy phase axis is restored naturally, the texture details of the deep weak signal are clear, and no obvious splicing traces are found. contrast Figure 11 (b) The reconstruction results based on the projection of convex sets (POCS) show obvious step-like false frequencies at near offsets, and the continuity of the phase axis is disrupted. This indicates that traditional interpolation algorithms are difficult to accurately predict rapidly changing wavefront curvature when facing non-stationary wave fields caused by complex structures. contrast Figure 11 (c) Although the reconstruction results based on conventional supervised learning methods restore the overall shape of the wave field, there is insufficient energy focusing at the shallow phase axis intersection, and some weak reflection signals are masked by smoothing noise, resulting in a lower signal-to-noise ratio than the method of this invention.

[0073] Please see the appendix Figure 12 To verify the fidelity of the reconstructed data in amplitude and phase from a microscopic perspective, this embodiment further extracted single-channel seismic records at different offset positions (300m, 1760m, 1790m, and 2000m) for detailed comparison. In the figure, the black solid line represents the original true seismic record, the red solid line represents the reconstructed record of each method, and the blue elliptical area marks the key comparison positions with drastic wavefield changes.

[0074] contrast Figure 12 (b) It can be seen that the reconstruction results based on the projection of convex sets (POCS) show obvious misalignment between the red curve and the black curve at the peak marked by the blue ellipse, which is manifested as obvious phase lag and amplitude energy attenuation, indicating that it has limitations in dealing with rapid wave field changes caused by complex structures. contrast Figure 12 (c) Although the reconstructed waveform based on conventional supervised learning methods has improved in phase, there is still amplitude fitting error at the peak of some strong reflection interfaces, which cannot completely restore the dynamic characteristics of the original signal. And observation Figure 12 (a) The red waveform reconstructed using the method of this invention closely matches the black true waveform at all selected offset locations, especially in the complex wavegroup region marked by the blue ellipse. This strongly demonstrates that the nonlinear mapping relationship of the focusing domain learned by the deep learning network in this invention can accurately capture and recover the details of seismic waveforms under complex geological conditions, achieving high-precision full-wavefield reconstruction.

[0075] The experimenter artificially removed data with near offset distances (0-200 meters) as the test set and compared the method of the present invention with the interpolation method based on convex set projection (POCS) and conventional supervised learning interpolation methods. Please see the appendix Figure 7 , Figure 7 The image visually demonstrates a local magnified comparison of the reconstruction effects of the three methods in the missing area of ​​the flat model (i.e., the area from -200 meters to 200 meters on the horizontal coordinate in the figure).

[0076] like Figure 7As shown in (a), the earthquake phase axis reconstructed using the method of the present invention exhibits a smooth and continuous hyperbolic shape when crossing the zero offset position of the center. The wave group characteristics are clear and the amplitude energy is consistent, accurately restoring the dynamic characteristics of the reflection of the underground medium. contrast Figure 7 (b) The reconstruction results based on the convex set projection method show obvious waveform distortion at the missing center. The phase axis exhibits non-physical broken line or sharp point protrusion at zero offset, indicating that the method cannot accurately fit the curvature change of the wave field when dealing with large-scale gaps. contrast Figure 7 (c) Although the reconstruction results based on conventional supervised learning methods are closer to the real wave field in terms of morphology, there are still blurring effects at the edge of the waveform, and the focusing degree of the phase axis is slightly inferior to that of the method of the present invention. Please see the appendix Figure 8 To further quantitatively analyze the amplitude fidelity and phase consistency of the reconstructed data, near-offset traces (30m) within the missing range and far-offset traces (520m, 540m, 800m) within the known range were extracted for waveform detail comparison. In the figure, the black solid line represents the original true earthquake record, and the red solid line represents the reconstructed or predicted record of each method.

[0077] Observing the left image (convex set projection method) and the middle image (conventional supervised learning method), it can be found that in the reconstruction of missing data at the 30m offset, there is a significant deviation between the red solid line and the black solid line. Specifically, the alignment of the wave crests is poor (phase lag) and the energy attenuation of the antinodes is reduced (amplitude distortion). This indicates that traditional methods are difficult to accurately recover the dynamic characteristics of the wave field when there is a lack of short-path constraints.

[0078] In contrast, observing the right figure (the method of this invention), whether at the missing near offset (30m) or the retained far offset (520m, 540m, 800m), the red reconstructed waveform and the black true waveform achieve a high degree of overlap, with almost no visible residual. This confirms that by introducing physical constraints on the focusing domain, this invention not only eliminates interpolation aliasing but also maintains extremely high signal fidelity at the microscopic waveform scale, verifying the robustness of the method.

[0079] Waveform consistency analysis: Please see the appendix Figure 13 This figure shows the comparison results of single-channel seismic waveforms at the near offset missing location (Offset=100m). The meaning and analysis of each element in the figure are as follows: Coordinate axis explanation: The horizontal axis in the figure represents time, that is, the time for the seismic wave to travel, in seconds (s), showing the reflection time from shallow to deep layers; the vertical axis represents the normalized amplitude, which indicates the energy strength of the seismic signal, with a value range between -1 and 1.

[0080] Actual acquired data (solid black line): This curve represents the actual seismic waveform that should have been acquired under ideal conditions, without any missing data. As a standard for comparison, its waveform is complete and its phase is accurate.

[0081] Existing technology (black dashed line): This curve represents the waveform recovered using traditional interpolation methods (such as supervised learning). By comparison, it is clear that the peak of this curve does not coincide with the black solid line, indicating a significant time travel deviation (phase lag) in complex wave fields (such as 1.5s-2.0s); and the peak height is relatively low, indicating amplitude attenuation and inaccurate energy recovery.

[0082] The method of this invention (gray line with a circle): This curve represents the waveform reconstructed using the method of this patent. As can be seen from the figure, it almost completely overlaps with the black solid line, and the center of the circle mark falls precisely on the actual data curve. This shows that the physical constraints introduced by the focusing transformation of this invention effectively guide the neural network, not only recovering the phase information very accurately, but also far surpassing traditional methods in amplitude fidelity.

[0083] Comparative analysis of quantitative indicators: Please see the appendix Figure 14 This figure uses a bar chart to compare the signal-to-noise ratio (SNR) and structural similarity (SSIM) metrics of the three methods after reconstruction. The meaning and analysis of each element in the figure are as follows: Horizontal axis (comparison method): The chart lists three different technical solutions: convex set projection method (basic interpolation), conventional supervised learning method, and the method of this invention (joint focusing transformation and deep learning).

[0084] The left vertical axis and the dark gray bars (signal-to-noise ratio): correspond to the height of the dark gray bars in the figure, and the unit is decibels (dB). The higher the signal-to-noise ratio value, the higher the proportion of effective signal in the reconstructed data and the lower the background noise. As shown in the figure, the signal-to-noise ratio of the method of this invention reaches 18.5dB, which is significantly higher than that of the convex set projection method (12.3dB) and the supervised learning method (14.1dB), proving that the signal purity generated by this method is the highest.

[0085] The right vertical axis and light gray bars (structural similarity): correspond to the height of the light gray bars in the figure, with values ​​ranging from 0 to 1. The closer the value is to 1, the more similar the geological texture structure of the reconstructed image is to the actual situation. As shown in the figure, the SSIM score of the method of this invention is as high as 0.94, which is closest to 1, indicating that the reconstructed seismic phase axis has good continuity and no obvious artificial artifacts or breaks appear; while the convex set projection method scores only 0.82, indicating that it has the problem of local structural blurring when dealing with large gaps.

[0086] Columnar Height Overview: Columnar height directly reflects the performance of various methods on specific indicators. The method of this invention ranks highest in both sets of columns, strongly demonstrating its superior overall performance in reconstructing seismic data from large gaps and complex structures.

[0087] in conclusion: In summary, the experimental results show that the joint focusing transformation and deep learning method proposed in this invention can recover missing near-offset seismic data more accurately than traditional mathematical transformation methods in deep-water towed cable scenarios with large gaps and complex geological structures, achieving high-fidelity physical reconstruction.

Claims

1. A method for reconstructing multi-wave driven seismic data based on deep learning, characterized in that, Includes the following steps: A focusing operator is constructed using the acquired data domain seismic shot gather, and a focusing transformation is performed to obtain the corresponding focusing domain data; Remove focal points and noise regions from the focal domain data, and construct a mapping relationship between the data domain and the focal domain based on a preset missing pattern; A sliding window is used to extract data blocks from the focal domain data and the data domain seismic shot gather corresponding to the seismic data to construct training samples. The training samples include the missing focal domain data as input and the corresponding missing original seismic data as labels. After data augmentation, a training set is formed. Construct a deep neural network model and train it using the training set to establish a nonlinear mapping relationship between the features of the focus domain and the signals of the data domain; The complete focal domain data to be reconstructed is input into the deep neural network model, block prediction is performed to obtain prediction data blocks, and Gaussian weighted smoothing is performed on the overlapping prediction data blocks to obtain the reconstructed missing seismic data. The reconstructed missing seismic data is filled into the missing positions of the original data domain seismic shot gather to obtain a complete seismic data profile.

2. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The step of constructing a focusing operator using the acquired data domain seismic shot gather and performing a focusing transformation operation on the seismic shot gather includes: The seismic shot gather is preprocessed to separate the primary wave field with main energy characteristics as the focusing operator; A sampling operator matrix is ​​constructed, which is used to identify the spatial existence state of the receiver point relative to the shot point in the seismic record; The focusing operator is used to perform correlation operations with the seismic shot gather containing multiple waves, and the focusing domain data is generated by combining the sampling operator matrix. The focusing domain data includes the focal point located at the origin and a background term containing wavefield information.

3. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The step of removing focal points and noise regions from the focal domain data includes: A spatiotemporal cutoff window is defined, which is a parabolic function that moves with the shot point, and the curvature of the parabola is determined by the phase axis of the seismic record. The sample values ​​of the focal domain data within the spatiotemporal cutoff window are set to zero or smoothed and attenuated to cut off the focal point and the surrounding noise region to obtain the focal domain background item. The background term of the focal region is normalized, and the processed background term of the focal region is set to the same near offset missing value as the original data and used as the input feature of the deep neural network.

4. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The step of constructing training samples by extracting data blocks from the focused domain data and the data domain seismic shot gather using a sliding window includes: Set the size of the two-dimensional sliding window and the movement step size in the time and space directions, wherein the movement step size is smaller than the size of the corresponding dimension of the two-dimensional sliding window; Simultaneous scanning and cropping are performed on the profile corresponding to the focused domain data and the profile corresponding to the data domain seismic shot gather. Samples are extracted only from the effective data area actually acquired by the seismic gun geophone in the data domain.

5. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The data augmentation operation performed on the training samples includes: Based on the physical characteristics of the seismic wave field, a spatial trace order flipping operation is performed on the training samples; Perform phase reversal and random gain adjustment operations on the training samples; Band-limited Gaussian noise is superimposed on the training samples.

6. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The construction of the deep neural network model includes: Construct a ResUNet network structure that includes an encoder, decoder, and skip connections; The encoder is used to extract multi-scale energy distribution features of the focal domain data through cascaded convolution operations; The decoder is used to perform upsampling operations to restore the spatiotemporal resolution of seismic data; The skip connection is used to pass features from the encoder to the corresponding layer of the decoder.

7. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 1, characterized in that, The step of training the model using the training set includes: A hybrid loss function is constructed, which is composed of a weighted sum of an L1 norm loss term and a structural similarity loss term; The structural similarity loss term is calculated based on the local mean, standard deviation, and covariance of the predicted data and the label data; The parameters of the deep neural network model are updated using an adaptive moment estimation optimizer based on the gradient of the hybrid loss function.

8. The method for reconstructing multi-wave driven seismic data based on deep learning according to claim 3, characterized in that, The step of performing Gaussian-weighted smoothing on the overlapping predicted data blocks includes: Construct a two-dimensional weight matrix with the same size as the predicted data block. The numerical distribution of the two-dimensional weight matrix follows a two-dimensional Gaussian distribution, with the weight value being the largest at the center and decreasing exponentially towards the edge. Establish a full-size signal buffer matrix and weight buffer matrix; Each predicted data block is multiplied by the two-dimensional weight matrix and then accumulated into the signal cache matrix, while the two-dimensional weight matrix is ​​accumulated into the weight cache matrix. Divide the elements in the signal buffer matrix by the elements at the corresponding positions in the weight buffer matrix, and perform the inverse normalization operation corresponding to the normalization process.

9. A method for reconstructing multi-wave driven seismic data based on deep learning according to claim 2, characterized in that, The construction of the sampling operator matrix includes: Set a threshold for the missing range on one side of the firing point; The absolute value of the difference between the spatial coordinate index of the receiver point and the spatial coordinate index of the shot point is used as the shot-receiver distance. When the shot-receiver distance is greater than the missing range threshold, the corresponding element of the sampling operator matrix is ​​marked as valid; When the shot-receiver distance is less than or equal to the missing range threshold, the corresponding element of the sampling operator matrix is ​​marked as missing.

10. A method for reconstructing multi-wave driven seismic data based on deep learning according to claim 9, characterized in that, The step of filling the missing locations in the original seismic data with the reconstructed missing seismic data includes: Verify the dimensionality consistency between the reconstructed missing seismic data and the seismic shot gather in the data domain; Traverse the spatial coordinates of the seismic shot gather in the data domain and query the values ​​at the corresponding positions in the sampling operator matrix; If the corresponding position in the sampling operator matrix is ​​marked as valid, the original values ​​in the data domain seismic shot set are retained. If the corresponding position in the sampling operator matrix is ​​marked as missing, the value of the missing seismic data after reconstruction is written in.

Images