A microseismic interference positioning method, system, device and medium based on time-shifted correlation coefficient

Abstract

The present invention discloses a microseismic interference positioning method, system, device and medium based on time-shifted correlation coefficients, including: constructing DTW time-shifted correlation coefficients by using the dynamic time warping algorithm, multiplying the constructed time-shifted correlation coefficients by the imaging operator, and superimposing all time-shifted correlation gathers at the same grid point to obtain an imaging profile of a single grid point. Further, superimposing all imaging profiles of single grid points to obtain the final imaging effect characterizing the source position, and the position with the maximum imaging value is the source, thereby realizing source positioning. The present invention can solve the phenomenon of reduced imaging resolution of the cross-correlation imaging operator caused by wave propagation effects and noise, etc., and improve the positioning effect of the method under the influence of noise pollution.

Description

Technical Field

[0001] The present invention belongs to the technical field of seismic exploration, and particularly relates to a microseismic interference positioning method, system, device and medium based on time-lapse correlation coefficients. Background Art

[0002] During the process of oil and gas exploitation, hydraulic fracturing technology is usually adopted to increase the production of oil and gas wells. The fracturing operation causes the underground rocks to break, and thus radiate energy outward, causing the rocks to move. These movements are similar to natural earthquakes along the fault plane and have a small magnitude, so they are called "microseisms". The occurrence of a large number of microseismic events can provide sufficient seismic data for seismic processing and interpretation.

[0003] Microseismic monitoring technology uses microseismic events induced by the rupture of underground rock layers for reservoir monitoring and underground fracture description. Different from traditional seismic exploration, the location, intensity and occurrence time of the seismic source are unknown in microseismic monitoring. Therefore, seismic source location is the basis of microseismic monitoring. Traditional microseismic location methods use the picked arrival times and a given velocity model to perform ray-based linear travel-time inversion to determine the seismic source location. However, due to the small magnitude and low signal-to-noise ratio of microseismic events, it is difficult to accurately pick the arrival times and locate the seismic source. Compared with traditional travel-time inversion methods, waveform-based location methods do not require phase picking operations and can detect and locate more events using low signal-to-noise ratio data.

[0004] Waveform-based location methods can be roughly divided into reverse time migration imaging methods and waveform stacking methods that mainly use the first arrival phases. Reverse time migration imaging methods require wavefield extrapolation, which is time-consuming and the imaging conditions used have a great influence on the seismic source location results. Waveform stacking methods initially mainly used diffraction stacking operators to focus or back-project the energy along the travel-time curve onto spatial grid points. In recent years, such methods have been successfully used to locate natural earthquakes and induced seismic activities related to mining operations, geothermal exploitation and oil and gas reservoirs. However, when using the diffraction stacking method for imaging, it is necessary to search for the occurrence time of the seismic source, and the error in the occurrence time will also affect the final location accuracy.

[0005] The cross - correlation superposition source location method belonging to the interference imaging method is a relatively new waveform superposition location method. This method uses the travel - time difference information extracted by the cross - correlation operator for superposition imaging, which can avoid searching for the earthquake - occurrence time. Moreover, by weighting different travel - time difference information according to the characteristics of different seismic phases, the location accuracy can be improved. Under the influence of factors such as the source mechanism, medium change, and environmental noise, different time differences will be generated at different times of the micro - seismic waveforms at different geophones. The overall time shift calculated by cross - correlation cannot achieve good alignment of the changing waveforms of micro - seismic events, thereby affecting the location accuracy. Currently, most methods use characteristic functions to transform the original waveforms. Characteristic functions can effectively avoid the influence of the change of the first - arrival waveform and reduce the sensitivity to velocity errors, but they will sacrifice the spatial resolution of imaging and are easily affected by noise. Therefore, it is very necessary to find a method that can achieve precise alignment between changing waveforms and use it to achieve accurate micro - seismic location. Summary of the Invention

[0006] The purpose of the present invention is to provide a micro - seismic interference location method, system, device, and medium based on time - shift correlation coefficients to solve the problems existing in the above - mentioned prior art.

[0007] To achieve the above purpose, the present invention provides a micro - seismic interference location method based on time - shift correlation coefficients, including:

[0008] Obtain the velocity model and micro - seismic records of the monitoring area;

[0009] Based on a preset grid size, divide the velocity model into several grid points, calculate the travel times for each grid point, and obtain a travel - time table containing theoretical travel - time data;

[0010] Perform dynamic time warping calculation on every two seismic records in the micro - seismic records to obtain the time shift of corresponding sampling points between every two micro - seismic records, pick up the first - arrival time of the seismic event, and based on the time shift at the first - arrival moment, construct a time - shift sequence of each sampling point relative to the first - arrival moment. Calculate the correlation coefficient between seismic traces based on the time - shift sequence to obtain a time - shift correlation trace set containing actual observed travel - time difference information;

[0011] Obtain the imaging profile corresponding to each grid point based on the time - shift correlation trace set and the travel - time table;

[0012] Construct a final imaging map representing the source location based on the imaging profiles corresponding to each grid point, and determine the source location based on the final imaging map.

[0013] Preferably, the calculation of the travel times for each grid point specifically includes:

[0014] The theoretical travel time data from each grid point to each detector is calculated, and a travel time table is constructed based on the calculated theoretical travel time data.

[0015] Preferably, the time shift sequence of each sampling point relative to the first arrival time is constructed, and the specific calculation formula is:

[0016] τ TC (t i )=τ t (t i )-τ t (t fb ),t i =t first ,t max

[0017] In the formula, i is the corresponding sampling point of the sequence, t i is the time corresponding to the sampling point, t first to t max is the time range corresponding to the sampling point, τ t (t i ) is the calculated time shift corresponding to each sampling point, τ t (t fb ) is the first arrival time of the microseismic signal t fb The corresponding time shift, τ TC The calculated time t for each sampling point i The relative time shift.

[0018] Preferably, the correlation coefficient between seismic traces is calculated based on the time-shift sequence to obtain a time-shift correlation gather containing actual observed travel time difference information. The specific calculation formula is:

[0019]

[0020] Where i and j are detector channel numbers, τ is the time shift between the corresponding sampling points of the two channels, and τ TC is the relative time shift, u(t i ,i) is the time series corresponding to detector i, u(t i +τ+τ TC (t i ), j) is the relative time shift sequence calculated by detector j relative to detector i, t i is the time corresponding to the sampling point, t first to t max is the time range corresponding to the sampling point, C TC is the correlation coefficient between the two traces based on DTW.

[0021] Preferably, obtaining the imaging profile corresponding to each grid point based on the time-shift correlation gather and the travel time table specifically includes:

[0022] Multiply the time-shifted associated gather by the interference imaging operator containing the theoretical travel time data to perform interference imaging, and stack all the time-shifted associated gathers at the same grid point to obtain the imaging profile of a single grid point.

[0023] Preferably, determining the source position based on the final imaging map specifically includes:

[0024] The source is obtained by determining the position corresponding to the maximum imaging value in the final imaging map, thereby realizing source positioning.

[0025] A microseismic interference positioning system based on time-shifted correlation coefficients, comprising:

[0026] A data acquisition module, configured to acquire the velocity model and microseismic records of the monitoring area;

[0027] A travel time table calculation module, configured to divide the velocity model into a plurality of grid points according to a preset grid size, perform travel time calculation on each of the grid points, and obtain a travel time table containing theoretical travel time data;

[0028] A time-shifted correlation coefficient calculation module, configured to perform dynamic time warping calculation on every two seismic records in the microseismic records to obtain the time shift amount of corresponding sampling points between every two microseismic records, pick up the arrival time of the seismic event, and construct a time shift sequence of each sampling point relative to the arrival time based on the time shift at the arrival time, and calculate the correlation coefficient between seismic traces based on the time shift sequence to obtain a time-shifted associated gather containing actual observed travel time difference information;

[0029] A positioning module, configured to obtain the imaging profile corresponding to each grid point according to the time-shifted associated gather and the travel time table; construct a final imaging map representing the source position based on the imaging profiles corresponding to each grid point, and determine the source position based on the final imaging map.

[0030] An electronic device, comprising a memory and a processor, where the memory is used to store a computer program, and the processor runs the computer program to enable the electronic device to execute the described microseismic interference positioning method based on time-shifted correlation coefficients.

[0031] A computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, it implements the described microseismic interference positioning method based on time-shifted correlation coefficients.

[0032] The technical effects of the present invention are:

[0033] The present invention uses the Dynamic Time Warping (DTW) algorithm to construct the DTW time-shifted correlation coefficient. The constructed time-shifted correlation coefficient is multiplied by the imaging operator, and all time-shifted correlation gathers are stacked at the same grid point to obtain the imaging profile of a single grid point. Further, all the imaging profiles of single grid points are stacked to obtain the final imaging effect representing the source position. The position with the maximum imaging value is the source, thus realizing source positioning. By using the DTW method to construct the time-shifted correlation coefficient and multiplying and stacking the time-shifted correlation coefficient with the imaging operator, the present invention realizes the imaging of the source position. The imaging effect is more focused and better than the result obtained by the traditional cross-correlation-based interference imaging method. The present invention can solve the phenomenon of reduced imaging resolution of the cross-correlation imaging operator caused by wave propagation effects and noise, etc., and improve the positioning effect of the method under the influence of noise pollution. BRIEF DESCRIPTION OF THE DRAWINGS

[0034] In order to more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the following will briefly introduce the drawings required for use in the embodiments. Obviously, the drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, other drawings can be obtained based on these drawings without creative efforts.

[0035] The drawings forming a part of this application are used to provide a further understanding of this application. The schematic embodiments of this application and their descriptions are used to explain this application and do not constitute an improper limitation of this application. In the drawings:

[0036] Figure 1 is the schematic diagram of the implementation process in the embodiment of the present invention;

[0037] Figure 2 is the schematic diagram of the constructed time-shifted correlation coefficient and the traditional cross-correlation coefficient proposed in the embodiment of the present invention;

[0038] Figure 3 is the schematic diagram of the forward velocity model used in the embodiment of the present invention;

[0039] Figure 4 is the schematic diagram of the noisy seismic data generated by adding -12 dB Gaussian noise to the original seismic data in the embodiment of the present invention;

[0040] Figure 5 is the filtered seismic data obtained by performing 5 - 80 Hz band-pass filtering on the data with -12 dB Gaussian noise in the embodiment of the present invention;

[0041] Figure 6 is the schematic diagram of the positioning results obtained by using the traditional cross-correlation interference positioning and the interference positioning method based on the time-shifted correlation coefficient proposed by the present invention respectively using the filtered data in the embodiment of the present invention. Detailed Embodiments

[0042] The various exemplary embodiments of the present invention will now be described in detail. This detailed description should not be construed as a limitation on the present invention, but rather as a more detailed description of certain aspects, features, and implementation schemes of the present invention.

[0043] It should be understood that the terms used in the present invention are only for describing specific embodiments and are not intended to limit the present invention. Additionally, for the numerical ranges in the present invention, it should be understood that each intermediate value between the upper and lower limits of the range is also specifically disclosed. Each intermediate value within any stated value or stated range, as well as each smaller range between any other stated value or intermediate value within the stated range, is also included in the present invention. The upper and lower limits of these smaller ranges may be independently included or excluded from the range.

[0044] Without departing from the scope or spirit of the present invention, various improvements and changes can be made to the specific embodiments of the present invention's specification, which will be obvious to those skilled in the art. Other embodiments obtained from the specification of the present invention will be obvious to those skilled in the art. The specification and embodiments of this application are merely exemplary.

[0045] Regarding the use of "comprising", "including", "having", "containing", etc. in this article, they are all open-ended terms, meaning including but not limited to.

[0046] It should be noted that, without conflict, the embodiments in this application and the features in the embodiments can be combined with each other. The following will refer to the accompanying drawings and combine with embodiments to detail this application.

[0047] Embodiment 1

