##// END OF EJS Templates
schain
RSS
Duplicate property eliminated
Miguel Valdez -
r630:142de39f73f1
parent child
Show More
Show file before | Show file after | Hide comments
@@ -1,1132 +1,1132
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
13 13 def getNumpyDtype(dataTypeCode):
14 14
15 15 if dataTypeCode == 0:
16 16 numpyDtype = numpy.dtype([('real','<i1'),('imag','<i1')])
17 17 elif dataTypeCode == 1:
18 18 numpyDtype = numpy.dtype([('real','<i2'),('imag','<i2')])
19 19 elif dataTypeCode == 2:
20 20 numpyDtype = numpy.dtype([('real','<i4'),('imag','<i4')])
21 21 elif dataTypeCode == 3:
22 22 numpyDtype = numpy.dtype([('real','<i8'),('imag','<i8')])
23 23 elif dataTypeCode == 4:
24 24 numpyDtype = numpy.dtype([('real','<f4'),('imag','<f4')])
25 25 elif dataTypeCode == 5:
26 26 numpyDtype = numpy.dtype([('real','<f8'),('imag','<f8')])
27 27 else:
28 28 raise ValueError, 'dataTypeCode was not defined'
29 29
30 30 return numpyDtype
31 31
32 32 def getDataTypeCode(numpyDtype):
33 33
34 34 if numpyDtype == numpy.dtype([('real','<i1'),('imag','<i1')]):
35 35 datatype = 0
36 36 elif numpyDtype == numpy.dtype([('real','<i2'),('imag','<i2')]):
37 37 datatype = 1
38 38 elif numpyDtype == numpy.dtype([('real','<i4'),('imag','<i4')]):
39 39 datatype = 2
40 40 elif numpyDtype == numpy.dtype([('real','<i8'),('imag','<i8')]):
41 41 datatype = 3
42 42 elif numpyDtype == numpy.dtype([('real','<f4'),('imag','<f4')]):
43 43 datatype = 4
44 44 elif numpyDtype == numpy.dtype([('real','<f8'),('imag','<f8')]):
45 45 datatype = 5
46 46 else:
47 47 datatype = None
48 48
49 49 return datatype
50 50
51 51 def hildebrand_sekhon(data, navg):
52 52 """
53 53 This method is for the objective determination of the noise level in Doppler spectra. This
54 54 implementation technique is based on the fact that the standard deviation of the spectral
55 55 densities is equal to the mean spectral density for white Gaussian noise
56 56
57 57 Inputs:
58 58 Data : heights
59 59 navg : numbers of averages
60 60
61 61 Return:
62 62 -1 : any error
63 63 anoise : noise's level
64 64 """
65 65
66 66 sortdata = numpy.sort(data,axis=None)
67 67 lenOfData = len(sortdata)
68 68 nums_min = lenOfData/10
69 69
70 70 if (lenOfData/10) > 2:
71 71 nums_min = lenOfData/10
72 72 else:
73 73 nums_min = 2
74 74
75 75 sump = 0.
76 76
77 77 sumq = 0.
78 78
79 79 j = 0
80 80
81 81 cont = 1
82 82
83 83 while((cont==1)and(j<lenOfData)):
84 84
85 85 sump += sortdata[j]
86 86
87 87 sumq += sortdata[j]**2
88 88
89 89 if j > nums_min:
90 90 rtest = float(j)/(j-1) + 1.0/navg
91 91 if ((sumq*j) > (rtest*sump**2)):
92 92 j = j - 1
93 93 sump = sump - sortdata[j]
94 94 sumq = sumq - sortdata[j]**2
95 95 cont = 0
96 96
97 97 j += 1
98 98
99 99 lnoise = sump /j
100 100 stdv = numpy.sqrt((sumq - lnoise**2)/(j - 1))
101 101 return lnoise
102 102
103 103 class Beam:
104 104 def __init__(self):
105 105 self.codeList = []
106 106 self.azimuthList = []
107 107 self.zenithList = []
108 108
109 109 class GenericData(object):
110 110
111 111 flagNoData = True
112 112
113 113 def __init__(self):
114 114
115 115 raise ValueError, "This class has not been implemented"
116 116
117 117 def copy(self, inputObj=None):
118 118
119 119 if inputObj == None:
120 120 return copy.deepcopy(self)
121 121
122 122 for key in inputObj.__dict__.keys():
123 123 self.__dict__[key] = inputObj.__dict__[key]
124 124
125 125 def deepcopy(self):
126 126
127 127 return copy.deepcopy(self)
128 128
129 129 def isEmpty(self):
130 130
131 131 return self.flagNoData
132 132
133 133 class JROData(GenericData):
134 134
135 135 # m_BasicHeader = BasicHeader()
136 136 # m_ProcessingHeader = ProcessingHeader()
137 137
138 138 systemHeaderObj = SystemHeader()
139 139
140 140 radarControllerHeaderObj = RadarControllerHeader()
141 141
142 142 # data = None
143 143
144 144 type = None
145 145
146 146 datatype = None #dtype but in string
147 147
148 148 # dtype = None
149 149
150 150 # nChannels = None
151 151
152 152 # nHeights = None
153 153
154 154 nProfiles = None
155 155
156 156 heightList = None
157 157
158 158 channelList = None
159 159
160 160 flagDiscontinuousBlock = False
161 161
162 162 useLocalTime = False
163 163
164 164 utctime = None
165 165
166 166 timeZone = None
167 167
168 168 dstFlag = None
169 169
170 170 errorCount = None
171 171
172 172 blocksize = None
173 173
174 174 # nCode = None
175 175 #
176 176 # nBaud = None
177 177 #
178 178 # code = None
179 179
180 180 flagDecodeData = False #asumo q la data no esta decodificada
181 181
182 182 flagDeflipData = False #asumo q la data no esta sin flip
183 183
184 184 flagShiftFFT = False
185 185
186 186 # ippSeconds = None
187 187
188 188 # timeInterval = None
189 189
190 190 nCohInt = None
191 191
192 192 # noise = None
193 193
194 194 windowOfFilter = 1
195 195
196 196 #Speed of ligth
197 197 C = 3e8
198 198
199 199 frequency = 49.92e6
200 200
201 201 realtime = False
202 202
203 203 beacon_heiIndexList = None
204 204
205 205 last_block = None
206 206
207 207 blocknow = None
208 208
209 209 azimuth = None
210 210
211 211 zenith = None
212 212
213 213 beam = Beam()
214 214
215 215 profileIndex = None
216 216
217 217 def __init__(self):
218 218
219 219 raise ValueError, "This class has not been implemented"
220 220
221 221 def getNoise(self):
222 222
223 223 raise ValueError, "Not implemented"
224 224
225 225 def getNChannels(self):
226 226
227 227 return len(self.channelList)
228 228
229 229 def getChannelIndexList(self):
230 230
231 231 return range(self.nChannels)
232 232
233 233 def getNHeights(self):
234 234
235 235 return len(self.heightList)
236 236
237 237 def getHeiRange(self, extrapoints=0):
238 238
239 239 heis = self.heightList
240 240 # deltah = self.heightList[1] - self.heightList[0]
241 241 #
242 242 # heis.append(self.heightList[-1])
243 243
244 244 return heis
245 245
246 246 def getltctime(self):
247 247
248 248 if self.useLocalTime:
249 249 return self.utctime - self.timeZone*60
250 250
251 251 return self.utctime
252 252
253 253 def getDatatime(self):
254 254
255 255 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
256 256 return datatimeValue
257 257
258 258 def getTimeRange(self):
259 259
260 260 datatime = []
261 261
262 262 datatime.append(self.ltctime)
263 263 datatime.append(self.ltctime + self.timeInterval+60)
264 264
265 265 datatime = numpy.array(datatime)
266 266
267 267 return datatime
268 268
269 269 def getFmax(self):
270 270
271 271 PRF = 1./(self.ippSeconds * self.nCohInt)
272 272
273 273 fmax = PRF/2.
274 274
275 275 return fmax
276 276
277 277 def getVmax(self):
278 278
279 279 _lambda = self.C/self.frequency
280 280
281 281 vmax = self.getFmax() * _lambda
282 282
283 283 return vmax
284 284
285 285 def get_ippSeconds(self):
286 286 '''
287 287 '''
288 288 return self.radarControllerHeaderObj.ippSeconds
289 289
290 290 def set_ippSeconds(self, ippSeconds):
291 291 '''
292 292 '''
293 293
294 294 self.radarControllerHeaderObj.ippSeconds = ippSeconds
295 295
296 296 return
297 297
298 298 def get_dtype(self):
299 299 '''
300 300 '''
301 301 return getNumpyDtype(self.datatype)
302 302
303 303 def set_dtype(self, numpyDtype):
304 304 '''
305 305 '''
306 306
307 307 self.datatype = getDataTypeCode(numpyDtype)
308 308
309 309 def get_code(self):
310 310 '''
311 311 '''
312 312 return self.radarControllerHeaderObj.code
313 313
314 314 def set_code(self, code):
315 315 '''
316 316 '''
317 317 self.radarControllerHeaderObj.code = code
318 318
319 319 return
320 320
321 321 def get_ncode(self):
322 322 '''
323 323 '''
324 324 return self.radarControllerHeaderObj.nCode
325 325
326 326 def set_ncode(self, nCode):
327 327 '''
328 328 '''
329 329 self.radarControllerHeaderObj.nCode = nCode
330 330
331 331 return
332 332
333 333 def get_nbaud(self):
334 334 '''
335 335 '''
336 336 return self.radarControllerHeaderObj.nBaud
337 337
338 338 def set_nbaud(self, nBaud):
339 339 '''
340 340 '''
341 341 self.radarControllerHeaderObj.nBaud = nBaud
342 342
343 343 return
344 344 # def getTimeInterval(self):
345 345 #
346 346 # raise IOError, "This method should be implemented inside each Class"
347 347
348 348 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
349 349 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
350 350 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
351 351 #noise = property(getNoise, "I'm the 'nHeights' property.")
352 352 datatime = property(getDatatime, "I'm the 'datatime' property")
353 353 ltctime = property(getltctime, "I'm the 'ltctime' property")
354 354 ippSeconds = property(get_ippSeconds, set_ippSeconds)
355 355 dtype = property(get_dtype, set_dtype)
356 356 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
357 357 code = property(get_code, set_code)
358 358 nCode = property(get_ncode, set_ncode)
359 359 nBaud = property(get_nbaud, set_nbaud)
360 360
361 361 class Voltage(JROData):
362 362
363 363 #data es un numpy array de 2 dmensiones (canales, alturas)
364 364 data = None
365 365
366 366 def __init__(self):
367 367 '''
368 368 Constructor
369 369 '''
370 370
371 371 self.useLocalTime = True
372 372
373 373 self.radarControllerHeaderObj = RadarControllerHeader()
374 374
375 375 self.systemHeaderObj = SystemHeader()
376 376
377 377 self.type = "Voltage"
378 378
379 379 self.data = None
380 380
381 381 # self.dtype = None
382 382
383 383 # self.nChannels = 0
384 384
385 385 # self.nHeights = 0
386 386
387 387 self.nProfiles = None
388 388
389 389 self.heightList = None
390 390
391 391 self.channelList = None
392 392
393 393 # self.channelIndexList = None
394 394
395 395 self.flagNoData = True
396 396
397 397 self.flagDiscontinuousBlock = False
398 398
399 399 self.utctime = None
400 400
401 401 self.timeZone = None
402 402
403 403 self.dstFlag = None
404 404
405 405 self.errorCount = None
406 406
407 407 self.nCohInt = None
408 408
409 409 self.blocksize = None
410 410
411 411 self.flagDecodeData = False #asumo q la data no esta decodificada
412 412
413 413 self.flagDeflipData = False #asumo q la data no esta sin flip
414 414
415 415 self.flagShiftFFT = False
416 416
417 417 self.flagDataAsBlock = False #Asumo que la data es leida perfil a perfil
418 418
419 419 self.profileIndex = 0
420 420
421 421 def getNoisebyHildebrand(self, channel = None):
422 422 """
423 423 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
424 424
425 425 Return:
426 426 noiselevel
427 427 """
428 428
429 429 if channel != None:
430 430 data = self.data[channel]
431 431 nChannels = 1
432 432 else:
433 433 data = self.data
434 434 nChannels = self.nChannels
435 435
436 436 noise = numpy.zeros(nChannels)
437 437 power = data * numpy.conjugate(data)
438 438
439 439 for thisChannel in range(nChannels):
440 440 if nChannels == 1:
441 441 daux = power[:].real
442 442 else:
443 443 daux = power[thisChannel,:].real
444 444 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
445 445
446 446 return noise
447 447
448 448 def getNoise(self, type = 1, channel = None):
449 449
450 450 if type == 1:
451 451 noise = self.getNoisebyHildebrand(channel)
452 452
453 453 return 10*numpy.log10(noise)
454 454
455 455 def getPower(self, channel = None):
456 456
457 457 if channel != None:
458 458 data = self.data[channel]
459 459 else:
460 460 data = self.data
461 461
462 462 power = data * numpy.conjugate(data)
463 463
464 464 return 10*numpy.log10(power.real)
465 465
466 466 def getTimeInterval(self):
467 467
468 468 timeInterval = self.ippSeconds * self.nCohInt
469 469
470 470 return timeInterval
471 471
472 472 noise = property(getNoise, "I'm the 'nHeights' property.")
473 473 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
474 474
475 475 class Spectra(JROData):
476 476
477 477 #data es un numpy array de 2 dmensiones (canales, perfiles, alturas)
478 478 data_spc = None
479 479
480 480 #data es un numpy array de 2 dmensiones (canales, pares, alturas)
481 481 data_cspc = None
482 482
483 483 #data es un numpy array de 2 dmensiones (canales, alturas)
484 484 data_dc = None
485 485
486 486 nFFTPoints = None
487 487
488 488 # nPairs = None
489 489
490 490 pairsList = None
491 491
492 492 nIncohInt = None
493 493
494 494 wavelength = None #Necesario para cacular el rango de velocidad desde la frecuencia
495 495
496 496 nCohInt = None #se requiere para determinar el valor de timeInterval
497 497
498 498 ippFactor = None
499 499
500 500 profileIndex = 0
501 501
502 502 def __init__(self):
503 503 '''
504 504 Constructor
505 505 '''
506 506
507 507 self.useLocalTime = True
508 508
509 509 self.radarControllerHeaderObj = RadarControllerHeader()
510 510
511 511 self.systemHeaderObj = SystemHeader()
512 512
513 513 self.type = "Spectra"
514 514
515 515 # self.data = None
516 516
517 517 # self.dtype = None
518 518
519 519 # self.nChannels = 0
520 520
521 521 # self.nHeights = 0
522 522
523 523 self.nProfiles = None
524 524
525 525 self.heightList = None
526 526
527 527 self.channelList = None
528 528
529 529 # self.channelIndexList = None
530 530
531 531 self.pairsList = None
532 532
533 533 self.flagNoData = True
534 534
535 535 self.flagDiscontinuousBlock = False
536 536
537 537 self.utctime = None
538 538
539 539 self.nCohInt = None
540 540
541 541 self.nIncohInt = None
542 542
543 543 self.blocksize = None
544 544
545 545 self.nFFTPoints = None
546 546
547 547 self.wavelength = None
548 548
549 549 self.flagDecodeData = False #asumo q la data no esta decodificada
550 550
551 551 self.flagDeflipData = False #asumo q la data no esta sin flip
552 552
553 553 self.flagShiftFFT = False
554 554
555 555 self.ippFactor = 1
556 556
557 557 #self.noise = None
558 558
559 559 self.beacon_heiIndexList = []
560 560
561 561 self.noise_estimation = None
562 562
563 563
564 564 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
565 565 """
566 566 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
567 567
568 568 Return:
569 569 noiselevel
570 570 """
571 571
572 572 noise = numpy.zeros(self.nChannels)
573 573
574 574 for channel in range(self.nChannels):
575 575 daux = self.data_spc[channel,xmin_index:xmax_index,ymin_index:ymax_index]
576 576 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
577 577
578 578 return noise
579 579
580 580 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
581 581
582 582 if self.noise_estimation != None:
583 583 return self.noise_estimation #this was estimated by getNoise Operation defined in jroproc_spectra.py
584 584 else:
585 585 noise = self.getNoisebyHildebrand(xmin_index, xmax_index, ymin_index, ymax_index)
586 586 return noise
587 587
588 588
589 589 def getFreqRange(self, extrapoints=0):
590 590
591 591 deltafreq = self.getFmax() / (self.nFFTPoints*self.ippFactor)
592 592 freqrange = deltafreq*(numpy.arange(self.nFFTPoints+extrapoints)-self.nFFTPoints/2.) - deltafreq/2
593 593
594 594 return freqrange
595 595
596 596 def getVelRange(self, extrapoints=0):
597 597
598 598 deltav = self.getVmax() / (self.nFFTPoints*self.ippFactor)
599 599 velrange = deltav*(numpy.arange(self.nFFTPoints+extrapoints)-self.nFFTPoints/2.) - deltav/2
600 600
601 601 return velrange
602 602
603 603 def getNPairs(self):
604 604
605 605 return len(self.pairsList)
606 606
607 607 def getPairsIndexList(self):
608 608
609 609 return range(self.nPairs)
610 610
611 611 def getNormFactor(self):
612 612
613 613 pwcode = 1
614 614
615 615 if self.flagDecodeData:
616 616 pwcode = numpy.sum(self.code[0]**2)
617 617 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
618 618 normFactor = self.nProfiles*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
619 619
620 620 return normFactor
621 621
622 622 def getFlagCspc(self):
623 623
624 624 if self.data_cspc is None:
625 625 return True
626 626
627 627 return False
628 628
629 629 def getFlagDc(self):
630 630
631 631 if self.data_dc is None:
632 632 return True
633 633
634 634 return False
635 635
636 636 def getTimeInterval(self):
637 637
638 638 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles
639 639
640 640 return timeInterval
641 641
642 642 def setValue(self, value):
643 643
644 644 print "This property should not be initialized"
645 645
646 646 return
647 647
648 648 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
649 649 pairsIndexList = property(getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
650 650 normFactor = property(getNormFactor, setValue, "I'm the 'getNormFactor' property.")
651 651 flag_cspc = property(getFlagCspc, setValue)
652 652 flag_dc = property(getFlagDc, setValue)
653 653 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
654 654 timeInterval = property(getTimeInterval, setValue, "I'm the 'timeInterval' property")
655 655
656 656 class SpectraHeis(Spectra):
657 657
658 658 data_spc = None
659 659
660 660 data_cspc = None
661 661
662 662 data_dc = None
663 663
664 664 nFFTPoints = None
665 665
666 666 # nPairs = None
667 667
668 668 pairsList = None
669 669
670 670 nCohInt = None
671 671
672 672 nIncohInt = None
673 673
674 674 def __init__(self):
675 675
676 676 self.radarControllerHeaderObj = RadarControllerHeader()
677 677
678 678 self.systemHeaderObj = SystemHeader()
679 679
680 680 self.type = "SpectraHeis"
681 681
682 682 # self.dtype = None
683 683
684 684 # self.nChannels = 0
685 685
686 686 # self.nHeights = 0
687 687
688 688 self.nProfiles = None
689 689
690 690 self.heightList = None
691 691
692 692 self.channelList = None
693 693
694 694 # self.channelIndexList = None
695 695
696 696 self.flagNoData = True
697 697
698 698 self.flagDiscontinuousBlock = False
699 699
700 700 # self.nPairs = 0
701 701
702 702 self.utctime = None
703 703
704 704 self.blocksize = None
705 705
706 706 self.profileIndex = 0
707 707
708 708 self.nCohInt = 1
709 709
710 710 self.nIncohInt = 1
711 711
712 712 def getNormFactor(self):
713 713 pwcode = 1
714 714 if self.flagDecodeData:
715 715 pwcode = numpy.sum(self.code[0]**2)
716 716
717 717 normFactor = self.nIncohInt*self.nCohInt*pwcode
718 718
719 719 return normFactor
720 720
721 721 def getTimeInterval(self):
722 722
723 723 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
724 724
725 725 return timeInterval
726 726
727 727 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
728 728 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
729 729
730 730 class Fits(JROData):
731 731
732 732 heightList = None
733 733
734 734 channelList = None
735 735
736 736 flagNoData = True
737 737
738 738 flagDiscontinuousBlock = False
739 739
740 740 useLocalTime = False
741 741
742 742 utctime = None
743 743
744 744 timeZone = None
745 745
746 746 # ippSeconds = None
747 747
748 748 # timeInterval = None
749 749
750 750 nCohInt = None
751 751
752 752 nIncohInt = None
753 753
754 754 noise = None
755 755
756 756 windowOfFilter = 1
757 757
758 758 #Speed of ligth
759 759 C = 3e8
760 760
761 761 frequency = 49.92e6
762 762
763 763 realtime = False
764 764
765 765
766 766 def __init__(self):
767 767
768 768 self.type = "Fits"
769 769
770 770 self.nProfiles = None
771 771
772 772 self.heightList = None
773 773
774 774 self.channelList = None
775 775
776 776 # self.channelIndexList = None
777 777
778 778 self.flagNoData = True
779 779
780 780 self.utctime = None
781 781
782 782 self.nCohInt = 1
783 783
784 784 self.nIncohInt = 1
785 785
786 786 self.useLocalTime = True
787 787
788 788 self.profileIndex = 0
789 789
790 790 # self.utctime = None
791 791 # self.timeZone = None
792 792 # self.ltctime = None
793 793 # self.timeInterval = None
794 794 # self.header = None
795 795 # self.data_header = None
796 796 # self.data = None
797 797 # self.datatime = None
798 798 # self.flagNoData = False
799 799 # self.expName = ''
800 800 # self.nChannels = None
801 801 # self.nSamples = None
802 802 # self.dataBlocksPerFile = None
803 803 # self.comments = ''
804 804 #
805 805
806 806
807 807 def getltctime(self):
808 808
809 809 if self.useLocalTime:
810 810 return self.utctime - self.timeZone*60
811 811
812 812 return self.utctime
813 813
814 814 def getDatatime(self):
815 815
816 816 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
817 817 return datatime
818 818
819 819 def getTimeRange(self):
820 820
821 821 datatime = []
822 822
823 823 datatime.append(self.ltctime)
824 824 datatime.append(self.ltctime + self.timeInterval)
825 825
826 826 datatime = numpy.array(datatime)
827 827
828 828 return datatime
829 829
830 830 def getHeiRange(self):
831 831
832 832 heis = self.heightList
833 833
834 834 return heis
835 835
836 836 def isEmpty(self):
837 837
838 838 return self.flagNoData
839 839
840 840 def getNHeights(self):
841 841
842 842 return len(self.heightList)
843 843
844 844 def getNChannels(self):
845 845
846 846 return len(self.channelList)
847 847
848 848 def getChannelIndexList(self):
849 849
850 850 return range(self.nChannels)
851 851
852 852 def getNoise(self, type = 1):
853 853
854 854 #noise = numpy.zeros(self.nChannels)
855 855
856 856 if type == 1:
857 857 noise = self.getNoisebyHildebrand()
858 858
859 859 if type == 2:
860 860 noise = self.getNoisebySort()
861 861
862 862 if type == 3:
863 863 noise = self.getNoisebyWindow()
864 864
865 865 return noise
866 866
867 867 def getTimeInterval(self):
868 868
869 869 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
870 870
871 871 return timeInterval
872 872
873 873 datatime = property(getDatatime, "I'm the 'datatime' property")
874 874 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
875 875 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
876 876 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
877 877 noise = property(getNoise, "I'm the 'nHeights' property.")
878 datatime = property(getDatatime, "I'm the 'datatime' property")
878
879 879 ltctime = property(getltctime, "I'm the 'ltctime' property")
880 880 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
881 881
882 882 class Correlation(JROData):
883 883
884 884 noise = None
885 885
886 886 SNR = None
887 887
888 888 pairsAutoCorr = None #Pairs of Autocorrelation
889 889
890 890 #--------------------------------------------------
891 891
892 892 data_corr = None
893 893
894 894 data_volt = None
895 895
896 896 lagT = None # each element value is a profileIndex
897 897
898 898 lagR = None # each element value is in km
899 899
900 900 pairsList = None
901 901
902 902 calculateVelocity = None
903 903
904 904 nPoints = None
905 905
906 906 nAvg = None
907 907
908 908 bufferSize = None
909 909
910 910 def __init__(self):
911 911 '''
912 912 Constructor
913 913 '''
914 914 self.radarControllerHeaderObj = RadarControllerHeader()
915 915
916 916 self.systemHeaderObj = SystemHeader()
917 917
918 918 self.type = "Correlation"
919 919
920 920 self.data = None
921 921
922 922 self.dtype = None
923 923
924 924 self.nProfiles = None
925 925
926 926 self.heightList = None
927 927
928 928 self.channelList = None
929 929
930 930 self.flagNoData = True
931 931
932 932 self.flagDiscontinuousBlock = False
933 933
934 934 self.utctime = None
935 935
936 936 self.timeZone = None
937 937
938 938 self.dstFlag = None
939 939
940 940 self.errorCount = None
941 941
942 942 self.blocksize = None
943 943
944 944 self.flagDecodeData = False #asumo q la data no esta decodificada
945 945
946 946 self.flagDeflipData = False #asumo q la data no esta sin flip
947 947
948 948 self.pairsList = None
949 949
950 950 self.nPoints = None
951 951
952 952 def getLagTRange(self, extrapoints=0):
953 953
954 954 lagTRange = self.lagT
955 955 diff = lagTRange[1] - lagTRange[0]
956 956 extra = numpy.arange(1,extrapoints + 1)*diff + lagTRange[-1]
957 957 lagTRange = numpy.hstack((lagTRange, extra))
958 958
959 959 return lagTRange
960 960
961 961 def getLagRRange(self, extrapoints=0):
962 962
963 963 return self.lagR
964 964
965 965 def getPairsList(self):
966 966
967 967 return self.pairsList
968 968
969 969 def getCalculateVelocity(self):
970 970
971 971 return self.calculateVelocity
972 972
973 973 def getNPoints(self):
974 974
975 975 return self.nPoints
976 976
977 977 def getNAvg(self):
978 978
979 979 return self.nAvg
980 980
981 981 def getBufferSize(self):
982 982
983 983 return self.bufferSize
984 984
985 985 def getPairsAutoCorr(self):
986 986 pairsList = self.pairsList
987 987 pairsAutoCorr = numpy.zeros(self.nChannels, dtype = 'int')*numpy.nan
988 988
989 989 for l in range(len(pairsList)):
990 990 firstChannel = pairsList[l][0]
991 991 secondChannel = pairsList[l][1]
992 992
993 993 #Obteniendo pares de Autocorrelacion
994 994 if firstChannel == secondChannel:
995 995 pairsAutoCorr[firstChannel] = int(l)
996 996
997 997 pairsAutoCorr = pairsAutoCorr.astype(int)
998 998
999 999 return pairsAutoCorr
1000 1000
1001 1001 def getNoise(self, mode = 2):
1002 1002
1003 1003 indR = numpy.where(self.lagR == 0)[0][0]
1004 1004 indT = numpy.where(self.lagT == 0)[0][0]
1005 1005
1006 1006 jspectra0 = self.data_corr[:,:,indR,:]
1007 1007 jspectra = copy.copy(jspectra0)
1008 1008
1009 1009 num_chan = jspectra.shape[0]
1010 1010 num_hei = jspectra.shape[2]
1011 1011
1012 1012 freq_dc = jspectra.shape[1]/2
1013 1013 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
1014 1014
1015 1015 if ind_vel[0]<0:
1016 1016 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
1017 1017
1018 1018 if mode == 1:
1019 1019 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
1020 1020
1021 1021 if mode == 2:
1022 1022
1023 1023 vel = numpy.array([-2,-1,1,2])
1024 1024 xx = numpy.zeros([4,4])
1025 1025
1026 1026 for fil in range(4):
1027 1027 xx[fil,:] = vel[fil]**numpy.asarray(range(4))
1028 1028
1029 1029 xx_inv = numpy.linalg.inv(xx)
1030 1030 xx_aux = xx_inv[0,:]
1031 1031
1032 1032 for ich in range(num_chan):
1033 1033 yy = jspectra[ich,ind_vel,:]
1034 1034 jspectra[ich,freq_dc,:] = numpy.dot(xx_aux,yy)
1035 1035
1036 1036 junkid = jspectra[ich,freq_dc,:]<=0
1037 1037 cjunkid = sum(junkid)
1038 1038
1039 1039 if cjunkid.any():
1040 1040 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
1041 1041
1042 1042 noise = jspectra0[:,freq_dc,:] - jspectra[:,freq_dc,:]
1043 1043
1044 1044 return noise
1045 1045
1046 1046 def getTimeInterval(self):
1047 1047
1048 1048 timeInterval = self.ippSeconds * self.nCohInt * self.nPoints
1049 1049
1050 1050 return timeInterval
1051 1051
1052 1052 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1053 1053 # pairsList = property(getPairsList, "I'm the 'pairsList' property.")
1054 1054 # nPoints = property(getNPoints, "I'm the 'nPoints' property.")
1055 1055 calculateVelocity = property(getCalculateVelocity, "I'm the 'calculateVelocity' property.")
1056 1056 nAvg = property(getNAvg, "I'm the 'nAvg' property.")
1057 1057 bufferSize = property(getBufferSize, "I'm the 'bufferSize' property.")
1058 1058
1059 1059
1060 1060 class Parameters(JROData):
1061 1061
1062 1062 #Information from previous data
1063 1063
1064 1064 inputUnit = None #Type of data to be processed
1065 1065
1066 1066 operation = None #Type of operation to parametrize
1067 1067
1068 1068 normFactor = None #Normalization Factor
1069 1069
1070 1070 groupList = None #List of Pairs, Groups, etc
1071 1071
1072 1072 #Parameters
1073 1073
1074 1074 data_param = None #Parameters obtained
1075 1075
1076 1076 data_pre = None #Data Pre Parametrization
1077 1077
1078 1078 data_SNR = None #Signal to Noise Ratio
1079 1079
1080 1080 # heightRange = None #Heights
1081 1081
1082 1082 abscissaList = None #Abscissa, can be velocities, lags or time
1083 1083
1084 1084 noise = None #Noise Potency
1085 1085
1086 1086 utctimeInit = None #Initial UTC time
1087 1087
1088 1088 paramInterval = None #Time interval to calculate Parameters in seconds
1089 1089
1090 1090 #Fitting
1091 1091
1092 1092 data_error = None #Error of the estimation
1093 1093
1094 1094 constants = None
1095 1095
1096 1096 library = None
1097 1097
1098 1098 #Output signal
1099 1099
1100 1100 outputInterval = None #Time interval to calculate output signal in seconds
1101 1101
1102 1102 data_output = None #Out signal
1103 1103
1104 1104
1105 1105
1106 1106 def __init__(self):
1107 1107 '''
1108 1108 Constructor
1109 1109 '''
1110 1110 self.radarControllerHeaderObj = RadarControllerHeader()
1111 1111
1112 1112 self.systemHeaderObj = SystemHeader()
1113 1113
1114 1114 self.type = "Parameters"
1115 1115
1116 1116 def getTimeRange1(self):
1117 1117
1118 1118 datatime = []
1119 1119
1120 1120 if self.useLocalTime:
1121 1121 time1 = self.utctimeInit - self.timeZone*60
1122 1122 else:
1123 1123 time1 = utctimeInit
1124 1124
1125 1125 # datatime.append(self.utctimeInit)
1126 1126 # datatime.append(self.utctimeInit + self.outputInterval)
1127 1127 datatime.append(time1)
1128 1128 datatime.append(time1 + self.outputInterval)
1129 1129
1130 1130 datatime = numpy.array(datatime)
1131 1131
1132 1132 return datatime
General Comments 0
You need to be logged in to leave comments. Login now