schain Commit - r1401:c018ae3094fa · Jicamarca Repository

Límite de bloques, hora, zona horaria

joabAM -

r1401:c018ae3094fa

parent child

Context file:

r1401:c018ae3094fa

Collapse all files

schainpy/model/io/jroIO_param.py +30 -16

             import os
             import time
             import datetime
             import numpy
             import h5py
             import schainpy.admin
             from schainpy.model.data.jrodata import *
             from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
             from schainpy.model.io.jroIO_base import *
             from schainpy.utils import log
             class HDFReader(Reader, ProcessingUnit):
                 """Processing unit to read HDF5 format files
                 This unit reads HDF5 files created with `HDFWriter` operation contains
                 by default two groups Data and Metadata all variables would be saved as `dataOut`
                 attributes.
                 It is possible to read any HDF5 file by given the structure in the `description`
                 parameter, also you can add extra values to metadata with the parameter `extras`.
                 Parameters:
                 -----------
                 path : str
                     Path where files are located.
                 startDate : date
                     Start date of the files
                 endDate : list
                     End date of the files
                 startTime : time
                     Start time of the files
                 endTime : time
                     End time of the files
                 description : dict, optional
                     Dictionary with the description of the HDF5 file
                 extras : dict, optional
                     Dictionary with extra metadata to be be added to `dataOut`
                 Examples
                 --------
                 desc = {
                     'Data': {
                         'data_output': ['u', 'v', 'w'],
                         'utctime': 'timestamps',
                     }  ,
                     'Metadata': {
                         'heightList': 'heights'
                     }
                 }
                 desc = {
                     'Data': {
                         'data_output': 'winds',
                         'utctime': 'timestamps'
                     },
                     'Metadata': {
                         'heightList': 'heights'
                     }
                 }
                 extras = {
                     'timeZone': 300
                 }
                 reader = project.addReadUnit(
                     name='HDFReader',
                     path='/path/to/files',
                     startDate='2019/01/01',
                     endDate='2019/01/31',
                     startTime='00:00:00',
                     endTime='23:59:59',
                     # description=json.dumps(desc),
                     # extras=json.dumps(extras),
                     )
                 """
                 __attrs__ = ['path', 'startDate', 'endDate', 'startTime', 'endTime', 'description', 'extras']
                 def __init__(self):
                     ProcessingUnit.__init__(self)
                     self.dataOut = Parameters()
                     self.ext = ".hdf5"
                     self.optchar = "D"
                     self.meta = {}
                     self.data = {}
                     self.open_file = h5py.File
                     self.open_mode = 'r'
                     self.description = {}
                     self.extras = {}
                     self.filefmt = "*%Y%j***"
                     self.folderfmt = "*%Y%j"
                     self.utcoffset = 0
                 def setup(self, **kwargs):
                     self.set_kwargs(**kwargs)
                     if not self.ext.startswith('.'):
                         self.ext = '.{}'.format(self.ext)
                     if self.online:
                         log.log("Searching files in online mode...", self.name)
                         for nTries in range(self.nTries):
                             fullpath = self.searchFilesOnLine(self.path, self.startDate,
                                 self.endDate, self.expLabel, self.ext, self.walk,
                                 self.filefmt, self.folderfmt)
                             pathname, filename = os.path.split(fullpath)
                             #print(pathname,filename)
                             try:
                                 fullpath = next(fullpath)
                             except:
                                 fullpath = None
                             if fullpath:
                                 break
                             log.warning(
                                 'Waiting {} sec for a valid file in {}: try {} ...'.format(
                                     self.delay, self.path, nTries + 1),
                                 self.name)
                             time.sleep(self.delay)
                         if not(fullpath):
                             raise schainpy.admin.SchainError(
                                 'There isn\'t any valid file in {}'.format(self.path))
                         pathname, filename = os.path.split(fullpath)
                         self.year = int(filename[1:5])
                         self.doy = int(filename[5:8])
                         self.set = int(filename[8:11]) - 1
                     else:
                         log.log("Searching files in {}".format(self.path), self.name)
                         self.filenameList = self.searchFilesOffLine(self.path, self.startDate,
                             self.endDate, self.expLabel, self.ext, self.walk, self.filefmt, self.folderfmt)
                     self.setNextFile()
                     return
                 def readFirstHeader(self):
                     '''Read metadata and data'''
                     self.__readMetadata()
                     self.__readData()
                     self.__setBlockList()
                     if 'type' in self.meta:
+                        ##print("Creting dataOut...")
                         self.dataOut = eval(self.meta['type'])()
+                        ##print(vars(self.dataOut))
                     for attr in self.meta:
-                        #print("attr: ", attr)
+                        ##print("attr: ", attr)
+                        ##print(type(self.dataOut).__name__)
                         setattr(self.dataOut, attr, self.meta[attr])
                     self.blockIndex = 0
                     return
                 def __setBlockList(self):
                     '''
                     Selects the data within the times defined
                     self.fp
                     self.startTime
                     self.endTime
                     self.blockList
                     self.blocksPerFile
                     '''
                     startTime = self.startTime
                     endTime = self.endTime
                     thisUtcTime = self.data['utctime'] + self.utcoffset
                     self.interval = numpy.min(thisUtcTime[1:] - thisUtcTime[:-1])
                     thisDatetime = datetime.datetime.utcfromtimestamp(thisUtcTime[0])
                     self.startFileDatetime = thisDatetime
                     thisDate = thisDatetime.date()
                     thisTime = thisDatetime.time()
                     startUtcTime = (datetime.datetime.combine(thisDate, startTime) - datetime.datetime(1970, 1, 1)).total_seconds()
                     endUtcTime = (datetime.datetime.combine(thisDate, endTime) - datetime.datetime(1970, 1, 1)).total_seconds()
                     ind = numpy.where(numpy.logical_and(thisUtcTime >= startUtcTime, thisUtcTime < endUtcTime))[0]
                     self.blockList = ind
                     self.blocksPerFile = len(ind)
                     self.blocksPerFile = len(thisUtcTime)
                     return
                 def __readMetadata(self):
                     '''
                     Reads Metadata
                     '''
                     meta = {}
                     if self.description:
                         for key, value in self.description['Metadata'].items():
                             meta[key] = self.fp[value][()]
                     else:
                         grp = self.fp['Metadata']
                         for name in grp:
                             meta[name] = grp[name][()]
                     if self.extras:
                         for key, value in self.extras.items():
                             meta[key] = value
                     self.meta = meta
                     return
                 def checkForRealPath(self, nextFile, nextDay):
                     # print("check FRP")
                     # dt = self.startFileDatetime + datetime.timedelta(1)
                     # filename = '{}.{}{}'.format(self.path, dt.strftime('%Y%m%d'), self.ext)
                     # fullfilename = os.path.join(self.path, filename)
                     # print("check Path ",fullfilename,filename)
                     # if os.path.exists(fullfilename):
                     #     return fullfilename, filename
                     # return None, filename
                     return None,None
                 def __readData(self):
                     data = {}
                     if self.description:
                         for key, value in self.description['Data'].items():
                             if isinstance(value, str):
                                 if isinstance(self.fp[value], h5py.Dataset):
                                     data[key] = self.fp[value][()]
                                 elif isinstance(self.fp[value], h5py.Group):
                                     array = []
                                     for ch in self.fp[value]:
                                         array.append(self.fp[value][ch][()])
                                     data[key] = numpy.array(array)
                             elif isinstance(value, list):
                                 array = []
                                 for ch in value:
                                     array.append(self.fp[ch][()])
                                 data[key] = numpy.array(array)
                     else:
                         grp = self.fp['Data']
                         for name in grp:
                             if isinstance(grp[name], h5py.Dataset):
                                 array = grp[name][()]
                             elif isinstance(grp[name], h5py.Group):
                                 array = []
                                 for ch in grp[name]:
                                     array.append(grp[name][ch][()])
                                 array = numpy.array(array)
                             else:
                                 log.warning('Unknown type: {}'.format(name))
                             if name in self.description:
                                 key = self.description[name]
                             else:
                                 key = name
                             data[key] = array
                     self.data = data
                     return
                 def getData(self):
                     if not self.isDateTimeInRange(self.startFileDatetime, self.startDate, self.endDate, self.startTime, self.endTime):
                         self.dataOut.flagNoData = True
                         self.blockIndex = self.blocksPerFile
                         #self.dataOut.error = True      TERMINA EL PROGRAMA, removido
                         return
                     for attr in self.data:
+                        #print("attr ",attr)
                         if self.data[attr].ndim == 1:
                             setattr(self.dataOut, attr, self.data[attr][self.blockIndex])
                         else:
                             setattr(self.dataOut, attr, self.data[attr][:, self.blockIndex])
-                    self.dataOut.flagNoData = False
                     self.blockIndex += 1
                     if self.blockIndex == 1:
                         log.log("Block No. {}/{} -> {}".format(
                             self.blockIndex,
                             self.blocksPerFile,
                             self.dataOut.datatime.ctime()), self.name)
                     else:
                         log.log("Block No. {}/{} ".format(
                             self.blockIndex,
                             self.blocksPerFile),self.name)
