Mine Ordovician limestone water pressure and top interface crack detection method based on micro-seismic technology

Abstract

The invention provides a mine Ordovician limestone water pressure and top interface crack detection method based on a micro-seismic technology. The method comprises the following steps: collecting actual micro-seismic data; de-noising processing is carried out on the collected actual micro-seismic data by using a wavelet threshold de-noising algorithm; a P-S wave combined positioning algorithm based on P wave and S wave first arrival is adopted to position the denoised micro-seismic event; seismic source parameter calculation and seismic source mechanism inversion are carried out based on the positioned micro-seismic event waveform and the initial polarity; ground stress inversion work is carried out, and the change condition of Ordovician limestone water pressure is deduced; calculating advanced attributes representing stratum stress or fracture characteristics on the basis of micro-seismic event source parameters, determining an advanced attribute threshold range causing Ordovician limestone water damage, and judging possible crack conditions; the Ordovician limestone top crack condition and Ordovician limestone water pressure are comprehensively analyzed, and whether the risk that Ordovician limestone water protrudes along the crack exists or not is judged. The method can monitor the water pressure change of the Ordovician limestone in real time, accurately detect the water flowing fracture of the Ordovician limestone top interface, and provide a scientific basis for mine water disaster prevention and control.

Description

Technical Field

[0001] This invention relates to the field of mine water hazard prevention technology, specifically to a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology. Background Technology

[0002] In the field of coal mining, especially in the North China coalfields, complex and variable geological conditions and hydrological environments are encountered. Among these, Ordovician limestone (referred to as "Ordovician limestone"), due to its intense karst development, high water content, and high water pressure, has become a key focus and challenge in mine water hazard prevention and control. Ordovician limestone water hazards not only affect normal mine production but can also trigger serious safety accidents, threatening the lives of mine workers. Therefore, effectively monitoring changes in Ordovician limestone water pressure and accurately detecting water conduits has become a crucial technical issue in mine water hazard prevention and control.

[0003] Traditional methods for preventing mine water hazards mainly rely on geological exploration, hydrological observation, and experience-based judgment. While these methods can identify some water-bearing channels to a certain extent, they have the following limitations: Insufficient precision: Traditional methods often fail to accurately determine the specific location and scale of water diversion channels, resulting in a lack of targeted prevention and control measures.

[0004] Poor real-time performance: Traditional methods rely on periodic or irregular observation data, making it difficult to achieve real-time monitoring and early warning of water pressure changes.

[0005] High cost: Some traditional methods require a large investment of manpower, material resources and financial resources, which increases the operating costs of the mine.

[0006] Currently, although some monitoring technologies have been applied to mine water hazard prevention and control, such as water level monitoring and flow rate monitoring, these technologies still have the following shortcomings in monitoring changes in Ordovician limestone water pressure and detecting water-conducting fractures: Single parameter monitoring: Existing technologies mostly focus on monitoring a single parameter, such as water level or flow rate, which is difficult to fully reflect the complex process of water disasters.

[0007] Low spatial resolution: Some monitoring technologies have low spatial resolution, making it difficult to accurately identify minute water-conducting fissures or changes.

[0008] Weak anti-interference ability: In the complex mining environment, existing monitoring technologies are easily affected by various interference factors, resulting in inaccurate or distorted monitoring data.

[0009] With increasing mining depth and more complex geological conditions, changes in Ordovician limestone water pressure have a significant impact on mine safety. Therefore, developing a technology capable of real-time monitoring of Ordovician limestone water pressure changes is crucial. This technology should possess high precision, real-time performance, and strong anti-interference capabilities to promptly detect water pressure anomalies and implement appropriate preventative measures. Summary of the Invention

[0010] This invention provides a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology. This method can monitor changes in Ordovician limestone water pressure in real time and accurately detect water-conducting cracks at the top interface of Ordovician limestone, providing a scientific basis for mine water hazard prevention and control.

[0011] The present invention achieves the above objectives through the following technical solutions: A method for detecting water pressure and top interface fractures in Ordovician limestone mines based on microseismic technology includes: Collect actual microseismic data and conduct preliminary analysis of the collected actual microseismic data; The wavelet threshold denoising algorithm is applied to denoise the collected actual microseismic data. The optimal parameter combination is determined by transforming the key filtering parameters to improve the recognition accuracy of microseismic signals. A PS-wave joint localization algorithm based on the first arrival of P-waves and S-waves is used to locate the denoised microseismic events in order to determine their spatial location. Based on the waveform and first arrival polarity of the microseismic event after location, source parameters are calculated and source mechanism is inverted; Using the source mechanism inversion results, geostress inversion work was carried out to infer the variation of Ordovician limestone water pressure. Specifically, the orientation and relative magnitude of triaxial geostress in the distribution area of ​​microseismic events were obtained, and the variation of Ordovician limestone water pressure was inferred by analyzing the changes in the relative magnitude of the maximum and minimum horizontal principal stresses in the geostress field. Based on the source parameters of microseismic events, high-level properties characterizing formation stress or fracture features are calculated, and their spatial distribution is analyzed to assess the incubation and development of water-conducting fractures at the top interface of Ordos limestone; the threshold range of high-level properties leading to water damage in Ordos limestone is determined, and the possible fracture conditions are identified. A comprehensive analysis of the cracks in the top of the Ordovician limestone and the water pressure in the Ordovician limestone was conducted to determine whether there was a risk of water protruding along the cracks.

[0012] According to the present invention, a method for detecting water pressure and top interface fractures in Ordovician limestone mines based on microseismic technology includes the following steps in the wavelet threshold denoising algorithm: The original waveform data, consisting of M samples of actual microseismic data, is first subjected to a 5-100Hz bandpass filter to initially remove some noise interference from the original microseismic data, thus obtaining the filtered microseismic data. For the microseismic data after bandpass filtering, wavelet basis functions and decomposition levels J are selected. The Mallat fast algorithm is used to transform the filtered microseismic data from the time domain to the wavelet domain, obtaining the wavelet coefficients for each decomposition level. and , is represented as:

[0013] in, D j ( n ), s ( n () represents the filtered microseismic data sequence. ; Based on the characteristics and noise level of the microseismic data, a threshold function and a threshold value are selected to evaluate the wavelet coefficients at each decomposition level. w dj ( j = J , 1) Perform thresholding to remove wavelet coefficients corresponding to noise, and retain or enhance wavelet coefficients corresponding to effective microseismic signals to obtain processed wavelet coefficients; The wavelet coefficients after thresholding are compared with the unprocessed approximation coefficients. The wavelet domain data is reconstructed using the following inverse wavelet transform formula to restore it to the time domain, yielding the denoised actual microseismic data:

[0014] in, A j (n) and D j (n) These are wavelet basis functions. A J ( n ) is the scaling function. D j ( n ) is a wavelet function.

[0015] According to the present invention, a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology is provided. This method uses a PS wave joint localization algorithm to locate denoised microseismic events, and includes the following steps: In actual microseismic monitoring, the first arrival information of P-waves and S-waves generated by microseismic events is acquired. When the signal-to-noise ratio of S-waves is low or invisible, the time difference of arrival method based on the first arrival information of P-waves is used first for positioning. If the signal-to-noise ratio of S-waves meets the requirements, the first arrival information of P-waves and S-waves is used simultaneously. The arrival time of P-waves and S-waves received by each sensor is recorded through the sensor network in the microseismic monitoring system. Under known geological conditions in a mine, the travel time between the seismic source and each detector is calculated based on the wave velocity model of the medium within the mine. For each assumed location of the seismic source, the theoretical travel time from the source to each sensor is calculated. The calculated time difference between the arrival times of the same seismic source at different sensors is set as follows:

[0016] in, t i and t j From the epicenter to the first i The and the first j The theoretical travel time of each sensor, while the actual measured time difference between detectors is... ; Construct an objective function with the travel time residuals of the seismic source points as variables. The objective function expression is as follows:

[0017] in, n The objective function, representing the number of sensors, measures the difference between the theoretical time difference and the actual measured time difference.

[0018] According to the present invention, a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology is provided. Based on the acquired first arrival information of P-waves and S-waves, the Geiger positioning method and grid search method are applied to perform microseismic source inversion and positioning. During the positioning process, ray tracing forward modeling is combined with constraints applied using the time residuals of theoretical and actual P-waves and S-waves, as well as the angular residuals of the source vector, to obtain the source positioning results. A network of microseismic monitoring sensors deployed in the mine continuously monitors and collects microseismic signals in real time. When a microseismic event occurs, the sensors capture waveform signals containing microseismic information, which at least include relevant information from the P-waves and S-waves. For each microseismic event that has been located, return the waveform data for inspection; if the residual of the microseismic event is large, adjust the initial arrival of the event, update the initial arrival information, and then re-perform the positioning operation. When the theoretical arrival curve is basically consistent with the actual arrival curve, it indicates that the positioning is relatively accurate; if the difference is large, it proves that the positioning accuracy is poor; for events with poor accuracy, return to the initial arrival picking interface to pick again and perform secondary positioning until the accurate positioning position is obtained. After all microseismic events are located, the results of all located microseismic events are reviewed. For microseismic events that do not meet the quality control requirements, return to the microseismic event picking window and adjust the initial arrival picking position to meet the quality control requirements. For microseismic events where the source location deviates from the main seam, it is necessary to examine the waveform information a second time. Specifically, this involves observing whether the waveform travel time characteristics match the detector arrangement characteristics and determining whether the detector closest to the source location receives the valid signal first.

[0019] According to the present invention, a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology includes the following steps when calculating source parameters: Based on the waveform data of microseismic events recorded by the microseismic monitoring system, the Brune model was selected to calculate the source parameters and moment magnitude. Calculate the seismic moment based on the Brune model. The calculation formula is:

[0020] in, For the density of the medium, The S-wave velocity is used when the recorded microseismic S-wave energy is dominant. If the recorded microseismic P-wave energy is dominant, it is replaced by the wave velocity corresponding to the P-wave. This refers to the displacement spectrum amplitude; r The distance from the seismic source to the geophone is determined based on the seismic source location results and the geophone location information. Rv As a correction factor; Determine the radius of the earthquake source r 0, when the S-wave energy of a microseismic event is dominant, is represented as:

[0021] in, For S-wave velocity, The corner frequency is the amplitude spectrum of the S-wave. When the recorded microseismic P-wave energy is dominant, the wave velocity and corner frequency can be changed to the corresponding values ​​of the P-wave. Based on the calculated seismic moment The moment magnitude is calculated using the following formula, expressed as:

[0022] Calculate stress drop as needed , is represented as:

[0023] Among them, stress drop reflects the stress changes caused by fault slip.

[0024] According to the present invention, a method for detecting water pressure and top interface fractures in Ordovician limestone mines based on microseismic technology includes the following steps in the source mechanism inversion: The sign information of the initial polarity of the P-wave generated by the microseismic event is obtained by using a microseismic monitoring system, which is used to reflect the initial directional characteristics of particle vibration during the propagation of the P-wave. The grid search algorithm in DeepListen software is used to find the pair of orthogonal nodal surfaces with the smallest discrepancy sign ratio Ψ, with the goal of finding the best fit between the theoretical P-wave initial motion polarity and the observed P-wave initial motion polarity sign. During the search process, the possible range of source mechanism parameters is divided into grids. By traversing all grid points, the discrepancy sign ratio between the theoretical P-wave initial motion polarity and the observed value is calculated for each grid point. A set of orthogonal nodal surfaces is obtained through a grid search algorithm, one of which is a potential fault plane and the other is an auxiliary surface; The optimal fault plane was determined by comparing the degree of matching between different fault plane combinations and the geostress state, as well as the influencing factors on the water pressure and top interface fracture formation in the mine. The determined optimal source mechanism results are fed back into the mine Ordovician limestone water pressure and top interface fracture detection system based on microseismic technology. This is used to analyze the correlation between microseismic events and changes in Ordovician limestone water pressure and the development of top interface fractures, so as to optimize the mine water hazard prevention and geological structure monitoring strategies.

[0025] According to the present invention, a method for detecting the water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology is provided. Given that the source mechanism inversion yields a set of orthogonal nodal surfaces and there are multiple solutions, the fault instability coefficient is used to determine the true crack surface. Project any nodal surface solution onto the normalized Mohr circle, and draw a straight line parallel to the fracture criterion line at the projection point. The perpendicular distance from the location of the normalized maximum normal stress to this line is the fault instability coefficient I of the nodal surface solution. The stress field inversion is based on the Wallace-B assumption, which states that the slip direction of the fault plane is parallel to the direction of the maximum shear stress along the fault plane. A stress field inversion model is constructed based on the Wallace-B assumption, clarifying the crack slip vector s and the normalized stress vector t (t = [τ]). 11 τ 12 τ 13 τ 22 τ 23 ] The relationship between them is At = s, where A is a 3×5 matrix, and the matrix elements are calculated from the normal vectors of the fault plane, expressed as:

[0026] in, n 1. n 2. n 3 represents the components of the normal vector of the fault plane.

[0027] When multiple composite focal mechanism solutions exist, the stress moment tensor is obtained by least squares inversion. The slip vectors s and matrices A corresponding to the multiple focal mechanism solutions are integrated to construct an overdetermined set of equations. The least squares method is used to solve the set of equations to obtain the optimal stress moment tensor. The stress moment tensor obtained by least squares inversion is decomposed into eigenvalues. The orientation and tilt angle of multiple principal stresses are calculated by the eigenvalue decomposition method. The order of magnitude of the principal stresses is determined according to the magnitude of the eigenvalues. The direction of the eigenvector corresponding to the eigenvalue is the direction of the principal stress, thereby determining the state of the stress field. The stress field results obtained by inversion were applied to the detection of Ordovician limestone water pressure and top interface fractures in the mine, and the intrinsic relationship between stress field and Ordovician limestone water pressure changes and fracture development was analyzed.

[0028] According to the present invention, a method for detecting hydraulic pressure and top interface fractures in Ordovician limestone mines based on microseismic technology is provided. When using cluster analysis for advanced attribute calculation, the cluster analysis is divided into static and dynamic modes. During microseismic monitoring, the calculation mode is determined based on the occurrence of microseismic events: if there are microseismic events that have not yet occurred, i.e., the monitoring process is still ongoing and new events may occur, the dynamic parameter calculation mode is used to describe the rupture process of the strata; when all microseismic events have occurred, the static parameter calculation mode is switched to. The following dynamic parameters are calculated according to the dynamic parameter calculation mode: Coseismic Deformation and Strain Rate: By analyzing the waveform data generated by microseismic events, and combining the mechanical properties of the formation medium and the geometric characteristics of the monitoring area, the coseismic deformation and corresponding strain rate of the formation at the time of each microseismic event are calculated to reflect the instantaneous deformation of the formation under microseismic action. Event time, distance interval, and event density: The occurrence time of each microseismic event is recorded, and the time interval between events is calculated; simultaneously, based on the source location results, the distance interval between events is calculated. Based on the temporal and spatial distribution information of events, the event density per unit time and unit spatial area is calculated to analyze the spatiotemporal evolution characteristics of microseismic activity; Radiated energy and seismic moment: By combining microseismic waveform information with seismological formulas, the radiated energy and corresponding seismic moment of each microseismic event are calculated to quantify the intensity and scale of the microseismic event.

[0029] According to the present invention, a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology is provided. This method calculates important dynamic clustering attributes based on dynamic parameters, namely, stress index, plasticity index, and diffusion index. The formula for calculating the stress index SI is:

[0030] Where ηc is the cumulative seismic efficiency of clustered events, c is the coefficient for different lithologies, times, and volumes, and DI is the diffusion index; when the SI value is low, it indicates that the stress transmission is unstable and confined to a very small area of ​​the reservoir; when the SI value is large, it represents that the microseismic events release a large amount of energy, the stress field is more stable, and the final rupture range is larger. The formula for calculating the plasticity index PI is:

[0031] Where μ is the dynamic shear modulus and ηc is the cumulative seismic efficiency of clustering events; a low PI value indicates that the reservoir is not easily deformed; a high PI value indicates that the reservoir is extremely easily deformed. The formula for calculating the diffusion index DI is:

[0032] Among them, X 2 denoted as the square of the average distance between clustered events, and t as the average time interval between clustered events; a lower diffusion index (DI) value indicates a large time interval between events but a small distance interval; a larger DI value indicates that the calculated stress index, plasticity index, and diffusion index are comprehensively analyzed to establish a multi-factor relationship model between them and the water pressure of the Ordovician limestone mine and the development of cracks at the top interface.

[0033] According to the present invention, a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology is provided. Static clustering attributes are obtained by calculating static parameters, namely, microseismic radiation energy, apparent stress, and apparent volume. Microseismic radiation energy: For each microseismic event, its displacement spectrum is integrated; using microseismic waveform data recorded by the microseismic monitoring system, displacement spectrum information is obtained through spectrum analysis, and the energy released by the microseismic event is calculated according to the following formula. :

[0034] in, ρ Density of the formation medium; The velocity of seismic waves; R The distance from the earthquake source to the monitoring point; F Ω is the attenuation factor; Ω is the displacement spectral function. fFor frequency; Apparent stress: used to obtain the shear modulus of reservoir rocks μ and the seismic moment of each microseismic event M 0. Calculate the apparent stress according to the apparent stress calculation formula. :

[0035] in, It is the shear modulus of the reservoir rock. The seismic moment of this microseismic event; Apparent volume: Apparent volume Characterizing the volume of inelastically deformable rock at the seismic source, the calculated energy released by the microseismic event is used. E Shear modulus of reservoir rocks μ Substituting the source parameters of the microseismic event into the following apparent volume calculation formula, calculate the apparent volume characterizing the inelastic deformation of the rock at the source. :

[0036] Therefore, compared with existing technologies, the method for detecting Ordovician limestone water pressure and its top interface cracks based on microseismic technology proposed in this invention can monitor the dynamic changes of Ordovician limestone water pressure in real time by continuously collecting and analyzing microseismic event data. Compared with traditional monitoring methods that rely on single parameters such as water level and flow rate observed periodically or irregularly, this invention achieves continuous and real-time tracking of water pressure changes, and has the following beneficial effects: 1. The microseismic monitoring system of this invention can operate 24 hours a day without interruption, promptly detecting abnormal fluctuations in water pressure, providing immediate early warning for mines, and effectively preventing water inrush accidents caused by sudden changes in water pressure. Microseismic technology can capture minute pressure changes, which are often precursors to water hazards. By monitoring these changes in real time, preventative measures can be taken in advance to reduce the risk of accidents.

