Approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS

Abstract

This invention discloses an approximate constant-Q viscous acoustic wave equation and a reverse-time migration imaging method based on GSLS, belonging to the field of exploration geophysics. First, a decoupled attenuation-compensated viscous acoustic wave equation is constructed to separate amplitude attenuation from phase dispersion effects, and the equation is modified to maintain accuracy at arbitrary reference frequencies. Second, the accuracy of the equation is verified through numerical simulation. In practical applications, the model and source parameters are input, and forward modeling is performed to obtain seismic records. Subsequently, forward and reverse extensions are performed on the source and seismic records respectively, the corresponding wavefields are preserved, and attenuation-compensated reverse-time migration imaging is achieved using imaging conditions. The beneficial effects of this invention are that this method not only effectively compensates for amplitude attenuation but also can be solved efficiently using the finite difference method, with computational efficiency superior to the traditional pseudospectral method, significantly improving the quality of the migration imaging profile.

Description

Technical Field

[0001] This invention relates to the field of exploration geophysics, and in particular to an approximate constant Q viscous acoustic wave equation and a reverse time migration imaging method based on GSLS. Background Technology

[0002] Seismic waves are affected by absorption effects when propagating in subsurface media, leading to wavefield energy attenuation and phase dispersion, which manifests as reduced amplitude and waveform distortion on the in-phase axis in seismic records. This phenomenon is particularly pronounced in strongly attenuating formations such as natural gas reservoirs. Without compensation for the attenuation effect, imaging results will exhibit positional shifts and amplitude distortions, significantly reducing seismic imaging quality. Therefore, to obtain high-precision subsurface structural imaging, accurate compensation for wavefield attenuation is essential, necessitating the development of stable and efficient attenuation-compensated reverse time migration imaging techniques capable of achieving this in viscous media.

[0003] Traditional viscous acoustic wave equations containing fractional-order Laplace operators, while having decoupled amplitude attenuation and phase dispersion terms, allow for separate compensation of amplitude attenuation and correction of phase dispersion during reverse-time migration. However, the fractional-order Laplace operator in these equations typically requires pseudospectral methods, involving numerous Fourier transforms, resulting in high computational complexity and low efficiency. In contrast, GSLS-based equations can be solved using various time-domain numerical methods such as the finite difference method, offering higher computational efficiency. However, most existing models use reference frequencies... Fixed as Accurate characterization of attenuation and dispersion effects can only be achieved at this specific reference frequency. Once the reference frequency is flexibly selected according to the actual exploration data, the accuracy of the equation will decrease significantly. At the same time, the amplitude attenuation and phase dispersion terms in this type of equation are essentially coupled and cannot handle these two effects independently. Therefore, it is not suitable for Q-compensated reverse time migration imaging.

[0004] How to develop a method that operates at arbitrary reference frequencies? A reliable, stable, and efficient Q-compensated reverse time migration method is one of the urgent problems to be solved in exploration geophysics. Summary of the Invention

[0005] To address the aforementioned technical problems, this invention discloses an approximate constant Q viscous acoustic wave equation and a reverse time migration imaging method based on GSLS.

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

[0007] The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS include the following steps:

[0008] S1. Constructing a viscous acoustic wave equation with attenuation compensation based on GSLS;

[0009] S2. Modify the equation to make it accurate at any reference frequency;

[0010] S3. Numerical simulation is used for accuracy analysis;

[0011] S4. Input the model and source parameters;

[0012] S5. Forward modeling yields earthquake records;

[0013] S6. Input the positive extension of the seismic source and save the wave field at the forward propagation time.

[0014] S7. Input the reverse extension of the seismic record and save the wavefield at the time of reverse propagation.

[0015] S8. Apply imaging conditions to perform attenuation-compensated reverse time migration imaging.

