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

decode(s)[source]

Return the Python representation of s (a str instance containing a JSON document).

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

exception qopen.core.ParseError[source]
exception qopen.core.QopenError[source]
exception qopen.core.SkipError[source]
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

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_arrivals(event)[source]

Arrivals of appropriate origin from event

qopen.core.get_eventid(event)[source]

Event id from event

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.get_origin(event)[source]

Preferred or first origin from event

qopen.core.get_pair(tr)[source]

Station and event id from trace

qopen.core.get_picks(arrivals, station)[source]

Picks for specific station from arrivals

qopen.core.get_station(seedid)[source]

Station name from seed id

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

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 ignore events 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

qopen.core.run_cmdline(args=None)[source]

Main entry point from the command line

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

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

qopen.imaging.plot_sds(freq, result, ax=None, annotate=False, va='bottom', annotate_label=None, seismic_moment_method=None, seismic_moment_options={}, cmap='viridis_r', vmin=None, vmax=None, max_nobs=None, color=None, **kwargs)[source]

Plot source displacement spectrum and fitted source model

qopen.imaging.plot_sites(result, mean=None, xlim=None, ylim=(0.01, 100.0), nx=None, cmap='viridis_r', vmin=None, vmax=None, xlabel='frequency (Hz)', ylabel='site amplification', show_excluded=True, sortkey=None, annotate=True, **kwargs)[source]

Plot site amplification factors

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

qopen.rt.plot_rt(c, g0, t=None, r=None, N=1000, log=False, include_direct=False, la=None, type='rt3d')[source]

Plot Green’s function as a function of distance

qopen.rt.plot_t(c, g0, r, t=None, N=1000, log=False, include_direct=False, la=None, type='rt3d', scale=False)[source]

Plot Green’s function as a function of time

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)
qopen.source.sds(W, f, v, rho)[source]

Calculate source displacement spectrum ωM from spectral source energy W

according to Sato & Fehler (2012, p.188)

Parameters:
  • W – spectral source energy (J/Hz)
  • f,v,rho – frequency, mean velocity, mean density
Returns:

source displacement spectrum (in Nm)

qopen.source.source_model(freq, M0, fc, n=2, gamma=1)[source]

Model for source displacement spectrum (Abercrombie 1995)

Parameters:
  • freq – frequencies
  • M0 – seismic moment (Nm)
  • fc – corner frequency
  • n – high frequency fall-of
  • gamma – corner sharpness

util Module

Some utility functions

qopen.util.gerr(data, **kwargs)[source]

Weighted or robust geometric mean and lower/upper error

qopen.util.gmean(data, **kwargs)[source]

Weighted or robust geometric mean

qopen.util.gmeanlist(*args, **kwargs)[source]

Weighted or robust geometric mean returned as list

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

qopen.util.smooth_func(f, t, window_len=None, window='flat')[source]

Smooth a function f at time samples t

qopen.util.weighted_stat(data, axis=None, weights=None, unbiased=True)[source]

Weighted mean and biased/unbiased standard deviation

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"}

}