Multi-dimensional quantitative evaluation method for crustal stability

Abstract

The invention belongs to the technical field of earth crust information, and aims to provide an earth crust stability multi-dimensional quantitative evaluation method which integrates earth surface deformation, underground activities, modern geodetic survey and long-term geological data. A crustal stability multi-dimensional monitoring and evaluation technology system composed of a plurality of first-level index and second-level index parameters is determined from different spatial and temporal scales, the objectivity and scientific rationality of safety quantitative evaluation of areas along urban activity fractures are remarkably improved, and guarantee is provided for engineering layout site selection and high-quality operation and maintenance.

Description

Technical Field

[0001] The invention belongs to the field of crust information technology, and specifically relates to a multi-dimensional quantitative evaluation method for crust stability. Background Art

[0002] The efficient and refined monitoring and evaluation of crustal stability is of great significance to the site selection of major national projects, regional planning and layout, and disaster prevention and mitigation. Especially in the period of rapid economic development, engineering construction is intensive, and earthquake activities and geological disasters are frequent. Accompanied by the over-intensive mining of coal mining areas and excessive groundwater extraction, stratum collapse and ground fissures are becoming increasingly prominent, which seriously affects the high-quality development of social economy and the safe operation of major national infrastructure projects. It is urgent to carry out precise monitoring and scientific evaluation of the possibility and scale of geological disasters to provide guarantees for project site selection and layout, regional development planning, and disaster prevention and mitigation. Traditional crustal stability monitoring and evaluation mainly focuses on geological surveys, empirical judgments, statistical analysis, and post-disaster emergency response. There are problems such as single technical indicators, large temporal and spatial scales, and insufficient quantitative evaluation.

[0003] With the development of modern satellite geodetic technology represented by GNSS, InSAR, and gravity, a new generation of data products with high precision, high temporal and spatial resolution, and high currentness have been formed, providing opportunities for refined monitoring of regional crustal stability. The time scale of modern geodetic observations can be high-resolution dynamic information such as every year, every day, every hour, and every second. This is different from seismology (based on seismic waves, with a time scale of seconds) and geology (average of millions of years in a large spatial range). Modern geodetic technology fills the time domain gap between the two, advancing the modern crustal movement research with a time scale of millions of years to a new stage of current crustal movement research with a time scale of several years and decades. It can monitor the deformation caused by crustal rebound after the ice age and earthquake events, which is difficult to do with geological methods. Therefore, there is an urgent need for a multi-dimensional quantitative evaluation method for crustal stability to solve the major and urgent needs of national geographic conditions monitoring and disaster prevention and mitigation. Summary of the invention

[0004] In order to address the shortcomings of the prior art, the present invention aims to provide a multi-dimensional quantitative evaluation method for crustal stability. The method integrates surface deformation and underground activities, modern geodetic measurements and long-term geological data, and establishes a multi-dimensional crustal stability monitoring and evaluation technology system composed of multiple primary indicators and secondary indicator parameters from different time and space scales. It significantly improves the objectivity and scientific rationality of the quantitative evaluation of regional safety along urban active faults, and provides a guarantee for engineering layout site selection and high-quality operation and maintenance.

[0005] In order to achieve the above object, the technical solution adopted by the present invention is:

[0006] A multi-dimensional quantitative evaluation method for crustal stability, the method comprising the following specific steps:

[0007] S1. To construct an evaluation index system, it is necessary to consider the comprehensive impact of internal and external dynamic geological effects;

[0008] S2. Quantification of evaluation indicators: Due to the different data formats of evaluation indicator materials and the different scale accuracy of the sources, it is necessary to unify the data format and quantify the evaluation indicators;

[0009] S3. Weight allocation of evaluation indicators: In order to ensure the rationality of weight allocation and reduce the correlation between evaluation indicators, each evaluation indicator is weighted according to the factor analysis method of multivariate statistics;

[0010] S4. Establish a multi-dimensional comprehensive evaluation model for crustal stability: establish a crustal stability evaluation grading standard based on the evaluation index system constructed in step S1, use multivariate statistical factor analysis to weight each evaluation index, and establish a multi-dimensional comprehensive evaluation model for crustal stability based on the quantified evaluation indexes on the ArcGIS platform;

[0011] S5. Comprehensive evaluation: Determine the comprehensive score of the crustal stability in the study area based on the multi-dimensional comprehensive evaluation model of crustal stability.

