Global phase correction and high-pass filtered seismic wave reflection matrix imaging method and system

Abstract

The application discloses a global phase correction and high-pass filtering seismic wave reflection matrix imaging method and system, belongs to the technical field of geophysical exploration and seismic signals, and comprises the following steps: firstly, seismic array data is acquired and a reflection matrix is constructed; wideband focusing reflection matrix is obtained through time Fourier transform and Green function propagation; then, the focusing matrix is converted into a wave number difference domain, global phase correction is realized by using a power method, low wave number components are inhibited by combining high-pass spectral filtering, and wave field aberration correction is completed through closed loop iteration; finally, the corrected matrix is converted back into a spatial domain and normalized, and a high-resolution reflectivity image is generated. The application can overcome the interference of strong scattering medium, improve imaging resolution and algorithm convergence speed, is suitable for magma channel detection in a volcanic area and high-precision depiction of a non-homogeneous reservoir and a fracture network in oil and gas exploration, and has scientific research and industrial application values.

Description

Technical Field

[0001] This invention relates to the fields of geophysical exploration and seismic signal technology, and in particular to a global phase correction and high-pass filtering seismic wave reflection matrix imaging method and system. Background Technology

[0002] Studying the subsurface structure of complex geological formations is of great significance for volcanic mechanism analysis, resource exploration, and potential geological risk assessment. The subsurface media in volcanic activity zones and oil and gas exploration zones exhibit strong heterogeneity, which triggers intense seismic wave scattering, resulting in a significant reduction in the quality of the coherent wavefield. This limits the imaging effectiveness of conventional seismic imaging techniques such as travel-time tomography and migration imaging, making it difficult to accurately depict the geometry of deep magma channels or oil and gas reservoirs.

[0003] Seismic reflection matrix imaging (SEM) extracts volume wave reflection information by cross-correlating seismic environmental noise or active source data to construct a focused reflection matrix, making it a promising alternative technology. However, traditional SEM methods have significant drawbacks: in subsurface environments dominated by low wavenumber components, traditional methods exhibit slow convergence speeds and may even fail to correct aberrations; furthermore, existing point-by-point phase correction strategies lack stability when processing complex scattered wave fields, failing to meet the high-precision imaging requirements of strongly scattering media.

[0004] To address the shortcomings of existing technologies, there is an urgent need for a seismic wave reflection matrix imaging method and system based on global phase correction and high-pass filtering. By introducing global phase correction and high-pass spectral filtering techniques, the strong scattering effect can be overcome, and the problems of traditional matrix imaging being prone to getting trapped in local minima and slow convergence can be solved, thereby achieving high-resolution structural imaging in strongly scattering crustal environments. Summary of the Invention

[0005] This invention provides a global phase correction and high-pass filtering seismic wave reflection matrix imaging method and system. It introduces global phase correction and high-pass spectral filtering technology to overcome strong scattering effects and solve the problems of traditional matrix imaging being prone to getting trapped in local minima and slow convergence, thereby achieving high-resolution structural imaging in strong scattering crustal environments.

[0006] To achieve the above objectives, the technical solution adopted by the present invention is as follows: This invention provides a global phase correction and high-pass filtering seismic wave reflection matrix imaging method, comprising: S1: Acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions and construct reflection matrices; S2: Perform a time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. Calculate the Green's function propagator using the seismic wave velocity model. Propagate the monochromatic matrix to the focusing depth of the underground target. Coherently sum the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix. S3: The broadband focusing reflection matrix is ​​transformed to the wavenumber domain and mapped to the wavenumber difference domain of the incident and reflected wavefields. A closed-loop single-scattering accumulation algorithm is used to correct wavefield aberrations, which includes global phase correction and high-pass filtering. S4: Convert the corrected reflection matrix back to the spatial domain, perform intensity normalization, and generate a high-resolution underground reflectance image.

[0007] Furthermore, the global phase correction in step S3 specifically involves: extracting the maximum singular value of the mapping matrix and its corresponding singular vector using the power method, reshaping the singular vector into a complex phase map, and performing global phase correction on the incident side; performing the same process on the transpose matrix to estimate the phase aberration on the reflecting side.

