Focus parameter prediction method based on physical constraint and array geometric adaptation

Abstract

The invention relates to a seismic source parameter prediction method based on physical constraint and array geometric adaptation, which comprises the following steps of: determining input, output and initial prediction models, sequentially obtaining station characteristics and predicted seismic source parameters, comparing the predicted seismic source parameters with a seismic source true value to obtain regression loss (if the seismic source true value does not exist, the regression loss is not recorded); obtaining phase loss according to the predicted seismic source parameter and the seismic phase label, and then obtaining total loss; carrying out a training process on the model according to the previous steps, and storing an optimal model; inputting related data into the optimal model, and obtaining predicted seismic source parameters; the phase loss only needs the arrival time data of the P wave and the S wave, and is the signal feature which is the most easily obtained signal feature, so that the application range is relatively wide; phase loss takes the physical law of seismic wave propagation as a hard constraint to be embedded into a model optimization process, so that the interpretability of a prediction result is improved; therefore, the method is wide in application range and high in interpretability when used for predicting the seismic source parameters.

Description

Technical Field

[0001] This invention relates to a method for predicting earthquake source parameters, belonging to the field of earthquake monitoring and seismology technology, and particularly to a method for predicting earthquake source parameters based on physical constraints and array geometry adaptation. Background Technology

[0002] Predicting earthquake source parameters, i.e. earthquake location, is a core task of earthquake monitoring and geophysics. It not only enables earthquake hazard analysis and disaster risk prevention, but also helps in the study of earthquake mechanisms and the internal structure of the Earth.

[0003] Chinese patent application number 202410952406.3, filed on July 16, 2024, discloses a method and apparatus for predicting earthquake source parameters based on deep learning. The method includes: acquiring earthquake event waveforms, P-wave arrival times, and epicentral distances from a target station; preprocessing the earthquake event waveforms based on the P-wave arrival times to obtain T-component earthquake event waveforms; constructing an auxiliary waveform with the same length as the T-component earthquake event waveform; wherein the auxiliary waveforms include waveforms corresponding to the epicentral distance, waveforms corresponding to the maximum amplitude of the T-component earthquake event waveforms, and waveforms corresponding to the station information of the target station; inputting the T-component earthquake event waveforms and auxiliary waveforms into a preset source parameter prediction model, and having the source parameter prediction model output the predicted source parameters. Although this design can obtain source parameters, it still has the following drawbacks:

[0004] In this design, the core input is the amplitude, phase, and spectral information of the complete waveform, and it is highly dependent on the accuracy of the waveform data. It is often difficult to obtain high-precision waveform data in microseismic monitoring, sparse seismic networks, or far-field recording. If the accuracy of the waveform data is not high, it will lead to low prediction accuracy of the source parameters. Therefore, the application scope of this design is relatively narrow.

[0005] The information disclosed in this background section is intended only to enhance understanding of the overall background of this application and should not be construed as an admission or in any way implying that the information constitutes prior art known to those skilled in the art. Summary of the Invention

[0006] The purpose of this invention is to overcome the shortcomings and problems of the narrow application range of the existing technology in predicting source parameters, and to provide a source parameter prediction method based on physical constraints and array geometry adaptation with a wider application range.

[0007] To achieve the above objectives, the technical solution of the present invention is:

[0008] A method for predicting seismic source parameters based on physical constraints and array geometry adaptation, the method comprising the following steps:

[0009] Step 1: First, determine the continuous waveform data and location information corresponding to multiple stations in the dataset as input, then determine the predicted source parameters as output, and then set an initial prediction model;

[0010] The second step is to extract features from the continuous waveform data and location information to obtain the predicted source parameters.

[0011] Step 3: First determine whether the dataset includes earthquake phase labels and the true value of the earthquake source. If the dataset includes earthquake phase labels and the true value of the earthquake source, proceed to step 4. If the dataset does not include earthquake phase labels and the true value of the earthquake source, proceed to step 5.

[0012] Step 4: First, compare the predicted source parameters with the actual source values ​​to obtain the regression loss. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, weight the regression loss and phase loss to obtain the total loss. Then, proceed to step 6.

[0013] Step 5: First, obtain the phase labels based on the dataset using the phase picking method. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, use the phase loss as the total loss and then proceed to step 6.

[0014] Step 6: First, train the initial prediction model, then select the initial prediction model with the minimum total loss as the optimal model.

[0015] Step 7: First, input the continuous waveform data and location information of multiple stations into the optimal model, and then obtain the predicted source parameters from the optimal model.

[0016] In the fourth and fifth steps, the actual arrival data refers to the actual P-wave arrival data and the actual S-wave arrival data, and the predicted arrival data refers to the predicted P-wave arrival data and the predicted S-wave arrival data.

[0017] The calculation steps for the predicted P-wave arrival time data and the predicted S-wave arrival time data are as follows:

[0018] First, calculate the epicentral distance:

[0019] ;

[0020] In the formula: The epicentral distance is... The average radius of the Earth This is the angular distance between the seismic source and the station;

[0021] Then calculate the distance from the epicenter:

[0022] ;

[0023] In the formula: The distance from the epicenter, The focal depth;

[0024] Then calculate the theoretical time-to-run time:

[0025] ; ;

[0026] In the formula: To predict the arrival time of the P-wave, To predict the arrival time of the S-wave, The distance from the epicenter, For P-wave velocity, This represents the S-wave velocity.

[0027] In steps four and five, obtaining the phase loss based on the predicted arrival time data and the actual arrival time data refers to:

[0028] First, calculate the phase loss of the predicted P-wave arrival data:

[0029] ;

[0030] In the formula: The phase loss of the P-wave arrival data. For all stations, To predict the arrival time of the P-wave, This is the actual arrival time data for the P-wave. The confidence level for picking up seismic phases;

