An algorithm of local earthquake detection from digital records

The problem of automatical detection of earthquake signals in seismograms and definition of first arrivals of p and s waves is considered. The algorithm is based on the analysis of t(A) function which represents the time of first appearence of a number of going one after another swings of amplitudes greather than A in seismic signals. It allows to explore such common features of seismograms of earthquakes as sudden first p-arrivals of amplitude greater than general amplitude of noise and after the definite interval of time before s-arrival the amplitude of which overcomes the amplitude of p-arrival. The method was applied to 3-channel recods of Friuli aftershocks, ?'-arrivals were defined correctly in all cases; p-arrivals were defined in most cases using strict criteria of detection. Any false signals were not detected. All p-arrivals were defined using soft criteria of detection but less reliability and two false events were obtained.

Il metodo è stato applicato alle registrazioni su tre componenti delle repliche del periodo sismico del Friuli (1976). Gli arrivi delle onde s furono definiti correttamente in tutti i casi; gli arrivi delle p furono definiti nella maggioranza dei casi usando criteri restrittivi per la detezione. Tutti gli arrivi delle p furono definiti adoperando criteri più larghi di detezione, ma si sono ottenuti due falsi allarmi e minore accuratezza. INTRODUCTION The problem of automatic seismic signal detection arises in literature from time to time. Usually, it is solved by the methods of stochastic process extrapolation. Here it is presented another approach to this problem based on the general features common for local seismic signals: sudden p and s arrivals, separated by definite time intervals and so on. It is supposed that amplitude of p-waves is larger than general level of background noise and s-amplitude is larger than p-amplitude. To eliminate random spikes in noise we consider minimal amplitude in a group of neighbour swings and construct an algorithm for automatic identification of prominent jumps of this characteristic as p and 5 arrivals. The algorithm was applied to 3-channel digital records of Friuli aftershocks supplied by the Istituto Nazionale di Geofisica, Roma.

DESCRIPTION OF SEISMIC SIGNALS BY t (A)-FUNCTION.
Let us denote by U, (i), j = 1, 2, 3 the displacements on the vertical and two horizontal instruments respectively in a time interval (T,, T 2 ), 0 ^ t ^ T 2 -7Y To exclude drift of zero of the instrument and make procedure independent of amplification factor, records were normalized as follows where JJj and a,-2 are the mean value and variance for channel j.
where {tk,j} is the sequence of zeroes of u, (i) in the interval CTi, Ti).
We assume that seismic signals are different from random spikes in noise by occurrence of a few prominent swings of positive and negative value following each other. For this reason we shall consider the minimum of n sequential maximums bk,j Now we define function A,-(t) as follows is a monotonous step-function, which in some sense is the left envelope of the seismic signal as it is shown in Fig. 1.
This function allows us to detect p and s arrivals under assumption that value minimax coresponding to p-waves is larger than all minimaxes in background noise before p arrival and minimax corresponding to s wave is larger than all minimaxes before 5 arrival including p-wave minimax.

As
Apm.x A P For the sake of convenience we switch now to the inverse function of Aj(t) and summarize results for three channels in function t (4) The representation of seismic signal by /(^4) function is sufficient for our algorithm of p and 5 detection.

RIVALS
For detection of seismic signals we need some volume of background noise placed in time before the real signal. Therefore we are looking for the maximum of seismic signal only in a second part of time interval T, T -Ti-T1, sliding with a step T/2. In other words, the first condition for the presence of the seismic singal in interval T is Detection of s first arrival is based on two assumptions: (ts-tp) should be rather big comparatively to period of 5 waves and/or ratio of amplitudes A smax /A pmax should be prominent. We combine this two conditions in one criterium defining:

where -¿(A)] is the jump of step function t(A) in the point A.
The largest value of A, which satisfies condition where x so is a parameter to be empirically determined, we assume as 4 pma x (it may be also the minimax amplitude of some other phase before s arrival). Correspondingly, we assume where t s is the time of s first arrival. Evidently, the amplitude of s first arrival is = A(t s ) . [11] Similar criterium we use for detection of p first arrival A and find t p and A p when To improve accuracy of s and p first arrival determinations we use as an estimate of t s or t p in the minimax satisfying condition [9] and [13] the following formula where ik + i.j and ik + 2"• are the second and the third zero in the definition of minimax [2-3], Some random disturbances may have maximum minimax amplitude in the second part of the interval T, so the condition [7] fails to reject disturbances in such case.
There is some probability that both conditions [9] and [13] will be satisfied by chance. Those cases we dicard with the help of two additional conditions which are naturally satisfied for all real signals. The first one is the condition of limited distance between the source and the seismic station, or in other words The second condition is the experimental observation that usually amplitude of Amax is much larger than amplitude of p first arrival A p