[0016] Furthermore, in step S1, during the process of separating the amplitude attenuation and phase dispersion terms, the traditional first-order viscous acoustic wave equation in the frequency domain is approximated to obtain a new first-order viscous acoustic wave equation in the frequency domain and its corresponding time domain expression. The new viscous acoustic wave equation contains factors that control the amplitude attenuation and phase dispersion, so equations that control the amplitude attenuation and phase dispersion effects can be obtained respectively. After a series of mathematical operations, a decoupled first-order viscous acoustic wave equation is constructed. The equation does not contain complex terms and can be solved directly using the finite difference method.

[0017] The traditional expression for the first-order viscous acoustic wave equation in the time domain based on GSLS is as follows:

[0018] ;

[0019] In the equation, For stress wave field, and For the density and velocity of the medium, and These represent the velocity components in the horizontal and vertical directions, respectively. For stress relaxation time, For strain relaxation time, both are used to describe the decay characteristics of viscous media, and R is a memory variable.

[0020] Taking the Fourier transform of the above first-order time-domain viscous acoustic wave equation, we can obtain the expression for the traditional first-order frequency-domain viscous acoustic wave equation as follows:

[0021] .

[0022] when When, the imaginary part in the above equation:

[0023] .

[0024] Given that Q is the ratio of the real part to the imaginary part, we can deduce that the real part is... .

[0025] Through the above derivation, the new first-order viscous acoustic wave equation in the frequency domain is obtained as follows:

[0026] ;

[0027] in, For the stress wave field in the frequency domain, and Let represent the spatial derivatives of the particle velocity, and K be the bulk modulus of the medium. Q is a dimensionless parameter calculated from stress relaxation time and strain relaxation time, and Q is the quality factor, the magnitude of which represents the degree of attenuation of seismic wave energy by the medium.

[0028] A new first-order time-domain viscous acoustic wave equation is obtained through inverse Fourier transform:

[0029] ;

[0030] In the newly introduced first-order viscous acoustic wave equation To control the phase dispersion operator, and To control the amplitude decay operator, we remove the amplitude decay operator from the complete equation, resulting in a new wave equation dominated by phase dispersion:

[0031] ;

[0032] Removing the operator controlling phase dispersion from the complete equation yields a new wave equation dominated by amplitude decay:

[0033] ;

[0034] However, the new wave equation dominated by amplitude decay contains complex terms, making direct solution very difficult. To address this issue, we subtract the phase dispersion operator from the initial first-order time-domain viscous acoustic wave equation to derive a new wave equation decoupled from phase dispersion and amplitude decay, without complex terms:

[0035] ;

[0036] Adding the expressions on the right-hand sides of the wave equations dominated by phase dispersion and amplitude decay, and then subtracting the operator controlling the propagation of the wave equation, yields a new wave equation decoupled from phase dispersion and amplitude decay, without complex terms:

[0037] ;

[0038] Furthermore, in step S2, the dispersion-dominant, attenuation-dominant, and decoupled viscous acoustic wave equations obtained in step S1 are modified so that these equations are not only valid at the reference frequency. It is accurate in time and accurate at any reference frequency;

[0039] First, the modified equation for viscous acoustic waves containing only phase dispersion is:

[0040] ;

[0041] Similarly, the equation for viscous acoustic waves with amplitude attenuation is modified:

[0042] ;

[0043] Adding the two equations above and subtracting the operator controlling the propagation of seismic waves, we obtain the modified equation for viscous acoustic waves with decoupled amplitude attenuation and phase dispersion:

[0044] .

[0045] Furthermore, in step S3, during the accuracy analysis of the numerical simulation, models with different Q values ​​are used for numerical simulation. The simulation results of the acoustic wave equation, the dispersion-dominated equation, the attenuation-dominated equation, and the complete viscous acoustic wave equation are compared to analyze the accuracy of the equations under different Q values.

[0046] The phase dispersion-dominated equation for viscous acoustic waves:

[0047] ;

[0048] The equation for viscous acoustic waves dominated by amplitude attenuation:

[0049] ;

[0050] Complete equation for viscous acoustic waves:

[0051] .

