Angle domain adaptive local optimization reverse time migration method, electronic device and medium

Abstract

The application discloses an angle domain adaptive local optimization reverse time migration method, an electronic device and a medium. The method can comprise: extrapolating a shot point and a receiver point reverse time wave field, constructing an energy target functional under gradient constraint on each time slice; iteratively solving a shot point end and a receiver point end wave field vector direction at a current time; constructing an adaptive threshold value through global amplitude scanning, screening to generate an effective point table at the current time, and calculating reflection opening angles or azimuth angles of effective points in the table; performing interpolation calculation of the reflection opening angles or the azimuth angles within a set range of each effective point in the table; and applying an azimuth angle domain cross-correlation imaging condition to determine a gather at an underground imaging point. The application solves the precision problem of the angle domain reverse time migration common imaging point gather, makes the angle domain reverse time migration common imaging point gather result more accurate and reliable, and is beneficial to subsequent interpretation processing.

Description

Technical Field

[0001] This invention relates to the field of seismic migration imaging technology, and more specifically, to an angle-domain adaptive local optimization reverse time migration method, electronic equipment, and medium. Background Technology

[0002] As the difficulty of oil and gas exploration and development in my country gradually increases, exploration targets are shifting towards complex exploration areas such as fractured and lithologic reservoirs. Researching high-precision and highly adaptable depth migration algorithms is crucial for subsequent seismic interpretation, attribute analysis, and reservoir prediction. Reverse-time migration (RTM) technology based on wave equations can meet the high-precision and highly adaptable imaging requirements of complex reservoirs, and the development of graphics card computing technology has led to its widespread application in practical production. Angle-domain common imaging point gathers are a high-dimensional representation of RTM imaging results, serving as a vital bridge between pre-stack seismic data and subsequent seismic interpretation and reservoir description. On the one hand, they provide more accurate residual time differences for tomographic inversion, adapting to fine-grained velocity model updates; on the other hand, they provide amplitude and phase information reflecting subsurface lithological characteristics. However, efficiently and robustly extracting angle-domain common imaging point gathers during reverse-time migration remains a challenge.

[0003] Both ray migration (Kirchhoff, ray beam, etc.) and wave equation migration can extract common imaging point gathers in the angle domain and represent incident wave field information at different angles at underground reflection points. However, due to the inherent limitations of ray theory, the imaging gathers extracted by ray migration have serious artifacts. In contrast, RTM based on the wave equation can image arbitrarily complex waves such as multiple waves and reflection waves, and the generated angle domain imaging gathers have fewer artifacts and better continuity.

[0004] Algorithm complexity is a key factor affecting the efficient extraction of angular domain common imaging point gathers (ADCIGs) by Real-Time Migration (RTM). Existing algorithms commonly used in RTM fall into two categories: local plane wave decomposition methods and wavefield vector direction estimation methods. Local plane wave decomposition methods directly acquire ADCIGs using RTM. This method determines the propagation direction of the wavefield at the shot-receiver point in the frequency-wavenumber domain, then calculates the reflection angle information and directly outputs high-quality ADCIGs. It also uses leakage-resistant Fourier transform to solve the angular resolution problem. However, this method requires storing and performing a four-dimensional Fourier transform on the wavefield at all times during the reverse-time migration process, resulting in significant computational and storage costs. Methods using analytical-time wavefield extrapolation and local wavefield sparse decomposition to generate ADCIGs can determine the wavefield propagation direction without performing a time-dimensional Fourier transform, solving the storage problem of the entire extrapolated wavefield in the shot-receiver direction. However, this method requires extrapolating the imaginary part of the analytical-time wavefield, and the additional computational cost is not negligible, especially in anisotropic media, where implementation will require a huge computational investment.

[0005] Therefore, it is necessary to develop an angle-domain adaptive local optimization reverse time offset method, electronic device, and medium.

[0006] The information disclosed in the background section of this invention is intended only to enhance the understanding of the general background of this invention, and should not be construed as an admission or in any way implying that such information constitutes prior art known to those skilled in the art. Summary of the Invention