[0012] Furthermore, step S1 specifically includes the following steps:

[0013] S1-1. Data collection: including satellite geodetic data GNSS CORS network observations, D-InSAR observation data and gravity data, as well as earthquake, geological data and geophysical data in the study area;

[0014] S1-2. Acquisition of surface deformation information: Acquisition of changes in velocity field, strain field, gravity field, and terrestrial water storage in the study area;

[0015] S1-3. Inversion of underground tectonic activity and establishment of a three-dimensional node grid structure on the fault plane:

[0016] The double-difference positioning method is used to obtain the precise positioning results of medium and small earthquakes in the study area since the seismic phase data was available; its velocity model uses the results of artificial seismic depth sounding and seismic tomography; and the spatial distribution morphology, deep occurrence, and depth of the seismogenic layer of each fault are characterized through the precise positioning profile of small earthquakes, so as to establish a three-dimensional node grid structure of the fault;

[0017] S1-4. Selection of evaluation index factors: Integrating surface deformation and underground activities, geodetic measurements and long-term geological data, a number of index parameters are established from different time and space scales to form a multi-dimensional monitoring and evaluation technology system for crustal stability.

[0018] Furthermore, the strain field in step S1-2 includes maximum shear strain and surface expansion. The shear deformation mainly reflects the tectonic movement and deformation characteristics of fault zone slip or shear, and the surface expansion strain mainly reflects the fault thrust or tension movement. The combination of the two can more comprehensively measure the surface deformation movement characteristics.

[0019] Furthermore, step S1-3 specifically includes:

[0020] S1-3-1. The nodes are distributed along and perpendicular to the isobaths, and each isobath contains the same number of nodes. The long-term fault slip rate obtained by geological methods and the fault slip vector data derived from the earthquake focal mechanism solution are used as prior constraints, and the GNSS-N, E, U deformation and InSAR-LOS deformation of the fault neighbors and surrounding areas are used as model inputs to invert the block rotation angular velocity, block internal strain, and the locking coefficient PHI and slip loss parameters at each node on the fault plane.

[0021] S1-3-2. The locking coefficient, locking depth and sliding loss distribution of the above inversion are comprehensively used as important indicators for determining the different locking coupling degrees and creep characteristics of the fault; the judgment basis is: when the locking coefficient PHI = 0, there is no sliding loss and the locking depth is shallow, indicating that it is in a shallow creep state; when PHI = 1, the sliding loss is large and the locking depth is deep, it means that it is in a strong locking state, indicating that the relative movement of the two plates of the fault is completely converted into strain energy, which is highly dangerous;

[0022] S1-3-3. The seismic moment accumulation rate is calculated by combining the fault locking depth and slip rate as the basis for evaluating the seismic potential of the fault; the seismic moment accumulation rate M is calculated by combining the deep locking and slip rate, that is,

[0023] M=∫μ·D·V·L

[0024] Where μ represents the rock shear modulus, which ranges from 3.0 to 3.3×10 10 N / m 2 ; D is the fault locking depth; V is the slip rate; L is the fault length;

[0025] The later the seismogenic period progresses, the greater the locked depth becomes. At this time, the surface deformation is small, but the calculated seismic moment accumulation rate is still large. This is an important basis for determining how fast the movement of the two plates of the fault is converted into strain energy and the earthquake hazard during the seismogenic period.

[0026] Furthermore, the indicator parameters established in step S1-4 are divided into 7 first-level indicators and 12 second-level indicators according to different factor types, and the factor types include geodetic data, seismic geological data and geophysical data; the first-level indicators in the geodetic data include strain characteristics, underground activity inversion and abnormal changes in terrestrial water reserves, and the corresponding second-level indicators include maximum shear strain rate, surface expansion rate, degree of closure, and equivalent water height; the first-level indicators in the seismic geological data include seismic activity, geological disasters and fault activity, and the corresponding second-level indicators include historical earthquakes, peak seismic acceleration, slope gradient, river system and average annual precipitation; the first-level indicators in the geophysical data include deep geophysics, and the corresponding second-level indicators include Bouguer gravity anomaly values.

[0027] Furthermore, the quantification process of the evaluation index in step S2 includes:

[0028] The maximum shear strain rate, surface expansion rate, equivalent water height, average annual precipitation, Bouguer gravity anomaly, peak seismic acceleration, and slope gradient data were converted into raster data based on the ArcGIS platform and reclassified;

[0029] The degree of closure and the distribution of sliding losses were obtained by extracting the depth values ​​corresponding to the strong closure state (0.8-1.0) to obtain the underground depth map, and then the data was screened, and the closure depth values ​​of the contour area were selected, which were expanded in the plane and interpolated by Kriging, and then converted into raster data based on the ArcGIS platform and reclassified;

[0030] The historical earthquake data was interpolated and analyzed based on ArcGIS spatial analysis functions to generate magnitude distribution maps and stability classification, which were then converted into raster data based on the ArcGIS platform and reclassified.

[0031] River systems and active faults were analyzed and processed based on the Euclidean distance buffer, and were classified and assigned values ​​according to the vertical distance from the fault. They were then converted into raster data and reclassified based on the ArcGIS platform.

[0032] Furthermore, the evaluation index weight distribution calculation steps in step S3 are as follows:

[0033] S3-1. First, classify the influencing factors and classify the influencing factors with high correlation into the same category; if the correlation between factors of different categories is low, each category of factors is regarded as a common factor; determine the weight of the indicator evaluation according to the solved common factor score coefficient;

[0034] Assume that all observed variables are represented by X = [X1 X2 … X p ], there are p variables in total, and there are

[0035] X i =μ i+a i1 f1+a i2 f2+…+a im f m +e i ,

[0036] In the formula, μ i Represents the observed value X i The mean value of j represents the common factor; e i represents a special factor; a ij Represents the observed value X i and the common factor f j The correlation coefficient of ; where i ranges from 1 to p, and j ranges from 1 to m;

[0037] S3-2. Data standardization; X i Converted to a form with a mean of 0 and a variance of 1, the standardized mathematical model is:

[0038] X i =α i1 f1+α i2 f2+…+α im f m +e i ;

[0039] S3-3. Establish a correlation coefficient matrix; let A = (α ij ) P×m , matrix A is the factor loading matrix, and the correlation coefficient matrix R is expressed as R = Α·A T ,but

[0040]

[0041] Among them, the matrix R represents the correlation between the observed variables, r 11 ,…,r pp represents the correlation coefficient of the observed variables, there are p variables in total;

[0042] S3-4. Solve the eigenvalues ​​and eigenvectors of the correlation coefficient matrix and determine the common factors; calculate the eigenvalues ​​of the correlation coefficient matrix R λ=[λ1 λ2 … λ p ], and its corresponding standard orthogonalized eigenvector l i ; where λ1, ..., λ p They represent the eigenvalues ​​corresponding to the eigenvectors respectively;

[0043] S3-5. Calculate the observed value X i and the common factor f j The correlation coefficient α ij , establish the load matrix; use principal component analysis to calculate the correlation coefficient αij , the calculation formula is:

[0044]

[0045] Find α ij , we can get the load matrix A=(α ij ) P×m ; where l ij Represents the element value of the eigenvector of the R matrix;

[0046] S3-6. Solve the factor score coefficient; use the regression method to calculate the factor score coefficient, substitute it into the load matrix A, and solve the common factor score coefficient matrix b j ;

[0047] S3-7. Calculate the common factor scores; obtain the score coefficients of each common factor and calculate the scores:

[0048] scor(m,j)=b j1 X 1m +b j2 X 2m +…+b jp X pm ,

[0049] Finally, the weight of the evaluation index is obtained through the score coefficient; where scor(m,j) represents the score of the common factor; b j1 ,...,b jp represents the element of the score coefficient matrix of the jth column; X 1m ,...,X pm Represents the element corresponding to the mth row of the observed variable.

[0050] Furthermore, in step S4, the weights determined by the factor analysis method of the multivariate analysis in step S3 are combined to establish a multi-dimensional comprehensive evaluation model of crustal stability according to the multi-factor weighted superposition method:

[0051]

[0052] Among them, K is the regional crust stability index, L i is each evaluation index, w i is the weight of each evaluation index.

[0053] The beneficial effects of the present invention are:

[0054] The present invention provides a multi-dimensional quantitative evaluation method for crustal stability. It uses the high-density observation data of the global navigation satellite system (GNSS) autonomous networking and station building as a constraint, combines underground fault activity with GNSS monitoring data of surface deformation, analyzes the crustal stability of the study area from seven aspects, namely, strain characteristics, underground activity inversion, abnormal changes in terrestrial water reserves, seismic activity, geological disasters, fault activity, and deep geophysical characteristics, and establishes a multi-dimensional monitoring and evaluation system for crustal stability. The weights of the evaluation indicators are established according to the factor analysis method of multivariate analysis, and a multi-dimensional quantitative evaluation model for crustal stability is constructed in combination with the multi-factor weighted superposition method, which avoids the subjectivity of the weight determination process and makes the evaluation model more scientific. BRIEF DESCRIPTION OF THE DRAWINGS

[0055] Figure 1 This is a flow chart of the multi-dimensional quantitative evaluation method of the present invention. DETAILED DESCRIPTION

[0056] The principles and features of the present invention are described below in conjunction with the accompanying drawings. The examples given are only used to explain the present invention and are not intended to limit the scope of use of the present invention.

[0057] like Figure 1 As shown, the present invention proposes a multi-dimensional quantitative evaluation method for crustal stability, which includes the following specific steps:

[0058] S1. Construct an evaluation index system; regional crust stability evaluation is a multi-factor, multi-angle comprehensive evaluation process, which needs to consider the comprehensive impact of internal and external dynamic geological effects. Specifically, it includes the following steps:

[0059] S1-1. Data collection: including satellite geodetic data GNSS CORS network observations, D-InSAR observation data and gravity data, as well as earthquake, geological data and geophysical data in the study area;

[0060] S1-2. Acquisition of surface deformation information: Obtain changes in velocity field, strain field, gravity field, and terrestrial water storage in the study area. The strain field includes maximum shear strain and surface expansion. Shear deformation mainly reflects the tectonic movement and deformation characteristics of fault zone slip or shear, and surface expansion strain mainly reflects fault thrust or tensile movement. The combination of the two can more comprehensively measure the characteristics of surface deformation movement.

[0061] S1-3. Inversion of underground tectonic activity and establishment of a three-dimensional node grid structure on the fault plane.

[0062] The double-difference positioning method is used to obtain the precise positioning results of medium and small earthquakes in the study area since the seismic phase data was available. Its velocity model uses the results of artificial seismic depth sounding and seismic tomography; and the spatial distribution morphology, deep occurrence, depth of the seismogenic layer, etc. of each fault are characterized through the precise positioning profile of small earthquakes, so as to establish a three-dimensional node grid structure of the fault. Specifically, it includes:

[0063] S1-3-1. The nodes are distributed along and perpendicular to the isobaths, and each isobath contains the same number of nodes. The long-term fault slip rate obtained by geological methods and the fault slip vector data derived from the earthquake focal mechanism solution are used as prior constraints, and the GNSS-N, E, U deformation and InSAR-LOS deformation of the fault neighbors and surrounding areas are used as model inputs to invert the block rotation angular velocity, internal strain of the block, and the locking coefficient PHI and slip loss parameters at each node on the fault plane.

[0064] S1-3-2. The locking coefficient, locking depth and sliding loss distribution of the above inversion are comprehensively used as important indicators for determining the different locking coupling degrees and creep characteristics of the fault. The judgment basis is: when PHI = 0, there is no sliding loss and the locking depth is shallow, indicating that it is in a shallow creep state; when PHI = 1, the sliding loss is large and the locking depth is deep, it means that it is in a strong locking state, indicating that the relative movement of the two plates of the fault is completely converted into strain energy, which is highly dangerous.

[0065] S1-3-3. Combine the fault locking depth and slip rate to calculate the seismic moment accumulation rate as a basis for evaluating the seismic potential of the fault. Combine the deep locking and slip rate to calculate the seismic moment accumulation rate M, that is,

[0066] M=∫μ·D·V·L

[0067] Where μ represents the rock shear modulus, which ranges from 3.0 to 3.3×10 10 N / m 2 ; D represents the fault locking depth; V represents the slip rate; L represents the fault length.

[0068] The later the seismogenic period is, the greater the locked depth is. Although the surface deformation is small at this time, the calculated seismic moment accumulation rate is still large. This is an important basis for determining the speed of the conversion of the two-plate movement of the fault into strain energy and the earthquake hazard during the seismogenic period.

