Utils

Codes to run basic utility functions for integration with seisan and to find peaks in noisy data.

Sfile_util

Part of the EQcorrscan module to read nordic format s-files and write them EQcorrscan is a python module designed to run match filter routines for seismology, within it are routines for integration to seisan and obspy. With obspy integration (which is necessary) all main waveform formats can be read in and output.

Code generated by Calum John Chamberlain of Victoria University of Wellington, 2015.

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

class Sfile_util.EVENTINFO(time, loc_mod_ind, dist_ind, ev_id, latitude, longitude, depth, depth_ind, loc_ind, agency, nsta, t_RMS, Mag_1, Mag_1_type, Mag_1_agency, Mag_2, Mag_2_type, Mag_2_agency, Mag_3, Mag_3_type, Mag_3_agency)[source]

Header information for seisan events

class Sfile_util.PICK(station, channel, impulsivity, phase, weight, polarity, time, coda, amplitude, peri, azimuth, velocity, AIN, SNR, azimuthres, timeres, finalweight, distance, CAZ)[source]

Pick information for seisan implimentation

Sfile_util.blanksfile(wavefile, evtype, userID, outdir, overwrite=False, evtime=False)[source]

Module to generate an empty s-file with a populated header for a given waveform.

Parameters:
  • wavefile (String) – Wavefile to associate with this S-file, the timing of the S-file will be taken from this file if evtime is not set
  • evtype (String) – L,R,D
  • userID (String) – 4-charectar SEISAN USER ID
  • outdir (String) – Location to write S-file
  • overwrite (Bool) – Overwrite an existing S-file, default=False
  • evtime (UTCDateTime) – If given this will set the timing of the S-file
Returns:

String, S-file name

Sfile_util.populateSfile(sfilename, picks)[source]

Module to populate a blank nordic format S-file with pick information, arguments required are the filename of the blank s-file and the picks where picks is a dictionary of picks including station, channel, impulsivity, phase, weight, polarity, time, coda, amplitude, peri, azimuth, velocity, SNR, azimuth residual, Time-residual, final weight, epicentral distance & azimuth from event to station.

This is a full pick line information from the seisan manual, P. 341

Parameters:
  • sfilename (str) – Path to S-file to populate, must have a header already
  • picks (List of :class: PICK) – List of the picks to be written out
Sfile_util.readheader(sfilename)[source]

Fucntion to read the header information from a seisan nordic format S-file.

Parameters:sfilename (str) – Path to the s-file
Returns:
class:EVENTINFO
Sfile_util.readpicks(sfilename)[source]

Function to read pick information from the s-file

Parameters:sfilename (String) – Path to sfile
Returns:List of :class: PICK
Sfile_util.readwavename(sfilename)[source]

Convenience function to extract the waveform filename from the s-file, returns a list of waveform names found in the s-file as multiples can be present.

Parameters:sfilename (str) – Path to the sfile
Returns:List of str
Sfile_util.test_rw()[source]

Function to test the functions herein.

findpeaks

Function to find peaks in data above a certain threshold as part of the EQcorr package written by Calum Chamberlain of Victoria University of Wellington in early 2015.

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

findpeaks.find_peaks2(arr, thresh, trig_int, debug=0, maxwidth=10, starttime=UTCDateTime(1970, 1, 1, 0, 0), samp_rate=1.0)[source]

Function to determine peaks in an array of data using scipy find_peaks_cwt, works fast in certain cases, but for match_filter cccsum peak finding, find_peaks2_short works better. Test it out and see which works best for your application.

Parameters:
  • arr (ndarray) – 1-D numpy array is required
  • thresh (float) – The threshold below which will be considered noise and peaks will not be found in.
  • trig_int (int) – The minimum difference in samples between triggers, if multiple peaks within this window this code will find the highest.
  • debug (int) – Optional, debug level 0-5
  • maxwidth (int) – Maximum peak width to look for in samples
Returns:

peaks, locs: Lists of peak values and locations.

findpeaks.find_peaks2_short(arr, thresh, trig_int, debug=0, starttime=UTCDateTime(1970, 1, 1, 0, 0), samp_rate=1.0)[source]

