1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6from __future__ import absolute_import
8import time
9import logging
10from collections import defaultdict
12import numpy as num
14from pyrocko.gui.snuffling import Snuffling, Param, Switch
15from pyrocko.trace import Trace, NoData
16from pyrocko import util
18logger = logging.getLogger('pyrocko.gui.snuffling.rms')
21class RootMeanSquareSnuffling(Snuffling):
23 '''
24 Create traces with blockwise root mean square values.
25 '''
27 def setup(self):
28 '''
29 Customization of the snuffling.
30 '''
32 self.set_name('RMS')
34 self.add_parameter(Param(
35 'Highpass [Hz]', 'highpass', None, 0.001, 1000.,
36 low_is_none=True))
38 self.add_parameter(Param(
39 'Lowpass [Hz]', 'lowpass', None, 0.001, 1000.,
40 high_is_none=True))
42 self.add_parameter(Param(
43 'Block Length [s]', 'block_length', 100., 0.1, 3600.))
45 self.add_parameter(Switch(
46 'Group channels', 'group_channels', True))
48 self.add_parameter(Switch(
49 'Log RMS', 'log', False))
51 self.add_trigger(
52 'Copy passband from Main', self.copy_passband)
54 self.set_live_update(False)
56 def copy_passband(self):
57 viewer = self.get_viewer()
58 self.set_parameter('lowpass', viewer.lowpass)
59 self.set_parameter('highpass', viewer.highpass)
61 def call(self):
62 '''
63 Main work routine of the snuffling.
64 '''
66 self.cleanup()
68 tinc = self.block_length
70 tpad = 0.0
71 for freq in (self.highpass, self.lowpass):
72 if freq is not None:
73 tpad = max(tpad, 1./freq)
75 targets = {}
76 tlast = time.time()
77 for batch in self.chopper_selected_traces(
78 tinc=tinc,
79 tpad=tpad,
80 want_incomplete=False,
81 style='batch',
82 mode='visible',
83 progress='Calculating RMS',
84 responsive=True,
85 fallback=True):
87 tcur = batch.tmin + 0.5 * tinc
89 if self.group_channels:
90 def grouping(nslc):
91 return nslc[:3] + (nslc[3][:-1],)
93 else:
94 def grouping(nslc):
95 return nslc
97 rms = defaultdict(list)
98 for tr in batch.traces:
100 # don't work traces produced by this (and other) snuffling
101 if tr.meta and 'tabu' in tr.meta and tr.meta['tabu']:
102 continue
104 if self.lowpass is not None:
105 tr.lowpass(4, self.lowpass, nyquist_exception=True)
107 if self.highpass is not None:
108 tr.highpass(4, self.highpass, nyquist_exception=True)
110 try:
111 tr.chop(batch.tmin, batch.tmax)
112 val = 0.0
113 if tr.ydata.size != 0:
114 y = tr.ydata.astype(float)
115 if num.any(num.isnan(y)):
116 logger.error(
117 'NaN value in trace %s (%s - %s)' % (
118 '.'.join(tr.nslc_id),
119 util.time_to_str(tr.tmin),
120 util.time_to_str(tr.tmax)))
122 val = num.sqrt(num.sum(y**2)/tr.ydata.size)
124 rms[grouping(tr.nslc_id)].append(val)
126 except NoData:
127 continue
129 tnow = time.time()
131 insert_now = False
132 if tnow > tlast + 0.8:
133 insert_now = True
134 tlast = tnow
136 to_add = []
137 for key, values in rms.items():
138 target = targets.get(key, None)
140 if not target \
141 or abs((target.tmax + tinc) - tcur) > 0.01 * tinc \
142 or insert_now:
144 if target:
145 to_add.append(target)
147 target = targets[key] = Trace(
148 network=key[0],
149 station=key[1],
150 location=key[2],
151 channel=key[3] + '-RMS',
152 tmin=tcur,
153 deltat=tinc,
154 ydata=num.zeros(0, dtype=float),
155 meta={'tabu': True})
157 value = num.atleast_1d(
158 num.sqrt(num.sum(num.array(values, dtype=float)**2)))
160 if self.log and value != 0.0:
161 value = num.log(value)
163 targets[key].append(value)
165 if to_add:
166 self.add_traces(to_add)
168 # add the newly created traces to the viewer
169 if targets.values():
170 self.add_traces(list(targets.values()))
173def __snufflings__():
174 '''
175 Returns a list of snufflings to be exported by this module.
176 '''
178 return [RootMeanSquareSnuffling()]