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.
Note that these functions do not provide full functionality between quakeML and seisan s-files. Currently (as of version 0.1.1) these only convert pick times and phase information, along with amplitude information for local magnitudes between seisan and quakeML. Location information including hypocentre, origin time and magnitudes are also handled.
A series of wrappers and conversions is included between the legacy PICK and EVENTINFO classes, however these will be depreciated along with these classes for version 0.1.0. Users should transition to using obspy.core.event classes as these have more support and functionality.
We have not implimented any handling of focal mechanism solutions between the two formats.
Code written by Calum John Chamberlain and Chet Hopp both of Victoria University of Wellington, 2015 & 2016.
Copyright 2015, 2016 the authors.
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=<mock.Mock object at 0x3505850>, loc_mod_ind=u' ', dist_ind=u' ', ev_id=u' ', latitude=nan, longitude=nan, depth=nan, depth_ind=u' ', loc_ind=u' ', agency=u' ', nsta=0, t_RMS=nan, Mag_1=nan, Mag_1_type=u' ', Mag_1_agency=u' ', Mag_2=nan, Mag_2_type=u' ', Mag_2_agency=u' ', Mag_3=nan, Mag_3_type=u' ', Mag_3_agency=u' ')[source]¶ Header information for seisan events, again all fields can be left blank for a default empty header. The print function for header will print important information, but not as seen in an S-file.
For more information on parameters see the seisan manual.
- Attributes:
type time: obspy.UTCDateTime param time: Event origin time type loc_mod_ind: str param loc_mod_ind: type dist_ind: str param dist_ind: Distance flag, usually ‘L’ for local, ‘R’ for regional and ‘D’ for distant type ev_id: str param ev_id: Often blank, ‘E’ denotes explosion and fixes depth to 0km type latitude: float param latitude: Hypocentre latitude in decimal degrees type longitude: float param lognitude: Hypocentre longitude in decimal degrees type depth: float param depth: hypocentre depth in km type depth_ind: str param depth_ind: type loc_ind: str param loc_ind: type agency: str param agency: Reporting agency, three letters type nsta: int param nsta: Number of stations recording type t_RMS: float param t_RMS: Root-mean-squared time residual type Mag_1: float param Mag_1: first magnitude type Mag_1_type: str param Mag_1_type: Type of magnitude for Mag_1 (‘L’, ‘C’, ‘W’) type Mag_1_agency: str param Mag_1_agency: Reporting agency for Mag_1 type Mag_2: float param Mag_2: second magnitude type Mag_2_type: str param Mag_2_type: Type of magnitude for Mag_2 (‘L’, ‘C’, ‘W’) type Mag_2_agency: str param Mag_2_agency: Reporting agency for Mag_2 type Mag_3: float param Mag_3: third magnitude type Mag_3_type: str param Mag_3_type: Type of magnitude for Mag_3 (‘L’, ‘C’, ‘W’) type Mag_3_agency: str param Mag_3_agency: Reporting agency for Mag_3
Note: Depreciated legacy function, use the obspy.core.event classes. This will be removed in future releases.
-
class
sfile_util.
PICK
(station=u' ', channel=u' ', impulsivity=u' ', phase=u' ', weight=999, polarity=u' ', time=<mock.Mock object at 0x3505850>, coda=999, amplitude=nan, peri=nan, azimuth=nan, velocity=nan, AIN=999, SNR=nan, azimuthres=999, timeres=nan, finalweight=999, distance=nan, CAZ=999)[source]¶ Pick information for seisan implimentation, note all fields can be left blank to obtain a default pick: picks have a print function which will print them as they would be seen in an S-file.
- Attributes:
type station: str param station: Station name, less than five charectars required as standard type channel: str param channel: Two or three charactar channel name, stored as two charactars in S-file type impulsivity: str param impulsivity: either ‘C’ or ‘D’ for compressive and dilatational type phase: str param phase: Any allowable phase name in two characters type weight: int param weight: 0-4 with 0=100%, 4=0%, use weight=9 for unknown timing type polarity: str type time: obspy.UTCDateTime() param time: Pick time as an obspy.UTCDateTime object type coda: int param coda: Length of coda in seconds type amplitude: float param amplitude: Amplitude (zero-peak), type is given in phase type peri: float param peri: Period of amplitude type azimuth: float param azimuth: Direction of approach in degrees type velocity: float param velocity: Phase velocity (km/s) type AIN: int param AIN: Angle of incidence. type SNR: float param SNR: Signal to noise ratio type azimuthres: int param azimuthres: Residual azimuth type timeres: float param timeres: Time residual in seconds type finalweight: int param finalweight: Final weight used in location type distance: float param distance: Source-reciever distance in km type CAZ: int param CAZ: Azimuth at source.
Note: Depreciated legacy function, use the obspy.core.event classes. This will be removed in future releases.
-
sfile_util.
_evmagtonor
(mag_type)[source]¶ Convenience tool to switch from obspy event magnitude types to seisan syntax
-
sfile_util.
_float_conv
(string)[source]¶ Convenience tool to convert from string to float, if empty string return NaN rather than an error
-
sfile_util.
_int_conv
(string)[source]¶ Convenience tool to convert from string to integer, if empty string return a 999 rather than an error.
-
sfile_util.
_nortoevmag
(mag_type)[source]¶ Convenience tool to switch from nordic type magnitude notation to obspy event magnitudes.
-
sfile_util.
_str_conv
(number, rounded=False)[source]¶ Convenience tool to convert a number, either float or into into a string, if the int is 999, or the float is NaN, returns empty string.
-
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.
eventtopick
(event)[source]¶ Wrapper function to convert from obspy.core.event to legacy PICK and EVENT classes.
Parameters: event (obspy.core.event.Event) – A single obspy event Returns: List of PICK(), and a single EVENTINFO() Note
This is a wrapper to simplify transition from PICK and EVENT classes to obspy.core.event classes. This will not be maintained beyond v 0.1.0.
New in version 0.1.0.
-
sfile_util.
eventtosfile
(event, userID, evtype, outdir, wavefiles, explosion=False, overwrite=False)[source]¶ Function to take an obspy.event and write the relevant information to a nordic formatted s-file
Parameters: - event (obspy.event.core.Event) – A single obspy event
- userID (str) – Up to 4 character user ID
- evtype (str) – Single character string to describe the event, either L, R or D.
- outdir (str) – Path to directory to write to
- wavefiles (list of str) – Waveforms to associate the sfile with
- explosion (bool) – Note if the event is an explosion, will be marked by an E.
- overwrite (bool) – force to overwrite old files, defaults to False
Returns: str: name of sfile written
Note
Seisan can find waveforms either by their relative or absolute path, or by looking for the file recursiuvely in directories within the WAV directory in your seisan install. Because all lines need to be less than 79 charecters long (fortran hangover) in the s-files, you will need to determine whether the full-path is okay or not.
-
sfile_util.
nordpick
(event)[source]¶ Function to print information from an obspy.event class to nordic format.
Parameters: event – A single obspy event. Returns: List of String Note
Currently finalweight is unsupported, nor is velocity, or angle of incidence. This is because obspy.event stores slowness in s/deg and takeoff angle, which would require computation from the values stored in seisan. Multiple weights are also not supported in Obspy.event.
New in version 0.1.0.
-
sfile_util.
picktoevent
(evinfo, picks)[source]¶ Wrapper function to convert from EVENTINFO and PICK classes to obspy.core.event.Event.
Parameters: - evinfo (EVENTINFO) – Event header info for a single event
- picks (List of PICK) – List of picks associated with the event
Returns: obspy.core.event.Event
Note
This is a legacy support function, users should avoid this as it will be removed for version 0.1.1. Written to aid transition from in-built classes to obspy.core.event classes.
New in version 0.1.0.
-
sfile_util.
populatesfile
(sfile, event)[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: - sfile (str) – Path to S-file to populate, must have a header already
- picks – A single event to be written to a single S-file.
-
sfile_util.
readheader
(sfile)[source]¶ Function to read the header information from a seisan nordic format S-file. Returns an obspy.core.event.Catalog type: note this changed for version 0.1.0 from the inbuilt class types.
Parameters: sfile (str) – Path to the s-file Returns: class: obspy.core.event.Event
-
sfile_util.
readpicks
(sfile)[source]¶ Function to read pick information from the s-file and store this in an obspy.event.Catalog type. This was changed for version 0.1.0 from using the inbuilt PICK class.
Parameters: sfile (String) – Path to sfile Returns: obspy.core.event.Event Warning
Currently finalweight is unsupported, nor is velocity, or angle of incidence. This is because obspy.event stores slowness in s/deg and takeoff angle, which would require computation from the values stored in seisan. Multiple weights are also not supported in Obspy.event.
-
sfile_util.
readwavename
(sfile)[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: sfile (str) – Path to the sfile Returns: List of str
-
sfile_util.
stationtoseisan
(station)[source]¶ Convert obspy station inventory to simple string to copy in to STATION0.HYP file for seisan locations.
Parameters: station (obspy.core.inventory.station.Station) – Inventory containing a single station. Returns: str Note
Only works to the low-precision level at the moment (see seisan manual for explanation).