Function to determine peaks in an array of data above a certain threshold. Uses a mask to remove data below threshold and finds peaks in what is left.

Parameters:
  • arr (ndarray) – 1-D numpy array is required
  • thresh (float) – The threshold below which will be considered noise and peaks will not be found in.
  • trig_int (int) – The minimum difference in samples between triggers, if multiple peaks within this window this code will find the highest.
  • debug (int) – Optional, debug level 0-5
Returns:

peaks, locs: Lists of peak values and locations.

findpeaks.find_peaks_dep(arr, thresh, trig_int, debug=0, starttime=UTCDateTime(1970, 1, 1, 0, 0), samp_rate=1.0)[source]

Function to determine peaks in an array of data above a certain threshold.

Depreciated peak-finding routine, very slow, but accurate. If all else fails this one should work.

Parameters:
  • arr (ndarray) – 1-D numpy array is required
  • thresh (float) – The threshold below which will be considered noise and peaks will not be found in.
  • trig_int (int) – The minimum difference in samples between triggers, if multiple peaks within this window this code will find the highest.
Returns:

peaks, locs: Lists of peak values and locations.

findpeaks.is_prime(number)[source]
Function to test primality of a number. Function lifted from online resource:
http://www.codeproject.com/Articles/691200/Primality-test-algorithms-Prime-test-The-fastest-w
This function is distributed under a seperate licence:
This article, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)
Parameters:number (int) – Integer to test for primality
Returns:bool

clustering

Code to compute the linkage between seismograms and cluster them accordingly

Written by Calum Chamberlain, in alpha stages of development as of 24/06/2015

Implimented to streamline templates after template detection in beamforming methods, employed by implimentation of Frank et al. code.

As such this code is designed to work only for templates with the same channels

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

clustering.SVD(templates)[source]

Function to compute the SVD of a number of templates and return the singular vectors and singular values of the templates.

Parameters:templates (List of Obspy.Stream) – List of the templates to be analysed
Returns:SVector(list of ndarray), SValues(list) for each channel, Uvalues(list of ndarray) for each channel, stachans, List of String (station.channel)

Note

It is recommended that you align the data before computing the SVD, e.g., the P-arrival on all templates for the same channel should appear at the same time in the trace.

clustering.SVD_2_stream(SVectors, stachans, k, sampling_rate)[source]

Function to convert the singular vectors output by SVD to streams, one for each singular vector level, for all channels.

Parameters:
  • SVectors (List of np.ndarray) – Singular vectors
  • stachans (List of Strings) – List of station.channel Strings
  • k (int) – Number of streams to return = number of SV’s to include
  • sampling_rate (float) – Sampling rate in Hz
Returns:

SVstreams, List of Obspy.Stream, with SVStreams[0] being composed of the highest rank singular vectors.

clustering.cluster(templates, show=True)[source]

Function to take a set of templates and cluster them, will return distance matrix of grouped templates

Parameters:
  • templates (List of Obspy.Stream) – List of templates to compute clustering for
  • show (bool) – plot linkage on screen if True, defaults to True
Returns:

List of cluster groups, array of length len(templates), with each number relating to a cluster

Note: Not fully feautured yet, returns the Z matrix, but doesn’t tell you what can be clustered.

clustering.corr_cluster(traces, thresh=0.9)[source]

Group traces based on correlations above threshold with the stack - will run twice, once with a lower threshold, then again with your threshold to remove large outliers

Parameters:
  • traces (List of :class:obspy.Trace) – Traces to compute similarity between
  • thrsh – Correlation threshold between -1-1
Returns:

np.ndarray of bool

clustering.cross_chan_coherence(st1, st2)[source]

Function to determine the cross-channel coherancy between two streams of multichannel seismic data.

Parameters:
  • st1 (obspy Stream) – Stream one
  • st2 (obspy Stream) – Stream two
Returns:

cross channel coherence, float - normalized by number of channels

clustering.distance_matrix(templates)[source]

Function to compute the distance matrix for all templates - will give distance as 1-abs(cccoh), e.g. a well correlated pair of templates will have small distances, and an equally well correlated reverse image will have the same distance as apositively correlated image - this is an issue

