In multicomponent seismic exploration, establishing depth domain interval velocity models for compressional waves (P-wave) and shear waves (S-wave) is a key step in data imaging processing and inversion. Compared to techniques such as tomography and migration velocity analysis that only use travel time to obtain velocity information, elastic wave full waveform inversion (EFWI) utilizes information such as travel time, phase, and amplitude of elastic waves to establish P- and S-wave velocity models by minimizing the residual between the observed multicomponent data and simulated multicomponent seismic data of a specific velocity model. Therefore, EFWI has the potential to reveal structural details and lithology in complex geological backgrounds, and theoretically can obtain high-resolution depth domain P- and S-wave velocity models.

The idea of full waveform inversion (FWI) of seismic waves was first proposed by Lailly (1983) and Tarantola (1984). Under the theoretical framework of generalized least squares, FWI calculates the gradient by correlating the simulated wave field and the reverse wave field of the residual record between simulated and observed seismic data, and then updates the parameter model through continuous iteration. Since FWI can make full use of the kinematic and dynamic information of seismic waves to estimate the propagation speed of seismic waves and theoretically can obtain a higher resolution velocity model, the idea of FWI has received extensive research interest since it was proposed and research progress has been made in the objective function optimization (Datta et al., 2016; Zhu et al., 2016), multi-scale inversion (Bunks et al., 1995; Boonyasiriwat et al., 2009; Xu et al., 2014), mixed domain inversion (Kim et al., 2013; Jun et al., 2014; Xu et al., 2014), envelope inversion (Wu et al., 2013; Wu et al., 2014; Ao et al., 2015), source wavelet inversion (Tarantola, 1984; Song et al., 1995; Hu et al., 2017), and inversion efficiency improvement (Krebs et al., 2009; Wang et al., 2011) of FWI. Currently, FWI has been widely used in field seismic data imaging and many application examples have significantly improved the imaging quality of seismic data (Sirgue et al., 2010; Lewis et al., 2014; Liu et al., 2014; Zhong et al., 2017).

Tarantola (1986) and Pratt (1990) extended the compressional wave FWI to EFWI. Although some ideas and techniques in longitudinal wave FWI can be extended or even directly applied in EFWI, the non-linearity of inversion is further exacerbated by the inclusion of more inversion parameters in EFWI, and due to the presence of various types of noise in seismic data, EFWI based on three-component waveform matching is very sensitive to these noises, and the inversion is easily affected by these noises and falls into local extrema. At the same time, if the minimum available frequency in the three-component seismic data is too high, it will greatly enhance the dependence of EFWI on the initial model, and it is difficult to obtain good inversion results when the accuracy of the initial model is low. In addition, when the source wavelet is inaccurate, there is a significant difference between the synthetic data obtained from the erroneous source wavelet and observed data, which often leads to the wrong direction of EFWI and increases the difficulty of its application.

The lack of low-frequency information in observed three-component seismic data is one of the key factors leading to a decrease in the accuracy of EFWI. Baeten et al. (2013) pointed out that the low-frequency information of 1.5∼2 Hz in seismic records is particularly important for alleviating the “cycle skipping” phenomenon in FWI. The key to solving this problem is how to provide an accurate initial model for FWI when low-frequency information is missing. Envelope inversion is often used to construct the initial model of FWI. Bozdağ et al. (2011) pointed out that the seismic envelope contains rich low-frequency information and using the envelope as input data for FWI inversion can establish a more accurate initial model. The research of Huang et al. (2015) indicates that envelope inversion can significantly improve the accuracy of the inversion of the compressional and shear wave velocity model when seismic data lack low-frequency information. Wu and Chen (2017, 2018, 2020) proposed the direct envelope Fréchet derivative and the direct envelope inversion (DEI) method, which can map the ultra-low frequency envelope data perturbation to the velocity perturbation directly and can invert the large scale strong-scattering velocity model without low-frequency information in original common shot gathers. Chen et al. (2018) combined the DEI method with the wavefield direction decomposition method and proposed a reflection DEI method, which can improve the inversion effects of the velocity structures in the strong-scattering shielding area. However, traditional envelope inversion is based on the accurate source wavelet assumption. When the source wavelet is inaccurate, envelope inversion cannot construct a reasonable initial model. Therefore, Ao et al. (2015) proposed a convolutional envelope objective function based on the compressional equation, which can eliminate the impact of wavelet inaccuracy on the accuracy of compressional envelope inversion.