[0008] Furthermore, the high-pass spectral filtering in step S3 specifically involves applying a high-pass spectral filter to the wavenumber difference domain matrix during each iteration to suppress low wavenumber components in the spatial wavenumber domain, prevent the algorithm from getting trapped in local minima, and ensure stable convergence of the iteration.

[0009] Furthermore, during the iterative process of the high-pass filtering, closed-loop iterative corrections are performed alternately on the incident and reflection sides until the normalized phase aberration reaches the convergence criterion.

[0010] Furthermore, step S4 includes: normalizing the cross-sectional intensity by using the average value of each depth to compensate for the intensity attenuation with depth, and obtaining a normalized seismic intensity map as a high-resolution subsurface reflectivity image.

[0011] Furthermore, in step S2, the diagonal elements of the broadband focusing reflection matrix correspond to confocal imaging, while the off-diagonal elements capture scattered energy caused by inaccurate velocity models.

[0012] Furthermore, the earthquake waveform data is natural earthquake background noise data, with a target frequency band of 1–5 Hz.

[0013] Furthermore, the seismic waveform data is artificially generated active source seismic data, with a target frequency band of 10–50 Hz and a focusing depth corresponding to the depth range of the oil and gas reservoir.

[0014] This invention also provides a seismic wave reflection matrix imaging system based on global phase correction and high-pass filtering, comprising: The data acquisition module is used to acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions and construct reflection matrices; The focusing matrix construction module is used to perform time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. It uses the seismic wave velocity model to calculate the Green's function propagator, propagates the monochromatic matrix to the focusing depth of the underground target, and coherently sums the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix. The aberration correction module is used to convert the broadband focusing reflection matrix to the wavenumber domain and map it to the wavenumber difference domain of the incident and reflected wavefields. It uses a closed-loop single-scattering accumulation algorithm to perform wavefield aberration correction, which includes global phase correction and high-pass filtering. The imaging output module is used to convert the corrected reflection matrix back to the spatial domain, perform intensity normalization processing, and generate a high-resolution underground reflectance image.

[0015] Furthermore, the aberration correction module uses the power method to extract the maximum singular vector to achieve global phase correction, and employs a high-pass spectral filter to suppress low wavenumber components during the iteration.

[0016] Compared with the prior art, the technical solution disclosed in this invention has the following beneficial effects: This invention constructs a reflection matrix and focuses it to the target depth, refocusing the scattered wavefield to effectively overcome imaging blurring caused by strong scattering and achieve high-resolution structural imaging. It employs a closed-loop single-scattering accumulation algorithm + global phase correction + high-pass filtering to replace unstable point-by-point correction, significantly improving convergence speed and avoiding aberration correction failure. High-pass filtering suppresses low wavenumber components, preventing iteration from getting trapped in local minima and ensuring reliable imaging results. This overcomes the strong scattering effect and solves the problems of traditional matrix imaging being prone to getting trapped in local minima and slow convergence, achieving high-resolution structural imaging in strongly scattering crustal environments. This method and system are not only applicable to deep monitoring of natural subduction zones and fault systems; its high-resolution impedance difference identification capability can also be used in oil and gas exploration to finely characterize reservoir microfracture networks and fluid accumulation interfaces, possessing extremely high engineering application value. Attached Figure Description

[0017] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0018] Figure 1 This is a schematic diagram of the seismic scattering attenuation background and spatial distribution of seismic stations in the application area of ​​this invention (taking the St. Helens Volcanic Area in the United States as an example); Figure 2 This is an example diagram of the normalized phase aberration convergence curve of the improved matrix imaging method in the iterative process according to an embodiment of the present invention; Figure 3 These are three-dimensional high-resolution horizontal slice images of reflectivity at different depths obtained by inversion using the method of the embodiments of the present invention, clearly showing strong reflectors (H1-H7) at different depths. Figure 4 The vertical cross-sectional reflectance image obtained by the method of the present invention reveals multi-layered deep fluid and magma transport channels. Figure 5 This is an explanatory diagram of the conceptual model of underground complex structures, fracture networks, and fluid channels established based on the high-resolution imaging results of embodiments of the present invention; Figure 6 This is a schematic diagram of the global phase correction and high-pass filtering seismic wave reflection matrix imaging method provided in an embodiment of the present invention. Detailed Implementation

