1
|
|
|
import os, sys, platform |
2
|
|
|
import pkg_resources as pr |
3
|
|
|
import time |
4
|
|
|
import math |
5
|
|
|
import numpy as np |
6
|
|
|
from datetime import datetime, timedelta |
7
|
|
|
import rsudp.raspberryshake as rs |
8
|
|
|
from rsudp import printM, printW, printE |
9
|
|
|
import rsudp |
10
|
|
|
import linecache |
11
|
|
|
sender = 'plot.py' |
12
|
|
|
QT = False |
13
|
|
|
QtGui = False |
14
|
|
|
PhotoImage = False |
15
|
|
|
try: # test for matplotlib and exit if import fails |
16
|
|
|
from matplotlib import use |
17
|
|
|
try: # no way to know what machines can handle what software, but Tk is more universal |
18
|
|
|
use('Qt5Agg') # try for Qt because it's better and has less threatening errors |
19
|
|
|
from PyQt5 import QtGui |
20
|
|
|
QT = True |
21
|
|
|
except Exception as e: |
22
|
|
|
printW('Qt import failed. Trying Tk...') |
23
|
|
|
printW('detail: %s' % e, spaces=True) |
24
|
|
|
try: # fail over to the more reliable Tk |
25
|
|
|
use('TkAgg') |
26
|
|
|
from tkinter import PhotoImage |
27
|
|
|
except Exception as e: |
28
|
|
|
printE('Could not import either Qt or Tk, and the plot module requires at least one of them to run.', sender) |
29
|
|
|
printE('Please make sure either PyQt5 or Tkinter is installed.', sender, spaces=True) |
30
|
|
|
printE('detail: %s'% e, sender, spaces=True) |
31
|
|
|
raise ImportError('Could not import either Qt or Tk, and the plot module requires at least one of them to run') |
32
|
|
|
import matplotlib.pyplot as plt |
33
|
|
|
import matplotlib.dates as mdates |
34
|
|
|
import matplotlib.image as mpimg |
35
|
|
|
from matplotlib import rcParams |
36
|
|
|
from matplotlib.ticker import EngFormatter |
37
|
|
|
rcParams['toolbar'] = 'None' |
38
|
|
|
plt.ion() |
39
|
|
|
MPL = True |
40
|
|
|
|
41
|
|
|
# avoiding a matplotlib user warning |
42
|
|
|
import warnings |
43
|
|
|
warnings.filterwarnings('ignore', category=UserWarning, module='rsudp') |
44
|
|
|
|
45
|
|
|
except Exception as e: |
46
|
|
|
printE('Could not import matplotlib, plotting will not be available.', sender) |
47
|
|
|
printE('detail: %s' % e, sender, spaces=True) |
48
|
|
|
MPL = False |
49
|
|
|
|
50
|
|
|
ICON = 'icon.ico' |
51
|
|
|
ICON2 = 'icon.png' |
52
|
|
|
|
53
|
|
|
class Plot: |
54
|
|
|
''' |
55
|
|
|
.. role:: json(code) |
56
|
|
|
:language: json |
57
|
|
|
|
58
|
|
|
GUI plotting algorithm, compatible with both Qt5 and TkAgg backends (see :py:func:`matplotlib.use`). |
59
|
|
|
This module can plot seismogram data from a list of 1-4 Shake channels, and calculate and display a spectrogram beneath each. |
60
|
|
|
|
61
|
|
|
By default the plotted :json:`"duration"` in seconds is :json:`30`. |
62
|
|
|
The plot will refresh at most once per second, but slower processors may take longer. |
63
|
|
|
The longer the duration, the more processor power it will take to refresh the plot, |
64
|
|
|
especially when the spectrogram is enabled. |
65
|
|
|
To disable the spectrogram, set :json:`"spectrogram"` to :json:`false` in the settings file. |
66
|
|
|
To put the plot into fullscreen window mode, set :json:`"fullscreen"` to :json:`true`. |
67
|
|
|
To put the plot into kiosk mode, set :json:`"kiosk"` to :json:`true`. |
68
|
|
|
|
69
|
|
|
:param cha: channels to plot. Defaults to "all" but can be passed a list of channel names as strings. |
70
|
|
|
:type cha: str or list |
71
|
|
|
:param int seconds: number of seconds to plot. Defaults to 30. |
72
|
|
|
:param bool spectrogram: whether to plot the spectrogram. Defaults to True. |
73
|
|
|
:param bool fullscreen: whether to plot in a fullscreen window. Defaults to False. |
74
|
|
|
:param bool kiosk: whether to plot in kiosk mode (true fullscreen). Defaults to False. |
75
|
|
|
:param deconv: whether to deconvolve the signal. Defaults to False. |
76
|
|
|
:type deconv: str or bool |
77
|
|
|
:param bool screencap: whether or not to save screenshots of events. Defaults to False. |
78
|
|
|
:param bool alert: whether to draw the number of events at startup. Defaults to True. |
79
|
|
|
:param queue.Queue q: queue of data and messages sent by :class:`rsudp.c_consumer.Consumer` |
80
|
|
|
:raise ImportError: if the module cannot import either of the Matplotlib Qt5 or TkAgg backends |
81
|
|
|
|
82
|
|
|
''' |
83
|
|
|
|
84
|
|
|
def __init__(self, cha='all', q=False, |
85
|
|
|
seconds=30, spectrogram=True, |
86
|
|
|
fullscreen=False, kiosk=False, |
87
|
|
|
deconv=False, screencap=False, |
88
|
|
|
alert=True): |
89
|
|
|
""" |
90
|
|
|
Initialize the plot process. |
91
|
|
|
|
92
|
|
|
""" |
93
|
|
|
super().__init__() |
94
|
|
|
self.sender = 'Plot' |
95
|
|
|
self.alive = True |
96
|
|
|
self.alarm = False # don't touch this |
97
|
|
|
self.alarm_reset = False # don't touch this |
98
|
|
|
|
99
|
|
|
if MPL == False: |
100
|
|
|
sys.stdout.flush() |
101
|
|
|
sys.exit() |
102
|
|
|
if QT == False: |
103
|
|
|
printW('Running on %s machine, using Tk instead of Qt' % (platform.machine()), self.sender) |
104
|
|
|
if q: |
105
|
|
|
self.queue = q |
106
|
|
|
else: |
107
|
|
|
printE('no queue passed to consumer! Thread will exit now!', self.sender) |
108
|
|
|
sys.stdout.flush() |
109
|
|
|
sys.exit() |
110
|
|
|
|
111
|
|
|
self.master_queue = None # careful with this, this goes directly to the master consumer. gets set by main thread. |
112
|
|
|
|
113
|
|
|
self.stream = rs.Stream() |
114
|
|
|
self.raw = rs.Stream() |
115
|
|
|
self.stn = rs.stn |
116
|
|
|
self.net = rs.net |
117
|
|
|
self.chans = [] |
118
|
|
|
cha = rs.chns if (cha == 'all') else cha |
119
|
|
|
cha = list(cha) if isinstance(cha, str) else cha |
120
|
|
|
l = rs.chns |
121
|
|
|
for c in l: |
122
|
|
|
n = 0 |
123
|
|
|
for uch in cha: |
124
|
|
|
if (uch.upper() in c) and (c not in str(self.chans)): |
125
|
|
|
self.chans.append(c) |
126
|
|
|
n += 1 |
127
|
|
|
if len(self.chans) < 1: |
128
|
|
|
self.chans = rs.chns |
129
|
|
|
printM('Plotting channels: %s' % self.chans, self.sender) |
130
|
|
|
self.totchns = rs.numchns |
131
|
|
|
self.seconds = seconds |
132
|
|
|
self.pkts_in_period = rs.tr * rs.numchns * self.seconds # theoretical number of packets received in self.seconds |
133
|
|
|
self.spectrogram = spectrogram |
134
|
|
|
|
135
|
|
|
self.deconv = deconv if (deconv in rs.UNITS) else False |
136
|
|
|
if self.deconv and rs.inv: |
137
|
|
|
deconv = deconv.upper() |
138
|
|
|
if self.deconv in rs.UNITS: |
139
|
|
|
self.units = rs.UNITS[self.deconv][0] |
140
|
|
|
self.unit = rs.UNITS[self.deconv][1] |
141
|
|
|
printM('Signal deconvolution set to %s' % (self.deconv), self.sender) |
142
|
|
|
else: |
143
|
|
|
self.units = rs.UNITS['CHAN'][0] |
144
|
|
|
self.unit = rs.UNITS['CHAN'][1] |
145
|
|
|
self.deconv = False |
146
|
|
|
printM('Seismogram units are %s' % (self.units), self.sender) |
147
|
|
|
|
148
|
|
|
self.per_lap = 0.9 |
149
|
|
|
self.fullscreen = fullscreen |
150
|
|
|
self.kiosk = kiosk |
151
|
|
|
self.num_chans = len(self.chans) |
152
|
|
|
self.delay = rs.tr if (self.spectrogram) else 1 |
153
|
|
|
self.delay = 0.5 if (self.chans == ['SHZ']) else self.delay |
154
|
|
|
|
155
|
|
|
self.screencap = screencap |
156
|
|
|
self.save_timer = 0 |
157
|
|
|
self.save_pct = 0.7 |
158
|
|
|
self.save = [] |
159
|
|
|
self.events = 0 |
160
|
|
|
self.event_text = ' - detected events: 0' if alert else '' |
161
|
|
|
self.last_event = [] |
162
|
|
|
self.last_event_str = False |
163
|
|
|
# plot stuff |
164
|
|
|
self.bgcolor = '#202530' # background |
165
|
|
|
self.fgcolor = '0.8' # axis and label color |
166
|
|
|
self.linecolor = '#c28285' # seismogram color |
167
|
|
|
|
168
|
|
|
printM('Starting.', self.sender) |
169
|
|
|
|
170
|
|
|
def deconvolve(self): |
171
|
|
|
''' |
172
|
|
|
Send the streams to the central library deconvolve function. |
173
|
|
|
''' |
174
|
|
|
rs.deconvolve(self) |
175
|
|
|
|
176
|
|
|
def getq(self): |
177
|
|
|
''' |
178
|
|
|
Get data from the queue and test for whether it has certain strings. |
179
|
|
|
ALARM and TERM both trigger specific behavior. |
180
|
|
|
ALARM messages cause the event counter to increment, and if |
181
|
|
|
:py:data:`screencap==True` then aplot image will be saved when the |
182
|
|
|
event is :py:data:`self.save_pct` of the way across the plot. |
183
|
|
|
''' |
184
|
|
|
d = self.queue.get() |
185
|
|
|
self.queue.task_done() |
186
|
|
|
if 'TERM' in str(d): |
187
|
|
|
plt.close() |
188
|
|
|
if 'SELF' in str(d): |
189
|
|
|
printM('Plot has been closed, plot thread will exit.', self.sender) |
190
|
|
|
self.alive = False |
191
|
|
|
rs.producer = False |
192
|
|
|
|
193
|
|
|
elif 'ALARM' in str(d): |
194
|
|
|
self.events += 1 # add event to count |
195
|
|
|
self.save_timer -= 1 # don't push the save time forward if there are a large number of alarm events |
196
|
|
|
event = [self.save_timer + int(self.save_pct*self.pkts_in_period), |
197
|
|
|
rs.fsec(rs.get_msg_time(d))] # event = [save after count, datetime] |
198
|
|
|
self.last_event_str = '%s UTC' % (event[1].strftime('%Y-%m-%d %H:%M:%S.%f')[:22]) |
199
|
|
|
printM('Event time: %s' % (self.last_event_str), sender=self.sender) # show event time in the logs |
200
|
|
|
if self.screencap: |
201
|
|
|
printM('Saving png in about %i seconds' % (self.save_pct * (self.seconds)), self.sender) |
202
|
|
|
self.save.append(event) # append |
203
|
|
|
self.fig.suptitle('%s.%s live output - detected events: %s' # title |
204
|
|
|
% (self.net, self.stn, self.events), |
205
|
|
|
fontsize=14, color=self.fgcolor, x=0.52) |
206
|
|
|
self.fig.canvas.set_window_title('(%s) %s.%s - Raspberry Shake Monitor' % (self.events, self.net, self.stn)) |
207
|
|
|
|
208
|
|
|
if rs.getCHN(d) in self.chans: |
209
|
|
|
self.raw = rs.update_stream( |
210
|
|
|
stream=self.raw, d=d, fill_value='latest') |
211
|
|
|
return True |
212
|
|
|
else: |
213
|
|
|
return False |
214
|
|
|
|
215
|
|
|
def set_sps(self): |
216
|
|
|
''' |
217
|
|
|
Get samples per second from the main library. |
218
|
|
|
''' |
219
|
|
|
self.sps = rs.sps |
220
|
|
|
|
221
|
|
|
# from https://docs.obspy.org/_modules/obspy/imaging/spectrogram.html#_nearest_pow_2: |
222
|
|
|
def _nearest_pow_2(self, x): |
223
|
|
|
""" |
224
|
|
|
Find power of two nearest to x |
225
|
|
|
|
226
|
|
|
>>> _nearest_pow_2(3) |
227
|
|
|
2.0 |
228
|
|
|
>>> _nearest_pow_2(15) |
229
|
|
|
16.0 |
230
|
|
|
|
231
|
|
|
:type x: float |
232
|
|
|
:param x: Number |
233
|
|
|
:rtype: Int |
234
|
|
|
:return: Nearest power of 2 to x |
235
|
|
|
|
236
|
|
|
Adapted from the `obspy <https://obspy.org>`_ library |
237
|
|
|
""" |
238
|
|
|
a = math.pow(2, math.ceil(np.log2(x))) |
239
|
|
|
b = math.pow(2, math.floor(np.log2(x))) |
240
|
|
|
if abs(a - x) < abs(b - x): |
241
|
|
|
return a |
242
|
|
|
else: |
243
|
|
|
return b |
244
|
|
|
|
245
|
|
|
def handle_close(self, evt): |
246
|
|
|
''' |
247
|
|
|
Handles a plot close event. |
248
|
|
|
This will trigger a full shutdown of all other processes related to rsudp. |
249
|
|
|
''' |
250
|
|
|
self.master_queue.put(rs.msg_term()) |
251
|
|
|
|
252
|
|
|
def handle_resize(self, evt=False): |
253
|
|
|
''' |
254
|
|
|
Handles a plot window resize event. |
255
|
|
|
This will allow the plot to resize dynamically. |
256
|
|
|
''' |
257
|
|
|
if evt: |
258
|
|
|
h = evt.height |
259
|
|
|
else: |
260
|
|
|
h = self.fig.get_size_inches()[1]*self.fig.dpi |
261
|
|
|
plt.tight_layout(pad=0, h_pad=0.1, w_pad=0, |
262
|
|
|
rect=[0.02, 0.01, 0.98, 0.90 + 0.045*(h/1080)]) # [left, bottom, right, top] |
263
|
|
|
|
264
|
|
|
def _eventsave(self): |
265
|
|
|
''' |
266
|
|
|
This function takes the next event in line and pops it out of the list, |
267
|
|
|
so that it can be saved and others preserved. |
268
|
|
|
Then, it sets the title to something having to do with the event, |
269
|
|
|
then calls the save figure function, and finally resets the title. |
270
|
|
|
''' |
271
|
|
|
self.save.reverse() |
272
|
|
|
event = self.save.pop() |
273
|
|
|
self.save.reverse() |
274
|
|
|
|
275
|
|
|
event_time_str = event[1].strftime('%Y-%m-%d-%H%M%S') # event time for filename |
276
|
|
|
title_time_str = event[1].strftime('%Y-%m-%d %H:%M:%S.%f')[:22] # pretty event time for plot |
277
|
|
|
|
278
|
|
|
# change title (just for a moment) |
279
|
|
|
self.fig.suptitle('%s.%s detected event - %s UTC' # title |
280
|
|
|
% (self.net, self.stn, title_time_str), |
281
|
|
|
fontsize=14, color=self.fgcolor, x=0.52) |
282
|
|
|
|
283
|
|
|
# save figure |
284
|
|
|
self.savefig(event_time=event[1], event_time_str=event_time_str) |
285
|
|
|
|
286
|
|
|
# reset title |
287
|
|
|
self._set_fig_title() |
288
|
|
|
|
289
|
|
|
|
290
|
|
|
def savefig(self, event_time=rs.UTCDateTime.now(), |
291
|
|
|
event_time_str=rs.UTCDateTime.now().strftime('%Y-%m-%d-%H%M%S')): |
292
|
|
|
''' |
293
|
|
|
Saves the figure and puts an IMGPATH message on the master queue. |
294
|
|
|
This message can be used to upload the image to various services. |
295
|
|
|
|
296
|
|
|
:param obspy.core.utcdatetime.UTCDateTime event_time: Event time as an obspy UTCDateTime object. |
297
|
|
|
:param str event_time_str: Event time as a string. This is used to set the filename. |
298
|
|
|
''' |
299
|
|
|
figname = os.path.join(rsudp.scap_dir, '%s-%s.png' % (self.stn, event_time_str)) |
300
|
|
|
elapsed = rs.UTCDateTime.now() - event_time |
301
|
|
|
if int(elapsed) > 0: |
302
|
|
|
printM('Saving png %i seconds after alarm' % (elapsed), sender=self.sender) |
303
|
|
|
plt.savefig(figname, facecolor=self.fig.get_facecolor(), edgecolor='none') |
304
|
|
|
printM('Saved %s' % (figname), sender=self.sender) |
305
|
|
|
printM('%s thread has saved an image, sending IMGPATH message to queues' % self.sender, sender=self.sender) |
306
|
|
|
# imgpath requires a UTCDateTime and a string figure path |
307
|
|
|
self.master_queue.put(rs.msg_imgpath(event_time, figname)) |
308
|
|
|
|
309
|
|
|
|
310
|
|
|
def _set_fig_title(self): |
311
|
|
|
''' |
312
|
|
|
Sets the figure title back to something that makes sense for the live viewer. |
313
|
|
|
''' |
314
|
|
|
self.fig.suptitle('%s.%s live output - detected events: %s' # title |
315
|
|
|
% (self.net, self.stn, self.events), |
316
|
|
|
fontsize=14, color=self.fgcolor, x=0.52) |
317
|
|
|
|
318
|
|
|
|
319
|
|
|
def _init_plot(self): |
320
|
|
|
''' |
321
|
|
|
Initialize plot elements and calculate parameters. |
322
|
|
|
''' |
323
|
|
|
self.fig = plt.figure(figsize=(11,3*self.num_chans)) |
324
|
|
|
self.fig.canvas.mpl_connect('close_event', self.handle_close) |
325
|
|
|
self.fig.canvas.mpl_connect('resize_event', self.handle_resize) |
326
|
|
|
|
327
|
|
|
if QT: |
328
|
|
|
self.fig.canvas.window().statusBar().setVisible(False) # remove bottom bar |
329
|
|
|
self.fig.canvas.set_window_title('%s.%s - Raspberry Shake Monitor' % (self.net, self.stn)) |
330
|
|
|
self.fig.patch.set_facecolor(self.bgcolor) # background color |
331
|
|
|
self.fig.suptitle('%s.%s live output%s' # title |
332
|
|
|
% (self.net, self.stn, self.event_text), |
333
|
|
|
fontsize=14, color=self.fgcolor,x=0.52) |
334
|
|
|
self.ax, self.lines = [], [] # list for subplot axes and lines artists |
335
|
|
|
self.mult = 1 # spectrogram selection multiplier |
336
|
|
|
if self.spectrogram: |
337
|
|
|
self.mult = 2 # 2 if user wants a spectrogram else 1 |
338
|
|
|
if self.seconds > 60: |
339
|
|
|
self.per_lap = 0.9 # if axis is long, spectrogram overlap can be shorter |
340
|
|
|
else: |
341
|
|
|
self.per_lap = 0.975 # if axis is short, increase resolution |
342
|
|
|
# set spectrogram parameters |
343
|
|
|
self.nfft1 = self._nearest_pow_2(self.sps) |
344
|
|
|
self.nlap1 = self.nfft1 * self.per_lap |
345
|
|
|
|
346
|
|
|
|
347
|
|
|
def _init_axes(self, i): |
348
|
|
|
''' |
349
|
|
|
Initialize plot axes. |
350
|
|
|
''' |
351
|
|
|
if i == 0: |
352
|
|
|
# set up first axes (axes added later will share these x axis limits) |
353
|
|
|
self.ax.append(self.fig.add_subplot(self.num_chans*self.mult, |
354
|
|
|
1, 1, label=str(1))) |
355
|
|
|
self.ax[0].set_facecolor(self.bgcolor) |
356
|
|
|
self.ax[0].tick_params(colors=self.fgcolor, labelcolor=self.fgcolor) |
357
|
|
|
self.ax[0].xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S')) |
358
|
|
|
self.ax[0].yaxis.set_major_formatter(EngFormatter(unit='%s' % self.unit.lower())) |
359
|
|
|
if self.spectrogram: |
360
|
|
|
self.ax.append(self.fig.add_subplot(self.num_chans*self.mult, |
361
|
|
|
1, 2, label=str(2)))#, sharex=ax[0])) |
362
|
|
|
self.ax[1].set_facecolor(self.bgcolor) |
363
|
|
|
self.ax[1].tick_params(colors=self.fgcolor, labelcolor=self.fgcolor) |
364
|
|
|
else: |
365
|
|
|
# add axes that share either lines or spectrogram axis limits |
366
|
|
|
s = i * self.mult # plot selector |
367
|
|
|
# add a subplot then set colors |
368
|
|
|
self.ax.append(self.fig.add_subplot(self.num_chans*self.mult, |
369
|
|
|
1, s+1, sharex=self.ax[0], label=str(s+1))) |
370
|
|
|
self.ax[s].set_facecolor(self.bgcolor) |
371
|
|
|
self.ax[s].tick_params(colors=self.fgcolor, labelcolor=self.fgcolor) |
372
|
|
|
self.ax[s].xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S')) |
373
|
|
|
self.ax[s].yaxis.set_major_formatter(EngFormatter(unit='%s' % self.unit.lower())) |
374
|
|
|
if self.spectrogram: |
375
|
|
|
# add a spectrogram and set colors |
376
|
|
|
self.ax.append(self.fig.add_subplot(self.num_chans*self.mult, |
377
|
|
|
1, s+2, sharex=self.ax[1], label=str(s+2))) |
378
|
|
|
self.ax[s+1].set_facecolor(self.bgcolor) |
379
|
|
|
self.ax[s+1].tick_params(colors=self.fgcolor, labelcolor=self.fgcolor) |
380
|
|
|
|
381
|
|
|
|
382
|
|
|
def _set_icon(self): |
383
|
|
|
''' |
384
|
|
|
Set RS plot icons. |
385
|
|
|
''' |
386
|
|
|
mgr = plt.get_current_fig_manager() |
387
|
|
|
ico = pr.resource_filename('rsudp', os.path.join('img', ICON)) |
388
|
|
|
if QT: |
389
|
|
|
mgr.window.setWindowIcon(QtGui.QIcon(ico)) |
390
|
|
|
else: |
391
|
|
|
try: |
392
|
|
|
ico = PhotoImage(file=ico) |
393
|
|
|
mgr.window.tk.call('wm', 'iconphoto', mgr.window._w, ico) |
394
|
|
|
except: |
395
|
|
|
printW('Failed to set PNG icon image, trying .ico instead', sender=self.sender) |
396
|
|
|
try: |
397
|
|
|
ico = pr.resource_filename('rsudp', os.path.join('img', ICON2)) |
398
|
|
|
ico = PhotoImage(file=ico) |
399
|
|
|
mgr.window.tk.call('wm', 'iconphoto', mgr.window._w, ico) |
400
|
|
|
except: |
401
|
|
|
printE('Failed to set window icon.') |
402
|
|
|
|
403
|
|
|
|
404
|
|
|
def _format_axes(self): |
405
|
|
|
''' |
406
|
|
|
Setting up axes and artists. |
407
|
|
|
''' |
408
|
|
|
# calculate times |
409
|
|
|
start = np.datetime64(self.stream[0].stats.endtime |
410
|
|
|
)-np.timedelta64(self.seconds, 's') # numpy time |
411
|
|
|
end = np.datetime64(self.stream[0].stats.endtime) # numpy time |
412
|
|
|
|
413
|
|
|
im = mpimg.imread(pr.resource_filename('rsudp', os.path.join('img', 'version1-01-small.png'))) |
414
|
|
|
self.imax = self.fig.add_axes([0.015, 0.944, 0.2, 0.056], anchor='NW') # [left, bottom, right, top] |
415
|
|
|
self.imax.imshow(im, aspect='equal', interpolation='sinc') |
416
|
|
|
self.imax.axis('off') |
417
|
|
|
# set up axes and artists |
418
|
|
|
for i in range(self.num_chans): # create lines objects and modify axes |
419
|
|
|
if len(self.stream[i].data) < int(self.sps*(1/self.per_lap)): |
420
|
|
|
comp = 0 # spectrogram offset compensation factor |
421
|
|
|
else: |
422
|
|
|
comp = (1/self.per_lap)**2 # spectrogram offset compensation factor |
423
|
|
|
r = np.arange(start, end, np.timedelta64(int(1000/self.sps), 'ms'))[-len( |
424
|
|
|
self.stream[i].data[int(-self.sps*(self.seconds-(comp/2))):-int(self.sps*(comp/2))]):] |
425
|
|
|
mean = int(round(np.mean(self.stream[i].data))) |
426
|
|
|
# add artist to lines list |
427
|
|
|
self.lines.append(self.ax[i*self.mult].plot(r, |
428
|
|
|
np.nan*(np.zeros(len(r))), |
429
|
|
|
label=self.stream[i].stats.channel, color=self.linecolor, |
430
|
|
|
lw=0.45)[0]) |
431
|
|
|
# set axis limits |
432
|
|
|
self.ax[i*self.mult].set_xlim(left=start.astype(datetime), |
433
|
|
|
right=end.astype(datetime)) |
434
|
|
|
self.ax[i*self.mult].set_ylim(bottom=np.min(self.stream[i].data-mean) |
435
|
|
|
-np.ptp(self.stream[i].data-mean)*0.1, |
436
|
|
|
top=np.max(self.stream[i].data-mean) |
437
|
|
|
+np.ptp(self.stream[i].data-mean)*0.1) |
438
|
|
|
# we can set line plot labels here, but not imshow labels |
439
|
|
|
ylabel = self.stream[i].stats.units.strip().capitalize() if (' ' in self.stream[i].stats.units) else self.stream[i].stats.units |
440
|
|
|
self.ax[i*self.mult].set_ylabel(ylabel, color=self.fgcolor) |
441
|
|
|
self.ax[i*self.mult].legend(loc='upper left') # legend and location |
442
|
|
|
if self.spectrogram: # if the user wants a spectrogram, plot it |
443
|
|
|
# add spectrogram to axes list |
444
|
|
|
sg = self.ax[1].specgram(self.stream[i].data, NFFT=8, pad_to=8, |
445
|
|
|
Fs=self.sps, noverlap=7, cmap='inferno', |
446
|
|
|
xextent=(self.seconds-0.5, self.seconds))[0] |
447
|
|
|
self.ax[1].set_xlim(0,self.seconds) |
448
|
|
|
self.ax[i*self.mult+1].set_ylim(0,int(self.sps/2)) |
449
|
|
|
self.ax[i*self.mult+1].imshow(np.flipud(sg**(1/float(10))), cmap='inferno', |
450
|
|
|
extent=(self.seconds-(1/(self.sps/float(len(self.stream[i].data)))), |
451
|
|
|
self.seconds,0,self.sps/2), aspect='auto') |
452
|
|
|
|
453
|
|
|
|
454
|
|
|
def _setup_fig_manager(self): |
455
|
|
|
''' |
456
|
|
|
Setting up figure manager and |
457
|
|
|
''' |
458
|
|
|
# update canvas and draw |
459
|
|
|
figManager = plt.get_current_fig_manager() |
460
|
|
|
if self.kiosk: |
461
|
|
|
figManager.full_screen_toggle() |
462
|
|
|
else: |
463
|
|
|
if self.fullscreen: # set fullscreen |
464
|
|
|
if QT: # maximizing in Qt |
465
|
|
|
figManager.window.showMaximized() |
466
|
|
|
else: # maximizing in Tk |
467
|
|
|
figManager.resize(*figManager.window.maxsize()) |
468
|
|
|
|
469
|
|
|
|
470
|
|
|
def setup_plot(self): |
471
|
|
|
""" |
472
|
|
|
Sets up the plot. Quite a lot of stuff happens in this function. |
473
|
|
|
Matplotlib backends are not threadsafe, so things are a little weird. |
474
|
|
|
See code comments for details. |
475
|
|
|
""" |
476
|
|
|
# instantiate a figure and set basic params |
477
|
|
|
self._init_plot() |
478
|
|
|
|
479
|
|
|
for i in range(self.num_chans): |
480
|
|
|
self._init_axes(i) |
481
|
|
|
|
482
|
|
|
for axis in self.ax: |
483
|
|
|
# set the rest of plot colors |
484
|
|
|
plt.setp(axis.spines.values(), color=self.fgcolor) |
485
|
|
|
plt.setp([axis.get_xticklines(), axis.get_yticklines()], color=self.fgcolor) |
486
|
|
|
|
487
|
|
|
# rs logos |
488
|
|
|
self._set_icon() |
489
|
|
|
|
490
|
|
|
# draw axes |
491
|
|
|
self._format_axes() |
492
|
|
|
|
493
|
|
|
self.handle_resize() |
494
|
|
|
|
495
|
|
|
# setup figure manager |
496
|
|
|
self._setup_fig_manager() |
497
|
|
|
|
498
|
|
|
# draw plot, loop, and resize the plot |
499
|
|
|
plt.draw() # draw the canvas |
500
|
|
|
self.fig.canvas.start_event_loop(0.005) # wait for canvas to update |
501
|
|
|
self.handle_resize() |
502
|
|
|
|
503
|
|
|
|
504
|
|
|
def _set_ch_specific_label(self, i): |
505
|
|
|
''' |
506
|
|
|
Set the formatter units if the deconvolution is channel-specific. |
507
|
|
|
''' |
508
|
|
|
if self.deconv: |
509
|
|
|
if (self.deconv in 'CHAN'): |
510
|
|
|
ch = self.stream[i].stats.channel |
511
|
|
|
if ('HZ' in ch) or ('HN' in ch) or ('HE' in ch): |
512
|
|
|
unit = rs.UNITS['VEL'][1] |
513
|
|
|
elif ('EN' in ch): |
514
|
|
|
unit = rs.UNITS['ACC'][1] |
515
|
|
|
else: |
516
|
|
|
unit = rs.UNITS['CHAN'][1] |
517
|
|
|
self.ax[i*self.mult].yaxis.set_major_formatter(EngFormatter(unit='%s' % unit.lower())) |
518
|
|
|
|
519
|
|
|
|
520
|
|
|
def _draw_lines(self, i, start, end, mean): |
521
|
|
|
''' |
522
|
|
|
Updates the line data in the plot. |
523
|
|
|
|
524
|
|
|
:param int i: the trace number |
525
|
|
|
:param numpy.datetime64 start: start time of the trace |
526
|
|
|
:param numpy.datetime64 end: end time of the trace |
527
|
|
|
:param float mean: the mean of data in the trace |
528
|
|
|
''' |
529
|
|
|
comp = 1/self.per_lap # spectrogram offset compensation factor |
530
|
|
|
r = np.arange(start, end, np.timedelta64(int(1000/self.sps), 'ms'))[-len( |
531
|
|
|
self.stream[i].data[int(-self.sps*(self.seconds-(comp/2))):-int(self.sps*(comp/2))]):] |
532
|
|
|
self.lines[i].set_ydata(self.stream[i].data[int(-self.sps*(self.seconds-(comp/2))):-int(self.sps*(comp/2))]-mean) |
533
|
|
|
self.lines[i].set_xdata(r) # (1/self.per_lap)/2 |
534
|
|
|
self.ax[i*self.mult].set_xlim(left=start.astype(datetime)+timedelta(seconds=comp*1.5), |
535
|
|
|
right=end.astype(datetime)) |
536
|
|
|
self.ax[i*self.mult].set_ylim(bottom=np.min(self.stream[i].data-mean) |
537
|
|
|
-np.ptp(self.stream[i].data-mean)*0.1, |
538
|
|
|
top=np.max(self.stream[i].data-mean) |
539
|
|
|
+np.ptp(self.stream[i].data-mean)*0.1) |
540
|
|
|
|
541
|
|
|
|
542
|
|
|
def _update_specgram(self, i, mean): |
543
|
|
|
''' |
544
|
|
|
Updates the spectrogram and its labels. |
545
|
|
|
|
546
|
|
|
:param int i: the trace number |
547
|
|
|
:param float mean: the mean of data in the trace |
548
|
|
|
''' |
549
|
|
|
self.nfft1 = self._nearest_pow_2(self.sps) # FFTs run much faster if the number of transforms is a power of 2 |
550
|
|
|
self.nlap1 = self.nfft1 * self.per_lap |
551
|
|
|
if len(self.stream[i].data) < self.nfft1: # when the number of data points is low, we just need to kind of fake it for a few fractions of a second |
552
|
|
|
self.nfft1 = 8 |
553
|
|
|
self.nlap1 = 6 |
554
|
|
|
sg = self.ax[i*self.mult+1].specgram(self.stream[i].data-mean, |
555
|
|
|
NFFT=self.nfft1, pad_to=int(self.nfft1*4), # previously self.sps*4), |
556
|
|
|
Fs=self.sps, noverlap=self.nlap1)[0] # meat & potatoes |
557
|
|
|
self.ax[i*self.mult+1].clear() # incredibly important, otherwise continues to draw over old images (gets exponentially slower) |
558
|
|
|
# cloogy way to shift the spectrogram to line up with the seismogram |
559
|
|
|
self.ax[i*self.mult+1].set_xlim(0.25,self.seconds-0.25) |
560
|
|
|
self.ax[i*self.mult+1].set_ylim(0,int(self.sps/2)) |
561
|
|
|
# imshow to update the spectrogram |
562
|
|
|
self.ax[i*self.mult+1].imshow(np.flipud(sg**(1/float(10))), cmap='inferno', |
563
|
|
|
extent=(self.seconds-(1/(self.sps/float(len(self.stream[i].data)))), |
564
|
|
|
self.seconds,0,self.sps/2), aspect='auto') |
565
|
|
|
# some things that unfortunately can't be in the setup function: |
566
|
|
|
self.ax[i*self.mult+1].tick_params(axis='x', which='both', |
567
|
|
|
bottom=False, top=False, labelbottom=False) |
568
|
|
|
self.ax[i*self.mult+1].set_ylabel('Frequency (Hz)', color=self.fgcolor) |
569
|
|
|
self.ax[i*self.mult+1].set_xlabel('Time (UTC)', color=self.fgcolor) |
570
|
|
|
|
571
|
|
|
|
572
|
|
|
def update_plot(self): |
573
|
|
|
''' |
574
|
|
|
Redraw the plot with new data. |
575
|
|
|
Called on every nth loop after the plot is set up, where n is |
576
|
|
|
the number of channels times the data packet arrival rate in Hz. |
577
|
|
|
This has the effect of making the plot update once per second. |
578
|
|
|
''' |
579
|
|
|
obstart = self.stream[0].stats.endtime - timedelta(seconds=self.seconds) # obspy time |
580
|
|
|
start = np.datetime64(self.stream[0].stats.endtime |
581
|
|
|
)-np.timedelta64(self.seconds, 's') # numpy time |
582
|
|
|
end = np.datetime64(self.stream[0].stats.endtime) # numpy time |
583
|
|
|
self.raw = self.raw.slice(starttime=obstart) # slice the stream to the specified length (seconds variable) |
584
|
|
|
self.stream = self.stream.slice(starttime=obstart) # slice the stream to the specified length (seconds variable) |
585
|
|
|
i = 0 |
586
|
|
|
for i in range(self.num_chans): # for each channel, update the plots |
587
|
|
|
mean = int(round(np.mean(self.stream[i].data))) |
588
|
|
|
self._draw_lines(i, start, end, mean) |
589
|
|
|
self._set_ch_specific_label(i) |
590
|
|
|
if self.spectrogram: |
591
|
|
|
self._update_specgram(i, mean) |
592
|
|
|
else: |
593
|
|
|
# also can't be in the setup function |
594
|
|
|
self.ax[i*self.mult].set_xlabel('Time (UTC)', color=self.fgcolor) |
595
|
|
|
|
596
|
|
|
|
597
|
|
|
def figloop(self): |
598
|
|
|
""" |
599
|
|
|
Let some time elapse in order for the plot canvas to draw properly. |
600
|
|
|
Must be separate from :py:func:`update_plot()` to avoid a broadcast error early in plotting. |
601
|
|
|
""" |
602
|
|
|
self.fig.canvas.start_event_loop(0.005) |
603
|
|
|
|
604
|
|
|
|
605
|
|
|
def mainloop(self, i, u): |
606
|
|
|
''' |
607
|
|
|
The main loop in the :py:func:`rsudp.c_plot.Plot.run`. |
608
|
|
|
|
609
|
|
|
:param int i: number of plot events without clearing the linecache |
610
|
|
|
:param int u: queue blocking counter |
611
|
|
|
:return: number of plot events without clearing the linecache and queue blocking counter |
612
|
|
|
:rtype: int, int |
613
|
|
|
''' |
614
|
|
|
if i > 10: |
615
|
|
|
linecache.clearcache() |
616
|
|
|
i = 0 |
617
|
|
|
else: |
618
|
|
|
i += 1 |
619
|
|
|
self.stream = rs.copy(self.stream) # essential, otherwise the stream has a memory leak |
620
|
|
|
self.raw = rs.copy(self.raw) # and could eventually crash the machine |
621
|
|
|
self.deconvolve() |
622
|
|
|
self.update_plot() |
623
|
|
|
if u >= 0: # avoiding a matplotlib broadcast error |
624
|
|
|
self.figloop() |
625
|
|
|
|
626
|
|
|
if self.save: |
627
|
|
|
# save the plot |
628
|
|
|
if (self.save_timer > self.save[0][0]): |
629
|
|
|
self._eventsave() |
630
|
|
|
u = 0 |
631
|
|
|
time.sleep(0.005) # wait a ms to see if another packet will arrive |
632
|
|
|
sys.stdout.flush() |
633
|
|
|
return i, u |
634
|
|
|
|
635
|
|
|
def qu(self, u): |
636
|
|
|
''' |
637
|
|
|
Get a queue object and increment the queue counter. |
638
|
|
|
This is a way to figure out how many channels have arrived in the queue. |
639
|
|
|
|
640
|
|
|
:param int u: queue blocking counter |
641
|
|
|
:return: queue blocking counter |
642
|
|
|
:rtype: int |
643
|
|
|
''' |
644
|
|
|
u += 1 if self.getq() else 0 |
645
|
|
|
return u |
646
|
|
|
|
647
|
|
|
|
648
|
|
|
def run(self): |
649
|
|
|
""" |
650
|
|
|
The heart of the plotting routine. |
651
|
|
|
|
652
|
|
|
Begins by updating the queue to populate a :py:class:`obspy.core.stream.Stream` object, then setting up the main plot. |
653
|
|
|
The first time through the main loop, the plot is not drawn. After that, the plot is drawn every time all channels are updated. |
654
|
|
|
Any plots containing a spectrogram and more than 1 channel are drawn at most every second (1000 ms). |
655
|
|
|
All other plots are drawn at most every quarter second (250 ms). |
656
|
|
|
""" |
657
|
|
|
self.getq() # block until data is flowing from the consumer |
658
|
|
|
for i in range((self.totchns)*2): # fill up a stream object |
659
|
|
|
self.getq() |
660
|
|
|
self.set_sps() |
661
|
|
|
self.deconvolve() |
662
|
|
|
self.setup_plot() |
663
|
|
|
|
664
|
|
|
n = 0 # number of iterations without plotting |
665
|
|
|
i = 0 # number of plot events without clearing the linecache |
666
|
|
|
u = -1 # number of blocked queue calls (must be -1 at startup) |
667
|
|
|
while True: # main loop |
668
|
|
|
while True: # sub loop |
669
|
|
|
if self.alive == False: # break if the user has closed the plot |
670
|
|
|
break |
671
|
|
|
n += 1 |
672
|
|
|
self.save_timer += 1 |
673
|
|
|
if self.queue.qsize() > 0: |
674
|
|
|
self.getq() |
675
|
|
|
time.sleep(0.009) # wait a ms to see if another packet will arrive |
676
|
|
|
else: |
677
|
|
|
u = self.qu(u) |
678
|
|
|
if n > (self.delay * rs.numchns): |
679
|
|
|
n = 0 |
680
|
|
|
break |
681
|
|
|
if self.alive == False: # break if the user has closed the plot |
682
|
|
|
printM('Exiting.', self.sender) |
683
|
|
|
break |
684
|
|
|
i, u = self.mainloop(i, u) |
685
|
|
|
return |
686
|
|
|
|