[0031] Next, calculate the phase loss of the data when the predicted S-wave arrives:

[0032] ;

[0033] In the formula: The phase loss of the data when the S-wave arrives. For all stations, To predict the arrival time of the S-wave, This is the actual arrival time data for the S-wave. The confidence level for picking up seismic phases;

[0034] The phase loss is then obtained by averaging the phase loss of the predicted P-wave arrival data and the phase loss of the predicted S-wave arrival data.

[0035] ;

[0036] In the formula, This is for phase loss.

[0037] In the fourth step, the regression loss obtained by comparing the predicted source parameters with the actual source values ​​refers to:

[0038] ;

[0039] In the formula, To regress the loss, To predict earthquake source parameters, This represents the true value of the earthquake source.

[0040] In the fourth step, the weighted average of regression loss and phase loss to obtain the total loss refers to:

[0041] ;

[0042] ; ;

[0043] In the formula, For the total loss, The weights for phase loss, The weights for the regression loss;

[0044] In the fifth step, referring to the phase loss as the total loss means:

[0045] ;

[0046] In the formula, This represents the total loss.

[0047] In the second step, feature extraction from continuous waveform data and location information to obtain predicted source parameters refers to: first, performing feature extraction for a single station. The feature extraction process for a single station involves: obtaining the local time-frequency features of the station based on its continuous waveform data; then obtaining the geometric distribution features of the station based on its location information; and finally fusing the geometric distribution features with the local time-frequency features to obtain the station's features. This process is repeated for each single station until the features of all stations are obtained.

[0048] Then, the characteristics of all stations are averaged to obtain the global event characteristics, and the predicted source parameters are obtained based on the global event characteristics.

[0049] In the second step, obtaining the local time-frequency features of the station based on the continuous waveform data of the station means: using a CNN convolutional neural network to obtain the local time-frequency features based on the continuous waveform data of the station;

[0050] In the second step, obtaining the geometric distribution features of a station based on its location information means: mapping the latitude and longitude coordinates of all stations to a multidimensional space through a linear layer and activation function to obtain the initial location features of all stations; calculating the distance information of all stations based on the initial location features of all stations; calculating a distance matrix based on all the distance information; defining the edge weights of the distance matrix using a Gaussian kernel; generating graph structure data based on the edge weights; and extracting the geometric distribution features of all stations from the graph structure data using a GNN graph neural network.

[0051] In the second step, calculating the distance information of all stations based on the initial location features of all stations means:

[0052] ;

[0053] ;

[0054] ; ;

[0055] In the formula, For the first The station and the first The surface distance between each station The average radius of the Earth No. The station and the first angular distance between stations It is the first The latitude (in arcs) of each station. It is the first Longitude (in radians) of each station. It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station;

[0056] In the second step, defining the edge weights of the distance matrix using a Gaussian kernel and then generating graph structure data based on those edge weights refers to:

[0057] ;

[0058] In the formula, For graph-structured data, Distance threshold For width parameters;

[0059] In the second step, fusing the geometric distribution characteristics of the station with local time-frequency characteristics to obtain the station characteristics refers to:

[0060] ;

[0061] ;

[0062] ;

[0063] In the formula, Features of the station It is a local time-frequency feature. For parameter items, It exhibits geometric distribution characteristics.

[0064] In the third step, obtaining the predicted source parameters based on global event features means: using a multilayer perceptron (MLP) to perform nonlinear transformation and mapping on the global event features to obtain high-dimensional information, and then outputting the high-dimensional information as the predicted source parameters.

[0065] In the sixth step, the phase picking method includes: PhaseNet method, STA / LTA method, or EQtransform method.

[0066] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0067] 1. The present invention discloses a method for predicting seismic source parameters based on physical constraints and array geometry adaptation. The method includes the following steps: Step 1: First, determine continuous waveform data and location information as input, then determine the predicted seismic source parameters as output, and then set an initial prediction model; Step 2: First, obtain the predicted seismic source parameters based on the continuous waveform data and location information; Step 3: First, determine whether the dataset includes phase labels and the true seismic source value. If the dataset includes phase labels and the true seismic source value, proceed to Step 4; if the dataset does not include phase labels and the true seismic source value, proceed to Step 5; Step 4: First, compare the predicted seismic source parameters with the true seismic source value to obtain the regression loss, and then obtain the true arrival time based on the phase labels. The phase loss is obtained by comparing the predicted arrival time data obtained based on the predicted source parameters, and then the total loss is obtained by weighting the regression loss and phase loss. Then, proceed to step six. Step five: First, obtain the phase labels based on the dataset using the phase picking method. Then, compare the actual arrival time data obtained based on the phase labels with the predicted arrival time data obtained based on the predicted source parameters to obtain the phase loss. Then, use the phase loss as the total loss and proceed to step six. Step six: Train the model according to steps two to five, and select the initial prediction model with the minimum total loss on the validation set as the optimal model. Step seven: Input continuous waveform data and location information into the optimal model to obtain the predicted source parameters. The advantages of this invention also include:

[0068] Firstly, this invention constructs a unified deep learning model that directly translates raw continuous waveform data and station location information into source parameters. Using this model, predicted source parameters can be obtained.

[0069] Secondly, phase loss calculation only requires the arrival time data of P-waves and S-waves, which are the most basic, reliable and easy-to-obtain signal characteristics in earthquake monitoring. Even when the signal-to-noise ratio is low and the waveform is incomplete, the phase arrival time is relatively reliable. Therefore, compared with existing technologies, it has a wider range of applications.

[0070] Thirdly, introducing phase loss into the training loss forces the model to predict the source location so that its calculated theoretical travel time matches the observed P-wave and S-wave arrival times as closely as possible. In essence, this embeds the physical laws of seismic wave propagation as hard constraints into the model optimization process, making the model's prediction no longer a pure "data fitting" but one with a clear physical causal relationship. This greatly improves the physical rationality and interpretability of the prediction results.