[0037] 2. This invention is not limited to real-time water pressure monitoring; it also utilizes advanced attribute analysis of microseismic events to accurately detect the distribution of water-conducting fractures at the top interface of Ordovician limestone. For example: Three-dimensional positioning capability: Combining the positioning information of microseismic events with the source mechanism inversion results, three-dimensional spatial positioning of water-conducting fractures can be achieved, accurately reflecting the fracture's orientation, dip angle, and scale. Quantitative assessment: By calculating advanced attributes such as stress index and energy index, the development degree and water-conducting capacity of fractures can be quantitatively assessed, providing a scientific basis for prevention and control measures such as grouting modification. Multi-source data fusion: This invention integrates and analyzes microseismic monitoring data with multi-source data such as geological exploration and hydrological observation, improving the accuracy and reliability of water-conducting fracture detection.

[0038] 3. Compared to traditional methods, this invention, through an automated and intelligent microseismic monitoring system, significantly shortens the data acquisition and processing cycle, improving work efficiency. Microseismic technology, characterized by high precision and high resolution, can detect minute cracks and pressure changes that are difficult to detect using traditional methods, thereby improving the accuracy of water hazard prediction. Through real-time monitoring and precise detection, this invention can promptly identify and eliminate potential water hazard risks, effectively reducing mine operational risks. Simultaneously, scientifically sound prevention and control measures also reduce unnecessary resource waste and production costs.

[0039] 4. The implementation of this invention also helps to promote the intelligent construction of mines and improve the overall safety level: As an important component of intelligent mine construction, the microseismic monitoring system can automatically collect, transmit, and process data, providing intelligent decision support for mine management. Through real-time monitoring and early warning, it enhances the safety awareness of mine workers and promotes the construction and development of a safety culture. The application of this invention helps to achieve safe, efficient, and green mining, and promotes the sustainable development of the coal industry.

[0040] The present invention will now be described in further detail with reference to the accompanying drawings and specific embodiments. Attached Figure Description

[0041] Figure 1 This is a flowchart of an embodiment of the present invention, which is a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology.

[0042] Figure 2 This is a schematic diagram of an embodiment of the present invention, which is a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology.

[0043] Figure 3 This is a schematic diagram of the overall hydrogeological columnar structure of the working face bottom plate in an embodiment of the present invention, which is based on microseismic technology for detecting water pressure and cracks at the top interface of Ordovician limestone mine.

[0044] Figure 4 This is an XY plane view of the location of a microseismic event in an embodiment of a method for detecting water pressure and cracks at the top interface of Ordovician limestone mines based on microseismic technology according to the present invention.

[0045] Figure 5 This is an XZ side view of the location of a microseismic event in an embodiment of a method for detecting water pressure and cracks at the top interface of Ordovician limestone mines based on microseismic technology according to the present invention.

[0046] Figure 6 This is a YZ side view of the location of a microseismic event in an embodiment of a method for detecting water pressure and cracks at the top interface of Ordovician limestone mines based on microseismic technology according to the present invention.

[0047] Figure 7 This is a schematic diagram illustrating the source mechanism inversion in an embodiment of a method for detecting hydraulic pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0048] Figure 8 This is a schematic diagram of the underground fault displacement and P-wave initial polarity in an embodiment of a method for detecting underground fault displacement and top interface fractures in Ordovician limestone mines based on microseismic technology according to the present invention.

[0049] Figure 9 This is a schematic diagram of the source mechanism of 20 microseismic events in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0050] Figure 10 This is a schematic diagram of the fault instability coefficient in an embodiment of a method for detecting water pressure and top interface fractures in Ordovician limestone mines based on microseismic technology according to the present invention.

[0051] Figure 11 This is a schematic diagram illustrating the principle of microseismic event clustering analysis in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0052] Figure 12 This is a planar distribution diagram of the dynamic parameters of the stress index of microseismic events in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0053] Figure 13 This is a planar distribution diagram of the dynamic parameters of the plasticity index of microseismic events in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0054] Figure 14 This is a planar distribution diagram of the dynamic parameters of the microseismic event diffusion index in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0055] Figure 15 This is a schematic diagram of the displacement spectrum of a typical microseismic event in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0056] Figure 16 This is a schematic diagram illustrating the energy distribution properties in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0057] Figure 17 This is a schematic diagram of the apparent stress properties in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention.

[0058] Figure 18This is a schematic diagram of the apparent volume properties in an embodiment of a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology according to the present invention. Detailed Implementation

[0059] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of this invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without creative effort are within the scope of protection of this invention.

[0060] In this document, the term "embodiment" means that a particular feature, structure, or characteristic described in connection with an embodiment may be included in at least one embodiment of this application. The appearance of this phrase in various places throughout the specification does not necessarily refer to the same embodiment, nor is it a separate or alternative embodiment mutually exclusive with other embodiments. It will be explicitly and implicitly understood by those skilled in the art that the embodiments described herein can be combined with other embodiments.

[0061] In this embodiment, a geostress inversion method based on the source mechanism inversion of microseismic events is used to obtain the orientation and relative magnitude of triaxial geostress in the distribution area of ​​microseismic events. The variation of Ordovician limestone water pressure is analyzed by examining changes in the geostress field, especially the relative magnitudes of the maximum and minimum horizontal principal stresses. Advanced properties of microseismic events, such as stress index and energy index, based on source parameters (seismic moment, stress drop, etc.), are studied to analyze the incubation and development of water-conducting fractures at the top interface of the Ordovician limestone. The relationship between these advanced properties and the occurrence of water hazards in the Ordovician limestone is statistically analyzed, and the threshold range of these advanced properties is determined. Ultimately, a method for preventing and controlling water hazards in the Ordovician limestone based on microseismic monitoring is developed.

[0062] See Figures 1 to 18 This embodiment provides a method for detecting water pressure and top interface cracks in Ordovician limestone mines based on microseismic technology, including: Step S1: Collect actual microseismic data and perform preliminary analysis on the collected data; Step S2: Apply wavelet threshold denoising algorithm to denoise the collected actual microseismic data. In this step, the optimal parameter combination is determined by transforming the key filtering parameters to improve the recognition accuracy of the microseismic signal. Step S3: The PS wave joint localization algorithm based on the first arrival of P wave and S wave is used to locate the denoised microseismic event in order to determine the spatial location of the microseismic event. Step S4: Calculate source parameters and invert source mechanism based on the waveform and first arrival polarity of the located microseismic event; Step S5: Use the source mechanism inversion results to carry out in-situ stress inversion work and infer the changes in Ordovician limestone water pressure; in this step, obtain the orientation and relative magnitude of the triaxial in-situ stress in the distribution area of ​​microseismic events, and infer the changes in Ordovician limestone water pressure by analyzing the changes in the relative magnitude of the maximum and minimum horizontal principal stresses in the in-situ stress field. Step S6: Based on the source parameters of the microseismic event, calculate the high-level properties characterizing the stress or fracture features of the formation, analyze their spatial distribution to assess the incubation and development of water-conducting fractures at the top interface of the Ordovician limestone; determine the threshold range of the high-level properties that lead to water damage in the Ordovician limestone, and determine the possible fracture conditions. Step S7: Analyze the top crack condition and water pressure of the Ordovician limestone to determine whether there is a risk of water protruding along the crack.

[0063] As can be seen, the method provided in this embodiment first analyzes the quality of the original microseismic waveform data and performs preprocessing such as removing bad channels and calculating noise levels; secondly, it uses wavelet threshold denoising method to denoise the preprocessed data, thereby improving the signal-to-noise ratio of microseismic events, especially weak events. The P-wave and S-wave first arrivals of microseismic events are picked up automatically or manually, and the polarity information of the P-wave first arrival is saved for focal mechanism inversion; then, the P-wave joint location method is used to locate the microseismic events with high precision. Using the waveforms and first arrival polarities of the located microseismic events, we can calculate the source parameters (moment magnitude, seismic moment, source radius, stress drop, corner frequency, etc.) and invert the source mechanism. Using the crack orientation, dip angle and slip angle information obtained by focal mechanism inversion, in-situ stress inversion was carried out to obtain the orientation and relative magnitude of the maximum and minimum horizontal principal stress and vertical stress in the distribution area of ​​microseismic events in Ordos limestone and its top interface. By comparing the changes in stress field over time with the changes in measured water pressure in Ordos limestone, the correspondence between the two was determined. Based on the source parameters, high-level attributes such as stress index and energy index are calculated, and high-level attribute planar thermal maps and cross-sectional thermal maps of the microseismic distribution area of ​​Ordovician limestone and its top interface are drawn. Areas with high values ​​of stress index and energy index reflect fracture development areas, and the distribution of fracture zones is characterized based on this. Finally, based on the results of the above steps, the detection of water pressure in Ordovician limestone and water-conducting fractures at its top interface in mines is realized using microseismic technology.

