5.2.8.1.5. eqcorrscan.utils.mag_calc.svd_moments

eqcorrscan.utils.mag_calc.svd_moments(u, s, v, stachans, event_list, n_svs=2)[source]

Calculate relative moments/amplitudes using singular-value decomposition.

Convert basis vectors calculated by singular value decomposition (see the SVD functions in clustering) into relative moments.

For more information see the paper by Rubinstein & Ellsworth (2010).

Parameters:
  • u (list) – List of the numpy.ndarray input basis vectors from the SVD, one array for each channel used.

  • s (list) – List of the numpy.ndarray of singular values, one array for each channel.

  • v (list) – List of numpy.ndarray of output basis vectors from SVD, one array per channel.

  • stachans (list) – List of station.channel input

  • event_list (list) – List of events for which you have data, such that event_list[i] corresponds to stachans[i], U[i] etc. and event_list[i][j] corresponds to event j in U[i]. These are a series of indexes that map the basis vectors to their relative events and channels - if you have every channel for every event generating these is trivial (see example).

  • n_svs (int) – Number of singular values to use, defaults to 4.

Returns:

M, array of relative moments

Return type:

numpy.ndarray

Returns:

events_out, list of events that relate to M (in order), does not include the magnitude information in the events, see note.

Return type:

obspy.core.event.event.Event

Note

M is an array of relative moments (or amplitudes), these cannot be directly compared to true moments without calibration.

Note

When comparing this method with the method used for creation of subspace detectors (Harris 2006) it is important to note that the input design set matrix in Harris contains waveforms as columns, whereas in Rubinstein & Ellsworth it contains waveforms as rows (i.e. the transpose of the Harris data matrix). The U and V matrices are therefore swapped between the two approaches. This is accounted for in EQcorrscan but may lead to confusion when reviewing the code. Here we use the Harris approach.

Example

>>> from eqcorrscan.utils.mag_calc import svd_moments
>>> from obspy import read
>>> import glob
>>> import os
>>> from eqcorrscan.utils.clustering import svd
>>> import numpy as np
>>> # Do the set-up
>>> testing_path = 'eqcorrscan/tests/test_data/similar_events_processed'
>>> stream_files = glob.glob(os.path.join(testing_path, '*'))
>>> stream_list = [read(stream_file) for stream_file in stream_files]
>>> event_list = []
>>> remove_list = [('WHAT2', 'SH1'), ('WV04', 'SHZ'), ('GCSZ', 'EHZ')]
>>> for i, stream in enumerate(stream_list):
...     st_list = []
...     for tr in stream:
...         if (tr.stats.station, tr.stats.channel) not in remove_list:
...             stream.remove(tr)
...             continue
...         st_list.append(i)
...     event_list.append(st_list) 
>>> event_list = np.asarray(event_list).T.tolist()
>>> SVec, SVal, U, stachans = svd(stream_list=stream_list) 
['GCSZ.EHZ', 'WV04.SHZ', 'WHAT2.SH1']
>>> M, events_out = svd_moments(u=U, s=SVal, v=SVec, stachans=stachans,
...                             event_list=event_list)