[0052] Furthermore, in step S6, during the step of inputting the earthquake source forward extension and saving the wavefield at the forward propagation time, the amplitude attenuation term in the newly obtained first-order viscous acoustic equation is changed to compensate for amplitude attenuation. A regularization factor is introduced to ensure the stability of the equation. The earthquake source is then input for forward modeling, as shown in the following formula:

[0053] ;

[0054] ;

[0055] In the formula, For any given time, the propagating wavefield data, Let Q be the velocity vector and Q be the quality factor. and In the 2D plane direction and Direction, t is time, The epicenter was located there.

[0056] Furthermore, in step S7, during the step of inputting the reverse continuation of the seismic record and saving the wavefield at the time of reverse propagation, the simulated reverse propagation wavefield continuation of the seismic record is input into the compensated first-order viscous acoustic wave equation, as shown in the following formula:

[0057] ;

[0058] ;

[0059] In the formula, For any given time, the propagating wavefield data, For earthquake records.

[0060] Further, in step S8, imaging conditions are applied to the forward propagation wavefield and reverse propagation wavefield saved in steps S6 and S7 to obtain the offset image, formula:

[0061] ;

[0062] In the formula, The offset image of the j-th shot. To correspond to the forward propagation wave field of each shot, This corresponds to the reverse propagation wave field for each shot.

[0063] The beneficial effects of this invention are as follows: Based on the Generalized Standard Linear Solid (GSLS) model, a novel decoupled first-order viscous acoustic wave equation is proposed. This equation maintains precise decoupling of amplitude attenuation and phase dispersion at any reference frequency. It is efficiently solved using the finite difference method (FDM), achieving independent compensation for amplitude attenuation and phase dispersion. This method significantly improves the imaging accuracy and computational efficiency of reverse time migration by separating attenuation and dispersion operators and combining regularization terms to suppress instability of high-frequency components. Verification using the BP gas cloud model and the Marmousi model shows that this invention can effectively recover the deep signal energy and phase consistency in strongly attenuated regions. At the same time, compared with the wave equation solved by the traditional pseudospectral method (PSM), the computation time is significantly reduced. In addition, this method can be extended to multiple SLS mechanisms, which is more in line with the constant Q theory of actual media.

[0064] This invention solves the problems of low compensation efficiency, poor stability, and applicability only at a fixed reference frequency in the traditional viscous acoustic wave equation in Q-RTM, providing efficient and reliable technical support for industrial-grade three-dimensional attenuation compensation imaging. Attached Figure Description

[0065] Figure 1 This is a flowchart of the present invention;

[0066] Figure 2 This is a comparison diagram of wavefields obtained from forward modeling of different equations in Example 1;

[0067] Figure 3 In Example 1 Comparison of amplitude over time obtained from simulations using different equations at different Q values;

[0068] Figure 4 In Example 1 A comparison of amplitude versus frequency obtained from simulations using different equations at different Q values;

[0069] Figure 5 In Example 1 Comparison of amplitude over time obtained from simulations using different equations at different Q values;

[0070] Figure 6 In Example 1 A comparison of amplitude versus frequency obtained from simulations using different equations at different Q values;

[0071] Figure 7 In Example 1 Comparison of amplitude over time obtained from simulations using different equations at different Q values;

[0072] Figure 8 In Example 1 A comparison of amplitude versus frequency obtained from simulations using different equations at different Q values;

[0073] Figure 9 a represents the BP cloud velocity model in Example 1; Figure 9 b represents the BP gas cloud Q model in Example 1;

[0074] Figure 10 a represents the seismic record of shot 45 obtained from the acoustic wave equation simulation in Example 1. Figure 10 b represents the seismic record of shot 45 obtained from the acoustic wave equation simulation in Example 1;

[0075] Figure 11 a is a conventional acoustic RTM image using unattenuated data from the BP gas cloud model in Example 1; Figure 11 b is a conventional acoustic RTM image using attenuation data in the BP gas cloud model of Example 1; Figure 11 c represents the RTM image compensated by attenuation data Q in Example 1 using the BP cloud model;

