An ultra-deep reflection imaging method based on passive source seismic data and related device

Abstract

The application discloses a kind of ultra-deep reflection imaging method and related device based on passive source seismic data, it is related to deep structure imaging technical field, the method is related to the following steps: first, for original continuous seismic waveform data, the seismic waveform data of P-wave first arrival and S-wave non-arrival are preprocessed, to obtain teleseismic event data section.Subsequently, high signal-to-noise ratio teleseismic event is identified using time-frequency-wave number domain analysis, and those data greatly affected by noise are removed, then, using multidimensional deconvolution technique, the reconstruction of reflection wave field is carried out to high-quality teleseismic data segment and deconvolution source wavelet, to obtain ultra-deep reflection wave field data.Finally, reflection imaging processing is carried out in combination with velocity analysis and stacking imaging technology, the analysis ability of deep reflection interface is optimized, and finally the deep reflection imaging profile of earth is obtained, the continuous structural interface extending from shallow crust to deep mantle can be clearly revealed, and the coupling characteristics between shallow and deep earth structure can be effectively captured.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of deep structure imaging, in particular to an ultra-deep reflection imaging method based on passive source seismic data and a related device. BACKGROUND

[0002] To image the deep structure of the lithosphere, the mantle transition zone (MTZ), the lower mantle and even the core, researchers have adopted a variety of passive source detection techniques, including SS / PP precursor waves, receiver functions, body wave travel time tomography and triple seismic phase methods. Each technique has its unique advantages: tomography can reveal the velocity structure of a large area, but its resolution is relatively low; SS / PP precursor waves can cover a global range, but their frequency is low (<0.05Hz), making it difficult to analyze small-scale structural details; receiver functions have a higher frequency (0.2-1Hz) and can capture the characteristics of interfaces (such as arrival time, amplitude, wave impedance strength), but their detection range is limited to the local area under the station, and they are easily disturbed by shallow structures, and are not sensitive enough to gradual interfaces; triple seismic phase methods are more sensitive to structures with gradual velocity changes, and are between global and regional scales, but their imaging relies on waveform fitting and the distribution of seismic events. At present, the imaging accuracy and detail resolution capability of these methods cannot meet the needs of accurately depicting the deep structure of the mantle transition zone and its dynamic mechanisms, which directly affects the understanding of the characteristics of the deep mantle. Therefore, it is urgent to develop innovative imaging techniques with high resolution to accurately depict the complex seismic wave velocity discontinuity structure within the mantle transition zone, analyze the distribution characteristics of its layered structure, and further reveal the material evolution mechanisms within the mantle transition zone. SUMMARY

[0003] The purpose of the present application is to provide an ultra-deep reflection imaging method based on passive source seismic data and a related device, which can provide clearer interface imaging from the shallow crust to the deep mantle than receiver functions and tomography.

[0004] To achieve the above-mentioned purpose, the present application provides the following solutions:

[0005] In a first aspect, the present application provides an ultra-deep reflection imaging method based on passive source seismic data, comprising:

[0006] Obtaining original continuous seismic waveform data, intercepting seismic waveform data with P-wave first arrival and S-wave non-arrival and pre-processing to obtain teleseismic event data segments; the pre-processing includes removing instrument response, normalization and band-pass filtering.

[0007] For teleseismic event data segment, high signal-to-noise ratio teleseismic events are identified by time-frequency-wave number domain analysis method, and data with signal-to-noise ratio lower than a preset threshold is removed to obtain high-quality teleseismic data segments.

[0008] The pure P-wave multidimensional deconvolution method is adopted to reconstruct the reflection wave field of the high-quality teleseismic data segments and deconvolve the source wavelet to obtain super-deep reflection wave field data.

[0009] Reflection imaging processing is performed on the super-deep reflection wave field data, and the analysis capability of the deep reflection interface is optimized by combining velocity analysis and stacking imaging to obtain the deep reflection imaging profile of the earth; in the stacking imaging, the common midpoint domain stacking is selected for regular observation system, and the common shot point domain stacking or common receiver point domain stacking is selected for irregular observation system.

[0010] Optionally, the original continuous seismic waveform data is obtained, the seismic waveform data with P-wave first arrival and S-wave non-arrival is intercepted and preprocessed to obtain the teleseismic event data segment, which specifically includes:

[0011] The original continuous seismic waveform data is obtained, and the seismic event data segment is extracted from the original continuous seismic waveform data according to the public earthquake occurrence time.

[0012] According to the P-wave first arrival time of the seismic event data segment, the seismic waveform data with P-wave first arrival and S-wave non-arrival is intercepted and preprocessed to obtain the teleseismic event data segment.

[0013] Optionally, high signal-to-noise ratio teleseismic events are identified according to the following formula:

[0014]