+                    self.dataOut.flagNoData = False
+                    self.dataOut.error = False
                     return
                 def run(self, **kwargs):
                     if not(self.isConfig):
                         self.setup(**kwargs)
                         self.isConfig = True
                     if self.blockIndex == self.blocksPerFile:
                         self.setNextFile()
                     self.getData()
                     return
             @MPDecorator
             class HDFWriter(Operation):
                 """Operation to write HDF5 files.
                 The HDF5 file contains by default two groups Data and Metadata where
                 you can save any `dataOut` attribute specified by `dataList` and `metadataList`
                 parameters, data attributes are normaly time dependent where the metadata
                 are not.
                 It is possible to customize the structure of the HDF5 file with the
                 optional description parameter see the examples.
                 Parameters:
                 -----------
                 path : str
                     Path where files will be saved.
                 blocksPerFile : int
                     Number of blocks per file
                 metadataList : list
                     List of the dataOut attributes that will be saved as metadata
                 dataList : int
                     List of the dataOut attributes that will be saved as data
                 setType : bool
                     If True the name of the files corresponds to the timestamp of the data
                 description : dict, optional
                     Dictionary with the desired description of the HDF5 file
                 Examples
                 --------
                 desc = {
                     'data_output': {'winds': ['z', 'w', 'v']},
                     'utctime': 'timestamps',
                     'heightList': 'heights'
                 }
                 desc = {
                     'data_output': ['z', 'w', 'v'],
                     'utctime': 'timestamps',
                     'heightList': 'heights'
                 }
                 desc = {
                     'Data': {
                         'data_output': 'winds',
                         'utctime': 'timestamps'
                     },
                     'Metadata': {
                         'heightList': 'heights'
                     }
                 }
                 writer = proc_unit.addOperation(name='HDFWriter')
                 writer.addParameter(name='path', value='/path/to/file')
                 writer.addParameter(name='blocksPerFile', value='32')
                 writer.addParameter(name='metadataList', value='heightList,timeZone')
                 writer.addParameter(name='dataList',value='data_output,utctime')
                 # writer.addParameter(name='description',value=json.dumps(desc))
                 """
                 ext = ".hdf5"
                 optchar = "D"
                 filename = None
                 path = None
                 setFile = None
                 fp = None
                 firsttime = True
                 #Configurations
                 blocksPerFile = None
                 blockIndex = None
                 dataOut = None
                 #Data Arrays
                 dataList = None
                 metadataList = None
                 currentDay = None
                 lastTime = None
+                typeTime = "ut"
+                hourLimit = 3
+                breakDays = True
                 def __init__(self):
                     Operation.__init__(self)
                     return
-                def setup(self, path=None, blocksPerFile=10, metadataList=None, dataList=None, setType=None, description=None):
+                def setup(self, path=None, blocksPerFile=10, metadataList=None, dataList=None, setType=None,
+                            description=None,typeTime = "ut",hourLimit = 3, breakDays=True):
                     self.path = path
                     self.blocksPerFile = blocksPerFile
                     self.metadataList = metadataList
                     self.dataList = [s.strip() for s in dataList]
                     self.setType = setType
                     self.description = description
+                    self.timeZone = typeTime
+                    self.hourLimit = hourLimit
+                    self.breakDays = breakDays
                     if self.metadataList is None:
                         self.metadataList = self.dataOut.metadata_list
                     tableList = []
                     dsList = []
                     for i in range(len(self.dataList)):
                         dsDict = {}
                         if hasattr(self.dataOut, self.dataList[i]):
                             dataAux = getattr(self.dataOut, self.dataList[i])
                             dsDict['variable'] = self.dataList[i]
                         else:
                             log.warning('Attribute {} not found in dataOut', self.name)
                             continue
                         if dataAux is None:
                             continue
                         elif isinstance(dataAux, (int, float, numpy.integer, numpy.float)):
                             dsDict['nDim'] = 0
                         else:
                             dsDict['nDim'] = len(dataAux.shape)
                             dsDict['shape'] = dataAux.shape
                             dsDict['dsNumber'] = dataAux.shape[0]
                             dsDict['dtype'] = dataAux.dtype
                         dsList.append(dsDict)
                     self.dsList = dsList
                     self.currentDay = self.dataOut.datatime.date()
                 def timeFlag(self):
                     currentTime = self.dataOut.utctime
-                    timeTuple = time.localtime(currentTime)
+                    if self.timeZone == "lt":
+                        timeTuple = time.localtime(currentTime)
+                    elif self.timeZone == "ut":
+                        timeTuple = time.gmtime(currentTime)
                     dataDay = timeTuple.tm_yday
                     #print("time UTC: ",currentTime, self.dataOut.datatime)
                     if self.lastTime is None:
                         self.lastTime = currentTime
                         self.currentDay = dataDay
                         return False
                     timeDiff = currentTime - self.lastTime
-                    #Si el dia es diferente o si la diferencia entre un dato y otro supera la hora
+                    #Si el dia es diferente o si la diferencia entre un dato y otro supera self.hourLimit
-                    if dataDay != self.currentDay:
+                    if (dataDay != self.currentDay) and self.breakDays:
                         self.currentDay = dataDay
                         return True
-                    elif timeDiff > 3*60*60:
+                    elif timeDiff > self.hourLimit*60*60:
                         self.lastTime = currentTime
                         return True
                     else:
                         self.lastTime = currentTime
                         return False
-                def run(self, dataOut, path, blocksPerFile=10, metadataList=None,
+                def run(self, dataOut,**kwargs):
-                        dataList=[], setType=None, description={}):
                     self.dataOut = dataOut
                     if not(self.isConfig):
-                        self.setup(path=path, blocksPerFile=blocksPerFile,
+                        self.setup(**kwargs)
-                                   metadataList=metadataList, dataList=dataList,
-                                   setType=setType, description=description)
                         self.isConfig = True
                         self.setNextFile()
                     self.putData()
                     return
                 def setNextFile(self):
                     ext = self.ext
                     path = self.path
                     setFile = self.setFile
+                    if self.timeZone == "lt":
-                    timeTuple = time.gmtime(self.dataOut.utctime)
+                        timeTuple = time.localtime(self.dataOut.utctime)
+                    elif self.timeZone == "ut":
+                        timeTuple = time.gmtime(self.dataOut.utctime)
                     #print("path: ",timeTuple)
                     subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
                     fullpath = os.path.join(path, subfolder)
                     if os.path.exists(fullpath):
                         filesList = os.listdir(fullpath)
                         filesList = [k for k in filesList if k.startswith(self.optchar)]
                         if len( filesList ) > 0:
                             filesList = sorted(filesList, key=str.lower)
                             filen = filesList[-1]
                             # el filename debera tener el siguiente formato
                             # 0 1234 567 89A BCDE (hex)
                             # x YYYY DDD SSS .ext
                             if isNumber(filen[8:11]):
                                 setFile = int(filen[8:11]) #inicializo mi contador de seteo al seteo del ultimo file
                             else:
                                 setFile = -1
                         else:
                             setFile = -1 #inicializo mi contador de seteo
                     else:
                         os.makedirs(fullpath)
                         setFile = -1 #inicializo mi contador de seteo
                     if self.setType is None:
                         setFile += 1
                         file = '%s%4.4d%3.3d%03d%s' % (self.optchar,
                                                        timeTuple.tm_year,
                                                        timeTuple.tm_yday,
                                                        setFile,
                                                        ext )
                     else:
                         setFile = timeTuple.tm_hour*60+timeTuple.tm_min
                         file = '%s%4.4d%3.3d%04d%s' % (self.optchar,
                                                        timeTuple.tm_year,
                                                        timeTuple.tm_yday,
                                                        setFile,
                                                        ext )
                     self.filename = os.path.join( path, subfolder, file )
                     #Setting HDF5 File
                     self.fp = h5py.File(self.filename, 'w')
                     #write metadata
                     self.writeMetadata(self.fp)
                     #Write data
                     self.writeData(self.fp)
                 def getLabel(self, name, x=None):
                     if x is None:
                         if 'Data' in self.description:
                             data = self.description['Data']
                             if 'Metadata' in self.description:
                                 data.update(self.description['Metadata'])
                         else:
                             data = self.description
                         if name in data:
                             if isinstance(data[name], str):
                                 return data[name]
                             elif isinstance(data[name], list):
                                 return None
                             elif isinstance(data[name], dict):
                                 for key, value in data[name].items():
                                     return key
                         return name
                     else:
                         if 'Metadata' in self.description:
                             meta = self.description['Metadata']
                         else:
                             meta = self.description
                         if name in meta:
                             if isinstance(meta[name], list):
                                 return meta[name][x]
                             elif isinstance(meta[name], dict):
                                 for key, value in meta[name].items():
                                     return value[x]
                         if 'cspc' in name:
                             return 'pair{:02d}'.format(x)
                         else:
                             return 'channel{:02d}'.format(x)
                 def writeMetadata(self, fp):
                     if self.description:
                         if 'Metadata' in self.description:
                             grp = fp.create_group('Metadata')
                         else:
                             grp = fp
                     else:
                         grp = fp.create_group('Metadata')
                     for i in range(len(self.metadataList)):
                         if not hasattr(self.dataOut, self.metadataList[i]):
                             log.warning('Metadata: `{}` not found'.format(self.metadataList[i]), self.name)
                             continue
                         value = getattr(self.dataOut, self.metadataList[i])
                         if isinstance(value, bool):
                             if value is True:
                                 value = 1
                             else:
                                 value = 0
                         grp.create_dataset(self.getLabel(self.metadataList[i]), data=value)
                     return
                 def writeData(self, fp):
                     if self.description:
                         if 'Data' in self.description:
                             grp = fp.create_group('Data')
                         else:
                             grp = fp
                     else:
                         grp = fp.create_group('Data')
                     dtsets = []
                     data = []
                     for dsInfo in self.dsList:
                         if dsInfo['nDim'] == 0:
                             ds = grp.create_dataset(
                                 self.getLabel(dsInfo['variable']),
                                 (self.blocksPerFile, ),
                                 chunks=True,
                                 dtype=numpy.float64)
                             dtsets.append(ds)
                             data.append((dsInfo['variable'], -1))
                         else:
                             label = self.getLabel(dsInfo['variable'])
                             if label is not None:
                                 sgrp = grp.create_group(label)
                             else:
                                 sgrp = grp
                             for i in range(dsInfo['dsNumber']):
                                 ds = sgrp.create_dataset(
                                     self.getLabel(dsInfo['variable'], i),
                                     (self.blocksPerFile, ) + dsInfo['shape'][1:],
                                     chunks=True,
                                     dtype=dsInfo['dtype'])
                                 dtsets.append(ds)
                                 data.append((dsInfo['variable'], i))
                     fp.flush()
                     log.log('Creating file: {}'.format(fp.filename), self.name)
                     self.ds = dtsets
                     self.data = data
                     self.firsttime = True
                     self.blockIndex = 0
                     return
                 def putData(self):
                     if (self.blockIndex == self.blocksPerFile) or self.timeFlag():
                         self.closeFile()
                         self.setNextFile()
                     for i, ds in enumerate(self.ds):
                         attr, ch = self.data[i]
                         if ch == -1:
                             ds[self.blockIndex] = getattr(self.dataOut, attr)
                         else:
                             ds[self.blockIndex] = getattr(self.dataOut, attr)[ch]
                     self.fp.flush()
                     self.blockIndex += 1
                     if self.blockIndex == 1:
                         log.log('Block No. {}/{} --> {}'.format(self.blockIndex, self.blocksPerFile,self.dataOut.datatime.ctime()), self.name)
                     else:
                         log.log('Block No. {}/{}'.format(self.blockIndex, self.blocksPerFile), self.name)
                     return
                 def closeFile(self):
                     if self.blockIndex != self.blocksPerFile:
                         for ds in self.ds:
                             ds.resize(self.blockIndex, axis=0)
                     if self.fp:
                         self.fp.flush()
                         self.fp.close()
                 def close(self):
                     self.closeFile()

