A complex geological body imaging method based on reverse time migration of reflected waves

Abstract

The present invention belongs to the technical field of oil and gas seismic exploration, and specifically relates to a method for imaging complex geological bodies based on reverse time migration of reflection waves. The method finally images complex geological bodies through reflection wave field feature analysis, reflection wave observation system design and optimization, reflection wave imaging, and reflection wave seismic acquisition test data imaging processing. The joint imaging profile of reflection waves and reflection waves effectively improves the imaging accuracy of structures. Reflection wave imaging alone is used to identify and interpret high-steep structures, and key parts such as turning points and breakpoints of complex structures can be more accurately determined. The reverse time migration imaging method based on two-way wave extension imaging is selected to effectively extract the imaging information of reflection waves. Angle decomposition is applied at local imaging points. By judging the relationship between the local stratum dip and the incident angle, the imaging results of the reflection waves are extracted, which significantly improves the imaging accuracy.

Description

Technical Field

[0001] The present invention belongs to the technical field of geophysical exploration, and in particular relates to a complex geological body imaging method based on reverse time migration of reflected waves. Background Art

[0002] In oil and gas seismic exploration, conventional reflection wave imaging methods have significant limitations when dealing with complex structures with steep dips, such as low imaging accuracy and poor illumination. While single-way wave equation migration and Kirchhoff integral migration can handle steep structures to a certain extent, they still struggle to achieve high-precision imaging results under complex geological conditions. Reflected waves, with their advantages such as deep propagation depth, excellent illumination of steep structures, and ability to provide low-wavenumber information for inversion, are key to addressing these challenges.

[0003] Research on reverberation waves primarily focuses on their propagation characteristics and basic imaging methods. When the velocity of the subsurface increases linearly with depth, the reverberation waves exhibit circular arc-shaped motion. However, existing technologies still have significant shortcomings when dealing with complex structures with steep dips. These limitations include the paraxial approximation limitations of the one-way wave equation migration method and the limited handling of caustics and multivalued traveltimes by the Kirchhoff integral migration method. These limitations make it difficult for existing technologies to meet the demands for high-precision imaging in complex geological conditions.

[0004] Currently, there is little research on seismic reflection waves at home and abroad. How to reasonably utilize the propagation characteristics of reflection waves in the imaging process to finely characterize steep-angle complex structures, rather than simply suppressing them through complex data preprocessing, is an urgent problem that needs to be solved. It is of great significance to the utilization of seismic reflection wave fields in complex structure imaging. Summary of the Invention

[0005] The present invention effectively solves the imaging problem of existing methods in complex structures with high steep dips by introducing a two-way wave equation wave field bidirectional continuation method and combining it with reverse time migration technology.

[0006] To achieve the above object, the present invention adopts the following technical solutions for a complex geological body imaging method based on reverse time migration of return waves:

[0007] S1. Analysis of the characteristics of the reflected wave field

[0008] Based on seismological principles, the formation mechanism of the back-reflection wave and the geological conditions required are analyzed. The complex geological bodies in the target work area are simulated and analyzed to determine their impact on the back-reflection wave propagation path and wave field distribution characteristics. Geological-velocity models of different types of high- and steep-angle structures are established, and forward modeling of the back-reflection wave is performed to obtain and analyze its wave field characteristics, providing a basis for observation system design and migration imaging.

[0009] S2. Design and optimization of the back-reflected wave observation system

[0010] In combination with the analysis results of step S1 and the geological objectives of the work area, the observation system parameters are optimized through forward simulation tests; using seismic record separation technology, the return wave containing reflection information (reflection return wave) and the return wave not containing reflection information (non-reflection return wave) are extracted from the original wave field respectively; the non-reflection return wave field is back-propagated and its wave field square is calculated to obtain the source end illumination representing the energy distribution of the earthquake source; the reflection return wave field is back-propagated and its wave field square is calculated to obtain the detector end illumination representing the energy distribution received by the detector; the source end illumination and the detector end illumination are used to evaluate and optimize the illumination effect of the observation system on the target geological body;

[0011] S3, Reverse Time Migration Imaging