[0071] Fourthly, a dual training path (i.e., the fifth and sixth steps) is designed, which allows for flexible selection of the fully supervised mode of "regression loss + phase loss" or the weakly supervised or self-supervised mode of "phase loss" depending on whether the data has a "true source value" label.

[0072] Therefore, this invention has a wide range of applications and strong interpretability when predicting earthquake source parameters.

[0073] 2. In the method for predicting source parameters based on physical constraints and array geometry adaptation of the present invention, the calculation steps for the predicted P-wave arrival time data and the predicted S-wave arrival time data in the fourth and fifth steps are as follows: first calculate the epicentral distance, then calculate the source distance, and then calculate the theoretical travel time, thus obtaining the predicted P-wave arrival time data and the predicted S-wave arrival time data; obtaining the phase loss based on the predicted arrival time data and the actual arrival time data means: first calculate the phase loss of the predicted P-wave arrival time data, then calculate the phase loss of the predicted S-wave arrival time data, and then average the phase loss of the predicted P-wave arrival time data and the phase loss of the predicted S-wave arrival time data to obtain the phase loss; in the fourth step, obtaining the regression loss by comparing the predicted source parameters with the actual source values ​​means: substituting the predicted source parameters and the actual source values ​​into a specified formula to calculate the regression loss. In the fourth step, the regression loss and phase loss are weighted to obtain the total loss. In application, the constraints of the fourth and fifth steps ensure that the predicted arrival times of the predicted source parameters match the actual observation times of a large number of stations. This conforms to the basic verification logic of seismology, making the results inherently physically reasonable. Since some real observation data lack precise "true source value" labels, this method automatically picks up the true source values ​​through phase picking, thereby transforming massive amounts of unlabeled data into effective training resources. This greatly reduces the model's dependence on expensive manually labeled data and effectively improves the model's generalization ability. The final total loss is a weighted sum of the regression loss (data fitting term) and the phase loss (physical constraint term). This design allows the model optimization process to simultaneously find a balance between "learning from labeled data" and "obeying physical laws." Therefore, this invention has good generalization ability.

[0074] 3. In the method for predicting seismic source parameters based on physical constraints and array geometry adaptation of the present invention, the second step of extracting features from continuous waveform data and location information to obtain predicted seismic source parameters refers to: firstly, using a CNN convolutional neural network to obtain local time-frequency features from the continuous waveform data of the stations; then, mapping the latitude and longitude coordinates of all stations to a multidimensional space through a linear layer and activation function to obtain the initial location features of all stations; then, calculating the distance information of all stations based on the initial location features of all stations; then, calculating a distance matrix based on all the distance information; then, defining the edge weights of the distance matrix using a Gaussian kernel; then, generating graph structure data based on the edge weights; then, extracting the geometric distribution features of all stations from the graph structure data using a GNN graph neural network; and then fusing the local time-frequency features with the geometric distribution features to obtain the station features; repeating the aforementioned steps until the station features of all stations are obtained. In application, the CNN is used directly... By learning local time-frequency features from the original continuous waveform, the model can automatically learn effective waveform patterns related to the seismic source, effectively capturing the local morphology and sequence dependencies of P-waves, S-waves, and other seismic phases, which is crucial for determining seismic events. The Geometric Node Network (GNN) can explicitly learn and aggregate spatial dependencies between stations, regardless of whether the stations are densely or sparsely distributed, or exhibit special geometries such as linear or network structures. In sparse or poorly covered network structures, traditional methods can lead to inaccurate positioning, while the GNN can more robustly integrate information from a limited number of stations through finite but well-defined connections in the graph, suppressing the resulting ambiguity and significantly improving positioning accuracy and stability in non-ideal network conditions. Local time-frequency features tell the model the "content" of the seismic signal, while the geometric distribution characteristics of the stations tell the model the spatial context in which the signal was observed. The combination of these two provides the model with the most comprehensive information necessary for determining the seismic source location. Therefore, this invention is adaptable to various network distributions.

[0075] 4. In the source parameter prediction method based on physical constraints and array geometry adaptation of the present invention, the second step of obtaining the predicted source parameters based on global event features refers to: using a multilayer perceptron (MLP) to perform nonlinear transformation and mapping on the global event features to obtain high-dimensional information, and then outputting the high-dimensional information as the predicted source parameters. In application, the input of the multilayer perceptron (MLP) is a single "global event feature" obtained after global average pooling. This feature has already fused the information of all stations. This approach forces the model to make a final judgment based on the overall consensus of information from all stations, rather than relying on the local features of one or a few stations. This effectively avoids the problem of serious deviations in prediction results due to poor data quality from individual stations (such as high noise or instrument failure), and significantly improves the robustness and anti-interference ability of the model. Therefore, the present invention has strong robustness.

[0076] 5. In the source parameter prediction method based on physical constraints and array geometry adaptation of the present invention, the sixth step, the phase picking method includes: PhaseNet method, STA / LTA method, or EQtransform method. In application, PhaseNet and EQtransform are advanced AI-based picking methods, while STA / LTA is a classic and reliable algorithm, widely verified in production practice, capable of providing phase arrival times with high accuracy and confidence. As a "physical supervision signal," it is relatively reliable, and the "phase loss" calculated using it can stably and effectively guide the model to learn the true wave propagation law, ensuring the effectiveness of the physical constraints and avoiding model training divergence or failure due to poor supervision signal quality. Therefore, the physical constraints of the present invention are stable. Attached Figure Description

[0077] Figure 1 This is a schematic diagram of the structure of the present invention.

[0078] Figure 2 This is a schematic diagram of the structure of Example 2.

[0079] Figure 3 This is a schematic diagram of the structure of Example 3. Detailed Implementation