[0019] 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 some embodiments of the present invention, and not all embodiments. 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.

[0020] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0021] This invention provides a global phase correction and high-pass filtering seismic wave reflection matrix imaging method and system. It introduces global phase correction and high-pass spectral filtering technology to overcome strong scattering effects and solve the problems of traditional matrix imaging being prone to getting trapped in local minima and slow convergence, thereby achieving high-resolution structural imaging in a strongly scattering crustal environment.

[0022] like Figures 1-6 As shown, this embodiment of the invention provides a global phase correction and high-pass filtering seismic wave reflection matrix imaging method, including: S1: Acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions, and construct reflection matrices.

[0023] S1 can be summarized as: acquiring seismic data and constructing a reflection matrix.

[0024] Specifically, this includes acquiring continuous seismic waveform data recorded by dense seismic arrays (or detectors) in the target area.

[0025] Ø Extract the vertical component data, perform cross-correlation analysis on all possible station pairs for environmental noise (or active sources), and generate cross-correlation functions (CCFs) characterizing the impulse response between stations.

[0026] These cross-correlation functions are combined into a reflection matrix, which physically represents the seismic wave field recorded by other receivers when a pulse source is excited at a given receiver location.

[0027] The earthquake waveform data is natural earthquake background noise data, with a target frequency band of 1–5 Hz.

[0028] The seismic waveform data is artificially generated active source seismic data, with a target frequency band of 10–50 Hz and a focusing depth corresponding to the depth range of the oil and gas reservoir.

[0029] S2: Perform a time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. Calculate the Green's function propagator using the seismic wave velocity model, propagate the monochromatic matrix to the focusing depth of the underground target, and coherently sum the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix.

[0030] S2 can be summarized as: constructing the focused reflection matrix.

[0031] Specifically: Perform a time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band (e.g., 1–5 Hz).

[0032] Green's function propagator is calculated using a one-dimensional or three-dimensional seismic wave velocity model, which propagates a monochromatic matrix in the transmission and reception directions to the target focusing depth underground.

[0033] By coherently summing the focusing matrices of all frequencies, a broadband focusing reflection matrix is ​​obtained. Its diagonal elements correspond to confocal imaging, while the off-diagonal elements capture scattered energy caused by inaccuracies in the velocity model.

[0034] In step S2, the diagonal elements of the broadband focusing reflection matrix correspond to confocal imaging, while the off-diagonal elements capture scattered energy caused by inaccurate velocity models.

[0035] S3: The broadband focusing reflection matrix is ​​transformed to the wavenumber domain and mapped to the wavenumber difference domain of the incident and reflected wavefields. A closed-loop single-scattering accumulation algorithm is used to correct wavefield aberrations, which includes global phase correction and high-pass filtering.

[0036] S3 can be summarized as: applying a closed-loop single-scattering accumulation algorithm for wavefield aberration correction.

[0037] Specifically, the broadband focusing reflection matrix is ​​transformed into the wavenumber domain and then mapped to the wavenumber differences domain of the incident and reflected wave fields to highlight the single scattering component.

[0038] Specifically, this can be divided into: Step 3.1 (Global Phase Correction): The maximum singular value of the mapping matrix and its corresponding singular vector are extracted using the power method, and it is reshaped into a complex phase map for global phase correction on the incident side. The same process is performed on the transpose matrix to estimate the phase aberration on the reflecting side.

[0039] Step 3.2 (High-pass spectral filtering): In each iteration, a high-pass spectral filter is applied to the wavenumber difference domain matrix to suppress low wavenumber components in the k-space, prevent the algorithm from getting trapped in local minima, and ensure stable convergence.

[0040] The closed-loop iterative correction is performed alternately on the incident and reflection sides until the calculated normalized phase aberration reaches the convergence criterion.

