An automated measurement method for local earthquake S-wave splitting based on multiple methods

Abstract

The present invention discloses a multi-method based near-seismic S-wave splitting automated measurement method, which relates to the technical field of seismic wave polarization analysis. The present invention combines artificial visualization measurement, cross-correlation function and minimum tangential energy algorithm to further eliminate abnormal points in search results, improve the accuracy and stability of the near-seismic S-wave splitting automated measurement method; select multiple filter ranges and analysis window lengths, eliminate the instability of search results caused by a single filter range and analysis window; limit the moving range of the analysis window by the length of the main period of the S-wave initial motion signal, and effectively avoid the interference of surface wave signals with larger energy on the measurement results.

Description

Technical Field

[0001] The present invention relates to the technical field of seismic wave polarization analysis, and in particular to a multi-method based automated measurement method for near-seismic S-wave splitting. Background Art

[0002] When propagating in an isotropic medium, the S wave will not split. Only radial and vertical S wave signals can be observed at the receiving point, and the energy of the S wave signal in the tangential direction is zero. On the contrary, when the S wave propagates in an anisotropic medium, the S wave will split and produce a pair of fast and slow waves with approximately orthogonal polarization directions. The polarization direction of the fast wave and the delay time of the slow wave are extremely sensitive to the direction and intensity of the deformation of the underground medium, and the fast and slow waves will show some regular changes in different propagation directions. For example, the energy of the tangential signal can be observed at the surface receiving point, and the waveforms of the fast and slow wave signals are very similar. Many scholars have used the above characteristics to develop a series of algorithms such as polarization analysis, cross-correlation function and minimum tangential energy to obtain the S wave splitting parameters. The measurement process of the artificial visualization measurement method is relatively simple, but its measurement results are easily affected by subjective consciousness and the calculation efficiency is low. The measurement results of the cross-correlation algorithm and the minimum tangential energy algorithm are easily affected by the length of the S wave analysis window, resulting in the surface wave signal after the S wave signal participating in the calculation process of the S wave splitting parameters, resulting in greater instability in the results. The MFAST Manual V2.2 program automated measurement and SAM system analysis method developed based on the cross-correlation function and minimum tangential energy algorithm does not simultaneously consider the impact of multiple filter ranges and S-wave analysis windows on the measurement results, thereby reducing the accuracy of the measurement results. Summary of the invention

[0003] In view of the above-mentioned deficiencies in the prior art, the present invention provides an automated measurement method for near-seismic S-wave splitting based on multiple methods, which solves the problem that a variety of filtering ranges and S-wave analysis windows are not taken into consideration, and the surface wave signal after the S-wave signal causes unstable measurement results and the existence of abnormal points.

[0004] In order to achieve the above-mentioned object of the invention, the technical solution adopted by the present invention is:

[0005] A multi-method based automated measurement method for near-earthquake S-wave splitting is provided, which comprises the following steps:

[0006] S1, obtaining the original waveform data of the earthquake event and setting the noise window;

[0007] S2, presetting a band-pass filter and screening the band-pass filter;

[0008] S3, divide the S wave analysis window;

[0009] S4, using the screened bandpass filters and S-wave analysis windows to form multiple search networks in pairs, and using the cross-correlation function and the minimum tangential energy algorithm to search for the initial S-wave splitting parameters in different search networks;

[0010] S5. Determine the initial S-wave splitting parameters obtained by each search network, obtain the final S-wave splitting parameters, and complete the automated measurement of the near-earthquake S-wave splitting.

[0011] Further, the original waveform data in step S1 is a three-component seismic record containing P-wave and S-wave initial motion signals of a near-seismic event, and its file format is SAC; wherein the three components are two horizontal components and a vertical component, namely, an N-direction component, an E-direction component, and a Z-direction component;

[0012] Step S1 includes the following steps:

[0013] S1-1, respectively picking up the arrival time of the initial motion signal corresponding to the S wave signal of the N-direction component and the P wave signal of the Z-direction component, and obtaining the arrival time of the corresponding S wave initial motion signal and the arrival time of the P wave initial motion signal;

[0014] S1-2, based on the arrival time of the P wave initial motion signal, selecting the time p of the P wave initial motion signal;

