##// END OF EJS Templates
schain
RSS
Block360_vRF4 and Weather_vRF_Plot Fixed and Tested
rflores -
r1447:73923d8d784b
parent child
Show More
Side by Side Unified
Context file:
r1447:73923d8d784b Loading diff...:

The requested changes are too big and content was truncated. Show full diff

Show file before | Show file after | Hide comments
@@ -1,2588 +1,2589
1 1 import os
2 2 import datetime
3 3 import numpy
4 4 from mpl_toolkits.axisartist.grid_finder import FixedLocator, DictFormatter
5 5
6 6 from schainpy.model.graphics.jroplot_base import Plot, plt
7 7 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot, SpectraCutPlot
8 8 from schainpy.utils import log
9 9 # libreria wradlib
10 10 import wradlib as wrl
11 11
12 12 EARTH_RADIUS = 6.3710e3
13 13
14 14
15 15 def ll2xy(lat1, lon1, lat2, lon2):
16 16
17 17 p = 0.017453292519943295
18 18 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
19 19 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
20 20 r = 12742 * numpy.arcsin(numpy.sqrt(a))
21 21 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
22 22 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
23 23 theta = -theta + numpy.pi/2
24 24 return r*numpy.cos(theta), r*numpy.sin(theta)
25 25
26 26
27 27 def km2deg(km):
28 28 '''
29 29 Convert distance in km to degrees
30 30 '''
31 31
32 32 return numpy.rad2deg(km/EARTH_RADIUS)
33 33
34 34
35 35
36 36 class SpectralMomentsPlot(SpectraPlot):
37 37 '''
38 38 Plot for Spectral Moments
39 39 '''
40 40 CODE = 'spc_moments'
41 41 # colormap = 'jet'
42 42 # plot_type = 'pcolor'
43 43
44 44 class DobleGaussianPlot(SpectraPlot):
45 45 '''
46 46 Plot for Double Gaussian Plot
47 47 '''
48 48 CODE = 'gaussian_fit'
49 49 # colormap = 'jet'
50 50 # plot_type = 'pcolor'
51 51
52 52 class DoubleGaussianSpectraCutPlot(SpectraCutPlot):
53 53 '''
54 54 Plot SpectraCut with Double Gaussian Fit
55 55 '''
56 56 CODE = 'cut_gaussian_fit'
57 57
58 58 class SnrPlot(RTIPlot):
59 59 '''
60 60 Plot for SNR Data
61 61 '''
62 62
63 63 CODE = 'snr'
64 64 colormap = 'jet'
65 65
66 66 def update(self, dataOut):
67 67
68 68 data = {
69 69 'snr': 10*numpy.log10(dataOut.data_snr)
70 70 }
71 71
72 72 return data, {}
73 73
74 74 class DopplerPlot(RTIPlot):
75 75 '''
76 76 Plot for DOPPLER Data (1st moment)
77 77 '''
78 78
79 79 CODE = 'dop'
80 80 colormap = 'jet'
81 81
82 82 def update(self, dataOut):
83 83
84 84 data = {
85 85 'dop': 10*numpy.log10(dataOut.data_dop)
86 86 }
87 87
88 88 return data, {}
89 89
90 90 class PowerPlot(RTIPlot):
91 91 '''
92 92 Plot for Power Data (0 moment)
93 93 '''
94 94
95 95 CODE = 'pow'
96 96 colormap = 'jet'
97 97
98 98 def update(self, dataOut):
99 99 data = {
100 100 'pow': 10*numpy.log10(dataOut.data_pow/dataOut.normFactor)
101 101 }
102 102 return data, {}
103 103
104 104 class SpectralWidthPlot(RTIPlot):
105 105 '''
106 106 Plot for Spectral Width Data (2nd moment)
107 107 '''
108 108
109 109 CODE = 'width'
110 110 colormap = 'jet'
111 111
112 112 def update(self, dataOut):
113 113
114 114 data = {
115 115 'width': dataOut.data_width
116 116 }
117 117
118 118 return data, {}
119 119
120 120 class SkyMapPlot(Plot):
121 121 '''
122 122 Plot for meteors detection data
123 123 '''
124 124
125 125 CODE = 'param'
126 126
127 127 def setup(self):
128 128
129 129 self.ncols = 1
130 130 self.nrows = 1
131 131 self.width = 7.2
132 132 self.height = 7.2
133 133 self.nplots = 1
134 134 self.xlabel = 'Zonal Zenith Angle (deg)'
135 135 self.ylabel = 'Meridional Zenith Angle (deg)'
136 136 self.polar = True
137 137 self.ymin = -180
138 138 self.ymax = 180
139 139 self.colorbar = False
140 140
141 141 def plot(self):
142 142
143 143 arrayParameters = numpy.concatenate(self.data['param'])
144 144 error = arrayParameters[:, -1]
145 145 indValid = numpy.where(error == 0)[0]
146 146 finalMeteor = arrayParameters[indValid, :]
147 147 finalAzimuth = finalMeteor[:, 3]
148 148 finalZenith = finalMeteor[:, 4]
149 149
150 150 x = finalAzimuth * numpy.pi / 180
151 151 y = finalZenith
152 152
153 153 ax = self.axes[0]
154 154
155 155 if ax.firsttime:
156 156 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
157 157 else:
158 158 ax.plot.set_data(x, y)
159 159
160 160 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
161 161 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
162 162 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
163 163 dt2,
164 164 len(x))
165 165 self.titles[0] = title
166 166
167 167
168 168 class GenericRTIPlot(Plot):
169 169 '''
170 170 Plot for data_xxxx object
171 171 '''
172 172
173 173 CODE = 'param'
174 174 colormap = 'viridis'
175 175 plot_type = 'pcolorbuffer'
176 176
177 177 def setup(self):
178 178 self.xaxis = 'time'
179 179 self.ncols = 1
180 180 self.nrows = self.data.shape('param')[0]
181 181 self.nplots = self.nrows
182 182 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
183 183
184 184 if not self.xlabel:
185 185 self.xlabel = 'Time'
186 186
187 187 self.ylabel = 'Range [km]'
188 188 if not self.titles:
189 189 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
190 190
191 191 def update(self, dataOut):
192 192
193 193 data = {
194 194 'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)
195 195 }
196 196
197 197 meta = {}
198 198
199 199 return data, meta
200 200
201 201 def plot(self):
202 202 # self.data.normalize_heights()
203 203 self.x = self.data.times
204 204 self.y = self.data.yrange
205 205 self.z = self.data['param']
206 206 self.z = 10*numpy.log10(self.z)
207 207 self.z = numpy.ma.masked_invalid(self.z)
208 208
209 209 if self.decimation is None:
210 210 x, y, z = self.fill_gaps(self.x, self.y, self.z)
211 211 else:
212 212 x, y, z = self.fill_gaps(*self.decimate())
213 213
214 214 for n, ax in enumerate(self.axes):
215 215
216 216 self.zmax = self.zmax if self.zmax is not None else numpy.max(
217 217 self.z[n])
218 218 self.zmin = self.zmin if self.zmin is not None else numpy.min(
219 219 self.z[n])
220 220
221 221 if ax.firsttime:
222 222 if self.zlimits is not None:
223 223 self.zmin, self.zmax = self.zlimits[n]
224 224
225 225 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
226 226 vmin=self.zmin,
227 227 vmax=self.zmax,
228 228 cmap=self.cmaps[n]
229 229 )
230 230 else:
231 231 if self.zlimits is not None:
232 232 self.zmin, self.zmax = self.zlimits[n]
233 233 ax.collections.remove(ax.collections[0])
234 234 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
235 235 vmin=self.zmin,
236 236 vmax=self.zmax,
237 237 cmap=self.cmaps[n]
238 238 )
239 239
240 240
241 241 class PolarMapPlot(Plot):
242 242 '''
243 243 Plot for weather radar
244 244 '''
245 245
246 246 CODE = 'param'
247 247 colormap = 'seismic'
248 248
249 249 def setup(self):
250 250 self.ncols = 1
251 251 self.nrows = 1
252 252 self.width = 9
253 253 self.height = 8
254 254 self.mode = self.data.meta['mode']
255 255 if self.channels is not None:
256 256 self.nplots = len(self.channels)
257 257 self.nrows = len(self.channels)
258 258 else:
259 259 self.nplots = self.data.shape(self.CODE)[0]
260 260 self.nrows = self.nplots
261 261 self.channels = list(range(self.nplots))
262 262 if self.mode == 'E':
263 263 self.xlabel = 'Longitude'
264 264 self.ylabel = 'Latitude'
265 265 else:
266 266 self.xlabel = 'Range (km)'
267 267 self.ylabel = 'Height (km)'
268 268 self.bgcolor = 'white'
269 269 self.cb_labels = self.data.meta['units']
270 270 self.lat = self.data.meta['latitude']
271 271 self.lon = self.data.meta['longitude']
272 272 self.xmin, self.xmax = float(
273 273 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
274 274 self.ymin, self.ymax = float(
275 275 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
276 276 # self.polar = True
277 277
278 278 def plot(self):
279 279
280 280 for n, ax in enumerate(self.axes):
281 281 data = self.data['param'][self.channels[n]]
282 282
283 283 zeniths = numpy.linspace(
284 284 0, self.data.meta['max_range'], data.shape[1])
285 285 if self.mode == 'E':
286 286 azimuths = -numpy.radians(self.data.yrange)+numpy.pi/2
287 287 r, theta = numpy.meshgrid(zeniths, azimuths)
288 288 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
289 289 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
290 290 x = km2deg(x) + self.lon
291 291 y = km2deg(y) + self.lat
292 292 else:
293 293 azimuths = numpy.radians(self.data.yrange)
294 294 r, theta = numpy.meshgrid(zeniths, azimuths)
295 295 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
296 296 self.y = zeniths
297 297
298 298 if ax.firsttime:
299 299 if self.zlimits is not None:
300 300 self.zmin, self.zmax = self.zlimits[n]
301 301 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
302 302 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
303 303 vmin=self.zmin,
304 304 vmax=self.zmax,
305 305 cmap=self.cmaps[n])
306 306 else:
307 307 if self.zlimits is not None:
308 308 self.zmin, self.zmax = self.zlimits[n]
309 309 ax.collections.remove(ax.collections[0])
310 310 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
311 311 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
312 312 vmin=self.zmin,
313 313 vmax=self.zmax,
314 314 cmap=self.cmaps[n])
315 315
316 316 if self.mode == 'A':
317 317 continue
318 318
319 319 # plot district names
320 320 f = open('/data/workspace/schain_scripts/distrito.csv')
321 321 for line in f:
322 322 label, lon, lat = [s.strip() for s in line.split(',') if s]
323 323 lat = float(lat)
324 324 lon = float(lon)
325 325 # ax.plot(lon, lat, '.b', ms=2)
326 326 ax.text(lon, lat, label.decode('utf8'), ha='center',
327 327 va='bottom', size='8', color='black')
328 328
329 329 # plot limites
330 330 limites = []
331 331 tmp = []
332 332 for line in open('/data/workspace/schain_scripts/lima.csv'):
333 333 if '#' in line:
334 334 if tmp:
335 335 limites.append(tmp)
336 336 tmp = []
337 337 continue
338 338 values = line.strip().split(',')
339 339 tmp.append((float(values[0]), float(values[1])))
340 340 for points in limites:
341 341 ax.add_patch(
342 342 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
343 343
344 344 # plot Cuencas
345 345 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
346 346 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
347 347 values = [line.strip().split(',') for line in f]
348 348 points = [(float(s[0]), float(s[1])) for s in values]
349 349 ax.add_patch(Polygon(points, ec='b', fc='none'))
350 350
351 351 # plot grid
352 352 for r in (15, 30, 45, 60):
353 353 ax.add_artist(plt.Circle((self.lon, self.lat),
354 354 km2deg(r), color='0.6', fill=False, lw=0.2))
355 355 ax.text(
356 356 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
357 357 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
358 358 '{}km'.format(r),
359 359 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
360 360
361 361 if self.mode == 'E':
362 362 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
363 363 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
364 364 else:
365 365 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
366 366 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
367 367
368 368 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
369 369 self.titles = ['{} {}'.format(
370 370 self.data.parameters[x], title) for x in self.channels]
371 371
372 372 class WeatherPlot(Plot):
373 373 CODE = 'weather'
374 374 plot_name = 'weather'
375 375 plot_type = 'ppistyle'
376 376 buffering = False
377 377
378 378 def setup(self):
379 379 self.ncols = 1
380 380 self.nrows = 1
381 381 self.width =8
382 382 self.height =8
383 383 self.nplots= 1
384 384 self.ylabel= 'Range [Km]'
385 385 self.titles= ['Weather']
386 386 self.colorbar=False
387 387 self.ini =0
388 388 self.len_azi =0
389 389 self.buffer_ini = None
390 390 self.buffer_azi = None
391 391 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
392 392 self.flag =0
393 393 self.indicador= 0
394 394 self.last_data_azi = None
395 395 self.val_mean = None
396 396
397 397 def update(self, dataOut):
398 398
399 399 data = {}
400 400 meta = {}
401 401 if hasattr(dataOut, 'dataPP_POWER'):
402 402 factor = 1
403 403 if hasattr(dataOut, 'nFFTPoints'):
404 404 factor = dataOut.normFactor
405 405 #print("DIME EL SHAPE PORFAVOR",dataOut.data_360.shape)
406 406 data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
407 407 data['azi'] = dataOut.data_azi
408 408 data['ele'] = dataOut.data_ele
409 409 return data, meta
410 410
411 411 def get2List(self,angulos):
412 412 list1=[]
413 413 list2=[]
414 414 for i in reversed(range(len(angulos))):
415 415 diff_ = angulos[i]-angulos[i-1]
416 416 if diff_ >1.5:
417 417 list1.append(i-1)
418 418 list2.append(diff_)
419 419 return list(reversed(list1)),list(reversed(list2))
420 420
421 421 def fixData360(self,list_,ang_):
422 422 if list_[0]==-1:
423 423 vec = numpy.where(ang_<ang_[0])
424 424 ang_[vec] = ang_[vec]+360
425 425 return ang_
426 426 return ang_
427 427
428 428 def fixData360HL(self,angulos):
429 429 vec = numpy.where(angulos>=360)
430 430 angulos[vec]=angulos[vec]-360
431 431 return angulos
432 432
433 433 def search_pos(self,pos,list_):
434 434 for i in range(len(list_)):
435 435 if pos == list_[i]:
436 436 return True,i
437 437 i=None
438 438 return False,i
439 439
440 440 def fixDataComp(self,ang_,list1_,list2_):
441 441 size = len(ang_)
442 442 size2 = 0
443 443 for i in range(len(list2_)):
444 444 size2=size2+round(list2_[i])-1
445 445 new_size= size+size2
446 446 ang_new = numpy.zeros(new_size)
447 447 ang_new2 = numpy.zeros(new_size)
448 448
449 449 tmp = 0
450 450 c = 0
451 451 for i in range(len(ang_)):
452 452 ang_new[tmp +c] = ang_[i]
453 453 ang_new2[tmp+c] = ang_[i]
454 454 condition , value = self.search_pos(i,list1_)
455 455 if condition:
456 456 pos = tmp + c + 1
457 457 for k in range(round(list2_[value])-1):
458 458 ang_new[pos+k] = ang_new[pos+k-1]+1
459 459 ang_new2[pos+k] = numpy.nan
460 460 tmp = pos +k
461 461 c = 0
462 462 c=c+1
463 463 return ang_new,ang_new2
464 464
465 465 def globalCheckPED(self,angulos):
466 466 l1,l2 = self.get2List(angulos)
467 467 if len(l1)>0:
468 468 angulos2 = self.fixData360(list_=l1,ang_=angulos)
469 469 l1,l2 = self.get2List(angulos2)
470 470
471 471 ang1_,ang2_ = self.fixDataComp(ang_=angulos2,list1_=l1,list2_=l2)
472 472 ang1_ = self.fixData360HL(ang1_)
473 473 ang2_ = self.fixData360HL(ang2_)
474 474 else:
475 475 ang1_= angulos
476 476 ang2_= angulos
477 477 return ang1_,ang2_
478 478
479 479 def analizeDATA(self,data_azi):
480 480 list1 = []
481 481 list2 = []
482 482 dat = data_azi
483 483 for i in reversed(range(1,len(dat))):
484 484 if dat[i]>dat[i-1]:
485 485 diff = int(dat[i])-int(dat[i-1])
486 486 else:
487 487 diff = 360+int(dat[i])-int(dat[i-1])
488 488 if diff > 1:
489 489 list1.append(i-1)
490 490 list2.append(diff-1)
491 491 return list1,list2
492 492
493 493 def fixDATANEW(self,data_azi,data_weather):
494 494 list1,list2 = self.analizeDATA(data_azi)
495 495 if len(list1)== 0:
496 496 return data_azi,data_weather
497 497 else:
498 498 resize = 0
499 499 for i in range(len(list2)):
500 500 resize= resize + list2[i]
501 501 new_data_azi = numpy.resize(data_azi,resize)
502 502 new_data_weather= numpy.resize(date_weather,resize)
503 503
504 504 for i in range(len(list2)):
505 505 j=0
506 506 position=list1[i]+1
507 507 for j in range(list2[i]):
508 508 new_data_azi[position+j]=new_data_azi[position+j-1]+1
509 509 return new_data_azi
510 510
511 511 def fixDATA(self,data_azi):
512 512 data=data_azi
513 513 for i in range(len(data)):
514 514 if numpy.isnan(data[i]):
515 515 data[i]=data[i-1]+1
516 516 return data
517 517
518 518 def replaceNAN(self,data_weather,data_azi,val):
519 519 data= data_azi
520 520 data_T= data_weather
521 521 if data.shape[0]> data_T.shape[0]:
522 522 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
523 523 c = 0
524 524 for i in range(len(data)):
525 525 if numpy.isnan(data[i]):
526 526 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
527 527 else:
528 528 data_N[i,:]=data_T[c,:]
529 529 c=c+1
530 530 return data_N
531 531 else:
532 532 for i in range(len(data)):
533 533 if numpy.isnan(data[i]):
534 534 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
535 535 return data_T
536 536
537 537 def const_ploteo(self,data_weather,data_azi,step,res):
538 538 if self.ini==0:
539 539 #-------
540 540 n = (360/res)-len(data_azi)
541 541 #--------------------- new -------------------------
542 542 data_azi_new ,data_azi_old= self.globalCheckPED(data_azi)
543 543 #------------------------
544 544 start = data_azi_new[-1] + res
545 545 end = data_azi_new[0] - res
546 546 #------ new
547 547 self.last_data_azi = end
548 548 if start>end:
549 549 end = end + 360
550 550 azi_vacia = numpy.linspace(start,end,int(n))
551 551 azi_vacia = numpy.where(azi_vacia>360,azi_vacia-360,azi_vacia)
552 552 data_azi = numpy.hstack((data_azi_new,azi_vacia))
553 553 # RADAR
554 554 val_mean = numpy.mean(data_weather[:,-1])
555 555 self.val_mean = val_mean
556 556 data_weather_cmp = numpy.ones([(360-data_weather.shape[0]),data_weather.shape[1]])*val_mean
557 557 data_weather = self.replaceNAN(data_weather=data_weather,data_azi=data_azi_old,val=self.val_mean)
558 558 data_weather = numpy.vstack((data_weather,data_weather_cmp))
559 559 else:
560 560 # azimuth
561 561 flag=0
562 562 start_azi = self.res_azi[0]
563 563 #-----------new------------
564 564 data_azi ,data_azi_old= self.globalCheckPED(data_azi)
565 565 data_weather = self.replaceNAN(data_weather=data_weather,data_azi=data_azi_old,val=self.val_mean)
566 566 #--------------------------
567 567 start = data_azi[0]
568 568 end = data_azi[-1]
569 569 self.last_data_azi= end
570 570 if start< start_azi:
571 571 start = start +360
572 572 if end <start_azi:
573 573 end = end +360
574 574
575 575 pos_ini = int((start-start_azi)/res)
576 576 len_azi = len(data_azi)
577 577 if (360-pos_ini)<len_azi:
578 578 if pos_ini+1==360:
579 579 pos_ini=0
580 580 else:
581 581 flag=1
582 582 dif= 360-pos_ini
583 583 comp= len_azi-dif
584 584 #-----------------
585 585 if flag==0:
586 586 # AZIMUTH
587 587 self.res_azi[pos_ini:pos_ini+len_azi] = data_azi
588 588 # RADAR
589 589 self.res_weather[pos_ini:pos_ini+len_azi,:] = data_weather
590 590 else:
591 591 # AZIMUTH
592 592 self.res_azi[pos_ini:pos_ini+dif] = data_azi[0:dif]
593 593 self.res_azi[0:comp] = data_azi[dif:]
594 594 # RADAR
595 595 self.res_weather[pos_ini:pos_ini+dif,:] = data_weather[0:dif,:]
596 596 self.res_weather[0:comp,:] = data_weather[dif:,:]
597 597 flag=0
598 598 data_azi = self.res_azi
599 599 data_weather = self.res_weather
600 600
601 601 return data_weather,data_azi
602 602
603 603 def plot(self):
604 604 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
605 605 data = self.data[-1]
606 606 r = self.data.yrange
607 607 delta_height = r[1]-r[0]
608 608 r_mask = numpy.where(r>=0)[0]
609 609 r = numpy.arange(len(r_mask))*delta_height
610 610 self.y = 2*r
611 611 # RADAR
612 612 #data_weather = data['weather']
613 613 # PEDESTAL
614 614 #data_azi = data['azi']
615 615 res = 1
616 616 # STEP
617 617 step = (360/(res*data['weather'].shape[0]))
618 618
619 619 self.res_weather, self.res_azi = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_azi=data['azi'],step=step,res=res)
620 620 self.res_ele = numpy.mean(data['ele'])
621 621 ################# PLOTEO ###################
622 622 for i,ax in enumerate(self.axes):
623 623 self.zmin = self.zmin if self.zmin else 20
624 624 self.zmax = self.zmax if self.zmax else 80
625 625 if ax.firsttime:
626 626 plt.clf()
627 627 cgax, pm = wrl.vis.plot_ppi(self.res_weather,r=r,az=self.res_azi,fig=self.figures[0], proj='cg', vmin=self.zmin, vmax=self.zmax)
628 628 else:
629 629 plt.clf()
630 630 cgax, pm = wrl.vis.plot_ppi(self.res_weather,r=r,az=self.res_azi,fig=self.figures[0], proj='cg', vmin=self.zmin, vmax=self.zmax)
631 631 caax = cgax.parasites[0]
632 632 paax = cgax.parasites[1]
633 633 cbar = plt.gcf().colorbar(pm, pad=0.075)
634 634 caax.set_xlabel('x_range [km]')
635 635 caax.set_ylabel('y_range [km]')
636 636 plt.text(1.0, 1.05, 'Azimuth '+str(thisDatetime)+" Step "+str(self.ini)+ " EL: "+str(round(self.res_ele, 1)), transform=caax.transAxes, va='bottom',ha='right')
637 637
638 638 self.ini= self.ini+1
639 639
640 640
641 641 class WeatherRHIPlot(Plot):
642 642 CODE = 'weather'
643 643 plot_name = 'weather'
644 644 plot_type = 'rhistyle'
645 645 buffering = False
646 646 data_ele_tmp = None
647 647
648 648 def setup(self):
649 649 print("********************")
650 650 print("********************")
651 651 print("********************")
652 652 print("SETUP WEATHER PLOT")
653 653 self.ncols = 1
654 654 self.nrows = 1
655 655 self.nplots= 1
656 656 self.ylabel= 'Range [Km]'
657 657 self.titles= ['Weather']
658 658 if self.channels is not None:
659 659 self.nplots = len(self.channels)
660 660 self.nrows = len(self.channels)
661 661 else:
662 662 self.nplots = self.data.shape(self.CODE)[0]
663 663 self.nrows = self.nplots
664 664 self.channels = list(range(self.nplots))
665 665 print("channels",self.channels)
666 666 print("que saldra", self.data.shape(self.CODE)[0])
667 667 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
668 668 print("self.titles",self.titles)
669 669 self.colorbar=False
670 670 self.width =12
671 671 self.height =8
672 672 self.ini =0
673 673 self.len_azi =0
674 674 self.buffer_ini = None
675 675 self.buffer_ele = None
676 676 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
677 677 self.flag =0
678 678 self.indicador= 0
679 679 self.last_data_ele = None
680 680 self.val_mean = None
681 681
682 682 def update(self, dataOut):
683 683
684 684 data = {}
685 685 meta = {}
686 686 if hasattr(dataOut, 'dataPP_POWER'):
687 687 factor = 1
688 688 if hasattr(dataOut, 'nFFTPoints'):
689 689 factor = dataOut.normFactor
690 690 print("dataOut",dataOut.data_360.shape)
691 691 #
692 692 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
693 693 #
694 694 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
695 695 data['azi'] = dataOut.data_azi
696 696 data['ele'] = dataOut.data_ele
697 697 #print("UPDATE")
698 698 #print("data[weather]",data['weather'].shape)
699 699 #print("data[azi]",data['azi'])
700 700 return data, meta
701 701
702 702 def get2List(self,angulos):
703 703 list1=[]
704 704 list2=[]
705 705 for i in reversed(range(len(angulos))):
706 706 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
707 707 diff_ = angulos[i]-angulos[i-1]
708 708 if abs(diff_) >1.5:
709 709 list1.append(i-1)
710 710 list2.append(diff_)
711 711 return list(reversed(list1)),list(reversed(list2))
712 712
713 713 def fixData90(self,list_,ang_):
714 714 if list_[0]==-1:
715 715 vec = numpy.where(ang_<ang_[0])
716 716 ang_[vec] = ang_[vec]+90
717 717 return ang_
718 718 return ang_
719 719
720 720 def fixData90HL(self,angulos):
721 721 vec = numpy.where(angulos>=90)
722 722 angulos[vec]=angulos[vec]-90
723 723 return angulos
724 724
725 725
726 726 def search_pos(self,pos,list_):
727 727 for i in range(len(list_)):
728 728 if pos == list_[i]:
729 729 return True,i
730 730 i=None
731 731 return False,i
732 732
733 733 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
734 734 size = len(ang_)
735 735 size2 = 0
736 736 for i in range(len(list2_)):
737 737 size2=size2+round(abs(list2_[i]))-1
738 738 new_size= size+size2
739 739 ang_new = numpy.zeros(new_size)
740 740 ang_new2 = numpy.zeros(new_size)
741 741
742 742 tmp = 0
743 743 c = 0
744 744 for i in range(len(ang_)):
745 745 ang_new[tmp +c] = ang_[i]
746 746 ang_new2[tmp+c] = ang_[i]
747 747 condition , value = self.search_pos(i,list1_)
748 748 if condition:
749 749 pos = tmp + c + 1
750 750 for k in range(round(abs(list2_[value]))-1):
751 751 if tipo_case==0 or tipo_case==3:#subida
752 752 ang_new[pos+k] = ang_new[pos+k-1]+1
753 753 ang_new2[pos+k] = numpy.nan
754 754 elif tipo_case==1 or tipo_case==2:#bajada
755 755 ang_new[pos+k] = ang_new[pos+k-1]-1
756 756 ang_new2[pos+k] = numpy.nan
757 757
758 758 tmp = pos +k
759 759 c = 0
760 760 c=c+1
761 761 return ang_new,ang_new2
762 762
763 763 def globalCheckPED(self,angulos,tipo_case):
764 764 l1,l2 = self.get2List(angulos)
765 765 ##print("l1",l1)
766 766 ##print("l2",l2)
767 767 if len(l1)>0:
768 768 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
769 769 #l1,l2 = self.get2List(angulos2)
770 770 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
771 771 #ang1_ = self.fixData90HL(ang1_)
772 772 #ang2_ = self.fixData90HL(ang2_)
773 773 else:
774 774 ang1_= angulos
775 775 ang2_= angulos
776 776 return ang1_,ang2_
777 777
778 778
779 779 def replaceNAN(self,data_weather,data_ele,val):
780 780 data= data_ele
781 781 data_T= data_weather
782 782 if data.shape[0]> data_T.shape[0]:
783 783 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
784 784 c = 0
785 785 for i in range(len(data)):
786 786 if numpy.isnan(data[i]):
787 787 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
788 788 else:
789 789 data_N[i,:]=data_T[c,:]
790 790 c=c+1
791 791 return data_N
792 792 else:
793 793 for i in range(len(data)):
794 794 if numpy.isnan(data[i]):
795 795 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
796 796 return data_T
797 797
798 798 def check_case(self,data_ele,ang_max,ang_min):
799 799 start = data_ele[0]
800 800 end = data_ele[-1]
801 801 number = (end-start)
802 802 len_ang=len(data_ele)
803 803 print("start",start)
804 804 print("end",end)
805 805 print("number",number)
806 806
807 807 print("len_ang",len_ang)
808 808
809 809 #exit(1)
810 810
811 811 if start<end and (round(abs(number)+1)>=len_ang or (numpy.argmin(data_ele)==0)):#caso subida
812 812 return 0
813 813 #elif start>end and (round(abs(number)+1)>=len_ang or(numpy.argmax(data_ele)==0)):#caso bajada
814 814 # return 1
815 815 elif round(abs(number)+1)>=len_ang and (start>end or(numpy.argmax(data_ele)==0)):#caso bajada
816 816 return 1
817 817 elif round(abs(number)+1)<len_ang and data_ele[-2]>data_ele[-1]:# caso BAJADA CAMBIO ANG MAX
818 818 return 2
819 819 elif round(abs(number)+1)<len_ang and data_ele[-2]<data_ele[-1] :# caso SUBIDA CAMBIO ANG MIN
820 820 return 3
821 821
822 822
823 823 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min):
824 824 ang_max= ang_max
825 825 ang_min= ang_min
826 826 data_weather=data_weather
827 827 val_ch=val_ch
828 828 ##print("*********************DATA WEATHER**************************************")
829 829 ##print(data_weather)
830 830 if self.ini==0:
831 831 '''
832 832 print("**********************************************")
833 833 print("**********************************************")
834 834 print("***************ini**************")
835 835 print("**********************************************")
836 836 print("**********************************************")
837 837 '''
838 838 #print("data_ele",data_ele)
839 839 #----------------------------------------------------------
840 840 tipo_case = self.check_case(data_ele,ang_max,ang_min)
841 841 print("check_case",tipo_case)
842 842 #exit(1)
843 843 #--------------------- new -------------------------
844 844 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
845 845
846 846 #-------------------------CAMBIOS RHI---------------------------------
847 847 start= ang_min
848 848 end = ang_max
849 849 n= (ang_max-ang_min)/res
850 850 #------ new
851 851 self.start_data_ele = data_ele_new[0]
852 852 self.end_data_ele = data_ele_new[-1]
853 853 if tipo_case==0 or tipo_case==3: # SUBIDA
854 854 n1= round(self.start_data_ele)- start
855 855 n2= end - round(self.end_data_ele)
856 856 print(self.start_data_ele)
857 857 print(self.end_data_ele)
858 858 if n1>0:
859 859 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
860 860 ele1_nan= numpy.ones(n1)*numpy.nan
861 861 data_ele = numpy.hstack((ele1,data_ele_new))
862 862 print("ele1_nan",ele1_nan.shape)
863 863 print("data_ele_old",data_ele_old.shape)
864 864 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
865 865 if n2>0:
866 866 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
867 867 ele2_nan= numpy.ones(n2)*numpy.nan
868 868 data_ele = numpy.hstack((data_ele,ele2))
869 869 print("ele2_nan",ele2_nan.shape)
870 870 print("data_ele_old",data_ele_old.shape)
871 871 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
872 872
873 873 if tipo_case==1 or tipo_case==2: # BAJADA
874 874 data_ele_new = data_ele_new[::-1] # reversa
875 875 data_ele_old = data_ele_old[::-1]# reversa
876 876 data_weather = data_weather[::-1,:]# reversa
877 877 vec= numpy.where(data_ele_new<ang_max)
878 878 data_ele_new = data_ele_new[vec]
879 879 data_ele_old = data_ele_old[vec]
880 880 data_weather = data_weather[vec[0]]
881 881 vec2= numpy.where(0<data_ele_new)
882 882 data_ele_new = data_ele_new[vec2]
883 883 data_ele_old = data_ele_old[vec2]
884 884 data_weather = data_weather[vec2[0]]
885 885 self.start_data_ele = data_ele_new[0]
886 886 self.end_data_ele = data_ele_new[-1]
887 887
888 888 n1= round(self.start_data_ele)- start
889 889 n2= end - round(self.end_data_ele)-1
890 890 print(self.start_data_ele)
891 891 print(self.end_data_ele)
892 892 if n1>0:
893 893 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
894 894 ele1_nan= numpy.ones(n1)*numpy.nan
895 895 data_ele = numpy.hstack((ele1,data_ele_new))
896 896 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
897 897 if n2>0:
898 898 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
899 899 ele2_nan= numpy.ones(n2)*numpy.nan
900 900 data_ele = numpy.hstack((data_ele,ele2))
901 901 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
902 902 # RADAR
903 903 # NOTA data_ele y data_weather es la variable que retorna
904 904 val_mean = numpy.mean(data_weather[:,-1])
905 905 self.val_mean = val_mean
906 906 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
907 907 self.data_ele_tmp[val_ch]= data_ele_old
908 908 else:
909 909 #print("**********************************************")
910 910 #print("****************VARIABLE**********************")
911 911 #-------------------------CAMBIOS RHI---------------------------------
912 912 #---------------------------------------------------------------------
913 913 ##print("INPUT data_ele",data_ele)
914 914 flag=0
915 915 start_ele = self.res_ele[0]
916 916 tipo_case = self.check_case(data_ele,ang_max,ang_min)
917 917 #print("TIPO DE DATA",tipo_case)
918 918 #-----------new------------
919 919 data_ele ,data_ele_old = self.globalCheckPED(data_ele,tipo_case)
920 920 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
921 921
922 922 #-------------------------------NEW RHI ITERATIVO-------------------------
923 923
924 924 if tipo_case==0 : # SUBIDA
925 925 vec = numpy.where(data_ele<ang_max)
926 926 data_ele = data_ele[vec]
927 927 data_ele_old = data_ele_old[vec]
928 928 data_weather = data_weather[vec[0]]
929 929
930 930 vec2 = numpy.where(0<data_ele)
931 931 data_ele= data_ele[vec2]
932 932 data_ele_old= data_ele_old[vec2]
933 933 ##print(data_ele_new)
934 934 data_weather= data_weather[vec2[0]]
935 935
936 936 new_i_ele = int(round(data_ele[0]))
937 937 new_f_ele = int(round(data_ele[-1]))
938 938 #print(new_i_ele)
939 939 #print(new_f_ele)
940 940 #print(data_ele,len(data_ele))
941 941 #print(data_ele_old,len(data_ele_old))
942 942 if new_i_ele< 2:
943 943 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
944 944 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
945 945 self.data_ele_tmp[val_ch][new_i_ele:new_i_ele+len(data_ele)]=data_ele_old
946 946 self.res_ele[new_i_ele:new_i_ele+len(data_ele)]= data_ele
947 947 self.res_weather[val_ch][new_i_ele:new_i_ele+len(data_ele),:]= data_weather
948 948 data_ele = self.res_ele
949 949 data_weather = self.res_weather[val_ch]
950 950
951 951 elif tipo_case==1 : #BAJADA
952 952 data_ele = data_ele[::-1] # reversa
953 953 data_ele_old = data_ele_old[::-1]# reversa
954 954 data_weather = data_weather[::-1,:]# reversa
955 955 vec= numpy.where(data_ele<ang_max)
956 956 data_ele = data_ele[vec]
957 957 data_ele_old = data_ele_old[vec]
958 958 data_weather = data_weather[vec[0]]
959 959 vec2= numpy.where(0<data_ele)
960 960 data_ele = data_ele[vec2]
961 961 data_ele_old = data_ele_old[vec2]
962 962 data_weather = data_weather[vec2[0]]
963 963
964 964
965 965 new_i_ele = int(round(data_ele[0]))
966 966 new_f_ele = int(round(data_ele[-1]))
967 967 #print(data_ele)
968 968 #print(ang_max)
969 969 #print(data_ele_old)
970 970 if new_i_ele <= 1:
971 971 new_i_ele = 1
972 972 if round(data_ele[-1])>=ang_max-1:
973 973 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
974 974 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
975 975 self.data_ele_tmp[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1]=data_ele_old
976 976 self.res_ele[new_i_ele-1:new_i_ele+len(data_ele)-1]= data_ele
977 977 self.res_weather[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1,:]= data_weather
978 978 data_ele = self.res_ele
979 979 data_weather = self.res_weather[val_ch]
980 980
981 981 elif tipo_case==2: #bajada
982 982 vec = numpy.where(data_ele<ang_max)
983 983 data_ele = data_ele[vec]
984 984 data_weather= data_weather[vec[0]]
985 985
986 986 len_vec = len(vec)
987 987 data_ele_new = data_ele[::-1] # reversa
988 988 data_weather = data_weather[::-1,:]
989 989 new_i_ele = int(data_ele_new[0])
990 990 new_f_ele = int(data_ele_new[-1])
991 991
992 992 n1= new_i_ele- ang_min
993 993 n2= ang_max - new_f_ele-1
994 994 if n1>0:
995 995 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
996 996 ele1_nan= numpy.ones(n1)*numpy.nan
997 997 data_ele = numpy.hstack((ele1,data_ele_new))
998 998 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
999 999 if n2>0:
1000 1000 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1001 1001 ele2_nan= numpy.ones(n2)*numpy.nan
1002 1002 data_ele = numpy.hstack((data_ele,ele2))
1003 1003 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1004 1004
1005 1005 self.data_ele_tmp[val_ch] = data_ele_old
1006 1006 self.res_ele = data_ele
1007 1007 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1008 1008 data_ele = self.res_ele
1009 1009 data_weather = self.res_weather[val_ch]
1010 1010
1011 1011 elif tipo_case==3:#subida
1012 1012 vec = numpy.where(0<data_ele)
1013 1013 data_ele= data_ele[vec]
1014 1014 data_ele_new = data_ele
1015 1015 data_ele_old= data_ele_old[vec]
1016 1016 data_weather= data_weather[vec[0]]
1017 1017 pos_ini = numpy.argmin(data_ele)
1018 1018 if pos_ini>0:
1019 1019 len_vec= len(data_ele)
1020 1020 vec3 = numpy.linspace(pos_ini,len_vec-1,len_vec-pos_ini).astype(int)
1021 1021 #print(vec3)
1022 1022 data_ele= data_ele[vec3]
1023 1023 data_ele_new = data_ele
1024 1024 data_ele_old= data_ele_old[vec3]
1025 1025 data_weather= data_weather[vec3]
1026 1026
1027 1027 new_i_ele = int(data_ele_new[0])
1028 1028 new_f_ele = int(data_ele_new[-1])
1029 1029 n1= new_i_ele- ang_min
1030 1030 n2= ang_max - new_f_ele-1
1031 1031 if n1>0:
1032 1032 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1033 1033 ele1_nan= numpy.ones(n1)*numpy.nan
1034 1034 data_ele = numpy.hstack((ele1,data_ele_new))
1035 1035 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1036 1036 if n2>0:
1037 1037 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1038 1038 ele2_nan= numpy.ones(n2)*numpy.nan
1039 1039 data_ele = numpy.hstack((data_ele,ele2))
1040 1040 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1041 1041
1042 1042 self.data_ele_tmp[val_ch] = data_ele_old
1043 1043 self.res_ele = data_ele
1044 1044 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1045 1045 data_ele = self.res_ele
1046 1046 data_weather = self.res_weather[val_ch]
1047 1047 #print("self.data_ele_tmp",self.data_ele_tmp)
1048 1048 return data_weather,data_ele
1049 1049
1050 1050
1051 1051 def plot(self):
1052 1052 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
1053 1053 data = self.data[-1]
1054 1054 r = self.data.yrange
1055 1055 delta_height = r[1]-r[0]
1056 1056 r_mask = numpy.where(r>=0)[0]
1057 1057 ##print("delta_height",delta_height)
1058 1058 #print("r_mask",r_mask,len(r_mask))
1059 1059 r = numpy.arange(len(r_mask))*delta_height
1060 1060 self.y = 2*r
1061 1061 res = 1
1062 1062 ###print("data['weather'].shape[0]",data['weather'].shape[0])
1063 1063 ang_max = self.ang_max
1064 1064 ang_min = self.ang_min
1065 1065 var_ang =ang_max - ang_min
1066 1066 step = (int(var_ang)/(res*data['weather'].shape[0]))
1067 1067 ###print("step",step)
1068 1068 #--------------------------------------------------------
1069 1069 ##print('weather',data['weather'].shape)
1070 1070 ##print('ele',data['ele'].shape)
1071 1071
1072 1072 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1073 1073 ###self.res_azi = numpy.mean(data['azi'])
1074 1074 ###print("self.res_ele",self.res_ele)
1075 1075 plt.clf()
1076 1076 subplots = [121, 122]
1077 1077 cg={'angular_spacing': 20.}
1078 1078 if self.ini==0:
1079 1079 self.data_ele_tmp = numpy.ones([self.nplots,int(var_ang)])*numpy.nan
1080 1080 self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
1081 1081 print("SHAPE",self.data_ele_tmp.shape)
1082 1082
1083 1083 for i,ax in enumerate(self.axes):
1084 1084 self.res_weather[i], self.res_ele = self.const_ploteo(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1085 1085 self.res_azi = numpy.mean(data['azi'])
1086 1086 if i==0:
1087 1087 print("*****************************************************************************to plot**************************",self.res_weather[i].shape)
1088 1088 self.zmin = self.zmin if self.zmin else 20
1089 1089 self.zmax = self.zmax if self.zmax else 80
1090 1090 if ax.firsttime:
1091 1091 #plt.clf()
1092 1092 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj=cg,vmin=self.zmin, vmax=self.zmax)
1093 1093 #fig=self.figures[0]
1094 1094 else:
1095 1095 #plt.clf()
1096 1096 if i==0:
1097 1097 print(self.res_weather[i])
1098 1098 print(self.res_ele)
1099 1099 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj=cg,vmin=self.zmin, vmax=self.zmax)
1100 1100 caax = cgax.parasites[0]
1101 1101 paax = cgax.parasites[1]
1102 1102 cbar = plt.gcf().colorbar(pm, pad=0.075)
1103 1103 caax.set_xlabel('x_range [km]')
1104 1104 caax.set_ylabel('y_range [km]')
1105 1105 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
1106 1106 print("***************************self.ini****************************",self.ini)
1107 1107 self.ini= self.ini+1
1108 1108
1109 1109 class Weather_vRF_Plot(Plot):
1110 1110 CODE = 'PPI'
1111 1111 plot_name = 'PPI'
1112 1112 #plot_type = 'ppistyle'
1113 1113 buffering = False
1114 1114
1115 1115 def setup(self):
1116 1116
1117 1117 self.ncols = 1
1118 1118 self.nrows = 1
1119 1119 self.width =8
1120 1120 self.height =8
1121 1121 self.nplots= 1
1122 1122 self.ylabel= 'Range [Km]'
1123 1123 self.titles= ['PPI']
1124 1124 self.polar = True
1125 1125 if self.channels is not None:
1126 1126 self.nplots = len(self.channels)
1127 1127 self.nrows = len(self.channels)
1128 1128 else:
1129 1129 self.nplots = self.data.shape(self.CODE)[0]
1130 1130 self.nrows = self.nplots
1131 1131 self.channels = list(range(self.nplots))
1132 1132
1133 1133 if self.CODE == 'POWER':
1134 1134 self.cb_label = r'Power (dB)'
1135 1135 elif self.CODE == 'DOPPLER':
1136 1136 self.cb_label = r'Velocity (m/s)'
1137 1137 self.colorbar=True
1138 1138 self.width =8
1139 1139 self.height =8
1140 1140 self.ini =0
1141 1141 self.len_azi =0
1142 1142 self.buffer_ini = None
1143 1143 self.buffer_ele = None
1144 1144 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1145 1145 self.flag =0
1146 1146 self.indicador= 0
1147 1147 self.last_data_ele = None
1148 1148 self.val_mean = None
1149 1149
1150 1150 def update(self, dataOut):
1151 1151
1152 1152 data = {}
1153 1153 meta = {}
1154 1154 if hasattr(dataOut, 'dataPP_POWER'):
1155 1155 factor = 1
1156 1156 if hasattr(dataOut, 'nFFTPoints'):
1157 1157 factor = dataOut.normFactor
1158 1158
1159 1159 if 'pow' in self.attr_data[0].lower():
1160 1160 data['data'] = 10*numpy.log10(getattr(dataOut, self.attr_data[0])/(factor))
1161 1161 else:
1162 1162 data['data'] = getattr(dataOut, self.attr_data[0])/(factor)
1163 1163
1164 1164 data['azi'] = dataOut.data_azi
1165 1165 data['ele'] = dataOut.data_ele
1166 1166
1167 1167 return data, meta
1168 1168
1169 1169 def plot(self):
1170 1170 data = self.data[-1]
1171 1171 r = self.data.yrange
1172 1172 delta_height = r[1]-r[0]
1173 1173 r_mask = numpy.where(r>=0)[0]
1174 1174 self.r_mask = r_mask
1175 1175 r = numpy.arange(len(r_mask))*delta_height
1176 1176 self.y = 2*r
1177 1177 res = 1
1178 1178
1179 var_ang =ang_max - ang_min
1180 step = (int(var_ang)/(res*data['data'].shape[0]))
1179 #var_ang = ang_max - ang_min
1180 #step = (int(var_ang)/(res*data['data'].shape[0]))
1181 1181
1182 1182 z = data['data'][self.channels[0]][:,r_mask]
1183 1183
1184 1184 self.titles = []
1185 1185
1186 1186 self.ymax = self.ymax if self.ymax else numpy.nanmax(r)
1187 1187 self.ymin = self.ymin if self.ymin else numpy.nanmin(r)
1188 1188 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
1189 1189 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
1190 1190 self.ang_min = self.ang_min if self.ang_min else 0
1191 self.ang_max = self.ang_max if self.ang_max else 2*numpy.pi
1191 self.ang_max = self.ang_max if self.ang_max else 360
1192 1192
1193 1193 subplots = [121, 122]
1194 1194
1195 1195 r, theta = numpy.meshgrid(r, numpy.radians(data['azi']) )
1196 1196
1197 1197 for i,ax in enumerate(self.axes):
1198 1198
1199 1199 if ax.firsttime:
1200 1200 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1201 1201 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1202 1202
1203 1203 else:
1204 1204 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1205 1205 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1206 1206
1207 1207 if len(self.channels) !=1:
1208 1208 self.titles = ['{} Ele: {} Channel {}'.format(self.CODE.upper(), str(round(numpy.mean(data['ele']),2)), x) for x in range(self.nrows)]
1209 1209 else:
1210 1210 self.titles = ['{} Ele: {} Channel {}'.format(self.CODE.upper(), str(round(numpy.mean(data['ele']),2)), self.channels[0])]
1211 1211
1212 1212 class WeatherRHI_vRF2_Plot(Plot):
1213 1213 CODE = 'weather'
1214 1214 plot_name = 'weather'
1215 1215 plot_type = 'rhistyle'
1216 1216 buffering = False
1217 1217 data_ele_tmp = None
1218 1218
1219 1219 def setup(self):
1220 1220 print("********************")
1221 1221 print("********************")
1222 1222 print("********************")
1223 1223 print("SETUP WEATHER PLOT")
1224 1224 self.ncols = 1
1225 1225 self.nrows = 1
1226 1226 self.nplots= 1
1227 1227 self.ylabel= 'Range [Km]'
1228 1228 self.titles= ['Weather']
1229 1229 if self.channels is not None:
1230 1230 self.nplots = len(self.channels)
1231 1231 self.nrows = len(self.channels)
1232 1232 else:
1233 1233 self.nplots = self.data.shape(self.CODE)[0]
1234 1234 self.nrows = self.nplots
1235 1235 self.channels = list(range(self.nplots))
1236 1236 print("channels",self.channels)
1237 1237 print("que saldra", self.data.shape(self.CODE)[0])
1238 1238 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
1239 1239 print("self.titles",self.titles)
1240 1240 self.colorbar=False
1241 1241 self.width =8
1242 1242 self.height =8
1243 1243 self.ini =0
1244 1244 self.len_azi =0
1245 1245 self.buffer_ini = None
1246 1246 self.buffer_ele = None
1247 1247 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1248 1248 self.flag =0
1249 1249 self.indicador= 0
1250 1250 self.last_data_ele = None
1251 1251 self.val_mean = None
1252 1252
1253 1253 def update(self, dataOut):
1254 1254
1255 1255 data = {}
1256 1256 meta = {}
1257 1257 if hasattr(dataOut, 'dataPP_POWER'):
1258 1258 factor = 1
1259 1259 if hasattr(dataOut, 'nFFTPoints'):
1260 1260 factor = dataOut.normFactor
1261 1261 print("dataOut",dataOut.data_360.shape)
1262 1262 #
1263 1263 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
1264 1264 #
1265 1265 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
1266 1266 data['azi'] = dataOut.data_azi
1267 1267 data['ele'] = dataOut.data_ele
1268 1268 data['case_flag'] = dataOut.case_flag
1269 1269 #print("UPDATE")
1270 1270 #print("data[weather]",data['weather'].shape)
1271 1271 #print("data[azi]",data['azi'])
1272 1272 return data, meta
1273 1273
1274 1274 def get2List(self,angulos):
1275 1275 list1=[]
1276 1276 list2=[]
1277 1277 for i in reversed(range(len(angulos))):
1278 1278 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
1279 1279 diff_ = angulos[i]-angulos[i-1]
1280 1280 if abs(diff_) >1.5:
1281 1281 list1.append(i-1)
1282 1282 list2.append(diff_)
1283 1283 return list(reversed(list1)),list(reversed(list2))
1284 1284
1285 1285 def fixData90(self,list_,ang_):
1286 1286 if list_[0]==-1:
1287 1287 vec = numpy.where(ang_<ang_[0])
1288 1288 ang_[vec] = ang_[vec]+90
1289 1289 return ang_
1290 1290 return ang_
1291 1291
1292 1292 def fixData90HL(self,angulos):
1293 1293 vec = numpy.where(angulos>=90)
1294 1294 angulos[vec]=angulos[vec]-90
1295 1295 return angulos
1296 1296
1297 1297
1298 1298 def search_pos(self,pos,list_):
1299 1299 for i in range(len(list_)):
1300 1300 if pos == list_[i]:
1301 1301 return True,i
1302 1302 i=None
1303 1303 return False,i
1304 1304
1305 1305 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
1306 1306 size = len(ang_)
1307 1307 size2 = 0
1308 1308 for i in range(len(list2_)):
1309 1309 size2=size2+round(abs(list2_[i]))-1
1310 1310 new_size= size+size2
1311 1311 ang_new = numpy.zeros(new_size)
1312 1312 ang_new2 = numpy.zeros(new_size)
1313 1313
1314 1314 tmp = 0
1315 1315 c = 0
1316 1316 for i in range(len(ang_)):
1317 1317 ang_new[tmp +c] = ang_[i]
1318 1318 ang_new2[tmp+c] = ang_[i]
1319 1319 condition , value = self.search_pos(i,list1_)
1320 1320 if condition:
1321 1321 pos = tmp + c + 1
1322 1322 for k in range(round(abs(list2_[value]))-1):
1323 1323 if tipo_case==0 or tipo_case==3:#subida
1324 1324 ang_new[pos+k] = ang_new[pos+k-1]+1
1325 1325 ang_new2[pos+k] = numpy.nan
1326 1326 elif tipo_case==1 or tipo_case==2:#bajada
1327 1327 ang_new[pos+k] = ang_new[pos+k-1]-1
1328 1328 ang_new2[pos+k] = numpy.nan
1329 1329
1330 1330 tmp = pos +k
1331 1331 c = 0
1332 1332 c=c+1
1333 1333 return ang_new,ang_new2
1334 1334
1335 1335 def globalCheckPED(self,angulos,tipo_case):
1336 1336 l1,l2 = self.get2List(angulos)
1337 1337 ##print("l1",l1)
1338 1338 ##print("l2",l2)
1339 1339 if len(l1)>0:
1340 1340 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
1341 1341 #l1,l2 = self.get2List(angulos2)
1342 1342 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
1343 1343 #ang1_ = self.fixData90HL(ang1_)
1344 1344 #ang2_ = self.fixData90HL(ang2_)
1345 1345 else:
1346 1346 ang1_= angulos
1347 1347 ang2_= angulos
1348 1348 return ang1_,ang2_
1349 1349
1350 1350
1351 1351 def replaceNAN(self,data_weather,data_ele,val):
1352 1352 data= data_ele
1353 1353 data_T= data_weather
1354 1354 if data.shape[0]> data_T.shape[0]:
1355 1355 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
1356 1356 c = 0
1357 1357 for i in range(len(data)):
1358 1358 if numpy.isnan(data[i]):
1359 1359 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1360 1360 else:
1361 1361 data_N[i,:]=data_T[c,:]
1362 1362 c=c+1
1363 1363 return data_N
1364 1364 else:
1365 1365 for i in range(len(data)):
1366 1366 if numpy.isnan(data[i]):
1367 1367 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1368 1368 return data_T
1369 1369
1370 1370 def check_case(self,data_ele,ang_max,ang_min):
1371 1371 start = data_ele[0]
1372 1372 end = data_ele[-1]
1373 1373 number = (end-start)
1374 1374 len_ang=len(data_ele)
1375 1375 print("start",start)
1376 1376 print("end",end)
1377 1377 print("number",number)
1378 1378
1379 1379 print("len_ang",len_ang)
1380 1380
1381 1381 #exit(1)
1382 1382
1383 1383 if start<end and (round(abs(number)+1)>=len_ang or (numpy.argmin(data_ele)==0)):#caso subida
1384 1384 return 0
1385 1385 #elif start>end and (round(abs(number)+1)>=len_ang or(numpy.argmax(data_ele)==0)):#caso bajada
1386 1386 # return 1
1387 1387 elif round(abs(number)+1)>=len_ang and (start>end or(numpy.argmax(data_ele)==0)):#caso bajada
1388 1388 return 1
1389 1389 elif round(abs(number)+1)<len_ang and data_ele[-2]>data_ele[-1]:# caso BAJADA CAMBIO ANG MAX
1390 1390 return 2
1391 1391 elif round(abs(number)+1)<len_ang and data_ele[-2]<data_ele[-1] :# caso SUBIDA CAMBIO ANG MIN
1392 1392 return 3
1393 1393
1394 1394
1395 1395 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min,case_flag):
1396 1396 ang_max= ang_max
1397 1397 ang_min= ang_min
1398 1398 data_weather=data_weather
1399 1399 val_ch=val_ch
1400 1400 ##print("*********************DATA WEATHER**************************************")
1401 1401 ##print(data_weather)
1402 1402 if self.ini==0:
1403 1403 '''
1404 1404 print("**********************************************")
1405 1405 print("**********************************************")
1406 1406 print("***************ini**************")
1407 1407 print("**********************************************")
1408 1408 print("**********************************************")
1409 1409 '''
1410 1410 #print("data_ele",data_ele)
1411 1411 #----------------------------------------------------------
1412 1412 tipo_case = case_flag[-1]
1413 1413 #tipo_case = self.check_case(data_ele,ang_max,ang_min)
1414 1414 print("check_case",tipo_case)
1415 1415 #exit(1)
1416 1416 #--------------------- new -------------------------
1417 1417 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
1418 1418
1419 1419 #-------------------------CAMBIOS RHI---------------------------------
1420 1420 start= ang_min
1421 1421 end = ang_max
1422 1422 n= (ang_max-ang_min)/res
1423 1423 #------ new
1424 1424 self.start_data_ele = data_ele_new[0]
1425 1425 self.end_data_ele = data_ele_new[-1]
1426 1426 if tipo_case==0 or tipo_case==3: # SUBIDA
1427 1427 n1= round(self.start_data_ele)- start
1428 1428 n2= end - round(self.end_data_ele)
1429 1429 print(self.start_data_ele)
1430 1430 print(self.end_data_ele)
1431 1431 if n1>0:
1432 1432 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
1433 1433 ele1_nan= numpy.ones(n1)*numpy.nan
1434 1434 data_ele = numpy.hstack((ele1,data_ele_new))
1435 1435 print("ele1_nan",ele1_nan.shape)
1436 1436 print("data_ele_old",data_ele_old.shape)
1437 1437 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
1438 1438 if n2>0:
1439 1439 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
1440 1440 ele2_nan= numpy.ones(n2)*numpy.nan
1441 1441 data_ele = numpy.hstack((data_ele,ele2))
1442 1442 print("ele2_nan",ele2_nan.shape)
1443 1443 print("data_ele_old",data_ele_old.shape)
1444 1444 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1445 1445
1446 1446 if tipo_case==1 or tipo_case==2: # BAJADA
1447 1447 data_ele_new = data_ele_new[::-1] # reversa
1448 1448 data_ele_old = data_ele_old[::-1]# reversa
1449 1449 data_weather = data_weather[::-1,:]# reversa
1450 1450 vec= numpy.where(data_ele_new<ang_max)
1451 1451 data_ele_new = data_ele_new[vec]
1452 1452 data_ele_old = data_ele_old[vec]
1453 1453 data_weather = data_weather[vec[0]]
1454 1454 vec2= numpy.where(0<data_ele_new)
1455 1455 data_ele_new = data_ele_new[vec2]
1456 1456 data_ele_old = data_ele_old[vec2]
1457 1457 data_weather = data_weather[vec2[0]]
1458 1458 self.start_data_ele = data_ele_new[0]
1459 1459 self.end_data_ele = data_ele_new[-1]
1460 1460
1461 1461 n1= round(self.start_data_ele)- start
1462 1462 n2= end - round(self.end_data_ele)-1
1463 1463 print(self.start_data_ele)
1464 1464 print(self.end_data_ele)
1465 1465 if n1>0:
1466 1466 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
1467 1467 ele1_nan= numpy.ones(n1)*numpy.nan
1468 1468 data_ele = numpy.hstack((ele1,data_ele_new))
1469 1469 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
1470 1470 if n2>0:
1471 1471 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
1472 1472 ele2_nan= numpy.ones(n2)*numpy.nan
1473 1473 data_ele = numpy.hstack((data_ele,ele2))
1474 1474 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1475 1475 # RADAR
1476 1476 # NOTA data_ele y data_weather es la variable que retorna
1477 1477 val_mean = numpy.mean(data_weather[:,-1])
1478 1478 self.val_mean = val_mean
1479 1479 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1480 1480 print("eleold",data_ele_old)
1481 1481 print(self.data_ele_tmp[val_ch])
1482 1482 print(data_ele_old.shape[0])
1483 1483 print(self.data_ele_tmp[val_ch].shape[0])
1484 1484 if (data_ele_old.shape[0]==91 or self.data_ele_tmp[val_ch].shape[0]==91):
1485 1485 import sys
1486 1486 print("EXIT",self.ini)
1487 1487
1488 1488 sys.exit(1)
1489 1489 self.data_ele_tmp[val_ch]= data_ele_old
1490 1490 else:
1491 1491 #print("**********************************************")
1492 1492 #print("****************VARIABLE**********************")
1493 1493 #-------------------------CAMBIOS RHI---------------------------------
1494 1494 #---------------------------------------------------------------------
1495 1495 ##print("INPUT data_ele",data_ele)
1496 1496 flag=0
1497 1497 start_ele = self.res_ele[0]
1498 1498 #tipo_case = self.check_case(data_ele,ang_max,ang_min)
1499 1499 tipo_case = case_flag[-1]
1500 1500 #print("TIPO DE DATA",tipo_case)
1501 1501 #-----------new------------
1502 1502 data_ele ,data_ele_old = self.globalCheckPED(data_ele,tipo_case)
1503 1503 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1504 1504
1505 1505 #-------------------------------NEW RHI ITERATIVO-------------------------
1506 1506
1507 1507 if tipo_case==0 : # SUBIDA
1508 1508 vec = numpy.where(data_ele<ang_max)
1509 1509 data_ele = data_ele[vec]
1510 1510 data_ele_old = data_ele_old[vec]
1511 1511 data_weather = data_weather[vec[0]]
1512 1512
1513 1513 vec2 = numpy.where(0<data_ele)
1514 1514 data_ele= data_ele[vec2]
1515 1515 data_ele_old= data_ele_old[vec2]
1516 1516 ##print(data_ele_new)
1517 1517 data_weather= data_weather[vec2[0]]
1518 1518
1519 1519 new_i_ele = int(round(data_ele[0]))
1520 1520 new_f_ele = int(round(data_ele[-1]))
1521 1521 #print(new_i_ele)
1522 1522 #print(new_f_ele)
1523 1523 #print(data_ele,len(data_ele))
1524 1524 #print(data_ele_old,len(data_ele_old))
1525 1525 if new_i_ele< 2:
1526 1526 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
1527 1527 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
1528 1528 self.data_ele_tmp[val_ch][new_i_ele:new_i_ele+len(data_ele)]=data_ele_old
1529 1529 self.res_ele[new_i_ele:new_i_ele+len(data_ele)]= data_ele
1530 1530 self.res_weather[val_ch][new_i_ele:new_i_ele+len(data_ele),:]= data_weather
1531 1531 data_ele = self.res_ele
1532 1532 data_weather = self.res_weather[val_ch]
1533 1533
1534 1534 elif tipo_case==1 : #BAJADA
1535 1535 data_ele = data_ele[::-1] # reversa
1536 1536 data_ele_old = data_ele_old[::-1]# reversa
1537 1537 data_weather = data_weather[::-1,:]# reversa
1538 1538 vec= numpy.where(data_ele<ang_max)
1539 1539 data_ele = data_ele[vec]
1540 1540 data_ele_old = data_ele_old[vec]
1541 1541 data_weather = data_weather[vec[0]]
1542 1542 vec2= numpy.where(0<data_ele)
1543 1543 data_ele = data_ele[vec2]
1544 1544 data_ele_old = data_ele_old[vec2]
1545 1545 data_weather = data_weather[vec2[0]]
1546 1546
1547 1547
1548 1548 new_i_ele = int(round(data_ele[0]))
1549 1549 new_f_ele = int(round(data_ele[-1]))
1550 1550 #print(data_ele)
1551 1551 #print(ang_max)
1552 1552 #print(data_ele_old)
1553 1553 if new_i_ele <= 1:
1554 1554 new_i_ele = 1
1555 1555 if round(data_ele[-1])>=ang_max-1:
1556 1556 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
1557 1557 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
1558 1558 self.data_ele_tmp[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1]=data_ele_old
1559 1559 self.res_ele[new_i_ele-1:new_i_ele+len(data_ele)-1]= data_ele
1560 1560 self.res_weather[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1,:]= data_weather
1561 1561 data_ele = self.res_ele
1562 1562 data_weather = self.res_weather[val_ch]
1563 1563
1564 1564 elif tipo_case==2: #bajada
1565 1565 vec = numpy.where(data_ele<ang_max)
1566 1566 data_ele = data_ele[vec]
1567 1567 data_weather= data_weather[vec[0]]
1568 1568
1569 1569 len_vec = len(vec)
1570 1570 data_ele_new = data_ele[::-1] # reversa
1571 1571 data_weather = data_weather[::-1,:]
1572 1572 new_i_ele = int(data_ele_new[0])
1573 1573 new_f_ele = int(data_ele_new[-1])
1574 1574
1575 1575 n1= new_i_ele- ang_min
1576 1576 n2= ang_max - new_f_ele-1
1577 1577 if n1>0:
1578 1578 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1579 1579 ele1_nan= numpy.ones(n1)*numpy.nan
1580 1580 data_ele = numpy.hstack((ele1,data_ele_new))
1581 1581 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1582 1582 if n2>0:
1583 1583 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1584 1584 ele2_nan= numpy.ones(n2)*numpy.nan
1585 1585 data_ele = numpy.hstack((data_ele,ele2))
1586 1586 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1587 1587
1588 1588 self.data_ele_tmp[val_ch] = data_ele_old
1589 1589 self.res_ele = data_ele
1590 1590 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1591 1591 data_ele = self.res_ele
1592 1592 data_weather = self.res_weather[val_ch]
1593 1593
1594 1594 elif tipo_case==3:#subida
1595 1595 vec = numpy.where(0<data_ele)
1596 1596 data_ele= data_ele[vec]
1597 1597 data_ele_new = data_ele
1598 1598 data_ele_old= data_ele_old[vec]
1599 1599 data_weather= data_weather[vec[0]]
1600 1600 pos_ini = numpy.argmin(data_ele)
1601 1601 if pos_ini>0:
1602 1602 len_vec= len(data_ele)
1603 1603 vec3 = numpy.linspace(pos_ini,len_vec-1,len_vec-pos_ini).astype(int)
1604 1604 #print(vec3)
1605 1605 data_ele= data_ele[vec3]
1606 1606 data_ele_new = data_ele
1607 1607 data_ele_old= data_ele_old[vec3]
1608 1608 data_weather= data_weather[vec3]
1609 1609
1610 1610 new_i_ele = int(data_ele_new[0])
1611 1611 new_f_ele = int(data_ele_new[-1])
1612 1612 n1= new_i_ele- ang_min
1613 1613 n2= ang_max - new_f_ele-1
1614 1614 if n1>0:
1615 1615 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1616 1616 ele1_nan= numpy.ones(n1)*numpy.nan
1617 1617 data_ele = numpy.hstack((ele1,data_ele_new))
1618 1618 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1619 1619 if n2>0:
1620 1620 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1621 1621 ele2_nan= numpy.ones(n2)*numpy.nan
1622 1622 data_ele = numpy.hstack((data_ele,ele2))
1623 1623 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1624 1624
1625 1625 self.data_ele_tmp[val_ch] = data_ele_old
1626 1626 self.res_ele = data_ele
1627 1627 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1628 1628 data_ele = self.res_ele
1629 1629 data_weather = self.res_weather[val_ch]
1630 1630 #print("self.data_ele_tmp",self.data_ele_tmp)
1631 1631 return data_weather,data_ele
1632 1632
1633 1633
1634 1634 def plot(self):
1635 1635 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
1636 1636 data = self.data[-1]
1637 1637 r = self.data.yrange
1638 1638 delta_height = r[1]-r[0]
1639 1639 r_mask = numpy.where(r>=0)[0]
1640 1640 ##print("delta_height",delta_height)
1641 1641 #print("r_mask",r_mask,len(r_mask))
1642 1642 r = numpy.arange(len(r_mask))*delta_height
1643 1643 self.y = 2*r
1644 1644 res = 1
1645 1645 ###print("data['weather'].shape[0]",data['weather'].shape[0])
1646 1646 ang_max = self.ang_max
1647 1647 ang_min = self.ang_min
1648 1648 var_ang =ang_max - ang_min
1649 1649 step = (int(var_ang)/(res*data['weather'].shape[0]))
1650 1650 ###print("step",step)
1651 1651 #--------------------------------------------------------
1652 1652 ##print('weather',data['weather'].shape)
1653 1653 ##print('ele',data['ele'].shape)
1654 1654
1655 1655 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1656 1656 ###self.res_azi = numpy.mean(data['azi'])
1657 1657 ###print("self.res_ele",self.res_ele)
1658 1658 plt.clf()
1659 1659 subplots = [121, 122]
1660 1660 try:
1661 1661 if self.data[-2]['ele'].max()<data['ele'].max():
1662 1662 self.ini=0
1663 1663 except:
1664 1664 pass
1665 1665 if self.ini==0:
1666 1666 self.data_ele_tmp = numpy.ones([self.nplots,int(var_ang)])*numpy.nan
1667 1667 self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
1668 1668 print("SHAPE",self.data_ele_tmp.shape)
1669 1669
1670 1670 for i,ax in enumerate(self.axes):
1671 1671 self.res_weather[i], self.res_ele = self.const_ploteo(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min,case_flag=self.data['case_flag'])
1672 1672 self.res_azi = numpy.mean(data['azi'])
1673 1673
1674 1674 if ax.firsttime:
1675 1675 #plt.clf()
1676 1676 print("Frist Plot")
1677 1677 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1678 1678 #fig=self.figures[0]
1679 1679 else:
1680 1680 #plt.clf()
1681 1681 print("ELSE PLOT")
1682 1682 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1683 1683 caax = cgax.parasites[0]
1684 1684 paax = cgax.parasites[1]
1685 1685 cbar = plt.gcf().colorbar(pm, pad=0.075)
1686 1686 caax.set_xlabel('x_range [km]')
1687 1687 caax.set_ylabel('y_range [km]')
1688 1688 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
1689 1689 print("***************************self.ini****************************",self.ini)
1690 1690 self.ini= self.ini+1
1691 1691
1692 1692 class WeatherRHI_vRF_Plot(Plot):
1693 1693 CODE = 'weather'
1694 1694 plot_name = 'weather'
1695 1695 plot_type = 'rhistyle'
1696 1696 buffering = False
1697 1697 data_ele_tmp = None
1698 1698
1699 1699 def setup(self):
1700 1700 print("********************")
1701 1701 print("********************")
1702 1702 print("********************")
1703 1703 print("SETUP WEATHER PLOT")
1704 1704 self.ncols = 1
1705 1705 self.nrows = 1
1706 1706 self.nplots= 1
1707 1707 self.ylabel= 'Range [Km]'
1708 1708 self.titles= ['Weather']
1709 1709 if self.channels is not None:
1710 1710 self.nplots = len(self.channels)
1711 1711 self.nrows = len(self.channels)
1712 1712 else:
1713 1713 self.nplots = self.data.shape(self.CODE)[0]
1714 1714 self.nrows = self.nplots
1715 1715 self.channels = list(range(self.nplots))
1716 1716 print("channels",self.channels)
1717 1717 print("que saldra", self.data.shape(self.CODE)[0])
1718 1718 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
1719 1719 print("self.titles",self.titles)
1720 1720 self.colorbar=False
1721 1721 self.width =8
1722 1722 self.height =8
1723 1723 self.ini =0
1724 1724 self.len_azi =0
1725 1725 self.buffer_ini = None
1726 1726 self.buffer_ele = None
1727 1727 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1728 1728 self.flag =0
1729 1729 self.indicador= 0
1730 1730 self.last_data_ele = None
1731 1731 self.val_mean = None
1732 1732
1733 1733 def update(self, dataOut):
1734 1734
1735 1735 data = {}
1736 1736 meta = {}
1737 1737 if hasattr(dataOut, 'dataPP_POWER'):
1738 1738 factor = 1
1739 1739 if hasattr(dataOut, 'nFFTPoints'):
1740 1740 factor = dataOut.normFactor
1741 1741 print("dataOut",dataOut.data_360.shape)
1742 1742 #
1743 1743 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
1744 1744 #
1745 1745 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
1746 1746 data['azi'] = dataOut.data_azi
1747 1747 data['ele'] = dataOut.data_ele
1748 1748 data['case_flag'] = dataOut.case_flag
1749 1749 #print("UPDATE")
1750 1750 #print("data[weather]",data['weather'].shape)
1751 1751 #print("data[azi]",data['azi'])
1752 1752 return data, meta
1753 1753
1754 1754 def get2List(self,angulos):
1755 1755 list1=[]
1756 1756 list2=[]
1757 1757 #print(angulos)
1758 1758 #exit(1)
1759 1759 for i in reversed(range(len(angulos))):
1760 1760 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
1761 1761 diff_ = angulos[i]-angulos[i-1]
1762 1762 if abs(diff_) >1.5:
1763 1763 list1.append(i-1)
1764 1764 list2.append(diff_)
1765 1765 return list(reversed(list1)),list(reversed(list2))
1766 1766
1767 1767 def fixData90(self,list_,ang_):
1768 1768 if list_[0]==-1:
1769 1769 vec = numpy.where(ang_<ang_[0])
1770 1770 ang_[vec] = ang_[vec]+90
1771 1771 return ang_
1772 1772 return ang_
1773 1773
1774 1774 def fixData90HL(self,angulos):
1775 1775 vec = numpy.where(angulos>=90)
1776 1776 angulos[vec]=angulos[vec]-90
1777 1777 return angulos
1778 1778
1779 1779
1780 1780 def search_pos(self,pos,list_):
1781 1781 for i in range(len(list_)):
1782 1782 if pos == list_[i]:
1783 1783 return True,i
1784 1784 i=None
1785 1785 return False,i
1786 1786
1787 1787 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
1788 1788 size = len(ang_)
1789 1789 size2 = 0
1790 1790 for i in range(len(list2_)):
1791 1791 size2=size2+round(abs(list2_[i]))-1
1792 1792 new_size= size+size2
1793 1793 ang_new = numpy.zeros(new_size)
1794 1794 ang_new2 = numpy.zeros(new_size)
1795 1795
1796 1796 tmp = 0
1797 1797 c = 0
1798 1798 for i in range(len(ang_)):
1799 1799 ang_new[tmp +c] = ang_[i]
1800 1800 ang_new2[tmp+c] = ang_[i]
1801 1801 condition , value = self.search_pos(i,list1_)
1802 1802 if condition:
1803 1803 pos = tmp + c + 1
1804 1804 for k in range(round(abs(list2_[value]))-1):
1805 1805 if tipo_case==0 or tipo_case==3:#subida
1806 1806 ang_new[pos+k] = ang_new[pos+k-1]+1
1807 1807 ang_new2[pos+k] = numpy.nan
1808 1808 elif tipo_case==1 or tipo_case==2:#bajada
1809 1809 ang_new[pos+k] = ang_new[pos+k-1]-1
1810 1810 ang_new2[pos+k] = numpy.nan
1811 1811
1812 1812 tmp = pos +k
1813 1813 c = 0
1814 1814 c=c+1
1815 1815 return ang_new,ang_new2
1816 1816
1817 1817 def globalCheckPED(self,angulos,tipo_case):
1818 1818 l1,l2 = self.get2List(angulos)
1819 1819 print("l1",l1)
1820 1820 print("l2",l2)
1821 1821 if len(l1)>0:
1822 1822 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
1823 1823 #l1,l2 = self.get2List(angulos2)
1824 1824 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
1825 1825 #ang1_ = self.fixData90HL(ang1_)
1826 1826 #ang2_ = self.fixData90HL(ang2_)
1827 1827 else:
1828 1828 ang1_= angulos
1829 1829 ang2_= angulos
1830 1830 return ang1_,ang2_
1831 1831
1832 1832
1833 1833 def replaceNAN(self,data_weather,data_ele,val):
1834 1834 data= data_ele
1835 1835 data_T= data_weather
1836 1836 #print(data.shape[0])
1837 1837 #print(data_T.shape[0])
1838 1838 #exit(1)
1839 1839 if data.shape[0]> data_T.shape[0]:
1840 1840 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
1841 1841 c = 0
1842 1842 for i in range(len(data)):
1843 1843 if numpy.isnan(data[i]):
1844 1844 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1845 1845 else:
1846 1846 data_N[i,:]=data_T[c,:]
1847 1847 c=c+1
1848 1848 return data_N
1849 1849 else:
1850 1850 for i in range(len(data)):
1851 1851 if numpy.isnan(data[i]):
1852 1852 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1853 1853 return data_T
1854 1854
1855 1855
1856 1856 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min,case_flag):
1857 1857 ang_max= ang_max
1858 1858 ang_min= ang_min
1859 1859 data_weather=data_weather
1860 1860 val_ch=val_ch
1861 1861 ##print("*********************DATA WEATHER**************************************")
1862 1862 ##print(data_weather)
1863 1863
1864 1864 '''
1865 1865 print("**********************************************")
1866 1866 print("**********************************************")
1867 1867 print("***************ini**************")
1868 1868 print("**********************************************")
1869 1869 print("**********************************************")
1870 1870 '''
1871 1871 #print("data_ele",data_ele)
1872 1872 #----------------------------------------------------------
1873 1873
1874 1874 #exit(1)
1875 1875 tipo_case = case_flag[-1]
1876 1876 print("tipo_case",tipo_case)
1877 1877 #--------------------- new -------------------------
1878 1878 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
1879 1879
1880 1880 #-------------------------CAMBIOS RHI---------------------------------
1881 1881
1882 1882 vec = numpy.where(data_ele<ang_max)
1883 1883 data_ele = data_ele[vec]
1884 1884 data_weather= data_weather[vec[0]]
1885 1885
1886 1886 len_vec = len(vec)
1887 1887 data_ele_new = data_ele[::-1] # reversa
1888 1888 data_weather = data_weather[::-1,:]
1889 1889 new_i_ele = int(data_ele_new[0])
1890 1890 new_f_ele = int(data_ele_new[-1])
1891 1891
1892 1892 n1= new_i_ele- ang_min
1893 1893 n2= ang_max - new_f_ele-1
1894 1894 if n1>0:
1895 1895 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1896 1896 ele1_nan= numpy.ones(n1)*numpy.nan
1897 1897 data_ele = numpy.hstack((ele1,data_ele_new))
1898 1898 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1899 1899 if n2>0:
1900 1900 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1901 1901 ele2_nan= numpy.ones(n2)*numpy.nan
1902 1902 data_ele = numpy.hstack((data_ele,ele2))
1903 1903 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1904 1904
1905 1905
1906 1906 print("ele shape",data_ele.shape)
1907 1907 print(data_ele)
1908 1908
1909 1909 #print("self.data_ele_tmp",self.data_ele_tmp)
1910 1910 val_mean = numpy.mean(data_weather[:,-1])
1911 1911 self.val_mean = val_mean
1912 1912 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1913 1913 self.data_ele_tmp[val_ch]= data_ele_old
1914 1914
1915 1915
1916 1916 print("data_weather shape",data_weather.shape)
1917 1917 print(data_weather)
1918 1918 #exit(1)
1919 1919 return data_weather,data_ele
1920 1920
1921 1921
1922 1922 def plot(self):
1923 1923 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
1924 1924 data = self.data[-1]
1925 1925 r = self.data.yrange
1926 1926 delta_height = r[1]-r[0]
1927 1927 r_mask = numpy.where(r>=0)[0]
1928 1928 ##print("delta_height",delta_height)
1929 1929 #print("r_mask",r_mask,len(r_mask))
1930 1930 r = numpy.arange(len(r_mask))*delta_height
1931 1931 self.y = 2*r
1932 1932 res = 1
1933 1933 ###print("data['weather'].shape[0]",data['weather'].shape[0])
1934 1934 ang_max = self.ang_max
1935 1935 ang_min = self.ang_min
1936 1936 var_ang =ang_max - ang_min
1937 1937 step = (int(var_ang)/(res*data['weather'].shape[0]))
1938 1938 ###print("step",step)
1939 1939 #--------------------------------------------------------
1940 1940 ##print('weather',data['weather'].shape)
1941 1941 ##print('ele',data['ele'].shape)
1942 1942
1943 1943 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1944 1944 ###self.res_azi = numpy.mean(data['azi'])
1945 1945 ###print("self.res_ele",self.res_ele)
1946 1946 plt.clf()
1947 1947 subplots = [121, 122]
1948 1948 if self.ini==0:
1949 1949 self.data_ele_tmp = numpy.ones([self.nplots,int(var_ang)])*numpy.nan
1950 1950 self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
1951 1951 print("SHAPE",self.data_ele_tmp.shape)
1952 1952
1953 1953 for i,ax in enumerate(self.axes):
1954 1954 self.res_weather[i], self.res_ele = self.const_ploteo(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min,case_flag=self.data['case_flag'])
1955 1955 self.res_azi = numpy.mean(data['azi'])
1956 1956
1957 1957 print(self.res_ele)
1958 1958 #exit(1)
1959 1959 if ax.firsttime:
1960 1960 #plt.clf()
1961 1961 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1962 1962 #fig=self.figures[0]
1963 1963 else:
1964 1964
1965 1965 #plt.clf()
1966 1966 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1967 1967 caax = cgax.parasites[0]
1968 1968 paax = cgax.parasites[1]
1969 1969 cbar = plt.gcf().colorbar(pm, pad=0.075)
1970 1970 caax.set_xlabel('x_range [km]')
1971 1971 caax.set_ylabel('y_range [km]')
1972 1972 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
1973 1973 print("***************************self.ini****************************",self.ini)
1974 1974 self.ini= self.ini+1
1975 1975
1976 1976 class WeatherRHI_vRF3_Plot(Plot):
1977 1977 CODE = 'weather'
1978 1978 plot_name = 'weather'
1979 1979 plot_type = 'rhistyle'
1980 1980 buffering = False
1981 1981 data_ele_tmp = None
1982 1982
1983 1983 def setup(self):
1984 1984 print("********************")
1985 1985 print("********************")
1986 1986 print("********************")
1987 1987 print("SETUP WEATHER PLOT")
1988 1988 self.ncols = 1
1989 1989 self.nrows = 1
1990 1990 self.nplots= 1
1991 1991 self.ylabel= 'Range [Km]'
1992 1992 self.titles= ['Weather']
1993 1993 if self.channels is not None:
1994 1994 self.nplots = len(self.channels)
1995 1995 self.nrows = len(self.channels)
1996 1996 else:
1997 1997 self.nplots = self.data.shape(self.CODE)[0]
1998 1998 self.nrows = self.nplots
1999 1999 self.channels = list(range(self.nplots))
2000 2000 print("channels",self.channels)
2001 2001 print("que saldra", self.data.shape(self.CODE)[0])
2002 2002 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
2003 2003 print("self.titles",self.titles)
2004 2004 self.colorbar=False
2005 2005 self.width =8
2006 2006 self.height =8
2007 2007 self.ini =0
2008 2008 self.len_azi =0
2009 2009 self.buffer_ini = None
2010 2010 self.buffer_ele = None
2011 2011 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
2012 2012 self.flag =0
2013 2013 self.indicador= 0
2014 2014 self.last_data_ele = None
2015 2015 self.val_mean = None
2016 2016
2017 2017 def update(self, dataOut):
2018 2018
2019 2019 data = {}
2020 2020 meta = {}
2021 2021 if hasattr(dataOut, 'dataPP_POWER'):
2022 2022 factor = 1
2023 2023 if hasattr(dataOut, 'nFFTPoints'):
2024 2024 factor = dataOut.normFactor
2025 2025 print("dataOut",dataOut.data_360.shape)
2026 2026 #
2027 2027 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
2028 2028 #
2029 2029 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
2030 2030 data['azi'] = dataOut.data_azi
2031 2031 data['ele'] = dataOut.data_ele
2032 2032 #data['case_flag'] = dataOut.case_flag
2033 2033 #print("UPDATE")
2034 2034 #print("data[weather]",data['weather'].shape)
2035 2035 #print("data[azi]",data['azi'])
2036 2036 return data, meta
2037 2037
2038 2038 def get2List(self,angulos):
2039 2039 list1=[]
2040 2040 list2=[]
2041 2041 for i in reversed(range(len(angulos))):
2042 2042 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
2043 2043 diff_ = angulos[i]-angulos[i-1]
2044 2044 if abs(diff_) >1.5:
2045 2045 list1.append(i-1)
2046 2046 list2.append(diff_)
2047 2047 return list(reversed(list1)),list(reversed(list2))
2048 2048
2049 2049 def fixData90(self,list_,ang_):
2050 2050 if list_[0]==-1:
2051 2051 vec = numpy.where(ang_<ang_[0])
2052 2052 ang_[vec] = ang_[vec]+90
2053 2053 return ang_
2054 2054 return ang_
2055 2055
2056 2056 def fixData90HL(self,angulos):
2057 2057 vec = numpy.where(angulos>=90)
2058 2058 angulos[vec]=angulos[vec]-90
2059 2059 return angulos
2060 2060
2061 2061
2062 2062 def search_pos(self,pos,list_):
2063 2063 for i in range(len(list_)):
2064 2064 if pos == list_[i]:
2065 2065 return True,i
2066 2066 i=None
2067 2067 return False,i
2068 2068
2069 2069 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
2070 2070 size = len(ang_)
2071 2071 size2 = 0
2072 2072 for i in range(len(list2_)):
2073 2073 size2=size2+round(abs(list2_[i]))-1
2074 2074 new_size= size+size2
2075 2075 ang_new = numpy.zeros(new_size)
2076 2076 ang_new2 = numpy.zeros(new_size)
2077 2077
2078 2078 tmp = 0
2079 2079 c = 0
2080 2080 for i in range(len(ang_)):
2081 2081 ang_new[tmp +c] = ang_[i]
2082 2082 ang_new2[tmp+c] = ang_[i]
2083 2083 condition , value = self.search_pos(i,list1_)
2084 2084 if condition:
2085 2085 pos = tmp + c + 1
2086 2086 for k in range(round(abs(list2_[value]))-1):
2087 2087 if tipo_case==0 or tipo_case==3:#subida
2088 2088 ang_new[pos+k] = ang_new[pos+k-1]+1
2089 2089 ang_new2[pos+k] = numpy.nan
2090 2090 elif tipo_case==1 or tipo_case==2:#bajada
2091 2091 ang_new[pos+k] = ang_new[pos+k-1]-1
2092 2092 ang_new2[pos+k] = numpy.nan
2093 2093
2094 2094 tmp = pos +k
2095 2095 c = 0
2096 2096 c=c+1
2097 2097 return ang_new,ang_new2
2098 2098
2099 2099 def globalCheckPED(self,angulos,tipo_case):
2100 2100 l1,l2 = self.get2List(angulos)
2101 2101 ##print("l1",l1)
2102 2102 ##print("l2",l2)
2103 2103 if len(l1)>0:
2104 2104 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
2105 2105 #l1,l2 = self.get2List(angulos2)
2106 2106 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
2107 2107 #ang1_ = self.fixData90HL(ang1_)
2108 2108 #ang2_ = self.fixData90HL(ang2_)
2109 2109 else:
2110 2110 ang1_= angulos
2111 2111 ang2_= angulos
2112 2112 return ang1_,ang2_
2113 2113
2114 2114
2115 2115 def replaceNAN(self,data_weather,data_ele,val):
2116 2116 data= data_ele
2117 2117 data_T= data_weather
2118 2118
2119 2119 if data.shape[0]> data_T.shape[0]:
2120 2120 print("IF")
2121 2121 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
2122 2122 c = 0
2123 2123 for i in range(len(data)):
2124 2124 if numpy.isnan(data[i]):
2125 2125 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
2126 2126 else:
2127 2127 data_N[i,:]=data_T[c,:]
2128 2128 c=c+1
2129 2129 return data_N
2130 2130 else:
2131 2131 print("else")
2132 2132 for i in range(len(data)):
2133 2133 if numpy.isnan(data[i]):
2134 2134 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
2135 2135 return data_T
2136 2136
2137 2137 def check_case(self,data_ele,ang_max,ang_min):
2138 2138 start = data_ele[0]
2139 2139 end = data_ele[-1]
2140 2140 number = (end-start)
2141 2141 len_ang=len(data_ele)
2142 2142 print("start",start)
2143 2143 print("end",end)
2144 2144 print("number",number)
2145 2145
2146 2146 print("len_ang",len_ang)
2147 2147
2148 2148 #exit(1)
2149 2149
2150 2150 if start<end and (round(abs(number)+1)>=len_ang or (numpy.argmin(data_ele)==0)):#caso subida
2151 2151 return 0
2152 2152 #elif start>end and (round(abs(number)+1)>=len_ang or(numpy.argmax(data_ele)==0)):#caso bajada
2153 2153 # return 1
2154 2154 elif round(abs(number)+1)>=len_ang and (start>end or(numpy.argmax(data_ele)==0)):#caso bajada
2155 2155 return 1
2156 2156 elif round(abs(number)+1)<len_ang and data_ele[-2]>data_ele[-1]:# caso BAJADA CAMBIO ANG MAX
2157 2157 return 2
2158 2158 elif round(abs(number)+1)<len_ang and data_ele[-2]<data_ele[-1] :# caso SUBIDA CAMBIO ANG MIN
2159 2159 return 3
2160 2160
2161 2161
2162 2162 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min,case_flag):
2163 2163 ang_max= ang_max
2164 2164 ang_min= ang_min
2165 2165 data_weather=data_weather
2166 2166 val_ch=val_ch
2167 2167 ##print("*********************DATA WEATHER**************************************")
2168 2168 ##print(data_weather)
2169 2169 if self.ini==0:
2170 2170
2171 2171 #--------------------- new -------------------------
2172 2172 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
2173 2173
2174 2174 #-------------------------CAMBIOS RHI---------------------------------
2175 2175 start= ang_min
2176 2176 end = ang_max
2177 2177 n= (ang_max-ang_min)/res
2178 2178 #------ new
2179 2179 self.start_data_ele = data_ele_new[0]
2180 2180 self.end_data_ele = data_ele_new[-1]
2181 2181 if tipo_case==0 or tipo_case==3: # SUBIDA
2182 2182 n1= round(self.start_data_ele)- start
2183 2183 n2= end - round(self.end_data_ele)
2184 2184 print(self.start_data_ele)
2185 2185 print(self.end_data_ele)
2186 2186 if n1>0:
2187 2187 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
2188 2188 ele1_nan= numpy.ones(n1)*numpy.nan
2189 2189 data_ele = numpy.hstack((ele1,data_ele_new))
2190 2190 print("ele1_nan",ele1_nan.shape)
2191 2191 print("data_ele_old",data_ele_old.shape)
2192 2192 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
2193 2193 if n2>0:
2194 2194 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
2195 2195 ele2_nan= numpy.ones(n2)*numpy.nan
2196 2196 data_ele = numpy.hstack((data_ele,ele2))
2197 2197 print("ele2_nan",ele2_nan.shape)
2198 2198 print("data_ele_old",data_ele_old.shape)
2199 2199 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
2200 2200
2201 2201 if tipo_case==1 or tipo_case==2: # BAJADA
2202 2202 data_ele_new = data_ele_new[::-1] # reversa
2203 2203 data_ele_old = data_ele_old[::-1]# reversa
2204 2204 data_weather = data_weather[::-1,:]# reversa
2205 2205 vec= numpy.where(data_ele_new<ang_max)
2206 2206 data_ele_new = data_ele_new[vec]
2207 2207 data_ele_old = data_ele_old[vec]
2208 2208 data_weather = data_weather[vec[0]]
2209 2209 vec2= numpy.where(0<data_ele_new)
2210 2210 data_ele_new = data_ele_new[vec2]
2211 2211 data_ele_old = data_ele_old[vec2]
2212 2212 data_weather = data_weather[vec2[0]]
2213 2213 self.start_data_ele = data_ele_new[0]
2214 2214 self.end_data_ele = data_ele_new[-1]
2215 2215
2216 2216 n1= round(self.start_data_ele)- start
2217 2217 n2= end - round(self.end_data_ele)-1
2218 2218 print(self.start_data_ele)
2219 2219 print(self.end_data_ele)
2220 2220 if n1>0:
2221 2221 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
2222 2222 ele1_nan= numpy.ones(n1)*numpy.nan
2223 2223 data_ele = numpy.hstack((ele1,data_ele_new))
2224 2224 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
2225 2225 if n2>0:
2226 2226 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
2227 2227 ele2_nan= numpy.ones(n2)*numpy.nan
2228 2228 data_ele = numpy.hstack((data_ele,ele2))
2229 2229 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
2230 2230 # RADAR
2231 2231 # NOTA data_ele y data_weather es la variable que retorna
2232 2232 val_mean = numpy.mean(data_weather[:,-1])
2233 2233 self.val_mean = val_mean
2234 2234 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
2235 2235 print("eleold",data_ele_old)
2236 2236 print(self.data_ele_tmp[val_ch])
2237 2237 print(data_ele_old.shape[0])
2238 2238 print(self.data_ele_tmp[val_ch].shape[0])
2239 2239 if (data_ele_old.shape[0]==91 or self.data_ele_tmp[val_ch].shape[0]==91):
2240 2240 import sys
2241 2241 print("EXIT",self.ini)
2242 2242
2243 2243 sys.exit(1)
2244 2244 self.data_ele_tmp[val_ch]= data_ele_old
2245 2245 else:
2246 2246 #print("**********************************************")
2247 2247 #print("****************VARIABLE**********************")
2248 2248 #-------------------------CAMBIOS RHI---------------------------------
2249 2249 #---------------------------------------------------------------------
2250 2250 ##print("INPUT data_ele",data_ele)
2251 2251 flag=0
2252 2252 start_ele = self.res_ele[0]
2253 2253 #tipo_case = self.check_case(data_ele,ang_max,ang_min)
2254 2254 tipo_case = case_flag[-1]
2255 2255 #print("TIPO DE DATA",tipo_case)
2256 2256 #-----------new------------
2257 2257 data_ele ,data_ele_old = self.globalCheckPED(data_ele,tipo_case)
2258 2258 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
2259 2259
2260 2260 #-------------------------------NEW RHI ITERATIVO-------------------------
2261 2261
2262 2262 if tipo_case==0 : # SUBIDA
2263 2263 vec = numpy.where(data_ele<ang_max)
2264 2264 data_ele = data_ele[vec]
2265 2265 data_ele_old = data_ele_old[vec]
2266 2266 data_weather = data_weather[vec[0]]
2267 2267
2268 2268 vec2 = numpy.where(0<data_ele)
2269 2269 data_ele= data_ele[vec2]
2270 2270 data_ele_old= data_ele_old[vec2]
2271 2271 ##print(data_ele_new)
2272 2272 data_weather= data_weather[vec2[0]]
2273 2273
2274 2274 new_i_ele = int(round(data_ele[0]))
2275 2275 new_f_ele = int(round(data_ele[-1]))
2276 2276 #print(new_i_ele)
2277 2277 #print(new_f_ele)
2278 2278 #print(data_ele,len(data_ele))
2279 2279 #print(data_ele_old,len(data_ele_old))
2280 2280 if new_i_ele< 2:
2281 2281 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
2282 2282 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
2283 2283 self.data_ele_tmp[val_ch][new_i_ele:new_i_ele+len(data_ele)]=data_ele_old
2284 2284 self.res_ele[new_i_ele:new_i_ele+len(data_ele)]= data_ele
2285 2285 self.res_weather[val_ch][new_i_ele:new_i_ele+len(data_ele),:]= data_weather
2286 2286 data_ele = self.res_ele
2287 2287 data_weather = self.res_weather[val_ch]
2288 2288
2289 2289 elif tipo_case==1 : #BAJADA
2290 2290 data_ele = data_ele[::-1] # reversa
2291 2291 data_ele_old = data_ele_old[::-1]# reversa
2292 2292 data_weather = data_weather[::-1,:]# reversa
2293 2293 vec= numpy.where(data_ele<ang_max)
2294 2294 data_ele = data_ele[vec]
2295 2295 data_ele_old = data_ele_old[vec]
2296 2296 data_weather = data_weather[vec[0]]
2297 2297 vec2= numpy.where(0<data_ele)
2298 2298 data_ele = data_ele[vec2]
2299 2299 data_ele_old = data_ele_old[vec2]
2300 2300 data_weather = data_weather[vec2[0]]
2301 2301
2302 2302
2303 2303 new_i_ele = int(round(data_ele[0]))
2304 2304 new_f_ele = int(round(data_ele[-1]))
2305 2305 #print(data_ele)
2306 2306 #print(ang_max)
2307 2307 #print(data_ele_old)
2308 2308 if new_i_ele <= 1:
2309 2309 new_i_ele = 1
2310 2310 if round(data_ele[-1])>=ang_max-1:
2311 2311 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
2312 2312 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
2313 2313 self.data_ele_tmp[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1]=data_ele_old
2314 2314 self.res_ele[new_i_ele-1:new_i_ele+len(data_ele)-1]= data_ele
2315 2315 self.res_weather[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1,:]= data_weather
2316 2316 data_ele = self.res_ele
2317 2317 data_weather = self.res_weather[val_ch]
2318 2318
2319 2319 elif tipo_case==2: #bajada
2320 2320 vec = numpy.where(data_ele<ang_max)
2321 2321 data_ele = data_ele[vec]
2322 2322 data_weather= data_weather[vec[0]]
2323 2323
2324 2324 len_vec = len(vec)
2325 2325 data_ele_new = data_ele[::-1] # reversa
2326 2326 data_weather = data_weather[::-1,:]
2327 2327 new_i_ele = int(data_ele_new[0])
2328 2328 new_f_ele = int(data_ele_new[-1])
2329 2329
2330 2330 n1= new_i_ele- ang_min
2331 2331 n2= ang_max - new_f_ele-1
2332 2332 if n1>0:
2333 2333 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
2334 2334 ele1_nan= numpy.ones(n1)*numpy.nan
2335 2335 data_ele = numpy.hstack((ele1,data_ele_new))
2336 2336 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
2337 2337 if n2>0:
2338 2338 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
2339 2339 ele2_nan= numpy.ones(n2)*numpy.nan
2340 2340 data_ele = numpy.hstack((data_ele,ele2))
2341 2341 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
2342 2342
2343 2343 self.data_ele_tmp[val_ch] = data_ele_old
2344 2344 self.res_ele = data_ele
2345 2345 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
2346 2346 data_ele = self.res_ele
2347 2347 data_weather = self.res_weather[val_ch]
2348 2348
2349 2349 elif tipo_case==3:#subida
2350 2350 vec = numpy.where(0<data_ele)
2351 2351 data_ele= data_ele[vec]
2352 2352 data_ele_new = data_ele
2353 2353 data_ele_old= data_ele_old[vec]
2354 2354 data_weather= data_weather[vec[0]]
2355 2355 pos_ini = numpy.argmin(data_ele)
2356 2356 if pos_ini>0:
2357 2357 len_vec= len(data_ele)
2358 2358 vec3 = numpy.linspace(pos_ini,len_vec-1,len_vec-pos_ini).astype(int)
2359 2359 #print(vec3)
2360 2360 data_ele= data_ele[vec3]
2361 2361 data_ele_new = data_ele
2362 2362 data_ele_old= data_ele_old[vec3]
2363 2363 data_weather= data_weather[vec3]
2364 2364
2365 2365 new_i_ele = int(data_ele_new[0])
2366 2366 new_f_ele = int(data_ele_new[-1])
2367 2367 n1= new_i_ele- ang_min
2368 2368 n2= ang_max - new_f_ele-1
2369 2369 if n1>0:
2370 2370 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
2371 2371 ele1_nan= numpy.ones(n1)*numpy.nan
2372 2372 data_ele = numpy.hstack((ele1,data_ele_new))
2373 2373 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
2374 2374 if n2>0:
2375 2375 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
2376 2376 ele2_nan= numpy.ones(n2)*numpy.nan
2377 2377 data_ele = numpy.hstack((data_ele,ele2))
2378 2378 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
2379 2379
2380 2380 self.data_ele_tmp[val_ch] = data_ele_old
2381 2381 self.res_ele = data_ele
2382 2382 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
2383 2383 data_ele = self.res_ele
2384 2384 data_weather = self.res_weather[val_ch]
2385 2385 #print("self.data_ele_tmp",self.data_ele_tmp)
2386 2386 return data_weather,data_ele
2387 2387
2388 2388 def const_ploteo_vRF(self,val_ch,data_weather,data_ele,res,ang_max,ang_min):
2389 2389
2390 2390 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,1)
2391 2391
2392 2392 data_ele = data_ele_old.copy()
2393 2393
2394 2394 diff_1 = ang_max - data_ele[0]
2395 2395 angles_1_nan = numpy.linspace(ang_max,data_ele[0]+1,int(diff_1)-1)#*numpy.nan
2396 2396
2397 2397 diff_2 = data_ele[-1]-ang_min
2398 2398 angles_2_nan = numpy.linspace(data_ele[-1]-1,ang_min,int(diff_2)-1)#*numpy.nan
2399 2399
2400 2400 angles_filled = numpy.concatenate((angles_1_nan,data_ele,angles_2_nan))
2401 2401
2402 2402 print(angles_filled)
2403 2403
2404 2404 data_1_nan = numpy.ones([angles_1_nan.shape[0],len(self.r_mask)])*numpy.nan
2405 2405 data_2_nan = numpy.ones([angles_2_nan.shape[0],len(self.r_mask)])*numpy.nan
2406 2406
2407 2407 data_filled = numpy.concatenate((data_1_nan,data_weather,data_2_nan),axis=0)
2408 2408 #val_mean = numpy.mean(data_weather[:,-1])
2409 2409 #self.val_mean = val_mean
2410 2410 print(data_filled)
2411 2411 data_filled = self.replaceNAN(data_weather=data_filled,data_ele=angles_filled,val=numpy.nan)
2412 2412
2413 2413 print(data_filled)
2414 2414 print(data_filled.shape)
2415 2415 print(angles_filled.shape)
2416 2416
2417 2417 return data_filled,angles_filled
2418 2418
2419 2419 def plot(self):
2420 2420 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
2421 2421 data = self.data[-1]
2422 2422 r = self.data.yrange
2423 2423 delta_height = r[1]-r[0]
2424 2424 r_mask = numpy.where(r>=0)[0]
2425 2425 self.r_mask =r_mask
2426 2426 ##print("delta_height",delta_height)
2427 2427 #print("r_mask",r_mask,len(r_mask))
2428 2428 r = numpy.arange(len(r_mask))*delta_height
2429 2429 self.y = 2*r
2430 2430 res = 1
2431 2431 ###print("data['weather'].shape[0]",data['weather'].shape[0])
2432 2432 ang_max = self.ang_max
2433 2433 ang_min = self.ang_min
2434 2434 var_ang =ang_max - ang_min
2435 2435 step = (int(var_ang)/(res*data['weather'].shape[0]))
2436 2436 ###print("step",step)
2437 2437 #--------------------------------------------------------
2438 2438 ##print('weather',data['weather'].shape)
2439 2439 ##print('ele',data['ele'].shape)
2440 2440
2441 2441 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
2442 2442 ###self.res_azi = numpy.mean(data['azi'])
2443 2443 ###print("self.res_ele",self.res_ele)
2444 2444
2445 2445 plt.clf()
2446 2446 subplots = [121, 122]
2447 2447 #if self.ini==0:
2448 2448 #self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
2449 2449 #print("SHAPE",self.data_ele_tmp.shape)
2450 2450
2451 2451 for i,ax in enumerate(self.axes):
2452 2452 res_weather, self.res_ele = self.const_ploteo_vRF(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],res=res,ang_max=ang_max,ang_min=ang_min)
2453 2453 self.res_azi = numpy.mean(data['azi'])
2454 2454
2455 2455 if ax.firsttime:
2456 2456 #plt.clf()
2457 2457 print("Frist Plot")
2458 2458 print(data['weather'][i][:,r_mask].shape)
2459 2459 print(data['ele'].shape)
2460 2460 cgax, pm = wrl.vis.plot_rhi(res_weather,r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
2461 2461 #cgax, pm = wrl.vis.plot_rhi(data['weather'][i][:,r_mask],r=r,th=data['ele'],ax=subplots[i], proj='cg',vmin=20, vmax=80)
2462 2462 gh = cgax.get_grid_helper()
2463 2463 locs = numpy.linspace(ang_min,ang_max,var_ang+1)
2464 2464 gh.grid_finder.grid_locator1 = FixedLocator(locs)
2465 2465 gh.grid_finder.tick_formatter1 = DictFormatter(dict([(i, r"${0:.0f}^\circ$".format(i)) for i in locs]))
2466 2466
2467 2467
2468 2468 #fig=self.figures[0]
2469 2469 else:
2470 2470 #plt.clf()
2471 2471 print("ELSE PLOT")
2472 2472 cgax, pm = wrl.vis.plot_rhi(res_weather,r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
2473 2473 #cgax, pm = wrl.vis.plot_rhi(data['weather'][i][:,r_mask],r=r,th=data['ele'],ax=subplots[i], proj='cg',vmin=20, vmax=80)
2474 2474 gh = cgax.get_grid_helper()
2475 2475 locs = numpy.linspace(ang_min,ang_max,var_ang+1)
2476 2476 gh.grid_finder.grid_locator1 = FixedLocator(locs)
2477 2477 gh.grid_finder.tick_formatter1 = DictFormatter(dict([(i, r"${0:.0f}^\circ$".format(i)) for i in locs]))
2478 2478
2479 2479 caax = cgax.parasites[0]
2480 2480 paax = cgax.parasites[1]
2481 2481 cbar = plt.gcf().colorbar(pm, pad=0.075)
2482 2482 caax.set_xlabel('x_range [km]')
2483 2483 caax.set_ylabel('y_range [km]')
2484 2484 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
2485 2485 print("***************************self.ini****************************",self.ini)
2486 2486 self.ini= self.ini+1
2487 2487
2488 2488 class WeatherRHI_vRF4_Plot(Plot):
2489 2489 CODE = 'RHI'
2490 2490 plot_name = 'RHI'
2491 2491 #plot_type = 'rhistyle'
2492 2492 buffering = False
2493 2493
2494 2494 def setup(self):
2495 2495
2496 2496 self.ncols = 1
2497 2497 self.nrows = 1
2498 2498 self.nplots= 1
2499 2499 self.ylabel= 'Range [Km]'
2500 2500 self.titles= ['RHI']
2501 2501 self.polar = True
2502 self.grid = True
2502 2503 if self.channels is not None:
2503 2504 self.nplots = len(self.channels)
2504 2505 self.nrows = len(self.channels)
2505 2506 else:
2506 2507 self.nplots = self.data.shape(self.CODE)[0]
2507 2508 self.nrows = self.nplots
2508 2509 self.channels = list(range(self.nplots))
2509 2510
2510 2511 if self.CODE == 'Power':
2511 2512 self.cb_label = r'Power (dB)'
2512 2513 elif self.CODE == 'Doppler':
2513 2514 self.cb_label = r'Velocity (m/s)'
2514 2515 self.colorbar=True
2515 2516 self.width =8
2516 2517 self.height =8
2517 2518 self.ini =0
2518 2519 self.len_azi =0
2519 2520 self.buffer_ini = None
2520 2521 self.buffer_ele = None
2521 2522 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
2522 2523 self.flag =0
2523 2524 self.indicador= 0
2524 2525 self.last_data_ele = None
2525 2526 self.val_mean = None
2526 2527
2527 2528 def update(self, dataOut):
2528 2529
2529 2530 data = {}
2530 2531 meta = {}
2531 2532 if hasattr(dataOut, 'dataPP_POWER'):
2532 2533 factor = 1
2533 2534 if hasattr(dataOut, 'nFFTPoints'):
2534 2535 factor = dataOut.normFactor
2535 2536
2536 2537 if 'pow' in self.attr_data[0].lower():
2537 2538 data['data'] = 10*numpy.log10(getattr(dataOut, self.attr_data[0])/(factor))
2538 2539 else:
2539 2540 data['data'] = getattr(dataOut, self.attr_data[0])/(factor)
2540 2541
2541 2542 data['azi'] = dataOut.data_azi
2542 2543 data['ele'] = dataOut.data_ele
2543 2544
2544 2545 return data, meta
2545 2546
2546 2547 def plot(self):
2547 2548 data = self.data[-1]
2548 2549 r = self.data.yrange
2549 2550 delta_height = r[1]-r[0]
2550 2551 r_mask = numpy.where(r>=0)[0]
2551 2552 self.r_mask =r_mask
2552 2553 r = numpy.arange(len(r_mask))*delta_height
2553 2554 self.y = 2*r
2554 2555 res = 1
2555 ang_max = self.ang_max
2556 ang_min = self.ang_min
2557 var_ang =ang_max - ang_min
2558 step = (int(var_ang)/(res*data['data'].shape[0]))
2556 #ang_max = self.ang_max
2557 #ang_min = self.ang_min
2558 #var_ang =ang_max - ang_min
2559 #step = (int(var_ang)/(res*data['data'].shape[0]))
2559 2560
2560 2561 z = data['data'][self.channels[0]][:,r_mask]
2561 2562
2562 2563 self.titles = []
2563 2564
2564 2565 self.ymax = self.ymax if self.ymax else numpy.nanmax(r)
2565 2566 self.ymin = self.ymin if self.ymin else numpy.nanmin(r)
2566 2567 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
2567 2568 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
2568 2569 self.ang_min = self.ang_min if self.ang_min else 0
2569 self.ang_max = self.ang_max if self.ang_max else 2*numpy.pi
2570 self.ang_max = self.ang_max if self.ang_max else 90
2570 2571
2571 2572 subplots = [121, 122]
2572 2573
2573 2574 r, theta = numpy.meshgrid(r, numpy.radians(data['ele']) )
2574 2575
2575 2576 for i,ax in enumerate(self.axes):
2576 2577
2577 2578 if ax.firsttime:
2578 2579 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
2579 2580 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
2580 2581
2581 2582 else:
2582 2583 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
2583 2584 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
2584 2585
2585 2586 if len(self.channels) !=1:
2586 2587 self.titles = ['RHI {} AZ: {} Channel {}'.format(self.labels[x], str(round(numpy.mean(data['azi']),1)), x) for x in range(self.nrows)]
2587 2588 else:
2588 2589 self.titles = ['RHI {} AZ: {} Channel {}'.format(self.labels[0], str(round(numpy.mean(data['azi']),1)), self.channels[0])]
Show file before | Show file after | Hide comments
@@ -1,838 +1,834
1 1 '''
2 2 Created on Jul 3, 2014
3 3
4 4 @author: roj-idl71
5 5 '''
6 6 # SUBCHANNELS EN VEZ DE CHANNELS
7 7 # BENCHMARKS -> PROBLEMAS CON ARCHIVOS GRANDES -> INCONSTANTE EN EL TIEMPO
8 8 # ACTUALIZACION DE VERSION
9 9 # HEADERS
10 10 # MODULO DE ESCRITURA
11 11 # METADATA
12 12
13 13 import os
14 14 import time
15 15 import datetime
16 16 import numpy
17 17 import timeit
18 18 from fractions import Fraction
19 19 from time import time
20 20 from time import sleep
21 21
22 22 import schainpy.admin
23 23 from schainpy.model.data.jroheaderIO import RadarControllerHeader, SystemHeader
24 24 from schainpy.model.data.jrodata import Voltage
25 25 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
26 26
27 27 import pickle
28 28 try:
29 29 import digital_rf
30 30 except:
31 31 pass
32 32
33 33
34 34 class DigitalRFReader(ProcessingUnit):
35 35 '''
36 36 classdocs
37 37 '''
38 38
39 39 def __init__(self):
40 40 '''
41 41 Constructor
42 42 '''
43 43
44 44 ProcessingUnit.__init__(self)
45 45
46 46 self.dataOut = Voltage()
47 47 self.__printInfo = True
48 48 self.__flagDiscontinuousBlock = False
49 49 self.__bufferIndex = 9999999
50 50 self.__codeType = 0
51 51 self.__ippKm = None
52 52 self.__nCode = None
53 53 self.__nBaud = None
54 54 self.__code = None
55 55 self.dtype = None
56 56 self.oldAverage = None
57 57 self.path = None
58 58
59 59 def close(self):
60 60 print('Average of writing to digital rf format is ', self.oldAverage * 1000)
61 61 return
62 62
63 63 def __getCurrentSecond(self):
64 64
65 65 return self.__thisUnixSample / self.__sample_rate
66 66
67 67 thisSecond = property(__getCurrentSecond, "I'm the 'thisSecond' property.")
68 68
69 69 def __setFileHeader(self):
70 70 '''
71 71 In this method will be initialized every parameter of dataOut object (header, no data)
72 72 '''
73 73 ippSeconds = 1.0 * self.__nSamples / self.__sample_rate
74 74 if not self.getByBlock:
75 75 nProfiles = 1.0 / ippSeconds # Number of profiles in one second
76 76 else:
77 77 nProfiles = self.nProfileBlocks # Number of profiles in one block
78 78
79 79 try:
80 80 self.dataOut.radarControllerHeaderObj = RadarControllerHeader(
81 81 self.__radarControllerHeader)
82 82 except:
83 83 self.dataOut.radarControllerHeaderObj = RadarControllerHeader(
84 84 txA=0,
85 85 txB=0,
86 86 nWindows=1,
87 87 nHeights=self.__nSamples,
88 88 firstHeight=self.__firstHeigth,
89 89 deltaHeight=self.__deltaHeigth,
90 90 codeType=self.__codeType,
91 91 nCode=self.__nCode, nBaud=self.__nBaud,
92 92 code=self.__code)
93 93
94 94 try:
95 95 self.dataOut.systemHeaderObj = SystemHeader(self.__systemHeader)
96 96 except:
97 97 self.dataOut.systemHeaderObj = SystemHeader(nSamples=self.__nSamples,
98 98 nProfiles=nProfiles,
99 99 nChannels=len(
100 100 self.__channelList),
101 101 adcResolution=14)
102 102 self.dataOut.type = "Voltage"
103 103
104 104 self.dataOut.data = None
105 105
106 106 self.dataOut.dtype = self.dtype
107 107
108 108 # self.dataOut.nChannels = 0
109 109
110 110 # self.dataOut.nHeights = 0
111 111
112 112 self.dataOut.nProfiles = int(nProfiles)
113 113
114 114 self.dataOut.heightList = self.__firstHeigth + \
115 115 numpy.arange(self.__nSamples, dtype=numpy.float) * \
116 116 self.__deltaHeigth
117 117
118 118 #self.dataOut.channelList = list(range(self.__num_subchannels))
119 119 self.dataOut.channelList = list(range(len(self.__channelList)))
120 120 if not self.getByBlock:
121 121
122 122 self.dataOut.blocksize = self.dataOut.nChannels * self.dataOut.nHeights
123 123 else:
124 124 self.dataOut.blocksize = self.dataOut.nChannels * self.dataOut.nHeights*self.nProfileBlocks
125 125
126 126 # self.dataOut.channelIndexList = None
127 127
128 128 self.dataOut.flagNoData = True
129 129 if not self.getByBlock:
130 130 self.dataOut.flagDataAsBlock = False
131 131 else:
132 132 self.dataOut.flagDataAsBlock = True
133 133 # Set to TRUE if the data is discontinuous
134 134 self.dataOut.flagDiscontinuousBlock = False
135 135
136 136 self.dataOut.utctime = None
137 137
138 138 # timezone like jroheader, difference in minutes between UTC and localtime
139 139 self.dataOut.timeZone = self.__timezone / 60
140 140
141 141 self.dataOut.dstFlag = 0
142 142
143 143 self.dataOut.errorCount = 0
144 144
145 145 try:
146 146 self.dataOut.nCohInt = self.fixed_metadata_dict.get(
147 147 'nCohInt', self.nCohInt)
148 148
149 149 # asumo que la data esta decodificada
150 150 self.dataOut.flagDecodeData = self.fixed_metadata_dict.get(
151 151 'flagDecodeData', self.flagDecodeData)
152 152
153 153 # asumo que la data esta sin flip
154 154 self.dataOut.flagDeflipData = self.fixed_metadata_dict['flagDeflipData']
155 155
156 156 self.dataOut.flagShiftFFT = self.fixed_metadata_dict['flagShiftFFT']
157 157
158 158 self.dataOut.useLocalTime = self.fixed_metadata_dict['useLocalTime']
159 159 except:
160 160 pass
161 161
162 162 self.dataOut.ippSeconds = ippSeconds
163 163
164 164 # Time interval between profiles
165 165 # self.dataOut.timeInterval = self.dataOut.ippSeconds * self.dataOut.nCohInt
166 166
167 167 self.dataOut.frequency = self.__frequency
168 168
169 169 self.dataOut.realtime = self.__online
170 170
171 171 def findDatafiles(self, path, startDate=None, endDate=None):
172 172
173 173 if not os.path.isdir(path):
174 174 return []
175 175
176 176 try:
177 177 digitalReadObj = digital_rf.DigitalRFReader(
178 178 path, load_all_metadata=True)
179 179 except:
180 180 digitalReadObj = digital_rf.DigitalRFReader(path)
181 181
182 182 channelNameList = digitalReadObj.get_channels()
183 183
184 184 if not channelNameList:
185 185 return []
186 186
187 187 metadata_dict = digitalReadObj.get_rf_file_metadata(channelNameList[0])
188 188
189 189 sample_rate = metadata_dict['sample_rate'][0]
190 190
191 191 this_metadata_file = digitalReadObj.get_metadata(channelNameList[0])
192 192
193 193 try:
194 194 timezone = this_metadata_file['timezone'].value
195 195 except:
196 196 timezone = 0
197 197
198 198 startUTCSecond, endUTCSecond = digitalReadObj.get_bounds(
199 199 channelNameList[0]) / sample_rate - timezone
200 200
201 201 startDatetime = datetime.datetime.utcfromtimestamp(startUTCSecond)
202 202 endDatatime = datetime.datetime.utcfromtimestamp(endUTCSecond)
203 203
204 204 if not startDate:
205 205 startDate = startDatetime.date()
206 206
207 207 if not endDate:
208 208 endDate = endDatatime.date()
209 209
210 210 dateList = []
211 211
212 212 thisDatetime = startDatetime
213 213
214 214 while(thisDatetime <= endDatatime):
215 215
216 216 thisDate = thisDatetime.date()
217 217
218 218 if thisDate < startDate:
219 219 continue
220 220
221 221 if thisDate > endDate:
222 222 break
223 223
224 224 dateList.append(thisDate)
225 225 thisDatetime += datetime.timedelta(1)
226 226
227 227 return dateList
228 228
229 229 def setup(self, path=None,
230 230 startDate=None,
231 231 endDate=None,
232 232 startTime=datetime.time(0, 0, 0),
233 233 endTime=datetime.time(23, 59, 59),
234 234 channelList=None,
235 235 nSamples=None,
236 236 online=False,
237 237 delay=60,
238 238 buffer_size=1024,
239 239 ippKm=None,
240 240 nCohInt=1,
241 241 nCode=1,
242 242 nBaud=1,
243 243 flagDecodeData=False,
244 244 code=numpy.ones((1, 1), dtype=numpy.int),
245 245 getByBlock=0,
246 246 nProfileBlocks=1,
247 247 **kwargs):
248 248 '''
249 249 In this method we should set all initial parameters.
250 250
251 251 Inputs:
252 252 path
253 253 startDate
254 254 endDate
255 255 startTime
256 256 endTime
257 257 set
258 258 expLabel
259 259 ext
260 260 online
261 261 delay
262 262 '''
263 263 self.path = path
264 264 self.nCohInt = nCohInt
265 265 self.flagDecodeData = flagDecodeData
266 266 self.i = 0
267 267
268 268 self.getByBlock = getByBlock
269 269 self.nProfileBlocks = nProfileBlocks
270 270 if not os.path.isdir(path):
271 271 raise ValueError("[Reading] Directory %s does not exist" % path)
272 272
273 273 try:
274 274 self.digitalReadObj = digital_rf.DigitalRFReader(
275 275 path, load_all_metadata=True)
276 276 except:
277 277 self.digitalReadObj = digital_rf.DigitalRFReader(path)
278 278
279 279 channelNameList = self.digitalReadObj.get_channels()
280 280
281 281 if not channelNameList:
282 282 raise ValueError("[Reading] Directory %s does not have any files" % path)
283 283
284 284 if not channelList:
285 285 channelList = list(range(len(channelNameList)))
286 286
287 287 ########## Reading metadata ######################
288 288
289 289 top_properties = self.digitalReadObj.get_properties(
290 290 channelNameList[channelList[0]])
291 291
292 292 self.__num_subchannels = top_properties['num_subchannels']
293 293 self.__sample_rate = 1.0 * \
294 294 top_properties['sample_rate_numerator'] / \
295 295 top_properties['sample_rate_denominator']
296 296 # self.__samples_per_file = top_properties['samples_per_file'][0]
297 297 self.__deltaHeigth = 1e6 * 0.15 / self.__sample_rate # why 0.15?
298 298
299 299 this_metadata_file = self.digitalReadObj.get_digital_metadata(
300 300 channelNameList[channelList[0]])
301 301 metadata_bounds = this_metadata_file.get_bounds()
302 302 self.fixed_metadata_dict = this_metadata_file.read(
303 303 metadata_bounds[0])[metadata_bounds[0]] # GET FIRST HEADER
304 304
305 305 try:
306 306 self.__processingHeader = self.fixed_metadata_dict['processingHeader']
307 307 self.__radarControllerHeader = self.fixed_metadata_dict['radarControllerHeader']
308 308 self.__systemHeader = self.fixed_metadata_dict['systemHeader']
309 309 self.dtype = pickle.loads(self.fixed_metadata_dict['dtype'])
310 310 except:
311 311 pass
312 312
313 313 self.__frequency = None
314 314
315 315 self.__frequency = self.fixed_metadata_dict.get('frequency', 1)
316 316
317 317 self.__timezone = self.fixed_metadata_dict.get('timezone', 18000)
318 318
319 319 try:
320 320 nSamples = self.fixed_metadata_dict['nSamples']
321 321 except:
322 322 nSamples = None
323 323
324 324 self.__firstHeigth = 0
325 325
326 326 try:
327 327 codeType = self.__radarControllerHeader['codeType']
328 328 except:
329 329 codeType = 0
330 330
331 331 try:
332 332 if codeType:
333 333 nCode = self.__radarControllerHeader['nCode']
334 334 nBaud = self.__radarControllerHeader['nBaud']
335 335 code = self.__radarControllerHeader['code']
336 336 except:
337 337 pass
338 338
339 339 if not ippKm:
340 340 try:
341 341 # seconds to km
342 342 ippKm = self.__radarControllerHeader['ipp']
343 343 except:
344 344 ippKm = None
345 345 ####################################################
346 346 self.__ippKm = ippKm
347 347 startUTCSecond = None
348 348 endUTCSecond = None
349 349
350 350 if startDate:
351 351 startDatetime = datetime.datetime.combine(startDate, startTime)
352 352 startUTCSecond = (
353 353 startDatetime - datetime.datetime(1970, 1, 1)).total_seconds() + self.__timezone
354 354
355 355 if endDate:
356 356 endDatetime = datetime.datetime.combine(endDate, endTime)
357 357 endUTCSecond = (endDatetime - datetime.datetime(1970,
358 358 1, 1)).total_seconds() + self.__timezone
359 359
360 360
361 print(startUTCSecond,endUTCSecond)
361 #print(startUTCSecond,endUTCSecond)
362 362 start_index, end_index = self.digitalReadObj.get_bounds(
363 363 channelNameList[channelList[0]])
364 364
365 print("*****",start_index,end_index)
366 print(metadata_bounds)
365 #print("*****",start_index,end_index)
367 366 if not startUTCSecond:
368 367 startUTCSecond = start_index / self.__sample_rate
369 368
370 369 if start_index > startUTCSecond * self.__sample_rate:
371 370 startUTCSecond = start_index / self.__sample_rate
372 371
373 372 if not endUTCSecond:
374 373 endUTCSecond = end_index / self.__sample_rate
375 print("1",endUTCSecond)
376 print(self.__sample_rate)
377 374 if end_index < endUTCSecond * self.__sample_rate:
378 375 endUTCSecond = end_index / self.__sample_rate #Check UTC and LT time
379 print("2",endUTCSecond)
380 376 if not nSamples:
381 377 if not ippKm:
382 378 raise ValueError("[Reading] nSamples or ippKm should be defined")
383 379 nSamples = int(ippKm / (1e6 * 0.15 / self.__sample_rate))
384 380
385 381 channelBoundList = []
386 382 channelNameListFiltered = []
387 383
388 384 for thisIndexChannel in channelList:
389 385 thisChannelName = channelNameList[thisIndexChannel]
390 386 start_index, end_index = self.digitalReadObj.get_bounds(
391 387 thisChannelName)
392 388 channelBoundList.append((start_index, end_index))
393 389 channelNameListFiltered.append(thisChannelName)
394 390
395 391 self.profileIndex = 0
396 392 self.i = 0
397 393 self.__delay = delay
398 394
399 395 self.__codeType = codeType
400 396 self.__nCode = nCode
401 397 self.__nBaud = nBaud
402 398 self.__code = code
403 399
404 400 self.__datapath = path
405 401 self.__online = online
406 402 self.__channelList = channelList
407 403 self.__channelNameList = channelNameListFiltered
408 404 self.__channelBoundList = channelBoundList
409 405 self.__nSamples = nSamples
410 406 if self.getByBlock:
411 407 nSamples = nSamples*nProfileBlocks
412 408
413 409
414 410 self.__samples_to_read = int(nSamples) # FIJO: AHORA 40
415 411 self.__nChannels = len(self.__channelList)
416 412 #print("------------------------------------------")
417 413 #print("self.__samples_to_read",self.__samples_to_read)
418 414 #print("self.__nSamples",self.__nSamples)
419 415 # son iguales y el buffer_index da 0
420 416 self.__startUTCSecond = startUTCSecond
421 417 self.__endUTCSecond = endUTCSecond
422 418
423 419 self.__timeInterval = 1.0 * self.__samples_to_read / \
424 420 self.__sample_rate # Time interval
425 421
426 422 if online:
427 423 # self.__thisUnixSample = int(endUTCSecond*self.__sample_rate - 4*self.__samples_to_read)
428 424 startUTCSecond = numpy.floor(endUTCSecond)
429 425
430 426 # por que en el otro metodo lo primero q se hace es sumar samplestoread
431 427 self.__thisUnixSample = int(startUTCSecond * self.__sample_rate) - self.__samples_to_read
432 428
433 429 #self.__data_buffer = numpy.zeros(
434 430 # (self.__num_subchannels, self.__samples_to_read), dtype=numpy.complex)
435 431 self.__data_buffer = numpy.zeros((int(len(channelList)), self.__samples_to_read), dtype=numpy.complex)
436 432
437 433
438 434 self.__setFileHeader()
439 435 self.isConfig = True
440 436
441 437 print("[Reading] Digital RF Data was found from %s to %s " % (
442 438 datetime.datetime.utcfromtimestamp(
443 439 self.__startUTCSecond - self.__timezone),
444 440 datetime.datetime.utcfromtimestamp(
445 441 self.__endUTCSecond - self.__timezone)
446 442 ))
447 443
448 444 print("[Reading] Starting process from %s to %s" % (datetime.datetime.utcfromtimestamp(startUTCSecond - self.__timezone),
449 445 datetime.datetime.utcfromtimestamp(
450 446 endUTCSecond - self.__timezone)
451 447 ))
452 448 self.oldAverage = None
453 449 self.count = 0
454 450 self.executionTime = 0
455 451
456 452 def __reload(self):
457 453 # print
458 454 # print "%s not in range [%s, %s]" %(
459 455 # datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
460 456 # datetime.datetime.utcfromtimestamp(self.__startUTCSecond - self.__timezone),
461 457 # datetime.datetime.utcfromtimestamp(self.__endUTCSecond - self.__timezone)
462 458 # )
463 459 print("[Reading] reloading metadata ...")
464 460
465 461 try:
466 462 self.digitalReadObj.reload(complete_update=True)
467 463 except:
468 464 self.digitalReadObj = digital_rf.DigitalRFReader(self.path)
469 465
470 466 start_index, end_index = self.digitalReadObj.get_bounds(
471 467 self.__channelNameList[self.__channelList[0]])
472 468
473 469 if start_index > self.__startUTCSecond * self.__sample_rate:
474 470 self.__startUTCSecond = 1.0 * start_index / self.__sample_rate
475 471
476 472 if end_index > self.__endUTCSecond * self.__sample_rate:
477 473 self.__endUTCSecond = 1.0 * end_index / self.__sample_rate
478 474 print()
479 475 print("[Reading] New timerange found [%s, %s] " % (
480 476 datetime.datetime.utcfromtimestamp(
481 477 self.__startUTCSecond - self.__timezone),
482 478 datetime.datetime.utcfromtimestamp(
483 479 self.__endUTCSecond - self.__timezone)
484 480 ))
485 481
486 482 return True
487 483
488 484 return False
489 485
490 486 def timeit(self, toExecute):
491 487 t0 = time.time()
492 488 toExecute()
493 489 self.executionTime = time.time() - t0
494 490 if self.oldAverage is None:
495 491 self.oldAverage = self.executionTime
496 492 self.oldAverage = (self.executionTime + self.count *
497 493 self.oldAverage) / (self.count + 1.0)
498 494 self.count = self.count + 1.0
499 495 return
500 496
501 497 def __readNextBlock(self, seconds=30, volt_scale=1):
502 498 '''
503 499 '''
504 500
505 501 # Set the next data
506 502 self.__flagDiscontinuousBlock = False
507 503 self.__thisUnixSample += self.__samples_to_read
508 504
509 505 if self.__thisUnixSample + 2 * self.__samples_to_read > self.__endUTCSecond * self.__sample_rate:
510 506 print ("[Reading] There are no more data into selected time-range")
511 507 if self.__online:
512 508 sleep(3)
513 509 self.__reload()
514 510 else:
515 511 return False
516 512
517 513 if self.__thisUnixSample + 2 * self.__samples_to_read > self.__endUTCSecond * self.__sample_rate:
518 514 return False
519 515 self.__thisUnixSample -= self.__samples_to_read
520 516
521 517 indexChannel = 0
522 518
523 519 dataOk = False
524 520
525 521 for thisChannelName in self.__channelNameList: # TODO VARIOS CHANNELS?
526 522 for indexSubchannel in range(self.__num_subchannels):
527 523 try:
528 524 t0 = time()
529 525 result = self.digitalReadObj.read_vector_c81d(self.__thisUnixSample,
530 526 self.__samples_to_read,
531 527 thisChannelName, sub_channel=indexSubchannel)
532 528 self.executionTime = time() - t0
533 529 if self.oldAverage is None:
534 530 self.oldAverage = self.executionTime
535 531 self.oldAverage = (
536 532 self.executionTime + self.count * self.oldAverage) / (self.count + 1.0)
537 533 self.count = self.count + 1.0
538 534
539 535 except IOError as e:
540 536 # read next profile
541 537 self.__flagDiscontinuousBlock = True
542 538 print("[Reading] %s" % datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone), e)
543 539 break
544 540
545 541 if result.shape[0] != self.__samples_to_read:
546 542 self.__flagDiscontinuousBlock = True
547 543 print("[Reading] %s: Too few samples were found, just %d/%d samples" % (datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
548 544 result.shape[0],
549 545 self.__samples_to_read))
550 546 break
551 547
552 548 self.__data_buffer[indexChannel, :] = result * volt_scale
553 549 indexChannel+=1
554 550
555 551 dataOk = True
556 552
557 553 self.__utctime = self.__thisUnixSample / self.__sample_rate
558 554
559 555 if not dataOk:
560 556 return False
561 557
562 558 print("[Reading] %s: %d samples <> %f sec" % (datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
563 559 self.__samples_to_read,
564 560 self.__timeInterval))
565 561
566 562 self.__bufferIndex = 0
567 563
568 564 return True
569 565
570 566 def __isBufferEmpty(self):
571 567
572 568 return self.__bufferIndex > self.__samples_to_read - self.__nSamples # 40960 - 40
573 569
574 570 def getData(self, seconds=30, nTries=5):
575 571 '''
576 572 This method gets the data from files and put the data into the dataOut object
577 573
578 574 In addition, increase el the buffer counter in one.
579 575
580 576 Return:
581 577 data : retorna un perfil de voltages (alturas * canales) copiados desde el
582 578 buffer. Si no hay mas archivos a leer retorna None.
583 579
584 580 Affected:
585 581 self.dataOut
586 582 self.profileIndex
587 583 self.flagDiscontinuousBlock
588 584 self.flagIsNewBlock
589 585 '''
590 586 #print("getdata")
591 587 err_counter = 0
592 588 self.dataOut.flagNoData = True
593 589
594 590
595 591 if self.__isBufferEmpty():
596 592 #print("hi")
597 593 self.__flagDiscontinuousBlock = False
598 594
599 595 while True:
600 596 if self.__readNextBlock():
601 597 break
602 598 if self.__thisUnixSample > self.__endUTCSecond * self.__sample_rate:
603 599 raise schainpy.admin.SchainError('Error')
604 600 return
605 601
606 602 if self.__flagDiscontinuousBlock:
607 603 raise schainpy.admin.SchainError('discontinuous block found')
608 604 return
609 605
610 606 if not self.__online:
611 607 raise schainpy.admin.SchainError('Online?')
612 608 return
613 609
614 610 err_counter += 1
615 611 if err_counter > nTries:
616 612 raise schainpy.admin.SchainError('Max retrys reach')
617 613 return
618 614
619 615 print('[Reading] waiting %d seconds to read a new block' % seconds)
620 616 sleep(seconds)
621 617
622 618
623 619 if not self.getByBlock:
624 620
625 621 #print("self.__bufferIndex",self.__bufferIndex)# este valor siempre es cero aparentemente
626 622 self.dataOut.data = self.__data_buffer[:, self.__bufferIndex:self.__bufferIndex + self.__nSamples]
627 623 self.dataOut.utctime = ( self.__thisUnixSample + self.__bufferIndex) / self.__sample_rate
628 624 self.dataOut.flagNoData = False
629 625 self.dataOut.flagDiscontinuousBlock = self.__flagDiscontinuousBlock
630 626 self.dataOut.profileIndex = self.profileIndex
631 627
632 628 self.__bufferIndex += self.__nSamples
633 629 self.profileIndex += 1
634 630
635 631 if self.profileIndex == self.dataOut.nProfiles:
636 632 self.profileIndex = 0
637 633 else:
638 634 # ojo debo anadir el readNextBLock y el __isBufferEmpty(
639 635 self.dataOut.flagNoData = False
640 636 buffer = self.__data_buffer[:,self.__bufferIndex:self.__bufferIndex + self.__samples_to_read]
641 637 buffer = buffer.reshape((self.__nChannels, self.nProfileBlocks, int(self.__samples_to_read/self.nProfileBlocks)))
642 638 self.dataOut.nProfileBlocks = self.nProfileBlocks
643 639 self.dataOut.data = buffer
644 640 self.dataOut.utctime = ( self.__thisUnixSample + self.__bufferIndex) / self.__sample_rate
645 641 self.profileIndex += self.__samples_to_read
646 642 self.__bufferIndex += self.__samples_to_read
647 643 self.dataOut.flagDiscontinuousBlock = self.__flagDiscontinuousBlock
648 644 return True
649 645
650 646
651 647 def printInfo(self):
652 648 '''
653 649 '''
654 650 if self.__printInfo == False:
655 651 return
656 652
657 653 # self.systemHeaderObj.printInfo()
658 654 # self.radarControllerHeaderObj.printInfo()
659 655
660 656 self.__printInfo = False
661 657
662 658 def printNumberOfBlock(self):
663 659 '''
664 660 '''
665 661 return
666 662 # print self.profileIndex
667 663
668 664 def run(self, **kwargs):
669 665 '''
670 666 This method will be called many times so here you should put all your code
671 667 '''
672 668
673 669 if not self.isConfig:
674 670 self.setup(**kwargs)
675 671
676 672 self.getData(seconds=self.__delay)
677 673
678 674 return
679 675
680 676 @MPDecorator
681 677 class DigitalRFWriter(Operation):
682 678 '''
683 679 classdocs
684 680 '''
685 681
686 682 def __init__(self, **kwargs):
687 683 '''
688 684 Constructor
689 685 '''
690 686 Operation.__init__(self, **kwargs)
691 687 self.metadata_dict = {}
692 688 self.dataOut = None
693 689 self.dtype = None
694 690 self.oldAverage = 0
695 691
696 692 def setHeader(self):
697 693
698 694 self.metadata_dict['frequency'] = self.dataOut.frequency
699 695 self.metadata_dict['timezone'] = self.dataOut.timeZone
700 696 self.metadata_dict['dtype'] = pickle.dumps(self.dataOut.dtype)
701 697 self.metadata_dict['nProfiles'] = self.dataOut.nProfiles
702 698 self.metadata_dict['heightList'] = self.dataOut.heightList
703 699 self.metadata_dict['channelList'] = self.dataOut.channelList
704 700 self.metadata_dict['flagDecodeData'] = self.dataOut.flagDecodeData
705 701 self.metadata_dict['flagDeflipData'] = self.dataOut.flagDeflipData
706 702 self.metadata_dict['flagShiftFFT'] = self.dataOut.flagShiftFFT
707 703 self.metadata_dict['useLocalTime'] = self.dataOut.useLocalTime
708 704 self.metadata_dict['nCohInt'] = self.dataOut.nCohInt
709 705 self.metadata_dict['type'] = self.dataOut.type
710 706 self.metadata_dict['flagDataAsBlock']= getattr(
711 707 self.dataOut, 'flagDataAsBlock', None) # chequear
712 708
713 709 def setup(self, dataOut, path, frequency, fileCadence, dirCadence, metadataCadence, set=0, metadataFile='metadata', ext='.h5'):
714 710 '''
715 711 In this method we should set all initial parameters.
716 712 Input:
717 713 dataOut: Input data will also be outputa data
718 714 '''
719 715 self.setHeader()
720 716 self.__ippSeconds = dataOut.ippSeconds
721 717 self.__deltaH = dataOut.getDeltaH()
722 718 self.__sample_rate = 1e6 * 0.15 / self.__deltaH
723 719 self.__dtype = dataOut.dtype
724 720 if len(dataOut.dtype) == 2:
725 721 self.__dtype = dataOut.dtype[0]
726 722 self.__nSamples = dataOut.systemHeaderObj.nSamples
727 723 self.__nProfiles = dataOut.nProfiles
728 724
729 725 if self.dataOut.type != 'Voltage':
730 726 raise 'Digital RF cannot be used with this data type'
731 727 self.arr_data = numpy.ones((1, dataOut.nFFTPoints * len(
732 728 self.dataOut.channelList)), dtype=[('r', self.__dtype), ('i', self.__dtype)])
733 729 else:
734 730 self.arr_data = numpy.ones((self.__nSamples, len(
735 731 self.dataOut.channelList)), dtype=[('r', self.__dtype), ('i', self.__dtype)])
736 732
737 733 file_cadence_millisecs = 1000
738 734
739 735 sample_rate_fraction = Fraction(self.__sample_rate).limit_denominator()
740 736 sample_rate_numerator = int(sample_rate_fraction.numerator)
741 737 sample_rate_denominator = int(sample_rate_fraction.denominator)
742 738 start_global_index = dataOut.utctime * self.__sample_rate
743 739
744 740 uuid = 'prueba'
745 741 compression_level = 0
746 742 checksum = False
747 743 is_complex = True
748 744 num_subchannels = len(dataOut.channelList)
749 745 is_continuous = True
750 746 marching_periods = False
751 747
752 748 self.digitalWriteObj = digital_rf.DigitalRFWriter(path, self.__dtype, dirCadence,
753 749 fileCadence, start_global_index,
754 750 sample_rate_numerator, sample_rate_denominator, uuid, compression_level, checksum,
755 751 is_complex, num_subchannels, is_continuous, marching_periods)
756 752 metadata_dir = os.path.join(path, 'metadata')
757 753 os.system('mkdir %s' % (metadata_dir))
758 754 self.digitalMetadataWriteObj = digital_rf.DigitalMetadataWriter(metadata_dir, dirCadence, 1, # 236, file_cadence_millisecs / 1000
759 755 sample_rate_numerator, sample_rate_denominator,
760 756 metadataFile)
761 757 self.isConfig = True
762 758 self.currentSample = 0
763 759 self.oldAverage = 0
764 760 self.count = 0
765 761 return
766 762
767 763 def writeMetadata(self):
768 764 start_idx = self.__sample_rate * self.dataOut.utctime
769 765
770 766 self.metadata_dict['processingHeader'] = self.dataOut.processingHeaderObj.getAsDict(
771 767 )
772 768 self.metadata_dict['radarControllerHeader'] = self.dataOut.radarControllerHeaderObj.getAsDict(
773 769 )
774 770 self.metadata_dict['systemHeader'] = self.dataOut.systemHeaderObj.getAsDict(
775 771 )
776 772 self.digitalMetadataWriteObj.write(start_idx, self.metadata_dict)
777 773 return
778 774
779 775 def timeit(self, toExecute):
780 776 t0 = time()
781 777 toExecute()
782 778 self.executionTime = time() - t0
783 779 if self.oldAverage is None:
784 780 self.oldAverage = self.executionTime
785 781 self.oldAverage = (self.executionTime + self.count *
786 782 self.oldAverage) / (self.count + 1.0)
787 783 self.count = self.count + 1.0
788 784 return
789 785
790 786 def writeData(self):
791 787 if self.dataOut.type != 'Voltage':
792 788 raise 'Digital RF cannot be used with this data type'
793 789 for channel in self.dataOut.channelList:
794 790 for i in range(self.dataOut.nFFTPoints):
795 791 self.arr_data[1][channel * self.dataOut.nFFTPoints +
796 792 i]['r'] = self.dataOut.data[channel][i].real
797 793 self.arr_data[1][channel * self.dataOut.nFFTPoints +
798 794 i]['i'] = self.dataOut.data[channel][i].imag
799 795 else:
800 796 for i in range(self.dataOut.systemHeaderObj.nSamples):
801 797 for channel in self.dataOut.channelList:
802 798 self.arr_data[i][channel]['r'] = self.dataOut.data[channel][i].real
803 799 self.arr_data[i][channel]['i'] = self.dataOut.data[channel][i].imag
804 800
805 801 def f(): return self.digitalWriteObj.rf_write(self.arr_data)
806 802 self.timeit(f)
807 803
808 804 return
809 805
810 806 def run(self, dataOut, frequency=49.92e6, path=None, fileCadence=1000, dirCadence=36000, metadataCadence=1, **kwargs):
811 807 '''
812 808 This method will be called many times so here you should put all your code
813 809 Inputs:
814 810 dataOut: object with the data
815 811 '''
816 812 # print dataOut.__dict__
817 813 self.dataOut = dataOut
818 814 if not self.isConfig:
819 815 self.setup(dataOut, path, frequency, fileCadence,
820 816 dirCadence, metadataCadence, **kwargs)
821 817 self.writeMetadata()
822 818
823 819 self.writeData()
824 820
825 821 ## self.currentSample += 1
826 822 # if self.dataOut.flagDataAsBlock or self.currentSample == 1:
827 823 # self.writeMetadata()
828 824 ## if self.currentSample == self.__nProfiles: self.currentSample = 0
829 825
830 826 return dataOut# en la version 2.7 no aparece este return
831 827
832 828 def close(self):
833 829 print('[Writing] - Closing files ')
834 830 print('Average of writing to digital rf format is ', self.oldAverage * 1000)
835 831 try:
836 832 self.digitalWriteObj.close()
837 833 except:
838 834 pass
Show file before | Show file after | Hide comments
1 NO CONTENT: modified file
The requested commit or file is too big and content was truncated. Show full diff
General Comments 0
You need to be logged in to leave comments. Login now