[0069] S1-4. Selection of evaluation index factors: Based on the surface deformation and underground activities, modern geodetic measurements and long-term geological data, a multi-dimensional monitoring and evaluation technology system for crustal stability consisting of 7 primary indexes and 12 secondary index parameters is established from different time and space scales, as shown in Table 1.

[0070] Table 1 Evaluation indexes of regional crustal stability

[0071]

[0072] S2. Quantification of evaluation indicators; Due to the different data formats of evaluation indicator materials and the different scale accuracy of the sources, it is necessary to unify the data format and quantify the evaluation indicators. The original data should be standardized before determining the weight to eliminate the differences in magnitude and dimension between different indicators. The KMO test (Kaiser-Meyer-Olkin Measure of Sampling Adequacy, a statistical test method for whether the sampling data is suitable for factor analysis) of the evaluation indicator must be greater than or equal to 0.6, otherwise it is not suitable for factor analysis. The specific implementation operations are as follows:

[0073] The data of maximum shear strain rate, surface expansion rate, equivalent water height, average annual precipitation, Bouguer gravity anomaly, peak seismic acceleration, and slope gradient were converted into raster data based on the ArcGIS platform and reclassified according to the crustal stability evaluation classification standard (as shown in Table 2).

[0074] The locking degree and sliding loss distribution were analyzed by extracting the depth values ​​corresponding to the strong locking state (0.8-1.0) to obtain the underground depth diagram, screening the data, selecting the locking depth values ​​of the contour area, unfolding them in the plane and performing Kriging interpolation. They were then converted into raster data based on the ArcGIS platform and reclassified according to the crustal stability evaluation classification standard.

[0075] The historical earthquake data were interpolated and analyzed based on the ArcGIS spatial analysis function to generate magnitude distribution maps and perform stability classification. The data were then converted into raster data based on the ArcGIS platform and reclassified according to the crustal stability evaluation classification standards.

[0076] River systems and active faults are analyzed and processed based on the Euclidean distance buffer, and are classified and assigned values ​​according to the vertical distance from the fault. They are then converted into raster data based on the ArcGIS platform and reclassified according to the crustal stability evaluation classification standard.

[0077] Table 2 Crustal stability evaluation classification standards

[0078]

[0079] The code column in the above table indicates the code corresponding to each evaluation indicator.

[0080] S3. Weight allocation of evaluation indicators: How to fully consider the influence of various factors to obtain reliable comprehensive evaluation results is the key to regional crustal stability evaluation; in order to ensure the rationality of weight allocation and reduce the correlation between evaluation indicators, each evaluation indicator is weighted according to the factor analysis method of multivariate statistics. The calculation steps of evaluation indicator weight allocation are as follows:

[0081] S3-1. First, classify the influencing factors and group the factors with high correlation into the same category. If the correlation between factors of different categories is low, each category of factors is regarded as a common factor. The weight of the indicator evaluation is determined according to the solved common factor score coefficient.

[0082] Assume that all observed variables are represented by X = [X1 X2…X p ], there are p variables in total, and there are

[0083] X i =μ i +a i1 f1+a i2 f2+…+a im f m +e i ,

[0084] In the formula, μ i Represents the observed value X i The mean value of j represents the common factor; e i represents a special factor; a ij Represents the observed value X i and the common factor f j The correlation coefficient of ; where i ranges from 1 to p, and j ranges from 1 to m.

[0085] S3-2. Data standardization; X i Converted to a form with a mean of 0 and a variance of 1, the standardized mathematical model is:

[0086] X i =α i1 f1+α i2 f2+…+α im f m +e i ;

[0087] S3-3. Establish a correlation coefficient matrix; let A = (α ij ) P×m , matrix A is the factor loading matrix, and the correlation coefficient matrix R is expressed as R = Α·A T ,but

[0088]

[0089] Among them, the matrix R represents the correlation between the observed variables, r 11 ,…,r pp represents the correlation coefficient of the observed variables, there are p variables in total;

[0090] S3-4. Solve the eigenvalues ​​and eigenvectors of the correlation coefficient matrix and determine the common factors; calculate the eigenvalues ​​of the correlation coefficient matrix R λ=[λ1λ2…λ p ], and its corresponding standard orthogonalized eigenvector l i ; where λ1, ..., λ p They represent the eigenvalues ​​corresponding to the eigenvectors respectively.