[0015] S1-3. The arrival time of the S wave initial motion signal is taken as time t1 and defined as the end position of the noise window. The time t0 before time t1 is selected as the starting position of the noise window to complete the setting of the noise window.

[0016] Further, step S2 includes the following steps:

[0017] S2-1, respectively analyzing the frequency spectra of the two horizontal components of the S wave signal to obtain the corresponding analysis results of the two horizontal components, that is, obtaining the main period of the S wave signal of the corresponding horizontal component; taking the average of the analysis results of the two horizontal components as the analysis result of the frequency spectrum of the S wave signal;

[0018] S2-2, based on the analysis results of the S wave signal spectrum, set M bandpass filters; where M is 10;

[0019] S2-3, using any bandpass filter to filter the S wave signals of the two horizontal components to obtain a filtered S wave signal;

[0020] S2-4, calculate the average energy of the filtered S-wave signal within a complete cycle;

[0021] S2-5, calculating the average energy of the filtered noise signal in the noise window;

[0022] S2-6, taking the ratio of the average energy of step S2-4 to the average energy of step S2-5 as the signal-to-noise ratio;

[0023] S2-7, determine whether the signal-to-noise ratio is not less than 3; if so, retain the bandpass filter and increase the value of m by 1; otherwise, discard the bandpass filter; proceed to step S2-8; wherein the initial value of m is 0;

[0024] S2-8, repeating steps S2-3 to S2-7 until all the band-pass filters are judged;

[0025] S2-9, determine whether the m value is 0; if so, determine that the original waveform data cannot meet the requirements of near-earthquake S-wave splitting and end the measurement; otherwise, enter step S3.

[0026] Furthermore, the specific process of step S3 is as follows:

[0027] The starting position of the noise window is used as the starting position of the S wave analysis window, and the 0.5 cycle of the S wave initial motion signal is used as the end position of the first S wave analysis window; the end position of the S wave analysis window is moved backward with a step size of 0.001s as the end position of the next S wave analysis window, until it ends at the 1.2 cycle of the S wave initial motion signal, and Q S wave analysis windows of different lengths are obtained.

[0028] Further, step S4 includes the following steps:

[0029] S4-1, combining the bandpass filter and the S-wave analysis window to obtain different search grids; wherein the number of search grids is the product of the number of the bandpass filter and the number of the S-wave analysis window;

[0030] S4-2, select any search grid, and rotate the S wave signals of the two horizontal components in the range of 0 to 180 degrees according to the step size of 1 degree to obtain fast wave signals and slow wave signals at different angles;

[0031] S4-3, in the interval of 0 to 0.1 s, time-correcting the fast wave signal and the slow wave signal obtained in step S4-2 according to a step length of 0.001 s, to obtain fast wave signals and slow wave signals after different time corrections;

[0032] S4-4, using the cross-correlation function to process the time-corrected fast wave signal and the slow wave signal to obtain the cross-correlation coefficients of the fast wave signal and the slow wave signal after different angles and time corrections; all the cross-correlation coefficients constitute a cross-correlation coefficient matrix;

[0033] S4-5, based on the cross-correlation coefficient matrix obtained in step S4-4, searching for the angles and correction times corresponding to the fast wave signal and the slow wave signal that make the waveforms most similar, and obtaining the fast wave polarization direction and the slow wave delay time, that is, the result of the cross-correlation function;

[0034] S4-6, obtaining a header file of a three-component seismic record and determining whether the header file contains a variable BAZ; if so, using the variable BAZ as the reverse azimuth of the seismic event; otherwise, searching for the reverse azimuth of the seismic event based on the waveform of the three-component seismic record near the time p of the P wave first motion signal;

[0035] S4-7, based on the back azimuth of the seismic event obtained in step S4-6, the time-corrected fast wave signal and slow wave signal are calculated using the minimum tangential energy algorithm to obtain the corresponding radial and tangential signals, and the average energy of the tangential signal in the S-wave analysis window is calculated; all tangential signal energies constitute a tangential signal energy matrix;

[0036] S4-8, searching according to the tangential signal energy matrix obtained in step S4-7, determining the fast wave polarization direction and slow wave delay time with the minimum average energy of the tangential signal, that is, the search result of the minimum tangential energy algorithm;