[0015] Wherein, TXFK represents high signal-to-noise ratio teleseismic events identified by time-frequency-wave number domain analysis method, X(t,x) is seismic waveform data in time-space domain, t is time, x is receiver position, f TX (X) is a feature extraction function in time-space domain, which captures the characteristics of seismic waveform data X, F(k,ω) is the amplitude distribution information extracted in frequency-wave number domain, k is wave number, ω is frequency, And are the preset screening criteria of high signal-to-noise ratio seismic waveform data, represents that the seismic waveform data contains high signal-to-noise ratio horizontal body wave signals, represents that the main energy of the seismic waveform data is distributed in the body wave area of the frequency-wave number spectrum.

[0016] Optionally, the reflection wave field reconstruction is performed according to the following formula:

[0017]

[0018] in, For effective reflection range, x represents R The virtual traction source at x A The z-component of the elastic dynamic reflection response of the P-wave induced at the receiver, V z (x R ,x B ,ω) represents x B The virtual traction source at x R The z-component of the complete wave field caused by the receiver at that location, d 2 x R V is the boundary area element. z (x A ,x B ,ω) represents x B The virtual traction source at x A The z-component of the complete wave field caused by the receiver at that location. x represents A The z-component of the P-wave field caused by the receiver at the location, where ω is the angular frequency.

[0019] Optionally, the frequency of the source wavelet is determined based on the maximum inter-station spacing; when the maximum inter-station spacing is 1000km, the source wavelet is selected as the Ricker wavelet with a main frequency of 0 to 0.4Hz.

[0020] Optionally, the reflection imaging processing performed on the ultra-deep reflection wavefield data includes denoising, static correction, and dynamic correction.

[0021] Secondly, this application provides an ultra-deep reflection imaging system based on passive source seismic data, comprising:

[0022] The seismic data acquisition and preprocessing module is used to acquire raw continuous seismic waveform data, extract seismic waveform data with P-wave arrival and S-wave arrival, and perform preprocessing to obtain teleseismic event data segments; the preprocessing includes instrument response removal, normalization, and bandpass filtering.

[0023] The high signal-to-noise ratio teleseismic data filtering module is used to identify high signal-to-noise ratio teleseismic events in teleseismic event data segments through time-frequency-wavenumber domain analysis and remove data with a signal-to-noise ratio lower than a preset threshold to obtain high-quality teleseismic data segments.

[0024] The ultra-deep reflection wavefield reconstruction module is used to reconstruct the reflection wavefield of high-quality teleseismic data segments and convolve the source wavelet using the pure P-wave multidimensional deconvolution method to obtain ultra-deep reflection wavefield data.

[0025] The earth deep profile reflection imaging module is used for reflection imaging processing of the super-deep reflection wave field data, and combines velocity analysis and stack imaging to optimize the analysis ability of the deep reflection interface, so as to obtain the earth deep reflection imaging profile; in the stack imaging, the common point domain stack is selected for the regular observation system, and the common shot point domain stack or the common receiver point domain stack is selected for the irregular observation system.

[0026] In a third aspect, the present application provides a computer device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the steps of the method for super-deep reflection imaging based on passive source seismic data.

[0027] In a fourth aspect, the present application provides a computer readable storage medium, which stores a computer program, wherein the computer program is executable on a processor to implement the steps of the method for super-deep reflection imaging based on passive source seismic data.

[0028] In a fifth aspect, the present application provides a computer program product, comprising a computer program, wherein the computer program is executable on a processor to implement the steps of the method for super-deep reflection imaging based on passive source seismic data.

[0029] According to the embodiments provided in the present application, the following technical effects are disclosed:

[0030] The present application provides a method for super-deep reflection imaging based on passive source seismic data and related devices, in which, for the obtained original continuous seismic waveform data, the seismic waveform data with P-wave first arrival and S-wave non-arrival is intercepted and preprocessed to obtain teleseismic event data segments, the physical separation of seismic phases is realized, and the reconstruction and imaging of pure P-wave super-deep reflection wave field are met; then, high signal-to-noise ratio teleseismic events are identified by the time-frequency-wave number domain analysis method, and data with large noise influence is removed to obtain high-quality teleseismic data segments, which effectively improves the reliability of seismic event screening and reduces the interference of noise and scattered wave energy on imaging; then, the pure P-wave multi-dimensional deconvolution method is used for reflection wave field reconstruction and deconvolution of source wavelet, which effectively removes the reconstruction artifacts introduced by source mechanism, distribution and uneven energy of the correlation method, and improves the signal-to-noise ratio; finally, in the reflection imaging, velocity analysis and stack imaging are combined, and for the irregular observation system, the common shot point domain stack or the common receiver point domain stack is selected, which also provides a reliable imaging method for reflection imaging of uneven receiver spacing and large offset distance. In summary, the present application can provide clearer mantle transition zone interface imaging than the receiver function and tomographic imaging, realize full-scale imaging from the surface to the mantle and the core, clearly reveal the continuous tectonic interface extending from the shallow crust to the deep mantle, and effectively capture the coupling characteristics between the shallow and deep earth structures. BRIEF DESCRIPTION OF DRAWINGS