[0091] S3-5. Calculate the observed value X i and the common factor f j The correlation coefficient α ij , establish the load matrix; use principal component analysis to calculate the correlation coefficient α ij , the calculation formula is:

[0092]

[0093] Find α ij , we can get the load matrix A=(α ij ) P×m ; where l ij Represents the element values ​​of the eigenvector of the R matrix.

[0094] S3-6. Solve the factor score coefficient; use the regression method to calculate the factor score coefficient, substitute it into the load matrix A, and solve the common factor score coefficient matrix b j ;

[0095] S3-7. Calculate the common factor scores; obtain the score coefficients of each common factor and calculate the scores:

[0096] scor(m,j)=b j1 X 1m +b j2 X 2m +…+b jp X pm ,

[0097] Finally, the weight of the evaluation index is obtained through the score coefficient. In the formula, scor(m,j) represents the score of the common factor; b j1 ,...,b jp represents the element of the score coefficient matrix of the jth column; X 1m ,...,X pm Represents the element corresponding to the mth row of the observed variable.

[0098] S4. Establish a multi-dimensional comprehensive evaluation model for crustal stability: establish a crustal stability evaluation grading standard based on the evaluation index system constructed in step S1, use multivariate statistical factor analysis to weight each evaluation index, and establish a multi-dimensional comprehensive evaluation model for crustal stability based on the quantified evaluation indexes on the ArcGIS platform;

[0099] Combined with the weights determined by the factor analysis method of multivariate analysis in step S3, a multi-dimensional comprehensive evaluation model for crustal stability is established based on the multi-factor weighted superposition method:

[0100]

[0101] Among them, K is the regional crust stability index, L i is each evaluation index, w i is the weight of each evaluation index.

[0102] S5. Comprehensive evaluation: Determine the comprehensive score of the crustal stability in the study area based on the multi-dimensional comprehensive evaluation model of crustal stability.

[0103] The multi-dimensional quantitative evaluation method for crustal stability provided by the present invention is based on the high-density observation data of the global navigation satellite system (GNSS) autonomous networking and station building, and combines the underground fault activity with the GNSS monitoring data of surface deformation. The crustal stability of the study area is analyzed from seven aspects, namely, strain characteristics, underground activity inversion, abnormal changes in terrestrial water reserves, seismic activity, geological disasters, fault activity, and deep geophysical characteristics, and a multi-dimensional monitoring and evaluation system for crustal stability is established. The weights of the evaluation indicators are established according to the factor analysis method of multivariate analysis, and a multi-dimensional quantitative evaluation model for crustal stability is constructed in combination with the multi-factor weighted superposition method, which avoids the subjectivity of the weight determination process and makes the evaluation model more scientific.

[0104] Obviously, the embodiments described above are only some embodiments of the present application, rather than all embodiments. The preferred embodiments of the present application are given in the accompanying drawings, but they do not limit the patent scope of the present application. The present application can be implemented in many different forms. On the contrary, the purpose of providing these embodiments is to make the understanding of the disclosure of the present application more thorough and comprehensive. Although the present application is described in detail with reference to the aforementioned embodiments, for those skilled in the art, it is still possible to modify the technical solutions recorded in the aforementioned specific implementation methods, or to perform equivalent replacement of some of the technical features therein. Any equivalent structure made using the contents of the specification and drawings of this application, directly or indirectly used in other related technical fields, is similarly within the scope of patent protection of this application.

Claims

1. A multi-dimensional quantitative evaluation method for crustal stability, characterized in that: The method comprises the following specific steps: S1. To construct an evaluation index system, it is necessary to consider the comprehensive impact of internal and external dynamic geological effects; S2. Quantification of evaluation indicators: Due to the different data formats of evaluation indicator materials and the different scale accuracy of the sources, it is necessary to unify the data format and quantify the evaluation indicators; S3. Weight allocation of evaluation indicators: In order to ensure the rationality of weight allocation and reduce the correlation between evaluation indicators, each evaluation indicator is weighted according to the factor analysis method of multivariate statistics; S4. Establish a multi-dimensional comprehensive evaluation model for crustal stability: establish a crustal stability evaluation grading standard based on the evaluation index system constructed in step S1, use multivariate statistical factor analysis to weight each evaluation index, and establish a multi-dimensional comprehensive evaluation model for crustal stability based on the quantified evaluation indexes on the ArcGIS platform; S5. Comprehensive evaluation: Determine the comprehensive score of the crustal stability in the study area based on the multi-dimensional comprehensive evaluation model of crustal stability.