REALIZATION OF THE ALGORITHM AND SELECTION OF PARAMETERS
The algorithm was applied to analog magnetic records of 54 Friuli aftershocks collected by the Istituto Nazionale di Geofisica a few days after the main shock in May 1976. The records were digitized on a special device (Lennartz Electronic) with preliminary selection of seismic events. Records were digitized when amplitudes overcome definite level in two points separated by interval of time Xi ^ x ^ x 2 . This hardware device was designed for catching, s-waves which usually continue rather long time comparatively with random spikes in the background noise.
The device has a delay network which allows to digitize 1200 points in advance of the signal which overcomes the threshold level (each interval between two points is equal to 4/150 sec.). After the last overcome of the threshold the device digitizes about 1800 points of the record.
To imitate the original, preselected situation in the problem of signal-noise discrimination and detection of 5 and p arrivals we joined these digitized pieces of the record corresponding to different aftershocks in one continuous record. To this record we applied our algorithm with sliding time intervals T of 1000 points.
Step of sliding was equal to 500 points.
Different intervals contained pure noise, or noise with the signal in the second part of the interval, or some disturbances, for example radio-impulses of a nearby station, or near parts of seismic signals already detected in previous intervals. The preselection of the record by the hardware device essentially enriches our record comparatively to the original one which containes enormous amount of pure noise.
For each interval T we calculate mean value and variance for each of these channels and normalize data by formula [1], Then we construct sequence {bk") of maximums of absolute values of signals between two neighbour zeroes tk, tk + i • Instead of zeroes we took as points of the sign change of the signal where tk is an integer number, measuring time in 4/150 seconds. This substitution does not affect the definition of the sequence bk,j • The sequence bk,j was then transformed by formulae [3-6] into function t (4) for discrete values of A with step A = 1 (is measured in a units of root mean square fluctuations of record in the interval T). We have tried different values of n = 2, 3 and 4 in formula [3] and found out that n = 3 is the most suitable value.
To find thresholds x so and x po , functions (A) and x p (A) were analysed (Fig. 2). From Figs. 2-3 it is evident that all open circles denoting first arrivals of s and p are situated to the right from solid points denoting ordinary discontinuities of t (4). Moreover, open circles are separated from solid points by rather big intervals which separate open circles from solid points for almost all of the earthquakes considered in this paper. The values Xso = 60 and x po = 150 were chosen from these common intervals as thresholds for criterium [9] and [13] respectively.
The maximum observed value of t s -t p = 250 equal to 6.6 sec) was fixed for the threshold x Q in criterium [14].
Usually the difference between maximum amplitude of s waves Ajmax and the amplitude of p first arrival A p is more than 3 a. We put the threshold A a = 1.5 in criterium [16] to discard random disturbances which are not increasing with time so much and have amplitude of « first arrival » very close to the maximum amplitude.  Tables I-II. When s arrival was detected but p was not found we made the second run of the algorithm with less severe conditions on p waves. We decreased the amount of noise analysel before s arrival. The new t (A) function was calculated in the interval of 200 point before the already detected s arrival. The algorithm of 200 point before the already detected 5 arrival. The algorithm value of x P o = 50.

RESULTS AND DISCUSSIONS
Algorithm described in section 4 was applied to the record of 158,000 points (about 70 minutes), which was divided into 315 overlapping intervals T of 1000 points (T = 26 sec). The record contained 54 seismograms of Friuli aftershocks of different magnitude.
In the first run of the algorithm 45 earthquakes were detected, of which two events were contained in the coda of earthquakes which triggerd the hardware device. No false alarms were detected in samples of random noise, the rear part of seismograms and any other kinds of disturbances.
For small earthquakes in some cases it is difficult to estimate whether the p determination was correct or not because their detection is also difficult for a geophysicist. The amplitude of such p-waves is lower than maximum amplitude in background noise and the algorithm get them only through the minimax criterium [2, 3].
S arrivals were correctly detected for all earthquakes but in some samples without earthquakes there were also detected 5 waves by mistake.
Such samples were discarded because the algorithm did not find an apropriate p arrival satisfying the conditions [15] and [16].
In the second run with less severe conditions on p waves all but one earthquakes were detected. The undetected earthquake had very large amplitude of s waves which was followed by the instrument saturation. That is why the mean values U, were calculated with significant errors and all the record before s arrival including noise became much larger than zero. Therefore the construction of bk,j and detection of p became impossible. Recalculating Uj in the second run over the part of the T interval which preceded the s arrival could probably allow to detect this earthquake.
P arrivals in 8 cases of the previous 45 essentially changed, that means lower reliability of the second run due to more soften criteria. Also two false events appeared in the second run.
Results are presented in Table I which contains maximum amplitudes of 5, p and the noise before p mesaured in units of root mean square fluctuations a calculated for the first 500 points of the interval T.
In figures 4-10 seismograms which were detected by our program as earthquakes, including false alarms, are plotted. In these figures amplitudes are normalized in units of root mean square of noise.
Finally, we may conclude that combination of the first and the second runs of the algorithm provide detection of earthquakes in rather a wide range of amplitudes, although we should notice that in the cases of disagreement between the two runs in p arrivals, the first run is more favourable and further analysis is needed.