API Documentation of Qopen¶
Qopen consists of the following modules:
core |
Qopen command line script and routines |
imaging |
Plotting functions |
rt |
Calculate or plot spectral energy densitiy Green’s function |
site |
Align site responses of different runs of Qopen |
source |
Fit source displacement spectrum with source model |
util |
Some utility functions |
core
Module¶
Qopen command line script and routines
run_cmdline()
is started by the qopen
command line script.
Import and call run()
if you want to use Qopen inside Python code:
>>> from qopen import run
>>> run(conf='conf.json')
Qopen will run the following functions top-down:
run_cmdline |
Main entry point from the command line |
run |
Main entry point for a direct call from Python |
invert_wrapper |
Qopen function for a list or Catalog of events |
invert |
Qopen function to invert events and stations simultaneously |
invert_fb |
Inverst streams in a specific frequency band for attenuation parameters |
-
class
qopen.core.
ConfigJSONDecoder
(*, object_hook=None, parse_float=None, parse_int=None, parse_constant=None, strict=True, object_pairs_hook=None)[source]¶ Decode JSON config with comments stripped
-
qopen.core.
Gsmooth
(G_func, r, t, v0, g0, smooth=None, smooth_window='flat')[source]¶ Return smoothed Green’s function as a function of time
-
qopen.core.
collect_results
(results, only=None, freqi=None)[source]¶ Collect g0, b, error, R, W, eventids and v0 from results of multiple events
Parameters: - results – result dictionary returned by
run()
- only – return only some of the above mentioned keys
- freqi – return only values at a single frequency
Returns: dictionary
- results – result dictionary returned by
-
qopen.core.
energy1c
(data, rho, df, fs=4)[source]¶ Spectral energy density of one channel
Parameters: - data – velocity data (m/s)
- rho – density (kg/m**3)
- df – filter width in Hz
- fs – free surface correction (default: 4)
Returns: energy density
-
qopen.core.
filter_width
(sr, freq=None, freqmin=None, freqmax=None, corners=2, zerophase=False, type='bandpass')[source]¶ Integrate over the squared filter response of a Butterworth filter
The result corresponds to the filter width, which equals approximately the difference of the corner frequencies. The energy density should be divided by the result to get the correct spectral energy density.
Parameters: - sr – sampling rate
- freq – corner frequencies of low- or highpass filter
- freqmin,freqmax – corner frequencies of bandpass filter
- corners – number of corners
- zerophase – if True number of corners are doubled
- type – ‘bandpass’, ‘highpass’ or ‘lowpass’
Returns: filter width
-
qopen.core.
get_freqs
(max=None, min=None, step=None, width=None, cfreqs=None, fbands=None)[source]¶ Determine frequency bands
Parameters: args – See example configuration file. Returns: ordered dictionary {central frequency: corner frequencies}
-
qopen.core.
get_magnitude
(event)[source]¶ Preferred or first magnitude of event
This is not the coda magnitude determined by the script, but the magnitude from the original event catalogue.
-
qopen.core.
init_data
(data, client_options=None, plugin=None, cache_waveforms=False, get_waveforms=None)[source]¶ Return appropriate get_waveforms function
See example configuration file for a description of the options
-
qopen.core.
invert
(events, inventory, get_waveforms, request_window, freqs, filter, rho0, vp=None, vs=None, remove_response=None, plot_remove_response=False, remove_response_options={}, skip=None, use_picks=False, correct_for_elevation=False, subtract_local_depth=False, njobs=None, seismic_moment_method=None, seismic_moment_options={}, plot_eventresult=False, plot_eventresult_options={}, plot_eventsites=False, plot_eventsites_options={}, plot_results=False, plot_results_options={}, plot_sites=False, plot_sites_options={}, plot_sds=False, plot_sds_options={}, plot_mags=False, plot_mags_options={}, cmd='go', input=None, coda_normalization=None, request_window_tolerance=0.5, **kwargs)[source]¶ Qopen function to invert events and stations simultaneously
Parameters: - events – is determined in
invert_wrapper()
- inventory,get_waveforms – are determined in
run()
All other options are described in the example configuration file.
Returns: result dictionary
- events – is determined in
-
qopen.core.
invert_fb
(freq_band, streams, filter, rho0, v0, coda_window, R0=1, free_surface=4, noise_windows=None, noise_windows_func='min', bulk_window=None, weight=None, optimize=None, g0_bounds=(1e-08, 0.001), b_bounds=(1e-05, 10), num_points_integration=1000, coda_normalization=None, smooth=None, smooth_window='flat', remove_noise=False, cut_coda=None, skip=None, adjust_sonset=None, adjust_sonset_options={}, plot_energies=False, plot_energies_options={}, plot_optimization=False, plot_optimization_options={}, plot_fits=False, plot_fits_options={}, ignore_network_code=False, borehole_stations=(), G_plugin='qopen.rt : G_rt3d', fix=False, fix_params=False, fix_sites=False, fix_sites_params=None, dump_optpkl=None, dump_fitpkl=None, **kwargs)[source]¶ Inverst streams in a specific frequency band for attenuation parameters
Parameters: freq_band, streams, borehole_stations, fix, fix_params , fix_sites, fix_sites_params – are determined in invert()
. All other options are described in the example configuration file.Returns: result tuple
-
qopen.core.
invert_wrapper
(events, plot_results=False, plot_results_options={}, plot_sites=False, plot_sites_options={}, plot_sds=False, plot_sds_options={}, plot_mags=False, plot_mags_options={}, invert_events_simultaneously=False, mean=None, noplots=False, **kwargs)[source]¶ Qopen function for a list or Catalog of events
Depending on ‘invert_events_simultaneously’ flag the function calls
invert()
for each event seperately or for all events once. In the first case mean results are calculated.Parameters: events – is determined in run()
. The rest of the options are described in the example configuration file.Returns: result dictionary
-
qopen.core.
observed_energy
(stream, rho, df, coda_normalization=None, fs=4, tolerance=1)[source]¶ Return trace with total spectral energy density of three component stream
Parameters: - stream – stream of a 3 component seismogram
- rho – density (kg/m**3)
- df – filter width in Hz
- fs – free surface correction (default: 4)
- tolerance – the number of samples the length of the traces in the 3 component stream may differ (default: 1)
Returns: trace with total energy density
-
qopen.core.
run
(cmd='go', conf=None, pdb=False, tutorial=False, eventid=None, get_waveforms=None, print_mag=False, plots=None, **args)[source]¶ Main entry point for a direct call from Python
Example usage:
>>> from qopen import run >>> run(conf='conf.json')
Parameters: - args –
All args correspond to the respective command line and configuration options. See the example configuration file for help and possible arguments. Options in args can overwrite the configuration from the file. E.g.
run(conf='conf.json', events=my_catalogue)
will ignoreevents
value in the configuration file.Exceptions from the description in configuration file:
- events – can be filename or ObsPy Catalog object
- inventory – can be filename or ObsPy Inventory object
- get_waveforms – function, if given the data option will be ignored. get_waveforms will be called as described in the example configuration file
Returns: result dictionary
- args –
-
qopen.core.
time2utc
(time, trace, starttime=None)[source]¶ Convert string with time information to UTCDateTime object
Parameters: - time –
can be one of:
”OT+-???s” seconds relative to origin time
”P+-???s” seconds relative to P-onset
”S+-???s” seconds relative to S-onset
”???Ptt” travel time relative to P-onset travel time
”???Stt” travel time relative to S-onset travel time
- ”???SNR” time after which SNR falls below this value
- (after time given in starttime)
- ”time>???SNR” time after which SNR falls below this value
- (after time given in front of expression)
- trace – Trace object with stats entries
- starttime – UTCDatetime object for SNR case.
- time –
-
qopen.core.
tw2utc
(tw, trace)[source]¶ Convert time window to UTC time window
Parameters: - tw – tuple of two values, both can be a string (see
time2utc()
) or a list of strings in which case the latest starttime and earliest endtime is taken. - trace – Trace object with stats entries
- tw – tuple of two values, both can be a string (see
imaging
Module¶
Plotting functions
Arguments supported by all plotting functions via its **kwargs are:
fname: | file name for the plot output (if not provided the figure will be left open) |
---|---|
title: | title of the plot |
figsize: | figure size (tuple of inches) |
dpi: | resolution of image file |
-
qopen.imaging.
calc_dependent
(quantity, value, freq=None, v0=None)[source]¶ Calculate dependent value (Qsc, Qi, lsc, li) from g0 and b
Parameters: - quantity (str) – one of Qsc, Qi, lsc, li
- value – value of g0 or b depending on requested quantity
- freq – frequency in Hz (needed for some calculations)
- v0 – velocity in m/s (needed for some calculations)
Returns: value of quantity
-
qopen.imaging.
plot_all_sds
(result, seismic_moment_method=None, seismic_moment_options=None, xlim=None, ylim=None, nx=None, annotate=None, va='top', annotate_label=None, annotate_evid=True, plot_only_ids=None, cmap='viridis_r', vmin=None, vmax=None, colors=None, xlabel='frequency (Hz)', ylabel='source displacement spectrum $\\omega$M (Nm)', **kwargs)[source]¶ Plot all source displacement spectra with fitted source models
-
qopen.imaging.
plot_energies
(energies, bulk_window=None, coda_window=None, downsample_to=None, xlim_lin=None, xlim_log=None, **kwargs)[source]¶ Plot observed spectral energy densities on different scales (linear, log)
-
qopen.imaging.
plot_eventresult
(result, v0=None, quantities=('g0', 'lsc', 'Qsc', 'b', 'li', 'Qi', 'error', 'W', 'sds'), seismic_moment_method=None, seismic_moment_options={}, xlabel='frequency (Hz)', **kwargs)[source]¶ Plot all results of one
invert()
call
-
qopen.imaging.
plot_eventsites
(result, xlabel='frequency (Hz)', ylabel='energy site amplification', **kwargs)[source]¶ Plot site amplification factors of one
invert()
call
-
qopen.imaging.
plot_fits
(energies, g0, b, W, R, v0, info, G_func, smooth=None, smooth_window='bartlett', xlim=None, ylim=None, xlabel='time (s)', ylabel='E (Jm$^{-3}$Hz$^{-1}$)', **kwargs)[source]¶ Plot fits of spectral energy densities
-
qopen.imaging.
plot_lstsq
(rec, event_station_pairs, ax=None, base=<Mock object>, **kwargs)[source]¶ Plot solution of weighted least squares inversion
-
qopen.imaging.
plot_mags
(result, xlim=None, ylim=None, plot_only_ids=None, xlabel='M from catalog', ylabel='Mw from inversion', secondary_yaxis=True, **kwargs)[source]¶ Plot Qopen moment magnitudes versus catalogue magnitudes
-
qopen.imaging.
plot_optimization
(record, record_g0, event_station_pairs, num=7, xlabel='g$_0$ ($\\mathrm{m}^{-1}$)', ylabel='misfit $\\epsilon$', **kwargs)[source]¶ Plot some steps of optimization
-
qopen.imaging.
plot_results
(result, v0=None, quantities=('g0', 'lsc', 'Qsc', 'b', 'li', 'Qi', 'error', 'nobs'), mean=None, llim=None, Qlim=None, xlabel='frequency (Hz)', **kwargs)[source]¶ Plot results
rt
Module¶
Calculate or plot spectral energy densitiy Green’s function
The Green’s function used by the qopen
command
line program can be switched with the “G_plugin” configuration option.
Available Green’s functions:
rt3d (default) ... Isotropic radiative transfer, 3D approximative
interpolation solution of Paasschens (1997)
rt2d/1d ... Isotropic radiative transfer, 2D or 1D solution
diff3d/2d/1d ... Diffusion, 3D, 2D or 1D
diffapprox3d/2d/1d ... Approximation of diffusion, 3D, 2D or 1D,
often used to determine coda Q
Used variables:
r ... Distance to source in m
t ... Time after source in s
c ... Velocity in m/s
l ... (transport) mean free path in m
la ... absorption length in m
g0 = 1/l ... (transport) scattering coefficient
Note
The formulas for radiative transfer are valid for the
scattering coefficient (g0) under assumption of isotropic scattering.
However, g0 is used by qopen.core
module as transport scattering
coefficient (g*) under the assumption of non-isotropic scattering.
g*=g0
is a reasonable assumption under these conditions.
-
qopen.rt.
G
(r, t, c, g0, type='rt3d', include_direct=True)[source]¶ Full Green’s function with direct wave term (optional)
-
qopen.rt.
G_diff3d
(r, t, c, g0)[source]¶ Green’s function for 3d diffusion
See Paaschens (1997), equation 1.
-
qopen.rt.
G_diffapprox3d
(r, t, c, g0)[source]¶ Green’s function for 3d diffusion, discarding second term.
Often used to determine Qcoda
See Paaschens (1997), equation 1.
-
qopen.rt.
G_rt1d
(r, t, c, g0)[source]¶ Full Green’s function for 1d radiative transfer
See Paaschens (1997), equation 2, originally from Hemmer (1961).
-
qopen.rt.
G_rt2d
(r, t, c, g0)[source]¶ Full Green’s function for 2d radiative transfer.
See Paaschens (1997), equation 26.
-
qopen.rt.
G_rt3d
(r, t, c, g0)[source]¶ Full Green’s function for 3d radiative transfer (approximation).
See Paaschens (1997), equation 36.
-
qopen.rt.
plot_r
(c, g0, t, r=None, N=1000, log=False, include_direct=False, la=None, type='rt3d', scale=False)[source]¶ Plot Green’s function as a function of distance
site
Module¶
Align site responses of different runs of Qopen
Consider two runs of Qopen with stations A, B in run 1 and stations
B, C in run 2. Qopen assumes a mean site amplification for each run of 1
with the default configuration.
Use align_site_responses
or the command line option qopen --align-sites
to correct source power and site amplification factors afterwards.
In the above case site responses and source powers will be adjusted such that
the site response of station B is the same for both runs.
-
qopen.site.
align_site_responses
(results, station=None, response=1.0, use_sparse=True, seismic_moment_method=None, seismic_moment_options=None, ignore_stations=None)[source]¶ Align station site responses (amplification) and correct source parameters
Determine best factor for each event so that site response is the same for each station and different events.
Parameters: results – original result dictionary. For the other options see the help for the corresponding command line options or configuration parameters. Returns: corrected result dictionary
source
Module¶
Fit source displacement spectrum with source model
(and some other functions dealing with the source)
If you want to fit source displacement spectra on the command line again
use the qopen --calc-source-params
option.
-
qopen.source.
calculate_source_properties
(results, rh0=None, v0=None, seismic_moment_method=None, seismic_moment_options=None)[source]¶ Calculate source porperties for results
-
qopen.source.
fit_sds
(freq, omM, method='mean', fc=None, n=2, gamma=1, fc_lim=None, n_lim=(0.5, 10), gamma_lim=(0.5, 10), fc0=10, n0=2, gamma0=1, fall_back=5, num_points=None, **opt_kw)[source]¶ Fit source displacement spectrum and calculate seismic moment
Parameters: - freq,omM – frequencies, source displacement spectrum (same length)
- method – ‘mean’ - take mean of sds of frequencies below fc, ‘fit’, ‘robust_fit’ - fit source model to obtain M0. If one or more of fc, n, gamma are None, M0 and these values are simultaneously determined. Robust version uses a robust linear model (which downweights outliers).
- fc,n,gamma – corner frequency and coefficients for source model
- fc_lim,gamma_lim – bounds for corner frequency and gamma (used for optimization if respective variable is set to None)
- fc0,gamma0 – starting values of fc and gamma for optimization (only used for optimization for fc and gamma)
- fall_back – use robust fit only if number of data points >= fall_back
- num_points – determine M0 only if number of data points >= num_points All other kwargs are passed to scipy.optimization, e.g.
- tol – tolerance for optimization
Returns: dictionary with M0 and optimized variables fc, n, and gamma if applicable. If M0 is not determined the function will return None
-
qopen.source.
insert_source_properties
(freq, evresult, v0, rho0, seismic_moment_method, seismic_moment_options, catmag=None)[source]¶ Insert sds, Mw and possibly Mcat in evresult dictionary
-
qopen.source.
moment_magnitude
(M0, inverse=False)[source]¶ Moment magnitude Mw from seismic moment M0
Based on Kanamori (1977), an alternative definition is based on Hanks and Kanamori (1979) with an offset of -6.03. Difference due to rounding errors.
Parameters: - M0 – seismic moment in Nm
- inverse – return the inverse relation ship M0(Mw)
util
Module¶
Some utility functions
-
qopen.util.
gstat
(data, axis=None, weights=None, unbiased=True, robust=False, fall_back=5)[source]¶ Weighted or robust geometric mean and error (log scale)
-
qopen.util.
linear_fit
(y, x, m=None, method='robust', **kw)[source]¶ Linear fit between x and y
Parameters: - y,x – data
- m – fix slope at specific value
- method – one of (‘least_squares’, ‘weighted’, ‘robust’)
Returns: slope a and intercept b of y = ax + b
-
qopen.util.
robust_stat
(data, axis=None, fall_back=5)[source]¶ Robust mean and mean absolute deviation
See also: statsmodels.sf.net/stable/rlm.html
-
qopen.util.
smooth
(x, window_len=None, window='flat', method='zeros')[source]¶ Smooth the data using a window with requested size.
This method is based on the convolution of a scaled window with the signal.
Parameters: - x – the input signal
- window_len – the dimension of the smoothing window; should be an odd integer
- window – the type of window from ‘flat’, ‘hanning’, ‘hamming’, ‘bartlett’, ‘blackman’ flat window will produce a moving average smoothing.
- method –
handling of border effects
’zeros’: zero padding on both ends (len(smooth(x)) = len(x))
’reflect’: pad reflected signal on both ends (same)
’clip’: pad signal on both ends with the last valid value (same)
None: no handling of border effects (len(smooth(x)) = len(x) - len(window_len) + 1)
See also: www.scipy.org/Cookbook/SignalSmooth
Template Configuration File¶
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 | ### Configuration file for qopen package in JSON format # Comments are indicated with "#" and ignored while parsing { ### Options for input and output ### # Filename of events file (format supported by ObsPy) "events": "example_events.xml", # Filter event catalog with Obspy method Catalog.filter # example: "filter_events": ["magnitude > 2.5"], "filter_events": null, # Filename of inventory of stations (format supported by ObsPy) "inventory": "example_inventory.xml", # Filter inventory with Obspy method Inventory.select # example: "filter_inventory": {"channel": "HH?"}, "filter_inventory": null, # Data can be # 1. a glob expression of files. All files are read at once. If you have # a huge dataset think about using option 3. # 2. one of the client modules supported by ObsPy (e.g "fdsn") # for getting the data from a webservice. # Option "client_options" is available. # 3. "plugin" e.g. for a heterogenous dataset. You have to tell Qopen # how to access the data. In this case the option "plugin" is # available. "data": "example_data.mseed", # Options for the webservices which are passed to Client.__init__. # See the documentation of the clients in ObsPy for available options. "client_options": {"user": "name@insitution.com"}, # Plugin "module : func" # 'module' has to be importable (located in current path or PYTHON_PATH # environment variable). # 'func' has to be the function which returns the data. # Kwargs passed to func are: network, station, location, channel, # starttime, endtime and event # Example using FSDN client: # # in module.py # from obspy.clients.fsdn import Client # client = Client() # def func(**kwargs): # kwargs.pop('event') # event kwarg is not needed by client # return client.get_waveforms(**kwargs) "plugin": "module : func", # Request window of data in s relative to origin time "request_window": [-10, 210], # Directory for optional waveform caching from webservice or plugin. # This works only if the waveforms are requested with the exact same # parameters (start and end time) i.e. if options # events and request_window are left untouched # This feature needs the joblib library "cache_waveforms": null, # Output results as json on stdout or in a file # possible options: null, "stdout", filename "output": "results.json", # Indent of json output, possible values: null, 0, 2, 4, ... "indent": 2, # Possibility to define your own logging configuration. This dictionary # will be passed to logging.config.dictConfig. # null will log to a file and stdout depending on 'logfile', 'loglevel' # and '-v' flag "logging": null, # Log file or null. Output to stdout can be configured with the -v flag "logfile": "example.log", # The higher the log level the more verbose is the logging to the file (1-3) "loglevel": 3, # Use multiple cores for calculating the results # 1 -> sequentuell processing # 2 -> use 2 cores # null -> use all available cores "njobs": null, ### Options for inversion ### # If true: all events are inverted together for one b and g0 # If false: events are inverted independently for b and g0 # You only want to set this option to true if only 1 or 2 stations are available "invert_events_simultaneously": false, # If previous option is false, the results are determined as a mean # of the inversions for individual events. "mean" can take the values # "robust" (downweights outliers), "weigthed" (uses errors), all other values # use the normal average. "mean": "robust", # Mean shear wave velocity in m/s for calculating the Green's function. # Has also impact on the calculation of the source displacement spectrum # ωM: ωM ~ v0^(5/2) # Furthermore, necessary for calculating Qsc from g0 and li from b # This parameter can also be calculated from picks. If you want to # use this feature set "v0": null and "use_picks": true. "v0": 3400, # Mean density in kg/m**3. # Important for calculating the correct observed energy and therefore # scales the spectral source energy W: W ~ rho0 and the # source displacement spectrum ωM ~ rho0 "rho0": 2700, # Geometric mean of station amplification factors for energy, # can be assumed as 1. # Important for calculating the correct modeled energy and therefore # scales W: W ~ 1/R0, ωM ~ 1/R0^(1/2) "R0": 1, # Correction factor for free surface. # According to Emoto (2010) the factor is around 4 for S-waves. # Specify a tuple like (4, 1) to use the first value for surface stations # and the second value for borehole stations. By default correction factors # for surface and borehole stations are the same to allow the direct comparison # of obtained site responses. "free_surface": 4, # Optimization options for the scattering coefficient. # Values are passed to scipy.optimize.minimize_scalar. # 'method' can be one of 'brent', 'golden', 'bounded' # 'tol' defines the relative tolerance for finding the optimal solution "optimize": {"method": "golden", "tol": 1e-3}, # Bounds for optimization. Optimal parameters near the bounds will be # treated as unphysical and corresponding results will be discarded # 'g0_bounds' bounds for (transport) scattering coefficient # 'b_bounds' bounds for intrinsic attenuation "g0_bounds": [1e-8, 1e-4], "b_bounds": [1e-3, 10], # Remove instrument response. One of null, "sensitivity", "full". # Instrument response has to be available in the inventory. # This option is important for the calculation of W and ωM "remove_response": "sensitivity", # only for full response removal, see help in ObsPy's Trace.remove_response "plot_remove_response": false, "remove_response_options" : { "fname": "plots/remove_response_{evid}_{tr.id}.png"}, # When calculating the distance, correct for the elevation of the station. # Useful if depth of events catalog is given relative to sea level. "correct_for_elevation": false, # Subtract the local sensor depth from the channel (if not present: station) elevation # If StationXML definition is honored use false. # If elevation is only specified on the station level use true. # Only relevant for setting "correct_for_elevation": true. "subtract_local_depth": false, # Filter options passed to obspy.core.Stream.filter for filtering the # different frequency bands (see obspy.core.Stream.filter). "filter": {"corners": 2, "zerophase": true}, # Define frequency bands: # step - difference between neighboring central frequencies in octave # (x octave corresponds to a factor of 2**x between the frequencies) # width - difference between higher and lower corner frequency in octave # max - maximal central frequency, # min - minimal possible central frequency. # cfreqs - list of central frequencies (min, max and step are ignored) # fbands - list frequency bands (all other entries are ignored) #"freqs": {"width": 1, "cfreqs": [1.5, 3, 6]}, #"freqs": {"fbands": [[1, 2], [2, 4], [4, 8]]}, "freqs": {"step": 1, "width": 1, "max": 6, "min": 0.3}, # Use picks given in events file to determine the P- and S-onset. # Picks can also be used to calculate v0. If a station does not have a # S pick (S or Sg) it is ignored. "use_picks": false, # Otherwise these velocities are used together with the distance between # preferred origin and station to determine the onsets "vp": 6000, "vs": 3400, # Determine the S-onset again. This option is for testing purposes when # no picks are available. It should be prefered to use picks. # Possible values: null, "maximum" (take maximum in given window as new # S-onset) # Example: # "adjust_sonset": "maximum", # "adjust_sonset_options": {"window": ["S-10s", "S+10s"]}, "adjust_sonset": null, # Remove noise level from data "remove_noise": true, # Definition of three different time windows # - List of time windows to determine the noise level (noise_windows) # (the minimum of the noise levels determined in the individual time # windows is used) # - Time window for calculation of the direct energy (bulk_window) # The energy is integrated in this time window # - Time window defining the coda (coda_window) # Time windows are a list of two values (start and end time). # The start time and end time can be a string, e.g.: # - 'P-10s', 'S+10s', 'OT+10s' (seconds relative to onset or origin time), # - '1.2Stt', '0.8Ptt' (travel time relative to onset travel time '1Stt'), # - '3SNR' (only for end time of codawindow, time after start time when # the energy hits this level relative to the noise level) or # - 'time>3SNR' (same as previous, but time after which the noise level # is checked is explicitly given, e.g. "S+10s>3SNR") # The start time and end time can also be a list, in this case the largest # start time and the smallest end time is chosen. "noise_windows": [["OT-10s", "OT-5s"], ["OT-5s", "OT+0s"]], "bulk_window": ["S-2s", "S+20s"], "coda_window": ["S+20s", ["S+150s", "3SNR"]], # Coda window ends 150s after S-onset or if SNR of 3 is reached. # The weight of the bulk window (list of value and unit) # unit is one of: # - codawindow (the same weight has the whole coda window) # - bulkwindow (the same weight as the length of the bulk window # relative to the coda window) # - samples (the same weight as one sample in the coda window) "weight" : [1, "bulkwindow"], # Smooth the coda over 'smooth' seconds (null for no smoothing) "smooth": 1, # The window used for smoothing. Can be one of # "flat", "hanning", "hamming", "bartlett" or "blackman". # "flat" corresponds to the moving average, for the rest see # http://docs.scipy.org/doc/numpy/reference/routines.window.html "smooth_window": "bartlett", # Cut coda if a local minimum is present in data (e.g. because of a second # earthquake with bulk waves arriving in the coda). Local maxima in the coda # are compared to previous local minima. If the ratio between a maximum and # the lowest minimum before the maximum is above a threshold the coda window # will end at time of the minimum. # Coda can optionally be smoothed before this operation. # Example: # "cut_coda": {"smooth": 5, "ratio": 5}, "cut_coda": null, # Skip station if one of these conditions is fulfilled. # Possible conditions: # - "coda_window" shorter than x seconds # - "distance" between station and event larger than x km # - "maximum" of envelope is not in given window # - "num_pairs" skip event if number of station event pairs smaller than x # Example: # "skip": {"coda_window": 5, # "distance": 500, # "maximum": ["0.5Stt", "S+10s"], # "num_pairs": 2}, "skip": null, # Calculate the seismic moment from the source displacement spectrum. # Possible values: # null: no calculation of the seismic moment # "mean" of the values for frequencies smaller than 'fc', # "fit" fit source displacement spectrum to source model # "robust_fit" robustly fit source displacement spectrum to source model # The model is defined in qopen.source.source_model() # The options passed to qopen.source.fit_sds() (see documentation there) are # defined in the seismic_moment_options dictionary. It is for example possible # to invert the source displacement spectrum for fc, n and gamma simultaneously # or to fix a subset of the parameters and invert for the others. #"seismic_moment_method": "mean", #"seismic_moment_options": {"fc": 1.5}, "seismic_moment_method": "robust_fit", "seismic_moment_options": {"fc": null, "n": null, "gamma": 2, "fc_lim": [0.1, 10], "num_points": 4}, # Optionally, a user defined alternative Green's function for scattering can be used. # The module has to define a function with the following arguments: # G(r, t, c, g0): # return Green's function # Some Green's functions are provided in qopen.rt module: # G_rt3d, G_rt2d, G_rt1d (radiative transfer 3d, 2d, 1d) # G_diff3d (diffusion 3d) # G_diffapprox3d (diffusion approximation 3d, can be used with # "optimize": null, "bulk_window": null # to invert envelopes for coda Q) #"G_plugin": "qopen.rt : G_rt3d", # Option for coda normalization of the envelope. In the case of coda normalization, # the source spectrum and site amplification are not inverted. The window for coda # normalization is specified in seconds after origin time. "coda_normalization": null, #"coda_normalization": [180, 200], ### Plot options ### # For all parameters which can be used in the plot_*_options dictionaries see the qopen.imaging module. # Plot observed energies # xlim_min: time limits for linear plot (relative to origin time) # xlim_log: time limits for log plot (relative to origin time) # fname: filename gets populated by eventid and freq_band "plot_energies": true, "plot_energies_options": {"fname": "plots/energies_{evid}_{f1:05.2f}Hz-{f2:05.2f}Hz.png", "xlim_lin": [0, 210], "xlim_log": [10, 210]}, # Plot optimization routine # num: Number of fit plots # fname: filename gets populated by eventid and freq_band "plot_optimization": true, "plot_optimization_options": {"num": 7, "fname": "plots/optimization_{evid}_{f1:05.2f}Hz-{f2:05.2f}Hz.png"}, # Plot fits for optimization result # fname: filename gets populated by eventid and freq_band "plot_fits": true, "plot_fits_options": {"fname": "plots/fits_{evid}_{f1:05.2f}Hz-{f2:05.2f}Hz.png"}, # Plot scattering parameters and magnitude estimation for one event # fname: filename gets populated by eventid "plot_eventresult": true, "plot_eventresult_options": {"fname": "plots/eventresult_{evid}.png"}, # Plot site corrections for one event # fname: filename gets populated by eventid "plot_eventsites": true, "plot_eventsites_options": {"fname": "plots/eventsites_{evid}.png"}, # Plot mean attenuation and scattering properties # fname: filename "plot_results": true, "plot_results_options": {"fname": "plots/results.pdf", "llim": [1e0, 1e3], "Qlim": [1e-5, 1e-1]}, # Plot mean site amplification factors # fname: filename "plot_sites": true, "plot_sites_options": {"fname": "plots/sites.pdf"}, # Plot source displacement spectra # fname: filename "plot_sds": true, "plot_sds_options": {"fname": "plots/sds.pdf"}, # Plot comparison between magnitude determined by this script and magnitude # from event file # fname: filename "plot_mags": true, "plot_mags_options": {"fname": "plots/mags.pdf"} } |