[0048] As Figure 1 - Figure 6 shown, in this embodiment, a microseismic interference positioning method based on time-shift correlation coefficient is provided, including: obtaining a velocity model and microseismic records of a monitoring area; dividing the velocity model into several grid points based on a preset grid size, calculating the travel times for each of the grid points to obtain a travel time table containing theoretical travel time data; performing dynamic time warping calculation on every two seismic records in the microseismic records to obtain the time shift amount of corresponding sampling points between every two microseismic records, picking up the first arrival time of the seismic event, and constructing a time shift sequence of each sampling point relative to the first arrival time based on the time shift at the first arrival moment, calculating the correlation coefficient between seismic traces based on the time shift sequence to obtain a time-shift correlation trace gather containing actual observed travel time difference information; obtaining an imaging profile corresponding to each grid point based on the time-shift correlation trace gather and the travel time table; constructing a final imaging map representing the source location based on the imaging profiles corresponding to each grid point, and determining the source location based on the final imaging map.

[0049] In this embodiment, the dynamic time warping (DTW) algorithm is used to construct the DTW time-shift correlation coefficient. The constructed time-shift correlation coefficient is multiplied by the imaging operator, and all time-shifted correlation gathers are stacked at the same grid point to obtain the imaging profile of a single grid point. Further, all the imaging profiles of single grid points are stacked to obtain the final imaging result representing the source position. The position with the maximum imaging value is the source, thus realizing source positioning. In this embodiment, by using the DTW method to construct the time-shift correlation coefficient and multiplying and stacking the time-shift correlation coefficient with the imaging operator, the imaging of the source position is realized. The imaging effect is more focused and better than the result obtained by the traditional cross-correlation-based interferometric imaging method. This embodiment can solve the phenomenon of reduced imaging resolution of the cross-correlation imaging operator caused by factors such as wave propagation effects and noise, and improve the positioning effect of the method under the influence of noise pollution.

[0050] The specific technical solution adopted in this embodiment is as follows:

[0051] Step 1: Calculate the travel time table. The velocity model of the monitoring area is divided into several grid points according to a certain grid size. Each grid point is regarded as a potential source, and the theoretical travel time from each grid point to each geophone is calculated to obtain the corresponding travel time table. Optionally, in order to calculate the travel time efficiently, this embodiment uses the open-source Python package fteikpy to perform fast solution of the Eikonal equation. Compared with the prior art, the use of fteikpy significantly reduces the computational complexity while ensuring the accuracy of the solution;

[0052] Step 2: Calculate the time-shift correlation coefficient. Perform DTW calculation between every two seismic traces in the seismic gather. The calculation result of the DTW algorithm is the time shift amount of the corresponding sampling points between every two seismic records, from which the time shift τ of each sampling point of the waveform corresponding to the microseismic event can be obtained. t Pick up the arrival time t of the corresponding seismic event first . Take the time shift corresponding to the picked arrival time t first as τ fb , and construct the relative time shift sequence of this waveform with respect to τ fb . Optionally, in order to pick up the arrival time of the corresponding seismic event, this embodiment uses the EventPicker Python package. Through the EventPicker package, the arrival time t can be accurately picked up. fb . On the basis of picking up the arrival time, take the time shift corresponding to the arrival time as the reference to construct the relative time shift sequence of this waveform with respect to time as follows:

[0053] τ TC (t i ) = τ t (ti ) - τ t (t fb ), t i = t first , t max

[0054] Wherein, i is the sampling point corresponding to the sequence, t i is the time corresponding to the sampling point, t first to t max is its belonging range, τ t (t i ) is the time shift amount corresponding to each sampling point time obtained by calculation, τ t (t fb ) is the initial arrival time t fb of the microseismic signal corresponding to the time shift amount, τ TC is the relative time shift amount of each sampling point time t i calculated.

[0055] Calculate the correlation coefficient between the seismic trace after relative time shift and another seismic trace to obtain a time shift associated trace gather containing the actual observed travel time difference information:

[0056]

[0057] Wherein, i and j are the detector trace numbers, τ is the time shift amount of the corresponding sampling points between the two traces, τ TC is the relative time shift amount, u(t i , i) is the time series corresponding to detector i, u(t i + τ + τ TC (t i ), j) is the relative time shift series calculated for detector j relative to detector i, t i is the time corresponding to the sampling point, t first to t max is the range to which the time corresponding to the sampling point belongs, C TC is the correlation coefficient based on DTW between the two traces.

