Imaging domain least square migration method and system based on psf prednet

Abstract

This invention discloses an imaging domain least-squares migration imaging method and system based on PsfPredNet, comprising: acquiring observational seismic shot gather data of the target area; obtaining the migration velocity field through conventional travel-time tomography inversion; and obtaining the reverse-time migration imaging result using a conventional reverse-time migration method; based on the observation system and the migration velocity field, performing conventional inverse migration and reverse-time migration on the reflectivity model of discrete scattering points to obtain the point spread function of the discrete points; constructing a training dataset for a prediction network using the migration velocity field, the reverse-time migration imaging result, and the point spread function of the discrete points; training a PsfPredNet based on an encoder-decoder structure based on the training dataset; and embedding the trained PsfPredNet into an imaging domain least-squares reverse-time migration inversion method with sparse constraints and total variational constraints to generate high-resolution seismic imaging results. This invention can improve the resolution and reliability of seismic imaging and is of great significance for the exploration of complex oil and gas reservoirs.

Description

Technical Field

[0001] This invention relates to the field of geophysical exploration technology, and in particular to an imaging domain least squares migration imaging method and system based on PsfPredNet. Background Technology

[0002] Seismic imaging technology is a crucial tool for acquiring information on subsurface geological structures and identifying the distribution of oil and gas reservoirs. Among traditional seismic imaging methods, Reverse Time Migration (RTM) has become one of the most mainstream imaging techniques due to its strong adaptability to complex structures and high-precision imaging capabilities for steeply dipping strata. However, the imaging operator used in RTM is essentially an adjoint operator of the forward modeling operator. Its imaging result is not the true stratigraphic reflectance, but rather the convolution of reflectance with the Hessian matrix. The Hessian matrix characterizes imaging crosstalk and blurring effects between spatial points. This convolution effect inevitably causes problems such as imaging amplitude distortion, resolution reduction, and imaging artifacts, thus severely restricting subsequent structural interpretation and quantitative reservoir characterization.

[0003] To correct the blurring effect of the Hessian matrix and improve imaging accuracy, Least-Square Reverse Time Migration (LSRTTM) has been proposed and extensively studied. LSRTM can be implemented through model-driven or data-driven approaches. Model-driven LSRTM formulates imaging as a least-squares inversion problem, iteratively recovering the optimal reflectivity model while suppressing imaging noise and compensating for uneven amplitudes. In this type of method, imaging-domain LSRTM approximates the Hessian matrix using the point spread function (PSF), avoiding repeated wavefield simulations during iteration and achieving high computational efficiency. However, the main bottleneck of imaging-domain LSRTM lies in the explicit calculation and storage of the PSF. Traditional methods typically calculate the PSF at sparse spatial locations and, based on the assumption of its gradual spatial variation, use bilinear or trilinear interpolation to obtain the remaining PSFs. This method achieves good results in regions with gentle velocity fields and simple geological structures because the spatial distribution of the PSF has strong continuity, resulting in smaller interpolation errors. However, in complex tectonic zones (such as areas containing faults, salt domes, and thrust-nappes), the velocity field exhibits drastic lateral variations, and wave field propagation is severely distorted, leading to significant spatial nonstationarity of the PSF. In such cases, traditional interpolation methods struggle to accurately characterize the rapid changes in PSF, thus affecting the imaging quality of LSRTM. In recent years, deep learning has demonstrated powerful nonlinear mapping and feature representation capabilities. In seismic imaging, deep learning has been used to generate high-quality images from conventional migration imaging results; this type of method is also classified as data-driven LSRTM. In this framework, neural networks learn an approximate representation of the inverse Hessian matrix through a training dataset, which can be considered an alternative model to imaging-domain LSRTM. However, in actual seismic exploration, it is impossible to obtain real stratigraphic reflectance for labeling. Such models typically rely on training based on synthetic datasets, which faces significant limitations in practical applications: synthetic data and measured data often have difficulty matching features, resulting in insufficient generalization ability of the trained model in real-world data processing.

[0004] In summary, the interpolation methods used by existing imaging domain LSRTM methods in calculating PSF cannot adapt to the nonlinear spatial distribution characteristics of complex structural regions, resulting in insufficient LSRTM imaging accuracy. Although data-driven methods can avoid explicit PSF calculation, they lack effective labels for actual data and rely heavily on synthetic data for training, resulting in poor model generalization ability and difficulty in stably applying them to actual work areas. Summary of the Invention

[0005] The purpose of this invention is to provide an imaging domain least squares migration imaging method and system based on PsfPredNet to solve the problem of insufficient interpolation accuracy of traditional interpolation methods for PSF interpolation in drastically spatially varied environments.

[0006] To achieve the above objectives, the present invention adopts the following technical solution:

[0007] In a first aspect, the present invention provides an imaging domain least-squares migration imaging method based on PsfPredNet, comprising:

[0008] The seismic shot gather data of the target work area were acquired, the migration velocity field was obtained through conventional travel-time tomography inversion, and the reverse-time migration imaging results were obtained using conventional reverse-time migration methods.

[0009] Based on the observation system and the migration velocity field, the reflectivity model of discrete scattering points is subjected to conventional inverse migration and inverse time migration to obtain the point spread function of discrete points.

[0010] Using the offset velocity field, reverse time offset imaging results, and point spread function of discrete points, a training dataset for the prediction network is constructed, and PsfPredNet based on the encoder-decoder structure is trained based on the training dataset.

[0011] The trained PsfPredNet is embedded into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variational constraints. The reverse time migration imaging results are iteratively optimized by the alternating direction multiplier method to generate high-resolution seismic imaging results.

[0012] Furthermore, the step of constructing a training dataset for the prediction network using the offset velocity field, reverse time migration imaging results, and point spread function of discrete points, and training PsfPredNet based on the encoder-decoder structure on the training dataset includes:

[0013] The input samples of the training dataset consist of a set of point spread functions of neighborhood scattering points, the reverse time migration imaging results of local sub-regions, and the migration velocity field. The output label is the point spread function of the target point. The set of point spread functions of neighborhood scattering points is weighted according to Euclidean distance.

[0014] Furthermore, training PsfPredNet based on the encoder-decoder structure using the training dataset includes:

[0015] The prediction network is a PsfPredNet neural network convolutional neural network built on the ResUNet residual encoder-decoder network architecture; its training process includes:

[0016] Initialize network parameters;

[0017] Input the training dataset into the network model;

[0018] The loss function is the weighted sum of the mean squared error loss between the predicted point spread function and the true point spread function and the frequency domain loss. The network parameters are iteratively optimized by gradient descent until the loss function converges or the set maximum number of iterations is reached, thus completing the model training.

[0019] Furthermore, during model training, the loss function consists of a time-domain mean absolute error term and a frequency-domain constraint term, defined as follows:

[0020]

[0021] in, The first real PSF One element, These are the predicted values ​​from the network model. Represents the Fast Fourier Transform. The total number of elements. and These are the weighting coefficients for the time-domain loss and the frequency-domain loss, respectively.

[0022] Furthermore, batch training is employed during the training process, with each training batch containing 2-8 samples and an initial learning rate of 10. -3 ~10 -4 .

[0023] Furthermore, the step of embedding the trained PsfPredNet into an imaging domain least-squares reverse time migration inversion method with sparse constraints and total variational constraints, and iteratively optimizing the reverse time migration imaging results through the alternating direction multiplier method to generate high-resolution seismic imaging results includes:

[0024] The least-squares reverse time migration (RTM) process integrated with PsfPredNet is as follows: Input the RTM imaging results, migration velocity field, and discrete point spread function; set the maximum number of iterations and the regularization parameter of the objective function; use the PsfPredNet network to predict the PSF of the entire model space, and construct a Hessian matrix with the PSF of each point as a column in the matrix; iteratively update the RTM imaging results using the alternating direction multiplier algorithm; when the number of iterations reaches the preset maximum number of iterations, output the high-resolution imaging result of the current iteration number.

[0025] The objective function is as follows:

[0026]

[0027] In the formula, For offset images, This represents the approximate Hessian matrix formed by the PSFs predicted by PsfPredNet. For PsfPredNet network parameters, Distribution of underground reflectance coefficients , and For regularization parameters, For L1 norm constraint terms, For various same-sex TV regularization terms;

[0028] in,

[0029]

[0030] In the formula, As auxiliary variables, where , , and These are the first-order difference operators along the horizontal and vertical directions, respectively. This represents the number of grid points in the reflectivity model.

[0031] Furthermore, the method of iteratively updating the RTM imaging results using the alternating direction multiplier method is described in the following specific way:

[0032] Construct the augmented Lagrangian function for the constrained least squares problem, which takes the following form:

[0033]

[0034] In the formula, For Lagrange multipliers, The augmented quadratic term serves as a penalty coefficient; it penalizes parts that violate equality rules, thus accelerating the convergence of the alternating direction multiplier algorithm. Perform a simple variable substitution A scaled version of the alternating direction multiplier algorithm is derived:

[0035]

[0036] During the ADMM algorithm update process, the first... The reflectance image at the next iteration is represented as follows:

[0037]

[0038] This problem is a least squares subproblem with L1 norm constraints, which is solved using a fast iterative soft thresholding algorithm.

[0039] After that, fix and For auxiliary variables Update:

[0040]

[0041] This subproblem can be solved by applying the generalized soft thresholding operator element by element to the auxiliary variable, resulting in a closed-form solution.

[0042]

[0043]

[0044] in, .

[0045] Finally, update the scaling dual variable as follows:

[0046]

[0047] Secondly, the present invention provides an imaging domain least-squares migration imaging system based on PsfPredNet, comprising:

[0048] The data acquisition module is used to acquire the observed seismic shot gather data of the target work area, obtain the migration velocity field through conventional travel-time tomography inversion, and obtain the reverse-time migration imaging results using conventional reverse-time migration methods.

[0049] The point spread function calculation module is used to perform conventional inverse migration and inverse time migration on the reflectivity model of discrete scattering points based on the observation system and the migration velocity field, so as to obtain the point spread function of discrete points.

[0050] The training module is used to construct a training dataset for the prediction network using the offset velocity field, inverse time offset imaging results, and point spread function of discrete points. Based on the training dataset, PsfPredNet with an encoder-decoder structure is trained.