[0007] This invention proposes an angle-domain adaptive local optimization reverse-time migration method, electronic device, and medium, which solves the accuracy problem of angle-domain reverse-time migration common imaging point gathers, making the angle-domain reverse-time migration common imaging point gather results more accurate and reliable, which is beneficial to subsequent interpretation and processing.

[0008] In a first aspect, embodiments of this disclosure provide an angle-domain adaptive local optimization inverse time migration method, including:

[0009] Extrapolate the reverse-time wave field at the shot point and receiver point, and construct the energy target functional under gradient constraints on each time slice;

[0010] Iteratively solve for the wave field vector directions at the current shot point and receiver point;

[0011] An adaptive threshold is constructed by global amplitude scanning, and a table of valid points at the current moment is generated by filtering. The reflection angle or azimuth of the valid points in the table is then calculated.

[0012] Interpolation calculations of the reflection angle or azimuth angle are performed within the set range of each valid point in the table;

[0013] By applying the cross-correlation imaging conditions in the azimuth and angle domain, the gather at the underground imaging point is determined.

[0014] Preferably, the energy objective functional under the gradient constraint is:

[0015]

[0016] Where x, y, z represent the three-dimensional Cartesian coordinate system, and P is the stress wave field. For spatial gradient operators, E represents energy, u, v, w are the gradient components in the x, y, z directions respectively, α is the weighting coefficient, and P... x ,P y ,P z P represents the spatial gradient of the wave field. t This represents the time gradient of the wave field.

[0017] Preferably, the energy target functional is minimized, and the wave field vector directions at the current shot point and receiver point are iteratively solved.

[0018] Preferably, the solution is obtained iteratively using formula (2):

[0019]

[0020] in, P x ,P y ,P z P represents the spatial gradient of the wave field. t The wave field's temporal gradient is represented by ; u, v, w represent the neighborhood mean of u, v, w; i, j, k represent the spatial exponents on x, y, z, respectively.

[0021] Preferably, the reflection angle is calculated using formula (3):

[0022]

[0023] Where θ is the reflection angle, S S S R These represent the wavefield direction vectors at the shot point and receiver point, respectively.

[0024] Preferably, the azimuth angle is calculated using formula (4):

[0025]

[0026] Where ψ is the azimuth angle, S S S R Let n represent the wavefield direction vectors at the shot point and receiver point, respectively. x = (1,0,0) and n z = (0,0,1) are the unit vectors in the x and z directions, respectively.

[0027] Preferably, interpolation is performed using formula (5):

[0028]

[0029] Where, x valid θ represents the current valid point position. valid This represents the reflection angle or azimuth of the valid point, where x is the position of the invalid point within the current valid point setting range, and θ is... invalid The reflection angle or azimuth of the invalid point is represented by N, the number of times the invalid point is interpolated, and σ is the degree of interpolation smoothness.

[0030] Preferably, the gather at the underground imaging point is:

[0031] I(x,θ,ψ)=∑ t u s(t,x,θ,ψ)u R (t,x,θ,ψ) (6)

[0032] Where u is the gather at the underground imaging point x. s (t,x,θ,ψ) represents the end-shot wavefield at imaging point x in the angular domain, u R (t,x,θ,ψ) represents the end wave field of the detector point x in the angular domain.

[0033] As one specific implementation of this disclosure,

[0034] Secondly, embodiments of this disclosure also provide an electronic device, the electronic device comprising:

[0035] Memory, which stores executable instructions;

[0036] A processor that executes the executable instructions in the memory to implement the angle-domain adaptive local optimization reverse-time offset method.

[0037] Thirdly, embodiments of this disclosure also provide a computer-readable storage medium storing a computer program that, when executed by a processor, implements the angle domain adaptive local optimization reverse time offset method.