We extend the source-independent FWI (Choi et al., 2005) to EFWI to improve the accuracy of envelope inversion when the source wavelet is inaccurate. On the basis of previous research, a misfit function for the elastic wave convolution envelope was established, and corresponding gradient and adjoint source formulas were derived. Applying these methods to mixed domain EFWI can improve the inversion accuracy of compressional and shear wave velocity models in cases of wavelet inaccuracy and missing low-frequency information.

The misfit function of EFWI in the time domain can be written as (Tarantola, 1986; Prat, 1990): E m = 1 2 ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T u i , j − d i , j 2 d t , (1) where m = v p , v s , v p is the compression wave velocity, v s is the shear wave velocity, N s and N r are the shot number and receiver number per shot, respectively, and u and d are the modeled and observed three-component vector wave fields, respectively.

For each shot gather, select u i , k and d i , k as the reference traces of the modeled and the observed data, respectively, and k as the channel number of the reference trace. Convolve the simulated traces with the reference traces of observed wavefields to obtain the modeled wavefields, and convolve the observed traces with the reference traces of simulated wavefields to obtain the modeled wavefields. So Eq. 1 becomes a new convolution misfit function: E c v m = 1 2 ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T u i , j *   d i , k − d i , j * u i , k 2 d t , (2) since the seismic traces of any components in the three-component seismic data can be regarded as Green’s function and source wavelet convolution, that is: d t = g t * s t , (3) where d is a certain component in the three-component data, g is the Green’s function, and s is the source wavelet.

Equation 2 can be expressed as: E c v m = 1 2 ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T g i , j u * g i , k d * s i d * s i u − g i , j d * g i , k u * s i d * s i u 2 d t , (4) where g u is the Green’s function of modeled three-component data, g d is the Green’s function of observed three-component data, and s u and s d are the source wavelet of modeled and observed three-component data, respectively.

From Eq. 4, we can see that for a shot gather, two terms in the time integration contain the same form of source s i d * s i u ; therefore, we can use the convolution misfit function in EFWI to eliminate the influence of the source wavelet. This is the basic principle of the time domain convolution method to eliminate the influence of wavelets.

For frequency domain EFWI, since the time domain convolution operation is equivalent to the frequency domain multiplication operation, the source free misfit function of EFWI in the frequency domain in the two-dimensional case can be expressed as: E c v 1 m = 1 2 R ∑ i = 1 N s ∑ j = 1 N r U i , k x D i , j x − D i , k x U i , j x U i , k z D i , j z − D i , k z U i , j z T U i , k x D i , j x − D i , k x U i , j x * U i , k z D i , j z − D i , k z U i , j z * , (5) where U i , j x , U i , j z , D i , j x , and D i , j z are the X and Z components of the modeled wavefields and observed wavefields in the frequency domain, respectively, U i , k x , U i , k z , D i , k x , and D i , k z are the X and Z components of reference traces of modeled wavefields and observed wavefields in the frequency domain, respectively, the superscript T is the transpose operator, the superscript * is the conjugate operator, and R is the real part operator.

The seismic wave in the frequency domain can also be regarded as the product of Green’s function and the source wavelet, so there are: D ω = G ω S ω , (6) where D is the frequency domain seismic wave, G is the frequency domain Green’s function, S is the source wavelet in the frequency domain, and ω is the angular frequency.

According to Eq. 6, Eq. 5 can be expressed as: E c v 1 m = 1 2 R ∑ i = 1 N s ∑ j = 1 N r G i , j U x G i , k D x S i D x S i U x − G i , j D x G i , k U x S i D x S i U x G i , j U z G i , k D z S i D z S i U z − G i , j D z G i , k U z S i D z S i U z T G i , j U x G i , k D x S i D x S i U x − G i , j D x G i , k U x S i D x S i U x * G i , j U z G i , k D z S i D z S i U z − G i , j D z G i , k U z S i D z S i U z * , (7) where G i , j U x and G i , j U z are the Green’s functions of modeled X and Z components in the frequency domain, respectively, G i , j D x and G i , j D z are the Green’s functions of observed X and Z components in the frequency domain, respectively, S i , j U x and S i , j U z are the source wavelets of modeled X and Z components in the frequency domain, respectively, and S i , j D x and S i , j D z are the source wavelets of observed X and Z components in the frequency domain, respectively.

From Eq. 7, it can be seen that both sides of each minus have the same form of source terms. By using this misfit function, the influence of wavelet differences between simulated and observed records is eliminated, which is the basic principle of eliminating the influence of wavelets in frequency domain EFWI.