Parameters:template – List of the streams to compute the distance matrix for
Returns:ndarray - distance matrix
clustering.empirical_SVD(templates, linear=True)[source]

Empirical subspace detector generation function. Takes a list of templates and computes the stack as the first order subspace detector, and the differential of this as the second order subspace detector following the emprical subspace method of Barrett & Beroza, 2014 - SRL.

Parameters:
  • templates (list of stream) – list of template streams to compute the subspace detectors from
  • linear (Bool) – Set to true by default to compute the linear stack as the first subspace vector, False will use the phase-weighted stack as the first subspace vector.
Returns:

list of two streams

clustering.extract_detections(detections, templates, extract_len=90.0, outdir=None, extract_Z=True, additional_stations=[])[source]

Function to extract the waveforms associated with each detection in a list of detections for the template, template. Waveforms will be returned as a list of obspy.Streams containing segments of extract_len. They will also be saved if outdir is set. The default is unset. The default extract_len is 90 seconds per channel.

Parameters:
  • detections (List tuple of of :class: datetime.datetime, string) – List of datetime objects, and their associated template name
  • templates (List of tuple of string and :class: obspy.Stream) – A list of the tuples of the template name and the template Stream used to detect detections.
  • extract_len (float) – Length to extract around the detection (will be equally cut around the detection time) in seconds. Default is 90.0.
  • outdir (Bool or String) – Default is None, with None set, no files will be saved, if set each detection will be saved into this directory with files named according to the detection time, NOT than the waveform start time. Detections will be saved into template subdirectories.
  • extract_Z (Bool) – Set to True to also extract Z channels for detections delays will be the same as horizontal channels, only applies if only horizontal channels were used in the template.
  • additional_stations (List of tuple) – List of stations, chanels and networks to also extract data for using an average delay.
Returns:

List of :class: obspy.Stream

clustering.group_delays(templates)[source]

Function to group template waveforms according to their delays

Parameters:templates (List of obspy.Stream) – List of the waveforms you want to group
Returns:List of List of obspy.Streams where each initial list is a group with the same delays

pre_processing

Utilities module for the EQcorrscan package written by Calum Chamberlain of Victoria University Wlelington. These functions are designed to do the basic processing of the data using obspy modules (which also rely on scipy and numpy).

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

pre_processing._check_daylong(tr)[source]

Function to check the data quality of the daylong file - check to see that the day isn’t just zeros, with large steps, if it is then the resampling will hate it.

Parameters:tr (obspy.Trace) – Trace to check if the data are daylong.
Return qual:bool
pre_processing.dayproc(tr, lowcut, highcut, filt_order, samp_rate, debug, starttime)[source]

Basic function to bandpass, downsample and check headers and length of trace to ensure files start at the start of a day and are daylong. Works in place on data. This is employed to ensure all parts of the data are processed in the same way.

Parameters:
  • tr (obspy.Trace) – Trace to process
  • highcut (float) – High cut in Hz for bandpass
  • filt_order (int) – Corners for bandpass
  • samp_rate (float) – Desired sampling rate in Hz
  • debug (int) – Debug output level from 0-5, higher numbers = more output
  • starttime (obspy.UTCDateTime) – Desired start of trace
Returns:

obspy.Stream

pre_processing.shortproc(st, lowcut, highcut, filt_order, samp_rate, debug=0)[source]

Basic function to bandpass, downsample. Works in place on data. This is employed to ensure all parts of the data are processed in the same way.

Parameters:
  • st (obspy.Stream) – Stream to process
  • highcut (float) – High cut for bandpass in Hz
  • lowcut (float) – Low cut for bandpass in Hz
  • filt_order (int) – Number of corners for bandpass filter
  • samp_rate (float) – Sampling rate desired in Hz
  • debug (int) – Debug flag from 0-5, higher numbers = more output
Returns:

obspy.Stream

EQcorrscan_plotting

Utility code for most of the plots used as part of the EQcorrscan package.

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