[0038] Its beneficial effects are as follows: Addressing the accuracy problem of angle-domain reverse-time migration common imaging point gathers, based on the conventional Poynting vector algorithm, an effective point table is generated through adaptive threshold global filtering at each time slice of the reverse-time wavefield extrapolation for the shot and receiver points. The reflection angle or azimuth is calculated only within this table, improving the stability of angle estimation and reducing the number of actual calculated subsurface grid points. Furthermore, local optimization is performed within a certain range of effective points, ensuring the continuity of the final RTM-generated ADCIGs and reducing imaging noise. Compared to the traditional Poynting vector algorithm, this invention solves the accuracy problem of angle-domain reverse-time migration common imaging point gathers, making the results more accurate and reliable, which is beneficial for subsequent interpretation and processing.

[0039] The methods and apparatus of the present invention have other features and advantages that will be apparent from or will be set forth in detail in the accompanying drawings and following detailed description, which together serve to explain the particular principles of the invention. Attached Figure Description

[0040] The above and other objects, features and advantages of the present invention will become more apparent from the more detailed description of exemplary embodiments of the invention in conjunction with the accompanying drawings, wherein the same reference numerals generally represent the same parts.

[0041] Figure 1a , Figure 1b , Figure 1c The diagrams show the cross-correlation imaging results of the angular domain imaging reflection angle according to an embodiment of the present invention, and the beam before and after local optimization of the effective point.

[0042] Figure 2 A flowchart illustrating the steps of an angle-domain adaptive local optimization reverse-time offset method according to an embodiment of the present invention is shown.

[0043] Figure 3 A schematic diagram of a work area offset velocity model according to an embodiment of the present invention is shown.

[0044] Figure 4 A schematic diagram of the angular domain RTM common imaging point gather at 780–810CDP along line 1100, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0045] Figure 5 A schematic diagram of an angular domain RTM common imaging point gather at 780–810 CDP along 1100 lines, generated using this method according to an embodiment of the present invention, is shown.

[0046] Figure 6 A schematic diagram of the angular domain RTM common imaging point gather at 840–870 CDP on the 1100 line, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0047] Figure 7 A schematic diagram of an angular domain RTM common imaging point gather at 840–870 CDP using the method according to an embodiment of the present invention is shown.

[0048] Figure 8 A schematic diagram of the angular domain RTM common imaging point gather at 900–930 CDP along the 1100 line, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0049] Figure 9 A schematic diagram of an angular domain RTM common imaging point gather at 900–930 CDP using the method according to an embodiment of the present invention is shown.

[0050] Figure 10a , Figure 10bSchematic diagrams of the final omnidirectional angle domain RTM common imaging point gather and overlay profile are shown respectively according to an embodiment of the present invention. Detailed Implementation

[0051] Preferred embodiments of the invention will now be described in more detail. While preferred embodiments of the invention are described below, it should be understood that the invention can be implemented in various forms and should not be limited to the embodiments set forth herein.

[0052] This invention provides an angle-domain adaptive local optimization reverse-time migration method, comprising:

[0053] Extrapolate the reverse-time wave field at the shot point and receiver point, and construct the energy target functional under gradient constraints on each time slice;

[0054] Iteratively solve for the wave field vector directions at the current shot point and receiver point;

[0055] An adaptive threshold is constructed by global amplitude scanning, and a table of valid points at the current moment is generated by filtering. The reflection angle or azimuth of the valid points in the table is then calculated.

[0056] Interpolation calculations of the reflection angle or azimuth angle are performed within the set range of each valid point in the table;

[0057] By applying the cross-correlation imaging conditions in the azimuth and angle domain, the gather at the underground imaging point is determined.

[0058] In one example, the energy objective functional under gradient constraints is:

[0059]

[0060] Where x, y, z represent the three-dimensional Cartesian coordinate system, and P is the stress wave field. For spatial gradient operators, E represents energy, u, v, and w are the gradient components in the x, y, and z directions, respectively, and α is a tradeoff coefficient that controls the balance between the data terms and the smoothing terms. It represents the reliability of the wavefield data and smoothing constraints, and plays a role in adjusting the wavefield energy conservation constraints and global smoothing. It is generally taken as α = 1. x ,P y ,P z P represents the spatial gradient of the wave field. t This represents the time gradient of the wave field.

