##// END OF EJS Templates
schain
RSS
Add MP150KmRTIPlot and several AverageDriftsPlot classes for testing the 150kM plotting.
imanay -
r1768:008ae31786c1
parent child
Show More
Side by Side Unified
Context file:
r1768:008ae31786c1 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,381 +1,612
1 1 import os
2 2 import datetime
3 3 import numpy
4 4
5 5 from schainpy.model.graphics.jroplot_base import Plot, plt
6 6 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot, SpectraCutPlot
7 7 from schainpy.utils import log
8 8
9 9 EARTH_RADIUS = 6.3710e3
10 10
11 11
12 12 def ll2xy(lat1, lon1, lat2, lon2):
13 13
14 14 p = 0.017453292519943295
15 15 a = 0.5 - numpy.cos((lat2 - lat1) * p) / 2 + numpy.cos(lat1 * p) * \
16 16 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
17 17 r = 12742 * numpy.arcsin(numpy.sqrt(a))
18 18 theta = numpy.arctan2(numpy.sin((lon2 - lon1) * p) * numpy.cos(lat2 * p), numpy.cos(lat1 * p)
19 19 * numpy.sin(lat2 * p) - numpy.sin(lat1 * p) * numpy.cos(lat2 * p) * numpy.cos((lon2 - lon1) * p))
20 20 theta = -theta + numpy.pi / 2
21 21 return r * numpy.cos(theta), r * numpy.sin(theta)
22 22
23 23
24 24 def km2deg(km):
25 25 '''
26 26 Convert distance in km to degrees
27 27 '''
28 28
29 29 return numpy.rad2deg(km / EARTH_RADIUS)
30 30
31 31
32 32
33 33 class SpectralMomentsPlot(SpectraPlot):
34 34 '''
35 35 Plot for Spectral Moments
36 36 '''
37 37 CODE = 'spc_moments'
38 38 # colormap = 'jet'
39 39 # plot_type = 'pcolor'
40 40
41 41 class DobleGaussianPlot(SpectraPlot):
42 42 '''
43 43 Plot for Double Gaussian Plot
44 44 '''
45 45 CODE = 'gaussian_fit'
46 46 # colormap = 'jet'
47 47 # plot_type = 'pcolor'
48 48
49 49
50 50 class DoubleGaussianSpectraCutPlot(SpectraCutPlot):
51 51 '''
52 52 Plot SpectraCut with Double Gaussian Fit
53 53 '''
54 54 CODE = 'cut_gaussian_fit'
55 55
56 56
57 57 class SpectralFitObliquePlot(SpectraPlot):
58 58 '''
59 59 Plot for Spectral Oblique
60 60 '''
61 61 CODE = 'spc_moments'
62 62 colormap = 'jet'
63 63 plot_type = 'pcolor'
64 64
65 65
66 66
67 67 class SnrPlot(RTIPlot):
68 68 '''
69 69 Plot for SNR Data
70 70 '''
71 71
72 72 CODE = 'snr'
73 73 colormap = 'jet'
74 74
75 75 def update(self, dataOut):
76 76
77 77 data = {
78 78 'snr': 10 * numpy.log10(dataOut.data_snr)
79 79 }
80 80
81 81 return data, {}
82 82
83 83 class DopplerPlot(RTIPlot):
84 84 '''
85 85 Plot for DOPPLER Data (1st moment)
86 86 '''
87 87
88 88 CODE = 'dop'
89 89 colormap = 'RdBu_r'
90 90
91 91 def update(self, dataOut):
92 92
93 93 data = {
94 94 'dop': dataOut.data_dop
95 95 }
96 96
97 97 return data, {}
98 98
99 99 class PowerPlot(RTIPlot):
100 100 '''
101 101 Plot for Power Data (0 moment)
102 102 '''
103 103
104 104 CODE = 'pow'
105 105 colormap = 'jet'
106 106
107 107 def update(self, dataOut):
108 108
109 109 data = {
110 110 'pow': 10 * numpy.log10(dataOut.data_pow / dataOut.normFactor)
111 111 }
112 112
113 113 return data, {}
114 114
115 115 class SpectralWidthPlot(RTIPlot):
116 116 '''
117 117 Plot for Spectral Width Data (2nd moment)
118 118 '''
119 119
120 120 CODE = 'width'
121 121 colormap = 'jet'
122 122
123 123 def update(self, dataOut):
124 124
125 125 data = {
126 126 'width': dataOut.data_width
127 127 }
128 128
129 129 return data, {}
130 130
131 131 class SkyMapPlot(Plot):
132 132 '''
133 133 Plot for meteors detection data
134 134 '''
135 135
136 136 CODE = 'param'
137 137
138 138 def setup(self):
139 139
140 140 self.ncols = 1
141 141 self.nrows = 1
142 142 self.width = 7.2
143 143 self.height = 7.2
144 144 self.nplots = 1
145 145 self.xlabel = 'Zonal Zenith Angle (deg)'
146 146 self.ylabel = 'Meridional Zenith Angle (deg)'
147 147 self.polar = True
148 148 self.ymin = -180
149 149 self.ymax = 180
150 150 self.colorbar = False
151 151
152 152 def plot(self):
153 153
154 154 arrayParameters = numpy.concatenate(self.data['param'])
155 155 error = arrayParameters[:, -1]
156 156 indValid = numpy.where(error == 0)[0]
157 157 finalMeteor = arrayParameters[indValid, :]
158 158 finalAzimuth = finalMeteor[:, 3]
159 159 finalZenith = finalMeteor[:, 4]
160 160
161 161 x = finalAzimuth * numpy.pi / 180
162 162 y = finalZenith
163 163
164 164 ax = self.axes[0]
165 165
166 166 if ax.firsttime:
167 167 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
168 168 else:
169 169 ax.plot.set_data(x, y)
170 170
171 171 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
172 172 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
173 173 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
174 174 dt2,
175 175 len(x))
176 176 self.titles[0] = title
177 177
178 178
179 179 class GenericRTIPlot(Plot):
180 180 '''
181 181 Plot for data_xxxx object
182 182 '''
183 183
184 184 CODE = 'param'
185 185 colormap = 'viridis'
186 186 plot_type = 'pcolorbuffer'
187 187
188 188 def setup(self):
189 189 self.xaxis = 'time'
190 190 self.ncols = 1
191 191 self.nrows = self.data.shape('param')[0]
192 192 self.nplots = self.nrows
193 193 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
194 194
195 195 if not self.xlabel:
196 196 self.xlabel = 'Time'
197 197
198 198 self.ylabel = 'Range [km]'
199 199 if not self.titles:
200 200 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
201 201
202 202 def update(self, dataOut):
203 203
204 204 data = {
205 205 'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)
206 206 }
207 207
208 208 meta = {}
209 209
210 210 return data, meta
211 211
212 212 def plot(self):
213 213 # self.data.normalize_heights()
214 214 self.x = self.data.times
215 215 self.y = self.data.yrange
216 216 self.z = self.data['param']
217 217
218 218 self.z = numpy.ma.masked_invalid(self.z)
219 219
220 220 if self.decimation is None:
221 221 x, y, z = self.fill_gaps(self.x, self.y, self.z)
222 222 else:
223 223 x, y, z = self.fill_gaps(*self.decimate())
224 224
225 225 for n, ax in enumerate(self.axes):
226 226
227 227 self.zmax = self.zmax if self.zmax is not None else numpy.max(
228 228 self.z[n])
229 229 self.zmin = self.zmin if self.zmin is not None else numpy.min(
230 230 self.z[n])
231 231
232 232 if ax.firsttime:
233 233 if self.zlimits is not None:
234 234 self.zmin, self.zmax = self.zlimits[n]
235 235
236 236 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
237 237 vmin=self.zmin,
238 238 vmax=self.zmax,
239 239 cmap=self.cmaps[n]
240 240 )
241 241 else:
242 242 if self.zlimits is not None:
243 243 self.zmin, self.zmax = self.zlimits[n]
244 244 ax.plt.remove()
245 245 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
246 246 vmin=self.zmin,
247 247 vmax=self.zmax,
248 248 cmap=self.cmaps[n]
249 249 )
250 250
251 251
252 252 class PolarMapPlot(Plot):
253 253 '''
254 254 Plot for weather radar
255 255 '''
256 256
257 257 CODE = 'param'
258 258 colormap = 'seismic'
259 259
260 260 def setup(self):
261 261 self.ncols = 1
262 262 self.nrows = 1
263 263 self.width = 9
264 264 self.height = 8
265 265 self.mode = self.data.meta['mode']
266 266 if self.channels is not None:
267 267 self.nplots = len(self.channels)
268 268 self.nrows = len(self.channels)
269 269 else:
270 270 self.nplots = self.data.shape(self.CODE)[0]
271 271 self.nrows = self.nplots
272 272 self.channels = list(range(self.nplots))
273 273 if self.mode == 'E':
274 274 self.xlabel = 'Longitude'
275 275 self.ylabel = 'Latitude'
276 276 else:
277 277 self.xlabel = 'Range (km)'
278 278 self.ylabel = 'Height (km)'
279 279 self.bgcolor = 'white'
280 280 self.cb_labels = self.data.meta['units']
281 281 self.lat = self.data.meta['latitude']
282 282 self.lon = self.data.meta['longitude']
283 283 self.xmin, self.xmax = float(
284 284 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
285 285 self.ymin, self.ymax = float(
286 286 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
287 287 # self.polar = True
288 288
289 289 def plot(self):
290 290
291 291 for n, ax in enumerate(self.axes):
292 292 data = self.data['param'][self.channels[n]]
293 293
294 294 zeniths = numpy.linspace(
295 295 0, self.data.meta['max_range'], data.shape[1])
296 296 if self.mode == 'E':
297 297 azimuths = -numpy.radians(self.data.yrange) + numpy.pi / 2
298 298 r, theta = numpy.meshgrid(zeniths, azimuths)
299 299 x, y = r * numpy.cos(theta) * numpy.cos(numpy.radians(self.data.meta['elevation'])), r * numpy.sin(
300 300 theta) * numpy.cos(numpy.radians(self.data.meta['elevation']))
301 301 x = km2deg(x) + self.lon
302 302 y = km2deg(y) + self.lat
303 303 else:
304 304 azimuths = numpy.radians(self.data.yrange)
305 305 r, theta = numpy.meshgrid(zeniths, azimuths)
306 306 x, y = r * numpy.cos(theta), r * numpy.sin(theta)
307 307 self.y = zeniths
308 308
309 309 if ax.firsttime:
310 310 if self.zlimits is not None:
311 311 self.zmin, self.zmax = self.zlimits[n]
312 312 ax.plt = ax.pcolormesh(# r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
313 313 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
314 314 vmin=self.zmin,
315 315 vmax=self.zmax,
316 316 cmap=self.cmaps[n])
317 317 else:
318 318 if self.zlimits is not None:
319 319 self.zmin, self.zmax = self.zlimits[n]
320 320 ax.plt.remove()
321 321 ax.plt = ax.pcolormesh(# r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
322 322 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
323 323 vmin=self.zmin,
324 324 vmax=self.zmax,
325 325 cmap=self.cmaps[n])
326 326
327 327 if self.mode == 'A':
328 328 continue
329 329
330 330 # plot district names
331 331 f = open('/data/workspace/schain_scripts/distrito.csv')
332 332 for line in f:
333 333 label, lon, lat = [s.strip() for s in line.split(',') if s]
334 334 lat = float(lat)
335 335 lon = float(lon)
336 336 # ax.plot(lon, lat, '.b', ms=2)
337 337 ax.text(lon, lat, label.decode('utf8'), ha='center',
338 338 va='bottom', size='8', color='black')
339 339
340 340 # plot limites
341 341 limites = []
342 342 tmp = []
343 343 for line in open('/data/workspace/schain_scripts/lima.csv'):
344 344 if '#' in line:
345 345 if tmp:
346 346 limites.append(tmp)
347 347 tmp = []
348 348 continue
349 349 values = line.strip().split(',')
350 350 tmp.append((float(values[0]), float(values[1])))
351 351 for points in limites:
352 352 ax.add_patch(
353 353 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
354 354
355 355 # plot Cuencas
356 356 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
357 357 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
358 358 values = [line.strip().split(',') for line in f]
359 359 points = [(float(s[0]), float(s[1])) for s in values]
360 360 ax.add_patch(Polygon(points, ec='b', fc='none'))
361 361
362 362 # plot grid
363 363 for r in (15, 30, 45, 60):
364 364 ax.add_artist(plt.Circle((self.lon, self.lat),
365 365 km2deg(r), color='0.6', fill=False, lw=0.2))
366 366 ax.text(
367 367 self.lon + (km2deg(r)) * numpy.cos(60 * numpy.pi / 180),
368 368 self.lat + (km2deg(r)) * numpy.sin(60 * numpy.pi / 180),
369 369 '{}km'.format(r),
370 370 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
371 371
372 372 if self.mode == 'E':
373 373 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
374 374 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
375 375 else:
376 376 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
377 377 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
378 378
379 379 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
380 380 self.titles = ['{} {}'.format(
381 381 self.data.parameters[x], title) for x in self.channels]
382
383 class MP150KmRTIPlot(Plot):
384 '''
385 Plot for data_xxxx object
386 '''
387
388 CODE = 'param'
389 colormap = 'viridis'
390 plot_type = 'pcolorbuffer'
391
392 def setup(self):
393 self.xaxis = 'time'
394 self.ncols = 1
395 self.nrows = self.data.shape('param')[0]
396 self.nplots = self.nrows
397 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
398
399 if not self.xlabel:
400 self.xlabel = 'Time'
401
402 self.ylabel = 'Range [km]'
403 if not self.titles:
404 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
405
406 def update(self, dataOut):
407 data = {
408 #'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)[0:3,:] # SNL, VERTICAL, ZONAL
409 'param' : dataOut.data_output[0:3,:] # SNL, VERTICAL, ZONAL
410 }
411
412 meta = {}
413
414 return data, meta
415
416 def plot(self):
417 # self.data.normalize_heights()
418 self.x = self.data.times
419 self.y = self.data.yrange
420 self.z = self.data['param']
421
422
423 self.z = numpy.ma.masked_invalid(self.z)
424
425 if self.decimation is None:
426 x, y, z = self.fill_gaps(self.x, self.y, self.z)
427 else:
428 x, y, z = self.fill_gaps(*self.decimate())
429
430 for n, ax in enumerate(self.axes):
431 self.zmax = self.zmax if self.zmax is not None else numpy.max(
432 self.z[n])
433 self.zmin = self.zmin if self.zmin is not None else numpy.min(
434 self.z[n])
435
436 if ax.firsttime:
437 if self.zlimits is not None:
438 self.zmin, self.zmax = self.zlimits[n]
439
440 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
441 vmin=self.zmin,
442 vmax=self.zmax,
443 cmap=self.cmaps[n]
444 )
445 else:
446 if self.zlimits is not None:
447 self.zmin, self.zmax = self.zlimits[n]
448 ax.plt.remove()
449 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
450 vmin=self.zmin,
451 vmax=self.zmax,
452 cmap=self.cmaps[n]
453 )
454
455 class AverageDriftsPlot_v2(Plot):
456 '''
457 Plot for average 150 Km echoes
458 '''
459
460 CODE = 'average'
461 plot_type = 'scatterbuffer'
462
463 def setup(self):
464 self.xaxis = 'time'
465 self.ncols = 1
466
467 self.nplots = 2
468 self.nrows = 2
469
470 self.ylabel = 'Velocity\nm/s'
471 self.xlabel = 'Local time'
472 #self.titles = ['VERTICAL VELOCITY: AVERAGE AND ERRORS', 'ZONAL VELOCITY: AVERAGE AND ERRORS']
473 self.titles = ['VERTICAL VELOCITY: AVERAGE', 'ZONAL VELOCITY: AVERAGE']
474
475 self.colorbar = False
476 self.plots_adjust.update({'hspace':0.5, 'left': 0.1, 'bottom': 0.1, 'right':0.95, 'top': 0.95 })
477
478
479 def update(self, dataOut):
480
481 data = {}
482 meta = {}
483
484 #data['average']= numpy.nanmean(dataOut.data_output[1:3,:], axis=1) # VERTICAL, ZONAL
485 data['average']= numpy.nanmean(dataOut.data_output[1:3,:], axis=1) # VERTICAL, ZONAL
486 data['error']= numpy.nanmean(dataOut.data_output[3:,:], axis=1) # ERROR VERTICAL, ERROR ZONAL
487 meta['yrange'] = numpy.array([])
488
489 return data, meta
490
491 def plot(self):
492
493 self.x = self.data.times
494 #self.xmin = self.data.min_time
495 #self.xmax = self.xmin + self.xrange * 60 * 60
496 self.y = self.data['average']
497 print('self.y:', self.y.shape)
498 self.y_error = self.data['error']
499 print('self.y_error:', self.y_error.shape)
500
501 for n, ax in enumerate(self.axes):
502 if ax.firsttime:
503 self.ymin = self.ymin if self.ymin is not None else -50
504 self.ymax = self.ymax if self.ymax is not None else 50
505 self.axes[n].plot(self.x, self.y[n], c='r', ls=':', lw=1)
506 else:
507 self.axes[n].lines[0].set_data(self.x, self.y[n])
508 '''
509 for n, ax in enumerate(self.axes):
510
511 if ax.firsttime:
512 self.ymin = self.ymin if self.ymin is not None else -50
513 self.ymax = self.ymax if self.ymax is not None else 50
514 ax.scatter(self.x, self.y[n], c='g', s=0.8)
515 #ax.errorbar(self.x, self.y[n], yerr = self.y_error[n,:], ecolor='r', elinewidth=0.2, fmt='|')
516 else:
517 ax.scatter(self.x, self.y[n], c='g', s=0.8)
518 #ax.errorbar(self.x, self.y[n], yerr = self.y_error[n,:], ecolor='r', elinewidth=0.2, fmt='|')
519 '''
520 class AverageDriftsPlot_bck(Plot):
521 '''
522 Plot for average 150 Km echoes
523 '''
524 CODE = 'average'
525 plot_type = 'scatterbuffer'
526
527 def setup(self):
528 self.xaxis = 'time'
529 self.ncols = 1
530 self.nplots = 2
531 self.nrows = 2
532 self.ylabel = 'Velocity\nm/s'
533 self.xlabel = 'Time'
534 self.titles = ['VERTICAL VELOCITY: AVERAGE', 'ZONAL VELOCITY: AVERAGE']
535 self.colorbar = False
536 self.plots_adjust.update({'hspace':0.5, 'left': 0.1, 'bottom': 0.1, 'right':0.95, 'top': 0.95 })
537
538
539
540 def update(self, dataOut):
541
542 data = {}
543 meta = {}
544 data['average']= numpy.nanmean(dataOut.data_output[1:3,:], axis=1) # VERTICAL, ZONAL
545 meta['yrange'] = numpy.array([])
546
547 return data, meta
548
549 def plot(self):
550
551 self.x = self.data.times
552 self.y = self.data['average']
553
554 for n, ax in enumerate(self.axes):
555 if ax.firsttime:
556
557 if self.zlimits is not None:
558 self.axes[n].set_ylim(self.zlimits[n])
559 self.axes[n].plot(self.x, self.y[n], c='r', ls='-', lw=1)
560 else:
561
562 if self.zlimits is not None:
563 ax.set_ylim((self.zlimits[n]))
564 self.axes[n].lines[0].set_data(self.x, self.y[n])
565
566 class AverageDriftsPlot(Plot):
567 '''
568 Plot for average 150 Km echoes
569 '''
570 CODE = 'average'
571 plot_type = 'scatterbuffer'
572
573 def setup(self):
574 self.xaxis = 'time'
575 self.ncols = 1
576 self.nplots = 2
577 self.nrows = 2
578 self.ylabel = 'Velocity\nm/s'
579 self.xlabel = 'Time'
580 self.titles = ['VERTICAL VELOCITY: AVERAGE', 'ZONAL VELOCITY: AVERAGE']
581 self.colorbar = False
582 self.plots_adjust.update({'hspace':0.5, 'left': 0.1, 'bottom': 0.1, 'right':0.95, 'top': 0.95 })
583
584
585
586 def update(self, dataOut):
587
588 data = {}
589 meta = {}
590 data['average']= dataOut.avg_output[0:2] # VERTICAL, ZONAL velocities
591 meta['yrange'] = numpy.array([])
592
593 return data, meta
594
595 def plot(self):
596
597 self.x = self.data.times
598 self.y = self.data['average']
599
600 for n, ax in enumerate(self.axes):
601
602 if ax.firsttime:
603
604 if self.zlimits is not None:
605 ax.set_ylim((self.zlimits[n]))
606 self.axes[n].plot(self.x, self.y[n], c='r', ls='-', lw=1)
607 else:
608
609 if self.zlimits is not None:
610 ax.set_ylim((self.zlimits[n]))
611 self.axes[n].lines[0].set_data(self.x, self.y[n])
612
General Comments 0
You need to be logged in to leave comments. Login now