[0080] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0081] Please see Figure 1 — Figure 3 A method for predicting seismic source parameters based on physical constraints and array geometry adaptation, the method comprising the following steps:

[0082] Step 1: First, determine the continuous waveform data and location information corresponding to multiple stations in the dataset as input, then determine the predicted source parameters as output, and then set an initial prediction model;

[0083] The second step is to extract features from the continuous waveform data and location information to obtain the predicted source parameters.

[0084] Step 3: First determine whether the dataset includes earthquake phase labels and the true value of the earthquake source. If the dataset includes earthquake phase labels and the true value of the earthquake source, proceed to step 4. If the dataset does not include earthquake phase labels and the true value of the earthquake source, proceed to step 5.

[0085] Step 4: First, compare the predicted source parameters with the actual source values ​​to obtain the regression loss. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, weight the regression loss and phase loss to obtain the total loss. Then, proceed to step 6.

[0086] Step 5: First, obtain the phase labels based on the dataset using the phase picking method. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, use the phase loss as the total loss and then proceed to step 6.

[0087] Step 6: First, train the initial prediction model, then select the initial prediction model with the minimum total loss as the optimal model.

[0088] Step 7: First, input the continuous waveform data and location information of multiple stations into the optimal model, and then obtain the predicted source parameters from the optimal model.

[0089] In the fourth and fifth steps, the actual arrival data refers to the actual P-wave arrival data and the actual S-wave arrival data, and the predicted arrival data refers to the predicted P-wave arrival data and the predicted S-wave arrival data.

[0090] The calculation steps for the predicted P-wave arrival time data and the predicted S-wave arrival time data are as follows:

[0091] First, calculate the epicentral distance:

[0092] ;

[0093] In the formula: The epicentral distance is... The average radius of the Earth This is the angular distance between the seismic source and the station;

[0094] Then calculate the distance from the epicenter:

[0095] ;

[0096] In the formula: The distance from the epicenter, The focal depth;

[0097] Then calculate the theoretical time-to-run time:

[0098] ; ;

[0099] In the formula: To predict the arrival time of the P-wave, To predict the arrival time of the S-wave, The distance from the epicenter, For P-wave velocity, This represents the S-wave velocity.

[0100] In steps four and five, obtaining the phase loss based on the predicted arrival time data and the actual arrival time data refers to:

[0101] First, calculate the phase loss of the predicted P-wave arrival data:

[0102] ;

[0103] In the formula: The phase loss of the P-wave arrival data. For all stations, To predict the arrival time of the P-wave, This is the actual arrival time data for the P-wave. The confidence level for picking up seismic phases;

[0104] Next, calculate the phase loss of the data when the predicted S-wave arrives:

[0105] ;

[0106] In the formula: The phase loss of the data when the S-wave arrives. For all stations, To predict the arrival time of the S-wave, This is the actual arrival time data for the S-wave. The confidence level for picking up seismic phases;

[0107] The phase loss is then obtained by averaging the phase loss of the predicted P-wave arrival data and the phase loss of the predicted S-wave arrival data.

[0108] ;

[0109] In the formula, This is for phase loss.

[0110] In the fourth step, the regression loss obtained by comparing the predicted source parameters with the actual source values ​​refers to:

[0111] ;

[0112] In the formula, To regress the loss, To predict earthquake source parameters, This represents the true value of the earthquake source.

[0113] In the fourth step, the weighted average of regression loss and phase loss to obtain the total loss refers to:

[0114] ;

[0115] ; ;

[0116] In the formula, For the total loss, The weights for phase loss, The weights for the regression loss;

[0117] In the fifth step, referring to the phase loss as the total loss means:

[0118] ;

[0119] In the formula, This represents the total loss.

[0120] In the second step, feature extraction from continuous waveform data and location information to obtain predicted source parameters refers to: first, performing feature extraction for a single station. The feature extraction process for a single station involves: obtaining the local time-frequency features of the station based on its continuous waveform data; then obtaining the geometric distribution features of the station based on its location information; and finally fusing the geometric distribution features with the local time-frequency features to obtain the station's features. This process is repeated for each single station until the features of all stations are obtained.

[0121] Then, the characteristics of all stations are averaged to obtain the global event characteristics, and the predicted source parameters are obtained based on the global event characteristics.

[0122] In the second step, obtaining the local time-frequency features of the station based on the continuous waveform data of the station means: using a CNN convolutional neural network to obtain the local time-frequency features based on the continuous waveform data of the station;

[0123] In the second step, obtaining the geometric distribution features of a station based on its location information means: mapping the latitude and longitude coordinates of all stations to a multidimensional space through a linear layer and activation function to obtain the initial location features of all stations; calculating the distance information of all stations based on the initial location features of all stations; calculating a distance matrix based on all the distance information; defining the edge weights of the distance matrix using a Gaussian kernel; generating graph structure data based on the edge weights; and extracting the geometric distribution features of all stations from the graph structure data using a GNN graph neural network.

[0124] In the second step, calculating the distance information of all stations based on the initial location features of all stations means:

[0125] ;

[0126] ;

[0127] ; ;

[0128] In the formula, For the first The station and the first The surface distance between each station The average radius of the Earth No. The station and the first angular distance between stations It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station. It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station;

[0129] In the second step, defining the edge weights of the distance matrix using a Gaussian kernel and then generating graph structure data based on those edge weights refers to:

[0130] ;

[0131] In the formula, For graph-structured data, Distance threshold For width parameters;

[0132] In the second step, fusing the geometric distribution characteristics of the station with local time-frequency characteristics to obtain the station characteristics refers to:

[0133] ;

[0134] ;

[0135] ;

[0136] In the formula, Features of the station It is a local time-frequency feature. For parameter items, It exhibits geometric distribution characteristics.