[0061] In one example, the energy objective functional is minimized, and the wave field vector directions at the current shot point and receiver point are iteratively solved.

[0062] In one example, the solution is obtained iteratively using formula (2):

[0063]

[0064] in, P x ,P y ,P z P represents the spatial gradient of the wave field. t Represents the time gradient of the wave field; represents the neighborhood mean of u, v, w; i, j, k represent the spatial indices on x, y, z, respectively.

[0065] In one example, the reflection angle is calculated using formula (3):

[0066]

[0067] Where θ is the reflection angle, S S S R These represent the wavefield direction vectors at the shot point and receiver point, respectively.

[0068] In one example, the azimuth angle is calculated using formula (4):

[0069]

[0070] Where ψ is the azimuth angle, S S S R Let n represent the wavefield direction vectors at the shot point and receiver point, respectively. x = (1,0,0) and n z = (0,0,1) are the unit vectors in the x and z directions, respectively.

[0071] In one example, interpolation is performed using formula (5):

[0072]

[0073] Where, x valid θ represents the current valid point position. valid This represents the reflection angle or azimuth of the valid point, where x is the position of the invalid point within the current valid point setting range, and θ is... invalid The reflection angle or azimuth of the invalid point is represented by N, which is the number of times the invalid point is interpolated and used to calculate the arithmetic mean. σ is the smoothness of the interpolation, which is generally taken as 1.

[0074] In one example, the gather at the underground imaging point is:

[0075] I(x,θ,ψ)=∑ t u s (t,x,θ,ψ)u R(t,x,θ,ψ) (6)

[0076] Where u is the gather at the underground imaging point x. s (t,x,θ,ψ) represents the end-shot wavefield at imaging point x in the angular domain, u R (t,x,θ,ψ) represents the end wave field of the detector point x in the angular domain.

[0077] Specifically, the offset velocity field and initial offset parameters are read in, and the inverse time wave field of the shot point and receiver point is extrapolated to construct the energy target functional under gradient constraints on each time slice.

[0078] The energy flux density vector, originally used to measure the electromagnetic field energy passing through a unit area perpendicular to the vector direction per unit time, was first introduced from the electromagnetic field into the field of seismic exploration in 2006 and applied to migration imaging to obtain the wave field propagation direction. Its calculation formula is as follows:

[0079]

[0080] Where x, y, z represent the three-dimensional Cartesian coordinate system, υ represents the displacement velocity vector, t is time, and P is the stress wave field. For the spatial gradient operator. As can be seen from formula (7), when applying the scalar acoustic wave equation, the propagation direction of the wave field is obtained by the energy flux density vector method. It is necessary to calculate the spatial derivative and time derivative of the wave field at the current moment. The calculation efficiency is high, it is easy to implement, and it has a high angular resolution. However, it is easily affected by noise, especially at the receiver end. When the surface received data is extrapolated in reverse time, the continuity of the calculated wave field vector direction is poor, and the accuracy is difficult to guarantee. In addition, if the wavefront crosses or reaches the interface, the reflection angle cannot be correctly estimated. At each time slice of the inverse time wave field extrapolation at the shot point and receiver point, the energy target functional under the gradient constraint is constructed as formula (1).

[0081] Minimize the energy target functional and iteratively solve the wave field vector directions at the current shot point and receiver point using formula (2);

[0082] In the process of extracting angular domain imaging gathers, it is necessary to convert the wavefield propagation direction calculated by the Poynting vector or optical flow vector into the reflection angle and azimuth angle. Here, S is used. S S R The wave field direction vectors at the shot point and receiver point are respectively represented. An adaptive threshold is constructed by global amplitude scanning, and a table of valid points at the current time is generated. The reflection angle of the valid points in the table is calculated by formula (3), and the azimuth angle is calculated by formula (4).

[0083] Not all imaging points within the current imaging angle range contribute to the final RTM imaging result; some may even introduce imaging noise. When calculating the reflection angle and azimuth angle, grid points with small or even zero denominators in formulas (3) and (4) will make division at that location extremely unstable, leading to incorrect angle estimation.

[0084] At each moment of the reverse-time wave field extrapolation at the shot point and receiver point, an adaptive threshold is constructed through global amplitude scanning to filter and generate a table of valid points for the current moment. Only the reflection angle or azimuth of the valid points in the table needs to be calculated, which improves the stability of angle estimation and also reduces the number of actual underground grid points to be calculated. The specific definition of the reflection angle is as follows:

[0085]

[0086] Where T represents the adaptive threshold at the current time t; ε is a human-controlled parameter that depends on the signal-to-noise ratio of the data; AMP is the denominator in formula (3), which is the product of the magnitudes of the wave field direction vector at the shot point and the wave field direction vector at the receiver point.

[0087] Within the set range of each valid point in the table, the reflection angle or azimuth angle is interpolated using formula (5);

[0088] Figure 1a , Figure 1b , Figure 1c The diagrams show the cross-correlation imaging results of the angular domain imaging reflection angle according to an embodiment of the present invention, and the beam before and after local optimization of the effective point.

[0089] To illustrate the effectiveness of the above optimization approach, a theoretical model example is given. A reverse-time migration process is selected at a certain moment, and the imaging results of the cross-correlation between the forward and reverse wave fields at that moment are as follows: Figure 1a As shown, observation reveals that the main imaging energy is concentrated on both sides (red boxes), which allows us to determine the threshold range, pinpoint the effective point location at the current moment, and identify its corresponding reflection angle. Figure 1b As shown, the number of calculation grid points has been reduced to 3.83% of the original. However, the angle result is biased towards high frequencies, which may affect the quality of the final generated ADCIGs. Therefore, it can be locally interpolated within a certain range using formula (5). The result of local optimization is shown in the figure. Figure 1c As shown. Unlike the distribution range of the reflection angle, the azimuth angle distribution is relatively smooth and does not require high resolution. Therefore, it can be defined in the same way, which ensures the continuity of the final azimuth angle domain common imaging point gather.

[0090] After forward modeling and reverse extrapolation of the reverse time migration, the shot end wavefield u of the imaging point x in the angular domain can be obtained respectively. s (t,x,θ,ψ) and the wavefield u at the receiver point R (t,x,θ,ψ), using the cross-correlation imaging conditions in the azimuth angle domain, the gather at the underground imaging point is determined by formula (6).

[0091] The present invention also provides an electronic device, comprising: a memory storing executable instructions; and a processor executing the executable instructions in the memory to implement the above-described angle domain adaptive local optimization reverse time offset method.

[0092] The present invention also provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the above-described angle-domain adaptive local optimization reverse-time offset method.

[0093] To facilitate understanding of the solutions and effects of the embodiments of the present invention, three specific application examples are given below. Those skilled in the art should understand that these examples are merely for the purpose of understanding the present invention, and any specific details therein are not intended to limit the present invention in any way.

[0094] Example 1

[0095] Figure 2 A flowchart illustrating the steps of an angle-domain adaptive local optimization reverse-time offset method according to an embodiment of the present invention is shown.

[0096] like Figure 2 As shown, the angle-domain adaptive local optimization reverse-time migration method includes: Step 101, extrapolating the reverse-time wavefield of the shot point and receiver point, and constructing an energy target functional under gradient constraints on each time slice; Step 102, iteratively solving the wavefield vector direction at the shot point end and receiver point end at the current time; Step 103, constructing an adaptive threshold through global amplitude scanning, filtering and generating a table of valid points at the current time, and calculating the reflection angle or azimuth of the valid points in the table; Step 104, performing interpolation calculation of the reflection angle or azimuth within a set range for each valid point in the table; Step 105, applying the azimuth-angle domain cross-correlation imaging conditions to determine the gather at the underground imaging point.