[0037] S4-9, repeating steps S4-2 to S4-8 until the search results of all search grids are obtained, i.e., the initial S-wave splitting parameters;

[0038] The initial S-wave splitting parameters include the search results of the cross-correlation function and the minimum tangential energy algorithm. Further, the calculation formula of step S4-3 is:

[0039] S F =S(t)·cosφ

[0040] S S =S(t-δt)·sinφ

[0041] Among them, S F , S S They represent the fast wave signal and slow wave signal after time correction, S(t) represents the composite signal of the S wave signals in the two horizontal components, cos and sin represent the cosine function and the sine function, t represents time, φ represents the angle between the fast wave polarization direction and the radial direction, δt represents the slow wave delay time, and S(t-δt) represents the S wave signal after a delay of δt.

[0042] Furthermore, the calculation formula of step S4-4 is:

[0043]

[0044] Here, COV represents the mutual correlation coefficient, and |·| represents the absolute value.

[0045] Furthermore, the calculation formula of step S4-7 is:

[0046] S R =S F ·cosφ+S S ·sinφ

[0047] S T =S S ·cosφ-S F ·sinφ

[0048] Among them, S R represents the radial signal, S T Indicates tangential signal.

[0049] Furthermore, the specific process of step S5 is:

[0050] S5-1, respectively normalizing the mutual correlation coefficient matrix and the tangential signal energy matrix obtained by any search network to obtain a normalized mutual correlation coefficient matrix and a normalized tangential signal energy matrix;

[0051] S5-2, searching the normalized mutual correlation coefficient matrix to obtain the fast wave polarization direction and slow wave delay time with a correlation coefficient greater than 95%, and calculating the corresponding standard deviations respectively;

[0052] S5-3, calculating the sum of the two standard deviations in step S5-2;

[0053] S5-4, searching the normalized tangential signal energy matrix to obtain the fast wave polarization direction and slow wave delay time with energy less than 5%, and calculating the corresponding standard deviations respectively;

[0054] S5-5, calculating the sum of the two standard deviations in step S5-4;

[0055] S5-6, summing the standard deviations of step S5-3 and step S5-5 to obtain the standard deviation of all search networks;

[0056] S5-7, compare the standard deviations of different search networks, and select the initial S-wave splitting parameters of the search network with the smallest standard deviation;

[0057] S5-8. Compare the initial S-wave splitting parameters obtained in step S5-7; if the search results of the cross-correlation function and the minimum tangential energy algorithm are consistent, the search result of the minimum tangential energy algorithm is selected as the final S-wave splitting parameter; if the results of the two algorithms are inconsistent, the two search results are inspected by a method based on manual visual inspection and visualization, and selected as the final S-wave splitting parameter.

[0058] The beneficial effects of the present invention are as follows: the method can further eliminate abnormal points in the search results by combining the cross-correlation function and the minimum tangential energy algorithm, thereby improving the accuracy and stability of the program; multiple filter ranges and analysis window lengths are selected to eliminate the instability of the search results caused by a single filter range and analysis window; the moving range of the analysis window is limited by the length of the main period of the S-wave initial motion signal, thereby effectively avoiding the interference of the surface wave signal with relatively large energy on the measurement results; and the measurement results of the cross-correlation function and the minimum tangential energy algorithm are further screened by an artificial visualization measurement method, thereby improving the accuracy of the final measurement results and having strong adaptability. BRIEF DESCRIPTION OF THE DRAWINGS

[0059] Figure 1 is a method flow chart;

[0060] Figure 2 It is a schematic diagram of three-component seismic recording data and noise window;

[0061] Figure 3 This is a schematic diagram of S wave splitting;

[0062] Figure 4 A schematic diagram of measurement results of an S-wave splitting measurement program according to an embodiment;

[0063] Figure 5 This is a schematic diagram of the artificial visualization measurement results;

[0064] Figure 6 This is a diagram showing the measurement results of the MFAST Manual V2.2 program. DETAILED DESCRIPTION