[0137] In the third step, obtaining the predicted source parameters based on global event features means: using a multilayer perceptron (MLP) to perform nonlinear transformation and mapping on the global event features to obtain high-dimensional information, and then outputting the high-dimensional information as the predicted source parameters.

[0138] In the sixth step, the phase picking method includes: PhaseNet method, STA / LTA method, or EQtransform method.

[0139] The following are supplementary descriptions of the present invention:

[0140] The continuous waveform data described in this invention refers to seismic waveform data containing three components: vertical, north-south, and east-west.

[0141] The location information of the station mentioned in this invention refers to the station's latitude, longitude, and elevation information.

[0142] Example 1:

[0143] Please see Figure 1 — Figure 3 A method for predicting seismic source parameters based on physical constraints and array geometry adaptation, the method comprising the following steps:

[0144] Step 1: First, determine the continuous waveform data and location information corresponding to multiple stations in the dataset as input, then determine the predicted source parameters as output, and then set an initial prediction model;

[0145] The second step is to extract features from the continuous waveform data and location information to obtain the predicted source parameters.

[0146] Step 3: First determine whether the dataset includes earthquake phase labels and the true value of the earthquake source. If the dataset includes earthquake phase labels and the true value of the earthquake source, proceed to step 4. If the dataset does not include earthquake phase labels and the true value of the earthquake source, proceed to step 5.

[0147] Step 4: First, compare the predicted source parameters with the actual source values ​​to obtain the regression loss. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, weight the regression loss and phase loss to obtain the total loss. Then, proceed to step 6.

[0148] Step 5: First, obtain the phase labels based on the dataset using the phase picking method. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, use the phase loss as the total loss and then proceed to step 6.

[0149] Step 6: First, train the initial prediction model, then select the initial prediction model with the minimum total loss as the optimal model.

[0150] Step 7: First, input the continuous waveform data and location information of multiple stations into the optimal model, and then obtain the predicted source parameters from the optimal model.

[0151] Preferably, in the sixth step, the training process is as follows: input a set of untrained datasets into the initial prediction model, then execute the second step until the total loss corresponding to the set of datasets is obtained, and at the same time, backpropagate to update the parameters of the initial prediction model.

[0152] Example 2:

[0153] The basic content is the same as in Example 1, except that:

[0154] Please see Figure 1 — Figure 2 In the fourth and fifth steps, the actual arrival data refers to the actual P-wave arrival data and the actual S-wave arrival data, and the predicted arrival data refers to the predicted P-wave arrival data and the predicted S-wave arrival data.

[0155] The calculation steps for the predicted P-wave arrival time data and the predicted S-wave arrival time data are as follows:

[0156] First, calculate the epicentral distance:

[0157] ;

[0158] In the formula: The epicentral distance is... The average radius of the Earth This is the angular distance between the seismic source and the station;

[0159] Then calculate the distance from the epicenter:

[0160] ;

[0161] In the formula: The distance from the epicenter, The focal depth;

[0162] Then calculate the theoretical time-to-run time:

[0163] ; ;

[0164] In the formula: To predict the arrival time of the P-wave, To predict the arrival time of the S-wave, The distance from the epicenter, For P-wave velocity, This represents the S-wave velocity.

[0165] In steps four and five, obtaining the phase loss based on the predicted arrival time data and the actual arrival time data refers to:

[0166] First, calculate the phase loss of the predicted P-wave arrival data:

[0167] ;

[0168] In the formula: The phase loss of the P-wave arrival data. For all stations, To predict the arrival time of the P-wave, This is the actual arrival time data for the P-wave. The confidence level for picking up seismic phases;

[0169] Next, calculate the phase loss of the data when the predicted S-wave arrives:

[0170] ;

[0171] In the formula: The phase loss of the data when the S-wave arrives. For all stations, To predict the arrival time of the S-wave, This is the actual arrival time data for the S-wave. The confidence level for picking up seismic phases;

[0172] The phase loss is then obtained by averaging the phase loss of the predicted P-wave arrival data and the phase loss of the predicted S-wave arrival data.

[0173] ;

[0174] In the formula, This is for phase loss.

[0175] In the fourth step, the regression loss obtained by comparing the predicted source parameters with the actual source values ​​refers to:

[0176] ;

[0177] In the formula, To regress the loss, To predict earthquake source parameters, This represents the true value of the earthquake source.

[0178] In the fourth step, the weighted average of regression loss and phase loss to obtain the total loss refers to:

[0179] ;

[0180] ; ;

[0181] In the formula, For the total loss, The weights for phase loss, The weights for the regression loss;

[0182] In the fifth step, referring to the phase loss as the total loss means:

[0183] ;

[0184] In the formula, This represents the total loss.

[0185] Preferably, the two-dimensional velocity model used to calculate the theoretical travel time can be replaced with a three-dimensional velocity model, which can improve the accuracy of the theoretical travel time calculation.

[0186] In application, this method embeds the physical model of seismic wave propagation (specifically, the spherical travel time calculation formula) as a differentiable module directly into the architecture of a deep neural network. First, it calculates the epicentral distance, i.e., the surface arc distance from the epicenter to the station; then, it calculates the focal distance, i.e., the straight-line spatial distance from the source to the station. Based on the straight-line spatial distance, P-wave velocity, S-wave velocity, and velocity model, it obtains predicted P-wave arrival time data and predicted S-wave arrival time data. The phase loss of the predicted P-wave arrival time data and the phase loss of the predicted S-wave arrival time data are then used to measure the difference between the predicted travel time and the actual arrival time, directly reflecting the physical constraints. Confidence weighting (i.e., only reliable phases with a confidence level higher than 0.5 are used in the loss calculation) automatically filters out unreliable phase picking results, making the training process tolerant of noise from the automatic picker. Finally, it compares the predicted source parameters with the actual source values ​​to obtain the regression loss, which measures the difference between the predicted and actual travel times. The method predicts the difference between the source parameters and the true labels. Then, depending on whether the training data includes phase labels, different loss integration methods are selected. If the training data includes phase labels, a weighted summation method is used, and the total loss is obtained by utilizing both data labels and physical laws. If the training data does not include phase labels, only the phase loss is used as the total loss, realizing weakly supervised learning and improving data utilization. This method can solve the pure data-driven defect of existing end-to-end methods, that is, existing models learn the statistical correlation between waveforms and parameters, and their prediction results may be physically inconsistent (for example, the predicted source location contradicts the arrival time information in the waveform). In contrast, during inference, the source parameters output by the model and the calculated theoretical travel time must meet the preset velocity model, avoiding the physical contradiction of "the predicted location is very far but the travel time is very short", that is, realizing physical constraints, which greatly improves the interpretability and credibility of the model.

