seismic package¶
Subpackages¶
- seismic.ASDFdatabase package
- Submodules
- seismic.ASDFdatabase.ClientUtils module
- seismic.ASDFdatabase.FederatedASDFDataSet module
- seismic.ASDFdatabase.asdf2mseed module
- seismic.ASDFdatabase.asdf_preprocess module
- seismic.ASDFdatabase.create_small_chunks module
- seismic.ASDFdatabase.make_field_notesDB module
- seismic.ASDFdatabase.plot_data_quality module
- seismic.ASDFdatabase.query_input_yes_no module
- seismic.ASDFdatabase.sc3toasdf module
- seismic.ASDFdatabase.seisds module
- seismic.ASDFdatabase.utils module
- Module contents
- seismic.amb_noise package
- seismic.gps_corrections package
- seismic.inventory package
- Subpackages
- Submodules
- seismic.inventory.add_time_corrections module
- seismic.inventory.engd2stxml module
- seismic.inventory.fdsnxml_convert module
- seismic.inventory.generate_network_plots module
- seismic.inventory.inventory_merge module
- seismic.inventory.inventory_split module
- seismic.inventory.inventory_util module
- seismic.inventory.iris_query module
- seismic.inventory.pdconvert module
- seismic.inventory.plotting module
- seismic.inventory.response module
- seismic.inventory.table_format module
- seismic.inventory.update_iris_inventory module
- Module contents
- seismic.inversion package
- Subpackages
- seismic.inversion.mcmc package
- seismic.inversion.wavefield_decomp package
- Submodules
- seismic.inversion.wavefield_decomp.call_count_decorator module
- seismic.inversion.wavefield_decomp.plot_nd_batch module
- seismic.inversion.wavefield_decomp.runners module
- seismic.inversion.wavefield_decomp.solvers module
- seismic.inversion.wavefield_decomp.wavefield_continuation_tao module
- seismic.inversion.wavefield_decomp.wfd_plot module
- Module contents
- Module contents
- Subpackages
- seismic.pick_harvester package
- seismic.receiver_fn package
- Submodules
- seismic.receiver_fn.bulk_rf_report module
- seismic.receiver_fn.generate_rf module
- seismic.receiver_fn.plot_ccp module
- seismic.receiver_fn.plot_ccp_batch module
- seismic.receiver_fn.plot_spatial_map module
- seismic.receiver_fn.pointsets2grid module
- seismic.receiver_fn.rf_3dmigrate module
- seismic.receiver_fn.rf_deconvolution module
- seismic.receiver_fn.rf_h5_file_event_iterator module
- seismic.receiver_fn.rf_h5_file_station_iterator module
- seismic.receiver_fn.rf_handpick_tool module
- seismic.receiver_fn.rf_inversion_export module
- seismic.receiver_fn.rf_network_dict module
- seismic.receiver_fn.rf_plot_utils module
- seismic.receiver_fn.rf_plot_vespagram module
- seismic.receiver_fn.rf_process_io module
- seismic.receiver_fn.rf_quality_filter module
- seismic.receiver_fn.rf_stacking module
- seismic.receiver_fn.rf_synthetic module
- seismic.receiver_fn.rf_util module
- Module contents
- seismic.synthetics package
- seismic.traveltime package
- Submodules
- seismic.traveltime.cluster_grid module
- seismic.traveltime.events_stations_rays_visualization2 module
- seismic.traveltime.gather_events module
- seismic.traveltime.mpiops module
- seismic.traveltime.parametrisation module
- seismic.traveltime.parse_param module
- seismic.traveltime.plotviews module
- seismic.traveltime.pslog module
- seismic.traveltime.sort_rays module
- seismic.traveltime.zone_rays module
- Module contents
- seismic.xcorqc package
- Submodules
- seismic.xcorqc.analytic_plot_utils module
- seismic.xcorqc.client_data module
- seismic.xcorqc.correlator module
- seismic.xcorqc.fft module
- seismic.xcorqc.generate_dispersion_curves module
- seismic.xcorqc.generate_test_data module
- seismic.xcorqc.utils module
- seismic.xcorqc.validate_xcorr_setup module
- seismic.xcorqc.xcorqc module
- seismic.xcorqc.xcorr_station_clock_analysis module
- Module contents
Submodules¶
seismic.analyze_station_orientations module¶
Analyze a data set of seismic arrival events on a per-station basis and try to detect and estimate any station orientation error.
The event waveform dataset provided as input to this script (or to NetworkEventDataset if calling directly into function analyze_station_orientations) is generated by the script seismic/extract_event_traces.py.
In future, consider moving this script to the inventory module and applying corrections to the station inventory xml (to the azimuth tag).
Reference: Wilde-Piórko, M., Grycuk, M., Polkowski, M. et al. On the rotation of teleseismic seismograms based on the receiver function technique. J Seismol 21, 857-868 (2017). https://doi.org/10.1007/s10950-017-9640-x
-
seismic.analyze_station_orientations.
analyze_station_orientations
(ned, curation_opts, config_filtering, config_processing, save_plots_path=None)[source]¶ Main processing function for analyzing station orientation using 3-channel event waveforms. Uses method of Wilde-Piorko https://doi.org/10.1007/s10950-017-9640-x
One should not worry about estimates that come back with error of less than about 20 degrees from zero, since this analysis provides only an estimate.
- Parameters
ned – NetworkEventDataset containing waveforms to analyze. Note: the data in this dataset will be modified by this function.
curation_opts – Seismogram curation options. Safe default to use is DEFAULT_CURATION_OPTS.
config_filtering – Seismogram filtering options for RF computation. Safe default to use is DEFAULT_CONFIG_FILTERING.
config_processing – Seismogram RF processing options. Safe default to use is DEFAULT_CONFIG_PROCESSING.
save_plots_path – Optional folder in which to save plot per station of mean arrival RF amplitude as function of correction angle
- Returns
Dict of estimated orientation error with net.sta code as the key.
-
seismic.analyze_station_orientations.
process_event_file
(src_h5_event_file, curation_opts=None, config_filtering=None, config_processing=None, dest_file=None, save_plots_path=None)[source]¶ Use event dataset from an HDF5 file to analyze station for orientation errors.
- Parameters
src_h5_event_file – HDF5 file to load. Typically one created by extract_event_traces.py script
curation_opts – Seismogram curation options. Safe default to use is DEFAULT_CURATION_OPTS.
config_filtering – Seismogram filtering options for RF computation. Safe default to use is DEFAULT_CONFIG_FILTERING.
config_processing – Seismogram RF processing options. Safe default to use is DEFAULT_CONFIG_PROCESSING.
dest_file – File in which to save results in JSON format
save_plots_path – Optional folder in which to save plot per station of mean arrival RF amplitude as function of correction angle
- Returns
None
seismic.extract_event_traces module¶
Use waveform database and station inventory to extract raw traces for all seismic events within a given magnitude and time range.
-
seismic.extract_event_traces.
asdf_get_waveforms
(asdf_dataset, network, station, location, channel, starttime, endtime)[source]¶ Custom waveform getter function to retrieve waveforms from FederatedASDFDataSet.
- Parameters
asdf_dataset (seismic.ASDFdatabase.FederatedASDFDataSet) – Instance of FederatedASDFDataSet to query
network (str) – Network code
station (str) – Station code
location (str) – Location code
channel (str) – Channel code
starttime (str in UTC datetime format) – Start time of the waveform query
endtime (str in UTC datetime format) – End time of the waveform query
- Returns
Stream containing channel traces
- Return type
obspy.Stream of obspy.Traces
-
seismic.extract_event_traces.
get_events
(lonlat, starttime, endtime, cat_file, distance_range, magnitude_range, early_exit=True)[source]¶ Load event catalog (if available) or create event catalog from FDSN server.
- Parameters
lonlat (tuple(float, float)) – (Longitude, latitude) of reference location for finding events
starttime (obspy.UTCDateTime or str in UTC datetime format) – Start time of period in which to query events
endtime (obspy.UTCDateTime or str in UTC datetime format) – End time of period in which to query events
cat_file (str or Path) – File containing event catalog, or file name in which to store event catalog
distance_range (tuple(float, float)) – Range of distances over which to query seismic events
magnitude_range (tuple(float, float)) – Range of event magnitudes over which to query seismic events.
early_exit (bool, optional) – If True, exit as soon as new catalog has been generated, defaults to True
- Returns
Event catalog
- Return type
-
seismic.extract_event_traces.
is_url
(resource_path)[source]¶ Convenience function to check if a given resource path is a valid URL
-
seismic.extract_event_traces.
timestamp_filename
(fname, t0, t1)[source]¶ - Append pair of timestamps (start and end time) to file name in format that is
compatible with filesystem file naming.
- Parameters
fname (str or path) – File name
t0 (obspy.UTCDateTime) – first timestamp
t1 (obspy.UTCDateTime) – second timestamp
seismic.model_properties module¶
seismic.network_event_dataset module¶
Class encapsulating a collection of event waveforms for stations of one network.
-
class
seismic.network_event_dataset.
NetworkEventDataset
(stream_src, network=None, station=None, location='', ordering='ZNE')[source]¶ Bases:
object
Collection of 3-channel ZNE streams with traces aligned to a fixed time window about seismic P-wave arrival events, for a given network.
Two indexes are provided. One indexes hierarchically by station code and event ID, yielding a 3-channel ZNE stream per event, so that you can easily gather all traces for a given station by iterating over events.
The other index indexes hierarchically by event ID and station code, yielding a 3-channel ZNE stream per station. Using this index you can easily gather all traces for a given event across multiple stations.
Preferably each input trace will already have an ‘event_id’ attribute in its stats. If not, an event ID will be invented based on station identifiers and time window.
-
apply
(_callable)[source]¶ Apply a callable across all streams. Use to apply uniform processing steps to the whole dataset.
- Parameters
_callable (Any callable compatible with the call signature.) – Callable object that takes an obspy Stream as input and applies itself to that Stream. Expect that stream may be mutated in-place by the callable.
- Returns
None
-
by_event
()[source]¶ Iterate over event sub-dictionaries
- Returns
Iterable over the discrete events, each element consisting of pair containing (event id, station dict).
-
by_station
()[source]¶ Iterate over station sub-dictionaries
- Returns
Iterable over the stations, each element consisting of pair containing (station code, event dict).
-
curate
(curator)[source]¶ Curate the dataset according to a callable curator. Curator call signature takes station code, event id and stream as input, and returns boolean whether to keep Stream or not.
- Parameters
curator (Callable) – Function or callable delegate to adjudicate whether to keep each given stream.
- Returns
None
-
event
(event_id)[source]¶ Accessor for stations for a given event. :param event_id: ID of event to look up :return: Station index for given event, if event ID is found :rtype: SortedDict
-
prune
(items, cull=True)[source]¶ Remove a given sequence of (station, event) pairs from the dataset.
- Parameters
items – Iterable of (station, event) pairs
cull – If True, then empty entries in the top level index will be removed.
- Returns
None
-
seismic.plot_network_event_dataset module¶
Bulk plotting helper functions based on NetworkEventDataset
-
seismic.plot_network_event_dataset.
plot_ned_seismograms
(ned, output_file, channel_order='ZNE')[source]¶ Plot seismograms in NetworkEventDataset to PDF file. If dataset is very large, this may take a long time to run.
- Parameters
ned – NetworkEventDataset containing waveforms.
output_file – Output file name
seismic.stream_io module¶
Helper functions for seismic stream IO.
-
seismic.stream_io.
iter_h5_stream
(src_file, headonly=False)[source]¶ Iterate over hdf5 file containing streams in obspyh5 format.
- Parameters
src_file – str or path to file to read
headonly – Only read trace stats, do not read actual time series data
- Yield
obspy.Stream containing traces for a single seismic event.
-
seismic.stream_io.
read_h5_stream
(src_file, network=None, station=None, loc='', root='/waveforms')[source]¶ Helper function to load stream data from hdf5 file saved by obspyh5 HDF5 file IO. Typically the source file is generated using extract_event_traces.py script. For faster loading time, a particular network and station may be specified.
- Parameters
- Returns
All the loaded data in an obspy Stream.
- Return type
obspy.Stream
-
seismic.stream_io.
sac2hdf5
(src_folder, basenames, channels, dest_h5_file, tt_model_id='iasp91')[source]¶ Convert collection of SAC files from a folder into a single HDF5 stream file.
- Parameters
src_folder (str or Path) – Path to folder containing SAC files
basenames (list of str) – List of base filenames (file name excluding extension) to load.
channels (List of str) – List of channels to load. For each base filename in basenames, there is expected to be a file with each channel as the filename extension.
dest_h5_file (str or Path) – Path to output file. Will be created, or overwritten if already exists.
tt_model_id (str) – Which travel time earth model to use for synthesizing trace metadata. Must be known to obspy.taup.TauPyModel
- Returns
None
-
seismic.stream_io.
write_h5_event_stream
(dest_h5_file, stream, mode='a', ignore=())[source]¶ Write stream to HDF5 file in event indexed format using obspy.
- Parameters
dest_h5_file (str or pathlib.Path) – File in which to write the stream.
stream (obspy.Stream) – The stream to write
mode (str) – Write mode, such as ‘w’ or ‘a’. Use ‘a’ to iteratively write multiple streams to one file.
ignore (Any iterable of str) – List of headers to ignore when writing attributes to group. Passed on directly to obspyh5.writeh5
seismic.stream_processing module¶
Utility stream processing functions.
-
seismic.stream_processing.
assert_homogenous_stream
(stream, funcname)[source]¶ Verify that the given stream does not contain mixture of stations or channels/components.
- Parameters
stream (obspy.core.stream.Stream or rf.RFStream) – Stream containing one or more traces
-
seismic.stream_processing.
back_azimuth_filter
(baz, baz_range)[source]¶ Check if back azimuth baz is within range. Inputs must be in the range [0, 360] degrees.
-
seismic.stream_processing.
correct_back_azimuth
(_event_id, stream, baz_correction)[source]¶ Apply modification to the back azimuth value in the stream stats
- Parameters
_event_id – Ignored
stream – Stream to which correction is applied
baz_correction – Any object with a registered scalarize function for generating an angle correction for a trace in degrees. E.g. could be a numeric value, a dictionary of correction values, or a file produced by script analyze_station_orientations.py
- Returns
Stream with modified back azimuth
-
seismic.stream_processing.
negate_channel
(_event_id, stream, channel)[source]¶ Negate the data in the given channel of the stream
- Parameters
_event_id – Ignored
stream – Stream containing channel to flip
channel – Single character string indicating which component to flip
- Returns
Stream with channel data negated
-
seismic.stream_processing.
scalarize
(_obj, _stats)[source]¶ -
seismic.stream_processing.
scalarize
(val: numbers.Number, _stats) -
seismic.stream_processing.
scalarize
(d: dict, stats) -
seismic.stream_processing.
scalarize
(filename: str, stats)
-
seismic.stream_processing.
sinc_resampling
(t, y, t_new)[source]¶ Resample signal y for known times t onto new times t_new. Sampling rates do not need to match and time windows do not need to overlap. t_new should not have a lower sampling rate than t.
- Parameters
t – numpy.array of times
y – numpy.array of sample values
t_new – numpy.array of new times to interpolate onto
- Returns
numpy.array of new interpolated sample values
-
seismic.stream_processing.
swap_ne_channels
(_event_id, stream)[source]¶ Swap N and E channels on a stream. Changes the input stream.
- Parameters
_event_id – Ignored
stream – Stream whose N and E channels are to be swapped
- Returns
Stream with channel swapping applied
seismic.stream_quality_filter module¶
Helper functions for curating and quality controlling stream objects.
-
seismic.stream_quality_filter.
curate_stream3c
(ev_id, stream3c, logger=None)[source]¶ Apply quality curation criteria to a stream. Modifies the stream in-place if required. Traces in stream must be in ZNE order. Each trace in the stream is expected to have metadata with starttime, endtime, channel and inclination.
The following checks are made. If any of these checks fails, the function returns False: * Inclination value is not NaN * The stream has 3 channels * Each trace in the stream has the same number of samples * None of the traces have any NaN values in the time series data * None of the traces have zero variance
The following cleanups are attempted on the stream: * All 3 traces have the same time range
seismic.units_utils module¶
Constants and utility functions for unit conversion.