[0041] The global phase correction in step S3 specifically involves: extracting the maximum singular value of the mapping matrix and its corresponding singular vector using the power method, reshaping the singular vector into a complex phase map, and performing global phase correction on the incident side; performing the same process on the transpose matrix to estimate the phase aberration on the reflecting side.

[0042] The high-pass spectral filtering in step S3 specifically involves applying a high-pass spectral filter to the wavenumber difference domain matrix during each iteration to suppress low wavenumber components in the spatial wavenumber domain, prevent the algorithm from getting trapped in local minima, and ensure stable convergence of the iteration.

[0043] During the iterative process of the high-pass filtering, closed-loop iterative corrections are performed alternately on the incident and reflection sides until the normalized phase aberration reaches the convergence criterion.

[0044] S4: Convert the corrected reflection matrix back to the spatial domain, perform intensity normalization, and generate a high-resolution underground reflectance image.

[0045] S4 can be summarized as: generating high-resolution reflectance images.

[0046] Specifically, after convergence, the corrected reflection matrix is ​​transformed back into the spatial domain.

[0047] The cross-sectional intensity is normalized by averaging the values ​​at various depths to compensate for the intensity attenuation with depth. The resulting normalized seismic intensity map is a high-resolution subsurface reflectivity image.

[0048] Step S4 specifically includes: normalizing the cross-sectional intensity by using the average value of each depth to compensate for the intensity attenuation with depth, and obtaining a normalized seismic intensity map as a high-resolution subsurface reflectivity image.

[0049] Based on the same idea, this invention also provides a seismic wave reflection matrix imaging system based on global phase correction and high-pass filtering, comprising: The data acquisition module is used to acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions and construct reflection matrices; The focusing matrix construction module is used to perform time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. It uses the seismic wave velocity model to calculate the Green's function propagator, propagates the monochromatic matrix to the focusing depth of the underground target, and coherently sums the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix. The aberration correction module is used to convert the broadband focusing reflection matrix to the wavenumber domain and map it to the wavenumber difference domain of the incident and reflected wavefields. It uses a closed-loop single-scattering accumulation algorithm to perform wavefield aberration correction, which includes global phase correction and high-pass filtering.

[0050] The aberration correction module uses the power method to extract the maximum singular vector to achieve global phase correction, and uses a high-pass spectral filter to suppress low wavenumber components during the iteration.

[0051] The imaging output module is used to convert the corrected reflection matrix back to the spatial domain, perform intensity normalization processing, and generate a high-resolution underground reflectance image.

[0052] The specific implementation method is as follows: Example 1: Volcano and fault zone detection based on natural earthquake background noise (taking the St. Helens Volcano area as an example) This embodiment details the specific steps for deep geological structure exploration using natural earthquake environmental noise data: 1. Data Acquisition and Preprocessing: Deploy a dense network of short-period nodes and broadband seismic arrays. Extract continuous noise data of the vertical component, and after mean removal, trend removal, instrument response correction, and frequency-time normalization (e.g., bandpass filtering to the 1-5 Hz band), calculate the cross-correlation functions (CCFs) of all station pairs.

[0053] 2. Reflection Matrix Focusing: Approximately 40,000 CCFs are constructed into a reflection matrix. The Green's function is calculated using a regional one-dimensional (or three-dimensional) Vp velocity model to "synthesize and reposition" the station from the surface to focal points at various underground depths, obtaining a broadband focused reflection matrix.

[0054] 3. Aberration Correction and Iteration: A closed-loop single-scattering accumulation algorithm is applied. The broadband focusing matrix is ​​transformed to the wavenumber difference domain, and the dominant singular vector is extracted using the power method to estimate the global phase aberration. Key step: In each iteration, a high-pass spectral filter is introduced to remove low wavenumber components, preventing the iteration from getting trapped in local minima. After approximately 15 alternating iterations (incident and reflection sides), the phase aberration achieves stable convergence.