[0064] In step S1 above, this embodiment relies on the microseismic monitoring project for the Ordovician limestone water pressure and top interface fracture detection of the 11916 working face of the Hebei Gequan Mine East Shaft. The Hebei Gequan Mine East Shaft is a major production mine of Jizhong Coal Mine Co., Ltd., with a designed production capacity of 90 × 10⁴ t·a⁻¹. The main coal seam is No. 9 coal, and the mining technology is strike longwall fully mechanized mining. The strike length of the 11916 working face is 1080 m, the dip length is 70 m, the coal seam dip angle is 7°–21°, the average thickness is 5.5 m, the elevation difference between the two roadways is about 20 m, and the 11916 material transport roadway is separated from the goaf of the 11915 working face by only a 4 m thick coal pillar. Furthermore, there is a sinkhole with a diameter exceeding 50 m in the middle of the working face. The hydrogeological columnar section of the working face floor is shown below. Figure 3 As shown, the area between the base plate and the Benxi limestone is a water-resistant layer with an average thickness of approximately 20.3 m and moderate water-blocking performance. The area between the base plate and the Ordovician limestone aquifer is a water-resistant layer with a thickness of 36.0–43.6 m, with an average thickness of 41.1 m. The Benxi karst fissure aquifer (referred to as "Benxi limestone") has moderate water-bearing capacity and is relatively thin, while the extremely thick Ordovician karst fissure aquifer (referred to as "Ordovician limestone") has strong water-bearing capacity. Benxi limestone and Ordovician limestone are the main aquifers in the working face, with Ordovician limestone karst water being the main target for mine water hazard prevention and control.

[0065] According to borehole measurements, the water pressure in the Ordovician limestone aquifer at the bottom of the working face is 1.71–2.21 MPa, and the water inrush coefficient is 0.047–0.061 MPa·m⁻¹, indicating a risk of karst water inrush at the working face. Before mining, grouting modification of the aquifer at the bottom had been implemented. Due to a water inrush at working face 11913, the aquifer was modified later for safety reasons. However, the bottom face has a large dip angle, significant undulations, and large variations in the thickness of the aquifer. Furthermore, a collapse column with a diameter exceeding 50 m exists in the adjacent working face 11915. Drilling confirmed that the collapse column contained water, with a flow rate of 60 m³·h⁻¹. After grouting modification, the inspection borehole still showed a flow rate of 2 m³·h⁻¹, indicating a risk of water inrush. Therefore, this embodiment utilizes microseismic monitoring technology to establish a water inrush monitoring system to monitor the water conduit.

[0066] Microseismic sensors can be installed in the material transport roadway, transport roadway, and east wing transport roadway near the 11916 working face. The east wing transport roadway will not collapse during mining. In order to protect the sensors and avoid damage to the cables, and to better monitor the fracture of the goaf floor, microseismic sensors are installed in the east wing transport roadway and material transport roadway respectively.

[0067] Therefore, considering the actual topography and features of the survey area, this project, based on the microseismic monitoring of the Ordovician limestone water pressure and its top interface crack detection at the 11916 working face of the Gequan Mine in Hebei Province, designed a total of 43 physical measuring points. Among them, 14 are located in each of the roadways on both sides of the 11916 working face, and 3 are arranged in the boreholes of each of the roadways. In order to better monitor the depth of damage to the floor, this embodiment selects the No. 2 return airway (81308) for monitoring, because this roadway will not collapse during the mining process of the 11916 working face, and the sensors deployed in this roadway can better monitor the damage to the floor in the goaf.

[0068] The microseismic data processing and interpretation in this embodiment adopts... Figure 2 The diagram illustrates the technical flow for microseismic data processing and interpretation. First, the raw microseismic waveform data undergoes raw data analysis and preprocessing, and the SEG-Y data is edited, converting the raw SEG-Y to standard SEG-Y. Second, wavelet threshold denoising technology is used to suppress noise and improve the signal-to-noise ratio. Then, the first arrival of microseismic events is manually picked to ensure accuracy. Next, a layered uniform velocity model is established, and the picked first arrivals are used to locate microseismic events. Using the located microseismic events and waveforms, source parameters, including moment magnitude, energy, seismic moment, and focal radius, are calculated. Source mechanism inversion is performed using the located microseismic events and first arrival polarity, followed by geostress inversion based on the source mechanism inversion. Finally, the high-level properties of the microseismic events are comprehensively analyzed by integrating the source parameters, source mechanism inversion, and geostress inversion results.

[0069] In step S2 above, a combined denoising method of bandpass filtering and wavelet thresholding is used to improve the signal-to-noise ratio of microseismic events. In addition to the more common bandpass filtering denoising, wavelet transform, a method developed based on Fourier analysis that can effectively handle non-stationary signals, is also selected for the identification and denoising of previous seismic signals.

[0070] Bandpass filtering is a common noise reduction technique for microseismic monitoring data. The bandpass filter used in this data preprocessing was a Butterworth bandpass filter with an 8th order, a minimum cutoff frequency of 5Hz, and a maximum cutoff frequency of 100Hz.

[0071] Compared to Fourier transform and short-time Fourier transform, wavelet transform is an adaptive time- and frequency-local transform with excellent time-frequency localization characteristics, making it suitable for extracting information from abruptly changing signals. Wavelet transform possesses multi-resolution analysis (also known as multi-scale analysis) characteristics, meaning it can decompose signals into subspaces with different resolutions, allowing for the extraction of different features at various resolutions depending on the specific problem. Wavelet transform not only offers the Mallat fast algorithm for signal tower-style multi-resolution decomposition and reconstruction but also provides a wealth of research results on threshold denoising for selection and reference.

[0072] Specifically, the wavelet thresholding denoising algorithm includes the following steps: The original waveform data, consisting of M samples of actual microseismic data, is first subjected to a 5-100Hz bandpass filter to initially remove some noise interference from the original microseismic data, thus obtaining the filtered microseismic data. For the microseismic data after bandpass filtering, a suitable wavelet basis function and decomposition level J are selected. The Mallat fast algorithm is used to transform the filtered microseismic data from the time domain to the wavelet domain, and the wavelet coefficients under each decomposition level are obtained. and , is represented as:

[0073] in, D j ( n ), s ( n () represents the filtered microseismic data sequence. ; Based on the characteristics and noise level of the microseismic data, appropriate threshold functions and thresholds are selected to evaluate the wavelet coefficients at each decomposition level. w dj ( j = J , 1) Perform thresholding to remove wavelet coefficients corresponding to noise, and retain or enhance wavelet coefficients corresponding to effective microseismic signals to obtain processed wavelet coefficients; The wavelet coefficients after thresholding are compared with the unprocessed approximation coefficients. The wavelet domain data is reconstructed using the following inverse wavelet transform formula to restore it to the time domain, yielding the denoised actual microseismic data:

[0074] in, A j (n) and D j (n) These are wavelet basis functions. A J ( n ) is the scaling function. D j ( n ) is a wavelet function.

[0075] For original records with a high signal-to-noise ratio, the wavelet denoising process described above requires that the true seismic signal be preserved to the greatest extent possible, and that no change in the phase of the first arrival wave is observed. For original records with a very low signal-to-noise ratio and where the effective signal is almost "buried" in the noise, the wavelet denoising process requires that the seismic signal be recovered to meet the requirements of subsequent first arrival picking and processing.

[0076] In step S3 above, the PS wave joint localization algorithm is used to locate the denoised microseismic events, including the following steps: In actual microseismic monitoring, the first arrival information of P-waves and S-waves generated by microseismic events is acquired. When the S-wave signal-to-noise ratio is low or invisible, the time difference of arrival method based on the first arrival information of P-waves is preferentially used for positioning. If the S-wave signal-to-noise ratio meets the requirements, the first arrival information of both P-waves and S-waves is used simultaneously. The arrival time of P-waves and S-waves received by each sensor is recorded through the sensor network in the microseismic monitoring system. Under known geological conditions in a mine, the travel time between the seismic source and each detector is calculated based on the wave velocity model of the medium within the mine. For each assumed source location, the theoretical travel time from the source to each sensor is calculated. The calculated time difference between the arrival times of the same seismic source at different sensors is set as follows:

[0077] in, t i and t j From the epicenter to the first i The and the first j The theoretical travel time of each sensor, while the actual measured time difference between detectors is... ; Construct an objective function with the travel time residuals of the seismic source points as variables. The objective function expression is as follows:

[0078] in, n The objective function, representing the number of sensors, measures the difference between the theoretical time difference and the actual measured time difference.

[0079] Since there are three unknowns (the three-dimensional coordinates of the seismic source), this is a three-dimensional localization problem, requiring data from at least four sensors. Optimization algorithms, such as least squares, genetic algorithms, or particle swarm optimization, are used to minimize the objective function to determine the location of the seismic source. Minimizing the objective function's value, i.e., minimizing the error between the theoretical and actual time difference, yields the precise location of the microseismic event. However, this is not included in the inversion parameters. 0, theoretically, has better inversion convergence.

[0080] Based on the acquired first arrival information of the P and S waves, the Geiger localization method and grid search method are applied to perform microseismic source inversion and localization. During the localization process, ray tracing forward modeling is combined with constraints based on the time residuals of the theoretical and actual P and S waves and the angular residuals of the source vectors, thereby improving the accuracy of the localization and ultimately obtaining the source localization results. In particular, a network of microseismic monitoring sensors deployed in the mine continuously monitors and collects microseismic signals in real time. When a microseismic event occurs, the sensors capture waveform signals containing microseismic information, which at least contain relevant information of the P and S waves. For each microseismic event that has been located, waveform data is returned for inspection. If the residual of a microseismic event is large, the initial arrival is finely adjusted for that event, the initial arrival information is updated, and the location operation is repeated to ensure that the location result is more accurate.

[0081] Accuracy analysis was conducted on the positioning results, with the residual error between the picked initial arrival and the theoretical initial arrival set at less than 5ms and the positioning uncertainty at less than 10m as accuracy standards. These standards were used to determine whether the positioning results met the accuracy requirements, providing a basis for subsequent operations.

[0082] The waveform window interface verifies the accuracy of the positioning result by comparing the difference between the theoretical arrival time and the actual arrival time curve. When the theoretical arrival time and the actual arrival time curve are basically consistent, it indicates that the positioning is relatively accurate; if the difference is large, it proves that the positioning accuracy is poor. For events with poor accuracy, return to the initial arrival picking interface to re-pick and perform secondary positioning until the accurate positioning position is obtained.

[0083] After all microseismic events have been located, review the results of all located microseismic events. For microseismic events that do not meet the quality control requirements, return to the microseismic event picking window and adjust the initial arrival picking position to meet the quality control requirements.

[0084] For microseismic events where the source location deviates from the main seam, it is necessary to examine the waveform information a second time. Specifically, this involves observing whether the waveform travel time characteristics match the detector arrangement characteristics and determining whether the detector closest to the source location receives the valid signal first. This ensures the reliability of the positioning results under special circumstances.

[0085] This embodiment processed a continuously acquired microseismic dataset from June 25, 2019 to October 28, 2019, successfully locating 2576 microseismic events. The locations of the microseismic events are as follows: Figure 4-6 As shown.

[0086] In step S4 above, source parameters were calculated, and the moment magnitude was obtained. Mw Parameters such as seismic moment, focal radius, stress drop, and corner frequency were obtained. After obtaining the location results, the moment magnitude, focal radius, seismic moment, stress drop, corner frequency, and radiated energy were recalculated using manual first arrival and waveform data. Due to the potential inaccuracy of detector sensitivity and gain information, the moment magnitude only provides a relative magnitude.

[0087] Source parameters are important parameters describing the energy of microseismic events, including moment magnitude. M w Seismic moment, focal radius, stress drop, and corner frequency are generally calculated during the moment magnitude calculation. Focal parameters provide a measure of the relative magnitude of microseismic events and offer fundamental parameters for magnitude-frequency relationship analysis. Moment magnitude M w Moment magnitude is currently the most ideal physical quantity for measuring earthquake magnitude. Compared to other traditional magnitude scales, moment magnitude does not saturate. It can measure moment magnitude for all earthquakes, regardless of size or depth, and regardless of whether far-field or near-field seismic wave data, or what type of geodetic or geological data is used. Furthermore, it can be compared with well-known magnitude scales such as surface wave magnitude. M S Connecting to each other. Moment magnitude is a uniform magnitude scale, suitable for statistics covering a wide range of magnitudes.

[0088] When calculating the source parameters, the Brune model is used and the corner frequency is calculated. The following steps are performed on the original waveform of the microseismic event, the first arrival waveform of the P-wave, the spectrum, and the fitting based on the Brune model: Based on the waveform data of microseismic events recorded by the microseismic monitoring system, the Brune model was selected for calculating source parameters and moment magnitude. This model is suitable for moment magnitude calculation, can be integrated with well-known magnitude scales such as surface wave magnitude, and uses a uniform magnitude scale, making it suitable for statistical analysis of a wide range of magnitudes.

[0089] Calculate the seismic moment based on the Brune model. The calculation formula is:

[0090] in, The density of the medium can be obtained from mine geological data; The S-wave velocity is used when the recorded microseismic S-wave energy is dominant. If the recorded microseismic P-wave energy is dominant, it is replaced by the wave velocity corresponding to the P-wave. The displacement spectrum amplitude is obtained through spectral analysis of microseismic waveform data; r The distance from the seismic source to the geophone is determined based on the seismic source location results and the geophone location information. RvThe correction factors for radiation factors and geometric diffusion effects related to instrument response and propagation path can be calculated based on detector characteristics and propagation path theory or determined through calibration experiments.

[0091] Determine the radius of the earthquake source r 0, when the S-wave energy of the microseismic event is dominant:

[0092] in, For S-wave velocity, The corner frequency is the amplitude spectrum of the S-wave. When the recorded microseismic P-wave energy is dominant, the wave velocity and corner frequency can be changed to the corresponding values ​​for the P-wave.

[0093] Based on the calculated seismic moment The moment magnitude is calculated using the following formula:

[0094] Based on the calculation of moment magnitude, stress drop can be calculated as needed. :

[0095] Among them, stress drop reflects the stress changes caused by fault slippage, which can provide a reference for analyzing the causes of microseismic events and the stress state of the geological structure of mines.

[0096] In step S4 above, the focal mechanism is a high-level interpretation result of microseismic monitoring, and much of the subsequent in-depth analysis relies on the focal mechanism. For example... Figure 7 As shown, focal mechanism information, besides providing two possible fault planes for the fault (hydraulic fracture), can also provide the strike, dip, and slip angle of the hydraulic fracture. Combined with the focal radius calculated from the focal parameters, it can quantitatively describe the strike, dip, and fracture length of the hydraulic fracture. Furthermore, focal mechanism can be used for in-situ stress inversion, providing the direction and relative magnitude of the maximum principal stress in three directions, which is fundamental for conducting microseismic geomechanical research. The purpose of focal mechanism inversion in this embodiment is primarily to obtain the strike and dip information of underground fractures corresponding to microseismic events, and to use this information to model the hydraulic fracture using a discrete fracture network based on the focal mechanism. Focal mechanism information not only helps to understand the mechanism of rock fractures and underground stress distribution, but also provides important parameters for generating discrete fracture networks (DFNs) and estimating reservoir stimulation volumes (SRVs). Focal mechanisms are generally described using fault plane solutions and focal moment tensors (MTs). For pure shear rupture events, the focal mechanism can generally be characterized using a dual-couple (DC) model. However, for events related to the propagation of hydraulic fractures, the source mechanism is complex and usually requires the use of moment tensors for characterization.

[0097] Among these methods, using the initial motion polarity of P-waves to solve for focal mechanism solutions is widely applied in earthquake event analysis, focal parameter determination, and stress field research. The theoretical basis for using P-wave initial motion polarity information to invert fracture focal mechanism solutions is that when a fault undergoes pure shear displacement underground, the initial motion polarity of the direct P-waves received at different locations on the ground differs (e.g., ...). Figure 8 As shown in the figure, its projection on the focal sphere exhibits a four-quadrant distribution. There are several methods for solving focal mechanism solutions using the P-wave initial motion polarity, which can be broadly classified into iterative methods and grid search methods. Because iterative algorithms are more dependent on data quality, the grid search method is currently the mainstream solution method.

[0098] Specifically, the focal mechanism inversion includes the following steps: The initial polarity sign information of P-waves generated by microseismic events is obtained using a microseismic monitoring system. This information is the basic data for subsequent grid search algorithm processing and reflects the initial directional characteristics of particle vibration during the propagation of P-waves.

[0099] The grid search algorithm in DeepListen software was employed, aiming to find the pair of orthogonal nodal planes with the optimal sign fit between the theoretical and observed P-wave initial motion polarities, specifically the pair with the smallest discrepancy sign ratio Ψ (Ψ = number of discrepancies / total number of symbols, 0 ≤ Ψ ≤ 1). During the search, the possible range of source mechanism parameters was divided into grids. By traversing all grid points, the discrepancy sign ratio between the theoretical and observed P-wave initial motion polarities for each grid point was calculated.

[0100] A set of orthogonal nodal surfaces was obtained through a grid search algorithm, one of which is a potential fault plane and the other is an auxiliary surface. However, for the case of a pure shear source (double couple source), these two surfaces are equivalent, and there are two possible combinations of strike and dip angles. The fault plane information cannot be uniquely determined by relying solely on the P-wave initial motion focal mechanism solution.

[0101] Due to the ambiguity of the two sets of fault planes given by the focal mechanism inversion results, geostress inversion was conducted. Combining the geostress characteristics, geological structure background, and other relevant auxiliary information of the mining area, a comprehensive analysis and verification of the two possible fault planes was performed. By comparing the degree of matching between different fault plane combinations and the geostress state, as well as the influence on the Ordovician limestone water pressure and the formation of fractures at the top interface of the mine, the optimal fault plane was determined, eliminating ambiguity and providing a reliable basis for accurately analyzing the causes of microseismic events and the stability of the mine's geological structure.

[0102] The determined optimal focal mechanism results are fed back into the mine Ordovician limestone water pressure and top interface fracture detection system based on microseismic technology. This allows for a more accurate analysis of the correlation between microseismic events and changes in Ordovician limestone water pressure and top interface fracture development, thereby optimizing mine water hazard prevention and geological structure monitoring strategies. Simultaneously, based on the actual application results, necessary adjustments and optimizations are made to the focal mechanism inversion method and parameters.

[0103] Given that the focal mechanism inversion yields a set of orthogonal nodal surfaces (one potential fault surface and one auxiliary surface) and exhibits multiple solutions, a fault instability coefficient is used to determine the true fracture surface. Any nodal surface solution is projected onto a normalized Mohr's circle, and a straight line parallel to the rupture criterion line is drawn at the projection point. The perpendicular distance from the location of the normalized maximum normal stress to this line is the fault instability coefficient (I) of that nodal surface solution. A larger instability coefficient indicates a more unstable fracture surface and a greater likelihood of slippage. Nodal surfaces with larger instability coefficients are identified as the true fracture surfaces to eliminate the influence of multiple solutions in the focal mechanism.

[0104] This embodiment performed focal mechanism inversion on all microseismic events. For example... Figure 9 As shown, Figure 9 The image shows beach balls representing the focal mechanisms of the 20 microseismic events with the highest signal-to-noise ratios. Different colored solid lines represent the optimal fault plane orientation. The above focal mechanism inversion results indicate that the region mainly develops normal faults (fractures), with most fractures trending northeast.

[0105] In this embodiment, the focal mechanism inversion yields a set of orthogonal nodal planes, one of which is the fault plane and the other is an auxiliary plane. However, for a pure shear source (double couple source), these two planes are equivalent. Therefore, the fault plane information cannot be determined solely by the focal mechanism solution of the P-wave initial motion; other auxiliary information is needed. For example, the formation stress state obtained from well logging data can be used to guide the selection of the fault plane; or, after obtaining the focal mechanism solutions of a large number of microseismic events, their focal mechanism spheres can be plotted on the same map according to the event location coordinates. The fracture attitudes of adjacent microseismic events should be relatively similar, thus determining which set of nodal plane solutions to use to represent the fracture attitude; or, the optimal set of fault plane solutions can be obtained through geostress inversion. Therefore, conducting reservoir-scale geostress inversion based on focal mechanism can obtain the distribution characteristics of local geostress field (relative magnitudes of maximum and minimum horizontal principal stresses and vertical stresses, as well as their azimuth and dip angles). On the one hand, it can eliminate the multiple solutions of focal mechanism and obtain the optimal set of fault / fracture surface strike and dip angle from the two sets of strike and dip angle information of shear fractures. On the other hand, the subsequent modeling of hydraulic fracture network uses this optimal set of strike and dip angle.

[0106] Focal mechanism inversion can provide two sets of nodal surface solutions. Without other data constraints, it's impossible to determine which one represents the true fracture surface or fault plane. Using an incorrect nodal surface solution for stress field inversion can significantly impact the accuracy of the solution. This project intends to use the fault instability coefficient to determine the true fracture surface and employ an iterative method to invert the stress field parameters based on the true fault surface solution (Vavry). (uk, 2014), and uniquely determine a set of optimal fault / fracture surface strikes and dips, eliminating the ambiguity of focal mechanisms. like Figure 10 As shown, any nodal surface solution can be projected onto the normalized Mohr circle, and a straight line parallel to the fracture criterion line can be drawn at the projection point. The perpendicular distance from the location of the normalized maximum normal stress to this line is the fault instability coefficient I of that nodal surface solution. The larger the instability coefficient, the more unstable the fracture surface, and the easier it is for slip to occur. Therefore, in the calculation, the nodal surface with a larger instability coefficient can be regarded as the real fracture surface.

[0107] The stress field inversion is primarily based on the Wallace-B assumption, which states that the slip direction of the fault plane is parallel to the direction of maximum shear stress along the fault plane. A stress field inversion model is constructed based on this assumption, clarifying the crack slip vector s and the normalized stress vector t (t = [τ]). 11 τ 12 τ 13 τ 22 τ 23 ] The relationship between them is At = s, where A is a 3×5 matrix, and the matrix elements are calculated from the normal vectors of the fault plane, specifically:

[0108] in n 1. n 2. n 3 represents the components of the normal vector of the fault plane.

[0109] When multiple composite focal mechanism solutions exist, the stress moment tensor is obtained through least squares inversion. The slip vectors s and matrices A corresponding to the multiple focal mechanism solutions are integrated to construct an overdetermined system of equations. The least squares method is used to solve this system of equations to obtain the optimal stress moment tensor, so as to comprehensively consider the influence of multiple microseismic events on the stress field.

[0110] The stress moment tensor obtained through least-squares inversion is decomposed into eigenvalues. The orientation and tilt angle of the three principal stresses are then calculated using this method. The magnitude order of the principal stresses is determined based on the eigenvalues, and the direction of the eigenvector corresponding to each eigenvalue is the direction of the principal stress, thus completely determining the state of the stress field.

[0111] The stress field results obtained from the inversion were applied to the detection of Ordovician limestone water pressure and top interface fractures in the mine, analyzing the intrinsic relationship between the stress field and changes in Ordovician limestone water pressure and fracture development. Simultaneously, the inverted stress field results were calibrated and verified by combining actual mine geological conditions, other geological exploration data, and subsequent monitored microseismic events, ensuring the accuracy and reliability of the stress field inversion results and providing a scientific basis for safe mine production and geological disaster prevention.

[0112] In this embodiment, geostress inversion was performed on the focal mechanism inversion information of the microseismic event. As shown in Table 1, after geostress inversion, information on the local geostress field of the area covered by the microseismic event in the upper boundary strata of the Ordovician limestone was obtained, including the azimuth, dip angle and relative magnitude of the maximum horizontal principal stress, minimum horizontal principal stress and vertical stress, and the optimal fault attitude was given.

[0113] Table 1: Results of geostress inversion

[0114] The stress field inversion results for this region indicate that: (1) Relative magnitude of reservoir stress in the study area: maximum horizontal principal stress S Hmax Minimum horizontal principal stress S hmin Vertical stress S v ; (2) Vertical stress Sv tilt angle 88.5°, maximum horizontal principal stress S Hmax The azimuth is 55.85° (N55.85°E), and the minimum horizontal principal stress is S. hmin The bearing is 325.83° (N34.17°W).

[0115] In step S5 above, analyzing only the spatial (xyz), temporal (t), magnitude (mw), and energy attributes of a single event is insufficient to accurately analyze the rupture and deformation process of the formation. Cluster analysis of multiple events, however, can accurately analyze the rupture response of the formation. When using cluster analysis for advanced attribute calculations, cluster analysis is divided into static and dynamic modes. During microseismic monitoring, the parameter calculation mode is determined based on the occurrence of microseismic events: if there are microseismic events that have not yet occurred, i.e., the monitoring process is still ongoing and new events may occur, the dynamic parameter calculation mode is used to describe the rupture process of the formation; once all microseismic events have occurred, the static parameter calculation mode is switched to. The principle of microseismic event cluster analysis is as follows: Figure 11 As shown. ; where, according to the dynamic parameter calculation mode, the following dynamic parameters are calculated: Coseismic deformation and strain rate: By analyzing the waveform data generated by microseismic events, and combining the mechanical properties of the formation medium and the geometric characteristics of the monitoring area, the coseismic deformation and corresponding strain rate of the formation at each microseismic event are calculated to reflect the instantaneous deformation of the formation under microseismic action. Event time, distance interval, and event density: The occurrence time of each microseismic event is recorded, and the time interval between events is calculated; simultaneously, based on the source location results, the distance interval between events is calculated. Based on the temporal and spatial distribution information of events, the event density per unit time and unit spatial area is calculated to analyze the spatiotemporal evolution characteristics of microseismic activity; Radiated energy and seismic moment: Using information such as the amplitude and frequency of microseismic waveforms, combined with relevant seismological formulas, the radiated energy and corresponding seismic moment of each microseismic event are calculated to quantify the intensity and scale of the microseismic event.

[0116] Important dynamic clustering properties were obtained based on dynamic parameter calculations, namely stress index, plasticity index, and diffusion index: The formula for calculating the stress index SI is:

[0117] Where ηc is the cumulative seismic efficiency of clustered events, c is the coefficient for different lithologies, times, and volumes, and DI is the diffusion index. The geological significance represented by the magnitude of the stress index is analyzed. When the SI value is low, it indicates that stress transmission is unstable and confined to a small area of ​​the reservoir, usually driven by fluids with increased pore pressure, resulting in microseismic events with relatively small energy release. When the SI value is high, it represents that the microseismic events release large amounts of energy, the stress field is more stable, and the final rupture range is larger, usually a fracture or fault pointing towards the direction of maximum principal stress, triggered by stress-induced microseismic events. The calculated stress index is applied to the detection of Ordovician limestone water pressure and top interface fractures in the mine, analyzing its relationship with changes in Ordovician limestone water pressure and fracture development, providing a basis for assessing formation stress state and fracture hazard. The formula for calculating the plasticity index PI is:

[0118] Where μ is the dynamic shear modulus and ηc is the cumulative seismic efficiency of clustered events. Understanding the reservoir deformation characteristics reflected by the plasticity index (PI) is crucial. A lower PI value indicates that the reservoir is not easily deformable; a higher PI value indicates that the reservoir is extremely easily deformable. By combining the plasticity index with orogenic hydraulic pressure and top interface fracture detection in the mine, its impact on reservoir deformation and fracture development is analyzed, providing a reference for predicting reservoir stability.

[0119] The formula for calculating the diffusion index DI is:

[0120] Among them, X 2 The diffusion index (DI) is the square of the average distance between clustered events, where t is the average time interval between clustered events. A lower DI value indicates a large time interval between events but a small distance interval; a higher DI value indicates a smaller time interval between events but a larger distance interval. The calculated stress index, plasticity index, and diffusion index are comprehensively analyzed to establish a multi-factor relationship model between them and the Ordovician limestone water pressure and top interface fracture development in the mine. Based on monitoring data and geological changes in practical applications, the calculation method and analysis model for clustering attributes are optimized through feedback, continuously improving the accuracy and reliability of Ordovician limestone water pressure and top interface fracture detection based on microseismic technology.