[0187] Example 3:

[0188] The basic content is the same as in Example 1, except that:

[0189] Please see Figure 1 — Figure 3 In the second step, feature extraction from continuous waveform data and location information to obtain predicted source parameters refers to: first, performing feature extraction for a single station, which involves: obtaining the local time-frequency features of the station based on its continuous waveform data; then obtaining the geometric distribution features of the station based on its location information; then fusing the geometric distribution features with the local time-frequency features to obtain the station features; repeating the feature extraction process for a single station until the station features of all stations are obtained; then averaging all the station features to obtain global event features; and finally obtaining the predicted source parameters based on the global event features. In the second step, obtaining the local time-frequency features of the station based on the continuous waveform data of the station means: using a CNN convolutional neural network to obtain the local time-frequency features based on the continuous waveform data of the station; in the second step, obtaining the geometric distribution features of the station based on the location information of the station means: mapping the latitude and longitude coordinates of all stations to a multi-dimensional space through a linear layer and activation function to obtain the initial location features of all stations, then calculating the distance information of all stations based on the initial location features of all stations, then calculating a distance matrix based on all the distance information, then using a Gaussian kernel to define the edge weights of the distance matrix, then generating graph structure data based on the edge weights, and then extracting the geometric distribution features of all stations from the graph structure data through a GNN graph neural network.

[0190] In the second step, calculating the distance information of all stations based on the initial location features of all stations means:

[0191] ;

[0192] ;

[0193] ; ;

[0194] In the formula, For the first The station and the first The surface distance between each station The average radius of the Earth No. The station and the first angular distance between stations It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station. It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station;

[0195] In the second step, defining the edge weights of the distance matrix using a Gaussian kernel and then generating graph structure data based on those edge weights refers to:

[0196] ;

[0197] In the formula, For graph-structured data, Distance threshold For width parameters;

[0198] Preferably, the step of calculating a distance matrix based on all distance information, defining the edge weights of the distance matrix using a Gaussian kernel, and then generating graph structure data based on the edge weights refers to: constructing the edges (connections) of the graph based on distance information and distance thresholds, assigning a weight reflecting spatial proximity to each edge using a Gaussian kernel function, and then using a multi-head graph attention network (GAT) to transmit information on the graph, allowing the features of each station to adaptively aggregate the information of its neighboring stations according to the attention weights; the entire graph structure (including node coordinates and edge information) is registered as a buffer for the model and deployed along with the model.

[0199] In application, CNN automatically and layer-by-layer learns and extracts useful local time-frequency features related to seismic events from continuous waveform data of the station, including the start position (arrival time) of P-wave and S-wave phases, the amplitude, period (or frequency components), duration, and waveform morphology of the phases. Traditional methods require manual algorithm design to pick up phases or calculate features, while this method achieves end-to-end automatic feature learning through CNN, significantly improving the picking efficiency. When mapping the coordinates of the station, the original latitude and longitude (in radians) are first mapped to a high-dimensional space through a linear layer to obtain the initial position features, and then the spherical trigonometric formula is used to accurately calculate any position. The surface arc distance between two stations was calculated by comprehensively considering the difference in latitude and longitude and the influence of latitude on the distance, ensuring the accuracy of distance calculation at any latitude (including near the equator and poles), which is the foundation for constructing reliable spatial relationships. Then, the discrete station network was abstracted into a graph, where stations are nodes and spatial proximity relationships between stations are edges. Gaussian kernels were used to assign weights to the edges, with weights ranging from 0 to 1. The closer the distance, the closer the weight is to 1, indicating a stronger connection and potentially more similar waveforms. Furthermore, a distance threshold was used to control the sparsity of the graph, ensuring that each station interacts only with a limited number of geographically nearby "neighbor" stations, consistent with seismic wavefields. Spatial-related physical intuition also improves computational efficiency; the width parameter and threshold can be adjusted according to the average station spacing of the network and the regional crustal characteristics, enabling the method to flexibly adapt to monitoring tasks at different scales (such as regional networks and dense arrays); this technique can compensate for the shortcomings of traditional CNNs or fully connected networks in handling variable and unevenly distributed station numbers, especially under geometric structures such as linear arrays or sparse networks where the performance of traditional CNNs or fully connected networks degrades sharply; this technique naturally adapts to any station layout through graph structure; and then uses a multi-head graph attention network (GAT) to dynamically calculate the attention coefficient between each station and all its "neighbor" stations. (Instead of simply treating all neighbors equally), this coefficient determines the contribution of each "neighbor" station's features when aggregating information, i.e., it can distinguish which "neighbor" stations' information is more critical, thus mathematically alleviating the ill-posed problem and significantly improving the prediction robustness in complex array environments; the final output geometric distribution features encode rich spatial structure information, registering the complete graph structure containing all node coordinates and edge information as a buffer for the model, saving and deploying it together with the model parameters. For fixed network stations, this means that there is no need to repeatedly calculate the graph structure during the application phase; it can be directly loaded for efficient inference, which is very suitable for real-time seismic monitoring systems.