[0058] Step 3: Obtain a single trace gather imaging profile. Multiply the time shift associated trace gather by the interference imaging operator containing the theoretical travel time difference information, perform interference imaging, and stack all the time shift associated trace gathers at the same grid point to obtain a single grid point imaging profile:

[0059]

[0060] Wherein, x is the source position vector, δ is the Dirac Delta function, τ is the time, N is the number of detectors, τ i,x , τ j,x are the theoretical travel times from the source x to detectors i and j, S TC(x, i) is the imaging profile of a single grid point. The interference imaging operator utilizes the travel time difference of the source-detector pair {i, j}.

[0061] Step 4: Superimpose all the imaging profiles of single grid points to obtain the final imaging result characterizing the source position. The position with the maximum imaging value is the source, thus achieving source localization.

[0062]

[0063] In the formula, N is the number of detectors, and S TC (x) is the final imaging result.

[0064] Compared with the traditional interference imaging method based on the cross-correlation coefficient, the time-shifted correlation coefficient is constructed in Equation 2 of this embodiment. The comparison between the time-shifted correlation coefficient and the cross-correlation coefficient is as Figure 2 shown. Figure 6 This is the comparison between the microseismic interference imaging localization effect based on the time-shifted correlation coefficient and the traditional interference imaging localization effect based on the cross-correlation coefficient in this embodiment. The imaging effect based on the time-shifted correlation coefficient has been greatly improved.

[0065] A microseismic interference positioning system based on the time-shifted correlation coefficient, comprising:

[0066] A data acquisition module, configured to obtain the velocity model and microseismic records of the monitoring area;

[0067] A travel time table calculation module, configured to divide the velocity model into several grid points according to a preset grid size, perform travel time calculation on each of the grid points, and obtain a travel time table containing theoretical travel time data;

[0068] A time-shifted correlation coefficient calculation module, configured to perform dynamic time warping calculation on every two seismic records in the microseismic records to obtain the time shift amount of corresponding sampling points between every two microseismic records, pick up the arrival time of the seismic event, construct a time shift sequence of each sampling point relative to the arrival time, and calculate the correlation coefficient between seismic traces based on the time shift sequence to obtain a time-shifted correlation trace gather containing actual observed travel time difference information;

[0069] A positioning module, configured to obtain the imaging profile corresponding to each grid point according to the time-shifted correlation trace gather and the travel time table; construct a final imaging map characterizing the source position based on the imaging profiles corresponding to each of the grid points, and determine the source position based on the final imaging map.

[0070] An electronic device, comprising a memory and a processor. The memory is used to store a computer program, and the processor runs the computer program to enable the electronic device to execute the described microseismic interference positioning method based on the time-shifted correlation coefficient.

[0071] A computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the described microseismic interference positioning method based on time-shifted correlation coefficients is implemented.

[0072] As described above, the above is only a preferred specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any changes or substitutions that can be easily thought of by those skilled in the technical field of the present application within the technical scope disclosed by the present application should be covered within the protection scope of the present application. Therefore, the protection scope of the present application should be subject to the protection scope of the claims.

Claims

1. A microseismic interference positioning method based on time-shifted correlation coefficient, characterized in that, Including: Obtain the velocity model and microseismic records of the monitoring area; Divide the velocity model into a number of grid points based on a preset grid size, perform travel time calculations on each of the grid points, and obtain a travel time table containing theoretical travel time data; Perform dynamic time warping calculations on every two seismic records in the microseismic records to obtain the time shift amount of corresponding sampling points between every two microseismic records, pick up the first arrival time of the seismic event, and construct a time shift sequence of each sampling point relative to the first arrival time based on the time shift at the first arrival moment. Calculate the correlation coefficient between seismic traces based on the time shift sequence to obtain a time shift correlation trace set containing actual observed travel time difference information; The specific calculation formula for calculating the correlation coefficient between seismic traces based on the time shift sequence to obtain a time shift correlation trace set containing actual observed travel time difference information is: Where i and j are detector channel numbers, τ is the time shift between the corresponding sampling points of the two channels, and τ TC is the relative time shift, u(t i ,i) is the time series corresponding to detector i, u(t i +τ+τ TC (t i ), j) is the relative time shift sequence calculated by detector j relative to detector i, t i is the time corresponding to the sampling point, t first to t max is the time range corresponding to the sampling point, C TC is the correlation coefficient between the two traces based on DTW; Obtain the imaging profile corresponding to each grid point based on the time shift correlation trace set and the travel time table; Construct a final imaging map representing the source location based on the imaging profiles corresponding to each grid point, and determine the source location based on the final imaging map.