[0055] 4. Imaging and Geological Interpretation: The corrected matrix was converted back to the spatial domain and intensity normalized. The imaging results not only identified fluid enrichment and degassing channels (H1) that corresponded to the earthquake swarm in the shallow layer (~3.5-7 km), but also precisely located the broad magma-fluid accumulation platform and mechanically sensitive valve at the top of the Siletz terrane (H2) at a depth of 10 km, revealing a multi-layered magma conduit system.

[0056] Example 2: Detection of Heterogeneous Reservoirs and Fracture Networks in Oil and Gas Basins Based on Artificial Seismic Sources This embodiment illustrates the specific application steps of extending this method to industrial oil and gas exploration scenarios: 1. Data Acquisition and Reflection Matrix Construction: In oil and gas exploration areas, a high-density surface geophone network is used to receive active source seismic wave signals excited by explosives or controlled seismic sources. Body wave reflection data is extracted, and a multi-source, multi-receiver spatial reflection matrix is ​​directly constructed based on the geometric relationship between the shot point and the geophone.

[0057] 2. Bandpass Filtering and Focusing of Target Layer: Bandpass filtering is applied to the reflection matrix based on the characteristic frequencies of the oil and gas reservoir (e.g., 10-50 Hz). Combined with existing 3D logging or pre-stack migration velocity models, the Green's function is calculated to focus and extend the reflected wave field downwards to the target oil and gas reservoir depth range (e.g., 2-4 km).

[0058] 3. Anti-scattering aberration correction: To address the high-frequency, strong scattering background caused by carbonate karst pores or shale fractures, this invention employs a closed-loop algorithm that includes global phase correction and high-pass spectral filtering. High-pass filtering at this stage significantly suppresses low-frequency, long-wavelength distortions caused by complex overlying strata (such as salt domes or shallow gas clouds), rapidly converging and restoring the true single-scattering high-frequency reflection information within the reservoir.

[0059] 4. High-Resolution Reservoir Characterization and Evaluation: Corrected and normalized high-frequency reflectivity images are used for reservoir evaluation. Strong reflectivity anomalies and local clumps in the imaging results accurately correspond to strong impedance difference boundaries (such as water / gas-bearing fluid contact surfaces) or highly developed microfracture networks. This super-resolution imaging overcomes the "beaded" blurring effect of conventional migration methods in fractured-vuggy reservoirs, providing high-precision geophysical evidence for horizontal well trajectory design and sweet spot selection in fracturing.

[0060] The basic principles of the present invention have been described above with reference to specific embodiments. However, it should be noted that the advantages, benefits, and effects mentioned in the present invention are merely examples and not limitations, and should not be considered as essential features of each embodiment of the present invention. Furthermore, the specific details disclosed above are for illustrative and facilitative purposes only, and are not limitations. These details do not limit the present invention to the necessity of employing the aforementioned specific details.

[0061] The block diagrams of devices, apparatuses, devices, and systems involved in this invention are merely illustrative examples and are not intended to require or imply that they must be connected, arranged, or configured in the manner shown in the block diagrams. As those skilled in the art will recognize, these devices, apparatuses, devices, and systems can be connected, arranged, and configured in any manner. Words such as “comprising,” “including,” “having,” etc., are open-ended terms meaning “including but not limited to,” and are used interchangeably with them. The terms “or” and “and” as used herein refer to the terms “and / or,” and are used interchangeably with them unless the context clearly indicates otherwise. The term “such as” as used herein refers to the phrase “such as but not limited to,” and is used interchangeably with it.

[0062] It should also be noted that in the apparatus, device, and method of the present invention, the components or steps can be disassembled and / or recombined. These disassemblies and / or recombinations should be considered as equivalent solutions of the present invention.

[0063] The above description of the disclosed aspects is provided to enable any person skilled in the art to make or use the invention. Various modifications to these aspects will be readily apparent to those skilled in the art, and the general principles defined herein can be applied to other aspects without departing from the scope of the invention. Therefore, the invention is not intended to be limited to the aspects shown herein, but rather to be carried out within the widest scope consistent with the principles and novel features disclosed herein.

[0064] It should be understood that the qualifying terms "first", "second", "third", "fourth", "fifth" and "sixth" used in the description of the embodiments of the present invention are only used to more clearly illustrate the technical solutions and are not intended to limit the scope of protection of the present invention.