EQcorrscan_plotting.NR_plot(stream, NR_stream, detections, false_detections=False, size=(18.5, 10), save=False, title=False)[source]

Function to plot the Network response alongside the streams used - highlights detection times in the network response

Parameters:
  • stream – Stream to plot
  • NR_stream – Stream for the network response
  • detections (List of datetime objects) – List of the detections
  • false_detections (List of datetime) – Either False (default) or list of false detection times
  • size (tuple) – Size of figure, default is (18.5,10)
  • save (bool) – Save figure or plot to screen, if not False, must be string of save path
  • title (str) – String for the title of the plot, set to False
EQcorrscan_plotting.Noise_plotting(station, channel, PAZ, datasource)[source]

Function to make use of obspy’s PPSD functionality to read in data from a single station and the poles-and-zeros for that station before plotting the PPSD for this station. See McNamara(2004) for more details.

Parameters:
  • station (String) – Station name as it is in the filenames in the database
  • channel (String) – Channel name as it is in the filenames in the database
  • PAZ (Dict) – Must contain, Poles, Zeros, Sensitivity, Gain :type Poles: List of Complex :type Zeros: List of Complex :type Sensitivity: Float :type Gain: Float
  • datasource (String) – The directory in which data can be found, can contain wildcards.
Returns:

PPSD object

EQcorrscan_plotting.SVD_plot(SVStreams, SValues, stachans, title=False)[source]

Function to plot the singular vectors from the clustering routines, one plot for each stachan

Parameters:
  • SVStreams (List of :class:Obspy.Stream) – See clustering.SVD_2_Stream - will assume these are ordered by power, e.g. first singular vector in the first stream
  • SValues (List of float) – List of the singular values corresponding to the SVStreams
  • stachans (List) – List of station.channel
EQcorrscan_plotting.cumulative_detections(dates, template_names, save=False, savefile='')[source]

Simple plotting function to take a list of datetime objects and plot a cumulative detections list. Can take dates as a list of lists and will plot each list seperately, e.g. if you have dates from more than one template it will overlay them in different colours.

Parameters:
  • dates (list of lists of datetime.datetime) – Must be a list of lists of datetime.datetime objects
  • template_names (list of strings) – List of the template names in order of the dates
  • save (Boolean, optional) – Save figure or show to screen
  • savefile (String, optional) – String to save to.
EQcorrscan_plotting.detection_multiplot(stream, template, times, streamcolour='k', templatecolour='r')[source]

Function to plot the stream of data that has been detected in, with the template on top of it timed according to a list of given times, just a pretty way to show a detection!

Parameters:
  • stream (obspy.Stream) – Stream of data to be plotted as the base (black)
  • template (obspy.Stream) – Template to be plotted on top of the base stream (red)
  • times (List of datetime.datetime) – list of times of detections in the order of the channels in template.
  • streamcolour (str) – String of matplotlib colour types for the stream
  • templatecolour (str) – Colour to plot the template in.
EQcorrscan_plotting.interev_mag(times, mags)[source]

Function to plot interevent times against magnitude for given times and magnitudes.

Parameters:
  • times (list of datetime) – list of the detection times, must be sorted the same as mags
  • mags (list of float) – list of magnitudes
EQcorrscan_plotting.interev_mag_sfiles(sfiles)[source]

Function to plot interevent-time versus magnitude for series of events.

Parameters:sfiles (List) – List of sfiles to read from
EQcorrscan_plotting.multi_event_singlechan(streams, picks, clip=10.0, pre_pick=2.0, freqmin=False, freqmax=False, realign=False, cut=(-3.0, 5.0), PWS=False, title=False)[source]

Function to plot data from a single channel at a single station for multiple events - data will be alligned by their pick-time given in the picks

Parameters:
  • streams (List of :class:obspy.stream) – List of the streams to use, can contain more traces than you plan on plotting
  • picks (List of :class:PICK) – List of picks, one for each stream
  • clip (float) – Length in seconds to plot, defaults to 10.0
  • pre_pick (Float) – Length in seconds to extract and plot before the pick, defaults to 2.0
  • freqmin (float) – Low cut for bandpass in Hz
  • freqmax (float) – High cut for bandpass in Hz
  • realign (Bool) – To compute best alignement based on correlation or not.
  • cut (tuple:) – tuple of start and end times for cut in seconds from the pick
  • PWS (bool) – compute Phase Weighted Stack, if False, will compute linear stack
  • title (str) – Plot title.