[0097] Figure 3 A schematic diagram of a work area offset velocity model according to an embodiment of the present invention is shown.

[0098] Figure 4 A schematic diagram of the angular domain RTM common imaging point gather at 780–810CDP along line 1100, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0099] Figure 5 A schematic diagram of an angular domain RTM common imaging point gather at 780–810 CDP along 1100 lines, generated using this method according to an embodiment of the present invention, is shown.

[0100] Figure 6 A schematic diagram of the angular domain RTM common imaging point gather at 840–870 CDP on the 1100 line, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0101] Figure 7 A schematic diagram of an angular domain RTM common imaging point gather at 840–870 CDP using the method according to an embodiment of the present invention is shown.

[0102] Figure 8 A schematic diagram of the angular domain RTM common imaging point gather at 900–930 CDP along the 1100 line, generated using the Poynting vector algorithm according to an embodiment of the present invention, is shown.

[0103] Figure 9 A schematic diagram of an angular domain RTM common imaging point gather at 900–930 CDP using the method according to an embodiment of the present invention is shown.

[0104] Figure 10a , Figure 10b Schematic diagrams of the final omnidirectional angle domain RTM common imaging point gather and overlay profile are shown respectively according to an embodiment of the present invention.

[0105] To test the robustness of the angle domain adaptive local optimization algorithm, data from a real land work area were processed. Figure 3 For the work area offset velocity model; Figure 4 This is a common imaging point gather in the angular domain at 780–810 CDP along the 1100 line, generated using the Poynting vector algorithm. Figure 5 This is the angular domain RTM common imaging point gather at 780-810CDP of 1100 lines generated using the method of this invention; Figure 6 This is a common imaging point gather for the angular domain RTM at 840–870°CDP, generated using the Poynting vector algorithm. Figure 7 This is a common imaging point gather of the 1100-line angular domain RTM at 840-870CDP, generated using the method of this invention; Figure 8 This is a common imaging point gather for the angular domain RTM at 900–930°CDP, generated using the Poynting vector algorithm. Figure 9The image shows the angular domain RTM common imaging point gather at 900–930 CDP along 1100 lines, generated using the method of this invention. Figure 10 shows the final omnidirectional angular domain RTM common imaging point gather and the overlay profile. A comparison reveals that the common imaging point gather generated by the adaptive local optimization algorithm significantly reduces large-angle noise, improves the signal-to-noise ratio in the middle and shallow layers, and maintains angular extensibility.

[0106] Example 2

[0107] This disclosure provides an electronic device comprising: a memory storing executable instructions; and a processor executing the executable instructions in the memory to implement the aforementioned angle-domain adaptive local optimization reverse-time offset method.

[0108] An electronic device according to an embodiment of the present disclosure includes a memory and a processor.

[0109] This memory is used to store non-transitory computer-readable instructions. Specifically, the memory may include one or more computer program products, which may include various forms of computer-readable storage media, such as volatile memory and / or non-volatile memory. The volatile memory may, for example, include random access memory (RAM) and / or cache memory. The non-volatile memory may, for example, include read-only memory (ROM), hard disk, flash memory, etc.

[0110] The processor may be a central processing unit (CPU) or other form of processing unit with data processing capabilities and / or instruction execution capabilities, and may control other components in the electronic device to perform desired functions. In one embodiment of this disclosure, the processor is used to execute computer-readable instructions stored in the memory.

[0111] Those skilled in the art will understand that, in order to solve the technical problem of how to achieve a good user experience, this embodiment may also include well-known structures such as communication buses and interfaces, and these well-known structures should also be included within the protection scope of this disclosure.

[0112] For a detailed description of this embodiment, please refer to the corresponding descriptions in the foregoing embodiments, which will not be repeated here.

[0113] Example 3

[0114] This disclosure provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the angle domain adaptive local optimization reverse time offset method.

[0115] A computer-readable storage medium according to embodiments of the present disclosure stores non-transitory computer-readable instructions. When these non-transitory computer-readable instructions are executed by a processor, all or part of the steps of the methods described in the foregoing embodiments of the present disclosure are performed.