[0031] In order to more clearly illustrate the technical solutions in the embodiments of the present application or the prior art, the accompanying drawings needed in the embodiments will be briefly introduced as follows. Obviously, the accompanying drawings in the following description only only some embodiments of the present application, and for those skilled in the art, other drawings can be obtained based on these drawings without creative labor.

[0032] Figure 1 A flow chart of a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0033] Figure 2 A flow chart of step S1 in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0034] Figure 3 A gray scale map of a two-dimensional waveform of a teleseismic event data segment arranged by receivers in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0035] Figure 4 An F-K spectrum map of a teleseismic event data segment in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0036] Figure 5 A gray scale map of a two-dimensional waveform of a selected P-wave first arrival wave arranged by receivers in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0037] Figure 6 An F-K spectrum map of a selected P-wave first arrival wave in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0038] Figure 7 A reflection imaging profile map of a deep mantle structure generated at different stacking velocities in a super-deep reflection imaging method based on passive source seismic data provided by an embodiment of the present application.

[0039] Figure 8 A functional module schematic diagram of a super-deep reflection imaging system based on passive source seismic data provided by an embodiment of the present application.

[0040] Figure 9 A structural schematic diagram of a computer device provided by an embodiment of the present application. DETAILED DESCRIPTION

[0041] With reference to the drawings of the embodiments of the present application, the technical solutions in the embodiments of the present application will be described clearly and completely. Obviously, the described embodiments are only a part of the embodiments of the present application, but not all the embodiments of the present application. Based on the embodiments of the present application, all other embodiments obtained by a person of ordinary skill in the art without creative effort belong to the scope of protection of the present application.

[0042] As described above, the structure imaging of the mantle transition zone mainly relies on the methods of SS / PP precursor wave, body wave travel time tomography, receiver function and triple seismic phase. These techniques have preliminarily revealed the existence of the 410-km and 660-km interfaces and their partial discontinuous characteristics, depicted the average depth of the two interfaces on a global scale, and found that the strength and width are affected by water content and chemical composition. However, for the complex tectonic region, there is still great uncertainty in the deep morphology and interface fluctuation. In recent years, the progress of multidimensional deconvolution method (MDD) and body wave surface wave separation technology provides a new path for fine imaging of the ultra-deep structure. The present application intends to expand the imaging research of the mantle transition zone and deeper part of the earth on this basis, which has the advantages of being cutting-edge and competitive.

[0043] The ultra-deep reflection imaging method based on passive source seismic data provided by the present application aims to realize high-precision reflection imaging of the deep structure of the lithosphere, the mantle transition zone, the lower mantle and even the earth core. The method reconstructs the ultra-deep reflection wave field of the passive source data (here, the ultra-deep reflection is defined relative to the traditional deep seismic reflection) by the pure P-wave (sound wave) multidimensional deconvolution MDD method, and performs high-precision imaging on the structure of the deep part of the earth.

[0044] In order to make the above-mentioned purposes, features and advantages of the present application more obvious and easy to understand, the present application will be further described in detail below with reference to the drawings and specific embodiments.

[0045] The ultra-deep reflection imaging method based on passive source seismic data provided by the present application, in an exemplary embodiment, as shown in Figure 1 The method comprises the following steps:

[0046] S1, obtaining original continuous seismic waveform data, cutting the seismic waveform data of P-wave first arrival and S-wave non-arrival and pre-processing to obtain teleseismic event data segment; the pre-processing includes removing instrument response, normalization and band-pass filtering. In the embodiment, as shown in Figure 2 The step S1 specifically comprises the following steps:

[0047] S11, obtaining original continuous seismic waveform data, and extracting a seismic event data segment from the original continuous seismic waveform data according to the public seismic occurrence time.

[0048] Specifically, the original continuous seismic waveform data taken can be seismic event data or continuous seismic waveform data of a fixed station of a seismic network, or can come from a mobile broadband station, or even a short-period nodal instrument (because the short-period instrument has limited ability to receive low-frequency signals, it is often used for imaging on the lithospheric scale).

[0049] The original seismic waveform data is obtained by downloading / exporting the miniseed or sac format seismic waveform from the instrument. Because the processing involves multiple seismic detectors, the continuous seismic records of these detectors are arranged in sequence according to the spatial distribution, and then the corresponding two-dimensional waveform data body can be obtained, such as Figure 3 In each column, the continuous seismic amplitude data measured by a detector is arranged, and each row is the data of a sampling point. Because the actual time corresponding to the sampling point in the continuous seismic record is known, the seismic event can be cut off according to the earthquake occurrence time published by the ISC or the seismic network.