[0065] The specific implementation modes of the present invention are described below so that those skilled in the art can understand the present invention. However, it should be clear that the present invention is not limited to the scope of the specific implementation modes. For those of ordinary skill in the art, as long as various changes are within the spirit and scope of the present invention as defined and determined by the attached claims, these changes are obvious, and all inventions and creations utilizing the concept of the present invention are protected.

[0066] like Figure 1 As shown, a multi-method based automated measurement method for near-earthquake S-wave splitting includes the following steps:

[0067] S1, obtaining the original waveform data of the earthquake event and setting the noise window;

[0068] The original waveform data in step S1 is a three-component seismic record containing P-wave and S-wave initial motion signals of a near-seismic event, and its file format is SAC; wherein the three components are two horizontal components and a vertical component, namely, an N-direction component, an E-direction component, and a Z-direction component;

[0069] like Figure 2 As shown, step S1 includes the following steps:

[0070] S1-1, respectively picking up the arrival time of the initial motion signal corresponding to the S wave signal of the N-direction component and the P wave signal of the Z-direction component, and obtaining the arrival time of the corresponding S wave initial motion signal and the arrival time of the P wave initial motion signal;

[0071] S1-2, based on the arrival time of the P wave initial motion signal, selecting the time p of the P wave initial motion signal;

[0072] S1-3. The arrival time of the S wave initial motion signal is taken as time t1 and defined as the end position of the noise window. The time t0 before time t1 is selected as the starting position of the noise window to complete the setting of the noise window.

[0073] In order to calculate the signal-to-noise ratio of the S-wave signal and the noise signal in a period of time before its initial motion, a suitable time is selected before the initial motion of the S-wave signal as the starting position t0 of the noise window, and t1 as the end position of the noise window. The selection of t0 should make the signal-to-noise ratio of the S-wave signal higher, that is, the energy of the noise signal in this period is weaker.

[0074] S2, presetting a band-pass filter and screening the band-pass filter;

[0075] Step S2 includes the following steps:

[0076] S2-1, respectively analyzing the frequency spectra of the two horizontal components of the S wave signal to obtain the corresponding analysis results of the two horizontal components, that is, obtaining the main period of the S wave signal of the corresponding horizontal component; taking the average of the analysis results of the two horizontal components as the analysis result of the frequency spectrum of the S wave signal;

[0077] S2-2, based on the analysis results of the S wave signal spectrum, set M bandpass filters; where M is 10;

[0078] S2-3, using any bandpass filter to filter the S wave signals of the two horizontal components to obtain a filtered S wave signal;

[0079] S2-4, calculate the average energy of the filtered S-wave signal within a complete cycle;

[0080] S2-5, calculating the average energy of the filtered noise signal in the noise window;

[0081] S2-6, taking the ratio of the average energy of step S2-4 to the average energy of step S2-5 as the signal-to-noise ratio;

[0082] S2-7, determine whether the signal-to-noise ratio is not less than 3; if so, retain the bandpass filter and increase the value of m by 1; otherwise, discard the bandpass filter; proceed to step S2-8; wherein the initial value of m is 0;

[0083] S2-8, repeating steps S2-3 to S2-7 until all the band-pass filters are judged;

[0084] S2-9, determine whether the m value is 0; if so, determine that the original waveform data cannot meet the requirements of near-earthquake S-wave splitting and end the measurement; otherwise, enter step S3.

[0085] S3, divide the S wave analysis window;

[0086] The specific process of step S3 is as follows:

[0087] The starting position of the noise window is used as the starting position of the S wave analysis window, and the 0.5 cycle of the S wave initial motion signal is used as the end position of the first S wave analysis window; the end position of the S wave analysis window is moved backward with a step size of 0.001s as the end position of the next S wave analysis window, until it ends at the 1.2 cycle of the S wave initial motion signal, and Q S wave analysis windows of different lengths are obtained.

[0088] S4, using the screened bandpass filters and S-wave analysis windows to form multiple search networks in pairs, and using the cross-correlation function and the minimum tangential energy algorithm to search for the initial S-wave splitting parameters in different search networks;

[0089] Step S4 includes the following steps:

[0090] S4-1, combining the bandpass filter and the S-wave analysis window to obtain different search grids; wherein the number of search grids is the product of the number of the bandpass filter and the number of the S-wave analysis window;

