-
Notifications
You must be signed in to change notification settings - Fork 10
Noise Toolkit PDF PSD package V.2
The PDF/PSD package provides three highly configurable Python 3 scripts to calculate waveform spectra. This package takes advantage of FDSN Web service client for ObsPy to retrieve necessary waveform data and it also allows users to process waveform data from their local files. This package provides PSD file collections similar to the PQLX package (McNamara and Boaz, 2005, https://www.usgs.gov/software/pqlx-a-software-tool-evaluate-seismic-station-performance).
The scripts included in this package are:
-
ntk_computePSD.py – requests waveforms and response data for given station(s) using the ObsPy FDSN client OR
reads user’s waveform data files in SAC, MSEED, CSS, etc. formats and computes PSDs and populates a file-based
PSD database
NOTE: ntk_computePSD.py first identifies the appropriate FDSN data provider for the requested station using the Fedcatalog service from IRIS (service.iris.edu/irisws/fedcatalog/1/ ) and then requests waveform/response data for that station using ObsPy’s FDSN client
-
ntk_extractPsdHour.py – extracts PSDs for a given channel and bounding parameters from the PSD database
(the output is similar to PQLX’s exPSDhour script, https://pubs.usgs.gov/of/2010/1292/ ) - ntk_binPsdDay.py – bins PSD’s to daily files for a given channel and bounding parameter
- install Python 3.8 or higher on your computer
- additional required Python modules are:
. obspy 1.2.2
. matplotlib 3.3.3
. numpy 6.2
. scipy 1.5.3
- download the package under the installation directory.
- under the root directory of the Noise Toolkit execute the following command to plot PSD of TA.O18A.—.BHZ for 12 pm to 1 pm on 2008-08-14
python bin/ntk_computePSD.py net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T13:00:00 xtype=period plot=1
- review all parameter files to get familiar with the script capabilities
- use examples given below to run/test the scripts and become familiar with them
- copy and edit the parameter files under the param directory to customize parameters
- for more information visit the IRIS DMC Noise Toolkit Data Product web page at:
ds.iris.edu/ds/products/noise-toolkit - examples below may turn on the verbose mode. Once you have configured the script, you can turn it off
- We welcome patches and enhancements to this software. When developing patches, please pay particular attention to ease of use and maintenance and also keep dependencies to a minimum.
- for issues, file a ticket under Issues
- 2020-11-16 V.2.0.0: Python 3, use of Fedcatalog, adoption of CSD changes and adoption of PEP 8 style guide.
- 2019-09-09 Robert Anthony (USGS, Albuquerque Seismological Laboratory):
Using CSD to compute the cross spectral density of two signals - 2017-01-18 V.0.9.5: support for reading data and metadata from files only with no Internet requirement
- 2016-11-01 V.0.9.0: added support for obtaining channel responses from local station XML response files by introducing the following two functions in tsLIB.py:
getResponseInventory – to build a list of response inventories under a given met directory
getResponseFromFile – to find response for a given Network, Station, Location and channel
added respDirectory to common.py parameter file to disable looking for response files on local drives, set this parameter to None. Otherwise, set it to the response directory path - 2016-01-25 V.0.8.2: added polarization parameters to common.par and also made changes to some libraries in support of the polarization package
added user and password parameters to common.par and ntk_computePSD.py in support of restricted data access - 2015-04-30 V.0.8.1: added powerDirectory and imageDirectory parameters to common.par in support of ME package.
- 2015-04-07 V.0.8.0: now produces clear messages for missing parameters in the parameter file. Addresses reported maximum period values. Minor bugs and enhancements
- 2015-02-24 V.0.7.0: introduced two new parameters (performInstrumentCorrection, applyScale) to allow user avoid instrument correction also now user can turn od decon. filter
- 2014-11-24 V.0.6.2: documentation update and expansion of the common.py file to support Microseism Energy package
- 2014-10-16 V.0.6.0: ntk_extractPsdHour.py output file name now includes the x-axis type
- 2014-10-01 V.0.5.0: initial Beta release. Made the smoothing configurable, reorganized parameters for easier maintenance
- 2014-05-20 modified the output format to:
1. compute hourly PSDs and store them under the psdDb directory (this will be similar to PQLX’s database) as basis for script outputs
2. provide the extractPsdHourly.py script to extract PSD’s for a given channel and bounding parameters. The output is similar to PQLX’s exPSDhour script
3. provide the extractPdf.py that will extract PSDs to create (hourly and daily) PDFs for a given channel and bounding parameters.
The output is similar to the current PDF output from the PDF-PSD data product (ds.iris.edu/ds/products/pdf-psd/)
- 2014-03-19: Added data file input option
- 2014-03-07: Initial Alpha release
Use below examples to run/test the scripts and become familiar with them and edit the parameter files under the param directory to change parameters based on your needs.
All parameters are under the param directory. The “shared.py” parameter file contains basic parameters shared by all parameter files. Please review ALL parameter to familiarize yourself with the script capabilities. For more information about this package, visit the IRIS DMC Noise Toolkit Data Product web page at:
http://ds.iris.edu/ds/products/noise-toolkit/A Python 3 script that calculates average power spectral densities for a given station. The script:
* identifies the FDSN data provider for the requested station using the Fedcatalog service
from IRIS (https://service.iris.edu/irisws/fedcatalog/1/)
* requests waveform and response data for the given station(s)/channels(s) using ObsPy’s FDSN client
OR
* reads user-supplied waveform data files in SAC, MSEED, CSS, etc. format from a local disk
Then
* computes PSDs for the waveform data and populates a file-based PSD database
python bin/ntk_computePSD.py to display the usage message OR python bin/ntk_computePSD.py param=FileName client=[FDSN|FILES] net=network sta=station loc=location chan=channel(s) start=YYYY-MM-DDTHH:MM:SS end=YYYY-MM-DDTHH:MM:SS xtype=[period|frequency] plot=[0|1] verbose=[0|1] timing=[0|1] to perform computations where: param [default: computePSD] the configuration file name client [default: FDSN] client to use to make data/metadata requests (FDSN or FILES) net [required] network code sta [required] station code loc [required] location ID chan [default: BH?] channel ID(s); separate multiple channel codes by comma (no space) xtype [period or frequency, default: period] X-axis type (period or frequency for outputs and plots) start [required] start date-time (UTC) of the interval for which PSDs to be computed (format YYYY-MM-DDTHH:MM:SS) end [required] end date-time (UTC) of the interval for which PSDs to be computed (format YYYY-MM-DDTHH:MM:SS) sw_width [default: 0.25] Smoothing window width in octave sw_shift [default: 0.125] Smoothing window shift in fraction of octave plotnm [0 or 1, default 1] plot the New High/Low Noise Models [0|1] plot [0 or 1, default: 0] to run in plot mode set to 1 timing [0 or 1, default: 0] to run in timing mode (set to 1 to output run times for different segments of the script) verbose [0 or 1, default: 0] to run in verbose mode set to 1Use the run in plot (plot=1) or verbose mode (verbose=1_ to tune the parameters before a production run (verbose=0 plot=0):
python bin/ntk_computePSD.py net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T13:00:00 xtype=period plot=1
The following example limits channel to BHZ, processes only one hour (default) and plot the PSD with more smoothing (sw_width=0.5 and sw_shift=0.25):
python bin/ntk_computePSD.py param=computePSD net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T13:00:00 xtype=period sw_width=0.5 sw_shift=0.25 plot=1
With even more smoothing (sw_width=1):
python bin/ntk_computePSD.py param=computePSD net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T13:30:00 xtype=period sw_width=1 plot=1
Your request could extend over multiple days:
python bin/ntk_computePSD.py param=computePSD net=NM sta=SLM loc=DASH chan=BHZ start=2009-03-01T00:00:00 end=2009-03-31T00:00:00 xtype=period plot=0 verbose=0
set
requestClient = "FILES"
to flag the script that the waveform data are coming from file and fileTag = "{IRIS_NTK_PSD}/data/TEST/SAC/W*.SAC"
to tell it which files to look at
python bin/ntk_computePSD.py net=TA sta=W5 loc=DASH start=2014-03-17T04:30:00 end=2014-03-17T05:30:00 type=period client=FILES plot=1
Requesting data from a FDSN data center, other than IRIS, is automatic and depends on the station you are requesting. For example requesting BHZ channel for GR.BFO:
python bin/ntk_computePSD.py param=computePSD net=GR sta=BFO loc=DASH chan=BHZ start=2020-10-01T00:00:00 end=2020-10-01T01:00:00 xtype=period plot=1
ntk_extractPsdHour.py to display the usage message OR python bin/ntk_extractPsdHour.py param=FileName net=network sta=station loc=location chan=channel(s) start=YYYY-MM-DDTHH:MM:SS end=YYYY-MM-DDTHH:MM:SS xtype=[period|frequency] verbose=[0|1] to perform binning: param [default: extractPsdHour] the configuration file name net [required] network code sta [required] station code loc [required] location ID chan [required] channel ID xtype [required] X-axis type (period or frequency for output) start [required] start date-time (UTC) of the interval for which PSDs to be computed (format YYYY-MM-DDTHH:MM:SS) end [required] end date-time (UTC) of the interval for which PSDs to be computed (format YYYY-MM-DDTHH:MM:SS) verbose [0 or 1, default: 0] to run in verbose mode set to 1The ntk_computePSD.py output files
Date Hour X-value (period/frequency) Power with values separated using the “separator” character specified in the parameter file.
Full path to the output data file is provided at the end of the run. The output file name has the form: net.sta.loc.chan.start.end.xtype.txt for example: TA.O18A.--.BHZ.2008-08-14.2008-08-14.period.txt
- usage:
python bin/ntk_extractPsdHour.py
- Assuming that you already have tried the following ntk_compute_PSD.py example:
python bin/ntk_computePSD.py net=TA sta=O18A loc=DASH start=2008-08-14T12:00:00 end=2008-08-14T13:30:00
- Then you can perform PSD extraction via:
python bin/ntk_extractPsdHour.py net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T12:30:00 xtype=period
A Python 3 script to bin PSDs to daily files for a given channel and bounding parameters.
ntk_binPsdDay.py to display the usage message OR ntk_binPsdDay.py param=FileName net=network sta=station loc=location chan=channel(s) start=YYYY-MM-DDTHH:MM:SS end=YYYY-MM-DDTHH:MM:SS xtype=[period|frequency] verbose=[0|1] where: param [default: binPsdDay] the configuration file name net [required] network code sta [required] station code loc [required] location ID chan [required] channel ID xtype [required] X-axis type for output (period or frequency) start [required] start date-time (UTC) of the binning interval (format YYYY-MM-DDTHH:MM:SS) end [required] end date-time (UTC) of the binning interval (format YYYY-MM-DDTHH:MM:SS) verbose [0 or 1, default: 0] to run in verbose mode set to 1The ntk_computePSD.py output files
hour X-value (period/frequency) Power with values separated using the “separator” character specified in the parameter file.
h3. OUTPUT FILE:
Full paths to the daily and hourly output data files are provided at the end of the run. You may turn off hourly file output via the configuration file.
The daily and hourly output file names have the following form respectively:
D?.bin and H?.bin (??? represents 3-digit day of the year)
for example:
D227.bin and H227.bin
* usage:
python bin/ntk_binPsdDay.py
* Assuming that you already have tried the following ntk_compute_PSD.py example:
python bin/ntk_computePSD.py net=TA sta=O18A loc=DASH start=2008-08-14T12:00:00 end=2008-08-14T13:30:00
* Then you can perform binning via:
python bin/ntk_binPsdDay.py net=TA sta=O18A loc=DASH chan=BHZ start=2008-08-14T12:00:00 end=2008-08-14T13:30:00 xtype=period
Previous release (V.1) docs: