Geological stress concentration area intelligent identification method

Abstract

The invention relates to the technical field of intelligent identification, and discloses an intelligent identification method for a geological stress concentration area. The method comprises the steps that a plurality of vector electromagnetic induction coils are arranged in a monitoring area in a latticed mode, and original electromagnetic pulse signals containing components in the X direction, the Y direction and the Z direction are collected; performing denoising and time-frequency analysis on the original electromagnetic pulse signal to obtain a multi-dimensional data cube; processing the multi-dimensional data cube by using an integer linear programming model and an adaptive fuzzy correction algorithm to obtain a corrected electromagnetic pulse signal data set; performing decomposition and mapping relation construction on the corrected electromagnetic pulse signal data set to obtain a geological stress distribution estimation value; and constructing a geological stress concentration area three-dimensional model according to the geological stress distribution estimation value. According to the method, the stress field boundary identification precision is enhanced, and the blind area is reduced, so that high-precision stress field reconstruction under complex geological conditions is realized.

Description

Technical Field

[0001] The present application relates to the field of intelligent identification technology, and in particular to a method for intelligent identification of geological stress concentration areas. Background Art

[0002] Geological stress concentration areas are high-risk areas for geological disasters, and their accurate identification is of great significance to mine safety, earthquake early warning, and geological disaster prevention and control. Traditional methods for identifying geological stress concentration areas mainly rely on borehole stress measurement, acoustic emission monitoring, or microseismic monitoring technology. These methods not only require a lot of field work and high equipment costs, but also have strong destructiveness and point limitations, making it difficult to achieve large-scale, real-time, and continuous geological stress field monitoring. With the development of geophysical exploration technology, methods based on detecting natural fields have gradually attracted attention. Among them, natural electromagnetic pulse signals have become one of the important means of identifying geological stress concentration areas due to their advantages of non-destructiveness, real-time, and large-scale monitoring.

[0003] Natural electromagnetic pulse signals are electromagnetic phenomena caused by changes in internal rock stress. The stress state of the geological body can be inferred by analyzing these signals. However, existing natural electromagnetic pulse signal monitoring technologies face many challenges, including weak signal strength, severe noise interference, difficulty in locating the signal source, and complex data interpretation. Especially in complex geological environments, traditional signal processing methods have difficulty in effectively distinguishing useful signals from interference signals, resulting in insufficient accuracy in the discrimination results. At the same time, most existing technologies use empirical formulas or simple statistical models for data analysis, and lack a systematic method for accurately identifying geological stress concentration areas. Summary of the invention

[0004] The present application provides a method for intelligently identifying geological stress concentration areas, which enhances the accuracy of stress field boundary recognition and reduces the blind area, thereby achieving high-precision stress field reconstruction under complex geological conditions.

[0005] In a first aspect, the present application provides a method for intelligently identifying geological stress concentration areas, the method comprising: Multiple vector electromagnetic induction coils are arranged in a grid shape in the monitoring area to collect original electromagnetic pulse signals containing three-directional components: X, Y, and Z. Performing denoising and time-frequency analysis on the original electromagnetic pulse signal to obtain a multidimensional data cube; Processing the multidimensional data cube using an integer linear programming model and an adaptive fuzzy correction algorithm to obtain a corrected electromagnetic pulse signal data set; Decomposing and mapping the corrected electromagnetic pulse signal data set to obtain an estimated value of geological stress distribution; A three-dimensional model of the geological stress concentration area is constructed according to the estimated value of the geological stress distribution.

[0006] In the technical solution provided by this application, the spatial resolution is improved through the grid-distributed multi-vector electromagnetic induction coil monitoring network, and the equilateral triangle arrangement mode and orthogonal three-component receiving method are adopted to realize the all-round three-dimensional monitoring of natural electromagnetic pulse signals, significantly enhance the accuracy of stress field boundary recognition, and reduce the blind area. The systematic signal preprocessing process includes power frequency filtering, wavelet denoising and adaptive bandpass filtering, combined with the time-frequency analysis of Hilbert-Huang transform, effectively extracts the effective components in weak electromagnetic pulse signals, and suppresses environmental interference and system noise. The method combining integer linear programming model with adaptive fuzzy correction algorithm has solved the signal baseline drift and electromagnetic signal ambiguity problems that have plagued this field for a long time. The multi-synchronous compression transform with time reallocation has broken through the limitations of traditional Fourier analysis and can effectively process non-stationary signals. The mapping relationship between signal characteristics and geological stress parameters constructed by the extreme gradient boosting algorithm has achieved high-precision stress field reconstruction under complex geological conditions. The stress concentration index evaluation system analyzes the correlation between stress field characteristics and geological structures, and performs spatial constraint optimization through the Markov random field model, which greatly reduces the false anomaly recognition rate. The three-dimensional model of the geological stress concentration area uses adaptive grid subdivision technology to finely depict the spatial geometric characteristics and internal stress distribution of the stress concentration area. BRIEF DESCRIPTION OF THE DRAWINGS

[0007] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings required for use in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are some embodiments of the present invention. For ordinary technicians in this field, other accompanying drawings can be obtained based on these accompanying drawings without paying creative work.

[0008] Figure 1 This is a schematic diagram of an embodiment of the method for intelligently identifying geological stress concentration areas in the embodiments of the present application. DETAILED DESCRIPTION

[0009] The embodiment of the present application provides a method for intelligently identifying geological stress concentration areas. The terms "first", "second", "third", "fourth", etc. (if any) in the specification and claims of the present application and the above-mentioned drawings are used to distinguish similar objects, and are not necessarily used to describe a specific order or sequence. It should be understood that the data used in this way can be interchangeable where appropriate, so that the embodiments described here can be implemented in an order other than that illustrated or described here. In addition, the terms "including" or "having" and any variations thereof are intended to cover non-exclusive inclusions, for example, a process, method, system, product or device that includes a series of steps or units is not necessarily limited to those steps or units that are clearly listed, but may include other steps or units that are not clearly listed or inherent to these processes, methods, products or devices.

[0010] For ease of understanding, the specific process of the embodiment of the present application is described below. Figure 1 , an embodiment of the method for intelligently identifying geological stress concentration areas in the embodiments of the present application includes: Step S101, arranging a plurality of vector electromagnetic induction coils in a grid pattern in a monitoring area to collect original electromagnetic pulse signals containing components in three directions of X, Y, and Z; It is understandable that the execution subject of the present application may be a geological stress concentration zone intelligent identification system, or a terminal or a server, which is not specifically limited here. The present application embodiment is described by taking a server as the execution subject as an example.

[0011] Specifically, the geological structure complexity analysis of the monitoring area is carried out to obtain key geological parameters such as the distribution characteristics, fault density, joint development degree and lithology differences of various structural units in the monitoring area. The analysis content includes a comprehensive analysis of existing geological maps, remote sensing image data and drilling results, supplemented by regional stress field simulation and geometric modeling of the fracture system, to obtain an analysis result reflecting the level of geological structure complexity. Based on the analysis results, it is determined which areas in the monitoring area have active geological activities and dense structural intersections, and which areas have single structures and slow changes. According to the complexity level of different areas, an adaptive strategy is used to determine the layout spacing of the receiving devices, and an equilateral triangle arrangement mode is uniformly adopted as the basic layout form. For areas with high structural complexity, the spacing between the receiving devices is set to 50 meters to improve the spatial resolution; while for areas with relatively simple structures, the spacing is relaxed to 200 meters to save resources and reduce the layout cost. Through the differentiated layout strategy, the layout plan of the receiving device is obtained. Each receiving device in the receiving device layout plan is equipped with an electromagnetic signal receiving unit containing three orthogonal vector electromagnetic induction coils. The coils arranged in three axes correspond to the three spatial directions of X, Y, and Z respectively, and the electromagnetic pulse vector components in the three directions are obtained at the same time to ensure the integrity and decoupling of the data in the spatial dimension. Each induction coil is connected to a highly sensitive fluxgate sensor with a sensor range of ±100nT to meet the observation requirements of weak natural electromagnetic signals and large frequency spans. At the same time, each electromagnetic signal receiving unit integrates a high-precision GPS clock module, which realizes unified time calibration between receiving devices, so that each site can synchronously collect data at the millisecond level, thereby ensuring that the responses of different points to the same electromagnetic event have time consistency. In order to improve the network flexibility and deployment efficiency of the overall system, each receiving device is embedded with a low-power wireless transmission module, through which a communication link is established with the surrounding devices to build a mesh network with self-organization and self-repair capabilities. This distributed mesh network can automatically find a way around when any node in the network fails, ensuring the continuity and robustness of data transmission. All electromagnetic signal receiving units form a distributed monitoring network with time synchronization and spatial coordination capabilities. While completing the layout of the main monitoring grid, several reference sites are selected outside the monitoring area. These sites are far away from the main geological tectonic activity zone. Their main function is to record electromagnetic interference signals outside the area, such as lightning activities, industrial electromagnetic waves, communication signals, etc. The interference signal is sometimes similar to the natural geological electromagnetic pulse in waveform, frequency and amplitude, which can easily cause misidentification. Therefore, it is independently recorded and feature extracted through the reference site. The reference interference data obtained will serve as an important reference standard for subsequent signal preprocessing. The entire distributed monitoring network collects electromagnetic pulse signals in real time 24 hours a day, and realizes multi-site synchronous recording through the GPS module.After the data is transmitted to the central processing system, the system automatically compares and analyzes the reference interference data with the electromagnetic signals received at each measuring point, and uses algorithms such as waveform matching, frequency band characteristics and phase difference discrimination to eliminate interference signals from non-geological sources, judge the validity of the remaining signals and mark them, and finally output the original electromagnetic pulse signal containing X, Y, and Z direction components.