[0051] The inversion module is used to embed the trained PsfPredNet into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variational constraints. It iteratively optimizes the reverse time migration imaging results through the alternating direction multiplier method to generate high-resolution seismic imaging results.

[0052] Thirdly, the present invention provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the steps of the PsfPredNet-based imaging domain least squares migration imaging method.

[0053] Fourthly, the present invention provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the steps of the PsfPredNet-based least squares offset imaging method for the imaging domain.

[0054] Compared with the prior art, the present invention has the following technical effects:

[0055] This invention utilizes deep learning to mine the nonlinear mapping between the PSF of neighboring scattering points, RTM imaging results, the migration velocity field, and the PSF of the target point. It can adaptively learn the drastic spatial variation characteristics of PSF caused by complex velocity fields. Compared to traditional imaging-domain LSRTM methods, this invention significantly reduces PSF prediction errors in complex structural regions such as fractures and salt domes, thereby effectively improving the imaging quality of subsequent LSRTM.

[0056] In this invention, the label data used to train PsfPredNet can be directly obtained from the observation system and migration velocity field of the target work area through the point spread function calculation module. Compared with traditional data-driven methods that rely on synthetic model data for training, it has stronger generalization ability and applicability in practical applications. Furthermore, in terms of the learning task, traditional methods require learning complex Hessian matrix inverse mappings, while PsfPredNet only needs to learn the relatively simple nonlinear relationship between the target point's PSF and its local influencing factors (neighborhood PSF, RTM imaging results, velocity field). This simpler learning objective reduces the requirement for a large amount of training data, making the network easier to train and converge. Attached Figure Description

[0057] To more clearly understand the technical solution of the present invention, the drawings used in the description of the embodiments or prior art will be briefly introduced below. The drawings in the following description are only schematic diagrams of some embodiments of the present invention and do not constitute a limitation of the present invention:

[0058] Figure 1 This is an overall flowchart of the imaging domain least squares migration imaging method based on PsfPredNet of the present invention.

[0059] Figure 2 This is a schematic diagram of the PsfPredNet deep neural network structure used in this invention;

[0060] Figure 3 This refers to the offset velocity field in this embodiment of the invention;

[0061] Figure 4 The results are RTM imaging results from the embodiments of the present invention;

[0062] Figure 5 This refers to the discrete scattering point PSF used for network training in this embodiment of the invention.

[0063] Figure 6 The target point PSF obtained by the method of the present invention in this embodiment of the invention;

[0064] Figure 7The target point PSF is obtained using the traditional bilinear interpolation method in this embodiment of the invention;

[0065] Figure 8 This is the true PSF of the target point obtained by inverse offset and RTM in the embodiments of the present invention;

[0066] Figure 9 The imaging domain LSRTM imaging results based on PsfPredNet in this embodiment of the invention are used.

[0067] Figure 10 This is the conventional imaging domain LSRTM imaging result in an embodiment of the present invention;

[0068] Figure 11 This is a true reflectance image from an embodiment of the present invention. Detailed Implementation

[0069] The present invention will be further described below with reference to the accompanying drawings:

[0070] Example 1, please refer to Figure 1 This invention provides an imaging domain least-squares migration imaging method based on PsfPredNet, comprising:

[0071] The seismic shot gather data of the target work area were acquired, the migration velocity field was obtained through conventional travel-time tomography inversion, and the reverse-time migration imaging results were obtained using conventional reverse-time migration methods.

[0072] Based on the observation system and the migration velocity field, the reflectivity model of discrete scattering points is subjected to conventional inverse migration and inverse time migration to obtain the point spread function of discrete points.

[0073] Using the offset velocity field, reverse time offset imaging results, and point spread function of discrete points, a training dataset for the prediction network is constructed, and PsfPredNet based on the encoder-decoder structure is trained based on the training dataset.

[0074] The trained PsfPredNet is embedded into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variational constraints. The reverse time migration imaging results are iteratively optimized by the alternating direction multiplier method to generate high-resolution seismic imaging results.

[0075] This method first utilizes PsfPredNet to learn the spatially varying characteristics of the PSF (Precipitated Seismic Array), overcoming the insufficient interpolation accuracy of traditional interpolation methods for PSFs with drastic spatial variations, thus achieving high-precision prediction of PSF in complex tectonic regions. Based on this, the trained PsfPredNet is integrated into the imaging domain LSRTM framework to construct a more accurate Hessian matrix. The RTM imaging results are then iteratively updated using an alternating direction multiplier algorithm, ultimately generating high-resolution seismic imaging results.

[0076] Example 2: This invention provides an imaging domain least-squares migration imaging method based on PsfPredNet, comprising:

[0077] Step S1: Obtain the observation seismic shot gather data of the target work area, obtain the migration velocity field through conventional travel-time tomography inversion method, and obtain the RTM imaging results through conventional RTM method;

[0078] Step S2: Based on the observation system and migration velocity field of the observed seismic data in S1, perform conventional inverse migration and RTM calculations on the reflectivity model of discrete scattering points to obtain the PSF of discrete points;

[0079] Step S3: Use the offset velocity field, RTM imaging results and discrete point PSF calculated in S1 and S2 to create the training dataset for PsfPredNet.

[0080] Step S4: Train PsfPredNet using the training dataset from S3.