[0121] The dynamic parameter calculation results are stored in the form of a two-dimensional dynamic parameter base map (heat map) generated at the time of each microseismic event, allowing viewing of the dynamic response of the formation at certain time intervals (every minute or at the time of each microseismic event). In this embodiment, a corresponding video is generated and displayed in the project results report multimedia (.pptx). Since the report can only display static images, only the dynamic parameter heat map corresponding to the time of the last microseismic event is shown. Figures 12-14 The figures show the planar distribution of stress index, plasticity index, and diffusion index, indicating that: the stress index is higher in the central and southwestern parts of the detector coverage area, and the strata fractures are more concentrated in these local areas; the plasticity index is relatively evenly distributed throughout the area, and the top Ordovician limestone strata are not prone to fracture; the diffusion index is higher in the central part of the detector coverage area, and the impact range is larger when the strata fracture in these local areas occurs.

[0122] Important static clustering properties were obtained based on static parameter calculations, namely microseismic radiation energy (referred to as "energy"), apparent stress, and apparent volume. Energy: For each microseismic event, its displacement spectrum is integrated; using microseismic waveform data recorded by the microseismic monitoring system, displacement spectrum information is obtained through spectrum analysis techniques, such as... Figure 15 As shown, the energy released by the microseismic event is then calculated using the following formula. :

[0123] in, ρ The density of the formation medium can be obtained from mine geological data; The velocity of the seismic wave is determined based on the lithology of the strata and the type of wave. R The distance from the earthquake source to the monitoring point is calculated based on the earthquake source location results and the monitoring point location information. FThe geometric diffusion and attenuation factors, which are related to instrument response and propagation path, can be determined through calibration experiments or theoretical calculations; Ω is the displacement spectral function. f For frequency.

[0124] Apparent stress: used to obtain the shear modulus of reservoir rocks μ and the seismic moment of each microseismic event M 0. Calculate the apparent stress according to the apparent stress calculation formula. :

[0125] in, It is the shear modulus of the reservoir rock. The seismic moment of this microseismic event is used to analyze the geological significance of the apparent stress value. A lower apparent stress value indicates the opening of natural fractures, while a higher apparent stress value indicates the formation of numerous new fractures. The calculated apparent stress is correlated with the mine's Ordovician limestone water pressure and the development of fractures at the top interface to assess the ease of fracture formation and the occurrence of new fractures, providing a reference for predicting fracture propagation and development. Energy and apparent stress are often used to describe the static response of mine formations.

[0126] Apparent volume: Apparent volume Characterizing the volume of inelastically deformable rock at the seismic source, the calculated energy released by the microseismic event is used. E Shear modulus of reservoir rocks μ And the source parameters of microseismic events (such as source intensity parameters, which can be obtained in conjunction with the source mechanism inversion results and used to calculate equivalent pressure) P The equivalent pressure already obtained P） Substitute the following formula for apparent volume to calculate the apparent volume representing the volume of the inelasticly deformed rock at the earthquake source. :

[0127] The apparent volume reflects the volumetric characteristics of inelastically deformable rock at the earthquake source. By combining mine hydraulic pressure and top interface fracture detection, the relationship between apparent volume and formation deformation and fracture spatial distribution is studied, providing a basis for assessing the deformation range and fracture development scale of the formation.

[0128] The calculated microseismic radiation energy, apparent stress, and apparent volume are comprehensively analyzed to establish their interrelationships and a comprehensive analytical model relating them to the Ordovician limestone water pressure and top interface fracture development in mines. Based on monitoring data and geological changes in practical applications, the static attribute analysis method and calculation model are optimized through feedback, continuously improving the accuracy and reliability of static attribute analysis in the detection of Ordovician limestone water pressure and top interface fractures based on microseismic technology. This provides more scientific decision support for safe mine production and geological disaster prevention.

[0129] Static property analysis results: Energy property map, apparent stress property map, and apparent volume property map are as follows: Figures 16-18 As shown, this indicates that: (1) The energy distribution is relatively uniform throughout the region, with no local abnormal areas; (2) The apparent stress distribution is relatively uniform throughout the area, with no local abnormal areas; The apparent volume values ​​are higher in the central and southwestern parts of the microseismic coverage area, indicating that the rocks in this area have undergone inelastic deformation.

[0130] In summary, taking the 11916 working face of the East Shaft of Gequan Mine in Hebei Province as an example, the method of this invention successfully monitored the changes in Ordovician limestone water pressure and detected potential water-conducting fracture zones within the working face. Specifically, 43 physical monitoring points were deployed, and microseismic data were continuously collected from June 25, 2019 to October 28, 2019. Through processing and analysis of this data, the changing trend of Ordovician limestone water pressure and the distribution of water-conducting fracture zones were obtained, providing strong support for safe mine production.

[0131] This invention proposes a method for monitoring water pressure and detecting water-conducting fractures in Ordovician limestone mines based on microseismic technology. This method achieves precise monitoring and detection of water pressure changes and water-conducting fractures in the Ordovician limestone mine through real-time monitoring of microseismic events, inversion of the geostress field, and analysis of advanced properties. This method has advantages such as high real-time performance, high accuracy, and high efficiency, and is of great significance for improving the level of mine water hazard prevention and control.

[0132] The technical features of the above embodiments can be combined in any way. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

[0133] The above embodiments are merely preferred embodiments of the present invention and should not be construed as limiting the scope of protection of the present invention. Any non-substantial changes and substitutions made by those skilled in the art based on the present invention shall fall within the scope of protection claimed by the present invention.

Claims

1. A method for detecting mine Ordovician limestone water pressure and its top interface fracture based on microseismic technology, characterized in that, The method comprises the following steps: Collecting actual microseismic data and performing preliminary analysis on the collected actual microseismic data; Applying a wavelet threshold denoising algorithm to the collected actual microseismic data for denoising processing, wherein the key filtering parameters are transformed to determine the optimal parameter combination, so as to improve the identification accuracy of the microseismic signal; Adopting a P-S wave joint positioning algorithm based on P wave and S wave first arrivals to position the denoised microseismic events, so as to determine the spatial position of the microseismic events; Performing source parameter calculation and source mechanism inversion based on the positioned microseismic event waveforms and first arrival polarities; Using the source mechanism inversion result to carry out geostress inversion work and deduce the change of the Ordovician limestone water pressure; Based on the source parameters of the microseismic events, calculating advanced attributes representing the stress or rupture characteristics of the stratum, analyzing the spatial distribution of the advanced attributes to evaluate the development of the water-conducting fracture of the Ordovician limestone top interface, determining the threshold range of the advanced attributes leading to the Ordovician limestone water disaster, and judging the possible fracture conditions; Comprehensively analyzing the Ordovician limestone top fracture conditions and the Ordovician limestone water pressure to determine whether there is a risk of Ordovician limestone water bursting along the fractures.

2. The method of claim 1, wherein, The wavelet threshold denoising algorithm comprises the following steps: The original waveform data composed of the collected actual microseismic data with M samples are first subjected to 5-100 Hz band-pass filtering processing to preliminarily remove part of the noise interference in the original microseismic data, so as to obtain filtered microseismic data; After the microseismic data is band-pass filtered, wavelet basis function and decomposition level J are selected, the filtered microseismic data is transformed from time domain to wavelet domain by using Mallat fast algorithm, and wavelet coefficients at each decomposition level are obtained and is expressed as: wherein D j ( n ), s ( n ) is a filtered microseismic data sequence, ; According to the characteristics and noise level of the microseismic data, a threshold function and a threshold are selected, the wavelet coefficients at each decomposition level are thresholded, the wavelet coefficients corresponding to the noise are removed, the wavelet coefficients corresponding to the effective microseismic signals are retained or enhanced, and the processed wavelet coefficients are obtained w dj ( j = J , ,1) are thresholded, the wavelet coefficients corresponding to the noise are removed, the wavelet coefficients corresponding to the effective microseismic signals are retained or enhanced, and the processed wavelet coefficients are obtained The wavelet coefficients after threshold processing and the unprocessed approximation coefficients The data in the wavelet domain is restored to the time domain by reconstructing through the following wavelet inverse transform formula to obtain the actual microseismic data after denoising: wherein A j (n) and D j (n) is a wavelet basis function, A J n is a scaling function, D j n is a wavelet function.​​ 3. The method of claim 1, wherein, The P-S wave joint positioning algorithm is used to position the denoised microseismic events, comprising the following steps: In the actual microseismic monitoring process, the P wave and S wave first arrival information generated by the microseismic events are acquired, the time difference method based on the P wave first arrival information is used for positioning when the S wave signal-to-noise ratio is low or invisible, and the P wave and S wave first arrival information is used simultaneously when the S wave signal-to-noise ratio meets the requirements; the arrival time of the P wave and S wave received by each sensor is recorded through the sensor network in the microseismic monitoring system; Under the known velocity of the mine geological conditions, the travel time between the source and each receiver is calculated according to the wave velocity model of the medium in the mine; The time difference between the same source reaching different sensors calculated is set as: wherein, t i and t j are the theoretical travel times from the source to the i thand j thsensor, respectively, while the actual measured time difference between the receivers is ; A target function with the travel time residual of the source point as the variable is constructed, and the expression of the target function is: wherein, n is the number of sensors, the objective function measures the degree of difference between the theoretical travel times and the actual measured travel times.