[0012] Step S102, performing denoising and time-frequency analysis on the original electromagnetic pulse signal to obtain a multidimensional data cube; Specifically, the original electromagnetic pulse signal is filtered for interference. Since there are a lot of background interferences caused by power systems, radio communications, and lightning activities in the actual environment, which affect the authenticity and identifiability of natural electromagnetic signals, the power frequency interference in the signal is suppressed by using a notch filter, especially 50Hz and its higher harmonic components. By designing the center frequency and bandwidth of the notch filter, the power system noise is effectively removed to obtain a preliminary filtered signal. Then the wavelet threshold denoising method is applied to the preliminary filtered signal to eliminate random noise and low-frequency drift. The db4 wavelet function in the Daubechies wavelet family is selected as the mother wavelet, and the signal is decomposed into 5 layers of wavelets. The detail coefficients are compressed by soft threshold processing to achieve effective noise removal and retention of useful signal features. After this processing process, a denoised signal that is closer to the true response of the geological event is output. The denoised signal is segmented in the time domain to enhance the resolution of time-frequency characteristics and adapt to the non-stationary characteristics of geological events. Each signal is divided into time domain segments of 60 seconds in length, and a 50% overlap rate is set between adjacent segments to achieve a good balance between time resolution and frequency resolution. The Hilbert-Huang transform is used to perform time-frequency analysis on each signal segment. This method decomposes the signal into several intrinsic mode functions through empirical mode decomposition, and then performs Hilbert transform on each mode function to extract its instantaneous frequency, instantaneous amplitude and phase information. In this way, the dynamic characteristics of frequency changes over time in natural electromagnetic pulses are effectively captured. The analysis results of all segmented signals are summarized to form a high-dimensional signal feature matrix, where each row corresponds to the analysis result within a time window, and each column represents a specific signal parameter indicator, such as instantaneous frequency or instantaneous energy. After obtaining the signal feature matrix, polarization analysis is performed on the three-component vector electromagnetic signal collected at each receiving point. By calculating the polarization ellipse parameters of the signal in each time period, the three key polarization parameters of its ellipticity, inclination and azimuth are extracted. These parameters reflect the spatial propagation characteristics and source mechanism of electromagnetic waves, which help to identify the geological source of the signal. At the same time, in order to remove potential interference signals from non-geological sources, a set of adaptive notch filters is constructed in combination with the regional external electromagnetic interference data recorded by the previous reference station. These filters dynamically attenuate the components in the signal spectrum that are consistent with the characteristic frequency of the interference source, and remove non-natural pulse interference such as lightning, electricity, and industry to the greatest extent without destroying the geological effective signal structure, so as to obtain pure geological source electromagnetic signals. The pure geological source signal is subjected to spatial correlation analysis. By calculating the signal correlation coefficients between each receiving device and the phase difference in the corresponding time window, the spatial propagation consistency and correlation of the signal are evaluated. High correlation and stable phase difference indicate that the signal originates from the same geological stress release event.All processed data are organized in a unified standard format to form a multidimensional data cube. Each data unit of the cube contains the acquisition time, measurement point location, amplitude, frequency and phase information in the three directions of X / Y / Z, constructing a standardized input data body that describes the spatiotemporal evolution characteristics of the geological electromagnetic field.

[0013] Step S103, using an integer linear programming model and an adaptive fuzzy correction algorithm to process the multidimensional data cube to obtain a corrected electromagnetic pulse signal data set; Specifically, a multidimensional data cube is used as input to construct and solve an integer linear programming model with sparsity constraints. The core goal of this model is to minimize the weighted residual sum of squares between the estimated value and the actual observed value of the natural electromagnetic pulse signal, so as to restore the true geological characteristics of the signal as much as possible. The introduction of the L1 norm sparsity constraint term in the optimization objective function helps to avoid the occurrence of overfitting, while enhancing the physical interpretability of the solution, ensuring that the final extracted signal is mainly concentrated on the components that are highly responsive to geological stress changes. Through this process, a set of preliminary optimized signal data is obtained. After obtaining the preliminary optimized signal data, the singular value decomposition method is used to perform feature analysis on the system matrix involved in the integer linear programming model, and the condition number of the system matrix is ​​calculated by decomposing its singular value spectrum. The condition number reflects the degree of pathological condition of the system matrix. If its value is significantly greater than 1000, it means that there is a problem of numerical instability or excessive stretching of the solution space during the system solution process. Therefore, the Tikhonov regularization method is introduced to correct it. Tikhonov regularization effectively reduces the volatility of the solution and improves its robustness by adding a penalty term related to the signal norm to the optimization objective function. In view of the baseline drift problem that is common in the actual measurement of electromagnetic signals, baseline correction is performed on the stabilized signal data. The signal baseline is modeled using a piecewise polynomial function, and the entire signal sequence is divided into several time periods. Each period is fitted with a polynomial curve to describe the background offset caused by non-geological factors. The polynomial coefficients are solved by minimizing the deviation between each signal segment and its fitting curve, thereby restoring the signal baseline. After baseline correction, the signal waveform can return to the real fluctuations, eliminating the influence of low-frequency trend terms caused by non-structural factors such as environmental disturbances and equipment drift. An adaptive fuzzy correction algorithm is introduced to eliminate the ambiguity of the corrected signal. The algorithm adopts the Mamdani fuzzy inference model built based on expert knowledge, and uses Gaussian membership functions to classify the input signal features. The input variables of the fuzzy system include the time scale, spatial correlation and environmental parameters of the signal during observation, while the output variable is the fuzzy correction coefficient for each signal segment. By reasoning based on the preset If-Then rules in the fuzzy rule library, the fuzzy control system adaptively generates the optimal correction strategy based on different input combinations, thereby classifying and clarifying the fuzzy information in the signal, and finally outputting the signal data with the ambiguity eliminated. In order to solve the problem of electromagnetic signal positioning error caused by complex transmission paths or ambiguous sources, the wave arrival time difference constraint condition is introduced on the basis of eliminating the ambiguity signal, and an overdetermined set of observation equations is constructed, which is solved by iterative reweighted least squares method. This method gradually converges to the global optimal solution by dynamically updating the weights of the signal components in each round, effectively dealing with the nonlinear characteristics in signal propagation. During the iteration process, the solution is considered to be the optimal solution when the Euclidean distance change is less than the set threshold (such as ) or reaching the maximum number of iterations (such as 100) as the stopping condition to ensure that the calculation results achieve the best balance between accuracy and convergence speed. A set of electromagnetic signal samples with high spatiotemporal consistency and ambiguity elimination is obtained. The quality of each signal segment is evaluated. By calculating its signal-to-noise ratio, amplitude stability, phase coherence and other indicators, the integrity, effectiveness and availability of the signal are quantitatively evaluated, and the signal segments that meet the requirements of engineering applications are screened out accordingly. All high-quality signal data that pass the screening are archived uniformly to form a corrected electromagnetic pulse signal data set.

[0014] Step S104, decomposing and mapping the corrected electromagnetic pulse signal data set to obtain an estimated value of geological stress distribution; Specifically, the high-quality electromagnetic pulse signal data is input into a time-redistributed multi-synchronous compression transformation algorithm, which dynamically analyzes the local frequency characteristics of the signal and adaptively adjusts the window length to achieve fine decomposition of the signal. In order to ensure the optimal time-frequency resolution in different frequency bands, the algorithm introduces an adaptive window length function, so that a longer window function is applied to the signal segment with slow frequency changes to improve the frequency resolution, and a shorter window function is used for the segment with rapid frequency changes to improve the time resolution. After the transformation process, the corrected electromagnetic pulse signal is decomposed into several single-component functions with narrow frequency bands, each of which represents a specific frequency component and has a clearer physical characteristic directionality. From the narrow-band single-component function, a set of high-dimensional signal feature parameter sets related to geological stress activities are extracted, which include the center frequency, bandwidth, distribution characteristics of signal energy in the time domain and frequency domain, duration, and instantaneous frequency change rate of each single component. These parameters comprehensively reflect the spectral structure, energy evolution, and dynamic propagation behavior of electromagnetic pulses. For the initial electromagnetic pulse wavefront information in the corrected electromagnetic pulse signal data set, the signal response at the corresponding time between multiple receiving points is selected as the basis, and the difference in these response times is used to measure the arrival time difference to form a set of arrival time difference data, which reflects the path characteristics of the electromagnetic pulse from the source point to each measuring point. Combined with the known spatial distribution coordinate information of the receiving device, the arrival time difference data and the spatial distance matrix are jointly modeled to construct an overdetermined geological stress location equation, whose basic expression is V·T=D, where V represents the velocity tensor of electromagnetic waves propagating in geological media, T is the arrival time difference matrix, and D is the spatial distance matrix between measuring points. Since this equation is a typical multivariate overdetermined problem, the singular value decomposition method is used for robust solution to obtain a spatially anisotropic propagation velocity tensor, which reflects the propagation ability of electromagnetic waves in different directions and indirectly reveals the elastic properties and stress distribution characteristics of the underground medium. After completing the propagation modeling, in order to establish a quantitative mapping relationship between signal characteristics and geological stress parameters, the extreme gradient boosting algorithm, namely the XGBoost model, is introduced to train and analyze the signal feature parameter set. As an enhanced tree model, XGBoost can efficiently capture nonlinear mapping relationships while taking into account overfitting control and feature selection capabilities. In this step, the known geological stress data in the historical monitoring data is used as a label to build a supervised learning model. The input features are the previously extracted spectrum, energy, and time-frequency dynamic parameters, and the model output is the stress size, principal stress direction, and stress gradient value of the target area. During model training, the number of decision trees, depth, learning rate and other hyperparameters are initially set, and the complex relationship between data is gradually fitted through the standard gradient boosting mechanism.In order to improve the generalization ability and prediction accuracy of the model, the importance of different input features in the model is sorted and screened. The feature importance analysis method is used in the screening process to identify the top M key features with the greatest contribution in the stress prediction task. This process effectively eliminates redundant features, reduces model complexity, and improves training efficiency. The Bayesian optimization method is combined to automatically adjust the hyperparameters of the XGBoost model, including the regularization term coefficient, the maximum depth of the tree, the learning rate, the minimum sample weight, etc., to maximize the performance on the validation set. The optimization process searches for the optimal solution point in the high-dimensional parameter space through the Bayesian framework to avoid the inefficiency of traditional grid search. The geological stress prediction model constructed by the above steps can accurately map the electromagnetic signals in any monitoring area and output structured geological stress distribution estimates, including the direction vector of the principal stress, the absolute value of the stress intensity, and the stress gradient distribution between different regions.

[0015] Fast Fourier transform is performed on the corrected electromagnetic pulse signal data set. The fast Fourier transform method is used to process the signal in the entire frequency band. With its advantages in computational efficiency and frequency resolution, a series of initial spectral characteristic parameters such as energy distribution, main frequency peak, and frequency band extension of the electromagnetic signal in the frequency domain are quickly obtained. The local eigenvalues ​​of the Hessian matrix are calculated for each frequency band based on the initial spectral characteristic parameters of the signal. The frequency change curvature model of the signal is constructed using the local Hessian matrix, and the main eigenvalues ​​of the Hessian matrix corresponding to each frequency point are calculated to obtain the intensity of the signal change in the frequency band where the point is located. The more drastic the frequency change interval, the larger the eigenvalue of the Hessian matrix, reflecting that the signal in this frequency band has a strong time-varying property, and a higher time resolution should be given in subsequent analysis. Based on this analysis result, the window length scaling factor is constructed in combination with the frequency change rate information. By designing an adaptive window length function, the time window length of different frequency intervals of the signal is dynamically allocated, so that the low-frequency stable segment obtains a longer window function to enhance the frequency resolution, while the high-frequency mutation segment is allocated a shorter window function to improve the time resolution, forming an adaptive window length function system that automatically adjusts with frequency changes. The short-time Fourier transform is used as the basic framework of time-frequency analysis, and the adaptive window length function is embedded in the short-time Fourier transform process. For each frequency band of the signal, the corresponding window function length is used for piecewise convolution calculation to generate a time-frequency distribution map with high-resolution time-frequency structure. The time-frequency distribution map is subjected to time redistribution operation. According to the instantaneous frequency and group velocity information of the signal, the energy distributed in the fuzzy area is refocused to its real time-frequency position, so as to improve the energy concentration and structural clarity of the time-frequency map, and form an enhanced time-frequency distribution feature map with stronger readability. The enhanced time-frequency distribution feature is input into the synchronous compression kernel function to perform ridge extraction, and the synchronous depth control strategy is applied to stabilize the ridge extraction process. The synchronous compression kernel function is used to identify the energy peak coherent trajectory and extract the time-frequency ridge sequence representing the main path of the signal. These ridges correspond to the most stable and geophysical frequency components of the signal in the time-frequency map. The synchronous depth control adjusts the threshold parameters in the extraction process according to the local consistency index of the signal at multiple scales to ensure that each extracted ridge has good energy coherence and physical consistency. Through this process, multiple original single-component functions are obtained, each of which corresponds to a signal sub-component with clear physical meaning and stable frequency characteristics, reflecting the independent electromagnetic response path generated during a specific geological activity. Energy contribution evaluation is performed on the original single-component function. This evaluation calculates the proportion of each single component in the total signal energy, and selects the parts with significant energy proportion and stable structure as the final analysis object. Only signal components with energy contribution exceeding the preset threshold are retained, and low-energy, dispersed structure, and difficult-to-interpret secondary components are eliminated, and finally multiple narrow-band single-component functions with high energy integration, spectral stability and time domain resolution are obtained.

[0016] Step S105: construct a three-dimensional model of the geological stress concentration area according to the estimated value of the geological stress distribution.

[0017] Specifically, the estimated value of geological stress distribution is used as the basis for modeling, and the Kriging interpolation method in spatial statistics is used to reconstruct the continuous field of stress parameters of all discrete sampling points. As the optimal unbiased estimation method in geological space interpolation, Kriging interpolation takes into account the spatial distance between measuring points and introduces the semivariance function to characterize the spatial correlation between measuring points. It shows a high degree of fitting accuracy when dealing with physical quantities with spatial heterogeneity such as geological stress. Through this method, the stress estimates distributed at different depths and locations are uniformly interpolated into a structurally coherent and continuously distributed stress field data to form a spatial stress tensor field expressed in the form of a three-dimensional grid, in which each grid unit contains the magnitude and direction information of the principal stress. Tensor analysis is carried out on the basis of the obtained stress field data to extract the direction and magnitude of the principal stress in each spatial unit, and calculate the three principal stress values. , and and its corresponding direction cosines are used to construct a complete principal stress ellipsoid expression. By statistically analyzing the principal stress distribution characteristics of the entire region, areas with high stress gradients and large principal stress deviations are identified, and these areas are candidates for the spatial distribution of potential stress concentration areas. In order to achieve intelligent discrimination, a set of clear criteria for distinguishing stress concentration areas is constructed in combination with the average stress state of the region. For example, if the maximum principal stress in a certain region is More than 1.5 times the average principal stress of the entire area, and the maximum and minimum principal stress difference If the difference is more than twice the regional average difference, the area is marked as a preliminarily identified stress concentration area. The Markov random field model is introduced to perform spatial optimization on the preliminary stress concentration area. The model establishes a statistical association between regional stress data and neighborhood relationships by setting local spatial adjacency and state transition probability, so that the optimization results take into account both data-driven and spatial structural laws. The model regards each spatial grid unit as a state node, and determines its final affiliation by comparing the stress state with its adjacent nodes, thereby eliminating isolated pseudo-outliers caused by abnormal sampling or interpolation errors in the identification results, and enhancing the spatial geometric coherence and geological structural consistency of the concentration area. After the optimization is completed, a stress concentration area with more reasonable space, clear boundaries and continuous structure is obtained. Based on the optimized stress concentration area, a stress concentration index (SCI) evaluation system is established to achieve quantitative grading of stress concentration levels in different regions. The SCI index comprehensively considers three factors: the relative intensity of the maximum principal stress, the principal stress difference and the fault distance. Its expression is defined as: ,in and Represent the average values ​​of the maximum principal stress and stress difference in the region, respectively, and D is the normalized distance index between the center point of the unit and the nearest active fault. The SCI value of each grid unit is calculated by this formula, and then the stress concentration area is divided into three levels: level I (SCI>4.0), level II (2.0≤SCI≤4.0) and level III (1.0≤SCI<2.0) according to the set threshold, so as to achieve the quantitative expression and classification of stress risk. Based on the classified stress concentration area data, three-dimensional modeling is carried out to integrate the spatial distribution of geological stress with the geological background information, and a comprehensive three-dimensional model of geological stress concentration area containing multiple types of information elements is constructed. The model is based on the stress tensor field, supplemented by lithology distribution, fault development degree and groundwater distribution data as geological structure constraints. In order to improve the accuracy and adaptability of the model, the adaptive meshing technology is used to automatically refine the mesh in areas with large stress gradients to more accurately describe the complex geological response characteristics. The model is presented as a multi-layer nested structure in three-dimensional space, and is visualized in various forms such as contour maps, slice maps, and vector field maps, fully presenting the internal structure of the stress concentration area, the stress evolution trend, and its coupling relationship with the surrounding tectonic system.

[0018] In the embodiment of the present application, the spatial resolution is improved by a multi-vector electromagnetic induction coil monitoring network with a grid distribution, and an equilateral triangle arrangement mode and an orthogonal three-component receiving method are adopted to achieve all-round three-dimensional monitoring of natural electromagnetic pulse signals, significantly enhance the accuracy of stress field boundary recognition, and reduce the blind area. The systematic signal preprocessing process includes power frequency filtering, wavelet denoising and adaptive bandpass filtering, combined with the time-frequency analysis of the Hilbert-Huang transform, effectively extracting the effective components in the weak electromagnetic pulse signal and suppressing environmental interference and system noise. The method combining integer linear programming model with adaptive fuzzy correction algorithm has solved the signal baseline drift and electromagnetic signal ambiguity problems that have plagued this field for a long time. The multi-synchronous compression transform with time reallocation has broken through the limitations of traditional Fourier analysis and can effectively process non-stationary signals. The mapping relationship between signal characteristics and geological stress parameters constructed by the extreme gradient boosting algorithm has achieved high-precision stress field reconstruction under complex geological conditions. The stress concentration index evaluation system analyzes the correlation between stress field characteristics and geological structures, and performs spatial constraint optimization through the Markov random field model, which greatly reduces the false anomaly recognition rate. The three-dimensional model of the geological stress concentration area uses adaptive grid subdivision technology to finely depict the spatial geometric characteristics and internal stress distribution of the stress concentration area.

[0019] In a specific embodiment, the process of executing step S101 may specifically include the following steps: Conducting geological structure complexity analysis on the monitoring area to obtain geological structure complexity analysis results, and determining an equilateral triangle arrangement mode based on the geological structure complexity analysis results to obtain a receiving device layout plan; Each receiving device in the receiving device arrangement scheme is equipped with three orthogonal vector electromagnetic induction coils to obtain an electromagnetic signal receiving unit; A GPS clock module is integrated into each electromagnetic signal receiving unit to obtain a time synchronization acquisition system, and the electromagnetic signal receiving units are connected to form a mesh network according to the time synchronization acquisition system and the low-power wireless transmission module to obtain a distributed monitoring network; Select reference sites in the monitoring area, and record the external electromagnetic interference signals in the area according to the reference sites to obtain reference interference data; The electromagnetic pulse signal of the distributed monitoring network is monitored, and the effective signal is marked in combination with the reference interference data to obtain the original electromagnetic pulse signal containing the X, Y, and Z direction components.

[0020] Specifically, the geological structure complexity analysis is carried out in the monitoring area to systematically identify the structural units, fault systems, lithology combinations and their distribution patterns in the area. The analysis process comprehensively utilizes geological maps, remote sensing images, drilling results and historical seismic data. By constructing a geological structure complexity evaluation model, the key factors such as joint density, fault interpenetration, number of structural intersections, number of stratigraphic unconformity interfaces in the region are quantitatively evaluated. Combined with the frequency of seismic activity, crustal movement speed and changes in geomorphic units, multi-factor weighted calculations are performed to obtain the geological structure complexity analysis results expressed in structural complexity levels. The analysis results are divided into multiple regional units according to space, and each unit is given a structural complexity level label and output in the form of a layer. According to the results of the geological structure complexity analysis, a layout strategy for the receiving device is formulated, and an equilateral triangle grid is used as the basic arrangement mode to ensure that the signal acquisition capability is balanced in all directions in space. In areas with high structural complexity, such as areas with dense fault intersections, dense folds or rock mutations, a higher density of receiving nodes is set, and the length of the triangle side is controlled between 50 meters and 100 meters to improve the local spatial resolution; in areas with relatively stable or single structures, the layout density is appropriately reduced, and the side length is expanded to 150 meters to 200 meters to reduce system resource consumption and simplify the layout engineering. Through the layout strategy driven by structural complexity, a receiving device layout plan is formed that has both macro coverage capabilities and local high-resolution acquisition. At each layout point, an electromagnetic signal receiving unit is installed according to a unified standard. Three orthogonally arranged vector electromagnetic induction coils are fixed in each receiving device, corresponding to the three spatial directions of X, Y, and Z, respectively, for synchronously collecting electromagnetic field change data in each direction. The inductor coils are all connected to a high-sensitivity fluxgate sensor, whose range is set to ±100nT and the sensitivity reaches 0.1pT / √Hz, to ensure that the extremely low amplitude electromagnetic signals generated during the release of natural geological stress can be captured. Magnetic shielding and non-magnetic shell material design are introduced into the inductive acquisition structure to prevent the device structure itself from interfering with or shielding the signal. At the same time, environmental sensors such as temperature, humidity, and vibration are equipped to record and collect environmental parameters in real time as important reference information for subsequent signal processing and interference identification. In order to achieve time synchronization between each receiving unit, a high-precision GPS clock module is integrated in each electromagnetic signal receiving unit. This module supports millisecond-level synchronization accuracy and provides a unified time reference for electromagnetic signal acquisition by receiving the standard time signal of the global navigation satellite system, ensuring that all devices have timestamp consistency when collecting the same natural electromagnetic event. After synchronization, each receiving unit compares the wave arrival time and inverts the propagation path of the electromagnetic signals from the same source point to support subsequent signal positioning and stress field modeling. On the basis of achieving time synchronization, a low-power wireless transmission module is added to each receiving device, and a mesh network architecture is constructed through a multi-hop communication mechanism to achieve automatic interconnection between nodes and coordinated data forwarding.The mesh network supports dynamic path reconstruction and self-repair mechanisms between devices. When a node fails or communication is interrupted, other nodes automatically rebuild the transmission path to ensure the continuous and stable operation of the entire network. Through this distributed mesh network, all receiving units upload the collected data to the central data processing server in real time to realize data centralized processing, remote monitoring and fault diagnosis functions, and build a data collection system with full coverage, high reliability and strong coordination. After the distributed collection system is deployed, several reference sites are selected at the outer edge of the monitoring area. These reference points try to avoid the main geological tectonic activity belt and strong artificial interference sources. Their main task is to record non-geological electromagnetic interference signals from outside the area or background sources. The data collected by the reference site is used to build a background interference library, including typical interference modes such as lightning activity, electromagnetic leakage of power lines, and radio frequency radiation. Its time, frequency, and waveform characteristics will be used as comparison templates to assist in the subsequent identification of effective geological electromagnetic events. After the entire distributed monitoring network is put into operation, each receiving unit will conduct uninterrupted real-time monitoring of the electromagnetic pulse signals that continue to appear in the natural environment. All collected signals are calibrated with GPS timestamps and then uploaded to the central server. The system automatically compares and analyzes them with the background interference data of the reference site. Multi-dimensional feature extraction and comparison methods such as frequency domain matching, waveform similarity, polarization parameters and wave arrival time are used to accurately identify and mark valid signals with clear sources and geological characteristics. The selected signals are the original electromagnetic pulse signals containing the three directional components of X, Y, and Z.

[0021] In a specific embodiment, the process of executing step S102 may specifically include the following steps: The original electromagnetic pulse signal is subjected to interference filtering to obtain a preliminary filtered signal, and the preliminary filtered signal is subjected to wavelet threshold denoising to obtain a denoised signal; Perform time domain segmentation processing on the denoised signal to obtain a segmented time domain signal, and perform time-frequency analysis on the segmented time domain signal to obtain a signal feature matrix; The polarization parameters including ellipticity, inclination and azimuth are calculated based on the three-component vector signal in the signal characteristic matrix, and the external interference signal recorded by the reference station is combined with adaptive notch filtering to obtain a pure geological source signal. The pure geological source signals are subjected to spatial correlation analysis, the mutual correlation coefficients and phase differences between adjacent measuring points are calculated, and they are organized in a unified format to form a multidimensional data cube containing time, space, amplitude, frequency, and phase.

