A real-time wave migration positioning method for microseismic signals

Abstract

The application provides a real-time wave fluctuation migration positioning method for microseismic signals and belongs to the technical field of seismic signal processing; solves the problems of complex data preprocessing, insufficient calculation efficiency and low timeliness of existing migration positioning processing technologies; and comprises the following steps: arranging multiple geophones as receiving points of microseismic data according to different depths in the range of a region to be monitored on a coal mine working face floor; counting the number N of triggered geophones, and distributing the data of the N geophones to different GPUs to carry out seismic wave field simulation of the N geophones in parallel; obtaining imaging results of longitudinal waves and transverse waves by using directional aperture constraints and storage-free imaging conditions, fusing the imaging results of the longitudinal waves and the transverse waves, and extracting imaging results with the most focused energy and consistent imaging point positions as microseismic information of microseismic events; and the application is applied to coal mine roof and floor monitoring.

Description

Technical Field

[0001] This invention provides a real-time fluctuation offset positioning method for microseismic signals, belonging to the field of microseismic signal processing technology. Background Technology

[0002] As coal mining in my country delves deeper, highly confined limestone water in the mine floor becomes a key target for mine water hazard prevention. Mining activities create fracture zones in the working face floor rock mass, which can easily communicate with the highly confined limestone water beneath, triggering water inrush accidents. Therefore, effectively assessing the integrity of the working face floor is a pressing issue. Rock fractures generate microseismic signals, which can be used to effectively locate fracture sites and calculate fracture strength. Therefore, microseismic monitoring technology can serve as a crucial means of evaluating the integrity of the working face floor rock mass. Locating fracture areas through seismic sources and analyzing the evolution of fracture zones can provide key criteria for water inrush early warning. Existing microseismic location technologies in coal mines, such as ray-based methods like time-difference and double-difference methods, suffer from insufficient accuracy because they ignore the propagation patterns of seismic waves in complex media. However, migration location technology based on wave theory can simulate the propagation patterns of seismic waves in complex media, and thus locate microseismic events through specific imaging conditions, possessing significant application potential.

[0003] However, current offset positioning processing techniques have the following shortcomings:

[0004] 1. Faced with complex data preprocessing schemes, such as when using a single type of P-wave positioning, it is necessary to use a time window to extract the first arrival waveform, while when choosing P-wave and S-wave positioning, it is inevitable to carry out P-wave and S-wave separation.

[0005] 2. The computational efficiency in offset positioning is insufficient, especially in actual 3D scenes, making it difficult to meet production requirements;

[0006] 3. The offset positioning conditions in the existing imaging conditions involve integration over time, which requires storing the three-dimensional spatial wave field for the entire simulation moment, resulting in high resource consumption.

[0007] 4. In sparse observation, offset positioning suffers from multiple solutions.

[0008] 5. The lack of real-time, coordinated processing solutions for multiple micro-seismic events makes it difficult to meet the high-efficiency calculation requirements of field applications.

[0009] Therefore, how to provide a method that can maximize computational efficiency and effectively apply microseismic wave offset positioning technology to field data without complex data preprocessing, so as to achieve real-time high-precision microseismic positioning, is a technical problem that urgently needs to be solved by those skilled in the art. Summary of the Invention

[0010] To address the problems of complex data preprocessing, insufficient computational efficiency, and low timeliness in existing offset positioning technologies, this invention proposes a real-time fluctuation offset positioning method for microseismic signals. This method features an adaptive P-wave and S-wave decoupling scheme, highly efficient three-dimensional wave equation simulation calculation, real-time linkage processing of microseismic events from the field, and provides the final positioning result.

[0011] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is: a real-time fluctuation offset positioning method for microseismic signals, comprising the following steps:

[0012] S1: Within the area that needs to be monitored on the floor of the coal mine working face, multiple geophones are arranged at different depths as receiving points for microseismic data. All geophones are connected to the acquisition host via the main line and acquire microseismic data after microseismic events occur in the coal seam working face and its roof.

[0013] S2: Count the number of triggered detectors N, and distribute the data of these N detectors to different GPUs. Perform seismic wave field simulation of N detectors in parallel and obtain wave field information at the imaging time in three-dimensional space. The seismic wave field simulation uses the three-dimensional elastic wave equation with P-wave and S-wave decoupling to adaptively decouple in the wave field solution. The differential solution of the wave equation is processed in parallel by GPU.

[0014] S3: The imaging results of P-waves and S-waves are obtained by using directional aperture constraint and non-storage migration imaging method. The imaging results of P-waves and S-waves are fused and the imaging results with the most focused energy and consistent imaging point positions are extracted as microseismic information of microseismic events. The energy information at the microseismic location is used to estimate the microseismic energy. Combined with the P-wave and S-wave energy relationship and the energy distribution within a certain range, the source mechanism is estimated and the location information of the current microseismic event is saved.

[0015] The number of triggered detectors N≥3, and all detectors have been accurately timed before use, with each detector equipped with a self-triggering system.

[0016] The form of the three-dimensional elastic wave equation after decoupling of longitudinal and transverse waves is as follows:

[0017]

[0018] In the above formula: ρ is density; α and β are the P-wave velocity and S-wave velocity, respectively; v x v y and v z τ represents the velocity components of the particle vibration in the x, y, and z directions, respectively; 111 τ 111x τ 111y τ 111z τ 12τ 13 and τ 23 These are the stress components of a particle vibrating along seven different directions; It is a first-order time differential operator; and The velocity components of the particles associated with the longitudinal wave vibrating along the x, y, and z directions. and The velocity components of the particles associated with the transverse wave vibrating along the x, y, and z directions. Separate operators for spatial directions associated with P-waves and S-waves, respectively; f x f y f z λ is the source term for the three-component data; λ and μ are a set of Lamé parameters, and their relationships with the P-wave and S-wave velocities are as follows:

[0019]

[0020]

[0021] Wavefield simulation involves transforming the known discretized P- and S-wave velocity models into Lamé parameters λ and μ within the decoupled three-dimensional elastic wave equations and substituting them into the equations. Using the coordinates of each geophone as the source location information, numerical simulations are performed with the data recorded by N geophones as source terms in the equations, thereby obtaining the seismic wavefields of these N geophones at imaging time T. and

[0022] In step S3, a specific imaging region in three-dimensional space is obtained using directional aperture constraints. The steps are as follows:

[0023] For the vector seismic wave characteristics of each receiving point, the real part Ω of the eigenvector corresponding to the largest eigenvalue of the covariance matrix in the time-frequency domain is calculated. x Ω y Ω z The azimuth angle θ = arctan(Ω) of the seismic wave propagation was estimated. x / Ω y and tilt angle Estimate the approximate incident direction of the seismic waves at the current receiving point. Similarly, calculate the incident directions of the seismic waves at other receiving points in three-dimensional space to obtain the spatial imaging region £ at the intersection. Within region £, use a migration imaging condition without wavefield storage to perform imaging in a specific area. The imaging results of P-waves and S-waves are obtained using the following imaging conditions:

[0024]

[0025]

[0026] In the above formula: i represents the wave field corresponding to different detector data, Π represents the product of wave fields, and I α Indicates the imaging result of the longitudinal wave, I β The image represents the imaging result of the transverse wave, and T represents the imaging time of the spatial wave field.

[0027] The detector is a three-component detector.

[0028] The advantages of this invention over the prior art are as follows:

[0029] (1) The present invention employs GPU thread parallelism, shared memory acceleration and multi-GPU parallel simulation, which can ensure that computational efficiency can be significantly improved.

[0030] (2) The wave equation of the present invention uses decoupled longitudinal and transverse waves, so there is no need to carry out data preprocessing scheme to separate longitudinal and transverse waves, which reduces errors and reduces workload.

[0031] (3) The vector characteristics of the seismic waves at the receiving point were used to determine the imaging domain of the wave equation, which greatly reduced the ambiguity in migration positioning and improved the interpretability of migration imaging.