[0050] S12, according to the P-wave first arrival time of the seismic event data segment, the P-wave first arrival and S-wave non-arrival seismic waveform data are cut off and preprocessed to obtain the teleseismic event data segment.

[0051] The occurrence time is initially calculated based on the theoretical model by using the commonly used software package TauP. Because there can be differences between the theoretical calculation and the actual observation data, after the initial data is cut off, an additional adjustment is needed. Specifically, according to the waveform information of the seismic event, the actual two-dimensional P-wave first arrival time can be identified by eye, and then a 300-second data window is cut off from this time point for subsequent wave field processing. It is worth noting that the reason for selecting 300 seconds is that according to the TauP theoretical arrival time calculation, the arrival time difference between the first P-wave and the first S-wave exceeds 300 seconds. Therefore, within 300 seconds after the P-wave first arrival, the S-wave phase is not included, thereby realizing the physical separation of the seismic phase. If the theoretical arrival time difference calculation is 500 seconds, then 500 seconds after the P-wave first arrival can be selected for reconstructing the two-way travel time of the reflected wave field.

[0052] Because the obtained seismic data is often affected by noise and instrument response, preprocessing including de-instrument response, normalization, and band-pass filtering is also needed. Figure 3 It is a gray-scale image of the two-dimensional waveform of the teleseismic event data segment after the spatial distribution of the detector is arranged and the preprocessing is performed.

[0053] S2, for the teleseismic event data segment, high signal-to-noise ratio teleseismic events are identified by time-frequency-wave number domain analysis method, and data with signal-to-noise ratio lower than a preset threshold are removed to obtain high-quality teleseismic data segments. Specifically, in this embodiment, high signal-to-noise ratio teleseismic events can be identified according to the following formula:

[0054]

[0055] where TXFK represents high signal-to-noise ratio teleseismic events identified by time- frequency-wavenumber domain analysis, X(t, x) is the time-space domain seismic waveform data, t is time, x is the receiver location, f TX (X) is the time-space domain feature extraction function that captures the characteristics of seismic waveform data X, F(k, ω) is the amplitude distribution information extracted in the frequency-wavenumber domain, k is the wavenumber, ω is the frequency, and respectively, are the preset screening criteria for high signal-to-noise ratio seismic waveform data, indicates that the seismic waveform data contains high signal-to-noise ratio horizontal body wave signals, indicates that the main energy of the seismic waveform data is distributed in the body wave region of the frequency-wavenumber spectrum.

[0056] Passive source data body wave and surface wave separation technology identifies high signal-to-noise ratio body wave signals (abbreviated as TXFK) by analyzing the characteristics of the time-distance (T-X) domain and the frequency-wavenumber (F-K) domain, which effectively improves the reliability of seismic event screening and reduces the interference of noise and scattered wave energy on imaging. Combined with the two-dimensional teleseismic waveform data shown in Figure 3 , the F-K spectrum of the teleseismic event data segment obtained by TXFK is shown in Figure 4 . From the waveform diagram and the F-K spectrum diagram, it can be seen that the body wave capacity (the body wave region below the two dashed lines) is mainly focused on the vicinity of k = 0, and k = 0 represents that these waves have very fast apparent velocity, and high speed usually corresponds to P wave.

[0057] In the process of super-deep reflection wave field reconstruction, the first arrival of P wave needs to be selected as the input of the sound wave case multi-dimensional deconvolution method. At this time, the TXFK method is also used for quality control of manual selection of the first arrival of P wave to prevent the introduction of artifacts in the virtual source super-deep reflection wave field. For example, by time-frequency-wavenumber domain analysis, the wave field energy of the selected first arrival can be evaluated to be mainly concentrated in the position with the fastest speed (wavenumber k = 0), i.e. the middle region of the F-K spectrum. As shown in Figure 5 and Figure 6 are the gray scale diagrams of the two-dimensional waveforms of the P wave first arrival arranged by receivers and the F-K spectrum diagrams thereof, respectively.

[0058] By ensuring that the arrival time difference between the P-wave first arrival and the S-wave first arrival is greater than 300s, irrelevant seismic phases are separated out. It is worth noting that this 300s corresponds to the depth that is intended to be imaged. According to the IASP91 velocity model, the two-way travel time of 300s corresponds to a depth of about 1500km. Therefore, if imaging of deeper parts of the earth is required, the arrival time difference between the P-wave first arrival and the S-wave first arrival needs to be greater, so as to meet the pure P-wave super-deep reflection wave field reconstruction and imaging.

