##// END OF EJS Templates
Rewrite Param reader/writer as HDF reader/writer
Juan C. Espinoza -
r1326:6584204c3602
parent child
Show More
@@ -1,1403 +1,1402
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 self.timeZone = None
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 self.timeZone = 0
486 487 # self.data = None
487 488 # self.dtype = None
488 489 # self.nChannels = 0
489 490 # self.nHeights = 0
490 491 self.nProfiles = None
491 492 self.heightList = None
492 493 self.channelList = None
493 494 # self.channelIndexList = None
494 495 self.pairsList = None
495 496 self.flagNoData = True
496 497 self.flagDiscontinuousBlock = False
497 498 self.utctime = None
498 499 self.nCohInt = None
499 500 self.nIncohInt = None
500 501 self.blocksize = None
501 502 self.nFFTPoints = None
502 503 self.wavelength = None
503 504 self.flagDecodeData = False # asumo q la data no esta decodificada
504 505 self.flagDeflipData = False # asumo q la data no esta sin flip
505 506 self.flagShiftFFT = False
506 507 self.ippFactor = 1
507 508 #self.noise = None
508 509 self.beacon_heiIndexList = []
509 510 self.noise_estimation = None
510 511
511 512 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
512 513 """
513 514 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
514 515
515 516 Return:
516 517 noiselevel
517 518 """
518 519
519 520 noise = numpy.zeros(self.nChannels)
520 521
521 522 for channel in range(self.nChannels):
522 523 daux = self.data_spc[channel,
523 524 xmin_index:xmax_index, ymin_index:ymax_index]
524 525 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
525 526
526 527 return noise
527 528
528 529 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
529 530
530 531 if self.noise_estimation is not None:
531 532 # this was estimated by getNoise Operation defined in jroproc_spectra.py
532 533 return self.noise_estimation
533 534 else:
534 535 noise = self.getNoisebyHildebrand(
535 536 xmin_index, xmax_index, ymin_index, ymax_index)
536 537 return noise
537 538
538 539 def getFreqRangeTimeResponse(self, extrapoints=0):
539 540
540 541 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
541 542 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.) - deltafreq / 2
542 543
543 544 return freqrange
544 545
545 546 def getAcfRange(self, extrapoints=0):
546 547
547 548 deltafreq = 10. / (self.getFmax() / (self.nFFTPoints * self.ippFactor))
548 549 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
549 550
550 551 return freqrange
551 552
552 553 def getFreqRange(self, extrapoints=0):
553 554
554 555 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
555 556 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
556 557
557 558 return freqrange
558 559
559 560 def getVelRange(self, extrapoints=0):
560 561
561 562 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
562 563 velrange = deltav * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.)
563 564
564 565 if self.nmodes:
565 566 return velrange/self.nmodes
566 567 else:
567 568 return velrange
568 569
569 570 def getNPairs(self):
570 571
571 572 return len(self.pairsList)
572 573
573 574 def getPairsIndexList(self):
574 575
575 576 return list(range(self.nPairs))
576 577
577 578 def getNormFactor(self):
578 579
579 580 pwcode = 1
580 581
581 582 if self.flagDecodeData:
582 583 pwcode = numpy.sum(self.code[0]**2)
583 584 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
584 585 normFactor = self.nProfiles * self.nIncohInt * self.nCohInt * pwcode * self.windowOfFilter
585 586
586 587 return normFactor
587 588
588 589 def getFlagCspc(self):
589 590
590 591 if self.data_cspc is None:
591 592 return True
592 593
593 594 return False
594 595
595 596 def getFlagDc(self):
596 597
597 598 if self.data_dc is None:
598 599 return True
599 600
600 601 return False
601 602
602 603 def getTimeInterval(self):
603 604
604 605 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
605 606 if self.nmodes:
606 607 return self.nmodes*timeInterval
607 608 else:
608 609 return timeInterval
609 610
610 611 def getPower(self):
611 612
612 613 factor = self.normFactor
613 614 z = self.data_spc / factor
614 615 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
615 616 avg = numpy.average(z, axis=1)
616 617
617 618 return 10 * numpy.log10(avg)
618 619
619 620 def getCoherence(self, pairsList=None, phase=False):
620 621
621 622 z = []
622 623 if pairsList is None:
623 624 pairsIndexList = self.pairsIndexList
624 625 else:
625 626 pairsIndexList = []
626 627 for pair in pairsList:
627 628 if pair not in self.pairsList:
628 629 raise ValueError("Pair %s is not in dataOut.pairsList" % (
629 630 pair))
630 631 pairsIndexList.append(self.pairsList.index(pair))
631 632 for i in range(len(pairsIndexList)):
632 633 pair = self.pairsList[pairsIndexList[i]]
633 634 ccf = numpy.average(self.data_cspc[pairsIndexList[i], :, :], axis=0)
634 635 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
635 636 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
636 637 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
637 638 if phase:
638 639 data = numpy.arctan2(avgcoherenceComplex.imag,
639 640 avgcoherenceComplex.real) * 180 / numpy.pi
640 641 else:
641 642 data = numpy.abs(avgcoherenceComplex)
642 643
643 644 z.append(data)
644 645
645 646 return numpy.array(z)
646 647
647 648 def setValue(self, value):
648 649
649 650 print("This property should not be initialized")
650 651
651 652 return
652 653
653 654 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
654 655 pairsIndexList = property(
655 656 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
656 657 normFactor = property(getNormFactor, setValue,
657 658 "I'm the 'getNormFactor' property.")
658 659 flag_cspc = property(getFlagCspc, setValue)
659 660 flag_dc = property(getFlagDc, setValue)
660 661 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
661 662 timeInterval = property(getTimeInterval, setValue,
662 663 "I'm the 'timeInterval' property")
663 664
664 665
665 666 class SpectraHeis(Spectra):
666 667
667 668 data_spc = None
668 669 data_cspc = None
669 670 data_dc = None
670 671 nFFTPoints = None
671 672 # nPairs = None
672 673 pairsList = None
673 674 nCohInt = None
674 675 nIncohInt = None
675 676
676 677 def __init__(self):
677 678
678 679 self.radarControllerHeaderObj = RadarControllerHeader()
679 680
680 681 self.systemHeaderObj = SystemHeader()
681 682
682 683 self.type = "SpectraHeis"
683 684
684 685 # self.dtype = None
685 686
686 687 # self.nChannels = 0
687 688
688 689 # self.nHeights = 0
689 690
690 691 self.nProfiles = None
691 692
692 693 self.heightList = None
693 694
694 695 self.channelList = None
695 696
696 697 # self.channelIndexList = None
697 698
698 699 self.flagNoData = True
699 700
700 701 self.flagDiscontinuousBlock = False
701 702
702 703 # self.nPairs = 0
703 704
704 705 self.utctime = None
705 706
706 707 self.blocksize = None
707 708
708 709 self.profileIndex = 0
709 710
710 711 self.nCohInt = 1
711 712
712 713 self.nIncohInt = 1
713 714
714 715 def getNormFactor(self):
715 716 pwcode = 1
716 717 if self.flagDecodeData:
717 718 pwcode = numpy.sum(self.code[0]**2)
718 719
719 720 normFactor = self.nIncohInt * self.nCohInt * pwcode
720 721
721 722 return normFactor
722 723
723 724 def getTimeInterval(self):
724 725
725 726 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
726 727
727 728 return timeInterval
728 729
729 730 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
730 731 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
731 732
732 733
733 734 class Fits(JROData):
734 735
735 736 heightList = None
736 737 channelList = None
737 738 flagNoData = True
738 739 flagDiscontinuousBlock = False
739 740 useLocalTime = False
740 741 utctime = None
741 timeZone = 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 # self.timeZone = None
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 self.timeZone = None
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
1070 1069 self.systemHeaderObj = SystemHeader()
1071
1072 1070 self.type = "Parameters"
1071 self.timeZone = 0
1073 1072
1074 1073 def getTimeRange1(self, interval):
1075 1074
1076 1075 datatime = []
1077 1076
1078 1077 if self.useLocalTime:
1079 1078 time1 = self.utctimeInit - self.timeZone * 60
1080 1079 else:
1081 1080 time1 = self.utctimeInit
1082 1081
1083 1082 datatime.append(time1)
1084 1083 datatime.append(time1 + interval)
1085 1084 datatime = numpy.array(datatime)
1086 1085
1087 1086 return datatime
1088 1087
1089 1088 def getTimeInterval(self):
1090 1089
1091 1090 if hasattr(self, 'timeInterval1'):
1092 1091 return self.timeInterval1
1093 1092 else:
1094 1093 return self.paramInterval
1095 1094
1096 1095 def setValue(self, value):
1097 1096
1098 1097 print("This property should not be initialized")
1099 1098
1100 1099 return
1101 1100
1102 1101 def getNoise(self):
1103 1102
1104 1103 return self.spc_noise
1105 1104
1106 1105 timeInterval = property(getTimeInterval)
1107 1106 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
1108 1107
1109 1108
1110 1109 class PlotterData(object):
1111 1110 '''
1112 1111 Object to hold data to be plotted
1113 1112 '''
1114 1113
1115 1114 MAXNUMX = 200
1116 1115 MAXNUMY = 200
1117 1116
1118 1117 def __init__(self, code, throttle_value, exp_code, buffering=True, snr=False):
1119 1118
1120 1119 self.key = code
1121 1120 self.throttle = throttle_value
1122 1121 self.exp_code = exp_code
1123 1122 self.buffering = buffering
1124 1123 self.ready = False
1125 1124 self.flagNoData = False
1126 1125 self.localtime = False
1127 1126 self.data = {}
1128 1127 self.meta = {}
1129 1128 self.__heights = []
1130 1129
1131 1130 if 'snr' in code:
1132 1131 self.plottypes = ['snr']
1133 1132 elif code == 'spc':
1134 1133 self.plottypes = ['spc', 'noise', 'rti']
1135 1134 elif code == 'cspc':
1136 1135 self.plottypes = ['cspc', 'spc', 'noise', 'rti']
1137 1136 elif code == 'rti':
1138 1137 self.plottypes = ['noise', 'rti']
1139 1138 else:
1140 1139 self.plottypes = [code]
1141 1140
1142 1141 if 'snr' not in self.plottypes and snr:
1143 1142 self.plottypes.append('snr')
1144 1143
1145 1144 for plot in self.plottypes:
1146 1145 self.data[plot] = {}
1147 1146
1148 1147
1149 1148 def __str__(self):
1150 1149 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1151 1150 return 'Data[{}][{}]'.format(';'.join(dum), len(self.times))
1152 1151
1153 1152 def __len__(self):
1154 1153 return len(self.data[self.key])
1155 1154
1156 1155 def __getitem__(self, key):
1157 1156
1158 1157 if key not in self.data:
1159 1158 raise KeyError(log.error('Missing key: {}'.format(key)))
1160 1159 if 'spc' in key or not self.buffering:
1161 1160 ret = self.data[key][self.tm]
1162 1161 elif 'scope' in key:
1163 1162 ret = numpy.array(self.data[key][float(self.tm)])
1164 1163 else:
1165 1164 ret = numpy.array([self.data[key][x] for x in self.times])
1166 1165 if ret.ndim > 1:
1167 1166 ret = numpy.swapaxes(ret, 0, 1)
1168 1167 return ret
1169 1168
1170 1169 def __contains__(self, key):
1171 1170 return key in self.data
1172 1171
1173 1172 def setup(self):
1174 1173 '''
1175 1174 Configure object
1176 1175 '''
1177 1176 self.type = ''
1178 1177 self.ready = False
1179 1178 del self.data
1180 1179 self.data = {}
1181 1180 self.__heights = []
1182 1181 self.__all_heights = set()
1183 1182 for plot in self.plottypes:
1184 1183 if 'snr' in plot:
1185 1184 plot = 'snr'
1186 1185 elif 'spc_moments' == plot:
1187 1186 plot = 'moments'
1188 1187 self.data[plot] = {}
1189 1188
1190 1189 if 'spc' in self.data or 'rti' in self.data or 'cspc' in self.data or 'moments' in self.data:
1191 1190 self.data['noise'] = {}
1192 1191 self.data['rti'] = {}
1193 1192 if 'noise' not in self.plottypes:
1194 1193 self.plottypes.append('noise')
1195 1194 if 'rti' not in self.plottypes:
1196 1195 self.plottypes.append('rti')
1197 1196
1198 1197 def shape(self, key):
1199 1198 '''
1200 1199 Get the shape of the one-element data for the given key
1201 1200 '''
1202 1201
1203 1202 if len(self.data[key]):
1204 1203 if 'spc' in key or not self.buffering:
1205 1204 return self.data[key].shape
1206 1205 return self.data[key][self.times[0]].shape
1207 1206 return (0,)
1208 1207
1209 1208 def update(self, dataOut, tm):
1210 1209 '''
1211 1210 Update data object with new dataOut
1212 1211 '''
1213 1212
1214 1213 self.profileIndex = dataOut.profileIndex
1215 1214 self.tm = tm
1216 1215 self.type = dataOut.type
1217 1216 self.parameters = getattr(dataOut, 'parameters', [])
1218 1217
1219 1218 if hasattr(dataOut, 'meta'):
1220 1219 self.meta.update(dataOut.meta)
1221 1220
1222 1221 if hasattr(dataOut, 'pairsList'):
1223 1222 self.pairs = dataOut.pairsList
1224 1223
1225 1224 self.interval = dataOut.getTimeInterval()
1226 1225 self.localtime = dataOut.useLocalTime
1227 1226 if True in ['spc' in ptype for ptype in self.plottypes]:
1228 1227 self.xrange = (dataOut.getFreqRange(1)/1000.,
1229 1228 dataOut.getAcfRange(1), dataOut.getVelRange(1))
1230 1229 self.__heights.append(dataOut.heightList)
1231 1230 self.__all_heights.update(dataOut.heightList)
1232 1231
1233 1232 for plot in self.plottypes:
1234 1233 if plot in ('spc', 'spc_moments', 'spc_cut'):
1235 1234 z = dataOut.data_spc/dataOut.normFactor
1236 1235 buffer = 10*numpy.log10(z)
1237 1236 if plot == 'cspc':
1238 1237 buffer = (dataOut.data_spc, dataOut.data_cspc)
1239 1238 if plot == 'noise':
1240 1239 buffer = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
1241 1240 if plot in ('rti', 'spcprofile'):
1242 1241 buffer = dataOut.getPower()
1243 1242 if plot == 'snr_db':
1244 1243 buffer = dataOut.data_SNR
1245 1244 if plot == 'snr':
1246 1245 buffer = 10*numpy.log10(dataOut.data_SNR)
1247 1246 if plot == 'dop':
1248 1247 buffer = dataOut.data_DOP
1249 1248 if plot == 'pow':
1250 1249 buffer = 10*numpy.log10(dataOut.data_POW)
1251 1250 if plot == 'width':
1252 1251 buffer = dataOut.data_WIDTH
1253 1252 if plot == 'coh':
1254 1253 buffer = dataOut.getCoherence()
1255 1254 if plot == 'phase':
1256 1255 buffer = dataOut.getCoherence(phase=True)
1257 1256 if plot == 'output':
1258 1257 buffer = dataOut.data_output
1259 1258 if plot == 'param':
1260 1259 buffer = dataOut.data_param
1261 1260 if plot == 'scope':
1262 1261 buffer = dataOut.data
1263 1262 self.flagDataAsBlock = dataOut.flagDataAsBlock
1264 1263 self.nProfiles = dataOut.nProfiles
1265 1264 if plot == 'pp_power':
1266 1265 buffer = dataOut.dataPP_POWER
1267 1266 self.flagDataAsBlock = dataOut.flagDataAsBlock
1268 1267 self.nProfiles = dataOut.nProfiles
1269 1268 if plot == 'pp_signal':
1270 1269 buffer = dataOut.dataPP_POW
1271 1270 self.flagDataAsBlock = dataOut.flagDataAsBlock
1272 1271 self.nProfiles = dataOut.nProfiles
1273 1272 if plot == 'pp_velocity':
1274 1273 buffer = dataOut.dataPP_DOP
1275 1274 self.flagDataAsBlock = dataOut.flagDataAsBlock
1276 1275 self.nProfiles = dataOut.nProfiles
1277 1276 if plot == 'pp_specwidth':
1278 1277 buffer = dataOut.dataPP_WIDTH
1279 1278 self.flagDataAsBlock = dataOut.flagDataAsBlock
1280 1279 self.nProfiles = dataOut.nProfiles
1281 1280
1282 1281 if plot == 'spc':
1283 1282 self.data['spc'][tm] = buffer
1284 1283 elif plot == 'cspc':
1285 1284 self.data['cspc'][tm] = buffer
1286 1285 elif plot == 'spc_moments':
1287 1286 self.data['spc'][tm] = buffer
1288 1287 self.data['moments'][tm] = dataOut.moments
1289 1288 else:
1290 1289 if self.buffering:
1291 1290 self.data[plot][tm] = buffer
1292 1291 else:
1293 1292 self.data[plot][tm] = buffer
1294 1293
1295 1294 if dataOut.channelList is None:
1296 1295 self.channels = range(buffer.shape[0])
1297 1296 else:
1298 1297 self.channels = dataOut.channelList
1299 1298
1300 1299 if buffer is None:
1301 1300 self.flagNoData = True
1302 1301 raise schainpy.admin.SchainWarning('Attribute data_{} is empty'.format(self.key))
1303 1302
1304 1303 def normalize_heights(self):
1305 1304 '''
1306 1305 Ensure same-dimension of the data for different heighList
1307 1306 '''
1308 1307
1309 1308 H = numpy.array(list(self.__all_heights))
1310 1309 H.sort()
1311 1310 for key in self.data:
1312 1311 shape = self.shape(key)[:-1] + H.shape
1313 1312 for tm, obj in list(self.data[key].items()):
1314 1313 h = self.__heights[self.times.tolist().index(tm)]
1315 1314 if H.size == h.size:
1316 1315 continue
1317 1316 index = numpy.where(numpy.in1d(H, h))[0]
1318 1317 dummy = numpy.zeros(shape) + numpy.nan
1319 1318 if len(shape) == 2:
1320 1319 dummy[:, index] = obj
1321 1320 else:
1322 1321 dummy[index] = obj
1323 1322 self.data[key][tm] = dummy
1324 1323
1325 1324 self.__heights = [H for tm in self.times]
1326 1325
1327 1326 def jsonify(self, tm, plot_name, plot_type, decimate=False):
1328 1327 '''
1329 1328 Convert data to json
1330 1329 '''
1331 1330
1332 1331 dy = int(self.heights.size/self.MAXNUMY) + 1
1333 1332 if self.key in ('spc', 'cspc'):
1334 1333 dx = int(self.data[self.key][tm].shape[1]/self.MAXNUMX) + 1
1335 1334 data = self.roundFloats(
1336 1335 self.data[self.key][tm][::, ::dx, ::dy].tolist())
1337 1336 else:
1338 1337 if self.key is 'noise':
1339 1338 data = [[x] for x in self.roundFloats(self.data[self.key][tm].tolist())]
1340 1339 else:
1341 1340 data = self.roundFloats(self.data[self.key][tm][::, ::dy].tolist())
1342 1341
1343 1342 meta = {}
1344 1343 ret = {
1345 1344 'plot': plot_name,
1346 1345 'code': self.exp_code,
1347 1346 'time': float(tm),
1348 1347 'data': data,
1349 1348 }
1350 1349 meta['type'] = plot_type
1351 1350 meta['interval'] = float(self.interval)
1352 1351 meta['localtime'] = self.localtime
1353 1352 meta['yrange'] = self.roundFloats(self.heights[::dy].tolist())
1354 1353 if 'spc' in self.data or 'cspc' in self.data:
1355 1354 meta['xrange'] = self.roundFloats(self.xrange[2][::dx].tolist())
1356 1355 else:
1357 1356 meta['xrange'] = []
1358 1357
1359 1358 meta.update(self.meta)
1360 1359 ret['metadata'] = meta
1361 1360 return json.dumps(ret)
1362 1361
1363 1362 @property
1364 1363 def times(self):
1365 1364 '''
1366 1365 Return the list of times of the current data
1367 1366 '''
1368 1367
1369 1368 ret = numpy.array([*self.data[self.key]])
1370 1369 if self:
1371 1370 ret.sort()
1372 1371 return ret
1373 1372
1374 1373 @property
1375 1374 def min_time(self):
1376 1375 '''
1377 1376 Return the minimun time value
1378 1377 '''
1379 1378
1380 1379 return self.times[0]
1381 1380
1382 1381 @property
1383 1382 def max_time(self):
1384 1383 '''
1385 1384 Return the maximun time value
1386 1385 '''
1387 1386
1388 1387 return self.times[-1]
1389 1388
1390 1389 @property
1391 1390 def heights(self):
1392 1391 '''
1393 1392 Return the list of heights of the current data
1394 1393 '''
1395 1394
1396 1395 return numpy.array(self.__heights[-1])
1397 1396
1398 1397 @staticmethod
1399 1398 def roundFloats(obj):
1400 1399 if isinstance(obj, list):
1401 1400 return list(map(PlotterData.roundFloats, obj))
1402 1401 elif isinstance(obj, float):
1403 1402 return round(obj, 2)
This diff has been collapsed as it changes many lines, (1329 lines changed) Show them Hide them
@@ -1,1435 +1,620
1 import numpy
2 import time
3 1 import os
4 import h5py
5 import re
2 import time
6 3 import datetime
7 4
5 import numpy
6 import h5py
7
8 8 import schainpy.admin
9 9 from schainpy.model.data.jrodata import *
10 10 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
11 11 from schainpy.model.io.jroIO_base import *
12 12 from schainpy.utils import log
13 13
14 14
15 class ParamReader(JRODataReader,ProcessingUnit):
16 '''
17 Reads HDF5 format files
18 path
19 startDate
20 endDate
21 startTime
22 endTime
23 '''
24
25 ext = ".hdf5"
26 optchar = "D"
27 timezone = None
28 startTime = None
29 endTime = None
30 fileIndex = None
31 utcList = None #To select data in the utctime list
32 blockList = None #List to blocks to be read from the file
33 blocksPerFile = None #Number of blocks to be read
34 blockIndex = None
35 path = None
36 #List of Files
37 filenameList = None
38 datetimeList = None
39 #Hdf5 File
40 listMetaname = None
41 listMeta = None
42 listDataname = None
43 listData = None
44 listShapes = None
45 fp = None
46 #dataOut reconstruction
47 dataOut = None
48
49 def __init__(self):#, **kwargs):
50 ProcessingUnit.__init__(self) #, **kwargs)
51 self.dataOut = Parameters()
52 return
53
54 def setup(self, **kwargs):
55
56 path = kwargs['path']
57 startDate = kwargs['startDate']
58 endDate = kwargs['endDate']
59 startTime = kwargs['startTime']
60 endTime = kwargs['endTime']
61 walk = kwargs['walk']
62 if 'ext' in kwargs:
63 ext = kwargs['ext']
64 else:
65 ext = '.hdf5'
66 if 'timezone' in kwargs:
67 self.timezone = kwargs['timezone']
68 else:
69 self.timezone = 'lt'
70
71 print("[Reading] Searching files in offline mode ...")
72 pathList, filenameList = self.searchFilesOffLine(path, startDate=startDate, endDate=endDate,
73 startTime=startTime, endTime=endTime,
74 ext=ext, walk=walk)
75
76 if not(filenameList):
77 print("There is no files into the folder: %s"%(path))
78 sys.exit(-1)
79
80 self.fileIndex = -1
81 self.startTime = startTime
82 self.endTime = endTime
83
84 self.__readMetadata()
85
86 self.__setNextFileOffline()
87
88 return
89
90 def searchFilesOffLine(self,
91 path,
92 startDate=None,
93 endDate=None,
94 startTime=datetime.time(0,0,0),
95 endTime=datetime.time(23,59,59),
96 ext='.hdf5',
97 walk=True):
98
99 expLabel = ''
100 self.filenameList = []
101 self.datetimeList = []
102
103 pathList = []
104
105 JRODataObj = JRODataReader()
106 dateList, pathList = JRODataObj.findDatafiles(path, startDate, endDate, expLabel, ext, walk, include_path=True)
107
108 if dateList == []:
109 print("[Reading] No *%s files in %s from %s to %s)"%(ext, path,
110 datetime.datetime.combine(startDate,startTime).ctime(),
111 datetime.datetime.combine(endDate,endTime).ctime()))
112
113 return None, None
114
115 if len(dateList) > 1:
116 print("[Reading] %d days were found in date range: %s - %s" %(len(dateList), startDate, endDate))
117 else:
118 print("[Reading] data was found for the date %s" %(dateList[0]))
119
120 filenameList = []
121 datetimeList = []
122
123 #----------------------------------------------------------------------------------
124
125 for thisPath in pathList:
126
127 fileList = glob.glob1(thisPath, "*%s" %ext)
128 fileList.sort()
129
130 for file in fileList:
131
132 filename = os.path.join(thisPath,file)
133
134 if not isFileInDateRange(filename, startDate, endDate):
135 continue
136
137 thisDatetime = self.__isFileInTimeRange(filename, startDate, endDate, startTime, endTime)
138
139 if not(thisDatetime):
140 continue
141
142 filenameList.append(filename)
143 datetimeList.append(thisDatetime)
144
145 if not(filenameList):
146 print("[Reading] Any file was found int time range %s - %s" %(datetime.datetime.combine(startDate,startTime).ctime(), datetime.datetime.combine(endDate,endTime).ctime()))
147 return None, None
148
149 print("[Reading] %d file(s) was(were) found in time range: %s - %s" %(len(filenameList), startTime, endTime))
150 print()
151
152 self.filenameList = filenameList
153 self.datetimeList = datetimeList
154
155 return pathList, filenameList
156
157 def __isFileInTimeRange(self,filename, startDate, endDate, startTime, endTime):
15 class HDFReader(Reader, ProcessingUnit):
16 """Processing unit to read HDF5 format files
17
18 This unit reads HDF5 files created with `HDFWriter` operation contains
19 by default two groups Data and Metadata all variables would be saved as `dataOut`
20 attributes.
21 It is possible to read any HDF5 file by given the structure in the `description`
22 parameter, also you can add extra values to metadata with the parameter `extras`.
23
24 Parameters:
25 -----------
26 path : str
27 Path where files are located.
28 startDate : date
29 Start date of the files
30 endDate : list
31 End date of the files
32 startTime : time
33 Start time of the files
34 endTime : time
35 End time of the files
36 description : dict, optional
37 Dictionary with the description of the HDF5 file
38 extras : dict, optional
39 Dictionary with extra metadata to be be added to `dataOut`
40
41 Examples
42 --------
43
44 desc = {
45 'Data': {
46 'data_output': ['u', 'v', 'w'],
47 'utctime': 'timestamps',
48 } ,
49 'Metadata': {
50 'heightList': 'heights'
51 }
52 }
53
54 desc = {
55 'Data': {
56 'data_output': 'winds',
57 'utctime': 'timestamps'
58 },
59 'Metadata': {
60 'heightList': 'heights'
61 }
62 }
63
64 extras = {
65 'timeZone': 300
66 }
67
68 reader = project.addReadUnit(
69 name='HDFReader',
70 path='/path/to/files',
71 startDate='2019/01/01',
72 endDate='2019/01/31',
73 startTime='00:00:00',
74 endTime='23:59:59',
75 # description=json.dumps(desc),
76 # extras=json.dumps(extras),
77 )
158 78
159 79 """
160 Retorna 1 si el archivo de datos se encuentra dentro del rango de horas especificado.
161
162 Inputs:
163 filename : nombre completo del archivo de datos en formato Jicamarca (.r)
164 startDate : fecha inicial del rango seleccionado en formato datetime.date
165 endDate : fecha final del rango seleccionado en formato datetime.date
166 startTime : tiempo inicial del rango seleccionado en formato datetime.time
167 endTime : tiempo final del rango seleccionado en formato datetime.time
168
169 Return:
170 Boolean : Retorna True si el archivo de datos contiene datos en el rango de
171 fecha especificado, de lo contrario retorna False.
172
173 Excepciones:
174 Si el archivo no existe o no puede ser abierto
175 Si la cabecera no puede ser leida.
176
177 """
178
179 try:
180 fp = h5py.File(filename,'r')
181 grp1 = fp['Data']
182
183 except IOError:
184 traceback.print_exc()
185 raise IOError("The file %s can't be opened" %(filename))
186
187 #In case has utctime attribute
188 grp2 = grp1['utctime']
189 # thisUtcTime = grp2.value[0] - 5*3600 #To convert to local time
190 thisUtcTime = grp2.value[0]
191
192 fp.close()
193
194 if self.timezone == 'lt':
195 thisUtcTime -= 5*3600
196
197 thisDatetime = datetime.datetime.fromtimestamp(thisUtcTime[0] + 5*3600)
198 thisDate = thisDatetime.date()
199 thisTime = thisDatetime.time()
200
201 startUtcTime = (datetime.datetime.combine(thisDate,startTime)- datetime.datetime(1970, 1, 1)).total_seconds()
202 endUtcTime = (datetime.datetime.combine(thisDate,endTime)- datetime.datetime(1970, 1, 1)).total_seconds()
203
204 #General case
205 # o>>>>>>>>>>>>>><<<<<<<<<<<<<<o
206 #-----------o----------------------------o-----------
207 # startTime endTime
208
209 if endTime >= startTime:
210 thisUtcLog = numpy.logical_and(thisUtcTime > startUtcTime, thisUtcTime < endUtcTime)
211 if numpy.any(thisUtcLog): #If there is one block between the hours mentioned
212 return thisDatetime
213 return None
214
215 #If endTime < startTime then endTime belongs to the next day
216 #<<<<<<<<<<<o o>>>>>>>>>>>
217 #-----------o----------------------------o-----------
218 # endTime startTime
219
220 if (thisDate == startDate) and numpy.all(thisUtcTime < startUtcTime):
221 return None
222
223 if (thisDate == endDate) and numpy.all(thisUtcTime > endUtcTime):
224 return None
225
226 if numpy.all(thisUtcTime < startUtcTime) and numpy.all(thisUtcTime > endUtcTime):
227 return None
228
229 return thisDatetime
230
231 def __setNextFileOffline(self):
232
233 self.fileIndex += 1
234 idFile = self.fileIndex
235
236 if not(idFile < len(self.filenameList)):
237 raise schainpy.admin.SchainError("No more Files")
238 return 0
239
240 filename = self.filenameList[idFile]
241 filePointer = h5py.File(filename,'r')
242 self.filename = filename
243 self.fp = filePointer
244
245 print("Setting the file: %s"%self.filename)
246
247 self.__setBlockList()
248 self.__readData()
249 self.blockIndex = 0
250 return 1
251
252 def __setBlockList(self):
253 '''
254 Selects the data within the times defined
255
256 self.fp
257 self.startTime
258 self.endTime
259
260 self.blockList
261 self.blocksPerFile
262
263 '''
264 fp = self.fp
265 startTime = self.startTime
266 endTime = self.endTime
267
268 grp = fp['Data']
269 thisUtcTime = grp['utctime'].value.astype(numpy.float)[0]
270
271 #ERROOOOR
272 if self.timezone == 'lt':
273 thisUtcTime -= 5*3600
274
275 thisDatetime = datetime.datetime.fromtimestamp(thisUtcTime[0] + 5*3600)
276
277 thisDate = thisDatetime.date()
278 thisTime = thisDatetime.time()
279
280 startUtcTime = (datetime.datetime.combine(thisDate,startTime) - datetime.datetime(1970, 1, 1)).total_seconds()
281 endUtcTime = (datetime.datetime.combine(thisDate,endTime) - datetime.datetime(1970, 1, 1)).total_seconds()
282
283 ind = numpy.where(numpy.logical_and(thisUtcTime >= startUtcTime, thisUtcTime < endUtcTime))[0]
284
285 self.blockList = ind
286 self.blocksPerFile = len(ind)
287
288 return
289
290 def __readMetadata(self):
291 '''
292 Reads Metadata
293
294 self.pathMeta
295 self.listShapes
296 self.listMetaname
297 self.listMeta
298
299 '''
300
301 filename = self.filenameList[0]
302 fp = h5py.File(filename,'r')
303 gp = fp['Metadata']
304
305 listMetaname = []
306 listMetadata = []
307 for item in list(gp.items()):
308 name = item[0]
309
310 if name=='array dimensions':
311 table = gp[name][:]
312 listShapes = {}
313 for shapes in table:
314 listShapes[shapes[0]] = numpy.array([shapes[1],shapes[2],shapes[3],shapes[4],shapes[5]])
315 else:
316 data = gp[name].value
317 listMetaname.append(name)
318 listMetadata.append(data)
319
320 self.listShapes = listShapes
321 self.listMetaname = listMetaname
322 self.listMeta = listMetadata
323
324 fp.close()
325 return
326
327 def __readData(self):
328 grp = self.fp['Data']
329 listdataname = []
330 listdata = []
331
332 for item in list(grp.items()):
333 name = item[0]
334 listdataname.append(name)
335
336 array = self.__setDataArray(grp[name],self.listShapes[name])
337 listdata.append(array)
338
339 self.listDataname = listdataname
340 self.listData = listdata
341 return
342
343 def __setDataArray(self, dataset, shapes):
344
345 nDims = shapes[0]
346 nDim2 = shapes[1] #Dimension 0
347 nDim1 = shapes[2] #Dimension 1, number of Points or Parameters
348 nDim0 = shapes[3] #Dimension 2, number of samples or ranges
349 mode = shapes[4] #Mode of storing
350 blockList = self.blockList
351 blocksPerFile = self.blocksPerFile
352
353 #Depending on what mode the data was stored
354 if mode == 0: #Divided in channels
355 arrayData = dataset.value.astype(numpy.float)[0][blockList]
356 if mode == 1: #Divided in parameter
357 strds = 'table'
358 nDatas = nDim1
359 newShapes = (blocksPerFile,nDim2,nDim0)
360 elif mode==2: #Concatenated in a table
361 strds = 'table0'
362 arrayData = dataset[strds].value
363 #Selecting part of the dataset
364 utctime = arrayData[:,0]
365 u, indices = numpy.unique(utctime, return_index=True)
366
367 if blockList.size != indices.size:
368 indMin = indices[blockList[0]]
369 if blockList[1] + 1 >= indices.size:
370 arrayData = arrayData[indMin:,:]
371 else:
372 indMax = indices[blockList[1] + 1]
373 arrayData = arrayData[indMin:indMax,:]
374 return arrayData
375
376 # One dimension
377 if nDims == 0:
378 arrayData = dataset.value.astype(numpy.float)[0][blockList]
379
380 # Two dimensions
381 elif nDims == 2:
382 arrayData = numpy.zeros((blocksPerFile,nDim1,nDim0))
383 newShapes = (blocksPerFile,nDim0)
384 nDatas = nDim1
385
386 for i in range(nDatas):
387 data = dataset[strds + str(i)].value
388 arrayData[:,i,:] = data[blockList,:]
389
390 # Three dimensions
391 else:
392 arrayData = numpy.zeros((blocksPerFile,nDim2,nDim1,nDim0))
393 for i in range(nDatas):
394
395 data = dataset[strds + str(i)].value
396
397 for b in range(blockList.size):
398 arrayData[b,:,i,:] = data[:,:,blockList[b]]
399
400 return arrayData
401
402 def __setDataOut(self):
403 listMeta = self.listMeta
404 listMetaname = self.listMetaname
405 listDataname = self.listDataname
406 listData = self.listData
407 listShapes = self.listShapes
408
409 blockIndex = self.blockIndex
410 # blockList = self.blockList
411
412 for i in range(len(listMeta)):
413 setattr(self.dataOut,listMetaname[i],listMeta[i])
414
415 for j in range(len(listData)):
416 nShapes = listShapes[listDataname[j]][0]
417 mode = listShapes[listDataname[j]][4]
418 if nShapes == 1:
419 setattr(self.dataOut,listDataname[j],listData[j][blockIndex])
420 elif nShapes > 1:
421 setattr(self.dataOut,listDataname[j],listData[j][blockIndex,:])
422 elif mode==0:
423 setattr(self.dataOut,listDataname[j],listData[j][blockIndex])
424 #Mode Meteors
425 elif mode ==2:
426 selectedData = self.__selectDataMode2(listData[j], blockIndex)
427 setattr(self.dataOut, listDataname[j], selectedData)
428 return
429
430 def __selectDataMode2(self, data, blockIndex):
431 utctime = data[:,0]
432 aux, indices = numpy.unique(utctime, return_inverse=True)
433 selInd = numpy.where(indices == blockIndex)[0]
434 selData = data[selInd,:]
435
436 return selData
437
438 def getData(self):
439
440 if self.blockIndex==self.blocksPerFile:
441 if not( self.__setNextFileOffline() ):
442 self.dataOut.flagNoData = True
443 return 0
444
445 self.__setDataOut()
446 self.dataOut.flagNoData = False
447
448 self.blockIndex += 1
449
450 return
451
452 def run(self, **kwargs):
453
454 if not(self.isConfig):
455 self.setup(**kwargs)
456 self.isConfig = True
457
458 self.getData()
459
460 return
461
462 @MPDecorator
463 class ParamWriter(Operation):
464 '''
465 HDF5 Writer, stores parameters data in HDF5 format files
466
467 path: path where the files will be stored
468 blocksPerFile: number of blocks that will be saved in per HDF5 format file
469 mode: selects the data stacking mode: '0' channels, '1' parameters, '3' table (for meteors)
470 metadataList: list of attributes that will be stored as metadata
471 dataList: list of attributes that will be stores as data
472 '''
473
474 ext = ".hdf5"
475 optchar = "D"
476 metaoptchar = "M"
477 metaFile = None
478 filename = None
479 path = None
480 setFile = None
481 fp = None
482 grp = None
483 ds = None
484 firsttime = True
485 #Configurations
486 blocksPerFile = None
487 blockIndex = None
488 dataOut = None
489 #Data Arrays
490 dataList = None
491 metadataList = None
492 dsList = None #List of dictionaries with dataset properties
493 tableDim = None
494 dtype = [('arrayName', 'S20'),('nDimensions', 'i'), ('dim2', 'i'), ('dim1', 'i'),('dim0', 'i'),('mode', 'b')]
495 currentDay = None
496 lastTime = None
497 setType = None
498
499 def __init__(self):
500
501 Operation.__init__(self)
502 return
503
504 def setup(self, dataOut, path=None, blocksPerFile=10, metadataList=None, dataList=None, mode=None, setType=None):
505 self.path = path
506 self.blocksPerFile = blocksPerFile
507 self.metadataList = metadataList
508 self.dataList = dataList
509 self.dataOut = dataOut
510 self.mode = mode
511 if self.mode is not None:
512 self.mode = numpy.zeros(len(self.dataList)) + mode
513 else:
514 self.mode = numpy.ones(len(self.dataList))
515
516 self.setType = setType
517
518 arrayDim = numpy.zeros((len(self.dataList),5))
519
520 #Table dimensions
521 dtype0 = self.dtype
522 tableList = []
523
524 #Dictionary and list of tables
525 dsList = []
526
527 for i in range(len(self.dataList)):
528 dsDict = {}
529 dataAux = getattr(self.dataOut, self.dataList[i])
530 dsDict['variable'] = self.dataList[i]
531 #--------------------- Conditionals ------------------------
532 #There is no data
533
534 if dataAux is None:
535
536 return 0
537
538 if isinstance(dataAux, (int, float, numpy.integer, numpy.float)):
539 dsDict['mode'] = 0
540 dsDict['nDim'] = 0
541 arrayDim[i,0] = 0
542 dsList.append(dsDict)
543
544 #Mode 2: meteors
545 elif self.mode[i] == 2:
546 dsDict['dsName'] = 'table0'
547 dsDict['mode'] = 2 # Mode meteors
548 dsDict['shape'] = dataAux.shape[-1]
549 dsDict['nDim'] = 0
550 dsDict['dsNumber'] = 1
551 arrayDim[i,3] = dataAux.shape[-1]
552 arrayDim[i,4] = self.mode[i] #Mode the data was stored
553 dsList.append(dsDict)
554
555 #Mode 1
556 else:
557 arrayDim0 = dataAux.shape #Data dimensions
558 arrayDim[i,0] = len(arrayDim0) #Number of array dimensions
559 arrayDim[i,4] = self.mode[i] #Mode the data was stored
560 strtable = 'table'
561 dsDict['mode'] = 1 # Mode parameters
562
563 # Three-dimension arrays
564 if len(arrayDim0) == 3:
565 arrayDim[i,1:-1] = numpy.array(arrayDim0)
566 nTables = int(arrayDim[i,2])
567 dsDict['dsNumber'] = nTables
568 dsDict['shape'] = arrayDim[i,2:4]
569 dsDict['nDim'] = 3
570
571 for j in range(nTables):
572 dsDict = dsDict.copy()
573 dsDict['dsName'] = strtable + str(j)
574 dsList.append(dsDict)
575
576 # Two-dimension arrays
577 elif len(arrayDim0) == 2:
578 arrayDim[i,2:-1] = numpy.array(arrayDim0)
579 nTables = int(arrayDim[i,2])
580 dsDict['dsNumber'] = nTables
581 dsDict['shape'] = arrayDim[i,3]
582 dsDict['nDim'] = 2
583
584 for j in range(nTables):
585 dsDict = dsDict.copy()
586 dsDict['dsName'] = strtable + str(j)
587 dsList.append(dsDict)
588
589 # One-dimension arrays
590 elif len(arrayDim0) == 1:
591 arrayDim[i,3] = arrayDim0[0]
592 dsDict['shape'] = arrayDim0[0]
593 dsDict['dsNumber'] = 1
594 dsDict['dsName'] = strtable + str(0)
595 dsDict['nDim'] = 1
596 dsList.append(dsDict)
597 80
598 table = numpy.array((self.dataList[i],) + tuple(arrayDim[i,:]),dtype = dtype0)
599 tableList.append(table)
600
601 self.dsList = dsList
602 self.tableDim = numpy.array(tableList, dtype = dtype0)
603 self.blockIndex = 0
604 timeTuple = time.localtime(dataOut.utctime)
605 self.currentDay = timeTuple.tm_yday
606
607 def putMetadata(self):
608
609 fp = self.createMetadataFile()
610 self.writeMetadata(fp)
611 fp.close()
612 return
613
614 def createMetadataFile(self):
615 ext = self.ext
616 path = self.path
617 setFile = self.setFile
618
619 timeTuple = time.localtime(self.dataOut.utctime)
620
621 subfolder = ''
622 fullpath = os.path.join( path, subfolder )
623
624 if not( os.path.exists(fullpath) ):
625 os.mkdir(fullpath)
626 setFile = -1 #inicializo mi contador de seteo
627
628 subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
629 fullpath = os.path.join( path, subfolder )
630
631 if not( os.path.exists(fullpath) ):
632 os.mkdir(fullpath)
633 setFile = -1 #inicializo mi contador de seteo
634
635 else:
636 filesList = os.listdir( fullpath )
637 filesList = sorted( filesList, key=str.lower )
638 if len( filesList ) > 0:
639 filesList = [k for k in filesList if k.startswith(self.metaoptchar)]
640 filen = filesList[-1]
641 # el filename debera tener el siguiente formato
642 # 0 1234 567 89A BCDE (hex)
643 # x YYYY DDD SSS .ext
644 if isNumber( filen[8:11] ):
645 setFile = int( filen[8:11] ) #inicializo mi contador de seteo al seteo del ultimo file
646 else:
647 setFile = -1
648 else:
649 setFile = -1 #inicializo mi contador de seteo
650
651 if self.setType is None:
652 setFile += 1
653 file = '%s%4.4d%3.3d%03d%s' % (self.metaoptchar,
654 timeTuple.tm_year,
655 timeTuple.tm_yday,
656 setFile,
657 ext )
658 else:
659 setFile = timeTuple.tm_hour*60+timeTuple.tm_min
660 file = '%s%4.4d%3.3d%04d%s' % (self.metaoptchar,
661 timeTuple.tm_year,
662 timeTuple.tm_yday,
663 setFile,
664 ext )
665
666 filename = os.path.join( path, subfolder, file )
667 self.metaFile = file
668 #Setting HDF5 File
669 fp = h5py.File(filename,'w')
670
671 return fp
672
673 def writeMetadata(self, fp):
674
675 grp = fp.create_group("Metadata")
676 grp.create_dataset('array dimensions', data = self.tableDim, dtype = self.dtype)
677
678 for i in range(len(self.metadataList)):
679 grp.create_dataset(self.metadataList[i], data=getattr(self.dataOut, self.metadataList[i]))
680 return
681
682 def timeFlag(self):
683 currentTime = self.dataOut.utctime
684
685 if self.lastTime is None:
686 self.lastTime = currentTime
687
688 #Day
689 timeTuple = time.localtime(currentTime)
690 dataDay = timeTuple.tm_yday
691
692 #Time
693 timeDiff = currentTime - self.lastTime
694
695 #Si el dia es diferente o si la diferencia entre un dato y otro supera la hora
696 if dataDay != self.currentDay:
697 self.currentDay = dataDay
698 return True
699 elif timeDiff > 3*60*60:
700 self.lastTime = currentTime
701 return True
702 else:
703 self.lastTime = currentTime
704 return False
705
706 def setNextFile(self):
707
708 ext = self.ext
709 path = self.path
710 setFile = self.setFile
711 mode = self.mode
712
713 timeTuple = time.localtime(self.dataOut.utctime)
714 subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
715
716 fullpath = os.path.join( path, subfolder )
717
718 if os.path.exists(fullpath):
719 filesList = os.listdir( fullpath )
720 filesList = [k for k in filesList if 'M' in k]
721 if len( filesList ) > 0:
722 filesList = sorted( filesList, key=str.lower )
723 filen = filesList[-1]
724 # el filename debera tener el siguiente formato
725 # 0 1234 567 89A BCDE (hex)
726 # x YYYY DDD SSS .ext
727 if isNumber( filen[8:11] ):
728 setFile = int( filen[8:11] ) #inicializo mi contador de seteo al seteo del ultimo file
729 else:
730 setFile = -1
731 else:
732 setFile = -1 #inicializo mi contador de seteo
733 else:
734 os.makedirs(fullpath)
735 setFile = -1 #inicializo mi contador de seteo
736
737 if self.setType is None:
738 setFile += 1
739 file = '%s%4.4d%3.3d%03d%s' % (self.optchar,
740 timeTuple.tm_year,
741 timeTuple.tm_yday,
742 setFile,
743 ext )
744 else:
745 setFile = timeTuple.tm_hour*60+timeTuple.tm_min
746 file = '%s%4.4d%3.3d%04d%s' % (self.optchar,
747 timeTuple.tm_year,
748 timeTuple.tm_yday,
749 setFile,
750 ext )
751
752 filename = os.path.join( path, subfolder, file )
753
754 #Setting HDF5 File
755 fp = h5py.File(filename,'w')
756 #write metadata
757 self.writeMetadata(fp)
758 #Write data
759 grp = fp.create_group("Data")
760 ds = []
761 data = []
762 dsList = self.dsList
763 i = 0
764 while i < len(dsList):
765 dsInfo = dsList[i]
766 #One-dimension data
767 if dsInfo['mode'] == 0:
768 ds0 = grp.create_dataset(dsInfo['variable'], (1,1), maxshape=(1,self.blocksPerFile) , chunks = True, dtype=numpy.float64)
769 ds.append(ds0)
770 data.append([])
771 i += 1
772 continue
773
774 elif dsInfo['mode'] == 2:
775 grp0 = grp.create_group(dsInfo['variable'])
776 ds0 = grp0.create_dataset(dsInfo['dsName'], (1,dsInfo['shape']), data = numpy.zeros((1,dsInfo['shape'])) , maxshape=(None,dsInfo['shape']), chunks=True)
777 ds.append(ds0)
778 data.append([])
779 i += 1
780 continue
781
782 elif dsInfo['mode'] == 1:
783 grp0 = grp.create_group(dsInfo['variable'])
784
785 for j in range(dsInfo['dsNumber']):
786 dsInfo = dsList[i]
787 tableName = dsInfo['dsName']
788
789
790 if dsInfo['nDim'] == 3:
791 shape = dsInfo['shape'].astype(int)
792 ds0 = grp0.create_dataset(tableName, (shape[0],shape[1],1) , data = numpy.zeros((shape[0],shape[1],1)), maxshape = (None,shape[1],None), chunks=True)
793 else:
794 shape = int(dsInfo['shape'])
795 ds0 = grp0.create_dataset(tableName, (1,shape), data = numpy.zeros((1,shape)) , maxshape=(None,shape), chunks=True)
796
797 ds.append(ds0)
798 data.append([])
799 i += 1
800
801 fp.flush()
802 fp.close()
803
804 log.log('creating file: {}'.format(filename), 'Writing')
805 self.filename = filename
806 self.ds = ds
807 self.data = data
808 self.firsttime = True
809 self.blockIndex = 0
810 return
811
812 def putData(self):
813
814 if self.blockIndex == self.blocksPerFile or self.timeFlag():
815 self.setNextFile()
816
817 self.readBlock()
818 self.setBlock() #Prepare data to be written
819 self.writeBlock() #Write data
820
821 return
822
823 def readBlock(self):
824
825 '''
826 data Array configured
827
828
829 self.data
830 '''
831 dsList = self.dsList
832 ds = self.ds
833 #Setting HDF5 File
834 fp = h5py.File(self.filename,'r+')
835 grp = fp["Data"]
836 ind = 0
837
838 while ind < len(dsList):
839 dsInfo = dsList[ind]
840
841 if dsInfo['mode'] == 0:
842 ds0 = grp[dsInfo['variable']]
843 ds[ind] = ds0
844 ind += 1
845 else:
846
847 grp0 = grp[dsInfo['variable']]
848
849 for j in range(dsInfo['dsNumber']):
850 dsInfo = dsList[ind]
851 ds0 = grp0[dsInfo['dsName']]
852 ds[ind] = ds0
853 ind += 1
854
855 self.fp = fp
856 self.grp = grp
857 self.ds = ds
858
859 return
860
861 def setBlock(self):
862 '''
863 data Array configured
864
865
866 self.data
867 '''
868 #Creating Arrays
869 dsList = self.dsList
870 data = self.data
871 ind = 0
872
873 while ind < len(dsList):
874 dsInfo = dsList[ind]
875 dataAux = getattr(self.dataOut, dsInfo['variable'])
876
877 mode = dsInfo['mode']
878 nDim = dsInfo['nDim']
879
880 if mode == 0 or mode == 2 or nDim == 1:
881 data[ind] = dataAux
882 ind += 1
883 # elif nDim == 1:
884 # data[ind] = numpy.reshape(dataAux,(numpy.size(dataAux),1))
885 # ind += 1
886 elif nDim == 2:
887 for j in range(dsInfo['dsNumber']):
888 data[ind] = dataAux[j,:]
889 ind += 1
890 elif nDim == 3:
891 for j in range(dsInfo['dsNumber']):
892 data[ind] = dataAux[:,j,:]
893 ind += 1
894
895 self.data = data
896 return
897
898 def writeBlock(self):
899 '''
900 Saves the block in the HDF5 file
901 '''
902 dsList = self.dsList
903
904 for i in range(len(self.ds)):
905 dsInfo = dsList[i]
906 nDim = dsInfo['nDim']
907 mode = dsInfo['mode']
908
909 # First time
910 if self.firsttime:
911 if type(self.data[i]) == numpy.ndarray:
912
913 if nDim == 3:
914 self.data[i] = self.data[i].reshape((self.data[i].shape[0],self.data[i].shape[1],1))
915 self.ds[i].resize(self.data[i].shape)
916 if mode == 2:
917 self.ds[i].resize(self.data[i].shape)
918 self.ds[i][:] = self.data[i]
919 else:
920
921 # From second time
922 # Meteors!
923 if mode == 2:
924 dataShape = self.data[i].shape
925 dsShape = self.ds[i].shape
926 self.ds[i].resize((self.ds[i].shape[0] + dataShape[0],self.ds[i].shape[1]))
927 self.ds[i][dsShape[0]:,:] = self.data[i]
928 # No dimension
929 elif mode == 0:
930 self.ds[i].resize((self.ds[i].shape[0], self.ds[i].shape[1] + 1))
931 self.ds[i][0,-1] = self.data[i]
932 # One dimension
933 elif nDim == 1:
934 self.ds[i].resize((self.ds[i].shape[0] + 1, self.ds[i].shape[1]))
935 self.ds[i][-1,:] = self.data[i]
936 # Two dimension
937 elif nDim == 2:
938 self.ds[i].resize((self.ds[i].shape[0] + 1,self.ds[i].shape[1]))
939 self.ds[i][self.blockIndex,:] = self.data[i]
940 # Three dimensions
941 elif nDim == 3:
942 self.ds[i].resize((self.ds[i].shape[0],self.ds[i].shape[1],self.ds[i].shape[2]+1))
943 self.ds[i][:,:,-1] = self.data[i]
944
945 self.firsttime = False
946 self.blockIndex += 1
947
948 #Close to save changes
949 self.fp.flush()
950 self.fp.close()
951 return
952
953 def run(self, dataOut, path, blocksPerFile=10, metadataList=None, dataList=None, mode=None, setType=None):
954
955 self.dataOut = dataOut
956 if not(self.isConfig):
957 self.setup(dataOut, path=path, blocksPerFile=blocksPerFile,
958 metadataList=metadataList, dataList=dataList, mode=mode,
959 setType=setType)
960
961 self.isConfig = True
962 self.setNextFile()
963
964 self.putData()
965 return
966
967
968
969 class ParameterReader(Reader, ProcessingUnit):
970 '''
971 Reads HDF5 format files
972 '''
81 __attrs__ = ['path', 'startDate', 'endDate', 'startTime', 'endTime', 'description', 'extras']
973 82
974 83 def __init__(self):
975 84 ProcessingUnit.__init__(self)
976 85 self.dataOut = Parameters()
977 86 self.ext = ".hdf5"
978 87 self.optchar = "D"
979 self.timezone = "lt"
980 self.listMetaname = []
981 self.listMeta = []
982 self.listDataname = []
983 self.listData = []
984 self.listShapes = []
88 self.meta = {}
89 self.data = {}
985 90 self.open_file = h5py.File
986 91 self.open_mode = 'r'
987 self.metadata = False
92 self.description = {}
93 self.extras = {}
988 94 self.filefmt = "*%Y%j***"
989 95 self.folderfmt = "*%Y%j"
990 96
991 97 def setup(self, **kwargs):
992 98
993 99 self.set_kwargs(**kwargs)
994 100 if not self.ext.startswith('.'):
995 101 self.ext = '.{}'.format(self.ext)
996 102
997 103 if self.online:
998 104 log.log("Searching files in online mode...", self.name)
999 105
1000 106 for nTries in range(self.nTries):
1001 107 fullpath = self.searchFilesOnLine(self.path, self.startDate,
1002 108 self.endDate, self.expLabel, self.ext, self.walk,
1003 109 self.filefmt, self.folderfmt)
1004
1005 110 try:
1006 111 fullpath = next(fullpath)
1007 112 except:
1008 113 fullpath = None
1009 114
1010 115 if fullpath:
1011 116 break
1012 117
1013 118 log.warning(
1014 119 'Waiting {} sec for a valid file in {}: try {} ...'.format(
1015 120 self.delay, self.path, nTries + 1),
1016 121 self.name)
1017 122 time.sleep(self.delay)
1018 123
1019 124 if not(fullpath):
1020 125 raise schainpy.admin.SchainError(
1021 126 'There isn\'t any valid file in {}'.format(self.path))
1022 127
1023 128 pathname, filename = os.path.split(fullpath)
1024 129 self.year = int(filename[1:5])
1025 130 self.doy = int(filename[5:8])
1026 131 self.set = int(filename[8:11]) - 1
1027 132 else:
1028 133 log.log("Searching files in {}".format(self.path), self.name)
1029 134 self.filenameList = self.searchFilesOffLine(self.path, self.startDate,
1030 135 self.endDate, self.expLabel, self.ext, self.walk, self.filefmt, self.folderfmt)
1031 136
1032 137 self.setNextFile()
1033 138
1034 139 return
1035 140
1036 141 def readFirstHeader(self):
1037 142 '''Read metadata and data'''
1038 143
1039 144 self.__readMetadata()
1040 145 self.__readData()
1041 146 self.__setBlockList()
147
148 if 'type' in self.meta:
149 self.dataOut = eval(self.meta['type'])()
150
151 for attr in self.meta:
152 setattr(self.dataOut, attr, self.meta[attr])
153
1042 154 self.blockIndex = 0
1043 155
1044 156 return
1045 157
1046 158 def __setBlockList(self):
1047 159 '''
1048 160 Selects the data within the times defined
1049 161
1050 162 self.fp
1051 163 self.startTime
1052 164 self.endTime
1053 165 self.blockList
1054 166 self.blocksPerFile
1055 167
1056 168 '''
1057 169
1058 170 startTime = self.startTime
1059 171 endTime = self.endTime
1060 172
1061 index = self.listDataname.index('utctime')
1062 thisUtcTime = self.listData[index]
173 thisUtcTime = self.data['utctime']
1063 174 self.interval = numpy.min(thisUtcTime[1:] - thisUtcTime[:-1])
1064 175
1065 if self.timezone == 'lt':
1066 thisUtcTime -= 5*3600
1067
1068 thisDatetime = datetime.datetime.fromtimestamp(thisUtcTime[0] + 5*3600)
176 thisDatetime = datetime.datetime.utcfromtimestamp(thisUtcTime[0])
1069 177
1070 178 thisDate = thisDatetime.date()
1071 179 thisTime = thisDatetime.time()
1072 180
1073 181 startUtcTime = (datetime.datetime.combine(thisDate,startTime) - datetime.datetime(1970, 1, 1)).total_seconds()
1074 182 endUtcTime = (datetime.datetime.combine(thisDate,endTime) - datetime.datetime(1970, 1, 1)).total_seconds()
1075 183
1076 184 ind = numpy.where(numpy.logical_and(thisUtcTime >= startUtcTime, thisUtcTime < endUtcTime))[0]
1077 185
1078 186 self.blockList = ind
1079 187 self.blocksPerFile = len(ind)
1080 188 return
1081 189
1082 190 def __readMetadata(self):
1083 191 '''
1084 192 Reads Metadata
1085 193 '''
1086 194
1087 listMetaname = []
1088 listMetadata = []
1089 if 'Metadata' in self.fp:
1090 gp = self.fp['Metadata']
1091 for item in list(gp.items()):
1092 name = item[0]
1093
1094 if name=='variables':
1095 table = gp[name][:]
1096 listShapes = {}
1097 for shapes in table:
1098 listShapes[shapes[0].decode()] = numpy.array([shapes[1]])
1099 else:
1100 data = gp[name].value
1101 listMetaname.append(name)
1102 listMetadata.append(data)
1103 elif self.metadata:
1104 metadata = json.loads(self.metadata)
1105 listShapes = {}
1106 for tup in metadata:
1107 name, values, dim = tup
1108 if dim == -1:
1109 listMetaname.append(name)
1110 listMetadata.append(self.fp[values].value)
1111 else:
1112 listShapes[name] = numpy.array([dim])
195 meta = {}
196
197 if self.description:
198 for key, value in self.description['Metadata'].items():
199 meta[key] = self.fp[value].value
1113 200 else:
1114 raise IOError('Missing Metadata group in file or metadata info')
201 grp = self.fp['Metadata']
202 for name in grp:
203 meta[name] = grp[name].value
1115 204
1116 self.listShapes = listShapes
1117 self.listMetaname = listMetaname
1118 self.listMeta = listMetadata
205 if self.extras:
206 for key, value in self.extras.items():
207 meta[key] = value
208 self.meta = meta
1119 209
1120 210 return
1121 211
1122 212 def __readData(self):
1123 213
1124 listdataname = []
1125 listdata = []
214 data = {}
1126 215
1127 if 'Data' in self.fp:
216 if self.description:
217 for key, value in self.description['Data'].items():
218 if isinstance(value, str):
219 if isinstance(self.fp[value], h5py.Dataset):
220 data[key] = self.fp[value].value
221 elif isinstance(self.fp[value], h5py.Group):
222 array = []
223 for ch in self.fp[value]:
224 array.append(self.fp[value][ch].value)
225 data[key] = numpy.array(array)
226 elif isinstance(value, list):
227 array = []
228 for ch in value:
229 array.append(self.fp[ch].value)
230 data[key] = numpy.array(array)
231 else:
1128 232 grp = self.fp['Data']
1129 for item in list(grp.items()):
1130 name = item[0]
1131 listdataname.append(name)
1132 dim = self.listShapes[name][0]
1133 if dim == 0:
233 for name in grp:
234 if isinstance(grp[name], h5py.Dataset):
1134 235 array = grp[name].value
1135 else:
236 elif isinstance(grp[name], h5py.Group):
1136 237 array = []
1137 for i in range(dim):
1138 array.append(grp[name]['table{:02d}'.format(i)].value)
238 for ch in grp[name]:
239 array.append(grp[name][ch].value)
1139 240 array = numpy.array(array)
1140
1141 listdata.append(array)
1142 elif self.metadata:
1143 metadata = json.loads(self.metadata)
1144 for tup in metadata:
1145 name, values, dim = tup
1146 listdataname.append(name)
1147 if dim == -1:
1148 continue
1149 elif dim == 0:
1150 array = self.fp[values].value
1151 241 else:
1152 array = []
1153 for var in values:
1154 array.append(self.fp[var].value)
1155 array = numpy.array(array)
1156 listdata.append(array)
242 log.warning('Unknown type: {}'.format(name))
243
244 if name in self.description:
245 key = self.description[name]
1157 246 else:
1158 raise IOError('Missing Data group in file or metadata info')
247 key = name
248 data[key] = array
1159 249
1160 self.listDataname = listdataname
1161 self.listData = listdata
250 self.data = data
1162 251 return
1163 252
1164 253 def getData(self):
1165 254
1166 for i in range(len(self.listMeta)):
1167 setattr(self.dataOut, self.listMetaname[i], self.listMeta[i])
1168
1169 for j in range(len(self.listData)):
1170 dim = self.listShapes[self.listDataname[j]][0]
1171 if dim == 0:
1172 setattr(self.dataOut, self.listDataname[j], self.listData[j][self.blockIndex])
255 for attr in self.data:
256 if self.data[attr].ndim == 1:
257 setattr(self.dataOut, attr, self.data[attr][self.blockIndex])
1173 258 else:
1174 setattr(self.dataOut, self.listDataname[j], self.listData[j][:,self.blockIndex])
259 setattr(self.dataOut, attr, self.data[attr][:, self.blockIndex])
1175 260
1176 self.dataOut.paramInterval = self.interval
1177 261 self.dataOut.flagNoData = False
1178 262 self.blockIndex += 1
1179 263
264 log.log("Block No. {}/{} -> {}".format(
265 self.blockIndex,
266 self.blocksPerFile,
267 self.dataOut.datatime.ctime()), self.name)
268
1180 269 return
1181 270
1182 271 def run(self, **kwargs):
1183 272
1184 273 if not(self.isConfig):
1185 274 self.setup(**kwargs)
1186 275 self.isConfig = True
1187 276
1188 277 if self.blockIndex == self.blocksPerFile:
1189 278 self.setNextFile()
1190 279
1191 280 self.getData()
1192 281
1193 282 return
1194 283
1195 284 @MPDecorator
1196 class ParameterWriter(Operation):
1197 '''
1198 HDF5 Writer, stores parameters data in HDF5 format files
1199
1200 path: path where the files will be stored
1201 blocksPerFile: number of blocks that will be saved in per HDF5 format file
1202 mode: selects the data stacking mode: '0' channels, '1' parameters, '3' table (for meteors)
1203 metadataList: list of attributes that will be stored as metadata
1204 dataList: list of attributes that will be stores as data
1205 '''
285 class HDFWriter(Operation):
286 """Operation to write HDF5 files.
287
288 The HDF5 file contains by default two groups Data and Metadata where
289 you can save any `dataOut` attribute specified by `dataList` and `metadataList`
290 parameters, data attributes are normaly time dependent where the metadata
291 are not.
292 It is possible to customize the structure of the HDF5 file with the
293 optional description parameter see the examples.
294
295 Parameters:
296 -----------
297 path : str
298 Path where files will be saved.
299 blocksPerFile : int
300 Number of blocks per file
301 metadataList : list
302 List of the dataOut attributes that will be saved as metadata
303 dataList : int
304 List of the dataOut attributes that will be saved as data
305 setType : bool
306 If True the name of the files corresponds to the timestamp of the data
307 description : dict, optional
308 Dictionary with the desired description of the HDF5 file
309
310 Examples
311 --------
312
313 desc = {
314 'data_output': {'winds': ['z', 'w', 'v']},
315 'utctime': 'timestamps',
316 'heightList': 'heights'
317 }
318 desc = {
319 'data_output': ['z', 'w', 'v'],
320 'utctime': 'timestamps',
321 'heightList': 'heights'
322 }
323 desc = {
324 'Data': {
325 'data_output': 'winds',
326 'utctime': 'timestamps'
327 },
328 'Metadata': {
329 'heightList': 'heights'
330 }
331 }
332
333 writer = proc_unit.addOperation(name='HDFWriter')
334 writer.addParameter(name='path', value='/path/to/file')
335 writer.addParameter(name='blocksPerFile', value='32')
336 writer.addParameter(name='metadataList', value='heightList,timeZone')
337 writer.addParameter(name='dataList',value='data_output,utctime')
338 # writer.addParameter(name='description',value=json.dumps(desc))
1206 339
340 """
1207 341
1208 342 ext = ".hdf5"
1209 343 optchar = "D"
1210 metaoptchar = "M"
1211 metaFile = None
1212 344 filename = None
1213 345 path = None
1214 346 setFile = None
1215 347 fp = None
1216 grp = None
1217 ds = None
1218 348 firsttime = True
1219 349 #Configurations
1220 350 blocksPerFile = None
1221 351 blockIndex = None
1222 352 dataOut = None
1223 353 #Data Arrays
1224 354 dataList = None
1225 355 metadataList = None
1226 dsList = None #List of dictionaries with dataset properties
1227 tableDim = None
1228 dtype = [('name', 'S20'),('nDim', 'i')]
1229 356 currentDay = None
1230 357 lastTime = None
1231 358
1232 359 def __init__(self):
1233 360
1234 361 Operation.__init__(self)
1235 362 return
1236 363
1237 def setup(self, path=None, blocksPerFile=10, metadataList=None, dataList=None, setType=None):
364 def setup(self, path=None, blocksPerFile=10, metadataList=None, dataList=None, setType=None, description=None):
1238 365 self.path = path
1239 366 self.blocksPerFile = blocksPerFile
1240 367 self.metadataList = metadataList
1241 368 self.dataList = dataList
1242 369 self.setType = setType
370 self.description = description
371
372 for s in ['type', 'timeZone', 'useLocalTime']:
373 if s not in self.metadataList:
374 self.metadataList.append(s)
1243 375
1244 376 tableList = []
1245 377 dsList = []
1246 378
1247 379 for i in range(len(self.dataList)):
1248 380 dsDict = {}
1249 381 dataAux = getattr(self.dataOut, self.dataList[i])
1250 382 dsDict['variable'] = self.dataList[i]
1251 383
1252 384 if dataAux is None:
1253 385 continue
1254 386 elif isinstance(dataAux, (int, float, numpy.integer, numpy.float)):
1255 387 dsDict['nDim'] = 0
1256 388 else:
1257 389 dsDict['nDim'] = len(dataAux.shape)
1258 390 dsDict['shape'] = dataAux.shape
1259 391 dsDict['dsNumber'] = dataAux.shape[0]
392 dsDict['dtype'] = dataAux.dtype
1260 393
1261 394 dsList.append(dsDict)
1262 tableList.append((self.dataList[i], dsDict['nDim']))
1263 395
1264 396 self.dsList = dsList
1265 self.tableDim = numpy.array(tableList, dtype=self.dtype)
1266 397 self.currentDay = self.dataOut.datatime.date()
1267 398
1268 399 def timeFlag(self):
1269 400 currentTime = self.dataOut.utctime
1270 401 timeTuple = time.localtime(currentTime)
1271 402 dataDay = timeTuple.tm_yday
1272 403
1273 404 if self.lastTime is None:
1274 405 self.lastTime = currentTime
1275 406 self.currentDay = dataDay
1276 407 return False
1277 408
1278 409 timeDiff = currentTime - self.lastTime
1279 410
1280 411 #Si el dia es diferente o si la diferencia entre un dato y otro supera la hora
1281 412 if dataDay != self.currentDay:
1282 413 self.currentDay = dataDay
1283 414 return True
1284 415 elif timeDiff > 3*60*60:
1285 416 self.lastTime = currentTime
1286 417 return True
1287 418 else:
1288 419 self.lastTime = currentTime
1289 420 return False
1290 421
1291 def run(self, dataOut, path, blocksPerFile=10, metadataList=None, dataList=None, setType=None):
422 def run(self, dataOut, path, blocksPerFile=10, metadataList=[],
423 dataList=[], setType=None, description={}):
1292 424
1293 425 self.dataOut = dataOut
1294 426 if not(self.isConfig):
1295 427 self.setup(path=path, blocksPerFile=blocksPerFile,
1296 428 metadataList=metadataList, dataList=dataList,
1297 setType=setType)
429 setType=setType, description=description)
1298 430
1299 431 self.isConfig = True
1300 432 self.setNextFile()
1301 433
1302 434 self.putData()
1303 435 return
1304 436
1305 437 def setNextFile(self):
1306 438
1307 439 ext = self.ext
1308 440 path = self.path
1309 441 setFile = self.setFile
1310 442
1311 443 timeTuple = time.localtime(self.dataOut.utctime)
1312 444 subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
1313 445 fullpath = os.path.join(path, subfolder)
1314 446
1315 447 if os.path.exists(fullpath):
1316 448 filesList = os.listdir(fullpath)
1317 449 filesList = [k for k in filesList if k.startswith(self.optchar)]
1318 450 if len( filesList ) > 0:
1319 451 filesList = sorted(filesList, key=str.lower)
1320 452 filen = filesList[-1]
1321 453 # el filename debera tener el siguiente formato
1322 454 # 0 1234 567 89A BCDE (hex)
1323 455 # x YYYY DDD SSS .ext
1324 456 if isNumber(filen[8:11]):
1325 457 setFile = int(filen[8:11]) #inicializo mi contador de seteo al seteo del ultimo file
1326 458 else:
1327 459 setFile = -1
1328 460 else:
1329 461 setFile = -1 #inicializo mi contador de seteo
1330 462 else:
1331 463 os.makedirs(fullpath)
1332 464 setFile = -1 #inicializo mi contador de seteo
1333 465
1334 466 if self.setType is None:
1335 467 setFile += 1
1336 468 file = '%s%4.4d%3.3d%03d%s' % (self.optchar,
1337 469 timeTuple.tm_year,
1338 470 timeTuple.tm_yday,
1339 471 setFile,
1340 472 ext )
1341 473 else:
1342 474 setFile = timeTuple.tm_hour*60+timeTuple.tm_min
1343 475 file = '%s%4.4d%3.3d%04d%s' % (self.optchar,
1344 476 timeTuple.tm_year,
1345 477 timeTuple.tm_yday,
1346 478 setFile,
1347 479 ext )
1348 480
1349 481 self.filename = os.path.join( path, subfolder, file )
1350 482
1351 483 #Setting HDF5 File
1352 484 self.fp = h5py.File(self.filename, 'w')
1353 485 #write metadata
1354 486 self.writeMetadata(self.fp)
1355 487 #Write data
1356 488 self.writeData(self.fp)
1357 489
490 def getLabel(self, name, x=None):
491
492 if x is None:
493 if 'Data' in self.description:
494 data = self.description['Data']
495 if 'Metadata' in self.description:
496 data.update(self.description['Metadata'])
497 else:
498 data = self.description
499 if name in data:
500 if isinstance(data[name], str):
501 return data[name]
502 elif isinstance(data[name], list):
503 return None
504 elif isinstance(data[name], dict):
505 for key, value in data[name].items():
506 return key
507 return name
508 else:
509 if 'Metadata' in self.description:
510 meta = self.description['Metadata']
511 else:
512 meta = self.description
513 if name in meta:
514 if isinstance(meta[name], list):
515 return meta[name][x]
516 elif isinstance(meta[name], dict):
517 for key, value in meta[name].items():
518 return value[x]
519 return 'channel{:02d}'.format(x)
520
1358 521 def writeMetadata(self, fp):
1359 522
1360 grp = fp.create_group("Metadata")
1361 grp.create_dataset('variables', data=self.tableDim, dtype=self.dtype)
523 if self.description:
524 if 'Metadata' in self.description:
525 grp = fp.create_group('Metadata')
526 else:
527 grp = fp
528 else:
529 grp = fp.create_group('Metadata')
1362 530
1363 531 for i in range(len(self.metadataList)):
1364 532 if not hasattr(self.dataOut, self.metadataList[i]):
1365 533 log.warning('Metadata: `{}` not found'.format(self.metadataList[i]), self.name)
1366 534 continue
1367 535 value = getattr(self.dataOut, self.metadataList[i])
1368 grp.create_dataset(self.metadataList[i], data=value)
536 if isinstance(value, bool):
537 if value is True:
538 value = 1
539 else:
540 value = 0
541 grp.create_dataset(self.getLabel(self.metadataList[i]), data=value)
1369 542 return
1370 543
1371 544 def writeData(self, fp):
1372 545
1373 grp = fp.create_group("Data")
546 if self.description:
547 if 'Data' in self.description:
548 grp = fp.create_group('Data')
549 else:
550 grp = fp
551 else:
552 grp = fp.create_group('Data')
553
1374 554 dtsets = []
1375 555 data = []
1376 556
1377 557 for dsInfo in self.dsList:
1378 558 if dsInfo['nDim'] == 0:
1379 559 ds = grp.create_dataset(
1380 dsInfo['variable'],
560 self.getLabel(dsInfo['variable']),
1381 561 (self.blocksPerFile, ),
1382 562 chunks=True,
1383 563 dtype=numpy.float64)
1384 564 dtsets.append(ds)
1385 565 data.append((dsInfo['variable'], -1))
1386 566 else:
1387 sgrp = grp.create_group(dsInfo['variable'])
567 label = self.getLabel(dsInfo['variable'])
568 if label is not None:
569 sgrp = grp.create_group(label)
570 else:
571 sgrp = grp
1388 572 for i in range(dsInfo['dsNumber']):
1389 573 ds = sgrp.create_dataset(
1390 'table{:02d}'.format(i),
574 self.getLabel(dsInfo['variable'], i),
1391 575 (self.blocksPerFile, ) + dsInfo['shape'][1:],
1392 chunks=True)
576 chunks=True,
577 dtype=dsInfo['dtype'])
1393 578 dtsets.append(ds)
1394 579 data.append((dsInfo['variable'], i))
1395 580 fp.flush()
1396 581
1397 582 log.log('Creating file: {}'.format(fp.filename), self.name)
1398 583
1399 584 self.ds = dtsets
1400 585 self.data = data
1401 586 self.firsttime = True
1402 587 self.blockIndex = 0
1403 588 return
1404 589
1405 590 def putData(self):
1406 591
1407 592 if (self.blockIndex == self.blocksPerFile) or self.timeFlag():
1408 593 self.closeFile()
1409 594 self.setNextFile()
1410 595
1411 596 for i, ds in enumerate(self.ds):
1412 597 attr, ch = self.data[i]
1413 598 if ch == -1:
1414 599 ds[self.blockIndex] = getattr(self.dataOut, attr)
1415 600 else:
1416 601 ds[self.blockIndex] = getattr(self.dataOut, attr)[ch]
1417 602
1418 603 self.fp.flush()
1419 604 self.blockIndex += 1
1420 605 log.log('Block No. {}/{}'.format(self.blockIndex, self.blocksPerFile), self.name)
1421 606
1422 607 return
1423 608
1424 609 def closeFile(self):
1425 610
1426 611 if self.blockIndex != self.blocksPerFile:
1427 612 for ds in self.ds:
1428 613 ds.resize(self.blockIndex, axis=0)
1429 614
1430 615 self.fp.flush()
1431 616 self.fp.close()
1432 617
1433 618 def close(self):
1434 619
1435 620 self.closeFile()
General Comments 0
You need to be logged in to leave comments. Login now