schainpy/model/proc/jroproc_spectra.py +1 -1

             # Copyright (c) 2012-2020 Jicamarca Radio Observatory
             # All rights reserved.
             #
             # Distributed under the terms of the BSD 3-clause license.
             """Spectra processing Unit and operations
             Here you will find the processing unit `SpectraProc` and several operations
             to work with Spectra data type
             """
             import time
             import itertools
             import numpy
             import math
             from schainpy.model.proc.jroproc_base import ProcessingUnit, MPDecorator, Operation
             from schainpy.model.data.jrodata import Spectra
             from schainpy.model.data.jrodata import hildebrand_sekhon
             from schainpy.utils import log
             from scipy.optimize import curve_fit
             class SpectraProc(ProcessingUnit):
                 def __init__(self):
                     ProcessingUnit.__init__(self)
                     self.buffer = None
                     self.firstdatatime = None
                     self.profIndex = 0
                     self.dataOut = Spectra()
                     self.id_min = None
                     self.id_max = None
                     self.setupReq = False #Agregar a todas las unidades de proc
                 def __updateSpecFromVoltage(self):
                     self.dataOut.timeZone = self.dataIn.timeZone
                     self.dataOut.dstFlag = self.dataIn.dstFlag
                     self.dataOut.errorCount = self.dataIn.errorCount
                     self.dataOut.useLocalTime = self.dataIn.useLocalTime
                     try:
                         self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
                     except:
                         pass
                     self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                     self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
                     self.dataOut.channelList = self.dataIn.channelList
                     self.dataOut.heightList = self.dataIn.heightList
                     self.dataOut.dtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
                     self.dataOut.nProfiles = self.dataOut.nFFTPoints
                     self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
                     self.dataOut.utctime = self.firstdatatime
                     self.dataOut.flagDecodeData = self.dataIn.flagDecodeData
                     self.dataOut.flagDeflipData = self.dataIn.flagDeflipData
                     self.dataOut.flagShiftFFT = False
                     self.dataOut.nCohInt = self.dataIn.nCohInt
                     self.dataOut.nIncohInt = 1
                     self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
                     self.dataOut.frequency = self.dataIn.frequency
                     self.dataOut.realtime = self.dataIn.realtime
                     self.dataOut.azimuth = self.dataIn.azimuth
                     self.dataOut.zenith = self.dataIn.zenith
                     self.dataOut.codeList = self.dataIn.codeList
                     self.dataOut.azimuthList = self.dataIn.azimuthList
                     self.dataOut.elevationList = self.dataIn.elevationList
                 def __getFft(self):
                     """
                     Convierte valores de Voltaje a Spectra
                     Affected:
                         self.dataOut.data_spc
                         self.dataOut.data_cspc
                         self.dataOut.data_dc
                         self.dataOut.heightList
                         self.profIndex
                         self.buffer
                         self.dataOut.flagNoData
                     """
                     fft_volt = numpy.fft.fft(
                         self.buffer, n=self.dataOut.nFFTPoints, axis=1)
                     fft_volt = fft_volt.astype(numpy.dtype('complex'))
                     dc = fft_volt[:, 0, :]
                     # calculo de self-spectra
                     fft_volt = numpy.fft.fftshift(fft_volt, axes=(1,))
                     spc = fft_volt * numpy.conjugate(fft_volt)
                     spc = spc.real
                     blocksize = 0
                     blocksize += dc.size
                     blocksize += spc.size
                     cspc = None
                     pairIndex = 0
                     if self.dataOut.pairsList != None:
                         # calculo de cross-spectra
                         cspc = numpy.zeros(
                             (self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
                         for pair in self.dataOut.pairsList:
                             if pair[0] not in self.dataOut.channelList:
                                 raise ValueError("Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" % (
                                     str(pair), str(self.dataOut.channelList)))
                             if pair[1] not in self.dataOut.channelList:
                                 raise ValueError("Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" % (
                                     str(pair), str(self.dataOut.channelList)))
                             cspc[pairIndex, :, :] = fft_volt[pair[0], :, :] * \
                                 numpy.conjugate(fft_volt[pair[1], :, :])
                             pairIndex += 1
                         blocksize += cspc.size
                     self.dataOut.data_spc = spc
                     self.dataOut.data_cspc = cspc
                     self.dataOut.data_dc = dc
                     self.dataOut.blockSize = blocksize
                     self.dataOut.flagShiftFFT = False
                 def run(self, nProfiles=None, nFFTPoints=None, pairsList=None, ippFactor=None, shift_fft=False):
                     if self.dataIn.type == "Spectra":
                         self.dataOut.copy(self.dataIn)
                         if shift_fft:
                             #desplaza a la derecha en el eje 2 determinadas posiciones
                             shift = int(self.dataOut.nFFTPoints/2)
                             self.dataOut.data_spc = numpy.roll(self.dataOut.data_spc, shift , axis=1)
                             if self.dataOut.data_cspc is not None:
                                 #desplaza a la derecha en el eje 2 determinadas posiciones
                                 self.dataOut.data_cspc = numpy.roll(self.dataOut.data_cspc, shift, axis=1)
                         if pairsList:
                             self.__selectPairs(pairsList)
                     elif self.dataIn.type == "Voltage":
                         self.dataOut.flagNoData = True
                         if nFFTPoints == None:
                             raise ValueError("This SpectraProc.run() need nFFTPoints input variable")
                         if nProfiles == None:
                             nProfiles = nFFTPoints
                         if ippFactor == None:
                             self.dataOut.ippFactor = 1
                         self.dataOut.nFFTPoints = nFFTPoints
                         if self.buffer is None:
                             self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nProfiles,
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                         if self.dataIn.flagDataAsBlock:
                             nVoltProfiles = self.dataIn.data.shape[1]
                             if nVoltProfiles == nProfiles:
                                 self.buffer = self.dataIn.data.copy()
                                 self.profIndex = nVoltProfiles
                             elif nVoltProfiles < nProfiles:
                                 if self.profIndex == 0:
                                     self.id_min = 0
                                     self.id_max = nVoltProfiles
                                 self.buffer[:, self.id_min:self.id_max,
                                             :] = self.dataIn.data
                                 self.profIndex += nVoltProfiles
                                 self.id_min += nVoltProfiles
                                 self.id_max += nVoltProfiles
                             else:
                                 raise ValueError("The type object %s has %d profiles, it should just has %d profiles" % (
                                     self.dataIn.type, self.dataIn.data.shape[1], nProfiles))
                                 self.dataOut.flagNoData = True
                         else:
                             self.buffer[:, self.profIndex, :] = self.dataIn.data.copy()
                             self.profIndex += 1
                         if self.firstdatatime == None:
                             self.firstdatatime = self.dataIn.utctime
                         if self.profIndex == nProfiles:
                             self.__updateSpecFromVoltage()
                             if pairsList == None:
                                 self.dataOut.pairsList = [pair for pair in itertools.combinations(self.dataOut.channelList, 2)]
                             else:
                                 self.dataOut.pairsList = pairsList
                             self.__getFft()
                             self.dataOut.flagNoData = False
                             self.firstdatatime = None
                             self.profIndex = 0
                     else:
                         raise ValueError("The type of input object '%s' is not valid".format(
                             self.dataIn.type))
                 def __selectPairs(self, pairsList):
                     if not pairsList:
                         return
                     pairs = []
                     pairsIndex = []
                     for pair in pairsList:
                         if pair[0] not in self.dataOut.channelList or pair[1] not in self.dataOut.channelList:
                             continue
                         pairs.append(pair)
                         pairsIndex.append(pairs.index(pair))
                     self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndex]
                     self.dataOut.pairsList = pairs
                     return
                 def selectFFTs(self, minFFT, maxFFT ):
                     """
                     Selecciona un bloque de datos en base a un grupo de valores de puntos FFTs segun el rango
                     minFFT<= FFT <= maxFFT
                     """
                     if (minFFT > maxFFT):
                         raise ValueError("Error selecting heights: Height range (%d,%d) is not valid" % (minFFT, maxFFT))
                     if (minFFT < self.dataOut.getFreqRange()[0]):
                         minFFT = self.dataOut.getFreqRange()[0]
                     if (maxFFT > self.dataOut.getFreqRange()[-1]):
                         maxFFT = self.dataOut.getFreqRange()[-1]
                     minIndex = 0
                     maxIndex = 0
                     FFTs = self.dataOut.getFreqRange()
                     inda = numpy.where(FFTs >= minFFT)
                     indb = numpy.where(FFTs <= maxFFT)
                     try:
                         minIndex = inda[0][0]
                     except:
                         minIndex = 0
                     try:
                         maxIndex = indb[0][-1]
                     except:
                         maxIndex = len(FFTs)
                     self.selectFFTsByIndex(minIndex, maxIndex)
                     return 1
                 def getBeaconSignal(self, tauindex=0, channelindex=0, hei_ref=None):
                     newheis = numpy.where(
                         self.dataOut.heightList > self.dataOut.radarControllerHeaderObj.Taus[tauindex])
                     if hei_ref != None:
                         newheis = numpy.where(self.dataOut.heightList > hei_ref)
                     minIndex = min(newheis[0])
                     maxIndex = max(newheis[0])
                     data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
                     heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
                     # determina indices
                     nheis = int(self.dataOut.radarControllerHeaderObj.txB /
                                 (self.dataOut.heightList[1] - self.dataOut.heightList[0]))
                     avg_dB = 10 * \
                         numpy.log10(numpy.sum(data_spc[channelindex, :, :], axis=0))
                     beacon_dB = numpy.sort(avg_dB)[-nheis:]
                     beacon_heiIndexList = []
                     for val in avg_dB.tolist():
                         if val >= beacon_dB[0]:
                             beacon_heiIndexList.append(avg_dB.tolist().index(val))
                     #data_spc = data_spc[:,:,beacon_heiIndexList]
                     data_cspc = None
                     if self.dataOut.data_cspc is not None:
                         data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
                         #data_cspc = data_cspc[:,:,beacon_heiIndexList]
                     data_dc = None
                     if self.dataOut.data_dc is not None:
                         data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
                         #data_dc = data_dc[:,beacon_heiIndexList]
                     self.dataOut.data_spc = data_spc
                     self.dataOut.data_cspc = data_cspc
                     self.dataOut.data_dc = data_dc
                     self.dataOut.heightList = heightList
                     self.dataOut.beacon_heiIndexList = beacon_heiIndexList
                     return 1
                 def selectFFTsByIndex(self, minIndex, maxIndex):
                     """
                     """
                     if (minIndex < 0) or (minIndex > maxIndex):
                         raise ValueError("Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex))
                     if (maxIndex >= self.dataOut.nProfiles):
                         maxIndex = self.dataOut.nProfiles-1
                     #Spectra
                     data_spc = self.dataOut.data_spc[:,minIndex:maxIndex+1,:]
                     data_cspc = None
                     if self.dataOut.data_cspc is not None:
                         data_cspc = self.dataOut.data_cspc[:,minIndex:maxIndex+1,:]
                     data_dc = None
                     if self.dataOut.data_dc is not None:
                         data_dc = self.dataOut.data_dc[minIndex:maxIndex+1,:]
                     self.dataOut.data_spc = data_spc
                     self.dataOut.data_cspc = data_cspc
                     self.dataOut.data_dc = data_dc
                     self.dataOut.ippSeconds = self.dataOut.ippSeconds*(self.dataOut.nFFTPoints / numpy.shape(data_cspc)[1])
                     self.dataOut.nFFTPoints = numpy.shape(data_cspc)[1]
                     self.dataOut.profilesPerBlock = numpy.shape(data_cspc)[1]
                     return 1
                 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
                     # validacion de rango
                     if minHei == None:
                         minHei = self.dataOut.heightList[0]
                     if maxHei == None:
                         maxHei = self.dataOut.heightList[-1]
                     if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
                         print('minHei: %.2f is out of the heights range' % (minHei))
                         print('minHei is setting to %.2f' % (self.dataOut.heightList[0]))
                         minHei = self.dataOut.heightList[0]
                     if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
                         print('maxHei: %.2f is out of the heights range' % (maxHei))
                         print('maxHei is setting to %.2f' % (self.dataOut.heightList[-1]))
                         maxHei = self.dataOut.heightList[-1]
                     # validacion de velocidades
                     velrange = self.dataOut.getVelRange(1)
                     if minVel == None:
                         minVel = velrange[0]
                     if maxVel == None:
                         maxVel = velrange[-1]
                     if (minVel < velrange[0]) or (minVel > maxVel):
                         print('minVel: %.2f is out of the velocity range' % (minVel))
                         print('minVel is setting to %.2f' % (velrange[0]))
                         minVel = velrange[0]
                     if (maxVel > velrange[-1]) or (maxVel < minVel):
                         print('maxVel: %.2f is out of the velocity range' % (maxVel))
                         print('maxVel is setting to %.2f' % (velrange[-1]))
                         maxVel = velrange[-1]
                     # seleccion de indices para rango
                     minIndex = 0
                     maxIndex = 0
                     heights = self.dataOut.heightList
                     inda = numpy.where(heights >= minHei)
                     indb = numpy.where(heights <= maxHei)
                     try:
                         minIndex = inda[0][0]
                     except:
                         minIndex = 0
                     try:
                         maxIndex = indb[0][-1]
                     except:
                         maxIndex = len(heights)
                     if (minIndex < 0) or (minIndex > maxIndex):
                         raise ValueError("some value in (%d,%d) is not valid" % (
                             minIndex, maxIndex))
                     if (maxIndex >= self.dataOut.nHeights):
                         maxIndex = self.dataOut.nHeights - 1
                     # seleccion de indices para velocidades
                     indminvel = numpy.where(velrange >= minVel)
                     indmaxvel = numpy.where(velrange <= maxVel)
                     try:
                         minIndexVel = indminvel[0][0]
                     except:
                         minIndexVel = 0
                     try:
                         maxIndexVel = indmaxvel[0][-1]
                     except:
                         maxIndexVel = len(velrange)
                     # seleccion del espectro
                     data_spc = self.dataOut.data_spc[:,
                                                      minIndexVel:maxIndexVel + 1, minIndex:maxIndex + 1]
                     # estimacion de ruido
                     noise = numpy.zeros(self.dataOut.nChannels)
                     for channel in range(self.dataOut.nChannels):
                         daux = data_spc[channel, :, :]
                         sortdata = numpy.sort(daux, axis=None)
                         noise[channel] = hildebrand_sekhon(sortdata, self.dataOut.nIncohInt)
                     self.dataOut.noise_estimation = noise.copy()
                     return 1
             class removeDC(Operation):
                 def run(self, dataOut, mode=2):
                     self.dataOut = dataOut
                     jspectra = self.dataOut.data_spc
                     jcspectra = self.dataOut.data_cspc
                     num_chan = jspectra.shape[0]
                     num_hei = jspectra.shape[2]
                     if jcspectra is not None:
                         jcspectraExist = True
                         num_pairs = jcspectra.shape[0]
                     else:
                         jcspectraExist = False
                     freq_dc = int(jspectra.shape[1] / 2)
                     ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
                     ind_vel = ind_vel.astype(int)
                     if ind_vel[0] < 0:
                         ind_vel[list(range(0, 1))] = ind_vel[list(range(0, 1))] + self.num_prof
                     if mode == 1:
                         jspectra[:, freq_dc, :] = (
                             jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2  # CORRECCION
                         if jcspectraExist:
                             jcspectra[:, freq_dc, :] = (
                                 jcspectra[:, ind_vel[1], :] + jcspectra[:, ind_vel[2], :]) / 2
                     if mode == 2:
                         vel = numpy.array([-2, -1, 1, 2])
                         xx = numpy.zeros([4, 4])
                         for fil in range(4):
                             xx[fil, :] = vel[fil]**numpy.asarray(list(range(4)))
                         xx_inv = numpy.linalg.inv(xx)
                         xx_aux = xx_inv[0, :]
                         for ich in range(num_chan):
                             yy = jspectra[ich, ind_vel, :]
                             jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
                             junkid = jspectra[ich, freq_dc, :] <= 0
                             cjunkid = sum(junkid)
                             if cjunkid.any():
                                 jspectra[ich, freq_dc, junkid.nonzero()] = (
                                     jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
                         if jcspectraExist:
                             for ip in range(num_pairs):
                                 yy = jcspectra[ip, ind_vel, :]
                                 jcspectra[ip, freq_dc, :] = numpy.dot(xx_aux, yy)
                     self.dataOut.data_spc = jspectra
                     self.dataOut.data_cspc = jcspectra
                     return self.dataOut
             # import matplotlib.pyplot as plt
             def fit_func( x, a0, a1, a2): #, a3, a4, a5):
                 z = (x - a1) / a2
                 y = a0 * numpy.exp(-z**2 / a2)  #+ a3 + a4 * x + a5 * x**2
                 return y
             class CleanRayleigh(Operation):
                 def __init__(self):
                     Operation.__init__(self)
                     self.i=0
                     self.isConfig = False
                     self.__dataReady = False
                     self.__profIndex = 0
                     self.byTime = False
                     self.byProfiles = False
                     self.bloques =  None
                     self.bloque0 = None
                     self.index = 0
                     self.buffer = 0
                     self.buffer2 = 0
                     self.buffer3 = 0
                 def setup(self,dataOut,min_hei,max_hei,n, timeInterval,factor_stdv):
                     self.nChannels = dataOut.nChannels
                     self.nProf = dataOut.nProfiles
                     self.nPairs = dataOut.data_cspc.shape[0]
                     self.pairsArray = numpy.array(dataOut.pairsList)
                     self.spectra = dataOut.data_spc
                     self.cspectra = dataOut.data_cspc
                     self.heights = dataOut.heightList #alturas totales
                     self.nHeights =  len(self.heights)
                     self.min_hei = min_hei
                     self.max_hei = max_hei
                     if (self.min_hei == None):
                         self.min_hei = 0
                     if (self.max_hei == None):
                         self.max_hei = dataOut.heightList[-1]
                     self.hval = ((self.max_hei>=self.heights) & (self.heights >= self.min_hei)).nonzero()
                     self.heightsClean = self.heights[self.hval] #alturas filtradas
                     self.hval = self.hval[0] # forma (N,), an solo N elementos -> Indices de alturas
                     self.nHeightsClean = len(self.heightsClean)
                     self.channels = dataOut.channelList
                     self.nChan = len(self.channels)
                     self.nIncohInt = dataOut.nIncohInt
                     self.__initime = dataOut.utctime
                     self.maxAltInd = self.hval[-1]+1
                     self.minAltInd = self.hval[0]
                     self.crosspairs = dataOut.pairsList
                     self.nPairs = len(self.crosspairs)
                     self.normFactor = dataOut.normFactor
                     self.nFFTPoints = dataOut.nFFTPoints
                     self.ippSeconds = dataOut.ippSeconds
                     self.currentTime = self.__initime
                     self.pairsArray = numpy.array(dataOut.pairsList)
                     self.factor_stdv = factor_stdv
                     #print("CHANNELS: ",[x for x in self.channels])
                     if n != None :
                         self.byProfiles = True
                         self.nIntProfiles = n
                     else:
                         self.__integrationtime = timeInterval
                     self.__dataReady = False
                     self.isConfig = True
                 def run(self, dataOut,min_hei=None,max_hei=None, n=None, timeInterval=10,factor_stdv=2.5):
                     #print (dataOut.utctime)
                     if not self.isConfig :
                         #print("Setting config")
                         self.setup(dataOut, min_hei,max_hei,n,timeInterval,factor_stdv)
                         #print("Config Done")
                     tini=dataOut.utctime
                     if self.byProfiles:
                         if self.__profIndex ==  self.nIntProfiles:
                             self.__dataReady = True
                     else:
                         if (tini - self.__initime) >= self.__integrationtime:
                             #print(tini - self.__initime,self.__profIndex)
                             self.__dataReady = True
                             self.__initime = tini
                     #if (tini.tm_min % 2) == 0 and (tini.tm_sec < 5 and self.fint==0):
                     if self.__dataReady:
                         #print("Data ready",self.__profIndex)
                         self.__profIndex = 0
                         jspc = self.buffer
                         jcspc = self.buffer2
                         #jnoise = self.buffer3
                         self.buffer = dataOut.data_spc
                         self.buffer2 = dataOut.data_cspc
                         #self.buffer3 = dataOut.noise
                         self.currentTime = dataOut.utctime
                         if numpy.any(jspc) :
                             #print( jspc.shape, jcspc.shape)
                             jspc = numpy.reshape(jspc,(int(len(jspc)/self.nChannels),self.nChannels,self.nFFTPoints,self.nHeights))
                             jcspc= numpy.reshape(jcspc,(int(len(jcspc)/self.nPairs),self.nPairs,self.nFFTPoints,self.nHeights))
                             self.__dataReady = False
                             #print( jspc.shape, jcspc.shape)
                             dataOut.flagNoData = False
                         else:
                             dataOut.flagNoData = True
                             self.__dataReady = False
                             return dataOut
                     else:
                         #print( len(self.buffer))
                         if numpy.any(self.buffer):
                             self.buffer = numpy.concatenate((self.buffer,dataOut.data_spc), axis=0)
                             self.buffer2 = numpy.concatenate((self.buffer2,dataOut.data_cspc), axis=0)
                             self.buffer3 += dataOut.data_dc
                         else:
                             self.buffer = dataOut.data_spc
                             self.buffer2 = dataOut.data_cspc
                             self.buffer3 = dataOut.data_dc
                         #print self.index, self.fint
                         #print self.buffer2.shape
                         dataOut.flagNoData = True                   ## NOTE: ?? revisar LUEGO
                         self.__profIndex += 1
                         return dataOut                             ## NOTE: REV
                     #index = tini.tm_hour*12+tini.tm_min/5
                     '''REVISAR'''
                     # jspc = jspc/self.nFFTPoints/self.normFactor
                     # jcspc = jcspc/self.nFFTPoints/self.normFactor
                     tmp_spectra,tmp_cspectra = self.cleanRayleigh(dataOut,jspc,jcspc,self.factor_stdv)
                     dataOut.data_spc = tmp_spectra
                     dataOut.data_cspc = tmp_cspectra
                     #dataOut.data_spc,dataOut.data_cspc = self.cleanRayleigh(dataOut,jspc,jcspc,self.factor_stdv)
                     dataOut.data_dc = self.buffer3
                     dataOut.nIncohInt *= self.nIntProfiles
                     dataOut.utctime = self.currentTime             #tiempo promediado
                     #print("Time: ",time.localtime(dataOut.utctime))
                     # dataOut.data_spc = sat_spectra
                     # dataOut.data_cspc = sat_cspectra
                     self.buffer = 0
                     self.buffer2 = 0
                     self.buffer3 = 0
                     return dataOut
                 def cleanRayleigh(self,dataOut,spectra,cspectra,factor_stdv):
                     #print("OP cleanRayleigh")
                     #import matplotlib.pyplot as plt
                     #for k in range(149):
                     #channelsProcssd = []
                     #channelA_ok = False
                     #rfunc = cspectra.copy() #self.bloques
                     rfunc = spectra.copy()
                     #rfunc = cspectra
                     #val_spc = spectra*0.0 #self.bloque0*0.0
                     #val_cspc = cspectra*0.0 #self.bloques*0.0
                     #in_sat_spectra = spectra.copy()  #self.bloque0
                     #in_sat_cspectra = cspectra.copy()  #self.bloques
                     ###ONLY FOR TEST:
                     raxs = math.ceil(math.sqrt(self.nPairs))
                     caxs = math.ceil(self.nPairs/raxs)
                     if self.nPairs <4:
                         raxs = 2
                         caxs = 2
                     #print(raxs, caxs)
                     fft_rev = 14 #nFFT to plot
                     hei_rev = ((self.heights >= 550) & (self.heights <= 551)).nonzero() #hei to plot
                     hei_rev = hei_rev[0]
                     #print(hei_rev)
                     #print numpy.absolute(rfunc[:,0,0,14])
                     gauss_fit, covariance = None, None
                     for ih in range(self.minAltInd,self.maxAltInd):
                         for ifreq in range(self.nFFTPoints):
                             '''
                             ###ONLY FOR TEST:
                             if ifreq ==fft_rev and ih==hei_rev: #TO VIEW A SIGNLE FREQUENCY
                                 fig, axs = plt.subplots(raxs, caxs)
                                 fig2, axs2 = plt.subplots(raxs, caxs)
                             col_ax = 0
                             row_ax = 0
                             '''
                             #print(self.nPairs)
                             for ii in range(self.nChan):  #PARES DE CANALES SELF y CROSS
                                 # if self.crosspairs[ii][1]-self.crosspairs[ii][0] > 1: # APLICAR SOLO EN PARES CONTIGUOS
                                 #     continue
                                 # if not self.crosspairs[ii][0] in channelsProcssd:
                                 #     channelA_ok = True
                                 #print("pair: ",self.crosspairs[ii])
                                 '''
                                 ###ONLY FOR TEST:
                                 if (col_ax%caxs==0 and col_ax!=0 and self.nPairs !=1):
                                     col_ax = 0
                                     row_ax += 1
                                 '''
                                 func2clean = 10*numpy.log10(numpy.absolute(rfunc[:,ii,ifreq,ih]))           #Potencia?
                                 #print(func2clean.shape)
                                 val = (numpy.isfinite(func2clean)==True).nonzero()
                                 if len(val)>0:                                                              #limitador
                                     min_val = numpy.around(numpy.amin(func2clean)-2) #> (-40)
                                     if min_val <= -40 :
                                         min_val = -40
                                     max_val = numpy.around(numpy.amax(func2clean)+2) #< 200
                                     if max_val >= 200 :
                                         max_val = 200
                                     #print min_val, max_val
                                     step = 1
                                     #print("Getting bins and the histogram")
                                     x_dist = min_val + numpy.arange(1 + ((max_val-(min_val))/step))*step
                                     y_dist,binstep = numpy.histogram(func2clean,bins=range(int(min_val),int(max_val+2),step))
                                     #print(len(y_dist),len(binstep[:-1]))
                                     #print(row_ax,col_ax, " ..")
                                     #print(self.pairsArray[ii][0],self.pairsArray[ii][1])
                                     mean = numpy.sum(x_dist * y_dist) / numpy.sum(y_dist)
                                     sigma = numpy.sqrt(numpy.sum(y_dist * (x_dist - mean)**2) / numpy.sum(y_dist))
                                     parg = [numpy.amax(y_dist),mean,sigma]
                                     newY = None
                                     try :
                                         gauss_fit, covariance = curve_fit(fit_func, x_dist, y_dist,p0=parg)
                                         mode = gauss_fit[1]
                                         stdv = gauss_fit[2]
                                         #print(" FIT OK",gauss_fit)
                                         '''
                                         ###ONLY FOR TEST:
                                         if ifreq ==fft_rev and ih==hei_rev: #TO VIEW A SIGNLE FREQUENCY
                                             newY = fit_func(x_dist,gauss_fit[0],gauss_fit[1],gauss_fit[2])
                                             axs[row_ax,col_ax].plot(binstep[:-1],y_dist,color='green')
                                             axs[row_ax,col_ax].plot(binstep[:-1],newY,color='red')
                                             axs[row_ax,col_ax].set_title("CH "+str(self.channels[ii]))
                                         '''
                                     except:
                                         mode = mean
                                         stdv = sigma
                                         #print("FIT FAIL")
                                         #continue
                                     #print(mode,stdv)
                                     #Removing echoes greater than mode + std_factor*stdv
                                     noval = (abs(func2clean - mode)>=(factor_stdv*stdv)).nonzero()
                                     #noval tiene los indices que se van a remover
                                     #print("Chan ",ii," novals: ",len(noval[0]))
                                     if len(noval[0]) > 0: #forma de array (N,) es igual a longitud (N)
                                         novall = ((func2clean - mode) >= (factor_stdv*stdv)).nonzero()
                                         #print(novall)
                                         #print("   ",self.pairsArray[ii])
                                         #cross_pairs = self.pairsArray[ii]
                                             #Getting coherent echoes which are removed.
                                         # if len(novall[0]) > 0:
                                         #
                                         #     val_spc[novall[0],cross_pairs[0],ifreq,ih] = 1
                                         #     val_spc[novall[0],cross_pairs[1],ifreq,ih] = 1
                                         #     val_cspc[novall[0],ii,ifreq,ih] = 1
                                             #print("OUT NOVALL 1")
                                         try:
                                             pair = (self.channels[ii],self.channels[ii + 1])
                                         except:
                                             pair = (99,99)
                                         #print("par ", pair)
                                         if ( pair in self.crosspairs):
                                             q = self.crosspairs.index(pair)
                                             #print("está aqui: ", q, (ii,ii + 1))
                                             new_a = numpy.delete(cspectra[:,q,ifreq,ih], noval[0])
                                             cspectra[noval,q,ifreq,ih] = numpy.mean(new_a) #mean CrossSpectra
                                         #if channelA_ok:
                                         #chA = self.channels.index(cross_pairs[0])
                                         new_b = numpy.delete(spectra[:,ii,ifreq,ih], noval[0])
                                         spectra[noval,ii,ifreq,ih] = numpy.mean(new_b) #mean Spectra Pair A
                                         #channelA_ok = False
                                         # chB = self.channels.index(cross_pairs[1])
                                         # new_c = numpy.delete(spectra[:,chB,ifreq,ih], noval[0])
                                         # spectra[noval,chB,ifreq,ih] = numpy.mean(new_c) #mean Spectra Pair B
                                         #
                                         # channelsProcssd.append(self.crosspairs[ii][0]) # save channel A
                                         # channelsProcssd.append(self.crosspairs[ii][1]) # save channel B
                                     '''
                                     ###ONLY FOR TEST:
                                     if ifreq ==fft_rev and ih==hei_rev: #TO VIEW A SIGNLE FREQUENCY
                                         func2clean = 10*numpy.log10(numpy.absolute(spectra[:,ii,ifreq,ih]))
                                         y_dist,binstep = numpy.histogram(func2clean,bins=range(int(min_val),int(max_val+2),step))
                                         axs2[row_ax,col_ax].plot(binstep[:-1],newY,color='red')
                                         axs2[row_ax,col_ax].plot(binstep[:-1],y_dist,color='green')
                                         axs2[row_ax,col_ax].set_title("CH "+str(self.channels[ii]))
                                     '''
                         '''
                             ###ONLY FOR TEST:
                                 col_ax += 1    #contador de ploteo columnas
                                 ##print(col_ax)
                             ###ONLY FOR TEST:
                             if ifreq ==fft_rev and ih==hei_rev: #TO VIEW A SIGNLE FREQUENCY
                                 title = str(dataOut.datatime)+" nFFT: "+str(ifreq)+"   Alt: "+str(self.heights[ih])+ " km"
                                 title2 = str(dataOut.datatime)+" nFFT: "+str(ifreq)+"   Alt: "+str(self.heights[ih])+ " km  CLEANED"
                                 fig.suptitle(title)
                                 fig2.suptitle(title2)
                                 plt.show()
                         '''
                         ##################################################################################################
                     #print("Getting average of the spectra and cross-spectra from incoherent echoes.")
                     out_spectra = numpy.zeros([self.nChan,self.nFFTPoints,self.nHeights], dtype=float) #+numpy.nan
                     out_cspectra = numpy.zeros([self.nPairs,self.nFFTPoints,self.nHeights], dtype=complex) #+numpy.nan
                     for ih in range(self.nHeights):
                         for ifreq in range(self.nFFTPoints):
                             for ich in range(self.nChan):
                                 tmp = spectra[:,ich,ifreq,ih]
                                 valid = (numpy.isfinite(tmp[:])==True).nonzero()
                                 if len(valid[0]) >0 :
                                     out_spectra[ich,ifreq,ih] = numpy.nansum(tmp)#/len(valid[0])
                             for icr in range(self.nPairs):
                                 tmp = numpy.squeeze(cspectra[:,icr,ifreq,ih])
                                 valid = (numpy.isfinite(tmp)==True).nonzero()
                                 if len(valid[0]) > 0:
                                     out_cspectra[icr,ifreq,ih] = numpy.nansum(tmp)#/len(valid[0])
                     return out_spectra, out_cspectra
                 def REM_ISOLATED_POINTS(self,array,rth):
             #         import matplotlib.pyplot as plt
                     if rth == None :
                         rth = 4
                     #print("REM ISO")
                     num_prof = len(array[0,:,0])
                     num_hei = len(array[0,0,:])
                     n2d = len(array[:,0,0])
                     for ii in range(n2d) :
                         #print ii,n2d
                         tmp = array[ii,:,:]
                         #print tmp.shape, array[ii,101,:],array[ii,102,:]
                       # fig = plt.figure(figsize=(6,5))
                       # left, bottom, width, height = 0.1, 0.1, 0.8, 0.8
                       # ax = fig.add_axes([left, bottom, width, height])
                       # x = range(num_prof)
                       # y = range(num_hei)
                       # cp = ax.contour(y,x,tmp)
                       # ax.clabel(cp, inline=True,fontsize=10)
                       # plt.show()
                         #indxs = WHERE(FINITE(tmp) AND tmp GT 0,cindxs)
                         tmp = numpy.reshape(tmp,num_prof*num_hei)
                         indxs1 = (numpy.isfinite(tmp)==True).nonzero()
                         indxs2 = (tmp > 0).nonzero()
                         indxs1 = (indxs1[0])
                         indxs2 = indxs2[0]
                         #indxs1 = numpy.array(indxs1[0])
                         #indxs2 = numpy.array(indxs2[0])
                         indxs = None
                         #print indxs1 , indxs2
                         for iv in range(len(indxs2)):
                             indv = numpy.array((indxs1 == indxs2[iv]).nonzero())
                             #print len(indxs2), indv
                             if len(indv[0]) > 0  :
                                 indxs =  numpy.concatenate((indxs,indxs2[iv]), axis=None)
               #           print indxs
                         indxs = indxs[1:]
                         #print(indxs, len(indxs))
                         if len(indxs) < 4 :
                             array[ii,:,:] = 0.
                             return
                         xpos = numpy.mod(indxs ,num_hei)
                         ypos = (indxs / num_hei)
                         sx = numpy.argsort(xpos) # Ordering respect to "x" (time)
                         #print sx
                         xpos = xpos[sx]
                         ypos = ypos[sx]
                  # *********************************** Cleaning isolated points **********************************
                         ic = 0
                         while True :
                             r = numpy.sqrt(list(numpy.power((xpos[ic]-xpos),2)+ numpy.power((ypos[ic]-ypos),2)))
                             #no_coh = WHERE(FINITE(r) AND (r LE rth),cno_coh)
                             #plt.plot(r)
                             #plt.show()
                             no_coh1 = (numpy.isfinite(r)==True).nonzero()
                             no_coh2 = (r <= rth).nonzero()
                             #print r, no_coh1, no_coh2
                             no_coh1 = numpy.array(no_coh1[0])
                             no_coh2 = numpy.array(no_coh2[0])
                             no_coh = None
                             #print valid1 , valid2
                             for iv in range(len(no_coh2)):
                                 indv = numpy.array((no_coh1 == no_coh2[iv]).nonzero())
                                 if len(indv[0]) > 0  :
                                     no_coh =  numpy.concatenate((no_coh,no_coh2[iv]), axis=None)
                             no_coh = no_coh[1:]
                             #print len(no_coh), no_coh
                             if len(no_coh) < 4 :
                                 #print xpos[ic], ypos[ic], ic
                 #                plt.plot(r)
                 #                plt.show()
                                 xpos[ic] = numpy.nan
                                 ypos[ic] = numpy.nan
                             ic = ic + 1
                             if  (ic == len(indxs)) :
                                 break
                             #print( xpos, ypos)
                         indxs = (numpy.isfinite(list(xpos))==True).nonzero()
                         #print indxs[0]
                         if len(indxs[0]) < 4 :
                             array[ii,:,:] = 0.
                             return
                         xpos = xpos[indxs[0]]
                         ypos = ypos[indxs[0]]
                         for i in range(0,len(ypos)):
                             ypos[i]=int(ypos[i])
                         junk = tmp
                         tmp = junk*0.0
                         tmp[list(xpos + (ypos*num_hei))] = junk[list(xpos + (ypos*num_hei))]
                         array[ii,:,:] = numpy.reshape(tmp,(num_prof,num_hei))
                         #print array.shape
                         #tmp = numpy.reshape(tmp,(num_prof,num_hei))
                         #print tmp.shape
                         # fig = plt.figure(figsize=(6,5))
                         # left, bottom, width, height = 0.1, 0.1, 0.8, 0.8
                         # ax = fig.add_axes([left, bottom, width, height])
                         # x = range(num_prof)
                         # y = range(num_hei)
                         # cp = ax.contour(y,x,array[ii,:,:])
                         # ax.clabel(cp, inline=True,fontsize=10)
                         # plt.show()
                     return array
             class IntegrationFaradaySpectra(Operation):
                 __profIndex = 0
                 __withOverapping = False
                 __byTime = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer_spc = None
                 __buffer_cspc = None
                 __buffer_dc = None
                 __dataReady = False
                 __timeInterval = None
                 n = None
                 def __init__(self):
                     Operation.__init__(self)
                 def setup(self, dataOut,n=None, timeInterval=None, overlapping=False, DPL=None):
                     """
                     Set the parameters of the integration class.
                     Inputs:
                         n        :    Number of coherent integrations
                         timeInterval   :    Time of integration. If the parameter "n" is selected this one does not work
                         overlapping    :
                     """
                     self.__initime = None
                     self.__lastdatatime = 0
                     self.__buffer_spc = []
                     self.__buffer_cspc = []
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     self.__dataReady = False
                     self.__byTime = False
                     #self.ByLags = dataOut.ByLags   ###REDEFINIR
                     self.ByLags = False
                     if DPL != None:
                         self.DPL=DPL
                     else:
                         #self.DPL=dataOut.DPL    ###REDEFINIR
                         self.DPL=0
                     if n is None and timeInterval is None:
                         raise ValueError("n or timeInterval should be specified ...")
                     if n is not None:
                         self.n = int(n)
                     else:
                         self.__integrationtime = int(timeInterval)
                         self.n = None
                         self.__byTime = True
                 def putData(self, data_spc, data_cspc, data_dc):
                     """
                     Add a profile to the __buffer_spc and increase in one the __profileIndex
                     """
                     self.__buffer_spc.append(data_spc)
                     if data_cspc is None:
                         self.__buffer_cspc = None
                     else:
                         self.__buffer_cspc.append(data_cspc)
                     if data_dc is None:
                         self.__buffer_dc = None
                     else:
                         self.__buffer_dc += data_dc
                     self.__profIndex += 1
                     return
                 def hildebrand_sekhon_Integration(self,data,navg):
                     sortdata = numpy.sort(data, axis=None)
                     sortID=data.argsort()
                     lenOfData = len(sortdata)
                     nums_min = lenOfData*0.75
                     if nums_min <= 5:
                         nums_min = 5
                     sump = 0.
                     sumq = 0.
                     j = 0
                     cont = 1
                     while((cont == 1)and(j < lenOfData)):
                         sump += sortdata[j]
                         sumq += sortdata[j]**2
                         if j > nums_min:
                             rtest = float(j)/(j-1) + 1.0/navg
                             if ((sumq*j) > (rtest*sump**2)):
                                 j = j - 1
                                 sump = sump - sortdata[j]
                                 sumq = sumq - sortdata[j]**2
                                 cont = 0
                         j += 1
                     #lnoise = sump / j
                     return j,sortID
                 def pushData(self):
                     """
                     Return the sum of the last profiles and the profiles used in the sum.
                     Affected:
                     self.__profileIndex
                     """
                     bufferH=None
                     buffer=None
                     buffer1=None
                     buffer_cspc=None
                     self.__buffer_spc=numpy.array(self.__buffer_spc)
                     self.__buffer_cspc=numpy.array(self.__buffer_cspc)
                     freq_dc = int(self.__buffer_spc.shape[2] / 2)
                     #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
                     for k in range(7,self.nHeights):
                         buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
                         outliers_IDs_cspc=[]
                         cspc_outliers_exist=False
                         #print("AQUIII")
                         for i in range(self.nChannels):#dataOut.nChannels):
                             buffer1=numpy.copy(self.__buffer_spc[:,i,:,k])
                             indexes=[]
                             #sortIDs=[]
                             outliers_IDs=[]
                             for j in range(self.nProfiles):
                                 # if i==0 and j==freq_dc: #NOT CONSIDERING DC PROFILE AT CHANNEL 0
                                 #     continue
                                 # if i==1 and j==0: #NOT CONSIDERING DC PROFILE AT CHANNEL 1
                                 #     continue
                                 buffer=buffer1[:,j]
                                 index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
                                 indexes.append(index)
                                 #sortIDs.append(sortID)
                                 outliers_IDs=numpy.append(outliers_IDs,sortID[index:])
                             outliers_IDs=numpy.array(outliers_IDs)
                             outliers_IDs=outliers_IDs.ravel()
                             outliers_IDs=numpy.unique(outliers_IDs)
                             outliers_IDs=outliers_IDs.astype(numpy.dtype('int64'))
                             indexes=numpy.array(indexes)
                             indexmin=numpy.min(indexes)
                             if indexmin != buffer1.shape[0]:
                                 cspc_outliers_exist=True
                                 ###sortdata=numpy.sort(buffer1,axis=0)
                                 ###avg2=numpy.mean(sortdata[:indexmin,:],axis=0)
                                 lt=outliers_IDs
                                 avg=numpy.mean(buffer1[[t for t in range(buffer1.shape[0]) if t not in lt],:],axis=0)
                                 for p in list(outliers_IDs):
                                     buffer1[p,:]=avg
                             self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
                             ###cspc IDs
                             #indexmin_cspc+=indexmin_cspc
                             outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
                         #if not breakFlag:
                         outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
                         if cspc_outliers_exist:
                             #sortdata=numpy.sort(buffer_cspc,axis=0)
                             #avg=numpy.mean(sortdata[:indexmin_cpsc,:],axis=0)
                             lt=outliers_IDs_cspc
                             avg=numpy.mean(buffer_cspc[[t for t in range(buffer_cspc.shape[0]) if t not in lt],:],axis=0)
                             for p in list(outliers_IDs_cspc):
                                 buffer_cspc[p,:]=avg
                         self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
                         #else:
                             #break
                     buffer=None
                     bufferH=None
                     buffer1=None
                     buffer_cspc=None
                     #print("cpsc",self.__buffer_cspc[:,0,0,0,0])
                     #print(self.__profIndex)
                     #exit()
                     buffer=None
                     #print(self.__buffer_spc[:,1,3,20,0])
                     #print(self.__buffer_spc[:,1,5,37,0])
                     data_spc = numpy.sum(self.__buffer_spc,axis=0)
                     data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
                     #print(numpy.shape(data_spc))
                     #data_spc[1,4,20,0]=numpy.nan
                     #data_cspc = self.__buffer_cspc
                     data_dc = self.__buffer_dc
                     n = self.__profIndex
                     self.__buffer_spc = []
                     self.__buffer_cspc = []
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     return data_spc, data_cspc, data_dc, n
                 def byProfiles(self, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if self.__profIndex == self.n:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def byTime(self, datatime, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if (datatime - self.__initime) >= self.__integrationtime:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def integrate(self, datatime, *args):
                     if self.__profIndex == 0:
                         self.__initime = datatime
                     if self.__byTime:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
                             datatime, *args)
                     else:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
                     if not self.__dataReady:
                         return None, None, None, None
                     return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
                 def run(self, dataOut, n=None, DPL = None,timeInterval=None, overlapping=False):
                     if n == 1:
                         return dataOut
                     dataOut.flagNoData = True
                     if not self.isConfig:
                         self.setup(dataOut, n, timeInterval, overlapping,DPL )
                         self.isConfig = True
                     if not self.ByLags:
                         self.nProfiles=dataOut.nProfiles
                         self.nChannels=dataOut.nChannels
                         self.nHeights=dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
                                                                                             dataOut.data_spc,
                                                                                             dataOut.data_cspc,
                                                                                             dataOut.data_dc)
                     else:
                         self.nProfiles=dataOut.nProfiles
                         self.nChannels=dataOut.nChannels
                         self.nHeights=dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
                                                                                             dataOut.dataLag_spc,
                                                                                             dataOut.dataLag_cspc,
                                                                                             dataOut.dataLag_dc)
                     if self.__dataReady:
                         if not self.ByLags:
                             dataOut.data_spc = numpy.squeeze(avgdata_spc)
                             dataOut.data_cspc = numpy.squeeze(avgdata_cspc)
                             dataOut.data_dc = avgdata_dc
                         else:
                             dataOut.dataLag_spc = avgdata_spc
                             dataOut.dataLag_cspc = avgdata_cspc
                             dataOut.dataLag_dc = avgdata_dc
                             dataOut.data_spc=dataOut.dataLag_spc[:,:,:,dataOut.LagPlot]
                             dataOut.data_cspc=dataOut.dataLag_cspc[:,:,:,dataOut.LagPlot]
                             dataOut.data_dc=dataOut.dataLag_dc[:,:,dataOut.LagPlot]
                         dataOut.nIncohInt *= self.n
                         dataOut.utctime = avgdatatime
                         dataOut.flagNoData = False
                     return dataOut
             class removeInterference(Operation):
                 def removeInterference2(self):
                     cspc = self.dataOut.data_cspc
                     spc = self.dataOut.data_spc
                     Heights = numpy.arange(cspc.shape[2])
                     realCspc = numpy.abs(cspc)
                     for i in range(cspc.shape[0]):
                         LinePower= numpy.sum(realCspc[i], axis=0)
                         Threshold = numpy.amax(LinePower)-numpy.sort(LinePower)[len(Heights)-int(len(Heights)*0.1)]
                         SelectedHeights = Heights[ numpy.where( LinePower < Threshold ) ]
                         InterferenceSum = numpy.sum( realCspc[i,:,SelectedHeights], axis=0 )
                         InterferenceThresholdMin = numpy.sort(InterferenceSum)[int(len(InterferenceSum)*0.98)]
                         InterferenceThresholdMax = numpy.sort(InterferenceSum)[int(len(InterferenceSum)*0.99)]
                         InterferenceRange = numpy.where( ([InterferenceSum > InterferenceThresholdMin]))# , InterferenceSum < InterferenceThresholdMax]) )
                         #InterferenceRange = numpy.where( ([InterferenceRange < InterferenceThresholdMax]))
                         if len(InterferenceRange)<int(cspc.shape[1]*0.3):
                             cspc[i,InterferenceRange,:] = numpy.NaN
                     self.dataOut.data_cspc = cspc
                 def removeInterference(self, interf = 2, hei_interf = None, nhei_interf = None, offhei_interf = None):
                     jspectra = self.dataOut.data_spc
                     jcspectra = self.dataOut.data_cspc
                     jnoise = self.dataOut.getNoise()
                     num_incoh = self.dataOut.nIncohInt
                     num_channel = jspectra.shape[0]
                     num_prof = jspectra.shape[1]
                     num_hei = jspectra.shape[2]
                     # hei_interf
                     if hei_interf is None:
                         count_hei = int(num_hei / 2)
                         hei_interf = numpy.asmatrix(list(range(count_hei))) + num_hei - count_hei
                         hei_interf = numpy.asarray(hei_interf)[0]
                     # nhei_interf
                     if (nhei_interf == None):
                         nhei_interf = 5
                     if (nhei_interf < 1):
                         nhei_interf = 1
                     if (nhei_interf > count_hei):
                         nhei_interf = count_hei
                     if (offhei_interf == None):
                         offhei_interf = 0
                     ind_hei = list(range(num_hei))
             #         mask_prof = numpy.asarray(range(num_prof - 2)) + 1
             #         mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
                     mask_prof = numpy.asarray(list(range(num_prof)))
                     num_mask_prof = mask_prof.size
                     comp_mask_prof = [0, num_prof / 2]
                     # noise_exist:    Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
                     if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
                         jnoise = numpy.nan
                     noise_exist = jnoise[0] < numpy.Inf
                     # Subrutina de Remocion de la Interferencia
                     for ich in range(num_channel):
                         # Se ordena los espectros segun su potencia (menor a mayor)
                         power = jspectra[ich, mask_prof, :]
                         power = power[:, hei_interf]
                         power = power.sum(axis=0)
                         psort = power.ravel().argsort()
                         # Se estima la interferencia promedio en los Espectros de Potencia empleando
                         junkspc_interf = jspectra[ich, :, hei_interf[psort[list(range(
                             offhei_interf, nhei_interf + offhei_interf))]]]
                         if noise_exist:
                             #    tmp_noise = jnoise[ich] / num_prof
                             tmp_noise = jnoise[ich]
                         junkspc_interf = junkspc_interf - tmp_noise
                         #junkspc_interf[:,comp_mask_prof] = 0
                         jspc_interf = junkspc_interf.sum(axis=0) / nhei_interf
                         jspc_interf = jspc_interf.transpose()
                         # Calculando el espectro de interferencia promedio
                         noiseid = numpy.where(
                             jspc_interf <= tmp_noise / numpy.sqrt(num_incoh))
                         noiseid = noiseid[0]
                         cnoiseid = noiseid.size
                         interfid = numpy.where(
                             jspc_interf > tmp_noise / numpy.sqrt(num_incoh))
                         interfid = interfid[0]
                         cinterfid = interfid.size
                         if (cnoiseid > 0):
                             jspc_interf[noiseid] = 0
                         # Expandiendo los perfiles a limpiar
                         if (cinterfid > 0):
                             new_interfid = (
                                 numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof) % num_prof
                             new_interfid = numpy.asarray(new_interfid)
                             new_interfid = {x for x in new_interfid}
                             new_interfid = numpy.array(list(new_interfid))
                             new_cinterfid = new_interfid.size
                         else:
                             new_cinterfid = 0
                         for ip in range(new_cinterfid):
                             ind = junkspc_interf[:, new_interfid[ip]].ravel().argsort()
                             jspc_interf[new_interfid[ip]
                                         ] = junkspc_interf[ind[nhei_interf // 2], new_interfid[ip]]
                         jspectra[ich, :, ind_hei] = jspectra[ich, :,
                                                              ind_hei] - jspc_interf  # Corregir indices
                         # Removiendo la interferencia del punto de mayor interferencia
                         ListAux = jspc_interf[mask_prof].tolist()
                         maxid = ListAux.index(max(ListAux))
                         if cinterfid > 0:
                             for ip in range(cinterfid * (interf == 2) - 1):
                                 ind = (jspectra[ich, interfid[ip], :] < tmp_noise *
                                        (1 + 1 / numpy.sqrt(num_incoh))).nonzero()
                                 cind = len(ind)
                                 if (cind > 0):
                                     jspectra[ich, interfid[ip], ind] = tmp_noise * \
                                         (1 + (numpy.random.uniform(cind) - 0.5) /
                                          numpy.sqrt(num_incoh))
                             ind = numpy.array([-2, -1, 1, 2])
                             xx = numpy.zeros([4, 4])
                             for id1 in range(4):
                                 xx[:, id1] = ind[id1]**numpy.asarray(list(range(4)))
                             xx_inv = numpy.linalg.inv(xx)
                             xx = xx_inv[:, 0]
                             ind = (ind + maxid + num_mask_prof) % num_mask_prof
                             yy = jspectra[ich, mask_prof[ind], :]
                             jspectra[ich, mask_prof[maxid], :] = numpy.dot(
                                 yy.transpose(), xx)
                         indAux = (jspectra[ich, :, :] < tmp_noise *
                                   (1 - 1 / numpy.sqrt(num_incoh))).nonzero()
                         jspectra[ich, indAux[0], indAux[1]] = tmp_noise * \
                             (1 - 1 / numpy.sqrt(num_incoh))
                     # Remocion de Interferencia en el Cross Spectra
                     if jcspectra is None:
                         return jspectra, jcspectra
                     num_pairs = int(jcspectra.size / (num_prof * num_hei))
                     jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
                     for ip in range(num_pairs):
                         #-------------------------------------------
                         cspower = numpy.abs(jcspectra[ip, mask_prof, :])
                         cspower = cspower[:, hei_interf]
                         cspower = cspower.sum(axis=0)
                         cspsort = cspower.ravel().argsort()
                         junkcspc_interf = jcspectra[ip, :, hei_interf[cspsort[list(range(
                             offhei_interf, nhei_interf + offhei_interf))]]]
                         junkcspc_interf = junkcspc_interf.transpose()
                         jcspc_interf = junkcspc_interf.sum(axis=1) / nhei_interf
                         ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
                         median_real = int(numpy.median(numpy.real(
                             junkcspc_interf[mask_prof[ind[list(range(3 * num_prof // 4))]], :])))
                         median_imag = int(numpy.median(numpy.imag(
                             junkcspc_interf[mask_prof[ind[list(range(3 * num_prof // 4))]], :])))
                         comp_mask_prof = [int(e) for e in comp_mask_prof]
                         junkcspc_interf[comp_mask_prof, :] = numpy.complex(
                             median_real, median_imag)
                         for iprof in range(num_prof):
                             ind = numpy.abs(junkcspc_interf[iprof, :]).ravel().argsort()
                             jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf // 2]]
                         # Removiendo la Interferencia
                         jcspectra[ip, :, ind_hei] = jcspectra[ip,
                                                               :, ind_hei] - jcspc_interf
                         ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
                         maxid = ListAux.index(max(ListAux))
                         ind = numpy.array([-2, -1, 1, 2])
                         xx = numpy.zeros([4, 4])
                         for id1 in range(4):
                             xx[:, id1] = ind[id1]**numpy.asarray(list(range(4)))
                         xx_inv = numpy.linalg.inv(xx)
                         xx = xx_inv[:, 0]
                         ind = (ind + maxid + num_mask_prof) % num_mask_prof
                         yy = jcspectra[ip, mask_prof[ind], :]
                         jcspectra[ip, mask_prof[maxid], :] = numpy.dot(yy.transpose(), xx)
                     # Guardar Resultados
                     self.dataOut.data_spc = jspectra
                     self.dataOut.data_cspc = jcspectra
                     return 1
                 def run(self, dataOut, interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None, mode=1):
                     self.dataOut = dataOut
                     if mode == 1:
                         self.removeInterference(interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None)
                     elif mode == 2:
                         self.removeInterference2()
                     return self.dataOut
             class IncohInt(Operation):
                 __profIndex = 0
                 __withOverapping = False
                 __byTime = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer_spc = None
                 __buffer_cspc = None
                 __buffer_dc = None
                 __dataReady = False
                 __timeInterval = None
                 n = None
                 def __init__(self):
                     Operation.__init__(self)
                 def setup(self, n=None, timeInterval=None, overlapping=False):
                     """
                     Set the parameters of the integration class.
                     Inputs:
                         n        :    Number of coherent integrations
                         timeInterval   :    Time of integration. If the parameter "n" is selected this one does not work
                         overlapping    :
                     """
                     self.__initime = None
                     self.__lastdatatime = 0
                     self.__buffer_spc = 0
                     self.__buffer_cspc = 0
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     self.__dataReady = False
                     self.__byTime = False
                     if n is None and timeInterval is None:
                         raise ValueError("n or timeInterval should be specified ...")
                     if n is not None:
                         self.n = int(n)
                     else:
                         self.__integrationtime = int(timeInterval)
                         self.n = None
                         self.__byTime = True
                 def putData(self, data_spc, data_cspc, data_dc):
                     """
                     Add a profile to the __buffer_spc and increase in one the __profileIndex
                     """
                     self.__buffer_spc += data_spc
                     if data_cspc is None:
                         self.__buffer_cspc = None
                     else:
                         self.__buffer_cspc += data_cspc
                     if data_dc is None:
                         self.__buffer_dc = None
                     else:
                         self.__buffer_dc += data_dc
                     self.__profIndex += 1
                     return
                 def pushData(self):
                     """
                     Return the sum of the last profiles and the profiles used in the sum.
                     Affected:
                     self.__profileIndex
                     """
                     data_spc = self.__buffer_spc
                     data_cspc = self.__buffer_cspc
                     data_dc = self.__buffer_dc
                     n = self.__profIndex
                     self.__buffer_spc = 0
                     self.__buffer_cspc = 0
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     return data_spc, data_cspc, data_dc, n
                 def byProfiles(self, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if self.__profIndex == self.n:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def byTime(self, datatime, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if (datatime - self.__initime) >= self.__integrationtime:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def integrate(self, datatime, *args):
                     if self.__profIndex == 0:
                         self.__initime = datatime
                     if self.__byTime:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
                             datatime, *args)
                     else:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
                     if not self.__dataReady:
                         return None, None, None, None
                     return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
                 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
                     if n == 1:
                         return dataOut
                     dataOut.flagNoData = True
                     if not self.isConfig:
                         self.setup(n, timeInterval, overlapping)
                         self.isConfig = True
                     avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
                                                                                         dataOut.data_spc,
                                                                                         dataOut.data_cspc,
                                                                                         dataOut.data_dc)
                     if self.__dataReady:
                         dataOut.data_spc = avgdata_spc
                         dataOut.data_cspc = avgdata_cspc
                         dataOut.data_dc = avgdata_dc
                         dataOut.nIncohInt *= self.n
                         dataOut.utctime = avgdatatime
                         dataOut.flagNoData = False
                     return dataOut
             class dopplerFlip(Operation):
                 def run(self, dataOut):
                     # arreglo 1: (num_chan, num_profiles, num_heights)
                     self.dataOut = dataOut
                     # JULIA-oblicua, indice 2
                     # arreglo 2: (num_profiles, num_heights)
                     jspectra = self.dataOut.data_spc[2]
                     jspectra_tmp = numpy.zeros(jspectra.shape)
                     num_profiles = jspectra.shape[0]
                     freq_dc = int(num_profiles / 2)
                     # Flip con for
                     for j in range(num_profiles):
                         jspectra_tmp[num_profiles-j-1]= jspectra[j]
                     # Intercambio perfil de DC con perfil inmediato anterior
                     jspectra_tmp[freq_dc-1]= jspectra[freq_dc-1]
                     jspectra_tmp[freq_dc]= jspectra[freq_dc]
                     # canal modificado es re-escrito en el arreglo de canales
                     self.dataOut.data_spc[2] = jspectra_tmp
                     return self.dataOut

General Comments 0

Write
Preview

You need to be logged in to leave comments. Login now

No TODOs yet

	Site-wide shortcuts
/	Use quick search box
g h	Goto home page
g g	Goto my private gists page
g G	Goto my public gists page
g 0-9	Goto bookmarked items from 0-9
n r	New repository page
n g	New gist page

	Repositories
g s	Goto summary page
g c	Goto changelog page
g f	Goto files page
g F	Goto files page with file search activated
g p	Goto pull requests page
g o	Goto repository settings
g O	Goto repository access permissions settings
t s	Toggle sidebar on some pages