[0076] Figure 12This is a comparison of single-channel images of different RTM models in Example 1 of the BP cloud model;

[0077] Figure 13 This is a comparison of wavenumber spectra of different RTM images of the BP cloud model in Example 1;

[0078] Figure 14 a represents the Marmousi velocity model in Example 2. Figure 14 b represents the Marmousi Q model in Example 2;

[0079] Figure 15 a represents the seismic record of shot 45 obtained from the acoustic wave equation simulation; Figure 15 b represents the seismic record of shot 45 obtained from the simulation of the viscoacoustic wave equation;

[0080] Figure 16 a is a conventional acoustic RTM image using unattenuated data from the Marmousi model in Example 2. Figure 16 b is a conventional acoustic RTM image using attenuation data in the Marmousi model of Example 2. Figure 16 c represents the RTM image compensated using attenuation data Q in the Marmousi model of Example 2;

[0081] Figure 17 This is a comparison of single-channel images of different RTMs in the Marmousi model in Example 2;

[0082] Figure 18 This is a comparison of wavenumber spectra of different RTM images of the Marmousi model in Example 2;

[0083] Figure 19 a represents the actual seismic record of shot 30 in Example 2. Figure 19 b represents the actual seismic record of shot 75 in Example 2;

[0084] Figure 20 'a' represents the actual data speed model in the application example. Figure 20 b represents the actual data Q model in the application example;

[0085] Figure 21 a is a conventional acoustic RTM image of actual data in an application example. Figure 21 b represents the actual Q-compensated RTM image in the application example;

[0086] Figure 22 a is Figure 21 a. The enlarged portion of the area within the red box. Figure 22 b is Figure 21 b. The enlarged portion of the area within the red box;

[0087] Figure 23 a is Figure 21a. The enlarged portion of the area within the blue box. Figure 23 b is Figure 21 b. The enlarged area within the blue box. Detailed Implementation

[0088] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of them. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0089] Example 1

[0090] This invention discloses an approximate constant Q viscous acoustic wave equation and a reverse time migration imaging method based on GSLS, such as... Figure 1 The specific steps are as shown:

[0091] S1. Construct the viscous acoustic wave equation based on GSLS.

[0092] Based on the GSLS derivation, a decoupled first-order viscous acoustic wave equation is constructed. Since the equation does not contain complex terms, it can be solved directly using the finite difference method.

[0093] .

[0094] S2. Modify the equation to make it accurate at any reference frequency.

[0095] First, the modified equation for viscous acoustic waves containing only phase dispersion is:

[0096] .

[0097] Similarly, the equation for viscous acoustic waves, which only involves amplitude attenuation, is modified as follows:

[0098] .

[0099] Adding the two equations above and then subtracting the operator controlling the propagation of seismic waves, we can obtain the modified equation for viscous acoustic waves with decoupled amplitude attenuation and phase dispersion.

[0100] .

[0101] Extending this to multiple SLSs yields:

[0102] .

[0103] S3. Numerical simulation is used for accuracy analysis.

[0104] A uniform model with a 601×601 grid size and a spatial sampling interval of 6.0 m was set up. The source was a Ricker wavelet with a dominant frequency of 20 Hz. The total simulation duration was 2 s, with a temporal sampling interval of 0.5 ms. Forward simulations were performed using the acoustic wave equation, the phase dispersion-dominated viscous acoustic wave equation, the amplitude attenuation-dominated viscous acoustic wave equation, and the complete viscous acoustic wave equation. The wavefield snapshots from each simulation were compared and analyzed.

[0105] The acoustic wave equation used in the simulation:

[0106] .

[0107] The phase dispersion-dominated equation for viscous acoustic waves:

[0108] .

[0109] The equation for viscous acoustic waves dominated by amplitude attenuation:

[0110] .

[0111] Complete equation for viscous acoustic waves:

[0112] .