[0091] S4-2, select any search grid, and rotate the S wave signals of the two horizontal components in the range of 0 to 180 degrees according to the step size of 1 degree to obtain fast wave signals and slow wave signals at different angles;

[0092] S4-3, in the interval of 0 to 0.1 s, time-correcting the fast wave signal and the slow wave signal obtained in step S4-2 according to a step length of 0.001 s, to obtain fast wave signals and slow wave signals after different time corrections;

[0093] like Figure 3As shown in the figure, N, F, R, E, S and T represent the due north direction, fast wave polarization direction, radial direction, due east direction, slow wave polarization direction and tangential direction respectively. In isotropic media, S waves do not split, and only radial signals can be observed in the horizontal direction. Assuming that the S wave signal is S(t), the radial component R(t) and the tangential component S(t) are obtained as follows:

[0094] R(t)=S(t)

[0095] T(t)=0

[0096] When the medium is anisotropic, the S wave will split, generating a pair of nearly orthogonal fast and slow waves with different velocities, so the calculation formula of step S4-3 is:

[0097] S F =S(t)·cosφ

[0098] S S =S(t-δt)·sinφ

[0099] Among them, S F , S S They represent the fast wave signal and slow wave signal after time correction, S(t) represents the composite signal of the S wave signals in the two horizontal components, cos and sin represent the cosine function and the sine function, t represents time, φ represents the angle between the fast wave polarization direction and the radial direction, δt represents the slow wave delay time, and S(t-δt) represents the S wave signal after a delay of δt.

[0100] S4-4, using the cross-correlation function to process the time-corrected fast wave signal and the slow wave signal to obtain the cross-correlation coefficients of the fast wave signal and the slow wave signal after different angles and time corrections; all the cross-correlation coefficients constitute a cross-correlation coefficient matrix;

[0101] S4-5, based on the cross-correlation coefficient matrix obtained in step S4-4, searching for the angles and correction times corresponding to the fast wave signal and the slow wave signal that make the waveforms most similar, and obtaining the fast wave polarization direction and the slow wave delay time, that is, the result of the cross-correlation function;

[0102] The calculation formula of step S4-5 is:

[0103]

[0104] Here, COV represents the mutual correlation coefficient, and |·| represents the absolute value.

[0105] S4-6, obtaining a header file of a three-component seismic record and determining whether the header file contains a variable BAZ; if so, using the variable BAZ as the reverse azimuth of the seismic event; otherwise, searching for the reverse azimuth of the seismic event based on the waveform of the three-component seismic record near the time p of the P wave first motion signal;

[0106] S4-7, based on the back azimuth of the seismic event obtained in step S4-6, the time-corrected fast wave signal and slow wave signal are calculated using the minimum tangential energy algorithm to obtain the corresponding radial and tangential signals, and the average energy of the tangential signal in the S-wave analysis window is calculated; all tangential signal energies constitute a tangential signal energy matrix;

[0107] The calculation formula of step S4-7 is:

[0108] S R =S F ·cosφ+S S ·sinφ

[0109] S T =S S ·cosφ-S F ·sinφ

[0110] Among them, S R represents the radial signal, φ represents the angle between the fast wave polarization direction and the radial direction, that is, the reverse azimuth, S T Indicates tangential signal.

[0111] S4-8, searching according to the tangential signal energy matrix obtained in step S4-7, determining the fast wave polarization direction and slow wave delay time with the minimum average energy of the tangential signal, that is, the search result of the minimum tangential energy algorithm;

[0112] S4-9, repeating steps S4-2 to S4-8 until the search results of all search grids are obtained, i.e., the initial S-wave splitting parameters;

[0113] The initial S-wave splitting parameters include the cross-correlation function and the search results of the minimum tangential energy algorithm.

[0114] S5. Determine the initial S-wave splitting parameters obtained by each search network, obtain the final S-wave splitting parameters, and complete the automated measurement of the near-earthquake S-wave splitting.

[0115] The specific process of step S5 is:

[0116] S5-1, respectively normalizing the mutual correlation coefficient matrix and the tangential signal energy matrix obtained by any search network to obtain a normalized mutual correlation coefficient matrix and a normalized tangential signal energy matrix;