[0012] The reverse time migration method based on the two-way wave equation is used for wave field extension. Specifically, the method includes: in the time domain, forward extension of the source wave field and reverse extension of the received reflected / reflected wave field (detection point wave field); in the wave field extension process, the local ray parameters of the source wave field and the detection point wave field at each point and time in space are calculated and extracted; based on the local ray parameters, the local formation dip at the imaging point is calculated; at the imaging point, different imaging conditions are applied according to the size of the local formation dip. Obtaining the back-reflected wave imaging value: When the local formation dip is greater than 90°, the first imaging condition is applied to calculate the back-reflected wave imaging value; when the local formation dip is less than 90°, the local imaging reflection angle is calculated. If the sum of the local imaging reflection angle and the local formation dip is greater than 90°, the second imaging condition is applied to calculate the back-reflected wave imaging value; the back-reflected wave imaging values at each imaging point that meets the conditions are accumulated to obtain the single-shot back-reflected wave imaging result; the back-reflected wave imaging results of all single shots are superimposed to obtain the final multi-shot back-reflected wave imaging volume;

[0013] S4. Imaging Processing and Effect Verification

[0014] The back-reflected wave imaging result obtained in step S3 is subjected to noise suppression post-processing; the imaging effect of the method is verified using numerical simulation data and actual seismic acquisition data of the target work area; and the quality of the imaging result is quantitatively evaluated by signal-to-noise ratio analysis, especially the effect of improving the imaging accuracy of key parts of turning points and breakpoints of steep structures.

[0015] The present invention has the following beneficial effects:

[0016] (1) The combined imaging profile of reflection waves and back-reflection waves effectively improves the imaging accuracy of structures. In addition, the use of back-reflection wave imaging alone to identify and interpret high-steep structures can more accurately determine key locations such as turning points and breakpoints in complex structures.

[0017] (2) The present invention uses a reverse time migration imaging method based on two-way wavelet continuation imaging, which can effectively extract imaging information of the return wave. By applying angle decomposition at the local imaging point and determining the relationship between the local formation dip and the incident angle, the imaging results of the return wave are extracted, significantly improving imaging accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 、 one Flowchart of a complex geological body imaging method based on reverse time migration of back-wave;

[0019] Figure 2 、 one Relationship between receiving distance and maximum diving depth of reflected wave at different surface velocities;

[0020] Figure 3 、 one Relationship between maximum diving depth and recording time at different surface velocities;

[0021] Figure 4 、 one Different groups Relationship between down-offset and maximum diving depth of reflected wave;

[0022] Figure 5 、 one Different groups Graph of the relationship between down-offset and recording time;

[0023] Figure 6 , velocity model diagram of the 82° strike-slip fault in Shunbei area according to an embodiment of the present invention;

[0024] Figure 7 , source end illumination (hanging wall excitation) diagram of the 82° strike-slip fault in Shunbei area according to an embodiment of the present invention;

[0025] Figure 8 , an illumination image (hanging wall excitation) of the 82° strike-slip fault in Shunbei area according to an embodiment of the present invention;

[0026] Figure 9 , source end illumination (footwall excitation) diagram of the 82° strike-slip fault in Shunbei area according to an embodiment of the present invention;

[0027] Figure 10 , schematic diagram of the propagation of the reflected wave;

[0028] Figure 11 , schematic diagram of underground high and steep dip propagation;

[0029] Figure 12 , a single shot imaging cross-section of the 82° strike-slip fault in Shunbei area according to an embodiment of the present invention;

[0030] Figure 13, FootHill model reflection wave imaging profile;

[0031] Figure 14 , FootHill model back-reflection wave imaging profile;

[0032] Figure 15 , FootHill model full wavefield imaging section;

[0033] Figure 16 , conventional imaging section of Shunbei working area in the embodiment of the present invention;

[0034] Figure 17 , the back-reflected wave imaging section of the Shunbei work area in the embodiment of the present invention;

[0035] Figure 18 , the back-reflected wave enhanced imaging section of the Shunbei work area in the embodiment of the present invention;

[0036] Figure 19 , conventional imaging section of Shunbei working area in the embodiment of the present invention;

