##// END OF EJS Templates
schain
RSS
update VOLT ACF PLOT
avaldez -
r1262:d63ca7e1ab2f
parent child
Show More
Side by Side Unified
Context file:
r1262:d63ca7e1ab2f Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,1264 +1,1281
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
11 11 from jroheaderIO import SystemHeader, RadarControllerHeader
12 12 from schainpy import cSchain
13 13
14 14
15 15 def getNumpyDtype(dataTypeCode):
16 16
17 17 if dataTypeCode == 0:
18 18 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
19 19 elif dataTypeCode == 1:
20 20 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
21 21 elif dataTypeCode == 2:
22 22 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
23 23 elif dataTypeCode == 3:
24 24 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
25 25 elif dataTypeCode == 4:
26 26 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
27 27 elif dataTypeCode == 5:
28 28 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
29 29 else:
30 30 raise ValueError, 'dataTypeCode was not defined'
31 31
32 32 return numpyDtype
33 33
34 34
35 35 def getDataTypeCode(numpyDtype):
36 36
37 37 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
38 38 datatype = 0
39 39 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
40 40 datatype = 1
41 41 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
42 42 datatype = 2
43 43 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
44 44 datatype = 3
45 45 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
46 46 datatype = 4
47 47 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
48 48 datatype = 5
49 49 else:
50 50 datatype = None
51 51
52 52 return datatype
53 53
54 54
55 55 def hildebrand_sekhon(data, navg):
56 56 """
57 57 This method is for the objective determination of the noise level in Doppler spectra. This
58 58 implementation technique is based on the fact that the standard deviation of the spectral
59 59 densities is equal to the mean spectral density for white Gaussian noise
60 60
61 61 Inputs:
62 62 Data : heights
63 63 navg : numbers of averages
64 64
65 65 Return:
66 66 -1 : any error
67 67 anoise : noise's level
68 68 """
69 69
70 70 sortdata = numpy.sort(data, axis=None)
71 71 # lenOfData = len(sortdata)
72 72 # nums_min = lenOfData*0.2
73 73 #
74 74 # if nums_min <= 5:
75 75 # nums_min = 5
76 76 #
77 77 # sump = 0.
78 78 #
79 79 # sumq = 0.
80 80 #
81 81 # j = 0
82 82 #
83 83 # cont = 1
84 84 #
85 85 # while((cont==1)and(j<lenOfData)):
86 86 #
87 87 # sump += sortdata[j]
88 88 #
89 89 # sumq += sortdata[j]**2
90 90 #
91 91 # if j > nums_min:
92 92 # rtest = float(j)/(j-1) + 1.0/navg
93 93 # if ((sumq*j) > (rtest*sump**2)):
94 94 # j = j - 1
95 95 # sump = sump - sortdata[j]
96 96 # sumq = sumq - sortdata[j]**2
97 97 # cont = 0
98 98 #
99 99 # j += 1
100 100 #
101 101 # lnoise = sump /j
102 102 #
103 103 # return lnoise
104 104
105 105 return cSchain.hildebrand_sekhon(sortdata, navg)
106 106
107 107
108 108 class Beam:
109 109
110 110 def __init__(self):
111 111 self.codeList = []
112 112 self.azimuthList = []
113 113 self.zenithList = []
114 114
115 115
116 116 class GenericData(object):
117 117
118 118 flagNoData = True
119 119
120 120 def copy(self, inputObj=None):
121 121
122 122 if inputObj == None:
123 123 return copy.deepcopy(self)
124 124
125 125 for key in inputObj.__dict__.keys():
126 126
127 127 attribute = inputObj.__dict__[key]
128 128
129 129 # If this attribute is a tuple or list
130 130 if type(inputObj.__dict__[key]) in (tuple, list):
131 131 self.__dict__[key] = attribute[:]
132 132 continue
133 133
134 134 # If this attribute is another object or instance
135 135 if hasattr(attribute, '__dict__'):
136 136 self.__dict__[key] = attribute.copy()
137 137 continue
138 138
139 139 self.__dict__[key] = inputObj.__dict__[key]
140 140
141 141 def deepcopy(self):
142 142
143 143 return copy.deepcopy(self)
144 144
145 145 def isEmpty(self):
146 146
147 147 return self.flagNoData
148 148
149 149
150 150 class JROData(GenericData):
151 151
152 152 # m_BasicHeader = BasicHeader()
153 153 # m_ProcessingHeader = ProcessingHeader()
154 154
155 155 systemHeaderObj = SystemHeader()
156 156
157 157 radarControllerHeaderObj = RadarControllerHeader()
158 158
159 159 # data = None
160 160
161 161 type = None
162 162
163 163 datatype = None # dtype but in string
164 164
165 165 # dtype = None
166 166
167 167 # nChannels = None
168 168
169 169 # nHeights = None
170 170
171 171 nProfiles = None
172 172
173 173 heightList = None
174 174
175 175 channelList = None
176 176
177 177 flagDiscontinuousBlock = False
178 178
179 179 useLocalTime = False
180 180
181 181 utctime = None
182 182
183 183 timeZone = None
184 184
185 185 dstFlag = None
186 186
187 187 errorCount = None
188 188
189 189 blocksize = None
190 190
191 191 # nCode = None
192 192 #
193 193 # nBaud = None
194 194 #
195 195 # code = None
196 196
197 197 flagDecodeData = False # asumo q la data no esta decodificada
198 198
199 199 flagDeflipData = False # asumo q la data no esta sin flip
200 200
201 201 flagShiftFFT = False
202 202
203 203 # ippSeconds = None
204 204
205 205 # timeInterval = None
206 206
207 207 nCohInt = None
208 208
209 209 # noise = None
210 210
211 211 windowOfFilter = 1
212 212
213 213 # Speed of ligth
214 214 C = 3e8
215 215
216 216 frequency = 49.92e6
217 217
218 218 realtime = False
219 219
220 220 beacon_heiIndexList = None
221 221
222 222 last_block = None
223 223
224 224 blocknow = None
225 225
226 226 azimuth = None
227 227
228 228 zenith = None
229 229
230 230 beam = Beam()
231 231
232 232 profileIndex = None
233 233
234 234 def getNoise(self):
235 235
236 236 raise NotImplementedError
237 237
238 238 def getNChannels(self):
239 239
240 240 return len(self.channelList)
241 241
242 242 def getChannelIndexList(self):
243 243
244 244 return range(self.nChannels)
245 245
246 246 def getNHeights(self):
247 247
248 248 return len(self.heightList)
249 249
250 250 def getHeiRange(self, extrapoints=0):
251 251
252 252 heis = self.heightList
253 253 # deltah = self.heightList[1] - self.heightList[0]
254 254 #
255 255 # heis.append(self.heightList[-1])
256 256
257 257 return heis
258 258
259 259 def getDeltaH(self):
260 260
261 261 delta = self.heightList[1] - self.heightList[0]
262 262
263 263 return delta
264 264
265 265 def getltctime(self):
266 266
267 267 if self.useLocalTime:
268 268 return self.utctime - self.timeZone * 60
269 269
270 270 return self.utctime
271 271
272 272 def getDatatime(self):
273 273
274 274 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
275 275 return datatimeValue
276 276
277 277 def getTimeRange(self):
278 278
279 279 datatime = []
280 280
281 281 datatime.append(self.ltctime)
282 282 datatime.append(self.ltctime + self.timeInterval + 1)
283 283
284 284 datatime = numpy.array(datatime)
285 285
286 286 return datatime
287 287
288 288 def getFmaxTimeResponse(self):
289 289
290 290 period = (10**-6) * self.getDeltaH() / (0.15)
291 291
292 292 PRF = 1. / (period * self.nCohInt)
293 293
294 294 fmax = PRF
295 295
296 296 return fmax
297 297
298 298 def getFmax(self):
299 299 PRF = 1. / (self.ippSeconds * self.nCohInt)
300 300 #print "ippSeconds",self.ippSeconds
301 301 #print "nCohInt",self.nIncohInt
302 302 #print PRF
303 303 #import time
304 304 #time.sleep(30)
305 305 fmax = PRF
306 306 return fmax
307 307
308 308 def getVmax(self):
309 309
310 310 _lambda = self.C / self.frequency
311 311
312 312 vmax = self.getFmax() * _lambda / 2
313 313
314 314 return vmax
315 315
316 316 def get_ippSeconds(self):
317 317 '''
318 318 '''
319 319 return self.radarControllerHeaderObj.ippSeconds
320 320
321 321 def set_ippSeconds(self, ippSeconds):
322 322 '''
323 323 '''
324 324
325 325 self.radarControllerHeaderObj.ippSeconds = ippSeconds
326 326
327 327 return
328 328
329 329 def get_dtype(self):
330 330 '''
331 331 '''
332 332 return getNumpyDtype(self.datatype)
333 333
334 334 def set_dtype(self, numpyDtype):
335 335 '''
336 336 '''
337 337
338 338 self.datatype = getDataTypeCode(numpyDtype)
339 339
340 340 def get_code(self):
341 341 '''
342 342 '''
343 343 return self.radarControllerHeaderObj.code
344 344
345 345 def set_code(self, code):
346 346 '''
347 347 '''
348 348 self.radarControllerHeaderObj.code = code
349 349
350 350 return
351 351
352 352 def get_ncode(self):
353 353 '''
354 354 '''
355 355 return self.radarControllerHeaderObj.nCode
356 356
357 357 def set_ncode(self, nCode):
358 358 '''
359 359 '''
360 360 self.radarControllerHeaderObj.nCode = nCode
361 361
362 362 return
363 363
364 364 def get_nbaud(self):
365 365 '''
366 366 '''
367 367 return self.radarControllerHeaderObj.nBaud
368 368
369 369 def set_nbaud(self, nBaud):
370 370 '''
371 371 '''
372 372 self.radarControllerHeaderObj.nBaud = nBaud
373 373
374 374 return
375 375
376 376 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
377 377 channelIndexList = property(
378 378 getChannelIndexList, "I'm the 'channelIndexList' property.")
379 379 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
380 380 #noise = property(getNoise, "I'm the 'nHeights' property.")
381 381 datatime = property(getDatatime, "I'm the 'datatime' property")
382 382 ltctime = property(getltctime, "I'm the 'ltctime' property")
383 383 ippSeconds = property(get_ippSeconds, set_ippSeconds)
384 384 dtype = property(get_dtype, set_dtype)
385 385 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
386 386 code = property(get_code, set_code)
387 387 nCode = property(get_ncode, set_ncode)
388 388 nBaud = property(get_nbaud, set_nbaud)
389 389
390 390
391 391 class Voltage(JROData):
392 392
393 393 # data es un numpy array de 2 dmensiones (canales, alturas)
394 394 data = None
395 395
396 396 def __init__(self):
397 397 '''
398 398 Constructor
399 399 '''
400 400
401 401 self.useLocalTime = True
402 402
403 403 self.radarControllerHeaderObj = RadarControllerHeader()
404 404
405 405 self.systemHeaderObj = SystemHeader()
406 406
407 407 self.type = "Voltage"
408 408
409 409 self.data = None
410 410
411 411 # self.dtype = None
412 412
413 413 # self.nChannels = 0
414 414
415 415 # self.nHeights = 0
416 416
417 417 self.nProfiles = None
418 418
419 419 self.heightList = None
420 420
421 421 self.channelList = None
422 422
423 423 # self.channelIndexList = None
424 424
425 425 self.flagNoData = True
426 426
427 427 self.flagDiscontinuousBlock = False
428 428
429 429 self.utctime = None
430 430
431 431 self.timeZone = None
432 432
433 433 self.dstFlag = None
434 434
435 435 self.errorCount = None
436 436
437 437 self.nCohInt = None
438 438
439 439 self.blocksize = None
440 440
441 441 self.flagDecodeData = False # asumo q la data no esta decodificada
442 442
443 443 self.flagDeflipData = False # asumo q la data no esta sin flip
444 444
445 445 self.flagShiftFFT = False
446 446
447 447 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
448 448
449 self.ippFactor = 1
450
449 451 self.profileIndex = 0
450 452
451 453 def getNoisebyHildebrand(self, channel=None):
452 454 """
453 455 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
454 456
455 457 Return:
456 458 noiselevel
457 459 """
458 460
459 461 if channel != None:
460 462 data = self.data[channel]
461 463 nChannels = 1
462 464 else:
463 465 data = self.data
464 466 nChannels = self.nChannels
465 467
466 468 noise = numpy.zeros(nChannels)
467 469 power = data * numpy.conjugate(data)
468 470
469 471 for thisChannel in range(nChannels):
470 472 if nChannels == 1:
471 473 daux = power[:].real
472 474 else:
473 475 daux = power[thisChannel, :].real
474 476 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
475 477
476 478 return noise
477 479
478 480 def getNoise(self, type=1, channel=None):
479 481
480 482 if type == 1:
481 483 noise = self.getNoisebyHildebrand(channel)
482 484
483 485 return noise
484 486
485 487 def getPower(self, channel=None):
486 488
487 489 if channel != None:
488 490 data = self.data[channel]
489 491 else:
490 492 data = self.data
491 493
492 494 power = data * numpy.conjugate(data)
493 495 powerdB = 10 * numpy.log10(power.real)
494 496 powerdB = numpy.squeeze(powerdB)
495 497
496 498 return powerdB
497 499
498 500 def getTimeInterval(self):
499 501
500 502 timeInterval = self.ippSeconds * self.nCohInt
501 503
502 504 return timeInterval
503 505
506 def getAcfRange(self, extrapoints=0):
507 #print "GET ACF RANGE"
508 #print "NFFTPoints",self.nFFTPoints
509 #print "IPPFactor", self.ippFactor
510 #deltafreq = 10. / ( self.getFmax() / (self.nFFTPoints * self.ippFactor) )
511 deltatime = 1./(self.getFmax()/ self.ippFactor)
512 #print "getFmax",self.getFmax()
513 #import time
514 #time.sleep(30)
515 timerange = deltatime * \
516 (numpy.arange(self.nProfiles + extrapoints))#- self.nProfiles / 2.)
517 #- deltafreq / 2
518 #print "timerange",timerange
519 return timerange
520
504 521 noise = property(getNoise, "I'm the 'nHeights' property.")
505 522 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
506 523
507 524
508 525 class Spectra(JROData):
509 526
510 527 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
511 528 data_spc = None
512 529
513 530 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
514 531 data_cspc = None
515 532
516 533 # data dc es un numpy array de 2 dmensiones (canales, alturas)
517 534 data_dc = None
518 535
519 536 # data power
520 537 data_pwr = None
521 538
522 539 nFFTPoints = None
523 540
524 541 # nPairs = None
525 542
526 543 pairsList = None
527 544
528 545 nIncohInt = None
529 546
530 547 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
531 548
532 549 nCohInt = None # se requiere para determinar el valor de timeInterval
533 550
534 551 ippFactor = None
535 552
536 553 profileIndex = 0
537 554
538 555 plotting = "spectra"
539 556
540 557 def __init__(self):
541 558 '''
542 559 Constructor
543 560 '''
544 561
545 562 self.useLocalTime = True
546 563
547 564 self.radarControllerHeaderObj = RadarControllerHeader()
548 565
549 566 self.systemHeaderObj = SystemHeader()
550 567
551 568 self.type = "Spectra"
552 569
553 570 # self.data = None
554 571
555 572 # self.dtype = None
556 573
557 574 # self.nChannels = 0
558 575
559 576 # self.nHeights = 0
560 577
561 578 self.nProfiles = None
562 579
563 580 self.heightList = None
564 581
565 582 self.channelList = None
566 583
567 584 # self.channelIndexList = None
568 585
569 586 self.pairsList = None
570 587
571 588 self.flagNoData = True
572 589
573 590 self.flagDiscontinuousBlock = False
574 591
575 592 self.utctime = None
576 593
577 594 self.nCohInt = None
578 595
579 596 self.nIncohInt = None
580 597
581 598 self.blocksize = None
582 599
583 600 self.nFFTPoints = None
584 601
585 602 self.wavelength = None
586 603
587 604 self.flagDecodeData = False # asumo q la data no esta decodificada
588 605
589 606 self.flagDeflipData = False # asumo q la data no esta sin flip
590 607
591 608 self.flagShiftFFT = False
592 609
593 610 self.ippFactor = 1
594 611
595 612 #self.noise = None
596 613
597 614 self.beacon_heiIndexList = []
598 615
599 616 self.noise_estimation = None
600 617
601 618 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
602 619 """
603 620 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
604 621
605 622 Return:
606 623 noiselevel
607 624 """
608 625
609 626 noise = numpy.zeros(self.nChannels)
610 627
611 628 for channel in range(self.nChannels):
612 629 daux = self.data_spc[channel,
613 630 xmin_index:xmax_index, ymin_index:ymax_index]
614 631
615 632 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
616 633
617 634 return noise
618 635
619 636 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
620 637
621 638 if self.noise_estimation is not None:
622 639 # this was estimated by getNoise Operation defined in jroproc_spectra.py
623 640 return self.noise_estimation
624 641 else:
625 642 noise = self.getNoisebyHildebrand(
626 643 xmin_index, xmax_index, ymin_index, ymax_index)
627 644 return noise
628 645
629 646 def getFreqRangeTimeResponse(self, extrapoints=0):
630 647
631 648 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
632 649 freqrange = deltafreq * \
633 650 (numpy.arange(self.nFFTPoints + extrapoints) -
634 651 self.nFFTPoints / 2.) - deltafreq / 2
635 652
636 653 return freqrange
637 654
638 655 def getAcfRange(self, extrapoints=0):
639 656 #print "GET ACF RANGE"
640 657 #print "NFFTPoints",self.nFFTPoints
641 658 #print "IPPFactor", self.ippFactor
642 659 #deltafreq = 10. / ( self.getFmax() / (self.nFFTPoints * self.ippFactor) )
643 660 deltatime = 1./(self.getFmax()/ self.ippFactor)
644 661 #print "getFmax",self.getFmax()
645 662 #import time
646 663 #time.sleep(30)
647 664 timerange = deltatime * \
648 665 (numpy.arange(self.nFFTPoints + extrapoints) -
649 666 self.nFFTPoints / 2.)
650 667 #- deltafreq / 2
651 668 #print "timerange",timerange
652 669 return timerange
653 670
654 671 def getFreqRange(self, extrapoints=0):
655 672
656 673 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
657 674 #print "deltafreq", deltafreq
658 675 freqrange = deltafreq * \
659 676 (numpy.arange(self.nFFTPoints + extrapoints) -
660 677 self.nFFTPoints / 2.) - deltafreq / 2
661 678 #print "freqrange",freqrange
662 679 return freqrange
663 680
664 681 def getVelRange(self, extrapoints=0):
665 682
666 683 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
667 684 velrange = deltav * (numpy.arange(self.nFFTPoints +
668 685 extrapoints) - self.nFFTPoints / 2.) # - deltav/2
669 686
670 687 return velrange
671 688
672 689 def getNPairs(self):
673 690
674 691 return len(self.pairsList)
675 692
676 693 def getPairsIndexList(self):
677 694
678 695 return range(self.nPairs)
679 696
680 697 def getNormFactor(self):
681 698
682 699 pwcode = 1
683 700
684 701 if self.flagDecodeData:
685 702 pwcode = numpy.sum(self.code[0]**2)
686 703 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
687 704 normFactor = self.nProfiles * self.nIncohInt * \
688 705 self.nCohInt * pwcode * self.windowOfFilter
689 706
690 707 return normFactor
691 708
692 709 def getFlagCspc(self):
693 710
694 711 if self.data_cspc is None:
695 712 return True
696 713
697 714 return False
698 715
699 716 def getFlagDc(self):
700 717
701 718 if self.data_dc is None:
702 719 return True
703 720
704 721 return False
705 722
706 723 def getTimeInterval(self):
707 724
708 725 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
709 726
710 727 return timeInterval
711 728
712 729 def getPower(self):
713 730
714 731 factor = self.normFactor
715 732 z = self.data_spc / factor
716 733 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
717 734 avg = numpy.average(z, axis=1)
718 735
719 736 return 10 * numpy.log10(avg)
720 737
721 738 def getCoherence(self, pairsList=None, phase=False):
722 739
723 740 z = []
724 741 if pairsList is None:
725 742 pairsIndexList = self.pairsIndexList
726 743 else:
727 744 pairsIndexList = []
728 745 for pair in pairsList:
729 746 if pair not in self.pairsList:
730 747 raise ValueError, "Pair %s is not in dataOut.pairsList" % (
731 748 pair)
732 749 pairsIndexList.append(self.pairsList.index(pair))
733 750 for i in range(len(pairsIndexList)):
734 751 pair = self.pairsList[pairsIndexList[i]]
735 752 ccf = numpy.average(
736 753 self.data_cspc[pairsIndexList[i], :, :], axis=0)
737 754 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
738 755 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
739 756 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
740 757 if phase:
741 758 data = numpy.arctan2(avgcoherenceComplex.imag,
742 759 avgcoherenceComplex.real) * 180 / numpy.pi
743 760 else:
744 761 data = numpy.abs(avgcoherenceComplex)
745 762
746 763 z.append(data)
747 764
748 765 return numpy.array(z)
749 766
750 767 def setValue(self, value):
751 768
752 769 print "This property should not be initialized"
753 770
754 771 return
755 772
756 773 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
757 774 pairsIndexList = property(
758 775 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
759 776 normFactor = property(getNormFactor, setValue,
760 777 "I'm the 'getNormFactor' property.")
761 778 flag_cspc = property(getFlagCspc, setValue)
762 779 flag_dc = property(getFlagDc, setValue)
763 780 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
764 781 timeInterval = property(getTimeInterval, setValue,
765 782 "I'm the 'timeInterval' property")
766 783
767 784
768 785 class SpectraHeis(Spectra):
769 786
770 787 data_spc = None
771 788
772 789 data_cspc = None
773 790
774 791 data_dc = None
775 792
776 793 nFFTPoints = None
777 794
778 795 # nPairs = None
779 796
780 797 pairsList = None
781 798
782 799 nCohInt = None
783 800
784 801 nIncohInt = None
785 802
786 803 def __init__(self):
787 804
788 805 self.radarControllerHeaderObj = RadarControllerHeader()
789 806
790 807 self.systemHeaderObj = SystemHeader()
791 808
792 809 self.type = "SpectraHeis"
793 810
794 811 # self.dtype = None
795 812
796 813 # self.nChannels = 0
797 814
798 815 # self.nHeights = 0
799 816
800 817 self.nProfiles = None
801 818
802 819 self.heightList = None
803 820
804 821 self.channelList = None
805 822
806 823 # self.channelIndexList = None
807 824
808 825 self.flagNoData = True
809 826
810 827 self.flagDiscontinuousBlock = False
811 828
812 829 # self.nPairs = 0
813 830
814 831 self.utctime = None
815 832
816 833 self.blocksize = None
817 834
818 835 self.profileIndex = 0
819 836
820 837 self.nCohInt = 1
821 838
822 839 self.nIncohInt = 1
823 840
824 841 def getNormFactor(self):
825 842 pwcode = 1
826 843 if self.flagDecodeData:
827 844 pwcode = numpy.sum(self.code[0]**2)
828 845
829 846 normFactor = self.nIncohInt * self.nCohInt * pwcode
830 847
831 848 return normFactor
832 849
833 850 def getTimeInterval(self):
834 851
835 852 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
836 853
837 854 return timeInterval
838 855
839 856 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
840 857 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
841 858
842 859
843 860 class Fits(JROData):
844 861
845 862 heightList = None
846 863
847 864 channelList = None
848 865
849 866 flagNoData = True
850 867
851 868 flagDiscontinuousBlock = False
852 869
853 870 useLocalTime = False
854 871
855 872 utctime = None
856 873
857 874 timeZone = None
858 875
859 876 # ippSeconds = None
860 877
861 878 # timeInterval = None
862 879
863 880 nCohInt = None
864 881
865 882 nIncohInt = None
866 883
867 884 noise = None
868 885
869 886 windowOfFilter = 1
870 887
871 888 # Speed of ligth
872 889 C = 3e8
873 890
874 891 frequency = 49.92e6
875 892
876 893 realtime = False
877 894
878 895 def __init__(self):
879 896
880 897 self.type = "Fits"
881 898
882 899 self.nProfiles = None
883 900
884 901 self.heightList = None
885 902
886 903 self.channelList = None
887 904
888 905 # self.channelIndexList = None
889 906
890 907 self.flagNoData = True
891 908
892 909 self.utctime = None
893 910
894 911 self.nCohInt = 1
895 912
896 913 self.nIncohInt = 1
897 914
898 915 self.useLocalTime = True
899 916
900 917 self.profileIndex = 0
901 918
902 919 # self.utctime = None
903 920 # self.timeZone = None
904 921 # self.ltctime = None
905 922 # self.timeInterval = None
906 923 # self.header = None
907 924 # self.data_header = None
908 925 # self.data = None
909 926 # self.datatime = None
910 927 # self.flagNoData = False
911 928 # self.expName = ''
912 929 # self.nChannels = None
913 930 # self.nSamples = None
914 931 # self.dataBlocksPerFile = None
915 932 # self.comments = ''
916 933 #
917 934
918 935 def getltctime(self):
919 936
920 937 if self.useLocalTime:
921 938 return self.utctime - self.timeZone * 60
922 939
923 940 return self.utctime
924 941
925 942 def getDatatime(self):
926 943
927 944 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
928 945 return datatime
929 946
930 947 def getTimeRange(self):
931 948
932 949 datatime = []
933 950
934 951 datatime.append(self.ltctime)
935 952 datatime.append(self.ltctime + self.timeInterval)
936 953
937 954 datatime = numpy.array(datatime)
938 955
939 956 return datatime
940 957
941 958 def getHeiRange(self):
942 959
943 960 heis = self.heightList
944 961
945 962 return heis
946 963
947 964 def getNHeights(self):
948 965
949 966 return len(self.heightList)
950 967
951 968 def getNChannels(self):
952 969
953 970 return len(self.channelList)
954 971
955 972 def getChannelIndexList(self):
956 973
957 974 return range(self.nChannels)
958 975
959 976 def getNoise(self, type=1):
960 977
961 978 #noise = numpy.zeros(self.nChannels)
962 979
963 980 if type == 1:
964 981 noise = self.getNoisebyHildebrand()
965 982
966 983 if type == 2:
967 984 noise = self.getNoisebySort()
968 985
969 986 if type == 3:
970 987 noise = self.getNoisebyWindow()
971 988
972 989 return noise
973 990
974 991 def getTimeInterval(self):
975 992
976 993 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
977 994
978 995 return timeInterval
979 996
980 997 datatime = property(getDatatime, "I'm the 'datatime' property")
981 998 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
982 999 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
983 1000 channelIndexList = property(
984 1001 getChannelIndexList, "I'm the 'channelIndexList' property.")
985 1002 noise = property(getNoise, "I'm the 'nHeights' property.")
986 1003
987 1004 ltctime = property(getltctime, "I'm the 'ltctime' property")
988 1005 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
989 1006
990 1007
991 1008 class Correlation(JROData):
992 1009
993 1010 noise = None
994 1011
995 1012 SNR = None
996 1013
997 1014 #--------------------------------------------------
998 1015
999 1016 mode = None
1000 1017
1001 1018 split = False
1002 1019
1003 1020 data_cf = None
1004 1021
1005 1022 lags = None
1006 1023
1007 1024 lagRange = None
1008 1025
1009 1026 pairsList = None
1010 1027
1011 1028 normFactor = None
1012 1029
1013 1030 #--------------------------------------------------
1014 1031
1015 1032 # calculateVelocity = None
1016 1033
1017 1034 nLags = None
1018 1035
1019 1036 nPairs = None
1020 1037
1021 1038 nAvg = None
1022 1039
1023 1040 def __init__(self):
1024 1041 '''
1025 1042 Constructor
1026 1043 '''
1027 1044 self.radarControllerHeaderObj = RadarControllerHeader()
1028 1045
1029 1046 self.systemHeaderObj = SystemHeader()
1030 1047
1031 1048 self.type = "Correlation"
1032 1049
1033 1050 self.data = None
1034 1051
1035 1052 self.dtype = None
1036 1053
1037 1054 self.nProfiles = None
1038 1055
1039 1056 self.heightList = None
1040 1057
1041 1058 self.channelList = None
1042 1059
1043 1060 self.flagNoData = True
1044 1061
1045 1062 self.flagDiscontinuousBlock = False
1046 1063
1047 1064 self.utctime = None
1048 1065
1049 1066 self.timeZone = None
1050 1067
1051 1068 self.dstFlag = None
1052 1069
1053 1070 self.errorCount = None
1054 1071
1055 1072 self.blocksize = None
1056 1073
1057 1074 self.flagDecodeData = False # asumo q la data no esta decodificada
1058 1075
1059 1076 self.flagDeflipData = False # asumo q la data no esta sin flip
1060 1077
1061 1078 self.pairsList = None
1062 1079
1063 1080 self.nPoints = None
1064 1081
1065 1082 def getPairsList(self):
1066 1083
1067 1084 return self.pairsList
1068 1085
1069 1086 def getNoise(self, mode=2):
1070 1087
1071 1088 indR = numpy.where(self.lagR == 0)[0][0]
1072 1089 indT = numpy.where(self.lagT == 0)[0][0]
1073 1090
1074 1091 jspectra0 = self.data_corr[:, :, indR, :]
1075 1092 jspectra = copy.copy(jspectra0)
1076 1093
1077 1094 num_chan = jspectra.shape[0]
1078 1095 num_hei = jspectra.shape[2]
1079 1096
1080 1097 freq_dc = jspectra.shape[1] / 2
1081 1098 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
1082 1099
1083 1100 if ind_vel[0] < 0:
1084 1101 ind_vel[range(0, 1)] = ind_vel[range(0, 1)] + self.num_prof
1085 1102
1086 1103 if mode == 1:
1087 1104 jspectra[:, freq_dc, :] = (
1088 1105 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
1089 1106
1090 1107 if mode == 2:
1091 1108
1092 1109 vel = numpy.array([-2, -1, 1, 2])
1093 1110 xx = numpy.zeros([4, 4])
1094 1111
1095 1112 for fil in range(4):
1096 1113 xx[fil, :] = vel[fil]**numpy.asarray(range(4))
1097 1114
1098 1115 xx_inv = numpy.linalg.inv(xx)
1099 1116 xx_aux = xx_inv[0, :]
1100 1117
1101 1118 for ich in range(num_chan):
1102 1119 yy = jspectra[ich, ind_vel, :]
1103 1120 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
1104 1121
1105 1122 junkid = jspectra[ich, freq_dc, :] <= 0
1106 1123 cjunkid = sum(junkid)
1107 1124
1108 1125 if cjunkid.any():
1109 1126 jspectra[ich, freq_dc, junkid.nonzero()] = (
1110 1127 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
1111 1128
1112 1129 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
1113 1130
1114 1131 return noise
1115 1132
1116 1133 def getTimeInterval(self):
1117 1134
1118 1135 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
1119 1136
1120 1137 return timeInterval
1121 1138
1122 1139 def splitFunctions(self):
1123 1140
1124 1141 pairsList = self.pairsList
1125 1142 ccf_pairs = []
1126 1143 acf_pairs = []
1127 1144 ccf_ind = []
1128 1145 acf_ind = []
1129 1146 for l in range(len(pairsList)):
1130 1147 chan0 = pairsList[l][0]
1131 1148 chan1 = pairsList[l][1]
1132 1149
1133 1150 # Obteniendo pares de Autocorrelacion
1134 1151 if chan0 == chan1:
1135 1152 acf_pairs.append(chan0)
1136 1153 acf_ind.append(l)
1137 1154 else:
1138 1155 ccf_pairs.append(pairsList[l])
1139 1156 ccf_ind.append(l)
1140 1157
1141 1158 data_acf = self.data_cf[acf_ind]
1142 1159 data_ccf = self.data_cf[ccf_ind]
1143 1160
1144 1161 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1145 1162
1146 1163 def getNormFactor(self):
1147 1164 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1148 1165 acf_pairs = numpy.array(acf_pairs)
1149 1166 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1150 1167
1151 1168 for p in range(self.nPairs):
1152 1169 pair = self.pairsList[p]
1153 1170
1154 1171 ch0 = pair[0]
1155 1172 ch1 = pair[1]
1156 1173
1157 1174 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1158 1175 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1159 1176 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1160 1177
1161 1178 return normFactor
1162 1179
1163 1180 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1164 1181 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1165 1182
1166 1183
1167 1184 class Parameters(Spectra):
1168 1185
1169 1186 experimentInfo = None # Information about the experiment
1170 1187
1171 1188 # Information from previous data
1172 1189
1173 1190 inputUnit = None # Type of data to be processed
1174 1191
1175 1192 operation = None # Type of operation to parametrize
1176 1193
1177 1194 # normFactor = None #Normalization Factor
1178 1195
1179 1196 groupList = None # List of Pairs, Groups, etc
1180 1197
1181 1198 # Parameters
1182 1199
1183 1200 data_param = None # Parameters obtained
1184 1201
1185 1202 data_pre = None # Data Pre Parametrization
1186 1203
1187 1204 data_SNR = None # Signal to Noise Ratio
1188 1205
1189 1206 # heightRange = None #Heights
1190 1207
1191 1208 abscissaList = None # Abscissa, can be velocities, lags or time
1192 1209
1193 1210 # noise = None #Noise Potency
1194 1211
1195 1212 utctimeInit = None # Initial UTC time
1196 1213
1197 1214 paramInterval = None # Time interval to calculate Parameters in seconds
1198 1215
1199 1216 useLocalTime = True
1200 1217
1201 1218 # Fitting
1202 1219
1203 1220 data_error = None # Error of the estimation
1204 1221
1205 1222 constants = None
1206 1223
1207 1224 library = None
1208 1225
1209 1226 # Output signal
1210 1227
1211 1228 outputInterval = None # Time interval to calculate output signal in seconds
1212 1229
1213 1230 data_output = None # Out signal
1214 1231
1215 1232 nAvg = None
1216 1233
1217 1234 noise_estimation = None
1218 1235
1219 1236 GauSPC = None # Fit gaussian SPC
1220 1237
1221 1238 def __init__(self):
1222 1239 '''
1223 1240 Constructor
1224 1241 '''
1225 1242 self.radarControllerHeaderObj = RadarControllerHeader()
1226 1243
1227 1244 self.systemHeaderObj = SystemHeader()
1228 1245
1229 1246 self.type = "Parameters"
1230 1247
1231 1248 def getTimeRange1(self, interval):
1232 1249
1233 1250 datatime = []
1234 1251
1235 1252 if self.useLocalTime:
1236 1253 time1 = self.utctimeInit - self.timeZone * 60
1237 1254 else:
1238 1255 time1 = self.utctimeInit
1239 1256
1240 1257 datatime.append(time1)
1241 1258 datatime.append(time1 + interval)
1242 1259 datatime = numpy.array(datatime)
1243 1260
1244 1261 return datatime
1245 1262
1246 1263 def getTimeInterval(self):
1247 1264
1248 1265 if hasattr(self, 'timeInterval1'):
1249 1266 return self.timeInterval1
1250 1267 else:
1251 1268 return self.paramInterval
1252 1269
1253 1270 def setValue(self, value):
1254 1271
1255 1272 print "This property should not be initialized"
1256 1273
1257 1274 return
1258 1275
1259 1276 def getNoise(self):
1260 1277
1261 1278 return self.spc_noise
1262 1279
1263 1280 timeInterval = property(getTimeInterval)
1264 1281 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
Show file before | Show file after | Hide comments
@@ -1,225 +1,426
1 1 '''
2 2 Created on Jul 9, 2014
3 3
4 4 @author: roj-idl71
5 5 '''
6 6 import os
7 7 import datetime
8 8 import numpy
9 9
10 10 from figure import Figure
11 from plotting_codes import *
11 12
12 13 class Scope(Figure):
13 14
14 15 isConfig = None
15 16
16 17 def __init__(self, **kwargs):
17 18 Figure.__init__(self, **kwargs)
18 19 self.isConfig = False
19 20 self.WIDTH = 300
20 21 self.HEIGHT = 200
21 22 self.counter_imagwr = 0
22 23
23 24 def getSubplots(self):
24 25
25 26 nrow = self.nplots
26 27 ncol = 3
27 28 return nrow, ncol
28 29
29 30 def setup(self, id, nplots, wintitle, show):
30 31
31 32 self.nplots = nplots
32 33
33 34 self.createFigure(id=id,
34 35 wintitle=wintitle,
35 36 show=show)
36 37
37 38 nrow,ncol = self.getSubplots()
38 39 colspan = 3
39 40 rowspan = 1
40 41
41 42 for i in range(nplots):
42 43 self.addAxes(nrow, ncol, i, 0, colspan, rowspan)
43 44
44 45 def plot_iq(self, x, y, id, channelIndexList, thisDatetime, wintitle, show, xmin, xmax, ymin, ymax):
45 46 yreal = y[channelIndexList,:].real
46 47 yimag = y[channelIndexList,:].imag
47 48
48 49 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
49 50 xlabel = "Range (Km)"
50 51 ylabel = "Intensity - IQ"
51 52
52 53 if not self.isConfig:
53 54 nplots = len(channelIndexList)
54 55
55 56 self.setup(id=id,
56 57 nplots=nplots,
57 58 wintitle='',
58 59 show=show)
59 60
60 61 if xmin == None: xmin = numpy.nanmin(x)
61 62 if xmax == None: xmax = numpy.nanmax(x)
62 63 if ymin == None: ymin = min(numpy.nanmin(yreal),numpy.nanmin(yimag))
63 64 if ymax == None: ymax = max(numpy.nanmax(yreal),numpy.nanmax(yimag))
64 65
65 66 self.isConfig = True
66 67
67 68 self.setWinTitle(title)
68 69
69 70 for i in range(len(self.axesList)):
70 71 title = "Channel %d" %(i)
71 72 axes = self.axesList[i]
72 73
73 74 axes.pline(x, yreal[i,:],
74 75 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
75 76 xlabel=xlabel, ylabel=ylabel, title=title)
76 77
77 78 axes.addpline(x, yimag[i,:], idline=1, color="red", linestyle="solid", lw=2)
78 79
79 80 def plot_power(self, x, y, id, channelIndexList, thisDatetime, wintitle, show, xmin, xmax, ymin, ymax):
80 81 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
81 82 yreal = y.real
82 83
83 84 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
84 85 xlabel = "Range (Km)"
85 86 ylabel = "Intensity"
86 87
87 88 if not self.isConfig:
88 89 nplots = len(channelIndexList)
89 90
90 91 self.setup(id=id,
91 92 nplots=nplots,
92 93 wintitle='',
93 94 show=show)
94 95
95 96 if xmin == None: xmin = numpy.nanmin(x)
96 97 if xmax == None: xmax = numpy.nanmax(x)
97 98 if ymin == None: ymin = numpy.nanmin(yreal)
98 99 if ymax == None: ymax = numpy.nanmax(yreal)
99 100
100 101 self.isConfig = True
101 102
102 103 self.setWinTitle(title)
103 104
104 105 for i in range(len(self.axesList)):
105 106 title = "Channel %d" %(i)
106 107 axes = self.axesList[i]
107 108 ychannel = yreal[i,:]
108 109 axes.pline(x, ychannel,
109 110 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
110 111 xlabel=xlabel, ylabel=ylabel, title=title)
111 112
112 113
113 114 def run(self, dataOut, id, wintitle="", channelList=None,
114 115 xmin=None, xmax=None, ymin=None, ymax=None, save=False,
115 116 figpath='./', figfile=None, show=True, wr_period=1,
116 117 ftp=False, server=None, folder=None, username=None, password=None, type='power', **kwargs):
117 118
118 119 """
119 120
120 121 Input:
121 122 dataOut :
122 123 id :
123 124 wintitle :
124 125 channelList :
125 126 xmin : None,
126 127 xmax : None,
127 128 ymin : None,
128 129 ymax : None,
129 130 """
130 131
131 132 if channelList == None:
132 133 channelIndexList = dataOut.channelIndexList
133 134 else:
134 135 channelIndexList = []
135 136 for channel in channelList:
136 137 if channel not in dataOut.channelList:
137 138 raise ValueError, "Channel %d is not in dataOut.channelList"
138 139 channelIndexList.append(dataOut.channelList.index(channel))
139 140
140 141 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
141 142
142 143 if dataOut.flagDataAsBlock:
143 144
144 145 for i in range(dataOut.nProfiles):
145 146
146 147 wintitle1 = wintitle + " [Profile = %d] " %i
147 148
148 149 if type == "power":
149 150 self.plot_power(dataOut.heightList,
150 151 dataOut.data[:,i,:],
151 152 id,
152 153 channelIndexList,
153 154 thisDatetime,
154 155 wintitle1,
155 156 show,
156 157 xmin,
157 158 xmax,
158 159 ymin,
159 160 ymax)
160 161
161 162 if type == "iq":
162 163 self.plot_iq(dataOut.heightList,
163 164 dataOut.data[:,i,:],
164 165 id,
165 166 channelIndexList,
166 167 thisDatetime,
167 168 wintitle1,
168 169 show,
169 170 xmin,
170 171 xmax,
171 172 ymin,
172 173 ymax)
173 174
174 175 self.draw()
175 176
176 177 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
177 178 figfile = self.getFilename(name = str_datetime) + "_" + str(i)
178 179
179 180 self.save(figpath=figpath,
180 181 figfile=figfile,
181 182 save=save,
182 183 ftp=ftp,
183 184 wr_period=wr_period,
184 185 thisDatetime=thisDatetime)
185 186
186 187 else:
187 188 wintitle += " [Profile = %d] " %dataOut.profileIndex
188 189
189 190 if type == "power":
190 191 self.plot_power(dataOut.heightList,
191 192 dataOut.data,
192 193 id,
193 194 channelIndexList,
194 195 thisDatetime,
195 196 wintitle,
196 197 show,
197 198 xmin,
198 199 xmax,
199 200 ymin,
200 201 ymax)
201 202
202 203 if type == "iq":
203 204 self.plot_iq(dataOut.heightList,
204 205 dataOut.data,
205 206 id,
206 207 channelIndexList,
207 208 thisDatetime,
208 209 wintitle,
209 210 show,
210 211 xmin,
211 212 xmax,
212 213 ymin,
213 214 ymax)
214 215
215 216 self.draw()
216 217
217 218 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S") + "_" + str(dataOut.profileIndex)
218 219 figfile = self.getFilename(name = str_datetime)
219 220
220 221 self.save(figpath=figpath,
221 222 figfile=figfile,
222 223 save=save,
223 224 ftp=ftp,
224 225 wr_period=wr_period,
225 226 thisDatetime=thisDatetime)
227
228
229 class VoltACFPLot(Figure):
230
231 isConfig = None
232 __nsubplots = None
233 PREFIX = 'voltacf'
234
235 def __init__(self, **kwargs):
236 Figure.__init__(self,**kwargs)
237 self.isConfig = False
238 self.__nsubplots = 1
239 self.PLOT_CODE = VOLT_ACF_CODE
240 self.WIDTH = 900
241 self.HEIGHT = 700
242 self.counter_imagwr= 0
243 self.FTP_WEI = None
244 self.EXP_CODE = None
245 self.SUB_EXP_CODE = None
246 self.PLOT_POS = None
247
248 def getSubplots(self) :
249 ncol = 1
250 nrow = 1
251 return nrow, ncol
252
253 def setup(self,id, nplots,wintitle,show):
254 self.nplots = nplots
255 ncolspan = 1
256 colspan = 1
257 self.createFigure(id=id,
258 wintitle = wintitle,
259 widthplot = self.WIDTH,
260 heightplot = self.HEIGHT,
261 show = show)
262 nrow, ncol = self.getSubplots()
263 counter = 0
264 for y in range(nrow):
265 for x in range(ncol):
266 self.addAxes(nrow, ncol*ncolspan,y, x*ncolspan, colspan, 1)
267
268 def run(self,dataOut, id, wintitle="",channelList = None , channel =None, nSamples = None,
269 nSampleList= None, resolutionFactor=None, xmin= None, xmax = None, ymin=None, ymax=None,
270 save= False, figpath='./', figfile= None,show= True, ftp= False, wr_period=1, server= None,
271 folder= None, username =None, password= None, ftp_wei=0 , exp_code=0,sub_exp_code=0,plot_pos=0,
272 xaxis="time"):
273
274 channel0 = channel
275 nSamples = nSamples
276 resFactor = resolutionFactor
277
278 if nSamples == None:
279 nSamples = 20
280
281 if resFactor == None:
282 resFactor = 5
283
284 if channel0 == None:
285 channel0 = 0
286 else:
287 if channel0 not in dataOut.channelList:
288 raise ValueError, "Channel %d is not in %s dataOut.channelList"%(channel0, dataOut.channelList)
289
290 if channelList == None:
291 channelIndexList = dataOut.channelIndexList
292 channelList = dataOut.channelList
293
294 else:
295 channelIndexList = []
296 for channel in channelList:
297 if channel not in dataOut.channelList:
298 raise ValueError, "Channel %d is not in dataOut.channelList"
299 channelIndexList.append(dataOut.channelList.index(channel))
300
301
302 #factor = dataOut.normFactor
303 y = dataOut.getHeiRange()
304 #print y, dataOut.heightList[0]
305 #print "pause"
306 #import time
307 #time.sleep(10)
308 deltaHeight = dataOut.heightList[1]-dataOut.heightList[0]
309 z = dataOut.data
310
311 shape = dataOut.data.shape
312 hei_index = numpy.arange(shape[2])
313 hei_plot = numpy.arange(nSamples)*resFactor
314
315 if nSampleList is not None:
316 for nsample in nSampleList:
317 if nsample not in dataOut.heightList/deltaHeight:
318 print "Lista available : %s "%(dataOut.heightList/deltaHeight)
319 raise ValueError, "nsample %d is not in %s dataOut.heightList"%(nsample,dataOut.heightList)
320
321 if nSampleList is not None:
322 hei_plot = numpy.array(nSampleList)*resFactor
323
324 if hei_plot[-1] >= hei_index[-1]:
325 print ("La cantidad de puntos en altura es %d y la resolucion es %f Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
326 raise ValueError, "resFactor %d multiplicado por el valor de %d nSamples es mayor a %d cantidad total de puntos"%(resFactor,nSamples,hei_index[-1])
327
328 #escalamiento -1 a 1 a resolucion (factor de resolucion en altura)* deltaHeight
329 #min = numpy.min(z[0,:,0])
330 #max =numpy.max(z[0,:,0])
331 for i in range(shape[0]):
332 for j in range(shape[2]):
333 min = numpy.min(z[i,:,j])
334 max = numpy.max(z[i,:,j])
335 z[i,:,j]= (((z[i,:,j]-min)/(max-min))*deltaHeight*resFactor + j*deltaHeight+dataOut.heightList[0])
336
337
338 if xaxis == "time":
339 x = dataOut.getAcfRange()*1000
340 zdB = z[channel0,:,hei_plot]
341 xlabel = "Time (ms)"
342 ylabel = "VOLT_ACF"
343
344
345 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
346 title = wintitle + "VOLT ACF Plot Ch %s %s" %(channel0,thisDatetime.strftime("%d-%b-%Y"))
347
348 if not self.isConfig:
349
350 nplots = 1
351
352 self.setup(id=id,
353 nplots=nplots,
354 wintitle=wintitle,
355 show=show)
356
357 if xmin == None: xmin = numpy.nanmin(x)#*0.9
358 if xmax == None: xmax = numpy.nanmax(x)#*1.1
359 if ymin == None: ymin = numpy.nanmin(zdB)
360 if ymax == None: ymax = numpy.nanmax(zdB)
361
362 print ("El parametro resFactor es %d y la resolucion en altura es %f"%(resFactor,deltaHeight ))
363 print ("La cantidad de puntos en altura es %d y la nueva resolucion es %f Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
364 print ("La altura maxima es %d Km"%(hei_plot[-1]*deltaHeight ))
365
366 self.FTP_WEI = ftp_wei
367 self.EXP_CODE = exp_code
368 self.SUB_EXP_CODE = sub_exp_code
369 self.PLOT_POS = plot_pos
370
371 self.isConfig = True
372
373 self.setWinTitle(title)
374
375 title = "VOLT ACF Plot: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
376 axes = self.axesList[0]
377
378 legendlabels = ["Range = %dKm" %y[i] for i in hei_plot]
379
380 axes.pmultilineyaxis( x, zdB,
381 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
382 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
383 ytick_visible=True, nxticks=5,
384 grid='x')
385
386 self.draw()
387
388 if figfile == None:
389 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
390 name = str_datetime
391 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
392 name = name + '_az' + '_%2.2f'%(dataOut.azimuth) + '_zn' + '_%2.2f'%(dataOut.zenith)
393 figfile = self.getFilename(name)
394
395 self.save(figpath=figpath,
396 figfile=figfile,
397 save=save,
398 ftp=ftp,
399 wr_period=wr_period,
400 thisDatetime=thisDatetime)
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
Show file before | Show file after | Hide comments
@@ -1,29 +1,33
1 1 '''
2 2 @author: roj-idl71
3 3 '''
4 #USED IN jroplot_voltage.py
5 VOLT_ACF_CODE =11 # Volt Autocorrelation function
6
4 7 #USED IN jroplot_spectra.py
5 8 RTI_CODE = 0 #Range time intensity (RTI).
6 9 SPEC_CODE = 1 #Spectra (and Cross-spectra) information.
7 10 CROSS_CODE = 2 #Cross-Correlation information.
8 11 COH_CODE = 3 #Coherence map.
9 12 BASE_CODE = 4 #Base lines graphic.
10 13 ROW_CODE = 5 #Row Spectra.
11 14 TOTAL_CODE = 6 #Total Power.
12 15 DRIFT_CODE = 7 #Drifts graphics.
13 16 HEIGHT_CODE = 8 #Height profile.
14 17 PHASE_CODE = 9 #Signal Phase.
15 18 ACF_CODE = 10 #Autocorrelation function.
16 19
20
17 21 POWER_CODE = 16
18 22 NOISE_CODE = 17
19 23 BEACON_CODE = 18
20 24
21 25 #USED IN jroplot_parameters.py
22 26 WIND_CODE = 22
23 27 MSKYMAP_CODE = 23
24 28 MPHASE_CODE = 24
25 29
26 30 MOMENTS_CODE = 25
27 31 PARMS_CODE = 26
28 32 SPECFIT_CODE = 27
29 33 EWDRIFT_CODE = 28
Show file before | Show file after | Hide comments
@@ -1,1666 +1,1667
1 1 import sys
2 2 import numpy
3 3 from scipy import interpolate
4 4 from schainpy import cSchain
5 5 from jroproc_base import ProcessingUnit, Operation
6 6 from schainpy.model.data.jrodata import Voltage
7 7 from time import time
8 8
9 9 import math
10 10
11 11 def rep_seq(x, rep=10):
12 12 L = len(x) * rep
13 13 res = numpy.zeros(L, dtype=x.dtype)
14 14 idx = numpy.arange(len(x)) * rep
15 15 for i in numpy.arange(rep):
16 16 res[idx + i] = x
17 17 return(res)
18 18
19 19
20 20 def create_pseudo_random_code(clen=10000, seed=0):
21 21 """
22 22 seed is a way of reproducing the random code without
23 23 having to store all actual codes. the seed can then
24 24 act as a sort of station_id.
25 25
26 26 """
27 27 numpy.random.seed(seed)
28 28 phases = numpy.array(
29 29 numpy.exp(1.0j * 2.0 * math.pi * numpy.random.random(clen)),
30 30 dtype=numpy.complex64,
31 31 )
32 32 return(phases)
33 33
34 34
35 35 def periodic_convolution_matrix(envelope, rmin=0, rmax=100):
36 36 """
37 37 we imply that the number of measurements is equal to the number of elements
38 38 in code
39 39
40 40 """
41 41 L = len(envelope)
42 42 ridx = numpy.arange(rmin, rmax)
43 43 A = numpy.zeros([L, rmax-rmin], dtype=numpy.complex64)
44 44 for i in numpy.arange(L):
45 45 A[i, :] = envelope[(i-ridx) % L]
46 46 result = {}
47 47 result['A'] = A
48 48 result['ridx'] = ridx
49 49 return(result)
50 50
51 51
52 52 B_cache = 0
53 53 r_cache = 0
54 54 B_cached = False
55 55 def create_estimation_matrix(code, rmin=0, rmax=1000, cache=True):
56 56 global B_cache
57 57 global r_cache
58 58 global B_cached
59 59
60 60 if not cache or not B_cached:
61 61 r_cache = periodic_convolution_matrix(
62 62 envelope=code, rmin=rmin, rmax=rmax,
63 63 )
64 64 A = r_cache['A']
65 65 Ah = numpy.transpose(numpy.conjugate(A))
66 66 B_cache = numpy.dot(numpy.linalg.inv(numpy.dot(Ah, A)), Ah)
67 67 r_cache['B'] = B_cache
68 68 B_cached = True
69 69 return(r_cache)
70 70 else:
71 71 return(r_cache)
72 72
73 73 class VoltageProc(ProcessingUnit):
74 74
75 75
76 76 def __init__(self, **kwargs):
77 77
78 78 ProcessingUnit.__init__(self, **kwargs)
79 79
80 80 # self.objectDict = {}
81 81 self.dataOut = Voltage()
82 82 self.flip = 1
83 83
84 84 def run(self):
85 85 if self.dataIn.type == 'AMISR':
86 86 self.__updateObjFromAmisrInput()
87 87
88 88 if self.dataIn.type == 'Voltage':
89 89 self.dataOut.copy(self.dataIn)
90 90
91 91 # self.dataOut.copy(self.dataIn)
92 92
93 93 def __updateObjFromAmisrInput(self):
94 94
95 95 self.dataOut.timeZone = self.dataIn.timeZone
96 96 self.dataOut.dstFlag = self.dataIn.dstFlag
97 97 self.dataOut.errorCount = self.dataIn.errorCount
98 98 self.dataOut.useLocalTime = self.dataIn.useLocalTime
99 99
100 100 self.dataOut.flagNoData = self.dataIn.flagNoData
101 101 self.dataOut.data = self.dataIn.data
102 102 self.dataOut.utctime = self.dataIn.utctime
103 103 self.dataOut.channelList = self.dataIn.channelList
104 104 # self.dataOut.timeInterval = self.dataIn.timeInterval
105 105 self.dataOut.heightList = self.dataIn.heightList
106 106 self.dataOut.nProfiles = self.dataIn.nProfiles
107 107
108 108 self.dataOut.nCohInt = self.dataIn.nCohInt
109 109 self.dataOut.ippSeconds = self.dataIn.ippSeconds
110 110 self.dataOut.frequency = self.dataIn.frequency
111 111
112 112 self.dataOut.azimuth = self.dataIn.azimuth
113 113 self.dataOut.zenith = self.dataIn.zenith
114 114
115 115 self.dataOut.beam.codeList = self.dataIn.beam.codeList
116 116 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
117 117 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
118 118 #
119 119 # pass#
120 120 #
121 121 # def init(self):
122 122 #
123 123 #
124 124 # if self.dataIn.type == 'AMISR':
125 125 # self.__updateObjFromAmisrInput()
126 126 #
127 127 # if self.dataIn.type == 'Voltage':
128 128 # self.dataOut.copy(self.dataIn)
129 129 # # No necesita copiar en cada init() los atributos de dataIn
130 130 # # la copia deberia hacerse por cada nuevo bloque de datos
131 131
132 132 def selectChannels(self, channelList):
133 133
134 134 channelIndexList = []
135 135
136 136 for channel in channelList:
137 137 if channel not in self.dataOut.channelList:
138 138 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
139 139
140 140 index = self.dataOut.channelList.index(channel)
141 141 channelIndexList.append(index)
142 142
143 143 self.selectChannelsByIndex(channelIndexList)
144 144
145 145 def selectChannelsByIndex(self, channelIndexList):
146 146 """
147 147 Selecciona un bloque de datos en base a canales segun el channelIndexList
148 148
149 149 Input:
150 150 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
151 151
152 152 Affected:
153 153 self.dataOut.data
154 154 self.dataOut.channelIndexList
155 155 self.dataOut.nChannels
156 156 self.dataOut.m_ProcessingHeader.totalSpectra
157 157 self.dataOut.systemHeaderObj.numChannels
158 158 self.dataOut.m_ProcessingHeader.blockSize
159 159
160 160 Return:
161 161 None
162 162 """
163 163
164 164 for channelIndex in channelIndexList:
165 165 if channelIndex not in self.dataOut.channelIndexList:
166 166 print channelIndexList
167 167 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
168 168
169 169 if self.dataOut.flagDataAsBlock:
170 170 """
171 171 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
172 172 """
173 173 data = self.dataOut.data[channelIndexList,:,:]
174 174 else:
175 175 data = self.dataOut.data[channelIndexList,:]
176 176
177 177 self.dataOut.data = data
178 178 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
179 179 # self.dataOut.nChannels = nChannels
180 180
181 181 return 1
182 182
183 183 def selectHeights(self, minHei=None, maxHei=None):
184 184 """
185 185 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
186 186 minHei <= height <= maxHei
187 187
188 188 Input:
189 189 minHei : valor minimo de altura a considerar
190 190 maxHei : valor maximo de altura a considerar
191 191
192 192 Affected:
193 193 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
194 194
195 195 Return:
196 196 1 si el metodo se ejecuto con exito caso contrario devuelve 0
197 197 """
198 198
199 199 if minHei == None:
200 200 minHei = self.dataOut.heightList[0]
201 201
202 202 if maxHei == None:
203 203 maxHei = self.dataOut.heightList[-1]
204 204
205 205 if (minHei < self.dataOut.heightList[0]):
206 206 minHei = self.dataOut.heightList[0]
207 207
208 208 if (maxHei > self.dataOut.heightList[-1]):
209 209 maxHei = self.dataOut.heightList[-1]
210 210
211 211 minIndex = 0
212 212 maxIndex = 0
213 213 heights = self.dataOut.heightList
214 214
215 215 inda = numpy.where(heights >= minHei)
216 216 indb = numpy.where(heights <= maxHei)
217 217
218 218 try:
219 219 minIndex = inda[0][0]
220 220 except:
221 221 minIndex = 0
222 222
223 223 try:
224 224 maxIndex = indb[0][-1]
225 225 except:
226 226 maxIndex = len(heights)
227 227
228 228 self.selectHeightsByIndex(minIndex, maxIndex)
229 229
230 230 return 1
231 231
232 232
233 233 def selectHeightsByIndex(self, minIndex, maxIndex):
234 234 """
235 235 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
236 236 minIndex <= index <= maxIndex
237 237
238 238 Input:
239 239 minIndex : valor de indice minimo de altura a considerar
240 240 maxIndex : valor de indice maximo de altura a considerar
241 241
242 242 Affected:
243 243 self.dataOut.data
244 244 self.dataOut.heightList
245 245
246 246 Return:
247 247 1 si el metodo se ejecuto con exito caso contrario devuelve 0
248 248 """
249 249
250 250 if (minIndex < 0) or (minIndex > maxIndex):
251 251 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
252 252
253 253 if (maxIndex >= self.dataOut.nHeights):
254 254 maxIndex = self.dataOut.nHeights
255 255
256 256 #voltage
257 257 if self.dataOut.flagDataAsBlock:
258 258 """
259 259 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
260 260 """
261 261 data = self.dataOut.data[:,:, minIndex:maxIndex]
262 262 else:
263 263 data = self.dataOut.data[:, minIndex:maxIndex]
264 264
265 265 # firstHeight = self.dataOut.heightList[minIndex]
266 266
267 267 self.dataOut.data = data
268 268 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
269 269
270 270 if self.dataOut.nHeights <= 1:
271 271 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
272 272
273 273 return 1
274 274
275 275
276 276 def filterByHeights(self, window):
277 277
278 278 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
279 279
280 280 if window == None:
281 281 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
282 282
283 283 newdelta = deltaHeight * window
284 284 r = self.dataOut.nHeights % window
285 285 newheights = (self.dataOut.nHeights-r)/window
286 286
287 287 if newheights <= 1:
288 288 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
289 289
290 290 if self.dataOut.flagDataAsBlock:
291 291 """
292 292 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
293 293 """
294 294 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
295 295 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
296 296 buffer = numpy.sum(buffer,3)
297 297
298 298 else:
299 299 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
300 300 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
301 301 buffer = numpy.sum(buffer,2)
302 302
303 303 self.dataOut.data = buffer
304 304 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
305 305 self.dataOut.windowOfFilter = window
306 306
307 307 def setH0(self, h0, deltaHeight = None):
308 308
309 309 if not deltaHeight:
310 310 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
311 311
312 312 nHeights = self.dataOut.nHeights
313 313
314 314 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
315 315
316 316 self.dataOut.heightList = newHeiRange
317 317
318 318 def deFlip(self, channelList = []):
319 319
320 320 data = self.dataOut.data.copy()
321 321
322 322 if self.dataOut.flagDataAsBlock:
323 323 flip = self.flip
324 324 profileList = range(self.dataOut.nProfiles)
325 325
326 326 if not channelList:
327 327 for thisProfile in profileList:
328 328 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
329 329 flip *= -1.0
330 330 else:
331 331 for thisChannel in channelList:
332 332 if thisChannel not in self.dataOut.channelList:
333 333 continue
334 334
335 335 for thisProfile in profileList:
336 336 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
337 337 flip *= -1.0
338 338
339 339 self.flip = flip
340 340
341 341 else:
342 342 if not channelList:
343 343 data[:,:] = data[:,:]*self.flip
344 344 else:
345 345 for thisChannel in channelList:
346 346 if thisChannel not in self.dataOut.channelList:
347 347 continue
348 348
349 349 data[thisChannel,:] = data[thisChannel,:]*self.flip
350 350
351 351 self.flip *= -1.
352 352
353 353 self.dataOut.data = data
354 354
355 355 def setRadarFrequency(self, frequency=None):
356 356
357 357 if frequency != None:
358 358 self.dataOut.frequency = frequency
359 359
360 360 return 1
361 361
362 362 def interpolateHeights(self, topLim, botLim):
363 363 #69 al 72 para julia
364 364 #82-84 para meteoros
365 365 if len(numpy.shape(self.dataOut.data))==2:
366 366 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
367 367 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
368 368 #self.dataOut.data[:,botLim:limSup+1] = sampInterp
369 369 self.dataOut.data[:,botLim:topLim+1] = sampInterp
370 370 else:
371 371 nHeights = self.dataOut.data.shape[2]
372 372 x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
373 373 y = self.dataOut.data[:,:,range(botLim)+range(topLim+1,nHeights)]
374 374 f = interpolate.interp1d(x, y, axis = 2)
375 375 xnew = numpy.arange(botLim,topLim+1)
376 376 ynew = f(xnew)
377 377
378 378 self.dataOut.data[:,:,botLim:topLim+1] = ynew
379 379
380 380 # import collections
381 381
382 382 class CohInt(Operation):
383 383
384 384 isConfig = False
385 385 __profIndex = 0
386 386 __byTime = False
387 387 __initime = None
388 388 __lastdatatime = None
389 389 __integrationtime = None
390 390 __buffer = None
391 391 __bufferStride = []
392 392 __dataReady = False
393 393 __profIndexStride = 0
394 394 __dataToPutStride = False
395 395 n = None
396 396
397 397 def __init__(self, **kwargs):
398 398
399 399 Operation.__init__(self, **kwargs)
400 400
401 401 # self.isConfig = False
402 402
403 403 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False):
404 404 """
405 405 Set the parameters of the integration class.
406 406
407 407 Inputs:
408 408
409 409 n : Number of coherent integrations
410 410 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
411 411 overlapping :
412 412 """
413 413
414 414 self.__initime = None
415 415 self.__lastdatatime = 0
416 416 self.__buffer = None
417 417 self.__dataReady = False
418 418 self.byblock = byblock
419 419 self.stride = stride
420 420
421 421 if n == None and timeInterval == None:
422 422 raise ValueError, "n or timeInterval should be specified ..."
423 423
424 424 if n != None:
425 425 self.n = n
426 426 self.__byTime = False
427 427 else:
428 428 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
429 429 self.n = 9999
430 430 self.__byTime = True
431 431
432 432 if overlapping:
433 433 self.__withOverlapping = True
434 434 self.__buffer = None
435 435 else:
436 436 self.__withOverlapping = False
437 437 self.__buffer = 0
438 438
439 439 self.__profIndex = 0
440 440
441 441 def putData(self, data):
442 442
443 443 """
444 444 Add a profile to the __buffer and increase in one the __profileIndex
445 445
446 446 """
447 447
448 448 if not self.__withOverlapping:
449 449 self.__buffer += data.copy()
450 450 self.__profIndex += 1
451 451 return
452 452
453 453 #Overlapping data
454 454 nChannels, nHeis = data.shape
455 455 data = numpy.reshape(data, (1, nChannels, nHeis))
456 456
457 457 #If the buffer is empty then it takes the data value
458 458 if self.__buffer is None:
459 459 self.__buffer = data
460 460 self.__profIndex += 1
461 461 return
462 462
463 463 #If the buffer length is lower than n then stakcing the data value
464 464 if self.__profIndex < self.n:
465 465 self.__buffer = numpy.vstack((self.__buffer, data))
466 466 self.__profIndex += 1
467 467 return
468 468
469 469 #If the buffer length is equal to n then replacing the last buffer value with the data value
470 470 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
471 471 self.__buffer[self.n-1] = data
472 472 self.__profIndex = self.n
473 473 return
474 474
475 475
476 476 def pushData(self):
477 477 """
478 478 Return the sum of the last profiles and the profiles used in the sum.
479 479
480 480 Affected:
481 481
482 482 self.__profileIndex
483 483
484 484 """
485 485
486 486 if not self.__withOverlapping:
487 487 data = self.__buffer
488 488 n = self.__profIndex
489 489
490 490 self.__buffer = 0
491 491 self.__profIndex = 0
492 492
493 493 return data, n
494 494
495 495 #Integration with Overlapping
496 496 data = numpy.sum(self.__buffer, axis=0)
497 497 # print data
498 498 # raise
499 499 n = self.__profIndex
500 500
501 501 return data, n
502 502
503 503 def byProfiles(self, data):
504 504
505 505 self.__dataReady = False
506 506 avgdata = None
507 507 # n = None
508 508 # print data
509 509 # raise
510 510 self.putData(data)
511 511
512 512 if self.__profIndex == self.n:
513 513 avgdata, n = self.pushData()
514 514 self.__dataReady = True
515 515
516 516 return avgdata
517 517
518 518 def byTime(self, data, datatime):
519 519
520 520 self.__dataReady = False
521 521 avgdata = None
522 522 n = None
523 523
524 524 self.putData(data)
525 525
526 526 if (datatime - self.__initime) >= self.__integrationtime:
527 527 avgdata, n = self.pushData()
528 528 self.n = n
529 529 self.__dataReady = True
530 530
531 531 return avgdata
532 532
533 533 def integrateByStride(self, data, datatime):
534 534 # print data
535 535 if self.__profIndex == 0:
536 536 self.__buffer = [[data.copy(), datatime]]
537 537 else:
538 538 self.__buffer.append([data.copy(),datatime])
539 539 self.__profIndex += 1
540 540 self.__dataReady = False
541 541
542 542 if self.__profIndex == self.n * self.stride :
543 543 self.__dataToPutStride = True
544 544 self.__profIndexStride = 0
545 545 self.__profIndex = 0
546 546 self.__bufferStride = []
547 547 for i in range(self.stride):
548 548 current = self.__buffer[i::self.stride]
549 549 data = numpy.sum([t[0] for t in current], axis=0)
550 550 avgdatatime = numpy.average([t[1] for t in current])
551 551 # print data
552 552 self.__bufferStride.append((data, avgdatatime))
553 553
554 554 if self.__dataToPutStride:
555 555 self.__dataReady = True
556 556 self.__profIndexStride += 1
557 557 if self.__profIndexStride == self.stride:
558 558 self.__dataToPutStride = False
559 559 # print self.__bufferStride[self.__profIndexStride - 1]
560 560 # raise
561 561 return self.__bufferStride[self.__profIndexStride - 1]
562 562
563 563
564 564 return None, None
565 565
566 566 def integrate(self, data, datatime=None):
567 567
568 568 if self.__initime == None:
569 569 self.__initime = datatime
570 570
571 571 if self.__byTime:
572 572 avgdata = self.byTime(data, datatime)
573 573 else:
574 574 avgdata = self.byProfiles(data)
575 575
576 576
577 577 self.__lastdatatime = datatime
578 578
579 579 if avgdata is None:
580 580 return None, None
581 581
582 582 avgdatatime = self.__initime
583 583
584 584 deltatime = datatime - self.__lastdatatime
585 585
586 586 if not self.__withOverlapping:
587 587 self.__initime = datatime
588 588 else:
589 589 self.__initime += deltatime
590 590
591 591 return avgdata, avgdatatime
592 592
593 593 def integrateByBlock(self, dataOut):
594 594
595 595 times = int(dataOut.data.shape[1]/self.n)
596 596 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
597 597
598 598 id_min = 0
599 599 id_max = self.n
600 600
601 601 for i in range(times):
602 602 junk = dataOut.data[:,id_min:id_max,:]
603 603 avgdata[:,i,:] = junk.sum(axis=1)
604 604 id_min += self.n
605 605 id_max += self.n
606 606
607 607 timeInterval = dataOut.ippSeconds*self.n
608 608 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
609 609 self.__dataReady = True
610 610 return avgdata, avgdatatime
611 611
612 612 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
613 613 if not self.isConfig:
614 614 self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
615 615 self.isConfig = True
616 616
617 617 if dataOut.flagDataAsBlock:
618 618 """
619 619 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
620 620 """
621 621 avgdata, avgdatatime = self.integrateByBlock(dataOut)
622 622 dataOut.nProfiles /= self.n
623 623 else:
624 624 if stride is None:
625 625 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
626 626 else:
627 627 avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
628 628
629 629
630 630 # dataOut.timeInterval *= n
631 631 dataOut.flagNoData = True
632 632
633 633 if self.__dataReady:
634 634 dataOut.data = avgdata
635 635 dataOut.nCohInt *= self.n
636 636 dataOut.utctime = avgdatatime
637 637 # print avgdata, avgdatatime
638 638 # raise
639 639 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
640 640 dataOut.flagNoData = False
641 641
642 642 class Decoder(Operation):
643 643
644 644 isConfig = False
645 645 __profIndex = 0
646 646
647 647 code = None
648 648
649 649 nCode = None
650 650 nBaud = None
651 651
652 652 def __init__(self, **kwargs):
653 653
654 654 Operation.__init__(self, **kwargs)
655 655
656 656 self.times = None
657 657 self.osamp = None
658 658 # self.__setValues = False
659 659 self.isConfig = False
660 660
661 661 def setup(self, code, osamp, dataOut):
662 662
663 663 self.__profIndex = 0
664 664
665 665 self.code = code
666 666
667 667 self.nCode = len(code)
668 668 self.nBaud = len(code[0])
669 669
670 670 if (osamp != None) and (osamp >1):
671 671 self.osamp = osamp
672 672 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
673 673 self.nBaud = self.nBaud*self.osamp
674 674
675 675 self.__nChannels = dataOut.nChannels
676 676 self.__nProfiles = dataOut.nProfiles
677 677 self.__nHeis = dataOut.nHeights
678 678
679 679 if self.__nHeis < self.nBaud:
680 680 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
681 681
682 682 #Frequency
683 683 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
684 684
685 685 __codeBuffer[:,0:self.nBaud] = self.code
686 686
687 687 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
688 688
689 689 if dataOut.flagDataAsBlock:
690 690
691 691 self.ndatadec = self.__nHeis #- self.nBaud + 1
692 692
693 693 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
694 694
695 695 else:
696 696
697 697 #Time
698 698 self.ndatadec = self.__nHeis #- self.nBaud + 1
699 699
700 700 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
701 701
702 702 def __convolutionInFreq(self, data):
703 703
704 704 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
705 705
706 706 fft_data = numpy.fft.fft(data, axis=1)
707 707
708 708 conv = fft_data*fft_code
709 709
710 710 data = numpy.fft.ifft(conv,axis=1)
711 711
712 712 return data
713 713
714 714 def __convolutionInFreqOpt(self, data):
715 715
716 716 raise NotImplementedError
717 717
718 718 def __convolutionInTime(self, data):
719 719
720 720 code = self.code[self.__profIndex]
721 721 for i in range(self.__nChannels):
722 722 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
723 723
724 724 return self.datadecTime
725 725
726 726 def __convolutionByBlockInTime(self, data):
727 727
728 728 repetitions = self.__nProfiles / self.nCode
729 729
730 730 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
731 731 junk = junk.flatten()
732 732 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
733 733 profilesList = xrange(self.__nProfiles)
734 734
735 735 for i in range(self.__nChannels):
736 736 for j in profilesList:
737 737 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
738 738 return self.datadecTime
739 739
740 740 def __convolutionByBlockInFreq(self, data):
741 741
742 742 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
743 743
744 744
745 745 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
746 746
747 747 fft_data = numpy.fft.fft(data, axis=2)
748 748
749 749 conv = fft_data*fft_code
750 750
751 751 data = numpy.fft.ifft(conv,axis=2)
752 752
753 753 return data
754 754
755 755
756 756 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
757 757
758 758 if dataOut.flagDecodeData:
759 759 print "This data is already decoded, recoding again ..."
760 760
761 761 if not self.isConfig:
762 762
763 763 if code is None:
764 764 if dataOut.code is None:
765 765 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
766 766
767 767 code = dataOut.code
768 768 else:
769 769 code = numpy.array(code).reshape(nCode,nBaud)
770 770 self.setup(code, osamp, dataOut)
771 771
772 772 self.isConfig = True
773 773
774 774 if mode == 3:
775 775 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
776 776
777 777 if times != None:
778 778 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
779 779
780 780 if self.code is None:
781 781 print "Fail decoding: Code is not defined."
782 782 return
783 783
784 784 self.__nProfiles = dataOut.nProfiles
785 785 datadec = None
786 786
787 787 if mode == 3:
788 788 mode = 0
789 789
790 790 if dataOut.flagDataAsBlock:
791 791 """
792 792 Decoding when data have been read as block,
793 793 """
794 794
795 795 if mode == 0:
796 796 datadec = self.__convolutionByBlockInTime(dataOut.data)
797 797 if mode == 1:
798 798 datadec = self.__convolutionByBlockInFreq(dataOut.data)
799 799 else:
800 800 """
801 801 Decoding when data have been read profile by profile
802 802 """
803 803 if mode == 0:
804 804 datadec = self.__convolutionInTime(dataOut.data)
805 805
806 806 if mode == 1:
807 807 datadec = self.__convolutionInFreq(dataOut.data)
808 808
809 809 if mode == 2:
810 810 datadec = self.__convolutionInFreqOpt(dataOut.data)
811 811
812 812 if datadec is None:
813 813 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
814 814
815 815 dataOut.code = self.code
816 816 dataOut.nCode = self.nCode
817 817 dataOut.nBaud = self.nBaud
818 818
819 819 dataOut.data = datadec
820 820
821 821 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
822 822
823 823 dataOut.flagDecodeData = True #asumo q la data esta decodificada
824 824
825 825 if self.__profIndex == self.nCode-1:
826 826 self.__profIndex = 0
827 827 return 1
828 828
829 829 self.__profIndex += 1
830 830
831 831 return 1
832 832 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
833 833
834 834
835 835 class ProfileConcat(Operation):
836 836
837 837 isConfig = False
838 838 buffer = None
839 839 concat_m =None
840 840
841 841 def __init__(self, **kwargs):
842 842
843 843 Operation.__init__(self, **kwargs)
844 844 self.profileIndex = 0
845 845
846 846 def reset(self):
847 847 self.buffer = numpy.zeros_like(self.buffer)
848 848 self.start_index = 0
849 849 self.times = 1
850 850
851 851 def setup(self, data, m, n=1):
852 852 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
853 853 self.nHeights = data.shape[1]#.nHeights
854 854 self.start_index = 0
855 855 self.times = 1
856 856
857 857 def concat(self, data):
858 858
859 859 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
860 860 self.start_index = self.start_index + self.nHeights
861 861
862 862 def run(self, dataOut, m):
863 863
864 864 self.concat_m= m
865 865 dataOut.flagNoData = True
866 866
867 867 if not self.isConfig:
868 868 self.setup(dataOut.data, m, 1)
869 869 self.isConfig = True
870 870
871 871 if dataOut.flagDataAsBlock:
872 872 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
873 873
874 874 else:
875 875 self.concat(dataOut.data)
876 876 self.times += 1
877 877 if self.times > m:
878 878 dataOut.data = self.buffer
879 879 self.reset()
880 880 dataOut.flagNoData = False
881 881 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
882 882 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
883 883 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
884 884 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
885 885 dataOut.ippSeconds *= m
886 886 dataOut.concat_m = int(m)
887 887
888 888 class ProfileSelector(Operation):
889 889
890 890 profileIndex = None
891 891 # Tamanho total de los perfiles
892 892 nProfiles = None
893 893
894 894 def __init__(self, **kwargs):
895 895
896 896 Operation.__init__(self, **kwargs)
897 897 self.profileIndex = 0
898 898
899 899 def incProfileIndex(self):
900 900
901 901 self.profileIndex += 1
902 902
903 903 if self.profileIndex >= self.nProfiles:
904 904 self.profileIndex = 0
905 905
906 906 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
907 907
908 908 if profileIndex < minIndex:
909 909 return False
910 910
911 911 if profileIndex > maxIndex:
912 912 return False
913 913
914 914 return True
915 915
916 916 def isThisProfileInList(self, profileIndex, profileList):
917 917
918 918 if profileIndex not in profileList:
919 919 return False
920 920
921 921 return True
922 922
923 923 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
924 924
925 925 """
926 926 ProfileSelector:
927 927
928 928 Inputs:
929 929 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
930 930
931 931 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
932 932
933 933 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
934 934
935 935 """
936 936
937 937 if rangeList is not None:
938 938 if type(rangeList[0]) not in (tuple, list):
939 939 rangeList = [rangeList]
940 940
941 941 dataOut.flagNoData = True
942 942
943 943 if dataOut.flagDataAsBlock:
944 944 """
945 945 data dimension = [nChannels, nProfiles, nHeis]
946 946 """
947 947 if profileList != None:
948 948 dataOut.data = dataOut.data[:,profileList,:]
949 949
950 950 if profileRangeList != None:
951 951 minIndex = profileRangeList[0]
952 952 maxIndex = profileRangeList[1]
953 953 profileList = range(minIndex, maxIndex+1)
954 954
955 955 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
956 956
957 957 if rangeList != None:
958 958
959 959 profileList = []
960 960
961 961 for thisRange in rangeList:
962 962 minIndex = thisRange[0]
963 963 maxIndex = thisRange[1]
964 964
965 965 profileList.extend(range(minIndex, maxIndex+1))
966 966
967 967 dataOut.data = dataOut.data[:,profileList,:]
968 968
969 969 dataOut.nProfiles = len(profileList)
970 970 dataOut.profileIndex = dataOut.nProfiles - 1
971 971 dataOut.flagNoData = False
972 972
973 973 return True
974 974
975 975 """
976 976 data dimension = [nChannels, nHeis]
977 977 """
978 978
979 979 if profileList != None:
980 980
981 981 if self.isThisProfileInList(dataOut.profileIndex, profileList):
982 982
983 983 self.nProfiles = len(profileList)
984 984 dataOut.nProfiles = self.nProfiles
985 985 dataOut.profileIndex = self.profileIndex
986 986 dataOut.flagNoData = False
987 987
988 988 self.incProfileIndex()
989 989 return True
990 990
991 991 if profileRangeList != None:
992 992
993 993 minIndex = profileRangeList[0]
994 994 maxIndex = profileRangeList[1]
995 995
996 996 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
997 997
998 998 self.nProfiles = maxIndex - minIndex + 1
999 999 dataOut.nProfiles = self.nProfiles
1000 1000 dataOut.profileIndex = self.profileIndex
1001 1001 dataOut.flagNoData = False
1002 1002
1003 1003 self.incProfileIndex()
1004 1004 return True
1005 1005
1006 1006 if rangeList != None:
1007 1007
1008 1008 nProfiles = 0
1009 1009
1010 1010 for thisRange in rangeList:
1011 1011 minIndex = thisRange[0]
1012 1012 maxIndex = thisRange[1]
1013 1013
1014 1014 nProfiles += maxIndex - minIndex + 1
1015 1015
1016 1016 for thisRange in rangeList:
1017 1017
1018 1018 minIndex = thisRange[0]
1019 1019 maxIndex = thisRange[1]
1020 1020
1021 1021 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
1022 1022
1023 1023 self.nProfiles = nProfiles
1024 1024 dataOut.nProfiles = self.nProfiles
1025 1025 dataOut.profileIndex = self.profileIndex
1026 1026 dataOut.flagNoData = False
1027 1027
1028 1028 self.incProfileIndex()
1029 1029
1030 1030 break
1031 1031
1032 1032 return True
1033 1033
1034 1034
1035 1035 if beam != None: #beam is only for AMISR data
1036 1036 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
1037 1037 dataOut.flagNoData = False
1038 1038 dataOut.profileIndex = self.profileIndex
1039 1039
1040 1040 self.incProfileIndex()
1041 1041
1042 1042 return True
1043 1043
1044 1044 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
1045 1045
1046 1046 return False
1047 1047
1048 1048 class Reshaper(Operation):
1049 1049
1050 1050 def __init__(self, **kwargs):
1051 1051
1052 1052 Operation.__init__(self, **kwargs)
1053 1053
1054 1054 self.__buffer = None
1055 1055 self.__nitems = 0
1056 1056
1057 1057 def __appendProfile(self, dataOut, nTxs):
1058 1058
1059 1059 if self.__buffer is None:
1060 1060 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
1061 1061 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
1062 1062
1063 1063 ini = dataOut.nHeights * self.__nitems
1064 1064 end = ini + dataOut.nHeights
1065 1065
1066 1066 self.__buffer[:, ini:end] = dataOut.data
1067 1067
1068 1068 self.__nitems += 1
1069 1069
1070 1070 return int(self.__nitems*nTxs)
1071 1071
1072 1072 def __getBuffer(self):
1073 1073
1074 1074 if self.__nitems == int(1./self.__nTxs):
1075 1075
1076 1076 self.__nitems = 0
1077 1077
1078 1078 return self.__buffer.copy()
1079 1079
1080 1080 return None
1081 1081
1082 1082 def __checkInputs(self, dataOut, shape, nTxs):
1083 1083
1084 1084 if shape is None and nTxs is None:
1085 1085 raise ValueError, "Reshaper: shape of factor should be defined"
1086 1086
1087 1087 if nTxs:
1088 1088 if nTxs < 0:
1089 1089 raise ValueError, "nTxs should be greater than 0"
1090 1090
1091 1091 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
1092 1092 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
1093 1093
1094 1094 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
1095 1095
1096 1096 return shape, nTxs
1097 1097
1098 1098 if len(shape) != 2 and len(shape) != 3:
1099 1099 raise ValueError, "shape dimension should be equal to 2 or 3. shape = (nProfiles, nHeis) or (nChannels, nProfiles, nHeis). Actually shape = (%d, %d, %d)" %(dataOut.nChannels, dataOut.nProfiles, dataOut.nHeights)
1100 1100
1101 1101 if len(shape) == 2:
1102 1102 shape_tuple = [dataOut.nChannels]
1103 1103 shape_tuple.extend(shape)
1104 1104 else:
1105 1105 shape_tuple = list(shape)
1106 1106
1107 1107 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
1108 1108
1109 1109 return shape_tuple, nTxs
1110 1110
1111 1111 def run(self, dataOut, shape=None, nTxs=None):
1112 1112
1113 1113 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
1114 1114
1115 1115 dataOut.flagNoData = True
1116 1116 profileIndex = None
1117 1117
1118 1118 if dataOut.flagDataAsBlock:
1119 1119
1120 1120 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
1121 1121 dataOut.flagNoData = False
1122 1122
1123 1123 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1124 1124
1125 1125 else:
1126 1126
1127 1127 if self.__nTxs < 1:
1128 1128
1129 1129 self.__appendProfile(dataOut, self.__nTxs)
1130 1130 new_data = self.__getBuffer()
1131 1131
1132 1132 if new_data is not None:
1133 1133 dataOut.data = new_data
1134 1134 dataOut.flagNoData = False
1135 1135
1136 1136 profileIndex = dataOut.profileIndex*nTxs
1137 1137
1138 1138 else:
1139 1139 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1140 1140
1141 1141 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1142 1142
1143 1143 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1144 1144
1145 1145 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1146 1146
1147 1147 dataOut.profileIndex = profileIndex
1148 1148
1149 1149 dataOut.ippSeconds /= self.__nTxs
1150 1150
1151 1151 class SplitProfiles(Operation):
1152 1152
1153 1153 def __init__(self, **kwargs):
1154 1154
1155 1155 Operation.__init__(self, **kwargs)
1156 1156
1157 1157 def run(self, dataOut, n):
1158 1158
1159 1159 dataOut.flagNoData = True
1160 1160 profileIndex = None
1161 1161
1162 1162 if dataOut.flagDataAsBlock:
1163 1163
1164 1164 #nchannels, nprofiles, nsamples
1165 1165 shape = dataOut.data.shape
1166 1166
1167 1167 if shape[2] % n != 0:
1168 1168 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1169 1169
1170 1170 new_shape = shape[0], shape[1]*n, shape[2]/n
1171 1171
1172 1172 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1173 1173 dataOut.flagNoData = False
1174 1174
1175 1175 profileIndex = int(dataOut.nProfiles/n) - 1
1176 1176
1177 1177 else:
1178 1178
1179 1179 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1180 1180
1181 1181 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1182 1182
1183 1183 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1184 1184
1185 1185 dataOut.nProfiles = int(dataOut.nProfiles*n)
1186 1186
1187 1187 dataOut.profileIndex = profileIndex
1188 1188
1189 1189 dataOut.ippSeconds /= n
1190 1190
1191 1191 class CombineProfiles(Operation):
1192 1192
1193 1193 def __init__(self, **kwargs):
1194 1194
1195 1195 Operation.__init__(self, **kwargs)
1196 1196
1197 1197 self.__remData = None
1198 1198 self.__profileIndex = 0
1199 1199
1200 1200 def run(self, dataOut, n):
1201 1201
1202 1202 dataOut.flagNoData = True
1203 1203 profileIndex = None
1204 1204
1205 1205 if dataOut.flagDataAsBlock:
1206 1206
1207 1207 #nchannels, nprofiles, nsamples
1208 1208 shape = dataOut.data.shape
1209 1209 new_shape = shape[0], shape[1]/n, shape[2]*n
1210 1210
1211 1211 if shape[1] % n != 0:
1212 1212 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1213 1213
1214 1214 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1215 1215 dataOut.flagNoData = False
1216 1216
1217 1217 profileIndex = int(dataOut.nProfiles*n) - 1
1218 1218
1219 1219 else:
1220 1220
1221 1221 #nchannels, nsamples
1222 1222 if self.__remData is None:
1223 1223 newData = dataOut.data
1224 1224 else:
1225 1225 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1226 1226
1227 1227 self.__profileIndex += 1
1228 1228
1229 1229 if self.__profileIndex < n:
1230 1230 self.__remData = newData
1231 1231 #continue
1232 1232 return
1233 1233
1234 1234 self.__profileIndex = 0
1235 1235 self.__remData = None
1236 1236
1237 1237 dataOut.data = newData
1238 1238 dataOut.flagNoData = False
1239 1239
1240 1240 profileIndex = dataOut.profileIndex/n
1241 1241
1242 1242
1243 1243 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1244 1244
1245 1245 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1246 1246
1247 1247 dataOut.nProfiles = int(dataOut.nProfiles/n)
1248 1248
1249 1249 dataOut.profileIndex = profileIndex
1250 1250
1251 1251 dataOut.ippSeconds *= n
1252 1252
1253 1253
1254 1254 class SSheightProfiles(Operation):
1255 1255
1256 1256 step = None
1257 1257 nsamples = None
1258 1258 bufferShape = None
1259 1259 profileShape = None
1260 1260 sshProfiles = None
1261 1261 profileIndex = None
1262 1262
1263 1263 def __init__(self, **kwargs):
1264 1264
1265 1265 Operation.__init__(self, **kwargs)
1266 1266 self.isConfig = False
1267 1267
1268 1268 def setup(self,dataOut ,step = None , nsamples = None):
1269 1269
1270 1270 if step == None and nsamples == None:
1271 1271 raise ValueError, "step or nheights should be specified ..."
1272 1272
1273 1273 self.step = step
1274 1274 self.nsamples = nsamples
1275 1275 self.__nChannels = dataOut.nChannels
1276 1276 self.__nProfiles = dataOut.nProfiles
1277 1277 self.__nHeis = dataOut.nHeights
1278 1278 shape = dataOut.data.shape #nchannels, nprofiles, nsamples
1279 1279 print "input nChannels",self.__nChannels
1280 1280 print "input nProfiles",self.__nProfiles
1281 1281 print "input nHeis",self.__nHeis
1282 1282 print "input Shape",shape
1283 1283
1284 1284
1285 1285
1286 1286 residue = (shape[1] - self.nsamples) % self.step
1287 1287 if residue != 0:
1288 1288 print "The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue)
1289 1289
1290 1290 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1291 1291 numberProfile = self.nsamples
1292 1292 numberSamples = (shape[1] - self.nsamples)/self.step
1293 1293
1294 1294 print "new numberProfile",numberProfile
1295 1295 print "new numberSamples",numberSamples
1296 1296
1297 1297 print "New number of profile: %d, number of height: %d, Resolution %f Km"%(numberProfile,numberSamples,deltaHeight*self.step)
1298 1298
1299 1299 self.bufferShape = shape[0], numberSamples, numberProfile # nchannels, nsamples , nprofiles
1300 1300 self.profileShape = shape[0], numberProfile, numberSamples # nchannels, nprofiles, nsamples
1301 1301
1302 1302 self.buffer = numpy.zeros(self.bufferShape , dtype=numpy.complex)
1303 1303 self.sshProfiles = numpy.zeros(self.profileShape, dtype=numpy.complex)
1304 1304
1305 1305 def run(self, dataOut, step, nsamples):
1306 1306
1307 1307 dataOut.flagNoData = True
1308 1308 dataOut.flagDataAsBlock = False
1309 1309 profileIndex = None
1310 1310
1311 1311
1312 1312 if not self.isConfig:
1313 1313 self.setup(dataOut, step=step , nsamples=nsamples)
1314 1314 self.isConfig = True
1315 1315
1316 1316 for i in range(self.buffer.shape[1]):
1317 1317 #self.buffer[:,i] = numpy.flip(dataOut.data[:,i*self.step:i*self.step + self.nsamples])
1318 1318 self.buffer[:,i] = dataOut.data[:,i*self.step:i*self.step + self.nsamples]
1319 1319 #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
1320 1320
1321 1321 for j in range(self.buffer.shape[0]):
1322 1322 self.sshProfiles[j] = numpy.transpose(self.buffer[j])
1323 1323
1324 1324 profileIndex = self.nsamples
1325 1325 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1326 1326 ippSeconds = (deltaHeight*1.0e-6)/(0.15)
1327 1327 #print "ippSeconds",ippSeconds
1328 1328 try:
1329 1329 if dataOut.concat_m is not None:
1330 1330 ippSeconds= ippSeconds/float(dataOut.concat_m)
1331 1331 #print "Profile concat %d"%dataOut.concat_m
1332 1332 except:
1333 1333 pass
1334 1334
1335 1335 dataOut.data = self.sshProfiles
1336 1336 dataOut.flagNoData = False
1337 1337 dataOut.heightList = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
1338 1338 dataOut.nProfiles = int(dataOut.nProfiles*self.nsamples)
1339 1339 dataOut.profileIndex = profileIndex
1340 1340 dataOut.flagDataAsBlock = True
1341 1341 dataOut.ippSeconds = ippSeconds
1342 1342 dataOut.step = self.step
1343 1343
1344 1344 class voltACFLags(Operation):
1345 1345
1346 1346 data_acf = None
1347 1347 lags = None
1348 1348 mode = None
1349 1349 fullBuffer = None
1350 1350 pairsList = None
1351 1351 tmp = None
1352 1352
1353 1353 def __init__(self, **kwargs):
1354 1354 Operation.__init__(self, **kwargs)
1355 1355 self.isConfig = False
1356 1356
1357 1357 def setup(self,dataOut ,lags = None,mode =None, fullBuffer= None ,pairsList = None,nAvg = 1):
1358 1358
1359 1359 self.lags = lags
1360 1360 self.mode = mode
1361 1361 self.fullBuffer= fullBuffer
1362 1362 self.nAvg = nAvg
1363 1363 nChannels = dataOut.nChannels
1364 1364 nProfiles = dataOut.nProfiles
1365 1365 nHeights = dataOut.nHeights
1366 1366 self.__nProfiles = dataOut.nProfiles
1367 1367 self.__nHeis = dataOut.nHeights
1368 1368
1369 1369 if mode == 'time':
1370 1370 print "Mode lags equal time for default."
1371 1371 else:
1372 1372 print "Mode lags equal height."
1373 1373
1374 1374 if pairsList == None:
1375 1375 print "Pairs list selected by default (1,0)"
1376 1376 pairsList = [(0,1)]
1377 1377 else:
1378 1378 pairList= pairsList
1379 1379
1380 1380 if lags == None:
1381 1381 if mode=='time':
1382 1382 self.lags = numpy.arange(0,nProfiles)# -nProfiles+1, nProfiles
1383 1383 print "self.lags", len(self.lags)
1384 1384 if mode=='height':
1385 1385 self.lags = numpy.arange(0,nHeights)# -nHeights+1, nHeights
1386 1386
1387 1387 if fullBuffer:
1388 1388 self.tmp = numpy.zeros((len(pairsList), len(self.lags), nProfiles, nHeights), dtype = 'complex')*numpy.nan
1389 1389 elif mode =='time':
1390 1390 self.tmp = numpy.zeros((len(pairsList), len(self.lags), nHeights),dtype='complex')
1391 1391 elif mode =='height':
1392 1392 self.tmp = numpy.zeros(len(pairsList), (len(self.lags), nProfiles),dtype='complex')
1393 1393
1394 1394 print "lags", len(self.lags)
1395 1395 print "mode",self.mode
1396 1396 print "nChannels", nChannels
1397 1397 print "nProfiles", nProfiles
1398 1398 print "nHeights" , nHeights
1399 1399 print "pairsList", pairsList
1400 1400 print "fullBuffer", fullBuffer
1401 1401 #print "type(pairsList)",type(pairsList)
1402 1402 print "tmp.shape",self.tmp.shape
1403 1403
1404 1404 def run(self, dataOut, lags =None,mode ='time',fullBuffer= False ,pairsList = None,nAvg = 1):
1405 1405
1406 1406 dataOut.flagNoData = True
1407 1407
1408 1408 if not self.isConfig:
1409 1409 self.setup(dataOut, lags = lags,mode = mode, fullBuffer= fullBuffer ,pairsList = pairsList,nAvg=nAvg)
1410 1410 self.isConfig = True
1411 1411
1412 1412 if dataOut.type == "Voltage":
1413 1413
1414 1414 data_pre = dataOut.data #data
1415 1415
1416 1416
1417 1417 # Here is the loop :D
1418 1418 for l in range(len(pairsList)):
1419 1419 ch0 = pairsList[l][0]
1420 1420 ch1 = pairsList[l][1]
1421 1421
1422 1422 for i in range(len(self.lags)):
1423 1423 idx = self.lags[i]
1424 1424 if self.mode == 'time':
1425 1425 acf0 = data_pre[ch0,:self.__nProfiles-idx,:]*numpy.conj(data_pre[ch1,idx:,:]) # pair,lag,height
1426 1426 else:
1427 1427 acf0 = data_pre[ch0,:,:self.__nHeights-idx]*numpy.conj(data_pre[ch1,:,idx:]) # pair,lag,profile
1428 1428
1429 1429 if self.fullBuffer:
1430 1430 self.tmp[l,i,:acf0.shape[0],:]= acf0
1431 1431 else:
1432 1432 self.tmp[l,i,:]= numpy.sum(acf0,axis=0)
1433 1433
1434 1434 if self.fullBuffer:
1435 1435
1436 1436 self.tmp = numpy.sum(numpy.reshape(self.tmp,(self.tmp.shape[0],self.tmp.shape[1],self.tmp.shape[2]/self.nAvg,self.nAvg,self.tmp.shape[3])),axis=3)
1437 1437 dataOut.nAvg = self.nAvg
1438 1438
1439 1439 if self.mode == 'time':
1440 1440 delta = dataOut.ippSeconds*dataOut.nCohInt
1441 1441 else:
1442 1442 delta = dataOut.heightList[1] - dataOut.heightList[0]
1443 1443
1444 1444 shape= self.tmp.shape # mode time
1445 1445 # Normalizando
1446 1446 for i in range(len(pairsList)):
1447 1447 for j in range(shape[2]):
1448 1448 self.tmp[i,:,j]= self.tmp[i,:,j].real / numpy.max(numpy.abs(self.tmp[i,:,j]))
1449 1449
1450 1450
1451 1451 #import matplotlib.pyplot as plt
1452 1452 #print "test",self.tmp.shape
1453 1453 #print self.tmp[0,0,0]
1454 1454 #print numpy.max(numpy.abs(self.tmp[0,:,0]))
1455 1455 #acf_tmp=self.tmp[0,:,100].real/numpy.max(numpy.abs(self.tmp[0,:,100]))
1456 1456 #print acf_tmp
1457 1457 #plt.plot(acf_tmp)
1458 1458 #plt.show()
1459 1459 #import time
1460 1460 #time.sleep(20)
1461 1461
1462 1462 dataOut.data = self.tmp
1463 1463 dataOut.mode = self.mode
1464 1464 dataOut.nLags = len(self.lags)
1465 dataOut.nProfiles = len(self.lags)
1465 1466 dataOut.pairsList = pairsList
1466 1467 dataOut.nPairs = len(pairsList)
1467 1468 dataOut.lagRange = numpy.array(self.lags)*delta
1468 1469 dataOut.flagDataAsBlock = True
1469 1470 dataOut.flagNoData = False
1470 1471
1471 1472 import time
1472 1473 #################################################
1473 1474
1474 1475 class decoPseudorandom(Operation):
1475 1476
1476 1477 nProfiles= 0
1477 1478 buffer= None
1478 1479 isConfig = False
1479 1480
1480 1481 def setup(self, clen= 10000,seed= 0,Nranges= 1000,oversample=1):
1481 1482 #code = create_pseudo_random_code(clen=clen, seed=seed)
1482 1483 code= rep_seq(create_pseudo_random_code(clen=clen, seed=seed),rep=oversample)
1483 1484 #print ("code_rx", code.shape)
1484 1485 #N = int(an_len/clen) # 100
1485 1486 B_cache = 0
1486 1487 r_cache = 0
1487 1488 B_cached = False
1488 1489 r = create_estimation_matrix(code=code, cache=True, rmax=Nranges)
1489 1490 #print ("code shape", code.shape)
1490 1491 #print ("seed",seed)
1491 1492 #print ("Code", code[0:10])
1492 1493 self.B = r['B']
1493 1494
1494 1495
1495 1496 def run (self,dataOut,length_code= 10000,seed= 0,Nranges= 1000,oversample=1):
1496 1497 #print((dataOut.data.shape))
1497 1498 if not self.isConfig:
1498 1499 self.setup(clen= length_code,seed= seed,Nranges= Nranges,oversample=oversample)
1499 1500 self.isConfig = True
1500 1501
1501 1502 dataOut.flagNoData = True
1502 1503 data =dataOut.data
1503 1504 #print "length_CODE",length_code
1504 1505 data_shape = (data.shape[1])
1505 1506 #print "data_shape",data_shape
1506 1507 n = (length_code /data_shape)
1507 1508 #print "we need this number of sample",n
1508 1509
1509 1510 if n>0 and self.buffer is None:
1510 1511 self.buffer = numpy.zeros([1, length_code], dtype=numpy.complex64)
1511 1512 self.buffer[0][0:data_shape] = data[0]
1512 1513 #print "FIRST CREATION",self.buffer.shape
1513 1514
1514 1515 else:
1515 1516 self.buffer[0][self.nProfiles*data_shape:(self.nProfiles+1)*data_shape]=data[0]
1516 1517
1517 1518 #print "buffer_shape",(self.buffer.shape)
1518 1519 self.nProfiles += 1
1519 1520 #print "count",self.nProfiles
1520 1521
1521 1522 if self.nProfiles== n:
1522 1523 temporal = numpy.dot(self.B, numpy.transpose(self.buffer))
1523 1524 #print temporal.shape
1524 1525 #import time
1525 1526 #time.sleep(40)
1526 1527 dataOut.data=numpy.transpose(temporal)
1527 1528
1528 1529 dataOut.flagNoData = False
1529 1530 self.buffer= None
1530 1531 self.nProfiles = 0
1531 1532
1532 1533 # import collections
1533 1534 # from scipy.stats import mode
1534 1535 #
1535 1536 # class Synchronize(Operation):
1536 1537 #
1537 1538 # isConfig = False
1538 1539 # __profIndex = 0
1539 1540 #
1540 1541 # def __init__(self, **kwargs):
1541 1542 #
1542 1543 # Operation.__init__(self, **kwargs)
1543 1544 # # self.isConfig = False
1544 1545 # self.__powBuffer = None
1545 1546 # self.__startIndex = 0
1546 1547 # self.__pulseFound = False
1547 1548 #
1548 1549 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1549 1550 #
1550 1551 # #Read data
1551 1552 #
1552 1553 # powerdB = dataOut.getPower(channel = channel)
1553 1554 # noisedB = dataOut.getNoise(channel = channel)[0]
1554 1555 #
1555 1556 # self.__powBuffer.extend(powerdB.flatten())
1556 1557 #
1557 1558 # dataArray = numpy.array(self.__powBuffer)
1558 1559 #
1559 1560 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1560 1561 #
1561 1562 # maxValue = numpy.nanmax(filteredPower)
1562 1563 #
1563 1564 # if maxValue < noisedB + 10:
1564 1565 # #No se encuentra ningun pulso de transmision
1565 1566 # return None
1566 1567 #
1567 1568 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1568 1569 #
1569 1570 # if len(maxValuesIndex) < 2:
1570 1571 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1571 1572 # return None
1572 1573 #
1573 1574 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1574 1575 #
1575 1576 # #Seleccionar solo valores con un espaciamiento de nSamples
1576 1577 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1577 1578 #
1578 1579 # if len(pulseIndex) < 2:
1579 1580 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1580 1581 # return None
1581 1582 #
1582 1583 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1583 1584 #
1584 1585 # #remover senales que se distancien menos de 10 unidades o muestras
1585 1586 # #(No deberian existir IPP menor a 10 unidades)
1586 1587 #
1587 1588 # realIndex = numpy.where(spacing > 10 )[0]
1588 1589 #
1589 1590 # if len(realIndex) < 2:
1590 1591 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1591 1592 # return None
1592 1593 #
1593 1594 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1594 1595 # realPulseIndex = pulseIndex[realIndex]
1595 1596 #
1596 1597 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1597 1598 #
1598 1599 # print "IPP = %d samples" %period
1599 1600 #
1600 1601 # self.__newNSamples = dataOut.nHeights #int(period)
1601 1602 # self.__startIndex = int(realPulseIndex[0])
1602 1603 #
1603 1604 # return 1
1604 1605 #
1605 1606 #
1606 1607 # def setup(self, nSamples, nChannels, buffer_size = 4):
1607 1608 #
1608 1609 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1609 1610 # maxlen = buffer_size*nSamples)
1610 1611 #
1611 1612 # bufferList = []
1612 1613 #
1613 1614 # for i in range(nChannels):
1614 1615 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1615 1616 # maxlen = buffer_size*nSamples)
1616 1617 #
1617 1618 # bufferList.append(bufferByChannel)
1618 1619 #
1619 1620 # self.__nSamples = nSamples
1620 1621 # self.__nChannels = nChannels
1621 1622 # self.__bufferList = bufferList
1622 1623 #
1623 1624 # def run(self, dataOut, channel = 0):
1624 1625 #
1625 1626 # if not self.isConfig:
1626 1627 # nSamples = dataOut.nHeights
1627 1628 # nChannels = dataOut.nChannels
1628 1629 # self.setup(nSamples, nChannels)
1629 1630 # self.isConfig = True
1630 1631 #
1631 1632 # #Append new data to internal buffer
1632 1633 # for thisChannel in range(self.__nChannels):
1633 1634 # bufferByChannel = self.__bufferList[thisChannel]
1634 1635 # bufferByChannel.extend(dataOut.data[thisChannel])
1635 1636 #
1636 1637 # if self.__pulseFound:
1637 1638 # self.__startIndex -= self.__nSamples
1638 1639 #
1639 1640 # #Finding Tx Pulse
1640 1641 # if not self.__pulseFound:
1641 1642 # indexFound = self.__findTxPulse(dataOut, channel)
1642 1643 #
1643 1644 # if indexFound == None:
1644 1645 # dataOut.flagNoData = True
1645 1646 # return
1646 1647 #
1647 1648 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1648 1649 # self.__pulseFound = True
1649 1650 # self.__startIndex = indexFound
1650 1651 #
1651 1652 # #If pulse was found ...
1652 1653 # for thisChannel in range(self.__nChannels):
1653 1654 # bufferByChannel = self.__bufferList[thisChannel]
1654 1655 # #print self.__startIndex
1655 1656 # x = numpy.array(bufferByChannel)
1656 1657 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1657 1658 #
1658 1659 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1659 1660 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1660 1661 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1661 1662 #
1662 1663 # dataOut.data = self.__arrayBuffer
1663 1664 #
1664 1665 # self.__startIndex += self.__newNSamples
1665 1666 #
1666 1667 # return
Show file before | Show file after | Hide comments
@@ -1,145 +1,145
1 1 #!/usr/bin/env python
2 2 '''
3 3 Created on Jul 7, 2014
4 4
5 5 @author: roj-idl71
6 6 '''
7 7 import os, sys
8 8
9 9 #path = os.path.dirname(os.getcwd())
10 10 #path = os.path.join(path, 'source')
11 11 #sys.path.insert(0, path)
12 12
13 13 from schainpy.controller import Project
14 14
15 15 if __name__ == '__main__':
16 16
17 17 desc = "High altitude experiment SHORT "
18 18 filename = "schain.xml"
19 dpath = '/home/soporte/test_avp'
19 dpath = '/media/soporte/APOLLO/hybrid'#'/home/soporte/test_avp'
20 20 figpath = "/home/soporte/pics"
21 21 remotefolder = "/home/wmaster/graficos"
22 22 t=['0','24']
23 23 db_range=['15','35']
24 24 db_range=['25','45']
25 25 db_range=['0','35']
26 26 period=60
27 27
28 28 controllerObj = Project()
29 29
30 30 controllerObj.setup(id = '191', name='test01', description=desc)
31 31
32 32 readUnitConfObj = controllerObj.addReadUnit(datatype='Voltage',
33 33 path=dpath,
34 startDate='2019/10/14',
34 startDate='2019/10/16',
35 35 # startDate='2018/06/18',
36 endDate='2019/10/14',
36 endDate='2019/10/16',
37 37 # endDate='2018/06/18',
38 startTime='14:00:00',
38 startTime='00:00:00',
39 39 endTime='23:59:59',
40 40 online=0,
41 41 walk=0,
42 42 expLabel='',
43 43 delay=20)
44 44
45 45
46 46 opObj00 = readUnitConfObj.addOperation(name='printInfo')
47 47 opObj00 = readUnitConfObj.addOperation(name='printNumberOfBlock')
48 48
49 49 procUnitConfObj1 = controllerObj.addProcUnit(datatype='Voltage', inputId=readUnitConfObj.getId())
50 50
51 51 opObj11 = procUnitConfObj1.addOperation(name='ProfileSelector', optype='other')
52 52 #opObj11.addParameter(name='profileList', value='(0, 1, 2, 3, 20, 21, 22, 23, 40, 41, 42, 43, 60, 61, 62, 63, 80, 81, 82, 83)', format='intlist')
53 53 opObj11.addParameter(name='rangeList', value='(0,127)', format='intlist') # (0,63),(127,149)
54 54
55 55 opObj11 = procUnitConfObj1.addOperation(name='selectHeights')
56 56 opObj11.addParameter(name='minHei', value='0', format='int')
57 57 opObj11.addParameter(name='maxHei', value='3000', format='int')
58 58
59 59
60 60 #opObj11 = procUnitConfObj1.addOperation(name='filterByHeights')
61 61 #opObj11.addParameter(name='window', value='4', format='int')
62 62
63 63 # opObj11 = procUnitConfObj1.addOperation(name='deFlip')
64 64 # opObj11.addParameter(name='channelList', value='1,3,5,7', format='intlist')
65 65
66 66 opObj11 = procUnitConfObj1.addOperation(name='SSheightProfiles', optype='other')
67 67 opObj11.addParameter(name='step', value='1', format='int')
68 68 opObj11.addParameter(name='nsamples', value='64', format='int')
69 69
70 70
71 71 procUnitConfObj1SPC = controllerObj.addProcUnit(datatype='Spectra', inputId=procUnitConfObj1.getId())
72 72 procUnitConfObj1SPC.addParameter(name='nFFTPoints', value='128', format='int')
73 73 procUnitConfObj1SPC.addParameter(name='nProfiles', value='64', format='int')
74 74 # procUnitConfObj1SPC.addParameter(name='pairsList', value='(1,0),(3,2),(5,4),(7,6)', format='pairsList')
75 75
76 76 opObj11 = procUnitConfObj1SPC.addOperation(name='IncohInt', optype='other')
77 77 opObj11.addParameter(name='timeInterval', value='10', format='int')
78 78
79 79 procUnitConfObj2SPC = controllerObj.addProcUnit(datatype='SpectraAFCProc', inputId=procUnitConfObj1SPC.getId())
80 80 procUnitConfObj2SPC.addParameter(name='real', value='1', format='int')
81 81
82 82 opObj11 = procUnitConfObj1SPC.addOperation(name='SpectraPlot', optype='other')
83 83 opObj11.addParameter(name='id', value='1', format='int')
84 84 opObj11.addParameter(name='wintitle', value='October 2019', format='str')
85 85 #opObj11.addParameter(name='xaxis', value='time', format='str')
86 86
87 87 #opObj11.addParameter(name='zmin', value=db_range[0], format='int')
88 88 #opObj11.addParameter(name='zmax', value=db_range[1], format='int')
89 89 #opObj11.addParameter(name='ymin', value='0', format='int')
90 90 #opObj11.addParameter(name='ymax', value='1500', format='int')
91 91 #opObj11.addParameter(name='xmin', value='-5', format='int')
92 92 #opObj11.addParameter(name='xmax', value='5', format='int')
93 93 opObj11.addParameter(name='showprofile', value='1', format='int')
94 94 opObj11.addParameter(name='save', value='1', format='int')
95 95 opObj11.addParameter(name='figpath', value=figpath)
96 96
97 97 opObj11 = procUnitConfObj1SPC.addOperation(name='RTIPlot', optype='other')
98 98 opObj11.addParameter(name='id', value='2', format='int')
99 99 opObj11.addParameter(name='wintitle', value='October 2019', format='str')
100 100 opObj11.addParameter(name='xmin', value=t[0], format='float')
101 101 opObj11.addParameter(name='xmax', value=t[1], format='float')
102 102 #opObj11.addParameter(name='ymin', value='0', format='int')
103 103 #opObj11.addParameter(name='ymax', value='1500', format='int')
104 104 #opObj11.addParameter(name='zmin', value=db_range[0], format='int')
105 105 #opObj11.addParameter(name='zmax', value=db_range[1], format='int')
106 106 opObj11.addParameter(name='showprofile', value='1', format='int')
107 107 opObj11.addParameter(name='save', value='1', format='int')
108 108 opObj11.addParameter(name='figpath', value=figpath)
109 109
110 110
111 111
112 112 opObj11 = procUnitConfObj2SPC.addOperation(name='ACFPlot', optype='other')
113 113 opObj11.addParameter(name='id', value='3', format='int')
114 114 opObj11.addParameter(name='wintitle', value = 'October 2019', format='str')
115 115 opObj11.addParameter(name='xaxis', value='time', format='str')
116 116 opObj11.addParameter(name='channel', value='0', format='int')
117 117 #opObj11.addParameter(name='nSampleList', value='(50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70)', format='intList')
118 118 opObj11.addParameter(name='resolutionFactor', value='15', format='int')
119 119 #opObj11.addParameter(name='zmin', value=0.5, format='int')
120 120 #opObj11.addParameter(name='zmax', value=-0.5, format='int')
121 121 #opObj11.addParameter(name='ymin', value='0', format='int')
122 122 #opObj11.addParameter(name='ymax', value='0.5', format='int')
123 123 #opObj11.addParameter(name='xmin', value='-1.2', format='int')
124 124 #opObj11.addParameter(name='xmax', value='1.2', format='int')
125 125 opObj11.addParameter(name='show', value='1', format='int')
126 126 opObj11.addParameter(name='save', value='1', format='int')
127 127 opObj11.addParameter(name='figpath', value=figpath)
128 128
129 129 '''
130 130 opObj11 = procUnitConfObj1SPC.addOperation(name='Noise', optype='other')
131 131 opObj11.addParameter(name='id', value='3', format='int')
132 132 opObj11.addParameter(name='wintitle', value='Short Valley August 2019', format='str')
133 133 opObj11.addParameter(name='xmin', value=t[0], format='float')
134 134 opObj11.addParameter(name='xmax', value=t[1], format='float')
135 135 opObj11.addParameter(name='ymin', value=db_range[0], format='int')
136 136 opObj11.addParameter(name='ymax', value=db_range[1], format='int')
137 137 opObj11.addParameter(name='showprofile', value='1', format='int')
138 138 opObj11.addParameter(name='save', value='1', format='int')
139 139 opObj11.addParameter(name='figpath', value=figpath)
140 140
141 141
142 142 '''
143 143
144 144
145 145 controllerObj.start()
General Comments 0
You need to be logged in to leave comments. Login now