[0022] Specifically, the three-component original electromagnetic pulse signals collected at each measuring point in the monitoring area are used as input data. The signal is strongly affected by power system interference, lightning electromagnetic waves, communication transmission signals and artificial electromagnetic sources, and interference filtering is performed on it. The notch filter is designed to suppress power frequency interference such as 50Hz and its higher harmonics. The accurate band-stop filter parameter setting is used to eliminate the noise brought by the power system to the maximum extent without affecting the effective components of the target frequency band; at the same time, for interference sources with wide frequency bands and fast time scale changes, a multi-channel filter group is introduced, and frequency division filtering is performed according to the frequency domain analysis results to obtain a preliminary filtered signal. The wavelet threshold denoising method is introduced for the preliminary filtered signal to perform multi-scale decomposition processing. The db4 in the Daubechies wavelet family is selected as the mother wavelet. The signal is divided into approximation coefficients and detail coefficients at multiple scales by performing 5-layer discrete wavelet decomposition. Then, the high-frequency noise in the detail coefficient is suppressed according to the soft threshold principle. The soft threshold function has good smoothing performance and can effectively remove high-frequency noise interference while retaining the mutation characteristics of effective signals. It is suitable for natural electromagnetic signals, which are non-stationary, nonlinear and contain mutation characteristics. Through the layer-by-layer denoising and reconstruction process, the denoised electromagnetic signal is obtained. The signal is characterized by continuous waveform, clear mutation points and stable baseline in the time domain, and has a concentrated and clear main frequency distribution in the frequency domain. The denoised signal is segmented in the time domain to achieve time-varying extraction of non-stationary signal characteristics. The sliding window mechanism is used to divide the entire signal sequence into multiple 60-second time periods, and the adjacent windows are set to overlap by 50% to improve the time resolution and enhance the ability to capture event continuity. Each time-domain segmented signal is input into the time-frequency analysis module, and the Hilbert-Huang transform is selected as the main analysis tool. This method decomposes each signal segment into a set of intrinsic mode functions through empirical mode decomposition, and then performs Hilbert transform on each mode function to extract time-varying characteristic parameters such as instantaneous frequency, instantaneous amplitude and instantaneous phase. All this information is summarized into a multidimensional signal feature matrix, where each row corresponds to a time period and each column represents a characteristic index, such as average frequency, frequency change rate, amplitude envelope, number of modes, etc., to form a characteristic description result with time series and frequency structure. The polarization parameter analysis of the three-component vector signal contained in the signal feature matrix is ​​performed to identify the polarization characteristics of electromagnetic waves during space propagation. By combining the signals in the three directions of X, Y, and Z to construct a spatial vector trajectory, the propagation ellipse is restored in three-dimensional space, and the ellipticity of the ellipse is calculated to determine the degree of deviation from circular polarization. At the same time, the signal inclination is derived according to the main axis direction to measure the inclination of the propagation path relative to the surface, while the azimuth reflects the main propagation direction of the horizontal component of the signal. The calculation of polarization parameters helps to identify the source type and propagation mode of the signal, especially in multi-source or multi-path propagation scenarios, where polarization characteristics are used to exclude secondary interference sources.An adaptive notch filter mechanism is implemented in combination with the background interference data recorded by the reference station. The reference station is located in an area far away from the main tectonic belt and engineering activities. The signals collected by the reference station are regarded as pure background interference models. By constructing the interference template spectrum and comparing it with the target signal in the frequency domain, the main frequency and energy distribution of the interference components are extracted. Then, a notch filter is designed in the target signal to dynamically attenuate these components, thereby eliminating the influence of external disturbances while retaining the geological effective components and outputting pure geological source electromagnetic signals. The pure geological source signal is subjected to spatial correlation analysis to evaluate the consistency and spatial coupling degree of the signal between each measuring point. By calculating the mutual correlation coefficient between adjacent receiving points, the similarity of the signal waveform in time is analyzed; at the same time, the instantaneous phase of each signal segment is extracted and the phase difference is calculated to reveal the propagation delay characteristics of the signal wavefront. All the above processing results are organized according to a unified data structure standard to construct a multidimensional data cube with time sequence, spatial distribution and frequency characteristic integrity. Each data unit of the cube contains parameters such as timestamp, measuring point spatial coordinates, normalized amplitude in three directions, main frequency, phase, etc.

[0023] In a specific embodiment, the process of executing step S103 may specifically include the following steps: The multidimensional data cube is input into the integer linear programming model, and the preliminary optimized signal data is obtained by minimizing the weighted residual square sum between the signal estimation value and the observation value and introducing the L1 norm sparsity constraint; Perform singular value decomposition on the initially optimized signal data to obtain the condition number, and apply the Tikhonov regularization method based on the condition number to obtain the stabilized signal data; performing baseline correction on the stabilized signal data using a piecewise polynomial function to obtain a baseline-corrected signal; Adaptive fuzzy correction algorithm is used to apply Mamdani fuzzy reasoning to the baseline-corrected signal, and the signal features are classified according to Gaussian membership function to obtain the signal data with ambiguity eliminated. The ambiguity-eliminating signal data is processed by iterative reweighted least square method using the time difference of arrival constraint condition to obtain the ambiguous signal; The signal-to-noise ratio, amplitude stability and phase coherence parameters are calculated based on the eliminated ambiguous signals, and the signal segments that meet the requirements are screened according to the quality indicators to obtain the corrected electromagnetic pulse signal data set.

[0024] Specifically, the multidimensional data cube is input into the integer linear programming model, in which an objective function is defined to minimize the weighted residual sum of squares between the electromagnetic signal estimation value and the observation value as the optimization goal, and by introducing the L1 norm sparsity constraint term, the solution space is made more compressible and physically interpretable. The optimization problem is expressed as minimizing , where A represents the system matrix, x is the signal vector to be solved, b is the actual observation vector, and λ is the sparse constraint coefficient. The weights are adaptively allocated according to the noise level and signal credibility of the signal, thereby enhancing the model's fitting accuracy in high signal-to-noise ratio areas and suppressing the interference contribution of noise areas. After solving the integer linear programming, a set of preliminary optimized signal data is obtained. Perform singular value decomposition on the preliminary optimized signal data and decompose the system matrix A into UΣV Tform, and calculate the condition number of the matrix by observing the changes in the singular value sequence. If the condition number is high, it means that the system matrix is ​​ill-conditioned during the solution process, which may easily lead to instability of the solution or be highly sensitive to the initial disturbance. At this time, the Tikhonov regularization method is introduced for robust processing. Tikhonov regularization controls the oscillation behavior of the solution by adding the L2 norm penalty term of the solution vector to the objective function, so that the solution shows stronger robustness in the face of redundant information, thereby obtaining a set of stabilized signal data after regularization. The stabilized signal data is baseline corrected to eliminate low-frequency trends or baseline drift problems caused by environmental changes, instrument zero drift and other factors. This step adopts a piecewise polynomial function modeling strategy to divide the entire signal into several time periods, and fit a polynomial baseline curve to each section. The common choice is a quadratic or cubic polynomial to accurately capture the nonlinear change characteristics of the baseline, and determine the polynomial coefficients by the least squares method. Then, the fitting result is subtracted point by point from the original signal with the polynomial baseline curve as a reference to obtain the baseline-corrected signal, so that the signal fluctuates around the zero value in value and does not contain additional drift components, which significantly improves the accuracy of subsequent ambiguity elimination and physical parameter extraction. On the basis of baseline correction, in order to enhance the clarity of signal structure and reduce interpretation ambiguity, an adaptive fuzzy correction algorithm is introduced to process the fuzziness of the signal, and the Mamdani fuzzy inference system is used as the core inference mechanism. The system uses Gaussian membership functions to represent the fuzzy set of input signal characteristic variables, sets the input dimensions including the instantaneous frequency change rate, amplitude fluctuation degree, phase continuity index, etc. of the signal, and uses the baseline correction factor as the output. Fuzzy reasoning is performed through the If-Then formal rules in the rule base, and the fuzzy correction value after adaptive adjustment is output to obtain a set of signal data with ambiguity removed. The ambiguity-eliminating signal data is processed by iterative reweighted least squares method using the time difference of arrival constraint condition. The signal's arrival time difference is calculated based on the GPS timestamp and the position of the spatial measurement point, and the propagation time error function is constructed. The signal estimation value is used as the optimization variable to minimize the propagation path error and signal residual. By assigning different weights to each round of calculation results, the outliers are punished, thereby suppressing the error amplification caused by inconsistent propagation models. During the iteration process, the weight matrix and residual function are continuously updated until the Euclidean distance between two consecutive iterative solutions is less than the preset threshold or the maximum number of iterations is reached. Finally, a set of optimal disambiguation signals in terms of time consistency and spatial propagation characteristics are output. In order to ensure the effectiveness of the final results in data analysis and model training, the disambiguated signals are evaluated in multiple dimensions, and quantitative scores are given by calculating indicators such as the signal-to-noise ratio, amplitude stability, and phase coherence of each signal segment. The signal-to-noise ratio reflects the proportion of effective components in the signal, the amplitude stability measures the degree of fluctuation of the signal in a continuous time window, and the phase coherence is used to determine whether the signal presents physical consistency of continuous propagation.By setting evaluation thresholds for these indicators, high-quality signal segments are screened out, and abnormal data that do not meet the standards are eliminated, a set of standardized, structured, and physically interpretable corrected electromagnetic pulse signal data sets is generated.

