@@ -9,7 +9,7 | |||
|
9 | 9 | import numpy |
|
10 | 10 | |
|
11 | 11 | from schainpy.model.graphics.jroplot_base import Plot, plt |
|
12 | ||
|
12 | import matplotlib.pyplot as plt | |
|
13 | 13 | |
|
14 | 14 | class SpectraHeisPlot(Plot): |
|
15 | 15 | |
@@ -33,17 +33,18 class SpectraHeisPlot(Plot): | |||
|
33 | 33 | meta = {} |
|
34 | 34 | spc = 10*numpy.log10(dataOut.data_spc / dataOut.normFactor) |
|
35 | 35 | data['spc_heis'] = spc |
|
36 | ||
|
37 |
return data, meta |
|
|
36 | ||
|
37 | return data, meta | |
|
38 | 38 | |
|
39 | 39 | def plot(self): |
|
40 | 40 | |
|
41 | 41 | c = 3E8 |
|
42 | deltaHeight = self.data.yrange[1] - self.data.yrange[0] | |
|
42 | deltaHeight = self.data.yrange[1] - self.data.yrange[0] # yrange = heightList | |
|
43 | 43 | x = numpy.arange(-1*len(self.data.yrange)/2., len(self.data.yrange)/2.)*(c/(2*deltaHeight*len(self.data.yrange)*1000)) |
|
44 | 44 | self.y = self.data[-1]['spc_heis'] |
|
45 | 45 | self.titles = [] |
|
46 | 46 | |
|
47 | Maintitle = "Range from %d km to %d km" %(int(self.data.yrange[0]),int(self.data.yrange[-1])) | |
|
47 | 48 | for n, ax in enumerate(self.axes): |
|
48 | 49 | ychannel = self.y[n,:] |
|
49 | 50 | if ax.firsttime: |
@@ -54,7 +55,7 class SpectraHeisPlot(Plot): | |||
|
54 | 55 | ax.plt.set_data(x, ychannel) |
|
55 | 56 | |
|
56 | 57 | self.titles.append("Channel {}: {:4.2f}dB".format(n, numpy.max(ychannel))) |
|
57 | ||
|
58 | plt.suptitle(Maintitle) | |
|
58 | 59 | |
|
59 | 60 | class RTIHeisPlot(Plot): |
|
60 | 61 | |
@@ -80,8 +81,8 class RTIHeisPlot(Plot): | |||
|
80 | 81 | spc = dataOut.data_spc / dataOut.normFactor |
|
81 | 82 | spc = 10*numpy.log10(numpy.average(spc, axis=1)) |
|
82 | 83 | data['rti_heis'] = spc |
|
83 | ||
|
84 |
return data, meta |
|
|
84 | ||
|
85 | return data, meta | |
|
85 | 86 | |
|
86 | 87 | def plot(self): |
|
87 | 88 |
@@ -22,11 +22,11 class ScopePlot(Plot): | |||
|
22 | 22 | def setup(self): |
|
23 | 23 | |
|
24 | 24 | self.xaxis = 'Range (Km)' |
|
25 | self.ncols = 1 | |
|
26 | self.nrows = 1 | |
|
27 | self.nplots = 1 | |
|
25 | self.nplots = len(self.data.channels) | |
|
26 | self.nrows = int(numpy.ceil(self.nplots/2)) | |
|
27 | self.ncols = int(numpy.ceil(self.nplots/self.nrows)) | |
|
28 | 28 | self.ylabel = 'Intensity [dB]' |
|
29 | self.titles = ['Scope'] | |
|
29 | self.titles = ['Channel '+str(self.data.channels[i])+" " for i in self.data.channels] | |
|
30 | 30 | self.colorbar = False |
|
31 | 31 | self.width = 6 |
|
32 | 32 | self.height = 4 |
@@ -56,15 +56,13 class ScopePlot(Plot): | |||
|
56 | 56 | |
|
57 | 57 | yreal = y[channelIndexList,:].real |
|
58 | 58 | yimag = y[channelIndexList,:].imag |
|
59 | title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y")) | |
|
59 | Maintitle = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y")) | |
|
60 | 60 | self.xlabel = "Range (Km)" |
|
61 | 61 | self.ylabel = "Intensity - IQ" |
|
62 | 62 | |
|
63 | 63 | self.y = yreal |
|
64 | 64 | self.x = x |
|
65 | 65 | |
|
66 | self.titles[0] = title | |
|
67 | ||
|
68 | 66 | for i,ax in enumerate(self.axes): |
|
69 | 67 | title = "Channel %d" %(i) |
|
70 | 68 | if ax.firsttime: |
@@ -75,6 +73,7 class ScopePlot(Plot): | |||
|
75 | 73 | else: |
|
76 | 74 | ax.plt_r.set_data(x, yreal[i,:]) |
|
77 | 75 | ax.plt_i.set_data(x, yimag[i,:]) |
|
76 | plt.suptitle(Maintitle) | |
|
78 | 77 | |
|
79 | 78 | def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle): |
|
80 | 79 | y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:]) |
@@ -99,6 +98,7 class ScopePlot(Plot): | |||
|
99 | 98 | else: |
|
100 | 99 | ax.plt_r.set_data(x, ychannel) |
|
101 | 100 | |
|
101 | ||
|
102 | 102 | def plot_weatherpower(self, x, y, channelIndexList, thisDatetime, wintitle): |
|
103 | 103 | |
|
104 | 104 |
@@ -384,7 +384,7 def isRadarFolder(folder): | |||
|
384 | 384 | |
|
385 | 385 | |
|
386 | 386 | def isRadarFile(file): |
|
387 |
try: |
|
|
387 | try: | |
|
388 | 388 | year = int(file[1:5]) |
|
389 | 389 | doy = int(file[5:8]) |
|
390 | 390 | set = int(file[8:11]) |
@@ -395,10 +395,10 def isRadarFile(file): | |||
|
395 | 395 | |
|
396 | 396 | |
|
397 | 397 | def getDateFromRadarFile(file): |
|
398 |
try: |
|
|
398 | try: | |
|
399 | 399 | year = int(file[1:5]) |
|
400 | 400 | doy = int(file[5:8]) |
|
401 |
set = int(file[8:11]) |
|
|
401 | set = int(file[8:11]) | |
|
402 | 402 | except: |
|
403 | 403 | return None |
|
404 | 404 | |
@@ -417,7 +417,7 def getDateFromRadarFolder(folder): | |||
|
417 | 417 | return thisDate |
|
418 | 418 | |
|
419 | 419 | def parse_format(s, fmt): |
|
420 | ||
|
420 | ||
|
421 | 421 | for i in range(fmt.count('%')): |
|
422 | 422 | x = fmt.index('%') |
|
423 | 423 | d = DT_DIRECTIVES[fmt[x:x+2]] |
@@ -484,7 +484,7 class Reader(object): | |||
|
484 | 484 | |
|
485 | 485 | def run(self): |
|
486 | 486 | |
|
487 |
raise NotImplementedError |
|
|
487 | raise NotImplementedError | |
|
488 | 488 | |
|
489 | 489 | def getAllowedArgs(self): |
|
490 | 490 | if hasattr(self, '__attrs__'): |
@@ -496,19 +496,19 class Reader(object): | |||
|
496 | 496 | |
|
497 | 497 | for key, value in kwargs.items(): |
|
498 | 498 | setattr(self, key, value) |
|
499 | ||
|
499 | ||
|
500 | 500 | def find_folders(self, path, startDate, endDate, folderfmt, last=False): |
|
501 | 501 | |
|
502 |
folders = [x for f in path.split(',') |
|
|
502 | folders = [x for f in path.split(',') | |
|
503 | 503 | for x in os.listdir(f) if os.path.isdir(os.path.join(f, x))] |
|
504 | 504 | folders.sort() |
|
505 | 505 | |
|
506 | 506 | if last: |
|
507 | 507 | folders = [folders[-1]] |
|
508 | 508 | |
|
509 |
for folder in folders: |
|
|
510 |
try: |
|
|
511 |
dt = datetime.datetime.strptime(parse_format(folder, folderfmt), folderfmt).date() |
|
|
509 | for folder in folders: | |
|
510 | try: | |
|
511 | dt = datetime.datetime.strptime(parse_format(folder, folderfmt), folderfmt).date() | |
|
512 | 512 | if dt >= startDate and dt <= endDate: |
|
513 | 513 | yield os.path.join(path, folder) |
|
514 | 514 | else: |
@@ -517,38 +517,38 class Reader(object): | |||
|
517 | 517 | log.log('Skiping folder {}'.format(folder), self.name) |
|
518 | 518 | continue |
|
519 | 519 | return |
|
520 | ||
|
521 |
def find_files(self, folders, ext, filefmt, startDate=None, endDate=None, |
|
|
520 | ||
|
521 | def find_files(self, folders, ext, filefmt, startDate=None, endDate=None, | |
|
522 | 522 | expLabel='', last=False): |
|
523 | ||
|
524 |
for path in folders: |
|
|
523 | ||
|
524 | for path in folders: | |
|
525 | 525 | files = glob.glob1(path, '*{}'.format(ext)) |
|
526 | 526 | files.sort() |
|
527 | 527 | if last: |
|
528 |
if files: |
|
|
528 | if files: | |
|
529 | 529 | fo = files[-1] |
|
530 |
try: |
|
|
530 | try: | |
|
531 | 531 | dt = datetime.datetime.strptime(parse_format(fo, filefmt), filefmt).date() |
|
532 |
yield os.path.join(path, expLabel, fo) |
|
|
533 |
except Exception as e: |
|
|
532 | yield os.path.join(path, expLabel, fo) | |
|
533 | except Exception as e: | |
|
534 | 534 | pass |
|
535 | 535 | return |
|
536 | 536 | else: |
|
537 | 537 | return |
|
538 | 538 | |
|
539 | 539 | for fo in files: |
|
540 |
try: |
|
|
541 |
dt = datetime.datetime.strptime(parse_format(fo, filefmt), filefmt).date() |
|
|
540 | try: | |
|
541 | dt = datetime.datetime.strptime(parse_format(fo, filefmt), filefmt).date() | |
|
542 | 542 | if dt >= startDate and dt <= endDate: |
|
543 | 543 | yield os.path.join(path, expLabel, fo) |
|
544 | 544 | else: |
|
545 | 545 | log.log('Skiping file {}'.format(fo), self.name) |
|
546 | 546 | except Exception as e: |
|
547 | 547 | log.log('Skiping file {}'.format(fo), self.name) |
|
548 |
continue |
|
|
548 | continue | |
|
549 | 549 | |
|
550 | 550 | def searchFilesOffLine(self, path, startDate, endDate, |
|
551 |
expLabel, ext, walk, |
|
|
551 | expLabel, ext, walk, | |
|
552 | 552 | filefmt, folderfmt): |
|
553 | 553 | """Search files in offline mode for the given arguments |
|
554 | 554 | |
@@ -561,12 +561,12 class Reader(object): | |||
|
561 | 561 | path, startDate, endDate, folderfmt) |
|
562 | 562 | else: |
|
563 | 563 | folders = path.split(',') |
|
564 | ||
|
564 | ||
|
565 | 565 | return self.find_files( |
|
566 |
folders, ext, filefmt, startDate, endDate, expLabel) |
|
|
566 | folders, ext, filefmt, startDate, endDate, expLabel) | |
|
567 | 567 | |
|
568 | 568 | def searchFilesOnLine(self, path, startDate, endDate, |
|
569 |
expLabel, ext, walk, |
|
|
569 | expLabel, ext, walk, | |
|
570 | 570 | filefmt, folderfmt): |
|
571 | 571 | """Search for the last file of the last folder |
|
572 | 572 | |
@@ -579,13 +579,13 class Reader(object): | |||
|
579 | 579 | Return: |
|
580 | 580 | generator with the full path of last filename |
|
581 | 581 | """ |
|
582 | ||
|
582 | ||
|
583 | 583 | if walk: |
|
584 | 584 | folders = self.find_folders( |
|
585 | 585 | path, startDate, endDate, folderfmt, last=True) |
|
586 | 586 | else: |
|
587 | 587 | folders = path.split(',') |
|
588 | ||
|
588 | ||
|
589 | 589 | return self.find_files( |
|
590 | 590 | folders, ext, filefmt, startDate, endDate, expLabel, last=True) |
|
591 | 591 | |
@@ -594,13 +594,13 class Reader(object): | |||
|
594 | 594 | |
|
595 | 595 | while True: |
|
596 | 596 | if self.fp != None: |
|
597 |
self.fp.close() |
|
|
597 | self.fp.close() | |
|
598 | 598 | |
|
599 | 599 | if self.online: |
|
600 | 600 | newFile = self.setNextFileOnline() |
|
601 | 601 | else: |
|
602 | 602 | newFile = self.setNextFileOffline() |
|
603 | ||
|
603 | ||
|
604 | 604 | if not(newFile): |
|
605 | 605 | if self.online: |
|
606 | 606 | raise schainpy.admin.SchainError('Time to wait for new files reach') |
@@ -609,10 +609,10 class Reader(object): | |||
|
609 | 609 | raise schainpy.admin.SchainWarning('No files found in the given path') |
|
610 | 610 | else: |
|
611 | 611 | raise schainpy.admin.SchainWarning('No more files to read') |
|
612 | ||
|
612 | ||
|
613 | 613 | if self.verifyFile(self.filename): |
|
614 | 614 | break |
|
615 | ||
|
615 | ||
|
616 | 616 | log.log('Opening file: %s' % self.filename, self.name) |
|
617 | 617 | |
|
618 | 618 | self.readFirstHeader() |
@@ -625,7 +625,7 class Reader(object): | |||
|
625 | 625 | self.filename |
|
626 | 626 | self.fp |
|
627 | 627 | self.filesize |
|
628 | ||
|
628 | ||
|
629 | 629 | Return: |
|
630 | 630 | boolean |
|
631 | 631 | |
@@ -633,7 +633,7 class Reader(object): | |||
|
633 | 633 | nextFile = True |
|
634 | 634 | nextDay = False |
|
635 | 635 | |
|
636 |
for nFiles in range(self.nFiles+1): |
|
|
636 | for nFiles in range(self.nFiles+1): | |
|
637 | 637 | for nTries in range(self.nTries): |
|
638 | 638 | fullfilename, filename = self.checkForRealPath(nextFile, nextDay) |
|
639 | 639 | if fullfilename is not None: |
@@ -643,18 +643,18 class Reader(object): | |||
|
643 | 643 | self.name) |
|
644 | 644 | time.sleep(self.delay) |
|
645 | 645 | nextFile = False |
|
646 |
continue |
|
|
647 | ||
|
646 | continue | |
|
647 | ||
|
648 | 648 | if fullfilename is not None: |
|
649 | 649 | break |
|
650 | ||
|
650 | ||
|
651 | 651 | self.nTries = 1 |
|
652 |
nextFile = True |
|
|
652 | nextFile = True | |
|
653 | 653 | |
|
654 | 654 | if nFiles == (self.nFiles - 1): |
|
655 | 655 | log.log('Trying with next day...', self.name) |
|
656 | 656 | nextDay = True |
|
657 |
self.nTries = 3 |
|
|
657 | self.nTries = 3 | |
|
658 | 658 | |
|
659 | 659 | if fullfilename: |
|
660 | 660 | self.fileSize = os.path.getsize(fullfilename) |
@@ -666,18 +666,18 class Reader(object): | |||
|
666 | 666 | self.flagNoMoreFiles = 0 |
|
667 | 667 | self.fileIndex += 1 |
|
668 | 668 | return 1 |
|
669 |
else: |
|
|
669 | else: | |
|
670 | 670 | return 0 |
|
671 | ||
|
671 | ||
|
672 | 672 | def setNextFileOffline(self): |
|
673 | 673 | """Open the next file to be readed in offline mode""" |
|
674 | ||
|
674 | ||
|
675 | 675 | try: |
|
676 | 676 | filename = next(self.filenameList) |
|
677 | 677 | self.fileIndex +=1 |
|
678 | 678 | except StopIteration: |
|
679 | 679 | self.flagNoMoreFiles = 1 |
|
680 |
return 0 |
|
|
680 | return 0 | |
|
681 | 681 | |
|
682 | 682 | self.filename = filename |
|
683 | 683 | self.fileSize = os.path.getsize(filename) |
@@ -685,22 +685,22 class Reader(object): | |||
|
685 | 685 | self.flagIsNewFile = 1 |
|
686 | 686 | |
|
687 | 687 | return 1 |
|
688 | ||
|
688 | ||
|
689 | 689 | @staticmethod |
|
690 | 690 | def isDateTimeInRange(dt, startDate, endDate, startTime, endTime): |
|
691 | 691 | """Check if the given datetime is in range""" |
|
692 | ||
|
693 | if startDate <= dt.date() <= endDate: | |
|
694 |
|
|
|
695 |
|
|
|
692 | startDateTime= datetime.datetime.combine(startDate,startTime) | |
|
693 | endDateTime = datetime.datetime.combine(endDate,endTime) | |
|
694 | if startDateTime <= dt <= endDateTime: | |
|
695 | return True | |
|
696 | 696 | return False |
|
697 | ||
|
697 | ||
|
698 | 698 | def verifyFile(self, filename): |
|
699 | 699 | """Check for a valid file |
|
700 | ||
|
700 | ||
|
701 | 701 | Arguments: |
|
702 | 702 | filename -- full path filename |
|
703 | ||
|
703 | ||
|
704 | 704 | Return: |
|
705 | 705 | boolean |
|
706 | 706 | """ |
@@ -711,7 +711,7 class Reader(object): | |||
|
711 | 711 | """Check if the next file to be readed exists""" |
|
712 | 712 | |
|
713 | 713 | raise NotImplementedError |
|
714 | ||
|
714 | ||
|
715 | 715 | def readFirstHeader(self): |
|
716 | 716 | """Parse the file header""" |
|
717 | 717 | |
@@ -783,8 +783,8 class JRODataReader(Reader): | |||
|
783 | 783 | Return: |
|
784 | 784 | str -- fullpath of the file |
|
785 | 785 | """ |
|
786 | ||
|
787 | ||
|
786 | ||
|
787 | ||
|
788 | 788 | if nextFile: |
|
789 | 789 | self.set += 1 |
|
790 | 790 | if nextDay: |
@@ -796,7 +796,7 class JRODataReader(Reader): | |||
|
796 | 796 | prefixFileList = ['d', 'D'] |
|
797 | 797 | elif self.ext.lower() == ".pdata": # spectra |
|
798 | 798 | prefixFileList = ['p', 'P'] |
|
799 | ||
|
799 | ||
|
800 | 800 | # barrido por las combinaciones posibles |
|
801 | 801 | for prefixDir in prefixDirList: |
|
802 | 802 | thispath = self.path |
@@ -816,9 +816,9 class JRODataReader(Reader): | |||
|
816 | 816 | |
|
817 | 817 | if os.path.exists(fullfilename): |
|
818 | 818 | return fullfilename, filename |
|
819 | ||
|
820 |
return None, filename |
|
|
821 | ||
|
819 | ||
|
820 | return None, filename | |
|
821 | ||
|
822 | 822 | def __waitNewBlock(self): |
|
823 | 823 | """ |
|
824 | 824 | Return 1 si se encontro un nuevo bloque de datos, 0 de otra forma. |
@@ -860,9 +860,9 class JRODataReader(Reader): | |||
|
860 | 860 | def __setNewBlock(self): |
|
861 | 861 | |
|
862 | 862 | if self.fp == None: |
|
863 |
return 0 |
|
|
864 | ||
|
865 |
if self.flagIsNewFile: |
|
|
863 | return 0 | |
|
864 | ||
|
865 | if self.flagIsNewFile: | |
|
866 | 866 | self.lastUTTime = self.basicHeaderObj.utc |
|
867 | 867 | return 1 |
|
868 | 868 | |
@@ -875,12 +875,12 class JRODataReader(Reader): | |||
|
875 | 875 | |
|
876 | 876 | currentSize = self.fileSize - self.fp.tell() |
|
877 | 877 | neededSize = self.processingHeaderObj.blockSize + self.basicHeaderSize |
|
878 | ||
|
878 | ||
|
879 | 879 | if (currentSize >= neededSize): |
|
880 | 880 | self.basicHeaderObj.read(self.fp) |
|
881 | 881 | self.lastUTTime = self.basicHeaderObj.utc |
|
882 | 882 | return 1 |
|
883 | ||
|
883 | ||
|
884 | 884 | if self.__waitNewBlock(): |
|
885 | 885 | self.lastUTTime = self.basicHeaderObj.utc |
|
886 | 886 | return 1 |
@@ -966,10 +966,10 class JRODataReader(Reader): | |||
|
966 | 966 | except IOError: |
|
967 | 967 | log.error("File {} can't be opened".format(filename), self.name) |
|
968 | 968 | return False |
|
969 | ||
|
969 | ||
|
970 | 970 | if self.online and self.waitDataBlock(0): |
|
971 | 971 | pass |
|
972 | ||
|
972 | ||
|
973 | 973 | basicHeaderObj = BasicHeader(LOCALTIME) |
|
974 | 974 | systemHeaderObj = SystemHeader() |
|
975 | 975 | radarControllerHeaderObj = RadarControllerHeader() |
@@ -996,7 +996,7 class JRODataReader(Reader): | |||
|
996 | 996 | dt2 = basicHeaderObj.datatime |
|
997 | 997 | if not self.isDateTimeInRange(dt1, self.startDate, self.endDate, self.startTime, self.endTime) and not \ |
|
998 | 998 | self.isDateTimeInRange(dt2, self.startDate, self.endDate, self.startTime, self.endTime): |
|
999 |
flag = False |
|
|
999 | flag = False | |
|
1000 | 1000 | |
|
1001 | 1001 | fp.close() |
|
1002 | 1002 | return flag |
@@ -1105,11 +1105,11 class JRODataReader(Reader): | |||
|
1105 | 1105 | return dateList |
|
1106 | 1106 | |
|
1107 | 1107 | def setup(self, **kwargs): |
|
1108 | ||
|
1108 | ||
|
1109 | 1109 | self.set_kwargs(**kwargs) |
|
1110 | 1110 | if not self.ext.startswith('.'): |
|
1111 | 1111 | self.ext = '.{}'.format(self.ext) |
|
1112 | ||
|
1112 | ||
|
1113 | 1113 | if self.server is not None: |
|
1114 | 1114 | if 'tcp://' in self.server: |
|
1115 | 1115 | address = server |
@@ -1131,36 +1131,36 class JRODataReader(Reader): | |||
|
1131 | 1131 | |
|
1132 | 1132 | for nTries in range(self.nTries): |
|
1133 | 1133 | fullpath = self.searchFilesOnLine(self.path, self.startDate, |
|
1134 |
self.endDate, self.expLabel, self.ext, self.walk, |
|
|
1134 | self.endDate, self.expLabel, self.ext, self.walk, | |
|
1135 | 1135 | self.filefmt, self.folderfmt) |
|
1136 | 1136 | |
|
1137 | 1137 | try: |
|
1138 | 1138 | fullpath = next(fullpath) |
|
1139 | 1139 | except: |
|
1140 | 1140 | fullpath = None |
|
1141 | ||
|
1141 | ||
|
1142 | 1142 | if fullpath: |
|
1143 | 1143 | break |
|
1144 | 1144 | |
|
1145 | 1145 | log.warning( |
|
1146 | 1146 | 'Waiting {} sec for a valid file in {}: try {} ...'.format( |
|
1147 |
self.delay, self.path, nTries + 1), |
|
|
1147 | self.delay, self.path, nTries + 1), | |
|
1148 | 1148 | self.name) |
|
1149 | 1149 | time.sleep(self.delay) |
|
1150 | 1150 | |
|
1151 | 1151 | if not(fullpath): |
|
1152 | 1152 | raise schainpy.admin.SchainError( |
|
1153 |
'There isn\'t any valid file in {}'.format(self.path)) |
|
|
1153 | 'There isn\'t any valid file in {}'.format(self.path)) | |
|
1154 | 1154 | |
|
1155 | 1155 | pathname, filename = os.path.split(fullpath) |
|
1156 | 1156 | self.year = int(filename[1:5]) |
|
1157 | 1157 | self.doy = int(filename[5:8]) |
|
1158 |
self.set = int(filename[8:11]) - 1 |
|
|
1158 | self.set = int(filename[8:11]) - 1 | |
|
1159 | 1159 | else: |
|
1160 | 1160 | log.log("Searching files in {}".format(self.path), self.name) |
|
1161 |
self.filenameList = self.searchFilesOffLine(self.path, self.startDate, |
|
|
1161 | self.filenameList = self.searchFilesOffLine(self.path, self.startDate, | |
|
1162 | 1162 | self.endDate, self.expLabel, self.ext, self.walk, self.filefmt, self.folderfmt) |
|
1163 | ||
|
1163 | ||
|
1164 | 1164 | self.setNextFile() |
|
1165 | 1165 | |
|
1166 | 1166 | return |
@@ -1181,7 +1181,7 class JRODataReader(Reader): | |||
|
1181 | 1181 | self.dataOut.useLocalTime = self.basicHeaderObj.useLocalTime |
|
1182 | 1182 | |
|
1183 | 1183 | self.dataOut.ippSeconds = self.radarControllerHeaderObj.ippSeconds / self.nTxs |
|
1184 | ||
|
1184 | ||
|
1185 | 1185 | def getFirstHeader(self): |
|
1186 | 1186 | |
|
1187 | 1187 | raise NotImplementedError |
@@ -1214,8 +1214,8 class JRODataReader(Reader): | |||
|
1214 | 1214 | """ |
|
1215 | 1215 | |
|
1216 | 1216 | Arguments: |
|
1217 |
path : |
|
|
1218 |
startDate : |
|
|
1217 | path : | |
|
1218 | startDate : | |
|
1219 | 1219 | endDate : |
|
1220 | 1220 | startTime : |
|
1221 | 1221 | endTime : |
@@ -1284,7 +1284,7 class JRODataWriter(Reader): | |||
|
1284 | 1284 | dtype_width = get_dtype_width(dtype_index) |
|
1285 | 1285 | |
|
1286 | 1286 | return dtype_width |
|
1287 | ||
|
1287 | ||
|
1288 | 1288 | def getProcessFlags(self): |
|
1289 | 1289 | |
|
1290 | 1290 | processFlags = 0 |
@@ -1322,9 +1322,9 class JRODataWriter(Reader): | |||
|
1322 | 1322 | |
|
1323 | 1323 | self.basicHeaderObj.size = self.basicHeaderSize # bytes |
|
1324 | 1324 | self.basicHeaderObj.version = self.versionFile |
|
1325 |
self.basicHeaderObj.dataBlock = self.nTotalBlocks |
|
|
1325 | self.basicHeaderObj.dataBlock = self.nTotalBlocks | |
|
1326 | 1326 | utc = numpy.floor(self.dataOut.utctime) |
|
1327 |
milisecond = (self.dataOut.utctime - utc) * 1000.0 |
|
|
1327 | milisecond = (self.dataOut.utctime - utc) * 1000.0 | |
|
1328 | 1328 | self.basicHeaderObj.utc = utc |
|
1329 | 1329 | self.basicHeaderObj.miliSecond = milisecond |
|
1330 | 1330 | self.basicHeaderObj.timeZone = self.dataOut.timeZone |
@@ -1465,9 +1465,9 class JRODataWriter(Reader): | |||
|
1465 | 1465 | if self.dataOut.datatime.date() > self.fileDate: |
|
1466 | 1466 | setFile = 0 |
|
1467 | 1467 | self.nTotalBlocks = 0 |
|
1468 | ||
|
1468 | ||
|
1469 | 1469 | filen = '{}{:04d}{:03d}{:03d}{}'.format( |
|
1470 |
self.optchar, timeTuple.tm_year, timeTuple.tm_yday, setFile, ext) |
|
|
1470 | self.optchar, timeTuple.tm_year, timeTuple.tm_yday, setFile, ext) | |
|
1471 | 1471 | |
|
1472 | 1472 | filename = os.path.join(path, subfolder, filen) |
|
1473 | 1473 | |
@@ -1515,11 +1515,11 class JRODataWriter(Reader): | |||
|
1515 | 1515 | self.ext = ext.lower() |
|
1516 | 1516 | |
|
1517 | 1517 | self.path = path |
|
1518 | ||
|
1518 | ||
|
1519 | 1519 | if set is None: |
|
1520 | 1520 | self.setFile = -1 |
|
1521 | 1521 | else: |
|
1522 |
self.setFile = set - 1 |
|
|
1522 | self.setFile = set - 1 | |
|
1523 | 1523 | |
|
1524 | 1524 | self.blocksPerFile = blocksPerFile |
|
1525 | 1525 | self.profilesPerBlock = profilesPerBlock |
@@ -347,7 +347,7 class AMISRReader(ProcessingUnit): | |||
|
347 | 347 | |
|
348 | 348 | else: |
|
349 | 349 | #get the last file - 1 |
|
350 |
self.filenameList = [self.filenameList[- |
|
|
350 | self.filenameList = [self.filenameList[-2]] | |
|
351 | 351 | new_dirnameList = [] |
|
352 | 352 | for dirname in self.dirnameList: |
|
353 | 353 | junk = numpy.array([dirname in x for x in self.filenameList]) |
@@ -421,7 +421,7 class AMISRReader(ProcessingUnit): | |||
|
421 | 421 | self.__selectDataForTimes(online=True) |
|
422 | 422 | filename = self.filenameList[0] |
|
423 | 423 | wait = 0 |
|
424 |
|
|
|
424 | self.__waitForNewFile=300 ## DEBUG: | |
|
425 | 425 | while self.__filename_online == filename: |
|
426 | 426 | print('waiting %d seconds to get a new file...'%(self.__waitForNewFile)) |
|
427 | 427 | if wait == 5: |
@@ -641,13 +641,14 class AMISRReader(ProcessingUnit): | |||
|
641 | 641 | #self.dataOut.utctime = self.timeset[self.profileIndex] + t_comp |
|
642 | 642 | #print(self.dataOut.utctime) |
|
643 | 643 | self.dataOut.profileIndex = self.profileIndex |
|
644 | #print("N profile:",self.profileIndex,self.newProfiles,self.nblocks,self.dataOut.utctime) | |
|
644 | 645 | self.dataOut.flagNoData = False |
|
645 | 646 | # if indexprof == 0: |
|
646 | 647 | # print self.dataOut.utctime |
|
647 | 648 | |
|
648 | 649 | self.profileIndex += 1 |
|
649 | 650 | |
|
650 | return self.dataOut.data | |
|
651 | #return self.dataOut.data | |
|
651 | 652 | |
|
652 | 653 | |
|
653 | 654 | def run(self, **kwargs): |
@@ -660,3 +661,4 class AMISRReader(ProcessingUnit): | |||
|
660 | 661 | self.isConfig = True |
|
661 | 662 | |
|
662 | 663 | self.getData() |
|
664 | return(self.dataOut.data) |
@@ -5,7 +5,6 from schainpy.model.data.jrodata import SpectraHeis | |||
|
5 | 5 | from schainpy.utils import log |
|
6 | 6 | |
|
7 | 7 | |
|
8 | ||
|
9 | 8 | class SpectraHeisProc(ProcessingUnit): |
|
10 | 9 | |
|
11 | 10 | def __init__(self):#, **kwargs): |
@@ -89,7 +88,7 class SpectraHeisProc(ProcessingUnit): | |||
|
89 | 88 | if self.dataIn.type == "Fits": |
|
90 | 89 | self.__updateObjFromFits() |
|
91 | 90 | self.dataOut.flagNoData = False |
|
92 |
return |
|
|
91 | return | |
|
93 | 92 | |
|
94 | 93 | if self.dataIn.type == "SpectraHeis": |
|
95 | 94 | self.dataOut.copy(self.dataIn) |
@@ -343,5 +342,5 class IncohInt4SpectraHeis(Operation): | |||
|
343 | 342 | # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nIncohInt |
|
344 | 343 | # dataOut.timeInterval = self.__timeInterval*self.n |
|
345 | 344 | dataOut.flagNoData = False |
|
346 | ||
|
347 | return dataOut No newline at end of file | |
|
345 | ||
|
346 | return dataOut |
This diff has been collapsed as it changes many lines, (559 lines changed) Show them Hide them | |||
@@ -12,12 +12,15 import time | |||
|
12 | 12 | import itertools |
|
13 | 13 | |
|
14 | 14 | import numpy |
|
15 | import math | |
|
15 | 16 | |
|
16 | 17 | from schainpy.model.proc.jroproc_base import ProcessingUnit, MPDecorator, Operation |
|
17 | 18 | from schainpy.model.data.jrodata import Spectra |
|
18 | 19 | from schainpy.model.data.jrodata import hildebrand_sekhon |
|
19 | 20 | from schainpy.utils import log |
|
20 | 21 | |
|
22 | from scipy.optimize import curve_fit | |
|
23 | ||
|
21 | 24 | |
|
22 | 25 | class SpectraProc(ProcessingUnit): |
|
23 | 26 | |
@@ -478,6 +481,562 class removeDC(Operation): | |||
|
478 | 481 | |
|
479 | 482 | return self.dataOut |
|
480 | 483 | |
|
484 | # import matplotlib.pyplot as plt | |
|
485 | ||
|
486 | def fit_func( x, a0, a1, a2): #, a3, a4, a5): | |
|
487 | z = (x - a1) / a2 | |
|
488 | y = a0 * numpy.exp(-z**2 / a2) #+ a3 + a4 * x + a5 * x**2 | |
|
489 | return y | |
|
490 | class CleanRayleigh(Operation): | |
|
491 | ||
|
492 | def __init__(self): | |
|
493 | ||
|
494 | Operation.__init__(self) | |
|
495 | self.i=0 | |
|
496 | self.isConfig = False | |
|
497 | self.__dataReady = False | |
|
498 | self.__profIndex = 0 | |
|
499 | self.byTime = False | |
|
500 | self.byProfiles = False | |
|
501 | ||
|
502 | self.bloques = None | |
|
503 | self.bloque0 = None | |
|
504 | ||
|
505 | self.index = 0 | |
|
506 | ||
|
507 | self.buffer = 0 | |
|
508 | self.buffer2 = 0 | |
|
509 | self.buffer3 = 0 | |
|
510 | #self.min_hei = None | |
|
511 | #self.max_hei = None | |
|
512 | ||
|
513 | def setup(self,dataOut,min_hei,max_hei,n, timeInterval,factor_stdv): | |
|
514 | ||
|
515 | self.nChannels = dataOut.nChannels | |
|
516 | self.nProf = dataOut.nProfiles | |
|
517 | self.nPairs = dataOut.data_cspc.shape[0] | |
|
518 | self.pairsArray = numpy.array(dataOut.pairsList) | |
|
519 | self.spectra = dataOut.data_spc | |
|
520 | self.cspectra = dataOut.data_cspc | |
|
521 | self.heights = dataOut.heightList #alturas totales | |
|
522 | self.nHeights = len(self.heights) | |
|
523 | self.min_hei = min_hei | |
|
524 | self.max_hei = max_hei | |
|
525 | if (self.min_hei == None): | |
|
526 | self.min_hei = 0 | |
|
527 | if (self.max_hei == None): | |
|
528 | self.max_hei = dataOut.heightList[-1] | |
|
529 | self.hval = ((self.max_hei>=self.heights) & (self.heights >= self.min_hei)).nonzero() | |
|
530 | self.heightsClean = self.heights[self.hval] #alturas filtradas | |
|
531 | self.hval = self.hval[0] # forma (N,), an solo N elementos -> Indices de alturas | |
|
532 | self.nHeightsClean = len(self.heightsClean) | |
|
533 | self.channels = dataOut.channelList | |
|
534 | self.nChan = len(self.channels) | |
|
535 | self.nIncohInt = dataOut.nIncohInt | |
|
536 | self.__initime = dataOut.utctime | |
|
537 | self.maxAltInd = self.hval[-1]+1 | |
|
538 | self.minAltInd = self.hval[0] | |
|
539 | ||
|
540 | self.crosspairs = dataOut.pairsList | |
|
541 | self.nPairs = len(self.crosspairs) | |
|
542 | self.normFactor = dataOut.normFactor | |
|
543 | self.nFFTPoints = dataOut.nFFTPoints | |
|
544 | self.ippSeconds = dataOut.ippSeconds | |
|
545 | self.currentTime = self.__initime | |
|
546 | self.pairsArray = numpy.array(dataOut.pairsList) | |
|
547 | self.factor_stdv = factor_stdv | |
|
548 | ||
|
549 | ||
|
550 | ||
|
551 | if n != None : | |
|
552 | self.byProfiles = True | |
|
553 | self.nIntProfiles = n | |
|
554 | else: | |
|
555 | self.__integrationtime = timeInterval | |
|
556 | ||
|
557 | self.__dataReady = False | |
|
558 | self.isConfig = True | |
|
559 | ||
|
560 | ||
|
561 | ||
|
562 | def run(self, dataOut,min_hei=None,max_hei=None, n=None, timeInterval=10,factor_stdv=2.5): | |
|
563 | #print (dataOut.utctime) | |
|
564 | if not self.isConfig : | |
|
565 | #print("Setting config") | |
|
566 | self.setup(dataOut, min_hei,max_hei,n,timeInterval,factor_stdv) | |
|
567 | #print("Config Done") | |
|
568 | tini=dataOut.utctime | |
|
569 | ||
|
570 | if self.byProfiles: | |
|
571 | if self.__profIndex == self.nIntProfiles: | |
|
572 | self.__dataReady = True | |
|
573 | else: | |
|
574 | if (tini - self.__initime) >= self.__integrationtime: | |
|
575 | #print(tini - self.__initime,self.__profIndex) | |
|
576 | self.__dataReady = True | |
|
577 | self.__initime = tini | |
|
578 | ||
|
579 | #if (tini.tm_min % 2) == 0 and (tini.tm_sec < 5 and self.fint==0): | |
|
580 | ||
|
581 | if self.__dataReady: | |
|
582 | print("Data ready",self.__profIndex) | |
|
583 | self.__profIndex = 0 | |
|
584 | jspc = self.buffer | |
|
585 | jcspc = self.buffer2 | |
|
586 | #jnoise = self.buffer3 | |
|
587 | self.buffer = dataOut.data_spc | |
|
588 | self.buffer2 = dataOut.data_cspc | |
|
589 | #self.buffer3 = dataOut.noise | |
|
590 | self.currentTime = dataOut.utctime | |
|
591 | if numpy.any(jspc) : | |
|
592 | #print( jspc.shape, jcspc.shape) | |
|
593 | jspc = numpy.reshape(jspc,(int(len(jspc)/self.nChannels),self.nChannels,self.nFFTPoints,self.nHeights)) | |
|
594 | jcspc= numpy.reshape(jcspc,(int(len(jcspc)/self.nPairs),self.nPairs,self.nFFTPoints,self.nHeights)) | |
|
595 | self.__dataReady = False | |
|
596 | #print( jspc.shape, jcspc.shape) | |
|
597 | dataOut.flagNoData = False | |
|
598 | else: | |
|
599 | dataOut.flagNoData = True | |
|
600 | self.__dataReady = False | |
|
601 | return dataOut | |
|
602 | else: | |
|
603 | #print( len(self.buffer)) | |
|
604 | if numpy.any(self.buffer): | |
|
605 | self.buffer = numpy.concatenate((self.buffer,dataOut.data_spc), axis=0) | |
|
606 | self.buffer2 = numpy.concatenate((self.buffer2,dataOut.data_cspc), axis=0) | |
|
607 | self.buffer3 += dataOut.data_dc | |
|
608 | else: | |
|
609 | self.buffer = dataOut.data_spc | |
|
610 | self.buffer2 = dataOut.data_cspc | |
|
611 | self.buffer3 = dataOut.data_dc | |
|
612 | #print self.index, self.fint | |
|
613 | #print self.buffer2.shape | |
|
614 | dataOut.flagNoData = True ## NOTE: ?? revisar LUEGO | |
|
615 | self.__profIndex += 1 | |
|
616 | return dataOut ## NOTE: REV | |
|
617 | ||
|
618 | # if self.index == 0 and self.fint == 1 : | |
|
619 | # if jspc != None: | |
|
620 | # print len(jspc), jspc.shape | |
|
621 | # jspc= numpy.reshape(jspc,(4,128,63,len(jspc)/4)) | |
|
622 | # print jspc.shape | |
|
623 | # dataOut.flagNoData = True | |
|
624 | # return dataOut | |
|
625 | # if path != None: | |
|
626 | # sys.path.append(path) | |
|
627 | # self.library = importlib.import_module(file) | |
|
628 | # | |
|
629 | # To be inserted as a parameter | |
|
630 | ||
|
631 | #Set constants | |
|
632 | #constants = self.library.setConstants(dataOut) | |
|
633 | #dataOut.constants = constants | |
|
634 | ||
|
635 | #snrth= 20 | |
|
636 | ||
|
637 | #crosspairs = dataOut.groupList | |
|
638 | #noise = dataOut.noise | |
|
639 | #print( nProf,heights) | |
|
640 | #print( jspc.shape, jspc.shape[0]) | |
|
641 | #print noise | |
|
642 | #print jnoise[len(jnoise)-1,:], numpy.nansum(jnoise,axis=0)/len(jnoise) | |
|
643 | #jnoise = jnoise/N | |
|
644 | #noise = numpy.nansum(jnoise,axis=0)#/len(jnoise) | |
|
645 | #print( noise) | |
|
646 | #power = numpy.sum(spectra, axis=1) | |
|
647 | #print power[0,:] | |
|
648 | #print("CROSSPAIRS",crosspairs) | |
|
649 | #nPairs = len(crosspairs) | |
|
650 | #print(numpy.shape(dataOut.data_spc)) | |
|
651 | ||
|
652 | #absc = dataOut.abscissaList[:-1] | |
|
653 | ||
|
654 | #print absc.shape | |
|
655 | #nBlocks=149 | |
|
656 | #print('spectra', spectra.shape) | |
|
657 | #print('noise print', crosspairs) | |
|
658 | #print('spectra', spectra.shape) | |
|
659 | #print('cspectra', cspectra.shape) | |
|
660 | #print numpy.array(dataOut.data_pre[1]).shape | |
|
661 | #spec, cspec = self.__DiffCoherent(snrth, spectra, cspectra, nProf, heights,nChan, nHei, nPairs, channels, noise*nProf, crosspairs) | |
|
662 | ||
|
663 | ||
|
664 | #index = tini.tm_hour*12+tini.tm_min/5 | |
|
665 | # jspc = jspc/self.nFFTPoints/self.normFactor | |
|
666 | # jcspc = jcspc/self.nFFTPoints/self.normFactor | |
|
667 | ||
|
668 | ||
|
669 | ||
|
670 | #dataOut.data_spc,dataOut.data_cspc = self.CleanRayleigh(dataOut,jspc,jcspc,crosspairs,heights,channels,nProf,nHei,nChan,nPairs,nIncohInt,nBlocks=nBlocks) | |
|
671 | #tmp_spectra,tmp_cspectra,sat_spectra,sat_cspectra = self.cleanRayleigh(dataOut,jspc,jcspc,self.min_hei,self.max_hei) | |
|
672 | tmp_spectra,tmp_cspectra = self.cleanRayleigh(dataOut,jspc,jcspc,self.factor_stdv) | |
|
673 | #jspectra = tmp_spectra*len(jspc[:,0,0,0]) | |
|
674 | #jcspectra = tmp_cspectra*len(jspc[:,0,0,0]) | |
|
675 | ||
|
676 | dataOut.data_spc = tmp_spectra | |
|
677 | dataOut.data_cspc = tmp_cspectra | |
|
678 | dataOut.data_dc = self.buffer3 | |
|
679 | dataOut.nIncohInt *= self.nIntProfiles | |
|
680 | dataOut.utctime = self.currentTime #tiempo promediado | |
|
681 | print("Time: ",time.localtime(dataOut.utctime)) | |
|
682 | # dataOut.data_spc = sat_spectra | |
|
683 | # dataOut.data_cspc = sat_cspectra | |
|
684 | self.buffer = 0 | |
|
685 | self.buffer2 = 0 | |
|
686 | self.buffer3 = 0 | |
|
687 | ||
|
688 | return dataOut | |
|
689 | ||
|
690 | def cleanRayleigh(self,dataOut,spectra,cspectra,factor_stdv): | |
|
691 | print("OP cleanRayleigh") | |
|
692 | #import matplotlib.pyplot as plt | |
|
693 | #for k in range(149): | |
|
694 | ||
|
695 | rfunc = cspectra.copy() #self.bloques | |
|
696 | val_spc = spectra*0.0 #self.bloque0*0.0 | |
|
697 | val_cspc = cspectra*0.0 #self.bloques*0.0 | |
|
698 | in_sat_spectra = spectra.copy() #self.bloque0 | |
|
699 | in_sat_cspectra = cspectra.copy() #self.bloques | |
|
700 | ||
|
701 | raxs = math.ceil(math.sqrt(self.nPairs)) | |
|
702 | caxs = math.ceil(self.nPairs/raxs) | |
|
703 | ||
|
704 | #print(self.hval) | |
|
705 | #print numpy.absolute(rfunc[:,0,0,14]) | |
|
706 | for ih in range(self.minAltInd,self.maxAltInd): | |
|
707 | for ifreq in range(self.nFFTPoints): | |
|
708 | # fig, axs = plt.subplots(raxs, caxs) | |
|
709 | # fig2, axs2 = plt.subplots(raxs, caxs) | |
|
710 | col_ax = 0 | |
|
711 | row_ax = 0 | |
|
712 | for ii in range(self.nPairs): #PARES DE CANALES SELF y CROSS | |
|
713 | #print("ii: ",ii) | |
|
714 | if (col_ax%caxs==0 and col_ax!=0): | |
|
715 | col_ax = 0 | |
|
716 | row_ax += 1 | |
|
717 | func2clean = 10*numpy.log10(numpy.absolute(rfunc[:,ii,ifreq,ih])) #Potencia? | |
|
718 | #print(func2clean.shape) | |
|
719 | val = (numpy.isfinite(func2clean)==True).nonzero() | |
|
720 | ||
|
721 | if len(val)>0: #limitador | |
|
722 | min_val = numpy.around(numpy.amin(func2clean)-2) #> (-40) | |
|
723 | if min_val <= -40 : | |
|
724 | min_val = -40 | |
|
725 | max_val = numpy.around(numpy.amax(func2clean)+2) #< 200 | |
|
726 | if max_val >= 200 : | |
|
727 | max_val = 200 | |
|
728 | #print min_val, max_val | |
|
729 | step = 1 | |
|
730 | #print("Getting bins and the histogram") | |
|
731 | x_dist = min_val + numpy.arange(1 + ((max_val-(min_val))/step))*step | |
|
732 | y_dist,binstep = numpy.histogram(func2clean,bins=range(int(min_val),int(max_val+2),step)) | |
|
733 | #print(len(y_dist),len(binstep[:-1])) | |
|
734 | #print(row_ax,col_ax, " ..") | |
|
735 | #print(self.pairsArray[ii][0],self.pairsArray[ii][1]) | |
|
736 | mean = numpy.sum(x_dist * y_dist) / numpy.sum(y_dist) | |
|
737 | sigma = numpy.sqrt(numpy.sum(y_dist * (x_dist - mean)**2) / numpy.sum(y_dist)) | |
|
738 | parg = [numpy.amax(y_dist),mean,sigma] | |
|
739 | gauss_fit, covariance = None, None | |
|
740 | newY = None | |
|
741 | try : | |
|
742 | gauss_fit, covariance = curve_fit(fit_func, x_dist, y_dist,p0=parg) | |
|
743 | mode = gauss_fit[1] | |
|
744 | stdv = gauss_fit[2] | |
|
745 | #print(" FIT OK",gauss_fit) | |
|
746 | ''' | |
|
747 | newY = fit_func(x_dist,gauss_fit[0],gauss_fit[1],gauss_fit[2]) | |
|
748 | axs[row_ax,col_ax].plot(binstep[:-1],y_dist,color='green') | |
|
749 | axs[row_ax,col_ax].plot(binstep[:-1],newY,color='red') | |
|
750 | axs[row_ax,col_ax].set_title("Pair "+str(self.crosspairs[ii]))''' | |
|
751 | except: | |
|
752 | mode = mean | |
|
753 | stdv = sigma | |
|
754 | #print("FIT FAIL") | |
|
755 | ||
|
756 | ||
|
757 | #print(mode,stdv) | |
|
758 | #Removing echoes greater than mode + 3*stdv | |
|
759 | #factor_stdv = 2 | |
|
760 | noval = (abs(func2clean - mode)>=(factor_stdv*stdv)).nonzero() | |
|
761 | #noval tiene los indices que se van a remover | |
|
762 | #print("Pair ",ii," novals: ",len(noval[0])) | |
|
763 | if len(noval[0]) > 0: #forma de array (N,) es igual a longitud (N) | |
|
764 | novall = ((func2clean - mode) >= (factor_stdv*stdv)).nonzero() | |
|
765 | #print(novall) | |
|
766 | #print(" ") | |
|
767 | cross_pairs = self.pairsArray[ii] | |
|
768 | #Getting coherent echoes which are removed. | |
|
769 | if len(novall[0]) > 0: | |
|
770 | ||
|
771 | val_spc[novall[0],cross_pairs[0],ifreq,ih] = 1 | |
|
772 | val_spc[novall[0],cross_pairs[1],ifreq,ih] = 1 | |
|
773 | val_cspc[novall[0],ii,ifreq,ih] = 1 | |
|
774 | #print("OUT NOVALL 1") | |
|
775 | #Removing coherent from ISR data | |
|
776 | ||
|
777 | #print(spectra[:,ii,ifreq,ih]) | |
|
778 | new_a = numpy.delete(cspectra[:,ii,ifreq,ih], noval[0]) | |
|
779 | mean_cspc = numpy.mean(new_a) | |
|
780 | new_b = numpy.delete(spectra[:,cross_pairs[0],ifreq,ih], noval[0]) | |
|
781 | mean_spc0 = numpy.mean(new_b) | |
|
782 | new_c = numpy.delete(spectra[:,cross_pairs[1],ifreq,ih], noval[0]) | |
|
783 | mean_spc1 = numpy.mean(new_c) | |
|
784 | spectra[noval,cross_pairs[0],ifreq,ih] = mean_spc0 | |
|
785 | spectra[noval,cross_pairs[1],ifreq,ih] = mean_spc1 | |
|
786 | cspectra[noval,ii,ifreq,ih] = mean_cspc | |
|
787 | ||
|
788 | ''' | |
|
789 | func2clean = 10*numpy.log10(numpy.absolute(cspectra[:,ii,ifreq,ih])) | |
|
790 | y_dist,binstep = numpy.histogram(func2clean,bins=range(int(min_val),int(max_val+2),step)) | |
|
791 | axs2[row_ax,col_ax].plot(binstep[:-1],newY,color='red') | |
|
792 | axs2[row_ax,col_ax].plot(binstep[:-1],y_dist,color='green') | |
|
793 | axs2[row_ax,col_ax].set_title("Pair "+str(self.crosspairs[ii])) | |
|
794 | ''' | |
|
795 | ||
|
796 | col_ax += 1 #contador de ploteo columnas | |
|
797 | ##print(col_ax) | |
|
798 | ''' | |
|
799 | title = str(dataOut.datatime)+" nFFT: "+str(ifreq)+" Alt: "+str(self.heights[ih])+ " km" | |
|
800 | title2 = str(dataOut.datatime)+" nFFT: "+str(ifreq)+" Alt: "+str(self.heights[ih])+ " km CLEANED" | |
|
801 | fig.suptitle(title) | |
|
802 | fig2.suptitle(title2) | |
|
803 | plt.show()''' | |
|
804 | ||
|
805 | ''' channels = channels | |
|
806 | cross_pairs = cross_pairs | |
|
807 | #print("OUT NOVALL 2") | |
|
808 | ||
|
809 | vcross0 = (cross_pairs[0] == channels[ii]).nonzero() | |
|
810 | vcross1 = (cross_pairs[1] == channels[ii]).nonzero() | |
|
811 | vcross = numpy.concatenate((vcross0,vcross1),axis=None) | |
|
812 | #print('vcros =', vcross) | |
|
813 | ||
|
814 | #Getting coherent echoes which are removed. | |
|
815 | if len(novall) > 0: | |
|
816 | #val_spc[novall,ii,ifreq,ih] = 1 | |
|
817 | val_spc[ii,ifreq,ih,novall] = 1 | |
|
818 | if len(vcross) > 0: | |
|
819 | val_cspc[vcross,ifreq,ih,novall] = 1 | |
|
820 | ||
|
821 | #Removing coherent from ISR data. | |
|
822 | self.bloque0[ii,ifreq,ih,noval] = numpy.nan | |
|
823 | if len(vcross) > 0: | |
|
824 | self.bloques[vcross,ifreq,ih,noval] = numpy.nan | |
|
825 | ''' | |
|
826 | ||
|
827 | print("Getting average of the spectra and cross-spectra from incoherent echoes.") | |
|
828 | out_spectra = numpy.zeros([self.nChan,self.nFFTPoints,self.nHeights], dtype=float) #+numpy.nan | |
|
829 | out_cspectra = numpy.zeros([self.nPairs,self.nFFTPoints,self.nHeights], dtype=complex) #+numpy.nan | |
|
830 | for ih in range(self.nHeights): | |
|
831 | for ifreq in range(self.nFFTPoints): | |
|
832 | for ich in range(self.nChan): | |
|
833 | tmp = spectra[:,ich,ifreq,ih] | |
|
834 | valid = (numpy.isfinite(tmp[:])==True).nonzero() | |
|
835 | # if ich == 0 and ifreq == 0 and ih == 17 : | |
|
836 | # print tmp | |
|
837 | # print valid | |
|
838 | # print len(valid[0]) | |
|
839 | #print('TMP',tmp) | |
|
840 | if len(valid[0]) >0 : | |
|
841 | out_spectra[ich,ifreq,ih] = numpy.nansum(tmp)#/len(valid[0]) | |
|
842 | #for icr in range(nPairs): | |
|
843 | for icr in range(self.nPairs): | |
|
844 | tmp = numpy.squeeze(cspectra[:,icr,ifreq,ih]) | |
|
845 | valid = (numpy.isfinite(tmp)==True).nonzero() | |
|
846 | if len(valid[0]) > 0: | |
|
847 | out_cspectra[icr,ifreq,ih] = numpy.nansum(tmp)#/len(valid[0]) | |
|
848 | ''' | |
|
849 | # print('##########################################################') | |
|
850 | print("Removing fake coherent echoes (at least 4 points around the point)") | |
|
851 | ||
|
852 | val_spectra = numpy.sum(val_spc,0) | |
|
853 | val_cspectra = numpy.sum(val_cspc,0) | |
|
854 | ||
|
855 | val_spectra = self.REM_ISOLATED_POINTS(val_spectra,4) | |
|
856 | val_cspectra = self.REM_ISOLATED_POINTS(val_cspectra,4) | |
|
857 | ||
|
858 | for i in range(nChan): | |
|
859 | for j in range(nProf): | |
|
860 | for k in range(nHeights): | |
|
861 | if numpy.isfinite(val_spectra[i,j,k]) and val_spectra[i,j,k] < 1 : | |
|
862 | val_spc[:,i,j,k] = 0.0 | |
|
863 | for i in range(nPairs): | |
|
864 | for j in range(nProf): | |
|
865 | for k in range(nHeights): | |
|
866 | if numpy.isfinite(val_cspectra[i,j,k]) and val_cspectra[i,j,k] < 1 : | |
|
867 | val_cspc[:,i,j,k] = 0.0 | |
|
868 | ||
|
869 | # val_spc = numpy.reshape(val_spc, (len(spectra[:,0,0,0]),nProf*nHeights*nChan)) | |
|
870 | # if numpy.isfinite(val_spectra)==str(True): | |
|
871 | # noval = (val_spectra<1).nonzero() | |
|
872 | # if len(noval) > 0: | |
|
873 | # val_spc[:,noval] = 0.0 | |
|
874 | # val_spc = numpy.reshape(val_spc, (149,nChan,nProf,nHeights)) | |
|
875 | ||
|
876 | #val_cspc = numpy.reshape(val_spc, (149,nChan*nHeights*nProf)) | |
|
877 | #if numpy.isfinite(val_cspectra)==str(True): | |
|
878 | # noval = (val_cspectra<1).nonzero() | |
|
879 | # if len(noval) > 0: | |
|
880 | # val_cspc[:,noval] = 0.0 | |
|
881 | # val_cspc = numpy.reshape(val_cspc, (149,nChan,nProf,nHeights)) | |
|
882 | tmp_sat_spectra = spectra.copy() | |
|
883 | tmp_sat_spectra = tmp_sat_spectra*numpy.nan | |
|
884 | tmp_sat_cspectra = cspectra.copy() | |
|
885 | tmp_sat_cspectra = tmp_sat_cspectra*numpy.nan | |
|
886 | ''' | |
|
887 | # fig = plt.figure(figsize=(6,5)) | |
|
888 | # left, bottom, width, height = 0.1, 0.1, 0.8, 0.8 | |
|
889 | # ax = fig.add_axes([left, bottom, width, height]) | |
|
890 | # cp = ax.contour(10*numpy.log10(numpy.absolute(spectra[0,0,:,:]))) | |
|
891 | # ax.clabel(cp, inline=True,fontsize=10) | |
|
892 | # plt.show() | |
|
893 | ''' | |
|
894 | val = (val_spc > 0).nonzero() | |
|
895 | if len(val[0]) > 0: | |
|
896 | tmp_sat_spectra[val] = in_sat_spectra[val] | |
|
897 | val = (val_cspc > 0).nonzero() | |
|
898 | if len(val[0]) > 0: | |
|
899 | tmp_sat_cspectra[val] = in_sat_cspectra[val] | |
|
900 | ||
|
901 | print("Getting average of the spectra and cross-spectra from incoherent echoes 2") | |
|
902 | sat_spectra = numpy.zeros((nChan,nProf,nHeights), dtype=float) | |
|
903 | sat_cspectra = numpy.zeros((nPairs,nProf,nHeights), dtype=complex) | |
|
904 | for ih in range(nHeights): | |
|
905 | for ifreq in range(nProf): | |
|
906 | for ich in range(nChan): | |
|
907 | tmp = numpy.squeeze(tmp_sat_spectra[:,ich,ifreq,ih]) | |
|
908 | valid = (numpy.isfinite(tmp)).nonzero() | |
|
909 | if len(valid[0]) > 0: | |
|
910 | sat_spectra[ich,ifreq,ih] = numpy.nansum(tmp)/len(valid[0]) | |
|
911 | ||
|
912 | for icr in range(nPairs): | |
|
913 | tmp = numpy.squeeze(tmp_sat_cspectra[:,icr,ifreq,ih]) | |
|
914 | valid = (numpy.isfinite(tmp)).nonzero() | |
|
915 | if len(valid[0]) > 0: | |
|
916 | sat_cspectra[icr,ifreq,ih] = numpy.nansum(tmp)/len(valid[0]) | |
|
917 | ''' | |
|
918 | #self.__dataReady= True | |
|
919 | #sat_spectra, sat_cspectra= sat_spectra, sat_cspectra | |
|
920 | #if not self.__dataReady: | |
|
921 | #return None, None | |
|
922 | #return out_spectra, out_cspectra ,sat_spectra,sat_cspectra | |
|
923 | return out_spectra, out_cspectra | |
|
924 | ||
|
925 | def REM_ISOLATED_POINTS(self,array,rth): | |
|
926 | # import matplotlib.pyplot as plt | |
|
927 | if rth == None : | |
|
928 | rth = 4 | |
|
929 | print("REM ISO") | |
|
930 | num_prof = len(array[0,:,0]) | |
|
931 | num_hei = len(array[0,0,:]) | |
|
932 | n2d = len(array[:,0,0]) | |
|
933 | ||
|
934 | for ii in range(n2d) : | |
|
935 | #print ii,n2d | |
|
936 | tmp = array[ii,:,:] | |
|
937 | #print tmp.shape, array[ii,101,:],array[ii,102,:] | |
|
938 | ||
|
939 | # fig = plt.figure(figsize=(6,5)) | |
|
940 | # left, bottom, width, height = 0.1, 0.1, 0.8, 0.8 | |
|
941 | # ax = fig.add_axes([left, bottom, width, height]) | |
|
942 | # x = range(num_prof) | |
|
943 | # y = range(num_hei) | |
|
944 | # cp = ax.contour(y,x,tmp) | |
|
945 | # ax.clabel(cp, inline=True,fontsize=10) | |
|
946 | # plt.show() | |
|
947 | ||
|
948 | #indxs = WHERE(FINITE(tmp) AND tmp GT 0,cindxs) | |
|
949 | tmp = numpy.reshape(tmp,num_prof*num_hei) | |
|
950 | indxs1 = (numpy.isfinite(tmp)==True).nonzero() | |
|
951 | indxs2 = (tmp > 0).nonzero() | |
|
952 | ||
|
953 | indxs1 = (indxs1[0]) | |
|
954 | indxs2 = indxs2[0] | |
|
955 | #indxs1 = numpy.array(indxs1[0]) | |
|
956 | #indxs2 = numpy.array(indxs2[0]) | |
|
957 | indxs = None | |
|
958 | #print indxs1 , indxs2 | |
|
959 | for iv in range(len(indxs2)): | |
|
960 | indv = numpy.array((indxs1 == indxs2[iv]).nonzero()) | |
|
961 | #print len(indxs2), indv | |
|
962 | if len(indv[0]) > 0 : | |
|
963 | indxs = numpy.concatenate((indxs,indxs2[iv]), axis=None) | |
|
964 | # print indxs | |
|
965 | indxs = indxs[1:] | |
|
966 | #print(indxs, len(indxs)) | |
|
967 | if len(indxs) < 4 : | |
|
968 | array[ii,:,:] = 0. | |
|
969 | return | |
|
970 | ||
|
971 | xpos = numpy.mod(indxs ,num_hei) | |
|
972 | ypos = (indxs / num_hei) | |
|
973 | sx = numpy.argsort(xpos) # Ordering respect to "x" (time) | |
|
974 | #print sx | |
|
975 | xpos = xpos[sx] | |
|
976 | ypos = ypos[sx] | |
|
977 | ||
|
978 | # *********************************** Cleaning isolated points ********************************** | |
|
979 | ic = 0 | |
|
980 | while True : | |
|
981 | r = numpy.sqrt(list(numpy.power((xpos[ic]-xpos),2)+ numpy.power((ypos[ic]-ypos),2))) | |
|
982 | #no_coh = WHERE(FINITE(r) AND (r LE rth),cno_coh) | |
|
983 | #plt.plot(r) | |
|
984 | #plt.show() | |
|
985 | no_coh1 = (numpy.isfinite(r)==True).nonzero() | |
|
986 | no_coh2 = (r <= rth).nonzero() | |
|
987 | #print r, no_coh1, no_coh2 | |
|
988 | no_coh1 = numpy.array(no_coh1[0]) | |
|
989 | no_coh2 = numpy.array(no_coh2[0]) | |
|
990 | no_coh = None | |
|
991 | #print valid1 , valid2 | |
|
992 | for iv in range(len(no_coh2)): | |
|
993 | indv = numpy.array((no_coh1 == no_coh2[iv]).nonzero()) | |
|
994 | if len(indv[0]) > 0 : | |
|
995 | no_coh = numpy.concatenate((no_coh,no_coh2[iv]), axis=None) | |
|
996 | no_coh = no_coh[1:] | |
|
997 | #print len(no_coh), no_coh | |
|
998 | if len(no_coh) < 4 : | |
|
999 | #print xpos[ic], ypos[ic], ic | |
|
1000 | # plt.plot(r) | |
|
1001 | # plt.show() | |
|
1002 | xpos[ic] = numpy.nan | |
|
1003 | ypos[ic] = numpy.nan | |
|
1004 | ||
|
1005 | ic = ic + 1 | |
|
1006 | if (ic == len(indxs)) : | |
|
1007 | break | |
|
1008 | #print( xpos, ypos) | |
|
1009 | ||
|
1010 | indxs = (numpy.isfinite(list(xpos))==True).nonzero() | |
|
1011 | #print indxs[0] | |
|
1012 | if len(indxs[0]) < 4 : | |
|
1013 | array[ii,:,:] = 0. | |
|
1014 | return | |
|
1015 | ||
|
1016 | xpos = xpos[indxs[0]] | |
|
1017 | ypos = ypos[indxs[0]] | |
|
1018 | for i in range(0,len(ypos)): | |
|
1019 | ypos[i]=int(ypos[i]) | |
|
1020 | junk = tmp | |
|
1021 | tmp = junk*0.0 | |
|
1022 | ||
|
1023 | tmp[list(xpos + (ypos*num_hei))] = junk[list(xpos + (ypos*num_hei))] | |
|
1024 | array[ii,:,:] = numpy.reshape(tmp,(num_prof,num_hei)) | |
|
1025 | ||
|
1026 | #print array.shape | |
|
1027 | #tmp = numpy.reshape(tmp,(num_prof,num_hei)) | |
|
1028 | #print tmp.shape | |
|
1029 | ||
|
1030 | # fig = plt.figure(figsize=(6,5)) | |
|
1031 | # left, bottom, width, height = 0.1, 0.1, 0.8, 0.8 | |
|
1032 | # ax = fig.add_axes([left, bottom, width, height]) | |
|
1033 | # x = range(num_prof) | |
|
1034 | # y = range(num_hei) | |
|
1035 | # cp = ax.contour(y,x,array[ii,:,:]) | |
|
1036 | # ax.clabel(cp, inline=True,fontsize=10) | |
|
1037 | # plt.show() | |
|
1038 | return array | |
|
1039 | ||
|
481 | 1040 | class removeInterference(Operation): |
|
482 | 1041 | |
|
483 | 1042 | def removeInterference2(self): |
General Comments 0
You need to be logged in to leave comments.
Login now