A method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field

Abstract

This invention discloses a method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field, relating to the field of geoelectric field observation technology. The method involves collecting geoelectric field observations and preprocessing them; separating the natural electric field, interference, and the geoelectric field from the preprocessed observations; extracting all 720 orders of tidal waves of the geoelectric field and calculating the amplitude and phase parameters of each harmonic; selecting 1440 daily geoelectric field data points from two channels (NS and EW) and performing trend fitting using the least squares method to obtain the fitting line parameters; and calculating the dominant azimuth angle α of the geoelectric field fractures in the rock mass based on the fitting line parameters. This invention can effectively reflect the dominant direction of rock mass fractures and may have higher sensitivity to suspected pre-earthquake anomalies, thus providing a useful supplement to the analysis of pre-earthquake geoelectric anomalies.

Description

Technical Field

[0001] This invention relates to the field of geoelectric field observation technology, and more specifically to a method for calculating the dominant azimuth angle of rock mass fractures in the geoelectric field based on trend fitting. Background Technology

[0002] Currently, the geoelectric field, as an important component of the Earth's electromagnetic environment, is closely related to the electrical structure and tectonic activity of the Earth's crust. According to the classification of geoelectric field origins, in addition to the stable natural electric fields generated by ore bodies, groundwater, and various water systems, the geoelectric field also includes rock fissure water seepage controlled by lunar and solar tidal forces and the geocell induced by ionospheric Sq current. Based on the understanding of the origin of the geoelectric field, researchers have constructed a rock fissure water (charge) seepage (movement) model using geoelectric field tidal waves. This model posits that dynamic changes in fissure structure (such as shearing and ordered arrangement) lead to abnormal fluctuations in the direction and intensity of the geoelectric field. Especially during the gestation of strong earthquakes, the accumulation of tectonic stress may cause rock fissures to evolve from disordered generation to directional arrangement or shear fracture, thus exhibiting a significant deflection of the dominant azimuth angle or directional arrangement of rock fissures before the earthquake.

[0003] Many scholars have conducted extensive research on the dominant azimuth of rock mass fractures in the geoelectric field, based on actual earthquake case data. Research results show that anomalies in the α-angle are well reflected in moderate-to-strong and shallow-focus earthquakes, and the stronger the earthquake, the earlier the anomaly appears; the closer the epicentral distance, the greater the likelihood of the anomaly. In terms of timing, the α-angle anomaly is a short-to-medium-term anomaly, with a particularly significant coseismic effect, and the anomaly exhibits quasi-synchronous and selective characteristics at different stations. Regarding the anomaly morphology, different stations may exhibit different anomaly morphologies. In terms of the anomaly mechanism, researchers have found that the α-angle anomaly is significantly influenced by tectonic effects, often located on the same tectonic line as the earthquake epicenter. The correlation between its dominant azimuth and the trend of the regional principal compressive stress P is basically consistent with rock physics theory, effectively capturing micro-changes in the structure of the seismogenic rock mass at the site.

[0004] In addition, the α-angle method has strong anti-interference capabilities and is less affected by typical interferences in geoelectric field observations such as high-voltage direct current transmission faults, geoelectric storms, artificial power supply for georesistivity observations, and subways. It is also currently the only widely applicable method for predicting earthquakes using geoelectric fields. However, the calculation of existing tidal geoelectric field rock mass fracture models still has certain limitations in terms of the dynamic evolution of the α-angle and the dependence on the quality of the original data, and there is still room for further optimization and improvement.

[0005] Therefore, how to provide a new method for calculating the α angle that can fully utilize the data of the Earth's electric field and improve the ability to identify pre-earthquake anomalies is a problem that urgently needs to be solved by those skilled in the art. Summary of the Invention

[0006] In view of this, the present invention provides a method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field, so as to solve the problems existing in the background art.