Returns:

Alligned and cut traces, and new picks

EQcorrscan_plotting.peaks_plot(data, starttime, samp_rate, save=False, peaks=[(0, 0)], savefile='')[source]

Simple utility code to plot the correlation peaks to check that the peak finding routine is running correctly, used in debugging for the EQcorrscan module.

Parameters:
  • data (numpy.array) – Numpy array of the data within which peaks have been found
  • starttime (obspy.UTCDateTime) – Start time for the data
  • samp_rate (float) – Sampling rate of data in Hz
  • save (Boolean, optional) – Save figure or plot to screen (False)
  • peaks (List of Tuple, optional) – List of peak locations and amplitudes (loc, amp)
  • savefile (String, optional) – Path to save to, only used if save=True
EQcorrscan_plotting.pretty_template_plot(template, size=(18.5, 10.5), save=False, title=False, background=False)[source]

Function to make a pretty plot of a single template, designed to work better than the default obspy plotting routine for short data lengths.

Parameters:
  • template – Template stream to plot
  • size (tuple) – tuple of plot size
  • save (Boolean) – if False will plot to screen, if True will save
  • title (Boolean) – String if set will be the plot title
  • background – Stream to plot the template within.
EQcorrscan_plotting.threeD_gridplot(nodes, save=False, savefile='')[source]

Function to plot in 3D a series of grid points.

Parameters:
  • nodes (List of tuples) – List of tuples of the form (lat, long, depth)
  • save (bool) – if True will save without plotting to screen, if False (default) will plot to screen but not save
  • savefile (str) – required if save=True, path to save figure to.
EQcorrscan_plotting.threeD_seismplot(stations, nodes)[source]

Function to plot seismicity and stations in a 3D, movable, zoomable space using matplotlibs Axes3D package.

Parameters:
  • stations (list of tuple) – list of one tuple per station of (lat, long, elevation), with up positive
  • nodes (list of tuple) – list of one tuple per event of (lat, long, depth) with down positive
EQcorrscan_plotting.triple_plot(cccsum, trace, threshold, save=False, savefile='')[source]

Main function to make a triple plot with a day-long seismogram, day-long correlation sum trace and histogram of the correlation sum to show normality

Parameters:
  • cccsum (numpy.array) – Array of the cross-channel cross-correlation sum
  • trace (obspy.Trace) – A sample trace from the same time as cccsum
  • threshold (float) – Detection threshold within cccsum
  • save (Bool, optional) – If True will svae and not plot to screen, vice-versa if False
  • savefile (String, optional) – Path to save figure to, only required if save=True

mag_calc

Functions to simulate Wood Anderson traces, pick maximum peak-to-peak amplitudes write these amplitudes and periods to SEISAN s-files and to calculate magnitudes from this and the informaiton within SEISAN s-files.

Written as part of the EQcorrscan package by Calum Chamberlain - first written to impliment magnitudes for the 2015 Wanaka aftershock sequence, written up by Warren-Smith [2014/15].

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

mag_calc.Amp_pick_sfile(sfile, datapath, respdir, chans=['Z'], var_wintype=True, winlen=0.9, pre_pick=0.2, pre_filt=True, lowcut=1.0, highcut=20.0, corners=4)[source]

Function to read information from a SEISAN s-file, load the data and the picks, cut the data for the channels given around the S-window, simulate a Wood Anderson seismometer, then pick the maximum peak-to-trough amplitude.

Output will be put into a mag_calc.out file which will be in full S-file format and can be copied to a REA database.