[0113] Compare the wave field snapshots simulated by the four equations, such as Figure 2 As shown, the wave field obtained by the dispersion-dominated viscous acoustic wave equation has a significantly earlier wavefront compared to the wave field obtained by the acoustic wave equation, but the amplitude is roughly the same. This indicates that the dispersion-dominated viscous acoustic wave equation causes a phase shift but does not result in amplitude loss. In contrast, the wave field obtained by the attenuation-dominated viscous acoustic wave equation has the same wavefront travel time as the wave field obtained by the acoustic wave equation, but the amplitude is significantly attenuated. The wave field obtained by forward modeling the complete viscous acoustic wave equation exhibits both amplitude attenuation and phase shift, demonstrating that the newly derived first-order viscous acoustic wave equation, decoupled from dispersion and attenuation, can accurately simulate the amplitude attenuation and phase dispersion phenomena of seismic wave propagation in subsurface media.

[0114] Next, a uniform model with a size of 601×601 was used, with a horizontal and vertical grid spacing of 8m. The model velocity was 3500m / s, and the density was 2g / cm^3. The time sampling interval was 0.5ms, and the duration was 2.2s. A Ricker wavelet with a center frequency of 25Hz was selected as the source and placed at the grid point (301, 301) of the model. Reference frequency. Select as By specifying different Q values, forward modeling simulations were performed using the acoustic wave equation, the dispersion-dominated viscous acoustic wave equation, the attenuation-dominated viscous acoustic wave equation, and the complete viscous acoustic wave equation, respectively. The simulation results are as follows: Figure 3 and Figure 4 As shown, the curves of the acoustic wave equation and the dispersion-dominant equation almost completely overlap, indicating that the dispersion-dominant equation only affects the phase of the seismic wave but does not change the magnitude of the seismic wave amplitude. Meanwhile, the curves corresponding to the complete viscous acoustic wave equation and the attenuation-dominant equation almost perfectly match, indicating that the attenuation-dominant viscous acoustic wave equation only simulates the amplitude attenuation effect of seismic waves in the subsurface medium and does not cause changes in the propagation time of the seismic wave.

[0115] To further verify the suitability of the proposed equation, the remaining model parameters and simulation conditions were kept constant, and only the value of the reference frequency was changed. Forward modeling was performed again using the acoustic wave equation, the dispersion-dominated viscosonic wave equation, the attenuation-dominated viscosonic wave equation, and the complete viscosonic wave equation constructed in this paper. The corresponding simulation results are shown in Figures 5–8. As can be seen from the figures, the waveforms of the acoustic wave equation and the dispersion-dominated equation still highly overlap, showing only a phase difference while the amplitude remains consistent. The waveforms of the complete viscosonic equation and the attenuation-dominated equation also match well, reflecting only amplitude attenuation while the propagation time shows no significant deviation. This is completely consistent with the simulation law under a fixed reference frequency, fully demonstrating that the viscosonic wave equation proposed in this paper has high simulation accuracy and stability under any reference frequency value.

[0116] S4. Input the model and source parameters.

[0117] This embodiment sets up a BP cloud model, such as Figure 9 As shown in a and 9b, the model has a grid size of 398×160, a spatial sampling interval of 10m, a total of 90 shots, each shot spaced 40m apart, a source with a main frequency of 25Hz Ricker wavelet, a total simulation duration of 2s, and a time sampling interval of 0.5ms.

[0118] S5. Seismic records are obtained through forward modeling.

[0119] This embodiment uses the acoustic wave equation and the viscous acoustic wave equation to perform forward modeling to obtain seismic records. Figure 10 (a) and Figure 10 (b) Represents the seismic record of shot 45 obtained by simulation of the acoustic wave equation and the viscous acoustic wave equation of the BP gas cloud model, respectively.

[0120] S6. Input the positive extension of the seismic source and save the wave field at the forward propagation time.

[0121] In the step of inputting the seismic source for forward continuation and saving the wavefield at the forward propagation time, the amplitude attenuation term in the newly obtained first-order viscous acoustic equation is modified to compensate for amplitude attenuation. A regularization factor is introduced to ensure the stability of the equation. The seismic source is then input for forward modeling, and the formula is as follows:

[0122] ;

[0123] ;

[0124] In the formula, For any given time, the propagating wavefield data, Let Q be the velocity vector and Q be the quality factor. and In the 2D plane direction and Direction, t is time, The epicenter;

[0125] S7. Input the reverse extension of the seismic record and save the wavefield at the time of reverse propagation.

[0126] In the step of inputting the attenuated seismic record for reverse continuation and saving the wavefield at the reverse propagation time, the seismic record is input into the compensated first-order viscous acoustic wave equation to simulate the reverse propagation wavefield continuation, as shown in the following formula:

[0127] ;

[0128] ;

[0129] In the formula, For any given time, the propagating wavefield data, For earthquake records.

[0130] S8. Apply imaging conditions to perform attenuation-compensated reverse time migration imaging.

[0131] Imaging conditions are applied to the forward and reverse propagation wavefields saved in steps S6 and S7 to obtain the offset image. The formula is as follows:

[0132] ;

[0133] In the formula The offset image of the j-th shot. To correspond to the forward propagation wave field of each shot, This corresponds to the reverse propagation wave field for each shot.

[0134] Figure 11 'a' represents the result of conventional acoustic RTM on the unattenuated data, which is selected as the reference profile for comparison. Figure 11 b represents the result of conventional acoustic RTM using attenuation data. Figure 11 c represents the reverse-time offset result using attenuation data for Q-compensation. At the same display scale, compared to the reference profile, the uncompensated profile... Figure 11 The amplitude of b is relatively weak, especially in deep regions where imaging is blurry, while the profile after attenuation compensation... Figure 11c. Amplitude information is recovered, and deep imaging results are relatively clear. Single-channel records are extracted from the three offset profiles, such as... Figure 12 To further analyze this, the corresponding wavenumber spectrum is shown below. Figure 13 As shown, the red, blue, and black curves represent the uncompensated offset profile, the compensated offset profile, and the single-channel curve on the acoustic offset profile, respectively. By comparing the red and black curves, there is obvious amplitude energy attenuation and phase misalignment on the uncompensated offset profile. After attenuation compensation, the blue and black curves have basically the same amplitude, and the phase misalignment is also corrected. This clearly proves that attenuation-compensated counter-time offset can compensate for energy attenuation and correct phase dispersion.

[0135] Example 2

[0136] Steps S1-S3 are the same as in Example 1.

[0137] Step S4: Input model and source parameters

[0138] This embodiment sets up the Marmousi model, such as Figure 14 As shown in a and 14b, the model has a grid size of 661×201, a spatial sampling interval of 6m, a total of 90 shots, each shot spaced 42m apart, a source with a main frequency of 25Hz Ricker wavelet, a total simulation duration of 2s, and a time sampling interval of 0.5ms.

[0139] S5. Seismic records obtained from forward modeling.

[0140] This embodiment uses the acoustic wave equation and the viscous acoustic wave equation to perform forward modeling to obtain seismic records, such as... Figure 15 (a) and Figure 10 (b) Represents the seismic record of shot 45 obtained by simulation of the Marmousi model acoustic wave equation and the viscous acoustic wave equation, respectively.

[0141] S6-S8 are the same as in Example 1.

[0142] Figure 16 The results of conventional acoustic RTM are used as a reference profile for comparison. Figure 16 b represents the result of acoustic RTM using attenuation data, while Figure 16 c represents the result of Q-compensated RTM using attenuation data. Comparative analysis shows that the uncompensated profile... Figure 16 b shows a significant decrease in amplitude relative to the reference profile. Single-channel data extracted from the three profiles show that... Figure 17 and obtain Figure 17 Corresponding wavenumber spectrum Figure 18The comparative results show that the uncompensated profile, compared with the acoustic imaging profile, not only exhibits significant amplitude differences but also significant phase displacement, and the deep reflected signal energy is weaker. In contrast, the attenuation-compensated profile shows a high degree of agreement with the curves extracted from the acoustic imaging profile, indicating that both amplitude and phase characteristics are effectively recovered. This result verifies that the first-order viscous acoustic wave equation proposed in this paper can effectively compensate for the amplitude attenuation and phase dispersion effects generated during the propagation of seismic waves in the subsurface medium, thereby significantly improving the quality of the migration imaging profile.

