##// END OF EJS Templates
schain
RSS
plot cities in SkyMap plot
jespinoza -
r1144:74565fa4d2f5
parent child
Show More
Side by Side Unified
Context file:
r1144:74565fa4d2f5 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,1085 +1,1088
1 1
2 2 import os
3 3 import time
4 4 import glob
5 5 import datetime
6 6 from multiprocessing import Process
7 7
8 8 import zmq
9 9 import numpy
10 10 import matplotlib
11 11 import matplotlib.pyplot as plt
12 12 from mpl_toolkits.axes_grid1 import make_axes_locatable
13 13 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
14 14
15 15 from schainpy.model.proc.jroproc_base import Operation
16 16 from schainpy.utils import log
17 17
18 18 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
19 19 blu_values = matplotlib.pyplot.get_cmap(
20 20 'seismic_r', 20)(numpy.arange(20))[10:15]
21 21 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
22 22 'jro', numpy.vstack((blu_values, jet_values)))
23 23 matplotlib.pyplot.register_cmap(cmap=ncmap)
24 24
25 25 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis', 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm', 'spectral')]
26 26
27 def ll2xy(lat1, lon1, lat2, lon2):
28
29 p = 0.017453292519943295
30 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
31 r = 12742 * numpy.arcsin(numpy.sqrt(a))
32 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)*numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
33 theta = -theta + numpy.pi/2
34 return r*numpy.cos(theta), r*numpy.sin(theta)
27 35
28 36 def figpause(interval):
29 37 backend = plt.rcParams['backend']
30 38 if backend in matplotlib.rcsetup.interactive_bk:
31 39 figManager = matplotlib._pylab_helpers.Gcf.get_active()
32 40 if figManager is not None:
33 41 canvas = figManager.canvas
34 42 if canvas.figure.stale:
35 43 canvas.draw()
36 44 canvas.start_event_loop(interval)
37 45 return
38 46
39 47
40 48 class PlotData(Operation, Process):
41 49 '''
42 50 Base class for Schain plotting operations
43 51 '''
44 52
45 53 CODE = 'Figure'
46 54 colormap = 'jro'
47 55 bgcolor = 'white'
48 56 CONFLATE = False
49 57 __missing = 1E30
50 58
51 59 __attrs__ = ['show', 'save', 'xmin', 'xmax', 'ymin', 'ymax', 'zmin', 'zmax',
52 60 'zlimits', 'xlabel', 'ylabel', 'xaxis','cb_label', 'title',
53 61 'colorbar', 'bgcolor', 'width', 'height', 'localtime', 'oneFigure',
54 62 'showprofile', 'decimation', 'ftp']
55 63
56 64 def __init__(self, **kwargs):
57 65
58 66 Operation.__init__(self, plot=True, **kwargs)
59 67 Process.__init__(self)
60 68
61 69 self.kwargs['code'] = self.CODE
62 70 self.mp = False
63 71 self.data = None
64 72 self.isConfig = False
65 73 self.figures = []
66 74 self.axes = []
67 75 self.cb_axes = []
68 76 self.localtime = kwargs.pop('localtime', True)
69 77 self.show = kwargs.get('show', True)
70 78 self.save = kwargs.get('save', False)
71 79 self.ftp = kwargs.get('ftp', False)
72 80 self.colormap = kwargs.get('colormap', self.colormap)
73 81 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
74 82 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
75 83 self.colormaps = kwargs.get('colormaps', None)
76 84 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
77 85 self.showprofile = kwargs.get('showprofile', False)
78 86 self.title = kwargs.get('wintitle', self.CODE.upper())
79 87 self.cb_label = kwargs.get('cb_label', None)
80 88 self.cb_labels = kwargs.get('cb_labels', None)
81 89 self.labels = kwargs.get('labels', None)
82 90 self.xaxis = kwargs.get('xaxis', 'frequency')
83 91 self.zmin = kwargs.get('zmin', None)
84 92 self.zmax = kwargs.get('zmax', None)
85 93 self.zlimits = kwargs.get('zlimits', None)
86 94 self.xmin = kwargs.get('xmin', None)
87 95 self.xmax = kwargs.get('xmax', None)
88 96 self.xrange = kwargs.get('xrange', 24)
89 97 self.xscale = kwargs.get('xscale', None)
90 98 self.ymin = kwargs.get('ymin', None)
91 99 self.ymax = kwargs.get('ymax', None)
92 100 self.yscale = kwargs.get('yscale', None)
93 101 self.xlabel = kwargs.get('xlabel', None)
94 102 self.decimation = kwargs.get('decimation', None)
95 103 self.showSNR = kwargs.get('showSNR', False)
96 104 self.oneFigure = kwargs.get('oneFigure', True)
97 105 self.width = kwargs.get('width', None)
98 106 self.height = kwargs.get('height', None)
99 107 self.colorbar = kwargs.get('colorbar', True)
100 108 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
101 109 self.channels = kwargs.get('channels', None)
102 110 self.titles = kwargs.get('titles', [])
103 111 self.polar = False
104 112
105 113 def __fmtTime(self, x, pos):
106 114 '''
107 115 '''
108 116
109 117 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
110 118
111 119 def __setup(self):
112 120 '''
113 121 Common setup for all figures, here figures and axes are created
114 122 '''
115 123
116 124 if self.CODE not in self.data:
117 125 raise ValueError(log.error('Missing data for {}'.format(self.CODE),
118 126 self.name))
119 127
120 128 self.setup()
121 129
122 130 self.time_label = 'LT' if self.localtime else 'UTC'
123 131 if self.data.localtime:
124 132 self.getDateTime = datetime.datetime.fromtimestamp
125 133 else:
126 134 self.getDateTime = datetime.datetime.utcfromtimestamp
127 135
128 136 if self.width is None:
129 137 self.width = 8
130 138
131 139 self.figures = []
132 140 self.axes = []
133 141 self.cb_axes = []
134 142 self.pf_axes = []
135 143 self.cmaps = []
136 144
137 145 size = '15%' if self.ncols == 1 else '30%'
138 146 pad = '4%' if self.ncols == 1 else '8%'
139 147
140 148 if self.oneFigure:
141 149 if self.height is None:
142 150 self.height = 1.4 * self.nrows + 1
143 151 fig = plt.figure(figsize=(self.width, self.height),
144 152 edgecolor='k',
145 153 facecolor='w')
146 154 self.figures.append(fig)
147 155 for n in range(self.nplots):
148 156 ax = fig.add_subplot(self.nrows, self.ncols,
149 157 n + 1, polar=self.polar)
150 158 ax.tick_params(labelsize=8)
151 159 ax.firsttime = True
152 160 ax.index = 0
153 161 ax.press = None
154 162 self.axes.append(ax)
155 163 if self.showprofile:
156 164 cax = self.__add_axes(ax, size=size, pad=pad)
157 165 cax.tick_params(labelsize=8)
158 166 self.pf_axes.append(cax)
159 167 else:
160 168 if self.height is None:
161 169 self.height = 3
162 170 for n in range(self.nplots):
163 171 fig = plt.figure(figsize=(self.width, self.height),
164 172 edgecolor='k',
165 173 facecolor='w')
166 174 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
167 175 ax.tick_params(labelsize=8)
168 176 ax.firsttime = True
169 177 ax.index = 0
170 178 ax.press = None
171 179 self.figures.append(fig)
172 180 self.axes.append(ax)
173 181 if self.showprofile:
174 182 cax = self.__add_axes(ax, size=size, pad=pad)
175 183 cax.tick_params(labelsize=8)
176 184 self.pf_axes.append(cax)
177 185
178 186 for n in range(self.nrows):
179 187 if self.colormaps is not None:
180 188 cmap = plt.get_cmap(self.colormaps[n])
181 189 else:
182 190 cmap = plt.get_cmap(self.colormap)
183 191 cmap.set_bad(self.bgcolor, 1.)
184 192 self.cmaps.append(cmap)
185 193
186 194 for fig in self.figures:
187 195 fig.canvas.mpl_connect('key_press_event', self.OnKeyPress)
188 196 fig.canvas.mpl_connect('scroll_event', self.OnBtnScroll)
189 197 fig.canvas.mpl_connect('button_press_event', self.onBtnPress)
190 198 fig.canvas.mpl_connect('motion_notify_event', self.onMotion)
191 199 fig.canvas.mpl_connect('button_release_event', self.onBtnRelease)
192 200 if self.show:
193 201 fig.show()
194 202
195 203 def OnKeyPress(self, event):
196 204 '''
197 205 Event for pressing keys (up, down) change colormap
198 206 '''
199 207 ax = event.inaxes
200 208 if ax in self.axes:
201 209 if event.key == 'down':
202 210 ax.index += 1
203 211 elif event.key == 'up':
204 212 ax.index -= 1
205 213 if ax.index < 0:
206 214 ax.index = len(CMAPS) - 1
207 215 elif ax.index == len(CMAPS):
208 216 ax.index = 0
209 217 cmap = CMAPS[ax.index]
210 218 ax.cbar.set_cmap(cmap)
211 219 ax.cbar.draw_all()
212 220 ax.plt.set_cmap(cmap)
213 221 ax.cbar.patch.figure.canvas.draw()
214 222 self.colormap = cmap.name
215 223
216 224 def OnBtnScroll(self, event):
217 225 '''
218 226 Event for scrolling, scale figure
219 227 '''
220 228 cb_ax = event.inaxes
221 229 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
222 230 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
223 231 pt = ax.cbar.ax.bbox.get_points()[:, 1]
224 232 nrm = ax.cbar.norm
225 233 vmin, vmax, p0, p1, pS = (
226 234 nrm.vmin, nrm.vmax, pt[0], pt[1], event.y)
227 235 scale = 2 if event.step == 1 else 0.5
228 236 point = vmin + (vmax - vmin) / (p1 - p0) * (pS - p0)
229 237 ax.cbar.norm.vmin = point - scale * (point - vmin)
230 238 ax.cbar.norm.vmax = point - scale * (point - vmax)
231 239 ax.plt.set_norm(ax.cbar.norm)
232 240 ax.cbar.draw_all()
233 241 ax.cbar.patch.figure.canvas.draw()
234 242
235 243 def onBtnPress(self, event):
236 244 '''
237 245 Event for mouse button press
238 246 '''
239 247 cb_ax = event.inaxes
240 248 if cb_ax is None:
241 249 return
242 250
243 251 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
244 252 cb_ax.press = event.x, event.y
245 253 else:
246 254 cb_ax.press = None
247 255
248 256 def onMotion(self, event):
249 257 '''
250 258 Event for move inside colorbar
251 259 '''
252 260 cb_ax = event.inaxes
253 261 if cb_ax is None:
254 262 return
255 263 if cb_ax not in [ax.cbar.ax for ax in self.axes if ax.cbar]:
256 264 return
257 265 if cb_ax.press is None:
258 266 return
259 267
260 268 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
261 269 xprev, yprev = cb_ax.press
262 270 dx = event.x - xprev
263 271 dy = event.y - yprev
264 272 cb_ax.press = event.x, event.y
265 273 scale = ax.cbar.norm.vmax - ax.cbar.norm.vmin
266 274 perc = 0.03
267 275
268 276 if event.button == 1:
269 277 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
270 278 ax.cbar.norm.vmax -= (perc * scale) * numpy.sign(dy)
271 279 elif event.button == 3:
272 280 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
273 281 ax.cbar.norm.vmax += (perc * scale) * numpy.sign(dy)
274 282
275 283 ax.cbar.draw_all()
276 284 ax.plt.set_norm(ax.cbar.norm)
277 285 ax.cbar.patch.figure.canvas.draw()
278 286
279 287 def onBtnRelease(self, event):
280 288 '''
281 289 Event for mouse button release
282 290 '''
283 291 cb_ax = event.inaxes
284 292 if cb_ax is not None:
285 293 cb_ax.press = None
286 294
287 295 def __add_axes(self, ax, size='30%', pad='8%'):
288 296 '''
289 297 Add new axes to the given figure
290 298 '''
291 299 divider = make_axes_locatable(ax)
292 300 nax = divider.new_horizontal(size=size, pad=pad)
293 301 ax.figure.add_axes(nax)
294 302 return nax
295 303
296 304 self.setup()
297 305
298 306 def setup(self):
299 307 '''
300 308 This method should be implemented in the child class, the following
301 309 attributes should be set:
302 310
303 311 self.nrows: number of rows
304 312 self.ncols: number of cols
305 313 self.nplots: number of plots (channels or pairs)
306 314 self.ylabel: label for Y axes
307 315 self.titles: list of axes title
308 316
309 317 '''
310 318 raise(NotImplementedError, 'Implement this method in child class')
311 319
312 320 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
313 321 '''
314 322 Create a masked array for missing data
315 323 '''
316 324 if x_buffer.shape[0] < 2:
317 325 return x_buffer, y_buffer, z_buffer
318 326
319 327 deltas = x_buffer[1:] - x_buffer[0:-1]
320 328 x_median = numpy.median(deltas)
321 329
322 330 index = numpy.where(deltas > 5 * x_median)
323 331
324 332 if len(index[0]) != 0:
325 333 z_buffer[::, index[0], ::] = self.__missing
326 334 z_buffer = numpy.ma.masked_inside(z_buffer,
327 335 0.99 * self.__missing,
328 336 1.01 * self.__missing)
329 337
330 338 return x_buffer, y_buffer, z_buffer
331 339
332 340 def decimate(self):
333 341
334 342 # dx = int(len(self.x)/self.__MAXNUMX) + 1
335 343 dy = int(len(self.y) / self.decimation) + 1
336 344
337 345 # x = self.x[::dx]
338 346 x = self.x
339 347 y = self.y[::dy]
340 348 z = self.z[::, ::, ::dy]
341 349
342 350 return x, y, z
343 351
344 352 def format(self):
345 353 '''
346 354 Set min and max values, labels, ticks and titles
347 355 '''
348 356
349 357 if self.xmin is None:
350 358 xmin = self.min_time
351 359 else:
352 360 if self.xaxis is 'time':
353 361 dt = self.getDateTime(self.min_time)
354 362 xmin = (dt.replace(hour=int(self.xmin), minute=0, second=0) -
355 363 datetime.datetime(1970, 1, 1)).total_seconds()
356 364 if self.data.localtime:
357 365 xmin += time.timezone
358 366 else:
359 367 xmin = self.xmin
360 368
361 369 if self.xmax is None:
362 370 xmax = xmin + self.xrange * 60 * 60
363 371 else:
364 372 if self.xaxis is 'time':
365 373 dt = self.getDateTime(self.max_time)
366 374 xmax = (dt.replace(hour=int(self.xmax), minute=59, second=59) -
367 375 datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds()
368 376 if self.data.localtime:
369 377 xmax += time.timezone
370 378 else:
371 379 xmax = self.xmax
372 380
373 381 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
374 382 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
375 383
376 384 Y = numpy.array([5, 10, 20, 50, 100, 200, 500, 1000, 2000, 5000])
377 385 i = 1 if numpy.where(abs(ymax-ymin) <= Y)[0][0] < 0 else numpy.where(abs(ymax-ymin) <= Y)[0][0]
378 386 ystep = Y[i] / 5.
379 387
380 388 for n, ax in enumerate(self.axes):
381 389 if ax.firsttime:
382 390 ax.set_facecolor(self.bgcolor)
383 391 ax.yaxis.set_major_locator(MultipleLocator(ystep))
384 392 ax.xaxis.set_major_locator(MultipleLocator(ystep))
385 393 if self.xscale:
386 394 ax.xaxis.set_major_formatter(FuncFormatter(lambda x, pos: '{0:g}'.format(x*self.xscale)))
387 395 if self.xscale:
388 396 ax.yaxis.set_major_formatter(FuncFormatter(lambda x, pos: '{0:g}'.format(x*self.yscale)))
389 397 if self.xaxis is 'time':
390 398 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
391 399 ax.xaxis.set_major_locator(LinearLocator(9))
392 400 if self.xlabel is not None:
393 401 ax.set_xlabel(self.xlabel)
394 402 ax.set_ylabel(self.ylabel)
395 403 ax.firsttime = False
396 404 if self.showprofile:
397 405 self.pf_axes[n].set_ylim(ymin, ymax)
398 406 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
399 407 self.pf_axes[n].set_xlabel('dB')
400 408 self.pf_axes[n].grid(b=True, axis='x')
401 409 [tick.set_visible(False)
402 410 for tick in self.pf_axes[n].get_yticklabels()]
403 411 if self.colorbar:
404 412 ax.cbar = plt.colorbar(
405 413 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
406 414 ax.cbar.ax.tick_params(labelsize=8)
407 415 ax.cbar.ax.press = None
408 416 if self.cb_label:
409 417 ax.cbar.set_label(self.cb_label, size=8)
410 418 elif self.cb_labels:
411 419 ax.cbar.set_label(self.cb_labels[n], size=8)
412 420 else:
413 421 ax.cbar = None
414 422
415 423 if not self.polar:
416 424 ax.set_xlim(xmin, xmax)
417 425 ax.set_ylim(ymin, ymax)
418 426 ax.set_title('{} {} {}'.format(
419 427 self.titles[n],
420 428 self.getDateTime(self.max_time).strftime('%Y-%m-%dT%H:%M:%S'),
421 429 self.time_label),
422 430 size=8)
423 431 else:
424 432 ax.set_title('{}'.format(self.titles[n]), size=8)
425 433 ax.set_ylim(0, 90)
426 434 ax.set_yticks(numpy.arange(0, 90, 20))
427 435 ax.yaxis.labelpad = 40
428 436
429 437 def __plot(self):
430 438 '''
431 439 '''
432 440 log.log('Plotting', self.name)
433 441
434 442 try:
435 443 self.plot()
436 444 self.format()
437 445 except Exception as e:
438 446 log.warning('{} Plot could not be updated... check data'.format(self.CODE), self.name)
439 447 log.error(str(e), '')
440 448 return
441 449
442 450 for n, fig in enumerate(self.figures):
443 451 if self.nrows == 0 or self.nplots == 0:
444 452 log.warning('No data', self.name)
445 453 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
446 454 fig.canvas.manager.set_window_title(self.CODE)
447 455 continue
448 456
449 457 fig.tight_layout()
450 458 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
451 459 self.getDateTime(self.max_time).strftime('%Y/%m/%d')))
452 460 fig.canvas.draw()
453 461
454 462 if self.save and (self.data.ended or not self.data.buffering):
455 463
456 464 if self.save_labels:
457 465 labels = self.save_labels
458 466 else:
459 467 labels = range(self.nrows)
460 468
461 469 if self.oneFigure:
462 470 label = ''
463 471 else:
464 472 label = '-{}'.format(labels[n])
465 473 figname = os.path.join(
466 474 self.save,
467 475 self.CODE,
468 476 '{}{}_{}.png'.format(
469 477 self.CODE,
470 478 label,
471 479 self.getDateTime(self.saveTime).strftime(
472 480 '%Y%m%d_%H%M%S'),
473 481 )
474 482 )
475 483 log.log('Saving figure: {}'.format(figname), self.name)
476 484 if not os.path.isdir(os.path.dirname(figname)):
477 485 os.makedirs(os.path.dirname(figname))
478 486 fig.savefig(figname)
479 487
480 488 def plot(self):
481 489 '''
482 490 '''
483 491 raise(NotImplementedError, 'Implement this method in child class')
484 492
485 493 def run(self):
486 494
487 495 log.log('Starting', self.name)
488 496
489 497 context = zmq.Context()
490 498 receiver = context.socket(zmq.SUB)
491 499 receiver.setsockopt(zmq.SUBSCRIBE, '')
492 500 receiver.setsockopt(zmq.CONFLATE, self.CONFLATE)
493 501
494 502 if 'server' in self.kwargs['parent']:
495 503 receiver.connect(
496 504 'ipc:///tmp/{}.plots'.format(self.kwargs['parent']['server']))
497 505 else:
498 506 receiver.connect("ipc:///tmp/zmq.plots")
499 507
500 508 while True:
501 509 try:
502 510 self.data = receiver.recv_pyobj(flags=zmq.NOBLOCK)
503 511 if self.data.localtime and self.localtime:
504 512 self.times = self.data.times
505 513 elif self.data.localtime and not self.localtime:
506 514 self.times = self.data.times + time.timezone
507 515 elif not self.data.localtime and self.localtime:
508 516 self.times = self.data.times - time.timezone
509 517 else:
510 518 self.times = self.data.times
511 519
512 520 self.min_time = self.times[0]
513 521 self.max_time = self.times[-1]
514 522
515 523 if self.isConfig is False:
516 524 self.__setup()
517 525 self.isConfig = True
518 526
519 527 self.__plot()
520 528
521 529 except zmq.Again as e:
522 530 # log.log('.', tag='', nl=False)
523 531 if self.data:
524 532 figpause(self.data.throttle)
525 533 else:
526 534 time.sleep(2)
527 535
528 536 def close(self):
529 537 if self.data:
530 538 self.__plot()
531 539
532 540
533 541 class PlotSpectraData(PlotData):
534 542 '''
535 543 Plot for Spectra data
536 544 '''
537 545
538 546 CODE = 'spc'
539 547 colormap = 'jro'
540 548
541 549 def setup(self):
542 550 self.nplots = len(self.data.channels)
543 551 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
544 552 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
545 553 self.width = 3.4 * self.ncols
546 554 self.height = 3 * self.nrows
547 555 self.cb_label = 'dB'
548 556 if self.showprofile:
549 557 self.width += 0.8 * self.ncols
550 558
551 559 self.ylabel = 'Range [km]'
552 560
553 561 def plot(self):
554 562 if self.xaxis == "frequency":
555 563 x = self.data.xrange[0]
556 564 self.xlabel = "Frequency (kHz)"
557 565 elif self.xaxis == "time":
558 566 x = self.data.xrange[1]
559 567 self.xlabel = "Time (ms)"
560 568 else:
561 569 x = self.data.xrange[2]
562 570 self.xlabel = "Velocity (m/s)"
563 571
564 572 if self.CODE == 'spc_mean':
565 573 x = self.data.xrange[2]
566 574 self.xlabel = "Velocity (m/s)"
567 575
568 576 self.titles = []
569 577
570 578 y = self.data.heights
571 579 self.y = y
572 580 z = self.data['spc']
573 581
574 582 for n, ax in enumerate(self.axes):
575 583 noise = self.data['noise'][n][-1]
576 584 if self.CODE == 'spc_mean':
577 585 mean = self.data['mean'][n][-1]
578 586 if ax.firsttime:
579 587 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
580 588 self.xmin = self.xmin if self.xmin else -self.xmax
581 589 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
582 590 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
583 591 ax.plt = ax.pcolormesh(x, y, z[n].T,
584 592 vmin=self.zmin,
585 593 vmax=self.zmax,
586 594 cmap=plt.get_cmap(self.colormap)
587 595 )
588 596
589 597 if self.showprofile:
590 598 ax.plt_profile = self.pf_axes[n].plot(
591 599 self.data['rti'][n][-1], y)[0]
592 600 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
593 601 color="k", linestyle="dashed", lw=1)[0]
594 602 if self.CODE == 'spc_mean':
595 603 ax.plt_mean = ax.plot(mean, y, color='k')[0]
596 604 else:
597 605 ax.plt.set_array(z[n].T.ravel())
598 606 if self.showprofile:
599 607 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
600 608 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
601 609 if self.CODE == 'spc_mean':
602 610 ax.plt_mean.set_data(mean, y)
603 611
604 612 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
605 613 self.saveTime = self.max_time
606 614
607 615
608 616 class PlotCrossSpectraData(PlotData):
609 617
610 618 CODE = 'cspc'
611 619 zmin_coh = None
612 620 zmax_coh = None
613 621 zmin_phase = None
614 622 zmax_phase = None
615 623
616 624 def setup(self):
617 625
618 626 self.ncols = 4
619 627 self.nrows = len(self.data.pairs)
620 628 self.nplots = self.nrows * 4
621 629 self.width = 3.4 * self.ncols
622 630 self.height = 3 * self.nrows
623 631 self.ylabel = 'Range [km]'
624 632 self.showprofile = False
625 633
626 634 def plot(self):
627 635
628 636 if self.xaxis == "frequency":
629 637 x = self.data.xrange[0]
630 638 self.xlabel = "Frequency (kHz)"
631 639 elif self.xaxis == "time":
632 640 x = self.data.xrange[1]
633 641 self.xlabel = "Time (ms)"
634 642 else:
635 643 x = self.data.xrange[2]
636 644 self.xlabel = "Velocity (m/s)"
637 645
638 646 self.titles = []
639 647
640 648 y = self.data.heights
641 649 self.y = y
642 650 spc = self.data['spc']
643 651 cspc = self.data['cspc']
644 652
645 653 for n in range(self.nrows):
646 654 noise = self.data['noise'][n][-1]
647 655 pair = self.data.pairs[n]
648 656 ax = self.axes[4 * n]
649 657 ax3 = self.axes[4 * n + 3]
650 658 if ax.firsttime:
651 659 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
652 660 self.xmin = self.xmin if self.xmin else -self.xmax
653 661 self.zmin = self.zmin if self.zmin else numpy.nanmin(spc)
654 662 self.zmax = self.zmax if self.zmax else numpy.nanmax(spc)
655 663 ax.plt = ax.pcolormesh(x, y, spc[pair[0]].T,
656 664 vmin=self.zmin,
657 665 vmax=self.zmax,
658 666 cmap=plt.get_cmap(self.colormap)
659 667 )
660 668 else:
661 669 ax.plt.set_array(spc[pair[0]].T.ravel())
662 670 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
663 671
664 672 ax = self.axes[4 * n + 1]
665 673 if ax.firsttime:
666 674 ax.plt = ax.pcolormesh(x, y, spc[pair[1]].T,
667 675 vmin=self.zmin,
668 676 vmax=self.zmax,
669 677 cmap=plt.get_cmap(self.colormap)
670 678 )
671 679 else:
672 680 ax.plt.set_array(spc[pair[1]].T.ravel())
673 681 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
674 682
675 683 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
676 684 coh = numpy.abs(out)
677 685 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
678 686
679 687 ax = self.axes[4 * n + 2]
680 688 if ax.firsttime:
681 689 ax.plt = ax.pcolormesh(x, y, coh.T,
682 690 vmin=0,
683 691 vmax=1,
684 692 cmap=plt.get_cmap(self.colormap_coh)
685 693 )
686 694 else:
687 695 ax.plt.set_array(coh.T.ravel())
688 696 self.titles.append(
689 697 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
690 698
691 699 ax = self.axes[4 * n + 3]
692 700 if ax.firsttime:
693 701 ax.plt = ax.pcolormesh(x, y, phase.T,
694 702 vmin=-180,
695 703 vmax=180,
696 704 cmap=plt.get_cmap(self.colormap_phase)
697 705 )
698 706 else:
699 707 ax.plt.set_array(phase.T.ravel())
700 708 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
701 709
702 710 self.saveTime = self.max_time
703 711
704 712
705 713 class PlotSpectraMeanData(PlotSpectraData):
706 714 '''
707 715 Plot for Spectra and Mean
708 716 '''
709 717 CODE = 'spc_mean'
710 718 colormap = 'jro'
711 719
712 720
713 721 class PlotRTIData(PlotData):
714 722 '''
715 723 Plot for RTI data
716 724 '''
717 725
718 726 CODE = 'rti'
719 727 colormap = 'jro'
720 728
721 729 def setup(self):
722 730 self.xaxis = 'time'
723 731 self.ncols = 1
724 732 self.nrows = len(self.data.channels)
725 733 self.nplots = len(self.data.channels)
726 734 self.ylabel = 'Range [km]'
727 735 self.cb_label = 'dB'
728 736 self.titles = ['{} Channel {}'.format(
729 737 self.CODE.upper(), x) for x in range(self.nrows)]
730 738
731 739 def plot(self):
732 740 self.x = self.times
733 741 self.y = self.data.heights
734 742 self.z = self.data[self.CODE]
735 743 self.z = numpy.ma.masked_invalid(self.z)
736 744
737 745 if self.decimation is None:
738 746 x, y, z = self.fill_gaps(self.x, self.y, self.z)
739 747 else:
740 748 x, y, z = self.fill_gaps(*self.decimate())
741 749
742 750 for n, ax in enumerate(self.axes):
743 751 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
744 752 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
745 753 if ax.firsttime:
746 754 ax.plt = ax.pcolormesh(x, y, z[n].T,
747 755 vmin=self.zmin,
748 756 vmax=self.zmax,
749 757 cmap=plt.get_cmap(self.colormap)
750 758 )
751 759 if self.showprofile:
752 760 ax.plot_profile = self.pf_axes[n].plot(
753 761 self.data['rti'][n][-1], self.y)[0]
754 762 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
755 763 color="k", linestyle="dashed", lw=1)[0]
756 764 else:
757 765 ax.collections.remove(ax.collections[0])
758 766 ax.plt = ax.pcolormesh(x, y, z[n].T,
759 767 vmin=self.zmin,
760 768 vmax=self.zmax,
761 769 cmap=plt.get_cmap(self.colormap)
762 770 )
763 771 if self.showprofile:
764 772 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
765 773 ax.plot_noise.set_data(numpy.repeat(
766 774 self.data['noise'][n][-1], len(self.y)), self.y)
767 775
768 776 self.saveTime = self.min_time
769 777
770 778
771 779 class PlotCOHData(PlotRTIData):
772 780 '''
773 781 Plot for Coherence data
774 782 '''
775 783
776 784 CODE = 'coh'
777 785
778 786 def setup(self):
779 787 self.xaxis = 'time'
780 788 self.ncols = 1
781 789 self.nrows = len(self.data.pairs)
782 790 self.nplots = len(self.data.pairs)
783 791 self.ylabel = 'Range [km]'
784 792 if self.CODE == 'coh':
785 793 self.cb_label = ''
786 794 self.titles = [
787 795 'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
788 796 else:
789 797 self.cb_label = 'Degrees'
790 798 self.titles = [
791 799 'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
792 800
793 801
794 802 class PlotPHASEData(PlotCOHData):
795 803 '''
796 804 Plot for Phase map data
797 805 '''
798 806
799 807 CODE = 'phase'
800 808 colormap = 'seismic'
801 809
802 810
803 811 class PlotNoiseData(PlotData):
804 812 '''
805 813 Plot for noise
806 814 '''
807 815
808 816 CODE = 'noise'
809 817
810 818 def setup(self):
811 819 self.xaxis = 'time'
812 820 self.ncols = 1
813 821 self.nrows = 1
814 822 self.nplots = 1
815 823 self.ylabel = 'Intensity [dB]'
816 824 self.titles = ['Noise']
817 825 self.colorbar = False
818 826
819 827 def plot(self):
820 828
821 829 x = self.times
822 830 xmin = self.min_time
823 831 xmax = xmin + self.xrange * 60 * 60
824 832 Y = self.data[self.CODE]
825 833
826 834 if self.axes[0].firsttime:
827 835 for ch in self.data.channels:
828 836 y = Y[ch]
829 837 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
830 838 plt.legend()
831 839 else:
832 840 for ch in self.data.channels:
833 841 y = Y[ch]
834 842 self.axes[0].lines[ch].set_data(x, y)
835 843
836 844 self.ymin = numpy.nanmin(Y) - 5
837 845 self.ymax = numpy.nanmax(Y) + 5
838 846 self.saveTime = self.min_time
839 847
840 848
841 849 class PlotSNRData(PlotRTIData):
842 850 '''
843 851 Plot for SNR Data
844 852 '''
845 853
846 854 CODE = 'snr'
847 855 colormap = 'jet'
848 856
849 857
850 858 class PlotDOPData(PlotRTIData):
851 859 '''
852 860 Plot for DOPPLER Data
853 861 '''
854 862
855 863 CODE = 'dop'
856 864 colormap = 'jet'
857 865
858 866
859 867 class PlotSkyMapData(PlotData):
860 868 '''
861 869 Plot for meteors detection data
862 870 '''
863 871
864 872 CODE = 'param'
865 873
866 874 def setup(self):
867 875
868 876 self.ncols = 1
869 877 self.nrows = 1
870 878 self.width = 7.2
871 879 self.height = 7.2
872 880 self.nplots = 1
873 881 self.xlabel = 'Zonal Zenith Angle (deg)'
874 882 self.ylabel = 'Meridional Zenith Angle (deg)'
875 883 self.polar = True
876 884 self.ymin = -180
877 885 self.ymax = 180
878 886 self.colorbar = False
879 887
880 888 def plot(self):
881 889
882 890 arrayParameters = numpy.concatenate(self.data['param'])
883 891 error = arrayParameters[:, -1]
884 892 indValid = numpy.where(error == 0)[0]
885 893 finalMeteor = arrayParameters[indValid, :]
886 894 finalAzimuth = finalMeteor[:, 3]
887 895 finalZenith = finalMeteor[:, 4]
888 896
889 897 x = finalAzimuth * numpy.pi / 180
890 898 y = finalZenith
891 899
892 900 ax = self.axes[0]
893 901
894 902 if ax.firsttime:
895 903 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
896 904 else:
897 905 ax.plot.set_data(x, y)
898 906
899 907 dt1 = self.getDateTime(self.min_time).strftime('%y/%m/%d %H:%M:%S')
900 908 dt2 = self.getDateTime(self.max_time).strftime('%y/%m/%d %H:%M:%S')
901 909 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
902 910 dt2,
903 911 len(x))
904 912 self.titles[0] = title
905 913 self.saveTime = self.max_time
906 914
907 915
908 916 class PlotParamData(PlotRTIData):
909 917 '''
910 918 Plot for data_param object
911 919 '''
912 920
913 921 CODE = 'param'
914 922 colormap = 'seismic'
915 923
916 924 def setup(self):
917 925 self.xaxis = 'time'
918 926 self.ncols = 1
919 927 self.nrows = self.data.shape(self.CODE)[0]
920 928 self.nplots = self.nrows
921 929 if self.showSNR:
922 930 self.nrows += 1
923 931 self.nplots += 1
924 932
925 933 self.ylabel = 'Height [km]'
926 934 if not self.titles:
927 935 self.titles = self.data.parameters \
928 936 if self.data.parameters else ['Param {}'.format(x) for x in xrange(self.nrows)]
929 937 if self.showSNR:
930 938 self.titles.append('SNR')
931 939
932 940 def plot(self):
933 941 self.data.normalize_heights()
934 942 self.x = self.times
935 943 self.y = self.data.heights
936 944 if self.showSNR:
937 945 self.z = numpy.concatenate(
938 946 (self.data[self.CODE], self.data['snr'])
939 947 )
940 948 else:
941 949 self.z = self.data[self.CODE]
942 950
943 951 self.z = numpy.ma.masked_invalid(self.z)
944 952
945 953 if self.decimation is None:
946 954 x, y, z = self.fill_gaps(self.x, self.y, self.z)
947 955 else:
948 956 x, y, z = self.fill_gaps(*self.decimate())
949 957
950 958 for n, ax in enumerate(self.axes):
951 959
952 960 self.zmax = self.zmax if self.zmax is not None else numpy.max(
953 961 self.z[n])
954 962 self.zmin = self.zmin if self.zmin is not None else numpy.min(
955 963 self.z[n])
956 964
957 965 if ax.firsttime:
958 966 if self.zlimits is not None:
959 967 self.zmin, self.zmax = self.zlimits[n]
960 968
961 969 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
962 970 vmin=self.zmin,
963 971 vmax=self.zmax,
964 972 cmap=self.cmaps[n]
965 973 )
966 974 else:
967 975 if self.zlimits is not None:
968 976 self.zmin, self.zmax = self.zlimits[n]
969 977 ax.collections.remove(ax.collections[0])
970 978 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
971 979 vmin=self.zmin,
972 980 vmax=self.zmax,
973 981 cmap=self.cmaps[n]
974 982 )
975 983
976 984 self.saveTime = self.min_time
977 985
978 986
979 987 class PlotOutputData(PlotParamData):
980 988 '''
981 989 Plot data_output object
982 990 '''
983 991
984 992 CODE = 'output'
985 993 colormap = 'seismic'
986 994
987 995
988 996 class PlotPolarMapData(PlotData):
989 997 '''
990 998 Plot for meteors detection data
991 999 '''
992 1000
993 1001 CODE = 'param'
994 1002 colormap = 'seismic'
995 1003
996 1004 def setup(self):
997 1005 self.ncols = 1
998 1006 self.nrows = 1
999 1007 self.width = 9
1000 1008 self.height = 8
1001 1009 if self.channels is not None:
1002 1010 self.nplots = len(self.channels)
1003 1011 self.nrows = len(self.channels)
1004 1012 else:
1005 1013 self.nplots = self.data.shape(self.CODE)[0]
1006 1014 self.nrows = self.nplots
1007 1015 self.channels = range(self.nplots)
1008 1016 if self.data.meta['mode'] == 'E':
1009 1017 self.xlabel = 'Zonal Distance (km)'
1010 1018 self.ylabel = 'Meridional Distance (km)'
1011 1019 else:
1012 1020 self.xlabel = 'Range (km)'
1013 1021 self.ylabel = 'Height (km)'
1014 1022 self.bgcolor = 'white'
1015 1023 self.cb_labels = self.data.meta['units']
1016 1024 # self.polar = True
1017 1025
1018 1026 def plot(self):
1019 1027
1020 1028 for n, ax in enumerate(self.axes):
1021 1029 data = self.data['param'][self.channels[n]]
1022 1030
1023 1031 zeniths = numpy.linspace(0, self.data.meta['max_range'], data.shape[1])
1024 1032 if self.data.meta['mode'] == 'E':
1025 1033 azimuths = -numpy.radians(self.data.heights)+numpy.pi/2
1026 1034 r, theta = numpy.meshgrid(zeniths, azimuths)
1027 1035 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
1028 1036 else:
1029 1037 azimuths = numpy.radians(self.data.heights)
1030 1038 r, theta = numpy.meshgrid(zeniths, azimuths)
1031 1039 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
1032 1040 self.y = zeniths
1033 1041
1034 1042 if ax.firsttime:
1035 1043 if self.zlimits is not None:
1036 1044 self.zmin, self.zmax = self.zlimits[n]
1037 1045 ax.plt = ax.pcolormesh(#r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
1038 1046 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
1039 1047 vmin=self.zmin,
1040 1048 vmax=self.zmax,
1041 1049 cmap=self.cmaps[n])
1042 1050 else:
1043 1051 if self.zlimits is not None:
1044 1052 self.zmin, self.zmax = self.zlimits[n]
1045 1053 ax.collections.remove(ax.collections[0])
1046 1054 ax.plt = ax.pcolormesh(# r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
1047 1055 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
1048 1056 vmin=self.zmin,
1049 1057 vmax=self.zmax,
1050 1058 cmap=self.cmaps[n])
1051 1059
1052 1060 if self.data.meta['mode'] == 'A':
1053 1061 continue
1054 '''
1055 f = open('/data/workspace/schain_scripts/map_lima.csv')
1062
1063 f = open('/home/jespinoza/workspace/schain_scripts/distrito.csv')
1056 1064
1057 1065 lat1 = -11.96
1058 1066 lon1 = -76.54
1059 1067
1060 1068 for line in f:
1061 label, x, y = [s.strip() for s in line.split(',') if s]
1062 lat2 = float(y)
1063 lon2 = float(x)
1064
1065 dx = (lon2-lon1)*40000*numpy.cos((lat1+lat2)*numpy.pi/360)/360
1066 dy = (lat1-lat2)*40000/360
1067 print label, dx, dy
1068 if label == 'map':
1069 print 'SDHSDHSDHSGHSDFHSDF'
1070 ax.plot([dx], [dy],'--k')
1071 else:
1072 ax.plot([dx], [dy],'.b', ms=2)
1073 '''
1069 label, lon, lat = [s.strip() for s in line.split(',') if s]
1070 lat = float(lat)
1071 lon = float(lon)
1072 x, y = ll2xy(lat1, lon1, lat, lon)
1073 ax.plot(x, y, '.b', ms=2)
1074 ax.text(x, y, label.decode('utf8'), ha='center', va='bottom', size='8', color='black')
1075
1076
1074 1077 if self.data.meta['mode'] == 'E':
1075 1078 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
1076 1079 label = 'E{:d}'.format(int(self.data.meta['elevation']))
1077 1080 else:
1078 1081 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
1079 1082 label = 'A{:d}'.format(int(self.data.meta['azimuth']))
1080 1083
1081 1084 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
1082 1085 self.titles = ['{} {}'.format(self.data.parameters[x], title) for x in self.channels]
1083 1086 self.saveTime = self.max_time
1084 1087
1085 1088
General Comments 0
You need to be logged in to leave comments. Login now