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