Parameters:
  • datapath (String) – Path to the waveform files - usually the path to the WAV directory
  • respdir (String) – Path to the response information directory
  • chans (List of strings) – List of the channels to pick on, defaults to [‘Z’] - should just be the orientations, e.g. Z,1,2,N,E
  • var_wintype (Bool) – If True, the winlen will be multiplied by the P-S time if both P and S picks are available, otherwise it will be multiplied by the hypocentral distance*0.34 - dervided using a p-s ratio of 1.68 and S-velocity of 1.5km/s to give a large window, defaults to True
  • winlen (Float) – Length of window, see above parameter, if var_wintype is False Then this will be in seconds, otherwise it is the multiplier to the p-s time, defaults to 0.5
  • pre_pick (Float) – Time before the s-pick to start the cut window, defaults to 0.2
  • pre_filt (Bool) – To apply a pre-filter or not, defaults to True
  • lowcut (Float) – Lowcut in Hz for the pre-filter, defaults to 1.0
  • highcut (Float) – Highcut in Hz for the pre-filter, defaults to 20.0
  • corners (Int) – Number of corners to use in the pre-filter
mag_calc._GSE2_PAZ_read(GSEfile)[source]

Function to read the instrument response information from a GSE Poles and Zeros file as generated by the SEISAN program RESP.

Format must be CAL2, not coded for any other format at the moment, contact the author to add others in.

Parameters:GSEfile (Str) – Path to GSE file
Returns:Dict of poles, zeros, gain and sensitivity
mag_calc._find_resp(station, channel, network, time, delta, directory)[source]

Helper function to find the response information for a given station and channel at a given time and return a dictionary of poles and zeros, gain and sensitivity.

Parameters:
  • station (String) – Station name (as in the response files)
  • channel (String) – Channel name (as in the response files)
  • network (String) – Network to scan for, can be a wildcard
  • time (datetime.datetime) – Date-time to look for repsonse information
  • delta (float) – Sample interval in seconds
  • directory (String) – Directory to scan for response information
Returns:

Dictionary

mag_calc._max_p2t(data, delta)[source]

Function to find the maximum peak to trough amplitude and period of this amplitude. Originally designed to be used to calculate magnitudes (by taking half of the peak-to-trough amplitude as the peak amplitude).

Parameters:
  • data (ndarray) – waveform trace to find the peak-to-trough in.
  • delta (float) – Sampling interval in seconds
Returns:

tuple of (amplitude, period, time) with amplitude in the same scale as given in the input data, and period in seconds, and time in seconds from the start of the data window.

mag_calc._sim_WA(trace, PAZ, seedresp, water_level)[source]

Function to remove the insturment response from a trace and return a de-meaned, de-trended, Wood Anderson simulated trace in it’s place.

Works in-place on data and will destroy your original data, copy the trace before giving it to this function!

Parameters:
  • trace (obspy.Trace) – A standard obspy trace, generally should be given without pre-filtering, if given with pre-filtering for use with amplitude determiniation for magnitudes you will need to worry about how you cope with the response of this filter yourself.
  • PAZ (dict) – Dictionary containing lists of poles and zeros, the gain and the sensitivity.
  • water_level (int) – Water level for the simulation.
Returns:

obspy.Trace

mag_calc.dist_calc(loc1, loc2)[source]

Function to calcualte the distance in km between two points, uses the flat Earth approximation

Parameters:
  • loc1 (Tuple) – Tuple of lat, lon, depth (in decimal degrees and km)
  • loc2 (Tuple) – Tuple of lat, lon, depth (in decimal degrees and km)

stacking

Utility module of the EQcorrscan package to allow for different methods of stacking of seismic signal in one place.

In alpha stages and only with linear stacking implimented thusfar

Calum Chamberlain 24/06/2015

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

stacking.PWS_stack(streams, weight=2)[source]

Function to compute the phase weighted stack of a series of streams. Recommend aligning the traces before stacking.

Parameters:
  • streams (list of obspy.Stream) – List of Stream to stack
  • weight (float) – Exponent to the phase stack used for weighting.
Returns:

obspy.Stream

stacking.align_traces(trace_list, shift_len)[source]

Function to allign traces relative to each other based on their cross-correlation value

Parameters:
  • trace_list (List of Traces) – List of traces to allign
  • shift_len (int) – Length to allow shifting within in samples
Returns:

list of shifts for best allignment in seconds