Make R i , j x = U i , k x D i , j x − D i , k x U i , j x and R i , j z = U i , k z D i , j z − D i , k z U i , j z , Eq. 7 can be expressed as: E c v 1 m = 1 2 R ∑ i = 1 N s ∑ j = 1 N r R i , j x R i , j z T ∙ R i , j x * R i , j z * . (8)

By taking the partial derivative of m on both sides of the above equation and using the elastic wave equation to obtain R i , j x and R i , j z , the gradient formula in the frequency domain can be obtained: ∂ E c v 1 m ∂ m = R ∑ i = 1 N s − U T ∂ A ∂ m T A − 1 α − β * , (9) where: A = v p 2 0   v s 2 0   0 0 v s 2 v p 2   0 v p 2 − v s 2 v p 2 − v s 2 0   ω 2 0   0 ω 2 , U = ∂ 2 U x ∂ x 2   ∂ 2 U x ∂ z 2   ∂ 2 U z ∂ x 2   ∂ 2 U z ∂ z 2   ∂ 2 U x ∂ x ∂ z   ∂ 2 U z ∂ x ∂ z   U x   U z T , α = d i , k x R i , 1 x * ⋯ d i , k x R i , j x * ⋯ d i , k x R i , N r x *   d i , k z R i , 1 z * ⋯ d i , k z R i , j z * ⋯ d i , k z R i , N r z * T , β = 0 ⋯ ∑ j = 1 N r d i , k x R i , j x * ⋯ 0   0 ⋯ ∑ j = 1 N r d i , k z R i , j z * ⋯ 0 T .

It can be seen that the gradient of source-independent frequency-domain EFWI is consistent with that of conventional frequency domain EFWI. For the calculation of the frequency domain wave field, DFT operation can be inserted while calculating the time domain wave field to convert the time domain wave field into the frequency domain, thereby achieving source-independent mixed domain EFWI.

The envelope inversion based on Hilbert transform continuously fits the Hilbert envelope of modeled data and the Hilbert envelope of observed data, and takes the velocity corresponding to the smallest fitting difference between the two as the optimal velocity model.

The analytic signal based on Hilbert transform can be expressed as: d ∼ t = d t + i H d t , (10) where d t is the seismic signal, d ∼ t is the analytical seismic signal, i = − 1 , and H is the Hilbert transform, which is defined as follows: H d t = − 1 π K ∫ − ∞ + ∞ d τ t − τ d τ , (11) where K is the Cauchy principal value.

The envelope of the signal can be expressed as: D t = d 2 t + H d t 2   . (12)

By obtaining the envelope of each trace, the envelope spectrum of the shot gather can be obtained. Figure 1 shows the X-component and its envelope spectrum of a shot gather; we can see that the waveform of the seismic record envelope spectrum is smoother, with fewer details in the waveform. The normalized spectrum is shown in Figure 2. It is found that the frequency of the original seismic data is concentrated near the main frequency, while the frequency of the envelope is mainly concentrated in the low-frequency parts below 5 Hz. Therefore, using envelopes containing rich low-frequency information for inversion is beneficial for establishing a better initial model and reducing the probability of EFWI falling into local extrema.

In order to verify the low-frequency information extraction ability of the seismic data envelope, low-frequency components below 3 Hz were filtered out from seismic data, and the normalized spectra of low-frequency missing seismic data and their envelopes were obtained, as shown in Figure 3.

From Figure 3, it can be seen that although low-frequency information is missing from the seismic data, their envelope still contains rich low-frequency components. Therefore, seismic data envelopes based on Hilbert transform have the ability to extract low-frequency information.

Furthermore, the misfit function of elastic wave envelope inversion based on Hilbert transform can be defined as: J e l m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T max D u i , j p − D d i , j p 2 d t , (13) where T max is the record length and D u i , j p and D d i , j p are the pth power of the vector envelope of the modeled and observed data, respectively.

In the two-dimensional case, there are: D u i , j n p = u i , j n 2 + H u i , j n 2 p 2 n = x , z , (14a) D d i , j n p = d i , j n 2 + H d i , j n 2 p 2 n = x , z , (14b) where D u i , j x and D u i , j z are the envelope of the X and Z components of the modeled data, respectively, D d i , j x and D d i , j z are the envelope of the X and Z components of the observed data, respectively, u i , j x and u i , j z are the X and Z components of the modeled data, respectively, d i , j x and d i , j z are the X and Z components of the observed data, respectively, i is the shot index number, and j is the channel index number.