[0116] The aforementioned computer-readable storage media include, but are not limited to: optical storage media (e.g., CD-ROM and DVD), magneto-optical storage media (e.g., MO), magnetic storage media (e.g., magnetic tape or portable hard drive), media with built-in rewritable non-volatile memory (e.g., memory card), and media with built-in ROM (e.g., ROM cartridge).

[0117] Those skilled in the art should understand that the above description of the embodiments of the present invention is only intended to illustrate the beneficial effects of the embodiments of the present invention, and is not intended to limit the embodiments of the present invention to any of the examples given.

[0118] The various embodiments of the present invention have been described above. These descriptions are exemplary and not exhaustive, nor are they limited to the disclosed embodiments. Many modifications and variations will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments.

Claims

1. An angle domain adaptive locally optimized reverse time migration method, characterized in that, The method comprises the following steps: Outward extrapolation of the shot point and the receiver point reverse-time wave field, construction of an energy objective function under gradient constraint on each time slice; Iterative solution of the shot point end and the receiver point end wave field vector direction at the current time; Construction of an adaptive threshold value through global amplitude scanning, screening to generate an effective point table at the current time, and calculation of the reflection opening angle or the azimuth angle of the effective points in the table; Interpolation calculation of the reflection opening angle or the azimuth angle within the set range of each effective point in the table; Application of an azimuth angle domain cross-correlation imaging condition to determine the gather at the underground imaging point.

2. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, The energy objective function under the gradient constraint is: wherein x, y, z represent a three-dimensional Cartesian coordinate system, P is a stress wave field, for a spatial gradient operator, E is energy, u, v, w are gradient components in vector x, y, z directions respectively, a is a weighting factor, P x ,P y ,P z represents a spatial gradient of the wave field, P t represents a time gradient of the wave field.

3. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, Minimization of the energy objective function, and iterative solution of the shot point end and the receiver point end wave field vector direction at the current time.

4. The angle-domain adaptive locally optimized reverse time migration method of claim 3, wherein, Iterative solution through formula (2): where P x ,P y ,P z denotes the spatial gradient of the wavefield, P t denotes the time gradient of the wavefield; denotes the neighborhood average of u, v, w; i, j, k denote the spatial indices in x, y, z, respectively.

5. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, Calculation of the reflection opening angle through formula (3): where θ is the reflection opening angle, S S , S R denote the wavefield direction vectors at the source and receiver ends, respectively.

6. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, Calculation of the azimuth angle through formula (4): where ψ is the azimuth angle, S S , S R denote the wavefield direction vectors at the source and receiver ends, respectively, n x = (1, 0, 0) and n z = (0, 0, 1) are the unit vectors in the x and z directions, respectively.

7. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, Interpolation calculation through formula (5): where x valid is the current valid point position, θ valid represents the reflection angle or azimuth of the valid point, x is the invalid point position within the current valid point setting range, θ invalid represents the reflection angle or azimuth of the invalid point, N is the calculation number of times that the invalid point is interpolated, and σ is the interpolation smoothness.

8. The angle-domain adaptive locally optimized reverse time migration method of claim 1, wherein, The gather at the underground imaging point is: I(x, θ, ψ) = ∑ t u s (t, x, θ, ψ)u R (t, x, θ, ψ) (6) where u(x, t) is the data set at the imaging point x, t s (t, x, θ, ψ) is the shot-side wavefield of the imaging point x in the angle domain, u R (t, x, θ, ψ) is the receiver-side wavefield of the imaging point x in the angle domain.

9. An electronic device, comprising: The electronic device comprises: A memory storing executable instructions; A processor running the executable instructions in the memory to implement the angle domain adaptive local optimization reverse-time migration method in any one of claims 1-8.

10. A computer-readable storage medium, characterized in that, The computer readable storage medium stores a computer program, which, when executed by a processor, implements the angle domain adaptive local optimization reverse-time migration method in any one of claims 1-8.

Images