[0007] To achieve the above objectives, the present invention adopts the following technical solution: A method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field includes: S1. Collect and preprocess the observed values ​​of the ground electric field; S2. Separate the natural electric field, interference, and earth electric field from the preprocessed geoelectric field observations; S3. Extract all 720 tidal waves of the earth's electric field and calculate the amplitude and phase parameters of each harmonic. S4. Select 1440 geoelectric field data points from the NS and EW survey channels on a single day, and perform trend fitting using the least squares method to obtain the parameters of the fitted straight line. S5. Calculate the dominant azimuth angle α of rock mass fractures in the geoelectric field based on the parameters of the fitted straight line.

[0008] Optionally, S2 specifically includes: Earth's electric field j Observations at time E j From natural electric field E spj Earth electric field E Tj and interference E rj Composition formula ; Take 1440 minutes of a single day E j Calculate the mean, which approximates the natural electric field. E spj Thus, the earth's electric field is obtained. E Tj = E j - E spj - E rj When the interference is negligible, E Tj = E j - E spj .

[0009] Optionally, the time series of the earth's electric field in S3 The expression is:

[0010] In the formula, , For the phase of the kth harmonic, according to Calculation of positive and negative values, k The order of the harmonics. t For the time points of the Earth's electric field sequence, This is the daily average value of the Earth's electric field. For the first k The amplitude of the first harmonic.

[0011] Optionally, the objective function for the least squares fitting line in S4 is:

[0012] In the formula, ( , ( ) is a series of points, a , b represents the parameters for fitting the straight line, and n represents the total number of data points for the Earth's electric field.

[0013] Optionally, the a The expression is: ; The expression for b is: .

[0014] Optionally, the formula for calculating the dominant azimuth angle α of the rock mass fracture in the geoelectric field is: when a When >0, ; when a When <0, .

[0015] Optionally, it also includes S6, using the sliding window method to calculate the distance standard deviation of the α angle and the maximum eigenvalue of the covariance matrix, and evaluating its polarization stability.

[0016] Optionally, the standard deviation of the distance is expressed as:

[0017] In the formula, For the first i points ( , ) to the center of mass ( u x , u y The Euclidean distance between the data points and the total number of data points. This represents the average distance from the data point to the centroid.

[0018] Optionally, the covariance matrix is ​​represented as:

[0019] In the formula, and Representing variables respectively x and y Its own variance It represents the covariance of two variables.

[0020] As can be seen from the above technical solution, compared with the prior art, this invention discloses a method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field. By separating the natural electric field and the geoelectric field, and combining the phase information of each order of electric field harmonics, the dominant azimuth angle of rock mass fractures is calculated using trend fitting of data from all time periods within a day. This invention can effectively reflect the dominant direction of rock mass fractures and may have higher sensitivity to suspected pre-earthquake anomalies, thus providing a useful supplement to the analysis of pre-earthquake geoelectric anomalies. Attached Figure Description

[0021] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on the provided drawings without creative effort.

[0022] Figure 1 A dominant orientation map of the geoelectric field in rock fractures based on trend fitting is provided for this invention. Figure 2 The main active faults and station distribution map of the research area provided for this invention; Figure 3 The flowchart of the evaluation method provided by the present invention; Figure 4 The present invention provides a plot showing the degree of concentration and dispersion of data around the central value under high and low standard deviations. Figure 5 The random noisy sine function graph provided for this invention (y(t)=sin(t)+ε(t), where the red solid line represents the principal component direction, the yellow solid line represents the local confidence ellipse, and the magnitude of the maximum eigenvalue λ measures the degree of concentration of the data in the principal direction). Figure 6a The α angle and dispersion index of the Pingliang Terrace calculated using traditional methods; Figure 6b To improve the α angle and dispersion index of Pingliangtai calculated based on trend fitting; Figure 6cThe α angle and dispersion index of the ancient Fengtai area calculated using traditional methods; Figure 6d To improve the α angle and dispersion index of Gufengtai calculated based on trend fitting; Figure 7a The angle α of Tianshuitai is calculated using traditional methods; Figure 7b The improved α angle of Ganzitai calculated based on trend fitting; Figure 7c The α angle of Chengdu Station calculated using traditional methods; Figure 7d The improved α angle of Chengdu Station calculated based on trend fitting; Figure 7e The α angle of Hanwangtai calculated using traditional methods; Figure 7f The improved α angle of Hanwangtai is calculated based on trend fitting. Detailed Implementation

[0023] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0024] This invention discloses a method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field, including: S1. Collect and preprocess the observed values ​​of the ground electric field; S2. Separate the natural electric field, interference, and earth electric field from the preprocessed geoelectric field observations; S3. Extract all 720 tidal waves of the earth's electric field and calculate the amplitude and phase parameters of each harmonic. S4. Select 1440 geoelectric field data points from the NS and EW survey channels on a single day, and perform trend fitting using the least squares method to obtain the parameters of the fitted straight line. S5. Calculate the dominant azimuth angle α of rock mass fractures in the geoelectric field based on the parameters of the fitted straight line.

[0025] The core steps for calculating the dominant azimuth angle of rock mass fractures based on trend fitting improvement are as follows: The specific principle is explained as follows: 1. Trend-fitting calculation method 1) Separating the Earth's electric field from the Earth's electric field Observations of the Earth's electric field at time j From natural electric field Earth's electric field and interference composition:

[0026] 1440 minutes a day The mean of the sums can basically eliminate The influence of this, if we choose data from a day when the natural electric field is stable, then in this case, the data for each day... The value and the day Since their means are the same, we have:

[0027] 2) Take all 720 orders of tidal waves of the Earth's electric field

[0028] For the time series of observed geoelectric field data Mathematically, this can be represented as:

[0029] In the formula: k The order of the harmonics. t For the time points of the Earth's electric field sequence, This is the daily average value of the Earth's electric field. For the first k The amplitude of the first harmonic. , and phase The calculation formulas are as follows:

[0030] The above harmonic decomposition process shows that the time series of the daily Earth electric field can be completely characterized by the first 720 harmonic components. Traditional methods usually only select the first 10 dominant harmonic components to reconstruct the signal and perform subsequent trend analysis. Although this can capture the main periodic terms, it inevitably loses high-frequency information. Therefore, this invention directly uses the time series of the Earth electric field, and its mathematical essence is equivalent to utilizing the complete physical information carried by all 720 harmonics.

[0031] 3) Fitting trend

[0032] Using two measurement channels, NS and EW, and combining the phase information of harmonics, an azimuth angle is obtained by fitting the trend change of 1440 data points per day (e.g., ...). Figure 1 The azimuth angle for the day is obtained by rotating counterclockwise from the direction of the rock fissure to due north, within the range of [0, 180) and not by directly obtaining the azimuth angle of the day using only the harmonic amplitude ratio of each order.

