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
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: 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: 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
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: 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:
-
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
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
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: 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: 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: 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