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