[0117] S5-2, searching the normalized mutual correlation coefficient matrix to obtain the fast wave polarization direction and slow wave delay time with a correlation coefficient greater than 95%, and calculating the corresponding standard deviations respectively;

[0118] S5-3, calculating the sum of the two standard deviations in step S5-2;

[0119] S5-4, searching the normalized tangential signal energy matrix to obtain the fast wave polarization direction and slow wave delay time with energy less than 5%, and calculating the corresponding standard deviations respectively;

[0120] S5-5, calculating the sum of the two standard deviations in step S5-4;

[0121] S5-6, summing the standard deviations of step S5-3 and step S5-5 to obtain the standard deviation of all search networks;

[0122] S5-7, compare the standard deviations of different search networks, and select the initial S-wave splitting parameters of the search network with the smallest standard deviation;

[0123] S5-8. Compare the initial S-wave splitting parameters obtained in step S5-7; if the search results of the cross-correlation function and the minimum tangential energy algorithm are consistent, the search result of the minimum tangential energy algorithm is selected as the final S-wave splitting parameter; if the results of the two algorithms are inconsistent, the two search results are inspected by a method based on manual visual inspection and visualization, and selected as the final S-wave splitting parameter.

[0124] In addition, when the measurement results of the three methods are less consistent, the measurement results of the earthquake events used are eliminated.

[0125] In one embodiment of the present invention, Figure 4 As shown in the figure, the present invention searches for an S-wave splitting parameter of a near-seismic event. After the S-wave splitting parameter is corrected, the energy of the tangential signal is reduced to the noise level, the fast and slow waves show strong similarity, and the particle motion trajectory of the two horizontal components changes from the original approximately elliptical shape to an approximately linear polarization state. Figure 4 In the figure, Figure a is the three-component seismic record of the earthquake event used; Figures b, c, d, e and f are the measurement results of the minimum tangential energy algorithm (110°, 0.044s); Figures g, h, i, j and k are the measurement results of the cross-correlation function algorithm (109°, 0.044s). When the search results of the two algorithms are displayed using contour distribution maps, in order to maintain the consistency of the contour lines, the normalized tangential signal energy is taken inversely, that is, the peak of the contour line is the position of the minimum tangential energy.

[0126] against Figure 4 The earthquake events given in Figure 5 and Figure 6 Some disadvantages of the existing technology are shown respectively. Figure 5As shown in the figure, the particle motion trajectory of the S wave initial motion signal shows two close linear polarization directions. The artificial visualization measurement method is affected by subjective consciousness, and it may be considered that the fast wave polarization direction is 102° and 105°. The two horizontal components are rotated according to the fast wave polarization directions of 102° and 105° respectively to obtain the waveform of the "fast and slow waves". When rotated to the 102° coordinate system, the waveform of the "fast and slow waves" shows an obvious positive correlation, and the "slow wave delay time" is measured to be 0.025s. When rotated to the 105° coordinate system, the waveform of the "fast and slow waves" shows an obvious negative correlation, and the "slow wave delay time" is measured to be 0.04s. The "fast and slow waves" obtained in the above two cases have a strong correlation, but the S wave splitting parameters have a large difference. Therefore, the results based on the artificial visualization measurement method and the cross-correlation function algorithm are relatively unstable. When using this method to measure the S wave splitting of a near earthquake, the results obtained may have a large deviation from the true value. Figure 6 The measurement results of the MFASTManual V2.2 program show that this method cannot control the length of the shear wave splitting analysis window. An overly long analysis window will cause irrelevant signals behind the S-wave initial motion signal to participate in the calculation process of the S-wave splitting parameters, which will cause the program to search for incorrect results (106°, 0.062s).

[0127] By comparing the measurement results of the present invention with those of the prior art, the present invention can further eliminate abnormal points in the search results in combination with the cross-correlation function and the minimum tangential energy algorithm, thereby improving the accuracy and stability of the program; select a variety of filter ranges and analysis window lengths, thereby eliminating the instability of the search results caused by a single filter range and analysis window; limit the moving range of the analysis window by the length of the main period of the S-wave initial motion signal, thereby effectively avoiding the interference of the surface wave signal with relatively large energy on the measurement results; further screen the measurement results of the cross-correlation function and the minimum tangential energy algorithm by an artificial visualization measurement method, thereby improving the accuracy of the final measurement results and having strong adaptability, and facilitating a detailed analysis of the direction and intensity of the anisotropy of the upper crust medium in the study area; and in combination with relevant geophysical data, more accurately discussing scientific issues closely related to regional tectonic deformation and stress accumulation, which is of great significance for exploring dynamic processes such as deep structures and strong earthquake-pregnant environments.

