APP下载

Inversion-based attenuation compensation with dip constraint

2022-06-02XiongLiLiHuoGuoLiHoLiQingLongMeng

Petroleum Science 2022年2期

Xiong M ,Li-Li Huo ,Guo-F Li ,*,Ho Li ,Qing-Long Meng

a School of Geophysics,China University of Petroleum-Beijing,State Key Lab of Petroleum Resources and Prospecting,Key Lab of Geophysical Exploration of CNPC,Changping,102249,Beijing,China

b Exploration and Development Research Institute,Jidong Oilfield,China National Petroleum Corporation,Tangshan,Hebei 063000,China

c Exploration and Development Research Institute,Dagang Oilfield,China National Petroleum Corporation,Tianjin,300280,China

Keywords:Attenuation compensation Instability Inverse scheme Seismic dip Seismic resolution Spatial continuity

ABSTRACT Instability is an inherent problem with the attenuation compensation methods and has been partially relieved by using the inverse scheme.However,the conventional inversion-based attenuation compensation approaches ignore the important prior information of the seismic dip.Thus,the compensated result appears to be distorted spatial continuity and has a low signal-to-noise ratio (S/N).To alleviate this issue,we have incorporated the seismic dip information into the inversion framework and have developed a dip-constrained attenuation compensation (DCAC) algorithm.The seismic dip information,calculated from the poststack seismic data,is the key to construct a dip constraint term.Benefiting from the introduction of the seismic dip constraint,the DCAC approach maintains the numerical stability and preserves the spatial continuity of the compensated result.Synthetic and field data examples demonstrate that the proposed method can not only improve seismic resolution,but also protect the continuity of seismic data.

1.Introduction

The attenuation phenomenon of seismic wave propagation in attenuating media has been acknowledged in theory for many years.Many attenuation models,including the Kolsky-Futterman model (Kolsky,1956;Futterman,1962),the power-law model(Strick,1967),Kjartansson's constant-Q model(Kjartansson,1979),Muller's model(Muller,1983),and the standard linear solid model(Zener,1948),have been established to describe these attenuationdispersion effects.The existence of these effects results in the distortion of a propagating wavelet and the decrease of the seismic resolution(Bickel and Natarajan 1985;Li and Wang,2007;Xie and Yang,2018).To compensate the seismic attenuation and enhance the seismogram resolution,a variety of inverse-Q filtering algorithms have been developed by many researchers(Hargreaves and Calvert,1991;Wang,2006;Oliveira and Lupinacci,2013;Chai et al.,2014;Wang and Chen,2014;Li et al.,2015;Yuan et al.,2016;Morozov et al.,2018).

Amplitude compensation performed in inverse Q filtering is inherently unstable and can cause undesirable noise amplification as well as bandwidth enhancement in the seismic data.Therefore,the current various inverse Q filtering methods are aimed to enhance the bandwidth while suppressing the noise.For example,Robinson (1979) proposes a phase-only inverse Q filtering algorithm which only corrects the phase distortion but neglects the amplitude compensation.Since the amplitude compensation is not considered,this method is unconditionally stable.Wang (2002)deals with the instability of inverse Q filtering by adding a regularization factor to the amplitude compensation operator.Furthermore,Zhang and Ulrych (2007) use the least-squares inversion with Cauchy-Gauss type regularization to stabilize the attenuation compensation solution.Braga and Moraes (2013) accomplish the inverse Q filtering in the wavelet domain and use the L2norm constraint to stabilize the compensation result.Wang et al.(2018)exploit a nonconvex L1-2norm constraint for the stabilized absorption compensation.Nevertheless,the above inversion-based compensation methods ignore an important prior information of the structural orientation or the seismic dip,thus,the invited results may appear an incorrect or distorted spatial continuity(Li and Oldenburg,2000;Lelievre and Oldenburg,2009;Zhang et al.,2013;Gholami,2015;Yuan et al.,2015;Karimi,2015;Wang and Wang,2017;Ma et al.,2020a).