[0081] Step S5: Integrate the PsfPredNet trained in S4 into the imaging domain LSRTM with sparse constraints and total variational regularization constraints, and iteratively update the RTM imaging results through the alternating direction multiplier algorithm to generate high-resolution seismic imaging results.

[0082] Preferably, step S3 specifically includes the following steps:

[0083] Step a: Normalize the discrete point PSF and RTM imaging results to the [-1,1] interval, normalize the migration velocity field to the [0,1] interval, and remove outliers;

[0084] Step b: Using the scattering points of the known PSF as the target points, extract the PSF set of the neighboring scattering points of the target point, the RTM imaging results of the local sub-region of the target point, and the migration velocity field;

[0085] Step c: Weight the PSF of the neighboring scattering points according to their Euclidean distance from the target point;

[0086] Step d: Perform image enhancement processing such as rotation, flipping, and noise addition, and combine them to form a training dataset with a unified format. At the same time, divide the training set and the validation set in a 7:3 ratio.

[0087] Preferably, in step b, the neighborhood scattering point PSF set refers to the four scattering points with the smallest Euclidean distance centered on the target point (excluding the target point itself). The spatial range of the local sub-region of the RTM imaging result and the migration velocity field is the same as that of the target point PSF, and its center point is the target point.

[0088] Preferably, in step c, the weighting process assigns weights based on the Euclidean distance between each scattering point and the target point. Scattering points closer to the target point are assigned greater weights, while scattering points farther away are assigned smaller weights. The calculation formula is as follows:

[0089]

[0090] in, Indicates the first The weights corresponding to each scattering point This represents the Euclidean distance between the scattering point and the target point. This represents the total number of scattering points in the neighborhood.

[0091] Preferably, in step d, the number of training samples is the same as the number of target points. The input data for each training sample includes: a set of neighborhood scattering point PSFs, the RTM imaging results of the local sub-region, and the offset velocity field, labeled as the target point PSF.

[0092] Preferably, step S4 specifically includes the following steps:

[0093] Step a: Build the PsfPredNet neural network based on the Residual Encoder-Decoder Network (ResUNet) architecture, which has strong feature extraction and multi-source information fusion capabilities, and initialize the convolutional layer parameters using the He normal distribution;

[0094] Step b: Input the training set into the neural network, using the mean square error between the predicted PSF and the true PSF and the frequency domain loss as the loss function, and use the gradient descent method (AdamW optimizer) to iteratively optimize the network parameters.

[0095] Step c: When the iteration termination condition is met, stop training and save the trained model.

[0096] Preferably, in step a, the specific network structure of PsfPredNet is as follows: Figure 1 As shown, the encoder is used to extract deep features from multi-source data, and the decoder is used to map the features to PSF output;

[0097] Preferably, in step b, the loss function is composed of a time-domain average absolute error term and a frequency-domain constraint term, and is defined as follows:

[0098]

[0099] in, The first real PSF One element, These are the predicted values ​​from the network model. Represents the Fast Fourier Transform. The total number of elements. and These are the weighting coefficients for the time-domain loss and the frequency-domain loss, respectively.

[0100] Preferably, in step c, the iteration termination condition is: the loss function value of the validation set no longer decreases, or the number of iterations reaches a set upper limit.

[0101] Preferably, step S5 specifically includes the following steps:

[0102] Step a: Input the RTM imaging results, the migration velocity field, and the PSF of discrete points; set the maximum number of iterations and the regularization parameter of the objective function.

[0103] Step b: Use the PsfPredNet network to predict the PSF of the entire model space, and construct a Hessian matrix with the PSF of each point as a column in the matrix.

[0104] Step c: Iteratively update the RTM imaging results using the alternating direction multiplier algorithm;

[0105] Step d: When the number of iterations reaches the preset maximum number of iterations, output the high-resolution imaging result of the current iteration number.

[0106] Preferably, in step a, the objective function is as follows:

[0107]

[0108] In the formula, For RTM imaging results, This represents the approximate Hessian matrix formed by the PSFs predicted by PsfPredNet. For PsfPredNet network parameters, For underground reflectivity, , and For regularization parameters, For sparse constraint terms, Let be the regularization term for the total variation with isotropic terms, and its formula is:

[0109]

[0110] In the formula, As auxiliary variables, where , , and These are the first-order difference operators along the horizontal and vertical directions, respectively. This represents the number of grid cells in the reflectivity model.

[0111] The objective function can be transformed into the following augmented Lagrangian function, in the following form:

[0112]

[0113] In the formula, For Lagrange multipliers, The penalty coefficient is used to augment the quadratic term.

[0114] Preferably, in step b, since directly predicting the PSF of all target points in the entire imaging space using PsfPredNet is time-consuming, a compromise strategy of "model prediction-interpolation completion" is adopted to predict the PSF of the entire model space, specifically including the following steps:

[0115] Within the entire model space, scattering points are selected as "model prediction target points" at intervals of 1 / 5 of the known PSF scattering point spacing.

[0116] Prepare multi-source input data (RTM imaging results, offset velocity field, and neighborhood scattering point PSF set) for each "model prediction target point" in the manner described in step S3, and input them into the trained PsfPredNet to predict the PSF corresponding to each "model prediction target point", forming a densely distributed scattering point PSF set in the entire model space.

