Passed
Pull Request — master (#25)
by Ian
06:01
created

build.rsudp.helpers.deconvolve()   B

Complexity

Conditions 7

Size

Total Lines 36
Code Lines 17

Duplication

Lines 0
Ratio 0 %

Code Coverage

Tests 13
CRAP Score 7.6383

Importance

Changes 0
Metric Value
eloc 17
dl 0
loc 36
ccs 13
cts 17
cp 0.7647
rs 8
c 0
b 0
f 0
cc 7
nop 1
crap 7.6383
1 1
import rsudp.raspberryshake as rs
2 1
from rsudp import COLOR, printM, printW
3 1
import rsudp
4 1
import os
5 1
import json
6
7
8 1
def dump_default(settings_loc, default_settings):
9
	'''
10
	Dumps a default settings file to a specified location.
11
12
	:param str settings_loc: The location to create the new settings JSON.
13
	:param str default_settings: The default settings to dump to file.
14
	'''
15
	print('Creating a default settings file at %s' % settings_loc)
16
	with open(settings_loc, 'w+') as f:
17
		f.write(default_settings)
18
		f.write('\n')
19
20
21 1
def default_settings(output_dir='%s/rsudp' % os.path.expanduser('~').replace('\\', '/'), verbose=True):
22
	'''
23
	Returns a formatted json string of default settings.
24
25
	:param str output_dir: the user's specified output location. defaults to ``~/rsudp``.
26
	:param bool verbose: if ``True``, displays some information as the string is created.
27
	:return: default settings string in formatted json
28
	:rtype: str
29
	'''
30 1
	def_settings = r"""{
31
"settings": {
32
    "port": 8888,
33
    "station": "Z0000",
34
    "output_dir": "%s",
35
    "debug": true},
36
"printdata": {
37
    "enabled": false},
38
"write": {
39
    "enabled": false,
40
    "channels": ["all"]},
41
"plot": {
42
    "enabled": true,
43
    "duration": 30,
44
    "spectrogram": true,
45
    "fullscreen": false,
46
    "kiosk": false,
47
    "eq_screenshots": false,
48
    "channels": ["HZ", "HDF"],
49
    "deconvolve": true,
50
    "units": "CHAN"},
51
"forward": {
52
    "enabled": false,
53
    "address": ["192.168.1.254"],
54
    "port": [8888],
55
    "channels": ["all"],
56
	"fwd_data": true,
57
	"fwd_alarms": false},
58
"alert": {
59
    "enabled": true,
60
    "channel": "HZ",
61
    "sta": 6,
62
    "lta": 30,
63
    "threshold": 1.7,
64
    "reset": 1.6,
65
    "highpass": 0,
66
    "lowpass": 50,
67
    "deconvolve": false,
68
    "units": "VEL"},
69
"alertsound": {
70
    "enabled": false,
71
    "mp3file": "doorbell"},
72
"custom": {
73
    "enabled": false,
74
    "codefile": "n/a",
75
    "win_override": false},
76
"tweets": {
77
    "enabled": false,
78
    "tweet_images": true,
79
    "api_key": "n/a",
80
    "api_secret": "n/a",
81
    "access_token": "n/a",
82
    "access_secret": "n/a",
83
	"extra_text": ""},
84
"telegram": {
85
    "enabled": false,
86
    "send_images": true,
87
    "token": "n/a",
88
    "chat_id": "n/a"},
89
"rsam": {
90
    "enabled": false,
91
    "quiet": true,
92
    "fwaddr": "192.168.1.254",
93
    "fwport": 8887,
94
    "fwformat": "LITE",
95
    "channel": "HZ",
96
    "interval": 10,
97
    "deconvolve": false,
98
    "units": "VEL"}
99
}
100
101
""" % (output_dir)
102 1
	if verbose:
103
		print('By default output_dir is set to %s' % output_dir)
104 1
	return def_settings
105
106
107 1
def read_settings(loc):
108
	'''
109
	Reads settings from a specific location.
110
111
	:param str loc: location on disk to read json settings file from
112
	:return: settings dictionary read from JSON, or ``None``
113
	:rtype: dict or NoneType
114
	'''
115
	settings_loc = os.path.abspath(os.path.expanduser(loc)).replace('\\', '/')
116
	settings = None
117
	with open(settings_loc, 'r') as f:
118
		try:
119
			data = f.read().replace('\\', '/')
120
			settings = json.loads(data)
121
		except Exception as e:
122
			print(COLOR['red'] + 'ERROR: Could not load settings file. Perhaps the JSON is malformed?' + COLOR['white'])
123
			print(COLOR['red'] + '       detail: %s' % e + COLOR['white'])
124
			print(COLOR['red'] + '       If you would like to overwrite and rebuild the file, you can enter the command below:' + COLOR['white'])
125
			print(COLOR['bold'] + '       shake_client -d %s' % loc + COLOR['white'])
126
			exit(2)
127
	return settings
128
129
130 1
def set_channels(self, cha):
131
	'''
132
	This function sets the channels available for plotting. Allowed units are as follows:
133
134
	- ``["SHZ", "EHZ", "EHN", "EHE"]`` - velocity channels
135
	- ``["ENZ", "ENN", "ENE"]`` - acceleration channels
136
	- ``["HDF"]`` - pressure transducer channel
137
	- ``["all"]`` - all available channels
138
139
	So for example, if you wanted to display the two vertical channels of a Shake 4D,
140
	(geophone and vertical accelerometer) you could specify:
141
142
	``["EHZ", "ENZ"]``
143
144
	You can also specify partial channel names.
145
	So for example, the following will display at least one channel from any
146
	Raspberry Shake instrument:
147
148
	``["HZ", "HDF"]``
149
150
	Or if you wanted to display only vertical channels from a RS4D,
151
	you could specify
152
153
	``["Z"]``
154
155
	which would match both ``"EHZ"`` and ``"ENZ"``.
156
157
	:param self self: self object of the class calling this function
158
	:param cha: the channel or list of channels to plot
159
	:type cha: list or str
160
	'''
161 1
	cha = rs.chns if ('all' in cha) else cha
162 1
	cha = list(cha) if isinstance(cha, str) else cha
163 1
	for c in rs.chns:
164 1
		n = 0
165 1
		for uch in cha:
166 1
			if (uch.upper() in c) and (c not in str(self.chans)):
167 1
				self.chans.append(c)
168 1
			n += 1
169 1
	if len(self.chans) < 1:
170
			self.chans = rs.chns
171
172
173 1
def fsec(ti):
174
	'''
175
	.. versionadded:: 0.4.3
176
177
	The Raspberry Shake records at hundredths-of-a-second precision.
178
	In order to report time at this precision, we need to do some time-fu.
179
180
	This function rounds the microsecond fraction of a
181
	:py:class:`obspy.core.utcdatetime.UTCDateTime`
182
	depending on its precision, so that it accurately reflects the Raspberry Shake's
183
	event measurement precision.
184
185
	This is necessary because datetime objects in Python are strange and confusing, and
186
	strftime doesn't support fractional returns, only the full integer microsecond field
187
	which is an integer right-padded with zeroes. This function uses the ``precision``
188
	of a datetime object.
189
190
	For example:
191
192
	.. code-block:: python
193
194
		>>> from obspy import UTCDateTime
195
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
196
		>>> fsec(ti)
197
		UTCDateTime(2020, 1, 1, 0, 0, 0, 600000)
198
199
	:param ti: time object to convert microseconds for
200
	:type ti: obspy.core.utcdatetime.UTCDateTime
201
	:return: the hundredth-of-a-second rounded version of the time object passed (precision is 0.01 second)
202
	:rtype: obspy.core.utcdatetime.UTCDateTime
203
	'''
204
	# time in python is weird and confusing, but luckily obspy is better than Python
205
	# at dealing with datetimes. all we need to do is tell it what precision we want
206
	# and it handles the rounding for us.
207 1
	return rs.UTCDateTime(ti, precision=2)
208
209
210 1
def conn_stats(TESTING=False):
211
	'''
212
	Print some stats about the connection.
213
214
	Example:
215
216
	.. code-block:: python
217
218
		>>> conn_stats()
219
		2020-03-25 01:35:04 [conn_stats] Initialization stats:
220
		2020-03-25 01:35:04 [conn_stats]                 Port: 18069
221
		2020-03-25 01:35:04 [conn_stats]   Sending IP address: 192.168.0.4
222
		2020-03-25 01:35:04 [conn_stats]     Set station name: R24FA
223
		2020-03-25 01:35:04 [conn_stats]   Number of channels: 4
224
		2020-03-25 01:35:04 [conn_stats]   Transmission freq.: 250 ms/packet
225
		2020-03-25 01:35:04 [conn_stats]    Transmission rate: 4 packets/sec
226
		2020-03-25 01:35:04 [conn_stats]   Samples per second: 100 sps
227
		2020-03-25 01:35:04 [conn_stats]            Inventory: AM.R24FA (Raspberry Shake Citizen Science Station)
228
229
	:param bool TESTING: if ``True``, text is printed to the console in yellow. if not, in white.
230
	'''
231 1
	s = 'conn_stats'
232 1
	pf = printW if TESTING else printM
233 1
	pf('Initialization stats:', sender=s, announce=False)
234 1
	pf('                Port: %s' % rs.port, sender=s, announce=False)
235 1
	pf('  Sending IP address: %s' % rs.firstaddr, sender=s, announce=False)
236 1
	pf('    Set station name: %s' % rs.stn, sender=s, announce=False)
237 1
	pf('  Number of channels: %s' % rs.numchns, sender=s, announce=False)
238 1
	pf('  Transmission freq.: %s ms/packet' % rs.tf, sender=s, announce=False)
239 1
	pf('   Transmission rate: %s packets/sec' % rs.tr, sender=s, announce=False)
240 1
	pf('  Samples per second: %s sps' % rs.sps, sender=s, announce=False)
241 1
	if rs.inv:
242 1
		pf('           Inventory: %s' % rs.inv.get_contents()['stations'][0],
243
			   sender=s, announce=False)
244
245
246 1
def msg_alarm(event_time):
247
	'''
248
	This function constructs the ``ALARM`` message as a bytes object.
249
	Currently this is only used by :py:class:`rsudp.p_producer.Producer`
250
	to construct alarm queue messages.
251
252
	For example:
253
254
	.. code-block:: python
255
256
		>>> from obspy import UTCDateTime
257
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
258
		>>> msg_alarm(ti)
259
		b'ALARM 2020-01-01T00:00:00.599Z'
260
261
	:param obspy.core.utcdatetime.UTCDateTime event_time: the datetime object to serialize and convert to bytes
262
	:rtype: bytes
263
	:return: the ``ALARM`` message, ready to be put on the queue
264
	'''
265 1
	return b'ALARM %s' % bytes(str(event_time), 'utf-8')
266
267
268 1
def msg_reset(reset_time):
269
	'''
270
	This function constructs the ``RESET`` message as a bytes object.
271
	Currently this is only used by :py:class:`rsudp.p_producer.Producer`
272
	to construct reset queue messages.
273
274
	For example:
275
276
	.. code-block:: python
277
278
		>>> from obspy import UTCDateTime
279
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
280
		>>> msg_reset(ti)
281
		b'RESET 2020-01-01T00:00:00.599Z'
282
283
	:param obspy.core.utcdatetime.UTCDateTime reset_time: the datetime object to serialize and convert to bytes
284
	:rtype: bytes
285
	:return: the ``RESET`` message, ready to be put on the queue
286
	'''
287 1
	return b'RESET %s' % bytes(str(reset_time), 'utf-8')
288
289
290 1
def msg_imgpath(event_time, figname):
291
	'''
292
	This function constructs the ``IMGPATH`` message as a bytes object.
293
	Currently this is only used by :py:class:`rsudp.c_plot.Plot`
294
	to construct queue messages containing timestamp and saved image path.
295
296
	For example:
297
298
	.. code-block:: python
299
300
		>>> from obspy import UTCDateTime
301
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
302
		>>> path = '/home/pi/rsudp/screenshots/test.png'
303
		>>> msg_imgpath(ti, path)
304
		b'IMGPATH 2020-01-01T00:00:00.599Z /home/pi/rsudp/screenshots/test.png'
305
306
	:param obspy.core.utcdatetime.UTCDateTime event_time: the datetime object to serialize and convert to bytes
307
	:param str figname: the figure path as a string
308
	:rtype: bytes
309
	:return: the ``IMGPATH`` message, ready to be put on the queue
310
	'''
311 1
	return b'IMGPATH %s %s' % (bytes(str(event_time), 'utf-8'), bytes(str(figname), 'utf-8'))
312
313
314 1
def msg_term():
315
	'''
316
	This function constructs the simple ``TERM`` message as a bytes object.
317
318
	.. code-block:: python
319
320
		>>> msg_term()
321
		b'TERM'
322
323
324
	:rtype: bytes
325
	:return: the ``TERM`` message
326
	'''
327 1
	return b'TERM'
328
329
330 1
def get_msg_time(msg):
331
	'''
332
	This function gets the time from ``ALARM``, ``RESET``,
333
	and ``IMGPATH`` messages as a UTCDateTime object.
334
335
	For example:
336
337
	.. code-block:: python
338
339
		>>> from obspy import UTCDateTime
340
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
341
		>>> path = '/home/pi/rsudp/screenshots/test.png'
342
		>>> msg = msg_imgpath(ti, path)
343
		>>> msg
344
		b'IMGPATH 2020-01-01T00:00:00.599Z /home/pi/rsudp/screenshots/test.png'
345
		>>> get_msg_time(msg)
346
		UTCDateTime(2020, 1, 1, 0, 0, 0, 599000)
347
348
	:param bytes msg: the bytes-formatted queue message to decode
349
	:rtype: obspy.core.utcdatetime.UTCDateTime
350
	:return: the time embedded in the message
351
	'''
352 1
	return rs.UTCDateTime.strptime(msg.decode('utf-8').split(' ')[1], '%Y-%m-%dT%H:%M:%S.%fZ')
353
354
355 1
def get_msg_path(msg):
356
	'''
357
	This function gets the path from ``IMGPATH`` messages as a string.
358
359
	For example:
360
361
	.. code-block:: python
362
363
		>>> from obspy import UTCDateTime
364
		>>> ti = UTCDateTime(2020, 1, 1, 0, 0, 0, 599000, precision=3)
365
		>>> path = '/home/pi/rsudp/screenshots/test.png'
366
		>>> msg = msg_imgpath(ti, path)
367
		>>> msg
368
		b'IMGPATH 2020-01-01T00:00:00.599Z /home/pi/rsudp/screenshots/test.png'
369
		>>> get_msg_path(msg)
370
		'/home/pi/rsudp/screenshots/test.png'
371
372
	:param bytes msg: the bytes-formatted queue message to decode
373
	:rtype: str
374
	:return: the path embedded in the message
375
	'''
376 1
	return msg.decode('utf-8').split(' ')[2]
377
378
379 1
def get_scap_dir():
380
	'''
381
	This function returns the screen capture directory from the init function.
382
	This allows the variable to be more threadsafe.
383
384
	.. code-block:: python
385
386
		>>> get_scap_dir()
387
		'/home/pi/rsudp/screenshots/'
388
389
	:return: the path of the screenshot directory
390
	'''
391 1
	return rsudp.scap_dir
392
393
394 1
def deconv_vel_inst(self, trace, output):
395
	'''
396
	.. role:: pycode(code)
397
		:language: python
398
	
399
	A helper function for :py:func:`rsudp.raspberryshake.deconvolve`
400
	for velocity channels.
401
402
	:param self self: The self object of the sub-consumer class calling this function.
403
	:param obspy.core.trace.Trace trace: the trace object instance to deconvolve
404
	'''
405 1
	if self.deconv not in 'CHAN':
406
		trace.remove_response(inventory=rs.inv, pre_filt=[0.1, 0.6, 0.95*self.sps, self.sps],
407
								output=output, water_level=4.5, taper=False)
408
	else:
409 1
		trace.remove_response(inventory=rs.inv, pre_filt=[0.1, 0.6, 0.95*self.sps, self.sps],
410
								output='VEL', water_level=4.5, taper=False)
411 1
	if 'ACC' in self.deconv:
412
		trace.data = rs.np.gradient(trace.data, 1)
413 1
	elif 'GRAV' in self.deconv:
414
		trace.data = rs.np.gradient(trace.data, 1) / rs.g
415
		trace.stats.units = 'Earth gravity'
416 1
	elif 'DISP' in self.deconv:
417
		trace.data = rs.np.cumsum(trace.data)
418
		trace.taper(max_percentage=0.1, side='left', max_length=1)
419
		trace.detrend(type='demean')
420
	else:
421 1
		trace.stats.units = 'Velocity'
422
423
424 1
def deconv_acc_inst(self, trace, output):
425
	'''
426
	.. role:: pycode(code)
427
		:language: python
428
	
429
	A helper function for :py:func:`rsudp.raspberryshake.deconvolve`
430
	for acceleration channels.
431
432
	:param self self: The self object of the sub-consumer class calling this function.
433
	:param obspy.core.trace.Trace trace: the trace object instance to deconvolve
434
	'''
435 1
	if self.deconv not in 'CHAN':
436
		trace.remove_response(inventory=rs.inv, pre_filt=[0.1, 0.6, 0.95*self.sps, self.sps],
437
								output=output, water_level=4.5, taper=False)
438
	else:
439 1
		trace.remove_response(inventory=rs.inv, pre_filt=[0.1, 0.6, 0.95*self.sps, self.sps],
440
								output='ACC', water_level=4.5, taper=False)
441 1
	if 'VEL' in self.deconv:
442
		trace.data = rs.np.cumsum(trace.data)
443
		trace.detrend(type='demean')
444 1
	elif 'DISP' in self.deconv:
445
		trace.data = rs.np.cumsum(rs.np.cumsum(trace.data))
446
		trace.detrend(type='linear')
447 1
	elif 'GRAV' in self.deconv:
448
		trace.data = trace.data / rs.g
449
		trace.stats.units = 'Earth gravity'
450
	else:
451 1
		trace.stats.units = 'Acceleration'
452 1
	if ('ACC' not in self.deconv) and ('CHAN' not in self.deconv):
453
		trace.taper(max_percentage=0.1, side='left', max_length=1)
454
455
456 1
def deconv_rbm_inst(self, trace, output):
457
	'''
458
	.. role:: pycode(code)
459
		:language: python
460
	
461
	A helper function for :py:func:`rsudp.raspberryshake.deconvolve`
462
	for Raspberry Boom pressure transducer channels.
463
464
	.. note::
465
466
		The Raspberry Boom pressure transducer does not currently have a
467
		deconvolution function. The Raspberry Shake team is working on a
468
		calibration for the Boom, but until then Boom units are given in
469
		counts.
470
471
	:param self self: The self object of the sub-consumer class calling this function.
472
	:param obspy.core.trace.Trace trace: the trace object instance to deconvolve
473
	'''
474
	trace.stats.units = ' counts'
475
476
477 1
def deconvolve(self):
478
	'''
479
	.. role:: pycode(code)
480
		:language: python
481
	
482
	A central helper function for sub-consumers (i.e. :py:class:`rsudp.c_plot.Plot` or :py:class:`rsudp.c_alert.Alert`)
483
	that need to deconvolve their raw data to metric units.
484
	Consumers with :py:class:`obspy.core.stream.Stream` objects in :pycode:`self.stream` can use this to deconvolve data
485
	if this library's :pycode:`rsudp.raspberryshake.inv` variable
486
	contains a valid :py:class:`obspy.core.inventory.inventory.Inventory` object.
487
488
	:param self self: The self object of the sub-consumer class calling this function. Must contain :pycode:`self.stream` as a :py:class:`obspy.core.stream.Stream` object.
489
	'''
490 1
	acc_channels = ['ENE', 'ENN', 'ENZ']
491 1
	vel_channels = ['EHE', 'EHN', 'EHZ', 'SHZ']
492 1
	rbm_channels = ['HDF']
493
494 1
	self.stream = self.raw.copy()
495 1
	for trace in self.stream:
496 1
		trace.stats.units = self.units
497 1
		output = 'ACC' if self.deconv == 'GRAV' else self.deconv	# if conversion is to gravity
498 1
		if self.deconv:
499 1
			if trace.stats.channel in vel_channels:
500 1
				deconv_vel_inst(self, trace, output)	# geophone channels
501
502 1
			elif trace.stats.channel in acc_channels:
503 1
				deconv_acc_inst(self, trace, output)	# accelerometer channels
504
505
			elif trace.stats.channel in rbm_channels:
506
				deconv_rbm_inst(self, trace, output)	# this is the Boom channel
507
508
			else:
509
				trace.stats.units = ' counts'	# this is a new one
510
511
		else:
512
			trace.stats.units = ' counts'		# this is not being deconvolved
513
514
515