The seismic dip,which reflects the orientations of the seismic events,is an essential attribute in a 2D seismic profile (Marfurt,2006).This information has been widely used in geophysical inversions to produce a result with improved S/N and better continuity (Yao et al.,2020).For example,Li and Oldenburg (2000)incorporate the geological dip information into resistivity inversion and induced-polarization inversion.Clapp et al.(2004)exploit the structural orientation information for the reflection tomography.Hamid and Pidlisecky (2016) use the structural angle information for the seismic impedance inversion.In addition,Zu et al.(2019) introduce the dip field information to the denoising algorithm and further propose a dip-oriented dictionary-learning approach for the random noise attenuation.

In this paper,we incorporate the seismic dip information into the attenuation compensation algorithm,and further propose a dip-constrained attenuation compensation (DCAC) algorithm.Benefiting from the introduction of the dip constraint,the proposed DCAC method will maintain the numerical stability and improve the spatial continuity of the compensation results.The compensation performance of this approach is verified by applying both synthetic and field data examples.

In the following of this paper,firstly,we briefly review the principle of seismic wave propagation in the attenuation medium and introduce the dip-unconstrained attenuation compensation(DUAC) algorithm.Then,we construct a dip constraint term by applying the seismic dip information.Next,we incorporate the dip constraint term into the inversion system and develop a DCAC algorithm.And then we use synthetic experiments and field data to testify the stability and superiority of the proposed DCAC method.Finally,we discuss the influence of the seismic dip estimation,the determination of the regularization parameters and Q analysis before drawing some conclusions.

2.Methodology

2.1.Seismic wave propagation in the attenuation medium

In the elastic medium,the seismogram can be synthesized by the well-known convolution model which can be expressed as,

where t and τ represent the record time,d(t)is the stationary seismic trace,w(t)is the source wavelet,and r(τ)is the reflectivity series.The schematic diagram of the convolution model is shown in Fig.1.As shown in Fig.1,the wavelet w(t)is constant when propagating in the non-attenuated medium.

Fig.1.The schematic diagram of the convolution model:(a) the stationary seismic trace,(b) the source wavelet matrix,(c) the reflectivity series.

When propagating in the attenuation medium,the propagation wavelet is time-varying and then the stationary convolution model(Equation (1)) can be modified as (Wang,2011):

where s(t)is the attenuated seismic trace,and(t,τ)is the timevarying wavelet which can be expressed as (Wang,2002):

According to the convolution theorems of Fourier transform,the frequency-domain expressions of Equation(1)and Equation(2)are,

where D(ω)and S(ω)are the frequency-domain non-attenuated trace and the frequency-domain attenuated trace respectively,W(ω)and W(ω,τ)are the frequency-domain source wavelet and the frequency-domain time-varying wavelet respectively,and R(ω)is the frequency-domain reflectivity sequences.

Substituting Equations (4) and (5) into Equation (6),we have,

Transforming Equation(7)into the time domain,we obtain the following formula,

Fig.2.The schematic diagram of the non-stationary convolution model:(a) the attenuated seismic trace,(b) the time-varying wavelet matrix,(c) the reflectivity series.

where N is the sampling points,s and d represent the discrete attenuated record and the discrete non-attenuated record respectively,and the matrix A stands for the time-domain attenuation matrix due to the earth's Q filtering effects.An illustration of this matrix-vector product is shown in Fig.3.

When considering a 2D seismic data with M traces,Equation(9)can be generalized to a multichannel system(Ma et al.,2020b):

where y and m are the concatenated attenuated data vector and the concatenated non-attenuated data vector respectively,and G is the block diagonal matrix standing for the multichannel attenuation effects.

Fig.3.An illustration of the matrix-vector product (Equation (6)):(a) the attenuated seismic trace,(b) the time-domain absorption matrix,(c) the non-attenuated seismic record.

The main objective of the attenuation compensation is to eliminate the seismic attenuation effects from the recorded seismic data,that is,we desire to obtain the non-attenuated data m by solving Equation (10) with the known (or estimated) Q value and the recorded attenuation data y.However,the attenuation function contains an exponential decay term,which means the block attenuation matrix G is ill-conditioned and the numerical solution of Equation (10) is often unstable.

2.2.Dip-unconstrained absorption compensation