[0117] Bilinear interpolation is performed on the spatially densely distributed set of scattering point PSFs to obtain the PSF of the entire model space. Since the data source of the interpolation is the high-density, high-precision scattering point PSFs predicted by the model, it avoids the problem of insufficient accuracy caused by the sparse data source of traditional interpolation. At the same time, it takes advantage of the spatial continuity of the interpolation method to ensure the smoothness of the PSF of the entire model space. Finally, a full model space PSF dataset covering all imaging points and with both high accuracy and spatial smoothness is formed.

[0118] Preferably, in step c, the formula for iteratively updating the RTM imaging results using the alternating direction multiplier algorithm is:

[0119]

[0120] This problem is a least-squares inversion problem with L1 norm constraints, which is solved using the traditional fast iterative soft thresholding algorithm.

[0121] After that, fix and For auxiliary variables Update:

[0122]

[0123] This subproblem can be solved by applying the generalized soft thresholding operator element by element to the auxiliary variable, resulting in a closed-form solution.

[0124]

[0125]

[0126] in, .

[0127] Finally, update the scaling dual variable as follows:

[0128]

[0129] Example 3: This invention provides an imaging domain least squares migration imaging method based on PsfPredNet, thereby solving the problem of insufficient resolution in imaging complex geological structures by existing imaging domain least squares migration imaging methods. (Refer to...) Figure 1 The embodiments of the present invention include the following steps:

[0130] Step S1: Acquisition of observation data and initial imaging results

[0131] Seismic shot gather data of the target work area were acquired, migration velocity field was obtained through conventional travel-time tomography inversion method, and reverse-time migration imaging results were obtained using conventional reverse-time migration method.

[0132] Step S2: Obtaining the point spread function of discrete scattering points

[0133] Based on the observation system and migration velocity field in S1, the reflectivity model of discrete scattering points is subjected to conventional inverse migration and reverse time migration to obtain the point spread function of discrete points.

[0134] Step S3: Training Dataset Construction

[0135] Using the scattering points of known PSFs as target points, we extract the set of PSFs of neighboring scattering points of the target point, the RTM imaging results of local sub-regions of the target point, and the migration velocity field. We then weight the PSFs of neighboring scattering points according to their Euclidean distance from the target point. Finally, we perform image enhancement processing such as rotation, flipping, and noise addition, and combine them to form a training dataset with a unified format. At the same time, we divide the training set and the validation set in a 7:3 ratio.

[0136] Among them, the neighborhood scattering point PSF set refers to the four scattering points with the smallest Euclidean distance centered on the target point (excluding the target point itself). The spatial range of the local sub-region of the RTM imaging result and the migration velocity field is the same as that of the target point PSF, and its center point is the target point.

[0137] The weighting process assigns weights based on the Euclidean distance between each scattering point and the target point. Scattering points closer to the target point are assigned greater weights, while those farther away are assigned smaller weights. The calculation formula is as follows:

[0138]

[0139] in, Indicates the first The weights corresponding to each scattering point This represents the Euclidean distance between the scattering point and the target point. This represents the total number of neighboring scattering points. In this embodiment of the invention, the number of neighboring scattering points is 4. Through the above normalization process, the sum of the weights of each scattering point can be guaranteed to be 1, avoiding the imbalance problem caused by scale differences, thereby more reasonably reflecting the influence of the neighborhood PSF on the target point PSF.

[0140] Step S4: PsfPredNet Network Training

[0141] The PsfPredNet neural network is built based on the ResUNet architecture, which has strong feature extraction and multi-source information fusion capabilities, and the convolutional layer parameters are initialized using the He normal distribution. The training set is input into the neural network, and the mean square error of the predicted PSF and the true PSF and the frequency domain loss are used as the loss function. The gradient descent method (AdamW optimizer) is used to iteratively optimize the network parameters. When the loss function value of the validation set no longer decreases, or when the number of iterations reaches the set upper limit, the training is stopped and the trained model is saved.

[0142] The loss function consists of a time-domain average absolute error term and a frequency-domain constraint term, and is defined as follows:

[0143]

[0144] in, The first real PSF One element, These are the predicted values ​​from the network model. Represents the Fast Fourier Transform. The total number of samples, and These are the weighting coefficients for the time-domain loss and the frequency-domain loss, respectively. In a preferred embodiment of the invention, the following is set: , .

[0145] Step S5: Iterative solution of LSRTM in the imaging domain based on PsfPredNet

[0146] Input the RTM imaging results, the migration velocity field, and the discrete point PSF, and set the maximum number of iterations and the regularization parameters of the objective function; use the PsfPredNet network to predict the PSF of the entire model space, and construct a Hessian matrix with the PSF of each point as a column in the matrix; iteratively update the RTM imaging results using the alternating direction multiplier algorithm; when the number of iterations reaches the preset maximum number of iterations, output the high-resolution imaging results of the current iteration number.

[0147] The objective function is as follows:

[0148]

[0149] In the formula, For offset images, This represents the approximate Hessian matrix formed by the PSFs predicted by PsfPredNet. For PsfPredNet network parameters, Distribution of underground reflectance coefficients , and For regularization parameters, For L1 norm constraint terms, For various same-sex TV regularization terms.