[0143] Application examples

[0144] The specific steps for applying the method of this invention to actual seismic data are as follows:

[0145] S1-S3 are the same as in Example 1.

[0146] S4. Input the model and source parameters.

[0147] The work area extends 15.61 km laterally and 6.01 km vertically, with the ground observation system comprising 82 non-uniformly distributed shot points. The velocity and Q-model of the actual seismic data are shown below. Figure 20 As shown in a and 20b, the velocity values ​​are distributed between 1600 and 6900 m / s, and the Q value varies from 47.2 to 235.7. The model has a grid size of 1561×601, a spatial sampling interval of 10m, a source with a main frequency of 20Hz Ricker wavelet, a total simulation duration of 2s, and a time sampling interval of 0.5ms.

[0148] S5. Seismic records are obtained through forward modeling.

[0149] The actual seismic records were obtained directly from the field. Figure 19 The images show the actual data from the seismic records of the 30th and 75th shots in the field.

[0150] S6-S8 are the same as in the embodiment.

[0151] Figure 21 a represents the result of conventional acoustic RTM imaging. Figure 21 b represents the Q-compensated RTM imaging result. For Figure 21 The areas marked by the red and blue boxes in the middle are enlarged for display, such as... Figure 22 and Figure 23 As shown in the arrow-indicated area, the stratigraphic contact relationship in the offset profile after attenuation compensation is clearer, indicating that the proposed attenuation compensation method can effectively restore the characteristics of deep strata and improve the quality of imaging profiles. The actual data processing results intuitively verify the effectiveness of the method.

[0152] This invention reconstructs the first-order viscous acoustic wave equation, decoupling and separating the amplitude attenuation and phase dispersion terms to construct an attenuation-compensated viscous acoustic wave equation. The constructed equation is then modified to maintain accuracy at any reference frequency. Simultaneously, it uses the finite difference method for direct solution, avoiding the complex pseudospectral calculations required for traditional fractional-order equations. The equation is extended to a form based on multiple SLSs to improve the approximation accuracy of constant Q over a wide frequency band, better reflecting the attenuation characteristics of actual strata. A regularization factor is introduced during the reverse-time migration process to suppress high-frequency components, ensuring imaging stability. This method not only overcomes the limitations of coupling amplitude attenuation and phase dispersion terms, easily achieving Q-compensated reverse-time migration imaging, but also uses the finite difference method for solution, resulting in higher computational efficiency than traditional equations solved using pseudospectral methods, and significantly improving the quality of the migrated imaging profile.

[0153] Of course, the above description is not intended to limit the present invention, and the present invention is not limited to the examples given above. Any changes, modifications, additions or substitutions made by those skilled in the art within the scope of the present invention should also fall within the protection scope of the present invention.

Claims

1. An approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS, characterized in that, Includes the following steps: S1. Constructing a viscous acoustic wave equation with attenuation compensation based on GSLS; S2. Modify the equation to make it accurate at any reference frequency; S3. Numerical simulation is used for accuracy analysis; S4. Input the model and source parameters; S5. Forward modeling yields earthquake records; S6. Input the positive extension of the seismic source and save the wave field at the forward propagation time. S7. Input the reverse extension of the seismic record and save the wavefield at the time of reverse propagation. S8. Apply imaging conditions to perform attenuation-compensated reverse time migration imaging.

2. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 1, characterized in that, In step S1, during the process of separating the amplitude attenuation and phase dispersion terms, an approximation is made to the traditional first-order viscous acoustic wave equation in the frequency domain to obtain a new first-order viscous acoustic wave equation in the frequency domain and its corresponding time domain expression. The traditional expression for the first-order time-domain viscous acoustic wave equation is: ； in, For stress wave field, , These represent the density and velocity of the medium, respectively. and These represent the velocity components in the horizontal and vertical directions, respectively. For stress relaxation time, R is the strain relaxation time, and R is the memory variable; By performing a Fourier transform on this traditional first-order viscous wave equation, we can approximate the equation in the frequency domain: in, Angular frequency; The new first-order viscous acoustic wave equation in the frequency domain is: ； in, For the stress wave field in the frequency domain, and Let represent the spatial derivatives of the particle velocity, and K be the bulk modulus of the medium. Let Q be a dimensionless parameter, and let Q be the quality factor. A new first-order time-domain viscous acoustic wave equation is obtained through inverse Fourier transform: ； in, To control the phase dispersion operator, Operators for controlling amplitude decay; Removing the operator controlling amplitude decay from the complete equation yields a new wave equation dominated by phase dispersion: ； Removing the operator controlling phase dispersion from the complete equation yields a new wave equation dominated by amplitude decay: ； By subtracting the phase dispersion operator from the traditional first-order time-domain viscous acoustic wave equation, a new wave equation decoupled from phase dispersion and amplitude attenuation without complex terms is derived as follows: ； Adding the expressions on the right-hand sides of the wave equations dominated by phase dispersion and amplitude decay, and then subtracting the operator controlling the propagation of the wave equation, yields a new wave equation decoupled from phase dispersion and amplitude decay, without complex terms: 。 3. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 2, characterized in that, In step S2, the dispersion-dominant, attenuation-dominant, and decoupled viscous acoustic wave equations obtained in step S1 are modified so that these equations are not only valid at the reference frequency. It is accurate in time and accurate at any reference frequency; The revised equation for viscous acoustic waves containing only phase dispersion is: ； in, For reference frequency, Main frequency; The equation for viscous acoustic waves, which only includes amplitude attenuation, is modified as follows: ； Adding the two equations above and subtracting the operator controlling the propagation of seismic waves, we obtain the modified equation for viscous acoustic waves with decoupled amplitude attenuation and phase dispersion: 。 4. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 3, characterized in that, In step S3, during the accuracy analysis of the numerical simulation, models with different Q values ​​are used for numerical simulation. The simulation results of the acoustic wave equation, the dispersion-dominated equation, the attenuation-dominated equation, and the complete viscous acoustic wave equation are compared to analyze the accuracy of the equations under different Q values. Simultaneously, different reference frequencies are set. The values ​​were used to verify that the proposed equations were accurate under different exponents; The phase dispersion-dominated equation for viscous acoustic waves: ； The equation for viscous acoustic waves dominated by amplitude attenuation: ； Complete equation for viscous acoustic waves: 。 5. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 4, characterized in that, In step S6, during the process of inputting the earthquake source forward continuation and saving the wavefield at the forward propagation time, the amplitude attenuation term in the newly obtained first-order viscous acoustic equation is modified, and a regularization factor is introduced to ensure the stability of the equation. The earthquake source is then input for forward modeling. The formula is: ； ； In the formula, For any given time, the propagating wavefield data, Let Q be the velocity vector and Q be the quality factor. and In the 2D plane direction and Direction, t is time, The epicenter was located there.

6. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 5, characterized in that, In step S7, during the step of inputting the reverse continuation of the seismic record and saving the wavefield at the time of reverse propagation, the simulated reverse propagation wavefield continuation of the seismic record is input into the compensated first-order viscous acoustic wave equation. The formula is: ； ； In the formula, For any given time, the propagating wavefield data, For earthquake records.

7. The approximate constant Q viscous acoustic wave equation and reverse time migration imaging method based on GSLS as described in claim 6, characterized in that, In step S8, imaging conditions are applied to the forward propagation wavefield and reverse propagation wavefield saved in steps S6 and S7 to obtain the offset image, using the formula: ； In the formula, The offset image of the j-th shot. To correspond to the forward propagation wave field of each shot, This corresponds to the reverse propagation wave field for each shot.

Images