[0025] In a specific embodiment, the process of executing step S104 may specifically include the following steps: The adaptive window length function is set by applying the multi-synchronous compression transform algorithm of time redistribution, and the local frequency characteristics of the signal are dynamically adjusted for the corrected electromagnetic pulse signal data set to obtain multiple narrow-band single-component functions; Extracting a signal characteristic parameter set including center frequency, bandwidth, energy distribution, duration and instantaneous frequency change rate from multiple narrowband single-component functions; Based on the initial electromagnetic pulse wavefront information in the corrected electromagnetic pulse signal data set, the arrival time difference between different measuring points is measured to obtain the arrival time difference data; Based on the time difference of arrival data and the spatial distribution coordinates of the receiving device, an overdetermined geological stress location equation is constructed, and the electromagnetic wave propagation velocity tensor is obtained by solving it. The extreme gradient boosting algorithm is used to establish the mapping relationship between the signal feature parameter set and the geological stress parameters, and the geological stress prediction model is obtained; The top M features that contribute most to stress prediction are selected through feature importance analysis. The extreme gradient boosting algorithm hyperparameters are automatically adjusted in combination with the Bayesian optimization method to optimize the geological stress prediction model and output the estimated geological stress distribution including the principal stress direction, stress magnitude and stress gradient.

[0026] Specifically, the corrected electromagnetic pulse signal data set is used as input, and the time-reassigned multi-synchronous compression transform algorithm is used to perform high-precision decomposition of its non-stationary time-frequency characteristics. The algorithm introduces a window length adjustment mechanism based on the calculation of the local frequency change rate of the signal, that is, by analyzing the frequency intensity of the signal in different time windows, the length of the analysis window function is dynamically set, so that a longer window is used in the frequency stable area to enhance the frequency resolution, and a shorter window is used in the area of ​​​​drastic frequency changes to enhance the time positioning ability, forming an adaptive window length function system, so that the signal shows a highly concentrated time-frequency distribution structure after transformation. This process effectively compresses the background noise and improves the accuracy of multi-component recognition while retaining the main physical characteristics of the signal. The signal processed by the multi-synchronous compression transform algorithm is decomposed into multiple narrow-band single-component functions with clear frequency structure and concentrated energy distribution. Each function corresponds to a frequency principal component with independent propagation path and geological significance. A set of characteristic parameters with strong representativeness and clear physical meaning are extracted from the above multiple narrow-band single-component functions to construct a signal characteristic parameter set. The parameter set includes the center frequency, which reflects the main frequency position where the signal energy is concentrated; the bandwidth, which indicates the frequency diffusion range and the stability of the propagation source; the energy distribution, which describes the change law of the signal strength in the frequency domain per unit time; the duration, which reflects the length and stability of the single component signal in the time dimension; and the instantaneous frequency change rate, which reveals the dynamic disturbance characteristics of the frequency component during the propagation process. The initial wavefront information in the electromagnetic signal is extracted, especially focusing on the first wave response time of each pulse event at different receiving points. By measuring the first wave arrival time difference between multiple measuring points, a wave arrival time difference data set is formed. This data contains the electromagnetic signal transmission delay between each pair of measuring points, which directly reflects the signal propagation speed and path differences. Combining these time difference data with the known three-dimensional spatial distribution coordinates of the receiving device, an electromagnetic signal propagation path model is established, and then an overdetermined geological stress location equation is constructed. Its basic form is V·T=D, where V is the unknown electromagnetic wave propagation velocity tensor, T is the wave arrival time difference matrix, and D is the spatial distance matrix between each measuring point. Since the number of sampling points is much larger than the tensor dimension, the equation is a typical overdetermined system. Robust solution algorithms such as singular value decomposition are used to infer the tensor characteristics of the propagation speed of electromagnetic waves in different directions while ensuring the stability of the solution. The propagation speed tensor represents the medium response capability of the electromagnetic signal, and also indirectly reflects the electrical anisotropy and stress state distribution of the underground rock mass. After obtaining the propagation tensor and signal characteristic parameter set, the extreme gradient boosting algorithm (XGBoost) is used to establish a nonlinear mapping relationship between the two and the actual geological stress parameters. In this supervised learning framework, the input of the model is the multi-dimensional variables in the signal characteristic parameter set, and the output is the geological stress index of the target area, including the principal stress magnitude, stress direction vector, and stress gradient value.As an integrated learning method, XGBoost constructs multiple weak classification tree models and uses weighted accumulation to achieve strong prediction performance, effectively identifying the nonlinear dependency structure between variables and outputs in high-dimensional complex feature space. During the model training process, data with known stress labels in historical monitoring data are used as training sets, and the generalization ability of the model is evaluated by cross-validation. Through this training mechanism, the XGBoost model can extract the potential relationship between multi-variable combinations and stress responses, and construct a geological stress prediction model with high robustness and high prediction accuracy. After the initial completion of the model construction, the feature importance analysis method is introduced to sort the importance of the input signal feature parameter set to identify the top M features that have the most significant impact on the prediction results. The feature importance is statistically evaluated based on indicators such as split gain, frequency or coverage within the model, which can effectively eliminate redundant or noise features, thereby simplifying the model structure and enhancing interpretability. The core hyperparameters of the XGBoost model are automatically adjusted by combining the Bayesian optimization method, including the maximum depth of the tree, learning rate, subsampling rate, regularization coefficient, etc. Bayesian optimization constructs a proxy model to probabilistically model the parameter space, and selects the optimal parameter combination based on the expected improvement function strategy, achieving global optimal performance while controlling the training time cost. After multiple rounds of training and optimization, the geological stress prediction model can quickly map any newly input electromagnetic signal feature parameter set into stress distribution estimates including principal stress direction, stress magnitude, and stress gradient.

[0027] In a specific embodiment, the execution step applies a time redistributed multi-synchronous compression transform algorithm to set an adaptive window length function, and dynamically adjusts the signal local frequency characteristics of the corrected electromagnetic pulse signal data set to obtain multiple narrowband single component functions. The process may specifically include the following steps: Performing fast Fourier transform on the corrected electromagnetic pulse signal data set to obtain initial spectrum characteristic parameters of the signal; Based on the initial spectrum characteristic parameters of the signal, the local eigenvalue of the Hessian matrix is ​​calculated for each frequency band, and the window length scaling factor is formulated by combining the frequency change rate with the time redistribution multi-synchronous compression transformation algorithm to obtain an adaptive window length function. Applying short-time Fourier transform to the corrected electromagnetic pulse signal data set, using adaptive window length function to set the time window length of different frequency bands, and obtaining the time-frequency distribution diagram; Performing a time reallocation operation on the time-frequency distribution graph to obtain an enhanced time-frequency distribution feature; The enhanced time-frequency distribution features are input into the synchronous compression kernel function to perform ridge extraction and synchronous depth control, and multiple original single-component functions are obtained; Energy contributions of multiple original single-component functions are evaluated to obtain multiple narrow-band single-component functions.