2. The multi-dimensional quantitative evaluation method for crustal stability according to claim 1, characterized in that: Step S1 specifically includes the following steps: S1-1. Data collection: including satellite geodetic data GNSS CORS network observations, D-InSAR observation data and gravity data, as well as earthquake, geological data and geophysical data in the study area; S1-2. Acquisition of surface deformation information: Acquisition of changes in velocity field, strain field, gravity field, and terrestrial water storage in the study area; S1-3. Inversion of underground tectonic activity and establishment of a three-dimensional node grid structure on the fault plane: The double-difference positioning method is used to obtain the precise positioning results of medium and small earthquakes in the study area since the seismic phase data was available; its velocity model uses the results of artificial seismic depth sounding and seismic tomography; and the spatial distribution morphology, deep occurrence, and depth of the seismogenic layer of each fault are characterized through the precise positioning profile of small earthquakes, so as to establish a three-dimensional node grid structure of the fault; S1-4. Selection of evaluation index factors: Integrating surface deformation and underground activities, geodetic measurements and long-term geological data, a number of index parameters are established from different time and space scales to form a multi-dimensional monitoring and evaluation technology system for crustal stability.

3. The multi-dimensional quantitative evaluation method for crustal stability according to claim 2, characterized in that: The strain field in step S1-2 includes the maximum shear strain and surface expansion. The shear deformation mainly reflects the tectonic movement and deformation characteristics of fault zone slip or shear, and the surface expansion strain mainly reflects the fault thrust or tension movement. The combination of the two can more comprehensively measure the surface deformation movement characteristics.

4. The multi-dimensional quantitative evaluation method for crustal stability according to claim 2, characterized in that: Step S1-3 specifically includes: S1-3-1. The nodes are distributed along and perpendicular to the isobaths, and each isobath contains the same number of nodes. The long-term fault slip rate obtained by geological methods and the fault slip vector data derived from the earthquake focal mechanism solution are used as prior constraints, and the GNSS-N, E, U deformation and InSAR-LOS deformation of the fault neighbors and surrounding areas are used as model inputs to invert the block rotation angular velocity, block internal strain, and the locking coefficient PHI and slip loss parameters at each node on the fault plane. S1-3-2. The locking coefficient, locking depth and sliding loss distribution of the above inversion are comprehensively used as important indicators for determining the different locking coupling degrees and creep characteristics of the fault; the judgment basis is: when the locking coefficient PHI = 0, there is no sliding loss and the locking depth is shallow, indicating that it is in a shallow creep state; when PHI = 1, the sliding loss is large and the locking depth is deep, it means that it is in a strong locking state, indicating that the relative movement of the two plates of the fault is completely converted into strain energy, which is highly dangerous; S1-3-3. The seismic moment accumulation rate is calculated by combining the fault locking depth and slip rate as the basis for evaluating the seismic potential of the fault; the seismic moment accumulation rate M is calculated by combining the deep locking and slip rate, that is, M=∫μ·D·V·L Where μ represents the rock shear modulus, which ranges from 3.0 to 3.3×10 10 N / m 2 ; D is the fault locking depth; V is the slip rate; L is the fault length; The later the seismogenic period progresses, the greater the locked depth becomes. At this time, the surface deformation is small, but the calculated seismic moment accumulation rate is still large. This is an important basis for determining how fast the movement of the two plates of the fault is converted into strain energy and the earthquake hazard during the seismogenic period.

5. The multi-dimensional quantitative evaluation method for crustal stability according to claim 2, characterized in that: The index parameters established in step S1-4 are divided into 7 primary indicators and 12 secondary indicators according to different factor types, and the factor types include geodetic data, seismic geological data and geophysical data; the primary indicators in geodetic data include strain characteristics, underground activity inversion and abnormal changes in terrestrial water reserves, and the corresponding secondary indicators include maximum shear strain rate, surface expansion rate, degree of closure, and equivalent water height; the primary indicators in seismic geological data include seismic activity, geological disasters and fault activity, and the corresponding secondary indicators include historical earthquakes, peak seismic acceleration, slope gradient, river system and average annual precipitation; the primary indicators in geophysical data include deep geophysics, and the corresponding secondary indicators include Bouguer gravity anomaly.