[0037] Figure 20 , the back-reflected wave enhanced imaging section of the Shunbei work area in the embodiment of the present invention;

[0038] Figure 21 、The present invention Figure 20 Enlarged view of point A in the middle;

[0039] Figure 22 , Signal-to-noise ratio analysis diagram at the strike-slip fault in Shunbei work area according to an embodiment of the present invention. DETAILED DESCRIPTION

[0040] To make the purpose, technical solutions, and advantages of this application more clear, the technical solutions of this application will be clearly and completely described below in conjunction with the specific embodiments of this application and the corresponding drawings. Obviously, the embodiments described are only part of the embodiments of this application, not all of them. Based on the embodiments in this specification, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of this application.

[0041] The Shunbei region is located on the northern edge of China's Tarim Basin, straddling the Aksu Prefecture of the Xinjiang Uyghur Autonomous Region and the Bayingolin Mongol Autonomous Prefecture. This area has unique geological structures within the basin and is a key area for deep oil and gas exploration. Located in the northern depression of the Tarim Basin, the Shunbei region boasts complex and diverse geological structures, including high-steep strike-slip faults. Seismic data with a 10km offset exhibit a high abundance of back-reflection waves. The proposed scheme is illustrated with the accompanying figures:

[0042] S1. Analysis of the characteristics of the reflected wave field

[0043] The formation mechanism of the reflex wave is analyzed based on the basic principles of seismology. The geological conditions for its formation are studied, and the influence of actual high-steep structures, igneous rocks, and other geological bodies on the propagation path and wavefield distribution of the reflex wave is analyzed. To achieve high-quality and high-precision imaging of high-steep structures, reflex wave simulations and wavefield characteristics of different types of high-steep structures are carried out, laying a solid foundation for subsequent observation system design and migration imaging.

[0044] The time-distance characteristics of the return wave under different underground medium parameters and ground receiving parameters are analyzed based on the time-distance curve formula of the return wave. The main factors affecting the time-distance characteristics of the return wave are: surface velocity , underground medium velocity change rate These factors together control the maximum propagation depth of the reflected wave received on the ground, as well as the recording time.

[0045] When the velocity change rate of underground medium When a certain time, given , the maximum propagation depth of the reflected waves with different surface velocities at different receiving distances, such as Figure 2 As shown in the figure, it can be concluded that when the velocity change rate is constant, the propagation depth of the reflected wave is the same for different surface velocities at a certain receiving distance. The propagation depth of the reflected wave is determined only by the surface receiving distance. The larger the receiving range, the deeper the reflected wave propagation depth. Therefore, to obtain deeper reflected wave information, a wide receiving range is required.

[0046] When the velocity change rate of underground medium When a certain time, given , different surface velocity reflected waves, receiving time at different receiving distances. Figure 3 It can be seen that when the velocity change rate is constant, at a given receiving distance, the greater the surface velocity, the shorter the return wave reception time. Since the velocity change rate is constant, the return wave propagation path at different receiving points is the same, so the greater the surface velocity, the shorter the propagation time. Consistent with the reflected wave, the greater the offset, the longer the reception time. Simultaneously, the greater the surface velocity, the shorter the reception time.