[0059] S3, using a pure P-wave multidimensional deconvolution method, reflection wave field reconstruction and deconvolution of the source wavelet are performed on a high-quality teleseismic data segment to obtain super-deep reflection wave field data.

[0060] Since the processing of conventional exploration seismic data relies on the effective separation of P-waves and S-waves, the traditional stacking or migration imaging technology is difficult to be directly applied to such reflection wave fields. However, for teleseismic events, the arrival of P-waves and S-waves has a significant travel time delay, and the delay amount changes greatly. Therefore, in the embodiment, the P-wave arrival and its subsequent coda can be selected for super-deep reflection wave field reconstruction and imaging, for example, 300 seconds after the arrival of the P-wave first arrival and before the arrival of the S-wave first arrival. Thus, in the case of pure P-waves, only the reflection wave field in the z direction dominated by P-waves needs to be reconstructed, and in the embodiment, the reflection wave field reconstruction is performed according to the following formula:

[0061]

[0062] wherein, is the effective reflection scope, represents the z component of the elastic dynamic reflection response of the P-wave caused by the virtual pulling source at x R at the receiver at x A , V z (x R ,x B ,ω) represents the z component of the complete wave field caused by the virtual pulling source at x B at the receiver at x R , d 2 x R is the boundary area element, V z (x A ,x B ,ω) represents the z component of the complete wave field caused by the virtual pulling source at x B at the receiver at x A , represents the z component of the P-wave field caused by the receiver at x A , and ω is the angular frequency.

[0063] After the super-deep reflection wave field is reconstructed, the source wavelet needs to be convolved. For the study of the mantle structure, such as the mantle transition zone, the main frequency of the wavelet should be 0.1-1 Hz. Specifically, the frequency of the source wavelet is determined according to the maximum inter-station distance; in the embodiment, when the maximum inter-station distance is 1000 km, a Ricker wavelet with a main frequency of 0-0.4 Hz is selected, and the super-deep reflection single-shot record data with high signal-to-noise ratio can be obtained by convolving the obtained super-deep reflection wave field with the source wavelet.

[0064] S4, reflection imaging processing is performed on the super-deep reflection wave field data, and the analysis ability of the deep reflection interface is optimized by combining velocity analysis and stack imaging to obtain a deep reflection imaging profile of the earth; in the stack imaging, the common midpoint domain stack is selected for a regular observation system, and the common shot point domain stack or the common receiver point domain stack is selected for a non-regular observation system.

[0065] Specifically, the reflection imaging processing performed on the super-deep reflection wave field data includes denoising, static correction and dynamic correction.

[0066] The denoising is used to filter the artifacts possibly introduced by the deconvolution and partially eliminate the influence of the converted wave phase, and the main parameters thereof are band-pass filtering (0.01-0.5 Hz) and F-K filtering (k value is selected between the maximum k and the minimum k, and the frequency range is 0-1.9 Hz). The static correction is used to eliminate the difference in wave propagation time caused by factors such as the change in surface elevation, the thickness or velocity change of the low-velocity zone (weathering layer), and the core function thereof is to correct the seismic record to a unified reference surface, ensure the alignment of the reflection wave events, and improve the accuracy of the subsequent processing and the imaging quality. Since the super-deep reflection imaging technology is used for imaging the super-deep part of the earth, the static correction mainly uses the elevation static correction, that is, the main parameter thereof is the elevation of the station.

[0067] The dynamic correction is to correct the travel time difference of the oblique propagation reflection wave (the difference in double-path travel time caused by the different propagation paths) to the travel time of the vertical incidence. The so-called velocity analysis is to analyze which velocity can best flatten the reflection events by trying a series of velocity values to determine the correct velocity model of the underground medium, and the velocity analysis and the dynamic correction are interdependent. In actual processing, the velocity scanning and the dynamic correction are an iterative process. For the imaging of the mantle scale, the commonly used velocity range is 8700-12500 m / s.

[0068] The superimposed imaging is to superimpose the gathers after the motion correction according to the CMP (common midpoint) to enhance the signal and suppress the noise. The common midpoint domain superposition is also commonly referred to as the common reflection point domain (CRP). For a regular observation system (linear, equally spaced arrangement of detectors), the common midpoint domain superposition is commonly used; for a non-regular observation system such as a wide frequency station and a dense array observation, the common shot point domain superposition or the common receiver point domain superposition can be used to solve the problem of failure of the CMP superposition.