[0150] in,

[0151]

[0152] In the formula, As auxiliary variables, where , , and These are the first-order difference operators along the horizontal and vertical directions, respectively. This represents the number of grid cells in the reflectivity model.

[0153] The method of updating the reflectance image using the alternating direction multiplier method is described in the following details:

[0154] Construct the augmented Lagrangian function for the constrained least squares problem, which takes the following form:

[0155]

[0156] In the formula, For Lagrange multipliers, This is the penalty coefficient for the augmented quadratic term. The augmented quadratic term can penalize parts that violate the equation, and can be used to accelerate the convergence of the alternating direction multiplier algorithm. Perform a simple variable substitution A scaled version of the alternating direction multiplier algorithm can be derived:

[0157]

[0158] During the update process of the alternating direction multiplier algorithm, the first... The reflectance image at the next iteration can be represented as:

[0159]

[0160] This problem is a least squares subproblem with L1 norm constraints, which is solved using a fast iterative soft thresholding algorithm.

[0161] After that, fix and For auxiliary variables Update:

[0162]

[0163] This subproblem can be solved by applying the generalized soft thresholding operator element by element to the auxiliary variable, resulting in a closed-form solution.

[0164]

[0165]

[0166] in, .

[0167] Finally, update the scaling dual variable as follows:

[0168]

[0169] To make the technical solution of the present invention more intuitive, a specific application example is given below for illustration:

[0170] In this embodiment, a Marmousi model with a size of 250×552 was selected as the experimental object. The point spread function size was set to 31×31, and the scattering points were uniformly distributed in the horizontal and vertical directions within the imaging domain, with an interval of 31 grid units. The observation system uniformly arranged 62 shot points at 90-meter intervals on the ground surface, and the seismic signals were received by 552 surface geophones with an interval of 10 meters.

[0171] When applying the method of the present invention to this Marmousi model, the offset velocity field (e.g.) is set. Figure 3 After (as shown), the RTM imaging results are obtained using conventional RTM methods (such as...). Figure 4 As shown); based on the observation system and the migration velocity field, the reflectivity model of discrete scattering points is calculated using conventional inverse migration and RTM to obtain the PSF of discrete points (e.g. Figure 5 (As shown).

[0172] Using the scattering points of known PSFs as target points, we extract the PSF set of neighboring scattering points of the target point, the RTM imaging results of the local sub-region of the target point, and the offset velocity field, perform image enhancement processing, and combine them to form a training dataset with a unified format.

[0173] The constructed multi-source sample dataset was divided into training and validation sets in a 7:3 ratio, and then input into the constructed PsfPredNet network model for training. The batch size during network training was set to 8, and the initial learning rate was set to 1×10. -3 During training, an adaptive learning rate strategy is employed (the learning rate is halved when the validation set loss function no longer decreases for 5 consecutive rounds), gradually decaying to 1×10⁻⁶. -6 To ensure smooth convergence during training; the AdamW method is selected as the optimizer; in the loss function, the weight coefficients for the time-domain loss and the frequency-domain loss are respectively... , Furthermore, a Dropout probability of 0.2 is introduced into the convolutional layers, and the upper limit of the training epochs is set to 300 to suppress model overfitting.

[0174] The trained PsfPredNet network model is used for dense prediction and interpolation completion to construct a set of point spread functions covering the entire imaging domain. Figure 6 The PSF (Point Spread Function) of some discrete target points is predicted by the PsfPredNet network. To verify the advantages of the proposed method in point spread function prediction in complex structural regions, the PSF of target points with similar positions (such as PsfPredNet) obtained using the bilinear interpolation method in the traditional imaging domain LSRTM is presented. Figure 7 As shown), and the true PSF calculated using inverse offset and RTM (as shown). Figure 8 As shown). PSF calculated by traditional methods (such as...) Figure 7 The PSF (as shown) has a smooth shape, but loses its nonlinear distortion features, resulting in a significant deviation from the true PSF; while the PSF predicted by the method of this invention (such as...) Figure 6 As shown in the figure, the shape is highly consistent with the real PSF, the amplitude distribution is accurate, and the nonlinear distortion characteristics of complex tectonic regions are reproduced. By comparing the three, the high-precision advantage of the method of this invention in predicting point spread function in complex tectonic regions can be intuitively demonstrated.

[0175] The Hessian matrix is ​​constructed using the point spread function of the full model space described above, and then applied to the LSRTM iterative inversion. , and The values ​​are 1, 0.001, and 0.1, respectively. After 100 iterations, the imaging result of the 100th iteration is output (e.g., ...). Figure 9(As shown). To illustrate the advantages of this invention in imaging quality, imaging results calculated using the traditional imaging domain LSRTM method after 100 iterations are presented (e.g., Figure 10 (as shown) and true reflectance image (as shown) Figure 11 As shown). Results of LSRTM imaging in the traditional imaging domain (e.g.) Figure 10 As shown, numerous false reflection layer artifacts appear in the deep left region, with uneven amplitude distribution and poor continuity of the reflection layer; while the imaging results using the technology of this invention (as shown) Figure 9 As shown, the fracture interface is continuous and clear, the reflection layer artifacts are well suppressed, the amplitude fidelity is good, and it is closer to the true reflectance image (as shown). Figure 11 As shown in the figure. The comparison shows that the method of the present invention significantly improves the imaging quality, especially in complex structural areas, achieving high resolution and high fidelity imaging effects, and providing reliable technical support for the fine depiction of complex underground structures.