[0033] For a series of points ( , Using the least squares method to fit the straight line y= a x+b, that is, minimizing the objective function:

[0034] Make the objective function S( a b) For a The partial derivative of b is 0:

[0035] Therefore, the parameters of the fitted line are:

[0036] Therefore, the fitting angle can be expressed as:

[0037] Therefore, for the first i The ground electric field at time EW and NS is and The formula for angle α is: when a When >0, ; when a When <0, .

[0038] Furthermore, this also includes the study area and validation design.

[0039] 2.1 Study Area and Data Sources

[0040] To verify the feasibility and advantages of the improved method, this invention selected 14 geoelectric stations within a radius of 560 km centered on the 2017 Jiuzhaigou Ms7.0 earthquake (103.82°E, 33.20°N), and collected geoelectric field monitoring data for three consecutive years from May 2015 to May 2018 (station distribution is shown in [reference]). Figure 2 The stability and anomaly detection capabilities of the calculation results from the traditional and improved methods were compared and analyzed.

[0041] 2.2 Methodology and Evaluation Criteria

[0042] 2.2.1 Method Flow

[0043] First, the data is preprocessed (error values ​​are removed, outliers are marked, etc.). Then, the trend is fitted based on the distribution of geodetic current values ​​in the NS and EW directions over a day. Finally, the standard deviation and covariance of the daily calculated α angle are calculated using the sliding window method to evaluate polarization stability. The method flow is as follows: Figure 3 .

[0044] 2.2.2 Evaluation Criteria

[0045] The quantification station has well-developed rock fractures ( To assess the temporal stability of the α-angle, this paper employs a dual-index evaluation system based on a sliding window: the first is the Distance Standard Deviation (DSD), reflecting the spatial clustering of data points within a local window; the second is the Maximum Eigenvalue of the Covariance Matrix (MECM), characterizing the dispersion of data along the principal component directions. Both indices reflect the polarization stability of the α-angle.

[0046] Definition of standard deviation of distance:

[0047] in, For the first i points ( , ) to the center of mass ( , The Euclidean distance between the data points and the total number of data points. This represents the average distance from the data point to the centroid.

[0048] DSD reflects the spatial consistency of the system within a local time window by calculating the radial oscillation of data points around the centroid. For example... Figure 4 As shown, a smaller σ value (lower standard deviation) indicates a narrower range of data fluctuations, closer to the central value, meaning that the α angle remains tightly clustered within this period, and its polarization direction consistency is better.

[0049] Definition of covariance matrix:

[0050] in Var ( x )and Var ( y ) represent variables respectively x and y Its own variance cov ( x,y The eigenvalue represents the covariance of two variables. Perform eigenvalue decomposition on matrix C and find the largest eigenvalue. λmax .

[0051] The largest eigenvalue of the covariance matrix represents the variance of the data along the principal component direction (i.e., the direction of greatest dispersion). A larger λmax indicates greater dispersion of the data along that direction; conversely, a smaller λmax indicates greater concentration. An example of this is a random, noisy sine function. Figure 5 The values ​​of λmax at the three yellow ellipses are 0.15, 0.18 and 0.36 respectively, indicating that the dispersion of the data in the three regions gradually increases in the main direction.

[0052] 3. Results

[0053] 3.1 Method Stability

[0054] For the orderly arrangement of rock fractures at the observation site ( At stations such as Pingliang Station and Gufeng Station in Gansu Province, both methods can maintain stable polarization during normal background changes, reflecting the correctness of this invention. Figure 6a , Figure 6b and Figure 6c , Figure 6d It can be seen that the calculation results of the two methods are very similar. Specific evaluation index parameter calculations show that the average SDS and MECM (averaging values ​​below the median SDS and MECM over three years) of the present invention during the calm period of the α angle at Pingliang Terrace decreased by 0.035% and 0.59%, respectively, indicating that the polarization stability of the present invention for calculating the α angle at Pingliang Terrace is slightly better than the traditional method. However, the average SDS and MECM of the present invention increased by 0.22% and 2.22%, respectively, during the calm period of the α angle at Gufeng Terrace, indicating that the polarization stability of the present invention for calculating the α angle at Gufeng Terrace is slightly worse than the traditional method. The difference in polarization stability between the two methods is not significant, and both can correctly reflect the changes in the α angle.

[0055] 3.2 Capability to detect anomalies in the α angle

[0056] The present invention exhibits a more significant effect on the difference before and after the α-angle anomaly, demonstrating stronger capture capability, as seen at the Tianshui Observatory in Gansu (…). Figure 7a , Figure 7b ), Chengdu TV Station, Sichuan Figure 7c , Figure 7d ) and Gansu Hanwangtai ( Figure 7e , Figure 7f ).

[0057] For Tianshuitai ( Figure 7a , Figure 7bThe α angle exhibited an anomalous jump before the 2017 Ms7.0 Jiuzhaigou earthquake. Between August 2016 and June 2017, and between August 2017 and February 2018, the α angle calculated by both methods mainly fluctuated around 90°. However, between June and August 2017 (the area within the black box in the figure), both methods showed a jump in the α angle. Because the α angle directions defined by the present invention and the traditional method differ (the α angle calculated by the present invention is equal to or complementary to that calculated by the traditional method), the anomaly morphology is slightly different.

[0058] For Chengdu TV ( Figure 7c , Figure 7d The α angle showed an anomaly of increased variation range before the 2017 Ms7.0 Jiuzhaigou earthquake. Traditional methods for the period from January to March 2017 before the earthquake... Increased from April to August (The part within the black box in the image) This invention was developed between January and March 2017, prior to the earthquake. Increased to April to August (The area within the black box in the image). Compared to traditional methods, the present invention shows a greater range of changes and more significant differences before an earthquake.

[0059] Similarly, for Hanwangtai ( Figure 7e , Figure 7f The α angle also showed an abnormally large range of variation before the 2017 Jiuzhaigou earthquake (Ms7.0). This invention focuses on the data from January to April 2017 before the earthquake. Increased to May to August (The part within the black box in the image), and traditional calculation results cannot reflect this feature.

[0060] The various embodiments in this specification are described in a progressive manner, with each embodiment focusing on its differences from other embodiments. Similar or identical parts between embodiments can be referred to interchangeably. For the apparatus disclosed in the embodiments, since they correspond to the methods disclosed in the embodiments, the description is relatively simple; relevant parts can be referred to the method section.

[0061] The above description of the disclosed embodiments enables those skilled in the art to make or use the invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the invention. Therefore, the invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features disclosed herein.

Claims

1. A method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field, characterized in that, include: S1. Collect and preprocess the observed values ​​of the ground electric field; S2. Separate the natural electric field, interference, and earth electric field from the preprocessed geoelectric field observations; S3. Extract all 720 tidal waves of the earth's electric field and calculate the amplitude and phase parameters of each harmonic. S4. Select 1440 geoelectric field data points from the NS and EW survey channels on a single day, and perform trend fitting using the least squares method to obtain the parameters of the fitted straight line. S5. Calculate the dominant azimuth angle α of rock mass fractures in the geoelectric field based on the parameters of the fitted straight line.

2. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field according to claim 1, characterized in that, S2 specifically includes: Earth's electric field j Observations at time E j From natural electric field E spj Earth electric field E Tj and interference E rj Composition formula ; Take 1440 minutes of a single day E j Calculate the mean, which approximates the natural electric field. E spj Thus, the earth's electric field is obtained. E Tj = E j - E spj - E rj When the interference is negligible, E Tj = E j - E spj .

3. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field according to claim 1, characterized in that, The time series of the earth electric field in S3 The expression is: In the formula, , For the phase of the kth harmonic, according to Calculation of positive and negative values, k The order of the harmonic is denoted by . t For the time points of the Earth's electric field sequence, This is the daily average value of the Earth's electric field. For the first k The amplitude of the first harmonic.

4. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field according to claim 1, characterized in that, The objective function for the least squares fitting line in S4 is: In the formula, ( , ( ) is a series of points, a , b represents the parameters for fitting the straight line, and n represents the total number of data points for the Earth's electric field.

5. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field according to claim 4, characterized in that, The a The expression is: ； The expression for b is: 。 6. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field according to claim 5, characterized in that, The formula for calculating the dominant azimuth angle α of the rock mass fracture in the geoelectric field is as follows: when a When >0, ; when a When <0, .

7. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of the geoelectric field according to claim 1, characterized in that, It also includes S6, which uses the sliding window method to calculate the distance standard deviation and the maximum eigenvalue of the covariance matrix of the α angle, and evaluates its polarization stability.

8. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field according to claim 7, characterized in that, The standard deviation of the distance is expressed as: In the formula, For the first i points ( , ) to the center of mass ( u x , u y The Euclidean distance between the data points and the total number of data points. This represents the average distance from the data point to the centroid.

9. The method for calculating the dominant azimuth angle of rock mass fractures based on trend fitting of geoelectric field according to claim 7, characterized in that, The covariance matrix is ​​expressed as follows: In the formula, and Representing variables respectively x and y Its own variance It represents the covariance of two variables.

Images