[0200] Example 4:

[0201] The basic content is the same as in Example 1, except that:

[0202] Please see Figure 1 — Figure 3In the second step, fusing the geometric distribution characteristics of the station with the local time-frequency characteristics to obtain the station characteristics refers to:

[0203] ;

[0204] ;

[0205] ;

[0206] In the formula, Features of the station It is a local time-frequency feature. For parameter items, It exhibits geometric distribution characteristics.

[0207] When applying, It is a local time-frequency feature extracted from the continuous waveform data of the station by a CNN convolutional neural network. It encodes the specific manifestation of the seismic signal in the waveform of the station, such as the phase characteristics, amplitude, and frequency components of P-wave and S-wave, and is the direct basis for judging "what was heard". It starts with the location information (latitude and longitude coordinates) of the station and is ultimately processed by a GNN (Graph Neural Network) to obtain geometric distribution features. It encodes the station's spatial location within the entire monitoring network and its relative relationship with other stations, serving as the direct basis for determining "where it was heard". It is an important adjustable parameter, which can be 0.1. It is used to control the proportion of geometric distribution features in the final fused features, to avoid the geometric distribution features drowning out the weak seismic phase features in the local time-frequency features, and at the same time to ensure that the geometric distribution features always have a gradient path and are stably optimized during training. It can be flexibly adapted to different situations. When the array layout is ideal (such as dense and uniform distribution), the waveform characteristics of each station can provide sufficient constraints, and smaller settings can be used. When the array has geometric defects (such as sparse or linear distribution), the model can be improved by increasing the [scale value]. Values ​​(e.g., 0.3 or higher) are used to increase the proportion of geometric distribution features, thereby utilizing limited spatial relationships to compensate for insufficient waveform orientation information and effectively improve positioning robustness under non-ideal network conditions; fused features It is a unified, high-dimensional representation that simultaneously includes the "physical signal" characteristics of the event and the "spatial" characteristics of the observation, providing a complete input for subsequent global information aggregation (average pooling) and final source parameter regression (MLP).

[0208] Example 5:

[0209] The basic content is the same as in Example 1, except that:

[0210] Please see Figure 1 — Figure 3 In the third step, obtaining the predicted source parameters based on global event features means: using a multilayer perceptron (MLP) to perform nonlinear transformation and mapping on the global event features to obtain high-dimensional information, and then outputting the high-dimensional information as the predicted source parameters.

[0211] In application, the input to the MLP is the global event feature, which is formed by aggregating (average pooling) the features from all stations in the second step. This feature contains comprehensive information about the earthquake event learned from the waveforms and spatial geometry of the entire station array. The MLP gradually transforms and maps the global event feature through multiple layers of nonlinear transformations, converting it into a higher-dimensional, more abstract representation. Ultimately, it understands the complex relationship between these features and the target source parameters, and the final output is the predicted source parameters, which typically include: source latitude and longitude, source depth, seismogenic event, and magnitude. For example... The CNN extracts features with a dimension of 256, and the GNN extracts location information features with a dimension of 256. These two features with the same dimension are concatenated to obtain a single-station concatenated feature (dimension 256). After feature fusion from all stations, the global features (256 dimensions) are passed through a 3-layer MLP (Multilayer Perceptron) (256→128→64→4) to output four-dimensional predicted source parameters. MLP can learn a highly complex and nonlinear mapping relationship between global features and source parameters, which is impossible for simple linear regression, thus ensuring high-precision prediction.

[0212] Example 6:

[0213] The basic content is the same as in Example 1, except that:

[0214] Please see Figure 1 — Figure 3 In the sixth step, the phase picking method includes: PhaseNet method, STA / LTA method or EQtransform method.

[0215] In application, the PhaseNet method is an advanced phase picking model based on deep learning, capable of automatically detecting the arrival times of P-waves and S-waves with high accuracy and recall from continuous waveforms; the STA / LTA method, or "short-time average / long-time average" ratio method, is a classic and reliable triggering and phase picking algorithm in seismology, identifying sudden signals by calculating the short-term to long-term energy ratio of waveforms; the EQtransform method is a transducer-based seismic signal processing method that uses an attention mechanism to process waveform sequences and achieve tasks such as phase picking; different picking methods have their own advantages in different scenarios (e.g., PhaseNet is sensitive to complex waveforms, while STA / LTA is computationally simple and robust), providing multiple options means that if one method fails or is not accurate enough, the results of other methods can be switched or fused, thereby ensuring the reliability of the phase labels on which the "phase loss" calculation in step six depends, and guaranteeing the effectiveness of physical constraint training; moreover, this design reserves an interface for future integration of more advanced phase picking algorithms, so that the performance of the entire method can be improved with the advancement of underlying picking technology.

[0216] Example 7:

[0217] The basic content is the same as in Example 1, except that:

[0218] Please see Figure 1 — Figure 3 A device and system for predicting source parameters using a method based on physical constraints and array geometry adaptation, comprising a feature extraction module, a global pooling and source parameter regression module, and a physical constraint module. The feature extraction module is used to perform feature extraction steps. The global pooling and source parameter regression module is used to average all station features to obtain global event features, and then obtain predicted source parameters based on the global event features. The physical constraint module is used to obtain regression loss, phase loss, and total loss based on regression loss and phase loss.

[0219] The above description is only a preferred embodiment of the present invention. The scope of protection of the present invention is not limited to the above embodiments. Any equivalent modifications or changes made by those skilled in the art based on the content disclosed in the present invention should be included within the scope of protection set forth in the claims.

Claims