To obtain a numerical stable result,the Tikhonov regularization is often used to establish the following objective function(Ma et al.,2020c):

where λ is a trade-off parameter which controls the relative strength of the constraint term (second term) to the data misfit(first term).Setting the derivative of Equation (11) with respect to m equals to zero,we obtain its solution which can be expressed as:

where the superscript T denotes the transpose.Since this attenuation compensation algorithm lacks of the dip constraint,we refer it to as DUAC method.

2.3.Constructing a dip constraint term

To improve the spatial continuity of the attenuation compensation algorithm,we incorporate the additional dip information into the inversion system.The dip orientation,estimated from the post-stack seismic data,plays an important role in constructing the dip constraint term.Therefore,in this section,we focus on estimating the dip orientation and constructing a dip constraint term.

According to Hamid and Pidlisecky (2016),the seismic dip is defined as the angle between the x-axis and the vector describing the minimum gradient direction (Fig.4).Based on the above definition,we calculate the seismic dip by the following expression:

where θ(t,x)is the seismic dip,dxand dtare the horizontal and vertical derivatives of the seismic data respectively,and Dx(·)and Dt(·)are the horizontal and vertical derivative operators respectively.

After estimating the seismic dip,we need to rotate the horizontal and vertical derivatives in the original coordinate system.And then we obtain two rotation operators,that is,the one parallel to the local dip and the other perpendicular to the local dip.According to Fig.4,the rotated and original derivative operators satisfy the following coordinate transformation formula:

where Dparlis the derivative operator parallel to the local dip direction and Dperpis the derivative operator perpendicular to the local dip direction.The matrix-vector form of Equation(15)can be expressed as,

The derivative parallel to the local dip can be used to enforce the smoothness and continuity of the seismic data in the dip direction,thus,we can construct the following dip constraint term:

Fig.4.An illustration for the seismic dip estimation.The blue curve represents a seismic event.dx and dt are the horizontal and vertical derivatives of the seismic event respectively.θ is the seismic dip.dparl is the derivative which is oriented parallel to the local dip direction and dperp is the derivative perpendicular to the local dip direction.

To better understand this constraint,we first consider a special case of the horizontally layered medium.In this case,the seismic dip is 0 and the operator Dparlwill degenerate to the horizontal derivative matrix Dx.As discussed by Hamid and Pidlisecky(2015),the effect of the horizontal derivative matrix Dxis to minimize the differences between the seismic data in the horizontal direction.This means the lateral constraintcan enforce the lateral continuity of the inverted result m (because the inverted result with the horizontal events can minimize the lateral constraint).Compared with the horizontal derivative matrix Dx,the dip-related derivative matrix Dparlcontains the extra dip information at each sampling point.Therefore,we can find that the effect of the diprelated derivative matrix Dparlis to minimize the differences between the seismic data in the local dip direction.Similarly,the dip constraintcan enforce the spatial continuity of the inverted result m along the dip direction (because the inverted result m with the similar seismic dip field to the dip-related derivative matrix Dparlcan minimize the dip constraint).It can be seen that the core issue of this dip constraint is the reliability of the seismic dip information.Theoretically,the most reliable seismic dip should be calculated from the true inverted result (or the nonattenuated data).In the field data,we do not have the true inverted result (or the non-attenuated data),but only have the attenuated data.Generally,the difference of the dip information between the attenuated and the compensated profiles will not be great.Thus,we can estimate the seismic dip information from the attenuated data.In summary,the role of the dip constraint is to ensure that the inverted results are structurally similar to the attenuated seismic data.

2.4.Dip-constrained attenuation compensation

The DUAC algorithm ignores the seismic dip information,thus,the inverted results may appear an incorrect or distorted spatial continuity when seismic data are contaminated by the random noise.To eliminate this issue,we incorporate the dip constraint terms (Equation (17)) into the inversion system and establish the following objective function:

where μ is the regularization parameter controlling the relative strength of the dip constraint (third term) to the data misfit term(first term).The least-squares solution of this problem is:

Since the proposed attenuation compensation algorithm includes a dip constraint term,we refer it to as DCAC method.It needs to be explained that although we use the matrix inversion to express the numerical results of the DCAC method,we do not directly use the inversion operation in programming(because the inversion operation is computationally expensive).Instead,we use the conjugate gradient method to solve the equivalent equation=GTy.Therefore,the algorithm has high computational efficiency.

The objective function of the DCAC method contains three terms,that is,the data misfit termthe vertical smoothness constraintand the dip constraintThe data misfit term contorls the inversion accuracy.The second term requires the inversion results to be smoothness (for suppressing the amplification of the high frequencies) and the parameter λ contorls the smoothing strength of the compensation results.The parameter λ has a similar action to the stabilization factor in Wang's method(Wang,2002).If the parameter λ is small,more high frequencies (both the effective signal and noise) can be recovered or amplified.In contrast,if the parameter λ is large,the high-frequency recovery or amplification(both the effective signals and noises) will be suppressed.For any frequency component,the magnification of the effective signals and noises is the same.In the proposed method,we have an additional dip constraint.The role of the dip constraint is to ensure that the compensation results are structurally similar to the seismic section before compensation.In other word,the dip constraint term is mainly to enhance the signals consistent with the seismic dip (or suppress the signals inconsistent with the seismic dip).This means the dip constraint can selectively enhance the effective signals while suppress the seismic noises(at the same bandwidth).In summary,the data misfit term is to ensure the accuracy of the inversion results.The second term is mainly to improve the seismic resolution and partially suppress the amplification of the high frequencies (similar to a low-pass filtering).The third term is to further selectively suppress the signals inconsistent with the seismic dip (similar to a dip-filtering).

3.Example

3.1.Synthetic data experiments

To testify the validity and superiority of the DCAC algorithm,we use Marmousi model shown in Fig.5a to conduct the attenuation compensation tests.Assuming that the model has a constant density,we obtain the reflectivity model displayed in Fig.5b.Convolving a 30 Hz Ricker wavelet with the above reflectivity model provides the non-attenuated seismic data(Fig.5c).The nonattenuated data can be served as the reference data to evaluate the compensation performance of Wang's method (Wang,2002),the DUAC method and the DCAC approach.Fig.5d shows the noisy attenuated seismogram with the quality factor Q=40 and contaminated by 20%random noise.This noisy data is exploited to conduct the compensation experiments and used to analyze the performance of Wang's,the DUAC and the DCAC approaches.We find from Fig.5d that the resolution of the attenuated seismic data is decreased due to the seismic attenuation,and the reflection events in the deep are almost drowned in the random noise.To compensate the seismic attenuation effects and enhance seismic resolution,we apply Wang's method,the DUAC approach and the DCAC algorithm to process the attenuated seismic data.

In order to implement the proposed DCAC algorithm,we first calculate the seismic dip by the following steps:1)Computing the vertical and horizontal derivatives of the attenuated seismic data.Note that the derivative calculation will amplify the seismic noise,thus,we estimate the derivatives from the smoothed seismic data.The obtained vertical and horizontal derivatives are shown in Fig.6a and b respectively.2) Estimating the seismic dip at every location using Equation(13)and depicting the result in Fig.6c.3)As we can see,the estimated dips are also effected by the seismic noise,thus,we use spatial Gaussian filter to smooth them and show the smoothed seismic dip section in Fig.6d.Based on the seismic dip information,we then construct the dip constraint and further implement the DCAC algorithm to compensate the seismic attenuation.For quantitatively assessing the compensation results,we use the average correlation coefficient (ACC) as an indicator:

where djandare the reference (non-attenuated) and the compensated traces respectively.

