Time-delay calculation

We can look at the filtered, tapered time series for a bunch of the stations for a single event (x-axis is time, y-axis is microns)




We can calculate the arrival time-delay between different stations using numpy.correlate. You first need to zero pad on either side so that you can slide one array next to the other. An example code snippet is below:

import numpy as np
arr1 = np.asarray([1,1,2,1])
arr2 = np.asarray([1,2,1,1])
times = np.arange(arr1.size)
ARR1 =\
    np.hstack((np.zeros(arr1.size), arr1, arr1.size))                      
corr = np.correlate(ARR1,arr2)
times = np.hstack((-times[::-1],times))

We get the time-delay between the two stations from the maximum of that correlation spectrum. 


Local Velocity

We want to get the velocity in our array by using the time-delay in arrival at different seismometers.




where tau is the measured time delay. We get phi using the azimuth spit out by geographiclib. 

In general, I use
as “North of East.” This means that 

.

It looks like α2 = 67o, which means φ=203o. The following pseudocode could do this (where stat1_data and stat2_data are seispy.trace.Trace objects):

import numpy as np
from seispy.trace import Trace

######## FETCH DATA ##########
# put something here
##############################

# do correlation
ARR1 =\
    np.hstack((np.zeros(stat1_data.size), stat1_data.value, stat1_data.size))                      
corr = np.correlate(ARR1,stat2_data.size)

# get times
times = stat1_data.times.value
times -= times[0]
times = np.hstack((-times[::-1],times))

# get locations from catalog in seispy
x_1 = stat1_data.get_location()
x_2 = stat2_data.get_location()
delta_x_vec = x_2 - x_1
omega = np.asarray([np.cos(np.radians(203)), np.sin(np.radians(203)), 0])
delta_x = np.dot(omega, delta_x_vec)



Method [Using only 1 event]
Frequency [Hz]
channel
velocity [m/s]
velocity low estimate [m/s]
velocity high estimate [m/s]
two channel [above]
0.4
HHZ
2982.7 +/- 53.3
2929.4
3036.0
arrival from location estimate
0.4
HHZ
1400



Some Issues

If the amplitude is low (which is true for HHR and HHT sometimes), then one doesn’t get an obvious maximum and often you end up with time-delays that are off by 1 period. Then you end up with three parallel lines in the above plot. Often, these amplitudes are dependent upon the depth of the station, so one will sometimes end up with these three parallel lines separating based on the vertical distance between stations (hence the color bar). I will sometimes make cuts based on this (insist time-lags must be less than 1 period). One could also envision just only considering time-lags modulo the period we’re looking at.