[0069] In addition, the converted wave (such as P-S, S-P) interferes with the reflection wave field reconstruction and high-precision imaging, especially superimposes with the main signal in the complex medium, and reduces the accuracy. In order to reduce the interference, the step S1 in the foregoing of the present application preferentially uses the vertical component record in the preprocessing, intercepts the signal before the S wave first arrival, suppresses the S wave and the surface wave, and ensures that the P wave is the main one; in the step S2, the P wave and the converted wave velocity are distinguished based on the frequency-wave number domain, the filter window is adjusted to suppress the low-velocity seismic phase, and the imaging accuracy of the 410-km, 520-km, 660-km or deeper mantle interface is improved.

[0070] As shown in Figure 7 The super-deep reflection imaging method provided by the above-mentioned embodiments of the present application directly and high-precision reflects the deep mantle structure for the first time; and the traditional methods of the receiver function, the tomography and the waveform inversion can only outline the rough imaging profile and are difficult to capture the fine reflection interface of the deep mantle. In Figure 7 , the first block corresponds to the reflection of the upper boundary of the mantle transition zone (about 410 km), also known as the discontinuity reflection; the second block corresponds to the reflection interface of the structure in the mantle transition zone (about 510 km); the third block corresponds to the reflection interface of the lower boundary of the mantle transition zone (about 660 km); and the fourth block corresponds to the reflection interface of the lower mantle structure (about 800 km).

[0071] Compared with the existing methods, the passive source data body wave and surface wave separation technology is used in the present application to improve the teleseismic data screening accuracy, and this technology has a good quality control effect on the first arrival wave selection; the pure P wave multi-dimensional deconvolution method is used to effectively remove the reconstruction artifacts introduced by the source mechanism, distribution and uneven energy of the mutual correlation method, improve the signal-to-noise ratio, and can provide clearer mantle transition zone interface imaging than the receiver function and the tomography. Moreover, the common shot point domain superposition or the common receiver point domain superposition is selected, which is suitable for a non-regular observation system, and provides a reliable imaging method for the reflection imaging of uneven detector spacing and large offset distance. The finally obtained imaging profile is a full-scale imaging from the surface to the mantle deep to the core of the earth, clearly reveals the continuous structural interface extending from the shallow crust to the deep mantle, and effectively captures the coupling characteristics between the shallow and deep earth structures.

[0072] Based on the same inventive concept, the embodiments of the present application also provide a system for implementing the above-mentioned passive-source seismic data-based super-deep reflection imaging method. The system provides a solution to the implementation similar to that described in the above method. In one exemplary embodiment, as shown in Figure 8 FIG. 1, a passive-source seismic data-based super-deep reflection imaging system is provided, comprising:

[0073] a seismic data acquisition and preprocessing module for acquiring original continuous seismic waveform data, intercepting P-wave first arrival and S-wave non-arrival seismic waveform data and pre-processing to obtain teleseismic event data segments; the preprocessing includes de-instrument response, normalization and band-pass filtering.

[0074] a high signal-to-noise ratio teleseismic data screening module for identifying high signal-to-noise ratio teleseismic events through time-frequency-wave number domain analysis method for teleseismic event data segments, and eliminating data with a signal-to-noise ratio lower than a preset threshold to obtain high-quality teleseismic data segments.

[0075] a super-deep reflection wave field reconstruction module for using a pure P-wave multi-dimensional deconvolution method to reconstruct and deconvolve source wavelets for high-quality teleseismic data segments to obtain super-deep reflection wave field data.

[0076] a deep earth profile reflection imaging module for performing reflection imaging processing on the super-deep reflection wave field data, and combining velocity analysis and stacking imaging to optimize the resolution capability of the deep reflection interface to obtain a deep earth reflection imaging profile; in stacking imaging, for regular observation systems, common midpoint domain stacking is selected, and for non-regular observation systems, common shot point domain stacking or common receiver point domain stacking is selected.

[0077] Of course, Figure 8 the architecture shown in FIG. 1 is only exemplary, and when implementing different functions, one or at least two components in the system shown in FIG. 1 can be omitted Figure 8 according to actual needs.

[0078] In one exemplary embodiment, a computer device is provided, which can be a server or a terminal, and an internal structure diagram thereof can be as shown in Figure 9As shown in the figure. The computer device includes a processor, a memory, an input / output interface (I / O for short) and a communication interface. Among them, the processor, the memory and the input / output interface are connected through the system bus, and the communication interface is connected to the system bus through the input / output interface. Among them, the processor of the computer device is used to provide computing and control capability. The memory of the computer device includes a non-volatile storage medium and an internal memory. The non-volatile storage medium stores an operating system, a computer program and a database. The internal memory provides an environment for the operation of the operating system and the computer program in the non-volatile storage medium. The input / output interface of the computer device is used to exchange information between the processor and external devices. The communication interface of the computer device is used to communicate with the terminal outside through network connection. The computer program is executed by the processor, and the passive source seismic data based on the super-deep reflection imaging method provided in the above embodiment can be realized.