In Wang's method,we choose σ2=0.01 and the inverse-Q filtering result is showed in Fig.7a.We can find that the seismic amplitude is enhanced especially in the mid-and deep-layers,and the spatial continuity of the compensated profile is slightly low.The ACC between it and the non-attenuated seismic data (Fig.5c) is 0.5975.In the DUAC method,we select the regularization parameter λ=0.007 and display the corresponding compensation result in Fig.7b.We see that the compensated result improves the seismic energy and seismic resolution,but the amplification of the high frequencies (both the effective signals and noises) is also evident.The ACC between it and the non-attenuated seismic data(Fig.5c)is 0.6124.To prove the superiority of the DCAC algorithm,we apply the proposed DCAC algorithm to process the attenuated data.In the DCAC algorithm,we fix the parameter λ=0.007 and determine the dip constraint parameter μ=0.1 by trial and error.The DCAC compensation result is shown in Fig.7c.Compared with Wang's and the DUAC compensation results,the DCAC compensation result shows better spatial continuity (see red arrows) and maintains a relatively higher S/N.Its ACC with the ideal data reaches 0.8602.

To clearly view the differences between the reference and the compensated data,we display the seismic traces extracted from the non-attenuated data(Fig.5c)and the compensated data(Fig.7) at CDP=101 and their spectra in Fig.8.Comparing the compensated results with the reference trace,we find that Wang's and the DUAC algorithms overcompensate the attenuated seismic trace due to the amplification of the high frequencies(both the effective signals and noises).The comparison of their spectra further confirms that the proposed method enhances the seismic bandwidth and maintains a relatively higher S/N.

We also use the attenuated data shown in Fig.5d to investigate the influence of the dip constraint parameter on the compensation results.In this experiment,we fix the parameters λ=0.005 and respectively choose μ as 1,0.1,0.01,and 0.001.The corresponding compensated data are displayed in Fig.9a-d.As shown in Fig.9a,when the dip constraint parameter μ is too large,the compensated profile appears to be over-smoothed and the resolution is slightly low.The ACC between it and the non-attenuated seismic data(Fig.5c)is 0.7874.However,when the parameter μ is too small,the amplification of the high frequencies is evident and the spatial continuity of the compensated result is poor (Fig.9d).In this case,the ACC between it and the non-attenuated seismic data (Fig.5c)becomes 0.8103.When the parameter μ is moderate(Fig.9b and c),the compensated results achieve a good balance between the noise suppression and the resolution enhancement.And the ACC between them and the non-attenuated seismic data (Fig.5c) are 0.8477 and 0.8298 respectively.This means the selection of the good regularization parameters is important to the compensation results and we will discuss the parameter selection strategy in the DISCUSSION.

Fig.5.Forward modeling for generating the synthetic data.(a)The velocity model,(b)the reflectivity model,(c)the synthetic noise-free non-attenuated seismic data serving as the reference,(d) the attenuated data with 20% Gaussian noise.

Fig.6.Estimation of the seismic dip.(a)The second derivative of the attenuated seismic data in the vertical direction,(b)the second derivative of the attenuated seismic data in the horizontal direction,(c) the estimated dip,(d) the smoothed seismic dip.

Fig.7.Attenuation compensation for the noisy attenuated seismic data.(a)Wang's method with σ2=0.01,(b)the DUAC algorithm with the regularization parameter λ=0.007,(c)the DCAC algorithm with the regularization parameters λ=0.007 and μ=0.1.

The accuracy and stability of the attenuation compensation results are affected by the noise level.Therefore,we compare the compensation results of Wang's method,the DUAC algorithm and the DCAC approach at different noise levels (5%,10%,15%,20% and 25% random noises).Different methods are used to process the noisy records and the parameters in these methods are fixed at different cases for an objective comparison.In Wang's method,we choose the stabilization factor σ2=0.01.In the DUAC method,the regularization parameter is λ=0.007.In the DCAC method,we fix the parameter λ=0.007 and select the dip constraint parameter μ=0.1.The ACCs of the compensation results are reported in Table 1.We can see that when the noise level is low(5%and 10%),all three methods provide a good result with high ACCs.With the noise level increases,the ACCs of the attenuation compensation results are decrease.Furthermore,for the relatively strong nosie cases(20%and 25%),the ACCs are decrease faster than that in the low noise levels.The main reason is that the high frequencies are partially over-compensated with above parameters for strong nosie cases.And in strong nosie cases,the accuracy of the dip estimation is decreased,which will further influence the accuracy of the proposed method.For the influence of the seismic dip estimation,we will discuss in the DISCUSSION.In addition,we find that the proposed method outperformes other methods in all cases,which further proves the superiority of the proposed method.