stacking.linstack(streams)[source]

Function to compute the linear stack of a series of seismic streams of multiplexed data

Parameters:stream – List of streams to stack
Returns:stack - Stream

catalogue2DD

Module written by Calum Chamberlain as part of the EQcorrscan package.

This module contains functions to convert a seisan catalogue to files ready for relocation in hypoDD - it will generate both a catalogue (dt.ct) file, event file (event.dat), station information file (station.dat), and a correlation oiutput file correlated every event in the catalogue with every other event to optimize the picks (dt.cc).

The correlation routine relies on obspy’s xcorrPickCorrection function from the obspy.signal.cross_correlation module. This function optimizes picks to better than sample accuracy by interpolating the correlation function and finding the maximum of this rather than the true maximum correlation value. The output from this function is stored in the dt.cc file.

Information for the station.dat file is read from SEISAN’s STATION0.HYP file

Earthquake picks and locations are taken from the catalogued s-files - these must be pre-located before entering this routine as origin times and hypocentre locations are needed for event.dat files.

Copyright 2015 Calum Chamberlain

This file is part of EQcorrscan.

EQcorrscan is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

EQcorrscan is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with EQcorrscan. If not, see <http://www.gnu.org/licenses/>.

catalogue2DD._av_weight(W1, W2)[source]

Function to convert from two seisan weights (0-4) to one hypoDD weight(0-1)

Parameters:
  • W1 (str) – Seisan input weight (0-4)
  • W2 (str) – Seisan input weight (0-4)
Returns:

str

catalogue2DD._cc_round(num, dp)[source]

Convenience function to take a float and round it to dp padding with zeros to return a string

Parameters:
  • num (float) – Number to round
  • dp (int) – Number of decimal places to round to.
Returns:

string

catalogue2DD._separation(loc1, loc2)[source]

Function to calculate the distance between two points in the earth

Parameters:
  • loc1 (tuple (float, float, float)) – First point location as lat, long, depth in deg, deg, km
  • loc2 (tuple (float, float, float)) – First point location as lat, long, depth in deg, deg, km
Returns:

distance in km (float)

catalogue2DD.readSTATION0(path, stations)[source]

Function to read the STATION0.HYP file on the path given. Outputs written in station.dat file.

Parameters:
  • path (String) – Path to the STATION0.HYP file
  • station (List) – Stations to look for
Returns:

List of tuples of station, lat, long, elevation

catalogue2DD.write_catalogue(event_list, max_sep=1, min_link=8)[source]

Function to write the dt.ct file needed by hypoDD - takes input event list from write_event as a list of tuples of event id and sfile. It will read the pick information from the seisan formated s-file using the Sfile_util utilities.

Parameters:
  • event_list (List of tuple) – List of tuples of event_id (int) and sfile (String)
  • max_sep (float) – Maximum seperation between event pairs in km
  • min_link (int) – Minimum links for an event to be paired
Returns:

List stations

catalogue2DD.write_correlations(event_list, wavbase, extract_len, pre_pick, shift_len, lowcut=1.0, highcut=10.0, max_sep=4, min_link=8, coh_thresh=0.0)[source]

Function to write a dt.cc file for hypoDD input - takes an input list of events and computes pick refienements by correlation.

Note that this is NOT fast.

Parameters:
  • event_list (List of tuple) – List of tuples of event_id (int) and sfile (String)
  • wavbase (string) – Path to the seisan wave directory that the wavefiles in the S-files are stored
  • extract_len (float) – Length in seconds to extract around the pick
  • pre_pick (float) – Time before the pick to start the correclation window
  • shift_len (float) – Time to allow pick to vary
  • lowcut (float) – Lowcut in Hz - default=1.0
  • highcut (float) – Highcut in Hz - deafult=10.0
  • max_sep (float) – Maximum seperation between event pairs in km
  • min_link (int) – Minimum links for an event to be paired
catalogue2DD.write_event(sfile_list)[source]

Function to write out an event.dat file of the events

Parameters:sfile_list (List) – List of s-files to sort and put into the database
Returns:List of tuples of event ID (int) and Sfile name