[0176] This invention establishes a nonlinear mapping between the PSF of neighboring scattering points, RTM imaging results, migration velocity field, and target point PSF using PsfPredNet mining. This mapping adaptively learns the drastic spatial variation characteristics of PSF caused by complex velocity fields. Compared to traditional imaging-domain LSRTM methods, this invention significantly reduces PSF prediction errors in complex structural areas such as fractures and salt domes, effectively improving the imaging quality of imaging-domain LSRTM, which is of great significance for the exploration of complex oil and gas reservoirs.

[0177] In another embodiment of the present invention, an imaging domain least squares migration imaging system based on PsfPredNet is provided, which can be used to implement the above-mentioned imaging domain least squares migration imaging method based on PsfPredNet. Specifically, the system includes:

[0178] The data acquisition module is used to acquire the observed seismic shot gather data of the target work area, obtain the migration velocity field through conventional travel-time tomography inversion, and obtain the reverse-time migration imaging results using conventional reverse-time migration methods.

[0179] The point spread function calculation module is used to perform conventional inverse migration and inverse time migration on the reflectivity model of discrete scattering points based on the observation system and the migration velocity field, so as to obtain the point spread function of discrete points.

[0180] The training module is used to construct a training dataset for the prediction network using the offset velocity field, inverse time offset imaging results, and point spread function of discrete points. Based on the training dataset, PsfPredNet with an encoder-decoder structure is trained.

[0181] The inversion module is used to embed the trained PsfPredNet into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variational constraints. It iteratively optimizes the reverse time migration imaging results through the alternating direction multiplier method to generate high-resolution seismic imaging results.

[0182] The module division in this embodiment of the invention is illustrative and represents only one logical functional division. In actual implementation, other division methods may be used. Furthermore, the functional modules in the various embodiments of the invention can be integrated into a single processor, exist as separate physical entities, or be integrated into a single module. The integrated modules described above can be implemented in hardware or as software functional modules.

[0183] In another embodiment of the present invention, a computer device is provided, comprising a processor and a memory. The memory stores a computer program, which includes program instructions. The processor executes the program instructions stored in the computer storage medium. The processor may be a Central Processing Unit (CPU), or other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. It is the computing and control core of the terminal, suitable for implementing one or more instructions, specifically suitable for loading and executing one or more instructions in the computer storage medium to achieve a corresponding method flow or corresponding function. The processor described in this embodiment of the present invention can be used for the operation of the least squares offset imaging method in the imaging domain based on PsfPredNet.

[0184] In another embodiment of the present invention, a storage medium is provided, specifically a computer-readable storage medium (Memory), which is a memory device in a computer device used to store programs and data. It is understood that the computer-readable storage medium here can include both the built-in storage medium in the computer device and extended storage media supported by the computer device. The computer-readable storage medium provides storage space that stores the terminal's operating system. Furthermore, the storage space also stores one or more instructions suitable for loading and execution by a processor. These instructions can be one or more computer programs (including program code). It should be noted that the computer-readable storage medium here can be high-speed RAM or non-volatile memory, such as at least one disk storage device. The processor can load and execute one or more instructions stored in the computer-readable storage medium to implement the corresponding steps of the PsfPredNet-based least-squares offset imaging method in the above embodiments.

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

[0186] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0187] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0188] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0189] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit it. Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the specific implementation of the present invention. Any modifications or equivalent substitutions that do not depart from the spirit and scope of the present invention should be covered within the scope of protection of the claims of the present invention.

Claims

1. A least-squares migration imaging method for the imaging domain based on PsfPredNet, characterized in that, include: The seismic shot gather data of the target work area were acquired, the migration velocity field was obtained through conventional travel-time tomography inversion, and the reverse-time migration imaging results were obtained using conventional reverse-time migration methods. Based on the observation system and the migration velocity field, the reflectivity model of discrete scattering points is subjected to conventional inverse migration and inverse time migration to obtain the point spread function of discrete points. Using the offset velocity field, reverse time offset imaging results, and point spread function of discrete points, a training dataset for the prediction network is constructed, and PsfPredNet based on the encoder-decoder structure is trained based on the training dataset. The trained PsfPredNet is embedded into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variation constraints. The reverse time migration imaging results are iteratively optimized by the alternating direction multiplier method to generate high-resolution seismic imaging results. The process involves constructing a training dataset for the prediction network using the offset velocity field, reverse time-shift imaging results, and the point spread function of discrete points. Training a PsfPredNet based on the encoder-decoder structure using this training dataset includes: The input samples of the training dataset consist of a set of point spread functions of neighborhood scattering points, the reverse time migration imaging results of local sub-regions, and the migration velocity field. The output label is the point spread function of the target point. The set of point spread functions of neighborhood scattering points is weighted according to Euclidean distance. The prediction network is a PsfPredNet neural network convolutional neural network built on the ResUNet architecture of a residual encoder-decoder network.