Table 1The ACC of the attenuation compensation results using different methods at different noise levels.

3.2.Field data tests

The field data shown in Fig.10a are used to further demonstrate the practicability of the DCAC algorithm.To compensate the energy attenuation and enhance seismic resolution,we exploit Wang's,the DUAC and the DCAC methods to process the original seismic profile.In Wang's method,we select σ2=0.02 and the inverse-Q filtering result is showed in Fig.10b.For the DUAC algorithm,we choose the regularization parameter λ=0.008.While in the DCAC algorithm,we fix the parameter λ=0.008,and then determine μ=0.2 by trial and error.The compensation results using the DUAC and DCAC methods are shown in Figs.10c and d respectively.The averaged amplitude spectra of the three compensation results are displayed in Fig.10e.From the five Figs.,we can draw the following conclusions:Firstly,we find that all three methods compensate the seismic attenuation and enhance the resolution of the raw data.Secondly,compared with Wang's and the DUAC compensation results,the DCAC result tend to have a higher S/R and smoother spatial continuity without losing evident vertical resolution,which indicates that the proposed DCAC algorithm enjoys better compensation performance.Thirdly,we see from the spectra that the frequency components of Wang's and the DUAC results are slightly higher than that of the DCAC result,which may cause by the suppression of the signals inconsistent with the seismic dip in the DCAC algorithm.

We also extract three reference traces from the original and the compensated data at CDP=50,CDP=150,and CDP=250 to further evaluate the compensation details of three compensation results.Note that we only show the seismic records from 2 s to 5 s and the original traces (green) are gained by a factor of 2 for a clearer comparison.As shown in Fig.11,all three methods improve the seismic energy in the deep reflections and they show similar performance in enhancing the vertical resolution.However,the DCAC result shows a relatively small amplitude compensation(see arrows) which may result from the suppression of the additional seismic noises inconsistent with the seismic dip.It indicates that the proposed DCAC method can not only recover the seismic events and improve the seismic resolution,but also keep the seismic noise at an acceptable level.

4.Discussions

Fig.8.The reference and the compensated traces and their spectra.(a)The reference and the compensated traces extracted from the reference data(Fig.5c)and the compensated data (Fig.7) at CDP=101,and (b) their spectra.

Fig.9.The investigation of the influence of the dip constraint parameter μ on the compensation results.We fix the parameters λ=0.005 and respectively choose μ as(a)1,(b)0.1,(c) 0.01,and (d) 0.001.

In the proposed DCAC method,there are two major factors affecting the compensation results,that is,the accuracy of the seismic dip estimation and the determination of the regularization parameters.Firstly,we discuss the influence of the seismic dip estimation.The structural dip information can be derived from many sources (e.g.geologic data,borehole data).In this paper,we derive the dip information directly from the seismic data.In the past few decades,many researchers have proposed or developed a large number of methods to estimate the seismic dip from the seismic data.For example,Bahorich and Farmer(1995)estimate the seismic dip based on the crosscorrelation of a set of windowed seismic data generated by using time lagging between nearby seismic traces.This algorithm is simple and efficient,but the stability of the seismic dip estimation is relatively poor.To improve the stability and accuracy of the dip estimation,many novel methods,including the semblance-based multiple window scanning methods(Marfurt et al.,1998),the structure-tensor-based methods(Fehmers and H¨ocker,2003;Wu,2017),and the plane-wave destruction (Fomel,2002),have been developed in recent years.Considering that the focus of this paper is not on developing a new method for the seismic dip estimation,therefore we select a relatively simple but efficient dip estimation algorithm proposed by Hamid and Pidlisecky (2016).For mid-and large-scale geologic structures,this method can provide a relatively reliable and stable result even if the seismic data are contaminated by some seismic noises.If the S/N of the seismic record is low,we have to increase the seismic noise attenuation before the seismic dip estimation to ensure the stability.In this case,the accuracy of the dip estimation will decrease.Actually,the seismic dip estimation is still an open and undergoing subject in geophysics and we will further research this problem in the future.In addition,we can incorporate any other dip estimation method into the proposed compensation algorithm by replacing the dip estimation method in this paper with the other selected method.This also shows that the proposed DCAC method has good flexibility.

