On Analysis of Seismic Vibrations Data Applying Doppler Effect Expression

In the paper, a possibility to develop the digital models of the seismic vibrations parameters is analyzed. To reach this goal, the observations at seismic station LUWI (Indonesia) were processed applying the statistical procedures. In fact, the biggest attention was given to the introduction of the Doppler effect expression and the employment of the theory of covariance functions. +e trend in vectors of vibrations intensities values was detected and estimated upon using the least-squares method and polynomial approximation. In addition, by this technique, the random errors were eliminated partially. +e self-developed computer programs based on Matlab programming package procedures were applied.


Introduction
Earthquake is one of the most costly, devastating, and deadly natural hazards. Every disaster damages thousands of buildings and displaces tens of thousands of people. e comprehensive knowledge of the earthquake nature and its behavior is extremely important. Here the main task for scientists is to constrain the suitable mathematical methods to analyse the earthquakes action, and most importantly to develop the earthquake model to predict its spread and to forecast its occurrence. e latest developments could be noted in [1][2][3][4][5][6][7], where the biggest efforts were taken for mathematical descriptions of wide earthquakes occurrence areas trying to construct the Ground Motion Prediction Equations. Deep analysis of different aspects of passive seismic methods like a horizontal to vertical spectral ratio, which is often used to describe the earthquake site, could be found in [8][9][10][11][12][13][14]. Some characteristics of concrete earthquake's sites from world's seismic zones are presented in [9,[15][16][17][18][19][20][21][22]. In some papers, the stress was done on significance of three-dimensional modelling of seismic waves propagation [18,[23][24][25][26][27].
What deals with scientific techniques and methods to investigate earthquakes application has InSAR technology, which enables detecting surface slips and Earth surface deformations [34,37,41]. For example, the 4-7 m surface slip in the area of Palu earthquake was detected [37] and the maximum horizontal deformation was from 1.8 m till 3.6 m [41], when ALOS-2 interferogram showed a peak slip of 6.5 m located at the south of Palu city [34]. GNSS plays a great role in the research of earthquakes giving very precise metrical parameters to improve the crustal deformation field and 3D geometric complexities of the faults in total [38,53,60,65,68]. Certainly, the main techniques to detect the technical parameters of earthquakes are seismograms and the combinations of some techniques as well [34,67]. So, from broadband regional seismograms, it was revealed that the 2018 Palu earthquake is a supershear rupture event from early on with an average rupture velocity of 4.1 km/s, and the total seismic moment of 2.64 × 1020 Nm (equivalent to Mw 7.55) was released within 40 s [34].
In this paper, we will show how the Doppler effect expression and the application of the theory of covariance functions could be employed for seismic waves modelling. e practical calculations were executed using the two fragments of the observations data of the intensity φ of the Earth's seismic field, which were chosen from LUWI seismic station (Sulawesi, Indonesia, latitude: −1.04180, longitude: 122.77170, elevation: 6.0 m): first on August 05, 2018, within one hour (11:30-12:30), and second on November 11, 2018, within two hours (5:00-7:00). At these periods, the seismic stations around the world have registered unusual vibrations of low frequencies. Wide basic information on Palu earthquake could be found in specialized portals [71,72]. e observations data were expressed by vectors N (North), E (East), and Z (Zenith). e time series views of the centered vectors N, E, and Z for both abovementioned periods are presented in Figures 1 and 2. In both figures, the time series views of the components E and N are similar. It looks like the influence of the unusual low frequencies vibrations in Figure 2 is possibly low. e systematic component of low frequency could be eliminated applying the 6degree polynomial approximation. It is presented in Figure 3. e accuracy of vectors N, E, and Z extracted from LUWI station data on August 05, 2018, is described by standard deviations S φ � (19480, 15926, 15810) cnt. ese numbers show that the accuracies of components of seismic vectors presented in Figure 1 are approximately the same. e accuracy of vectors N, E, and Z extracted from LUWI station on November 11, 2018 (if the systematic component is not eliminated), is described by the vector of standard deviations S φ � (1260, 559, 41) cnt. It shows that the accuracy of observations is slightly higher at this period. e accuracy of vectors N, E, and Z extracted on November 11, 2018 (if the systematic component is eliminated), is described by the vector of standard deviations S φ � (215, 257, 41) cnt. In this case, the obtained accuracy of processed observation data is considerably higher.
Mathematical-statistical methods are widely applied for data processing in geophysics, geodesy, and other Earth sciences [73][74][75]. To predict and develop the model of the spread of seismic vibrations, first of all, we assume that seismic waves from the quake hypocenter spread as harmonic vibrations of decreasing amplitudes in all the directions. So, we can assume also that the core structures of seismic observations at the tracking stations mounted in short distances from the hypocenter and at those more distant are possibly very similar. For mathematical treatment of the seismic observations, the covariance functions and the theory of Doppler effect were applied. e correlations between changes of intensities of seismic waves spreading in time and space were detected by introducing the variations of covariations of the seismic vibrations intensities vectors. Some equations were derived to obtain the estimates of covariation matrixes and autocovariances and cross-covariances of seismic field intensities vectors based on seismic observations data. e accuracies of corresponding calculated parameters were obtained also. e background of the mathematical model of observations data treatment is concept of a stationary random function and especially paying attention to statement that the errors of seismic vibrations observations are random and possibly are near the same precision. So we assume that the mathematical average of random errors MΔ � constant ⟶ 0, its dispersion DΔ � constant, and the covariances of the observations depend on the difference of the arguments only, so practically from the quantised intervals on the time scale.

Modelling of Seismic Vibrations
e observation data registered by the seismic station had been previously examined and processed upon reaching a goal to eliminate both random and possibly systematic errors. e most reliable values of the trend in the seismic vibrations arrays were detected employing the least-squares method. Application of least-squares technique gives a possibility to eliminate the random errors partially. While treating the big volumes of observations, the least-squares technique produces the asymptotically efficient values of the derived parameters also in case when a statistical distribution of the observations errors is not normal.
Any vector of seismic vibrations intensities could be treated as a random function, which involves the random errors of observations. By employing a least-squares technique to treat the vector of intensities φ, we can detect the most reliable value φ of the trend. A parametric equation of a single vector's element φ i will look like the following: where ε i is a random error of the vector's element, φ i is the value of the vector's element, and φ is vector's trend. e expression in matrix form of equation (1) will be as follows: where ε is vector of random errors, φ � (φ 1 , φ 2 , . . . , φ n ) T is vector of seismic field intensities, and e is vector of units (n × 1). e most reliable value of vector φ trend could be calculated by introducing the general condition of the leastsquare method: where P is diagonal matrix (n × n) of weights p i of the values φ i . Weights p i could be detected according to simple formula: where σ 0 is the standard deviation of the observation φ 0 , the weight of which is supposed to be equal to unit, that is, p 0 � 1. Furthermore, we can write the following equation: and we further obtain From formula (6), we can see that the value of σ φi pertains from the value of φ i . So, the components, which have the bigger values, are of a lower accuracy just because φ i ≫ σ ui .
Upon applying formula (4), we write where the accepted average value is σ 2 0 /σ 2 ui � 5 · 10 4 . To find the extremum of function (3), let us calculate its partial derivatives according to trend φ. We can write and solve the equation: en we will obtain −e T Pε � 0, e T Peφ − e T Pφ � 0.
us, we will get the following solution: where N � (e T Pe), ω � e T Pφ. e accuracy of the trend could be detected by calculating its covariance matrix K φ : where σ 0 ′ is the estimate of the standard deviation σ 0 . It is assessed by formula: e considerably high systematic component of the vibrations of the data of seismic station LUWI was eliminated upon applying a 6-degree polynomial approximation. Now, it is possible to calculate cross-covariance and autocovariance functions of the seismic vibrations as well as the shifts of the seismic vibrations, respectively, to each other by introducing the Doppler effect expression. Let us take the formula for the parameter z [76][77][78]: where f e is frequency of emitted vibrations and f o is frequency of observed vibrations. We accept that the changes of vibrations phases observed at tracking stations possibly correspond to the changes of the seismic vibrations intensities. Consequently, sum of the seismic vibrations intensities is proportional to the algebraic sum of the frequencies phases of vibrations accordingly; that is, δB ∼ δω; (14) here δB, δω are changes of vibrations intensities and vibrations frequencies phases, respectively. We can write the expressions for changes of vibrations intensities and the sum of them as follows: where ω e � 2πf e , 2πf o , and the initial phases φ 0 are supposed to be equal to zero; δa ⟶ δB.
By employing the parameter z of the Doppler effect formula, we can express the strength of the seismic vibrations at the moment in time t i : where B ei is the intensity of emitted seismic vibrations, B oi is the intensity of observed seismic vibrations, B ei ∼ ω ei , and B oi ∼ ω oi . In further developments, we employ the theory of covariance functions to detect the value of the argument z from the Doppler effect expression. Mathematical derivations are grounded on the conception of a stationary random function considering that errors of observations of seismic vibrations are random and possibly have similar precision.
It is possible to express a cross-covariance function of the straight algebraic sum ΔB ei � B ei − B oi of the two intensities B ei and B oi (emitted and observed) at a moment in time t i and a separate intensity B oi as follows: where By employing the theory of the covariance functions, it is possible to express the cross-covariance functions of the corresponding seismic vibrations vectors taking into account the fact that every vector of vibrations intensities could be treated as a random function as follows [77,79,80]: where u is argument of any seismic vibrations vector, τ � s · Δ is quantised interval, which is variable, s is number of quantised intervals, Δ is the value of the accepted unit of observations, and T is the diapason of the fluctuations of seismic vectors elements.
By using the vectors of observations data, an estimation K z ′ (τ) of the cross-covariance function could be calculated according to the following formula: where n is number of vector elements. Now, using formula (18) in the vector form, we get the formula to detect the mathematical average of the argument z of the Doppler effect expression: where σ ′2 B ⟶ K z ′ (0) is the estimate of the dispersion and m is number of cross-covariance values.

Analysis of the Experimental Results
e estimates of autocovariance and cross-covariance functions of the seismic vibrations intensities could be calculated employing formula (20). e values of the quantised intervals were assigned from 1 to n/2. Here, n � 144000 is the number of seismic vibrations vector components. e graphical images of autocovariance and cross-covariance functions were generated also. Some graphical images of covariance functions are shown in

Advances in Civil Engineering
Let us derive the estimates of the standard deviation of the argument z. It could be done in two ways. Firstly, z values could be detected upon applying the frequencies of vibrations emitted from earthquake source and observed frequencies of vibrations at seismic station. Secondly, z values could be detected using the vibrations intensities. e formula to calculate the estimate of standard deviation of the argument z using formula (13) could be written as follows: In formula (22), the estimate σ z ′ of the standard deviation of the argument z was detected assuming σ fe � σ fo , z � 1.0, and σ fe /σ fo � 1.4 · 10 − 8 . So, the ratio of the mathematical average of the argument z is σ z ′ /z � 3.0 · 10 − 8 .
Let us derive the accuracy of the argument z of the Doppler effect expression using the seismic vibrations intensities. Upon applying formulas (15), (16), and (22), we have and the above was calculated upon considering and z � 1.0. e ratio of z B will be σ z B ′ /z B � 3 · 10 − 8 . e results of the calculations demonstrate that the detected accuracies of the argument z of the Doppler effect expression are nearly the same in both cases: when registered phases of vibrations frequencies or intensities (strengths) of them are used in the calculation procedures.
Let us calculate the accuracy of the estimates of the motion speed v � z · c. Upon using the equation we can write the equation of the ratio as follows: and (σ v /v) ≈ 3.0 · 10 − 9 , when (σ c /c) � 3.0 · 10 − 9 .

Data Availability
All data used during the research are available in a repository online in accordance with funder data retention policies. Practically, data for this research were taken from the EIDA and GEOFON Data Archives (http://eida.gfz-potsdam.de/ webdc3/).

Disclosure
e research was done as a part of employment at Vilnius Gediminas Technical University, Lithuania.

Conflicts of Interest
e authors declare that there are no conflicts of interest.