1. A method for predicting source parameters based on physical constraints and array geometry adaptation, characterized in that: The method includes the following steps: Step 1: First, determine the continuous waveform data and location information corresponding to multiple stations in the dataset as input, then determine the predicted source parameters as output, and then set an initial prediction model; The second step is to extract features from the continuous waveform data and location information to obtain the predicted source parameters. Step 3: First determine whether the dataset includes earthquake phase labels and the true value of the earthquake source. If the dataset includes earthquake phase labels and the true value of the earthquake source, proceed to step 4. If the dataset does not include earthquake phase labels and the true value of the earthquake source, proceed to step 5. Step 4: First, compare the predicted source parameters with the actual source values ​​to obtain the regression loss. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, weight the regression loss and phase loss to obtain the total loss. Then, proceed to step 6. Step 5: First, obtain the phase labels based on the dataset using the phase picking method. Then, obtain the actual arrival time data based on the phase labels. Next, obtain the predicted arrival time data based on the predicted source parameters. Then, obtain the phase loss based on the predicted arrival time data and the actual arrival time data. Finally, use the phase loss as the total loss and then proceed to step 6. Step 6: First, train the initial prediction model, then select the initial prediction model with the minimum total loss as the optimal model. Step 7: First, input the continuous waveform data and location information of multiple stations into the optimal model, and then obtain the predicted source parameters from the optimal model.

2. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 1, characterized in that: In the fourth and fifth steps, the actual arrival data refers to the actual P-wave arrival data and the actual S-wave arrival data, and the predicted arrival data refers to the predicted P-wave arrival data and the predicted S-wave arrival data. The calculation steps for the predicted P-wave arrival time data and the predicted S-wave arrival time data are as follows: First, calculate the epicentral distance: ； In the formula: The epicentral distance is... The average radius of the Earth This is the angular distance between the seismic source and the station; Then calculate the distance from the epicenter: ； In the formula: The distance from the epicenter, The focal depth; Then calculate the theoretical time-to-run time: ； ； In the formula: To predict the arrival time of the P-wave, To predict the arrival time of the S-wave, The distance from the epicenter, For P-wave velocity, This represents the S-wave velocity.

3. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 2, characterized in that: In steps four and five, obtaining the phase loss based on the predicted arrival time data and the actual arrival time data refers to: First, calculate the phase loss of the predicted P-wave arrival data: ； In the formula: The phase loss of the P-wave arrival data. For all stations, To predict the arrival time of the P-wave, This is the actual arrival time data for the P-wave. The confidence level for picking up seismic phases; Next, calculate the phase loss of the data when the predicted S-wave arrives: ； In the formula: The phase loss of the data when the S-wave arrives. For all stations, To predict the arrival time of the S-wave, This is the actual arrival time data for the S-wave. The confidence level for picking up seismic phases; The phase loss is then obtained by averaging the phase loss of the predicted P-wave arrival data and the phase loss of the predicted S-wave arrival data. ； In the formula, This is for phase loss.

4. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 3, characterized in that: In the fourth step, the regression loss obtained by comparing the predicted source parameters with the actual source values ​​refers to: ； In the formula, To regress the loss, To predict earthquake source parameters, This represents the true value of the earthquake source.

5. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 4, characterized in that: In the fourth step, the weighted average of regression loss and phase loss to obtain the total loss refers to: ； ； ； In the formula, For the total loss, The weights for phase loss, The weights for the regression loss; In the fifth step, referring to the phase loss as the total loss means: ； In the formula, This represents the total loss.

6. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 1, characterized in that: In the second step, feature extraction of continuous waveform data and location information to obtain predicted source parameters refers to: first, performing feature extraction for a single station. The feature extraction process for a single station is as follows: obtaining the local time-frequency features of the station based on the continuous waveform data of the station, obtaining the geometric distribution features of the station based on the location information of the station, and then fusing the geometric distribution features of the station with the local time-frequency features to obtain the station features. Repeat the feature extraction process for each individual station until the station features for all stations are obtained. Then, the characteristics of all stations are averaged to obtain the global event characteristics, and the predicted source parameters are obtained based on the global event characteristics.

7. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 6, characterized in that: In the second step, obtaining the local time-frequency features of the station based on the continuous waveform data of the station means: using a CNN convolutional neural network to obtain the local time-frequency features based on the continuous waveform data of the station; In the second step, obtaining the geometric distribution features of a station based on its location information means: mapping the latitude and longitude coordinates of all stations to a multidimensional space through a linear layer and activation function to obtain the initial location features of all stations; calculating the distance information of all stations based on the initial location features of all stations; calculating a distance matrix based on all the distance information; defining the edge weights of the distance matrix using a Gaussian kernel; generating graph structure data based on the edge weights; and extracting the geometric distribution features of all stations from the graph structure data using a GNN graph neural network.

8. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 7, characterized in that: In the second step, calculating the distance information of all stations based on the initial location features of all stations means: ； ； ； ； In the formula, For the first The station and the first The surface distance between each station The average radius of the Earth No. The station and the first angular distance between stations It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station. It is the first The latitude (in radians) of each station. It is the first Longitude (in radians) of each station; In the second step, defining the edge weights of the distance matrix using a Gaussian kernel and then generating graph structure data based on those edge weights refers to: ； In the formula, For graph-structured data, Distance threshold For width parameters; In the second step, fusing the geometric distribution characteristics of the station with local time-frequency characteristics to obtain the station characteristics refers to: ； ； ； In the formula, Features of the station It is a local time-frequency feature. For parameter items, It exhibits geometric distribution characteristics.

9. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 6, characterized in that: In the third step, obtaining the predicted source parameters based on global event features means: using a multilayer perceptron (MLP) to perform nonlinear transformation and mapping on the global event features to obtain high-dimensional information, and then outputting the high-dimensional information as the predicted source parameters.

10. The method for predicting source parameters based on physical constraints and array geometry adaptation according to claim 1, characterized in that: In the sixth step, the phase picking method includes: PhaseNet method, STA / LTA method, or EQtransform method.

Images