[0028] Specifically, a fast Fourier transform is performed on the corrected electromagnetic pulse signal data set, and a full spectrum conversion operation is performed on each signal segment. With its efficient computing performance and global frequency perspective, the time domain signal is mapped to the frequency domain, and the energy distribution characteristics of the signal at each frequency are extracted, thereby obtaining the initial spectrum characteristic parameters containing key information such as the main frequency component, spectrum energy density, and spectrum width. The local eigenvalues ​​of the Hessian matrix are calculated for each frequency band based on the initial spectrum characteristic parameters of the signal. A Hessian matrix based on the second-order derivative is constructed for the local area of ​​each frequency band in the spectrum, and its eigenvalue is calculated to characterize the second-order curvature of the signal energy density change in the frequency band. The larger the eigenvalue, the more drastic the energy change and the more non-stationary the frequency band is, and a shorter time window function needs to be used to improve the time positioning capability; on the contrary, if the eigenvalue is small, it means that the energy change of the frequency band is slow and the frequency is stable, and a longer window function is used to improve the frequency resolution. According to the statistical results of the Hessian matrix eigenvalues ​​and the frequency change rate index of the initial spectrum, the window length scaling factor is constructed, and the analysis window length parameters exclusive to each frequency band are defined accordingly to form a frequency-dependent and adaptively changing window length function. Based on the constructed adaptive window length function, the short-time Fourier transform is performed on the corrected electromagnetic pulse signal data set. The window length function is bound to the frequency region, and the corresponding window function length is selected for different frequency bands, so as to dynamically adjust the balance weight between time and frequency in the analysis. In the signal frequency mutation section, the window length is shortened to improve the time resolution; while in the frequency stable region, the window length is appropriately lengthened to improve the frequency resolution. This strategy improves the ability to analyze complex signal structures and generates a nonlinear time-frequency distribution diagram with a time-varying window structure. The time-frequency distribution diagram is time-redistributed. The fuzzy energy distributed in the non-central area is concentrated and compressed back to its real physical position to obtain a more structured time-frequency representation. By utilizing the instantaneous frequency and group delay information of the signal, each point of energy is remapped to its real time-frequency center, removing the frequency leakage effect and fuzzy diffusion caused by the window function, and obtaining an enhanced time-frequency distribution feature diagram. The enhanced time-frequency distribution features are input into the synchronous compression kernel function to perform ridge extraction and synchronous depth control operations. The essence of ridge extraction is to identify the main energy path in the time-frequency diagram, which corresponds to the main response trajectory generated by different physical components in the electromagnetic signal during the propagation process; and the synchronous compression kernel function combines local consistency evaluation with global curvature fitting to ensure that the extracted ridges have significant characteristics in energy intensity, frequency coherence and time consistency. At the same time, in order to prevent the misidentification of non-principal components during the extraction process, a synchronous depth control strategy is introduced. According to the embedding depth of the signal ridge in the time-frequency diagram, the ridge density and the relative weight relationship with the adjacent ridges, only the main ridge structure with stable structure and physical rationality is retained, and multiple original single-component functions are output. Each function is a structural unit extracted from the enhanced time-frequency diagram, with a clear time-frequency distribution range and stable physical propagation directionality.In order to compress the data size and improve the analysis focus, it is necessary to evaluate the energy contribution of these original single-component functions. The evaluation process calculates the proportion of each single component in the total energy and comprehensively scores its energy concentration, time-frequency stability, and propagation consistency. In the evaluation results, components with high energy contribution, clear signals, and complete structures are selected and retained as the final results, and redundant, low-intensity, or morphologically unstable components are eliminated. The multiple narrowband single-component functions finally obtained have a high degree of physical interpretability and good mathematical stability and information decoupling capabilities.

[0029] In a specific embodiment, the process of executing step S105 may specifically include the following steps: Based on the estimated value of geological stress distribution, the stress parameters of discrete measuring points are interpolated into stress field data using Kriging interpolation method; Perform tensor analysis on stress field data to obtain the principal stress distribution characteristics, and formulate stress concentration area identification criteria based on the principal stress distribution characteristics to obtain the preliminarily identified stress concentration areas; The Markov random field model is used to perform spatial constraint optimization on the initially identified stress concentration area to obtain the optimized stress concentration area. A stress concentration index evaluation system is established based on the optimized stress concentration area, and the graded stress concentration area is calculated according to the stress concentration index evaluation system; Three-dimensional modeling is performed based on hierarchical stress sets to generate a three-dimensional model of the geological stress concentration area containing stress tensor field, lithology distribution, fracture development degree and groundwater distribution information.

[0030] Specifically, based on the estimated value of geological stress distribution, the stress parameters of discrete measuring points are interpolated into stress field data using the Kriging interpolation method. Kriging interpolation fits the variation function between measuring points, fully considers the spatial autocorrelation and the weight relationship between measuring point spacing, and performs interpolation calculations on the entire area to achieve the optimal estimation of stress values ​​at unobserved points, minimize the estimation error and maintain the rationality of the overall spatial structure, and generate a three-dimensional continuous stress field data with a geological interpretation basis and mathematical stability. Based on the stress field data, tensor analysis is carried out to extract the principal values ​​of the stress tensor and their principal axis directions in each grid cell, including the maximum principal stress. , intermediate principal stress , minimum principal stress and its corresponding direction cosines, construct the stress ellipsoid expression of each point. By clustering analysis and gradient operation on the spatial distribution of the principal stress values ​​of the regional stress tensor field, the area with highly concentrated and drastic changes in local stress, namely, the stress abnormality concentration zone, is identified. In order to standardize the identification mechanism, a clear criterion for distinguishing stress concentration areas is formulated, such as stipulating that the maximum principal stress in a certain area is More than 1.5 times the average principal stress of the whole area, and and If the difference is more than twice the average stress difference, the area is marked as a potential stress concentration area that is preliminarily identified. The Markov random field model is applied to the preliminarily identified stress concentration area for spatial constraint optimization. The model defines the adjacency relationship and state transition probability between grid cells, and establishes a joint probability distribution between the stress state and the surrounding environment, so that the spatial consistency constraint is strengthened while retaining the local identification results, thereby eliminating pseudo-outliers, filling spatial fault zones, and optimizing the structural markers of the entire stress field through maximum a posteriori estimation. A stress concentration index evaluation system is established based on the optimized stress concentration area. The index comprehensively considers the maximum principal stress intensity, principal stress difference, and the relative spatial relationship with the geological structure (especially faults). The formula is: ; in and They represent the average values ​​of the maximum principal stress and stress difference in the region, respectively, and D is the normalized distance index between the center point of the unit and the nearest active fault.

[0031] The stress concentration areas are objectively scored and ranked by the SCI index, and the classification threshold is set according to the score. For example, SCI>4.0 is a first-level high-risk area, 2.0≤SCI≤4.0 is a second-level medium-risk area, and 1.0≤SCI<2.0 is a third-level concern area. After the classification is completed, the stress concentration areas of different levels are treated separately, and different simulation weights and response priorities are assigned in the subsequent modeling. Based on the above-mentioned graded stress concentration areas, a three-dimensional geological stress concentration area model integrating multi-source information is constructed. In this model, the stress tensor field constitutes the main framework, which accurately expresses the spatial distribution law of the principal stress; the lithology distribution is used as the structural background of the model to reflect the mechanical response differences of different stratigraphic units; the degree of fracture development is constructed by superimposing the geological structural characteristic data such as faults, fracture zones, and shear zones to construct a spatial weak surface system, providing a basis for mechanical boundary and stress reconstruction; the groundwater distribution is embedded in the model as an important factor of fluid-induced stress field disturbance to evaluate the triple coupling relationship of water-rock-stress. In order to improve the model resolution, the adaptive grid generation method is applied in the area with high SCI values ​​to locally refine the grid units and enhance the modeling accuracy. The model is presented in a three-dimensional grid data structure, and the output forms include vector field diagrams, principal stress ellipsoid visualization, cross-section diagrams, isosurface diagrams and other expressions, which have spatial operability and analytical scalability.

[0032] Among them, the Markov random field model is applied to the initially identified stress concentration area for spatial constraint optimization to obtain the optimized stress concentration area, including: dividing the initially identified stress concentration area into regular grid units, establishing the connection relationship between each unit and its adjacent units, constructing a Markov random field space network with a four-neighborhood or eight-neighborhood structure, and obtaining a spatial constraint network structure; defining a single-point potential function and a paired potential function based on the spatial constraint network structure, wherein the single-point potential function characterizes the degree of conformity between the stress value of the grid unit and the observed data, and the paired potential function characterizes the spatial continuity constraint strength between adjacent units, and obtaining a Markov random field energy function; using the maximum likelihood estimation or quasi-maximum likelihood estimation method, combined with the stress distribution characteristics of the known geological structure area, the Markov random field energy function is estimated. The weight parameters in the function are trained and calculated to obtain the optimal model parameters; according to the boundaries of geological structural units, fault locations and lithological boundaries, hard boundary constraints and soft transition conditions are set for the Markov random field model to obtain constraints that take into account geological structural characteristics; the iterative conditional pattern algorithm or graph cut algorithm is applied to the Markov random field model with constraints, and the marking optimization is performed by minimizing the global energy function. The iterative calculation is performed until convergence or the preset maximum number of iterations is reached to obtain the optimized marking result; the optimized marking result is post-processed, including removing isolated areas with an area smaller than a threshold, merging similar areas with too small spacing, and smoothing regional boundaries. The results are verified according to the principle of geomechanical equilibrium to obtain the optimized stress concentration area that eliminates pseudo-anomalies and maintains spatial continuity.

[0033] Those skilled in the art can clearly understand that, for the convenience and brevity of description, the specific working processes of the above-described systems, systems and units can refer to the corresponding processes in the aforementioned method embodiments and will not be repeated here.

[0034] If the integrated unit is implemented in the form of a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention is essentially or partly contributed to the prior art or all or part of the technical solution can be embodied in the form of a software product. The computer software product is stored in a storage medium, including several instructions to enable a geological stress concentration zone intelligent identification device (which can be a personal computer, server, or network device, etc.) to perform all or part of the steps of the method described in each embodiment of the present invention. The aforementioned storage medium includes: U disk, mobile hard disk, read-only memory (ROM), random access memory (RAM), disk or optical disk and other media that can store program codes.