It should be noted that the value of p has a significant impact on the accuracy of envelope inversion. Chi et al. (2014) compared and analyzed the results of envelope inversion under different p-values, and their analysis results showed that a decrease in p-values would weaken the energy of deep reflected waves, thereby increasing the envelope inversion error of deep stratum. We study the first and second envelope inversion methods for elastic waves when p is taken as 1 and 2, respectively.

The gradient of the misfit function J e l over the model parameter m can be expressed as: ∂ J e l ∂ m = p ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T D u i , j p − D d i , j p D u i , j p − 2 u i , j − H D u i , j p − D d i , j p D u i , j p − 2 H u i , j ∂ u i , j ∂ m d t . (15)

In order to compare the formal similarities and differences between the gradient formula of elastic wave envelope inversion and the gradient formula of full waveform elastic wave inversion in the conventional time domain, the partial derivative of Eq. 1 with respect to model parameter m can be obtained as follows: ∂ E ∂ m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T u i , j − d i , j ∂ u i , j ∂ m d t . (16)

Comparing Eq. 15 and Eq. 16, it can be seen that the gradient of elastic wave envelope inversion is consistent in form with the gradient of conventional time domain elastic wave full waveform inversion. Therefore, envelope inversion can be carried out according to the process of time domain EFWI. The difference between the two is that the accompanying source of elastic wave envelope inversion is a function related to the seismic data envelope (Eq. 17), Therefore, it is only necessary to replace the accompanying source of the FWI of elastic waves in the conventional time domain to achieve elastic wave envelope inversion. f i , j e l = p D u i , j p − D d i , j p D u i , j p − 2 u i , j − H D u i , j p − D d i , j p D u i , j p − 2 H u i , j , (17) where f i , j e l is the accompanying source vector of the ith shot and the jth channel. In the two-dimensional case, there are: f i , j e l x = p { D u i , j x p − D d i , j x p D u i , j x p − 2 u i , j x − H D u i , j x p − D d i , j x p D u i , j x p − 2 H u i , j x } , (18a) f i , j e l z = p { D u i , j z p − D d i , j z p D u i , j z p − 2 u i , j z − H D u i , j z p − D d i , j z p D u i , j z p − 2 H u i , j z } , (18b) where f i , j e l x and f i , j e l z are the X and Z components of f i , j e l .

The elastic wave envelope inversion provided above is derived based on the assumption of accurate source wavelets, without considering the different effects of wavelets used in modeled wave fields and observed wave fields. In order to eliminate the influence of source wavelets on elastic wave envelope inversion, we extend the source-independent envelope inversion method to the field of elastic waves, so that the elastic wave envelope inversion can still establish a reliable initial model when the source wavelet is not accurate.

The misfit function of source-independent elastic wave envelope inversion can be expressed as: J c v 1 m = 1 2 ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T E c v i , j u − E c v i , j d 2 d t , (19) where E c v i , j u and E c v i , j d are the envelope vectors of the modeled convolutional wave field and the observed convolutional wave field of the ith shot and the jth channel, respectively.

In the two-dimensional case, Eq. 19 can be written as: E c v i , j u x = u i , j x * d i , k x 2 + H u i , j x * d i , k x 2 , (20a) E c v i , j u z = u i , j z * d i , k z 2 + H u i , j z * d i , k z 2 , (20b) E c v i , j d x = d i , j x * u i , k x 2 + H d i , j x * u i , k x 2 , (20c) E c v i , j d z = d i , j z * u i , k z 2 + H d i , j z * u i , k z 2 . (20d)

It can be seen that the misfit function is composed of the sum of squares of the envelope residuals of the convolution wave fields in the x and z directions. Taking Eq. 20a and Eq. 20c as examples, the seismic wave field in the formula can be expressed in the form of Green’s function and source wavelet convolution, so Eq. 20a and Eq. 20c can be expressed as: E c v i , j u x = u i , j x * d i , k x 2 + H u i , j x * d i , k x 2 = g i , j u x * g i , k d x * s i d x * s i u x 2 + H g i , j u x * g i , k d x * s i d x * s i u x 2 , (20e) E c v i , j d x = d i , j x * u i , k x 2 + H d i , j x * u i , k x 2 = g i , j d x * g i , k u x * s i d x * s i u x 2 + H g i , j d x * g i , k u x * s i d x * s i u x 2 . (20f)