[0079] Those skilled in the art can understand that, Figure 9 The structure shown in the figure is only a block diagram of part of the structure related to the scheme of the present application, and does not constitute a limitation on the computer device to which the scheme of the present application is applied. The specific computer device can include more or fewer components than those shown in the figure, or combine certain components, or have a different component arrangement.

[0080] In an exemplary embodiment, a computer device is also provided, including a memory and a processor, the memory storing a computer program, and the processor executing the computer program to realize the steps in the above method embodiments.

[0081] In an exemplary embodiment, a computer readable storage medium is provided, storing a computer program, which is executed by a processor to realize the steps in the above method embodiments.

[0082] In an exemplary embodiment, a computer program product is provided, including a computer program, which is executed by a processor to realize the steps in the above method embodiments.

[0083] It should be noted that the user information (including but not limited to user equipment information, user personal information, etc.) and data (including but not limited to data for analysis, stored data, displayed data, etc.) involved in the present application are all information and data authorized by the user or authorized by all parties, and the collection, use and processing of related data need to comply with relevant regulations.

[0084] Those skilled in the art can understand that all or part of the processes in the above-mentioned embodiment methods can be completed by instructing the relevant hardware through a computer program. The computer program can be stored in a non-volatile computer readable storage medium, and when the computer program is executed, the processes of the above-mentioned embodiments of the methods can be included. Any reference to memory, database or other medium used in the embodiments provided in the present application can include at least one of non-volatile and volatile memory. Non-volatile memory can include read-only memory (ROM), magnetic tape, floppy disk, flash memory, optical storage, high-density embedded non-volatile memory, resistive memory (ReRAM), magnetoresistive random access memory (MRAM), ferroelectric memory (FRAM), phase change memory (PCM), graphene memory, etc. Volatile memory can include random access memory (RAM) or external cache memory, etc. As an illustration but not limitation, RAM can be in various forms, such as static random access memory (SRAM) or dynamic random access memory (DRAM), etc.

[0085] The database involved in the embodiments provided in the present application can include at least one of a relational database and a non-relational database. The non-relational database can include a distributed database based on a blockchain, etc., without being limited thereto. The processor involved in the embodiments provided in the present application can be a general-purpose processor, a central processing unit, a graphics processing unit, a digital signal processor, a programmable logic device, a data processing logic device based on quantum computing, etc., without being limited thereto.

[0086] The technical features of the above embodiments can be combined arbitrarily. In order to make the description concise, all possible combinations of the technical features in the above embodiments are not described, but as long as the combinations of the technical features do not exist, they should be considered as the scope of the present application.

[0087] The principles and implementation modes of the present application are described by applying specific examples herein. The above description of the embodiments is only used to help understand the method and its core idea of the present application; meanwhile, for those skilled in the art, according to the idea of the present application, the specific implementation mode and application range can be changed. In conclusion, the content of the present application should not be understood as a limitation.

Claims

1. A method of ultra-deep reflection imaging based on passive source seismic data, characterized in that, The method comprises the following steps: obtaining original continuous seismic waveform data, cutting and preprocessing seismic waveform data with P-wave first arrival and S-wave non-arrival to obtain teleseismic event data segment; the preprocessing comprises de-instrument response, normalization and band-pass filtering; for the teleseismic event data segment, high signal-to-noise ratio teleseismic events are identified by time-frequency-wave number domain analysis method, and data with signal-to-noise ratio lower than a preset threshold is removed to obtain high-quality teleseismic data segment; a pure P-wave multidimensional deconvolution method is used to reconstruct the reflection wave field of the high-quality teleseismic data segment and deconvolve the source wavelet to obtain super-deep reflection wave field data; reflection imaging processing is performed on the super-deep reflection wave field data, and the resolving power of deep reflection interface is optimized by combining velocity analysis and stacking imaging to obtain deep reflection imaging profile of the earth; in the stacking imaging, common midpoint domain stacking is selected for regular observation system, and common shot point domain stacking or common receiver point domain stacking is selected for irregular observation system; high signal-to-noise ratio teleseismic events are identified according to the following formula: wherein TXFK represents a high signal-to-noise ratio teleseismic event identified by a time-frequency-wave number domain analysis method, X(t,x) is seismic waveform data in a time-space domain, t is time, x is a geophone position, f TX (X) is a feature extraction function in a time-space domain, capturing the features of seismic waveform data X, F(k,ω) is amplitude distribution information extracted in a frequency-wave number domain, k is a wave number, ω is a frequency, and respectively represent a preset screening criterion for high signal-to-noise ratio seismic waveform data, indicates that the seismic waveform data contains a high signal-to-noise ratio horizontal body wave signal, indicates that the main energy of the seismic waveform data is distributed in a body wave region of a frequency-wave number spectrum; the time-frequency-wave number domain analysis method is to identify a high signal-to-noise ratio body wave signal by analyzing the features of a time-distance domain and a frequency-wave number domain.