[0032] (4) The use of storage-free positioning and imaging conditions greatly reduces the hard disk storage of massive wave fields in conventional offset positioning and reduces the overhead of computing resources.

[0033] (5) The positioning method proposed in this invention can adapt to the positioning of microseismic signals recorded by different numbers of detectors, and can process microseismic signals in real time and save positioning result information. Attached Figure Description

[0034] The present invention will be further described below with reference to the accompanying drawings:

[0035] Figure 1 This is a diagram showing the wave field simulation results used in this invention;

[0036] Figure 2 This is a schematic diagram of GPU differential calculation in this invention;

[0037] Figure 3 This is a real-time positioning result diagram of the actual data collected in this invention. Detailed Implementation

[0038] like Figures 1 to 3 As shown, the present invention provides a real-time fluctuation offset localization method for microseismic signals, comprising the following steps:

[0039] S1: Within the monitoring area, multiple geophones are deployed at different depths as receiving points for microseismic data. All geophones are connected to the main acquisition unit via a main line, forming a three-dimensional observation system composed of multiple planes. All geophones are precisely timed and equipped with self-triggering systems. Due to microseismic events such as rock fracturing occurring at different locations, only geophones within the effective range can sense and trigger the recording signal. When a microseismic event occurs in the coal seam working face and its roof and floor, several geophones near the seismic source acquire microseismic data.

[0040] S2: Count the number of triggered detectors N (N≥3), and distribute this N across different GPUs to conduct high-performance seismic wavefield simulations. Specifically, in the P-wave and S-wave decoupling equations, input the known discretized P-wave and S-wave velocity models and the coordinate information of each detector, and use the data recorded by the N detectors as source terms in the equations to conduct numerical simulations, thereby obtaining the seismic wavefields of these N detectors at imaging time T. and

[0041]

[0042] Using the P-wave and S-wave decoupled three-dimensional elastic wave equations for adaptive decoupling in wavefield solving avoids the need for preprocessing methods that separate P-waves and S-waves. The form of the P-wave and S-wave decoupled three-dimensional elastic wave equations is as follows:

[0043]

[0044] In the above formula: ρ is density; α and β are the P-wave velocity and S-wave velocity, respectively; v x v y and v z τ represents the velocity components of the particle vibration in the x, y, and z directions, respectively; 111 τ 111x τ 111y τ 111z τ 12 τ 13 and τ 23 These are the stress components of a particle vibrating along seven different directions; It is a first-order time differential operator; and The velocity components of the particles associated with the longitudinal wave vibrating along the x, y, and z directions. and The velocity components of the particles associated with the transverse wave vibrating along the x, y, and z directions. Separate operators for spatial directions associated with P-waves and S-waves, respectively; f x f y f zλ is the source term for the three-component data; λ and μ are a set of Lamé parameters, and their relationships with the P-wave and S-wave velocities are as follows:

[0045] By solving the aforementioned wave equations, wavefield information related to and independent of both P-waves and S-waves from any seismic source can be obtained. The efficiency of solving these wave equations determines the efficiency of migration imaging. In this invention, the key parts of the differential calculation are mapped to a GPU mesh, and a thread-parallel computing mode is adopted to improve computational efficiency. At the same time, in order to reduce frequent data access in differential calculations, the data used for differential calculations is pre-copied to shared memory, further improving computational efficiency.

[0046] S3: To avoid the ambiguity of full-space imaging, a directional aperture-constrained, storageless migration imaging method is employed. Specifically, for each receiver point's vector seismic wave characteristics, the real part Ω of the eigenvector corresponding to the largest eigenvalue of the covariance matrix in the time-frequency domain is calculated. x Ω y Ω z The azimuth angle θ = arctan(Ω) can be estimated. x / Ω y and tilt angle The approximate incident direction of the seismic wave at the current receiving point can be estimated. Similarly, the incident direction of the seismic wave at other receiving points in three-dimensional space can be calculated to obtain the spatial imaging region £ at the intersection.