In equations Eq. 20e and Eq. 20f, the simulated convolutional wave field envelope E c v u x and the observed convolutional wave field envelope E c v d x have the same source wavelet s i d x * s i u x . Similarly, both E c v u z and E c v d z have the same source wavelet s i d z * s i u z . Therefore, inversion based on the objective function J c v 1 can eliminate the influence of source wavelets on elastic wave envelope inversion.

The derivative of the misfit function J c v 1 over the model parameter m can be expressed as: ∂ J c v 1 ∂ m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T E c v i , j u − E c v i , j d ∂ E c v i , j u ∂ m − ∂ E c v i , j d ∂ m d t . (21)

Since the observed convolution wave field contains the modeled reference trace wave field, its partial derivative to the model parameters is not zero, so Eq. 21 can be expanded to: ∂ J c v 1 ∂ m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T E c v i , j u − E c v i , j d 1 E c v i , j u u i , j * d i , k d i , k * ∂ u i , j ∂ m + H u i , j * d i , k H * d i , k * ∂ u i , j ∂ m } } d t − ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T E c v i , j u − E c v i , j d 1 E c v i , j d d i , j * u i , k d i , j * ∂ u i , k ∂ m + H d i , j * u i , k H * d i , j * ∂ u i , k ∂ m } } d t . (22)

According to Born’s approximation and convolution theorem (Choi et al., 2011), Eq. 22 can be transformed into: ∂ J c v 1 ∂ m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T W i , x u t ∫ 0 T g x , j u τ − t ∙ f i , j τ d τ d t + ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T W i , x u t ∫ 0 T g x , k u τ − t ∙ f i , k τ d τ d t , (23) where subscript x is the spatial position vector, W i , x u is the simulated multi-component wavefield of the ith shot, g x , j u is the Green’s function of the up-going wave, f i , j = f i , j x , f i , j z is the adjoint source vector of the ith shot and the jth channel, and f i , k = f i , k x , f i , k z is the adjoint source vector of the ith shot and the kth channel. Their calculation formula is as follows: f i , j x = d i , k x ⊗ E c v i , j u x − E c v i , j d x 1 E c v i , j u x u i , j x * d i , k x − H E c v i , j u x − E c v i , j d x 1 E c v i , j u x H u i , j x * d i , k x , f i , j z = d i , k z ⊗ E c v i , j u z − E c v i , j d z 1 E c v i , j u z u i , j z * d i , k z − H E c v i , j u z − E c v i , j d z 1 E c v i , j u z H u i , j z * d i , k z , f i , k x = − d i , j x ⊗ E c v i , j u x − E c v i , j d x 1 E c v i , j d x d i , j x * u i , k x − H E c v i , j u x − E c v i , j d x 1 E c v i , j d x H d i , j x * u i , k x , f i , k z = − d i , j z ⊗ E c v i , j u z − E c v i , j d z 1 E c v i , j d z d i , j z * u i , k z − H E c v i , j u z − E c v i , j d z 1 E c v i , j d z H d i , j z * u i , k z , (24) where ⊗ is the cross correlation operator, ∫ 0 T g x , j u τ − t ∙ f i , j τ d τ and ∫ 0 T g x , k u τ − t ∙ f i , k τ d τ are the adjoint wavefield of the jth channel obtained through f i , j backpropagation and the adjoint wavefield of the reference channel obtained through f i , k backpropagation, respectively.

Eq. 23 and Eq. 24 show that the partial derivative of the misfit function J c v 1 with respect to the model parameter m is still the cross correlation between the forward wavefield and the adjoint wavefield, which is consistent with the gradient formula of the EFWI in the time domain in form. Therefore, in actual inversion, only the adjoint source of the time domain EFWI needs to be replaced by Eq. 24 to obtain the reverse-time propagation wave field, then a zero delay cross correlation can be performed between the reverse-time propagation wavefield and the forward wavefield according to Eq. 19 to obtain the corresponding gradient.

Similarly, the misfit function of source-independent second-order elastic wave envelope inversion can be defined as: J c v 2 m = 1 2 ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T E c v i , j u 2 − E c v i , j d 2 2 d t . (25)