Claims

1. A multi-method based automated measurement method for local earthquake S-wave splitting, characterized by: The following steps are involved: S1, obtaining the original waveform data of the earthquake event and setting the noise window; S2, presetting a band-pass filter and screening the band-pass filter; S3, divide the S wave analysis window; S4, using the screened bandpass filters and S-wave analysis windows to form multiple search networks in pairs, and using the cross-correlation function and the minimum tangential energy algorithm to search for the initial S-wave splitting parameters in different search networks; S5. Determine the initial S-wave splitting parameters obtained by each search network, obtain the final S-wave splitting parameters, and complete the automated measurement of the S-wave splitting of the near-earthquake; The original waveform data in step S1 is a three-component seismic record containing P-wave and S-wave initial motion signals of a near-seismic event, and its file format is SAC; wherein the three components are two horizontal components and a vertical component, namely, an N-direction component, an E-direction component, and a Z-direction component; The step S1 comprises the following steps: S1-1, respectively picking up the arrival time of the initial motion signal corresponding to the S wave signal of the N-direction component and the P wave signal of the Z-direction component, and obtaining the arrival time of the corresponding S wave initial motion signal and the arrival time of the P wave initial motion signal; S1-2, based on the arrival time of the P wave initial motion signal, selecting the time p of the P wave initial motion signal; S1-3, the arrival time of the S wave initial motion signal is taken as time t1, and defined as the end position of the noise window, and the time t0 before time t1 is selected as the starting position of the noise window, and the setting of the noise window is completed; The step S2 comprises the following steps: S2-1, respectively analyzing the frequency spectra of the two horizontal components of the S wave signal to obtain the corresponding analysis results of the two horizontal components, that is, obtaining the main period of the S wave signal of the corresponding horizontal component; taking the average of the analysis results of the two horizontal components as the analysis result of the frequency spectrum of the S wave signal; S2-2, based on the analysis results of the S wave signal spectrum, set M bandpass filters; where M is 10; S2-3, using any bandpass filter to filter the S wave signals of the two horizontal components to obtain a filtered S wave signal; S2-4, calculate the average energy of the filtered S-wave signal within a complete cycle; S2-5, calculating the average energy of the filtered noise signal in the noise window; S2-6, taking the ratio of the average energy of step S2-4 to the average energy of step S2-5 as the signal-to-noise ratio; S2-7, determine whether the signal-to-noise ratio is not less than 3; if so, retain the bandpass filter and increase the value of m by 1; otherwise, discard the bandpass filter; proceed to step S2-8; wherein the initial value of m is 0; S2-8, repeating steps S2-3 to S2-7 until all the band-pass filters are judged; S2-9, determine whether the m value is 0; if so, determine that the original waveform data cannot meet the requirements of near-earthquake S-wave splitting and end the measurement; otherwise, enter step S3.

2. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 1, characterized in that: The specific process of step S3 is as follows: The starting position of the noise window is used as the starting position of the S wave analysis window, and the 0.5 cycle of the S wave initial motion signal is used as the end position of the first S wave analysis window; the end position of the S wave analysis window is moved backward with a step size of 0.001s as the end position of the next S wave analysis window, until it ends at the 1.2 cycle of the S wave initial motion signal, and Q S wave analysis windows of different lengths are obtained.

3. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 1, characterized in that: The step S4 comprises the following steps: S4-1, combining the bandpass filter and the S-wave analysis window to obtain different search grids; wherein the number of search grids is the product of the number of the bandpass filter and the number of the S-wave analysis window; S4-2, select any search grid, and rotate the S wave signals of the two horizontal components in the range of 0 to 180 degrees according to the step size of 1 degree to obtain fast wave signals and slow wave signals at different angles; S4-3, in the interval of 0 to 0.1 s, time-correcting the fast wave signal and the slow wave signal obtained in step S4-2 according to a step length of 0.001 s, to obtain fast wave signals and slow wave signals after different time corrections; S4-4, using the cross-correlation function to process the time-corrected fast wave signal and the slow wave signal to obtain the cross-correlation coefficients of the fast wave signal and the slow wave signal after different angles and time corrections; all the cross-correlation coefficients constitute a cross-correlation coefficient matrix; S4-5, based on the cross-correlation coefficient matrix obtained in step S4-4, searching for the angles and correction times corresponding to the fast wave signal and the slow wave signal that make the waveforms most similar, and obtaining the fast wave polarization direction and the slow wave delay time, that is, the result of the cross-correlation function; S4-6, obtaining a header file of a three-component seismic record and determining whether the header file contains a variable BAZ; if so, using the variable BAZ as the reverse azimuth of the seismic event; otherwise, searching for the reverse azimuth of the seismic event based on the waveform of the three-component seismic record near the time p of the P wave first motion signal; S4-7, based on the back azimuth of the seismic event obtained in step S4-6, the time-corrected fast wave signal and slow wave signal are calculated using the minimum tangential energy algorithm to obtain the corresponding radial and tangential signals, and the average energy of the tangential signal in the S-wave analysis window is calculated; all tangential signal energies constitute a tangential signal energy matrix; S4-8, searching according to the tangential signal energy matrix obtained in step S4-7, determining the fast wave polarization direction and slow wave delay time with the minimum average energy of the tangential signal, that is, the search result of the minimum tangential energy algorithm; S4-9, repeating steps S4-2 to S4-8 until the search results of all search grids are obtained, i.e., the initial S-wave splitting parameters; The initial S-wave splitting parameters include the cross-correlation function and the search results of the minimum tangential energy algorithm.

4. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 3 is characterized by: The calculation formula of step S4-3 is: S F =S(t)·cosφ S S =S(t-δt)·sinφ Among them, S F , S S They represent the fast wave signal and slow wave signal after time correction, S(t) represents the composite signal of the S wave signals in the two horizontal components, cos and sin represent the cosine function and the sine function, t represents time, φ represents the angle between the fast wave polarization direction and the radial direction, δt represents the slow wave delay time, and S(t-δt) represents the S wave signal after a delay of δt.

5. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 4 is characterized by: The calculation formula of step S4-4 is: Here, COV represents the mutual correlation coefficient, and |·| represents the absolute value.

6. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 4, characterized in that: The calculation formula of step S4-7 is: S R =S F ·cosφ+S S ·sinφ S T =S S ·cosφ-S F ·sinφ Among them, S R represents the radial signal, S T Indicates tangential signal.

7. The multi-method based automated measurement method for near-earthquake S-wave splitting according to claim 3 is characterized by: The specific process of step S5 is as follows: S5-1, respectively normalizing the mutual correlation coefficient matrix and the tangential signal energy matrix obtained by any search network to obtain a normalized mutual correlation coefficient matrix and a normalized tangential signal energy matrix; S5-2, searching the normalized mutual correlation coefficient matrix to obtain the fast wave polarization direction and slow wave delay time with a correlation coefficient greater than 95%, and calculating the corresponding standard deviations respectively; S5-3, calculating the sum of the two standard deviations in step S5-2; S5-4, searching the normalized tangential signal energy matrix to obtain the fast wave polarization direction and slow wave delay time with energy less than 5%, and calculating the corresponding standard deviations respectively; S5-5, calculating the sum of the two standard deviations in step S5-4; S5-6, summing the standard deviations of step S5-3 and step S5-5 to obtain the standard deviation of all search networks; S5-7, compare the standard deviations of different search networks, and select the initial S-wave splitting parameters of the search network with the smallest standard deviation; S5-8. Compare the initial S-wave splitting parameters obtained in step S5-7; if the search results of the cross-correlation function and the minimum tangential energy algorithm are consistent, the search result of the minimum tangential energy algorithm is selected as the final S-wave splitting parameter; if the results of the two algorithms are inconsistent, the two search results are inspected by a method based on manual visual inspection and visualization, and selected as the final S-wave splitting parameter.

Images