[0047] When the surface velocity is constant, the surface velocity is given Velocity, analysis of velocity change rates of different underground media Corresponding to the reflected wave, the relationship between different receiving distances and the propagation depth of the reflected wave, such as Figure 4 As shown. It can be seen that when the surface velocity is constant, at a certain receiving distance, the faster the velocity changes ( The larger the receiving distance, the greater the maximum propagation depth of the reflected wave. When the receiving distance increases, the maximum propagation depth of the reflected wave also increases. When it increases to several kilometers, the maximum propagation depth of the reflected wave is nearly linearly related to the receiving distance.

[0048] When the surface velocity is constant, the surface velocity is given Velocity, analysis of velocity change rates of different underground media Corresponding to the reflected wave, the relationship between different receiving distances and the reflected wave receiving time, such as Figure 5 As shown. Figure 5 This is the back-reflection wave curve in the seismic record. When the surface velocity is constant, at a certain receiving distance, the faster the velocity changes ( The larger the value, the shorter the return wave reception time. When the receiving distance increases, the return wave reception time also increases, and the apparent speed increases.

[0049] S2. Design and optimization of the back-reflected wave observation system

[0050] The observation system was optimized through forward modeling and observation system parameter testing. For the Shunbei geological model, the reflected and non-reflected waves were extracted by simple excision. The non-reflected wave was propagated back and the wave field was squared to obtain the source end illumination formula (1); the reflected wave was propagated back and the wave field was squared to obtain the detector end illumination formula (2). Figure 6 Velocity model of the 82° strike-slip fault in Shunbei area, Figure 7-Figure 9 They are source-end lighting and detector-end lighting at the best excitation points on the upper and lower walls of the strike-slip fault.

[0051] ;

[0052] in, is the illumination at the source end; is the illumination intensity at the detection end; is the non-reflected back-reflected wave field, is the reflected wave field, x is the offset, z is the wave field depth; is the corresponding frequency.

[0053] S3, Reverse Time Migration Imaging

[0054] A two-way wave equation bidirectional continuation method is used for back-reflected wave imaging. Conventional one-way wave equation wave field continuation methods primarily consider the downward propagation of the wave field, and are limited in propagation angle. Large-angle imaging errors are large, making it difficult to image steeply dipped structures. A two-way wave equation bidirectional continuation method for steeply dipped structures, known as reverse time migration, is developed. This method fully utilizes the amplitude-preserving properties of the two-way wave equation during wave field propagation and utilizes a combination of bidirectional wave field extrapolation and the two-way wave equation wave field continuation method to address seismic migration imaging of steeply dipped and complex structures.

[0055] Using the Poynting vector to generate angle gathers within the framework of the two-way wave equation, reverse time migration (RTM) of acoustic media backscatter waves can achieve high-precision imaging of structures with dips near or greater than 90 degrees compared to conventional reflection-wave RTM. This method overcomes the poor imaging of steeply dipped structures by conventional reflection waves and improves the imaging accuracy of complex, steep structures such as salt domes, igneous rocks, and overthrust nappes. The specific steps are as follows.

[0056] ① Extrapolate the source wave field and the receiver wave field to extract the local ray parameters of the wave field;

[0057] The wave equation is solved using the finite difference time domain method to achieve forward extrapolation of the source wave field in the time direction and reverse extrapolation of the detector wave field in the time direction. The wave field (source end or detector end) is , the velocity model is , the wave equation can be expressed as:

[0058] ;

[0059] in, They are the second-order partial derivatives of the wave field in directions j, k, and l, as well as the second-order partial derivatives with respect to time t.

[0060] The wave equation (3) is solved using the time-domain finite difference to obtain the wave field values (source end and detector end) at each time and each spatial position. Figure 10 As shown in Figure 2, for the reflected wave, the main imaging structure is a high-steep-angle structure. Figure 10 In the two-dimensional case, assume that the propagation path of a shot S1 on the left side of the steep-angle structure is as shown in the figure, reflected at point R1 and transmitted back to point G1. Similarly, a shot S2 on the right side of the steep-angle structure propagates to the reflection point R2 of the steep-angle structure, and then reflects and transmits back to point G2. The time derivative and spatial derivative of the wave field are obtained to obtain the local propagation direction of the wave field (ray parameter). Figure 11 for Figure 10 The local amplification of the source and detector wave field propagation directions (ray parameters) can be expressed as follows:

[0061] ;

[0062] in are the local ray parameters of the source and detector wave fields, respectively; x is the offset; z is the wave field depth; t is the time; are the first-order partial derivatives of the wave field x, z, t and the first-order partial derivative with respect to time t, are the first-order partial derivatives of the wave field x, z, t and the first-order partial derivatives with respect to x and z respectively.

[0063] like Figure 11As shown in Figure 1, it is an enlarged schematic diagram of the incident reflection at the local imaging point. The thick line in the figure is the local reflection interface, the dotted line AB is the normal of the reflection interface passing through the reflection point, and ∠S1RA and ∠S2RB are the local incident angles of the incident wave field of the left and right shot points, respectively (angle range 0~90°).

[0064] Let the local inclination be Angle range 0~180°), to form a folded wave, at least: or , that is, at the imaging point, the local structural inclination angle + the incident angle is greater than 90°.

[0065] For the same reflection point, incident waves from the left and right have opposite reflection coefficients, and direct superposition will weaken the imaging effect of steep structures. Therefore, when performing back-reflected wave imaging, it is necessary to further distinguish whether the wave field is incident from the left or right. Figure 11 In the example, the dip angle obtained from the left is greater than 90°, while the dip angle from the right is less than 90°. This means that we only need to determine the local dip angle to easily distinguish the incident direction of the wavefield. Using the reflection coefficient obtained for incidence angles less than 90° as a criterion, we perform phase inversion processing on the imaging results for incidence angles greater than 90°, resulting in phase-aligned imaging.

[0066] ②Calculate the local formation dip;

[0067] According to the calculated wave field local ray parameters Calculate local formation dip , the cosine of the angle can be expressed as:

[0068] ;

[0069] in, are the z-direction components of the ray parameters at the source and detection ends respectively.

[0070] ③ If the formation dip is greater than 90°, the back-reflected wave single-shot imaging result is calculated;

[0071] If the local formation dip is greater than 90°, the back-reflection wave imaging condition is met and the imaging point can be imaged by the back-reflection wave. Otherwise, the local imaging reflection angle needs to be calculated and determined. The back-reflection wave imaging value of the point is obtained using the following formula:

[0072] ;

[0073] Where I1 is the reflected wave imaging body when the local structural inclination is greater than 90°. is a sign function, which is used to ensure that the imaging phase is consistent and in-phase superposition is achieved when the reflected wave is incident from the left and right sides of the structure; is the local formation dip; is the reflected wave field; is the reflected wave field; x is the offset; z is the wave field depth; t is the time; dt is the time derivative; is the wave field propagation time.

[0074] ④ If the formation dip angle is less than 90°, the imaging reflection angle imaging result is calculated;

[0075] If the local formation dip is less than 90°, it is necessary to calculate the local imaging reflection angle and the local structural dip calculated in step 2 to comprehensively determine whether the back-reflection wave is imaged. If the local imaging reflection angle plus the local structural dip is greater than 90°, the back-reflection wave imaging condition is met and the imaging point can be imaged by the back-reflection wave, otherwise it cannot. The formula for calculating the cosine of the local imaging reflection angle is as follows:

[0076] ;

[0077] in, is the local imaging reflection angle, are the local ray parameters of the source and detector wave fields respectively; when When the angle is greater than 90°, the following imaging conditions apply:

[0078] ;

[0079] Wherein, I2 is the image of the reflected wave when the local structural dip angle is less than 90°; is the local imaging reflection angle; is the local formation dip; is the reflected wave field; is the reflected wave field; x is the shot detection distance; z is the wave field depth; t is time.

[0080] ⑤ Calculate the back-reflected wave imaging results;

[0081] The single shot reflection wave imaging is obtained by summing the reflection waves calculated in ③ and ④:

[0082] ;

[0083] Among them, I is the single-shot back-reflection wave imaging body, I1 is the back-reflection wave imaging body when the local structural inclination angle is greater than 90°, and I2 is the back-reflection wave imaging body when the local structural inclination angle is less than 90°.

[0084] ⑥ Accumulate the single shot results into the final imaging result;

[0085] The following formula is used to calculate the multi-shot back-reflected wave imaging results:

[0086] ;

[0087] in, is a multi-shot back-reflected wave imaging body, I is a single-shot back-reflected wave imaging body, is the minimum value of multiple shots, is the maximum value of multiple guns, It is the multi-gun value.

[0088] Repeat steps ① to ⑥ until all guns are calculated.

[0089] S4. Imaging Processing and Effect Verification

[0090] The imaging effect is verified through numerical simulation and actual data testing, and the imaging effects of horizontal structures and high-steep structures are analyzed.

[0091] ① Numerical simulation test

[0092] The single shot imaging results are shown in the following figure: Figure 12 As shown, it can be seen that this imaging condition can image the reflected wave very well.

[0093] The FootHill model was then used to verify the effectiveness of the imaging method in complex models. Figure 13 The FootHill model reflection wave imaging profile shows that the reflection wave imaging quality is poor in the high and steep structure area. Figure 14 This is the FootHill model back-reflection wave imaging section, and its main imaging area is the high and steep structure area. Figure 15 This is the full-wavefield imaging section of the FootHill model, which can well image both horizontal structures and high-steep structures.

[0094] ② Actual work area test

[0095] First, the observation system of Shunbei Work Area was analyzed. Its maximum offset was 8178m, and the offset range was 4~8km. Actual data test was conducted on Shunbei Work Area, and back-reflection wave data was collected. The data was imaged using conventional reverse time migration, back-reflection wave imaging alone, and back-reflection wave enhanced imaging. Figures 16-20 As shown in Figure 2, it is obvious that the back-reflection wave can well image the high and steep strike-slip fault.

[0096] Figure 21 , Figure 22 As shown in the signal-to-noise ratio analysis diagram of a strike-slip fault of the present invention, the strike-slip fault area was selected for spectrum analysis, and the signal-to-noise ratio was calculated using the spectral value method. The signal-to-noise ratio of the image after the reflected wave enhancement (29.76) was increased by 14.7% compared with the image of the conventional two-way wave extension (25.96).

[0097] Those skilled in the art will understand that the discussion of the above embodiments is merely illustrative and is not intended to limit the scope of the present invention to these examples. Within the spirit and principles of the present invention, the technical features of the above embodiments or different embodiments may be combined, the steps may be implemented in any order, and many other variations exist for the various aspects of the present invention described above, which are not provided in detail for the sake of clarity. Any omissions, modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention shall be included within the scope of protection of the present invention.

Claims

1. A complex geological body imaging method based on reverse time migration of return waves, characterized in that: The following steps are involved: S1. Analysis of the characteristics of the reflected wave field Based on seismological principles, the formation mechanism of the back-reflection wave and the geological conditions required are analyzed. The complex geological bodies in the target work area are simulated and analyzed to determine their impact on the back-reflection wave propagation path and wave field distribution characteristics. Geological-velocity models of different types of high- and steep-angle structures are established, and forward modeling of the back-reflection wave is performed to obtain and analyze its wave field characteristics, providing a basis for observation system design and migration imaging. S2. Design and optimization of the back-reflected wave observation system In combination with the analysis results of step S1 and the geological objectives of the work area, the observation system parameters are optimized through forward simulation tests; using seismic record separation technology, the reflected wave containing reflection information, i.e., the reflected reflected wave, and the reflected wave not containing reflection information, i.e., the non-reflected reflected wave, are extracted from the original wave field respectively; the non-reflected reflected wave field is back-propagated and its wave field square is calculated to obtain the source end illumination characterizing the energy distribution of the earthquake source; the reflected reflected wave field is back-propagated and its wave field square is calculated to obtain the detector end illumination characterizing the energy distribution received by the detector; the source end illumination and the detector end illumination are used to evaluate and optimize the illumination effect of the observation system on the target geological body; S3, Reverse Time Migration Imaging A reverse time migration method based on a two-way wave equation is used for wavefield extension. In the time domain, the source wavefield is forward extended, and the received receiver wavefield is reversely extended. During the wavefield extension process, local ray parameters of the source wavefield and the receiver wavefield at each point and time in space are calculated and extracted. Based on the local ray parameters, the local formation dip at the imaging point is calculated. At the imaging point, different imaging conditions are applied according to the magnitude of the local formation dip to extract the back-reflected wave imaging value. The back-reflected wave imaging values at each imaging point that meets the conditions are accumulated to obtain a single-shot back-reflected wave imaging result. The back-reflected wave imaging results of all single shots are superimposed to obtain the final multi-shot back-reflected wave imaging volume. S4. Imaging Processing and Effect Verification The back-reflected wave imaging result obtained in step S3 is subjected to noise suppression post-processing; the imaging effect of the method is verified using numerical simulation data and actual seismic acquisition data of the target work area; and the quality of the imaging result is quantitatively evaluated by signal-to-noise ratio analysis.

2. The complex geological body imaging method based on reverse time migration of return waves according to claim 1, characterized in that: In step S2, the calculation formulas for the source end illumination and the detection end illumination are: ； in, is the illumination at the source end; is the illumination intensity at the detection end; is the non-reflected back-reflected wave field, is the reflected wave field, x is the offset, z is the wave field depth; is the corresponding frequency.

3. The complex geological body imaging method based on reverse time migration of return waves according to claim 1, characterized in that: The specific steps in step S3 are: Extrapolate the source wave field and the receiver wave field to extract the local ray parameters of the wave field; The wave equation is solved using the finite-difference time-domain method to achieve forward extrapolation of the source wave field in the time direction and reverse extrapolation of the receiver wave field in the time direction. The wave field is , the velocity model is , the wave equation is expressed as: ； in are the second-order partial derivatives of the wave field in directions j, k, and l, and the second-order partial derivatives with respect to time t; The wave equation (3) is solved using the finite difference time domain to obtain the wave field values at the source and detector at each spatial position at each time. The propagation direction of the wave field at the source and detector is expressed as follows: ； in 、 are the local ray parameters of the source and detector wave fields, respectively; x is the offset; z is the wave field depth; t is the time; are the first-order partial derivatives of the wave field x, z, t and the first-order partial derivative with respect to time t, are the first-order partial derivatives of the wave field x, z, t and the first-order partial derivatives with respect to x and z respectively; According to the calculated wave field local ray parameters 、 Calculate local formation dip , the cosine of the angle is expressed as: ； in, 、 are the z-direction components of the ray parameters at the source end and the detection end respectively; If the formation dip is greater than 90°, the back-reflected wave single-shot imaging result is calculated; If the local formation dip is greater than 90°, the back-reflection wave imaging condition is met and the imaging point is imaged by the back-reflection wave. Otherwise, the local imaging reflection angle needs to be calculated and determined. The back-reflection wave imaging value of the point is obtained using the following formula: ； Where I1 is the reflected wave imaging body when the local structural inclination is greater than 90°. is a sign function, which is used to ensure that the imaging phase is consistent and in-phase superposition is achieved when the reflected wave is incident from the left and right sides of the structure; is the local formation dip; is the reflected wave field; is the reflected wave field; x is the offset; z is the wave field depth; t is the time; dt is the time derivative; is the wave field propagation time; If the formation dip angle is less than 90°, the imaging reflection angle imaging result is calculated; If the local formation dip is less than 90°, it is necessary to calculate the local imaging reflection angle and the local structural dip to comprehensively determine whether the back-reflection wave is imaged. If the local imaging reflection angle plus the local structural dip is greater than 90°, the back-reflection wave imaging condition is met and the imaging point is imaged by the back-reflection wave. Otherwise, it cannot be imaged. The formula for calculating the cosine of the local imaging reflection angle is as follows: ； in, is the local imaging reflection angle, are the local ray parameters of the source and detector wave fields respectively; when When the angle is greater than 90°, the following imaging conditions apply: ； Wherein, I2 is the image of the reflected wave when the local structural dip angle is less than 90°; is the local imaging reflection angle; is the local formation dip; is the reflected wave field; is the reflected wave field; x is the shot detection distance; z is the wave field depth; t is the time; The single-shot reflected wave imaging is obtained by accumulating the reflected waves calculated by formula (7) and formula (9): ； Where I is the single shot back-reflection wave imaging volume; I1 is the back-reflection wave imaging volume when the local structural inclination angle is greater than 90°; I2 is the back-reflection wave imaging volume when the local structural inclination angle is less than 90°; The single-shot results are accumulated into the final imaging result, and the multi-shot back-reflected wave imaging result is calculated using the following formula: ； in, is a multi-shot back-reflected wave imaging body, I is a single-shot back-reflected wave imaging body, is the minimum value of multiple shots, is the maximum value of multiple guns, It is the multi-gun value.

Images