4. The method according to claim 3, characterized in that: According to the picked P and S wave first arrival information, the Geiger positioning method and the grid search method are applied to perform microseismic source inversion positioning; in the positioning process, the ray tracing forward is combined, and the time residual of the theoretical and actual P and S waves and the angle residual of the source vector are applied to perform constraint, so that the source positioning result is acquired; wherein the microseismic signal in the mine is continuously monitored and data is collected through the microseismic monitoring sensor network arranged in the mine; when the microseismic event occurs, the waveform signal containing the microseismic information is captured through the sensor, and the waveform signal at least contains the related information of the P and S waves; Check the waveform data of each microseismic event after positioning; if the residual of the microseismic event is large, adjust the first arrival of the event, update the first arrival information, and then reposition the event; When the theoretical arrival time curve is basically consistent with the actual arrival time curve, it indicates that the positioning is relatively accurate; if the difference is large, it proves that the positioning accuracy is poor; for events with poor accuracy, return to the first arrival picking interface for re-picking and secondary positioning until the accurate positioning position is obtained; After all the microseismic events are positioned, check the results of all the positioned microseismic events, and return to the microseismic event picking window for adjusting the first arrival picking position to meet the quality control requirements; For microseismic events with source position deviating from the main suture, the waveform information needs to be checked again, specifically whether the observed waveform travel time characteristics and the detector arrangement characteristics are consistent, and whether the detector with the shortest straight-line distance from the source position receives the effective signal first.

5. The method of claim 1, wherein, When calculating the source parameters, the following steps are performed: Based on the microseismic event waveform data recorded by the microseismic monitoring system, the Brune model is selected for source parameter calculation and moment magnitude calculation. The seismic moment is calculated according to the Brune model The calculation formula is: wherein, is the medium density, is the S-wave velocity, which is used when the recorded microseismic S-wave energy is dominant, and is replaced by the P-wave corresponding wave velocity if the recorded microseismic P-wave energy is dominant; is the displacement spectrum amplitude; r is the source-to-detector distance, which is determined according to the source location result and the detector position information; Rv is the correction factor; Determining source radius r 0, when microseismic S-wave energy dominates, is given by: where, is the S-wave velocity, is the corner frequency of the S-wave amplitude spectrum, when the recorded microseismic P-wave energy is dominant, the wave velocity and the corner frequency are changed to the corresponding quantities of the P-wave; The moment of the earthquake is calculated The moment magnitude is calculated using the following formula, which is expressed as: Stress reduction is calculated as needed is expressed as: Among them, the stress drop reflects the stress change caused by fault slip.

6. The method of claim 1, wherein, The source mechanism inversion includes the following steps: Use the microseismic monitoring system to obtain the sign information of the P-wave initial polarity of the microseismic event, which reflects the initial direction characteristics of the particle vibration of the P-wave in the propagation process; Use the grid search algorithm in the DeepListen software to seek the pair of orthogonal joint surfaces with the smallest contradiction symbol ratio Ψ, that is, the best fitting of the theoretical P-wave initial polarity and the observed P-wave initial polarity symbol; in the search process, the possible value range of the source mechanism parameters is grid divided, and the contradiction symbol ratio of the theoretical P-wave initial polarity and the observed value corresponding to each grid point is calculated by traversing all grid points; A set of orthogonal joint surfaces is obtained through the grid search algorithm, one of which is the potential fault plane and the other is the auxiliary plane; Determine the optimal fault plane by comparing the matching degree of different fault plane combinations with the in-situ stress state and the influencing factors of the Ordovician limestone water pressure and the top interface fracture; Feed back the determined optimal source mechanism result to the mine Ordovician limestone water pressure and top interface fracture detection system based on microseismic technology, and analyze the correlation between microseismic events and Ordovician limestone water pressure changes and top interface fracture development to optimize mine water disaster prevention and geological structure monitoring strategies.

7. The method of claim 6, characterized in that: In view of the multiple solutions of the set of orthogonal joint surfaces obtained by the source mechanism inversion, the fault instability coefficient is used to determine the true crack plane; Project any joint surface solution onto the normalized Mohr circle, draw a straight line parallel to the fracture criterion line at the projection point, and the vertical distance from the position of the normalized maximum normal stress to the line is the fault instability coefficient I of the joint surface solution. The inversion of stress field is based on the Wallace-B assumption that the slip direction of fault plane is parallel to the direction of maximum shear stress along the fault plane. According to the Wallace-B assumption, the stress field inversion model is constructed, and the relationship between the slip vector s of the fracture and the normalized stress vector t (t = [τ 11 τ 12 τ 13 τ 22 τ 23 ] ) is At = s, wherein A is a 3x5 matrix, and the matrix elements are calculated from the normal vector of the fault plane, which is expressed as: wherein n 1、 n 2、 n 3 is a component of a fault plane normal vector. When there are multiple composite focal mechanism solutions, the stress moment tensor is obtained by the least square inversion method; the slip vectors s and matrix A corresponding to multiple focal mechanism solutions are integrated to construct an over-determined equation group, and the least square method is used to solve the equation group to obtain the optimal stress moment tensor; The stress moment tensor obtained by the least square inversion is subjected to eigenvalue decomposition, the azimuth and inclination of multiple principal stresses are calculated by the eigenvalue decomposition method, the size order of the principal stresses is determined according to the size of the eigenvalue, and the direction of the characteristic vector corresponding to the eigenvalue is the direction of the principal stress, so as to determine the state of the stress field; The stress field results obtained by inversion are applied to the detection of Ordovician limestone water pressure and top interface fracture in a mine, and the internal relationship between the stress field and the change of Ordovician limestone water pressure and the development of fracture is analyzed.

8. The method of any one of claims 1 to 7, characterized in that: In the use of cluster analysis for advanced attribute calculation, cluster analysis is divided into static and dynamic two modes, and in the microseismic monitoring process, the calculation mode is determined according to the occurrence of microseismic events: if there are microseismic events that have not yet occurred, that is, the monitoring process is still in progress and new events may occur, the dynamic parameter calculation mode is adopted to describe the rupture process of the stratum; when all microseismic events have occurred, the static parameter calculation mode is switched; wherein, according to the dynamic parameter calculation mode, the following dynamic parameters are calculated: Coseismic deformation and strain rate: by analyzing the waveform data generated by the microseismic event, combining the mechanical properties of the stratum medium and the geometric characteristics of the monitoring area, the coseismic deformation of the stratum and the corresponding strain rate when each microseismic event occurs are calculated to reflect the instantaneous deformation of the stratum under the action of microseismic; Event time, distance interval and event density: record the occurrence time of each microseismic event, and calculate the time interval between events; at the same time, according to the focal location results, the distance interval between events is calculated. Based on the time and spatial distribution information of the events, the event density in unit time and unit space range is calculated, which is used to analyze the temporal and spatial evolution characteristics of microseismic activity; Radiated energy and seismic moment: using microseismic waveform information combined with seismological formula, the energy radiated by each microseismic event and the corresponding seismic moment are calculated to quantify the intensity and scale of the microseismic event.

9. The method of claim 8, characterized in that: Based on the dynamic parameter calculation, important dynamic cluster attributes are obtained, which are stress index, plasticity index and diffusion index respectively: The calculation formula of stress index SI is: Wherein, ηc is the cumulative seismic efficiency of cluster events, c is the coefficient of different lithology, time and volume, and DI is the diffusion index; when the SI value is low, it indicates that the stress conduction is unstable and limited to a very small range of the reservoir; when the SI value is large, it represents that the microseismic event releases a large amount of energy, the stress field is more stable, and the final fracture range is larger; The calculation formula of plasticity index PI is: Wherein, μ is the dynamic shear modulus, and ηc is the cumulative seismic efficiency of cluster events; when the PI value is low, it represents that the reservoir is not easy to deform; when the PI value is large, it represents that the reservoir is extremely easy to deform; The calculation formula of diffusion index DI is: where X 2 is the square of the average distance of the clustered events, and t is the average time interval of the clustered events; a lower value of the diffusion index DI represents a large time interval of event occurrence but a small distance interval; a larger value of DI represents a small time interval of event occurrence but a large distance interval. The stress index, plastic index, and diffusion index are comprehensively analyzed, and a multi-factor relationship model between the indices and the pressure of Ordovician limestone water and the development of the top boundary fracture is established.

10. The method of claim 8, wherein: a static clustering attribute is calculated based on static parameters, which are microseismic radiated energy, apparent stress and apparent volume, respectively; Microseismic radiation energy: for each microseismic event, the displacement spectrum is integrated; using the microseismic waveform data recorded by the microseismic monitoring system, the displacement spectrum information is obtained through spectrum analysis technology, and the energy released by the microseismic event is calculated according to the following formula : wherein, p is the formation medium density; is the seismic wave velocity; R is the source-to-monitor distance; F is the attenuation factor; Ω is the displacement spectrum function, f is the frequency; apparent stress: shear modulus of reservoir rock m and seismic moment of each microseismic event M 0, apparent stress is calculated from the apparent stress formula : wherein, is the shear modulus of the reservoir rock, is the seismic moment of the microseismic event; apparent volume characterizing the volume of the rock deformed inelastically by the source of the microseismic event E , the shear modulus of the reservoir rock m and the source parameters of the microseismic event into the following apparent volume calculation formula to calculate the apparent volume characterizing the volume of the rock deformed inelastically by the source of the microseismic event : 。

Images