2. The microseismic interference positioning method based on time-shifted correlation coefficient according to claim 1, wherein The performing travel time calculations on each of the grid points specifically includes: Calculate the theoretical travel time data of each grid point to each geophone, and construct a travel time table based on the calculated theoretical travel time data.

3. The microseismic interference positioning method based on time-shifted correlation coefficient according to claim 1, characterized in that The specific calculation formula for constructing a time shift sequence of each sampling point relative to the first arrival time is: τ TC (t i ) = τ t (t i ) - τ t (t fb ), t i = [t first ,t max ​ where \(i\) is the sampling point corresponding to the sequence, \(t\) i is the time corresponding to the sampling point, \(t\) first to \(t\) max is the time range to which the sampling point belongs, \(\tau\) t (\(t\) i ) is the time shift amount calculated corresponding to the time of each sampling point, \(\tau\) t (\(t\) fb ) is the time shift amount corresponding to the arrival time \(t\) fb of the microseismic signal, \(\tau\) TC is the relative time shift amount of the time \(t\) i of each sampling point calculated.

4. The microseismic interference positioning method based on time-shifted correlation coefficient according to claim 1, wherein, The obtaining the imaging profile corresponding to each grid point based on the time shift correlation trace set and the travel time table specifically includes: Multiply the time shift correlation trace set by an interference imaging operator containing theoretical travel time data, perform interference imaging, and stack all time shift correlation trace sets at the same grid point to obtain an imaging profile of a single grid point.

5. A microseismic interference positioning method based on time-shifted correlation coefficient according to claim 1, characterized in that The determining the source location based on the final imaging map specifically includes: Locate the source by determining the position corresponding to the maximum imaging value in the final imaging map to achieve source location.

6. A microseismic interference positioning system based on time-shifted correlation coefficient, characterized in that, Including: A data acquisition module for obtaining the velocity model and microseismic records of the monitoring area; A travel time table calculation module for dividing the velocity model into a number of grid points according to a preset grid size, performing travel time calculations on each of the grid points, and obtaining a travel time table containing theoretical travel time data; A time shift correlation coefficient calculation module for performing dynamic time warping calculations on every two seismic records in the microseismic records to obtain the time shift amount of corresponding sampling points between every two microseismic records, picking up the first arrival time of the seismic event, and constructing a time shift sequence of each sampling point relative to the first arrival time based on the time shift at the first arrival moment. Calculate the correlation coefficient between seismic traces based on the time shift sequence to obtain a time shift correlation trace set containing actual observed travel time difference information; The specific calculation formula for calculating the correlation coefficient between seismic traces based on the time shift sequence to obtain a time shift correlation trace set containing actual observed travel time difference information is: where i and j are the detector channel numbers, τ is the time shift of the corresponding sampling points between the two channels, and τ TC is the relative time shift, u(t i , i) is the time series corresponding to detector i, u(t i +τ + τ TC (t i ), j) is the relative time shift sequence calculated for detector j relative to detector i, t i is the time corresponding to the sampling point, and t first to t max is the range to which the time corresponding to the sampling point belongs. C TC is the correlation coefficient based on DTW between the two channels; A positioning module for obtaining the imaging profile corresponding to each grid point according to the time shift correlation trace set and the travel time table; constructing a final imaging map representing the source location based on the imaging profiles corresponding to each grid point, and determining the source location based on the final imaging map.

7. An electronic device, characterized in that, It includes a memory and a processor. The memory is used to store a computer program, and the processor runs the computer program to enable the electronic device to execute a microseismic interference positioning method based on time-shifted correlation coefficients according to any one of claims 1-5.

8. A computer-readable storage medium, characterized in that, It stores a computer program, and when the computer program is executed by a processor, it implements a microseismic interference positioning method based on time-shifted correlation coefficients according to any one of claims 1-5.

Images