[0035] As described above, the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit the same. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that the technical solutions described in the aforementioned embodiments may still be modified, or some of the technical features thereof may be replaced by equivalents. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for intelligently identifying geological stress concentration areas, characterized in that: include: Multiple vector electromagnetic induction coils are arranged in a grid shape in the monitoring area to collect original electromagnetic pulse signals containing three-directional components: X, Y, and Z. Performing denoising and time-frequency analysis on the original electromagnetic pulse signal to obtain a multidimensional data cube; Use the integer linear programming model to process the multidimensional data cube to obtain the preliminary optimized signal data; perform singular value decomposition and apply the Tikhonov regularization method to obtain the stabilized signal data; use the piecewise polynomial function to perform baseline correction; use the adaptive fuzzy correction algorithm to apply the Mamdani fuzzy reasoning to process the signal; use the time difference of arrival constraint condition to perform iterative reweighted least squares processing to obtain the disambiguated signal, and select the signal segments that meet the requirements according to the calculation parameters of the disambiguated signal to obtain the corrected electromagnetic pulse signal data set; Decomposing and mapping the corrected electromagnetic pulse signal data set to obtain an estimated value of geological stress distribution; A three-dimensional model of the geological stress concentration area is constructed according to the estimated value of the geological stress distribution.

2. The intelligent identification method of geological stress concentration area according to claim 1 is characterized in that: The method of arranging a plurality of vector electromagnetic induction coils in a grid-like manner in the monitoring area to collect original electromagnetic pulse signals in three directions of X, Y and Z includes: Performing a geological structure complexity analysis on the monitoring area to obtain a geological structure complexity analysis result, and determining an equilateral triangle arrangement mode according to the geological structure complexity analysis result to obtain a receiving device arrangement plan; Each receiving device in the receiving device arrangement scheme is equipped with three orthogonal vector electromagnetic induction coils to obtain an electromagnetic signal receiving unit; Integrate a GPS clock module into each electromagnetic signal receiving unit to obtain a time synchronization acquisition system, and connect the electromagnetic signal receiving units to form a mesh network according to the time synchronization acquisition system and the low-power wireless transmission module to obtain a distributed monitoring network; Selecting a reference site in the monitoring area, and recording the external electromagnetic interference signal of the area according to the reference site to obtain reference interference data; The electromagnetic pulse signal of the distributed monitoring network is monitored, and the effective signal is marked in combination with the reference interference data to obtain the original electromagnetic pulse signal containing the components in three directions of X, Y and Z.

3. The intelligent identification method of geological stress concentration area according to claim 1 is characterized in that: The performing of denoising and time-frequency analysis on the original electromagnetic pulse signal to obtain a multidimensional data cube includes: Performing interference filtering on the original electromagnetic pulse signal to obtain a preliminary filtered signal, and performing wavelet threshold denoising on the preliminary filtered signal to obtain a denoised signal; Performing time-domain segmentation processing on the noise-reduced signal to obtain a segmented time-domain signal, and performing time-frequency analysis on the segmented time-domain signal to obtain a signal feature matrix; Calculating polarization parameters including ellipticity, inclination and azimuth according to the three-component vector signal in the signal characteristic matrix, and performing adaptive notch filtering in combination with the external interference signal recorded at the reference station to obtain a pure geological source signal; The pure geological source signal is subjected to spatial correlation analysis, the mutual correlation coefficient and phase difference between adjacent measuring points are calculated, and the signals are organized in a unified format to form a multi-dimensional data cube including time, space, amplitude, frequency and phase.

4. The intelligent identification method of geological stress concentration area according to claim 1 is characterized in that: The method uses an integer linear programming model to process a multidimensional data cube to obtain preliminary optimized signal data; performs singular value decomposition and applies the Tikhonov regularization method to obtain stabilized signal data; performs baseline correction using a piecewise polynomial function; uses an adaptive fuzzy correction algorithm and applies Mamdani fuzzy reasoning to process the signal; uses a time difference of arrival constraint condition to perform iterative reweighted least squares processing to obtain an ambiguous signal, and selects signal segments that meet the requirements according to the calculation-related parameters of the ambiguous signal to obtain a corrected electromagnetic pulse signal data set, including: Inputting the multidimensional data cube into an integer linear programming model, obtaining preliminary optimized signal data by minimizing the weighted residual square sum between signal estimation values ​​and observation values ​​and introducing L1 norm sparsity constraints; Performing singular value decomposition on the initially optimized signal data to obtain a condition number, and applying a Tikhonov regularization method according to the condition number to obtain stabilized signal data; Performing baseline correction on the stabilized signal data using a piecewise polynomial function to obtain a baseline-corrected signal; Adopting an adaptive fuzzy correction algorithm to apply Mamdani fuzzy reasoning to the baseline-corrected signal, classifying the signal features according to the Gaussian membership function, and obtaining signal data with ambiguity eliminated; Using the time difference of arrival constraint condition, the signal data to be deambiguated is processed by iterative reweighted least square method to obtain a deambiguated signal; The signal-to-noise ratio, amplitude stability and phase coherence parameters are calculated according to the eliminated ambiguous signal, and the signal segments meeting the requirements are screened according to the quality index to obtain a corrected electromagnetic pulse signal data set.

5. The method for intelligently identifying geological stress concentration areas according to claim 1, characterized in that: Decomposing and mapping the corrected electromagnetic pulse signal data set to obtain an estimated value of geological stress distribution includes: Applying a time-redistributed multi-synchronous compression transformation algorithm to set an adaptive window length function, and dynamically adjusting the signal local frequency characteristics of the corrected electromagnetic pulse signal data set to obtain multiple narrowband single-component functions; Extracting a signal characteristic parameter set including center frequency, bandwidth, energy distribution, duration and instantaneous frequency change rate from the plurality of narrowband single-component functions; For the initial electromagnetic pulse wavefront information in the corrected electromagnetic pulse signal data set, measuring the time difference of arrival between different measuring points to obtain the time difference of arrival data; Based on the time difference of arrival data and the spatial distribution coordinates of the receiving device, an overdetermined geological stress location equation is constructed, and the electromagnetic wave propagation velocity tensor is obtained by solving the equation; An extreme gradient boosting algorithm is used to establish a mapping relationship between the signal characteristic parameter set and the geological stress parameter to obtain a geological stress prediction model; The top M features that contribute most to stress prediction are selected through feature importance analysis, and the extreme gradient boosting algorithm hyperparameters are automatically adjusted in combination with the Bayesian optimization method to optimize the geological stress prediction model and output an estimated value of geological stress distribution including principal stress direction, stress magnitude and stress gradient.

6. The intelligent identification method of geological stress concentration area according to claim 5 is characterized in that: The multi-synchronous compression transformation algorithm using time redistribution is used to set an adaptive window length function, and the local frequency characteristics of the signal are dynamically adjusted for the corrected electromagnetic pulse signal data set to obtain multiple narrowband single-component functions, including: Performing a fast Fourier transform on the corrected electromagnetic pulse signal data set to obtain initial signal spectrum characteristic parameters; Based on the initial spectrum characteristic parameters of the signal, the local eigenvalues ​​of the Hessian matrix are calculated for each frequency band, and a time-redistributed multi-synchronous compression transformation algorithm is applied in combination with the frequency change rate to formulate a window length scaling factor to obtain an adaptive window length function; Applying short-time Fourier transform to the corrected electromagnetic pulse signal data set, using the adaptive window length function to set the time window lengths of different frequency bands, and obtaining a time-frequency distribution diagram; Performing a time reallocation operation on the time-frequency distribution graph to obtain an enhanced time-frequency distribution feature; Inputting the enhanced time-frequency distribution feature into a synchronous compression kernel function to perform ridge extraction and synchronous depth control to obtain a plurality of original single-component functions; Energy contribution evaluation is performed on the multiple original single-component functions to obtain multiple narrow-band single-component functions.

7. The method for intelligently identifying geological stress concentration areas according to claim 1, characterized in that: The step of constructing a three-dimensional model of a geological stress concentration area according to the estimated value of geological stress distribution comprises: Based on the estimated value of geological stress distribution, the stress parameters of discrete measuring points are interpolated into stress field data using Kriging interpolation method; Performing tensor analysis on the stress field data to obtain principal stress distribution characteristics, and formulating a stress concentration area discrimination criterion based on the principal stress distribution characteristics to obtain a preliminarily identified stress concentration area; Applying a Markov random field model to perform spatial constraint optimization on the initially identified stress concentration area to obtain an optimized stress concentration area; Establishing a stress concentration index evaluation system based on the optimized stress concentration area, and calculating the graded stress concentration area according to the stress concentration index evaluation system; Based on the hierarchical stress set, three-dimensional modeling is performed to generate a three-dimensional model of the geological stress concentration area including stress tensor field, lithology distribution, fracture development degree and groundwater distribution information.

Images