[0065] The above description has been given for purposes of illustration and description. Furthermore, this description is not intended to limit the embodiments of the invention to the forms disclosed herein. Although numerous exemplary aspects and embodiments have been discussed above, those skilled in the art will recognize certain variations, modifications, alterations, additions, and sub-combinations therein.

Claims

1. A global phase correction and high-pass filtering seismic wave reflection matrix imaging method, characterized in that, include: S1: Acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions and construct reflection matrices; S2: Perform a time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. Calculate the Green's function propagator using the seismic wave velocity model. Propagate the monochromatic matrix to the focusing depth of the underground target. Coherently sum the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix. S3: The broadband focusing reflection matrix is ​​transformed to the wavenumber domain and mapped to the wavenumber difference domain of the incident and reflected wavefields. A closed-loop single-scattering accumulation algorithm is used to correct wavefield aberrations, which includes global phase correction and high-pass filtering. S4: Convert the corrected reflection matrix back to the spatial domain, perform intensity normalization, and generate a high-resolution underground reflectance image.

2. The method according to claim 1, characterized in that, The global phase correction in step S3 specifically involves: extracting the maximum singular value of the mapping matrix and its corresponding singular vector using the power method, reshaping the singular vector into a complex phase map, and performing global phase correction on the incident side; performing the same process on the transpose matrix to estimate the phase aberration on the reflecting side.

3. The method according to claim 1, characterized in that, The high-pass spectral filtering in step S3 specifically involves applying a high-pass spectral filter to the wavenumber difference domain matrix during each iteration to suppress low wavenumber components in the spatial wavenumber domain, prevent the algorithm from getting trapped in local minima, and ensure stable convergence of the iteration.

4. The method according to claim 3, characterized in that, During the iterative process of the high-pass filtering, closed-loop iterative corrections are performed alternately on the incident and reflection sides until the normalized phase aberration reaches the convergence criterion.

5. The method according to claim 1, characterized in that, Step S4 includes: normalizing the cross-sectional intensity by using the average value of each depth to compensate for the intensity attenuation with depth, and obtaining a normalized seismic intensity map as a high-resolution subsurface reflectivity image.

6. The method according to claim 1, characterized in that, In step S2, the diagonal elements of the broadband focusing reflection matrix correspond to confocal imaging, while the off-diagonal elements capture scattered energy caused by inaccurate velocity models.

7. The method according to claim 1, characterized in that, The earthquake waveform data is natural earthquake background noise data, with a target frequency band of 1–5 Hz.

8. The method according to claim 1, characterized in that, The seismic waveform data is artificially generated active source seismic data, with a target frequency band of 10–50 Hz and a focusing depth corresponding to the depth range of the oil and gas reservoir.

9. A seismic wave reflection matrix imaging system based on global phase correction and high-pass filtering, characterized in that, include: The data acquisition module is used to acquire continuous seismic waveform data recorded by seismic arrays or detectors in the target area, extract vertical component data, perform environmental noise or active source cross-correlation analysis on all stations, generate cross-correlation functions and construct reflection matrices; The focusing matrix construction module is used to perform time Fourier transform on the reflection matrix to generate a single-frequency reflection matrix within the target frequency band. It uses the seismic wave velocity model to calculate the Green's function propagator, propagates the monochromatic matrix to the focusing depth of the underground target, and coherently sums the focusing matrices of all frequencies to obtain the broadband focusing reflection matrix. The aberration correction module is used to convert the broadband focusing reflection matrix to the wavenumber domain and map it to the wavenumber difference domain of the incident and reflected wavefields. It uses a closed-loop single-scattering accumulation algorithm to perform wavefield aberration correction, which includes global phase correction and high-pass filtering. The imaging output module is used to convert the corrected reflection matrix back to the spatial domain, perform intensity normalization processing, and generate a high-resolution underground reflectance image.

10. The system according to claim 9, characterized in that, The aberration correction module uses the power method to extract the maximum singular vector to achieve global phase correction, and uses a high-pass spectral filter to suppress low wavenumber components during the iteration.

Images