Fig.10.Seismic attenuation compensation for the field data.(a) The field data,(b) the compensation result from Wang's method,(c) the compensation result from the DUAC method,(d) the compensation result from the DCAC approach,(e) the amplitude spectra of the original data,Wang's result,the DUAC result and the DCAC result.

Fig.11.Comparison of the compensation performance using three reference traces extracted from Figs.10a-d at(a)CDP=50,(b)CDP=150 and(c)CDP=250.The original traces(green) are displayed with a gain factor of 2.

Then we discuss the strategy of selecting the good regularization parameters.In the DUAC algorithm,there is a regularization parameter λ to be determined.As discussed by Wang(2011),when the S/N of the seismic data is high,we can choose a small parameter λ to maximally enhance the bandwidth of the seismic data.While the S/N of the seismic data is lower,we should select a relatively larger value for suppressing the amplitude amplification of the high frequencies.This is the general principle for selecting a good regularization parameter λ.To quantitatively determine a suitable regularization parameter λ,Hansen and O'Leary(1993)propose a Lcurve technique,which applies a cross-plot of the data error versus the solution length as a function of λ.A good value for the parameter λ is the one located at the corner of the L-curve.For the proposed DCAC method,there are two regularization parameters,λ and μ,to be determined.The role of the parameter λ is the same as that in the DUAC algorithm,while the role of the parameter μ is to enhance the signals consistent with the seismic dip (or suppress the signals inconsistent with the seismic dip).Until now,the optimization for the hyperparameters functional remains a complicated problem(Du et al.,2018).In this paper,we select them by trial and error,but we apply a relatively elegant strategy.Firstly,we set the dip regularization parameter μ=0 (that means ignoring the dip constraint term)and use L-curve technique to choose a proper parameter λ.Secondly,we fix the parameter λ calculated in the first step and then search a rough value of the parameter μ from a wide range using a large search-step.Finally,we minorly adjust the value of λ and μ by trial and error to determine the most suitable parameters.

In addition,the Q estimation (or Q analysis) is a basic work for the attenuation compensation.Wang (2004) has proposed two methods(the attenuation-based Q analysis and the compensationbased Q analysis)for Q estimation from the reflection seismic data.These methods are the single-trace Q eatimation methods.In the proposed DCAC algorithm,we require a 2D Q model as the input.But the Q analysis processing is not necessarily based on the multichannel model.We can use the attenuation-based Q analysis or the compensation-based Q analysis to estimate Q value trace-bytrace,and then form the single-trace results to a 2D Q model.As showed in Equation(10),the multichannel model can be expressed asActually,this expression can be rewritten asFor the first trace,the Q estimation result will be included in the matrix A1.Similarly,for the second trace,the Q estimation result will be included in the matrix A2.And for the M-th trace,the Q estimation result will be included in the matrix AM.Therefore,the Q analysis is not necessarily on the multichannel base and the compensation results would not be affected while using the single-trace Q eatimation method.

5.Conclusions

In this paper,we propose a novel DCAC method by incorporating the dip constraint into the inversion-based compensation algorithm.The dip regularization term,calculated from the rotated horizontal and vertical derivatives,plays an important role in improving the spatial continuity of the compensated result.Compared with Wang's and DUAC approaches,the proposed DCAC method has advantages in protecting the spatial continuity and maintaining a relatively higher S/R of the compensation data.The superior performance of the proposed approach is demonstrated by using both synthetic tests and field data.The synthetic tests indicate that the proposed method can not only enhance the vertical resolution but also improve the spatial continuity.The field data experiment further confirms the practicability and stability of the DCAC algorithm.

Acknowledgments

We thank the financial support provided by National Natural Science Foundation of China(42074141),the Strategic Cooperation Technology Projects of CNPC and CUPB (ZLZX2020-03),and National Key R &D Program of China (2018YFA0702504).We also gratefully acknowledge the helpful comments from the editors and anonymous reviewers,which greatly improved this manuscript.