A seismic wave field reconstruction method and an elastic wave reverse time migration method
By employing a source wavefield reconstruction method and utilizing Gaussian filtering and dictionary learning algorithms, the high computational resource consumption and stability issues in reverse time migration were resolved. This resulted in efficient and stable source wavefield reconstruction, reduced storage requirements, and improved computational efficiency.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-04-12
- Publication Date
- 2026-03-24
AI Technical Summary
Existing reverse time migration methods in seismic exploration are computationally expensive, unstable in attenuating media, have excessive storage requirements, and low computational efficiency.
The source wavefield reconstruction method is adopted, including velocity model preprocessing, source wavefield forward modeling, dictionary learning dataset construction, wavefield feature dictionary learning, sparse representation and reconstruction of source wavefield. The source wavefield is reconstructed by using Gaussian filtering, staggered grid finite difference method and dictionary learning algorithm through sparse coding.
It effectively reduces storage requirements, improves computing efficiency and stability, ensures stable operation in decaying media, and achieves a balance between reconstruction accuracy and compression ratio.
Smart Images
Figure CN116400410B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of seismic migration imaging, specifically relating to a source wavefield reconstruction method and an elastic wave reverse time migration method. Background Technology
[0002] With the advancement of seismic exploration technology, the demand for imaging techniques for complex subsurface media in exploration projects is increasing. Reverse-time migration (RTM) offers an effective solution, but its associated state method requires simultaneous access to both the source and detector wavefields for high-resolution subsurface imaging and inversion. Therefore, whether applying imaging conditions in RTM or calculating gradient equations in full-waveform inversion, it involves storing source wavefields with large time steps, placing an unbearable storage burden on computing equipment. Consequently, this problem remains a pressing issue to be addressed in seismic imaging and inversion applications.
[0003] Currently, approaches to solving storage problems mainly fall into two categories: checkpoint-based methods and boundary value-based methods. Checkpoint-based methods significantly reduce storage requirements, but the additional recalculation greatly increases computational resource consumption; essentially, they utilize ample computational resources to alleviate the storage burden. Unlike checkpoint methods, boundary value-based methods use attenuated wavefield values in artificial boundary regions, requiring reconstruction of the source wavefield, and offer higher accuracy. However, this boundary value-based approach requires re-performing a forward modeling simulation when reconstructing the source wavefield, resulting in low computational efficiency when dealing with multi-wave, multi-component problems, and stability issues when facing attenuating media. Summary of the Invention
[0004] The technical problem to be solved by the present invention is to provide a source wavefield reconstruction method and an elastic wave reverse time migration method, which solves the problems of high computational resource consumption and stability of existing methods.
[0005] This invention is implemented as follows:
[0006] A source wavefield reconstruction method includes: velocity model preprocessing, source wavefield forward modeling, dictionary learning dataset construction, wavefield feature dictionary learning, sparse representation of source wavefield, and source wavefield reconstruction.
[0007] The velocity model preprocessing includes performing Gaussian filtering on the original velocity model to obtain a smoothed initial velocity model;
[0008] The forward modeling of the source wavefield includes setting the Ricker wavelet as the explosion source to excite the wavefield in forward time sequence, inputting a smoothed initial velocity model, using CPML absorbing boundary conditions to attenuate the false reflected waves generated at the boundary, and using the staggered grid finite difference method to solve the decoupled elastic wave equation.
[0009] The construction of the dictionary learning dataset includes filtering wavefield images from the source wavefield, separating the wavefield images using a sliding window technique to obtain a series of wavefield image patches, performing a normalization step on the wavefield image patches, and reconstructing them into a dictionary learning dataset.
[0010] The wavefield feature dictionary learning includes training a wavefield feature dictionary from the aforementioned dictionary learning dataset using the K-SVD dictionary learning algorithm;
[0011] The sparse representation of the source wavefield includes using the batch-OMP algorithm and a wavefield feature dictionary to sparsely encode the source wavefield, traversing the source wavefields corresponding to all time steps, and reassembling them into a joint sparse matrix.
[0012] The source wavefield reconstruction includes partitioning the joint sparse matrix according to the recombination order, reconstructing the wavefield image corresponding to each time step, and arranging all wavefield images in positive time order to complete the reconstruction of the source wavefield.
[0013] Furthermore, a Gaussian filter is applied to the original velocity model to obtain a smoothed initial velocity model, specifically including:
[0014] Given a velocity model V(x, z) and a Gaussian filter G, applying a Gaussian filter to the velocity model yields:
[0015] V f (x,z)=V(x,z)*G (1)
[0016] In the formula, V f (x, z) represents the initial velocity model after filtering the true velocity model.
[0017] Furthermore, solving the decoupled elastic wave equation includes:
[0018] Set V f (x, z) represents the longitudinal wave model, and the transverse wave velocity model V s Given (x, z), density model ρ(x, z), and source wavelet f(t), solve the following decoupled elastic wave equation to obtain the vector source wavefield.
[0019]
[0020]
[0021]
[0022]
[0023]
[0024]
[0025]
[0026] v xs =v x -v xp
[0027] v zs =v z -v zp (2)
[0028] Among them, (v x v z ) is the velocity vector of the particle, (σ) xx , σ zz , σ xz ) represents the proton stress vector, ρ represents the density, and S represents the density. xx and S zz For the source terms loaded in the 'xx' and 'zz' directions, λ and μ are elastic Lamé parameters, (v xp v zp ) and (v xs v zs ) are the longitudinal wave and transverse wave velocity vectors, respectively.
[0029] Furthermore, the construction of the dictionary learning dataset specifically includes:
[0030] Vector source wavefield images at h time points were selected. For each vector source wavefield image, a sliding window technique is used to separate it into a series of 8×8 image blocks. These separated image blocks are then normalized into corresponding column vectors and reassembled into a vector matrix in the positive order of separation.
[0031] Traversing all vector source wavefields, arranging all recombined matrices in ascending order of m, forms the joint vector matrix W. k , with W k As a dataset.
[0032] Furthermore, the sparse representation of the source wavefield specifically includes:
[0033] Input source wavelet f(t), P-wave velocity model V p (x, z), shear wave velocity model V s (x, z), density model ρ(x, z), and the vector source wave field at each time step is obtained by solving formula (2) using the staggered grid finite difference method of time order 2 and space order 10.
[0034] Using the sliding window technique to extract vector source wave fields The image is separated into a series of 8×8 blocks, and then the image blocks are normalized into corresponding column vectors. l is Includes an index of image patches;
[0035] Based on the wavefield feature dictionary, the column vector is obtained by solving formula (3) using batch-OMP. Corresponding sparse representation
[0036]
[0037] Where D k A dictionary of wavefield features;
[0038] Iterate through all column vectors and reorganize all sparse coefficients into a vector sparse matrix θ. k (t); Traverse all vector source wave fields and obtain all vector sparse matrices θ k (t) is reorganized into a joint vector sparse matrix θ in positive time order. k .
[0039] A method for performing elastic wave reverse time migration includes:
[0040] The multi-component seismic record {v} is simulated by forward modeling of the vector source wavefield. z v x Using the initial value of the receiver point, the finite difference method with a staggered grid of second-order time and tenth-order space is used to solve formula (2) to obtain the receiver point vector wave field at each time step.
[0041] During the reverse-time extrapolation of the detector points, the sparse matrix θ k The data is divided according to time, and the sparse coefficients in the sparse matrix of each division are... The corresponding column vector is reconstructed using the multi-component wavefield feature dictionary.
[0042]
[0043] Traverse all column vectors The image is normalized into corresponding image patches and placed in their original separated positions to obtain the reconstructed source wavefield.
[0044] The receiver point vector wave field obtained through calculation and reconstructed longitudinal wave vector wave field Single-shot PP and PS imaging results under the following imaging conditions:
[0045]
[0046] Among them, I ij (x, z) is a single-shot elastic wave reverse-time migration imaging profile, containing I pp (x, z) and I ps (x, z), where · is the scalar product. as well as These are the reconstructed source vector wavefield and vector receiver wavefield, respectively, where z and x are position coordinates and t is the time node. max The maximum value is the observation time, and δ is a local minimum value;
[0047] All imaging profiles are superimposed according to the number of shot points, and finally, a Laplace filter is used to filter out low-frequency noise from the superimposed profiles.
[0048] Compared with the prior art, the beneficial effects of this invention are as follows:
[0049] 1. This invention proposes a source wavefield reconstruction method based on sparse representation and dictionary learning. The sparse representation method is used to compress the source wavefield, and the resulting sparse domain wavefield file has a significantly reduced memory footprint compared to wavefield data stored in the third-order tensor mode, thus solving the problem of high storage requirements faced by the adjoint cases of wave equations.
[0050] 2. The technical approach of this invention is not synchronized with the attenuation compensation process, which ensures the stable operation of this invention in attenuation medium conditions and solves the instability problem of existing boundary value-based methods in this situation.
[0051] 3. In this invention, by adjusting the sparse approximation error and the size of the wavefield feature dictionary, the running time can be reasonably controlled while ensuring reconstruction accuracy and compression ratio, thus ensuring that this invention can achieve a balance between reconstruction accuracy, compression ratio, and computational efficiency.
[0052] 4. This invention uses source wavefield data as the dataset and employs the K-SVD dictionary learning algorithm to learn a wavefield feature dictionary from the dataset. Compared with a fixed-type dictionary model, this invention improves the adaptability and stability of the sparse representation process, further enhances the sparse representation performance of the dictionary on the source wavefield signal, and ensures that sparse coding can stably compress the source wavefield.
[0053] 5. This invention treats the source wavefield corresponding to each time step as a two-dimensional image and introduces image processing methods to divide the original source wavefield into blocks and process it in time. Compared with sparse representation of the entire wavefield, this reduces the reconstruction error of the source wavefield and improves the computational efficiency. Attached Figure Description
[0054] Figure 1This is a schematic diagram of a source wavefield reconstruction method based on sparse representation and dictionary learning in this embodiment.
[0055] Figure 2 This is a schematic diagram illustrating the application of a source wavefield reconstruction method based on sparse representation and dictionary learning provided in this embodiment to the reverse time migration of elastic waves.
[0056] Figure 3 This is a schematic diagram of wavefield feature dictionary learning provided in this embodiment.
[0057] Figure 4 The P-wave velocity model V provided in this embodiment p A schematic diagram of (x, z).
[0058] Figure 5 The S-wave velocity model V provided in this embodiment s A schematic diagram of (x, z).
[0059] Figure 6 This is a schematic diagram of the P-wave X-component wavefield of the seismic source provided in this embodiment (t = 1500 ms).
[0060] Figure 7 This is a schematic diagram of the P-wave X-component wavefield characteristic dictionary provided in this embodiment.
[0061] Figure 8 This is a schematic diagram of the reconstruction of the P-wave X-component wave field provided in this embodiment.
[0062] Figure 9 This is a schematic diagram of the residual wavefield between the original P-wave X component wavefield and the reconstructed P-wave X component provided in this embodiment.
[0063] Figure 10 This is a schematic diagram of the conventional elastic wave reverse time migration PP wave imaging results provided in this embodiment.
[0064] Figure 11 This is a schematic diagram of the conventional elastic wave reverse time migration PS wave imaging results provided in this embodiment.
[0065] Figure 12 This is a schematic diagram of the source wavefield reconstruction method based on sparse representation and dictionary learning provided in this embodiment, applied to the PP imaging results of elastic wave reverse time migration.
[0066] Figure 13 This is a schematic diagram of the PS imaging results of a source wavefield reconstruction method based on sparse representation and dictionary learning provided in this embodiment, applied to the reverse time migration of elastic waves.
[0067] Figure 14 This is a schematic diagram of the residual wave field between the original PP wave imaging result provided in this embodiment and the PP wave imaging result of the method described in this invention.
[0068] Figure 15 This is a schematic diagram of the residual wave field between the original PS wave imaging result provided in this embodiment and the PS wave imaging result of the method described in this invention. Detailed Implementation
[0069] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.
[0070] This embodiment uses a salt dome model as an example. This model has a transverse length of 1000m and a longitudinal length of 256m. The longitudinal wave model and transverse wave model are as follows: Figure 4 and Figure 5 As shown, the spatial step size is set to 10m, the number of shot points is set to 100, and shots are fired on the model surface. Detector points are evenly distributed on the model surface grid at 10m intervals. The observation and recording time is set to 5s, the time step size is 0.1ms, the density model is set to a constant of 2000, the number of time nodes is set to 10, the step size of the sliding window is set to 0, the sparse approximation error of the batch-OMP algorithm is set to 0.0001, and the number of iterations of the K-SVD algorithm is set to 60. The present invention will be further described in detail below with reference to the schematic diagram provided in this embodiment.
[0071] like Figure 1 As shown, the source wavefield reconstruction method provided in this embodiment includes the following steps: velocity model preprocessing, source wavefield forward modeling, dictionary learning dataset construction, wavefield feature dictionary learning, sparse representation of source wavefield, and source wavefield reconstruction.
[0072] The velocity model preprocessing involves applying a Gaussian filter to the original velocity model to obtain a smoothed initial velocity model. Specifically, this includes:
[0073] Given a velocity model V(x, z) and a Gaussian filter G, applying a Gaussian filter to the velocity model yields:
[0074] V f (x,z)=V(x,z)*G (1)
[0075] In the formula, V f (x, z) represents the initial velocity model after filtering the true velocity model.
[0076] Forward modeling of the source wavefield is performed by setting the Ricker wavelet as the source of the explosion and exciting the wavefield in forward time sequence. The smoothed initial velocity model is input, and the CPML absorbing boundary conditions are used to attenuate the spurious reflected waves generated at the boundary. The wave equation is solved using the staggered grid finite difference method.
[0077] Set V f (x, z) represents the longitudinal wave model, and the transverse wave velocity model V s Given (x, z), density model ρ(x, z), and source wavelet f(t), solve the following decoupled elastic wave equation to obtain the vector source wavefield.
[0078]
[0079]
[0080]
[0081]
[0082]
[0083]
[0084]
[0085] v xs =v x -v xp
[0086] v zs =v z -v zp (2)
[0087] Among them, (v x v z ) is the velocity vector of the particle, (σ) xx , σ zz , σ xz ) represents the proton stress vector, ρ represents the density, and S represents the density. xx and S zz For the source terms loaded in the 'xx' and 'zz' directions, λ and μ are elastic Lamé parameters, (v xp v zp ) and (v xs v zs ) are the longitudinal wave and transverse wave velocity vectors, respectively.
[0088] The construction of the dictionary learning dataset requires selecting wavefield images with sufficient information from the source wavefield, separating the wavefield images using the sliding window technique to obtain a series of wavefield image patches, performing a normalization step on the wavefield image patches, and reconstructing them into a dictionary learning dataset.
[0089] Vector source wavefield images at h time points were selected. For each vector source wavefield image, a sliding window technique is used to separate it into a series of 8×8 image blocks. These separated image blocks are then normalized into corresponding column vectors and reassembled into a vector matrix in the positive order of separation.
[0090] Traversing all vector source wavefields, arranging all recombined matrices in ascending order of m, forms the joint vector matrix W. k , with W k As a dataset.
[0091] The wavefield feature dictionary is trained using the K-SVD dictionary learning algorithm from the aforementioned dictionary learning dataset.
[0092] The sparse representation of the source wavefield utilizes the batch-OMP algorithm and a wavefield feature dictionary to sparsely encode the source wavefield, traversing the source wavefield corresponding to all time steps and reconstructing it into a joint sparse matrix. This includes: the input source wavelet f(t), the P-wave velocity model V... p (x, z), shear wave velocity model V s (x, z), density model ρ(x, z), and the vector source wave field at each time step is obtained by solving formula (2) using the staggered grid finite difference method of time order 2 and space order 10.
[0093] Using the sliding window technique to extract vector source wave fields The image is separated into a series of 8×8 blocks, and then the image blocks are normalized into corresponding column vectors. l is Includes an index of image patches;
[0094] Based on the wavefield feature dictionary, the column vector is obtained by solving formula (3) using batch-OMP. Corresponding sparse representation
[0095]
[0096] Where D k A dictionary of wavefield features;
[0097] Iterate through all column vectors and reorganize all sparse coefficients into a vector sparse matrix θ. k (t); Traverse all vector source wave fields and obtain all vector sparse matrices θ k (t) is reorganized into a joint vector sparse matrix θ in positive time order. k .
[0098] The source wavefield reconstruction is performed by partitioning the joint sparse matrix according to the recombination order, accurately reconstructing the wavefield image corresponding to each time step, and arranging all wavefield images in positive time order to complete the reconstruction of the source wavefield.
[0099] like Figure 2 As shown in Figure 3, the specific implementation steps of the method for elastic wave reverse time migration using the source wavefield reconstruction method described in this invention are as follows:
[0100] (1) Given a velocity model V(x, z) and a Gaussian filter G, applying a Gaussian filter to the velocity model yields:
[0101] V f (x, z) = V(x, z) * G (1)
[0102] In the formula, V f (x, z) represents the initial velocity model after filtering the true velocity model;
[0103] (2) Input longitudinal wave velocity model V f (x, z), shear wave velocity model V s Given (x, z), density model ρ(x, z), and source wavelet f(t), the following decoupled elastic wave equation is solved using the second-order time and tenth-order spatial staggered grid finite difference method to obtain the P-wave vector wavefield P at each time step. k (x, z, t)(k = x, z):
[0104]
[0105]
[0106]
[0107]
[0108]
[0109]
[0110]
[0111] v xs =v x -v xp
[0112] v zs =v z -v zp (2)
[0113] Among them, (v x v z) is the velocity vector of the particle, (σ) xx , σ zz , σ xz ) represents the proton stress vector, ρ represents the density, and S represents the density. xx and S zz For the source terms loaded in the 'xx' and 'zz' directions, λ and μ are elastic Lamé parameters, (v xp v zp ) and (v xs v zs ) represents the velocity vectors of P-waves and S-waves;
[0114] (2) Filter out the longitudinal wave vector wave field images at h time points. For each P-wave vector wavefield image, a sliding window technique is used to separate it into a series of 8×8 image blocks. These separated image blocks are then normalized into corresponding column vectors and reassembled into a vector matrix in the positive order of separation.
[0115] (3) Traverse all vector P-wave field images, arrange all recombined matrices in ascending order of m, and form a joint vector matrix W. k , with W k As a dataset, the P-wave field feature dictionary D was learned using the K-SVD algorithm. k .
[0116] (4) Input source wavelet f(t), P-wave velocity model V p (x, z), shear wave velocity model V s (x, z), density model ρ(x, z), and the longitudinal wave vector field P at each time step are obtained by solving formula (2) using the staggered grid finite difference method of time order 2 and space order 10. k (x, z, t)(k = x, z).
[0117] (5) In particular, further utilizing the sliding window technique to P k (x, z, t) is separated into a series of 8×8 image patches, and then the image patches are normalized into corresponding column vectors. For wave field P k (x, z, t) contains the image patch index;
[0118] (6) The wave field feature dictionary D obtained through step (3) k The column vector is obtained by solving the following formula using batch-OMP. Corresponding sparse representation
[0119]
[0120] (7) Traverse all column vectors and reorganize all sparse coefficients into a vector sparse matrix θ. k (t); Traverse all longitudinal wave vector wave fields and obtain all vector sparse matrices θ k (t) is reorganized into a joint vector sparse matrix θ in positive time order. k ;
[0121] (8) After the forward modeling simulation is completed, save the multi-component seismic records from the forward modeling simulation {v z v x}
[0122] (9) Input multi-component seismic records {v z v x Using this as the initial value for the receiver, the finite difference method with a staggered grid of second-order time and tenth-order space is further used to solve formula (2) to obtain the receiver vector wave field at each time step.
[0123] (10) During the reverse time extrapolation of the detector point, the sparse matrix θ i (t) is divided according to time, and the sparse coefficients in the sparse matrix of each division are... The corresponding column vector is reconstructed using the multi-component wavefield feature dictionary.
[0124]
[0125] (11) Traverse all The image is normalized to its corresponding patch and placed in its original separated position, thus obtaining the reconstructed longitudinal wave vector field.
[0126] (12) Calculated receiver point vector wave field and reconstructed longitudinal wave vector wave field Single-shot PP and PS imaging results under the following imaging conditions:
[0127]
[0128] Among them, I pj (z, x) is a single-shot elastic wave reverse-time migration imaging profile, containing I pp (z, x) and I ps (z, x), where · is the scalar product. as well as The reconstructed P-wave vector wavefield and receiver point vector wavefield are given, where z and x are position coordinates, and t is the time node. max The maximum observation time is represented by δ, which is a minimum value to ensure the stability of the profile calculation.
[0129] (13) All imaging profiles are superimposed according to the number of shot points, and finally the low-frequency noise of the superimposed profiles is filtered out using a Laplace filter.
[0130] The following describes the actual effects of the present invention in this embodiment:
[0131] The method provided in this invention was applied to salt dome model imaging, achieving satisfactory PP and PS imaging results. First, the X component wavefield of the source P-wave at 1500 ms was extracted, and then compared... Figure 6 , Figure 7 , Figure 8 and Figure 9 It can be seen that the reconstructed P-wave X-component wavefield is visually almost identical to the original P-wave X-component wavefield, while the amplitude of its residual wavefield image is only one percent of the original wavefield, with no loss of wavefield energy within the visible range. To compare the imaging effects, PP and PS imaging results obtained using the conventional elastic wave reverse-time migration method are shown below. Figure 10 , Figure 11 As shown, the source wavefield reconstruction method based on sparse representation and dictionary learning provided by this invention is applied to the PP (Positioning over Transient Time) imaging of elastic waves. The PS imaging results are as follows: Figure 12 and Figure 13 As shown, the residual image between the reconstructed result and the original imaging result is as follows: Figure 14 and Figure 15 As shown, the results of the method described in this invention are almost identical to those of conventional elastic wave reverse time migration methods, and the residual image reveals that the amplitude loss is relatively small. Compared to the huge storage requirements of the source wavefield, this invention, by reasonably adjusting the sparse approximation error and dictionary size, significantly reduces the memory requirements of the joint sparse matrix file compared to the original wavefield. When the sparse approximation error is 0.0001, the storage requirement is only 13.4% of the original, while the additional computational cost from sparse coding and reconstruction is 38.5% of the original; when the sparse approximation error is 0.001, the storage requirement is only 1.8% of the original, while the additional computational cost from sparse coding and reconstruction is 20.9% of the original. Compared to traditional checkpoint-based methods, this is still acceptable.
[0132] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the protection scope of the present invention.
Claims
1. A method of seismic wavefield reconstruction, characterized by, The method comprises: velocity model preprocessing, source wave field forward modeling, dictionary learning data set construction, wave field feature dictionary learning, source wave field sparse representation, and source wave field reconstruction; The velocity model preprocessing comprises Gaussian filtering of the original velocity model to obtain a smoothed initial velocity model; The source wave field forward modeling comprises setting a Ricker wavelet as an explosion source to excite a wave field in a time forward order, inputting the smoothed initial velocity model, using a CPML absorbing boundary condition to attenuate false reflection waves generated at a boundary, and using a staggered grid finite difference method to solve a decoupled elastic wave equation; The dictionary learning data set construction comprises screening wave field images from the source wave field, separating the wave field images through a sliding window technique to obtain a series of wave field image blocks, performing a normalization step on the wave field image blocks, and recombining the wave field image blocks into a dictionary learning data set; The wave field feature dictionary learning comprises using a K-SVD dictionary learning algorithm to train a wave field feature dictionary from the dictionary learning data set; The source wave field sparse representation comprises using a batch-OMP algorithm and the wave field feature dictionary to perform sparse coding on the source wave field, traversing all time steps of the source wave field to recombine a joint sparse matrix; The source wave field reconstruction comprises splitting the joint sparse matrix according to a recombination order, reconstructing a wave field image corresponding to each time step, and returning all wave field images in a time forward order to complete reconstruction of the source wave field.
2. The seismic wavefield reconstruction method according to claim 1, wherein, The Gaussian filtering of the original velocity model to obtain the smoothed initial velocity model specifically comprises: Given velocity model Gaussian filter Gaussian filtering the velocity model gives: （1）， In the formula, denotes the initial velocity model filtered with the true velocity model.
3. The seismic wavefield reconstruction method of claim 2, wherein, The solving of the decoupled elastic wave equation specifically comprises: Setting for the P-wave model, the S-wave velocity model , the density model , the source wavelet , solve the decoupled elastic wave equation as follows to obtain the vector source wave field : (2)， where is the particle velocity vector, is the proton stress vector, is the density, and is the source term in and directions, and are the elastic Lame parameters, and are the longitudinal and transverse velocity vectors, respectively.
4. The seismic wavefield reconstruction method of claim 1, wherein, The dictionary learning data set construction specifically comprises: Screening out Vectorial source wavefield images at individual time nodes , for each vectorial source wavefield image, separate it into a series of image tiles of size using a sliding window technique, and then normalize the separated image tiles into corresponding column vectors, reorganize them into a vectorial matrix in the order of separation; Traversing all vector source wave fields, according to Arrange all the recombined matrices in ascending order to form a joint vector matrix. ,by As a dataset.
5. The seismic wavefield reconstruction method of claim 3, wherein, The source wave field sparse representation specifically comprises: Input seismic source wavelet , P-wave velocity model , S-wave velocity model , Density model , the formula (2) is solved by using the staggered grid finite difference method with time 2 order and space 10 order to obtain the vector seismic source wave field at each time ; Separating the vector source wavefield into a series of image patches of size containing the image patch index; According to the wave field feature dictionary, the column vector is obtained by using batch-OMP to solve formula (3) The corresponding sparse representation : （3）， wherein is a wavefield signature dictionary; reorganize all sparse coefficients into vector sparse matrices by traversing all column vectors ; reorganize all vector sparse matrices into joint vector sparse matrix by traversing all vector source wavefields in time forward order in time forward order .
6. An elastic wave reverse time migration method using the seismic wavefield reconstruction method of any one of claims 1-5, characterized in that, The source wave field reconstruction specifically comprises: Multi-component seismic records from vector source wavefield forward modeling As the initial value of the receiver, the finite difference method of staggered grid of time 2 order and space 10 order is used to solve formula (2), and the vector wave field of each time of the receiver is obtained ; In performing the receiver point reverse-time extrapolation process, the sparse matrix is partitioned by time, and the sparse coefficients in each partitioned sparse matrix are reconstructed into a corresponding column vector using the multi-component wavefield feature dictionary （4）， traversing all column vectors , which are normalized to corresponding image blocks and placed in originally separate locations, to finally obtain the reconstructed source wave field ; The calculated receiver vector wavefield and the reconstructed P-wave vector wavefield The single-shot PP and PS imaging results are obtained by using formula (5) imaging conditions. （5）， wherein, is a single-shot elastic wave reverse-time migration image section, comprising and , is a scalar product, and are respectively a reconstructed source vector wavefield and a vector receiver wavefield, and are position coordinates, is a time node, is an observation time maximum, is a minimum. The stacking of all imaging sections according to the number of shots, and finally the filtering of low-frequency noise of the stacked sections using a Laplace filter.