[0047] Imaging is performed within a specific region using a wavefield-free migration imaging condition, and the imaging results for both P-waves and S-waves are obtained using the following imaging conditions:

[0048]

[0049]

[0050] In the above formula: i represents the wave field corresponding to different detector data, Π represents the product of wave fields, and I α Indicates the imaging result of the longitudinal wave, I β The image represents the imaging result of the transverse wave, and T represents the imaging time of the spatial wave field. The imaging results show that when the above imaging conditions are used, it is not necessary to store the wave field in space along time; only the wavelength value at the corresponding imaging time needs to be recorded in memory.

[0051] For the P-wave and S-wave imaging results of this microseismic event I α and I β Imaging fusion is performed to extract the imaging results with the most focused energy and consistent imaging point locations, which are used as microseismic information for the microseismic event. Furthermore, the energy information at the microseismic location is used to estimate the microseismic energy J = A.α +A β A α and A β Let A and B represent the energy at the source microseismic point on the imaging profile, respectively. α A represents the expansion and contraction properties of the earthquake source. β It represents the shear properties of the earthquake source, and combined with the relationship between the P-wave and S-wave energy and the energy distribution within a certain range, the ratio of the two energies can be used to determine and estimate the source mechanism, and save the location information of the current microseismic event.

[0052] The source intensity can be obtained from the P-wave and S-wave imaging results. By combining the detection direction of each detector receiving point, the location of the source can be estimated.

[0053] Following the above procedure, subsequent microseismic events will be handled using a similar approach, processing the microseismic events recorded by the observation system in real time. The only difference is that each time a microseismic event is recorded, the records and coordinate information of the corresponding triggered geophones must be loaded.

[0054] The detector in this invention can be a three-component detector.

[0055] The method of the present invention will be further described below with reference to the accompanying drawings.

[0056] In this embodiment, a spherical anomalous medium is used for simulation. The P-wave and S-wave velocities of the anomalous body are 3600 m / s and 1800 m / s, respectively, while the P-wave and S-wave velocities in the background are constant at 3000 m / s and 2160 m / s, respectively, with a constant density. The model is discretized into 301×301×301 grid nodes with a spatial interval of 15 m and a time step of 2.0 ms. The first derivative of the 18 Hz Gaussian function is used as the source, excited in the center of the model. The calculated decoupled P-wave and S-wave fields are as follows: Figure 1 As shown in the figure. It can be seen that by using the decoupling equation of the present invention, effective separation of longitudinal and transverse wave fields can be achieved without mutual interference.

[0057] like Figure 2 The diagram shows a high-performance computing scheme used in this invention, which copies CPU data to the GPU and uses gridded threads to perform differential calculations in batches in parallel, thereby improving computing efficiency.

[0058] Figure 3This study focuses on the location calculation of microseismic events collected during on-site monitoring at a coal mine. Sixteen geophones were deployed at different depths within the area. Source location and energy calculation were performed on the microseismic signals recorded over a period of time. Approximately 300 microseismic events were triggered in total. The program automatically selected data from those recorded by at least three geophones for each microseismic event, and valid location calculations were performed for a total of 245 microseismic events. The location method followed the steps outlined in this invention. The final location results are as follows... Figure 3 As shown, the blue dot represents the detector's micro-vibration, the triangle represents the micro-vibration event being located, and the size of the triangle represents the energy intensity of the micro-vibration event. The positioning results show that the real-time positioning scheme proposed in this invention is completely effective in practical field positioning applications.

[0059] This invention employs P-wave and S-wave decoupling equations to perform adaptive P-wave and S-wave decoupling and separation. GPU parallel processing is used in the differential solution of the wave equations, and GPU shared memory is used to accelerate and optimize the processing of frequently accessed data. This scheme proposes a flexible processing approach for a non-fixed number of detectors. For data from multiple detectors, multi-GPU parallel processing is used to obtain the wavefield information at the imaging time. Finally, the P-wave and S-wave images are obtained through multiplication operations between the wavefield data from N detectors, enabling the comprehensive acquisition of information such as source location, energy, and source mechanism.

[0060] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A real-time wave migration positioning method for microseismic signals, characterized in that: The method comprises the following steps: S1: In the area range of the coal mine working face floor needing to be monitored, a plurality of geophones are arranged as microseismic data receiving points according to different depths, the geophones are connected to the acquisition host through the main line, and microseismic data are collected after microseismic occurs in the coal seam working face and its roof; S2: The number N of triggered geophones is counted, the data of the N geophones are distributed to different GPUs, seismic wave field simulation of the N geophones is carried out in parallel, and wave field information at the imaging moment in the three-dimensional space is obtained, wherein the three-dimensional elastic wave wave equation of P-S wave decoupling is used for self-adaptive decoupling in wave field solving, and GPU parallel processing is used in the difference solving of the wave equation; The three-dimensional elastic wave wave equation of P-S wave decoupling has the following form: ； where: p is the density; a and b are P-wave and S-wave velocities, respectively; v x , y , and v z are velocity components of a particle's vibration in x, y, and z directions, respectively; t 111 , t 111x , t 111y , t 111z , t 12 , t 13 , and t 23 are stress components of a particle's vibration along seven different directions; t is a first-order time differential operator; , , and are velocity components of a particle's vibration in x, y, and z directions, respectively, associated with a longitudinal wave, , , and are velocity components of a particle's vibration in x, y, and z directions, respectively, associated with a transverse wave, are spatial directional derivative operators associated with longitudinal and transverse waves, respectively; f x , f y , and f z are source terms of three-component data; and l and m are a set of Lame parameters, and their relations with longitudinal and transverse wave velocities are: a = l + 2m and b = l - m. ； ； Wavefield simulation is to transform the known discrete P-S wave velocity model into Lame parameters λ and μ in the three-dimensional elastic wave equation decoupled in P-S wave, and bring it into the equation, and take the coordinate information of each geophone as the source position information, and carry out numerical simulation with the data recorded by N geophones as the source term of the equation, and then obtain the seismic wave field of the N geophones at the imaging time T and ; S3: The imaging results of P waves and S waves are obtained by using the directional aperture constraint and the storage-free migration imaging method, the imaging results of P waves and S waves are fused, the imaging result with the most focused energy and the consistent imaging point position is extracted as the microseismic information of the microseismic event, the energy of the microseismic event is estimated using the energy information of the microseismic, the focal mechanism is estimated combining the P-S wave energy relationship and the energy distribution within a certain range, and the positioning information of the current microseismic event is saved; In step S3, the specific imaging area in the three-dimensional space is obtained by using the directional aperture constraint, and the steps are as follows: For the vector seismic wave characteristics of each receiving point, the real part Ω of the eigenvector corresponding to the largest eigenvalue of the covariance matrix in the time-frequency domain is calculated. x Ω y Ω z The azimuth angle θ = arctan(Ω) of the seismic wave propagation was estimated. x / Ω y ) and tilt angle φ = arctan(Ω) z / The approximate incident direction of the seismic waves at the current receiving point is estimated. Similarly, the incident directions of the seismic waves at other receiving points in three-dimensional space are calculated to obtain the spatial imaging region £ at the intersection. Imaging is carried out in a specific area within region £ using a migration imaging condition without wavefield storage. The imaging results of P-waves and S-waves are obtained using the following imaging conditions: ； ； In the above formula: i represents the wave field corresponding to different detector data, Π represents the multiplication between wave fields, I α represents the imaging result of the longitudinal wave, I β represents the imaging result of the transverse wave, T represents the imaging time of the spatial wave field.

2. The real-time wave migration method for microseismic signal positioning according to claim 1, characterized in that: The number N of triggered geophones is greater than or equal to 3, and all the geophones are accurately time-synchronized before use, and each geophone is provided with a self-triggering system.

3. The real-time wave migration method for microseismic signal positioning according to claim 1 or 2, characterized in that: The geophone is a three-component geophone.

Images