6. The multi-dimensional quantitative evaluation method for crustal stability according to claim 5, characterized in that: The evaluation index quantification process in step S2 includes: The maximum shear strain rate, surface expansion rate, equivalent water height, average annual precipitation, Bouguer gravity anomaly, peak seismic acceleration, and slope gradient data were converted into raster data based on the ArcGIS platform and reclassified; The degree of closure and the distribution of sliding losses were obtained by extracting the depth values ​​corresponding to the strong closure state (0.8-1.0) to obtain the underground depth map, and then the data was screened, and the closure depth values ​​of the contour area were selected, which were expanded in the plane and interpolated by Kriging, and then converted into raster data based on the ArcGIS platform and reclassified; The historical earthquake data was interpolated and analyzed based on ArcGIS spatial analysis functions to generate magnitude distribution maps and stability classification, which were then converted into raster data based on the ArcGIS platform and reclassified. River systems and active faults were analyzed and processed based on the Euclidean distance buffer, and were classified and assigned values ​​according to the vertical distance from the fault. They were then converted into raster data and reclassified based on the ArcGIS platform.

7. The multi-dimensional quantitative evaluation method for crustal stability according to claim 1, characterized in that: The calculation steps for evaluation index weight distribution in step S3 are as follows: S3-1. First, classify the influencing factors and classify the influencing factors with high correlation into the same category; if the correlation between factors of different categories is low, each category of factors is regarded as a common factor; determine the weight of the indicator evaluation according to the solved common factor score coefficient; Assume that all observed variables are represented by X = [X1 X2 … X p ], there are p variables in total, and there are X i =μ i +a i1 f1+a i2 f2+…+a im f m +e i , In the formula, μ i Represents the observed value X i The mean value of j represents the common factor; e i represents a special factor; a ij Represents the observed value X i and the common factor f j The correlation coefficient of ; where i ranges from 1 to p, and j ranges from 1 to m; S3-2. Data standardization; X i Converted to a form with a mean of 0 and a variance of 1, the standardized mathematical model is: X i =a i1 f1+a i2 f2+…+a im f m +e i ; S3-3. Establish a correlation coefficient matrix; let A = (α ij ) P×m , matrix A is the factor loading matrix, and the correlation coefficient matrix R is expressed as R = Α·A T ,but Among them, the matrix R represents the correlation between the observed variables, r 11 ,…,r pp represents the correlation coefficient of the observed variables, there are p variables in total; S3-4. Solve the eigenvalues ​​and eigenvectors of the correlation coefficient matrix and determine the common factors; calculate the eigenvalues ​​of the correlation coefficient matrix R λ=[λ1 λ2 … λ p ], and its corresponding standard orthogonalized eigenvector l i ; where λ1, ..., λ p They represent the eigenvalues ​​corresponding to the eigenvectors respectively; S3-5. Calculate the observed value X i and the common factor f j The correlation coefficient α ij , establish the load matrix; use principal component analysis to calculate the correlation coefficient α ij , the calculation formula is: Find α ij , we can get the load matrix A=(α ij ) P×m ; where l ij Represents the element value of the eigenvector of the R matrix; S3-6. Solve the factor score coefficient; use the regression method to calculate the factor score coefficient, substitute it into the load matrix A, and solve the common factor score coefficient matrix b j ; S3-7. Calculate the common factor scores; obtain the score coefficients of each common factor and calculate the scores: scor(m,j)=b j1 X 1m +b j2 X 2m +…+b jp X pm , Finally, the weight of the evaluation index is obtained through the score coefficient; where scor(m,j) represents the score of the common factor; b j1 ,...,b jp represents the j-th column score coefficient matrix element; X 1m ,...,X pm Represents the element corresponding to the mth row of the observed variable.

8. The multi-dimensional quantitative evaluation method for crustal stability according to claim 7, characterized in that: In step S4, the weights determined by the factor analysis method of the multivariate analysis in step S3 are combined to establish a multi-dimensional comprehensive evaluation model for crustal stability according to the multi-factor weighted superposition method: Among them, K is the regional crust stability index, L i is each evaluation index, w i is the weight of each evaluation index.

Images