2. The passive source seismic data based ultra-deep reflectivity imaging method of claim 1, wherein, obtaining original continuous seismic waveform data, cutting and preprocessing seismic waveform data with P-wave first arrival and S-wave non-arrival to obtain teleseismic event data segment, specifically comprising: obtaining original continuous seismic waveform data, and extracting seismic event data segment from the original continuous seismic waveform data according to the public earthquake occurrence time; cutting and preprocessing seismic waveform data with P-wave first arrival and S-wave non-arrival according to the P-wave first arrival time of the seismic event data segment to obtain teleseismic event data segment.

3. The passive source seismic data based ultra-deep reflection imaging method of claim 1, wherein, reflection wave field reconstruction is performed according to the following formula: in, For effective reflection range, x represents R The virtual traction source at x A The z-component of the elastic dynamic reflection response of the P-wave induced at the receiver, V(x) R ,x B ,ω) represents x B The virtual traction source at x R The z-component of the complete wave field caused by the receiver at that location, d 2 x R V is the boundary area element. z (x A ,x B ,ω) represents x B The virtual traction source at x A The z-component of the complete wave field caused by the receiver at that location. x represents A The z-component of the P-wave field caused by the receiver at the location, where ω is the angular frequency.

4. The passive source seismic data based ultra-deep reflection imaging method of claim 1, wherein, the frequency of the source wavelet is determined according to the maximum station spacing; when the maximum station spacing is 1000km, the source wavelet selects Ricker wavelet with 0-0.4Hz as the main frequency.

5. The passive source seismic data based ultra-deep reflection imaging method of claim 1, wherein, The reflection imaging processing performed on the super-deep reflection wave field data includes denoising, static correction and dynamic correction.

6. An ultra-deep reflection imaging system based on passive source seismic data, characterized in that, The method comprises the following steps: a seismic data acquisition and preprocessing module is used to obtain original continuous seismic waveform data, cut and preprocess seismic waveform data with P-wave first arrival and S-wave non-arrival to obtain teleseismic event data segment; the preprocessing comprises de-instrument response, normalization and band-pass filtering; a high signal-to-noise ratio teleseismic data screening module is used to identify high signal-to-noise ratio teleseismic events by time-frequency-wave number domain analysis method for the teleseismic event data segment, and remove data with signal-to-noise ratio lower than a preset threshold to obtain high-quality teleseismic data segment; a super-deep reflection wave field reconstruction module is used to reconstruct the reflection wave field of the high-quality teleseismic data segment by a pure P-wave multidimensional deconvolution method and deconvolve the source wavelet to obtain super-deep reflection wave field data; a deep profile reflection imaging module is used to perform reflection imaging processing on the super-deep reflection wave field data, and optimize the resolving power of deep reflection interface by combining velocity analysis and stacking imaging to obtain deep reflection imaging profile of the earth; in the stacking imaging, common midpoint domain stacking is selected for regular observation system, and common shot point domain stacking or common receiver point domain stacking is selected for irregular observation system; in the high signal-to-noise ratio teleseismic data screening module, high signal-to-noise ratio teleseismic events are identified according to the following formula: wherein TXFK represents a high signal-to-noise ratio teleseismic event identified by a time-frequency-wave number domain analysis method, X(t,x) is seismic waveform data in a time-space domain, t is time, x is a geophone position, f TX (X) is a feature extraction function in a time-space domain, capturing the features of seismic waveform data X, F(k,ω) is amplitude distribution information extracted in a frequency-wave number domain, k is a wave number, ω is a frequency, and respectively are preset screening criteria for high signal-to-noise ratio seismic waveform data, indicates that the seismic waveform data contains a horizontal body wave signal with high signal-to-noise ratio, indicates that the main energy of the seismic waveform data is distributed in a body wave region of a frequency-wave number spectrum; the time-frequency-wave number domain analysis method is to identify a high signal-to-noise ratio body wave signal by analyzing the features in a time-distance domain and a frequency-wave number domain.

7. A computer device comprising: A memory, a processor, and a computer program stored on the memory and executable on the processor, characterized in that the processor executes the computer program to implement the passive source seismic data based ultra-deep reflection imaging method of any one of claims 1-5.

8. A computer-readable storage medium having stored thereon a computer program, characterized in that, The computer program is executed by the processor to implement the passive source seismic data based ultra-deep reflection imaging method of any one of claims 1-5.

9. A computer program product comprising a computer program, characterized in that, The computer program is executed by the processor to implement the passive source seismic data based ultra-deep reflection imaging method of any one of claims 1-5.

Images