2. The imaging domain least squares migration imaging method based on PsfPredNet according to claim 1, characterized in that, The training of PsfPredNet based on the encoder-decoder structure using the training dataset includes: Its training process includes: Initialize network parameters; Input the training dataset into the network model; The loss function is the weighted sum of the mean squared error loss between the predicted point spread function and the true point spread function and the frequency domain loss. The network parameters are iteratively optimized by gradient descent until the loss function converges or the set maximum number of iterations is reached, thus completing the model training.

3. The imaging domain least squares migration imaging method based on PsfPredNet according to claim 2, characterized in that, During model training, the loss function consists of a time-domain mean absolute error term and a frequency-domain constraint term, defined as follows: in, The first of the real PSF One element, These are the predicted values ​​from the network model. Represents the Fast Fourier Transform. The total number of elements, and These are the weighting coefficients for the time-domain loss and the frequency-domain loss, respectively.

4. The imaging domain least squares migration imaging method based on PsfPredNet according to claim 3, characterized in that, The training process uses batch training, with each batch containing 2-8 samples and an initial learning rate of 10. -3 ~10 -4 .

5. The imaging domain least squares migration imaging method based on PsfPredNet according to claim 1, characterized in that, The method involves embedding the trained PsfPredNet into an imaging domain least-squares reverse time migration inversion method with sparse constraints and total variational constraints. The reverse time migration imaging results are iteratively optimized using the alternating direction multiplier method to generate high-resolution seismic imaging results, including: The least-squares reverse time migration (RTM) process integrated with PsfPredNet is as follows: Input the RTM imaging results, migration velocity field, and discrete point spread function; set the maximum number of iterations and the regularization parameter of the objective function; use the PsfPredNet network to predict the PSF of the entire model space, and construct a Hessian matrix with the PSF of each point as a column in the matrix; iteratively update the RTM imaging results using the alternating direction multiplier algorithm; when the number of iterations reaches the preset maximum number of iterations, output the high-resolution imaging result of the current iteration number. The objective function is as follows: In the formula, For offset images, This represents the approximate Hessian matrix formed by the PSFs predicted by PsfPredNet. For PsfPredNet network parameters, Distribution of underground reflectance coefficients, , and For regularization parameters, For L1 norm constraint terms, For various same-sex TV regularization terms; in, In the formula, As auxiliary variables, where , , and These are the first-order difference operators along the horizontal and vertical directions, respectively. This represents the number of grid points in the reflectivity model.

6. The imaging domain least squares migration imaging method based on PsfPredNet according to claim 5, characterized in that, The method of iteratively updating the RTM imaging results using the alternating direction multiplier method is as follows: Construct the augmented Lagrangian function for the constrained least squares problem, which takes the following form: In the formula, For Lagrange multipliers, The augmented quadratic term serves as a penalty coefficient; it penalizes parts that violate equality rules, thus accelerating the convergence of the alternating direction multiplier algorithm. Perform a simple variable substitution A scaled version of the alternating direction multiplier algorithm is derived: During the ADMM algorithm update process, the first... The reflectance image at the next iteration is represented as follows: This problem is a least squares subproblem with L1 norm constraints, which is solved using a fast iterative soft thresholding algorithm. After that, fix and For auxiliary variables Update: This subproblem can be solved by applying the generalized soft thresholding operator element by element to the auxiliary variable, resulting in a closed-form solution. in, ; Finally, update the scaling dual variable as follows: 。 7. A least-squares migration imaging system for the imaging domain based on PsfPredNet, characterized in that, include: The data acquisition module is used to acquire the observed seismic shot gather data of the target work area, obtain the migration velocity field through conventional travel-time tomography inversion, and obtain the reverse-time migration imaging results using conventional reverse-time migration methods. The point spread function calculation module is used to perform conventional inverse migration and inverse time migration on the reflectivity model of discrete scattering points based on the observation system and the migration velocity field, so as to obtain the point spread function of discrete points. The training module is used to construct a training dataset for the prediction network using the offset velocity field, inverse time offset imaging results, and point spread function of discrete points. Based on the training dataset, PsfPredNet with an encoder-decoder structure is trained. The inversion module is used to embed the trained PsfPredNet into the imaging domain least squares reverse time migration inversion method with sparse constraints and total variational constraints. It iteratively optimizes the reverse time migration imaging results through the alternating direction multiplier method to generate high-resolution seismic imaging results. The process involves constructing a training dataset for the prediction network using the offset velocity field, reverse time-shift imaging results, and the point spread function of discrete points. Training a PsfPredNet based on the encoder-decoder structure using this training dataset includes: The input samples of the training dataset consist of a set of point spread functions of neighborhood scattering points, the reverse time migration imaging results of local sub-regions, and the migration velocity field. The output label is the point spread function of the target point. The set of point spread functions of neighborhood scattering points is weighted according to Euclidean distance. The prediction network is a PsfPredNet neural network convolutional neural network built on the ResUNet architecture of a residual encoder-decoder network.

8. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the steps of the imaging domain least squares offset imaging method based on PsfPredNet as described in any one of claims 1 to 6.

9. A computer-readable storage medium storing a computer program, characterized in that, When the computer program is executed by the processor, it implements the steps of the imaging domain least squares offset imaging method based on PsfPredNet as described in any one of claims 1 to 6.

Images