The derivative of J c v 2 with respect to m can be expressed as: ∂ J c v 2 ∂ m = ∑ i = 1 N s ∑ j = 1 N r ∫ 0 T 2 E c v u 2 − E c v d 2 u i , j * d i , k d i , k * ∂ u i , j ∂ m + 2 H u i , j * d i , k H * d i , k * ∂ u i , j ∂ m − 2 E c v u 2 − E c v d 2 d i , j * u i , k d i , j * ∂ u i , k ∂ m + 2 H d i , j * u i , k H * d i , j * ∂ u i , k ∂ m d t , (26) assuming that the adjoint source vectors of misfit function J c v 2 at receiver point j and reference channel k are f i , j ′ and f i , k ′ , respectively. The adjoint source expression in vector form is as follows: f i , j ′ = d i , k ⊗ 2 E c v u 2 − E c v d 2 ∙ u i , j * d i , k − 2 H E c v u 2 − E c v d 2 ∙ H u i , j * d i , k } , (27a) f i , k ′ = − d i , j ⊗ 2 E c v u 2 − E c v d 2 ∙ d i , j * u i , k − 2 H E c v u 2 − E c v d 2 ∙ H d i , j * u i , k } . (27b)

So far, the gradient and adjoint source formulas for the first and second-order envelope inversion of elastic waves independent of wavelets have been obtained.

In terms of inversion strategy, Wang et al. (2016) proposed a multi-step multi-scale strategy for the reconstruction of high-accuracy P- and S-wave velocities, which uses envelope-based EFWI as the first step to obtain the long wavelength components of P- and S-wave velocities, and then uses the above velocities as initial models to obtain the final inversion results through multi-scale inversion. Due to the more natural multi-scale framework of mixed domain EFWI, we combine source-independent elastic wave envelope inversion with source-independent mixed domain EFWI to form a step-by-step multiscale inversion strategy when the source wavelet is inaccurate. The specific process is as follows: firstly, the low-frequency model is obtained by using the first or second-order source-independent elastic wave envelope inversion; then, the above low-frequency mode is used to perform a mixed domain source-independent multi-scale EFWI, and the final inversion results of P- and S-wave velocities are obtained. Theoretically, the strategy can not only avoid the influence of wavelet inaccuracy on the conventional mixed domain EFWI but, on the other hand, using source-independent elastic wave envelope inversion to provide an initial model is beneficial for reducing the probability of EFWI falling into local extrema.

Our goal is to develop a depth domain P- and S-wave interval velocity models building method utilizing common shot gathers of multi-component seismic data. However, it is important to note that our current tests were conducted solely on synthetic data. Given the substantial disparities between synthetic and field data, it is crucial to exercise caution when applying this algorithm to real-world scenarios to address practical problems. The following data processing should be considered when working with field data:

(1) Preprocessing of multi-component seismic data (including denoising and surface consistency processing). The assumption in the algorithm we proposed assumes that the input multi-component seismic data do not contain noise and satisfy the assumption of surface consistency, while field data often do not meet above assumptions. Therefore, before using these algorithms for inversion of field data, it is necessary to perform noise suppression and surface consistency processing on the multi-component shot gathers to ensure that the input data are as close as possible to the above assumptions.

Usually, preprocessed multi-component seismic gathers often lack low-frequency information, which is one of the reasons for the failure of existing EFWI. Therefore, we use Hilbert transform to compensate for the low-frequency information in multi-component seismic data. But we emphasize that the Hilbert transform does not create low-frequency information out of thin air, but rather utilizes existing high-frequency information to predict and supplement low-frequency components. Through Hilbert envelope transformation, we can extract phase information from the observed high-frequency signals and use this phase information to synthesize a complex signal that can fully describe the spectral characteristics of the original signal (including the missing low-frequency components in the original signal). This supplementary low-frequency information is not completely new; it is inferred and restored from existing high-frequency signals through mathematical methods such as frequency domain analysis and phase correction. Therefore, Hilbert envelope transform is a mathematical method that effectively utilizes existing information to infer missing signal components.

When low-frequency time-domain FWI of information is missing in the multi-component seismic data, the EFWI strongly relies on an accurate initial model. When the initial model is poor, the inversion is easily trapped in local extrema. Elastic wave envelope inversion can establish a good initial model when low-frequency information is missing, but when the wavelet is not accurate, elastic wave envelope inversion cannot perform normal velocity model construction. We introduce the convolution-based source-independent method in the field of acoustic FWI into the elastic wave envelope inversion. We derived the gradient and adjoint source formulas for the source-independent elastic wave envelope inversion, and formed a step-by-step inversion strategy by combining the source-independent elastic wave envelope inversion and the source-independent mixed domain EFWI. The model testing results indicate that the source-independent elastic wave envelope inversion method we proposed can establish a good low-frequency model for P- and S-waves when seismic data lacks low-frequency information and the source wavelet is not accurate. Using the low-frequency models as the initial models for source-independent EFWI can effectively alleviate the impact of the “cycle skipping” phenomenon on the accuracy of EFWI.