##// END OF EJS Templates
schain
RSS
Fix python 2 compatibility
Juan C. Espinoza -
r1333:ac63616fef28 v3.0.0b2
parent child
Show More
Side by Side Unified
Context file:
r1333:ac63616fef28 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,8 +1,8
1 1 """Signal chain python package"""
2 2
3 3 try:
4 4 from .controller import Project
5 5 except:
6 6 pass
7 7
8 __version__ = '3.0.0b1'
8 __version__ = '3.0.0b2'
Show file before | Show file after | Hide comments
@@ -1,1400 +1,1400
1 1 '''
2 2
3 3 $Author: murco $
4 4 $Id: JROData.py 173 2012-11-20 15:06:21Z murco $
5 5 '''
6 6
7 7 import copy
8 8 import numpy
9 9 import datetime
10 10 import json
11 11
12 12 import schainpy.admin
13 13 from schainpy.utils import log
14 14 from .jroheaderIO import SystemHeader, RadarControllerHeader
15 15 from schainpy.model.data import _noise
16 16
17 17
18 18 def getNumpyDtype(dataTypeCode):
19 19
20 20 if dataTypeCode == 0:
21 21 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
22 22 elif dataTypeCode == 1:
23 23 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
24 24 elif dataTypeCode == 2:
25 25 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
26 26 elif dataTypeCode == 3:
27 27 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
28 28 elif dataTypeCode == 4:
29 29 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
30 30 elif dataTypeCode == 5:
31 31 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
32 32 else:
33 33 raise ValueError('dataTypeCode was not defined')
34 34
35 35 return numpyDtype
36 36
37 37
38 38 def getDataTypeCode(numpyDtype):
39 39
40 40 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
41 41 datatype = 0
42 42 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
43 43 datatype = 1
44 44 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
45 45 datatype = 2
46 46 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
47 47 datatype = 3
48 48 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
49 49 datatype = 4
50 50 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
51 51 datatype = 5
52 52 else:
53 53 datatype = None
54 54
55 55 return datatype
56 56
57 57
58 58 def hildebrand_sekhon(data, navg):
59 59 """
60 60 This method is for the objective determination of the noise level in Doppler spectra. This
61 61 implementation technique is based on the fact that the standard deviation of the spectral
62 62 densities is equal to the mean spectral density for white Gaussian noise
63 63
64 64 Inputs:
65 65 Data : heights
66 66 navg : numbers of averages
67 67
68 68 Return:
69 69 mean : noise's level
70 70 """
71 71
72 72 sortdata = numpy.sort(data, axis=None)
73 73 '''
74 74 lenOfData = len(sortdata)
75 75 nums_min = lenOfData*0.2
76 76
77 77 if nums_min <= 5:
78 78
79 79 nums_min = 5
80 80
81 81 sump = 0.
82 82 sumq = 0.
83 83
84 84 j = 0
85 85 cont = 1
86 86
87 87 while((cont == 1)and(j < lenOfData)):
88 88
89 89 sump += sortdata[j]
90 90 sumq += sortdata[j]**2
91 91
92 92 if j > nums_min:
93 93 rtest = float(j)/(j-1) + 1.0/navg
94 94 if ((sumq*j) > (rtest*sump**2)):
95 95 j = j - 1
96 96 sump = sump - sortdata[j]
97 97 sumq = sumq - sortdata[j]**2
98 98 cont = 0
99 99
100 100 j += 1
101 101
102 102 lnoise = sump / j
103 103 '''
104 104 return _noise.hildebrand_sekhon(sortdata, navg)
105 105
106 106
107 107 class Beam:
108 108
109 109 def __init__(self):
110 110 self.codeList = []
111 111 self.azimuthList = []
112 112 self.zenithList = []
113 113
114 114
115 115 class GenericData(object):
116 116
117 117 flagNoData = True
118 118
119 119 def copy(self, inputObj=None):
120 120
121 121 if inputObj == None:
122 122 return copy.deepcopy(self)
123 123
124 124 for key in list(inputObj.__dict__.keys()):
125 125
126 126 attribute = inputObj.__dict__[key]
127 127
128 128 # If this attribute is a tuple or list
129 129 if type(inputObj.__dict__[key]) in (tuple, list):
130 130 self.__dict__[key] = attribute[:]
131 131 continue
132 132
133 133 # If this attribute is another object or instance
134 134 if hasattr(attribute, '__dict__'):
135 135 self.__dict__[key] = attribute.copy()
136 136 continue
137 137
138 138 self.__dict__[key] = inputObj.__dict__[key]
139 139
140 140 def deepcopy(self):
141 141
142 142 return copy.deepcopy(self)
143 143
144 144 def isEmpty(self):
145 145
146 146 return self.flagNoData
147 147
148 148 def isReady(self):
149 149
150 150 return not self.flagNoData
151 151
152 152
153 153 class JROData(GenericData):
154 154
155 155 # m_BasicHeader = BasicHeader()
156 156 # m_ProcessingHeader = ProcessingHeader()
157 157
158 158 systemHeaderObj = SystemHeader()
159 159 radarControllerHeaderObj = RadarControllerHeader()
160 160 # data = None
161 161 type = None
162 162 datatype = None # dtype but in string
163 163 # dtype = None
164 164 # nChannels = None
165 165 # nHeights = None
166 166 nProfiles = None
167 167 heightList = None
168 168 channelList = None
169 169 flagDiscontinuousBlock = False
170 170 useLocalTime = False
171 171 utctime = None
172 172 timeZone = None
173 173 dstFlag = None
174 174 errorCount = None
175 175 blocksize = None
176 176 # nCode = None
177 177 # nBaud = None
178 178 # code = None
179 179 flagDecodeData = False # asumo q la data no esta decodificada
180 180 flagDeflipData = False # asumo q la data no esta sin flip
181 181 flagShiftFFT = False
182 182 # ippSeconds = None
183 183 # timeInterval = None
184 184 nCohInt = None
185 185 # noise = None
186 186 windowOfFilter = 1
187 187 # Speed of ligth
188 188 C = 3e8
189 189 frequency = 49.92e6
190 190 realtime = False
191 191 beacon_heiIndexList = None
192 192 last_block = None
193 193 blocknow = None
194 194 azimuth = None
195 195 zenith = None
196 196 beam = Beam()
197 197 profileIndex = None
198 198 error = None
199 199 data = None
200 200 nmodes = None
201 201
202 202 def __str__(self):
203 203
204 204 return '{} - {}'.format(self.type, self.getDatatime())
205 205
206 206 def getNoise(self):
207 207
208 208 raise NotImplementedError
209 209
210 210 def getNChannels(self):
211 211
212 212 return len(self.channelList)
213 213
214 214 def getChannelIndexList(self):
215 215
216 216 return list(range(self.nChannels))
217 217
218 218 def getNHeights(self):
219 219
220 220 return len(self.heightList)
221 221
222 222 def getHeiRange(self, extrapoints=0):
223 223
224 224 heis = self.heightList
225 225 # deltah = self.heightList[1] - self.heightList[0]
226 226 #
227 227 # heis.append(self.heightList[-1])
228 228
229 229 return heis
230 230
231 231 def getDeltaH(self):
232 232
233 233 delta = self.heightList[1] - self.heightList[0]
234 234
235 235 return delta
236 236
237 237 def getltctime(self):
238 238
239 239 if self.useLocalTime:
240 240 return self.utctime - self.timeZone * 60
241 241
242 242 return self.utctime
243 243
244 244 def getDatatime(self):
245 245
246 246 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
247 247 return datatimeValue
248 248
249 249 def getTimeRange(self):
250 250
251 251 datatime = []
252 252
253 253 datatime.append(self.ltctime)
254 254 datatime.append(self.ltctime + self.timeInterval + 1)
255 255
256 256 datatime = numpy.array(datatime)
257 257
258 258 return datatime
259 259
260 260 def getFmaxTimeResponse(self):
261 261
262 262 period = (10**-6) * self.getDeltaH() / (0.15)
263 263
264 264 PRF = 1. / (period * self.nCohInt)
265 265
266 266 fmax = PRF
267 267
268 268 return fmax
269 269
270 270 def getFmax(self):
271 271 PRF = 1. / (self.ippSeconds * self.nCohInt)
272 272
273 273 fmax = PRF
274 274 return fmax
275 275
276 276 def getVmax(self):
277 277
278 278 _lambda = self.C / self.frequency
279 279
280 280 vmax = self.getFmax() * _lambda / 2
281 281
282 282 return vmax
283 283
284 284 def get_ippSeconds(self):
285 285 '''
286 286 '''
287 287 return self.radarControllerHeaderObj.ippSeconds
288 288
289 289 def set_ippSeconds(self, ippSeconds):
290 290 '''
291 291 '''
292 292
293 293 self.radarControllerHeaderObj.ippSeconds = ippSeconds
294 294
295 295 return
296 296
297 297 def get_dtype(self):
298 298 '''
299 299 '''
300 300 return getNumpyDtype(self.datatype)
301 301
302 302 def set_dtype(self, numpyDtype):
303 303 '''
304 304 '''
305 305
306 306 self.datatype = getDataTypeCode(numpyDtype)
307 307
308 308 def get_code(self):
309 309 '''
310 310 '''
311 311 return self.radarControllerHeaderObj.code
312 312
313 313 def set_code(self, code):
314 314 '''
315 315 '''
316 316 self.radarControllerHeaderObj.code = code
317 317
318 318 return
319 319
320 320 def get_ncode(self):
321 321 '''
322 322 '''
323 323 return self.radarControllerHeaderObj.nCode
324 324
325 325 def set_ncode(self, nCode):
326 326 '''
327 327 '''
328 328 self.radarControllerHeaderObj.nCode = nCode
329 329
330 330 return
331 331
332 332 def get_nbaud(self):
333 333 '''
334 334 '''
335 335 return self.radarControllerHeaderObj.nBaud
336 336
337 337 def set_nbaud(self, nBaud):
338 338 '''
339 339 '''
340 340 self.radarControllerHeaderObj.nBaud = nBaud
341 341
342 342 return
343 343
344 344 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
345 345 channelIndexList = property(
346 346 getChannelIndexList, "I'm the 'channelIndexList' property.")
347 347 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
348 348 #noise = property(getNoise, "I'm the 'nHeights' property.")
349 349 datatime = property(getDatatime, "I'm the 'datatime' property")
350 350 ltctime = property(getltctime, "I'm the 'ltctime' property")
351 351 ippSeconds = property(get_ippSeconds, set_ippSeconds)
352 352 dtype = property(get_dtype, set_dtype)
353 353 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
354 354 code = property(get_code, set_code)
355 355 nCode = property(get_ncode, set_ncode)
356 356 nBaud = property(get_nbaud, set_nbaud)
357 357
358 358
359 359 class Voltage(JROData):
360 360
361 361 # data es un numpy array de 2 dmensiones (canales, alturas)
362 362 data = None
363 363 dataPP_POW = None
364 364 dataPP_DOP = None
365 365 dataPP_WIDTH = None
366 366 dataPP_SNR = None
367 367
368 368 def __init__(self):
369 369 '''
370 370 Constructor
371 371 '''
372 372
373 373 self.useLocalTime = True
374 374 self.radarControllerHeaderObj = RadarControllerHeader()
375 375 self.systemHeaderObj = SystemHeader()
376 376 self.type = "Voltage"
377 377 self.data = None
378 378 # self.dtype = None
379 379 # self.nChannels = 0
380 380 # self.nHeights = 0
381 381 self.nProfiles = None
382 382 self.heightList = None
383 383 self.channelList = None
384 384 # self.channelIndexList = None
385 385 self.flagNoData = True
386 386 self.flagDiscontinuousBlock = False
387 387 self.utctime = None
388 388 self.timeZone = 0
389 389 self.dstFlag = None
390 390 self.errorCount = None
391 391 self.nCohInt = None
392 392 self.blocksize = None
393 393 self.flagCohInt = False
394 394 self.flagDecodeData = False # asumo q la data no esta decodificada
395 395 self.flagDeflipData = False # asumo q la data no esta sin flip
396 396 self.flagShiftFFT = False
397 397 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
398 398 self.profileIndex = 0
399 399
400 400 def getNoisebyHildebrand(self, channel=None):
401 401 """
402 402 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
403 403
404 404 Return:
405 405 noiselevel
406 406 """
407 407
408 408 if channel != None:
409 409 data = self.data[channel]
410 410 nChannels = 1
411 411 else:
412 412 data = self.data
413 413 nChannels = self.nChannels
414 414
415 415 noise = numpy.zeros(nChannels)
416 416 power = data * numpy.conjugate(data)
417 417
418 418 for thisChannel in range(nChannels):
419 419 if nChannels == 1:
420 420 daux = power[:].real
421 421 else:
422 422 daux = power[thisChannel, :].real
423 423 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
424 424
425 425 return noise
426 426
427 427 def getNoise(self, type=1, channel=None):
428 428
429 429 if type == 1:
430 430 noise = self.getNoisebyHildebrand(channel)
431 431
432 432 return noise
433 433
434 434 def getPower(self, channel=None):
435 435
436 436 if channel != None:
437 437 data = self.data[channel]
438 438 else:
439 439 data = self.data
440 440
441 441 power = data * numpy.conjugate(data)
442 442 powerdB = 10 * numpy.log10(power.real)
443 443 powerdB = numpy.squeeze(powerdB)
444 444
445 445 return powerdB
446 446
447 447 def getTimeInterval(self):
448 448
449 449 timeInterval = self.ippSeconds * self.nCohInt
450 450
451 451 return timeInterval
452 452
453 453 noise = property(getNoise, "I'm the 'nHeights' property.")
454 454 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
455 455
456 456
457 457 class Spectra(JROData):
458 458
459 459 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
460 460 data_spc = None
461 461 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
462 462 data_cspc = None
463 463 # data dc es un numpy array de 2 dmensiones (canales, alturas)
464 464 data_dc = None
465 465 # data power
466 466 data_pwr = None
467 467 nFFTPoints = None
468 468 # nPairs = None
469 469 pairsList = None
470 470 nIncohInt = None
471 471 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
472 472 nCohInt = None # se requiere para determinar el valor de timeInterval
473 473 ippFactor = None
474 474 profileIndex = 0
475 475 plotting = "spectra"
476 476
477 477 def __init__(self):
478 478 '''
479 479 Constructor
480 480 '''
481 481
482 482 self.useLocalTime = True
483 483 self.radarControllerHeaderObj = RadarControllerHeader()
484 484 self.systemHeaderObj = SystemHeader()
485 485 self.type = "Spectra"
486 486 self.timeZone = 0
487 487 # self.data = None
488 488 # self.dtype = None
489 489 # self.nChannels = 0
490 490 # self.nHeights = 0
491 491 self.nProfiles = None
492 492 self.heightList = None
493 493 self.channelList = None
494 494 # self.channelIndexList = None
495 495 self.pairsList = None
496 496 self.flagNoData = True
497 497 self.flagDiscontinuousBlock = False
498 498 self.utctime = None
499 499 self.nCohInt = None
500 500 self.nIncohInt = None
501 501 self.blocksize = None
502 502 self.nFFTPoints = None
503 503 self.wavelength = None
504 504 self.flagDecodeData = False # asumo q la data no esta decodificada
505 505 self.flagDeflipData = False # asumo q la data no esta sin flip
506 506 self.flagShiftFFT = False
507 507 self.ippFactor = 1
508 508 #self.noise = None
509 509 self.beacon_heiIndexList = []
510 510 self.noise_estimation = None
511 511
512 512 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
513 513 """
514 514 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
515 515
516 516 Return:
517 517 noiselevel
518 518 """
519 519
520 520 noise = numpy.zeros(self.nChannels)
521 521
522 522 for channel in range(self.nChannels):
523 523 daux = self.data_spc[channel,
524 524 xmin_index:xmax_index, ymin_index:ymax_index]
525 525 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
526 526
527 527 return noise
528 528
529 529 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
530 530
531 531 if self.noise_estimation is not None:
532 532 # this was estimated by getNoise Operation defined in jroproc_spectra.py
533 533 return self.noise_estimation
534 534 else:
535 535 noise = self.getNoisebyHildebrand(
536 536 xmin_index, xmax_index, ymin_index, ymax_index)
537 537 return noise
538 538
539 539 def getFreqRangeTimeResponse(self, extrapoints=0):
540 540
541 541 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
542 542 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.) - deltafreq / 2
543 543
544 544 return freqrange
545 545
546 546 def getAcfRange(self, extrapoints=0):
547 547
548 548 deltafreq = 10. / (self.getFmax() / (self.nFFTPoints * self.ippFactor))
549 549 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
550 550
551 551 return freqrange
552 552
553 553 def getFreqRange(self, extrapoints=0):
554 554
555 555 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
556 556 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
557 557
558 558 return freqrange
559 559
560 560 def getVelRange(self, extrapoints=0):
561 561
562 562 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
563 563 velrange = deltav * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.)
564 564
565 565 if self.nmodes:
566 566 return velrange/self.nmodes
567 567 else:
568 568 return velrange
569 569
570 570 def getNPairs(self):
571 571
572 572 return len(self.pairsList)
573 573
574 574 def getPairsIndexList(self):
575 575
576 576 return list(range(self.nPairs))
577 577
578 578 def getNormFactor(self):
579 579
580 580 pwcode = 1
581 581
582 582 if self.flagDecodeData:
583 583 pwcode = numpy.sum(self.code[0]**2)
584 584 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
585 585 normFactor = self.nProfiles * self.nIncohInt * self.nCohInt * pwcode * self.windowOfFilter
586 586
587 587 return normFactor
588 588
589 589 def getFlagCspc(self):
590 590
591 591 if self.data_cspc is None:
592 592 return True
593 593
594 594 return False
595 595
596 596 def getFlagDc(self):
597 597
598 598 if self.data_dc is None:
599 599 return True
600 600
601 601 return False
602 602
603 603 def getTimeInterval(self):
604 604
605 605 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
606 606 if self.nmodes:
607 607 return self.nmodes*timeInterval
608 608 else:
609 609 return timeInterval
610 610
611 611 def getPower(self):
612 612
613 613 factor = self.normFactor
614 614 z = self.data_spc / factor
615 615 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
616 616 avg = numpy.average(z, axis=1)
617 617
618 618 return 10 * numpy.log10(avg)
619 619
620 620 def getCoherence(self, pairsList=None, phase=False):
621 621
622 622 z = []
623 623 if pairsList is None:
624 624 pairsIndexList = self.pairsIndexList
625 625 else:
626 626 pairsIndexList = []
627 627 for pair in pairsList:
628 628 if pair not in self.pairsList:
629 629 raise ValueError("Pair %s is not in dataOut.pairsList" % (
630 630 pair))
631 631 pairsIndexList.append(self.pairsList.index(pair))
632 632 for i in range(len(pairsIndexList)):
633 633 pair = self.pairsList[pairsIndexList[i]]
634 634 ccf = numpy.average(self.data_cspc[pairsIndexList[i], :, :], axis=0)
635 635 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
636 636 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
637 637 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
638 638 if phase:
639 639 data = numpy.arctan2(avgcoherenceComplex.imag,
640 640 avgcoherenceComplex.real) * 180 / numpy.pi
641 641 else:
642 642 data = numpy.abs(avgcoherenceComplex)
643 643
644 644 z.append(data)
645 645
646 646 return numpy.array(z)
647 647
648 648 def setValue(self, value):
649 649
650 650 print("This property should not be initialized")
651 651
652 652 return
653 653
654 654 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
655 655 pairsIndexList = property(
656 656 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
657 657 normFactor = property(getNormFactor, setValue,
658 658 "I'm the 'getNormFactor' property.")
659 659 flag_cspc = property(getFlagCspc, setValue)
660 660 flag_dc = property(getFlagDc, setValue)
661 661 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
662 662 timeInterval = property(getTimeInterval, setValue,
663 663 "I'm the 'timeInterval' property")
664 664
665 665
666 666 class SpectraHeis(Spectra):
667 667
668 668 data_spc = None
669 669 data_cspc = None
670 670 data_dc = None
671 671 nFFTPoints = None
672 672 # nPairs = None
673 673 pairsList = None
674 674 nCohInt = None
675 675 nIncohInt = None
676 676
677 677 def __init__(self):
678 678
679 679 self.radarControllerHeaderObj = RadarControllerHeader()
680 680
681 681 self.systemHeaderObj = SystemHeader()
682 682
683 683 self.type = "SpectraHeis"
684 684
685 685 # self.dtype = None
686 686
687 687 # self.nChannels = 0
688 688
689 689 # self.nHeights = 0
690 690
691 691 self.nProfiles = None
692 692
693 693 self.heightList = None
694 694
695 695 self.channelList = None
696 696
697 697 # self.channelIndexList = None
698 698
699 699 self.flagNoData = True
700 700
701 701 self.flagDiscontinuousBlock = False
702 702
703 703 # self.nPairs = 0
704 704
705 705 self.utctime = None
706 706
707 707 self.blocksize = None
708 708
709 709 self.profileIndex = 0
710 710
711 711 self.nCohInt = 1
712 712
713 713 self.nIncohInt = 1
714 714
715 715 def getNormFactor(self):
716 716 pwcode = 1
717 717 if self.flagDecodeData:
718 718 pwcode = numpy.sum(self.code[0]**2)
719 719
720 720 normFactor = self.nIncohInt * self.nCohInt * pwcode
721 721
722 722 return normFactor
723 723
724 724 def getTimeInterval(self):
725 725
726 726 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
727 727
728 728 return timeInterval
729 729
730 730 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
731 731 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
732 732
733 733
734 734 class Fits(JROData):
735 735
736 736 heightList = None
737 737 channelList = None
738 738 flagNoData = True
739 739 flagDiscontinuousBlock = False
740 740 useLocalTime = False
741 741 utctime = None
742 742 # ippSeconds = None
743 743 # timeInterval = None
744 744 nCohInt = None
745 745 nIncohInt = None
746 746 noise = None
747 747 windowOfFilter = 1
748 748 # Speed of ligth
749 749 C = 3e8
750 750 frequency = 49.92e6
751 751 realtime = False
752 752
753 753 def __init__(self):
754 754
755 755 self.type = "Fits"
756 756
757 757 self.nProfiles = None
758 758
759 759 self.heightList = None
760 760
761 761 self.channelList = None
762 762
763 763 # self.channelIndexList = None
764 764
765 765 self.flagNoData = True
766 766
767 767 self.utctime = None
768 768
769 769 self.nCohInt = 1
770 770
771 771 self.nIncohInt = 1
772 772
773 773 self.useLocalTime = True
774 774
775 775 self.profileIndex = 0
776 776
777 777 # self.utctime = None
778 778 self.timeZone = 0
779 779 # self.ltctime = None
780 780 # self.timeInterval = None
781 781 # self.header = None
782 782 # self.data_header = None
783 783 # self.data = None
784 784 # self.datatime = None
785 785 # self.flagNoData = False
786 786 # self.expName = ''
787 787 # self.nChannels = None
788 788 # self.nSamples = None
789 789 # self.dataBlocksPerFile = None
790 790 # self.comments = ''
791 791 #
792 792
793 793 def getltctime(self):
794 794
795 795 if self.useLocalTime:
796 796 return self.utctime - self.timeZone * 60
797 797
798 798 return self.utctime
799 799
800 800 def getDatatime(self):
801 801
802 802 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
803 803 return datatime
804 804
805 805 def getTimeRange(self):
806 806
807 807 datatime = []
808 808
809 809 datatime.append(self.ltctime)
810 810 datatime.append(self.ltctime + self.timeInterval)
811 811
812 812 datatime = numpy.array(datatime)
813 813
814 814 return datatime
815 815
816 816 def getHeiRange(self):
817 817
818 818 heis = self.heightList
819 819
820 820 return heis
821 821
822 822 def getNHeights(self):
823 823
824 824 return len(self.heightList)
825 825
826 826 def getNChannels(self):
827 827
828 828 return len(self.channelList)
829 829
830 830 def getChannelIndexList(self):
831 831
832 832 return list(range(self.nChannels))
833 833
834 834 def getNoise(self, type=1):
835 835
836 836 #noise = numpy.zeros(self.nChannels)
837 837
838 838 if type == 1:
839 839 noise = self.getNoisebyHildebrand()
840 840
841 841 if type == 2:
842 842 noise = self.getNoisebySort()
843 843
844 844 if type == 3:
845 845 noise = self.getNoisebyWindow()
846 846
847 847 return noise
848 848
849 849 def getTimeInterval(self):
850 850
851 851 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
852 852
853 853 return timeInterval
854 854
855 855 def get_ippSeconds(self):
856 856 '''
857 857 '''
858 858 return self.ipp_sec
859 859
860 860
861 861 datatime = property(getDatatime, "I'm the 'datatime' property")
862 862 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
863 863 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
864 864 channelIndexList = property(
865 865 getChannelIndexList, "I'm the 'channelIndexList' property.")
866 866 noise = property(getNoise, "I'm the 'nHeights' property.")
867 867
868 868 ltctime = property(getltctime, "I'm the 'ltctime' property")
869 869 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
870 870 ippSeconds = property(get_ippSeconds, '')
871 871
872 872 class Correlation(JROData):
873 873
874 874 noise = None
875 875 SNR = None
876 876 #--------------------------------------------------
877 877 mode = None
878 878 split = False
879 879 data_cf = None
880 880 lags = None
881 881 lagRange = None
882 882 pairsList = None
883 883 normFactor = None
884 884 #--------------------------------------------------
885 885 # calculateVelocity = None
886 886 nLags = None
887 887 nPairs = None
888 888 nAvg = None
889 889
890 890 def __init__(self):
891 891 '''
892 892 Constructor
893 893 '''
894 894 self.radarControllerHeaderObj = RadarControllerHeader()
895 895
896 896 self.systemHeaderObj = SystemHeader()
897 897
898 898 self.type = "Correlation"
899 899
900 900 self.data = None
901 901
902 902 self.dtype = None
903 903
904 904 self.nProfiles = None
905 905
906 906 self.heightList = None
907 907
908 908 self.channelList = None
909 909
910 910 self.flagNoData = True
911 911
912 912 self.flagDiscontinuousBlock = False
913 913
914 914 self.utctime = None
915 915
916 916 self.timeZone = 0
917 917
918 918 self.dstFlag = None
919 919
920 920 self.errorCount = None
921 921
922 922 self.blocksize = None
923 923
924 924 self.flagDecodeData = False # asumo q la data no esta decodificada
925 925
926 926 self.flagDeflipData = False # asumo q la data no esta sin flip
927 927
928 928 self.pairsList = None
929 929
930 930 self.nPoints = None
931 931
932 932 def getPairsList(self):
933 933
934 934 return self.pairsList
935 935
936 936 def getNoise(self, mode=2):
937 937
938 938 indR = numpy.where(self.lagR == 0)[0][0]
939 939 indT = numpy.where(self.lagT == 0)[0][0]
940 940
941 941 jspectra0 = self.data_corr[:, :, indR, :]
942 942 jspectra = copy.copy(jspectra0)
943 943
944 944 num_chan = jspectra.shape[0]
945 945 num_hei = jspectra.shape[2]
946 946
947 947 freq_dc = jspectra.shape[1] / 2
948 948 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
949 949
950 950 if ind_vel[0] < 0:
951 951 ind_vel[list(range(0, 1))] = ind_vel[list(
952 952 range(0, 1))] + self.num_prof
953 953
954 954 if mode == 1:
955 955 jspectra[:, freq_dc, :] = (
956 956 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
957 957
958 958 if mode == 2:
959 959
960 960 vel = numpy.array([-2, -1, 1, 2])
961 961 xx = numpy.zeros([4, 4])
962 962
963 963 for fil in range(4):
964 964 xx[fil, :] = vel[fil]**numpy.asarray(list(range(4)))
965 965
966 966 xx_inv = numpy.linalg.inv(xx)
967 967 xx_aux = xx_inv[0, :]
968 968
969 969 for ich in range(num_chan):
970 970 yy = jspectra[ich, ind_vel, :]
971 971 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
972 972
973 973 junkid = jspectra[ich, freq_dc, :] <= 0
974 974 cjunkid = sum(junkid)
975 975
976 976 if cjunkid.any():
977 977 jspectra[ich, freq_dc, junkid.nonzero()] = (
978 978 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
979 979
980 980 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
981 981
982 982 return noise
983 983
984 984 def getTimeInterval(self):
985 985
986 986 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
987 987
988 988 return timeInterval
989 989
990 990 def splitFunctions(self):
991 991
992 992 pairsList = self.pairsList
993 993 ccf_pairs = []
994 994 acf_pairs = []
995 995 ccf_ind = []
996 996 acf_ind = []
997 997 for l in range(len(pairsList)):
998 998 chan0 = pairsList[l][0]
999 999 chan1 = pairsList[l][1]
1000 1000
1001 1001 # Obteniendo pares de Autocorrelacion
1002 1002 if chan0 == chan1:
1003 1003 acf_pairs.append(chan0)
1004 1004 acf_ind.append(l)
1005 1005 else:
1006 1006 ccf_pairs.append(pairsList[l])
1007 1007 ccf_ind.append(l)
1008 1008
1009 1009 data_acf = self.data_cf[acf_ind]
1010 1010 data_ccf = self.data_cf[ccf_ind]
1011 1011
1012 1012 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1013 1013
1014 1014 def getNormFactor(self):
1015 1015 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1016 1016 acf_pairs = numpy.array(acf_pairs)
1017 1017 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1018 1018
1019 1019 for p in range(self.nPairs):
1020 1020 pair = self.pairsList[p]
1021 1021
1022 1022 ch0 = pair[0]
1023 1023 ch1 = pair[1]
1024 1024
1025 1025 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1026 1026 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1027 1027 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1028 1028
1029 1029 return normFactor
1030 1030
1031 1031 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1032 1032 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1033 1033
1034 1034
1035 1035 class Parameters(Spectra):
1036 1036
1037 1037 experimentInfo = None # Information about the experiment
1038 1038 # Information from previous data
1039 1039 inputUnit = None # Type of data to be processed
1040 1040 operation = None # Type of operation to parametrize
1041 1041 # normFactor = None #Normalization Factor
1042 1042 groupList = None # List of Pairs, Groups, etc
1043 1043 # Parameters
1044 1044 data_param = None # Parameters obtained
1045 1045 data_pre = None # Data Pre Parametrization
1046 1046 data_SNR = None # Signal to Noise Ratio
1047 1047 # heightRange = None #Heights
1048 1048 abscissaList = None # Abscissa, can be velocities, lags or time
1049 1049 # noise = None #Noise Potency
1050 1050 utctimeInit = None # Initial UTC time
1051 1051 paramInterval = None # Time interval to calculate Parameters in seconds
1052 1052 useLocalTime = True
1053 1053 # Fitting
1054 1054 data_error = None # Error of the estimation
1055 1055 constants = None
1056 1056 library = None
1057 1057 # Output signal
1058 1058 outputInterval = None # Time interval to calculate output signal in seconds
1059 1059 data_output = None # Out signal
1060 1060 nAvg = None
1061 1061 noise_estimation = None
1062 1062 GauSPC = None # Fit gaussian SPC
1063 1063
1064 1064 def __init__(self):
1065 1065 '''
1066 1066 Constructor
1067 1067 '''
1068 1068 self.radarControllerHeaderObj = RadarControllerHeader()
1069 1069 self.systemHeaderObj = SystemHeader()
1070 1070 self.type = "Parameters"
1071 1071 self.timeZone = 0
1072 1072
1073 1073 def getTimeRange1(self, interval):
1074 1074
1075 1075 datatime = []
1076 1076
1077 1077 if self.useLocalTime:
1078 1078 time1 = self.utctimeInit - self.timeZone * 60
1079 1079 else:
1080 1080 time1 = self.utctimeInit
1081 1081
1082 1082 datatime.append(time1)
1083 1083 datatime.append(time1 + interval)
1084 1084 datatime = numpy.array(datatime)
1085 1085
1086 1086 return datatime
1087 1087
1088 1088 def getTimeInterval(self):
1089 1089
1090 1090 if hasattr(self, 'timeInterval1'):
1091 1091 return self.timeInterval1
1092 1092 else:
1093 1093 return self.paramInterval
1094 1094
1095 1095 def setValue(self, value):
1096 1096
1097 1097 print("This property should not be initialized")
1098 1098
1099 1099 return
1100 1100
1101 1101 def getNoise(self):
1102 1102
1103 1103 return self.spc_noise
1104 1104
1105 1105 timeInterval = property(getTimeInterval)
1106 1106 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
1107 1107
1108 1108
1109 1109 class PlotterData(object):
1110 1110 '''
1111 1111 Object to hold data to be plotted
1112 1112 '''
1113 1113
1114 1114 MAXNUMX = 200
1115 1115 MAXNUMY = 200
1116 1116
1117 1117 def __init__(self, code, throttle_value, exp_code, localtime=True, buffering=True, snr=False):
1118 1118
1119 1119 self.key = code
1120 1120 self.throttle = throttle_value
1121 1121 self.exp_code = exp_code
1122 1122 self.buffering = buffering
1123 1123 self.ready = False
1124 1124 self.flagNoData = False
1125 1125 self.localtime = localtime
1126 1126 self.data = {}
1127 1127 self.meta = {}
1128 1128 self.__heights = []
1129 1129
1130 1130 if 'snr' in code:
1131 1131 self.plottypes = ['snr']
1132 1132 elif code == 'spc':
1133 1133 self.plottypes = ['spc', 'noise', 'rti']
1134 1134 elif code == 'cspc':
1135 1135 self.plottypes = ['cspc', 'spc', 'noise', 'rti']
1136 1136 elif code == 'rti':
1137 1137 self.plottypes = ['noise', 'rti']
1138 1138 else:
1139 1139 self.plottypes = [code]
1140 1140
1141 1141 if 'snr' not in self.plottypes and snr:
1142 1142 self.plottypes.append('snr')
1143 1143
1144 1144 for plot in self.plottypes:
1145 1145 self.data[plot] = {}
1146 1146
1147 1147 def __str__(self):
1148 1148 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1149 1149 return 'Data[{}][{}]'.format(';'.join(dum), len(self.times))
1150 1150
1151 1151 def __len__(self):
1152 1152 return len(self.data[self.key])
1153 1153
1154 1154 def __getitem__(self, key):
1155 1155
1156 1156 if key not in self.data:
1157 1157 raise KeyError(log.error('Missing key: {}'.format(key)))
1158 1158 if 'spc' in key or not self.buffering:
1159 1159 ret = self.data[key][self.tm]
1160 1160 elif 'scope' in key:
1161 1161 ret = numpy.array(self.data[key][float(self.tm)])
1162 1162 else:
1163 1163 ret = numpy.array([self.data[key][x] for x in self.times])
1164 1164 if ret.ndim > 1:
1165 1165 ret = numpy.swapaxes(ret, 0, 1)
1166 1166 return ret
1167 1167
1168 1168 def __contains__(self, key):
1169 1169 return key in self.data
1170 1170
1171 1171 def setup(self):
1172 1172 '''
1173 1173 Configure object
1174 1174 '''
1175 1175 self.type = ''
1176 1176 self.ready = False
1177 1177 del self.data
1178 1178 self.data = {}
1179 1179 self.__heights = []
1180 1180 self.__all_heights = set()
1181 1181 for plot in self.plottypes:
1182 1182 if 'snr' in plot:
1183 1183 plot = 'snr'
1184 1184 elif 'spc_moments' == plot:
1185 1185 plot = 'moments'
1186 1186 self.data[plot] = {}
1187 1187
1188 1188 if 'spc' in self.data or 'rti' in self.data or 'cspc' in self.data or 'moments' in self.data:
1189 1189 self.data['noise'] = {}
1190 1190 self.data['rti'] = {}
1191 1191 if 'noise' not in self.plottypes:
1192 1192 self.plottypes.append('noise')
1193 1193 if 'rti' not in self.plottypes:
1194 1194 self.plottypes.append('rti')
1195 1195
1196 1196 def shape(self, key):
1197 1197 '''
1198 1198 Get the shape of the one-element data for the given key
1199 1199 '''
1200 1200
1201 1201 if len(self.data[key]):
1202 1202 if 'spc' in key or not self.buffering:
1203 1203 return self.data[key].shape
1204 1204 return self.data[key][self.times[0]].shape
1205 1205 return (0,)
1206 1206
1207 1207 def update(self, dataOut, tm):
1208 1208 '''
1209 1209 Update data object with new dataOut
1210 1210 '''
1211 1211
1212 1212 self.profileIndex = dataOut.profileIndex
1213 1213 self.tm = tm
1214 1214 self.type = dataOut.type
1215 1215 self.parameters = getattr(dataOut, 'parameters', [])
1216 1216
1217 1217 if hasattr(dataOut, 'meta'):
1218 1218 self.meta.update(dataOut.meta)
1219 1219
1220 1220 if hasattr(dataOut, 'pairsList'):
1221 1221 self.pairs = dataOut.pairsList
1222 1222
1223 1223 self.interval = dataOut.getTimeInterval()
1224 1224 if True in ['spc' in ptype for ptype in self.plottypes]:
1225 1225 self.xrange = (dataOut.getFreqRange(1)/1000.,
1226 1226 dataOut.getAcfRange(1), dataOut.getVelRange(1))
1227 1227 self.__heights.append(dataOut.heightList)
1228 1228 self.__all_heights.update(dataOut.heightList)
1229 1229
1230 1230 for plot in self.plottypes:
1231 1231 if plot in ('spc', 'spc_moments', 'spc_cut'):
1232 1232 z = dataOut.data_spc/dataOut.normFactor
1233 1233 buffer = 10*numpy.log10(z)
1234 1234 if plot == 'cspc':
1235 1235 buffer = (dataOut.data_spc, dataOut.data_cspc)
1236 1236 if plot == 'noise':
1237 1237 buffer = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
1238 1238 if plot in ('rti', 'spcprofile'):
1239 1239 buffer = dataOut.getPower()
1240 1240 if plot == 'snr_db':
1241 1241 buffer = dataOut.data_SNR
1242 1242 if plot == 'snr':
1243 1243 buffer = 10*numpy.log10(dataOut.data_SNR)
1244 1244 if plot == 'dop':
1245 1245 buffer = dataOut.data_DOP
1246 1246 if plot == 'pow':
1247 1247 buffer = 10*numpy.log10(dataOut.data_POW)
1248 1248 if plot == 'width':
1249 1249 buffer = dataOut.data_WIDTH
1250 1250 if plot == 'coh':
1251 1251 buffer = dataOut.getCoherence()
1252 1252 if plot == 'phase':
1253 1253 buffer = dataOut.getCoherence(phase=True)
1254 1254 if plot == 'output':
1255 1255 buffer = dataOut.data_output
1256 1256 if plot == 'param':
1257 1257 buffer = dataOut.data_param
1258 1258 if plot == 'scope':
1259 1259 buffer = dataOut.data
1260 1260 self.flagDataAsBlock = dataOut.flagDataAsBlock
1261 1261 self.nProfiles = dataOut.nProfiles
1262 1262 if plot == 'pp_power':
1263 1263 buffer = dataOut.dataPP_POWER
1264 1264 self.flagDataAsBlock = dataOut.flagDataAsBlock
1265 1265 self.nProfiles = dataOut.nProfiles
1266 1266 if plot == 'pp_signal':
1267 1267 buffer = dataOut.dataPP_POW
1268 1268 self.flagDataAsBlock = dataOut.flagDataAsBlock
1269 1269 self.nProfiles = dataOut.nProfiles
1270 1270 if plot == 'pp_velocity':
1271 1271 buffer = dataOut.dataPP_DOP
1272 1272 self.flagDataAsBlock = dataOut.flagDataAsBlock
1273 1273 self.nProfiles = dataOut.nProfiles
1274 1274 if plot == 'pp_specwidth':
1275 1275 buffer = dataOut.dataPP_WIDTH
1276 1276 self.flagDataAsBlock = dataOut.flagDataAsBlock
1277 1277 self.nProfiles = dataOut.nProfiles
1278 1278
1279 1279 if plot == 'spc':
1280 1280 self.data['spc'][tm] = buffer
1281 1281 elif plot == 'cspc':
1282 1282 self.data['cspc'][tm] = buffer
1283 1283 elif plot == 'spc_moments':
1284 1284 self.data['spc'][tm] = buffer
1285 1285 self.data['moments'][tm] = dataOut.moments
1286 1286 else:
1287 1287 if self.buffering:
1288 1288 self.data[plot][tm] = buffer
1289 1289 else:
1290 1290 self.data[plot][tm] = buffer
1291 1291
1292 1292 if dataOut.channelList is None:
1293 1293 self.channels = range(buffer.shape[0])
1294 1294 else:
1295 1295 self.channels = dataOut.channelList
1296 1296
1297 1297 if buffer is None:
1298 1298 self.flagNoData = True
1299 1299 raise schainpy.admin.SchainWarning('Attribute data_{} is empty'.format(self.key))
1300 1300
1301 1301 def normalize_heights(self):
1302 1302 '''
1303 1303 Ensure same-dimension of the data for different heighList
1304 1304 '''
1305 1305
1306 1306 H = numpy.array(list(self.__all_heights))
1307 1307 H.sort()
1308 1308 for key in self.data:
1309 1309 shape = self.shape(key)[:-1] + H.shape
1310 1310 for tm, obj in list(self.data[key].items()):
1311 1311 h = self.__heights[self.times.tolist().index(tm)]
1312 1312 if H.size == h.size:
1313 1313 continue
1314 1314 index = numpy.where(numpy.in1d(H, h))[0]
1315 1315 dummy = numpy.zeros(shape) + numpy.nan
1316 1316 if len(shape) == 2:
1317 1317 dummy[:, index] = obj
1318 1318 else:
1319 1319 dummy[index] = obj
1320 1320 self.data[key][tm] = dummy
1321 1321
1322 1322 self.__heights = [H for tm in self.times]
1323 1323
1324 1324 def jsonify(self, tm, plot_name, plot_type, decimate=False):
1325 1325 '''
1326 1326 Convert data to json
1327 1327 '''
1328 1328
1329 1329 dy = int(self.heights.size/self.MAXNUMY) + 1
1330 1330 if self.key in ('spc', 'cspc'):
1331 1331 dx = int(self.data[self.key][tm].shape[1]/self.MAXNUMX) + 1
1332 1332 data = self.roundFloats(
1333 1333 self.data[self.key][tm][::, ::dx, ::dy].tolist())
1334 1334 else:
1335 1335 if self.key is 'noise':
1336 1336 data = [[x] for x in self.roundFloats(self.data[self.key][tm].tolist())]
1337 1337 else:
1338 1338 data = self.roundFloats(self.data[self.key][tm][::, ::dy].tolist())
1339 1339
1340 1340 meta = {}
1341 1341 ret = {
1342 1342 'plot': plot_name,
1343 1343 'code': self.exp_code,
1344 1344 'time': float(tm),
1345 1345 'data': data,
1346 1346 }
1347 1347 meta['type'] = plot_type
1348 1348 meta['interval'] = float(self.interval)
1349 1349 meta['localtime'] = self.localtime
1350 1350 meta['yrange'] = self.roundFloats(self.heights[::dy].tolist())
1351 1351 if 'spc' in self.data or 'cspc' in self.data:
1352 1352 meta['xrange'] = self.roundFloats(self.xrange[2][::dx].tolist())
1353 1353 else:
1354 1354 meta['xrange'] = []
1355 1355
1356 1356 meta.update(self.meta)
1357 1357 ret['metadata'] = meta
1358 1358 return json.dumps(ret)
1359 1359
1360 1360 @property
1361 1361 def times(self):
1362 1362 '''
1363 1363 Return the list of times of the current data
1364 1364 '''
1365 1365
1366 ret = numpy.array([*self.data[self.key]])
1366 ret = numpy.array([t for t in self.data[self.key]])
1367 1367 if self:
1368 1368 ret.sort()
1369 1369 return ret
1370 1370
1371 1371 @property
1372 1372 def min_time(self):
1373 1373 '''
1374 1374 Return the minimun time value
1375 1375 '''
1376 1376
1377 1377 return self.times[0]
1378 1378
1379 1379 @property
1380 1380 def max_time(self):
1381 1381 '''
1382 1382 Return the maximun time value
1383 1383 '''
1384 1384
1385 1385 return self.times[-1]
1386 1386
1387 1387 @property
1388 1388 def heights(self):
1389 1389 '''
1390 1390 Return the list of heights of the current data
1391 1391 '''
1392 1392
1393 1393 return numpy.array(self.__heights[-1])
1394 1394
1395 1395 @staticmethod
1396 1396 def roundFloats(obj):
1397 1397 if isinstance(obj, list):
1398 1398 return list(map(PlotterData.roundFloats, obj))
1399 1399 elif isinstance(obj, float):
1400 1400 return round(obj, 2)
Show file before | Show file after | Hide comments
@@ -1,713 +1,716
1 1
2 2 import os
3 3 import sys
4 4 import zmq
5 5 import time
6 6 import numpy
7 7 import datetime
8 from queue import Queue
8 try:
9 from queue import Queue
10 except:
11 from Queue import Queue
9 12 from functools import wraps
10 13 from threading import Thread
11 14 import matplotlib
12 15
13 16 if 'BACKEND' in os.environ:
14 17 matplotlib.use(os.environ['BACKEND'])
15 18 elif 'linux' in sys.platform:
16 19 matplotlib.use("TkAgg")
17 20 elif 'darwin' in sys.platform:
18 21 matplotlib.use('WxAgg')
19 22 else:
20 23 from schainpy.utils import log
21 24 log.warning('Using default Backend="Agg"', 'INFO')
22 25 matplotlib.use('Agg')
23 26
24 27 import matplotlib.pyplot as plt
25 28 from matplotlib.patches import Polygon
26 29 from mpl_toolkits.axes_grid1 import make_axes_locatable
27 30 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
28 31
29 32 from schainpy.model.data.jrodata import PlotterData
30 33 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
31 34 from schainpy.utils import log
32 35
33 36 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
34 37 blu_values = matplotlib.pyplot.get_cmap(
35 38 'seismic_r', 20)(numpy.arange(20))[10:15]
36 39 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
37 40 'jro', numpy.vstack((blu_values, jet_values)))
38 41 matplotlib.pyplot.register_cmap(cmap=ncmap)
39 42
40 43 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis',
41 44 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm')]
42 45
43 46 EARTH_RADIUS = 6.3710e3
44 47
45 48 def ll2xy(lat1, lon1, lat2, lon2):
46 49
47 50 p = 0.017453292519943295
48 51 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
49 52 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
50 53 r = 12742 * numpy.arcsin(numpy.sqrt(a))
51 54 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
52 55 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
53 56 theta = -theta + numpy.pi/2
54 57 return r*numpy.cos(theta), r*numpy.sin(theta)
55 58
56 59
57 60 def km2deg(km):
58 61 '''
59 62 Convert distance in km to degrees
60 63 '''
61 64
62 65 return numpy.rad2deg(km/EARTH_RADIUS)
63 66
64 67
65 68 def figpause(interval):
66 69 backend = plt.rcParams['backend']
67 70 if backend in matplotlib.rcsetup.interactive_bk:
68 71 figManager = matplotlib._pylab_helpers.Gcf.get_active()
69 72 if figManager is not None:
70 73 canvas = figManager.canvas
71 74 if canvas.figure.stale:
72 75 canvas.draw()
73 76 try:
74 77 canvas.start_event_loop(interval)
75 78 except:
76 79 pass
77 80 return
78 81
79 82
80 83 def popup(message):
81 84 '''
82 85 '''
83 86
84 87 fig = plt.figure(figsize=(12, 8), facecolor='r')
85 88 text = '\n'.join([s.strip() for s in message.split(':')])
86 89 fig.text(0.01, 0.5, text, ha='left', va='center',
87 90 size='20', weight='heavy', color='w')
88 91 fig.show()
89 92 figpause(1000)
90 93
91 94
92 95 class Throttle(object):
93 96 '''
94 97 Decorator that prevents a function from being called more than once every
95 98 time period.
96 99 To create a function that cannot be called more than once a minute, but
97 100 will sleep until it can be called:
98 101 @Throttle(minutes=1)
99 102 def foo():
100 103 pass
101 104
102 105 for i in range(10):
103 106 foo()
104 107 print "This function has run %s times." % i
105 108 '''
106 109
107 110 def __init__(self, seconds=0, minutes=0, hours=0):
108 111 self.throttle_period = datetime.timedelta(
109 112 seconds=seconds, minutes=minutes, hours=hours
110 113 )
111 114
112 115 self.time_of_last_call = datetime.datetime.min
113 116
114 117 def __call__(self, fn):
115 118 @wraps(fn)
116 119 def wrapper(*args, **kwargs):
117 120 coerce = kwargs.pop('coerce', None)
118 121 if coerce:
119 122 self.time_of_last_call = datetime.datetime.now()
120 123 return fn(*args, **kwargs)
121 124 else:
122 125 now = datetime.datetime.now()
123 126 time_since_last_call = now - self.time_of_last_call
124 127 time_left = self.throttle_period - time_since_last_call
125 128
126 129 if time_left > datetime.timedelta(seconds=0):
127 130 return
128 131
129 132 self.time_of_last_call = datetime.datetime.now()
130 133 return fn(*args, **kwargs)
131 134
132 135 return wrapper
133 136
134 137 def apply_throttle(value):
135 138
136 139 @Throttle(seconds=value)
137 140 def fnThrottled(fn):
138 141 fn()
139 142
140 143 return fnThrottled
141 144
142 145
143 146 @MPDecorator
144 147 class Plot(Operation):
145 148 '''
146 149 Base class for Schain plotting operations
147 150 '''
148 151
149 152 CODE = 'Figure'
150 153 colormap = 'jet'
151 154 bgcolor = 'white'
152 155 buffering = True
153 156 __missing = 1E30
154 157
155 158 __attrs__ = ['show', 'save', 'ymin', 'ymax', 'zmin', 'zmax', 'title',
156 159 'showprofile']
157 160
158 161 def __init__(self):
159 162
160 163 Operation.__init__(self)
161 164 self.isConfig = False
162 165 self.isPlotConfig = False
163 166 self.save_counter = 1
164 167 self.sender_time = 0
165 168 self.data = None
166 169 self.firsttime = True
167 170 self.sender_queue = Queue(maxsize=60)
168 171 self.plots_adjust = {'left': 0.125, 'right': 0.9, 'bottom': 0.15, 'top': 0.9, 'wspace': 0.2, 'hspace': 0.2}
169 172
170 173 def __fmtTime(self, x, pos):
171 174 '''
172 175 '''
173 176
174 177 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
175 178
176 179 def __setup(self, **kwargs):
177 180 '''
178 181 Initialize variables
179 182 '''
180 183
181 184 self.figures = []
182 185 self.axes = []
183 186 self.cb_axes = []
184 187 self.localtime = kwargs.pop('localtime', True)
185 188 self.show = kwargs.get('show', True)
186 189 self.save = kwargs.get('save', False)
187 190 self.save_period = kwargs.get('save_period', 1)
188 191 self.colormap = kwargs.get('colormap', self.colormap)
189 192 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
190 193 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
191 194 self.colormaps = kwargs.get('colormaps', None)
192 195 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
193 196 self.showprofile = kwargs.get('showprofile', False)
194 197 self.title = kwargs.get('wintitle', self.CODE.upper())
195 198 self.cb_label = kwargs.get('cb_label', None)
196 199 self.cb_labels = kwargs.get('cb_labels', None)
197 200 self.labels = kwargs.get('labels', None)
198 201 self.xaxis = kwargs.get('xaxis', 'frequency')
199 202 self.zmin = kwargs.get('zmin', None)
200 203 self.zmax = kwargs.get('zmax', None)
201 204 self.zlimits = kwargs.get('zlimits', None)
202 205 self.xmin = kwargs.get('xmin', None)
203 206 self.xmax = kwargs.get('xmax', None)
204 207 self.xrange = kwargs.get('xrange', 12)
205 208 self.xscale = kwargs.get('xscale', None)
206 209 self.ymin = kwargs.get('ymin', None)
207 210 self.ymax = kwargs.get('ymax', None)
208 211 self.yscale = kwargs.get('yscale', None)
209 212 self.xlabel = kwargs.get('xlabel', None)
210 213 self.attr_time = kwargs.get('attr_time', 'utctime')
211 214 self.decimation = kwargs.get('decimation', None)
212 215 self.showSNR = kwargs.get('showSNR', False)
213 216 self.oneFigure = kwargs.get('oneFigure', True)
214 217 self.width = kwargs.get('width', None)
215 218 self.height = kwargs.get('height', None)
216 219 self.colorbar = kwargs.get('colorbar', True)
217 220 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
218 221 self.channels = kwargs.get('channels', None)
219 222 self.titles = kwargs.get('titles', [])
220 223 self.polar = False
221 224 self.type = kwargs.get('type', 'iq')
222 225 self.grid = kwargs.get('grid', False)
223 226 self.pause = kwargs.get('pause', False)
224 227 self.save_code = kwargs.get('save_code', None)
225 228 self.throttle = kwargs.get('throttle', 0)
226 229 self.exp_code = kwargs.get('exp_code', None)
227 230 self.plot_server = kwargs.get('plot_server', False)
228 231 self.sender_period = kwargs.get('sender_period', 60)
229 232 self.tag = kwargs.get('tag', '')
230 233 self.height_index = kwargs.get('height_index', None)
231 234 self.__throttle_plot = apply_throttle(self.throttle)
232 235 self.data = PlotterData(
233 236 self.CODE, self.throttle, self.exp_code, self.localtime, self.buffering, snr=self.showSNR)
234 237
235 238 if self.plot_server:
236 239 if not self.plot_server.startswith('tcp://'):
237 240 self.plot_server = 'tcp://{}'.format(self.plot_server)
238 241 log.success(
239 242 'Sending to server: {}'.format(self.plot_server),
240 243 self.name
241 244 )
242 245 if 'plot_name' in kwargs:
243 246 self.plot_name = kwargs['plot_name']
244 247
245 248 def __setup_plot(self):
246 249 '''
247 250 Common setup for all figures, here figures and axes are created
248 251 '''
249 252
250 253 self.setup()
251 254
252 255 self.time_label = 'LT' if self.localtime else 'UTC'
253 256
254 257 if self.width is None:
255 258 self.width = 8
256 259
257 260 self.figures = []
258 261 self.axes = []
259 262 self.cb_axes = []
260 263 self.pf_axes = []
261 264 self.cmaps = []
262 265
263 266 size = '15%' if self.ncols == 1 else '30%'
264 267 pad = '4%' if self.ncols == 1 else '8%'
265 268
266 269 if self.oneFigure:
267 270 if self.height is None:
268 271 self.height = 1.4 * self.nrows + 1
269 272 fig = plt.figure(figsize=(self.width, self.height),
270 273 edgecolor='k',
271 274 facecolor='w')
272 275 self.figures.append(fig)
273 276 for n in range(self.nplots):
274 277 ax = fig.add_subplot(self.nrows, self.ncols,
275 278 n + 1, polar=self.polar)
276 279 ax.tick_params(labelsize=8)
277 280 ax.firsttime = True
278 281 ax.index = 0
279 282 ax.press = None
280 283 self.axes.append(ax)
281 284 if self.showprofile:
282 285 cax = self.__add_axes(ax, size=size, pad=pad)
283 286 cax.tick_params(labelsize=8)
284 287 self.pf_axes.append(cax)
285 288 else:
286 289 if self.height is None:
287 290 self.height = 3
288 291 for n in range(self.nplots):
289 292 fig = plt.figure(figsize=(self.width, self.height),
290 293 edgecolor='k',
291 294 facecolor='w')
292 295 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
293 296 ax.tick_params(labelsize=8)
294 297 ax.firsttime = True
295 298 ax.index = 0
296 299 ax.press = None
297 300 self.figures.append(fig)
298 301 self.axes.append(ax)
299 302 if self.showprofile:
300 303 cax = self.__add_axes(ax, size=size, pad=pad)
301 304 cax.tick_params(labelsize=8)
302 305 self.pf_axes.append(cax)
303 306
304 307 for n in range(self.nrows):
305 308 if self.colormaps is not None:
306 309 cmap = plt.get_cmap(self.colormaps[n])
307 310 else:
308 311 cmap = plt.get_cmap(self.colormap)
309 312 cmap.set_bad(self.bgcolor, 1.)
310 313 self.cmaps.append(cmap)
311 314
312 315 def __add_axes(self, ax, size='30%', pad='8%'):
313 316 '''
314 317 Add new axes to the given figure
315 318 '''
316 319 divider = make_axes_locatable(ax)
317 320 nax = divider.new_horizontal(size=size, pad=pad)
318 321 ax.figure.add_axes(nax)
319 322 return nax
320 323
321 324 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
322 325 '''
323 326 Create a masked array for missing data
324 327 '''
325 328 if x_buffer.shape[0] < 2:
326 329 return x_buffer, y_buffer, z_buffer
327 330
328 331 deltas = x_buffer[1:] - x_buffer[0:-1]
329 332 x_median = numpy.median(deltas)
330 333
331 334 index = numpy.where(deltas > 5 * x_median)
332 335
333 336 if len(index[0]) != 0:
334 337 z_buffer[::, index[0], ::] = self.__missing
335 338 z_buffer = numpy.ma.masked_inside(z_buffer,
336 339 0.99 * self.__missing,
337 340 1.01 * self.__missing)
338 341
339 342 return x_buffer, y_buffer, z_buffer
340 343
341 344 def decimate(self):
342 345
343 346 # dx = int(len(self.x)/self.__MAXNUMX) + 1
344 347 dy = int(len(self.y) / self.decimation) + 1
345 348
346 349 # x = self.x[::dx]
347 350 x = self.x
348 351 y = self.y[::dy]
349 352 z = self.z[::, ::, ::dy]
350 353
351 354 return x, y, z
352 355
353 356 def format(self):
354 357 '''
355 358 Set min and max values, labels, ticks and titles
356 359 '''
357 360
358 361 if self.xmin is None:
359 362 xmin = self.data.min_time
360 363 else:
361 364 if self.xaxis is 'time':
362 365 dt = self.getDateTime(self.data.min_time)
363 366 xmin = (dt.replace(hour=int(self.xmin), minute=0, second=0) -
364 367 datetime.datetime(1970, 1, 1)).total_seconds()
365 368 if self.data.localtime:
366 369 xmin += time.timezone
367 370 else:
368 371 xmin = self.xmin
369 372
370 373 if self.xmax is None:
371 374 xmax = xmin + self.xrange * 60 * 60
372 375 else:
373 376 if self.xaxis is 'time':
374 377 dt = self.getDateTime(self.data.max_time)
375 378 xmax = self.xmax - 1
376 379 xmax = (dt.replace(hour=int(xmax), minute=59, second=59) -
377 380 datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds()
378 381 if self.data.localtime:
379 382 xmax += time.timezone
380 383 else:
381 384 xmax = self.xmax
382 385
383 386 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
384 387 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
385 388
386 389 for n, ax in enumerate(self.axes):
387 390 if ax.firsttime:
388 391
389 392 dig = int(numpy.log10(ymax))
390 393 if dig == 0:
391 394 digD = len(str(ymax)) - 2
392 395 ydec = ymax*(10**digD)
393 396
394 397 dig = int(numpy.log10(ydec))
395 398 ystep = ((ydec + (10**(dig)))//10**(dig))*(10**(dig))
396 399 ystep = ystep/5
397 400 ystep = ystep/(10**digD)
398 401
399 402 else:
400 403 ystep = ((ymax + (10**(dig)))//10**(dig))*(10**(dig))
401 404 ystep = ystep/5
402 405
403 406 if self.xaxis is not 'time':
404 407
405 408 dig = int(numpy.log10(xmax))
406 409
407 410 if dig <= 0:
408 411 digD = len(str(xmax)) - 2
409 412 xdec = xmax*(10**digD)
410 413
411 414 dig = int(numpy.log10(xdec))
412 415 xstep = ((xdec + (10**(dig)))//10**(dig))*(10**(dig))
413 416 xstep = xstep*0.5
414 417 xstep = xstep/(10**digD)
415 418
416 419 else:
417 420 xstep = ((xmax + (10**(dig)))//10**(dig))*(10**(dig))
418 421 xstep = xstep/5
419 422
420 423 ax.set_facecolor(self.bgcolor)
421 424 ax.yaxis.set_major_locator(MultipleLocator(ystep))
422 425 if self.xscale:
423 426 ax.xaxis.set_major_formatter(FuncFormatter(
424 427 lambda x, pos: '{0:g}'.format(x*self.xscale)))
425 428 if self.xscale:
426 429 ax.yaxis.set_major_formatter(FuncFormatter(
427 430 lambda x, pos: '{0:g}'.format(x*self.yscale)))
428 431 if self.xaxis is 'time':
429 432 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
430 433 ax.xaxis.set_major_locator(LinearLocator(9))
431 434 else:
432 435 ax.xaxis.set_major_locator(MultipleLocator(xstep))
433 436 if self.xlabel is not None:
434 437 ax.set_xlabel(self.xlabel)
435 438 ax.set_ylabel(self.ylabel)
436 439 ax.firsttime = False
437 440 if self.showprofile:
438 441 self.pf_axes[n].set_ylim(ymin, ymax)
439 442 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
440 443 self.pf_axes[n].set_xlabel('dB')
441 444 self.pf_axes[n].grid(b=True, axis='x')
442 445 [tick.set_visible(False)
443 446 for tick in self.pf_axes[n].get_yticklabels()]
444 447 if self.colorbar:
445 448 ax.cbar = plt.colorbar(
446 449 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
447 450 ax.cbar.ax.tick_params(labelsize=8)
448 451 ax.cbar.ax.press = None
449 452 if self.cb_label:
450 453 ax.cbar.set_label(self.cb_label, size=8)
451 454 elif self.cb_labels:
452 455 ax.cbar.set_label(self.cb_labels[n], size=8)
453 456 else:
454 457 ax.cbar = None
455 458 if self.grid:
456 459 ax.grid(True)
457 460
458 461 if not self.polar:
459 462 ax.set_xlim(xmin, xmax)
460 463 ax.set_ylim(ymin, ymax)
461 464 ax.set_title('{} {} {}'.format(
462 465 self.titles[n],
463 466 self.getDateTime(self.data.max_time).strftime(
464 467 '%Y-%m-%d %H:%M:%S'),
465 468 self.time_label),
466 469 size=8)
467 470 else:
468 471 ax.set_title('{}'.format(self.titles[n]), size=8)
469 472 ax.set_ylim(0, 90)
470 473 ax.set_yticks(numpy.arange(0, 90, 20))
471 474 ax.yaxis.labelpad = 40
472 475
473 476 if self.firsttime:
474 477 for n, fig in enumerate(self.figures):
475 478 fig.subplots_adjust(**self.plots_adjust)
476 479 self.firsttime = False
477 480
478 481 def clear_figures(self):
479 482 '''
480 483 Reset axes for redraw plots
481 484 '''
482 485
483 486 for ax in self.axes+self.pf_axes+self.cb_axes:
484 487 ax.clear()
485 488 ax.firsttime = True
486 489 if hasattr(ax, 'cbar') and ax.cbar:
487 490 ax.cbar.remove()
488 491
489 492 def __plot(self):
490 493 '''
491 494 Main function to plot, format and save figures
492 495 '''
493 496
494 497 self.plot()
495 498 self.format()
496 499
497 500 for n, fig in enumerate(self.figures):
498 501 if self.nrows == 0 or self.nplots == 0:
499 502 log.warning('No data', self.name)
500 503 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
501 504 fig.canvas.manager.set_window_title(self.CODE)
502 505 continue
503 506
504 507 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
505 508 self.getDateTime(self.data.max_time).strftime('%Y/%m/%d')))
506 509 fig.canvas.draw()
507 510 if self.show:
508 511 fig.show()
509 512 figpause(0.01)
510 513
511 514 if self.save:
512 515 self.save_figure(n)
513 516
514 517 if self.plot_server:
515 518 self.send_to_server()
516 519
517 520 def save_figure(self, n):
518 521 '''
519 522 '''
520 523
521 524 if self.save_counter < self.save_period:
522 525 self.save_counter += 1
523 526 return
524 527
525 528 self.save_counter = 1
526 529
527 530 fig = self.figures[n]
528 531
529 532 if self.save_code:
530 533 if isinstance(self.save_code, str):
531 534 labels = [self.save_code for x in self.figures]
532 535 else:
533 536 labels = self.save_code
534 537 else:
535 538 labels = [self.CODE for x in self.figures]
536 539
537 540 figname = os.path.join(
538 541 self.save,
539 542 labels[n],
540 543 '{}_{}.png'.format(
541 544 labels[n],
542 545 self.getDateTime(self.data.max_time).strftime(
543 546 '%Y%m%d_%H%M%S'
544 547 ),
545 548 )
546 549 )
547 550 log.log('Saving figure: {}'.format(figname), self.name)
548 551 if not os.path.isdir(os.path.dirname(figname)):
549 552 os.makedirs(os.path.dirname(figname))
550 553 fig.savefig(figname)
551 554
552 555 if self.throttle == 0:
553 556 figname = os.path.join(
554 557 self.save,
555 558 '{}_{}.png'.format(
556 559 labels[n],
557 560 self.getDateTime(self.data.min_time).strftime(
558 561 '%Y%m%d'
559 562 ),
560 563 )
561 564 )
562 565 fig.savefig(figname)
563 566
564 567 def send_to_server(self):
565 568 '''
566 569 '''
567 570
568 571 interval = self.data.tm - self.sender_time
569 572 if interval < self.sender_period:
570 573 return
571 574
572 575 self.sender_time = self.data.tm
573 576
574 577 attrs = ['titles', 'zmin', 'zmax', 'tag', 'ymin', 'ymax']
575 578 for attr in attrs:
576 579 value = getattr(self, attr)
577 580 if value:
578 581 if isinstance(value, (numpy.float32, numpy.float64)):
579 582 value = round(float(value), 2)
580 583 self.data.meta[attr] = value
581 584 if self.colormap == 'jet':
582 585 self.data.meta['colormap'] = 'Jet'
583 586 elif 'RdBu' in self.colormap:
584 587 self.data.meta['colormap'] = 'RdBu'
585 588 else:
586 589 self.data.meta['colormap'] = 'Viridis'
587 590 self.data.meta['interval'] = int(interval)
588 591 # msg = self.data.jsonify(self.data.tm, self.plot_name, self.plot_type)
589 592 try:
590 593 self.sender_queue.put(self.data.tm, block=False)
591 594 except:
592 595 tm = self.sender_queue.get()
593 596 self.sender_queue.put(self.data.tm)
594 597
595 598 while True:
596 599 if self.sender_queue.empty():
597 600 break
598 601 tm = self.sender_queue.get()
599 602 try:
600 603 msg = self.data.jsonify(tm, self.plot_name, self.plot_type)
601 604 except:
602 605 continue
603 606 self.socket.send_string(msg)
604 607 socks = dict(self.poll.poll(5000))
605 608 if socks.get(self.socket) == zmq.POLLIN:
606 609 reply = self.socket.recv_string()
607 610 if reply == 'ok':
608 611 log.log("Response from server ok", self.name)
609 612 time.sleep(0.2)
610 613 continue
611 614 else:
612 615 log.warning(
613 616 "Malformed reply from server: {}".format(reply), self.name)
614 617 else:
615 618 log.warning(
616 619 "No response from server, retrying...", self.name)
617 620 self.sender_queue.put(self.data.tm)
618 621 self.socket.setsockopt(zmq.LINGER, 0)
619 622 self.socket.close()
620 623 self.poll.unregister(self.socket)
621 624 time.sleep(0.1)
622 625 self.socket = self.context.socket(zmq.REQ)
623 626 self.socket.connect(self.plot_server)
624 627 self.poll.register(self.socket, zmq.POLLIN)
625 628 break
626 629
627 630 def setup(self):
628 631 '''
629 632 This method should be implemented in the child class, the following
630 633 attributes should be set:
631 634
632 635 self.nrows: number of rows
633 636 self.ncols: number of cols
634 637 self.nplots: number of plots (channels or pairs)
635 638 self.ylabel: label for Y axes
636 639 self.titles: list of axes title
637 640
638 641 '''
639 642 raise NotImplementedError
640 643
641 644 def plot(self):
642 645 '''
643 646 Must be defined in the child class
644 647 '''
645 648 raise NotImplementedError
646 649
647 650 def run(self, dataOut, **kwargs):
648 651 '''
649 652 Main plotting routine
650 653 '''
651 654
652 655 if self.isConfig is False:
653 656 self.__setup(**kwargs)
654 657
655 658 t = getattr(dataOut, self.attr_time)
656 659
657 660 if self.localtime:
658 661 self.getDateTime = datetime.datetime.fromtimestamp
659 662 else:
660 663 self.getDateTime = datetime.datetime.utcfromtimestamp
661 664
662 665 if self.xmin is None:
663 666 self.tmin = t
664 667 if 'buffer' in self.plot_type:
665 668 self.xmin = self.getDateTime(t).hour
666 669 else:
667 670 self.tmin = (
668 671 self.getDateTime(t).replace(
669 672 hour=int(self.xmin),
670 673 minute=0,
671 674 second=0) - self.getDateTime(0)).total_seconds()
672 675
673 676 self.data.setup()
674 677 self.isConfig = True
675 678 if self.plot_server:
676 679 self.context = zmq.Context()
677 680 self.socket = self.context.socket(zmq.REQ)
678 681 self.socket.connect(self.plot_server)
679 682 self.poll = zmq.Poller()
680 683 self.poll.register(self.socket, zmq.POLLIN)
681 684
682 685 tm = getattr(dataOut, self.attr_time)
683 686
684 687 if self.data and (tm - self.tmin) >= self.xrange*60*60:
685 688 self.save_counter = self.save_period
686 689 self.__plot()
687 690 if 'time' in self.xaxis:
688 691 self.xmin += self.xrange
689 692 if self.xmin >= 24:
690 693 self.xmin -= 24
691 694 self.tmin += self.xrange*60*60
692 695 self.data.setup()
693 696 self.clear_figures()
694 697
695 698 self.data.update(dataOut, tm)
696 699
697 700 if self.isPlotConfig is False:
698 701 self.__setup_plot()
699 702 self.isPlotConfig = True
700 703
701 704 if self.throttle == 0:
702 705 self.__plot()
703 706 else:
704 707 self.__throttle_plot(self.__plot)#, coerce=coerce)
705 708
706 709 def close(self):
707 710
708 711 if self.data and not self.data.flagNoData:
709 712 self.save_counter = self.save_period
710 713 self.__plot()
711 714 if self.data and not self.data.flagNoData and self.pause:
712 715 figpause(10)
713 716
Show file before | Show file after | Hide comments
@@ -1,207 +1,203
1 1 '''
2 2 Base clases to create Processing units and operations, the MPDecorator
3 3 must be used in plotting and writing operations to allow to run as an
4 4 external process.
5 5 '''
6 6
7 7 import inspect
8 8 import zmq
9 9 import time
10 10 import pickle
11 11 import traceback
12 try:
13 from queue import Queue
14 except:
15 from Queue import Queue
16 12 from threading import Thread
17 13 from multiprocessing import Process, Queue
18 14 from schainpy.utils import log
19 15
20 16
21 17 class ProcessingUnit(object):
22 18 '''
23 19 Base class to create Signal Chain Units
24 20 '''
25 21
26 22 proc_type = 'processing'
27 23
28 24 def __init__(self):
29 25
30 26 self.dataIn = None
31 27 self.dataOut = None
32 28 self.isConfig = False
33 29 self.operations = []
34 30
35 31 def setInput(self, unit):
36 32
37 33 self.dataIn = unit.dataOut
38 34
39 35 def getAllowedArgs(self):
40 36 if hasattr(self, '__attrs__'):
41 37 return self.__attrs__
42 38 else:
43 39 return inspect.getargspec(self.run).args
44 40
45 41 def addOperation(self, conf, operation):
46 42 '''
47 43 '''
48 44
49 45 self.operations.append((operation, conf.type, conf.getKwargs()))
50 46
51 47 def getOperationObj(self, objId):
52 48
53 49 if objId not in list(self.operations.keys()):
54 50 return None
55 51
56 52 return self.operations[objId]
57 53
58 54 def call(self, **kwargs):
59 55 '''
60 56 '''
61 57
62 58 try:
63 59 if self.dataIn is not None and self.dataIn.flagNoData and not self.dataIn.error:
64 60 return self.dataIn.isReady()
65 61 elif self.dataIn is None or not self.dataIn.error:
66 62 self.run(**kwargs)
67 63 elif self.dataIn.error:
68 64 self.dataOut.error = self.dataIn.error
69 65 self.dataOut.flagNoData = True
70 66 except:
71 67 err = traceback.format_exc()
72 68 if 'SchainWarning' in err:
73 69 log.warning(err.split('SchainWarning:')[-1].split('\n')[0].strip(), self.name)
74 70 elif 'SchainError' in err:
75 71 log.error(err.split('SchainError:')[-1].split('\n')[0].strip(), self.name)
76 72 else:
77 73 log.error(err, self.name)
78 74 self.dataOut.error = True
79 75
80 76 for op, optype, opkwargs in self.operations:
81 77 if optype == 'other' and not self.dataOut.flagNoData:
82 78 self.dataOut = op.run(self.dataOut, **opkwargs)
83 79 elif optype == 'external' and not self.dataOut.flagNoData:
84 80 op.queue.put(self.dataOut)
85 81 elif optype == 'external' and self.dataOut.error:
86 82 op.queue.put(self.dataOut)
87 83
88 84 return 'Error' if self.dataOut.error else self.dataOut.isReady()
89 85
90 86 def setup(self):
91 87
92 88 raise NotImplementedError
93 89
94 90 def run(self):
95 91
96 92 raise NotImplementedError
97 93
98 94 def close(self):
99 95
100 96 return
101 97
102 98
103 99 class Operation(object):
104 100
105 101 '''
106 102 '''
107 103
108 104 proc_type = 'operation'
109 105
110 106 def __init__(self):
111 107
112 108 self.id = None
113 109 self.isConfig = False
114 110
115 111 if not hasattr(self, 'name'):
116 112 self.name = self.__class__.__name__
117 113
118 114 def getAllowedArgs(self):
119 115 if hasattr(self, '__attrs__'):
120 116 return self.__attrs__
121 117 else:
122 118 return inspect.getargspec(self.run).args
123 119
124 120 def setup(self):
125 121
126 122 self.isConfig = True
127 123
128 124 raise NotImplementedError
129 125
130 126 def run(self, dataIn, **kwargs):
131 127 """
132 128 Realiza las operaciones necesarias sobre la dataIn.data y actualiza los
133 129 atributos del objeto dataIn.
134 130
135 131 Input:
136 132
137 133 dataIn : objeto del tipo JROData
138 134
139 135 Return:
140 136
141 137 None
142 138
143 139 Affected:
144 140 __buffer : buffer de recepcion de datos.
145 141
146 142 """
147 143 if not self.isConfig:
148 144 self.setup(**kwargs)
149 145
150 146 raise NotImplementedError
151 147
152 148 def close(self):
153 149
154 150 return
155 151
156 152
157 153 def MPDecorator(BaseClass):
158 154 """
159 155 Multiprocessing class decorator
160 156
161 157 This function add multiprocessing features to a BaseClass.
162 158 """
163 159
164 160 class MPClass(BaseClass, Process):
165 161
166 162 def __init__(self, *args, **kwargs):
167 163 super(MPClass, self).__init__()
168 164 Process.__init__(self)
169 165
170 166 self.args = args
171 167 self.kwargs = kwargs
172 168 self.t = time.time()
173 169 self.op_type = 'external'
174 170 self.name = BaseClass.__name__
175 171 self.__doc__ = BaseClass.__doc__
176 172
177 173 if 'plot' in self.name.lower() and not self.name.endswith('_'):
178 174 self.name = '{}{}'.format(self.CODE.upper(), 'Plot')
179 175
180 176 self.start_time = time.time()
181 177 self.err_queue = args[3]
182 178 self.queue = Queue(maxsize=1)
183 179 self.myrun = BaseClass.run
184 180
185 181 def run(self):
186 182
187 183 while True:
188 184
189 185 dataOut = self.queue.get()
190 186
191 187 if not dataOut.error:
192 188 try:
193 189 BaseClass.run(self, dataOut, **self.kwargs)
194 190 except:
195 191 err = traceback.format_exc()
196 192 log.error(err, self.name)
197 193 else:
198 194 break
199 195
200 196 self.close()
201 197
202 198 def close(self):
203 199
204 200 BaseClass.close(self)
205 201 log.success('Done...(Time:{:4.2f} secs)'.format(time.time()-self.start_time), self.name)
206 202
207 203 return MPClass
General Comments 0
You need to be logged in to leave comments. Login now